text:
# Testing the Collective Properties of Small-World Networks through Roughness Scaling

## Abstract

Motivated by a fundamental synchronization problem in scalable parallel computing and by a recent criterion for “mean-field” synchronizability in interacting systems, we study the Edwards-Wilkinson model on two variations of a small-world network. In the first version each site has exactly one random link of strength $`p`$, while in the second one each site on average has $`p`$ links of unit strength. We construct a perturbative description for the width of the stationary-state surface (a measure of synchronization), in the weak- and sparse-coupling limits, respectively, and verify the results by performing exact numerical diagonalization. The width remains finite in both cases, but exhibits anomalous scaling with $`p`$ in the latter for $`d2`$.

Crossovers from low-dimensional to mean-field-like behavior have been studied and found in various interacting systems BARRAT ; GITTERMAN ; xy\_sw ; ising\_sw ; phase\_sw  when the (original) regular and short-range interaction topology is modified to a small-world (SW) network WATTS98 ; SW\_book ; NEWMAN . Mean-field behavior is commonly observed, even when the random links are added to a one-dimensional original “substrate” xy\_sw ; ising\_sw ; phase\_sw . Recently, a general criterion for the crossover to mean-field behavior was given by Hastings HASTINGS03 .

The equilibrium Edwards-Wilkinson (EW) model EW  (considering it as the $`\varphi ^2`$ Gaussian model at the continuous phase transition) is a particularly interesting case in that it does not satisfy the mean-field criterion HASTINGS03  when each pair of sites is connected with probability $`p/N`$ with a link of unit strength on top of a one-dimensional regular substrate. Thus, the EW model on this SW network is expected to display drastically different scaling properties (as a function of $`p`$) than the mean-field version of the model in one dimension. In this paper we address this question.

The scaling properties of the EW model on regular or random networks have direct relevance to the scalability and synchronizability of parallel discrete-event simulations (PDES)KORNISS00 ; KORNISS03a . In PDES schemes the individual processing elements (PE) generate their own time streams for update attempts (local simulated times) FUJI ; LUBA . The synchronization/communication between neighboring PEs (following the interaction topology of an underlying short-range interacting system with asynchronous dynamics) leads to Kardar-Parisi-Zhang (KPZ)-like kinetic roughening KPZ ; BARABASI  of the simulated time horizon KORNISS00 ; KORNISS\_MRS . For a one-dimensional chain or ring of PEs, the steady-state performance of the PDES scheme is governed by the EW Hamiltonian KORNISS00 . In particular, for $`N`$ PEs, the width of the simulated time horizon diverges as $`N^{1/2}`$, seriously hindering efficient data collection and state saving KORNISS03a ; KORNISS\_ACM . Since scalable data management crucially depends on the finiteness of the width of the time horizon (the spread of the progress of the individual PEs), one must suppress the diverging fluctuations of the simulated time horizon. As an alternative to costly and frequent global synchronizations among the PEs, an autonomous small-world synchronization scheme was demonstrated to work KORNISS03a . This finding provided another concrete example for synchronizability in generalized multi-agent systems facilitated by a small-world network WATTS98 ; Strogatz\_review .

Motivated by both the generic scaling properties of SW networks, and also the specific applications to scalable PDES synchronization schemes, we will consider two different variations of a the EW model on a SW network. In the “soft” and frequently studied version of the SW network, random links of unit strength are added to the one-dimensional substrate with probability $`p/N`$ to each pair of sites NEWMAN ; NEWMAN\_WATTS ; MONA . In the “hard” version, each site has exactly one random link (in addition to the nearest-neighbors) and the strength of the interaction through the random links is $`p`$. That is, pairs of sites are selected at random, and once they are linked, they cannot be selected again. The latter construction originates from scalable PDES schemes, where fluctuations in the individual connectivity of the PEs are to be avoided KORNISS03a . We introduced the above terminology for the two versions of the network in regard to the eigenvalue spectrum of the respective coupling matrix MONA , discussed at the end.

We consider the equation (for a single realization of the small world)

$$_th_i=(2h_ih_{i+1}h_{i1})\underset{j=1}{\overset{N}{}}J_{ij}(h_ih_j)+\eta _i(t),$$
(1)
where $`h_i`$ is the surface height, $`\eta _i(t)`$ is a delta-correlated Gaussian noise with variance $`2`$ (without loss of generality), and we have dropped the $`t`$-dependence from the argument of $`h_i`$ for brevity. The symmetric matrix $`J_{ij}`$ represents the quenched random links on top of a one-dimensional lattice of length $`N`$ with periodic boundary conditions. For the hard version of the small world $`J_{ij}`$ has exactly one non-zero element (being equal to $`p`$) in each row and column. This construction results in $`_lJ_{il}=p`$ for all $`i`$. For the soft version of the small world, each element of $`J_{ij}`$ (e.g., above the diagonal) is 1 with probability $`p/N`$ and zero otherwise. In this case $`[_lJ_{il}]=p`$, where $`[\mathrm{}]`$ denotes the average over the network disorder, i.e., the average coordination number is $`p`$. We write Eq. (1) as $`_th_i=_j\mathrm{\Gamma }_{ij}h_j+\eta _i`$, where $`\mathrm{\Gamma }=\mathrm{\Gamma }^o+V`$. Here, $`\mathrm{\Gamma }^o`$ is the Laplacian of the original one-dimensional ring, while $`V_{ij}=J_{ij}+\delta _{ij}_lJ_{il}`$ is the Laplacian on the random part of the network.

In this Letter we focus on the behavior of the width, which probes the generic collective properties of the underlying networks by providing a sensitive measure of synchronization KORNISS03a . For a given realization of the small-world network the average surface width characterizing the roughness is equal to

$$w^2_N\frac{1}{N}\underset{i=1}{\overset{N}{}}(h_i\overline{h})^2=\frac{1}{N}\underset{k=1}{\overset{N1}{}}\frac{1}{\lambda _k},$$
(2)
where $`\overline{h}`$$`=`$$`(1/N)_{i=1}^Nh_i`$ is the mean height, $`\mathrm{}`$ denotes an ensemble average over the noise in Eq. (1), and $`\lambda _k`$ are the non-zero eigenvalues of the real symmetric coupling matrix $`\mathrm{\Gamma }`$. Our numerical scheme relied on the exact diagonalization NUM\_REC  of the coupling matrix $`\mathrm{\Gamma }`$ for a given realization of the network, and then we exploited the right side of Eq. (2) to obtain the width. Finally, we obtained the disorder-averaged (denoted by $`[\mathrm{}]`$) width $`[w^2_N]`$ by averaging over a large number of (ranging from $`100`$ to $`1000`$) realizations.

To motivate the results that follow, we consider the role of dimension $`d`$ for the EW model on the soft version of the network, along the lines of Ref. HASTINGS03 . As an approximation to the soft network, consider instead an alternate model, where in addition to the nearest-neighbor connections, each site is coupled to all others with strength $`p/N`$. This yields the trivial (inhomogeneous Landau) mean-field behavior with the mass $`p`$ meaf  and the corresponding correlation length $`\xi p^{1/2}`$, and a correlation volume $`\xi ^dp^{d/2}`$. Now, let us check if this model is a valid description of the soft network. For the soft network, within this volume one would have on average $`p\xi ^dp^{1d/2}`$ links. For $`d`$$`>`$$`2`$ this number diverges as $`p`$$``$$`0`$; then for small $`p`$ there are a large number of links leaving the volume and the sample-to-sample fluctuation in the mass in a correlation volume is negligible compared to the mass itself so that the trivial mean-field behavior is expected to be valid. For $`d`$$``$$`2`$, the trivial behavior breaks down as we will see in a perturbative calculation. In contrast, for the hard version of the network the density of random links is unity, and the trivial mean-field scaling is expected to hold for all $`d`$ in the $`p`$$``$$`0`$ limit.

The same result can be obtained by using the general criterion HASTINGS03  in terms of the critical exponents for the correlation length and for the susceptibility, $`\nu `$ and $`\gamma `$, respectively. When viewed as a model near a continuous phase transition, the EW model has $`\nu `$$`=`$$`1/2`$, $`\gamma `$$`=`$$`1`$. The criterion $`\nu d/2+\gamma /2>1`$ for the validity of mean-field behavior is thus violated for $`d`$$`=`$$`1`$, while $`d`$$`=`$$`2`$ is marginal.

We define the propagator $`G`$ to be equal to $`\mathrm{\Gamma }^1`$ in the space of non-zero eigenvalues of $`\mathrm{\Gamma }`$, while $`G`$ vanishes when acting on the zero mode of $`\mathrm{\Gamma }`$. Thus, $`G=P(\mathrm{\Gamma }+iϵ)^1`$, where $`P_{ij}=\delta _{ij}1/N`$ is the projector onto the vector space orthogonal to this zero mode. We now present the results of a perturbative approach for the disorder-averaged propagator $`[G]`$, using techniques of impurity averaged perturbation theory iapt , leaving details to be published elsewhere future . The perturbative expansion of $`[G]`$ can be obtained by $`[G]=G^o[G^oVG]=G^o[G^oVG^o]+[G^oVG^oVG^o]\mathrm{}`$, where $`G^o`$$`=`$$`P(\mathrm{\Gamma }^o+iϵ)^1`$ is the propagator of the Laplacian on the original one-dimensional lattice. To obtain $`[G]`$ it is necessary to average this expansion over the network disorder in $`V`$.

To deal with only one-particle irreducible disorder-averaged diagrams, we calculate the self-energy $`\mathrm{\Sigma }=([G])^1(G^o)^1`$ perturbatively for both the soft and hard versions of the small-world. In Fig. 1, a single line denotes the propagator $`G^o`$ while a double line denotes $`[G]`$. The relation between $`[G],G^o`$, and $`\mathrm{\Sigma }`$ is shown at the top of the figure. In these calculations, a cross with no dashed lines attached is used to denote the average over different realizations of the network of a diagonal term $`V_{ii}`$ in $`V`$. A pair of crosses connected by a dashed line is used to denote an average $`[V_{ii}V_{jj}][V_{ii}][V_{jj}]`$, while three or more crosses connected by dashed lines are used to denote higher cumulants. Similarly, circles connected by dashed lines are used to denote averages of off-diagonal terms $`J_{ij}`$ in $`V`$. Dashed lines can connect both circles and crosses. Terms in the expansion with a circle, not connected by dashed lines to other circles, vanish as $`1/N`$ for large $`N`$ and so may be neglected. For the hard network, terms in the expansion with two or more crosses connected by dashed lines vanish since there are no random fluctuations in $`\delta _{ij}_lJ_{il}`$ in this case, while for the soft network these terms do appear. In the hard network, each power of $`V`$ comes with a power of $`p`$, while in soft network, each set of circles or crosses connected by dashed lines comes with a power of $`p`$.

