Magnetic Mirror

Examples - Electromagnetism

Last edited: February 18th 2020


Introduction

A charged particle moving in a magnetic field will experience a magnetic force governed by Lorentz force law. This force is always perpendicular to the direction of motion of the particle. If the particle moves in uniform magnetic field and the velocity of the particle is normal to the direction of the field, the particle end up moving in a circle. Now, if the direction of the field is not normal to the velocity, the particle starts moving in a helix, unless its motion is completely parallel to the field in which the there is no Lorentz force and the particle moves in a perfectly straight line. Additionally, with a non-uniform field, you may be able to trap the helical moving particle in a so called magnetic mirror. This may happen if the particle is traveling from a region of weak magnetic field into a region with strong magnetic field.

In this example, we will construct a non-uniform magnetic field and trap the motion of the charged particle in a magnetic mirror. If you are familiar with our notebook Uniform magnetic field, you may skip some of the theory section, as its the same in both notebooks.

Theory

The Lorentz force law

A particle with charge $Q$ and mass $m_p$ moving with a velocity $\mathbf{v}$ in a magnetic field $\mathbf{B}$ experiences a force

\begin{equation} \mathbf{F}_{mag} = Q(\mathbf{v} \times \mathbf{B}) \label{Lorentz} \end{equation}

known as the Lorentz force [1].

Discretisation

The Lorentz force can be written as a second order ordinary differential equation (ODE)

\begin{equation} m_p\frac{d^2\mathbf{x}}{dt^2} = Q\big(\frac{d\mathbf{x}}{dt} \times \mathbf{B}\big), \label{LorentzODE} \end{equation}

where the velocity components are written as the time derivative of the positional components $v_i = dx_i/dt$. By decomposing the vectors, \eqref{LorentzODE} can be reduced into three sets of two first order ODE's each.

\begin{equation} \begin{aligned} &\frac{dx}{dt} = v_x\\[4pt] &\frac{dv_x}{dt} = \frac{Q}{m_{p}}(v_yB_z - v_zB_y)\\[10pt] &\frac{dy}{dt} = v_y\\[4pt] &\frac{dv_y}{dt}= \frac{Q}{m_{p}}(v_zB_x - v_xB_z)\\[10pt] &\frac{dz}{dt} = v_z\\[4pt] &\frac{dv_z}{dt} = \frac{Q}{m_{p}}(v_xB_y - v_yB_x) . \label{LorentzComponents} \end{aligned} \end{equation}

Using a numerical method such as the Euler method, Runge-Kutta fourth order, or an embedded Runge-Kutta pair, these equations can be numerically integrated using a specified time step, $\Delta t$, and the motion of a charged particle in a magnetic field can be found.

A non-uniform magnetic field

To to trap a charged particle in a mirror we want a magnetic field which mainly points, and whose magnitude varies, in the $z$-direction. A quadratic function of the direction of the field, plus a radial component is typical [3]

\begin{equation} \mathbf{B}(x, y, z) = B_r\hat{\mathbf{r}} + B_0 \big(1 + \frac{z^2}{L^2} \big)\hat{\mathbf{z}} \; . \label{field} \end{equation}

Here the mangetic field is on cylindrical form and $L$ is a parameter controlling the mirror length in $z$-direction. This field is axisymmetric (cylindrical symmetry), i.e. the azimuthal component $B_\theta = 0$ and $\partial/\partial\theta = 0$. From Maxwell's equations $\nabla \cdot \mathbf{B} = 0$ we derive (using the nabla operator corresponding to cylindrical coordinate frames)

\begin{equation} \begin{aligned} \nabla \cdot \mathbf{B} &= \frac{1}{r}\frac{\partial(rB_r)}{\partial r} + \frac{1}{r}\frac{\partial B_\theta}{\partial \theta} + \frac{\partial B_z}{\partial z} \\[4pt] &= \frac{1}{r}\frac{\partial(rB_r)}{\partial r} + \frac{\partial B_z}{\partial z} = 0 \end{aligned} \end{equation}

Carrying out the partial derivation of $B_z$ we get

\begin{equation} \begin{aligned} \frac{1}{r}\frac{\partial(rB_r)}{\partial r} &= -2B_0\frac{z}{L^2} \\[4pt] rB_r &= -\int_0^r 2B_0\frac{z}{L^2} r dr \\[4pt] \implies B_r &= -rB_0\frac{z}{L^2} , \end{aligned} \end{equation}

leading to the magnetic field

\begin{equation} \mathbf{B}(x, y, z) = -xB_0\frac{z}{L^2}\hat{\mathbf{x}} -yB_0\frac{z}{L^2}\hat{\mathbf{y}} + B_0 \big(1 + \frac{z^2}{L^2} \big)\hat{\mathbf{z}} , \label{BCartesian} \end{equation}

on cartesian form.

Magnetic moment

The usual definition for the magnetic moment of a current loop with a current $I$ and area $A$ is [1]

$$ \mu = IA.$$

