Categories: Electromagnetism, Physics, Plasma physics.

Magnetohydrodynamics

Magnetohydrodynamics (MHD) describes the dynamics of fluids that are electrically conductive. Notably, it is often suitable to describe plasmas, and can be regarded as a special case of the two-fluid model; we will derive it as such, but the results are not specific to plasmas.

In the two-fluid model, we described the plasma as two separate fluids, but in MHD we treat it as a single conductive fluid. The macroscopic pressure pp and electric current density J\vb{J} are:

p=pi+peJ=qiniui+qeneue\begin{aligned} p = p_i + p_e \qquad \quad \vb{J} = q_i n_i \vb{u}_i + q_e n_e \vb{u}_e \end{aligned}

Meanwhile, the macroscopic mass density ρ\rho and center-of-mass flow velocity u\vb{u} are as follows, although the ions dominate due to their large mass:

ρ=mini+meneminiu=1ρ(miniui+meneue)ui\begin{aligned} \rho = m_i n_i + m_e n_e \approx m_i n_i \qquad \quad \vb{u} = \frac{1}{\rho} \Big( m_i n_i \vb{u}_i + m_e n_e \vb{u}_e \Big) \approx \vb{u}_i \end{aligned}

With these quantities in mind, we add up the two-fluid continuity equations, multiplied by their respective particles’ masses:

0=minit+menet+mi(niui)+me(neue)\begin{aligned} 0 &= m_i \pdv{n_i}{t} + m_e \pdv{n_e}{t} + m_i \nabla \cdot (n_i \vb{u}_i) + m_e \nabla \cdot (n_e \vb{u}_e) \end{aligned}

After some straightforward rearranging, we arrive at the single-fluid continuity relation:

ρt+(ρu)=0\begin{aligned} \boxed{ \pdv{\rho}{t} + \nabla \cdot (\rho \vb{u}) = 0 } \end{aligned}

Next, consider the two-fluid momentum equations for the ions and electrons, respectively:

miniDuiDt=qini(E+ui×B)pifiemini(uiue)meneDueDt=qene(E+ue×B)pefeimene(ueui)\begin{aligned} m_i n_i \frac{\mathrm{D} \vb{u}_i}{\mathrm{D} t} &= q_i n_i (\vb{E} + \vb{u}_i \cross \vb{B}) - \nabla p_i - f_{ie} m_i n_i (\vb{u}_i - \vb{u}_e) \\ m_e n_e \frac{\mathrm{D} \vb{u}_e}{\mathrm{D} t} &= q_e n_e (\vb{E} + \vb{u}_e \cross \vb{B}) - \nabla p_e - f_{ei} m_e n_e (\vb{u}_e - \vb{u}_i) \end{aligned}

We will assume that electrons’ inertia is negligible compared to the Lorentz force. Let τchar\tau_\mathrm{char} be the characteristic timescale of the plasma’s dynamics, i.e. nothing noticable happens in times shorter than τchar\tau_\mathrm{char}, then this assumption can be written as:

1meneDue/Dtqeneue×Bmeneue/τcharqeneueB=meqeBτchar=1ωceτchar1\begin{aligned} 1 \gg \frac{\big| m_e n_e \mathrm{D} \vb{u}_e / \mathrm{D} t \big|}{\big| q_e n_e \vb{u}_e \cross \vb{B} \big|} \sim \frac{m_e n_e |\vb{u}_e| / \tau_\mathrm{char}}{q_e n_e |\vb{u}_e| |\vb{B}|} = \frac{m_e}{q_e |\vb{B}| \tau_\mathrm{char}} = \frac{1}{\omega_{ce} \tau_\mathrm{char}} \ll 1 \end{aligned}

Where we have recognized the cyclotron frequency ωc\omega_c (see Lorentz force article). In other words, our assumption is equivalent to the electron gyration period 2π/ωce2 \pi / \omega_{ce} being small compared to the macroscopic dynamics’ timescale τchar\tau_\mathrm{char}. By construction, we can thus ignore the left-hand side of the electron momentum equation, leaving:

miniDuiDt=qini(E+ui×B)pifiemini(uiue)0=qene(E+ue×B)pefeimene(ueui)\begin{aligned} m_i n_i \frac{\mathrm{D} \vb{u}_i}{\mathrm{D} t} &= q_i n_i (\vb{E} + \vb{u}_i \cross \vb{B}) - \nabla p_i - f_{ie} m_i n_i (\vb{u}_i - \vb{u}_e) \\ 0 &= q_e n_e (\vb{E} + \vb{u}_e \cross \vb{B}) - \nabla p_e - f_{ei} m_e n_e (\vb{u}_e - \vb{u}_i) \end{aligned}

We add up these momentum equations, recognizing the pressure pp and current J\vb{J}:

miniDuiDt=(qini+qene)E+J×Bpfiemini(ui ⁣ ⁣ue)feimene(ue ⁣ ⁣ui)=(qini+qene)E+J×Bp\begin{aligned} m_i n_i \frac{\mathrm{D} \vb{u}_i}{\mathrm{D} t} &= (q_i n_i + q_e n_e) \vb{E} + \vb{J} \cross \vb{B} - \nabla p - f_{ie} m_i n_i (\vb{u}_i \!-\! \vb{u}_e) - f_{ei} m_e n_e (\vb{u}_e \!-\! \vb{u}_i) \\ &= (q_i n_i + q_e n_e) \vb{E} + \vb{J} \cross \vb{B} - \nabla p \end{aligned}

Where we have used fiemini=feimenef_{ie} m_i n_i = f_{ei} m_e n_e because momentum is conserved by the underlying Rutherford scattering process, which is elastic. In other words, the momentum given by ions to electrons is equal to the momentum received by electrons from ions.

Since the two-fluid model assumes that the Debye length λD\lambda_D is small compared to a “blob” dV\dd{V} of the fluid, we can invoke quasi-neutrality qini+qene=0q_i n_i + q_e n_e = 0. Using that ρmini\rho \approx m_i n_i and uui\vb{u} \approx \vb{u}_i, we thus arrive at the momentum equation:

ρDuDt=J×Bp\begin{aligned} \boxed{ \rho \frac{\mathrm{D} \vb{u}}{\mathrm{D} t} = \vb{J} \cross \vb{B} - \nabla p } \end{aligned}

However, we found this by combining two equations into one, so some information was implicitly lost; we need a second momentum equation. Therefore, we return to the electrons’ momentum equation, after a bit of rearranging:

E+ue×Bpeqene=feimeqe(ueui)\begin{aligned} \vb{E} + \vb{u}_e \cross \vb{B} - \frac{\nabla p_e}{q_e n_e} = \frac{f_{ei} m_e}{q_e} (\vb{u}_e - \vb{u}_i) \end{aligned}

Again using quasi-neutrality qini=qeneq_i n_i = - q_e n_e, the current density J=qene(ue ⁣ ⁣ui)\vb{J} = q_e n_e (\vb{u}_e \!-\! \vb{u}_i), so:

E+ue×Bpeqene=ηJηfeimeneqe2\begin{aligned} \vb{E} + \vb{u}_e \cross \vb{B} - \frac{\nabla p_e}{q_e n_e} = \eta \vb{J} \qquad \quad \eta \equiv \frac{f_{ei} m_e}{n_e q_e^2} \end{aligned}

Where η\eta is the electrical resistivity of the plasma, see Spitzer resistivity for more information, and a rough estimate of this quantity for a plasma.

Now, using that uui\vb{u} \approx \vb{u}_i, we add (u ⁣ ⁣ui)×B0(\vb{u} \!-\! \vb{u}_i) \cross \vb{B} \approx 0 to the equation, and insert J\vb{J} again:

ηJ=E+u×B+(ueui)×Bpeqene=E+u×B+J×Bqenepeqene\begin{aligned} \eta \vb{J} &= \vb{E} + \vb{u} \cross \vb{B} + (\vb{u}_e - \vb{u}_i) \cross \vb{B} - \frac{\nabla p_e}{q_e n_e} \\ &= \vb{E} + \vb{u} \cross \vb{B} + \frac{\vb{J} \cross \vb{B}}{q_e n_e} - \frac{\nabla p_e}{q_e n_e} \end{aligned}

Next, we want to get rid of the pressure term. To do so, we take the curl of the equation:

×(ηJ)=Bt+×(u×B)+×J×Bqene×peqene\begin{aligned} \nabla \cross (\eta \vb{J}) = - \pdv{\vb{B}}{t} + \nabla \cross (\vb{u} \cross \vb{B}) + \nabla \cross \frac{\vb{J} \cross \vb{B}}{q_e n_e} - \nabla \cross \frac{\nabla p_e}{q_e n_e} \end{aligned}

Where we have used Faraday’s law. This is the induction equation, and is used to compute B\vb{B}. The pressure term can be rewritten using the ideal gas law pe=kBTenep_e = k_B T_e n_e:

×peqene=kBqe×(neTe)ne=kBqe×(Te+Tenene)\begin{aligned} \nabla \cross \frac{\nabla p_e}{q_e n_e} = \frac{k_B}{q_e} \nabla \cross \frac{\nabla (n_e T_e)}{n_e} = \frac{k_B}{q_e} \nabla \cross \Big( \nabla T_e + T_e \frac{\nabla n_e}{n_e} \Big) \end{aligned}