In the hard version of the network, we can proceed by expanding $`[G]`$ in powers of $`V`$. The first diagram in Fig. 1(a) yields the lowest order result $`\mathrm{\Sigma }=p`$. Adding the second diagram yields $`\mathrm{\Sigma }=pp^2[G]_{ii}=pp^2/(2\sqrt{\mathrm{\Sigma }})`$, where we use the fact that, for large $`N`$, $`[G]_{ii}=1/(2\sqrt{\mathrm{\Sigma }})`$ plus terms of order $`\sqrt{\mathrm{\Sigma }}`$ (see, e.g., TORO ). Inserting the leading order result for $`\mathrm{\Sigma }`$ into the second diagram we find

$$\mathrm{\Sigma }_{\mathrm{hard}}=p(1/2)p^{3/2}+\mathrm{},$$
(3)
so that the higher powers of $`V`$ lead to corrections to $`\mathrm{\Sigma }`$ which are higher order in $`p`$ as required.

In the soft version of the network, the above procedure does not work. The first diagram in Fig. 1(b) yields $`\mathrm{\Sigma }=p`$. Inserting this result for $`\mathrm{\Sigma }`$ into the next two diagrams leads to $`\mathrm{\Sigma }=pp^{1/2}`$, so that the expansion in $`V`$ does not lead to an expansion in $`p`$. To correct this, at “leading order” we instead sum up all diagrams involving a single link. The first terms in this sum are shown in the first pair of parenthesis in Fig. 1(b). The infinite sum yields $`\mathrm{\Sigma }=p(12[G]_{ii}+4[G]_{ii}^2\mathrm{})=p/(1+2[G]_{ii})`$. Solving this equation self-consistently to lowest order in $`p`$ yields $`\mathrm{\Sigma }=p^2`$. Physically, since the density of links is small, of order $`p`$, this result consists of exactly solving the interaction with a single link.

Having done this infinite summation, we can now consider another infinite series of diagrams (containing two links), starting with the last one shown in Fig. 1(b), and adding additional diagrams where each single interaction with the given link is replaced by two, three, or more interactions with the link, as well as diagrams in which crosses are replaced by circles. Let us determine the order of this summation in $`p`$. The infinite summation replaces the interaction with a given link with the sum: $`12[G]_{ii}+4[G]_{ii}^2\mathrm{}=1/(1+2[G]_{ii})p`$. Thus, the diagram has four such “scatterings”, leading to a result of order $`p^4`$. There is an additional factor of $`p^2`$, due to the appearance of two sets of circles/crosses connected by dashed lines. There is a further factor of $`[G]_{ii}^3p^3`$ from the three Green’s functions. Finally, there is a summation over the spacing between the two impurities; this spacing is of order $`p^1`$. As a result, this sum yields a result which is again of order $`p^2`$. Note that these diagrams lead to a self-energy $`\mathrm{\Sigma }_{ij}`$ which is off-diagonal in real-space, and hence is momentum dependent; we consider here the zero momentum part of $`\mathrm{\Sigma }`$. More complicated diagrams continue to yield results of order $`p^2`$. As a result, after this resummation, we are able to determine only that $`\mathrm{\Sigma }_{\mathrm{soft}}`$ scales as $`p^2`$, but not the exact coefficient, and not the higher order corrections:

$$\mathrm{\Sigma }_{\mathrm{soft}}p^2$$
(4)
Then, the disorder-averaged width can be found from $`\mathrm{\Sigma }`$ by $`\left[w^2\right]=\left[G\right]_{ii}\frac{1}{2\sqrt{\mathrm{\Sigma }}}`$. Thus, the asymptotic small-$`p`$ behavior of the width in the thermodynamic limit is $`[w^2]_{\mathrm{hard}}1/(2\sqrt{p})`$ and $`[w^2]_{\mathrm{soft}}1/p`$, for the respective versions of the small-world networks. These asymptotic small-$`p`$, infinite system-size behaviors are indicated with the two slopes next to the numerical data in Fig. 2. In an attempt to match the behavior of the width for finite systems in the hard network, we used the finite-system version of the propagator (see, e.g., TORO ) with the effective mass from Eq. (3). Note that we ignored various finite-size corrections to $`\mathrm{\Sigma }`$ itself, as being too complicated to be worth calculating. The results lead to good agreement and are systematically improved by including the next-to-leading order corrections for the hard network (Fig. 3).

These results indicate that in both SW versions, the width approaches a finite value for any non-zero value of $`p`$ as $`N`$$``$$`\mathrm{}`$. In the hard network, the scaling of the width asymptotically approaches that of the mean-field version of Eq. (1) meaf , while the soft network does not exhibit this scaling, as expected from the general criterion HASTINGS03 . Similarly, for $`d`$$`=`$$`2`$, one finds $`\mathrm{\Sigma }_{\mathrm{hard}}=p\mathrm{}`$, while for the soft case the resummation of the expansion enables one to find a logarithmic correction to the trivial mean-field behavior future :

$$\mathrm{\Sigma }_{\mathrm{soft}}=2\pi p/|\mathrm{log}(p)|\mathrm{}.$$
(5)
For $`d`$$`>`$$`2`$, $`\mathrm{\Sigma }`$ is asymptotically of order $`p`$ for both soft and hard cases, i.e, the systems are effectively mean-field.

We now interpret the effective mass obtained perturbatively, Eqs. (3) and (4), in light of earlier results. The density of states $`\rho (\lambda )`$ (the eigenvalue spectrum of the the coupling matrix $`\mathrm{\Gamma }`$) has been recently investigated in the context of diffusion on the “soft” version of the small-world network MONA . The soft construction of the network allows for the existence of arbitrarily long “pure” chain segments when the system size $`N`$ goes to infinity MONA . Although the probability of these quasilinear chain segments of length $`l`$ is exponentially small ($`e^{pl}`$), they can contribute to eigenvalues of order $`1/l^2`$ LIFS ; BRAY . Summing up over large $`l`$ values with the exponential weight above yields $`\rho (\lambda )(1/\sqrt{\lambda })e^{cp/\sqrt{\lambda }}`$ for small $`\lambda `$ MONA , where $`c`$ is a constant. Thus, there is no gap in the spectrum. The density of states, however, vanishes exponentially fast for small $`\lambda `$ as a result of the essential singularity in the exponent, and one can refer to $`p^2`$ as a pseudo-gap. From the right side of Eq. (2) it follows that in the $`N`$$``$$`\mathrm{}`$ limit, the disorder-averaged width can be expressed in terms of the density of states as $`[w^2]=(1/\lambda )\rho (\lambda )𝑑\lambda `$. The small-$`\lambda `$ behavior of $`\rho (\lambda )`$ determines whether the width remains finite or diverges in the thermodynamic limit. The exponentially small $`\rho (\lambda )`$ above more than compensates for the term $`\lambda `$ in the denominator yielding $`[w^2]_{\mathrm{soft}}1/p`$.

Our perturbative scheme can be extended to calculate $`\rho (\lambda )`$ by defining $`[G(\lambda )]=P[(\mathrm{\Gamma }\lambda +iϵ)^1]`$ so that $`\rho (\lambda )=\mathrm{Im}[G(\lambda )]/\pi `$. For the hard version of the network, we expect that there is indeed a true gap of order $`p`$. The self-energy is only weakly $`\lambda `$-dependent. Each site has one link and for any fixed segment of the chain, as $`N`$$``$$`\mathrm{}`$, there is a vanishing probability that any of these links connect two sites in the given segment, thus preventing the construction of quasilinear segments as above.

We thank Z. Toroczkai, H. Guclu, Z. Rácz, M.A. Novotny, P.A. Rikvold, and G. Györgyi for discussions. G.K. thanks CNLS, Los Alamos National Laboratory for their hospitality, where this manuscript was completed. We acknowledge the support of NSF through Grant No. DMR-0113049 and the support of the Research Corporation Grant No. RI0761. MBH was supported by US DOE W-7405-ENG-36.


format: text

text:
# Cooling of X-ray Emitting Gas by Heat Conduction in the Center of Cooling Flow Clusters

## 1 Introduction

Recent Chandra and XMM-Newton observations appear to be in contradiction with the simplest cluster cooling flow models. While hot, $`T18\times 10^7`$ K, X-ray emitting gas seems to cool at rates of $`10500M_{}`$ yr-1, there is little indication for gas cooling below a temperature of $`10^7`$ K (e.g., Molendi & Pizzolato 2001; Kaastra et al. 2001; Peterson et al. 2001, 2002, 2003; Fabian et al. 2002a). Namely, the limit on the cooling rate below $`10^7`$ K inferred from X-ray observations is $`<20\%`$ of the mass cooling rates cited in the past (Fabian 2002), but compatible with the moderate cooling flow model of Soker et al. (2001). On the other hand, X-ray observations reveal properties that seem to support the cooling of gas. First, there is no indication for violent ongoing heating (e.g., shocks or cluster mergers) to offset radiative cooling. Hot bubbles of radio plasma may supply heat, but they also compress the gas, which reduces its cooling time (Soker et al. 2001). High spatial resolution X-ray images from Chandra reveal complex X-ray structures in the inner $`30`$ kpc of cooling flow clusters, (e.g., McNamara 2002 and references therein), including X-ray deficient bubbles that rise buoyantly transporting energy from the AGN to the intracluster medium (ICM).

