text:
# Radiation from the first forming stars

## 1 Introduction

As the temperature of the cosmic bath decreases, atoms start to recombine and therefore decouple from the CMB radiation at redshift $`1100`$. The baryonic Jeans mass immediately after this event is given by

$$M_J4.5\times 10^4\left(\frac{\mathrm{\Omega }_b}{\mathrm{\Omega }}\right)(\mathrm{\Omega }h^2)^{1/2}M_{},$$
(1)
where $`\mathrm{\Omega }`$ and $`\mathrm{\Omega }_b`$ are the matter and baryon density parameters, respectively. Masses larger than $`M_J`$ are gravitationally unstable and should, in principle, collapse. However, in order for the actual collapse to occur a more severe condition must be satisfied, i.e. that the cooling time of the gas is shorter than the Hubble time at that epoch. In fact, radiative losses provide the only way for the gas to lose pressure and settle down in the potential well of the host dark matter halo. Since the virial temperature corresponding to the masses of the first objects (PopIII) is typically $`\stackrel{<}{}8000`$ K, cooling by hydrogen Ly$`\alpha `$ excitation is strongly quenched, and the only viable coolant in a primordial H-He plasma is molecular hydrogen. H2 is produced during the recombination phase, but its relic abundance is very small ($`f_{\mathrm{H}_2}10^6`$). This primordial fraction is usually lower than the one required to fulfill the above time-scale constraint, but during the collapse phase the molecular hydrogen abundance reaches values that are high enough to allow the collapse to continue. It follows that the fate of a virialized lump depends crucially on its ability to rapidly increase its H2 content during the collapse phase.

Haiman, Thoul & Loeb (1996) and Tegmark et al. (1997) have addressed this question by calculating the evolution of the H2 abundance for different halo masses and initial conditions for a standard CDM cosmology. They conclude that if the prevailing conditions are such that a molecular hydrogen fraction of order of $`f_{\mathrm{H}_2}5\times 10^4`$ can be produced during the collapse, then the lump will cool, fragment and eventually form stars. This criterion is met only by the largest PopIIIs, implying that at each virialization redshift one can define a critical mass, $`M_{crit}`$, such that only protogalaxies with total mass $`M>M_{crit}`$ will be able to eventually form stars. As an example, a 3$`\sigma `$ fluctuation of a Cold Dark Matter primordial spectrum, has $`M_{crit}10^6M_{}`$ and collapses at $`z30`$. Other studies have revised the value of $`M_{crit}`$ as a result of a more accurate treatment of the micro-physics (Nishi & Susa 1999; Fuller & Couchman 2000; Machacek, Bryan & Abel 2001); however, these results do not qualitatively alter the basic argument given above.

Recently, some works have started to tackle the formation and collapse of PopIII objects through numerical simulations (Abel et al. 1998; Nakamura & Umemura 1999; Nakamura & Umemura 2001; Bromm, Coppi & Larson 1999, 2001; Abel, Bryan & Norman 2000 - ABN00; Bromm, Coppi & Larson 2001 - BCL01; Bromm et al. 2001) based on hierarchical scenarios of structure formation. These studies have shown that gravitational collapse induces fragmentation of the first pregalactic objects with initial baryonic mass $`10^5M_{}`$ into smaller clumps of typical mass of about $`10^{23}M_{}`$, which corresponds to the Jeans mass set by molecular hydrogen cooling.

Tracking the subsequent evolution of these clumps is a very challenging problem as it requires the simultaneous solution of the hydrodynamic equations and of radiative transfer in lines. At present, due to the tremendous computational demand, this as only been possible in one dimensional formulations of the problem. Star formation in the early universe naively appears easier to understand as compared to the present, because several complicating effects can be neglected to a first approximation: among these are magnetic fields, dust grains and metal enrichment.

In their pioneering work, Omukai & Nishi (1998) (ON98) showed that at a certain stage of the collapse a quasi-hydrostatic central core is formed whose specific mass is almost independent on the details of the problem, with a typical value of $`5\times 10^3M_{}`$. However, as the accretion rate of the infalling gas is found to be rather high, around $`10^2M_{}`$ yr-1, the core mass at its formation might be unimportant in setting the final mass of the star. In the absence of any effect quenching accretion, a large fraction of the initial object can become soon part of the protostar. Pushing this conclusion a bit further, one might predict a top heavy IMF for this first generation of sources.

If the infall can instead be stopped, different scenarios can be envisaged. It has yet to be firmly established if standard ways proposed to halt the infall continue to work under primordial conditions: radiation force might be opacity limited; bipolar flows need some MHD acceleration process and therefore seem to be excluded by the weak magnetic field present at the main formation epoch (Gnedin, Ferrara & Zweibel 2000), although the magnetic field could be amplified during the protostar formation process. Angular momentum barrier related to disk formation and competitive accretion among different clumps are slightly more promising, but much more difficult to model accurately. The estimate of the first mechanism (radiation reaction) requires to accurately derive the properties of the radiation field produced during the collapse and to evaluate the effects of the radiation force on the infalling envelope. This is one of the main aims of the present paper.

The second motivation for our study consists in assessing the detection chances of distant PopIII objects via the (infrared) radiation they release during the collapse phase. Ciardi & Ferrara (2001) have already pointed out that IR H2 molecular lines produced by the cooling of supernova shocked PopIII gas will be very likely observed by the Next Generation Space Telescope. However, as H2 roto-vibrational lines are essentially the only carrier of the entire gravitational energy release, it seems worthwhile to compute the intensity, shape and evolution of the emitted spectrum in detail. Recently, Kamaya & Silk (2001) (KS01) have estimated the H2 luminosity from primordial collapsing stars using a very approximated treatment.

To investigate the described problem we have performed 1D numerical simulations including a detailed treatment of the molecular line radiative transfer. The study is germane to the one of ON98 but it includes the radiative force term in the dynamical equations adding a number of imporvement, as we outline throughout the paper.

The paper is organized as follows. In §2, we present our numerical methodology. The choice and justification for the initial conditions, and the results of the simulation are discussed in §3. In §4 we discuss the luminosities and the detectability, and in the final section we give a brief summary.

## 2 Numerical methods: the code

We developed a 1-D Lagrangian hydrodynamical code, including full treatment of radiative transfer and chemical evolution (with careful treatment of the H2 molecule, which is the main coolant in the initial phases of the collapse). Such code is similar to the one described by ON98. Here we only summarize the relevant differences between the two; a more detailed description is given in the Appendix.

### 2.1 Dynamics: radiative force

We included the radiative force term into the momentum equation, which then becomes

$$\frac{dv}{dt}=4\pi r^2\frac{dp}{dM_r}\frac{GM_r}{r^2}+f_{\mathrm{rad}},$$
(2)
where we are exploiting the spherical symmetry, and where

$$f_{\mathrm{rad}}=\frac{1}{c}\kappa _\nu F_\nu 𝑑\nu $$
(3)
is the radiative force per unit mass, $`F_\nu `$ is the specific energy flux at frequency $`\nu `$ (both $`f_{\mathrm{rad}}`$ and $`F_\nu `$ are considered in the radial outward direction), and $`\kappa _\nu `$ is the opacity at the frequency $`\nu `$.

### 2.2 Chemistry and thermodynamics: non ideal effects

As discussed by several authors (e.g.  Hummer & Mihalas, 1988; Mihalas, Däppen & Hummer, 1988; Saumon, Chabrier & Van Horn, 1995; hereafter SCVH95), at high densities ($`\rho \stackrel{>}{}10^3\mathrm{g}\mathrm{cm}^3`$) particles are jammed closely together, and the bound electron orbitals filling too large a volume fail to survive so that electrons migrate from atom to atom. Such an effect is known as *pressure ionization* and generally it is not taken into account in the Saha equations. As an example, a pure hydrogen gas with $`\rho =1`$ g cm-3 and $`T=3\times 10^4`$ K is largely ionized (the mean distance between nuclei is only twice the Bohr radius of an H atom on the fundamental level, and the electrons are free in a Fermi gas, interacting simultaneously with several ions), while solving the Saha equations results in a completely molecular gas, which is clearly unphysical. In such regimes, the gas is also non-ideal, and the perfect gas equations are poor approximations.

For these reasons, at high densities and temperatures we find the chemical abundances and thermodynamical properties (such as the pressure $`p`$, the internal energy per unit mass $`u`$ and the entropy per unit mass $`s`$) by interpolation of the values tabulated by SCVH95, which are calculated using the *free energy minimization* technique. SCVH95 calculations include several non-ideal effects, such as, (a) the partial electron degeneracy, (b) the Coulomb interactions among charged particles, (c) the finite extension of “hard-sphere” particles, (d) weak attractive interactions (e.g. van der Waals), and (e) modifications of the internal partition functions due to the shift of the energy levels.

## 3 The simulations

### 3.1 Initial conditions

The main parameters specifying the initial condition of our runs are: the total cloud mass $`M_{\mathrm{tot}}`$, the temperature, density, and velocity profiles $`T_0(M_r)`$, $`\rho _0(M_r)`$, $`v_0(M_r)`$, respectively, and the hydrogen molecular fraction $`f_{\mathrm{H}_2}`$.

Several previous investigations set constraints on such parameters. Theoretical predictions of the mass of fragments leading to the birth of the first stars are quite uncertain, ranging from a minimum of $`0.1M_{}`$ (Palla, Salpeter & Stahler 1983), to $`\stackrel{>}{}200M_{}`$ (Hutchins 1976). On the other hand, there is agreement about the value of $`f_{\mathrm{H}_2}5\times 10^4`$, and about the value of $`T_0(M_r)`$ which is of order of few hundreds degrees (e.g. Tegmark et al. 1997). Detailed numerical simulations (BCL01, ABN00) suggest that the runaway collapse starts soon after the (number) density reaches a value of $`10^410^5`$ cm-3. Such value is determined by the transition from non-LTE cooling (with $`\mathrm{\Lambda }n^2`$) to LTE cooling (with $`\mathrm{\Lambda }n`$), slowing down the contraction, and allowing the gas to reach a quasi-hydrostatic state. This transition occurs when $`T300400`$ K and $`f_{\mathrm{H}_2}510\times 10^4`$, inside a cloud with a mass of a few thousand $`M_{}`$, a value close to the Jeans mass of the cloud.

