Overview of Magnetohydrodynamics Theory in Toroidal Plasma Confinement

The biggest challenges mankind ever faces are falling energy sources and food shortages. If controlled nuclear fusion were achieved with net energy yield, the energy source problem would be solved. If natural photosynthesis were reproduced, food shortage concern would be addressed. Though both nuclear fusion and photosynthesis are universal, the difficulties to achieve them are disproportionally great. Instead, those discoveries harmful to nature, though naturally unpopular, are invented relatively easily. This tendency reminds us of a bible verse (Genesis 3:19): "In the sweat of thy face shalt thou eat bread". This verse basically sketches the dependence of efforts (sweat) demanded for scientific discoveries on the usefulness (bread) of the discoveries to mankind (see Fig. 1). The more the discovery is relevant to mankind, the more the sweat is needed for that discovery. This may explain why controlled nuclear fusion is so difficult and its underlying plasma physics is so complicated.


Introduction
In this chapter we address magnetohydrodynamics (MHD) theory for magnetically confined fusion plasmas. To be specific we focus on toroidal confinement of fusion plasmas, especially tokamak physics.
The biggest challenges mankind ever faces are falling energy sources and food shortages. If controlled nuclear fusion were achieved with net energy yield, the energy source problem would be solved. If natural photosynthesis were reproduced, food shortage concern would be addressed. Though both nuclear fusion and photosynthesis are universal, the difficulties to achieve them are disproportionally great. Instead, those discoveries harmful to nature, though naturally unpopular, are invented relatively easily. This tendency reminds us of a bible verse (Genesis 3:19): "In the sweat of thy face shalt thou eat bread". This verse basically sketches the dependence of efforts (sweat) demanded for scientific discoveries on the usefulness (bread) of the discoveries to mankind (see Fig. 1). The more the discovery is relevant to mankind, the more the sweat is needed for that discovery. This may explain why controlled nuclear fusion is so difficult and its underlying plasma physics is so complicated. Fig. 1. Schematic interpretation of Genesis 3:19: the dependence of efforts (sweat) demanded for scientific discoveries on the usefulness (bread) of the discoveries to mankind.
When God created universe each day, He always claimed "it was so" and "it was good" (Genesis 1). The hard aspect of fusion plasma physics lies in that we often miss simplicity (it was so) and beauty (it was good) in theoretical formalism. MHD theory seems to be unique in plasma physics. Though many charged particle system, with long mean free path and long range correlation, is intrinsically complicated, MHD theory is relatively simple and, nonetheless, gives rise to rather relevant theoretical predictions for experiments: Tokamaks are designed according to MHD equilibrium theory and nowadays none would expect that a magnetic confinement of fusion plasmas could survive if MHD theory predicted major instabilities. As one will see even with MHD description the theoretical formulation of magnetically confined plasmas in toroidal geometry can still be hard to deal with. Thanks to decades-long efforts many beautiful MHD theoretical formulations for toroidal confinement of fusion plasmas have been laid out. In this chapter we try to give a comprehensive review of these prominent theories. Four key types of modes: interchange/peeling modes Mercier (1962) Greene & Johnson (1962) Glasser et al. (1975) Lortz (1975) Wesson (1978), ballooning modes Connor et al. (1979) Chance et al. (1979), toroidal Alfvén egenmodes (TAEs) Cheng et al. (1985) Rosenbluth et al. (1992) Betti & Freidberg (1992) Zheng & Chen (1998), and kinetically driven modes, such as kinetic ballooning modes (KBMs) Tsai & Chen (1993) Chen (1994) and energetic particle driven modes (EPMs) Zheng et al. (2000), are addressed. Besides, we also describe an advanced numerical method (AEGIS Zheng & Kotschenreuther (2006)) for systematically investigating MHD stability of toroidally confined fusion plasmas. Description of global formulation used for numerical computation can also provide an overall picture of MHD eigen mode structure for toroidal plasmas.
MHD theory is a single fluid description of plasmas. Fluid approach is based on the assumption that particle movements are spatially localized so that a local thermal equilibrium can be established. In the conventional fluid theory particle collision is the ingredient for particle localization. However, for magnetically confined fusion plasmas collision frequency usually is low. One cannot expect particle collisions to play the role for localizing particles spatially. The relevance of partial fluid description of magnetically confined fusion plasmas relies on the presence of strong magnetic field. Charged particles are tied to magnetic field lines due to gyro-motions. Therefore, in the direction perpendicular to magnetic field lines magnetic field can play the role of localization, so that MHD description becomes relevant at least in lowest order. One can expect that perpendicular MHD description needs modification only when finite Larmor radius effects become significant.
In the direction parallel to magnetic field, however, particles can move rather freely. Collision frequency is not strong enough hold charged particles together to establish local thermal equilibrium. One cannot define local thermal parameters, such as fluid density, velocity, temperature, etc. The trapped particle effect, wave-particle resonances, and parallel electric field effects need to be included. Plasma behavior in parallel direction is intrinsically non-fuild and needs kinetic description. Surprisingly, even under this circumstance MHD description still yields valuable and relevant theoretical predictions without major modifications in the concerned low (ω ≪ ω si ) and intermediate (ω si ≪ ω ≪ ω se ) frequency regimes, where ω is mode frequency, ω si and ω se represent respectively ion and electron acoustic frequencies. In the low frequency regime coupling of parallel motion results only in an enhanced apparent mass effect; while in the intermediate frequency regime kinetic effects only gives rise to a new phenomenological ratio of special heats in leading order. This chapter is arranged as follows: In Sec. 2 the basic set of ideal MHD equations is described; In Sec. 3 MHD equilibrium is discussed; In Sec. 4 analytical or semi-analytical theories for four types of major MHD modes are presented; In Sec. 5 the formulation of global numerical analyses of MHD modes are given; In the last section the results are summarized. Gyrokinetic and resistive effects are also discussed in this last section.