Morris & Fabian (2003) suggest that the weakness of the X-ray emission lines below $`10^7`$ K can be at least partially reconciled with larger cooling rates if most of the metals are concentrated in a small fraction ($`10`$%) of the cooling gas. Soker et al. (2001) suggest that the cooling flow model can be put on a more solid ground if some assumptions are relaxed, e.g., the cooling flow age is taken to be much below the age of the clusters (Allen et al. 2001; Bayer-Kim et al. 2002). These relaxed assumptions lead to much lower mass cooling rates, in what was called the “moderate cooling flow model”. In the present paper, we continue along these lines, and show that moderate cooling flow models can be compatible with the non-detection of line emission from cooler gas if some of the cooling hot X-ray gas in moderate cooling flow models is cooling from $`10^7`$ K to $`10^4`$ K by heat conduction. Heat conduction has similar effects to those of mixing, as proposed by, e.g., Fabian et al. (2001, 2002a) and Bayer-Kim et al. (2002). Observations show that magnetic fields become strong at the center of cooling flow clusters (Eilek & Owen 2002, and references therein). Magnetic fields inhibit heat conduction across field lines. We propose that reconnection of magnetic field lines occurs between hot X-ray emitting gas and cooler gas at $`10^4`$ K. We argue that the conduction fronts are narrow, making heat conduction very efficient.

We note here that for strong magnetic fields, the self-inhibiting of the heat-flux along magnetic field lines due to scattering of conducting electrons off plasma waves is not efficient (Pistinner & Eichler 1998). When the magnetic field is weak, heat conduction is suppressed along field lines (Levinson & Eichler 1992), while when the magnetic field is strong, the suppression is small (Pistinner & Eichler 1998). Pistinner, Levinson, & Eichler (1996) consider this effect in cooling flow clusters. They show that in the outer parts where the thermal pressure is much higher than the magnetic pressure, the heat conduction is suppressed by this mechanism by one to three orders of magnitude, while in the inner parts it is not suppressed much if the magnetic field is strong (see their Figs. 1-3; eq. 69 of Pistinner & Eichler 1998 shows that modest suppression may occur). In addition to that effect, immediately after reconnection, heat conduction may be quite efficient in the turbulent region (Lazarian & Cho 2003). We therefore assume here that heat conduction is efficient along magnetic field lines in regions where the magnetic field is close to equipartition.

Our magnetic reconnection model is different from papers dealing with large scale ($`100`$ kpc) magnetic fields and conduction between exterior very hot intracluster gas and hot gas within the cooling flow region (e.g., Bregman & David 1988; Pistinner & Shaviv 1996; Norman & Meiksin 1996; Narayan & Medvedev 2001; Loeb 2002; Voigt et al. 2002; Fabian, Voigt, & Morris 2002b; Ruszkowski & Begelman 2002; Markevitch, Vikhlinin, & Forman 2002). In many of these models, heat conduction heats the gas in the cooling flow region, and may reduce the cooling rate or even eliminate it completely. Other studies argued that heat conduction is suppressed by a large factor, and hence plays no role in the cooling flow process. We instead examine local heat conduction at lower temperatures, after the gas has lost most of its thermal energy. The gas continues to cool to lower temperatures, but does so conductively, rather than radiatively. Since it cools rapidly due to conduction, we do not see strong low ionization X-ray lines (e.g., Peterson et al. 2001) or evidence for lower temperature X-ray multiphase gas (e.g., Molendi & Pizzolato 2001) in the X-ray spectra of cooling flow regions in clusters.

## 2 The Proposed Scenario

The basic scenario is as follows. Gas cools radiatively from the ambient cluster temperature to $`10^7`$ K, and slowly streams inward, becoming gravitationally bound to the cD galaxy. As a result of the cD galaxy motion through the ICM and radio bubble activity, there is a constant mixing of gas in the inner regions ($`r1030`$ kpc) of cooling flow clusters. Different parcels of gas collide with one another. In particular, we will assume that cold ($`T10^4`$ K) clouds move rapidly through the hotter ICM. Some processes involving cold clouds in cooling flows, among them heat conduction, were studied before (Loewenstein & Fabian 1990; Ferland, Fabian, & Johnstone 2002). We note that colder material is known to exist in many cooling flow clusters, as inferred from optical filaments (see § 4 below), and from recent CO observations conducted by Edge (2001; see also Edge & Frayer 2003, and Salome & Combes 2003) who found the molecular mass to be in the range of $`10^910^{11.5}M_{}`$. These cold clouds may result in large part from stripping of the ISM from galaxies and subclusters. However, our scenario is in accord with the presence of a moderate rate of cooling of intracluster gas, and this may produce all or part of the cold gas. Evidence for cooling flows still exists in the Chandra and XMM/Newton observations, but with much lower cooling rates than previously suggested (e.g., Bayer-Kim et al. 2002; David et al. 2001; Blanton et al. 2003). The cooling gas may cool to a temperature of $`10^4`$ K and below. In particular, the heat conduction process itself cools the hot gas, forming more cold clouds, although it can also evaporate the cold clouds. The sign of the conductive mass flux (evaporation or cooling?) is discussed further in § 3.

The observed line widths of the optical emission line clouds in cooling flows are typically a few hundred km s-1 (e.g., Heckman et al. 1989). Similar line widths are seen for the molecular hydrogen clouds seen in C0 (e.g., Edge 2001). Since we identify the cold clouds in our model with these observed materials (see § 4 below), we will assume that the cold clouds move through the ambient hot gas at speeds of a few hundred km s-1. These motions are mildly subsonic in the ambient gas, and are comparable or slightly greater than the velocity dispersions of stars in the same central regions of the cD galaxies. Thus, it is likely that these motions are orbital or infall velocities, perhaps supplemented by other dynamical effects (e.g., the expansion of radio bubbles).

We note that in the moderate cooling flow model, the mass cooling rate evolves with time, being zero at the the beginning (Soker et al. 2001; their discussion following eq. \[4.17\]). Therefore, for young clusters the number of clouds from cooling gas is small (cold clouds originating from striped galaxies may still be present), but so is the mass cooling rate. Therefore, even in young clusters the X-ray emission from gas at $`T10^7`$ K should be weak on average.

As a result of cooling and inflow, the density of the ICM increases by a factor of $`10`$. An isotropically-tangled magnetic field in the gas will have its pressure increased by a factor of $`10^{4/3}`$. Radial inward motion (Soker & Sarazin 1990), stretching of ICM magnetic fields by expanding jets and bubbles (Soker 1997), as well as a possible cluster dynamo (Godon, Soker, & White 1998), may result in a much larger increase in the magnetic pressure. A dynamically important magnetic field in the inner regions of cluster cooling flows is inferred from observations (Eilek & Owen 2002; Taylor, Fabian, & Allen 2002). The enhanced magnetic fields and the gas motions in the inner cooling flow region will most likely result in reconnection of magnetic field lines. The energy released by reconnection can play an important role in the energetics of the optical filaments (Sabra, Shields, & Filippenko 2000), via heating produced directly by reconnection (Jafelice & Friaca 1996; Godon, Soker & White 1998) and/or by Alfvén waves (Friaca et al. 1997), and/or by facilitating heat conduction (Böhringer & Fabian 1989).

At present, we are mainly interested in reconnection between the magnetic fields in cold clouds ($`10^4`$ K) and in the hot ($`10^7`$ K) ICM. Such reconnection can enhance heat conduction from the hot gas into the cold clouds, and this may result in the hot gas cooling conductively rather than radiatively. That is, the thermal energy released in cooling the hot gas is delivered to the cold clouds, and is radiated as optical/UV line emission at $`10^4`$ K. Replacing radiative cooling of the hot gas with conductive cooling will reduce the emission at low but X-ray emitting temperatures, and could account for the weak soft X-ray luminosity from gas at intermediate temperatures $`10^6T10^7`$ K. In addition to allowing conductive heat flow, reconnection may also cause mixing of the gas via the fast flow emanating from the reconnection site (Lazarian & Cho 2003). Mixing was also suggested as a mechanism to suppress X-ray emission lines below $`1`$ keV (e.g., Oegerle et al. 2001; Fabian et al. 2001, 2002a; Johnstone et al. 2002). Fabian et al. (2002a) briefly mention that heat conduction to cold gas plays a similar role to that of mixing. Recently, Cho et al. (2003) have argued that magnetic fields are relatively ineffective at suppressing turbulent heat conductions. It is unclear how turbulent in the regions around the clouds will be; the existence of large-scale abundance gradients in cooling flows indicates that global turbulence is not very efficient (Böhringer et al. 2004).

In principle, reconnection might break the field lines connecting hot and cold gas, rather than reconnect them, and this would suppress thermal conduction. Here, we assume that the field lines in the cold clouds are initially disconnected from the field lines in the hot ICM, as much as the magnetic field lines of the Earth are disconnected from the magnetic field lines of the Solar wind. We assume that reconnection only occurs near the stagnation point at the front of the cloud where ram pressure is maximum, and in the tail behind the motion of the cloud. A schematic description of the evolution of a single ICM magnetic field line is given in Figure 1. This is similar to what is observed in the Earth-Solar wind interaction (e.g., Oieroset et al. 2001). For this configuration, heat conduction substantially increases after the reconnection event.

One worry with this scenario might be that reconnection near the stagnation point at the front of the cloud would be too slow to occur before the magnetic fields in the ambient hot gas were advected downstream (in the frame of the cloud). Exactly at the stagnation point, the relative velocity is zero, but the relative velocity returns to its full value as the ambient gas moves around the cloud. Thus, magnetic fields lines might slip around the cloud before they reconnect. In strongly magnetized plasmas with magnetic field reversals, reconnection can occur at a fraction of the Alfven speed (e.g., Petschek and Thorne 1967). As noted above, the speed of the cloud relative to the ambient gas is a fraction of the sound speed in the ambient gas. At the front of the cloud near the stagnation point, this motion will be somewhat slower. Thus, there is (just) enough time for magnetic fields to reconnect at the front of the cloud. It is also possible that viscosity or Kelvin-Helmholtz instabilities will help to keep magnetic fields in the ambient gas close to the clouds for a somewhat longer period.

As drawn, Figure 1 shows a smooth laminar reconnection of the magnetic field lines behind the cloud. If there is a turbulent wake, the geometry will be more complex, but the same basic process may occur. What is crucial is that the field lines from the clouds are eventually disconnected from the ICM. A turbulent wake may even enhance the rate of reconnection. In any case, the stretched field lines will overcome the turbulence pressure and reconnect.

