## Magnetohydrodynamics

The dynamics of conducting fluids in an electromagnetic field

### Introduction

Magnetohydrodynamics, or MHD, is a branch of the science of the dynamics of matter moving in an electromagnetic field, especially where currents established in the matter by induction modify the field, so that the field and dynamics equations are coupled. It treats, in particular, conducting fluids, whether liquid or gaseous, in which certain simplifying postulates are accepted. These are, generally, that the Maxwell displacement current is neglected, and the fluid may be treated as a continuum, without mean-free-path effects. It is distinguished from the closely related plasma dynamics in which these postulates are relaxed, but there is still a large intermediate area in which similar treatment is possible.

Solid matter is generally excluded from MHD, but it should be realized that the same principles apply. Electrical conduction in metals, and the Hall Effect, are two examples. In an electric motor, the magnetic field produced by the armature current affects the operation of the motor in an important way, so that the mechanical and electrical analyses are coupled, just as in MHD. Electromagnetic forces are an essential part of motors and generators, though they generally do not produce significant elastic deformations, and the motions occur with the help of rotating and sliding contacts. Homopolar generators (ones that produce DC currents) are, indeed, closely related to MHD analogues.

MHD was originally applied to astrophysical and geophysical problems, where it is still very important, but more recently to the problem of fusion power, where the application is the creation and containment of hot plasmas by electromagnetic forces, since material walls would be destroyed. Astrophysical problems include solar structure, especially in the outer layers, the solar wind bathing the earth and other planets, and interstellar magnetic fields. The primary geophysical problem is planetary magnetism, produced by currents deep in the planet, a problem that has not been solved to any degree of satisfaction.