To check how our results depend on the choice of the initial conditions, we run a set of models with different initial mass $`M_{\mathrm{tot}}`$, ranging from $`10`$ to $`1500M_{}`$. We use a constant profile for the H2 abundance and zero velocity (the cloud is at rest). The density and temperature profiles inside the cloud are derived imposing hydrostatic equilibrium and a polytropic index equal to 5/3 or 1 (isothermal). In addition we have considered ‘numerical’ profiles as given by ABN00 from cosmological simulations. In runs P10, due to the high value of the initial central density, we used a higher value for $`f_{\mathrm{H}_2}`$. The initial electron fraction adopted is $`f_e^{}=10^{10}`$. With the exception of run P100/2, we always divide our simulated cloud into 100 Lagrangian shells, with massive external shells and much lighter internal ones (e.g. , in run P100 each external shell comprises $`\stackrel{>}{}1M_{}`$ of gas, while the central one has a mass of only $`3\times 10^7M_{}`$). Such choice is aimed to follow with good resolution the fast evolution of the centre of the cloud, which accounts only for a tiny fraction of the total mass. Details of the initial conditions and general set up are given in Table 1. Note that runs P100 and P10 are directly comparable to those of ON98. Run P100/2 is identical to run P100, but the object is divided into 140 (rather than 100) shells; the “extra” shells were used to improve the resolution in the particularly important region where the shock forms ($`3\stackrel{<}{}\mathrm{log}M_r\stackrel{<}{}1`$).

### 3.2 Dynamical evolution

The dynamics of the collapse is illustrated in Figure 1 that shows the evolution of the cloud in run P100/2 at 10 different stages (see Table 2). In Figures 2 and 3 we show the evolution of the temperature and of the main hydrogen species in the most central shell as a function of density. In the rest of this section we will describe the results obtained in run P100/2; results of other runs are quite similar. In runs P100, P100/2 and P10 we find good agreement with the results of ON98.

In the beginning, the evolution is very well described by a Larson-Penston self-similar solution (Larson 1969 - L69; Penston 1969) with a density profile scaling as $`\rho r^\alpha `$ (with $`\alpha 2.2`$). The self-similar evolution breaks down between stages 6 and 7, when a quasi-hydrostatic core forms in the centre and starts growing in mass.

#### 3.2.1 The collapse phase

During collapse the central regions of the cloud undergo several important chemical and thermal transformations. When the number density reaches $`\stackrel{>}{}10^8`$ cm-3(stage 1), 3-body reactions become efficient, converting the bulk of H into H2 (Fig. 1, stages 2-3; Fig. 3). The growing H2 fraction increases the gas cooling rate, only partially counterbalanced by radiative transfer effects (see next section 3.3). The evolution of the thermal properties of the gas is shown in Figure 4, where the adiabatic heating $`\mathrm{\Gamma }_{\mathrm{ad}}`$, the radiative cooling rate $`\mathrm{\Lambda }_{\mathrm{rad}}`$ and the effective heating $`\mathrm{\Gamma }_{\mathrm{eff}}=\mathrm{\Gamma }_{\mathrm{ad}}+\mathrm{\Gamma }_{\mathrm{ch}}`$ ($`\mathrm{\Gamma }_{\mathrm{ch}}`$ is the heating - or cooling - due to chemical reactions) are plotted vs. the gas temperature. All these quantities refer to the central shell. As $`T\stackrel{<}{}1550`$ K, the chemical heating due to the H2 formation dominates over $`\mathrm{\Gamma }_{\mathrm{ad}}`$ and heating and cooling terms almost balance, leading to a very slow temperature growth. For $`T\stackrel{>}{}1550`$ K H2 is destroyed rather than formed, and the chemical term becomes a net cooling (this is the reason why the curves $`\mathrm{\Gamma }_{\mathrm{eff}}`$ and $`\mathrm{\Gamma }_{\mathrm{ad}}`$ cross at $`T1550`$ K). At this stage the chemical and radiative cooling terms are of the same order. Radiative cooling is dominated by line radiation for $`T\stackrel{<}{}1800`$ K, while continuum emission is important above. However continuum radiative cooling soon becomes highly inefficient because of the large optical depth: for $`T\stackrel{>}{}2350`$ K the innermost regions become almost adiabatic, and compressional heating is only matched by chemical cooling.

We just mentioned that H2 starts to dissociate when $`T\stackrel{>}{}1550`$ K. This process is slow because of the density increase and the chemical cooling associated to the dissociation process (4.48 eV per dissociated molecule are lost by the gas thermal energy). As a result, the H2 fraction at the end of the self-similar phase ($`T_\mathrm{c}15000`$ K, stage 6) is still $`5\%`$, even in the hottest central regions.

#### 3.2.2 Core formation

During the last stages of the self-similar collapse (stages 5-6; evolution from point A to point B in Fig. 2), when the centre of the cloud is essentially adiabatic, H2 dissociation acts as a “thermostate”, preventing a faster increase in temperature and pressure, which would halt the collapse. In fact, dissociation absorbs $`\stackrel{>}{}50`$% of the heat generated by gas compression, which (due to the opaqueness of the cloud’s centre) could not have been dissipated otherwise. When nearly all H2 molecules are dissociated, the thermostatic effect is no longer effective, and the central temperature and pressure rise rapidly stopping the collapse of the central regions. At the end of H2 dissociation ($`T_\mathrm{c}20000`$ K) the chemical energy and the thermal energy per H0 atom are similar ($`2.2`$ eV and $`2.6`$ eV, respectively). This fact reflects the virial theorem: half of the cloud gravitational energy is converted into thermal energy after compression, while the rest is tapped into chemical energy through dissociation.

Outside the core, the gas is still collapsing with a velocity close to free-fall. This results in a shock front forming at the edge of the core, where the infalling material is suddenly halted (see Fig. 1.b). Unfortunately, the numerical method (i.e. the introduction of an artificial viscosity) used for the treatment of hydrodynamical effects is inaccurate when dealing with strong shocks, since they are spread upon several shells. However, though the treatment of the shock physics is only approximate, the general picture of the cloud evolution is essentially correct. In run P100/2 we also improve the situation by refining the shells in the region where the shock forms.

The core evolution is driven by two main processes. First, as can be seen in Figures 1.b, and 6, the core mass (defined as the mass inside the radius where the infall velocity is $`<0.1v_{\mathrm{ff}}`$, $`v_{\mathrm{ff}}\sqrt{2GM_r/r}`$) grows, due to the accretion of infalling material. Second, even if the infall velocity is drastically reduced, the core is compressed by the infalling material, and then heats up (as can be seen in Fig. 2, beyond point B).

#### 3.2.3 Cloud structure after core formation

After core formation and during its subsequent evolution, the H2 fraction in the centre never rises again to important levels111The central value of $`f_{\mathrm{H}_2}`$ at late stages is quite uncertain: using SCVH95 results we obtain $`f_{\mathrm{H}_2}0.03`$; had we used results such as those of Mihalas, Däppen & Hummer 1988, we would have obtained $`f_{\mathrm{H}_2}<0.01`$; ON98 find $`f_{\mathrm{H}_2}0.4`$, but this is likely to be too high since they do not consider pressure ionization effects.. H2 tends to form a layer of several solar masses between the small core (where H is progressively ionized) and the outer layers (where H is essentially atomic, see Figure 5). After stage 6 we can recognize three different regions, i.e., the core, the molecular envelope, and the outer region (see Figure 5).

* The central core is a region whose infall velocity is much smaller than the free-fall value. At core formation, the gas here is basically neutral, but at stage 9 it has become mostly ionized (but $`n_{\mathrm{H}^+}/n_\mathrm{H}`$ is still only $`2`$) because of the pressure ionization process. The outer edge of the core is marked by a shock front, where material accretes nearly with free-fall velocity. Since the artificial viscosity technique “smooths” the shock, we have a thin transition region, where neutral H is the dominant species; that may be an numerical artefact.
* Outside the shock defining the core outer edge, there is an envelope of $`3M_{}`$ composed by molecular hydrogen (the shock also marks a sharp transition in chemical composition, from H2 to H0 or H+). Such envelope is essentially in free fall ($`v/v_{\mathrm{ff}}`$ is in the range 0.3-0.7, decreasing outward). Its outer edge is roughly at a radius $`\stackrel{>}{}10^{15}`$ cm.
* The remaining outer gas, accounting for most of the cloud mass, has a chemical composition only slightly different from the initial one, since reactions leading to H2 formation are slow (mainly because of the low density), and involve only a minor fraction of H atoms. This layer is the extension to the inner molecular envelope, and has a velocity $`v/v_{\mathrm{ff}}`$ between 0.1 and 0.3.

At stage 9 the core temperature and density are $`10^5`$ K and $`0.3`$ g cm-3, respectively. In these conditions the Courant time-scale is extremely short ($`\stackrel{<}{}300`$ sec.), and it becomes more and more difficult (and costly, in terms of CPU time) to follow the evolution further. For this reason, we stopped our calculations when the central temperature reached $`10^5`$ K.

#### 3.2.4 The accretion phase

After core formation, the evolution continues because of the accretion of the outer envelope onto the core.

Despite our approximate treatment of the shock, we find that the accretion rate, as inferred from the full calculation, is extremely high: the core grows from $`0.3\times 10^2M_{}`$ when it forms, to $`4.5\times 10^2M_{}`$ in $`0.30`$ yr, which means $`\dot{M}0.14M_{}`$ yr-1. Such accretion rate is much larger than typical estimates for present-day star formation ($`\dot{M}\stackrel{<}{}10^4M_{}`$ yr-1).

In order to extend our results to a more advanced state of accretion, we modified our code “freezing” the evolution of the shells deeply embedded in the core. In practice, when a shell has i) a small infall velocity ($`v/v_{\mathrm{ff}}<10^3`$), ii) a high temperature ($`>5\times 10^4`$ K), and iii) is located at least 5 shells inside the shell where the infall velocity is maximum, we completely stop its evolution (i.e., we set velocity at 0 and keep all the other quantities constant to the last computed values).

This greatly speeds-up our calculations because there is no more need to follow the shells with the shortest time-scales, drastically reducing the number of integration steps needed. Such “freezing” is appropriate, since these shells do evolve on time-scales much longer than those imposed by the Courant condition and have little influence upon the layers where accretion occurs. There are two major concerns:

* the neglect of the secular changes in the core radius: for example, in the full calculation the radius of the 10 innermost shells undergoes a reduction of about $`40`$% (from $`1.6\times 10^{10}`$ cm to $`1.0\times 10^{10}`$ cm) as $`T_\mathrm{c}`$ goes from $`50000`$ K (stage 8) to $`10^5`$ K (stage 9), while in the simplified calculation performed between stages 8 and 9, we keep this radius constant. This change is less important in more external shells: the difference in the shock position between the “full” and the “frozen” run at stage 9 is found to be $``$4%.
* the radiation produced upon the shock surface could ionize a thin shell outside the shock itself, as reported in Staher, Palla & Salpeter (1986) (SPS86), modifying the dynamics of accretion (see next sections).

Instead, the influence of the lowered temperatures in central shells is completely negligible, since the optical depth from the “frozen” shells is very high ($`\tau _\mathrm{c}\stackrel{>}{}10^{12}`$).

Figure 6 shows the result of the accretion calculations performed with the modified code: we stop $`\stackrel{>}{}10`$ years after core formation. In the first 0.34 years this results are almost identical to those obtained with the “complete” code; instead, points at later times should be considered indicative.

ON98 investigated the evolution of the mass accretion rate on the newly formed hydrostatic core using the Larson-Penston self-similar solution (L69, Penston 1969) and its extensions after core formation (Hunter 1977, Yahil 1983, Suto & Silk 1988). For a model close to our model P100 (and P100/2), they found that the core mass grows following a power-law

$$M_{}(t)=0.11M_{}((tt_{\mathrm{cf}})/\mathrm{yr})^{0.73}.$$
(4)
We also find an empirical model which could better fit the numerical results, based on the fact that in our calculations we find that the mass flux $`M_{\mathrm{flux}}(M_r)=4\pi r^2\rho v_{\mathrm{in}}`$ is well approximated by a power-law. Such empirical model leads us to predict that

$$M_{}(t)=9.2\times 10^2M_{}\left(\frac{t}{\mathrm{yr}}+0.0063\right)^{0.657}.$$
(5)
Both these analytical expressions are compared with the numerical results in Figure 6.

#### 3.2.5 Importance of initial conditions

The above description essentially applies to all the models we computed: the evolutionary sequence is always the same, even if there are small quantitative differences. Models lose memory of the initial conditions almost completely during the self-similar phase, as can be seen in Fig. 7, where tracks referring to the different models are very close for $`n_\mathrm{c}>10^{14}`$ cm-3, i.e. well before the cloud becomes optically thick (point A) or forms the central hydrostatic core (point B).

There remain only relatively small differences (within a factor of 2 or 3), in particular:

* Models with higher total mass ($`M_{\mathrm{tot}}`$), have higher temperatures than low $`M_{\mathrm{tot}}`$ models, for any given central density $`n_\mathrm{c}`$;
* Models with high $`M_{\mathrm{tot}}`$ have cores that grow faster than low $`M_{\mathrm{tot}}`$ models (see Fig. 6).

Both these trends can be easily explained: when we impose that in the initial stage the clouds are close to equilibrium, we implicitly force high-mass models to have higher temperatures at any given density; even if quantitative differences are strongly reduced in the self-similar phase, high-mass models remain warmer than low-mass models at all stages. Furthermore, the accretion rate depends on temperature, since $`\dot{M}M_{\mathrm{Jeans}}(T,\rho )/t_{\mathrm{ff}}(\rho )T^{3/2}`$ (see Stahler, Shu & Taam, 1980) and temperature differences explain also those in $`\dot{M}`$.

The behaviour of the N1500 model (the only one for which we do not impose quasi-equilibrium initial conditions) is particularly relevant: after some “wandering” in the $`n_\mathrm{c}T_\mathrm{c}`$ plane, its central conditions evolve almost exactly along the evolutionary track of model P100, indicating that the collapse does not preserve the memory of initial conditions and the cloud evolves along an equilibrium track.

### 3.3 The radiation field

The evolution of the radiation field is intimately coupled with dynamics and our physical description of forming primordial stars fully accounts for their emission properties. Here we give details on the evolution of the optical depth, luminosity and spectra during the two phases of collapse and core formation, and at the very beginning of the accretion phase.

#### 3.3.1 Luminosity from the collapsing cloud

At the onset of dynamical evolution, radiation results from the cooling of the optically thin H2 present in the cloud. At stage 1, as 3-body reactions keep converting H atoms into H2 molecules, the line optical depth rapidly grows, reaching significant values. At stage 2 there are already 14 lines with $`\tau >0.1`$ (although none of them with $`\tau =1`$), and this number rises to 42 lines (half with $`\tau >1`$) by the time the central temperature reaches $`1000`$ K (between stage 2 and stage 3). In Figure 8a we show, as an illustrative example, the evolution of the optical depth profile for a selected line during the self-similar phase.

The increase of the optical depth obviously affects the gas cooling properties (as can be seen from the flattening of the radiative cooling curve in Fig. 4), but it does not influence the dynamics in a fundamental way. The increase in H2 abundance is associated with an increase of the gas temperature and of the number of relevant transitions of the H2 molecule. As the temperature grows, new, high quantum number lines are initially optically thin, providing a minimum level of cooling. Though line cooling is severely quenched by opacity (e.g., at $`T_\mathrm{c}2000`$ K only a fraction $`0.02`$ of the line flux is not re-absorbed locally) its absolute value *increases* anyhow as $`T_\mathrm{c}`$ increases; this prevents the temperature from rising too fast, possibly stopping the collapse in its early phases.

As already mentioned, when $`T_\mathrm{c}=1000`$ K we start to include cooling from the continuum, which rapidly gains importance, becoming the dominant process in the central $`0.1M_{}`$ already when $`T_\mathrm{c}1800`$ K. The evolution of the continuum opacity profile is shown in Figure 8b, where the formation of the core sharp edge is apparent.

The corresponding emitted luminosities vs. time are plotted in Figure 9. Note that even when continuum radiation dominates the central emission, the H2 lines still have a comparable or slightly larger overall luminosity. This is due to the fact that line radiation comes from the outer, massive H2 envelope, while continuum comes from an inner, much smaller region (only slightly larger than the core itself); this can be seen in fig. 10, where we show the radial dependence of the emitted luminosity; in particular, note that the drop in the central line luminosity has negligible influence upon the total line luminosity, and that the growth of continuum luminosity stops at $`M_r10^210^1M_{}`$, while line luminosity is still growing at $`M_r1M_{}`$.

#### 3.3.2 Spectra

In Fig. 11, we show the emitted spectra at 3 different stages. In the upper panel the spectrum comes from line emission only, while in the lower panels we assumed that the continuum spectrum is a black-body at the effective temperature of the central shell (1500 K - panel b) and of the continuum photosphere (i.e.  of the region where the radial outward optical depth is $`1`$). Such temperature remains roughly constant at $`23002400`$ K from the birth of the continuum photosphere to the end of our calculations.

In appendix B we give a full listing of the luminosities of all the calculated lines, plus the integrated continuum luminosity at 4 stages (the 3 depicted in fig. 11 plus the stage when $`T_\mathrm{c}=1000\mathrm{K}`$). We note that the durations of the various stages are very different: line luminosities resemble those listed in the $`L_{650}`$ column (corresponding to fig. 11.a) for $`10^4\mathrm{yr}`$, while they resemble those given in columns $`L_{1000}`$, $`L_{1500}`$ and $`L_{\tau _\mathrm{c}=1}`$ only for $`300\mathrm{yr}`$, $`40\mathrm{yr}`$ and $`2\mathrm{yr}`$, respectively222The last duration ($`2\mathrm{yr}`$) is actually only a lower limit, since line emission is still at high levels when our calculations stop; however, the H2 envelope where line emission originates is not likely to last for $`\stackrel{>}{}100\mathrm{yr}`$ after core formation, since it could be accreted onto the core, or photodissociated because of the outgoing radiation..

We note that the spectra are produced only in the central collapsed region: line emission comes essentially from the molecular envelope, while continuum emission is produced close to the surface of the hydrostatic core. In practice, the whole protostellar luminosity comes from the central $`3M_{}`$. The structure of this central region is quite similar in all the models, since it is largely determined by the self-similar collapse, rather than by initial conditions. For this reason, the spectra emitted by different clouds during these phases should be similar.

#### 3.3.3 Radiation pressure effects

We have found that the role of radiative force is negligible during all the evolution we followed, e.g.  at least until the first phases of accretion (our stage 9). In a test run we have removed the radiative term in the momentum equation and re-calculated evolution. We then compared the results with those of the “full” code, and found that the difference is completely negligible.

Returning to the early phase of accretion, despite the very high luminosity ($`10^{36}\mathrm{erg}\mathrm{s}^1`$ in the continuum, plus other $`10^{36}\mathrm{erg}\mathrm{s}^1`$ in H2 lines), which is comparable to the Eddington luminosity $`L_{\mathrm{Edd}}1.25\times 10^{36}[M/(0.01M_{})]\mathrm{erg}\mathrm{s}^1`$, the opacity of the infalling gas is much smaller than the electron scattering opacity ($`\kappa _{\mathrm{es}}=\sigma _\mathrm{T}/m_\mathrm{H}0.4\mathrm{cm}^2\mathrm{g}^1`$), since the infalling material is essentially molecular, with a very low ionized fraction, and other opacity sources (Rayleigh scattering, bound-free and free-free transitions, collision-induced absorption of H2 ) give small contributions. We find that the opacity is everywhere less than $`0.05\kappa _{\mathrm{es}}`$ at all times before the formation of the hydrostatic core; typical values are actually several order of magnitude less than this upper limits. The same is true also after the formation of the hydrostatic core, if we consider only the accreting layers.

#### 3.3.4 Bridging collapse models with accretion models

SPS86 and recently Omukai & Palla (2001) studied the subsequent phase of accretion by modeling the hydrostatic core as a star of increasing mass, surrounded by an infalling envelope whose structure was determined by an assumed steady mass accretion rate $`\dot{M_{}}=4.4\times 10^3M_{}\mathrm{yr}^1`$. They found that the shock surface is inside an opticaly thick ionized region, whose outer edge is roughly coincident with the photosphere. They also found that between the shock front and the photosphere the radiative force term is important, since $`\kappa \kappa _{\mathrm{es}}`$ and the luminosity is not far from the the Eddington value. Despite this fact, both works find that the radiative force is not able to stop the accretion before the core mass grows to very high values ($`>100M_{}`$).