Rayleigh-Taylor and Kelvin-Helmholtz instability modes are expected to occur in the front interface between the cloud and the ICM. It is likely that the Kelvin-Helmholtz instabilities will be partly suppressed by the strong magnetic fields (e.g., Vikhlinin, Markevitch, & Murray 2001). In any case, the breakup of clouds by instabilities may not be a problem for the present scenario, since we require small clouds. In the next section the size range which is optimal for the proposed scenario is derived.

One concern with this argument is raised by the Chandra observations of “cold fronts” in clusters (see the review by Markevitch et al. 2002). A cold front is a contact discontinuity observed in many clusters of galaxies between a large region of relatively cool gas ($`10^7`$ K) which is moving rapidly through the hot ICM. In some cases, the cold fronts are the cooling cores associated with a merging subcluster. The sharp temperature jump observed in cold fronts implies that the magnetic field there practically inhibits heat conduction. Since the geometry is essentially identical to the geometry of cold clouds moving through hot ICM, why doesn’t magnetic field line reconnection occur in cold fronts and lead to rapid thermal conduction? The answer may lie in the relative strength of the magnetic field at the stagnation point at the front of the cold front or cold cloud. In the proposed scenario, the magnetic field in the cold clouds is more or less in pressure equipartition with the thermal gas, and it is quite strong, reaching $`10\%`$ of equipartition, in the ICM. The magnetic field in either sides of cold fronts may be increased by stretching and reach a value of $`10\mu \mathrm{G}`$, as inferred from the requirement to suppress instabilities (Markevitch et al. 2002), but it is still below equipartition values. Note that we are examining reconnection in the ambient post-shock and compressed gas which has a slow speed relative to the cloud. In this case it is the magnetic pressure that drives the reconnection, and not a direct fast collision of magnetic field lines. In this case the condition for fast reconnection to occur is that the magnetic pressure be comparable or larger than the thermal pressure (Petschek & Thorne 1967; see discussion in section IIb of Soker & Sarazin 1990).

## 3 Magnetic Field Reconnection and Heat Conduction

In this section we quantify some of the arguments, and estimate the importance of some of the processes, discussed in the previous section. We start by considering the rate of heat conduction across the interfaces between hot ICM and cold cloud gas in the inner regions of cooling flows. Our treatment of the evolution of these conduction fronts follows Borkowski, Balbus, & Fristrom (1990), which should be consulted for a more detailed treatment. Heat conduction proceeds on a time scale of $`\tau _{\mathrm{con}}(5n\mathrm{k}T/2)(\chi T^{7/2}L^2)^1`$, where the first parenthesis is the enthalpy per unit volume, and the second parenthesis is the heat loss per unit volume per unit time because of heat conduction. Here, $`\chi `$ is (nearly) a constant defined such that $`\chi T^{5/2}`$ is the heat conduction coefficient, $`n`$ is the total number density, and $`L`$ the width of the conduction front between the cold cloud and the ICM in its surrounding. Substituting typical values, we obtain

$$\tau _{\mathrm{con}}4\times 10^5\left(\frac{n_e}{0.04\mathrm{cm}^3}\right)\left(\frac{T}{10^7\mathrm{K}}\right)^{5/2}\left(\frac{L}{0.1\mathrm{kpc}}\right)^2\mathrm{yr},$$
(1)
where $`n_e`$ is the electron density. Basically, $`\tau _{\mathrm{con}}`$ is the time required to cool the hot phase gas located inside the conduction front via heat conduction. Our scenario requires that this time be much shorter than the radiative cooling time (see eq. 5 below). For the conduction time to be short, the conduction front must remain narrow. However, the heat conduction front thickness does increase with time. After a time $`t_{\mathrm{br}}`$, the conduction front thickness is (Balbus 1986)

$$L(t_{\mathrm{br}})0.15\left(\frac{t_{\mathrm{br}}}{10^6\mathrm{yr}}\right)^{1/2}\left(\frac{n_e}{0.04\mathrm{cm}^3}\right)^{1/2}\left(\frac{T}{10^7\mathrm{K}}\right)^{5/4}\mathrm{kpc}.$$
(2)
The isobaric radiative cooling time in the temperature range $`2\times 10^5T4\times 10^7`$ K is (see cooling curve in, e.g., Gaetz, Edgar, & Chevalier 1988)

$$\tau _{\mathrm{rad}}2\times 10^8\left(\frac{n_e}{0.04\mathrm{cm}^3}\right)^1\left(\frac{T}{10^7\mathrm{K}}\right)^{3/2}\mathrm{yr}.$$
(3)
To explain the low luminosity of the X-ray emission lines for gas in the temperature range $`10^6T5\times 10^6`$ K by heat conduction, we require $`\tau _{\mathrm{con}}<\tau _{\mathrm{rad}}`$ at these temperatures. This translates into a condition on the width of the conduction front. Let us first consider the non-magnetic case. Assuming the gas is isobaric and using equations (1) and (3), we find the condition to be

$$L<0.28\left(\frac{T_{\mathrm{min}}}{5\times 10^6}\right)^3\mathrm{kpc}.$$
(4)
Here, $`T_{\mathrm{min}}`$ is the temperature below which the gas must cool conductively to avoid detection by X-ray line emission. Thus, heat conduction dominates at those lower temperatures as long as heat conduction occurs in very narrow fronts, with a scale much shorter than the size of the cooling flow region.

Equation (4) is a necessary condition for conduction to dominate for radiative cooling, but it is not sufficient. This condition guarantees that conduction dominates over radiative cooling within the narrow conduction fronts. However, although conductive cooling only occurs within these narrow fronts, radiative cooling goes on everywhere within the cooling flow region. Thus, while conduction might dominate locally, radiative cooling could be more important overall, and the resulting soft X-ray lines might still exceed the observational limits. Let $`V_{\mathrm{con}}`$ be the total volume of the hot gas within conduction fronts (the total volume cooling via heat conduction) at any given time. The conductive cooling time of the hot gas in this volume is $`\tau _{\mathrm{con}}`$ (eq. 1). Let $`V_{\mathrm{CF}}`$ be the total volume of the whole cooling flow region, where the isobaric radiative cooling time is $`\tau _{\mathrm{rad}}`$ (eq. 3). We assume that $`V_{\mathrm{con}}V_{\mathrm{CF}}`$. The condition that conduction cooling dominate radiative cooling for the gas at temperatures of $`T_{\mathrm{min}}`$ is that on average at all times

$$\frac{V_{\mathrm{con}}}{\tau _{\mathrm{con}}}>\frac{V_{\mathrm{CF}}}{\tau _{\mathrm{rad}}}.$$
(5)
For the last equation to hold the following should be true in the environment of cooling flow clusters. First, the conduction at these low temperatures must be fast, which requires narrow conduction fronts (eq. 4). We note here that since we assume the presence of many clouds, the width of the conduction front around any cold cloud can be much smaller than the overall size of the inner cooling flow region, $`L20`$ kpc. Second, the clouds should be small, so that the cooling fronts are narrow, and there should be a large number of them so that a large enough volume is covered, according to equation (5). Third, in order to cool the hot gas via conduction, the magnetic field lines of the cold clouds should be reconnected to those of the ICM. This requires, as stated earlier, reconnection at the head of the cloud.

Usually it is expected that heat conduction is more important at high temperatures than at low temperatures. Thus, one might not expect conduction to play an important role in the coolest X-ray emitting gas detected in clusters, as we are suggesting. It will, however, occur for cool gas if the conduction fronts are narrow. Also, if the heat conduction is regulated by magnetic reconnection, then it may occur preferentially in cooler gas. This might occur because the magnetic pressure can become dynamically important in the cooler inner region. If the magnetic field is isotropically turbulent, then the magnetic pressure grows with density as $`P_Bn_e^{4/3}`$. The ratio of the magnetic to thermal pressure then varies as $`P_B/P_{\mathrm{th}}n_e^{1/3}T^1`$. The lowest temperature inferred from X-ray emission in CF clusters is typically $`T_{\mathrm{min}}T_c/3`$, where $`T_c`$ is the cluster temperature outside the cooling flow region. The density in the inner region where the temperature is $`T_{\mathrm{min}}`$ is $`10`$ times higher than in the region where the temperature is $`T_c`$. Since in an inflow the magnetic pressure can increase much faster than as $`n_e^{4/3}`$ (Soker & Sarazin 1990), the ratio $`P_B/P_{\mathrm{th}}`$ is probably increased by a factor of $`10`$ in going from the general ICM to the inner cooling flow region. Observations also support the presence of strong magnetic fields in the centers of cooling flow clusters (Eilek & Owen 2002). The stretching of ICM magnetic field lines around the cold clouds may increase the ratio $`P_B/P_{\mathrm{th}}`$ by an additional factor. This may be sufficient to allow reconnection between the magnetic field in the ICM and that inside the fast moving clouds. We argue that the requirement for reconnection sets the lower temperature in CF clusters. The magnetic pressure in the ICM is still slightly below the equipartition value, and ram pressure and stretching due to the motions of cold clouds may amplify the field to achieve equipartition near the fronts of the clouds. In this scenario, the negligible heat conduction across cold fronts is because the magnetic field is below its equipartition value, as discussed earlier.

From the discussion above (eq. 4), we conclude that the typical cloud sizes should be $`L0.1`$ kpc. We consider such small cold clouds moving through the ICM. Reconnection is likely to occur in two places, as is seen in the interaction of the Earth’s magnetic field with the Solar wind (e.g., Oieroset et al. 2001). At the leading edge of the cloud, the magnetopause, the magnetic field of the ICM and the cloud’s magnetic field are compressed together, leading to reconnection if the field lines are not aligned. As the cloud moves relative to the ICM, the field lines are stretched behind the cloud, and eventually reconnect at a distance of several cloud radii behind the cloud. For example, in numerical magnetohydrodynamic simulations of the interaction of a cloud with a moving medium, Gregori et al. (2000) find that the reconnection region extends to 5-10 cloud radii. Although the results of Gregori et al. (2000) do not correspond one-to-one with the conditions in clusters (e.g., the clouds are not magnetized in their simulations), we note that they show that magnetic effects tend to be important only if the initial Mach Alfvenic number is sufficiently small. This applies to the situation studied here. Once reconnection occurs and leads to conductive cooling, part of the ICM may be accreted onto the cold cloud or may form another cold cloud. In either case, the frozen-in magnetic fields will be accreted as well, and this will lead to more field lines connecting the hot ICM with cold clouds, and to more conduction. Moreover, the compression of the ICM fields during cooling will lead to very strong magnetic fields within the cold clouds.