Basic set of ideal MHD equations
The basic set of ideal MHD equations are derived from single fluid and Maxwell's equations. They are given as follows where ρ m is mass density, v denotes fluid velocity, P is plasma pressure, Γ represents the ratio of specific heats, E and B represents respectively electric and magnetic fields, J is current density, μ 0 is vacuum permeability, and bold faces denote vectors.
The MHD equations (1)-(6) can be linearized. For brevity we will use the same symbols for both full and equilibrium quantities. Perturbed quantities will be tagged with δ, unless specified. Equilibrium equations are The linearized perturbed MHD equations become where ξ = v/(−iω) represents plasma displacement, and the time dependence of perturbed quantities is assumed to be of exponential type exp{−iωt}. Inserting Eqs. (11)-(13) into Eq. (10), one obtains a single equation for ξ: We first outline the derivation of Grad-Shafranov equation Grad & Rubin (1958) Shafranov (1966. The cylindrical coordinate system (X, Z, φ) is introduced, where X is radius from axi-symmetry axis of plasma torus, Z denotes vertical coordinate, and φ is toroidal axi-symmetric angle. We introduce the vector potential A to represent magnetic field B = ∇×A. Due to toroidal symmetry φ is an ignorable coordinate. Using curl expression in cylinder coordinates and noting that ∂A X /∂φ = ∂A Z /∂φ = 0, one can prove that the vector potential A in X and Z directions (A XZ ) can be expressed through the single toroidal component: A φ . Without losing generality one can express A φ = −χ∇φ. Therefore, total equilibrium magnetic field can be expressed as, by adding (X, Z) components and toroidal component, where B φ is toroidal component of magnetic field and g = XB φ . From Eq. (15) one can prove that B ·∇χ = 0 and therefore χ = const. labels magnetic surfaces. Equation (15) can be used to show that 2πχ is poloidal magnetic flux. One can also define the toroidal flux 2πψ T (χ). The safety factor is then defined as q = dψ T /dχ, which characterizes the field line winding on a magnetic surface.
Inserting Eqs. (15) and (16) into force balance equation (7) one obtains From Eq. (7) one can prove that B ·∇P = 0. Therefore, one can conclude that plasma pressure is a surface faction, i.e., P(χ). From Eq. (17) one can further determine that g is a surface function as well, through projecting Eq. (17) on ∇φ. Therefore, Eq. (17) can be reduced to the so-called Grad-Shafranov equation Here and later on we use prime to denote derivative with respect to flux coordinate chosen. This is a nonlinear equation for χ for given functions P(χ) and g(χ). It generally needs numerical solution. Since it is a two dimensional problem, one needs to introduce a poloidal angle coordinate θ eq around magnetic axis of plasma torus in addition to radial coordinate χ.
Instead of physical cylinder coordinates (X, Z, φ)o r( χ, θ eq , φ), magnetic flux coordinates are often used in theoretical analyses, which is characterized by that the magnetic field line is straight in the covariant representation of coordinate system. Note that the coordinate system (ψ, θ eq , φ) usually is not a flux coordinate system. In most equilibrium codes θ eq is just an equal-arc length poloidal coordinate. One of flux coordinate systems is the so-called PEST coordinate system Grimm et al. (1976) (χ, θ pest , φ), where θ pest is generalized poloidal coordinate, such that the equilibrium magnetic field can be represented as By equating Eqs. (19) and (15) one can find that the Jacobian of PEST coordinates should be In PEST coordinate system the flux coordinate is chosen as where X 0 is major radius at magnetic axis and v dτ denotes volume integration over entire plasma domain. The PEST poloidal angle θ pest can be related to physical angle coordinate θ eq as follows. Using Eq. (20), one can determine poloidal angle in PEST coordinate where J eq = 1/∇χ ×∇θ eq ·∇φ, which can be computed from equilibrium solution. Here, the integration is along the path of constant χ and φ.
Next, we discuss construction of general flux coordinates. The covariant type of representation as in Eq. (19) is not unique. It is preserved under the following coordinate transforms such that Here, θ and ζ are referred to as generalized poloidal and toroidal angles, respectively. PEST coordinates are characterized by its toroidal angle coordinate being axisymmetric toroidal angle. In this general case, by equating Eqs. (22) and (15) in ∇φ projection one can find that ∂ν ∂θ eq ψ,φ 6 Will-be-set-by- IN-TECH where J = 1/∇ψ ×∇θ ·∇ζ. Using J eq and J definitions, one can prove that One can solve Eq. (23), yielding where Eq. (24) has been used.
Equations (21)-(25) can be used to construct various types of flux coordinate systems. There are two classes of them: One is by specifying Jacobian (e.g., Hamada coordinates Hamada (1962) and Boozer coordinates Boozer (1982)) and the other by directly choosing generalized poloidal angle (e.g., equal arc-length coordinate). In the Hamada coordinates the volume inside a magnetic surface is used to label magnetic surfaces, i.e., ψ = V, and Jacobian J h = 1/∇V ·∇θ h ×∇ζ h is set to be unity. With Jacobian specified, Eq. (24) can be used to solve for θ h at given (V,φ). With ν determined by Eq. (25) the definition Eq. (21) can be used to specify ζ h . In the Boozer coordinates Jacobian is chosen to be J B = V ′ B 2 s /(4π 2 B 2 ), where · s represents surface average. The procedure for specifying Boozer poloidal and toroidal coordinates θ B and ζ B is similar to that for Hamada coordinates. In the equal-arc-length coordinates poloidal angle is directly specified as equal-arc-length coordinate θ e . In this case, Jacobian J e can be computed through Eq. (24). With ν determined by Eq. (25) the definition Eq. (21) can be used to specify ζ e .
We can also express current density vector in covariant representation with generalized flux coordinates. Using Ampere's law in Eq. (8) for determining J ·∇θ and Eq. (7) for J ·∇ζ, one can also express equilibrium current density in covariant representation This general coordinate expression for J can be alternatively obtained from PEST representation in Eq. (16) and Grad-Shafranov equation (18) where I(V) and J(V) are toroidal and poloidal current fluxes, The force balance equation (7) can be simply expressed as It is also interesting to discuss diamagnetic current and Pfirsch-Schlüter current in plasma torus. Due to the existence of plasma pressure there is diamagnetic current in tokamak system. The diamagnetic current alone is not divergence-free and is always accompanied by a return current in the parallel direction, i.e., the so-called Pfirsch-Schlüter current. The total equilibrium current is therefore can be expressed as where the second term is diamagnetic current and the first term denotes the Pfirsch-Schlüter current. We can determine the Pfirsch-Schlüter current from ∇·J = 0, where λ 0 is the integration constant and can be determined by Ohm's law in the parallel direction.

Linear MHD instabilities
In this subsection we overview the linear MHD stability theories in toroidal geometry. We will detail major analytical techniques developed in this field in the past decades, such as interchange, ballooning, TAE, and EPM/KBM theories. Due to space limitation, we focus ourselves on ideal MHD theory.

Decomposition of linearized MHD equations, three basic MHD waves
There are three fundamental waves in magnetic confined plasmas. The compressional Alfvén mode characterizes the oscillation due to compression and restoration of magnetic field. It mainly propagates in the derection perpendicular to magnetic field. Since plasmas are frozen in magnetic field, such a magnetic field compression also induces plasma compression. Note that the ratio of plasma pressure to magnetic pressure (referred to as plasma beta β) usually is low. The compression and restoration forces mainly result from magnetic field energy. The shear Alfvén wave describes the oscillation due to magnetic field line bending and restoration. It mainly propagates along the magnetic field lines. Since long wave length is allowed for shear Alfvén wave, shear Alfvén wave frequency (or restoration force) is usually lower than that of compressional Alfvén wave. Therefore, shear Alfvén wave is often coupled to plasma instabilities. Another fundamental wave in magnetic confined plasmas is parallel acoustic wave (sound wave). Since plasma can move freely along magnetic field lines without being affected by Lorentz's force. Parallel acoustic wave can prevail in plasmas. The various types of electrostatic drift waves are related to it. Due to low beta assumption, the frequency of ion sound wave is lower than that of shear Alfvén wave by oder β. The behaviors of these three waves in simplified geometry have been widely studied in many MHD books. Here, we focus on toroidal geometry theories. MHD equation (14) in toroidal geometry can be hard to deal with. One usually needs to separate the time scales for three fundamental waves to reduce the problem. This scale separation is realized through proper projections and reduction of MHD equation (14).
There are three projections for MHD equation, Eq. (14). We introduce three unit vectors: e b = B/B, e 1 = ∇ψ/|∇ψ|, and e 2 = e b × e 1 for projections. The e 2 projection of the MHD equation (14) gives Similarly, the e 1 projection of the MHD equation (14) yields The e b projection of MHD equation (14) can be reduced to, using ∇·ξ as an independent unknown, Noting that δJ and δB are determined completely by ξ ⊥ , one can see that the set of equations (31) -(33) is complete to determine two components of ξ ⊥ and scalar unknown ∇·ξ.

Singular layer equation: interchange and peeling modes
Interchange modes are most fundamental phenomena in magnetically confined plasmas. It resembles to the so-called Rayleigh-Taylor instability in conventional fluid theory. Through interchange of plasma flux tubes plasma thermal energy can be released, so that instability develops. Perturbation of magnetic energy from field line bending is minimized for interchange instability. In slab or cylinder configurations such an interchange happens due to the existence of bad curvature region. In toroidal geometry with finite q value, however, the curvature directions with respect to plasma pressure gradient are different on high and low field sides of plasma torus. Good and bad curvature regions appear alternately along magnetic field line. Therefore, one needs to consider toroidal average in evaluating the change of plasma and magnetic energies. This makes interchange mode theory in plasma torus become complicated. The interchange mode theory is the first successful toroidal theory in this field. It includes the derivations of the so-called singular later equation and interchange stability criterion, i.e., the so-called Mercier criterion Mercier (1962) Greene & Johnson (1962).
Early derivation of singular layer equation relies on the assumption that the modes are somewhat localized poloidally. This assumption was released in a later paper by Glasser et al Glasser et al. (1975). However, the details have been omitted in this paper and direct projection method, alternative to the original vorticity equation approach, is used. Here, we detail the derivation of singular layer equation by vorticity equation approach. These derivation can tell analytical techniques to separate the compressional Alfvén wave from low frequency interchange mode and to minimize field line bending effects due to shear Alfvén mode. The singular equation will be used to derive stability criterion for interchange and peeling modes.
In order to investigate the modes which localize around a particular rational (or singular) magnetic surface V 0 , we specialize the Hamada coordinates to the neighborhood of mode rational surface V 0 and introduce the localized Hamada coordinates x, u, θ as usual, where x = V − V 0 and u = mθ − nζ. In this coordinate system the parallel derivative becomes Using the coordinates (x, u, θ), we find that, in an axisymmetric torus, equilibrium scalars are independent of u, and therefore perturbations can be assumed to vary as exp{ik u u} with k u = 2πn/χ ′ . As in Refs. Johnson & Greene (1967) and Glasser et al. (1975), ξ and δB are projected in three directions as follows: We consider only singular modes whose wavelength across the magnetic surface λ ⊥ is much smaller than that on the surface and perpendicular to magnetic field line λ ∧ . This leads us to choose following ordering scheme as in Ref. Glasser et al. (1975): where ǫ ≪ 1, being a small parameter. Furthermore, we consider only the low-frequency regime where ω si is parallel ion acoustic frequency.
Since the modes vary on a slow time scale, they are decoupled from compressional Alfvén wave. It can be verified a posteriori that we can make following ordering assumptions: where δP (2) = −ΓP∇·ξ. These ordering assumptions are the same as those in Ref. Glasser et al. (1975), except that we use δP (2) as unknown to replace ν. With these ordering assumptions we can proceed to analyze the basic set of linearized MHD equations. As usual, perturbed quantities are separated into constant and oscillatory parts along the field lines: ξ =ξ +ξ, whereξ = ξ ≡ ξ dl/B/ dl/B, l is arc length of magnetic field line, and ξ = ξ − ξ .  The condition that δB be divergence free, as required by Eq. (11), yields After surface average it gives The two significant orders of induction equation, Eq. (11), in the ∇V-direction are The component of Eq. (11) in the ∇u-direction, in lowest order, yields To satisfy the component of Eq. (11) along the magnetic field line, one must set where Eq. (41) and ∇·(B ×∇V/B 2 )=2B × κ ·∇V/B 2 have been used, and κ = b ·∇b is magnetic field line curvature.
Next, we turn to momentum equation (14). The two components perpendicular to B of the momentum equation (14) both lead, in lowest order, to This is consistent to Eq. (42). Since both components yield the same information, we can directly work on the vorticity equation Eq. (34) and obtain 10 Topics in Magnetohydrodynamics www.intechopen.com

Overview of Magnetohydrodynamics Theory in Toroidal Plasma Confinement 11
We will derive the singular layer equation by averaging this equation. Therefore, it is needed to express unknowns in this equation in terms of ξ (1) .
From Eqs. (40) and (37) We need also to solve the equation of parallel motion, Eq. (33). Taking into consideration of low frequency assumption in Eq. (36) and the result in Eq. (42), the equation of parallel motion can be reduced to Noting that B ×∇V · κ/B 2 = χ ′ (∂σ/∂θ), equation (47) can be solved to yield Inserting these results into Eq. (45) and averaging over l, one obtains the singular layer equation where the total mass parameter M = M c + M t , Here, the mass factor M c results from perpendicular motion and M t from parallel motion due to toroidal coupling. M t is often referred to as apparent mass. In the kinetic description the apparent mass is enhanced by the so-called small parallel ion speed effect. In the large aspect ratio configurations this enhancement factor is of order √ R/a, where R and a are respectively major and minor radii Mikhailovsky (1974) Zheng & Tessarotto (1994b).
From Eq. (48) one can derive the Mercier criterion, i.e., the stability criterion for localized interchange modes in toroidal geometry. In the marginal stability ω 2 = 0, Eq. (48) becomes the Euler differential equation. Its solution is The system stability can be determined by Newcomb's theorem 5 Newcomb (1960): system is unstable, if and only if the solution of Eq. (48) vanishes two or more points. From the solution in Eq. (49) one can see that if −D I < 0 ξ becomes oscillated. Therefore, interchange mode stability criterion is simply −D I > 0.
Interchange modes are internal modes. When internal modes are stable, it is still possible to develop unstable external modes. For external modes one needs to consider the matching condition between plasma and vacuum solutions. As discussed in conventional MHD books, these matching conditions are: (1) the tangential magnetic perturbation (δB t ) should be continuous; and (2) total magnetic and thermal force (B · δB + δP) should balance across plasma-vacuum interface in the case without plasma surface current. It can be proved that for localized modes the vacuum contribution is of order ǫ 2 and therefore can be neglected Lortz (1975). Consequently, the boundary condition becomes that total magnetic and thermal forces on the plasma side of the plasma-vacuum interface should vanish. This gives the necessary and sufficient stability condition for peeling modes where b is the coordinate of plasma-vacuum interface, relative to the rational surface, and Note that the stability condition Eq. (50) can be alternatively obtained by the approach of minimization of plasma energy Lortz (1975) Wesson (1978. One can derive the peeling mode stability criterion by inserting Eq. (49) into Eq. (50) Wesson (1978). In the derivation of peeling stability criterion we assume system to be Mercier stable, i.e., −D I > 0. For the case with ∆ < 0 we assume that rational surface resides inside plasma region, so that b > 0. In this case the stability condition becomes For the case with ∆ > 0 we assume that rational surface resides outside plasma region, so that b < 0. In this case the stability condition becomes Therefore, both cases, Eqs. (51) and (52), give rise to the same stability criterion for peeling mode: Λ P < 0. This is more stringent than the Mercier criterion.

Ballooning modes
In this section we review high-n ballooning mode theory. The stability criterion for interchange modes takes into account only average magnetic well effect. As it is well-known tokamak plasmas have good and bad curvature regions, referring to whether pressure gradient and magnetic field line curvature point in same direction or not. Usually bad curvature region lies on low field side of plasma torus; good curvature region on high field side. Although tokamaks are usually designed to have average good curvature, i.e., Mercier stable, the ballooning modes can still develop as soon as the release of plasma thermal energy on bad curvature region is sufficient to counter the magnetic energy resulting from field line bending Connor et al. (1979) Chance et al. (1979). In difference from interchange modes ballooning modes have high toroidal mode number n, while interchange modes can be either low and high n. Also, ballooning modes allow normal and geodesic wave lengths to be of same order λ ⊥ ∼ λ ∧ , but both of them are much smaller than parallel wave length λ .
We first derive ballooning mode equation. In high n limit, both components of perpendicular momentum equation, Eqs. (31) and (32), give the same result In lowest order, one has ∇·ξ ⊥ ∼ ξ/R. This allows to introduce the so-called stream function δϕ Chance et al. (1979): Equation (33), meanwhile, can be reduced to where Eq. (53) has been used.
The key formalism to ballooning mode theory is the so-called ballooning representation Lee & Van Dam (1977) Connor et al. (1979. Here, we outline its physics basis and derivation,

13
Overview of Magnetohydrodynamics Theory in Toroidal Plasma Confinement www.intechopen.com especially to explain the equivalence of two kinds of representations in Refs. Lee & Van Dam (1977) and Connor et al. (1979). In tokamak geometry one can introduce the following Fourier decomposition: For simply to describe ballooning mode representation we have used nq as flux surface label. This is allowed for systems with finite magnetic shear, in which the ballooning representation applies. For high n modes the distance of mode rational surfaces is of order 1/n, which is much smaller than equilibrium scale length. Therefore, in lowest order we can neglect the spatial variance of equilibrium quantities and require mode Fourier harmonics to have the so-called ballooning invariance: so that the Fourier decomposistion in Eq. (57) can be expressed as We can further introduce the Laplace tranform Using this transform Eq. (58) can be written as Noting that 1 2π This indicates that we can represent high n modes at a reference surface as without concern of periodicity requirement. Here, β ≡ ζ − qθ. The periodic eigenfunction can always be formed through the summation in Eq. (61). This representation characterizes the most important feature of ballooning modes in a plasma torus that perpendicular wave number is much larger than parallel one: k ⊥ ≫ k . This reduction shows the equivalence of two kinds of representations in Eqs. (58) and (61)  Uniqueness and inversion of ballooning mode representation were proved in Ref. Hazeltine et al. (1981).
With ballooning mode representation described, we can proceed to derive ballooning mode equation. It is convenient to use the so-called Celbsch coordinates (ψ, β, θ) to construct ballooning mode equations. In this coordinates ∇→−in∇β and B ·∇= χ ′ (∂/∂θ). Applying Eq. (62) to Eqs. (54) and (55) and employ the high n ordering, one can obtain following coupled ballooning mode equations where δΞ = i∇·ξ/n. These two equations are coupled second order differential equations. The derivatives here are along a reference magnetic field line labeled by ψ and β. The boundary conditions are δϕ, δΞ → 0atθ →±∞ to guarantee the convergence of the Laplace transform in Eq. (59).
In studying ballooning stability at finite beta equilibrium, the so called steep-pressure-gradient equilibrium model is often used Connor et al. (1978) Greene & Chance (1981). In this model, finite beta modification is only taken into account for magnetic shear, while others remain to their low beta values. This model has been proved to be successful for ballooning mode studies. Here, we outline the formulation in Ref. Berk et al. (1983). Noting that β = ζ − qθ, one can see that the magnetic shear effect resides at the quantity ∇β in the ballooning mode equations (63) and (64). From Eq. (22) one can prove that where Λ s is the so-called shear parameter and can be obtained by applying operator B ×∇χ ·∇× ··· on Eq. (65), We need to determine finite beta modification to Λ. We assume that χ = χ 0 + χ 1 and β = β 0 + β 1 , where χ 0 and β 0 are low beta values and χ 1 and β 1 represent finite beta modifications. The linearized Ampere's law can be written as follows: Noting that in the curl operation on left hand side only the gradient component in ∇χ direction needs to be taken, i.e., ∇× →∇χ 0 ∂/∂χ×, equation (67) can be solved Now taking the dot product of Eq. (68) with ∇β 0 gives We can remove subscript 0 afterward for brevity. Now substituting Eq. (69) for ∂Q/∂χ and noting that ∂β 1 /∂χ ≡ Λ s1 , one finds that Therefore, the shear parameter can be evaluated as follows The second term here gives rise to the finite beta modification to shear parameter Λ s in steep pressure gradient model. The rest parameters here and in ballooning equations (63) and (64) can be evaluated with low beta values.
We now consider tokamak model equilibrium with circular cross section, low beta, and large aspect ratio (i.e., 1/ǫ = R/a ≫ 1). The magnetic field in this model can be expressed as B = B φ (r)/(1 + ǫ cos θ)e φ + B θ (r)e θ . The shear parameter can be expressed as Λ s = s(θ − θ k ) − α sin θ. Here, α = −(2Rq 2 /B 2 )(dP/dr), s = d ln q/d ln r,andθ k is integration constant. Therefore, ballooning equations (63) and (64) can be reduced to To further analyze this set of equations it is interesting to consider two limits: the low frequency (ω ≪ ω si ) and intermediate frequency limit (ω si ≪ ω ≪ ω se ). In the low frequency limit the second term on left hand side of Eq. (72) can be neglected and inertia term is only important in the outer region θ → ∞. Equation (72) can be solved to yield δΞ = 2ρ m q 2 ω 2 ΓPB θ sθ sin θδϕ. Here, it has been considered that in sound wave scale the slow variable sθ can be regarded as constant. Inserting Eq. (73) into Eq. (71) yields that Here, we see that the sound wave coupling results in the so-called apparent mass effect: i.e., the inertia term is enhanced by a factor (1 + 2q 2 ) Greene & Johnson (1962). In the kinetic description the 2q 2 term is further boosted by the so-called small particle speed effect to become of order 2q 2 / √ r/R for large aspect ratio case Mikhailovsky (1974) Zheng & Tessarotto (1994b). In the marginal stability ω 2 = 0 the ballooning stability can be determined by Newcomb's theorem 5 Newcomb (1960): system is unstable, if and only if the solution of Eq. (48) vanishes two or more points. Refs. Connor et al. (1978) and Lortz & Nührenberg (1978). have obtained the stability boundaries for ballooning modes.
In the intermediate frequency regime the first term in Eq. (72) can be neglected and therefore one obtains Inserting Eq. (75) into Eq. (71) yields Tang et al. (1980) The sound wave coupling term (3rd term) results in the so-called second harmonic TAE in the circular cross section case Zheng et al. (1999).

Toroidal Alfvén eigen modes
In this subsection we review TAE theory. In the last two subsections we see that interchange and ballooning modes are characterized by having only single dominant or resonant mode at resonance surfaces. In particular their resonance surfaces locates at mode rational surface where m − nq = 0. TAEs are different from them. TAEs involve two mode coupling. In particular, the first TAEs are centered at the surface where q =( m 0 + 1/2)/n.T w o neighboring Fourier modes (m 0 and m 0 + 1) propagate roughly with same speed v A /2Rq but in opposite directions. They form a standing wave. The toroidal geometry can induce the first frequency gap so that the standing wave becomes an eigen mode, i.e., TAEs Cheng et al. (1985) Rosenbluth et al. (1992. In the second TAE case, although they have same mode resonance surfaces as interchange and ballooning modes, m 0 ± 1 mode coupling is involved to form standing 2nd TAEs. The frequency gap for second TAEs in circular cross section case is due to plasma compressibility effect Zheng & Chen (1998).
To explain two mode coupling picture, we consider tokamak model equilibrium with circular cross section, low beta, and large aspect ratio (i.e., 1/ǫ = R/a ≫ 1). There is a review paper on TAEs Vlad et al. (1999). Here, we describe the local dispersion relation for even and odd modes and explain the 2nd TAEs together with the 1st TAEs. The magnetic field in this model can simply be expressed as B = B φ (r)/(1 + ǫ cos θ)e φ + B θ (r)e θ . The general case will be addressed in Sec. 5 with AEGIS code formalism. Since their frequency is much larger than shear Alfvén mode frequency, the compressional Alfvén modes are decoupled. Therefore, we can use Eqs. (54) and (55) as starting equations for TAE investigation. Noting that Alfvén frequency is much larger than sound wave frequency, the first term in Eq. (55) can be dropped. Adopting the Fourier decomposition in Eq. (56), the sound wave equation (55) becomes Using this solution for ∇·ξ, Eq. (54) can be reduced to Zheng et al. (1999) where E m = ϕ m /r, v 2 A = B 2 0 /ρ m , B 0 denotes magnetic field at magnetic axis, and w represents the rest magnetic well terms.
We first examine singular layer physics. In this layer only terms contains second order derivative in r need to be taken into consideration. From the first six terms in Eq. (77) one can see that the 2nd TAEs (coupling of E m−1 and E m+1 ) have structure similarity to the 1st TAEs (coupling of E m and E m+1 ). The 1st TAE coupling is due to finite value of aspect ratio; while the 2nd TAE coupling is due to finite beta value. For brevity we focus ourselves to discuss the 1st TAE case. Denoting ω 0 = ω A /2, q 0 =( m + 1/2)/n, δω = ω − ω 0 , and δq = q − q 0 , the singular layer equations describing the coupling of m and m + 1 modes becomes Introducing even and odd modes: δφ ± = δφ m ± δφ m+1 , these two equations become ∂ ∂δq Integrating once one obtains where D is 2 × 2 matrix and A ± are integration constants. Integration of Eq. (78) across singular layer (i.e., from δq − to δq + ) one obtains the dispersion relation Here, D ij are D matrix elements and two parameters ∆ ± ≡ δφ ± | δq + δq − /A ± are determined by the outer solutions to the left and right of singular layer. As soon as ∆ ± are computed from outer regions, Eq. (79) can be used to determine the frequency. In general this frequency can be complex.
One can obtain the gap width ∆ω = δω + − δω − = ǫω 0 1 − 1/(2m + 1) 2 . The 1st TAEs are Alfvén eigen modes with frequency inside this gap. They are marginally stable and tend to be excited by resonances with energetic particles. Note that the gap width is proportional to ǫ. In cylinder limit the gap vanishes. Therefore, existence of TAEs is due to toroidal effects. Also, we note that the dispersion relation, Eq. (79), allows two types of TAEs: even and odd types (ϕ ± ), depending on the values of ∆ ± .
We have discussed the 1st TEA theory through coupling of neighboring modes. In similar way one can also develop the 2nd TAE theory through coupling of m ± 1 modes Zheng & Chen (1998). If FLR effects are taken into consideration, the Alfvén types of singularities can be resolved, so that discrete modes can emerge in the continuum. This types of modes are referred to as kinetic TAEs (i.e., KTAEs). Due to correction of gyrokinetic theory Zheng et al. (2007), several missing FLR effects are recovered. Consequently, KTAE theories by far need to be reevaluated. In this subsection we describe the kinetically driven modes (KDMs), such as KBMs, EPMs, etc. The frequencies of these modes usually reside in continuum spectrum. Therefore, they are generally damped without driving effects. Unlike KTAEs, for which FLR effects are taken into account to resolve singularity, for KDMs strong kinetic effects from wave-particle-resonances are included to overcome continuum damping. That is why they are referred to as kinetically driven modes. Energetic particle drives to marginal stable TAEs can instantly lead unstable TAEs, but the drives to KDMs need to accumulate sufficient energy to overcome continuum damping for unstable KDMs to develop Tsai & Chen (1993) Zheng et al. (2000. In Secs. 4.3 and 4.4 one has seen that there are two types of modes: ballooning and TAEs. Therefore, KDMs also have two types. Those related to ballooning modes are referred to as KBMs, while EPMs are related to TAEs and usually driven by wave-energetic-particle resonances. We employ ballooning representation formalism to discuss them. We start with the ballooning mode equation in intermediate frequency regime, Eq. (76), with energetic particle effects included. Introducing the transformation ζ = ϕp 1/2 , Eq. (76) becomes where p = 1 + Λ 2 s , g = cos θ + Λ sin θ, δg h is perturbed distribution functions for hot ions, κ r and κ θ are respectively radial and poloidal components of magnetic field line curvature κ, Ω = ω/ω A , ω d is magnetic drift frequency, v is particle speed, the subscripts ⊥ and represent respectively perpendicular and parallel components to the equilibrium magnetic field line, ε = v 2 /2 is particle energy, μ = v 2 ⊥ /2B is magnetic moment, and F 0h is equilibrium distribution function for hot ions. For simplicity we have neglected the finite Larmor radius effects and only take into account the kinetic effects from energetic ions.
To study KDMs one need to investigate singular layer behavior. In ballooning representation space, singular layer corresponds to θ → ∞ limit. Again, we exclude the 2nd TAE from discussion (i.e., assuming Γ = 0). Equation (81) in θ → ∞ limit becomes: This is the well-known Mathieu equation. According to Floquet's Theorem, its solution takes following form where P(θ + 2π)=P(θ). Since modes with longer parallel-to-B wavelengths tend to be more unstable, we shall examine solutions corresponding to the two lowest periodicities. The first one is related to KBMs Tsai & Chen (1993) and the second one is related to EPMs Zheng et al. (2000).

Overview of Magnetohydrodynamics Theory in Toroidal Plasma Confinement 21
We first discuss KBMs. The KBM-type solution is given by Inserting Eq. (84) into Eq. (83), one obtains, noting ǫ ≪ 1, Therefore, at leading order, one has where ℑm{Ω} > 0 for causality. Note here that Eq. (85) is valid for general Ω, so that frequency at continuum is allowed, as soon as causality condition is satisfied.
Inserting Eq. (86) into Eq. (83), one obtains, for Ω 2 − 1/4 ∼O(ǫ) and ǫ ≪ 1, and Ω 2 +,− = 1/4 ± ǫΩ 2 . The leading order solution can, therefore, be expressed as The causality condition is ℑm (Ω 2 − Ω 2 + )(Ω 2 − Ω 2 − ) 1/2 < 0. Equation (88) can describe both TAEs and KDMs of TAE type (e.g., EPMs). Existence of TAE solution requires mode frequency to fall in the gap: Ω − < Ω < Ω + , as shown by the TAE theory in configuration space in Sec. 4.4. For KDMs mode frequency can be in the continuum, i.e., outside the gap as soon as causality condition is satisfied. For TAEs ζ T contains an O(1) back scattering and, hence, the continuum damping is either suppressed or much reduced. On the other hand, for KDMs ζ K contains no back scattering from periodic potential in Eq. (83), and, consequently, there is significant mount of continuum damping. Note that in principle both types of solutions can co-exist at |Ω|≈1/2 . However, the TAE solution tends to be more unstable in this case than KDMs, since its continuum damping is much less or absent while the instability drives are generally comparable.
With outer solutions given by Eq. (85) or Eq. (88), one can obtain the corresponding dispersion relation by matching outer and inner solutions. For KBMs Eq. (81) can be used to construct the following quadratic form in inner region: where Here, the superscript * represents complex conjugate. Matching inner (Eq. (89)) and outer (Eq. (85)) solutions one obtains the dispersion relation Tsai & Chen (1993) − iΩ + δW f + δW K = 0.
Here, we note that kinetic effects from core plasma should also be taken into account in outer region. As proved in Ref. Zheng & Tessarotto (1994a) this results in the so-called apparent mass effect and leads Ω in the first term of Eq. (90) to become complicated function of actual mode frequency.
Similarly, for KDMs of TAE type, for example EPMs, one need to consider even and odd modes. For even modes the dispersion relation is given by Refs. Zheng et al. (2000) and Tsai & Chen (1993) − where δT f represents MHD fluid contribution and δT K is energetic-particle contribution to the quadratic form in inner region.
The dispersion relations in Eqs. (90) and (91) extend respectively MHD ballooning modes in diamagnetic gap and TAEs in Alfvén gap to respective continua. Kinetic drives are the causes to make causality conditions satisfied.

Global numerical analyses of MHD modes: AEGIS code formalism
In Sec. 2.3 analytical or semi-analytical theories are presented to describe four types of MHD modes in toroidal geometry. Due to the developments of modern numerical method and computer hardware, conventional asymptotic expansion methods for global modes have become outdated and been substituted by direct numerical computation. Several excellent numerical codes have been developed in the past to study linear MHD stability of toroidally confined plasmas, such as such as PEST Grimm et al. (1976, GATO Bernard et al. (1981), DCON Glasser (1997), AEGIS Zheng & Kotschenreuther (2006), etc. In this section we focus on description of AEGIS code, in view of that AEGIS is an adaptive MHD shooting code capable to study MHD continuum Zheng et al. (2005). Through describing AEGIS formalism, we can further explain the general features of MHD eigen modes in toroidally confined plasmas.
Let us first describe the toroidal system to be investigated. The core part is plasma torus, which is surround by a resistive wall; Between plasma torus and resistive wall there is inter vacuum region and outside the resistive wall there is outer vacuum region, which extends to infinity. For simplicity, it is assumed that the wall is thin. We denote respectively the interfaces between plasma torus and inner vacuum region, inner vacuum region and wall, and wall and outer vacuum region as ψ a , ψ b− , and ψ b+ .

22
Topics in Magnetohydrodynamics www.intechopen.com Overview of Magnetohydrodynamics Theory in Toroidal Plasma Confinement 23

MHD equations and numerical solution method for plasma region
In this subsection we describe how to reduce MHD equations for global mode analyses. The starting equation is the single fluid MHD equation (14). This is a vector equation and can be projected onto three directions to get scalar equations. The parallel projection has be derived in Eq. (33). From parallel equation one can solve for ∇·ξ, which is the only unknown needed for two-perpendicular equations to become a complete set of equations. In principle the parallel motion can not be described by MHD model, since particles are not localized along magnetic field line. There are wave-partcle resonance, trapped particle, and so-called small parallel particle speed effects, etc. Nevertheless, from analyses in Sec. 4.3 one can see that in low frequency regime the parallel coupling results only in the so-called apparent mass effect, while in intermediate regime the parallel coupling mainly gives rise to the 2nd TAEs. Note that apparent mass effect can be absorbed by rescaling mode frequency and inclusion of the 2nd TAE effect is straightforward as discussed in Sec. 4.4. For brevity we limit ourselves to treat only two perpendicular components of Eq. (14) with Γ set to zero. AEGIS-K code has been developed to include parallel dynamics using kinetic description .
Using general flux coordnates in Eq. (22), the magnetic field line displacement is decomposed as follows Since we deal with linear problem, the Fourier transform can be used to decompose perturbed quantities in poloidal and toroidal directions, with ξ m = π −π dθξ exp{−imθ}/ √ 2π. With the toroidal symmetry assumed, only a single toroidal Fourier component needs to be considered. As usual, equilibrium quantities can be decomposed as matrices in poloidal Fourier space, for example In the poloidal Fourier decomposition, the Fourier components are cut off both from lower and upper sides respectively by m min and m max . Therefore, the total Fourier component under consideration is M = m max − m min + 1. We also use bold face (or alternatively [[···]]) to represent Fourier space vectors, and calligraphic capital letters (or alternatively ··· )t o represent the corresponding equilibrium matrices (e.g., J for J) in poloidal Fourier space.
To get scalar equations, we project Eq. (14) respectively onto two directions J 2 ∇θ ×∇ ζ · B × [··· ×B]/B 2 and (1/qχ ′ )J 2 ∇ζ ×∇ ψ · B × [··· ×B]/B 2 , and then introduce the Fourier transformation in Eq. (93) to two projected equations. These procedures lead to the following set of differential equations in matrices Here, the equilibrium matrices contain two parts: plasma/field force and inertia contributions, e.g., denotes the dimensionless growth rate normalized by the Alfvén frequency at magnetic axis, ρ N is the dimensionless mass density normalized by the mass density at magnetic axis, subscript "0" refers to quantities at magnetic axis, and We can reduce the set of equations (94) and (95) into a set of first order differential equations as in the DCON formalism Glasser (1997). By solving Eq. (95), one obtains Inserting this solution into Eq. (94), we get where F = D−B † A −1 B, K = E−B † A −1 C, and G = H−C † A −1 C. These matrices can be further simplified as follows Glasser (1997) Introducing the expanded 2M unknowns u = ξ u 2 , where u 2 = F ξ ′ + Kξ, Eq. (96) is reduced to the set of 2M first order equations where 2M × 2M matrix We note that ξ and u 2 in plasma region are related to the magnetic field and pressure as follows The set of eigen mode equations in Eq. (99) can be solved numerically by independent solution method together with multiple region matching technique as described in Ref. Zheng & Kotschenreuther (2006). With M boundary conditions imposed at magnetic axis, there remain only M independent solutions: where the superscripts are used to label independent solutions. We use the cylinder limit to describe the boundary condition at magnetic axis, i.e., ξ ψ,m ∝ r m . The general solution can be then obtained as a combination of the M independent solutions,

The solution of vacuum region
For completeness, in this subsection we briefly review the vacuum solutions in Ref. Zheng & Kotschenreuther (2006). The vacuum regions are described by the Laplace equation where u is the magnetic scalar potential and is related to the perturbed magnetic field by δB = −∇u. Here, we note that this representation of vacuum magnetic field, although being simple, excludes the consideration of n = 0 modes. To study n = 0 modes, one more scalar is needed to represent the vacuum magnetic field. For the sake of conciseness, we outline the general solutions for inner and outer vacuum regions simultaneously.
As in the plasma region, Fourier decompositions are introduced for both poloidal and toroidal directions to solve Eq. (103). Then Eq. (103) becomes a set of second-order differential equations of number M for u. This set of second-order differential equations can be transformed into a set of first-order differential equations of number 2M, by introducing an additional field v = −[[J∇ψ · δB]], which is related to the magnetic scalar potential in Fourier space as follows: v = J|∇ψ| 2 ∂u ∂ψ + iJ∇ψ ·∇θ Mu.
There are 2M independent solutions for Eq. (103), which can be used to construct the following independent solution matrices: The general solutions in the vacuum regions can be expressed as a linear combination of the independent solutions: where c v and d v are constant vectors in the independent solution space. To distinguish the inner and outer vacuum solutions, we let c v1 and d v1 denote the constants for inner vacuum region and c v2 and d v2 for outer vacuum region.
In the outer vacuum region, the scalar potential u is subjected to M boundary conditions at infinite ψ. With these M boundary conditions imposed, there are only M independent solutions left. Without loss of generality, we can set c v2 to be zero in this case. Consequently, eliminating d v2 in Eq. (104), we obtain where the M × M matrix T is given by T = U 2 V −1 2 | ψ b+ . The matrix T can be computed by means of the Green function method Chance (1997).

Overview of Magnetohydrodynamics Theory in Toroidal Plasma Confinement 27
In the inner vacuum region, the independent solutions can be constructed, for example, with the use of an inward numerical shooting Zheng & Kotschenreuther (2006), with the following boundary conditions imposed at ψ b− : where O is M × M zero matrix. Since the boundary conditions in Eq. (105) give δB ·∇ψ = 0 at wall, these conditions correspond to a set of solutions that corresponds to the perfectly conducting wall type. On the other hand, since the boundary conditions in Eq. (106) guarantee that the independent solutions to be continuous with outer vacuum solutions, these conditions correspond to a set of solutions that corresponds to the no-wall type. Using the general expression for the solutions in Eq. (104), we can express the normal and parallel magnetic fields at the plasma-vacuum interface as follows:

Eigenvalue problem
The solutions in the plasma and vacuum regions described in the last two subsections can be used to construct the eigen value problem Zheng & Kotschenreuther (2006 (107) and (108)] solutions at the interface ψ a gives where Note that δW ∞ and δW b correspond to the energy matrices without a wall and with a perfectly conducting wall at ψ b , respectively, as can be seen from the boundary conditions in Eqs. (105) and (106).
We now consider the matching across the thin resistive wall. For the radial magnetic field, the Maxwell equation ∇·δB = 0 and the thin wall assumption lead to The current in the resistive wall causes a jump in the scalar magnetic potential. This can be obtained from the Ampére law where σ is the wall conductivity. Equation (110) can be reduced to where τ w = μ 0 σdb/τ A , d is the wall thickness, b is the average wall minor radius, and V = M J|∇ψ||∇θ|−J|∇ψ ·∇θ| 2 /(|∇ψ||∇θ|) M +n 2 J|∇φ| 2 |∇ψ|/|∇θ| .
Since c v2 = 0, we find that Eqs. (104) -(106) yield From Eqs. (109), (111), and (112) we find the eigen mode equations The dispersion relation for this eigen value problem is given by the determinant equation det |D 0 (γ N )| = 0. In general the Nyquist diagram can be used to determine the roots of this dispersion relation. For RWMs, however, the growth rate is much smaller than the Alfvén frequency. Therefore, the growth rate dependence of δW ∞ δW −1 b can be neglected for determining the stability condition. Consequently, one can use the reduced eigen value problem with the RWM mode growth rate γ N on the right hand side of this equation used as the eigen value to determine the stability.

Discussion
Now let us discuss the connection of current global theory with localized analytical theories described in Sec. 4. The singular layer equation in Eq. (48)

Summary and discussion
In this chapter we have given an overview of MHD theory in toroidal confinement of fusion plasmas. Four types of fundamental MHD modes in toroidal geometry: interchange, ballooning, TAEs, and KDMs, are discussed. In describing these modes we detail some fundamental analytical treatments of MHD modes in toroidal geometry, such as the average technique for singular layer modes, ballooning representation, mode coupling treatment in TAEs/KDMs theories. Note that analytical approach is often limited for toroidal plasma physics. Global numerical treatment of MHD modes is also reviewed in this chapter, especially the AEGIS code formalism. These theories are reviewed in ideal MHD framework.
Here, we briefly discuss kinetic and resistive modifications to ideal MHD, as well as the connection of MHD instabilities to transport.
Let us first discuss kinetic effects. Since strong magnetic field is used to contain plasmas in magnetically confined fusion experiments, MHD theory can be rather good to describe fusion plasmas in the direction perpendicular to magnetic field. This is because strong magnetic field can hold plasmas together in perpendicular movement. Therefore, MHD is a very good model to describe perpendicular physics, if FLR effects are insignificant. However, in parallel direction the Lorentz force vanishes and particle collisions are insufficient to keep particles to move collectively. Consequently, kinetic description in parallel direction is generally necessary. Kinetic effect is especially important when wave-particle resonance effect prevails in the comparable frequency regime ω ∼ ω si Zheng & Tessarotto (1994a). In the low frequency regime ω ≪ ω si , wave-particle resonances can be so small that kinetic description results only in an enhancement of apparent mass effect. Kinetic effect in this case can be included by introducing enhanced apparent mass. Another non-resonance case is the intermediate frequency regime ω si ≪ ω ≪ ω se . In this regime kinetic description results in a modification of ratio of special heats. By introducing proper Γ MHD can still be a good approximation. Recovery of perpendicular MHD from gyrokinetics has been studied in details in Ref. Zheng et al. (2007). Perturbed current at rational surface due to interchange modes leads to field line reconnection and formation of magnetic islands.

29
Overview of Magnetohydrodynamics Theory in Toroidal Plasma Confinement www.intechopen.com Next, let us discuss resistivity effects. Resistivity usually is small in magnetically confined fusion plasmas. Due to its smallness resistivity effects are only important in the singular layer region. With ideal MHD singular layer theory detailed in Sec. 4.2 one can rederive resistive singular layer equations given in Ref. Glasser et al. (1975). However, it should be pointed out that, when kinetic enhancement of apparent mass effect is taken into account, the ratio of resistivity and inertia layer widths changes. This leads kinetic description of resistive MHD modes to become substantially different from fluid description Zheng & Tessarotto (1996) Zheng & Tessarotto (1995. Kinetic analysis of low frequency resistive MHD modes becomes necessary. The driving force for ideal MHD instabilities is related to pressure gradient. Resistivity can instead cause field line reconnection and induce the so-called tearing modes. It is important to note that if current gradient is taken into account pressure driven modes and tearing modes are coupled to each other. The underlying driving mechanism for pressure driven modes is the release of plasma thermal energy from the interchange of magnetic flux tubes. Actually, interchange-type modes exchange not only thermal and magnetic energies between flux tubes, but also current. In a plasma with a current (or resistivity) gradient, such an interchange can create a current sheet at a mode resonance surface and result in the excitation of current interchange tearing modes (CITMs) as shown in Fig. 2 Zheng & Furukawa (2010).
Instabilities of interchange type have been widely used to explain anomalous transport in tokamaks in terms of the formation of turbulent eddies through nonlinear coupling. However, the explanation for experimental observations that the electron energy transport is much larger than what one would expect from diffusive process due to Coulomb collisions is still unsatisfactory. The electron Larmor radius is much smaller than ion one. Nonetheless, the electron thermal transport often is stronger than ion transport. In Ref. Rechester & Rosenbluth (1978), the broken magnetic surfaces due to formation of magnetic island and stochastic field lines are used to explain the enhanced electron transport. But, how magnetic islands are formed in axisymmetric tokamak plasmas has not been given. CITM theory shows that interchange-type instabilities can directly convert to current interchange tearing modes. This helps to clarify the source of electron transport in tokamaks.
Another transport issue we need to discuss is the so-called flow shear de-correlation of turbulences. This concept has been widely used for explaining suppression of plasma turbulences. In fact, this picture is not right for systems with magnetic shear. We use Fig.  3 to explain it (L. J. Zheng and M. Tessarotto, private communication). In Fig. 3, the dashed long arrow represents a magnetic field line on a given magnetic surface ψ 0 , and two solid long arrows denote the magnetic field lines respectively at two time sequences t 0 and t 0 + ∆t on an adjacent magnetic surface ψ 1 . Let us examine the correlation pattern in the local frame moving together with equilibrium velocity of the dashed long arrow on surface ψ 0 . The modes are supposed to locate around the point "O" initially at t = t 0 . After a time interval ∆t, the field line on surface ψ 1 moves relatively to the dashed long arrow on the surface ψ 0 due to flow shear. From Fig. 3 one can see that the fixed pattern has not been de-correlated by flow shear, instead the pattern just propagates from point "O" at time t = t 0 to point "O ′ " at subsequent time t = t 0 + ∆t. This indicates that flow shear does not de-correlate turbulence eddies. Only flow curvature can result in the de-correlation. This resembles to ballooning mode behavior in rotating plasmas with Cooper representation Waelbroeck & Chen (1991).  In conclusion significant progresses have been made for linear ideal MHD theories and numerical codes in the past dacades. However, the kinetic effects on MHD remains considerably open. Although correction of gyrokinetics theory has been made recently Zheng et al. (2007), the applications of the new gyrokinetics theory remain to be worked out. The theories for FLR effects on ballooning modes, KTAEs, energetic particle effects, etc. need to be modified with newly corrected gyrokinetics theory. The extension of toroidal resistive MHD theory Glasser et al. (1975) to take into account the small parallel ion speed effect Zheng & Tessarotto (1996) and current interchange effects Zheng & Furukawa (2010) is under consideration. To understand plasma physics intuitively one need to master the MHD behaviors. As sciences advance, gap between published textbooks and cutting-edge researches gradually develops. Connection from textbook knowledge to up-to-dated research results can often be tough. Review articles can help. This book contains eight topical review papers on MHD. For magnetically confined fusion one can find toroidal MHD theory for tokamaks, magnetic relaxation process in spheromaks, and the formation and stability of field-reversed configuration. In space plasma physics one can get solar spicules and X-ray jets physics, as well as general sub-fluid theory. For numerical methods one can find the implicit numerical methods for resistive MHD and the boundary control formalism. For low temperature plasma physics one can read theory for Newtonian and non-Newtonian fluids etc.