{'fraction_non_alphanumeric': 0.033992486108299025, 'fraction_numerical': 0.01370047917932105, 'mean_word_length': 4.661953891224048, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Protostar Formation in Magnetic Molecular Clouds beyond Ion Detachment: I. Formulation of the Problem and Method of Solution

We formulate the problem of the formation of magnetically supercritical cores in magnetically subcritical parent molecular clouds, and the subsequent collapse of the cores to high densities, past the detachment of ions from magnetic field lines and into the opaque regime. We employ the six-fluid MHD equations, accounting for the effects of grains (negative, positive and neutral) including their inelastic collisions with other species. We do not assume that the magnetic flux is frozen in any of the charged species. We derive a generalized Ohm's law that explicitly distinguishes between flux advection (and the associated process of ambipolar diffusion) and Ohmic dissipation, in order to assess the contribution of each mechanism to the increase of the mass-to-flux ratio of the central parts of a collapsing core and possibly to the resolution of the magnetic flux problem of star formation. We show how our formulation is related to and can be transformed into the traditional, directional formulation of the generalized Ohm's law, and we derive formulae for the perpendicular, parallel and Hall conductivities entering the latter, which include, for the first time, the effect of inelastic collisions between grains. In addition, we present a general (valid in any geometry) solution for the velocities of charged species as functions of the velocity of the neutrals and of the effective flux velocity (which can in turn be calculated from the dynamics of the system and Faraday's law). The last two sets of formulae can be adapted for use in any general non-ideal MHD code to study phenomena beyond star formation in magnetic clouds. The results, including a detailed parameter study, are presented in two accompanying papers.

In the Milky Way, stars are observed to form inside molecular clouds, located primarily along spiral arms. Molecular clouds do not exhibit significant ordered expansion or contraction velocities. They are however observed to have supersonic but sub-Alfvénic linewidths. Polarization (Hildebrand et al. 1999;Rao et al. 1998;Lai et al. 2002) and Zeeman (Crutcher, Kazés, & Troland 1987;Crutcher et al. 1994) observations have revealed that such clouds are threaded by ordered magnetic fields, of magnitude sufficient to support them against their self-gravity.
In such magnetically supported (magnetically subcritical) clouds, magnetically supercritical fragments (or cores) form by the process of ambipolar diffusion, which allows neutral molecules to diffuse through charged particles and magnetic field lines, on a timescale 1 -10 Myr (depending on the local physical conditions). When selfgravity becomes strong enough to overwhelm the magnetic forces, the fragments begin to contract dynamically. The initial stages of this contraction are isothermal. However, the column density eventually increases enough that radiation is trapped, the temperature rises, and hydrostatic protostellar cores form in the fragments. These late stages of contraction and the formation of hydrostatic protostellar cores are the subject matter of this paper.
The leftover magnetic flux at the end of the early, isothermal phase of contraction (n 10 10 cm −3 ) exceeds observed stellar magnetic fluxes by two to four orders of magnitude (Mouschovias 1987) and must therefore be redistributed or dissipated during the pre-main-sequence evolution.
Approximate theoretical estimates concluded that Ohmic dissipa-tion is effective in destroying the magnetic flux at densities above 10 10 cm −3 (Pneuman & Mitchell 1965;Nishi, Nakano, & Umebayashi 1991). Pneuman & Mitchell (1965) considered both ambipolar diffusion and Ohmic dissipation as magnetic flux-loss mechanisms and performed detailed calculations of the charged species abundances, accounting for ionization due to cosmic rays, radioactivity and thermal ionization. They estimated the relative importance of the flux-reduction mechanisms by comparing their timescales to the free-fall timescale, and concluded that Ohmic dissipation is more effective than ambipolar diffusion at densities 10 10 cm −3 . The field is recoupled to the matter when, at high enough temperatures, thermal ionization raises the degree of ionization sufficiently. They estimated the magnetic field of a newly born solar-mass star to be ∼ 40 Gauss.
Two decades later Nakano & Umebayashi (1986a,b) repeated these calculations using a more refined chemical model.
Through timescale comparisons for a series of different geometries (disk-like or spherical), they also found that the magnetic field decouples from the matter because of Ohmic dissipation. This decoupling took place at densities between 10 11 and 10 12 cm −3 .
In this series of papers, we pursue the self-initiated, ambipolar-diffusion-controlled star formation to densities ≈ 10 15 cm −3 , by which opaque hydrostatic protostellar cores form. It is in this density regime (after the detachment of ions and before thermal ionization becomes important and the increased degree of ionization reattaches the matter to the field) that the field may become completely detached from the matter. A significant contribution to the resolution of the magnetic flux problem of star formation is expected at this stage. Whether ambipolar diffusion or Ohmic dissipation is more effective in increasing the mass-to-flux ratio can only be determined through such detailed evolutionary calculations, not through comparisons of timescales.
In the present paper we formulate the problem in terms of the six-fluid MHD equations (neutral molecules; atomic and molecular ions; electrons; negative, positive and neutral grains). We do not assume that the magnetic flux is frozen in any fluid (in order to account for the possibility that all charges become detached from the magnetic field). We derive a new expression for the generalized Ohm's law, which lies at the heart of our approach to the problem at hand. This new formulation (1) includes the effects of inelastic coupling between grain species, and (2) explicitly distinguishes between flux advection, ambipolar diffusion, and Ohmic dissipation. We show how our formulation is related to and can be transformed into the traditional, "directional" formulation of the generalized Ohm's law, and we derive formulae for the perpendicular, parallel and Hall conductivities entering the latter, which include, for the first time, the effect of inelastic collisions between grains. In addition, we present a general (valid in any geometry) solution for the velocities of charged species as functions of the velocity of the neutrals and of the effective flux velocity (which can in turn be calculated from the dynamics of the system and Faraday's law). The last two sets of formulae can be adapted for use in any general MHD code to study phenomena beyond star formation in magnetic clouds.
Our results for a fiducial model cloud and a detailed parameter study are presented in two companion papers (Tassis & Mouschovias 2006b, c; hereafter papers II and III, respectively). As shown in paper III, at the high densities treated here memory of the mass-to-flux ratio of the parent cloud is lost, and the evolution of the contracting core is essentially independent of the initial conditions (for values within the parameter space determined by observations). This point is of particular importance, especially in light of the recent debate concerning the formation processes of molecular cloud cores. (For reviews and detailed discussion of the arguments on both sides of this debate, see, for example, MacLow & Klessen 2004;Tassis & Mouschovias 2004;Ballesteros-Paredes & Hartmann 2006;and Mouschovias, Tassis & Kunz 2006 and references therein). The need for a nonideal MHD treatment of the high-density regime of core evolution, when the magnetic flux problem of star formation has to be resolved, and hence the relevance of this work, is independent of the adopted theory of core formation.

MODEL CLOUD
We consider a self-gravitating, magnetic, weaklyionized model molecular cloud consisting of neutral particles (H 2 with 20% He by number), ions (both molecular HCO + and atomic Na + , Mg + ), electrons, neutral, singly negatively-charged and positively-charged grains. The abundances of charged species are determined from the chemical reaction network discussed in Tassis & Mouschovias (2005a, Appendix A). For a detailed discussion of the chemistry, see Ciolek & Mouschovias (1993, 1995, and Desch & Mouschovias (2001). Cosmic rays of energy 100 MeV, able to penetrate deep inside the cloud, drive the ionization for a wide range of densities. However, cosmic rays are appreciably attenuated for column densities 100 g cm −2 . At even higher densities, radioactive decays become the dominant source of ionization (Pneuman & Mitchell 1965;Umebayashi & Nakano 1981;Umebayashi 1983). The ionization rate ζ is the sum of the cosmic-ray ionization rate (ζ CR ) and the rate of ionization due to natural radioactivity (ζ RD ). The cosmic-ray ionization rate can be described (Umebayashi & Nakano 1980) by the relation: where σ is the column density of the gas, and ζ CR,0 is the unattenuated cosmic-ray ionization rate. Values of ζ CR,0 are estimated to lie in the range 10 −17 -10 −15 s −1 , with 5×10 −17 s −1 accepted as a canonical value (Spitzer 1978). The ionization rate due to natural radioactivity is dominated by the the isotope 40 K because of its long half-life of 1.25 Gyr. In the present calculations ζ RD = 6.9 × 10 −23 s −1 (Umebayashi & Nakano 1980). When the isotope 26 Al is considered, ζ RD = 1.94 × 10 −19 s −1 (Consolmagno & Jokipii 1978). UV radiation provides an additional ionization mechanism, but it only affects the envelope of molecular clouds (Ciolek & Mouschovias 1995). In this paper, UV radiation is not taken into consideration, as we concentrate our attention to the very dense, central parts of molecular clouds. For simplicity, we consider spherical grains of uniform size, since the chemical model used by Desch & Mouschovias (2001) accounted for grains of different sizes and showed that the effect of the grain size distribution on the evolution of the cloud is minimal. In the case of ion neutralization through collisions of ions (molecular or atomic) with grains, we assume that the ions do not get attached to the grains, but that they get neutralized and the resulting neutral particle escapes into the gas phase. Thus the total abundance of metals as well as the HCO abundance remain constant. Grain growth and depletion of elements by attachment onto grains are not considered here.
The axisymmetric model cloud is initially in an exact equilibrium state with self-gravity and external pressure being balanced by internal magnetic and thermalpressure forces. We follow the ambipolar-diffusioninitiated evolution from mean molecular cloud densities to a stage beyond the formation of a hydrostatic proto-stellar core. Isothermality is an excellent approximation for the early stages of star formation while the density is below 10 10 cm −3 ( Gaustad 1963;Hayashi 1966). The density at which isothermality breaks down depends on grain opacities. With a value of 0.039 cm 2 g −1 for spherical grains (Draine & Lee 1984), and 0.056 cm 2 g −1 for fluffy grains (Ossenkopf & Henning 1994), the required dust column density for the trapping of radiation is reached at densities ≥ 10 11 cm −3 . However, the optical depth is not the determining factor for the breakdown of isothermality (Masunaga & Inutsuka 1999). Instead, the rate of compressional heating must be compared to the rate of escape of thermal radiation. With these considerations, Desch & Mouschovias (2001) found that deviations from T = 10 K should occur at densities ≥ 3 × 10 10 cm −3 but that the temperature is not expected to rise rapidly, since a small temperature gradient suffices to allow the energy generated by the collapse to escape. To avoid the complicated and computationally taxing radiative transfer calculation, we switch from an isothermal to an adiabatic equation of state at a critical density which we treat as a free parameter.

STAR FORMATION
We extend the multifluid formalism described by Tassis & Mouschovias (2005) to relax the assumption of flux-freezing in the electrons and to account for the breakdown of isothermality at high densities. The equations governing the behavior of the six-fluid system (neutrals, electrons, ions, negative, positive, and neutral grains) are 3 : 3 A set of equations that comprise a logical entity, or a system of equations that need to be solved simultaneously, are numbered as a block; e.g., Na, Nb, Nc, etc..., where N is an integer.
The quantities ρ s , n s and v s denote the mass density, number density, and velocity of species s. The subscripts n, i, e, g − , g + and g 0 refer to the neutrals, ions, electrons, negatively-charged grains, positively-charged grains, and neutral grains, respectively. The symbols g, E and B denote the gravitational, electric and magnetic fields, respectively. The magnetic field satisfies the condition ∇ · B = 0 everywhere at all times. The definitions and derivations of F sn (the per unit volume frictional force on species s due to elastic collisions with neutrals) and F γδ,inel (the per unit volume force on grain fluid γ due to the conversion of dust particles of fluid δ into dust particles of fluid γ) are discussed in detail in §3 of Tassis & Mouschovias (2005).
In the force equations of the electrons, ions, and grains, the acceleration terms have been neglected because of the small inertia of these species; it was shown by Mouschovias, Paleologou, & Fiedler (1985) that the plasma reaches a terminal drift velocity very fast. Similarly, the thermal-pressure and gravitational forces (per unit volume) have been dropped from the force equations of all species other than the neutrals because they are negligible compared to the electromagnetic and collisional forces. The inelastic momentum transfer by the electron and ion fluids due to attachment onto grains and neutralization are negligible compared to the momentum transfer due to elastic collisions with neutrals, and they have been omitted from the force equations (5) and (6).
The equation of state is either isothermal or adiabatic depending on the density. For densities smaller than a critical density ρ opq , the isothermal equation of state is used P n = ρ n C 2 , where C = (k B T iso /µm H ) 1/2 is the isothermal speed of sound in the neutrals, k B the Boltzmann constant, µ the mean mass per particle in units of the atomic-hydrogen mass m H , and T iso = 10 K the temperature of the isothermal stages. For densities greater than ρ opq , adiabaticity is assumed and the equation of state becomes where γ is the adiabatic exponent. For simplicity we set to account for the rotational degrees of freedom of H 2 which only become accessible at high temperatures. The temperature T is calculated using Equation (39). For the early stages of star formation and for the isothermally infalling gas outside the optically thick inner region, the temperatures are much smaller than the grain sublimation temperature (≈ 1500 K), and a mass continuity equation (3) for the grains (charged + neutral) is appropriate.

DERIVATION OF A GENERALIZED OHM'S LAW
The evolution of the magnetic field in a fluid depends on the degree of coupling between the field and the fluid. If the magnetic field changes in time, an electric field is generated in the fluid, as described by Faraday's law of induction 4 (eq. [12]). This electric field is then capable of driving electric currents through the fluid, which in turn generate magnetic fields that tend to oppose the original change. If the conductivity of the fluid is infinite, there are no dissipative losses of the electric currents, and the magnetic field is said to be frozen in the fluid. The current generated in response to an electric field is described by Ohm's law. A generalized Ohm's law can be derived from the force equations of the plasma.
The force equations of all charged species and neutral grains, with the elastic and inelastic collision terms written out, are where are the relevant timescales for inelastic collisions between grains. In particular, τ − is the timescale for a neutral grain to be involved in any inelastic reaction involving the conversion of negative grains to neutral grains or vice versa, and τ + is the timescale for a neutral grain to participate in any inelastic reaction involving conversion between positive and neutral grains. Dividing each of the first four of the above equations by the mass of the corresponding species, multiplying by eτ sn (where s denotes the species), multiplying the ion and positive grain equations by -1, and adding up the 4 Note that the sequence of events described here does not imply causality. An equivalent picture, consistent with Maxwell's equations, would be that the magnetic field changes in response to a ∇ × E. resulting equations we obtain In equation (25), the third term vanishes because of charge neutrality, while the fourth term is the current density, j. We define the plasma conductivity σ p by where the individual contributions to σ p are Then equation (25) is rewritten as where (33) The quantity j 0 has units of current density and the corresponding term in the generalized Ohm's law (eq. [32]) represents the modification of this law because of the conversion of neutral grains to charged grains and vice versa through inelastic collisions.
Equation (32) is a different form of the generalized Ohm's law from the usual expression (Parks 1991), where E n is the electric field in the frame of the neutrals, and ⊥ denote directions parallel and perpendicular to the magnetic field,ê B is the unit vector parallel to the magnetic field, and η , η ⊥ and η H are the parallel, perpendicular, and Hall resistivities, respectively. Equation (32) can be transformed into equation (34), with the resistivities appropriately modified to account for the effects of the inelastic couplings between grains (see Appendix A). However, equation (32) offers several advantages in explicitly expressing the physical significance of each process, distinguishing between different mechanisms responsible for increasing the massto-flux ratio (through Faraday's law of induction) and simplifying the numerical implementation of the problem at hand. First, equation (32) calculates the electric field in the "lab" frame rather than in a frame comoving with the neutrals. Although from a mathematical point of view the two formulations are equivalent (the electric field in the lab frame is simply obtained by subtracting v n ×B/c from the electric field in the frame of the neutrals), the interpretation of the results is aided greatly by using the lab frame. In this frame, the neutrals are falling toward a center of gravity and the magnetic flux follows, at an infall speed equal to that of charged species attached to the field lines. Flux freezing would be a good approximation if the coupling between attached species and the neutrals were perfect, in which case the magnetic flux would fall in at the same rate as the neutrals. However, under the conditions present in the dense cores of molecular clouds, this coupling is imperfect, and the flux is "left behind", in a process known as ambipolar diffusion. After all species except the underabundant electrons (which interact weakly with the rest of the species) have detached, the inward dragging of the magnetic flux with the neutral fluid becomes ineffective and its advection speed decreases almost to zero, as shown in Paper II.
If this process is viewed in the frame comoving with the neutrals, it appears that the flux has to escape the fluid element by moving outwards. And at late times, when ambipolar diffusion is resurrected (Desch & Mouschovias 2001), it would be the responsibility of the extremely few electrons to transport this flux outwards -a process that appears improbable and counter-intuitive to some workers, who then suggest alternative mechanisms, such as grain drifts (Nakano, Nishi, & Umebayashi 2002). 5 Such drifts, however, cannot possibly be relevant because all charged species except the electrons have detached from the field lines by the densities of interest ( 10 11 cm −3 ), and the contribution of detached species to the drift velocity is insignificant.
Second, equation (34) expresses the components of the electric field parallel and perpendicular to the field lines, while equation (32) splits the electric field in a more physical way; when it is substituted in Faraday's law of induction (eq. [12]), it yields The first term on the right-hand side represents flux advection 6 (including the effect of collisions, i.e., ambipo-  (28) to derive their results. They used equations (7)-(13), which do account for the effects of all grain species (although not the effects of inelastic couplings between grain species, which were also not considered by Nakano, Nishi, & Umebayashi (2002), but are explicitly treated here. For the generalized Ohm's law including inelastic coupling between grains, see Appendix A). Equation (28) in Desch & Mouschovias (2001) is a simplification presented to assist the physical interpretation of the results and is not used in any of the numerical simulations. 6 Here we distinguish between flux advection due to the finite ve-lar diffusion), and the other two terms represent flux dissipation due to disruption of currents by collisions (i.e., resistivity). The last two terms can be combined into a single term if one defines an "effective" current j eff = j + j 0 . In an infinitely conducting medium, flux is simply transported inwards with a velocity v p , which is the conductivity-weighted velocity of the charged species. Third, if one uses the electron velocity to advect the magnetic flux rather than the ion velocity, then the Hall term in the generalized Ohm's law is absorbed in the ambipolar diffusion term (see Mouschovias 1987, eqs. [72a] and [76]). Here, we use a similar approach, which aims at simplifying the generalized Ohm's law by taking the advection of the magnetic flux to occur with the conductivity-weighted plasma velocity. In this form, the components which are responsible for flux dissipation and flux advection are clearly separated. Ohmic losses are represented by the j eff term, while flux advection is represented by the v p term 7 . The effect of ambipolar diffusion is then quantified by the difference between v p and v n . When all species detach (in which case v p becomes equal to v n ), ambipolar diffusion ceases to operate.

APPROXIMATION
Numerical simulations of the ambipolar-diffusioninduced evolution of axisymmetric, isothermal, molecular clouds have shown that force balance is rapidly established along field lines and is maintained throughout the evolution of the model cloud, even during the dynamic phase of contraction perpendicular to the field lines that follows the formation of a thermally and magnetically supercritical core (Fiedler & Mouschovias 1992Desch & Mouschovias 2001). Following the approach of Ciolek & Mouschovias (1993), we take advantage of these results and model the cloud as an axisymmetric thin (but not infinitesimally thin) disk, with its axis of symmetry aligned with the z-axis of a cylindrical polar coordinate system (r, φ, z) and with force balance along the field lines maintained at all times. The sense in which the thin-disk approximation is used is that the radial variation of any physical quantity is small over a distance comparable to the local disk half-thickness. The upper and lower surfaces of the disk are at z = +Z(r, t) and z = −Z(r, t), while the radius R of the disk satisfies the condition R ≫ Z. The magnetic field is poloidal and has the configuration described in Tassis & Mouschovias (2005, eqs. [17a] and [17b]). The model cloud is embedded in an external medium of constant thermal pressure P ext .
For computational simplicity, we reduce the mathematical dimensionality of the problem by integrating the equations analytically over z, assuming z-independence of all physical quantities inside the model cloud. This is an excellent approximation since thermal-pressure forces smooth out density gradients over lengthscales smaller locities of charged species and any flux redistribution due to Ohmic dissipation. The first term on the right hand side of equation (35) refers only to the former effect, not to Ohmic dissipation. 7 An effective flux advection velocity, representing the apparent velocity of magnetic field lines when Ohmic losses have also been accounted for, can also be defined and is equal to (c/B 2 )E × B. This effective flux velocity is our v f (see eq. [45]) and should not be confused with vp. than the critical thermal lengthscale (Mouschovias 1991), which is comparable to Z(r, t). A detailed derivation of the thin-disk equations under the assumption of isothermality and flux freezing in the ion fluid is given in Ciolek & Mouschovias (1993). Here we extend the formalism to account for breakdown of isothermality and of flux-freezing in the electrons.
The continuity equations for the neutrals and the grains in the disk are obtained through integration of equations (2) and (3) over z, and have the form given in Tassis & Mouschovias (2005, eqs. [18a-c]). In the isothermal regime, integration of the z-component of the force equation of the neutrals (eq. [4]), with the requirement of force balance along field lines, yields where P ext is the constant external pressure. In the adiabatic regime, the pressure becomes so the force balance along the z−axis yields in this case In the isothermal regime, equation (36) is used to calculate the density ρ n (r, t) from the column density σ n (r, t) and the external pressure. In the adiabatic regime, though, equation (38) must be solved simultaneously with the ρ n − T relation to obtain the density ρ n (r, t) and the temperature T (r, t) in the disk.
Integration over z of the r-component of the neutral force (eq. [4]), using the total (thermal-pressure plus Maxwell) stress tensor (Ciolek & Mouschovias 1993) to combine the thermal-pressure and magnetic-force terms, yields The expression for the total radial magnetic force (per unit area), includes the magnetic tension force (the B r terms) acting on the surfaces of the disk, and is formally identical to Tassis & Mouschovias (2005), equation (18g). The radial component of the thermal-pressure force in the isothermal regime is where and accounts for the contribution of the radial force exerted by the external pressure P ext on the upper and lower surfaces of the disk, since the surfaces are not horizontal. In the adiabatic regime, the thermal-pressure force becomes (see Appendix B for derivation) where The r-component of the gravitational field calculated from Poisson's equation for a thin disk has the same form as Tassis & Mouschovias (2005), equation (18h).
The assumption that the model cloud is embedded in a hot, tenuous, current-free medium implies that the magnetic field in the external medium is irrotational (∇ × B = 0). This allows the use of a scalar magnetic potential, Ψ, from which the radial component of the magnetic field at the top surface of the disk can be obtained (see Tassis & Mouschovias 2005 eq. [18i]).
In the thin-disk approximation the magnetic field is purely poloidal and from Ampère's law (eq. [10]) the current density must be purely toroidal, and thus they are orthogonal to each other. Consequently, the E-field has no component parallel to the magnetic field, as required by the generalized Ohm's law (eq. [32]). We may then define a new quantity, v f , which has units of velocity and can be used to describe the electric field through Substitution in Faraday's law of induction (eq. [12]) yields the evolution equation for the magnetic field In the plane of the disk (z = 0), equation (46) becomes Note that the velocity v f defines a frame of reference in which the magnetic flux is conserved. The velocities in the r and φ directions of all species, except the neutrals, and v f can be obtained from the solution of the steady-state force equations of the species and Ampère's law (eqs. [5]-[9] and [10]). The derivation is given in Appendix C. All the radial velocities have the where Θ s is the indirect attachment parameter of species s: for Θ s ≫ 1 the r-velocity of species s is approximately equal to that of the field lines, and the species is well attached to the magnetic field, while for Θ s ≪ 1 the species is detached and is falling in with the neutrals (v s ≈ v n ).

THE DIMENSIONLESS PROBLEM AND METHOD OF SOLUTION
We put the equations in dimensionless form by choosing units natural to the model molecular clouds. The unit of velocity is the isothermal sound speed in the neutral fluid, C; the units of column density σ c,ref and acceleration 2πGσ c,ref are those of the neutral column density and gravitational acceleration at the surface on the axis of symmetry of the reference state from which the initial equilibrium state is calculated. The unit of magnetic field is its uniform value in the reference state, B ref . Although not necessary, we define for computational convenience the unit of temperature as the temperature of the isothermal lower density parts of the model cloud T iso = 10 K. The temperature above which the rotational degrees of freedom of H 2 are excited is denoted by T rot .
We denote the dimensionless form of the cylindrical polar coordinates (r, z) by (ξ, ζ) and the dimensionless time by τ , so that the dimensionless thin-disk equations become Equations (49-n) and the equations for the velocities of all species (Appendix C) comprise the dimensionless system to be solved numerically. The dependent variables are now dimensionless but in this section only, for ease in identification, we have retained the same symbols as for their dimensional counterparts. The free parameters of the problem are: the initial central mass-to-flux ratio µ d,c0 = (2πG 1/2 σ c,ref )/B ref in units of its critical value for gravitational collapse in disk geometry, 1/(2πG 1/2 ) (Nakano & Nakamura 1978); the constant external pressure relative to the central gravitational "pressure" along the axis of symmetry of the reference stateP ext ≡ P ext /( π 2 Gσ c,ref ); and the density n opq , above which the equation of state becomes adiabatic.
The reaction rates of the chemical network and the value of the radius of the dust grains are given in Appendix A of Tassis & Mouschovias (2005). [The thermal (Jeans) length does not appear explicitly in the equations because of the choice of the units we are using: The unit of length C 2 /2πGσ c,ref is proportional to the Jeans length.] A control-volume formulation is employed for solving numerically the dimensionless equations that govern the evolution of the model cloud. The numerical method consists of an implicit time integrator; a conservative upwind advective difference scheme that incorporates the second-order accurate monotonicity algorithm of van Leer (1979); a second-order difference approximation of forces inside a mass zone; an integral approximation of the gravitational and magnetic fields; and an adaptive mesh capable of resolving the smallest natural lengthscale (i.e., thermal lengthscale λ T,crit ). Using convergence tests, we found that 85 zones in the radial direction are adequate to follow the evolution of the core to the densities of interest with an accuracy better than 1%. The algebraic equations for the equilibrium abundances of charged species are solved at each time step. A detailed description of the numerical methods is given by Morton, Mouschovias, & Ciolek (1994).
In order to follow shocks that might arise during the evolution, we change the originally central-difference scheme for the calculation of spatial derivatives to a second-order accurate one-sided difference scheme (Paleologou & Mouschovias 1983). In addition, tensor artificial viscosity (Tscharnuter & Winkler 1979) is introduced. The expression for the artificial viscosity term that is added to the RHS of the radial force equation (54) is given in Tassis & Mouschovias (2005).
Our initial and boundary conditions are the same as those discussed in detail in §5 of Tassis & Mouschovias (2005).

DISCUSSION
We have presented the formulation of the problem of the ambipolar-diffusion-controlled formation and evolution of magnetically supercritical molecular cloud cores inside initially magnetically subcritical parent clouds, and the core collapse to high densities, past the detachment of ions from the magnetic field and into the opaque regime, in which we approximated the treatment of radiation trapping using an adiabatic equation of state. We have used the six-fluid MHD equations to describe the model cloud. The magnetic flux was not assumed to be frozen in any of the charged species, in order to allow for the study of the possible detachment of the lightest charged species, the electrons, which would signify the complete decoupling of the magnetic field from the matter.
We have reduced the numerical complexity of the prob-lem by integrating the equations analytically along field lines, and we have used force-balance along the same direction to calculate the thickness of the formed magnetic disk.
We have derived a new expression for the generalized Ohm's law. It includes the effects of inelastic couplings between grains, which not only modify the expression for the contribution of each charged species to the total conductivity, but also account for the contribution of neutral grains to the conductivity (because neutral grains become charged part of the time). Moreover, this form of Ohm's law, when substituted in Faraday's law of induction, explicitly distinguishes between advection of the magnetic flux by charged species and Ohmic dissipation, which involves disruption of currents by collisions.
In two companion papers, we present the solution of the problem for a model cloud with fiducial parameters (Paper II) and a detailed parameter study (Paper III).
KT would like to thank Arieh Königl, Matt Kunz, and Vasiliki Pavlidou for useful discussions, and Glenn Ciolek for providing the base version of the code used in this work. The work of KT was supported in part by the University of Illinois through a Dissertation Completion Fellowship. This research was partially supported by a grant from the American Astronomical Society and NSF grants AST 02-06216 and AST 02-39759.

GENERALIZED OHM'S LAW: EQUIVALENCE OF DIFFERENT FORMULATIONS
In this Appendix we present the derivation of a generalized Ohm's law (in its usual directional formulation) from the six-fluid MHD force equations, accounting explicitly for the inelastic coupling between grain species. The inelastic grain coupling not only modifies the perpendicular, parallel, and Hall resistivities, but also introduces a second dissipative term associated with the neutral grain fluid.
In § 4 we presented an alternative formulation of the generalized Ohm's law, which is written in a maximally simplified functional form by determining the appropriate, conductivity-weighted plasma velocity, with which the magnetic flux is being advected. In this Appendix we discuss how the directional formulation of the generalized Ohm's law is related to this maximally simplified formulation. In addition, we present a general solution of the velocities of all charged species as a function of the velocity of the neutrals and an effective magnetic-flux velocity.

General Solution for the Velocities of Charged Species and Neutral Grains
We adopt a local cartesian coordinate system such that the x−axis is parallel to j × B, the y−axis is parallel to (j × B) × B, and the z−axis parallel to B. We can then write whereê x ,ê y andê z are the three orthonormal unit vectors. In addition, we define an effective magnetic-flux velocity v f by The physical meaning of this velocity is that it represents the velocity of the boundary of a constant-flux surface. The effects responsible for the motion of this boundary can be motion of the fluid element, flux redistribution between fluid elements, or flux dissipation due to Ohmic losses. With these definitions, we can write the force equations (5)-(9) in component form in the x − y plane as Once the x− and y− components of the velocity of the neutrals and the magnetic-flux velocity have been specified, the x− and y− velocities of the other species can be found by solving the linear, algebraic system of the force equations of the species (in steady state). In matrix form, the system of equations (A3) -(A12) is where A is the matrix of the coefficients of the unknowns We use Cramer's method to solve the above system. We define In addition, we use the notation D nx sx to represent the determinant of A where the column corresponding to the x−component of species s has been substituted by C nx , the notation D nx sy for the determinant derived from D if the y−component of species s column is replaced by C nx , and similarly for the other C vectors. Then, the solution of the system is, for the x−components of the velocities, and v sy = D nx for the y−components. The determinants defined above satisfy the relations where i and j can take one of the values x, y and δ ij is the Kronecker delta. These properties can be easily proved as follows. The determinant D ni sj + D fi sj is the same as the determinant which we get if we substitute the sj-column of D by C f i + C ni . If from this column we then subtract the 3 other j-columns, we get either the determinant of the coefficients D (if i = j), or a determinant which is identically zero (if i = j). We define With these definitions and the properties (A22), the solution can be re-written as for the x−components and v sy = 1 Θ sy + 1 v ny + Θ sy for the y−components.
Drift of Magnetic Flux Relative to the Neutrals as a Function of the Current Density Let us define w s = v s − v n and w f = v f − v n . Adding up the force equations for all species (including the neutral grains) we obtain the plasma force equation, Taking the x− and y−components of equation (A29) we find and s 1 (ω s τ sn ) 2 σ s σ p w sy = 0 .
Hence, equations (A30) and (A31) can be rewritten as and w fy For convenience and brevity, we define the coefficients of the above equations to be with i again being either x or y. Then, the solution of equations (A34) and (A35) with respect to w fx and w fy is Equivalence of Equation (32) and (34) for Ohm's Law Adding (v n /c) × B to both sides of equation (32) we obtain where E n is the electric field in a frame comoving with the neutrals. We rewrite the first term in the right-hand side of equation (A38) in terms of w fx and w fy , using the definition of v p (eq. [31]) to introduce the drift velocities w and equations (A32) and (A33) to eliminate w cx and w cy in favor of w fx and w fy : Note that k counts over charged species only, while s counts over all species (including neutral grains). Now substituting w fx and w fy using the solution (A37), and defining the sums above to be or, defining and we can write the advection term as The other two components of equation (A38) can be analyzed as and, with the definition σ g0 = e 2 n g0 τ gn /m g , If we now define and we can rewrite the j 0 term as Substituting equations (A44), (A45), and (A49) in equation (A38), we find that Note that the inelastic collisions between grains enter this expression both explicitly (through the new term j 0 and the new contributions to the perpendicular and Hall conductivities η 0⊥ and η 0H ) as well as implicitly, by altering the expressions for η adv,⊥ and η adv,H . In this Appendix we showed how a generalized Ohm's law can be derived from the six-fluid MHD force equations with the inelastic coupling between grains taken into account. We discussed two different formulations of the generalized Ohm's law. In the first case, the electric field is separated directionally, into a component parallel to the magnetic field, a component parallel to the component of the current which is perpendicular to the magnetic field, and a third component perpendicular to both magnetic field and current. In the second case, the contributions to the electric field are separated into an advection term and a term quantifying the contribution of each charge species to the Ohmic losses when substituted into Faraday's law of induction. The formulations are mathematically equivalent (and we have explicitly showed how one transforms into the other). Each offers different advantages for physical understanding and interpretation of the formalism and for application in numerical simulations.
We have presented new formule for all relevant resistivities for the case that inelastic grain couplings are taken into account. We have shown that the inelastic grain coupling not only modifies the perpendicular, parallel and Hall resistivities, but also introduces an extra dissipative term associated with the resistance to the motion of the neutral grain fluid. Finally, we have presented a general solution of the velocities of all charged species as a function of the velocity of the neutrals and an effective magnetic-flux velocity. The expressions presented in this paper can readily be used in non-ideal MHD codes for the study of weakly ionized astrophysical fluids. Depending on the problem and the code, one or the other form of the generalized Ohm's law can be used.

DERIVATION OF THE THERMAL -PRESSURE FORCE IN THE THIN-DISK APPROXIMATION
We define a volume element V of the disk as the cylinder-like volume that an arbitrary surface A in the midplane (z = 0) sweeps when it is translated vertically until it reaches the top z = Z(x, y) and bottom z = −Z(x, y) surfaces. The total force exerted on that volume element V in the i direction is or, using the divergence theorem where T is the thermal-pressure stress tensor, S is the bounding surface of the volume V , andn the local unit vector normal to the surface. This force can be expressed in terms of the area A on the midplane, perpendicular to the field lines: The integrand F A is the force per unit area on the midplane of the disk. We need to determine F A . We split the bounding surface S into three surfaces: top, bottom, and side. Then the contribution of each one to F can be calculated separately, and The where ∇ is the derivative in the direction parallel to the midplane, perpendicular to the field lines inside the disk. The surface element of the top/bottom surface is related to the surface element of the midplane area A as dS top/bottom = 1 + (∇ Z) 2 dA.
The intersection of the side surface with the midplane is a closed curve l. Thus the surface element of the side surface is and the normal unit vector is always on the midplane. So Similarly, and The contribution of each surface to the thermal-pressure force is The total thermal-pressure force in the adiabatic regime is Force balance in the z−direction implies so that Z = σ n 2ρ n = σ n C 2 T Tiso 2ρ n C 2 T Tiso = σ n C 2 T Tiso 2P ext + πGσ 2 n .
Then the pressure force is written in a compact form as In the r−direction, ∇ σ n = ∂σ n ∂r and ∇ T T iso = ∂ ∂r T T iso , so that the radial thermal-pressure force per unit area in the adiabatic regime is Clearly, if T = T iso the second term vanishes, and equation (B21) reduces to the expression (28a) of Ciolek & Mouschovias (1993) for the thermal-pressure force. If T > T iso , the temperature gradient contributes to the thermal-pressure force.
where σ p,g0 has units of conductivity and is defined as σ p,g0 = e 2 n g0 τ gn m g , and τ pn is a properly averaged collisional timescale between the plasma (charged species) and the neutrals and is defined by 1 τ pn = ρ e ρ n (ω e τ en ) 2 τ en + ρ i ρ n (ω i τ in ) 2 τ in + ρ g+ ρ n (ω g τ gn ) 2 τ gn + ρ g− ρ n (ω g τ gn ) 2 τ gn . (C20) Note that in equation (C18), σ n is the column density of the neutrals, not a conductivity. We solve the system of equations (C7) -(C17) in terms of v n and v f , and then substitute all velocities in equation (C18) to find v f in terms of v n . In matrix form, the system of equations (C7) -(C17) is where A is the matrix of the coefficients of the unknowns and