The curl of a gradient is always zero, and we notice that ne/ne= ⁣ln(ne)\nabla n_e / n_e = \nabla\! \ln(n_e). Then we use the vector identity ×(fg)=f×g\nabla \cross (f \nabla g) = \nabla f \cross \nabla g, leading to:

×peqene=kBqe×(Te ⁣ln(ne))=kBqe(Te× ⁣ln(ne))=kBqene(Te×ne)\begin{aligned} \nabla \cross \frac{\nabla p_e}{q_e n_e} = \frac{k_B}{q_e} \nabla \cross \big( T_e \: \nabla\! \ln(n_e) \big) = \frac{k_B}{q_e} \big( \nabla T_e \cross \nabla\! \ln(n_e) \big) = \frac{k_B}{q_e n_e} \big( \nabla T_e \cross \nabla n_e \big) \end{aligned}

It is reasonable to assume that Te\nabla T_e and ne\nabla n_e point in roughly the same direction, in which case the pressure term can be neglected. Consequently, pep_e has no effect on the dynamics of B\vb{B}, so we argue that it can be dropped from the original (non-curled) equation too, leaving:

E+u×B+J×Bqene=ηJ\begin{aligned} \boxed{ \vb{E} + \vb{u} \cross \vb{B} + \frac{\vb{J} \cross \vb{B}}{q_e n_e} = \eta \vb{J} } \end{aligned}

This is known as the generalized Ohm’s law, since it contains the relation E=ηJ\vb{E} = \eta \vb{J}.

Next, consider Ampère’s law, where we would like to neglect the last term:

×B=μ0J+1c2Et\begin{aligned} \nabla \cross \vb{B} = \mu_0 \vb{J} + \frac{1}{c^2} \pdv{\vb{E}}{t} \end{aligned}

From Faraday’s law, we can obtain a scale estimate for E\vb{E}. Recall that τchar\tau_\mathrm{char} is the characteristic timescale of the plasma, and let λcharλD\lambda_\mathrm{char} \gg \lambda_D be its characteristic lengthscale:

×E=Bt    EλcharτcharB\begin{aligned} \nabla \cross \vb{E} = - \pdv{\vb{B}}{t} \quad \implies \quad |\vb{E}| \sim \frac{\lambda_\mathrm{char}}{\tau_\mathrm{char}} |\vb{B}| \end{aligned}

From this, we find when we can neglect the last term in Ampère’s law: the characteristic velocity vcharv_\mathrm{char} must be tiny compared to cc, i.e. the plasma must be non-relativistic:

1(E/t)/c2×BE/τcharBc2/λcharBλchar2/τchar2Bc2=vchar2c21\begin{aligned} 1 \gg \frac{\big| (\ipdv{\vb{E}}{t}) / c^2 \big|}{\big| \nabla \cross \vb{B} \big|} \sim \frac{|\vb{E}| / \tau_\mathrm{char}}{|\vb{B}| c^2 / \lambda_\mathrm{char}} \sim \frac{|\vb{B}| \lambda_\mathrm{char}^2 / \tau_\mathrm{char}^2}{|\vb{B}| c^2} = \frac{v_\mathrm{char}^2}{c^2} \ll 1 \end{aligned}

We thus have the following reduced form of Ampère’s law, in addition to Faraday’s law:

×B=μ0J×E=Bt\begin{aligned} \boxed{ \nabla \cross \vb{B} = \mu_0 \vb{J} } \qquad \quad \boxed{ \nabla \cross \vb{E} = - \pdv{\vb{B}}{t} } \end{aligned}

Finally, we revisit the thermodynamic equation of state, for a single fluid this time. Using the product rule of differentiation yields:

0=DDt(pργ)=DpDtργpγργ1DρDt\begin{aligned} 0 &= \frac{\mathrm{D}}{\mathrm{D} t} \Big( \frac{p}{\rho^\gamma} \Big) = \frac{\mathrm{D} p}{\mathrm{D} t} \rho^{-\gamma} - p \gamma \rho^{-\gamma - 1} \frac{\mathrm{D} \rho}{\mathrm{D} t} \end{aligned}

The continuity equation allows us to rewrite the material derivative Dρ/Dt\mathrm{D} \rho / \mathrm{D} t as follows:

ρt+(ρu)=ρt+ρu+uρ=ρu+DρDt=0\begin{aligned} \pdv{\rho}{t} + \nabla \cdot (\rho \vb{u}) = \pdv{\rho}{t} + \rho \nabla \cdot \vb{u} + \vb{u} \cdot \nabla \rho = \rho \nabla \cdot \vb{u} + \frac{\mathrm{D} \rho}{\mathrm{D} t} = 0 \end{aligned}

Inserting this into the equation of state leads us to a differential equation for pp:

0=DpDt+pγ1ρρu    DpDt=pγu\begin{aligned} 0 = \frac{\mathrm{D} p}{\mathrm{D} t} + p \gamma \frac{1}{\rho} \rho \nabla \cdot \vb{u} \quad \implies \quad \boxed{ \frac{\mathrm{D} p}{\mathrm{D} t} = - p \gamma \nabla \cdot \vb{u} } \end{aligned}

This closes the set of 14 MHD equations for 14 unknowns. Originally, the two-fluid model had 16 of each, but we have merged nin_i and nen_e into ρ\rho, and pip_i and pip_i into pp.

Ohm’s law variants

It is worth discussing the generalized Ohm’s law in more detail. Its full form was:

E+u×B+J×Bqene=ηJ\begin{aligned} \vb{E} + \vb{u} \cross \vb{B} + \frac{\vb{J} \cross \vb{B}}{q_e n_e} = \eta \vb{J} \end{aligned}

However, most authors neglect some of its terms: this form is used for Hall MHD, where J×B\vb{J} \cross \vb{B} is called the Hall term. This term can be dropped in any of the following cases:

1J×B/qeneu×Bρvchar/τcharvcharBqiniminiBqiniτchar=1ωciτchar11J×B/qeneηJJBqe2nefeimeJqene=Bqefeime=ωcefei1\begin{gathered} 1 \gg \frac{\big| \vb{J} \cross \vb{B} / q_e n_e \big|}{\big| \vb{u} \cross \vb{B} \big|} \sim \frac{\rho v_\mathrm{char} / \tau_\mathrm{char}}{v_\mathrm{char} |\vb{B}| q_i n_i} \approx \frac{m_i n_i}{|\vb{B}| q_i n_i \tau_\mathrm{char}} = \frac{1}{\omega_{ci} \tau_\mathrm{char}} \ll 1 \\ 1 \gg \frac{\big| \vb{J} \cross \vb{B} / q_e n_e \big|}{\big| \eta \vb{J} \big|} \sim \frac{|\vb{J}| |\vb{B}| q_e^2 n_e}{f_{ei} m_e |\vb{J}| q_e n_e} = \frac{|\vb{B}| q_e}{f_{ei} m_e} = \frac{\omega_{ce}}{f_{ei}} \ll 1 \end{gathered}

Where we have used the MHD momentum equation with p0\nabla p \approx 0 to obtain the scale estimate J×Bρvchar/τchar\vb{J} \cross \vb{B} \sim \rho v_\mathrm{char} / \tau_\mathrm{char}. In other words, if the ion gyration period is short τcharωci\tau_\mathrm{char} \gg \omega_{ci}, and/or if the electron gyration period is long compared to the electron-ion collision period ωcefei\omega_{ce} \ll f_{ei}, then we are left with this form of Ohm’s law, used in resistive MHD:

E+u×B=ηJ\begin{aligned} \vb{E} + \vb{u} \cross \vb{B} = \eta \vb{J} \end{aligned}

Finally, we can neglect the resisitive term ηJ\eta \vb{J} if the Lorentz force is much larger. We formalize this condition as follows, where we have used Ampère’s law to find JB/μ0λchar\vb{J} \sim \vb{B} / \mu_0 \lambda_\mathrm{char}:

1u×BηJvcharBηJvcharBηB/μ0λchar=Rm1\begin{aligned} 1 \ll \frac{\big| \vb{u} \cross \vb{B} \big|}{\big| \eta \vb{J} \big|} \sim \frac{v_\mathrm{char} |\vb{B}|}{\eta \vb{J}} \sim \frac{v_\mathrm{char} |\vb{B}|}{\eta |\vb{B}| / \mu_0 \lambda_\mathrm{char}} = \mathrm{R_m} \gg 1 \end{aligned}

Where we have defined the magnetic Reynolds number Rm\mathrm{R_m} as follows, which is analogous to the fluid Reynolds number Re\mathrm{Re}:

Rmvcharλcharη/μ0\begin{aligned} \boxed{ \mathrm{R_m} \equiv \frac{v_\mathrm{char} \lambda_\mathrm{char}}{\eta / \mu_0} } \end{aligned}

If Rm1\mathrm{R_m} \ll 1, the plasma is “electrically viscous”, such that resistivity needs to be accounted for, whereas if Rm1\mathrm{R_m} \gg 1, the resistivity is negligible, in which case we have ideal MHD:

E+u×B=0\begin{aligned} \vb{E} + \vb{u} \cross \vb{B} = 0 \end{aligned}

References

  1. P.M. Bellan, Fundamentals of plasma physics, 1st edition, Cambridge.
  2. M. Salewski, A.H. Nielsen, Plasma physics: lecture notes, 2021, unpublished.