Notes on single particle motion II

                Equipo de Propulsión Espacial y Plasmas, Mario Merino, 2021

Variable fields

We will consider only cases with:

ℓ≪L;v∥/L≪ωc;∣∂E/∂t∣/E≪ωc;∣∂B/∂t∣/B≪ωc\ell \ll L; \quad v_\parallel / L \ll \omega_c; \quad |\partial \bm E/\partial t|/E \ll \omega_c; \quad |\partial \bm B/\partial t|/B \ll \omega_c

where LL is the gradient length of changes in E\bm E and B\bm B (L∼E/∣∇E∣∼B/∣∇B∣L\sim E/|\nabla \bm E|\sim B/|\nabla \bm B|).
This ensures that the particle does not "see" large changes in E\bm E, B\bm B within an orbit.

We introduce the small parameter ε=ℓ/L\varepsilon = \ell / L and rewrite m/q→εm^/q^m/q \to \varepsilon \hat m/\hat q, where m^,q^\hat m,\hat q are arbitrary normalization constants:

εdvdt=±q^m^(E+v×B);drdt=v\varepsilon\frac{d \bm v}{dt} = \pm \frac{\hat q}{\hat m}(\bm E + \bm v \times \bm B); \quad \frac{d \bm r}{d t} = \bm v

  • We will study the limit ε→0+\varepsilon \to 0^+ and drop the hat symbols for simplicity

Two-scale averaging method

There are two distinct timescales in the problem: τc=1/ωc\tau_c = 1/\omega_c for the fast changes of the gyrophase φ\varphi, and τ\tau for the slow changes of the position of the gyrocenter:

τc≪τ⇒τcτ∼ε≪1\tau_c \ll \tau \quad\Rightarrow\quad \frac{\tau_c}{\tau} \sim \varepsilon \ll 1

  • While the gyrophase is a function of time, φ=φ(t)\varphi = \varphi(t), we will write all variables as functions of tt and φ\varphi as if they were two independent variables:

r(t)→r(t,φ(t))→r(t,φ)\bm r(t) \quad \rightarrow \quad \bm r(t,\varphi(t)) \quad \rightarrow \quad \bm r(t,\varphi)

  • We will compute derivatives as: drdt=∂r∂t+φ˙∂r∂φ\dfrac{d \bm r}{dt} = \dfrac{\partial \bm r}{\partial t} + \dot\varphi \dfrac{\partial \bm r}{\partial \varphi}, with φ˙≈ωc=O(1/ε)\dot\varphi \approx \omega_c = O(1/\varepsilon)

Next, we express our solution in the following form:

r(t,φ)=R(t)+ερ(t,φ);v(t,φ)=U(t)+u(t,φ)\bm r(t,\varphi) = \bm R(t) + \varepsilon \bm \rho(t,\varphi); \quad \bm v(t,\varphi) = \bm U(t)+ \bm u(t,\varphi)

  • R\bm R and U=U⊥+U∥1b\bm U = \bm U_\perp + U_\parallel \bm 1_b describe the gyrocenter, with U=dR/dt\bm U = d\bm R/dt
  • ρ,u\bm \rho,\bm u represent the fast gyromotion:
    • They are perpendicular to B(R,t)\bm B(\bm R,t), i.e., ρ⋅1b=0\bm \rho \cdot \bm 1_b =0, u⋅1b=0\bm u \cdot \bm 1_b =0
    • They are 2π2\pi-periodic in φ\varphi for fixed tt
    • They have zero φ\varphi-average:

<ρ>=12π∮ρdφ=0;<u>=12π∮udφ=0\left<\bm \rho\right> = \frac{1}{2\pi}\oint \bm \rho d\varphi = \bm 0; \quad \left<\bm u\right> = \frac{1}{2\pi}\oint \bm u d\varphi = \bm 0

The vector basis {1b,1r,1φ}\{\bm 1_b, \bm 1_r, \bm 1_\varphi\} is defined at (R,t)(\bm R,t).

Finally, we expand:

R=R0+εR1+ε2R2+…U=U0+εU1+ε2U2+…ρ=ρ0+ερ1+ε2ρ2+…u=u0+εu1+ε2u2+…\bm R = \bm R_0 +\varepsilon\bm R_1 +\varepsilon^2\bm R_2 + \ldots \\ \bm U = \bm U_0 +\varepsilon\bm U_1 +\varepsilon^2\bm U_2 + \ldots \\ \bm \rho = \bm \rho_0 +\varepsilon\bm \rho_1 +\varepsilon^2\bm \rho_2 + \ldots \\ \bm u = \bm u_0 +\varepsilon\bm u_1 +\varepsilon^2\bm u_2 + \ldots