In this geometry, and based on the numerical simulations of Gregori et al. (2000), the heat conduction from the ICM to the cloud occurs along the field lines which enclose an ICM volume only $`\mathrm{\Gamma }10`$ times larger than the cloud volume at each moment. For the proposed mechanism to work, we must, however, postulate a much more efficient reconnection, with an ICM volume of $`\mathrm{\Gamma }100`$ times larger than the cold cloud volume. With a temperature ratio of $`10^3`$ and pressure equilibrium, this means that the ICM mass is $`10^3\mathrm{\Gamma }`$ times the cloud mass. The radiative coefficient $`\mathrm{\Lambda }`$ at $`10^7`$ K is about 10 times that for ionized gas at $`10^4`$ K, so the ratio of the luminosity of the cold cloud to that of the ICM connected to it is $`L_{\mathrm{cold}}/L_{\mathrm{hot}}10^5/\mathrm{\Gamma }10^3`$. Most of the luminosity comes in the optical and UV, rather than in the X-ray band. Without magnetic field lines enclosing the interaction region, clouds with a size smaller than the Field length will evaporate (e.g., McKee & Begelman 1990). In the present magnetic field geometry, on the other hand, the ICM is being cooled rather than the cloud being evaporated. This is because at any given time the total mass of the gas in the hot phase which is connected to the cloud is $`0.1`$ times the mass of the cold cloud. (Recall our estimate above, based on the simulations of Gregori et al. 2000, that only a small volume of the surrounding gas is magnetically connected to the cold cloud). Thus, even without any radiative cooling, thermal equilibrium will bring the total mass enclosed within the magnetic field lines to a temperature of only $`0.05`$ times the temperature of the hot phase. Thus, there is not enough thermal energy in the hot phase in contact with the clouds to evaporate them. Again, some non-negligible fraction of the hot-phase energy is still expected to be radiated in the X-ray band. From equation (1), the cold cloud will cool the region connected to it in $`5\times 10^5`$ yr. In a time $`t`$ years, the clouds will cool a mass of $`10^3(t/5\times 10^5)\mathrm{\Gamma }M_{\mathrm{cold}}`$, where $`M_{\mathrm{cold}}`$ is the total mass of the cold clouds in the inner region. The cooling rate via heat conduction is then

$$\dot{M}(\mathrm{conduction})20\left(\frac{\mathrm{\Gamma }}{100}\right)\left(\frac{M_{\mathrm{cloud}}}{10^8M_{}}\right)M_{}\mathrm{yr}^1,$$
(6)
with strong dependence on the width and geometry of the conduction front. For example, parcels of hot-phase gas cooling first, will cool via heat conduction parcels of gas still hot and on the same field lines, i.e., the freshly cooling gas will serve as cold clouds. If large masses of cold dense gas are present in cooling flow clusters (Ferland et al. 2002), then a significant amount of gas can cool via heat conduction. As evident from the last equation, however, we still expect the conductive cooling rates to be much smaller than the radiative cooling rates suggested by pre-Chandra and XMM-Newton X-ray observations, and to be consistent with the amounts of cold matter and rates of star formation observed in the centers of cooling flows (e.g., Soker et al. 2001). The heat conduction cooling rate may become a significant, and even dominant, fraction of the cooling rate of the hot phase in the moderate cooling flow model (Soker et al. 2001).

The implication is that the cooling ICM mainly loses its thermal energy below $`1`$ keV via heat conduction rather than line emission. Some line emission is still expected. Mixing can also lead to a reduction in the radiative cooling rate (Bayer-Kim et al. 2002), and was discussed in a turbulent environment by Begelman & Fabian (1990).

To summarize, combining the moderate cooling flow model and the heat conduction process proposed here, the total line emission expected for gas cooling from $`10^7`$ K to $`10^6`$ K might be factor of $`30100`$ lower than that expected in older versions of the cooling flow model, where only line emission is considered for the cooling process, and where the cooling rates are much higher. This factor comes from a $`10`$ times lower mass cooling rate of the hot gas (Soker et al. 2001), and from another factor of $`310`$ for cooling via heat conduction. This much lower expected radiative cooling rate at low X-ray emitting temperatures is compatible with limits inferred from observations. This limit is $`<20\%`$ of the expected cooling rate in old versions of the cooling flow model (Fabian 2002). Of course, energy conservation implies that eventually the gas cools by radiation. This occurs at much lower temperatures, $`10^4`$ K. Therefore, the large ratio of the number of $`H\alpha `$ photons to the number of cooling hydrogen atoms is explained by this scenario, as well as by mixing (Bayer-Kim et al. 2002, their section 8). From energy conservation, the total luminosity of optical/UV/IR emission from the cold gas is

$`L_{\mathrm{opt}}`$ $``$ $`{\displaystyle \frac{5}{2}}{\displaystyle \frac{kT_{\mathrm{min}}}{\mu m_p}}\dot{M_X}`$ (7)
$``$ $`7\times 10^{42}\left({\displaystyle \frac{T_{\mathrm{min}}}{10^7\mathrm{K}}}\right)\left({\displaystyle \frac{\dot{M_X}}{30M_{}\mathrm{yr}^1}}\right)\mathrm{erg}\mathrm{s}^1.`$
The fraction of this luminosity which comes out in any one line depends on the abundances, temperatures, densities, and excitation conditions in the gas, but roughly one expects a few percent of the luminosity to emerge in H$`\alpha `$ (e.g., Heckman et al. 1989). Thus, the expected H$`\alpha `$ luminosities due to conductive heating are

$$L(\mathrm{H}\alpha )2\times 10^{41}\left(\frac{T_{\mathrm{min}}}{10^7\mathrm{K}}\right)\left(\frac{\dot{M_X}}{30M_{}\mathrm{yr}^1}\right)\mathrm{erg}\mathrm{s}^1.$$
(8)
This is crudely comparable to the optical line luminosities observed from cooling flows (e.g., Heckman et al. 1989).

## 4 Correlation between Optical and X-ray Emitting Gas in Cooling Flow Clusters

A strong spatial correlation between bright regions of X-ray emitting gas and optical ($`H\alpha `$) filaments has been found with Chandra observations of several cooling flow clusters including Hydra A (McNamara et al. 2000), Abell 1795 (Fabian et al. 2001), Abell 2052 (Blanton et al. 2001), and Virgo/M87 (Young et al. 2002). One of the clearest correlations of bright X-ray and $`H\alpha `$ emission is in Abell 2052, where the $`H\alpha `$ emission corresponds with the brightest parts of the shells surrounding bubbles evacuated by the cluster’s central radio source. A temperature map of Abell 2052 (Blanton et al. 2003) shows that these regions have cooled to about 0.8 keV ($`9\times 10^6`$ K), or about 1/4 of the cluster temperature outside of the cooling flow region. The $`H\alpha `$ emission represents gas at a temperature of approximately $`10^4`$ K, and it is possible that cooling between $`10^7`$ and $`10^4`$ K has occurred by conduction.

Currently, the best measurements of mass-deposition rates in cooling flow clusters can be made with high-resolution grating spectroscopy using Chandra and XMM-Newton. These observations are able to resolve individual emission lines expected from gas emitting at a range of temperatures. Peterson et al. (2003) present results from 14 cooling flow clusters observed with the Reflection Grating Spectrometer on XMM-Newton. The mass-deposition rates are measured from $`kT`$ to $`kT/2`$, $`kT/2`$ to $`kT/4`$, $`kT/4`$ to $`kT/8`$, and $`kT/8`$ to $`kT/16`$, where $`kT`$ is the temperature of the cluster outside of the cooling flow region. Mass-deposition rates ($`\dot{M}`$)’s were successfully measured for 12 of the 14 clusters, and of these 12, 9 of them had well-constrained $`\dot{M}`$’s in the $`kT`$ to $`kT/2`$ and $`kT/2`$ to $`kT/4`$ ranges. For the remaining three, only upper limits were possible in the $`kT/2`$ to $`kT/4`$ range. We have used the $`\dot{M}`$’s in the $`kT/2`$ to $`kT/4`$ range from Table 5 in Peterson et al. (2003), and taken $`kT/4`$ as the limiting temperature for each cluster. For all clusters that have measured $`H\alpha `$ luminosities (Heckman et al. 1989, Allen 1995), we have calculated the expected $`H\alpha `$ luminosity ($`L_{H\alpha }`$) from cooling by conduction using equation 8. We were able to make this comparison using ten of the clusters. A comparison of the predicted vs. the measured $`H\alpha `$ luminosities is presented in Table 1, and plotted in Figure 2. The error bars represent 1-$`\sigma `$ errors calculated by the propagation of errors (of the temperature and mass-deposition values). The dotted line represents $`L_{H\alpha }`$(measured) $`=L_{H\alpha }`$(predicted). Points fall on both sides of the line, in other words, equation 8 sometimes over-predicts and sometimes under-predicts the amount of $`H\alpha `$ emission. Still, the points generally follow the line with higher predicted $`L_{H\alpha }`$ values corresponding with higher measured $`L_{H\alpha }`$ values. This is consistent with the correlation of $`\dot{M}`$ and measured $`L_{H\alpha }`$ described in Heckman et al. (1989).

## 5 Summary

We proposed that the upper limit on the mass cooling rate of gas cooling below a temperature of $`10^7`$ K in cooling flow clusters (e.g., Kaastra et al. 2001; Peterson et al. 2001, 2002, 2003) can be made compatible with the existence of cooling flows in the following way. First, one must consider the so called “moderate cooling flow model” (Soker et al. 2001), where the effective age of the ICM is much lower than the age of the cluster. This implies a much lower mass cooling rate, by a factor of $`10`$, than that in the old version of the cooling flow model. Second, the possibility that a large fraction of the gas at $`10^7`$ K cools by heat conduction to lower temperatures, rather than by radiative cooling, further reduces the X-ray luminosity expected from gas below $`T10^7`$ K. With these two processes, the limit on the cooling rate below $`10^7`$ K inferred from X-ray observations, which is $`<20\%`$ of the mass cooling rates cited in the past (Fabian 2002), is easily met.

