Hydromagnetic Waves in Weakly Ionised Media. I. Basic Theory, and Application to Interstellar Molecular Clouds

We present a comprehensive study of MHD waves and instabilities in a weakly ionised system, e.g., an interstellar molecular cloud. We determine all the critical wavelengths of perturbations across which the sustainable wave modes can change radically (and so can their decay rates), and various instabilities are present or absent. These critical wavelengths are essential for understanding the effects of MHD waves (or turbulence) on the structure and evolution of molecular clouds. Depending on the angle of propagation relative to the magnetic field and the physical parameters of a cloud, there are wavelength ranges in which no wave can be sustained as such. Yet, for other directions of propagation or different properties of a model cloud, there may exist some wave mode(s) at all wavelengths. For a typical cloud, magnetically-driven ambipolar diffusion leads to removal of any support against gravity that most short-wavelength waves (or turbulence) may have had, and gravitationally-driven ambipolar diffusion sets in and leads to cloud fragmentation into stellar-size masses, as first suggested by Mouschovias three decades ago -- a single-stage fragmentation theory of star formation, distinct from the hierarchical fragmentation picture. Phase velocities, decay times, and eigenvectors (the densities and velocities of neutral particles and the plasma, and components of the magnetic field) are determined as functions of the wavelength of the disturbances and are explained physically. Comparison of the results with nonlinear analytical or numerical calculations is also presented, excellent agreement is found, and confidence in the analytical, linear approach is gained to explore phenomena difficult to study through numerical simulations. Mode splitting and mode merging, which are impossible in single-fluid systems for linear perturbations, occur naturally in multifluid systems.