We find that also in our models the shock front is inside the optically thick region, but the gas between the shock and the photosphere is not ionized, mainly because the density in the region is relatively high ($`10^6\mathrm{g}\mathrm{cm}^3`$), about $`1000`$ times larger than the value that can be estimated from Fig. 12 of SPS86; instead, we find that the temperature is only slightly lower than in the SPS86 model. The ionization fraction inside such region is low and remains such ($`\stackrel{<}{}0.05`$) even if we make the hypotesis that radiation coming from the shock heats the gas inside the photospheric shells, so that gas and radiation temperature are the same ($`12000\mathrm{K}`$, the post-shock temperature found both by ours and SPS86 calculations). So, this gas still has a low opacity, and the radiative force is small.

Although numerical effects could play some role in this discrepancy, the most likely explanation is that the models of SPS86 are admittedly indicative for core masses $`\stackrel{<}{}0.1M_{}`$, i.e.  during the “decay of transients” phase, when the arbitrarily chosen initial conditions are still important; instead in our models the core mass is always in this range (here we do not consider the simplified calculations aimed to the study the growth of the core mass). It is quite possible that both results are correct and that a transition occurs after the last stage of our calculations. For instance, we note that the difference in the density of the photospheric layers (which keeps the ionization fraction at low levels) is mainly due to the large difference (a factor $`100`$) between the accretion rate we find and that assumed by SPS86; since in our models we also find that the mass accretion rate is decreasing, such difference will reduce.

### 3.4 Comparison with metal-rich protostars

In this section we discuss the similarities and differences between the PopIII (high redshift, zero metal) star formation process we addressed and the PopI-PopII (low redshift, solar metallicity) star formation.

#### 3.4.1 Hydrostatic core

In our work (and also in that of ON98) the birth of the quasi-hydrostatic core is quite different to that described by L69 (and more recently by Masunaga, Miyama & Inutsuka, 1998 and by Masunaga & Inutsuka, 2000) for a present-day protostar because of the very different conditions. In fact, those authors found that a molecular core forms first, and only later on, H2 dissociation leads to the formation of a second more internal core, which is now ionized and at much higher temperature and density, as can be seen in table 3. The two cores coexist for a short time, while the second core “absorbs” the first one.

Instead in our models we find a single core, completely different from the first (molecular) core of metal-rich models, but closely resembling the second one.

The reason of both the differences and the similarities between the two cases is essentially explained by an argument due to Omukai (2000). He found that the locus of points in the $`nT`$ plane where $`\tau 1`$, $`MM_{\mathrm{Jeans}}(T,n)`$ and $`\mathrm{\Lambda }_{\mathrm{rad}}\mathrm{\Gamma }_{\mathrm{ad}}`$ (cfr. Silk, 1977) is approximately described by a line (L) with $`Tn^{2/5}`$, which is shown in Fig. 2; on the right of L, a core with mass $`M_{\mathrm{Jeans}}(T,n)`$ is optically thick.

Our metal free models become optically thick at point A (remarkably close to L), where $`T_\mathrm{c}`$ is high and H2 dissociation is immediately effective, preventing the formation of a hydrostatic core. Instead, metallic protostars ($`Z\stackrel{>}{}10^6Z_{}`$) reach L at temperatures too low for H2 dissociation, and the first (molecular) core forms. Such core evolves adiabatically, i.e. along a $`Tn^{\gamma _{\mathrm{ad}}1}`$ line, which remains close to L, because $`\gamma _{\mathrm{ad}}\gamma _{\mathrm{ad},\mathrm{H}_2}7/5`$. In this way cores of all metallicities converge towards the vicinity of point A. After point A, the evolution of the central regions is essentially determined by the balance between compressional heating and H2 dissociation cooling and no more depends on metal content; thus, metal rich and metal free models form a core when H2 is exhausted, close to point B.

The initial core mass $`3\times 10^3M_{}`$ (only slightly dependent upon initial conditions) agrees with $`M_{\mathrm{Jeans}}(T_\mathrm{B},n_\mathrm{B})`$ (where subscript $`B`$ denotes the values at point B). Such mass is also close to the minimum value reached by $`M_{\mathrm{Jeans}}(T_\mathrm{c},n_\mathrm{c})`$ during the cloud evolution (see Fig. 12). This minimum Jeans mass is generally interpreted as the minimum fragment mass (i.e. the mass below which further fragmentation is impossible) and we find it depends very weakly upon initial conditions and on metallicity. An important corollary is that the minimum fragment mass can not be simply found imposing $`\tau \stackrel{>}{}1`$, since we have seen that such condition, although necessary, is not sufficient for the formation of a stable hydrostatic core.

It is also interesting to note that in our models there is no sign of instability inside the hydrostatic core due, for example, to H ionization, even if the core is mostly atomic at stage 7 and mostlyly ionized at stage 9.

#### 3.4.2 Accretion rate

Even if important core properties, such as the initial mass and the evolutionary path of central regions in the $`nT`$ plane are almost independent of metal content, the mass accretion rates are completely different, with metal-free objects accreting at a rate $`\stackrel{>}{}10^3`$ times larger than their metal enriched counterparts. For example, in our P100 model core mass reaches $`0.3M_{}`$ in $`6`$ years (see fig. 6) while in the model of Masunaga & Inutsuka (2000) about $`20000`$ years are required for the accretion of the same mass (see their Fig. 7).

As already discussed by several authors (e.g. Stahler et al. 1980; SPS86), the reason for this difference is found in the initial conditions, in particular in the higher temperature, since for an object with $`MM_{\mathrm{Jeans}}`$ we have

$$\dot{M}(k/m_H)^{3/2}T^{3/2}/G1.8\times 10^7T^{3/2}M_{}\mathrm{yr}^1.$$
(6)
If we take a typical initial temperature of our models $`300`$ K we obtain $`\dot{M}10^3M_{}\mathrm{yr}^1`$, but if we take $`T5`$ K, as appropriate for present-day protostellar clouds, we get a value $`500`$ times smaller for $`\dot{M}`$.

We note that (6) gives the *mean* accretion rate for the *complete* accretion of the protostellar cloud onto the core; should we consider accretion of only a fraction of the cloud, equation (6) would give higher values of $`\dot{M}`$ because the appropriate reference temperature would be higher.

#### 3.4.3 Radiation feedback

In the previous sections we showed that radiation pressure is not effective in stopping the accretion flow upon a zero-metallicity protostellar core. The case of metal rich protostars is different: luminosities are smaller because of smaller accretion rates, but the presence of dust brings the opacity at much higher levels. For example, Stahler, Shu & Taam (1981) estimate that at the start of accretion radiative deceleration retards the inflow of 42% from the free fall value, even if it can not completely stop accretion (see also the models of Masunaga et al. 1998 and Masunaga & Inutsuka 2000, which include radiation pressure effects in their momentum equation). In a different context, radiative forces upon dust were also found to act decisively in halting the accretion inflow upon massive stars (Wolfire & Cassinelli, 1987; Larson & Starrfield 1971).

For sufficiently high density and temperature, ignition of nuclear reactions occurs in the core interior. To check if the effect can, in principle, stop the accretion phase, we compare the accretion time-scale with the diffusion time-scale of the energy produced.

As can be seen in Fig. 8, the interior of the core has optical depth $`\tau _{}\stackrel{>}{}10^{12}`$. If energy is transported radiatively, the time-scale for its diffusion through the core is $`t_{\mathrm{diff}}\tau _{}R_{}/c10^5\mathrm{yr}[\tau _{}/10^{12}][R_{}/(10^{11}\mathrm{cm})]`$, where the core radius $`R_{}`$ has a value of the order of $`5\times 10^{11}`$ cm. We note that, as $`T_\mathrm{c}`$ rises, the opacity $`\kappa `$ inside the core decreases; anyway, $`\kappa \stackrel{>}{}\kappa _{\mathrm{es}}`$, which leads to a strong lower limit $`\tau _{}\stackrel{>}{}10^{10}`$.

An identical calculation for present-day protostellar cores leads to a time-scale of the same order of magnitude. Instead, the accretion times differ enormously. For example, if we extrapolate equation (5), we find that our core would take only $`t_140`$ years after formation to reach a mass of $`M_{}1M_{}`$ (but we note that $`t_1250\mathrm{yr}`$ also in the models of SPS86), while in the model by Masunaga & Inutsuka (2000) more than $`t_11.4\times 10^5`$ years are required for the core to reach the same mass.

So, $`t_1/t_{\mathrm{diff}}1`$ in metal-rich cores, but $`t_1/t_{\mathrm{diff}}1`$ in the metal-free case. Accordingly, in the metal-free case the energy generated in the core centre (for instance, because of ignition of nuclear reactions) cannot affect the accretion, while in the metal-rich case such energy could be able to reach the surface and change the conditions of the infalling envelope.

The arguments discussed in the present section are subject to important uncertainties. For example, we did not investigate any energy transport mechanisms, such as convection, different from radiative diffusion, which could be much more efficient (e.g. , see the “luminosity wave” described in SPS86).

## 4 Line flux and detectability

The work by KS01 has addresses the important problem of detectability of the emission from collapsing primordial protostar. They employed approximate analytical calculations (partially derived from the ON98 work) in order to estimate order-of-magnitude luminosities of collapsing objects, assuming a central temperature $`T1000\mathrm{K}`$.

In table 4 we list the luminosities of 5 important lines, at all the four stages we tabulated in table 6, plus the total line luminosity, the continuum luminosity, and their sum. They are compared to the values predicted by KS01.

It can be seen that KS01 predict a total luminosity comparable to our results for the initial phases of collapse; instead, in the last phases our total luminosities are $`\stackrel{>}{}10`$ times larger than the KS01 value. This is quite reasonable, since the approximations they used are much better at the beginning of the collapse rather than at the end. If we consider the luminosities of single lines, a more important discrepancy emerges. KS01 correctly predict that the rotational transition at $`\lambda =9.69\mu \mathrm{m}`$ is the most important one at the beginning, but they estimate a very high line luminosity, that no lines actually achieve, even in advanced stages. The reason is that they attribute too large a fraction of the total luminosity to single lines. Instead, the emission is actually split into a growing number of lines: in the early phases, the two rotational lines listed in table 4 accounts for $`50\%`$ of total emission (first tabulated stage), but the luminosities of single lines “saturate” when their central optical depth becomes $`\stackrel{>}{}1`$. As the total line emission rises, we have mainly an increase in the number of important lines, rather than a growth of the luminosity of single lines: at the last stage no line accounts for more than $`5\%`$ of total line emission.