The purpose of this article is to explain the fundamental principles of MHD as best I understand them, and not to treat any of the applications in detail. This requires an acquaintance with several areas of theoretical physics, including electromagnetism and fluid mechanics, which are discussed in other areas of this website. I will use the index notation freely, in addition to the usual vector notation. Index notation is explained in Euclidean Tensors. The electromagnetic units used will be Gaussian. Much MHD literature is in emu (which is like Gaussian without the c's), while modern references will probably sink to MKSA. I shall relate all quantities to practical units, anyway.

### Conduction

Electrons have a charge of -e = -4.803 x 10-10 esu and a mass of 9.11 x 10-28 g. Because they are light, they are very mobile and are responsible for most of the electric current in all the materials we shall study. They also have an angular momentum, the spin, of (1/2)(h/2π) g-cm2/s and an associated magnetic moment. However, unless a macroscopic number of the magnetic moments are aligned, this magnetic moment, and the angular momentum, can be neglected here. Also because they are light, they take up quite a bit of room, meaning that the density of states they can occupy is not large. Furthermore, electrons, as spin-1/2 particles, obey the Fermi exclusion principle that a given spatial state can be occupied by at most two electrons with opposite spins. Electrons begin filling the states from the lowest ones up to an energy called the Fermi level when not thermally agitated. If this energy level is much higher than the average thermal energy, then thermal agitation modifies the situation very little, and the electrons are described as degenerate. In momentum space, the electrons occupy a sphere called the Fermi sphere of radius pF = √2mEF, while in physical space they are uniformly distributed with density n.

The electrons in metals are very highly degenerate, since the electron density is high. In copper, there are n = 8.5 x 1022 electrons per cm3, so EF = 7.0 eV, equivalent to a temperature of 82,000K. The electrons at the top of the Fermi sphere are, therefore, darting about at 1.56 x 108 cm/s, so any external fields make only a small impression. In many materials, semiconductors as well as metals, the allowable energies are restricted to bands. A full band results in no conduction, since none of the electrons can change its state. In an almost empty band, the electrons will act almost as if they were free, except that their effective mass may be different. In an almost full band, the "holes" will act like positive charges, again with a modified mass. In gases and plasmas, the electron density may be so low that many states are available and the electrons act as free particles.

The quantities of positive and negative charge in any volume are very close to equal, since even a small unbalance results in very strong forces that act to correct the situation. Conduction is the motion of charge carriers caused by an electric field. These charge carriers can be light, negatively charged electrons, or heavy, positively charged ions when they are considered as free and independent. The electronic and ionic components of the current are both in the direction of the applied electric field, although the velocities of the charge carriers are oppositely directed.

If F = qE is the force on an electron of charge q = -e, the electron receives an increase in momentum of dp = (h/2π)dk = qEτ, where τ is the average length of time between collisions of the rapidly-moving electrons with the lattice ions, in which it transfers all of the momentum gained to the lattice. For copper, this is τ = 2 x 10-14 s. This steady momentum increase represents a bodily displacement of the Fermi sphere in momentum space by dk = mdv, where dv is the average velocity of the electrons in the direction of the field. The current density J = nqdv = (nq2τ/m)E, so the conductivity is, by definition, σ = nq2τ/m in (statohm-cm)-1. Note that the square of the charge enters, so the sign of the charge is not significant. Charges of different sign, however, produce currents in opposite directions for the same electric field. We see a net current because the electrons move, while the ions do not. Conductivities in Gaussian units are large numbers. 1 (ohm-cm)-1, a familiar practical unit, is 9 x 1011 (statohm-cm)-1. The conductivity of copper is 5.98 x 105 (ohm-cm)-1, or 5.38 x 1017 esu. Fom this simple theory we can estimate the mean free path of an electron in copper as 30 nm. Direct measurement gives 45 nm, so the theory is pretty good. The resistivity is the reciprocal of the conductivity. For copper, it is 1.673 μΩ-cm or 1.86 x 10-18 statohm-cm. The conductivity of liquid mercury is 9.15 x 1015 esu at normal temperatures, and its density is 13.56 g/cm3.

Conduction shows us a frictional force mv/τ acting on the motion of the electrons under the electric field, so we have Ohm's Law V = IR and all its consequences, including power dissipation at a rate P = I2R, the energy appearing as heat. The relaxation time τ gives an upper bound on frequencies for which the continuum postulate of MHD is valid. The transition to plasma dynamics comes at frequencies ≈ 1/τ. Also, we may find charge separation effects, in which the electrons and the ions move with respect to one another, phenomena called plasma oscillations, in which a "two-fluid" model is required.

The Hall effect is another aspect of charge-carrier motion. Consider a current density Jx moving in a conductor bathed in a magnetic field Bz, as shown in the diagram. A charge +q experiences a force Fy = qvxBz/c, where v is the average velocity v = Jx/nq, in which n is the number of charge carriers per unit volume. This force is the familar Lorentz force on a charged particle. In equilibrium, this force is balanced by the force exerted by an electric field Ey, so Ey = (1/nqc)BzJx. The coefficient RH = 1/nqc is called the Hall Coefficient, and depends on the sign of the charge. For electrons, it is negative. It tells us directly what the net carrier charge density is in a material. For copper, it is indeed negative, but for cadmium it is positive, showing that the charge carriers are holes in the conduction band in that metal.

### Diffusion

We shall find that under some circumstances the magnetic field in a conducting fluid is described by the diffusion equation, ∂f/∂t = D ∂iif, where f is the diffusing quantity and D is the diffusivity, in cm2/s. This equation also applies to heat, to vorticity in a viscous fluid, to neutrons in a reactor, to molecular mixing of gases, and to other important physical phenomena. Therefore, it is worth considering here, where the example of thermal diffusion is used.

The thermal conductivity k is defined by F = -k (T2 - T1)/L, where F is the heat flux in erg/cm2-s, from temperature T1 to temperature T2 (in °C or K) across a uniform slab of thickness L. We use T to avoid confusion with the time t. The conductivity k has dimensions erg/cm-s-K. To convert to other units, use 1 J = 107 erg, 1 cal = 4.184 x 107 erg, 1 Btu = 1.054 x 1010 erg. The conductivity of copper is 4.01 J/s-cm-K, its density is 8.89 g/cm3, and its specific heat is 0.3846 J/g-K. For granite, these figures are 6.5 x 10-6 cal/cm-s-K, 2.65 g/cm3 and 2.28 x 10-4 cal/g-K. For mercury, they are 0.0830 J/s-cm-K, 13.56 g/cm3 and 0.1389 J/g-K. These figures can be used for calculations to get a feel for the magnitudes involved.

In tensor notation, the simple definition can be written Fi = -k∂iT. The condition that the heat content of a volume V must decrease at the rate that heat is conducted outwards away from it through its surface S is expressed differentially by ∂(c'T)/∂t = k∂iiT, where we have used the divergence theorem to change the surface integral to a volume integral. Here, c' is the heat capacity per unit volume, which is equal to cd, where c is the usual heat capacity per unit mass, and d is the density. Using c, we have ∂T/∂t = (k/cd)∂iiT = D ∂iiT, where D is the diffusivity, a kinematic coefficient with dimensions cm2/s. For some of the materials mentioned in the last paragraph, thermal diffusivities are 1.160 for copper, 0.01076 for granite, and 0.04407 for mercury, all in cm2/s. Diffusivities are directly comparable between substances and whatever is diffusing. We are probably directly familiar with thermal diffusion, so it makes a good comparison.

It will be good for us to have the solution of a popular problem in heat conduction, the cooling of a uniform sphere, as given by Churchill (see References). Let us consider a sphere of radius a, with a given radially symmetric temperature distribution T = f(r) at t = 0, and the surface of the sphere held at T = 0. It is easy to solve the diffusion equation in spherical coordinates, and to use the method of Fourier series to construct a solution satisfying the boundary conditions. We find a series of solutions depending on time through factors of exp(-n2π2Dt/a2), corresponding to radial temperature distributions proportional to (a/πr)sin(nπr/a). The n=1 term is T(r,t) = T0(a/πr)sin(πr/a)e2Dt/a2. The temperature will fall by a factor of 1/e in a time t' = a22D, which is the main result we desire.

If we apply this to a copper sphere of radius 1 cm, we find that t' = 0.0873 s. Even if this seems a bit rapid, there is nothing wrong with our analysis, simply that the condition of T = 0 on the boundary is impossible to arrange in practice. The actual time of cooling is dependent on the surface rate of heat transfer, which will govern in the practical case. If we consider instead a granite sphere the size of the earth, we find t' = 1.213 x 1011 years, about 1/28th of the time the earth has existed. Even a copper sphere the size of the earth would take 1.125 x 109 years to cool by a factor 1/e. These results will make results for magnetic fields a little less unbelievable.

If N particles are initially at the origin, and begin diffusing with diffusivity D at t = 0, the subsequent density of particles n(r,t) = N(4πDt)3/2exp(-r2/4Dt). The shape of this function is shown at the right. The average r at time t is ∫(0,∞)n(r,t)r(4πr2)dr. This integral can be found in tables, and the result is rav = (4/√π)(DT)1/2. The mean square r is ∫(0,∞)n(r,t)r2(4πr2)dr, with the result (r2)av = 6Dt. Either of these parameters can be used to describe the progress of the diffusion. The average r is 2.26(Dt)1/2, while the rms r is 2.45(Dt)1/2. As shown in the diagram, the density is down to 0.223 of the peak value at the origin at the rms r. Diffusion is a process of evening-out and spreading, and progresses proportionally to the square root of the time. Any quantity that satisfies the diffusion equation behaves in this manner, whether heat, electrons or magnetic field.

The average distance a particle will diffuse in its lifetime τ is the diffusion length L = √(Dτ), familiar from semiconductor electronics (where an electron may recombine) and reactor physics (where a neutron may be absorbed). This is another way to visualize the diffusivity.

### Electrodynamic Equations

In this section, we shall presume that the fluid velocity v is given, and study the consequent behavior of the electromagnetic fields.

The electromagnetic field is described by Maxwell's equations. We will postulate that the permeability μ = 1, so that B and H are the same thing. We will arbitrarly choose to denote the field by B, and measure it in gauss. It satisfies div B = 0, and curl B = (4π/c)J, where J is the current density in esu/s-cm2 and c = 3 x 1010 cm/s, approximately. The displacement current (1/c)∂D/∂t is neglected here, so that conduction current is the only source of B. The displacement current can be considered when necessary, but this is very seldom required, and it is usually quite negligible. The curl equation gives div J = 0, which is consistent with no free electric charge ρ.

The electric field E, measured in statvolt/cm, will be assumed to satisfy div E = 0, so there will be no free charge density, consistently with the preceding paragraph. However, curl E = -(1/c)∂B/∂t, so the changing magnetic field will be the source of the electrical field. The conduction current density J = σ[E + (1/c)v x B]. The second term is the apparent electric field E' in a coordinate system moving with the material at velocity v cm/s. This is a first-order relativistic effect, and velocities are low enough that it is a good approximation. We will seldom need the displacement D = εE, where ε is the dielectric constant. If we do need the displacement, we must be aware that ε may be significantly greater than 1.

Now the magnetic curl equation becomes curl B = (4π/c)E + (4π/c2)v x B, or E = (c/4πσ)[curl B - (1/c)v x B]. The electric curl equation then gives us the time derivative of B: ∂B/∂t = (c2/4πσ)div grad B + curl (v x B), where we have expanded curl curl B and used div B = 0. The first term gives us the diffusion equation for B, with diffusivity D = c2/4πσ, while the second term describes the dragging of the field lines by the fluid, as we shall see.

If σ is much greater than ω, a frequency characteristic of the time variation of the motions (conductivity in esu has dimensions 1/s) then the required J can be obtained with only a very small E'. The first term in the equation we just derived will be much smaller than the second, which will take over the duty of supplying the time rate of change of B. In this case, we have approximately that E = -(1/c)v x B, so the electric field is determined jointly by the fluid velocity and the magnetic field. At the opposite limit, of low conductivity, J ≈ 0 and the fields are almost unaffected by the fluid motion. The dimensionless magnetic Reynolds number, RM = VL/D = VL√(4πρσ) is a criterion for the relative importance of magnetic and inertial effects. V is a velocity, and L a length, characteristic of the situation as with the usual Reynolds number. If it is large, magnetic effects may be expected to be prominent. We may say that the conductivity is "infinite" when RM is large.

This result has a very important consequence. Consider the time rate of change of the magnetic flux through a surface S with boundary curve C that moves with the fluid. This rate, (d/dt)∫B·ndS, has two components. One is the change due to the change of B through S in dt, considering C fixed, and the other is the net flux added when each point of C moves a distance vdt, divided by dt. The first is ∫(∂B/∂t)·ndS, while the second is ∫(v x dl·B. Interchanging dot and cross, it becomes ∫(B x v)·dl. Now we can use Stokes's Theorem, and combine the two contributions as one volume integral, which is ∫[∂B/∂t - curl (v x B)]·ndS = 0! The flux passing through a curve C moving with the fluid is constant, or the curve never "cuts" magnetic flux: the magnetic field is rigidly attached to the fluid, "frozen in," and moves with it.

In this case, we can use Faraday's Law to obtain an expression for the time rate of change of the magnetic field. From the expression for E we obtained above, we find ∂B/∂t = curl v x B. When we use this expression, we imply that we are dealing with the case of high magnetic Reynolds numbers, or "infinite" conductivity, when the magnetic field is frozen in. The fluid can always move freely along the magnetic lines of force. These are not real lines, of course, but are merely a vivid representation of the direction and intensity of the magnetic field. When the magnetic Reynolds number is smaller, the lines of force only act as a drag of greater or less strength on the fluid that moves perpendicularly to them. Also, note that this always depends on the time scale. The magnetic field may be frozen in for milliseconds, but not for hours.

The only liquid metal that can be used easily in laboratory experiments is mercury. From its conductivity (given above), the diffusivity of the magnetic field is 7827 cm2/s, a rather high value. The field will diffuse out of a mercury sphere 1 cm in radius to 1/e in a time of only 13 μs, which makes the effect rather difficult to observe. Liquid sodium is better, with higher conductivity and smaller density, but is very nasty to work with. Measurements have been made on both these conducting liquids.

The ponderomotive force exerted per unit volume on a current density J esu/cm2 in a magnetic field B gauss is f = (1/c)J x B dyne/cm3. The force is perpendicular to both the current density and the magnetic field. A wire carrying a current of 1 A is acted upon by a force of only 0.1 dyne/cm in a field of 1 gauss. To get an idea of this magnitude, a centimeter of #22 AWG copper wire weighs 28.4 dyne, so it would require 284 gauss to support the wire's own weight. Nevertheless, these forces can be increased to practically useful values by increasing the magnetic field and the number of turns of wire. Magnetic fields of up to about 20,000 gauss are available with ferromagnetic materials to help. The direction of the force can be visualized most easily by considering on which side of the wire the sum of the applied fields and the fields caused by the current is larger, as shown in the diagram at the left. This increased flux pushes the wire in the direction of the weaker field.

If we move a wire in a magnetic field, and the motional emf causes a current to flow, then the force on this current opposes our moving of the wire. That is, we do mechanical work Fv per unit time, and this work appears as electrical work EI per unit time, conserving energy. One never occurs without the other, and the directions of the emf and the mechanical force are such as to support this, as shown in the diagram. Some time ago, NASA stupidly tried to generate power by moving a wire attached to a satellite in the earth's magnetic field, and thought they would get it for free! The mechanical reaction was unexpected and violent, and we have heard nothing more of this scheme.

The motional emf looks a lot like the magnetic term in the Lorentz force on a point charge, and there is a connection between them. However, the Lorentz expression is exact (for E and B in a given inertial system), while the formula for motional emf is an approximation for velocities much less than c.

It is appropriate to consider the Faraday homopolar disc generator at this time, since MHD analogies of it may exist. This is a conducting disk with sliding contacts at the axle and rim, which is rotated in a magnetic field normal to the plane of the disc. If the magnetic field is into the page, and the disc rotates clockwise, then the rim becomes positive relative to the axle. The radial electric field at a radius r is E = ωrB/c, so the generated emf is ωBa2/2c, where a is the radius of the disc and ω is the angular velocity of rotation. A table model 6" in diameter rotating at 1000 rpm in 1000 gauss will generate only 29.45 mV, but a 5 ft. model rotating in 10,000 gauss at 1000 rpm will generate a respectable 29.45 V. No practical homopolar generator of this kind has ever been designed, but there is nothing wrong with the idea.

### Hydrodynamic Equations

Now we shall assume that the electromagnetic fields are given, and study their effects on the fluid motion.

The decrease in the mass of fluid contained in a volume V is equal to the mass flux across the surface S bounding V. We take the positive direction of the normal to S as outward, so the outward mass flux of fluid of density ρ across S is ∫(S)ρv·ndS. By the divergence theorem, this is equal to the volume integral ∫(V)div (ρ)vdV, which must equal (d/dt)∫(V)ρdV = ∫(V)(∂ρ/∂t)dV. Hence, we have the differential equation of continuity ∂ρ/∂t + div(ρv) = 0. Any flow must necessarily satisfy this equation.

The density ρ may be a function of temperature and pressure: ρ = ρ(T,p), which is the equation of state. To avoid the complexities of thermodynamics, we assume that the flow processes are either isothermal (slow) or adiabatic (fast). Unless it is significant for the problem at hand, we shall not explicitly indicate the thermodynamic nature of the process. The bulk modulus k is defined by dρ/ρ = -dp/k, where k has the dimensions of the pressure p. The fractional change in density must be small for this linear approximation to be valid. For an isothermal process in an ideal gas, p = ρRT/M, where M is the molecular weight of the gas. In this case, k = p, and holds over a wide range of p. For an adiabatic process, pρ = constant, where γ is the ratio of the specific heat at constant pressure to the specific heat at constant volume, and is greater than 1. For adiabatic processes, k = γp.

For liquids, it is often satisfactory to set ρ = constant, since the density is nearly independent of pressure and temperature. The small changes that then occur will not influence the flow significantly. In this case, the equation of continuity becomes div v = 0.

The dynamics follows from Newton's Second Law, which is applied to an infinitesimal element of the fluid of mass ρdV moving with velocity v. The mass of this element (ρdV) times its acceleration, (dv/dt), will be equal to the net force acting upon it exerted by neighboring elements, or by external influences (fdV). The total time derivative must take into account the convection of the element, so d/dt = ∂/∂t + v·grad, which will be familiar to every student of hydrodynamics. This gives the equation of motion ρdv/dt = f. This is a nonlinear equation because of the appearance of the velocity in d/dt, which gives a term proportional to the square of the velocity. For slow motions, it is usual to neglect this term, which makes the equations more easily soluble.

Now we consider the various force terms. An important force in any fluid is the pressure force -grad p, where the isotropic pressure p is a scalar function of position. To derive this equation, consider the force on the surface S of a volume V. The total force is -∫(S)pndS (the normal points outward, the pressure presses inward!). By a theorem related to the divergence theorem (see the cover of Jackson), this is equal to -∫(V)(grad p)dV, which gives the desired result. The viscous forces, proportional to the rate of shear, are given by ηdiv grad v, where η is the dynamic viscosity in poise (dyne-s/cm2). These are the two forces arising from neighboring fluid particles.

A common external force is gravity, giving a force density of ρg, where the vector g points in the direction of the force. At the surface of the earth this is downward, of course. Another is the Coriolis force 2ω x v, where ω is the angular velocity of the rotating coordinate system to which v is referred. For the earth, ω = 7.29211 x 10-5 rad/s. There is also the centrifugal force in this case, but on the earth it is usually included in the gravity.

The external force of most interest to us is the ponderomotive force (1/c)J x B, which we discussed in connection with the electromagnetic equations. This is the term that couples the fields and the matter in MHD. It is worth noting that this force is always perpendicular to the magnetic field. Let us write the equation of motion, using indices, as ρdvi /dt = fi + (1/c)εijkJjBk, where fi represents all the non-electromagnetic forces. Now substitute Jj = σ(Ej + (1/c)εjlmvlBm) in this equation and simplify, using the formula for the contraction of two antisymmetric densities.

The result is ρdvi/dt = fi + (σ/c)εijkEjBk - (σ/c2)BkBkvi + (σ/c2)vkBkBi. Note the summation indices carefully in this equation. Let BkBk = B2, the square of the magnitude of the magnetic field. The second term can be written σ(B/c)2wi, where wi = εijkEjBk/B2 is a velocity in the direction of E x B. The final two terms can be written -σ(B/c)2[vi - vkBkBi/B2]. This can be recognized as the component of the fluid velocity perpendicular to the magnetic field, vpi. Finally, then, the equation of motion becomes ρdvi/dt = fi - Q(vpi - wi), where Q is the dragging coefficient ω(B/c)2.

This shows that the effect of the magnetic field on the fluid motion is to exert a force perpendicular to the magnetic field that tends to make the normal velocity equal to w, the "E x B drift." The higher the conductivity, the stronger is this force, and the more closely the magnetic field is dragged by the fluid (or vice-versa). However, the motion of the fluid along the magnetic field is unaffected.

We can also substitute for Ji in the ponderomotive force its expression in terms of the magnetic curl equation, Ji = (c/4π)εijkjBk. Note that we are making no assumption concerning the conductivity. Easy algebra gives us (1/4π)(BkkBi - BkiBk). It is important to keep straight what the partials operate upon in expressions like this. Now, ∂iB2 = ∂i(BkBk) = 2BiBk, so we can write the force density as -∂i(B2/8π) + (1/4π)BkkBi.

The first term is the force due to a magnetic pressure pm = B2/8π. This stress is familiar from the Maxwell stresses as the compressive stress normal to the direction of the field. It should be remembered that the Maxwell stresses do not give a pressure field, but that there is a tensile stress of the same magnitude in the direction of the field instead of a pressure. The interpretation of the second term is less obvious. Both Jackson and Cowling waffle on this, mentioning things that are not even true, and wishing the term would go away. However, it will not, and does play an important role. Since we made no assumption about the conductivity, the equation holds even with zero conductivity, when there is no magnetic force at all. In this case, the second term exactly cancels the magnetic pressure term, since now ∂iBk = ∂kBi.

If the conductivity is nonzero, then this term cancels the part of the magnetic pressure force that is along the lines of induction, since the field cannot cause any force in this direction on the fluid, as well as the force that allows the lines of induction to slip through the fluid. To see this more clearly, consider a magnetic field in the 3-direction. The only nonzero contribution to this term will be (1/4π)B33B3, which clearly cancels the force in the 3-direction from the magnetic pressure term. Other variations are necessary to correct the transverse forces, and these are not easy to illustrate in general. We can now see quite clearly the purpose of the correction term (1/4π)B·grad B in the equation of motion.

Where the flow is nonviscous, and the correction term is zero, and the magnetic field is frozen in the fluid, the equation of motion becomes simply ρdvi/dt = -∂i(p + pm + ψ), where ψ is the potential of forces like gravity. From this equation we can study a kind of magnetohydrostatics in the steady state, adding the magnetic pressure to the ordinary pressure. In any case, however, we must take care that the correction term is indeed zero, or the results will not be valid. We must also have high conductivity, or the transverse forces will also cancel in the steady state.

Suppose we have pipe flow with a constant axial magnetic field B(r). We presume that the axial current is zero, and in this case, the correction term is zero. In the steady state, p + pm = constant, so if p = po on the axis, where B = 0, then p will drop to zero at the radius where pm = po, or B = √(8πpo). This is the principle of magnetic confinement. For a pressure of 1 atm, a magnetic field of 5044 gauss is required. A solenoid would have to have 4000 A-t/cm to create this axial field, which would not be easy to do.

The usual "pinch" effect is confinement by the magnetic field produced by the longitudinal current flow itself. In this case, in the steady state -dp/dr = d(B2/8π)/dr + B2/4πr. The correction term appears quite prominently. The magnetic field is in rings around the current. Jackson shows that under certain reasonable assumptions, if the fluid is confined within some radius R, then the average pressure p = (1/2π)(I/Rc)2, where I is the total current. The magnetic field at radius R will be 2I/cR. This again gives p = B2/8π. If R = 1 cm, Jackson calculates that some 90,000 A would be required to pinch a typical 14-atm plasma.

As noted in Jackson, it is easy to see why magnetic confinement is unstable. If we consider the current restricted to a smaller cross-section, a "pinch," then the smaller radius means a larger magnetic field B = 2I/cr. This means a larger magnetic pressure, which further pinches the pinch. If the current deviates from a straight path in a "kink," then the lines of induction are brought closer together (B is increased) on the inside of the kink, and the opposite occurs on the outside. This further kinks the kink. These instabilities dismayed the people trying to confine a plasma in a possible fusion reactor.

### Hartmann Flow

A simple soluble problem in the flow of a conducting liquid in the presence of a magnetic field is solved in most references on MHD. This is the flow of an incompressible fluid between parallel planes with a magnetic field Bo normal to them, and is called Hartmann Flow after the investigator who first studied it theoretically and experimentally. The geometry and coordinates are shown in the figure. There can also be an electric field Eo in the y-direction, which produces a E x B drift velocity vxo = cEo/Bo of the fluid through the magnetic field. We will assume Eo = 0 so we do not have to carry it through the analysis. It can be added at the end, if desired.

This is a steady-state flow in the x-direction, supported by a pressure gradient dp/dx = -h. The flow velocity must be zero at z = 0 and z = a. The equation of motion for this problem is ∂ip = (1/c)εijkJjBk + η∂kkvi, and the equation of continuity is ∂kvk = 0. The current density is J2 = σ[Eo + (1/c)Bov1], where we will use only the second term to keep things simple. The velocity can be found from the x-equation only: -(h/η) = (M/a)2v + d2v/dz2, where the constant (M/a) = σBo2/ηc.

The dimensionless constant M is called the Hartmann Number, which indicates the ratio between the magnetic viscosity and the ordinary viscosity in the flow. In the equation of motion, the magnetic force is σ(B/c)2, while the viscous force is η/a2, the a factor standing in for the Laplacian operator with the same dimensions. The ratio of these is just M2. The case M = 0 corresponds to no magnetic field and normal viscous flow (provided the Reynolds number is low enough!). The equation is easily solved in this case, and we find the parabolic velocity profile v = (h/2η)z(a - z), as expected.

When M > 0, the equations are a little more difficult to solve, but it is just a linear equation of the second order with constant coefficients and a constant term. I will leave this as an exercise for the reader. Just find a particular solution with two adjustable constants, and apply the boundary condtions of v = 0 at z = 0 and z = a. After the use of several identities for the hyperbolic functions, the result can be cast into the convenient form v = (hc2/σBo2)sinh(Mz/a)sinh[((M/2)(1 - z/a)]/cosh(Mz/2a)cosh(M/2). For M→0, this gives the parabolic velocity profile, while for M→∞ it approaches a constant equal to the coefficient of the hyperbolic functions in the above expression, except in the immediate vicinity of the surfaces at z = 0 and z = a. As the Hartmann number increases, the velocity profile changes from a parabola to a constant value, showing the effect of dragging the magnetic field. For σ = ∞, we see that there is no flow at all!

It is now easy to find the x-component of the magnetic field and see how the lines of force are dragged by the flow, and to find the pressure gradient in the z-direction that adds to the magnetic pressure.

### Magnetohydrodynamic Waves

There is an interesting variety of magnetohydrodynamic waves, but the subject is one of complexity, so I will only try to give a good idea of the types of magnetohydrodynamic waves in the simplest cases, leaving all the complexity of finite conductivity, fluid viscosity and other factor aside. Let's assume that the fluid is in a constant applied magnetic field Bo, that the conductivity is high enough that the magnetic field is frozen in, that the fluid is compressible, but nonviscous, that the frequency is sufficiently low to eliminate displacement current or charge separation (plasma) effects, and, of course, that all the variable quantities in the wave are small enough to make linearization of the equations a satisfactory approximation. Even with all these simplifications, the algebra is still formidable. I shall use index notation, which is much more convenient than vector notation in this problem. One of the expressions would have a quadruple cross product, which would be unpleasant.

We begin with the case of a fluid in the absence of a magnetic field to get a taste of what has to be done. The continuity equation ∂ρ/∂t + ∂i(ρvi) = 0, the equation of motion ρ(∂vi/∂t + vkkvi) = -∂ip, and the equation of state p=p(ρ) are to be satisfied. We take the pressure as po + p, and the density as ρo + ρ, where now p << poand ρ << ρo are the small differences from the undisturbed state. Also, vi = ∂ξi/∂t, where ξi is the displacement in the wave, are small enough that their squares can be neglected.

We now have the linearized equations ∂ρ/∂t + ρoivi = 0, ∂vi/∂t = -∂ip/ρo, and p = -κρ/ρo = -κ∂kξk, where κ is the bulk modulus (we have a sufficient number of k's already). The equation of motion becomes ∂2ξi/∂t2 = c2ikξk, where c2 = κ/ρo. This is almost the wave equation, but not quite.

If we take the curl of the equation of motion, we find that (∂/∂t) curl v = 0. Integrating with respect to t, we find that curl v = constant. Since it is initially zero in the undisturbed fluid, it will continue to be zero. Integrating again, we find that curl ξ = 0 as well. This is an important property of waves in the medium: they are irrotational waves. We will see that this is not, in general, true for magnetohydrodynamic waves.

Now we introduce plane-wave solutions and find the conditions under which they satisfy the equations. We presume that the variable quantities depend on space and time through a factor ej(ωt - kixi). Then the partial derivative with respect to t can be replaced by jω, and ∂i by -jki. This makes the continuity equation ωρ - ρokivi, and the equation of motion ω2ξi = c2kikjξj. The vanishing of curl ξ becomes k x ξ = 0, which means that the displacement ξ is in the direction of k: the wave is longitudinal. Since then kiξj = kjξi, the right-hand side of the equation of motion becomes c2kjkjξi = ck2ξi. Therefore, the equation of motion is satisfied provided ω = ck, or c = ω/k, showing that c is indeed the phase velocity.

We now can express the velocity vi = jωξi, the density change ρ = jρokξ, and the pressure p = jc2kξ, in terms of the displacement, where k and ξ are the magnitudes of the vectors. We see that p = (κ/ρo)ρ, as it should. We shall call such a wave a sonic wave.

In the presence of a constant magnetic field Bo and high conductivity, the equation of continuity is the same, but the ponderomotive force is added to the right-hand side of the equation of motion. For the current in this expression, we use Ampére's Law. The ponderomotive force is then -(1/4π)εijkεklmBjlBm. Doing the contraction on the tensor density gives -∂i(B2/8π) + (1/4π)BjjBi. The third equation we need expresses the rate of change of B for large conductivity. This is ∂Bi/∂t = εijkεklmj(vlBm. Performing the contraction, we find ∂Bi/∂t = Bmmvi - vjjBi - Bijvj, where we have used ∂mBm = 0. It is necessary to be careful to note what the partials operate on.

These three equations are now linearized. Bi = Bio + bi, where b is the change in the magnetic field. We define V = Bo/√(4πρo), and call it the Alfvén velocity. The equation of motion now becomes ∂2ξi/∂t2 = c2ikξk - Vk{∂i[Vmmvk - Vkmξm] - ∂k[Vmmξi - Vimξm]. We now assume a plane wave solution, and write this equation in terms of ω and ki.

We now consider waves propagating in the direction of V, and perpendicular to it. In the latter case, we have kiVi = 0. Making use of this relation, we find that ω2ξi = (c2 + V2)kikjξj. This condition can be satisfied if k and ξ are in the same direction. In that case, kiξj = kjξi, and we find that ω2 = (c2 + V2)k2. This is a longitudinal wave whose phase velocity is the square root of the sum of the squares of c and V. This sort of wave is called a magnetosonic wave.

If k is in the direction of V, we write ki = (k/V)Vi, where (k/V) is the ratio of the magnitudes. Then the dispersion relation becomes [(ω/k)2 - V2i = (Viξi)[(c/V)2 - 1]Vi. Now we have two choices. Either ξ can be in the direction of V or normal to it. In the second case, the right-hand side vanishes, and the phase velocity becomes V. This transverse wave is called an Alfvén wave, and travels at the Alfvén velocity in the direction of the magnetic field. In the first case, ξ is in the direction of the magnetic field, and we have a longitudinal wave with phase velocity c. This is a normal sonic wave, and shows that the magnetic field does not affect motion parallel to it.

These are, then, the three pure types of magnetohydrodynamic waves: sonic, magnetosonic and Alfvén. In the magnetosonic waves, the magnetic field lines are pushed sideways, and the pressure between them aids the fluid pressure. In Alfvén waves, the magnetic field is strained from side to side and the tension in them gives the restoring force. There is no change in volume, so the fluid elasticity plays no part. In both cases, the magnetic field must be frozen into the field by high conductivity. The third wave type is a sonic wave in which the magnetic field plays no part. These waves all have definite polarizations with respect to the direction of propagation and the direction of the constant magnetic field.

### MHD and Flow Stability

A slow movement of a fluid may be characterized by streamlines, indicative of regular, orderly motion, called laminar flow. As the flow velocity is increased, or the space available made larger, there is a point where any small disturbance is amplified and the flow breaks up into turbulence. Then the velocity profile, and all the other flow characteristics, change radically. This is especially familiar in pipe flow, in the transition from laminar Poiseuille flow to turbulent Fanning-Hagen flow. The criterion is the dimensionless Reynolds number R = VLρ/μ, which indicates the relative importance of inertial and viscous forces. For flow in rough pipes, the transition occurs around R = 2500. We have already mentioned the analogous magnetic Reynolds number RM above.

A magnetic field may be expected to stabilize a flow against the transition to turbulent flow. The cases of the magnetic field in the direction of flow, and normal to the flow, have been studied, and the stabilizing effect observed. The critical Reynolds number is increased as the magnetic Reynolds number is increased. Experiments have been made that show this general result, but they are difficult and quantitatively not very conclusive.

Another kind of instability occurs when there is a temperature gradient and a pressure gradient in the same direction, such as occurs in the lower atmosphere or in a pot of liquid heated on a stove. If we imagine a small parcel of fluid from one level suddenly moved to a level where the pressure is different, it will expand or contract adiabatically (to a good approximation) until the pressures are equalized, and its temperature will change by a definite amount. The change in temperature as a function of difference of position defines the adiabatic lapse rate, to use meteorological terminology. If the actual lapse rate (temperature gradient) is greater, then if the parcel has risen, it will find itself hotter and less dense than its surroundings. Buoyancy will then encourage it to rise even higher. Should the parcel have sunk, then things are exactly reversed, and it will find itself cooler and denser, so it will be impelled to sink further. This is illustrated in the diagram. A parcel at 1 just arrived from 0 is hotter than its environment 1', while a parcel at 2 is cooler than its environment at 2'. If, on the other hand, the lapse rate is less than adiabatic, then a displaced parcel will be encouraged to return to its original position. In the former case, we have convective instability, and in the latter, convective stability. Indeed, the lower part of the atmosphere, in which the decrease in temperature with height is usually more rapid than the adiabatic decrease, and which is, therefore, convectively unstable, is called the troposphere ("sphere of turning"). The unusual opposite case is called an inversion. Convection is important not only in the troposphere and in pots, but also deep in the earth, and in the sun, where energy is struggling to get out.

A magnetic field in the direction of the temperature and pressure gradient will hinder the transverse motion essential to convection, and make the convective "cells" narrower and less efficient, reducing the rate of energy transfer. In the Sun, this is seen by the relative darkness of sunspots, where a vertical magnetic field of thousands of gauss reduces the efficiency of convection of heat from below, so the surface cools below the general level of the photosphere.

### MHD and the Sun

It was noted above that the only conducting fluids available for laboratory experiments are mercury and liquid sodium, both inconvenient for different reasons. It is very difficult to reach large magnetic Reynolds numbers in laboratory experiments, and so to verify important theoretical results with any accuracy. There is another MHD laboratory available, however, and that is the Sun. The surface of the sun is a hot, relatively dense plasma where the magnetic Reynolds number is very large, so the magnetic field is well and truly frozen into this fluid. However, we cannot change the experimental parameters, and we do not know what is going on below the level that we can see, so it is a less than perfect laboratory. Nevertheless, there are many interesting and varied phenomena that show the influence of MHD very well.

The visible surface of the sun is the photosphere, a fuzzy surface that we can see a little ways into, as into a cloud. It is composed of 90% hydrogen, 10% helium, with all the other elements as traces. Here energy rising from below encounters a steep declining temperature gradient as the surface layers cool by radiation into space. Convective instability produces active convection, with cells ("granules") averaging 700 miles in diameter, hot plasma rising in the center and sinking at the periphery after cooling a little. We observe a temperature of about 6000K looking straight down, closer to 5000K at the limb where we look obliquely and not so deep. The plasma emits a black-body spectrum, acting as an almost perfect emitter. The plasma cools rapidly above this level, becoming less dense and more transparent. In its lower temperature of about 4500K, neutral atoms absorb their characteristic lines, creating the Fraunhofer lines in the solar spectrum. This region is called the reversing layer, and extends some hundreds of miles vertically.

In solar eclipses, the moon covers the photosphere and the lower part of the reversing layer, but the upper part of this region is seen in its bright red Hα radiation at 656.3 nm. which dominates the flash spectrum observed at this instant, when all the solar Fraunhofer lines become bright emission lines. For this reason, this region is called the chromosphere, the lower part of the solar atmosphere. Above the reversing layer, the chromosphere becomes rarer and hotter, merging with the whitish corona above 13,000 miles. This is a rare plasma, in which the particles have very high kinetic energies, characteristic of 106K, but there is no thermal equilibrium here, so temperature is a doubtful concept. The magnetic Reynolds number is high throughout the solar atmosphere, and matter clings to the lines of force, sliding along them freely.

The Sun is almost perfectly spherical, subtending an angle of about 32' at mean distance. The equator is inclined 7° to the ecliptic, and we see the north pole in September, the south pole in March, from our position 1.495977 x 1013 cm distant. At the top of the atmosphere, we receive 1.37 x 106 erg/s-cm2, which corresponds to a total luminosity of 3.8 x 1033 erg/s. The average density of the sun is 1.4 g/cm2, about the same as Jupiter's. Its total mass is 2 x 1033 g, which makes its surface gravity about 28 times stronger than earth's. It is curious that its equatorial period of rotation is 25 days (27 days, as seen from earth that is revolving about the sun), but it rotates more and more slowly towards the poles. At a latitude of 45°, the rotational period is 28 days, and has been said to approach 33 days at the poles.

This difference in rotational periods is noted without much comment in astronomy texts, but it is really an extraordinary thing. In an ordinary liquid body, viscous forces would transfer angular momentum to the polar regions until they rotated at the same speed. In the sun, this momentum transfer also must occur, but the polar regions do not speed up. This means that some mechanism transfers angular momentum back to the equatorial regions at the same rate, and this mechanism can only be some magnetohydrodynamic effect, perhaps only in the outer layers of the sun. One side effect of this process may even be the appearance of sunspots when the necessary magnetic fields reach the sun's surface.

Sunspots, indeed, are evidence of solar magnetic fields. Magnetic flux enters or leaves the sun vertically at a sunspot, and the field magnitude can be thousands of gauss. We have already mentioned how this hinders convection in the region of high field, making the area cooler and darker than surrounding areas. The magnetic field of sunspots was discovered by Hale in 1908, and has been an interesting field of speculation since then. Sunspots occur in pairs of opposite polarity at roughly the same latitudes, so it seems that the field pierces the surface, bends over, and descends again in the companion spot. Polarities of the leading and following spots are reversed in the opposite hemisphere, and the polarities reverse in each successive 11.2-year sunspot cycle. The cycle is actually 22.4 years long, with two maxima in each cycle. No explanation is known for this cycle, but it is almost certainly magnetohydrodynamic. Sunspots do not occur in polar or in equatiorial regions, but are restricted to mid-latitudes. Each maximum begins at high latitudes, and approaches the equator as the cycle progresses. The bigger the sunspot, the longer it lives. Large sunspots can live for weeks, and appear to cross the sun's disc more than once, but most do not. An average sunspot is 1000 miles in diameter, while a really large one could hold the earth.

There is a number of phenomena surrounding sunspots and associated with the chromosphere and corona that make magnetic fields manifest. There are spicules and faculae and flocculae, coronal streamers and rays, all apparently supported by magnetic fields. The most impressive and long-lived are the filaments or prominences, best seen in elevation at the limb of the sun in their reddish Hα light, great arches of luminous matter supported by the magnetic field. The field itself seems to be the longest-lived of all the spot phenomena. The form of the coronal streamers suggested a dipole field, but the sun has no strong general field, which cannot be stronger than a gauss or two. Originally, the sun was suspected of having a general field as large as 25 gauss, but that has proved erroneous.

Solar flares are probably the result of some process of MHD acceleration of plasma in the active regions that breed sunspots. They are very hot, beginning with a blast of ultraviolet radiation that reaches the earth in a few minutes, creating extra ionization in the upper atmosphere. Then a jet of plasma is ejected that reaches the earth in about a day, causing ionospheric currents and ionization with magnetic storms, strong earth currents, disturbances to radio communication and the aurora borealis.

### References

J. D. Jackson, Classical Electrodynamics, 2nd ed. (New York: John Wiley & Sons, 1975). Chapter 10.

T. G. Cowling, Magnetohydrodynamics (New York: Interscience, 1957).

C. Kittel, Introduction to Solid State Physics, 3rd ed. (New York: John Wiley & Sons, 1966). Chapter 7.

R. V. Churchill, Fourier Series and Boundary Value Problems (New York: McGraw-Hill, 1941). pp 112-113.