INTRODUCTION -BACKGROUND
A typical molecular cloud which has not yet given birth to stars is a cold (T ≃ 10 K) but complex, partially ionised system, in which self-gravitational and magnetic forces are of comparable magnitude, with thermal-pressure forces becoming important at high densities ( > ∼ 3 × 10 8 cm −3 ) or along magnetic field lines.  showed that, barring external disturbances, if the magnetic field were to be frozen in the matter, interstellar clouds that have not yet given birth to stars would remain in magnetohydrostatic (MHS) equilibrium states. However, ambipolar diffusion (the relative motion of neutral particles and charged particles attached to magnetic field lines) is an unavoidable process in partially ionised media. It reveals itself in two distinct ways, depending on whether it is magnetically or gravitationally driven (see discussion in § 4). The two kinds of ambipolar diffusion acting together initiate fragmentation and star formation in molecular clouds (Mouschovias 1987a).
In this fragmentation theory, the evolutionary (or fragmentation, or core formation) timescale is the gravitationallydriven ambipolar-diffusion timescale, τ AD . This does not mean, however, that it takes a time equal to τ AD to form stars. The star-formation timescale can be a fraction or a multiple of τ AD , depending on the mass-to-flux ratio of the parent cloud and the degree to which hydromagnetic waves (HM) contribute to the support of the cloud: the closer to its critical value the mass-to-flux ratio is and/or the greater the contribution of HM waves to cloud support, the faster the evolution and the shorter the star-formation timescale (e.g., see Mouschovias 1987a;Fiedler & Mouschovias 1993, Fig. 9a;Ciolek & Basu 2001;Tassis & Mouschovias 2004, Fig. 4). 1 The collapse retardation factor, ν ff = τ ff /τni (where τ ff is the free-fall time and τni the neutral-ion collision time), is the factor by which magnetic forces slow down the contraction relative to free-fall. This was a new theory of fragmentation (or core formation), initiated by the decay, due to ambipolar diffusion, of relatively small-wavelength perturbations (Mouschovias 1987a(Mouschovias , 1991a. Magnetic braking operates on a timescale shorter than the ambipolar-diffusion timescale and even the free-fall timescale, and keeps a cloud (or fragment) essentially corotating with the background up to densities ≃ 10 4 − 10 6 cm −3 and thus resolves the angular momentum problem of star formation. More specifically, the entire range of periods of binary stars from 10 hr to 100 yr was shown to be accounted for by this self-initiated mode of star formation (Mouschovias 1977). Even single stars and planetary systems become dynamically possible (Mouschovias 1978(Mouschovias , 1983).
Star formation, whether self-initiated or triggerred (e.g., by a spiral density wave, or the expansion of an Hii region or a supernova remnant; see review by Woodward 1978) is an inherently nonlinear process. Ambipolar-diffusion−initiated star formation has been studied analytically (Mouschovias 1979(Mouschovias , 1991a and numerically using adaptive grid techniques in axisymmetric geometry up to densities ∼ 10 10 cm −3 , by which isothermality begins to break down (Fiedler & Mouschovias 1992Ciolek & Mouschovias 1993Basu & Mouschovias 1994, 1995a. More recently, these calculations were extended into the opaque phases of star formation (Desch & Mouschovias 2001;Tassis & Mouschovias 2007a, b, c;Kunz & Mouschovias 2009, 2010. The key conclusions of the earlier analytical calculations have been verified and numerous new, specific, quantitative predictions have been made, many of which have been confirmed by observations (e.g., see Crutcher et al. 1994;Ciolek & Basu 2000;Chiang et al. 2008; reviews by Mouschovias , 1996. One may nevertheless take a step back from the relatively complicated numerical calculations and ask which results, if any, of the nonlinear simulations can be recovered with a linear analysis. With one's confidence increased in the validity of the linear approach (within some selfevident limits), one may then make new predictions concerning phenomena that cannot be or have not been included yet in the nonlinear calculations.
The propagation, dissipation, and growth of perturbations in a physical system depends both on the nature of the perturbations and the properties of the system. Hydromagnetic waves (or MHD turbulence) seem to play a significant role in molecular clouds on lengthscales typically greater than ∼ 0.1 pc. They have been shown to account quantitatively for the observed supersonic but subAlfvénic spectral linewidths: an observational almost-scatter diagram of linewidth versus size is converted into an almost perfect straight line if plotted in accordance with a theoretical prediction by Mouschovias (1987a), which relates the linewidth, the size, and the magnetic field strength of an observed object (see Mouschovias & Psaltis 1995, Figs. 1 and 2, and update by Mouschovias et al. 2006, Figs. 1 and 2).
Using a linear analysis as a first step in understanding nonlinear phenomena is not, of course, a new idea. Jeans (1928) used it to obtain his famous instability criterion for the collapse of a cloud against thermal-pressure forces. Hardly any astrophysical system exists whose stability with respect to small-amplitude disturbances has not yet been studied by using at least an idealized, mathematically tractable model of the physical system. A magnetically supported molecular cloud, however, defies a simple linear analysis. First, no realistic equilibrium states have been obtained by analytical means. Second, to study the role of ambipolar diffusion in star formation, one must use at least the two-fluid magnetohydrodynamic (MHD) equations governing the motions of the neutral particles and the plasma (ions and electrons). In fact, as shown by Ciolek & Mouschovias (1993), charged (and even neutral) grains play a very significant role in the ambipolar-diffusion−initiated protostar formation. One then has to use at least the four-fluid (neutral molecules, plasma, negatively-charged and neutral grains) MHD equations even for a linear analysis to be realistic and relevant to typical molecular clouds. In this paper we use the two-fluid MHD equations to study the propagation, dissipation, and growth of HM waves in an idealized model molecular cloud. In a subsequent paper we consider the effects of the grain fluid(s). Langer (1978) studied the stability of a model molecular cloud (infinite in extent and uniform in density and magnetic field) with respect to small-amplitude, adiabatic perturbations in the presence of ambipolar diffusion. For propagation along the magnetic field lines, he recovered, as one would expect intuitively, the Jeans dispersion relation and instability in the absence of the magnetic field. He then investigated the wave propagation perpendicular to the field lines. He showed that the Jeans instability is still present, that the critical wavenumber for instability is independent of the magnetic field strength, but that the growth rate depends on both the field strength and the degree of ionisation. Aside from two spurious curves in his Figure 1, which exhibits the growth rate and decay time of some modes, our results for propagation of the low-frequency modes perpendicular to the magnetic field are in agreement with Langer's -he ignored the high-frequency ion modes. Yet even in this, previously studied case, we offer new analytical expressions and new physical insight and interpretation of the results. Moreover, we present not only the eigenvalues (frequencies or, equivalently, phase velocities) as functions of wavelength but the eigenvectors as well (i.e., material velocities, densities, magnetic-field components, etc.). We also study propagation at arbitrary angles with respect to the magnetic field, and we offer a thorough discussion of the wave modes, not just the ambipolar-diffusion-induced instability. Pudritz (1990) revisited Langer's problem (with the minor difference of considering isothermal perturbations) but introduced a new effect: he assumed that there exists a power-law spectrum of small-amplitude waves, and then he studied the effect that this spectrum has on the ambipolar-diffusion-induced, Jeans-like instability. 2 He concluded that the slope of the spectrum (considered as a function of wavelength) has an important effect on the growth rate of the instability; the steeper the spectrum, the greater the growth rate. (The growth rate of gravitationally-driven ambipolar diffusion, however, cannot possibly exceed the free-fall rate.) Several other papers have appeared in print since 1990, studying different aspects of weakly ionised systems, focusing usually on the stability of certain MHD modes or shocks, especially as it may relate to the formation of structures in molecular clouds and/or on the effect of the grain fluid(s) on the allowable wave modes or shocks (e.g., Wardle 1990;Balsara 1996;Zweibel 1998;Kamaya & Nishi 1998Mamun & Shukla 2001;Cramer et al. 2001;Tytarenko et al. 2002;Zweibel 2002;Ciolek et al. 2004;Lim et al. 2005;Oishi & Mac Low 2006;Roberge & Ciolek 2007;van Loo et al. 2008;Li & Houde 2008). In this paper we present a general theory of the propagation, dissipation and growth of MHD waves in partially ionised media in three dimensions, with emphasis on mathematical transparency of the formulation and analytical solution of the problem, the physical understanding and interpretation of all modes, including their eigenvectors, the many critical wavelengths that exist and which separate regimes dominated by different waves or instabilities, and on specific features relevant to the evolution of molecular clouds. As mentioned above, even when a particular result agrees with previous work, we offer new insight into its physical understanding.
In § 2 we present the equations governing the behaviour of a weakly ionised, magnetic, self-gravitating interstellar cloud. The equations are linearised, Fourier-analyzed, and put in dimensionless form. The free parameters of the problem are identified, their physical meaning explained, and their typical values given. The different hydromagnetic modes and their dependence on wavelength for different directions of propagation relative to the magnetic field are calculated and explained physically in § 3. Analytical expressions for the phase velocities, damping timescales, growth timescales, including critical or cutoff wavelengths, are also obtained. A physical discussion of the eigenvectors is an integral part of this presentation. Section 4 summarizes some of the results and their relevance to the formation of protostellar fragments (or cores) and to other observable phenomena. It also gives in two Tables all the critical wavelengths and the ranges of wavelengths in which different modes can exist in molecular clouds, for propagation parallel, perpendicular, and at arbitrary angles with respect to the magnetic field.

Basic Equations
We consider a weakly ionised medium (e.g., an interstellar molecular cloud) consisting of neutral particles (H2 with a 20% helium abundance by number; subscript n), electrons, and singly-charged positive ions (subscript i). For specificity we assume that the ions are molecular ions (such as HCO + ); for the densities of interest in this paper (∼ 10 3 cm −3 ), this is sufficient since atomic ions (such as Na + or Mg + ) are less abundant (a result of depletion of metals in dense clouds) and, in any case, they have masses comparable to that of HCO + (for more detailed treatments of the chemistry, see Ciolek & Mouschovias 1995, 1998, or the appendix of Mouschovias & Ciolek 1999). Interstellar grains, which have been shown to have significant effects on the formation and contraction of protostellar cores (Ciolek & Mouschovias 1993 and in the opaque phase of star formation (Tassis & Mouschovias 2007a, b, c;Kunz & Mouschovias 2009, 2010 are neglected in this analysis; they are accounted for in a subsequent paper. The magnetohydrodynamic (MHD) equations governing the evolution of the above two-fluid system are where ρα and vα are, respectively, the density and velocity of species α, is the time-derivative comoving with species α, ψ the gravitational potential, Pn the neutral pressure, B the magnetic field, T the temperature, Γn and Λn the heating and cooling rates (per unit volume) of the neutral gas, and the neutral-ion and ion-neutral mean collision (i.e., momentum-exchange) times. The quantity G is the universal gravitational constant, and kB is Boltzmann's constant; µ is the mean mass of a neutral particle in units of the atomic-hydrogen mass and is equal to 2.33 for a H2 gas with a 20% helium abundance by number. The quantities ζ CR and α dr in the ion mass continuity equation (1b) are, respectively, the cosmic-ray ionisation rate and the coefficient for dissociative recombination of molecular ions and electrons (in cm 3 s −1 ). In writing equation (1b), we have used the condition of local charge neutrality e(ni − ne) = 0, where ni and ne are, respectively, the number densities of ions and electrons. We may use the assumption of local charge neutrality because the various HM modes of interest here have frequencies much smaller than the electron plasma frequency ωp,e = (4πnee 2 /me) 1/2 = (4πxenne 2 /me) 1/2 = 5.64 × 10 2 (xe/10 −7 ) 1/2 (nn/10 3 cm −3 ) 1/2 s −1 (where xe = ne/nn is the abundance of electrons relative to the neutrals, and is equal to the degree of ionisation for weakly-ionised systems); hence, any excess charge density is quickly shielded by the mobile electrons, so that ne = ni for timescales > ∼ ω −1 p,e . Since we neglect the effects of grains in this paper, capture of ions onto grains is not included on the right-hand side of equation (1b) as a sink term for ions.
Heat conduction and viscosity are not important for the densities and lengthscales of interest (nn ≃ 10 2 − 10 5 cm −3 , L ≃ 10 −3 − 10 1 pc), and are therefore ignored as possible sources of heating/cooling in the model clouds. ( Because ρn ≫ ρe, ρi, we include only the neutral density as a source term in Poisson's equation (1h). Similarly, the gravitational and thermal-pressure forces (per unit volume) on the plasma (ions and electrons) have been neglected in the plasma force equation (1d). One can easily show that, for the physical conditions in typical molecular clouds, they are completely negligible in comparison to the magnetic force exerted on the plasma, except in a direction almost exactly parallel to the magnetic field. Ignoring them parallel to the magnetic field implies that we are neglecting the ion acoustic waves and the (extremely long-wavelength) Jeans instability in the ions.

Linear System
To investigate the propagation, dissipation, and growth of HM waves in molecular clouds, we follow the original analysis by Jeans (1928;see also Spitzer 1978, § 13.3a;and Binney & Tremaine 1987, § 5.1), and assume that the zeroth-order state is uniform, static (i.e., vα = 0), and in equilibrium. 3 We consider only adiabatic perturbations; therefore, the net heating rate (Γn − Λn) on the right-hand side of equation (1e) vanishes.
We write any scalar quantity or component of a vector qtot(r, t) in the form qtot(r, t) = q0 + q(r, t), where q0 refers to the zeroth-order state, and the first-order quantity q satisfies the condition |q| ≪ |q0|. We thus obtain from equations (1a) -(1i) the linearised system Equation (4b) has been simplified by using the relation which expresses equilibrium of the ion density in the zeroth-order state, as a result of balance between the rate of creation of ions from ionisation of neutral matter by high-energy (E > ∼ 100 MeV) cosmic rays and the rate of destruction of ions by electron−molecular-ion dissociative recombinations. This relation allowed us to replace ζ CR by α dr (ρi,0/mi) 2 (µmH/ρn,0) in equation (4b). The quantity xi,0 ≡ ni,0/nn,0 in equation (4b) is the degree of ionisation (where ni,0 and nn,0 are the number densities of ions and neutrals in the unperturbed state). For an ideal gas (with only translational degrees of freedom), γ = 5/3 in equation (4e).
We seek plane-wave solutions of the form q(r, t) = q exp(ik · r − iωt), where k is the propagation vector, ω the frequency, and q the amplitude (in general, complex) of the perturbation. Equations (4a) -(4i) reduce to where Ca,0 ≡ γPn,0 ρn,0 is the adiabatic speed of sound in the neutrals, and is the (one-dimensional) neutral free-fall timescale. The quantities T , Pn, and ψ have been eliminated by using equations (4e), (4g), and (4h), respectively.
We note thatτni,0 = 1/ν ff,0 , where ν ff,0 is the collapse retardation factor, which is a parameter that measures the effectiveness with which magnetic forces are transmitted to the neutrals via neutral-ion collisions (Mouschovias 1982), and appears naturally in the timescale for the formation of protostellar cores by ambipolar diffusion (e.g., see reviews by Mouschovias 1987a, § 2.2.5;1987b, § 3.4;1991b, § 2.3.1;anddiscussions in Fiedler &Mouschovias 1992, 1993;Ciolek & Mouschovias 1993Basu & Mouschovias 1994. It is essentially the factor by which ambipolar diffusion in a magnetically supported cloud retards the formation and contraction of a protostellar fragment (or core) relative to free fall up to the stage at which the mass-to-flux ratio exceeds the critical value for collapse. It is discussed further in § 3.2.1.
Equations (10a) -(10k) govern the behaviour of small-amplitude disturbances in a weakly ionised cloud; they (without eq. [10e]) form a 10 × 10 homogeneous system. In general, the dispersion relationω(k) can be obtained by setting the determinant of the coefficients equal to zero. To each root (eigenvalue) of the dispersion relation there corresponds an eigenvector (or "mode"), whose components are the dependent variables appearing in equations (10a) -(10k). (Note that, once the dependent variables in eqs.
[4e], [4g], and [4h] to solve for the perturbed quantities T , Pn, and ψ, respectively.) Since, in general,ω is complex, modes with Im{ω} < 0 decay and those with Im{ω} > 0 grow exponentially in time. In what follows we investigate the propagation, dissipation, and growth of the allowable HM modes in typical interstellar molecular clouds.
The second mode is one in which density enhancements in the ions rapidly decay by dissociative recombinations of molecular ions and electrons; because the degree of ionisation is so small (xi,0 = 1.58 × 10 −7 ), this mode does not involve any motion of the neutrals and leaves the neutral density essentially unchanged (see Figs. 3a and 3c, lines labeled "i,rec"). Solving equation (12b) withρn = 0 andṽi,x = 0, one finds thatṽ φ = 0 and which is equal to 1.28 × 10 −4 for the model described here. It is again the case that this mode is independent ofk (see Figs. 2b and 3a -3d, lines labeled "i,rec"). The remaining two longitudinal modes are low-frequency modes, with |ω| ≪ 1/τin,0. Because the inertia of the ions along the magnetic field is small (ρi,0 ≪ ρn,0), the neutrals are able to sweep up the ions, and, as a result,ṽi,x ≃ṽn,x (see Figs. 3a and 3b). For these conditions, equation (15) yields the thermal Jeans modes (e.g., see Chandrasekhar 1961, Ch. XIII;or Spitzer 1978, § 13.3a). Therefore, fork > 1, i.e., the two acoustic waves have the same phase velocity, modified by gravity, but propagate in opposite directions (along the field lines). In the limitλ ≪ 1,ṽ φ = ±1 andτ d = ∞; i.e., these modes are undamped sound waves (recall that the unit of speed is Ca,0). At longer wavelengths, gravitational forces become increasingly more important and the phase velocity becomes less than unity (see Fig. 2a). The waves are gravitationally suppressed (i.e.,ṽ φ = 0) at wavelengths greater than the thermal Jeans wavelength (= 2πCa,0τ ff,0 = 0.611 pc, dimensionally). Forλ >λ J,th (i.e,k < 1), it follows from equation (18) that each of the Jeans modes splits into two separate, conjugate modes. One is a gravitational growth (or fragmentation) mode, with timescalẽ (see Fig. 2c). This is the classical Jeans instability. Asλ → ∞,τ gr → 1; dimensionally, this is just the free-fall timescale, τ ff,0 . The corresponding eigenvector is labeled as "ff+" in Figures 3a -3d. The other mode is one of exponential decay, with damping timescaleτ d also given by equation (21) (see eq. [18]); it is the curve labeled by "ff−" in Figure 2b. The eigenvector, also labeled by "ff−", is shown in Figures 3a -3d. This mode is one in which an initial density enhancement causes expansive motion, opposed by gravity, at such a rate that the density enhancement decreases to zero at the same time that the velocity vanishes. Hence, this is a monotonically decaying mode; no wave motion is involved. It is similar to the well-known classical cosmological problem of an expanding "flat" universe. Note that, asλ → ∞,τ d → 1.
For the convenience of the reader, Table 1 contains a list of all abbreviations (and their meaning) used to label the curves in all the figures of this paper.

Transverse Modes for k B0: Eigenfrequencies and Eigenvectors
Comparing the systems of equations for the modes with motions only in the y and only in the z directions, (13a) -(13c) and (14a) -(14c), we note that they are identical. Hence, for propagation along the field (θ = 0 • ),ω(k) is degenerate for these, transverse modes. Figure 4a displays the absolute value of their phase velocity, obtained by solving the dispersion relatioñ (Note that, because of the degeneracy, there are four different modes.) Damping timescales are shown in Figure 4b as functions ofλ. None of the modes are unstable. The absolute values of the eigenvectors are exhibited in Figure 5: |ṽn,y| and |ṽn,z| in Figure 5a, |ṽi,y| and |ṽi,z| in Figure 5b, and |By| and |Bz| in Figures 5c and 5d. Note that the |By| (and |Bz|) axis in Figure  5d is logarithmic in order to show the behaviour of the modes at small wavelengths. From the dispersion relation (23) and Figures 4a and 5a, it is evident that small-wavelength, high-frequency (|ω| > ∼ 1/τin,0) ion modes propagate withṽn,y,ṽn,z ≃ 0. Solving equations (13a) -(13c) (or, equivalently, eqs. [14a] -[14c]) in these limits yields (In deriving eq.
[11a] and [11b].) Forλ less than the ion Alfvén cutoff wavelength λA,i = 1.57 × 10 −2 (i.e., λA,i = 1.53 × 10 −3 pc) for the model cloud parameters specified at the beginning of § 3. Hence, forλ ≪λA,i, the waves are Alfvén waves, withṽ φ = ±ṽA,i,0 = 2.81 × 10 3 (see Fig. 4a, curve labeled "i,A"). There are four   waves in all. The two polarised in the y-direction are normal, shear Alfvén waves, and the two polarised in the z-direction are modified Alfvén waves. (At θ > 0 the latter set of waves are fast waves.) All the waves are damped on the timescaleτ d = 2τin,0 (see eq. [24] and Fig. 4b, curve labeled "i,A") because of collisions with the neutrals. It is noteworthy that the damping time is longer by a factor of 2 than that (τin,0) referring to the dissipation (momentum exchange) of ion streaming motion relative to the neutrals. This is so because, although a typical ion indeed loses memory of the collective (wave) motion on a timescalẽ τin,0, half of the wave energy is stored as potential energy in the magnetic field. Therefore, it takes twice as long for collisions to damp the wave than it takes them to damp ion streaming. Atλ =λA,i the ion Alfvén waves are critically damped. Forλ ≥λA,i, the ion-neutral collision frequency 1/τin,0 is greater than the wave (angular) frequencyω, and the waves can no longer propagate (ṽ φ = 0). Atλ =λA,i there is a bifurcation in the ion modes (see Fig. 4b). Two of the modes (one polarised in the y-direction and the other in the z-direction), corresponding to the negative root in equation (24), become ion collisional-decay modes (discussed earlier in § 3.1.1; curves labeled "i,coll" in Figs. 4b and 5a -5d), with In the limitλ ≫λA,i,By,Bz → 0 (see Figs. 5c and 5d), andτ d →τin,0, just as in equation (16). Thus, asλ increases, magnetic restoring forces on the ions become negligible, and the motion of the ions simply decays on a timescaleτin,0 because of collisions with the neutrals (see Fig. 4b, curve labeled "i,coll"). The remaining two ion modes (again, one polarised in the y-direction and the other in the z-direction), corresponding to the positive roots of the dispersion relation (eq.
[24]), are magnetically-driven ambipolar-diffusion modes (see Fig. 4b, curve labeled "AD"), in which the ions (and electrons) diffuse quasistatically (i.e., with negligible acceleration; this is equivalent to having |ω| ≪ 1/τin,0 in eqs. [13b] and [14b]) relative to the stationary neutrals. The damping timescale for these modes is which is the curve labeled as "AD" in Figure 4b. In the limitλ ≫λA,i,τ d becomes equal to the ambipolar-diffusion timescale, wherẽ is the ion ambipolar-diffusion coefficient. For the values of the free parameters cited at the beginning of § 3,Da,i = 3.52. In Figures 4b and 5 it is also evident that there exists two small-wavelength, low-frequency (|ω| ≪ 1/τin,0) neutral modes. In these modes the plasma and magnetic field lines are essentially stationary (i.e.,ṽi,y ≃ 0,ṽi,z ≃ 0 andBy ≃ 0,Bz ≃ 0). Under these constraints, equations (13a) -(13c) (and, similarly, eqs. These modes are the neutral collisional-decay modes (see curves labeled by "n,coll" in Figs. 4 and 5), in which the motion of the neutrals in the y and z directions decays due to collisions with ions that are held fixed in space by the magnetic field; v φ = 0, andτ d =τni,0 = 0.226 for these modes (see Fig. 4b).
[34b] and [29b]). As explained in the case of the ion Alfvén waves, it takes twice as long to damp a wave than it takes to damp streaming motion (or diffusion). The existence of λA,n = πvA,nτni for Alfvén waves in the neutrals was first shown by Kulsrud & Pearce (1969), who studied the excitation and propagation of HM waves in the intercloud medium due to cosmic-ray streaming (see, also, Parker 1967). Mouschovias (1987aMouschovias ( , 1991a discussed the importance of the lengthscalẽ λA,n and the thermal Jeans lengthscaleλ J,th in the formation of protostellar fragments (or cores) in self-gravitating molecular clouds. He proposed that fragmentation is initiated by the decay of HM waves due to magnetically-driven ambipolar diffusion and the almost simultaneous onset of a Jeans-like instability, due to gravitationally-driven ambipolar diffusion (see discussion in § 4 below).

Transverse Modes for k B0: Further Discussion of Eigenvectors
More insight in the physics of the transverse modes can be gained by understanding analytically certain key features of the eigenvectors shown in Figures 5a -5d. In the case of the ion Alfvén mode, we may ignore the motion of and the dissipation due to the neutrals at short wavelengths and we may use equations (13b), (13c), and the dispersion relationω ≃ ∓ṽA,i,0k (see eq. [23]), to find that vi,ỹ By ≃ ∓ṽA,i,0.
It is also clear from Figures 5a -5d that the eigenvectors for both the neutral and the ion collisional-decay modes behave exactly as expected on the basis of our discussion of these modes in relation to Figure 4b.
Motions in the y-direction are governed bỹ ωṽi,y = ĩ τin,0ṽ n,y − ĩ τin,0ṽ i,y , Similarly, the equations for the z-components of the neutral and ion velocities arẽ

Longitudinal Modes for k ⊥ B0: Eigenfrequencies
The dispersion relation is easily obtained from equations (40a) -(40e): (ω + 2iρi,0α m,dr ) ω 4 + i 1 τin,0 Because it is a fifth-order polynomial, there are five longitudinal modes in all. The phase velocities, damping timescales, and growth timescales for the longitudinal modes are displayed in Figures 6a, 6b, and 6c, respectively; eigenvectors are shown as functions ofλ in Figure 7a -7e. At small wavelengths there are again two high-frequency (|ω| ∼ 1/τin,0) ion modes, withṽn,x ≃ 0 (they are degenerate with respect to the direction of propagation). In this limit, the solution of the dispersion relation is identical with that given by equation (24). The modes it represents in this case are ion magnetosonic modes. Because the speed of sound of the ions (< Ca,0, since mi > µmH) is negligible compared to the ion Alfvén speed (see eq. [11c]), the phase velocity of the waves for λ <λA,i is again given by equation (26) and the characteristic decay time isτ d = 2τin,0 (see Figs. 6a and 6b, curves labeled "i,ms"). Forλ ≥λA,i (= 1.57 × 10 −2 for the typical model cloud) the ion magnetosonic waves cannot propagate because of frequent ion-neutral collisions. Instead, each mode bifurcates (see Fig. 6b), just as in the case of ion Alfvén waves described in § 3.1.2. One of the two resulting modes is an ion collisional decay mode (τ d =τin,0; curve labeled "i,coll" in Fig. 6b) and the other is a magnetically-driven ion ambipolar-diffusion mode, curve labeled "AD" in Figure 6b (τ d =λ 2 /4π 2D a forλ ≫λA,i; see discussion preceding eq. [29b]).
The third mode is an ion dissociative-recombination decay mode, as discussed in § 3.1.1, in which density enhancements in the ions decay rapidly [τ d = (2ρi,0α m,dr ) −1 = 1.28 × 10 −4 for the typical model cloud parameters] because of dissociative recombinations between molecular ions and electrons (see Fig. 6b, curve labeled "i,rec"). The only nonvanishing component of the eigenvector for this mode isρi,0 (see Figs. 7e and 7a -7d).
The pressure-driven diffusion mode corresponds to the positive root of equation (44). The neutrals are diffusing quasistatically (i.e., with |ω| ≪ 1/τni,0) through a background of effectively stationary ions and magnetic field. The eigenvector for this mode is labeled "PD" in Figs. 7a -7e. Forλs,n ≤λ ≤λ J,th , the damping timescale is   wherẽ DP ≡τni,0 (48) (= C 2 a,0 τni,0 dimensionally) is the neutral pressure-driven diffusion coefficient. In obtaining equation (47b) we have used the fact thatτ 2 ni,0 ≪ 1. For the representative model cloud used in this paper,DP = 0.226 (≃Da/16). We note that, at λ =λ J,th ,τ d = ∞ (see Fig. 6b, curve labeled "PD"). At this wavelength, the restoring pressure forces in this mode are exactly balanced by self-gravitational forces, and the system is on the verge of gravitational instability (ω = 0). Thus the Jeans instability manifests itself atλ J,th even in the presence of a magnetic field, as originally recognized by Langer (1979). Ambipolar diffusion allows this to happen, but, because of neutral-ion collisions, the growth time of the instability is longer than that of the nonmagnetic Jeans instability (compare eq. [49] below with eq. [21]). Forλ >λ J,th , density perturbations in the neutrals grow exponentially in time as a result of gravitational contraction of neutrals through essentially stationary ions attached to rigid magnetic field lines. The growth time for this instability can be obtained easily from equation (44) under the conditionsλ >λ J,th and 4τ 2 ni,0 1 − λ J,th /λ 2 ≪ 1: where ν ff,0 = 1/τni,0 (= τ ff,0 /τni,0 = 1/0.226 = 4.42) is the collapse retardation factor, discussed in the penultimate paragraph of § 2.3. It follows from equation (49) that, in the limitλ ≫λ J,th ,τ gr → ν ff,0 . The growth time of this ambipolar-diffusioninduced fragmentation is shown in Figure 6c (the part of the curve labeled by "AD,fr"). Dimensionally, the growth time for this mode is ν ff,0 τ ff,0 = τ 2 ff,0 /τni,0. It is the same as the nonlinear solution found analytically by Mouschovias (1979;see also 1983, 1987a, b, 1989 for the timescale of formation and evolution of protostellar cores (due to gravitationally-driven ambipolar diffusion). Numerical simulations (including the effects of grains, UV ionisation, and magnetic braking) of the formation of protostellar cores in magnetically supported molecular clouds have also found that the evolution occurs on this timescale (Fiedler & Mouschovias 1993;Basu & Mouschovias 1994, 1995a. The same timescale was found in the one-dimensional similarity solution of Scott (1984). It is clear from Figure 6c that, as predicted,τ gr would tend to 1/τni,0 = 4.42 forλ ≫λ J,th ; see the inflection point in the curve (labeled "AD,fr"). However,τ gr falls below this would-be asymptotic value because, forλ > ∼λ J,mag , whereλJ,mag is the magnetic Jeans wavelength (= 25.6; see eq. [54] below), gravitational forces on the neutrals, transmitted to the ions by neutral-ion collisions, exceed the restoring magnetic forces on the ions, and the ions and magnetic field are no longer able to remain stationary; the mode becomes a gravitational (Jeans) instability against the magnetic field, as originally found by Chandrasekhar & Fermi (1953); 4 see Figure 6c, part of curve beyond the inflecton point, labeled "ff+". The approximationsṽi,x ≃ 0,Bz ≃ 0 (which were used in deriving eq. [44]) are no longer valid, and equation (49) no longer describes this mode. The proper expressions are derived below.
The conjugate Jeans (or classical cosmological) mode was discussed earlier in connection with equation (18). In this mode the neutrals, ions, and magnetic field lines are well coupled, but the gravitational forces prevent them from oscillating; motions are damped monotonically. The damping time, corresponding to the negative root of equation (51), forλ ≥λJ,mag is given bỹ It is clear that, in the limitλ → ∞,τ d → 1 (i.e., τ d = τ ff,0 ); see Figure 6b, curve labeled "ff−". Finally, the gravitational instability mode (see eq. [49]) also changes behaviour atλJ,mag, as discussed above. Forλ > λJ,mag, the neutrals, plasma, and magnetic field lines are well coupled and behave like a single fluid; self-gravity overwhelms restoring magnetic and thermal-pressure forces, and the mode behaves as a classical Jeans instability, in which density perturbations grow exponentially in time with a timescalẽ which is equal to the damping time of the conjugate Jeans mode (see eq. [56]). From this equation we see thatτ gr → 1 (i.e., τgr = τ ff,0 ) forλ → ∞, in agreement with the long-wavelength behaviour ofτ gr exhibited in Figure 6c.

Longitudinal Modes for k ⊥ B0: Eigenvectors
The main features of the eigenvectors are as follows.
(a) The dominant component of the ion magnetosonic and the ion ambipolar-diffusion modes isṽi,x (see Fig. 7b, curves labeled "i,ms" and "AD"). It does not vanish atλ =λA,i; in fact, it hardly changes from 1, because the ion-AD mode maintainsṽi,x large beyondλA,i. Only when the ion-AD mode induces motions in the neutrals doesṽi,x begin to decrease, as the magnitude ofṽn,x increases asλ →λms,n (see Figs. 7b and 7a, curves labeled "i,ms" and "AD"). The z-component of the magnetic fieldBz in Figure 7c actually does not vanish. It is equal toṽi,x/ṽi,ms (see eq. [35]), butṽi,ms ≃ vA,i,0 = 2.81×10 3 .
(b) The velocityṽi,x is also the dominant component of the eigenvector of the ion collisional-decay mode at all wavelengths (see Fig. 7b, curve labeled "i,coll").
(c) The eigenvector of the neutral acoustic mode has significant componentsṽn,x andρn. Atλs,n, beyond which neutral sound waves do not exist and at which the acoustic mode bifurcates into the neutral collisional-decay mode and the pressure-diffusion mode (see Fig. 6b), the collisional-decay mode is responsible for the increase in |ṽn,x| as |ρn| decreases to zero. Because motion is induced in the ions asλ approachesλs,n (see Fig. 7b, curve labeled "n,coll"), |ṽn,x| does not reach unity.
(d) The most significant component of the PD mode isρn (see Fig. 7d);ρn increases asλ increases fromλs,n toλ J,th , at which wavelengthρn reaches a maximum. Beyondλ J,th , ambipolar-diffusion-induced fragmentation sets in. At exactlyλ J,th , ρn is the only nonvanishing component of the eigenvector of the AD,fr mode (see Figs. 7a -7d). Asλ increases beyondλ J,th , the field lines begin to be compressed as the neutrals begin to couple to the ions, andBz increases (see Fig. 7c) -at the expense ofρn -whileṽi,x andṽn,x are negligible. Asλ approachesλJ,mag, the ambipolar-diffusion-fragmentation mode induces significant velocitiesṽn,x andṽi,x, whileρn remains large (see Figs. 7a, 7b and 7d). As discussed in relation to Figure 6c, the AD,fr mode turns into the Jeans free-fall mode (ff+) beyondλ ≃λJ,mag, as is clearly shown in Figures 7a and 7b (see curves labeled "AD,fr" and "ff+"). The bifurcation of the neutral magnetosonic mode into the AD and conjugate Jeans (or, cosmological, ff-) modes beyondλ ≃λJ,mag, discussed in relation to Figure 6b, is also seen clearly in Figures 7a and 7b. [Note: If we had plotted only the real part of the eigenvector, we would have found, for the neutral magnetosonic mode, that Re{ρn} decreases discontinuously to zero, then increases smoothly, reaches a maximum, and then vanishes again as the magnetosonic waves are suppressed by gravity.] (e) Slightly beyondλms,n, magnetosonic waves exist in the neutrals but, in addition to maintainingṽn,x large, they cause significant motion in the ions as well (see Figs. 7a and 7b, curves labeled "n,ms"). The quantitiesBz andρn are also nonnegligible. Note that atλ J,th ,ṽn,x andρn decrease whileṽi,x increases.

Transverse Modes for k ⊥ B0
Examination of equation (41c) and both of which are independent of wavelength. (Note: there are four modes in all, two corresponding to the first solution, and two corresponding to the second.) In the first mode the ions and neutrals move together withṽn,y =ṽi,y (orṽn,z =ṽi,z); hence, there are no frictional forces between the two species, andτ d = ∞. The second mode consists of oppositely flowing streams of ions and neutrals. The momentum of each species decays by collisions with the other; hence, ion-neutral and neutral-ion collisions occur in parallel, and the net damping time for this mode is the harmonic mean ofτin,0 andτni,0 (see eq. [58b]). For the typical model cloud, the ion fluid has a much smaller inertia than the neutral fluid and, therefore, the streaming velocity of the ions is much greater than that of the neutrals. As a result, the motion of the ions is collisionally damped on the timescaleτ d ≃τin,0.

Propagation at Angles
Equations (10a) -(10k) reveal that, for 0 • ≤ θ ≤ 90 • , the equations for the modes with motions in the y-direction are uncoupled from those with motions in the (x, z)-plane. Motions in the x-and z-direction, however, are coupled.

Transverse Modes
Figures 8a and 8b display the magnitude of the phase velocities and the damping timescales of the three transverse modes for propagation at an angle θ = 45 • with respect to B0. Because they are incompressible, none of these modes can become gravitationally unstable. Eigenvectors are displayed in Figure 9. The dispersion relation for ion Alfvén waves, obtained from equations (10f) -(10h), has solutions given by equation (24) Figure 11. Eigenvalues of transverse modes as functions of wavelength, at an angle of propagation θ = 80 • with respect to B 0 . All normalisations, labels, and symbols are as in Fig. 8  ± sign denotes degeneracy with respect to the direction of propagation. The damping timescale for these waves is 2τin,0 (see Fig. 8b, curve labeled "i,A").
Collisions with the neutrals cut off the propagation of these waves forλ ≥λA,i cos θ = 1.11 × 10 −2 (see Fig 8a). For wavelengths greater than this value, each mode bifurcates into an ion collisional-decay mode and an (ion) ambipolar-diffusion mode. The damping timescale of the former mode is given by equation (27) and is shown in Figure 8b (curve labeled "i,A"); it is identical withτ d of the normal Alfvén waves propagating along B0, except for the fact thatλA,i cos θ replacesλA,i in that expression (see discussion following eq. [24]). For very largeλ,τ d =τin,0 (see Fig. 8b, curve labeled "i,coll"). For the ambipolar-diffusion mode (Fig. 8b, curve labeled "AD"), the damping timescale is the same as that given by equation (28), but againλA,i cos θ replacesλA,i; hence, for long wavelengths, it follows from equation (29a) that At small wavelengths the third mode is a low-frequency neutral collisional-decay mode, withṽi,y ≃ 0 (see Fig. 9b, curve labeled "n,coll"). The purely imaginary frequency of this mode is given by equation (30),ω = −i/τni,0; hence the damping timescale isτ d =τni,0 (= 0.226 for the typical model cloud; see Fig. 8b, curve labeled "n,coll"). At longer wavelengths, the motions of the neutrals and the ions become better coupled, and the neutral collisional-decay mode merges with the ion ambipolar-diffusion mode (see Fig. 8b), as discussed in § § 3.1.2 and § 3.2.1. For longer wavelengths, Alfvén waves in the neutrals can be sustained, as seen in Figure 8a. The frequencies of these modes are given by equations (31a) and (31b), but withṽA,i,0 cos θ andṽA,n,0 cos θ replacingṽA,i,0 andṽA,n,0, respectively. The lower cutoff wavelength, below which these waves cannot propagate, isλA,n cos θ = 1.98. They have phase velocitỹ v φ = ±ṽA,n,0 cos θ 1 − λ A,n cos θ and damping timescalẽ Asλ → ∞, |ṽ φ | →ṽA,n,0 cos θ = 2.79. The components of the eigenvectors of these modes displayed in Figures 9a -9d are qualitatively similar to those shown in Figures 5a -5d (for transverse modes propagating along B0), but with (small) quantitative decrease in the values of the critical wavelengths. We therefore do not discuss them further for economy of space.
In Figures 10a and 10b we display |ṽ φ | andτ d , respectively, for the three transverse modes propagating at an angle θ = 10 • with respect to the magnetic field. The qualitative solutions of the dispersion relations for these modes are the same as for those described above for propagation at an angle θ = 45 • . Quantitatively, the differences in Figures 10a -10b and 8a -8b stem from the different values of the numerical factors cos θ and sin θ. Similarly, the magnitude of the phase velocities and the damping timescales for the transverse modes propagating at an angle θ = 80 • with respect to B0 are shown in Figures 11a  and 11b. They again differ from those for θ = 45 • because of the factors cos θ and sin θ. Although qualitatively the modes at θ = 10 • and θ = 80 • are the same, comparison of Figures 10 and 11 reveals that the quantitative differences in the magnitude of the phase velocities and the cutoff wavelengths are substantial. Thus the critical wavelengths that determine which modes can or cannot propagate in molecular clouds are expected to be very different depending on their direction of propagation with respect to the mean magnetic field B0.

Modes with Motions in the (x, z)-Plane
Figures 12a, 12b, and 12c exhibit, respectively, the magnitude of the phase velocityṽ φ , the damping timescaleτ d , and the growth timeτ gr as functions ofλ for the specific case of propagation at θ = 45 • with respect to B0 for modes with motions in the (x, z)-plane. Eigenvectors are displayed in Figures 13a -13g. There are 7 modes in all displayed in these figures.
As in the preceding sections, one of the ion modes is a collisional-decay mode (see Figs. 12 and 13, curves labeled "i,coll"), with the ions streaming through a fixed background of neutral particles (ṽn,x =ṽn,z = 0); for this caseBz = 0, and the ions move withṽi,x =ṽi,z. The frequency is given by equation (16); hence, the damping timescale isτ d =τin,0 (= 4.45 × 10 −7 for the typical model cloud). This is the horizontal line labeled by "i,coll" in Figure 12b.
There also exist two high-frequency ion wave modes. These waves are ion fast modes. For these modes,ṽi,x/ṽi,z = − sinθ/ cosθ = −1 at θ = 45 • (see Figs. 13c and 13d, curves labeled "i,fast"; in these Figures the "i,fast" curves coincide with the "i,coll" curves). Since the ion Alfvén speed and the magnetosonic speed (in this typical case) are essentially the same, the dispersion relation describing them is the same as equation (24); it is again the case that the waves are cut off atλ ≥λA,i (= 1.57 × 10 −2 ; see Fig. 12a). For wavelengths >λA,i, each ion mode bifurcates into an ion collisional-decay mode and an ion ambipolar-diffusion mode (see Fig. 12b). The decay timescales for these two modes are the same as those previously discussed for the cases θ = 0 • and 90 • (see § 3.1.2 and § 3.2.1). The fourth mode is the ion dissociative-decay mode, discussed previously in § § 3.1.1 and 3.2.1. In this mode, nonpropagating density perturbations in the ions decay by dissociative recombinations of molecular ions and electrons. The decay time isτ d = (2ρi,0α m,dr ) −1 (= 1.28 × 10 −4 for the typical model cloud); this is the horizontal line labeled by "i,rec" in Figure 12b.
For |ω| > ∼ 1/τni,0, the term in brackets in equation (63) is essentially unity, andṽi,x ≃ṽn,x cos 2 θ. The factor cos 2 θ is a measure of the opposition presented to the neutrals by the ions, which are attached to magnetic field lines. If θ = 0 • the ions are effectively inertialess, and the neutrals sweep them up easily, so thatṽi,x =ṽn,x. However, if θ = 90 • ,ṽi,x ≃ 0, i.e., the ions are held in place by the magnetic field as the neutrals move through them. In this case, the ions present the stiffest opposition to the neutral motion. For |ω| > ∼ 1/τni,0, we insert equation (63) in equation (10c) to find the solutioñ [Note that for θ = 0 • or 90 • , equation (64)  As,n(θ) , where As,n(θ) ≡ sin 4 θ + (2τni,0) 2 1 + (2τni,0) 2 these modes are sound waves, with v φ = ± 1 − λ As,n(θ) λs,n 2 1/2 (66) (see Fig. 12a) andτ d = 2τni,0/ sin 2 θ = 0.904 for θ = 45 • (see horizontal line labeled "acoustic" in Fig. 12b). As,n(θ) is the acoustic-wave angular (or, stiffness) parameter: for θ = 90 • equations (65a) and (66) yield our earlier result, that sound waves propagate only forλ ≤λs,n, while for θ = 0 • we recover our other earlier result thatλ J,th is the upper cutoff wavelength for these waves. The sound-wave upper cutoff wavelength (65a) is relevant only if the waves do not transition into neutral slow modes at larger wavelengths (see discussion below). It turns out that this depends on the angle θ of propagation: the upper cutoff equation is applicable only to sound waves propagating at θ > θmax (for the typical model cloud, θmax = 62.9 • -see below). For sound waves propagating at such angles, there is a mode bifurcation at the cutoff (65a). The dispersion relation (64) reveals that at largerλ one of the modes becomes a neutral collisional-decay mode, with damping timescaleτ d →τni,0/ sin 2 θ, and the other mode becomes a pressure-driven diffusion mode with damping timescalẽ This diffusion timescale is the same as that of equation (47a), except that it is multiplied by sin 2 θ, which reflects the reduced effectiveness of neutral-ion collisions in slowing down the neutrals at angles θ < 90 • with respect to the magnetic field. Atλ =λ J,th the timescale (67) becomes infinite. Forλ >λ J,th there is a neutral ambipolar-diffusion-induced gravitational fragmentation mode at the angle θ (> θmax), having a growth time τ gr = sin 2 θ/τni,0 The fragmentation time (68) is equal to the growth time (49) multiplied by sin 2 θ, again reflecting the reduced collisional resistance on the neutrals by magnetically-coupled ions when the angle of propagation is less than 90 • with respect to the field direction. Forλ ≫λ J,th , equation (68) yieldsτ gr → sin 2 θ/τni,0 = ν ff,0 sin 2 θ. However, similar to what occurs in the case for θ = 90 • (see § 3.2.1 and Fig. 5c), this limit will not be attained because the approximation of stationary field lines breaks down whenλ becomes > ∼λ J,mag . When this happens, magnetic forces are overwhelmed by self-gravitational forces and τ gr → 1 at these larger wavelengths.
From equation (64) one finds that collisional damping of the motion of the neutrals by ions causes |ṽ φ | to become less than unity forλ near the value given by the right-hand side of equation (65a). However, for 0 < θ < 90 • , the sound waves are not always cut off at this wavelength. This is due to the fact that the bracketed term in the expression for the x-component of the ion velocity (eq. [63]) is no longer essentially unity (because |ω| ∼ 1/τni,0) at these wavelengths; thus, the dispersion relation (64) is no longer valid. Physically, this is a result of the fact that the ions move readily with the neutrals in the x-direction at these frequencies; this means that the waves suffer less damping, because the frictional force on the neutrals is reduced when the ions and the neutrals move together. (We note that, for θ < 90 • , equation [63] yieldsṽi,x ≃ṽn,x in the limit |ω| ≪ 1/τni,0.) Forλ sufficiently large, such that |ω| < ∼ cos 2 θ/τni,0, the frequency of the neutral waves is Hence, for these conditions, the phase velocity is and These waves are neutral slow modes, modified by gravity, and cannot propagate at wavelengthsλ ≥λ J,th (= 2π for the typical model cloud); at θ = 0 • the slow mode dispersion relation (69) is identical to the relation for the Jeans mode (eq. [18]). From the low-frequency condition used to derive equation (69) it is found that slow modes exist only for wavelengths where Sn(θ) ≡ 1 + (2τni,0) 2 1/2 2| cos θ| , provided that As,n(θ)Sn(θ) < ∼ 1 .
The quantity Sn(θ) is the angular slow-mode factor. The relation (73) is derived from the requirement that the slow modes arise from the acoustic wave modes. For this mode conversion to occur, the minimum slow mode wavelength (72a) must be less than or equal to the acoustic-wave upper cutoff wavelength (65a); thus, the inequality (73) follows. This is equivalent to having θ ≤ θmax, where θmax is defined by the condition cos 2 θmax sin 4 θmax + (2τni,0) 2 = 1 4 .
If θ ≤ θmax, slow modes emerge from the acoustic modes (without a bifurcation) and propagate forλs,nSn(θ) < ∼λ ≤λ J,th . Otherwise, when θ > θmax, the acoustic waves have an upper cutoff and mode bifurcation occurs at the wavelength (65a); in that case, there are no slow modes. 5 For the typical model cloud, withτni,0 = 0.226, θmax = 62.9 • . At θ = 45 • the slow mode minimum wavelength λs,nSn(45 • ) = 2.01, and the transition from sound waves to slow modes can be seen in Figures 12a and 12b (curves labeled "n,slow") to occur at this wavelength; Figures 13a and 13c show that, for these modes,ṽn,x ≃ṽi,x at wavelengths greater than or equal to the transition wavelength.
Note that, asλ →λ J,th from below, equation (71) shows that the slow mode damping timescaleτ d → ∞ (see curve labeled by "n,slow" in Fig. 12b). Forλ >λ J,th there are again the two conjugate modes: the gravitational instability (or Jeans) mode (see Fig. 12c, curve labeled "ff+") and the "cosmological" mode (see curve labeled "ff−" in Fig. 12b; although this curve "crosses" the curve labeled "n,fast" in Fig. 12b, the two modes do not actually interact); asλ → ∞, both the growth timescale of the unstable mode and the damping timescale of the cosmological mode go to unity (i.e., in dimensional form, τgr = τ d = τ ff,0 ), as seen in Figures 12b and 12c, respectively.
The other mode affecting the neutrals at short wavelengths is a neutral collisional-decay mode ("n,coll"). The velocity of the neutrals at smallλ is predominantly in the z-direction for this mode (see Figs. 13a and 13b); the frequency is purely imaginary and given bỹ Thus,τ d =τni,0/ cos 2 θ = 0.452 (see Fig. 12b, curve labeled "n,coll"). Forλ greater than the value given by the right-hand side of equation (65a), it is again the case that the motion of the ions and magnetic field becomes better coupled to that of the neutrals. In this wavelength regime |ṽn,x| ≃ |ṽn,z| ≃ |ṽi,x| ≃ |ṽi,z| (see . This mode merges with the ion ambipolar-diffusion mode atλ =λms,n (see Fig. 12b), and, forλ >λms,n, fast modes are excited in the neutrals (curves labeled "n,fast" in Figs. 12a, 12b, and 13a -13d), degenerate with respect to the direction of propagation. In these modes, the polarisation is given byṽn,x/ṽn,z = − sinθ/ cosθ (= −1 for θ = 45 • ). Hence,ṽn · B0 = 0, i.e., the fast modes are polarised perpendicular to the magnetic field. The dispersion relation for these modes is essentially the same as equation (51) because of the fact thatṽ 2 A,n,0 ≫ 1. They decay, as they propagate, on the ambipolar-diffusion timescale (see Fig. 12b). Forλ ≥λJ,mag, these waves tend to get suppressed by gravity.
Forλ >λJ,mag thermal and magnetic restoring forces in the x-direction are overwhelmed by gravitational forces, making vn,x ≃ 0 (see Fig. 13a). Hence the modes become essentially incompressible. Waves are still able to propagate for longer wavelengths, however, because of the transverse restoring magnetic tension force (i.e.,Bz = 0; see Fig. 13g). Solving equations (10a) -(10d) and (10i) -(10k) with the conditionsṽn,x = 0 =ρn in the limit |ω| ≪ 1/τin,0, we find ω = ±ṽA,n,0k cos θ 1 − ṽA,n,0τni,0k 2 cosθ Hence, these modes are modified Alfvén waves in the neutrals, with v φ = ±ṽA,n,0 cos θ 1 − λ A,ñ λ cos θ 2 1/2 (77) In the limitλ → ∞, |ṽ φ | ≃ṽA,n,0 cos θ = 2.79, in agreement with the long-wavelength behaviour of these modes shown in Figure 12a (curves labeled "n,fast" and "n,A"). This is yet another example of a transition between wave modes without bifurcation. Figures 14a, 14b, and 14c show |ṽ φ |,τ d , andτ gr as functions ofλ for the seven different modes with motions in the (x, z)plane propagating at an angle θ = 10 • with respect to B0. Comparison with Figures 12a -12c reveals that the qualitative behaviour of the various modes as functions of wavelength is the same as in the case of propagation at θ = 45 • . The quantitative differences stem from the numerical factors cos θ and sin θ, which become substantial for θ approaching 0 • or 90 • . As Figure  14a shows clearly, wave modes in the neutrals exist at all wavelengths and their decay times are very long (see Fig. 14b, curves labeled "acoustic", "n,slow", and "n,fast"). Figures 15a, -15c show the same quantities as Figures 14a -14c but for propagation at θ = 80 • with respect to the unperturbed magnetic field B0. There are no slow modes in Figure 15 (unlike the cases in Figs. 12 and 14) because θ > θmax for that angle of propagation in the typical model. Instead, the sound waves are cut off at the maximum wavelengthλ =λs,n/As,n(80 • ) = 2.66, where there is a bifurcation. At wavelengths greater than this maximum, the modes are a pressure-driven diffusion mode ("PD" curve) and a neutral collisional-decay mode ("n,coll"). There is also an ambipolar diffusion-induced fragmentation mode seen in Figure 15c ("AD,fr" curve), which approaches the predicted limiting value (see eq. [68]) of sin 2 80 • /τni,0 = 4.29 atλ just belowλJ,mag (= 25.6 in the typical model cloud).

SUMMARY AND DISCUSSION
We have obtained and analyzed the dispersion relations for MHD wave modes and instabilities for different directions of propagation with respect to the zeroth-order magnetic field B0 in a two-fluid weakly ionised system, and we have applied the results to a typical interstellar molecular cloud. The system of equations has four dimensionless free parameters,τni,0, τin,0,ṽA,i,0, andα m,dr . They represent, respectively, the neutral-ion (momentum-exchange) collision time and the ion-neutral collision time in units of the free-fall time of the zeroth-order state, the Alfvén speed in the ions in units of the adiabatic speed of sound in the neutrals, and the dissociative recombination coefficient (see eq. [11d]). (Because of ionisation equilibrium in the zeroth-order state, the dimensionless cosmic-ray ionisation rateζ CR is expressible in terms ofα m,dr ,τin,0 andτni,0.) There are two distinct kinds of ambipolar diffusion, whose combined effect is unavoidable in typical molecular clouds and has crucial consequences on their evolution: (a) In the presence of hydromagnetic (HM) waves or turbulence, the tension of field lines (or the outward pressure due to compressed field lines) drives the motion of charged particles relative to the neutrals, with the tendency/consequence to straighten out the bent or tangled magnetic field lines (or to move compressed field lines apart, toward a more uniform configuration). The timescale of this process is proportional to the square of the wavelength of the HM waves (or the characteristic length of the field-line tangling, or the magnitude of the field gradient) -see eqs. (34a) and (80). For lengthscales typical of molecular cloud cores ( < ∼ 0.1 pc), it is much smaller than the free-fall time. This is the magnetically-driven ambipolar diffusion. It is this kind of ambipolar diffusion which is responsible for the observed large-scale ordering of polarisation vectors, indicating large-scale ordering of the magnetic field lines in molecular clouds.
(b) Gravitationally-driven ambipolar diffusion sets in with a short enough timescale, but longer than the free-fall time, in the deep interiors of self-gravitating clouds, where the degree of ionisation is xi < ∼ 10 −7 . Its onset can be spontaneous or initiated as a result of magnetically-driven ambipolar diffusion depriving a self-gravitating cloud of any support that most short-wavelength HM waves (or turbulence) may have provided against the cloud's self-gravity (Mouschovias 1987a). It (i) allows the clouds to fragment as neutral particles contract through almost stationary field lines (and the attached charged particles); and, consequently, (ii) increases the mass-to-flux ratio of the forming fragments (or cores). Once the critical massto-flux ratio (Mouschovias & Spitzer 1976) of fragments is exceeded, dynamical contraction ensues, while the cloud envelopes remain magnetically supported, as found by numerical simulations starting with Fiedler & Mouschovias (1993) and confirmed by numerous observations. Hydromagnetic (HM) waves with phase velocity |v φ | ≃ vA,n,0 = 0.96 B0 30 µG 2 × 10 3 cm −3 nn,0 can propagate in all directions with respect to B0, provided that λ > ∼ λA,n, where λA,n = πvA,n,0τni,0 = 0.22 B0 30 µG 2 × 10 3 cm −3 nn,0 The long-wavelength waves are long-lived; the decay time is the magnetically-driven ambipolar-diffusion timescale which is to be distinguished from the growth time of gravitationally-driven ambipolar-diffusion, relevant to fragmentation of molecular clouds into self-gravitating cores; namely, τ AD,fr = τ ff,0 2 τni,0 ≃ 1.1 × 10 6 xi,0 10 −7 yr.
The (one-dimensional) free-fall timescale at the density nn,0 = 2 × 10 3 cm −3 as λ → ∞ is τ ff,0 = (4πGρn,0) −1/2 = 3.9 × 10 5 yr. The nonlinear counterparts of these modes have been shown to explain quantitatively the observed highly supersonic but subAlfvénic linewidths in molecular clouds, their cores, and even in OH and H2O masers in which the strength of the magnetic field has been measured (Mouschovias & Psaltis 1995;Mouschovias et al. 2006). Most HM waves with λ < λA,n cannot propagate in the neutrals because they are damped rapidly by ambipolar diffusion. This means that there cannot be any contribution from this wavelength regime to the spectrum of HM "turbulence" in molecular clouds (Mouschovias & Psaltis 2011, in preparation), which may provide clouds with a source of nonthermal pressure. This led Mouschovias (1987a) to argue that the decay of HM waves by ambipolar diffusion on wavelength scales < ∼ 0.1−0.3 pc can initiate the formation of protostellar cores in otherwise magnetically supported clouds (see also Mouschovias 1991a). Damping of short-wavelength HM waves by ambipolar diffusion has also been proposed (Mouschovias 1987a, § 2.2.5) as the cause of the observed narrowing and thermalization of linewidths with increasing column density, as observed, for example, by Baudry et al. (1981) in the cloud TMC 2. Depending on the angle of propagation with respect to the unperturbed magnetic field B0, however, we find that certain long-lived modes in the neutrals exist at all wavelengths, while ion modes usually damp very rapidly even at short wavelengths. This may explain the observations by Li & Houde (2008) in the M17 molecular cloud, which show neutral motions at small lengthscales but much smaller ion motions on the same scales.
Gravitational instability is found to set in at λ = λ J,th (the thermal Jeans wavelength) for all θ, in agreement with Chandrasekhar & Fermi (1953). For θ ≃ 90 • , though, the timescale for the instability is ≃ ν ff,0 τ ff,0 = τ 2 ff,0 /τni,0, which is the gravitationally-driven ambipolar-diffusion timescale. This is in agreement with the result for the ambipolar-diffusioninitiated formation and contraction of protostellar fragments (or cores) obtained from nonlinear calculations analytically by Mouschovias (1979Mouschovias ( , 1989 and Mouschovias & Paleologou (1986), and numerically by Fiedler & Mouschovias (1993) for a two-fluid system (neutral molecules and plasma),  for a four-fluid system (i.e., neutral molecules, plasma, negatively-charged, and neutral grains) , and Basu & Mouschovias (1994, 1995a for a two-fluid system including rotation and magnetic braking. It is also in agreement with the results of simulations by Tassis & Mouschovias (2007a, b, c) and Kunz & Mouschovias (2009, 2010 for a six-fluid system (neutral molecules, atomic and molecular ions, electrons, negatively-charged, positively-charged, and neutral grains).
Another recent success of the linear theory (applied to model molecular clouds flattened along the magnetic field) in predicting inherently nonlinear phenomena was demonstrated by Kunz & Mouschovias (2010). They obtained analytically the Core Mass Function (CMF) resulting from gravitationally-driven ambipolar-diffusion-induced fragmentation of molecular clouds. The predicted CMF is in excellent agreement with observations of more than 300 cores in Orion (Nutter & Ward-Thompson 2007), not only at the high-mass end (for which many alternative models can obtain such agreement), but also at the turnover mass and low-mass end. ion Alfvén-wave upper cutoffλ A,i = 4πṽ A,i,0τin,0 0.0157 λ A,i = 4πv A,i,0 τ in,0 ion magnetosonic-wave upper cutoffλ A,i 0.0157 λ A,i acoustic-wave upper cutoffλs,n = 4πτ ni,0 /[1 + (2τ ni,0 ) 2 ] 1/2 2.59 λs,n = 4πC a,0 τ ni,0 /[1 + (2τ ni,0 ) 2 ] 1/2 neutral Alfvén-wave lower cutoffλ A,n = πṽ A,n,0τni,0 2.80 λ A,n = πv A,n,0 τ ni,0 neutral magnetosonic-wave lower cutoffλms,n = (ṽ A,n,0 /ṽ ms,n,0 )λ A,n 2.72 λms,n = (v A,n,0 /v ms,n,0 )λ A,n thermal Jeansλ J,th = 2π 2π λ J,th = 2πC a,0 τ ff,0 magnetic Jeansλ J,mag = 2πṽ ms,n,0 25.6 λ J,mag = 2πv ms,n,0 τ ff,0 † To convert to parsecs, multiply by the unit of length for the typical model cloud, C a,0 τ ff,0 = 9.72 × 10 −2 pc. Table 2 lists all the critical wavelengths present in a two-fluid system, such as a typical molecular cloud. The name of each critical wavelength, its defining dimensionless expression and typical value, and the corresponding dimensional expression are listed in columns 1 -4, respectively. Table 3 summarizes conveniently all the modes that can exist in a two-fluid system (typical molecular cloud) for propagation parallel, perpendicular, and at an arbitrary angle with respect to the unperturbed magnetic field. The name of each mode is shown in the first column, and the wavelength range in which the mode exists is shown in the second column.
Although in a single-fluid system linear modes are independent of one another, in a multifluid system (such as a molecular cloud) this is not the case. A mode in one fluid can bifurcate due to interaction with the other fluid and give rise to two distinct daughter modes. Such interaction is also responsible for the reverse phenomenon of mode merging. Moreover, in a multifluid system, a wave mode in one fluid/species (e.g. a sound wave in the neutrals) can transition into another wave mode (e.g., a slow MHD wave in the neutrals) without bifurcation. In other words, linear waves in a multifluid system exhibit behaviour that only nonlinear waves exhibit in a single-fluid system.
In a following paper we present a study of the free parameters, spanning the range of observationally allowed values, to examine their effect on HM waves and instabilities in molecular clouds.
This work was carried out in 1987 -1989, when the bulk of this paper was also written. The paper was updated in 1994 -1996, but never submitted for publication due to spatial dispersal of its three authors and their involvement in code development for nonlinear calculations. Its Introduction and some other parts have been updated to refer to relevant work published since then. The work was supported in part by the National Science Foundation under grant AST-93-20250, and also by the New York Center for Astrobiology under grant # NNA09DA80A. TChM acknowledges the hospitality of the Alexander von Humboldt Foundation and of the Max-Planck-Institut für Extraterrestrische Physik during the writing of part of this paper in 1994. GC is supported by the New York Center for Astrobiology (a member of the NASA Astrobiology Institute) under grant # NNA09DA80A.