The gyrophase evolves as:

dφdt=1εΩ−1+Ω0+εΩ1+…\frac{d \varphi}{dt} = \frac{1}{\varepsilon}\Omega_{-1} + \Omega_0 + \varepsilon \Omega_1 + \ldots

Magnetic field gradient:

Consider B(R,t)=B01z\bm B (\bm R,t) = B_0 \bm 1_z. The gradient of B\bm B has different terms:

∇B=[∂Bx/∂x∂By/∂x∂Bz/∂x∂Bx/∂y∂By/∂y∂Bz/∂y∂Bx/∂z∂By/∂z∂Bz/∂z]\nabla \bm B = \left[ \begin{array}{cc|c} \partial B_x/\partial x & \partial B_y/\partial x & \partial B_z/\partial x \\ \partial B_x/\partial y & \partial B_y/\partial y & \partial B_z/\partial y \\ \hline \partial B_x/\partial z & \partial B_y/\partial z & \partial B_z/\partial z \\ \end{array} \right]

  • Diagonal terms are the divergence terms. Observe that always ∇⋅B=0\nabla \cdot \bm B = 0
  • The variation of BzB_z with xx and yy are the gradient terms
  • Variations of BxB_x, ByB_y with zz are the curvature terms
  • Terms ∂Bx/∂y\partial B_x/\partial y and ∂By/∂x\partial B_y/\partial x are shear terms

We write:

B(r,t)≃B(R,t)+ε(ρ⋅∇)B(R,t)+…E(r,t)≃E(R,t)+ε(ρ⋅∇)E(R,t)+…\bm B(\bm r,t)\simeq \bm B(\bm R,t) + \varepsilon(\bm \rho\cdot \nabla) \bm B (\bm R,t) + \ldots \\ \bm E(\bm r,t)\simeq \bm E(\bm R,t) + \varepsilon(\bm \rho\cdot \nabla) \bm E (\bm R,t) + \ldots

and expand the magnetic force as:

±qv×B=±qv×B(R,t)±εqv×(ρ⋅∇)B(R,t)\pm q\bm v \times \bm B = \pm q\bm v \times \bm B(\bm R,t) \pm \varepsilon q\bm v \times (\bm \rho\cdot\nabla)\bm B(\bm R,t)

  • The second term on the RHS can be thought of as a small "extra force" F\bm F

Furthermore, we require that E∥E_\parallel be small:

E(R,t)=E⊥(R,t)+εE∥(R,t)1b\bm E(\bm R,t) = \bm E_\perp(\bm R,t) + \varepsilon E_\parallel(\bm R,t)\bm 1_b

Leading order solution

To order O(1/ε)O(1/\varepsilon), the equation of motion reads:

Ω−1du0dφ∓qmu0×B(R,t)=±qmE⊥(R,t)±qmU⊥0×B(R,t)\Omega_{-1}\frac{d \bm u_0}{d\varphi} \mp \frac{q}{m}\bm u_0 \times \bm B(\bm R,t) = \pm\frac{q}{m}\bm E_\perp(\bm R,t) \pm \frac{q}{m}\bm U_{\perp 0} \times \bm B(\bm R,t)

Averaging in φ\varphi:

0=E⊥(R,t)+U⊥0×B(R,t)\bm 0 = \bm E_\perp(\bm R,t) + \bm U_{\perp 0} \times \bm B(\bm R,t)

This yields

U⊥0=E⊥(R,t)×B(R,t)B2(R,t)\bm U_{\perp 0} = \frac{\bm E_\perp(\bm R,t)\times \bm B(\bm R,t)}{B^2(\bm R,t)}

Subtracting this from the equation:

Ω−1du0dφ∓qmu0×B(R,t)=0;Ω−1dρ0dφ=u0\Omega_{-1}\frac{d \bm u_0}{d\varphi} \mp \frac{q}{m}\bm u_0 \times \bm B(\bm R,t) = \bm 0; \quad \Omega_{-1}\frac{d \bm \rho_0}{d \varphi} = \bm u_0

In the local cylindrical vector basis, integrating and imposing periodicity in φ\varphi:

u0=∓u⊥01φ;Ω−1=∓qB(R,t)m=∓ωc(R,t);ρ0=ℓ01r\bm u_0 = \mp u_{\perp 0} \bm 1_\varphi; \quad \Omega_{-1} = \mp \frac{qB(\bm R,t)}{m} = \mp\omega_c(\bm R,t); \quad \bm \rho_0 = \ell_0 \bm 1_r

Thus, we recover the solution we already had for uniform, constant fields

Guiding center drifts

To order O(1)O(1) the equation reads:

Ω−1∂u1∂φ∓qmu1×B(R,t)+dU0dt+∂u0∂t+Ω0∂u0∂φ=±qmE∥(R,t)1b±qm(ρ0⋅∇)E(R,t)±qmU1⊥×B(R,t)±qmU⊥0×(ρ0⋅∇)B(R,t)±qmu0×(ρ0⋅∇)B(R,t)\Omega_{-1}\frac{\partial \bm u_1}{\partial \varphi} \mp \frac{q}{m}\bm u_1 \times \bm B(\bm R,t) + \frac{d \bm U_0}{dt} + \frac{\partial \bm u_0}{\partial t} + \Omega_{0}\frac{\partial \bm u_0}{\partial \varphi} \\ = \pm \frac{q}{m}E_\parallel(\bm R,t)\bm 1_b \pm \frac{q}{m}(\bm\rho_0\cdot\nabla)\bm E(\bm R,t) \pm \frac{q}{m}\bm U_{1\perp} \times \bm B(\bm R,t) \\ \pm \frac{q}{m}\bm U_{\perp 0} \times (\bm\rho_0\cdot\nabla)\bm B(\bm R,t) \pm \frac{q}{m}\bm u_0 \times (\bm\rho_0\cdot\nabla)\bm B(\bm R,t)

Averaging in φ\varphi eliminates most of the terms:

dU0dt=±qmE∥(R,t)1b±qmU1⊥×B(R,t)±qm<u0×(ρ0⋅∇)B(R,t)>\frac{d \bm U_0}{dt} = \pm \frac{q}{m}E_\parallel(\bm R,t)\bm 1_b \pm \frac{q}{m}\bm U_{1\perp} \times \bm B(\bm R,t) \pm \frac{q}{m}\left<\bm u_0 \times (\bm\rho_0\cdot\nabla)\bm B(\bm R,t)\right>

±qm<u0×(ρ0⋅∇)B(R,t)>=−qmu⊥0ℓ0<1φ×(1r⋅∇)B(R,t)>=qmu⊥0ℓ012π∮dφ[1b(1r⋅∇)B(R,t)⋅1r−1r(1r⋅∇)B(R,t)⋅1b]\pm \frac{q}{m}\left<\bm u_0 \times (\bm\rho_0\cdot\nabla)\bm B(\bm R,t)\right> = - \frac{q}{m} u_{\perp 0} \ell_0\left<\bm 1_\varphi \times (\bm 1_r\cdot\nabla)\bm B(\bm R,t)\right> \\ = \frac{q}{m} u_{\perp 0} \ell_0 \frac{1}{2\pi} \oint d\varphi \left[ \bm 1_b (\bm 1_r\cdot\nabla)\bm B(\bm R,t) \cdot \bm 1_r - \bm 1_r (\bm 1_r\cdot\nabla)\bm B(\bm R,t) \cdot \bm 1_b \right]

To work out this expression use a local Cartesian vector basis {1x,1y,1z}\{\bm 1_x,\bm 1_y, \bm 1_z\}, with

1r=cos⁡φ1x+sin⁡φ1y;1b=1z\bm 1_r = \cos\varphi\bm 1_x + \sin\varphi\bm 1_y; \quad \bm 1_b = \bm 1_z

and use:

∇⋅B=∂Bx∂x+∂By∂y+∂Bz∂z=0\nabla\cdot\bm B = \frac{\partial B_x}{\partial x} + \frac{\partial B_y}{\partial y} + \frac{\partial B_z}{\partial z} = 0

After integrating, and noting that B(R,t)=Bz(R,t)1z=B(R,t)1z\bm B(\bm R,t) = B_z(\bm R,t) \bm 1_z = B(\bm R,t) \bm 1_z:

±qm<u0×(ρ0⋅∇)B(R,t)>=−qmu⊥0ℓ012∇B(R,t)=−u⊥022B(R,t)∇B(R,t)\pm \frac{q}{m}\left<\bm u_0 \times (\bm\rho_0\cdot\nabla)\bm B(\bm R,t)\right> = - \frac{q}{m} u_{\perp 0} \ell_0 \frac{1}{2} \nabla B(\bm R,t) = - \frac{u_{\perp 0}^2}{2B(\bm R,t)} \nabla B(\bm R,t)

Introducing the (unsigned) magnetic moment at order O(1)O(1):

μ0(R,t)=mu⊥022B(R,t)\mu_0(\bm R,t) = \frac{m u_{\perp 0}^ 2}{2B(\bm R,t)}

Now, we can rewrite the averaged equation as

mdU0dt=±qE∥(R,t)1b±qU1⊥×B(R,t)−μ0(R,t)∇B(R,t)m\frac{d \bm U_0}{dt} = \pm q E_\parallel(\bm R,t)\bm 1_b \pm q\bm U_{1\perp} \times \bm B(\bm R,t) -\mu_0(\bm R,t)\nabla B(\bm R,t)

The parallel projection of this equation gives U∥0U_{\parallel 0}:

mdU∥0dt=±qE∥(R,t)−μ0(R,t)(1b⋅∇)B(R,t)−m1b⋅dU⊥0dtm\frac{d U_{\parallel 0}}{dt} = \pm q E_\parallel(\bm R,t) -\mu_0(\bm R,t)(\bm 1_b\cdot\nabla) B(\bm R,t) - m\bm 1_b \cdot\frac{d \bm U_{\perp 0}}{dt}

Dropping the (R,t)(\bm R,t) dependency for the sake of notation:

mdU∥0dt=±qE∥−μ0(1b⋅∇)B−m1b⋅dU⊥0dtm\frac{d U_{\parallel 0}}{dt} = \pm q E_\parallel -\mu_0(\bm 1_b\cdot\nabla) B - m\bm 1_b \cdot\frac{d \bm U_{\perp 0}}{dt}

  • The second term in the RHS is the magnetic mirror force
  • It pushes particles in the direction of decreasing BB regardless of their sign

The perpendicular projection gives:

U1⊥=±1bωc(R,t)×[dU0dt+μ0(R,t)∇B(R,t)]\bm U_{1\perp} = \pm \frac{\bm 1_b}{\omega_c(\bm R,t)} \times \left[ \frac{d \bm U_0}{dt} +\mu_0(\bm R,t)\nabla B(\bm R,t) \right]

Dropping the (R,t)(\bm R,t) dependency for the sake of notation:

U1⊥=±1bωc×[dU0dt+μ0m∇B]\bm U_{1\perp} = \pm \frac{\bm 1_b}{\omega_c} \times \left[ \frac{d \bm U_0}{dt} +\frac{\mu_0}{m}\nabla B \right]

  • The first term is known as the inertial drift
  • The second term is the magnetic drift

Magnetic drift

±μ0mωc1b×∇B\pm \frac{\mu_0}{m\omega_c} \bm 1_b \times \nabla B

  • It displaces the gyrocenter at right angles with ∂B/∂1⊥\partial B/\partial \bm 1_\perp
  • The effect is the opposite for positive and negative particles
  • It depends on mm and on u⊥0u_{\perp 0} through μ0\mu_0

Inertial drift

We split it in two parts. The derivative of U∥01bU_{\parallel 0}\bm 1_b gives:

±1bωc×dU∥01bdt=±U∥0ωc1b×d1bdt=±U∥0ωc1b×[∂1b∂t+(U⊥0⋅∇)1b+U∥0(1b⋅∇)1b]\pm \frac{\bm 1_b}{\omega_c} \times \frac{d \bm U_{\parallel 0} \bm 1_b}{dt} = \pm \frac{U_{\parallel 0}}{\omega_c} \bm 1_b \times \frac{d \bm 1_b}{dt} \\ = \pm \frac{U_{\parallel 0}}{\omega_c} \bm 1_b \times \left[ \frac{\partial\bm 1_b}{\partial t} + (\bm U_{\perp 0} \cdot \nabla) \bm 1_b + U_{\parallel 0} (\bm 1_b \cdot \nabla) \bm 1_b \right]

Using the curvature vector κ=(1b⋅∇)1b\bm \kappa = (\bm 1_b \cdot \nabla) \bm 1_b, the last part is the curvature drift:

±U∥02ωc1b×κ\pm \frac{U_{\parallel 0}^2}{\omega_c} \bm 1_b \times \bm \kappa

  • It displaces the gyrocenter at right angles with κ\bm \kappa
  • The effect is the opposite for positive and negative particles
  • It depends on mm and on U∥02U_{\parallel 0}^2

The derivative of U⊥0\bm U_{\perp 0} gives:

±1bωc×dU⊥0dt\pm \frac{\bm 1_b}{\omega_c} \times \frac{d \bm U_{\perp 0}}{dt}

  • This is known as the polarization drift
  • The effect is the opposite for positive and negative particles

Conservation of μ0\mu_0 to O(1)O(1)

The equation of the φ\varphi-averaged kinetic energy of the particle reads:

m2ddt(u⊥02+U02)=±qU∥0E∥±qU⊥1⋅E⊥±q<(ρ0⋅∇)E⋅u0>\frac{m}{2} \frac{d}{dt} (u_{\perp 0}^2 + U_0^2) = \pm q U_{\parallel 0}E_\parallel \pm q \bm U_{\perp 1}\cdot \bm E_\perp \pm q \left<(\bm \rho_0\cdot \nabla)\bm E \cdot \bm u_0 \right>

But:

mU0⋅dU0dt=±qU∥0E∥±qU⊥0⋅(U1⊥×B)−μ0U0⋅∇Bm\bm U_0\cdot\frac{d \bm U_0}{dt} = \pm q U_{\parallel 0}E_\parallel \pm q\bm U_{\perp 0}\cdot (\bm U_{1\perp} \times \bm B) -\mu_0\bm U_0\cdot\nabla B

And:

U⊥0=E⊥×BB2⇒±qU⊥0⋅(U1⊥×B)=±qU⊥1⋅E⊥U_{\perp 0} = \frac{\bm E_\perp\times\bm B}{B^2} \quad\Rightarrow\quad \pm q\bm U_{\perp 0}\cdot (\bm U_{1\perp} \times \bm B)= \pm q \bm U_{\perp 1}\cdot \bm E_\perp

Additionally:

±q<(ρ0⋅∇)E⋅u0>=mu⊥022πB∮dφ(1r⋅∇)E⋅1φ=−μ01b⋅∇×E=μ0∂B∂t\pm q \left<(\bm \rho_0\cdot \nabla)\bm E \cdot \bm u_0 \right> = \frac{m u_{\perp 0}^2}{2\pi B}\oint d\varphi (\bm 1_r \cdot \nabla)\bm E \cdot \bm 1_\varphi = - \mu_0 \bm 1_b \cdot \nabla \times \bm E = \mu_0 \frac{\partial B}{\partial t}

Therefore:

m2ddtu⊥02−μ0∂B∂t−μ0U0⋅∇B=0\frac{m}{2} \frac{d}{dt} u_{\perp 0}^2 - \mu_0 \frac{\partial B}{\partial t} -\mu_0 \bm U_0\cdot \nabla B = 0

Since μ0=u⊥02/(2B)\mu_0 = u_{\perp 0}^2/(2B) and dB/dt=∂B/∂t+U0⋅∇BdB/dt = \partial B/\partial t + \bm U_0 \cdot \nabla B in the averaged equation:

ddtμ0=0\frac{d}{dt} \mu_0 = 0

  • μ0\mu_0 is conserved to lowest order and is therefore an adiabatic invariant

See HAZE18 p 29 to learn more about Poincaré invariants and Adiabatic invariants

Magnetic mirror

In stationary fields,

H=12mv∥2+μB=constH = \frac{1}{2} mv_\parallel^2 + \mu B = \text{const}

μ=mv⊥22B=const\mu = \frac{mv_\perp^2}{2B} = \text{const}

  • Particles trade parallel to perpendicular kinetic energy as they move toward a region of maximum BB
  • If BmaxB_{max} at the maximum, only particles with v⊥2/v2<B/Bmaxv_\perp^2/v^2 < B/B_{max} can escape. This defines a loss cone in phase space

Further reading

  • I have followed mainly HAZE18 (with adaptations)
  • BITT04 gives an ad-hoc derivation for each drift, but easier to understand
  • litt83 is a beautiful derivation based on Hamiltonian dynamics—maybe a topic for another time

I follow mainly HAZE17, with adaptations