In table 4 we also give the expected specific fluxes coming from a single protostellar object at a redshift $`z=19`$ in an SCDM cosmological model at the first and the, for the first and the last listed stage. It can be seen that the most promising lines are those produced by purely rotational transitions, such as those at $`\lambda =9.69\mu \mathrm{m}`$ and $`\lambda =8.06\mu \mathrm{m}`$ (which are redshifted to $`\lambda _{\mathrm{obs}}190`$ and $`\lambda _{\mathrm{obs}}160\mu \mathrm{m}`$, respectively), whose expected specific flux is $`\stackrel{>}{}10\mathrm{nJy}`$ for at least $`10^4\mathrm{yr}`$); vibrational transitions, such as those at $`\lambda =2.41\mu \mathrm{m}`$ and $`\lambda =2.81\mu \mathrm{m}`$ (redshifted to $`\lambda _{\mathrm{obs}}48\mu \mathrm{m}`$ and $`\lambda _{\mathrm{obs}}56\mu \mathrm{m}`$) reach specific fluxes of the same order only in the last phases and their emission at such levels should last only for $`\stackrel{<}{}10^3\mathrm{yr}`$. Finally, we predict that the continuum specific flux should be of the order of $`0.002\mathrm{nJy}`$ at an observed wavelength of $`\lambda _{obs}40\mu \mathrm{m}`$ (from a rest-frame wavelength of $`2\mu \mathrm{m}`$); we cannot predict the duration of continuum emission, since it is deeply affected by the evolution of the accretion and by the evolution of the hydrostatic core at the centre.

Such fluxes are largely below the detection threshold of any observational facility that will be available in the near future. NGST sensitivity for continuum radiation should be of the order of $`1\mu \mathrm{Jy}`$ at $`\lambda =20\mu \mathrm{m}`$, while ALMA should reach a sensitivity of $`10\mathrm{mJy}`$ for a line $`25\mathrm{km}\mathrm{s}^1`$ wide at $`\lambda 350\mu \mathrm{m}`$. However, the possibility that primordial protostars form in large groups must not to be excluded, since a galactic-scale ($`10^{11}M_{}`$) fluctuation could host a large number of $`10^6M_{}`$ halos, and each of them could host several protostellar objects. For instance KS01 roughly estimate that a galactic-scale fluctuation could have a luminosity $`10^4`$ times larger than a single object, and that such luminosity could last for $`10^7\mathrm{yr}`$. If this is correct, specific fluxes from galactic-scale objects should be of the order of $`0.01\mu \mathrm{Jy}`$ and $`0.1\mathrm{mJy}`$ for continuum and line emission, respectively.

## 5 Conclusions

In this work we have investigated the evolution of the luminosity and spectra of the radiation emitted by primordial collapsing protostars. We find that initially radiation comes mostly from H2 rotational lines. Later, the contributions from H2 vibrational lines and also from continuum radiation (due to the accretion upon a centrally hydrostatic core) become important. At early stages, only two rotational lines account for about 50% of the total emission, but the subsequent increase in total line luminosity results in an increase of the number of prominent lines, rather than a large increase in the luminosity of the most important ones. As a result, the expected line specific fluxes are always $`\stackrel{<}{}50\mathrm{nJy}`$ (with line widths $`\stackrel{<}{}30\mathrm{km}\mathrm{s}^1`$), even for the most important lines, which should fall at observed wavelengths in the range $`150200\mu \mathrm{m}`$. Continuum radiation should reach a specific flux of the order of $`0.002\mathrm{nJ}`$, peaking at an observed wavelength between $`20`$ and $`40\mu \mathrm{m}`$. These fluxes should be quite universal, provided that a molecular central region of $`3M_{}`$ develops at the centre of each protostellar cloud. The predicted fluxes, both for line and continuum radiation, are largely below the detection threshold of existing or planned instruments (such as NGST and ALMA), but the possibility that primordial protostars form large “clusters” remains to be investigated, even if rough estimates do not look very promising (KS01).

In our work we also found that, different from the present-day case, radiation forces play virtually no role in the initial evolution of a protostellar object (i.e. , during the self-similar collapse, the formation of the core and the onset of accretion). Although we were not able to completely bridge the gap with accretion models such as those of SPS86 and Omukai & Palla 2001, we have shown that their initial assumptions are quite reasonable. Even if there are differences in the ionization structure, there is indication (in particular, the decrease in the mass accretion rate) that these will vanish as the hydrostatic core mass grows. Thus, undisturbed spherical collapse should proceed until the initial protostellar cloud has been largely accreted, or until radiation forces halt the infall at $`M_{}300M_{}`$ (Omukai & Palla 2001).

What can stop mass growth is the presence of HD molecules (even in a tiny fraction) that may result from shocks (Flower & Pineau des Forêts 2001, Uehara & Inutsuka 2000); in this case the fragmentation of the PopIII object might end up in clumps of roughly solar rather than thousand solar masses, thus drastically changing the final mass of the star.

If this process is not important and angular momentum and magnetic fields do not play a significant role either, the final mass of a primordial protostar is primarily set by the amount of gas available in the parent cloud. Then, it is conceivable that the primordial IMF is skewed towards high masses, and that very massive stars (VMS), with masses in excess of $`M100M_{}`$, can form. Such stars may end-up as a population of intermediate mass black holes (Portinari et al. 1998, Heger & Langer 2000) that represent the end product of such an episode of pregalactic star formation. Then, these seeds could be incorporated into the larger units as seen in the bulges of present-day galaxies (Madau & Rees 2001).

We thank an anonymous referee for a critical reading and of the manuscript and for stimulating remarks. We also thank B. Ciardi for discussion and for making her chemistry routines available to us, and K. Omukai, R. Nishi and F. Palla for discussion and clarifying remarks about their works.

## Appendix A The code

Hereafter we describe the main features of the code we developed.

### A.1 Generalities

#### A.1.1 Gas hydrodynamics

Despite the cosmological context of our calculations, the collapsing baryon clouds we are studying are assumed to be self-gravitating, and we can neglect the gravitational effects of host dark matter halo.

The code integrates the following 1-D Lagrangian hydrodynamical equations (conservation of mass, momentum and energy) using an explicit finite-difference scheme:

$$\frac{dM_r}{dr}=4\pi r^2\rho $$
(7)
$$\frac{dv}{dt}=4\pi r^2\frac{dp}{dM_r}\frac{GM_r}{r^2}+f_{\mathrm{rad}}$$
(8)
$$\frac{du}{dt}=p\frac{d}{dt}\left(\frac{1}{\rho }\right)+\frac{\mathrm{\Gamma }\mathrm{\Lambda }}{\rho },$$
(9)
where the Lagrangian coordinate $`M_r`$ is the mass inside a radius $`r`$. $`\rho `$, $`v`$, $`p`$, $`f_{\mathrm{rad}}`$, $`u`$, represent the density, velocity, pressure, radiative force per unit mass and internal energy per unit mass, respectively. $`\mathrm{\Gamma }`$ and $`\mathrm{\Lambda }`$ are the heating and cooling rate per unit volume, respectively.

Hydrodynamical effects are included by adding to the pressure an artificial viscosity term of the form

$$q\rho \mathrm{\Delta }v^2,$$
(10)
in regions where contiguous gas shells would compenetrate (e.g. Bowers & Wilson 1991; Thoul & Weinberg 1995; hereafter, TW95).

#### A.1.2 Chemistry

Because of the extremely low metallicity of the primordial material as predicted by standard big-bang nucleosynthesis ($`Z10^{12}`$; Kolb & Turner 1990), we can safely limit our treatment of chemistry to the various ionic species of atomic hydrogen (H0, H+, H-), molecular hydrogen (H2 , H$`{}_{}{}^{+}{}_{2}{}^{}`$), and helium (He0, He+, He++), plus free electrons (e-).

The He number abundance is taken to be 1/12 of the total H number abundance; this translates into mass fractions of X=0.75, Y=0.25 and Z=0.

Some authors (Flower 2000; Flower & Pineau de Forêts 2001; and also, in a different context, Uehara & Inutsuka 2000) suggested that deuterium is an important coolant in a metal-free gas, because of the permanent dipole moment of HD molecules and their low excitation temperature ($`160`$ K for HD instead of $`510`$ K for H2 ). However HD cooling is likely to be important only at temperatures below some hundred degrees, i.e. during the formation of protostellar clouds, rather than during their collapse. For this reason, we neglect the effects of deuterium. Note that BCL01 included HD cooling in their investigation of primordial cloud fragmentation, and found that (for a standard Deuterium abundance) the thermal evolution of the gas would not change significatively.

Initial conditions were chosen with chemical abundances far from equilibrium, and we follow their evolution by integrating a network of reactions (see Table 1). In particular, we include the 3-body reactions described by Palla, Salpeter & Stahler (1983).

As the gas gets denser ($`\rho \stackrel{>}{}\rho _{\mathrm{eq}}5\times 10^{11}`$ g cm-3), the abundances of the relevant species (H2 , H0 and He0, which then represent more than 99.9% of the total number abundance) are close to their equilibrium values, as reactions. We then switch to the equilibrium chemistry.

When the gas is almost ideal, i.e.  when $`p10^4\mathrm{dyne}\mathrm{cm}^2`$ ($`p`$ is the gas pressure) we find the equilibrium abundances by solving the Saha equations for the relevant species. However these equations become inappropriate at high temperatures and densities, so when $`p10^6\mathrm{dyne}\mathrm{cm}^2`$ we interpolate the abundance values from the tables given by SCVH95. In the intermediate regime ($`10^4\mathrm{dyne}\mathrm{cm}^3p10^6\mathrm{dyne}\mathrm{cm}^2`$), we evaluate the chemical composition both with Saha equations and using the SCVH95 tables, and then we appropriately interpolate the two results; this is done in order to ensure a smooth transition between the two regimes.

#### A.1.3 Thermodynamical properties

For the solution of the hydrodynamical equations, we need to specify the equation of state and to calculate the internal energy as a function of temperature, density and chemical composition.

At low densities and temperatures (in practice, we use $`p10^4\mathrm{dyne}\mathrm{cm}^2`$ as a criterion), the gas is ideal to a very good level of approximation and we can use the usual equation

$$p=nkT$$
(11)
where $`n`$ is the number density of the gas.

In this ideal regime we can also write the total internal energy per unit mass as

$$u=u_{\mathrm{th}}+u_{\mathrm{chem}}$$
(12)
with