In the present paper we showed that reconnection between the magnetic field lines in cold ($`10^4`$ K) clouds and the field lines in the ICM, if it occurs, forms a narrow conduction front. Despite the relatively low temperature, the narrow conduction front allows efficient heat conduction from the hot ICM to the cold clouds. We also showed that (eq. 6) a large fraction of the ICM at $`10^7`$ K will cool by heat conduction. The reconnection between the field lines in cold clouds and those in the ICM occurs only when the magnetic field in the ICM is strong enough, i.e., close to the equipartition value. This occurs only in the very inner regions of cooling flow cluster, at $`r1030`$ kpc.

Energy conservation implies that eventually the gas cools by radiation. As with mixing (Bayer-Kim et al. 2002, their section 8), the large ratio of the number of $`H\alpha `$ photons to the number of cooling hydrogen atoms is explained by this scenario. As evident from Figure 2, there is a reasonable agreement between the observed H$`\alpha `$ luminosities and the values we estimated based on our proposed scenario. To a certain extent, the model we have proposed reconciles the standard cooling flow model with the recent observations of weak X-ray line emission from gas below $`10^7`$ K. If conductive cooling to cold clouds is important, these observations might still be in accord with a moderate cooling flow model, and cooling at such lowered rates could explain the existence of cooler gas and young stars in the same central regions of cooling flows, and provide a source of fuel for the central AGNs often found there.

We benefited from discussions with Ray White at early stages of this research. While at the University of Virginia, N. S. was supported by a Celerity Foundation Distinguished Visiting Scholar grant. N. S. was supported in part by grants from the US-Israel Binational Science Foundation. E. L. B. and C. L. S. were supported in part by the National Aeronautics and Space Administration through Chandra Award Numbers GO0-1158X and GO1-2133X, issued by the Chandra X-ray Center, which is operated by the Smithsonian Astrophysical Observatory for and on behalf of NASA under contract NAS8-39073. Support for E. L. B. was provided by NASA through the Chandra Fellowship Program, grant award number PF1-20017, under NASA contract number NAS8-39073.


text:
# Structural and computational depth of diffusion limited aggregation

## Abstract

Diffusion limited aggregation is studied from the perspective of computational complexity. A parallel algorithm is exhibited that requires a number of steps that scales as the depth of the tree defined by the cluster. The existence of this algorithm suggests a connection between a fundamental computational and structural property of DLA.

Diffusion limited aggregation (DLA), introduced some 20 years ago by Witten and Sander Witten and Sander (1981), is a model of pattern formation and an example of self-organized criticality. The dynamical rules for DLA are simple but the patterns produced are complex and have thus far defied full theoretical understanding though there has been recent progress for two-dimensional DLA Barra et al. (2002); Ball and Somfai (2003). DLA models a number of physical systems including electrodeposition, fluid flow in porous media and the growth of bacterial colonies Vicsek (1992).

The stochastic growth rules for DLA can be couched in terms of random walk dynamics. In the present work we consider two-dimensional diffusion limited aggregates composed of circular particles of unit diameter. The initial condition is a single seed particle at the origin. The growth of the aggregate occurs one particle at a time and proceeds by starting the first particle at a random position on the “birth circle,” $`r_1(0)=r_b>1`$. This particle does a random walk, $`𝐫_1(t)`$, until it either drifts out to the “death circle,” $`r_1>r_dr_b`$, or comes in contact with the particle at the origin, $`r_1=1`$. If it contacts the particle at the origin it sticks and the aggregate grows. If it reaches the death circle it is reborn at a random point on the birth circle and the process is repeated until the particle sticks. After $`n`$ steps of this construction, the aggregate consists of $`n`$ connected particles just contained in a circle of radius $`R_n`$. On the $`n^{\mathrm{th}}`$ step, a particle is launched on the birth circle, $`r_b>R_n`$, and follows a random walk $`𝐫_n(t)`$ until it sticks on one of the existing particles in the aggregate at sticking position $`𝐫_n^{}`$ such that there is a $`k<n`$ and $`|𝐫_n^{}𝐫_k^{}|=1`$.

DLA is believed to form fractal clusters such that $`R_NN^{1/d_f}`$ where $`R_N`$ is the average radius of a cluster of $`N`$ particles and $`d_f`$ is the fractal dimension, estimated by numerical simulation to be $`1.715\pm .004`$ Tolman and Meakin (1989). DLA clusters define a graph where the nodes are the particles and the edges are the contacts between particles. Since each particle sticks to a single predecessor, this graph is a tree rooted at the origin. One of the quantities that we will be interested in is the structural depth, $`𝒟_s`$ of this tree, defined as the length of the path from the root to the outermost leaf. Numerical simulations Meakin et al. (1984) support the conjecture that the structural depth grows linearly in the radius of the cluster, $`𝒟_sR`$. The growth mechanism of DLA strongly favors adding particles at the outermost tips of the cluster and thus exerts a tension that tends to make the growth radial. Figure 1A color codes particles in a cluster according to their distance to the seed particle along the tree, sometimes referred to as the chemical distance to the origin. This figure reveals the close relation between chemical distance and Euclidean distance from the origin.

In addition to its structural properties DLA has computational properties that help elucidate its complexity. The growth rules for DLA and the resulting patterns suggest that DLA has a great deal of history dependence. The construction of a DLA cluster requires a long sequence of steps and some random choices made early in the growth process are frozen in and have a large impact on the structure of the aggregate thereafter. Indeed the stochastic history dependence present in DLA can be considered a marker for complexity. Biological evolution, the ultimate process for generating complexity, also has the feature that the present state is emerges from a long sequence of stochastic steps and that accidents occurring in early epochs are frozen in and create the ground rules for later epochs.

The intuitive notion of history dependence in stochastic models such as DLA can be formalized in the context of computational complexity theory. The idea that the length of a history can be measured in the framework of computation was introduced by Bennett Bennett (1988, 1990) though our treatment differs by emphasizing parallel rather than sequential computation and ensembles rather than individuals. We define the computational depth, $`𝒟_c`$ of a statistical physics models such as DLA as the minimum number of parallel computational steps needed to generate a typical system state of the model. Both statistical physics and computational complexity theory are concerned with scaling properties, for example, in the case of DLA, we are interested in the leading large $`N`$ behavior of $`R`$, $`𝒟_s`$ and $`𝒟_c`$.

The notion of computational depth requires that we specify the parallel computer used to simulate the system. Computational complexity theory Papadimitriou (1994); Greenlaw et al. (1995) provides a standard, idealized model of parallel computation called the parallel random access machine or PRAM. The PRAM consists of many simple processors all connected to a single global memory. All processor run the same program synchronously though each processor has a distinct label so that the program can directs different processors to do different calculations. In a single computational step, each processor carries out an elementary logical or arithmetic operation and exchanges information with the global memory. Since we are considering a stochastic system random numbers are needed and we assume a sufficient supply are stored in the global memory before the simulation begins. The PRAM is a massively parallel model of computation where the number of processors and memory cells is allowed to grow polynomially (i.e. as a power) in the size of the system, in our case the number of particles $`N`$. Since any processor in the PRAM can communicate with any memory element in a single time step, the PRAM is not a realistic, scalable model of parallel computation and the considerations of this paper are not intended to provide a practical means of simulating DLA. Instead, the goal is to elucidate whether assembling a DLA cluster requires a long sequence of steps. The PRAM is one of many equivalent models of parallel computation and the scaling behavior for $`𝒟_c`$ reflects an intrinsic property of DLA rather than a particular strength or weakness of the PRAM. With this preamble we can now define computational depth as the average number of steps needed to construct a system state on a PRAM using the fastest parallel algorithm.

The primary objective of this paper is to show that $`𝒟_c𝒟_s`$ and to motivate the conjecture that in fact $`𝒟_c𝒟_s`$ 111The symbols $``$ and $``$ are to be interpreted, respectively, as asymptotic inequality and equality up to logarithmic factors, e.g. $`N^3\mathrm{log}^2NN^3`$.. Previous work has placed upper bounds on the computational depth of various pattern formation processes in statistical physics. For example, it has been shown that the clusters formed by the Bak-Sneppen model Machta and Li (2001), internal DLA Moore and Machta (2000), invasion percolation, the Eden model and ballistic deposition Machta and Greenlaw (1994) all have polylog depth $`𝒟_c\mathrm{log}^{𝒪(1)}N`$, where $`N`$ is the number of degrees of freedom of the system. These results show that these models do not have a strong intrinsic history dependence since constructing typical states can be carried out in a small number of parallel steps. These models generate trees whose depth greatly exceeds the computational depth, $`𝒟_c𝒟_s`$. The situation is more complicated for DLA. It has been shown Machta and Greenlaw (1996) that random walk dynamics for DLA defines a P-complete problem, which strongly suggests that there is no polylog depth parallel construction for DLA, though it does not rule out a better power law than the linear scaling with $`N`$ of the conventional, one walk at a time growth rule. The question posed here is what is the best power law scaling for a parallel simulation of DLA. Previous work Moriarty et al. (1997); Barra et al. (2002) showed that $`𝒟_cN^{0.74}R^{1.26}`$. Here we improve that bound by exhibiting a faster parallel algorithm than the one described in Moriarty et al. (1997).

The parallel algorithm for DLA proposed here assembles the cluster iteratively. On each step, every particle is moved to a tentative position so that the true cluster is approximated with increasing fidelity. Before the iteration is begun, an ordered list of $`N`$ sufficiently long random walk trajectories $`\{𝐫_i(t)|i=1,\mathrm{},N\}`$ are generated and stored in memory 222Random walk trajectories may include several restarts from the birth circle each time the death circle is reached.. The final cluster generated by the parallel algorithm is the same as the cluster that would result using these random walks trajectories and the standard sequential rules for DLA.