An electron moving around in a circle of radius $r$, completing a loop with a frequency $f = \omega/2\pi$, generates a current $I = e\omega/2\pi$. Here $\omega$ is the angular frequency of the movement. With the area $A = \pi r^2 = \pi v^2\omega^2$, where $v$ is the velocity of the electron, the magnetic moment can be written as

\begin{equation} \mu = \frac{1}{2}\frac{v^2e}{\omega} \label{magneticMoment} \end{equation}

Larmor radius

As stated in the introduction, a charged particle will move in a circle if it has a velocity component perpendicular to an external magnetic field. This is called Larmor gyration, and the gyroradius (or Larmor radius), the radius of the circle, is given by [2]

\begin{equation} r_L = \frac{mv_{\perp}}{|Q| B} , \label{larmorRadius} \end{equation}

where $v_{\perp}$ is the velocity of the particle perpendicular to the magnetic field with magnitude $B$.

Plotting magnetic field

In [21]:
import numpy as np
import matplotlib.pyplot as plt
import sys
from mpl_toolkits.mplot3d import Axes3D
plt.style.use("bmh")
%matplotlib inline
import progressbar
import sys
import time
figsize = (6, 6)
dpi = 600

To begin with, we plot the magnetic field itself. To plot the streamlines of the field in the different planes, a function which returns the magnetic field vector at a given point is needed.

In [22]:
def magneticMirror(x, y, z, B_0, L):
    """Calculates and returns the magnetic field components
    given the cartesian position of the particle and B_0 and L."""
    # Component field strength in cartesian coordinates
    B_x = -x*B_0*z/L**2
    B_y = -y*B_0*z/L**2
    B_z = B_0*(1 + z**2/L**2)
    
    # Magnetic vector field matrix
    B = np.array([B_x, B_y, B_z])
    
    return B

Using a meshgrid, this function will return the magnetic vector field in any plane of interest.

In [23]:
B_0 = 30     # The magnetic field magnitude at origin
L = 0.4      # Mirror length parameter
a = 2*L      # To be used for grid creation

# Define meshgrid limits
[xmin, xmax, ymin, ymax, zmin, zmax] = [-a, a, -a, a, -a, a]
# Grid resolution
n = 100j

# 2D meshgrid for creating streamplot
y, x = np.mgrid[xmin:xmax:n, ymin:ymax:n]

def plotMagneticFieldInPlane(subplot, plane, xlabel, ylabel, abscissa, ordinate, B_a, B_o):
    """Function for streamplotting (in a specified subplot)the magnetic 
    field in a given plane with given grid and field components.
    Parameters:
        subplot        int; the position of the subplot
        plane          string; which plane is plottet? to be used in
                       in the title.
        xlabel         string; label of the subplot's abscissa.
        ylabel         string; label of the subplot's ordinate.
        abscissa       nxn array; the "x-values" of the mesh grid points.
        ordinate       nxn array; the "y-values" of the mesh grid points.
        B_a            nxn array; the magnetic field components in the direction
                       of the abscissa.
        B_o            nxn array; the magnetic field components in the direction
                       of the ordinate."""
    
    plt.subplot(subplot)
    plt.title("Magnetic field in the %s-plane"%(plane))
    # Adding linewidth
    # Thicker streamline represent higher magnetic field magnitude
    magnitude = np.sqrt(B_a**2 + B_o**2) # Absolute magnitude of the magnetic field
    lw = 5*magnitude/np.amax(magnitude)
    plt.streamplot(abscissa, ordinate, B_a, B_o, linewidth=lw)
    plt.plot(0, 0, 'o', color="blue")
    plt.xlabel(xlabel)
    plt.ylabel(ylabel)

plt.figure(figsize=(12,10), dpi=200)
# xy-plane
B = magneticMirror(x, y, x*0, B_0, L)
plotMagneticFieldInPlane(221, 'xy', r"$x$", r"$y$", x, y, B[0], B[1])
# yz-plane
# In case of unequal axis limits, re-create meshgrid
z, y = np.mgrid[ymin:ymax:n, zmin:zmax:n]
B    = magneticMirror(y*0, y, z, B_0, L)
plotMagneticFieldInPlane(222, 'yz', r"$y$", r"$z$", y, z, B[1], B[2])
# xz-plane
z, x = np.mgrid[xmin:xmax:n, zmin:zmax:n]
B    = magneticMirror(x, x*0, z, B_0, L)
plotMagneticFieldInPlane(223, 'xz', r"$x$", r"$z$", x, z, B[0], B[2])
# xy-plane with z=1
y, x = np.mgrid[xmin:xmax:n, ymin:ymax:n]
z    = np.ones(np.shape(x))
B    = magneticMirror(x, y, z, B_0, L)
plotMagneticFieldInPlane(224, 'xy', r"$x$", r"$y$", x, y, B[0], B[1])

plt.show()
/Users/eilifso/anaconda3/envs/numfys/lib/python3.7/site-packages/ipykernel_launcher.py:34: RuntimeWarning: invalid value encountered in true_divide