$$u_{\mathrm{th}}=\frac{1}{\gamma _{\mathrm{ad}}1}N_{\mathrm{tot}}kT,$$
(13)
$$u_{\mathrm{chem}}=N_{\mathrm{tot}}\underset{i=1}{\overset{9}{}}y_i\chi _i$$
(14)
where $`u_{\mathrm{th}}`$ and $`u_{\mathrm{chem}}`$ are the thermal and the chemical energy per unit mass, while $`N_{\mathrm{tot}}`$, $`y_i`$ and $`\chi _i`$ are the total number of particles per unit mass, the number fraction of particles of species $`i`$ ($`y_i=N_i/N_{\mathrm{tot}}`$), and the binding energy of a particle of species $`i`$, respectively. Last, $`\gamma _{\mathrm{ad}}`$ is the mean adiabatic index of the gas, obtained by a weighted average of the adiabatic indexes of the various chemical species: for the seven included mono-atomic species we have $`\gamma _{\mathrm{ad},\mathrm{X}}=5/3`$, while for the two diatomic species we use $`\gamma _{\mathrm{ad},\mathrm{H}_2^+}=7/5`$ and the expression (from ON98) is

$$\frac{1}{\gamma _{ad,H_2}1}=\frac{5}{2}+x^2\frac{e^x}{(e^x1)^2},x=\frac{6100K}{T}.$$
(15)
The value of $`\gamma `$ for H$`{}_{}{}^{+}{}_{2}{}^{}`$ neglects the vibrational degrees of freedom, but the abundance of this ion is always very small and this does not introduce any significant error.

When the gas is no more ideal ($`p10^6\mathrm{dyne}\mathrm{cm}^2`$), we find $`p(\rho ,T)`$ and $`u(\rho ,T)`$ by interpolating the values tabulated by SCVH95, which have been calculated jeeping non-ideal effects into account. In the transition between the two regimes, we interpolate the results of the ideal and of the SCVH95 values in order to obtain a smooth transition.

#### A.1.4 Radiative transfer and cooling

In the initial stages of collapse, gas cooling is dominated by the H2 roto-vibrational transitions, which soon become optically thick. A careful treatment of radiative transfer effects is then crucial in the study of the collapse.

H2 line cooling is important only when $`T\stackrel{<}{}2500`$ K, as at higher temperatures H2 molecules slowly dissociate, and continuum emission starts dominating the local cooling. The energy gaps between H2 rotational levels typically correspond to excitation temperatures of few$`\times 100`$ K - $`10^3`$ K, while the gap between vibrational levels corresponds to $`6000`$ K. Under such conditions, it is adequate to consider only H2 roto-vibrational levels with rotational quantum number $`0J20`$, and vibrational quantum number $`0v2`$. In this way, all the states whose excitation temperatures (relative to the ground state) is $`\stackrel{<}{}20000`$ K are considered. We treat the radiative transfer of all the lines (between the considered levels) permitted by the quadrupole selection rules.

Line opacities are computed using the Einstein A-coefficients from Turner, Kirby-Docken & Dalgarno (1977), and assuming local thermodynamic equilibrium (LTE) for estimating the level populations. This is correct since the density at which a line gives a relevant contribution to the total cooling is always greater than the critical density at which collisional de-excitations dominate upon radiative transitions.

Finally, the effects of thermal broadening and Doppler shifts due to bulk motions arising in the collapse phase have been included in the line transfer calculations.

When the temperature rises to a sufficiently high value ($`T>1000`$ K), we include the cooling by continuum radiation; continuum radiative transfer is treated as a gray problem, using the mean Planck opacity for $`Z=0`$ gas from Lenzuni, Chernoff & Salpeter (1991) when $`1000K<T<7000`$ K, and the Rosseland mean opacity for metal-free gas from Rogers & Iglesias (1992) when $`T>7000`$ K.

### A.2 Numerical details

#### A.2.1 Implementation of the basic equations

The numerical treatment of the dynamical evolution is made along the lines described by Bowers & Wilson (1991), and TW95.

With the exception of the energy equation, the dynamical equations listed in §2.1 were discretized as in TW95 (their equations 9, 10, 11, 12, 13 and 15), assuming a value of 4 for the artificial viscosity parameter $`c_q`$.

Such Lagrangian finite difference scheme is second order accurate in space and time.

For the energy equation (9), our treatment starts from equation (14) of TW95. The chemical term in the internal energy (see §2.3) introduces a complex dependence upon the final temperature. For this reason the energy equation is solved iteratively, using the bisection method. In this case the result is only first order accurate in time.

#### A.2.2 Time-steps

The integration time-step is chosen evaluating the minimum among:

* twice the previous integration time-step;
* the minimum dynamical time-scale over all shells $`i`$ ($`\mathrm{min}_i\{C_{\mathrm{dyn}}(Gr_i^3/M_{r,i})^{1/2}\}`$);
* the minimum Courant time-scale over all shells $`i`$ ($`\mathrm{min}_i\{C_{\mathrm{Cour}}\mathrm{\Delta }r_i/c_{\mathrm{s},i}\}`$);
* the minimum cooling time-scale over all shells $`i`$ ($`\mathrm{min}_i\{C_{\mathrm{cool}}u_i\rho /(\mathrm{\Gamma }_i\mathrm{\Lambda }_i)\}`$);
* the minimum crossing time-scale over all shells $`i`$ ($`\mathrm{min}_i\{C_{\mathrm{cross}}\mathrm{\Delta }r_i/\mathrm{\Delta }v_i\}`$);

In addition, we check *a posteriori* that the temperature does not change by more than 0.5% in a single time-step. The factors $`C_{\mathrm{dyn}}=1`$, $`C_{\mathrm{Cour}}=0.2`$, $`C_{\mathrm{cool}}=0.1`$ and $`C_{\mathrm{cross}}=0.05`$ are safety constants.

#### A.2.3 Boundary conditions

As internal boundary condition we simply require that the centre of our object is a geometrical point ($`r=0`$, $`M_r=M_0=0`$), with zero velocity ($`v_{\mathrm{in}}=0`$).

The external boundary condition is set by taking a constant pressure outside the most external shell.

#### A.2.4 Implementation of radiative transfer

General scheme

The solution of the radiative transfer equation in spherical symmetry requires to find the specific intensity $`I_\nu (r,\mu )`$, where $`\mu `$ is the cosine of the angle $`\theta `$ between the radial outward direction and the direction of the considered ray.

This is done by means of the standard tangent ray method (Hummer & Rybicki 1971; Bowers & Wilson 1991; ON98). In this scheme, the specific intensity $`I_\nu (r,\mu )`$ is evaluated at the external radius $`r_i`$ of each shell for $`1+2i`$ values of $`\mu `$.

Once $`I(r_i,\mu _j)`$ is known, we compute the monochromatic flux and the net outward bolometric luminosity $`L(r_i)`$ at each shell radius.

The cooling rate (per unit time per unit volume) of the shell $`i`$ is then given by

$$\mathrm{\Lambda }_{rad,i}=\frac{L(r_i)L(r_{i1})}{(4\pi /3)(r_i^3r_{i1}^3)}$$
(16)
H2 lines

To calculate the source function and the absorption coefficient, we assume LTE. This is a fair assumption, as the density is always much greater than the relevant critical values.

The absorption coefficient is given by

$$\alpha _\nu =\frac{h\nu }{4\pi }\varphi (\nu )\left(\frac{2h\nu ^3}{c^2}A_{21}\right)\left(n_1\frac{g_2}{g_1}n_2\right)$$
(17)
where $`\varphi (\nu )`$ is the line profile, $`A_{21}`$ is the Einstein coefficient for the relevant transition (taken from Turner et al. 1977), $`g_2`$ and $`g_1`$ are the statistical weights of the two levels involved in the transition, and $`n_2`$ and $`n_1`$ are number densities of molecules in the two considered levels, computed assuming LTE.

As noted by ON98, since the temperature typical of H2 line emission is quite low, the line profiles are narrow, and the Doppler shift caused by the gas bulk infall can be important. We assume a Gaussian profile for each line, sampled in 30 frequency bins, and we properly take into account Doppler shift in the line transfer. If the line absorption depth is $`<0.1`$, we use the optically thin approximation.
Continuum

The continuum transfer is treated as a gray problem. Again, we assume LTE to estimate the source function. The absorption coefficient is calculated from the (Planck) opacity tables of Lenzuni et al. 1991 when $`T<7000`$ K, and from the (Rosseland) opacity tables of Rogers & Iglesias 1992 for $`T>7000`$ K.

As the continuum optical depth can reach values $`1`$ (see fig. 8), the numerical treatment above can become inaccurate. For this reason, when $`\tau >4`$ in a single shell (see the conditions listed in ON98), we use the diffusion approximation, i.e.

$$F=\frac{16\sigma T^3}{3\alpha }\frac{T}{r}\frac{16\sigma T_i^3}{3\alpha _i}\frac{T_{i+1}T_{i1}}{r_{i+1}r_{i1}}$$
(18)
where $`\sigma `$ is the Stefan-Boltzmann constant.

## Appendix B List of emitted lines

In the following table we list the luminosity at infinity of each H2 line we included into our code. Lines are sorted according to their wavelength.

The columns give, respectively:

* the quantum levels describing the transition;
* the wavelength (in $`\mu \mathrm{m}`$) of the line;
* the radial optical depth from the centre of the simulated object to infinity, as evaluated at line centre, and at the epoch when the continumm optical depth is $`\tau _\mathrm{c}=1`$;
* the line total luminosity at infinity, evaluated at four different epochs, i.e.  when the central temperature is 650 K ($`L_{650}`$), 1000 K ($`L_{1000}`$) and 1500 K ($`L_{1500}`$), and when $`\tau _\mathrm{c}=1`$ ($`L_{\tau _\mathrm{c}=1}`$); such luminosities are given in $`\mathrm{erg}\mathrm{s}^1`$;
* the specific luminosity (in $`\mathrm{erg}\mathrm{s}^1\mathrm{Hz}^1`$) at the stage when $`\tau _\mathrm{c}=1`$, estimated using a line width equal to the Doppler broadening $`\mathrm{\Delta }\nu =(\nu _0/c)(2kT/\mu m_\mathrm{H})^{1/2}`$, where $`\nu _0`$ is the line central frequency, $`\mu 2`$ is the H2 molecular weight and $`T`$ is the temperature of the shell where the radial optical depth to infinity at line centre is equal to 1 (or $`T2338\mathrm{K}`$ for the few lines with $`\tau <1`$).