In each step, all particles are moved along their trajectories to tentative sticking points defined by a temporary cluster called the semi-secure cluster. On the first step, the semi-secure cluster consists only of the seed particle and each particle sticks at the point where it first contacts the seed particle, ignoring overlaps with other particles. The semi-secure cluster for the second step consists of the seed particle and any particle whose path to the seed did not cross the sticking point of any lower-numbered particle. On the second step every particle is moved independently along its trajectory until it first contacts this semi-secure cluster. Note that it is possible that a given particle on a given step may be in the semi-secure cluster in one location but stick in a different location.

At the beginning of the $`m^{\mathrm{th}}`$ step in the algorithm, we have a semi-secure cluster defined by $`N^{(m)}`$ particles at locations $`\{𝐬_j^{(m)}|jS^{(m)}\}`$ where $`S^{(m)}`$ is the set of semi-secure particle indices at step $`m`$. Each particle is then moved to a location based on the template provided by this semi-secure cluster. The sticking point of a particle is determined by independently moving it along its trajectory until it first contacts a semi-secure particle location with an index lower than its own. $`𝐫_i^{(m)}`$ is the sticking point of particle $`i`$ at step $`m`$ if it is the first time $`\tau _i^m`$ along the trajectory such that there is contact with a semi-secure particle $`j<i`$ and $`|𝐫_i(\tau _i^m)𝐬_j^{(m)}|=1`$.

After particles are moved to their sticking points the current semi-secure cluster is discarded and the locations of the sticking points used to determine the next semi-secure cluster. Each particle may be categorized as “semi-secure” or “not semi-secure.” A particle $`i`$ is semi-secure at the beginning of step $`m+1`$ if the path to its sticking point in step $`m`$ does not intersect the sticking point of any predecessor particle on step $`m`$. More formally, $`iS^{(m+1)}`$ and $`𝐬_i^{(m+1)}=𝐫_i^{(m)}`$, if for all $`j<i`$ and all $`t<\tau _i^m`$, $`|𝐫_i(t)𝐫_j^{(m)}|>1`$. A semi-secure particle such that all its predecessors are also semi-secure is called secure. It is easy to see that a secure particle is in the correct location in the cluster. The parallel algorithm has constructed a correct DLA cluster of size $`M`$ for the given trajectories as soon as particle $`M`$ is secure.

There is one additional rule for constructing the semi-secure cluster. It may happen that the semi-secure cluster defined by the above rules is a multiply connected “forest” composed of several “trees.” In this case, only the tree rooted at the seed particle is retained and the remaining particles, not connected to the seed, are removed from the semi-secure cluster.

As discussed in detail in Moriarty et al. (1997), random walk trajectories can be generated, sticking points of these trajectories to an existing cluster found and interferences among particles identified in polylog parallel time using polynomially many processors. Using a parallel algorithms for connected components Gibbons and Rytter (1988), one can identify the the tree connected to the origin in polylog time. Thus, the setup stage and each step of the algorithm require polylog time so that, up to logarithmic factors, the parallel time used by the algorithm scales as the number of steps.

How well does this algorithm perform? It is easy to see that at least one new particle becomes secure on each iteration so that the parallel algorithm requires no more than $`N`$ iterations. In fact, the performance is much better than this weak bound. We have simulated the parallel algorithm on a sequential computer for $`N`$ up to 20,000. Let $`T`$ be the number of iterations required by the algorithm and let $`\kappa =𝒟_s/T`$ be the ratio of the number of iterations divided by the structural depth. Figure 2 shows $`\kappa =𝒟_s/T`$ plotted against $`T`$. The average is over several thousand clusters for the small sizes and 63 cluster for $`N=20,000`$. Although $`\kappa `$ decreases slightly as $`T`$ and $`N`$ increase, the data strongly suggests that $`\kappa `$ approaches an asymptote near $`0.9`$. Thus, up to logarithmic factors, the computational depth of DLA is no greater than the structural depth. In terms of a dynamic exponent for the algorithm, defined by $`TR^z`$, the data suggests that $`z=1`$.

The algorithm of Moriarty et al. (1997) is similar to the present algorithm except that the template used in step $`m+1`$ is the secure cluster from step $`m`$ rather than the semi-secure cluster. The semi-secure cluster is larger than the secure cluster but, unlike the secure cluster, it contains particles at incorrect positions which must later be moved. Our numerical results show that the larger size of the semi-secure cluster more than compensates for the errors in the semi-secure cluster. The dynamic exponent of the “secure” algorithm of Moriarty et al. (1997) is $`1.26`$ compared to $`1`$ for the “semi-secure” algorithm described here.

If $`\kappa `$ is indeed asymptotically near one, it implies that on most steps of the parallel algorithm, one new level is added to the tree defined by the cluster. It is obvious that the algorithm cannot add more than one level to the tree in one parallel step so $`\kappa 1`$. It is instructive to consider the way that $`T`$ becomes larger than the structural depth. After each parallel step, every semi-secure particle can be assigned a deficit. A particle’s deficit is defined as the difference between the step on which it most recently became semi-secure at it’s current location and its chemical distance from the origin. Figure 3 shows how the deficit of a particle can increase by an arbitrary amount. The two panels show two successive parallel steps. Particle $`b`$ has just become semi-secure in the step preceding panel 1. Trajectories and sticking points for particles $`c`$, $`d`$ and $`e`$ are shown. The ordering for the particles is assumed to be $`b<c<d<e`$ and $`a<d`$. Before step 1, all particles except $`c`$ are semi-secure. Suppose that $`e`$ had become semi-secure at its current location $`k`$ steps before step 1. In panel 1, particle $`c`$ becomes semi-secure and interferes with $`d`$ so $`d`$ is no longer semi-secure. In panel 2, both $`d`$ and $`e`$ become semi-secure in new positions. The deficit of $`e`$ is increased by $`k+2`$ over what it was before panel 1 since its chemical distance is decreased by one and it has become semi-secure $`k+1`$ steps later. When the algorithm is finished and all particles are secure, deficits along branches of the tree are a non-decreasing function of chemical distance. The maximum, taken over leaves of the tree, of the sum of the deficit and chemical distance to the origin gives the running time of the algorithm in parallel steps.

Figure 1 shows three images of the same cluster colored to reveal different properties. This cluster consists of 20,000 particles, it has a structural depth of 439 and requires 479 steps to assemble. In 1A the particles are colored according to their structural depth (chemical distance from the origin). Note that the contours of equal structural depth are nearly circular. In 1B particles are colored according to the parallel step on which they first become semi-secure at their final locations. This image reveals that growth generated by the parallel algorithm conforms more to the shape of the branches of the cluster. In 1C particles are colored according to their deficits, 1C shows the difference between the quantities in 1B and 1A and it reveals that the dominant branches have small deficits but that large deficits can develop in branches which are less robust and screened by the dominant branches. The reason for this tendency can be seen in the example of Fig. 3 where a theft from one branch to another increases the deficit on the branch that loses the particle.

The algorithm proposed here assembles DLA clusters by adding nearly one level of structural depth in each parallel step. Our intuition is that it is not possible to generate DLA clusters in substantially fewer steps than the structural depth of the cluster so that $`𝒟_c𝒟_s`$. This conjecture links fundamental structural and computational properties of DLA. Like other lower bounds in computational complexity theory, this conjecture is likely to be difficult to prove.

###### Acknowledgements.

This work was supported by NSF grants DMR-9978233 and DMR-0242402.


text:
# Non-parabolicity of the conduction band of wurtzite GaN

## Abstract

Using cyclotron resonance, we measure the effective mass, $`m`$\*, of electrons in AlGaN/GaN heterostructures with densities, $`n_{2D}16\times 10^{12}`$cm-2. From our extensive data, we extrapolate a band edge mass of $`(0.208\pm 0.002)m_e`$. By comparing our $`m`$* data with the results of a multi-band k.p calculation we infer that the effect of remote bands is essential in explaining the observed conduction band non-parabolicity (NP). Our calculation of polaron mass corrections – including finite width and screening - suggests those to be negligible. It implies that the behavior of $`m`$*$`(n_{2D})`$ can be understood solely in terms of NP. Finally, using our NP and polaron corrections, we are able to reduce the large scatter in the published band edge mass values.

The magnitude of the conduction band non-parabolicity (NP) in wurtzite GaN currently remains controversial. NP of a band can be probed by measuring the carrier effective mass $`m`$* as a function of energy. Such experiments have been performed in the past in both bulk and in two-dimensional electron systems (2DES). The deduced band edge mass values, $`m_0`$\*, however, exhibit considerable scatter. Using cyclotron resonance (CR), Drechsler et al. Drechsler et. al. (1995) determined $`m_0`$* in bulk wurtzite GaN to be 0.20$`m_e`$, where $`m_e`$ is the free electron mass. Other methods such as infrared reflectivity on electron plasma Perlin et. al. (1996) and spectroscopy on shallow donors Moore et. al. (1997); Witowskia et. al. (1999); Meyer et. al. (1995) in bulk GaN have yielded 0.220$`m_e<m_0^{}<`$ 0.236$`m_e`$. An even wider range of values for the band edge mass, 0.185$`m_e<m_0^{}<`$ 0.231$`m_e`$, emerge from experiments in AlGaN/GaN heterostructures. From the temperature dependence of Shubnikov-de Haas (SdH) oscillations of 2DES, Lin et al. deduced $`m_0`$* = 0.22$`m_e`$ Lin et. al. (1998) while Hang et al. reported $`m_0`$* = 0.185$`m_e`$ Hang et. al. (2001). Cyclotron resonance experiments on heterostructures have revealed 0.223$`m_e<m_0^{}<`$ 0.231$`m_e`$ Knap et. al. (1996, 1997). This spread in $`m_0`$* suggests that the “extrapolation” from the various experiment values $`m`$*$`(E)`$ to the band edge remains poorly controlled.

Using CR we have measured $`m`$* in a series of high mobility ($``$20,000 cm2/Vsec) AlGaN/GaN structures. Our heterostructures are all grown by molecular beam epitaxy (MBE) on GaN templates prepared by hydride vapor phase epitaxy (HVPE). The specimens are described in detail elsewhere Syed et. al. (2003); Manfra et. al. (2002). Our data cover a density range of $`16\times 10^{12}`$ cm-2. In these 2D systems this implies energies from $``$27meV to $``$120meV above the band edge due to electron confinement and band filling. Therefore, our mass data probe the NP of the GaN conduction band in small steps over a wide energy range. We also perform extensive k.p calculations and determine that instead of the commonly used two-band model, a multi-band model is required to explain our experimental results. Additionally, our calculations on polaron correction of the effective mass of 2D electrons in GaN show them to be at the 1$`\%`$ level, considerably less than previously thought Knap et. al. (1996). A detailed comparison between our data and k.p calculation sets the band mass value to $`m_0^{}=(0.208\pm 0.002)m_e`$ and, when applied to other investigators’ results, considerably reduces the spread in band edge mass values.