The last lines of the table give the total line luminosities, the continuum luminosity and the sum of the two. Continuum luminosity is not calculated (being actually very small) when the central temperature is $`1000\mathrm{K}`$; continuum specific luminosity is only a rough estimate, since we took $`\mathrm{\Delta }\nu =\nu _{\mathrm{peak}}`$, where $`\nu _{\mathrm{peak}}1.4\times 10^{14}\mathrm{Hz}`$ is the frequency of the peak of blackbody emission at $`T=2338\mathrm{K}`$.


format: text

text:
# Electronic conduction in multi-walled carbon nanotubes: Role of intershell coupling and incommensurability

## Abstract

Geometry incommensurability between weakly coupled shells in multi-walled carbon nanotubes is shown to be the origin of unconventional electronic conduction mechanism, power-law scaling of the conductance, and remarkable magnetotransport and low temperature dependent conductivity when the dephasing mechanism is dominated by weak electron-electron coupling.

Single-walled carbon nanotubes (SWNTs) can be either metallic of semiconducting depending on their helicities, i.e. how the graphene sheet is rolled up. For weak uniform disorder, the application of the Fermi golden rule for metallic SWNTs have demonstrated $`\mu m`$ long mean free paths  in agreement with experiments, clearly pointing towards ballistic transport.

Multi-walled carbon nanotubes (MWNTs) typically consist of several tenths of coaxial shells with intershell coupling between neighboring layers. The understanding of intrinsic properties of MWNTs remain elusive, since different transport scenario have been reported from ballistic, diffusive and low resistive, to insulating behaviors. An underestimated feature of MWNT’s geometry is that depending on their helical vector $`(n,m)`$ (expressed in the basis of the honeycomb lattice), two consecutive shells are commensurate (resp. incommensurate), if the ratio between their respective unit cell lengths $`T_{(n,m)}`$ along the tube axis is a rational (resp. irrational) number. An illustrative case is given on Fig.1 where the two double-walled nanotubes $`(9,0)\mathrm{@}(18,0)`$ and $`(9,0)\mathrm{@}(10,10)`$ are shown. Apparently, these defect-free structures should conduct electrons in a similar fashion since they are constituted of metallic shells with nearly identical diameters. Notwithstanding, by introducing coupling between shells, the $`(9,0)\mathrm{@}(18,0)`$ remains a periodic system, with a unit cell of length $`T_{(9,0)}=T_{(18,0)}2.15\AA `$, whereas the $`(9,0)\mathrm{@}(10,10)`$ can not be described by a single cell, since $`T_{(10,10)}/T_{(9,0)}=\sqrt{3}`$, so that coupling strengths between nearest-neighbors carbon atoms located at different shells, are not translationally invariant along the tube axis. The system is therefore not-periodic and electronic wave-coherence is expected to be qualitatively different. Incommensurate structures have an increasingly higher probability to be synthetized than fully commensurate structures, when increasing the number of inner shells.

Incommensurability was recently shown to yield anomalous corrugation properties . Non-ballistic conduction was reported and implications on magnetotransport discussed. Hereafter after defining the hamiltonian together with the computational method used, several consequences of the non-ballistic conduction are analyzed for transport coefficients in MWNTs.

Hamiltonian parameters for the MWNTs.-Our tight-binding model features one $`p_{}`$-orbital per carbon atom, zero onsite energies, and constant nearest-neighbor hopping on each layer n (n.n.), and hopping between neighboring layers (n.l.) are defined by :

$$=\gamma _0\left[\underset{i,j\text{ }n.n.}{}|p_{}^jp_{}^i|\right]\beta \left[\underset{i,jn.l.}{}\mathrm{cos}(\theta _{ij})e^{\frac{d_{ij}a}{\delta }}|p_{}^jp_{}^i|\right]$$

where $`\theta _{ij}`$ is the angle between the $`p_{}^i`$ and $`p_{}^j`$ orbitals, and $`d_{ij}`$ denotes their relative distance. The parameters used here are : $`\gamma _0=2.9eV`$, $`a=3.34\AA `$, $`\delta =0.45\AA `$. Estimate based on ab-initio calculations for $`\beta `$ gives $`\beta \gamma _0/8`$, but in order to get insight in the effect of $`\beta `$ on transport properties (since small change in the interlayer distance yield large variation of the coupling strength), several values of $`\beta `$ ($`0\beta \gamma _0`$) have been considered. Synthetized MWNTs contain typically a few tenth of inner layers, but here the study is restricted to 2 and 3-walled nanotubes, taking the intershell distance of $`3.4\AA `$ as in graphite.

Quantum diffusion in MWNTs.-We investigate quantum dynamics through the propagation of wavepackets in long MWNTs with about $`1000`$ unit cells ($``$ half million of carbon atoms). The time-dependent Schrödinger equation is solved numerically by using a polynomial expansion of the evolution operator. It gives us access to the diffusion coefficient of a wavepacket $`|\psi `$. If $`|\psi `$ is initially localized at the center (x=0) of the nanotube, then its diffusivity is simply $`D_\psi (t)=L_\psi (t)^2/t`$ with $`L_\psi (t)=\sqrt{\psi |(\widehat{X}(t)\widehat{X}(0))^2|\psi }`$ the spatial spreading of the wavepacket ($`\widehat{X}`$ is the position operator along the tube axis), and the system size is taken sufficiently large to avoid boundary effects. For random-phase or energy-filtered initial states extented to the whole system, periodic boundary conditions are used for the tube of length $``$, and the diffusivity is approximated as $`D_\psi (t)4\pi ^2/^2I_\psi (t)/t`$, where $`I_\psi (t)=\psi |\widehat{A}^+(t)\widehat{A}(t)|\psi `$ and $`\widehat{A}(t)=\mathrm{exp}(2i\pi \widehat{X}(t)/)\mathrm{exp}(2i\pi \widehat{X}(0)/)11`$. This estimation of $`D_\psi (t)`$ is accurate only for diffusion lengthes smaller than $``$, as checked by finite-size scaling analysis. The average spreading $`L(t)`$ and the average diffusion coefficient $`D(t)`$ are defined by $`L(t)=\sqrt{L_\psi ^2(t)}=\sqrt{tD(t)}`$, where $``$ denotes an average over many initial wavepackets. $`D(\tau _\phi )`$ is the average diffusivity along the nanotube axis, if at $`\tau _\phi `$ the electronic wavefunction looses its phase memory due to some inelastic scattering. The diffusion coefficient at $`\tau _\phi `$ are thus connected to the Kubo conductivity $`\sigma =e^2\rho D(\tau _\phi )`$, where $`\rho `$ is the density of states. This approach provides a good qualitative picture of wavepacket propagation in MWNTs constituted of conducting shells, with energies away from the charge neutrality point. Experimental studies demonstrate that chemical potential away from the charge neutrality point can be obtained upon doping. Krüger et al. obtained variations of the Fermi energy of the order of $`0.3eV0.5eV`$, corresponding to $`1015`$ conducting channels instead of $`2`$ per metallic layer.

Anomalous conduction in incommensurate MWNTs.-For disorder-free commensurate systems, a wavepacket initially localized in the outermost shell of the MWNTs, quickly transfer its weight to inner shells, due to intershell coupling. By following its spreading properties, we found that the wavepacket propagation is ballistic along the tube-axis over large distances ($`\mu `$m), i.e. $`L(t)vt`$ similarly to the case of perfect metallic single-walled nanotube.

However, for the incommensurate case, the transfer of the electronic wavepacket is followed by a non ballistic propagation given by $`L(t)𝒜t^\eta `$. The coefficient $`\eta `$ is found to decrease from $`1`$ to $`1/2`$ by increasing the number of coupled incommensurate shells. On Fig.2, ballistic versus non ballistic conduction mechanisms are illustrated in the defect-free incommensurate double-walled nanotubes $`(9,0)\mathrm{@}(18,0)`$ and $`(9,0)\mathrm{@}(10,10)`$. The case $`(6,4)\mathrm{@}(10,10)\mathrm{@}(17,13)`$ is also shown to demonstrate the effect of increasing the ”strength” of incommensurability on the wavepacket propagation. The interlayer coupling parameter taken here is given by $`\beta =\gamma _0/3`$ which is slightly larger than the expected ab-initio value. By increasing $`\beta `$, a reduction of power-law exponent $`\eta `$ is obtained, similarly to what is found by increasing the number of coupled incommensurate shells, for a fixed value of $`\beta `$.

Magnetotransport properties.-The effect of magnetic field and disorder when charges propagate predominantly on the outermost shell has been analyzed in . Here, we focus on the magnetotransport in incommensurate MWNTs if the magnetic field is applied parallel to the tube axis. Magnetic field dependent diffusion coefficients of electronic wavepackets are given on Fig.3., for incommensurate 2-wall $`(9,0)\mathrm{@}(10,10)`$ (main frame) in the power-law diffusion regime (with two values of $`\beta =\gamma _0/3,\gamma _0`$). $`\mathrm{\Phi }_0`$-periodic oscillations and positive magnetoresistance at low field are observed, similar to what is predicted for the ballistic regime in the SWNT case. Note that the Aharonov-Bohm oscillations for the commensurate $`(9,0)\mathrm{@}(18,0)`$ are still $`\mathrm{\Phi }_0`$-periodic, but with much smaller amplitude (not shown here). By increasing $`\beta `$, the coupling strength between walls, the magnetoconductance is reduced without modification of the oscillatory behavior.

Differently, for the incommensurate nanotube $`(6,4)\mathrm{@}(10,10)\mathrm{@}(17,13)`$ (inset of Fig.3), in which diffusive-like propagation takes place, there is evidence for negative magnetoresistance at low field and $`\mathrm{\Phi }_0/2`$-periodicity of $`D(\tau _\varphi ,\mathrm{\Phi })`$. From the saturation of diffusivity, we deduce some effective elastic mean free path $`\stackrel{~}{l}_e`$, that turns out to be slightly smaller that the outer shell circumference (for $`\beta =\gamma _0`$). By considering a reduced value of the coupling parameter ($`\beta =\gamma _0/3`$-not shown here), despite the fact that the effective mean free path becomes larger than the outershell circumference, still negative magnetoresistance at low field and $`\mathrm{\Phi }_0/2`$-periodic oscillation persist. These results confirm that magnetotransport properties of MWNTs are very sensitive to the geometry, the number of shells carrying current, and the hamiltonian parameters. In disordered systems, the basic scheme is that of ballistic electrons scattering on random impurities. The $`\mathrm{\Phi }_0/2`$ periodicity showed here, also results from quantum interferences of the electronic pathways around the cylinder wrist. However the usual scheme of interfering backscattered electronic trajectories is not strictly applicable for incommensurate disorder-free MWNTs, as confirmed by a study of the probability of return to the origin of wavepackets. Fujiwara et al. observed a surprising $`\mathrm{\Phi }_0/3`$-oscillation of magnetoconductance and gave an interpretation in terms of superposition of phase shift due to spectral effects, assuming three internal tubes with different helicity (same situation as in $`(6,4)\mathrm{@}(10,10)\mathrm{@}(17,13)`$). The effect of magnetic field on density of states (DoS) has been discussed theoretically and should also be considered to fully understand the magnetotransport in SWNTs and moderately large diameter MWNTs.