A Fourier transform spectrometer with light pipe optics and a composite Si bolometer was used for the detection of far-infrared transmission. Magnetic field was applied normal to the 2D electron layer. The carrier density of each sample was determined in situ from the Shubnikov-de Haas (SdH) oscillations of the 2DES. All CR and SdH experiments were conducted at 4.2K.

Fig.1 shows the cyclotron resonance energies versus magnetic field, B of a sample with $`n_{2D}=2.3\times 10^{12}`$ cm-2. All data are taken with a resolution of 0.24 meV. The inset shows high field transmission spectra normalized to the spectrum taken at B=0T. The solid line in Fig.1 is a fit to the data at high fields, B$`>`$27T, and low fields, B$`<`$12T, resulting in an effective mass of $`m`$* = 0.228$`m_e`$. For 15T$`<`$B$`<`$25T, there is a pronounced deviation from this straight line. This anomaly in the CR represents a recent discovery, which is being analyzed and published elsewhere Syed et. al. (2003). Here, we observe that this anomaly is limited to a finite field region outside which all CR data can be fit by a straight line.

We have measured the effective mass in eleven samples with carrier density, $`n_{2D}`$, ranging from $`16\times 10^{12}`$ cm-2. In all cases, we observed either a broadening or a splitting of the CR line at intermediate fields but could fit the data away from this regime as well as the data seen in Fig. 1. Fig.2 shows the dependence of $`m`$* on $`n_{2D}`$. For comparison, we also plot data from Ref. Knap et. al. (1996, 1997); Wang et. al. (1996). The mass data from the different references are located in the general vicinity of our results but, due to their considerable error bar or sparsity, are difficult to extrapolate to zero density. The combined data of Fig.2 show an increase in $`m`$* by $``$17$`\%`$ as $`n_{2D}`$ changes from $`19\times 10^{12}`$ cm-2. The rise in $`m`$* with $`n_{2D}`$ reflects the non-parabolicity of the conduction band of the GaN host. Simply extrapolating our closely spaced data linearly to vanishing $`n_{2D}`$, we arrive at $`m_0^{}`$ =0.214$`m_e`$. This value is about $``$8$`\%`$ lower than previously published CR data Knap et. al. (1996, 1997). Since our data contain a small error bar and extend to very low $`n_{2D}`$, such a simple extrapolation should already be quite reliable.

Other groups have previously addressed NP in GaN. For example, Knap et al. explored NP using CR with different $`n_{2D}`$ (see our Fig.2) and accounted for the magnitude of NP using a simple two-band approximation Knap et. al. (1996, 1997). In such an approximation, which includes only coupling between the lowest conduction and highest valence bands, the effective mass varies as

$$m^{}(E)=m_0^{}(1+2KE/E_g)$$
(1)
with $`K`$=1 and $`E_g`$=3.5 eV. In a 2DES the energy, $`E=E_k+E_F`$, above the band minimum is composed of the average kinetic energy, $`E_k`$, of the electrons in the confining potential well and $`E_F=\pi \mathrm{}^2n_{2D}/m_0^{}`$ mas , the Fermi energy of the 2D electron system. The value of $`E_k`$ is dependent on the form of the wave function of the confined electrons. In a simple triangular potential approximation, the average kinetic energy is $`E_k=E_c/3`$, with $`E_c`$ being the confinement energy of the lowest subband Ando et. al. (1982). Using the more accurate Fang-Howard variational wave function gives $`E_k=\mathrm{}^2b^2/8m_0^{}`$, where $`b^3=48\pi m_0^{}e^2(N_{dep}+\frac{11}{32}n_{2D})/ϵ\mathrm{}^2`$ Ando et. al. (1982). In our samples, since both the MBE and the HVPE GaN are n-type, the depletion layer density, $`N_{dep}`$, can be set to zero. We use this two-band model with only one adjustable parameter, $`m_0^{}`$, to fit the high-density data of Ref. Knap et. al. (1997) (dashed line in Fig.2). $`E_k`$ was computed in the triangular approximation following the authors of Ref. Knap et. al. (1997). Clearly, although this procedure can describe the original density dependence of $`m`$* of Ref. Knap et. al. (1997) due to the rather large error bars, it fails to account for our data. Neither a simple vertical shift of the line (a different $`m_0^{}`$) nor the usage of the Fang-Howard wavefunction can resolve this discrepancy. What is required, is a much stronger dependence of $`m`$* on $`n_{2D}`$.

Empirically, we pursue an approach taken by Singleton et al. Singleton et. al. (1988), who used a *modified* two-band model to describe their NP data in GaAs. The authors considered $`K`$ in Eq.1 to be a second fitting parameter. The inclusion of a variable $`K>1`$ into the analytic expression incorporates the influence of higher conduction bands, simulating the results of a more elaborate, multi-band k.p calculation Pfeffer and Zawadzki (1990). Working with the Fang-Howard model to determine $`E_k`$ we find a very good fit to our data for $`K`$ = 2.5 and $`m_0^{}=(0.208\pm 0.002)m_e`$. Even if we use the less reliable triangular approximation, the required $`K`$ = 1.9. The effectiveness of the modified expression in fitting the density dependence of $`m`$* over a wide range of $`n_{2D}`$ demonstrates that the conduction band of GaN is more non-parabolic than was previously assumed Lin et. al. (1998); Knap et. al. (1996, 1997); Hang et. al. (2001).

Before concluding that the proposed NP model appropriately describes our CR data, we need to assure ourselves that polaron effects are a negligible contributor to $`m`$\*. This mass enhancement factor results from electron-LO phonon coupling in polar semiconductors, such as GaN. A determination of a polaron mass enhancement in heterostructures requires inclusion of screening in 2D and the finite width of the electronic wave function Pau , since both greatly reduce interaction between 2D carriers and LO phonons DasSarma (1983). Following Ref. DasSarma (1983), we calculated the polaron effective mass in AlGaN/GaN heterostructures using a a Frohlich constant, $`\alpha `$ = 0.49 Drechsler et. al. (1995), Fang-Howard variational wave function, and a static Thomas-Fermi screening model. We find a polaron enhancement of less than 1$`\%`$ for $`m`$* for $`n_{2D}=9\times 10^{12}`$cm-2. Since the effect decreases with decreasing density, corrections for lower density specimens are smaller yet. This mass enhancement is considerably smaller than the 10$`\%`$ estimated previously for $`n_{2D}=3.1\times 10^{12}`$cm-2 Knap et. al. (1996), where screening and finite width had been neglected. A 1$`\%`$ mass enhancement due to polaronic coupling lies within the error bars of our CR data.

The large spread in the value of the band edge mass $`m_0^{}`$ in the literature is mostly due to two reasons: the underestimation of the NP and the overestimation of polaronic corrections. Applying our NP and polaron corrections to the available effective mass data (Ref. Drechsler et. al. (1995); Knap et. al. (1996, 1997); Hang et. al. (2001); Lin et. al. (1998); Perlin et. al. (1996); Moore et. al. (1997); Witowskia et. al. (1999); Meyer et. al. (1995); Wang et. al. (1996)), we reach a much more coherent picture for the band edge effective mass, $`m_0^{}`$, in GaN, as shown in the inset to Fig. 2. We observe that the majority of the values for $`m_0^{}`$ are very close to an average $`m_0^{}`$ = 0.204$`m_e`$. However, most of the $`m_0^{}`$ data from SdH (displayed in stars) remain at variance from the CR, infra-red reflectivity and donor spectroscopy data, a fact that remains unexplained.

Since $`m`$*$`(E)`$ cannot be accurately represented with a two-band model, a five-band k.p calculation Pfeffer and Zawadzki (1990) in the zincblende approximation was performed to model our data. Figure 3 shows the results. The material parameters employed using the Koster notation are the spin-orbit splitting, $`\mathrm{\Delta }_0^{}`$, of the $`\mathrm{\Gamma }_5`$ conduction band and the momentum matrix elements, $`P^2`$ and $`\lambda ^2P^2`$, coupling the $`\mathrm{\Gamma }_1`$ conduction band with the $`\mathrm{\Gamma }_5`$ valence and conduction bands. They are not precisely known. The formalism proposed by Carlos et al. Carlos et. al. (1993) and used by Bayerl et al. Bayerl et. al. (2001) in a 5-band model to relate the parameters was utilized to calculate $`m`$*$`(E)`$. However, as seen in Fig.3, the parameters chosen by Bayerl et al. lead to a non-parabolicity even lower than in a two-level model. We find that the experimental data points can only be matched if an extra parameter $`C`$, taking remote bands into account at $`k`$=0, is included Hermann and Weisbuch (1977). If one assumes $`C`$= -1.5 (compare with $`C`$=-2 for GaAs Pfeffer and Zawadzki (1990)), $`P^2`$ is in the range 25-27 eV with $`\lambda ^2`$0.33-0.48 and $`\mathrm{\Delta }_0^{}`$120-180 meV. This set of parameters is close to the “best” set of values chosen by Kennedy et al. Kennedy et. al. (1999) and confirm that $`P^2`$ is relatively large for GaN ($`P^2`$26 eV).

In conclusion, our CR experiments set the conduction band edge mass in GaN to $`m_0^{}=(0.208\pm 0.002)m_e`$. Using our determination of the nonparabolicity and reviewing the polaron correction we reach much better agreement between several published data for $`m_0^{}`$.

###### Acknowledgements.

We thank M.S. Brandt, L.N. Pfeiffer, K.W. West, and W.I. Wang for helpful discussions. A portion of the work was performed at the National High Magnetic Field Laboratory, which is supported by NSF Cooperative Agreement No.DMR-0084173 and by the State of Florida. Financial support from the W. M. Keck Foundation is gratefully acknowledged.