Electron-electron interactions and temperature dependence of inelastic scattering times and conductivity-.We focus now on the combined effect of anomalous propagation of wavepackets (ASWP) together with a dephasing mechanism dominated by electron-electron interactions. Recently, the influence of ASWP on frequency and temperature dependences of inelastic scattering times and conductivities, in systems close to a metal-insulator transition, in aperiodic incommensurate and quasiperiodic structures has been demonstrated.

The effect of electron-electron interaction in transport properties of MWNTs can be addressed within two limits. The strong coupling case is related to the existence of a Luttinger liquid phase, which has been discussed theoretically recently for metallic commensurate MWNTs . The effect of incommensurability on Luttinger liquid was shown in other context to develop new phase diagram, result that could be valuable in the context of incommensurate carbon nanotubes. In the weak coupling limit, the electron-electron interaction acts as the main dephasing mechanism and monitor the temperature dependences of inelastic scattering times and conductivity.

The interference pattern between two waves $`\mathrm{\Psi }_1(t)`$ and $`\mathrm{\Psi }_2(t)`$ which propagate coherently within the system is obtained from the evaluation of the relative phase difference accumulated in time. In the limit of low temperature, the inelastic coupling with other electrons, which produces a loss of phase coherence, can be expressed by an external time-dependent potential $`V(L_\mathrm{\Psi }(t),t)`$, whose effect is to superimpose new phase factors

$$e^{i\varphi }=e^{i({\scriptscriptstyle V(L_\mathrm{\Psi }(t),t)𝑑t/\mathrm{}})}=𝒫(\varphi )e^{i\varphi }𝑑\varphi $$

with $`\varphi `$ defined by a probability distribution $`𝒫(\varphi )`$. As the potential is time-dependent, the initial phase interference pattern is reduced in amplitude since the phase factors are modified by uncorrelated phase shifts. The loss of phase coherence has been shown to be driven by the increase of phase uncertainty between the superimposed phases accumulated by the two non-interacting waves $`\mathrm{\Psi }_1(t)`$ and $`\mathrm{\Psi }_2(t)`$ and related (under certain approximation, see ) to $`\delta \varphi _{12}^2=\varphi _{12}^2\varphi _{12}^2`$ ($`\varphi _{12}=𝑑tV(L_{\mathrm{\Psi }_1}(t),t)V(L_{\mathrm{\Psi }_2}(t),t)`$) with $`V(L_\mathrm{\Psi }(t),t)=\frac{e}{c}\frac{dL_\mathrm{\Psi }(t)}{dt}𝒜(L_\mathrm{\Psi }(t),t)`$, where $`𝒜(L_\mathrm{\Psi }(t),t)`$ is the potential vector. The inelastic phase coherence time $`\tau _\phi `$ defines the limit of coherent phase interference and thus comes out from

$$\delta \varphi _{12}^2e^2k_BT_0^{\tau _\phi }𝑑tL_{\mathrm{\Psi }_1}(t)L_{\mathrm{\Psi }_2}(t)^{2d}1$$

In a diffusive medium, $`L_{\mathrm{\Psi }_1}(t)L_{\mathrm{\Psi }_2}(t)t^{1/2}`$, which triggers $`\tau _\phi ^1T^{2/3}`$, behavior clearly identified experimentally at low temperature ($`T<4K`$) in metallic wires.

In defect-free incommensurate MWNTs, ASWP implies a power-law dependent diffusivity $`D(\tau _\phi )\tau _\phi ^{2\eta 1}`$, which thus results in a $`\tau _\phi ^1T^{1/(1+\eta )}`$ power-law for the inelastic scattering times. Subsequently, an anomalous temperature dependence of the conductivity follows

$$\sigma (T)T^{{\displaystyle \frac{(12\eta )}{(1+\eta )}}}$$

In the main frame of Fig.4 such behavior is exemplified for different exponents $`\eta `$. One notes that a temperature independent conductivity is found for the limit $`\eta 1/2`$ (diffusive-like motion), and the decrease of conductivity with temperature may be completely opposite for subdiffusive regime ($`\eta <1/2`$). Although in our calculations $`\eta `$ remains above $`1/2`$, such possibility, often found in quasiperiodic systems, can not be fully discarded for more complex situations. Recently Natelson et al.  have reported a continuous change of the power law exponent of the inelastic scattering times $`\tau _\phi `$ by varying the geometry of a quantum wire, i.e. by tuning quantum interferences. Here, we have ascribed similar behavior to the anomalous propagation steming from incommensurability in MWNTs. Along the same lines, one also notes that a study of the frequency dependent conductivity could be an interesting way to unveil anomalous quantum diffusion. Frequency dependent study of electronic impedance has been recently demonstrated to be relevant for carbon nanotubes.

Anomalous scaling of the conductance-.The relation between the electronic conductivity and the conductance allows to anticipate an anomalous length dependence of the Landauer conductance. Indeed, in a non-ballistic conduction, electrons propagate according $`L(t)𝒜t^\eta `$, with $`𝒜`$ and $`\eta `$ depending on the MWNT geometry but also the number of shells participating in transport ($`N_c`$). The conductance is related to the conductivity through $`G(L)\sigma (L)/L`$. By defining $`D_{nb}(L)`$ (resp. $`D_b(L)`$) the value of the non-ballistic diffusion coefficient (resp. ballistic diffusion coefficient for a fully commensurate metallic MWNT with $`N_c`$ shells), then the conductance at scale $`L`$ writes

$`{\displaystyle \frac{G_{nb}(L)}{2N_ce^2/h}}`$ $``$ $`{\displaystyle \frac{D_{nb}(L)}{D_b(L)}}`$ (1)
$`G_{nb}(L)`$ $``$ $`{\displaystyle \frac{2N_ce^2}{h}}\times \alpha \times L^{(\eta 1)/\eta }`$ (2)
with a factor $`\alpha (𝒜/v_F)^{1/\eta }`$ that will depend on the particular system and chemical potential considered. Indeed, if the chemical potential is close to the charge neutrality point, fewer conducting channels will be available because of the presence of semiconducting shells weakly conducting at such energy. The ballistic case $`\eta =1`$ corresponds to $`G\frac{N_ce^2}{h}`$, that is a length-independent conductance, whereas in the diffusive limit, the conductance scales as $`1/L`$ (physically related with the distance between electrodes). To test such physical effect, experiments should be carried out on highly purified MWNTs with typical diameter of $`1020nm`$, since it will allow to tune the participation in transport of more than 10 shells, under reasonably weak chemical doping ($`1\%`$), or directly by varying the applied voltage $`V_{bias}`$ (the number of shells $`N_C`$, participating in conduction being a function of $`V_{bias}`$). Furthermore, the study of the conductance with the distance between leads probing the sample requires perfect ohmic contacts between electrodes and the nanotube.

Conclusions-.Incommensurability in MWNTs is a geometrical degree of freedom at the origin of non ballistic conduction. The redistribution of the wavepacket in the whole object is followed by a specific multiple wave scattering that jeopardize ballistic conduction. All transport coefficients should be sensitive to such effects, especially the temperature dependence of inelastic scattering times and conductivity driven by electron-electron interactions at low temperature, as well as the length-scaling of the Landauer conductance. Experiments should be performed in conjunction with theoretical calculations of the relevant microscopic parameters ($`𝒜,\eta `$) that depend on the geometry of the MWNT as well as the number of shells participating in conduction.

Acknowledgments: Financial support from NAMITECH \[ERBFMRX-CT96-0067(DG12-MITH)\], DGES (PB98-0345), COMELCAN( HPRN-CT-2000-00128), JCyL (VA28/99) and C4 are acknowledged. S.R. is indebted to Riichiro Saito and Hélène Bouchiat for stimulating discussions.

Figure captions:

Figure 1: Geometrical representation of $`(9,0)\mathrm{@}(18,0)`$(right) and $`(9,0)\mathrm{@}(10,10)`$(left) double-walled carbon nanotubes ASK me by e-mail ! at sroche@cea.fr.

Figure 2: Main Frame: Averaged diffusion coefficient (arb. unit) for $`(9,0)\mathrm{@}(18,0)`$ and $`(9,0)\mathrm{@}(10,10)`$ with $`\beta =\gamma _0/3`$. Inset: Avergared diffusion coefficient for the incommensurate MWNT $`(6,4)\mathrm{@}(10,10)\mathrm{@}(17,13)`$. One notes that $`T_{(6,4)}=3\sqrt{19}a_{cc}`$, $`T_{(10,10)}=\sqrt{3}a_{cc}`$, and $`T_{(17,13)}=3\sqrt{679}a_{cc}`$ where $`a_{cc}=1.42\AA `$ is the interatomic distance between carbon atoms.

Figure 3: Main frame : Magnetic field dependence of the diffusion coefficient for the incommensurate double-walled $`(9,0)\mathrm{@}(10,10)`$ for two values of $`\beta =\gamma _0/3`$ (solid line, right y-axis) and $`\beta =\gamma _0`$ (dashed line, left y-axis) evaluated at the same time $`\tau _\varphi `$. Inset : Case for $`(6,4)\mathrm{@}(10,10)\mathrm{@}(17,13)`$ at two different times $`\tau _\varphi `$($`\mathrm{}/\gamma _0`$ unit) for $`\beta =\gamma _0`$.

Figure 4:Main Frame: Scaling of the temperature dependence of the conductivity with decreasing exponent $`\eta `$. Inset: Length dependence of the conductance for the same parameters. Bold solid lines are for $`\eta 1`$, and bold dashed lines for $`\eta 1/2`$, whereas other solid lines interpolate between 1 and 1/2.