text:
# FEEDBACK HEATING IN CLUSTER AND GALACTIC COOLING FLOWS

## 1 Introduction

Recent Chandra and XMM observations of hot gas in clusters and galaxies have radically altered our previous models for cooling flows and presented new possibilities for understanding them in a different way. X-ray spectra taken with the XMM RGS (Reflection Grating Spectrometer) fail to show line emission from ions having intermediate or low temperatures, implying that the cooling gas is somehow hidden from view or that the cooling rate is at least 5 or 10 times less than previously assumed (e.g. Peterson et al. 2001; Tamura et al. 2001; Kaastra et al. 2001; Xu et al. 2001; Sakelliou, et al. 2002). XMM-EPIC observations also fail to detect cooling gas (Molendi & Pizzolato 2001; Böhringer et al. 2002). However, X-ray images reveal that the hot gas inside the E or cD galaxies located at the centers of cooling flows is often highly disturbed (e.g. Böhringer et al. 1993; Fabian et al 2000; McNamara et al. 2000; Loewenstein et al. 2000; Blanton et al 2001; Jones et al. 2002; Forman et al. 2001; Buote et al. 2002; Trinchieri & Goudfrooij 2002). Evidently, massive black holes, thought to lie at the centers of all stellar bulges, become intermittently active possibly stimulated by the inflow of hot gas. Of particular interest are the cavities in the hot gas that are often (but not always) associated with radio lobes as in Perseus/NGC 1275 (Fabian et al 2000) and HydraA/3C295 (McNamara et al. 2000; David et al. 2001; Allen et al 2001). Cavities occur when ultrahot or relativistic gas displaces the hot thermal gas. Since the rims surrounding the X-ray cavities are typically no hotter than other ambient gas (Fabian 2001; McNamara 2001), strong shocks do not seem to be involved in producing the holes so the innerhole gas is in approximate pressure equilibrium with gas in the rims (Nulsen et al 2002; Soker, Blanton & Sarazin 2002; Brighenti & Mathews 2002c). Obviously the holes must be buoyant (Churazov et al. 2001) and their formation and subsequent evolution must feed energy into the cooling flow gas.

Can the heating visible in Chandra images explain the absence of cooling in XMM spectra? This question can be answered only from studies of the global effect of many heating episodes over many Gyrs. Recently we presented a number of evolutionary models of so-called cooling flows that were heated by a variety of scenarios initiated by gas flowing into the central black hole (Brighenti & Mathews 2002b). None of the flow models we considered were satisfactory; whenever the heating was sufficient to reduce the rate that gas cools near the center, the flow was highly disrupted and spatially distributed cooling occurred at larger radii at essentially the same rate as in the unheated flows.

In the light of several recent developments – including the significant continued interest in heated flows – we are encouraged to reconsider such models again here. Our first motivation is the argument by Narayan & Medvedev (2001) that thermal conductivity is reduced only by $`0.2`$ in the presence of tangled magnetic fields which are expected in cooling flows. Secondly, in a recent paper Ruszkowski & Begelman (2002) demonstrated that heating and thermal conduction can conspire to produce successful cluster cooling flows that have the characteristic positive temperature gradient at small radii, a common observed feature in all these flows. Finally, in our models we follow the prescription of Ruszkowski & Begelman and other authors (e.g. Binney & Tabor 1995; Tucker & David 1997; Ciotti & Ostriker 2001) by considering cluster and galactic flows in which the mechanical energy released by the central AGN is exactly proportional to the mass of gas that flows into the center of the flow, optimizing the effect of the energy feedback.

In the following models we consider the long term evolution of hot gas on both cluster and galactic scales but with different approaches. Our models are intended to be generic, but as a guide to their success we compare them with recent observations of Abell 1795 and NGC 4472 for cluster and galactic flows respectively. Since cluster flows are occasionally disturbed by merging events, we begin with a quiescent post-merger flow initially in hydrostatic equilibrium and follow its evolution for several Gyrs under the influence of radiative cooling, AGN feedback heating in the core, and thermal conduction. For galactic flows, we begin with a model elliptical having an additional complement of circumgalactic hot gas and follow its evolution for 12 Gyrs under the influence of radiative cooling, AGN heating and conduction as well as the additional effects of Type Ia supernovae and stellar mass loss appropriate to galactic flows. We find that it is possible to get marginally satisfactory solutions for cluster flows provided the parameters are chosen within rather narrow limits. However, our models of galactic cooling flows with the same choice of parameters are inadequate or require additional possibly unsatisfactory assumptions.

## 2 Computational Methods

The calculations described below are based on the same flow equations described in detail in our previous work on heated cooling flows (Brighenti & Mathews 2002b) and the reader is referred to that paper for details. We assume no distributed radiative cooling, but this may occur naturally in convective regions when mass elements collide and compress. In such instances we remove the gas from the flow when the temperature has dropped to $`T_{cool}=5\times 10^5`$ K as described in Brighenti & Mathews (2002b). Results of 1D and 2D are in essential agreement when applied to the same flow (see Brighenti & Mathews 2002b) although 2D flows are essential to describe flows with off-center heating.

The grid spacing we use varies with the overall scale of the flow and the number of dimensions considered. 2D Cluster: For these flows we use $`250\times 250`$ cylindrical zones. The grid is uniform at 2.5 kpc spacing for the inner $`150\times 150`$ zones extending to $`375\times 375`$ kpc. Beyond this region the zone size increases geometrically to 2 Mpc. We find very similar results with a $`490\times 490`$ grid where the inner $`300\times 300`$ zones are 0.4 kpc wide. 1D Cluster: These calculations are done with 360 zones increasing in size from 0.5 kpc at the center to an outer boundary at 3 Mpc. 2D Galaxy/Group: These flows are calculated with $`280\times 280`$ cylindrical zones. The grid is uniform at 0.25 kpc spacing for the inner $`200\times 200`$ zones that extend to $`50\times 50`$ kpc. Beyond this region the zone size increases geometrically to 1.7 Mpc. 1D Galaxy/Group: Here we use 350 zones increasing in size from 0.05 kpc at the center to an outer boundary at 1.6 Mpc.

## 3 Cluster Flows with AGN Heating, Cooling and Conduction

Our cluster flow calculations begin in hydrostatic equilibrium with temperature and density profiles based on the well-observed cluster Abell 1795 assumed to be at a distance of 243 Mpc. Abell 1795 is a typical rich cluster with a central cD galaxy and reasonably relaxed overall structure (Boute & Tsai 1996). Abell 1795 has the usual attributes of normal cooling flows: strong central peak in X-ray surface brightness (e.g. Tamura et al. 2001), a central radiative cooling time $`3\times 10^8`$ yrs that is much less than the cluster age (e.g. Edge et al. 1992; Fabian et al. 2001), optical line emission near the central cD (Cowie et al. 1983), an excess of blue and ultraviolet light from massive young stars (Johnstone, Fabian & Nulsen 1987; Cardiel, Gorgas & Aragon-Salamanca 1998; Mittaz et a. 2001) and a central radio source 4C 26.42 (McNamara et al. 1996a,b). Near the center of Abell 1795 Chandra images reveal an X-ray emission feature aligned with an optical filament (Fabian et al. 2001). This and the central total mass profile, $`Mr^{0.6}`$ inside 40 kpc, which is somewhat flatter than NFW, suggest a possible deviation from hydrostatic equilibrium. Deprojections of ROSAT images of Abell 1795 led Allen et al. (2000) to conclude that the total cooling rate is $`\dot{M}500`$ $`M_{}`$ yr-1. However, the XMM RGS spectrum showed no evidence of gas with temperatures less than $`2`$ keV, corresponding to a much smaller upper limit to the cooling rate, $`\dot{M}<150`$ $`M_{}`$ yr-1. This is consistent with a rate $`\dot{M}100`$ $`M_{}`$ yr-1 estimated from Chandra observations (Ettori et al. 2002).

In the upper two rows of Figure 1 we show electron density and temperature profiles of Abell 1795 determined from observations with XMM (filled triangles; Tamura et al. 2001) and Chandra (open triangles; Ettori et al. 2002). The Chandra temperatures have been spatially deprojected and are therefore slightly lower. We inserted the observed density and temperature profiles into the equation of hydrostatic equilibrium to determine the mass distribution. In this process we ignore the possible deviations from equilibrium suggested by the Chandra observations. Nevertheless, our results described below are independent of the particular mass profile assumed which could have been strictly NFW for example. Our adopted fits $`n(r)`$ and $`T(r)`$ to the observed profiles are shown as dot-dashed lines in the first column of Figure 1. The central density and temperature in these fits are $`n(0)=0.1`$ cm-3 and $`T(0)=2.5\times 10^7`$ K.

### 3.1 Unheated Cluster Flow

In the first column of Figure 1 we illustrate the evolution of Abell 1795 as determined by radiative losses only with no AGN or conductive heating. Our calculations are done in 2D and the results in Figure 1 are azimuthally averaged. For standard cosmological parameters ($`\mathrm{\Omega }=0.3`$; $`\mathrm{\Lambda }=0.7`$; $`H=70`$ km s-1 Mpc-1) large clusters like Abell 1795 were formed only recently, so we consider their evolution for only 7 Gyrs. As gas within $`30`$ kpc cools, a subsonic inflow develops to reestablish approximate equilibrium, forcing the central density to rise above the observations. In the bottom row of Figure 1 we plot the evolution of the cooling rate $`\dot{M}(t)`$ which increases from $`\dot{M}100`$ $`M_{}`$ yr-1 at the beginning of the calculation toward an asymptotic value of $`300`$ $`M_{}`$ yr-1 after $`4`$ Gyrs when a quasi-steady inflow is established. Not only does this $`\dot{M}(t)`$ exceed the cooling rate allowed by XMM spectra, all the cooling occurs at the very center where $`3\times 10^{11}`$ $`M_{}`$ accumulates every Gyr. Such an unacceptably huge mass would double the mass of the central cD in a few Gyrs. The results of each evolutionary calculation are summarized in Table 1.

### 3.2 Cluster Flows with Heating

In the next series of models we consider the effect of central and off-central AGN heating. The heating power $`L_h=\epsilon \dot{M}(0)c^2`$ erg s-1 is proportional to the rate $`\dot{M}(0)`$ that gas flows into the central grid zone. The efficiency $`\epsilon `$ is a combined measure of both the AGN power generated and the fraction that is delivered to the hot cluster gas. As in our previous paper (Brighenti & Mathews 2002b) we assume that the heating is instantaneously transferred to the hot gas near the flow center and that the degree of heating has a Gaussian profile $`\mathrm{exp}[(r/r_h)^2]`$. The spatial scale of the heated region $`r_h`$ has only a secondary influence on our results; we restrict our discussion here to $`r_h=25`$ kpc. Our two dimensional (2D) calculations of heated flows create buoyant bubbles that transport energy to larger radii where they eventually dissolve into the background. The merging of the hot bubbles with the background flow is probably an artifact of numerical diffusion in our code, but a similar dissipation may be naturally accomplished by thermal conduction at the bubble surfaces, where the surface area is greatly enlarged due to Rayleigh-Taylor instabilities. We do not attempt to follow the detailed bubble physics here, but only ensure that the injected energy is globally conserved.

The second and third columns in Figure 1 shows the consequences of Gaussian feedback heating on the 2D flow in the first column using $`\epsilon =10^4`$ and $`10^3`$ respectively. When $`\epsilon <10^4`$, the flow evolves in a manner similar to the unheated flow in the first column. Irregularities in the density and temperature profiles are introduced by the heating intermittancy. With $`\epsilon =10^4`$ the density and temperature profiles after 4 and 7 Gyrs are almost acceptable, but $`\dot{M}`$ (third row of Fig. 1) is too large, $`\begin{array}{c}>\hfill \\ \hfill \end{array}150`$ $`M_{}`$ yr-1 beyond 4.5 Gyrs. When the heating efficiency is increased to $`\epsilon =10^3`$ (third column of Fig. 1), the central cooling rate, $`\dot{M}<20`$ $`M_{}`$ yr-1, is acceptable. However, for this higher efficiency the gas is much too hot within the heated region and resembles no known cluster.

The results of similar evolutionary calculations with off-center heating (similar to a one-sided jet) are shown in the final two columns of Figure 1. The Gaussian heating is displaced 50 kpc from the flow center along the $`z`$-axis, $`(z,R)_h=(50\mathrm{kpc},0\mathrm{kpc})`$, but all other parameters are identical to the centrally heated flows. Figure 1 shows azimuthally averaged profiles of 2D calculations. The $`\epsilon =10^4`$ case has a marginally acceptable central cooling rate, $`\dot{M}\begin{array}{c}<\hfill \\ \hfill \end{array}30`$ $`M_{}`$ yr-1, but the gas temperature profile within $`100`$ kpc is flatter than normal. As before $`\epsilon =10^3`$ lowers $`\dot{M}`$ to acceptable values but clearly overheats the gas.

In summary, both heated and unheated flows deviate from typical cluster profiles after a few Gyrs. When the heating is large enough to satisfy the upper limits on $`\dot{M}`$ imposed by XMM spectra, the flow profiles are too distorted to be acceptable. There is only a very small range of heating efficiencies for which these problems are minimized (but not removed). These general conclusions also apply to similar calculations with different heating scales $`r_h`$. Finally, all the gas cools at or near the origin in these heated flows (Table 1).

### 3.3 Cluster Flows with Conduction and Heating

Since the temperature gradient is negative throughout the observed region in Abell 1795, the central parts of the flow are conductively heated from the outside. However, just after the center of the flow is heated with AGN energy, thermal conduction rapidly transfers this energy outward toward larger radii. We first discuss the influence of conduction alone and then the combined effect of conduction plus heating on the evolution of Abell 1795. Because of the shorter computational time steps when conduction is included, we describe the evolution of these flows in 1D. Based on examples discussed in our previous paper (Brighenti & Mathews 2002b), we are confident that the 1D calculations are trustworthy since 1D and 2D models are in essential agreement in every case that we have tested.

Thermal conductivity in a hot plasma, $`\kappa =1.84\times 10^5(\mathrm{ln}\mathrm{\Lambda })^1T^{5/2}`$ erg/sec cm K can be important at high temperatures, but is reduced by magnetic fields. The effect of magnetic suppression on the thermal flux $`F_{cond}=f\kappa (dT/dr)`$ is represented with a coefficient $`f1`$ which is assumed to be uniform in the flow. It is well known that thermal conduction perpendicular to a uniform magnetic field is lowered by many orders of magnitude, so it has generally been assumed that $`f1`$ applied to any plasma with a reasonably tangled magnetic field. However, Narayan & Medvedev (2001) have recently shown that $`f0.2`$ is appropriate for thermal conduction in a hot plasma with chaotic magnetic field fluctuations and this opens up new possibilities for cooling flows.

The evolution of Abell 1795 with thermal conduction but no AGN heating is illustrated in the first column of Figure 2. Each curve in Figure 2 is labeled with three parameters, ($`\epsilon `$,$`f`$,$`t`$ in Gyrs). The influence of conduction on unheated cluster flows is quite sensitive to $`f`$. When $`f=0.5`$ the conductive flux is so high that the solution begins to deviate from the observations after only 1 Gyr and by 4 Gyrs the temperature profile is almost isothermal. This $`T(r)`$ is clearly discordant with observations although $`\dot{M}`$ is nicely suppressed (Table 1). The flow with $`f=0.25`$ is similar to that of the unheated flow in Figure 1. The value of $`f`$ required to maintain the observed profiles in Abell 1795 is $`0.35\pm 0.10`$ (see Table 1), but this optimal value probably varies from cluster to cluster. The fine-tuning of the $`f`$ parameter as well as the difficulty in balancing radiative losses near the centers of cooling flows with the inward conduction of heat from larger radii are well known from previous studies (Stewart, Canizares, Fabian & Nulsen 1984; Bertschinger & Meiksin 1986; Meiksin 1988; Bregman & David 1988; Loewenstein, Zweibel & Begelman 1991; Brighenti & Mathews 2002b).

However, the effect of thermal conduction in heated flows helps to spread the AGN heating into the flow, preserving smooth temperature profiles resembling the observations. Flows with AGN heating and conduction are shown in the second and third columns of Figure 2. For $`f=0.25`$ and $`10^4\begin{array}{c}<\hfill \\ \hfill \end{array}\epsilon \begin{array}{c}<\hfill \\ \hfill \end{array}10^3`$ (second column of Fig. 2) both the density and temperature profiles are maintained for up to 4 Gyrs. These results are similar to the model described by Ruszkowski & Begelman (2002) in which a cluster flow is heated by conduction from the outside and the entire flow is instantaneously heated from the inside when gas flows into the central AGN. We have duplicated the Ruszkowski-Begelman cluster flow with our own program using their type of AGN heating. In general a combination of central AGN heating and thermal conduction can produce cluster flows with many desirable features. While the details of the flows depend on the spatial scale $`r_h`$ of the centrally heated regions, satisfactory gasdynamic flows in general agreement with cluster observations can be found for a wide range of feedback heating scales $`r_h`$. However, unlike Ruszkowski & Begelman, we find that the results are sensitive to the value of the magnetic suppression coefficient $`f`$, as illustrated in the third column of Figure 2. For a given heating efficiency $`\epsilon `$ the flow evolves toward isothermality if $`f=0.5`$ and cools similar to the unheated flow (column 1, Fig. 1) if $`f=0.1`$. The $`\dot{M}`$ for the $`f=0.1`$ solution also exceeds the observed limits for $`t>4`$ Gyrs. Happily, these bad outcomes bracket the value $`f=0.2`$ proposed by Narayan & Medvedev (2001). All the heated flows shown in Figure 2 vary with time within $`rr_h=25`$ kpc. But the negative thermal gradients produced by intermittent heating in this region are short-lived because of the efficiency of the outward heat flux due to thermal conduction. Since none of the times at which the heated flows are plotted in Figure 2 directly follows a heating episode, the temperature gradients shown there are generally positive.

The combination of heating and thermal conduction also has a dramatic influence on the cooling rate $`\dot{M}`$ shown in the bottom row of Figure 2 and in Table 1. Of the four combinations of $`\epsilon `$ and $`f`$ shown in the last two columns of Figure 2, only $`\epsilon =10^4`$ and $`f=0.1`$ has a cooling rate $`\dot{M}`$ that exceeds the XMM upper limits for NGC 1795 after time 4 Gyrs. Since we don’t know how old Abell 1795 is or when it may have been thermally upset by a major merger, the age of the cluster flow must be regarded as an additional parameter in addition to $`\epsilon `$ and $`f`$. The density and temperature profiles plotted in the second column of Figure 2 are in satisfactory agreement with Abell 1795 observations for $`t\begin{array}{c}<\hfill \\ \hfill \end{array}4`$ Gyrs. At later times both flows continue to agree well with the observed profiles until $`t89`$ Gyrs when $`\dot{M}`$ exceeds the XMM upper limits and cooling features develop just beyond the region of AGN heating at $`2550`$ kpc. Within the uncertainties, the two solutions in column 2 successfully match the observations of Abell 1795 for $`t\begin{array}{c}<\hfill \\ \hfill \end{array}8`$ Gyrs.

Temperature profiles in the outer regions of clusters can also be used to constrain the thermal conduction parameter $`f`$, but the observations are controversial. Temperature profiles observed in cluster samples by Markevitch et al. (1998) with ASCA and De Grandi & Molendi (2002) with BeppoSAX characteristically decline by a factor of $`2`$ between 0.1 and 0.5 of the virial radius $`r_{180}`$. Loeb (2002) points out that the outward conductive flux in this region would erase the negative temperature gradients unless $`f\begin{array}{c}<\hfill \\ \hfill \end{array}f_{max}=0.15(t_{cl}/10\mathrm{Gyr})^1(\overline{T}/10\mathrm{keV})^{3/2}`$ where $`t_{cl}`$ is the cluster age and the mean temperature $`\overline{T}`$ is weighted toward larger radii. When applied to Abell 1795 (assuming $`T_{max}/2<\overline{T}<T_{max}=6.85`$ keV) we find $`f_{max}(0.260.75)(t_{cl}/10\mathrm{Gyr})^1`$, which is comparable to the maximum $`f=0.5`$ that we consider in our models. However, the thermal profiles $`T(r)`$ determined by various detectors and X-ray data reduction procedures are inconsistent. For example, the XMM-PN gas temperature profile of Abell 1795 from Tamura et al. (2001) that we plot in Figures 1 and 2 continues to rise slowly at $`r0.4r_{180}`$ whereas the original combined XMM data is nearly isothermal from 0.1 to 0.4 $`r_{180}`$ (Arnaud et al. 2001). Both XMM data sets differ greatly from the mean cluster profile of Markevitch et al. for which $`dT/dr<0`$ throughout this region. Several additional clusters observed more recently with XMM also have nearly isothermal thermal profiles (Arnaud et al. 2002; Pratt & Arnaud 2002). Since the original cluster temperature is likely to decrease with radius (e.g. Loken et al. 2002), these nearly isothermal clusters clusters suggest that $`f`$ may exceed the maximum set by Loeb.

In many respects a combination of AGN heating with efficiency $`\epsilon \begin{array}{c}>\hfill \\ \hfill \end{array}10^4`$ and conductive suppression $`f0.30\pm 0.05`$ gives excellent results for the evolution of rich clusters like Abell 1795. However, Abell 1795 is far away and details of the flow near its central cD galaxy are not available. In central E or cD galaxies in nearby clusters – M87 in Virgo, NGC 4874 in Coma and the cD in Abell 1060 – the gas temperature drops to $`1`$ keV inside the optical galaxy (M87 with XMM: Böhringer et al. 2001; Molendi & Gastaldello 2001; NGC 4874 with Chandra: Vikhlinin et al. 2001; Abell 1060 with Chandra: Yamasaki, Ohashi & Furusho 2002). It is clear from the second column of Figure 2 that the temperature at a few kpc exceeds 1 keV, particularly if $`\epsilon >10^4`$. All our calculations for Abell 1795 included a massive central E galaxy with stellar mass loss etc. so the temperature comparison can be made rather accurately. In our recent paper on M87 and NGC 4874 (Brighenti & Mathews 2002a) we showed that the sharp central temperature drop in the central E galaxy is possible only if $`f\begin{array}{c}<\hfill \\ \hfill \end{array}0.1`$. Perhaps our results in Figure 2 would improve if $`f`$ decreased at small galactic radius.

## 4 Galaxy Group Flows with Conduction and Heating

We now apply the same combination of thermal conduction and heating to flows on galactic scales. For a representative X-ray galaxy we choose NGC 4472, a well-observed massive E1 galaxy and the brightest galaxy in the Virgo cluster, assumed to be at a distance of d = 17 Mpc. NGC 4472 is also a typical group-centered E galaxy since it is surrounded by an extended cloud of hotter gas at $`kT1.3`$ keV, similar to the virial temperature of the galaxy group from which NGC 4472 formed. Although the hot gas distribution in NGC 4472 is asymmetric at large radii, perhaps due to an interaction with the Virgo cluster, its azimuthally averaged density, shown in the top row of Figure 3, is similar to many other group-dominant E galaxies with extended X-ray emission.

Our evolutionary calculation for the galaxy/group differs somewhat from that of for the cluster. We assume that the stellar configuration in NGC 4472 was assembled at some early time $`t_{in}=1`$ Gyr in an NFW dark halo of total mass $`M_h=4\times 10^{13}`$ $`M_{}`$. The stellar density distribution $`\rho _{}(r)`$ has a constant de Vaucouleurs profile (total mass: $`M_t=7.26\times 10^{11}`$ $`M_{}`$; effective radius: $`R_e=1.733^{}=8.57`$ kpc) with a core $`\rho _{,core}(r)=\rho _{,deV}(r_b)(r/r_b)^{0.90}`$ within the break radius $`r_b=2.41^{\prime \prime }=200`$ pc (Gebhardt 1996; Faber et al. 1997). As the (essentially) single burst stellar population evolves, it expels mass at a rate $`\alpha _{}\rho _{}`$ gm cm-3 s-1 where $`\alpha _{}=4.7\times 10^{20}(t/t_n)^{1.3}`$ s-1 and $`t_n=13`$ Gyrs is the current cosmic time. We assume that gas ejected from stars rapidly merges with the hot gas (in $`\begin{array}{c}<\hfill \\ \hfill \end{array}10^5`$ yrs). A modest amount of additional heating is provided by Type Ia supernovae at a rate SNu$`(t)=`$SNu$`(t_n)(t/t_n)^s`$ where SNu$`(t_n)=0.06`$ in SNu-units (supernovae per $`10^{10}L_B`$ per 100 years) gives good agreement with the currently observed iron abundance in the hot gas of NGC 4472. The circumgalactic hot gas component, required to reproduce the X-ray image far beyond $`R_e`$, is also added at time $`t_{in}`$ in hydrostatic equilibrium with the NFW halo; we find that a variety of different scenarios for creating the circumgalactic gas gives similar results at time $`t_n`$.

### 4.1 Unheated Galactic Flow

As before we begin with a 2D cooling flow model for NGC 4472 without AGN or conductive heating, but with radiative cooling. The evolved state of this model at time $`t_n=13`$ Gyrs is shown in the first column of Figure 3. The temperature distribution is in good agreement with the observed profile, the gas cools by radiation losses as it flows inward. However, the computed density clearly exceeds the observed density within about 10 kpc and fails to develop the small central core characteristic of E galaxies observed with Chandra (Loewenstein et al. 2001). This density discrepancy is a generic difficulty with all simple cooling flows of this type. The current cooling rate at the center of NGC 4472 is $`\dot{M}2.4`$ $`M_{}`$ yr-1, but it was larger in the past. Since time $`t_{in}=1`$ Gyr, a total mass $`M_{cold}5.9\times 10^{10}`$ $`M_{}`$ has cooled, but $`M_{cold}`$ can be reduced somewhat if 4472 is assumed to form at a later time. Aside from the unrealistic rise in density, the problems with this model are that $`M_{cold}`$ is too large (for any reasonable $`t_{in}`$) and $`\dot{M}`$ is probably about 5 or 10 times too large, based on XMM RGS observations of NGC 4636, an E galaxy similar to NGC 4472 (Xu et al. 2001).

### 4.2 Galactic Flows with Heating

The second and third columns in Figure 3 show the effect of feedback heating for three efficiencies $`\epsilon `$ using a Gaussian heating profile with $`r_h=2`$ kpc. These are 2D flows, azimuthally averaged in the Figure. Centrally heated flows with $`\epsilon =10^5`$ and $`10^4`$ (in column 2) differ little from the unheated flow except for fluctuations that accompany quasi-periodic heating episodes captured in the Figure at time $`t_n`$ during a short-lived excursion. The global properties of these two models in Table 1 are comparable. The total cooling rate $`\dot{M}`$ is not appreciably reduced until $`\epsilon `$ is increased to $`10^3`$, but at this level of cooling the computed density and (especially the) temperature profiles strongly disagree with observations. The flows shown in the third column of Figure 3 are heated off-center with the Gaussian heating shifted to half the effective radius, $`4.29`$ kpc. The azimuthally averaged profiles of these 2D flows shown at time $`t_n`$ in Figure 3 have almost nothing in common with the observations and must be rejected. Finally, we note that the gas cools in a larger volume ($`r\begin{array}{c}<\hfill \\ \hfill \end{array}30`$ kpc) in the heated flows (Table 1). This distant cooling arises when convective elements collide at times prior to a heating episode when the mean gas density is high; see Brighenti & Mathews (2002b) for a more detailed discussion of this.

### 4.3 Galactic Flows Including Thermal Conduction

The influences of thermal conductivity and AGN heating on 1D flows are shown in the fourth and fifth columns of Figure 3.

In the fourth column of Figure 3 we show the influence of adding thermal conduction but with no additional AGN heating. Each flow is labeled in the Figure with values of ($`\epsilon `$, $`f`$) also listed in Table 1. The 1D ($`0`$, $`f`$) solutions in column four fit the observed inner density profile at time $`t_n`$ but this is a momentary fluctuation. These solutions undergo quasi-cyclic transitions in $`r\begin{array}{c}<\hfill \\ \hfill \end{array}1`$ kpc between high density states (as in column 1 of Fig. 3) and the low density states shown in column 4, Figure 3. The low density phase of the cyclic variation occurs just after an episode of enhanced cooling very near the center of the flow. Then the density slowly increases toward a high state that exceeds observed densities in $`r\begin{array}{c}<\hfill \\ \hfill \end{array}1`$ kpc. These cycles do not appear in similar flows without thermal conduction. The sharp rise in the gas temperature within $`5`$ kpc (column 4, Figure 3) is due to gravitational compression by the large centrally concentrated mass of cooled gas, which was not included in the previously discussed 2D calculations. This feature is of debatable relevance and would be much less pronounced if the galaxy formation time $`t_{in}`$ were later (see below). Although the cooling rate $`\dot{M}`$ shown in column 4 of Figure 3 is clearly reduced by conduction, all of the cooling occurs at the center of the flow where the cold gas mass at $`t_n`$ is $`M_{cold}1.52\times 10^{10}`$ $`M_{}`$ for all $`f`$ considered.

These centrally cooled masses are unacceptably large, but most of this cooling occurred at very early times ($`\alpha _{}t^{1.3}`$) when the galactic evolution is uncertain. For a cosmology with $`H_0=70`$ km s-1 Mpc-1, $`\mathrm{\Omega }=0.3`$ and $`\mathrm{\Lambda }=0.7`$ our initial time $`t_{in}=1`$ Gyr corresponds to a redshift $`z5`$, somewhat larger than may be plausible. In the fourth column of Figure 3 we plot a third solution with $`\epsilon =0`$ and $`f=1`$ but with the calculation starting at time $`t_{in}=6`$ Gyrs corresponding to $`z1.2`$. In this case less gas is expelled from the stars and a much smaller mass of gas cools at the center (Table 1). In this flow the central density profile exceeds observed values, similar to the unheated flow without conduction in the first column.

In the final column of Figure 3 we show two flows with both AGN heating and thermal conduction. The flow with $`\epsilon =10^3`$ and $`f=0.25`$ is an excellent fit to the density and temperature profiles, disregarding the central spike in the temperature created as gas is compressed by the point-like gravitational attraction of centrally cooled gas. But the global cooling rate for this solution, $`\dot{M}1.6`$ $`M_{}`$ yr-1, is too large to be consistent with typical XMM spectral limits. Another solution of this type with $`\epsilon =5\times 10^3`$ and $`f=1`$ is seen to have a good density profile and a very low $`\dot{M}`$, but the temperature is clearly too high for $`r\begin{array}{c}<\hfill \\ \hfill \end{array}10`$ kpc. This temperature rise is an important defect in the flow and one that cannot be removed by fine-tuning $`\epsilon `$ and $`f`$. Because of the lower temperature and conductivity $`\kappa T^{5/2}`$ in galaxy/group flows, the conduction needs to be much closer to the Spitzer value in order to achieve results similar to the cluster flows.

Another quite different problem with these marginally acceptable galaxy/group flows is the iron abundance in the hot gas. Typically the iron abundance is $``$solar near the cores of group-centered E galaxies (e.g. Buote et al. 2002) where the gas is enriched by Type Ia supernovae, each producing $`0.7`$ $`M_{}`$ of iron. At least some of this iron must go into the hot gas phase since the iron abundance in the gas increases toward the center of the stellar system. However in otherwise acceptable flows such as those in the last two columns of Figure 3, where the gas density approximately agrees with observations, our computed iron abundance is much larger than observed. The iron abundance at $`t_n`$ is shown in Figure 4 for the flow with $`\epsilon =5\times 10^3`$ and $`f=1`$. For example, the iron abundance in the computed flow greatly exceeds that observed by Buote (2000). One way to avoid this problem, and also retain some of the good features of the heating plus conduction solutions, may be to carefully regulate (i.e., fine-tune) the fraction of Type Ia iron that ultimately goes into the hot gas. Alternatively, the gas may be cooling in some fashion that cannot be observed with XMM.

## 5 Final Remarks and Conclusions

Cluster flows computed with both AGN heating and thermal conduction reduced by a factor $`0.35\pm 0.10`$ from the Spitzer value can provide acceptable fits to cluster flow observations. Within broad limits, successful flows are insensitive to the (Gaussian) width $`r_h`$ of the feedback-heated region. The observed density and temperature profiles can be maintained with lower cooling rates as required by XMM spectra. Nevertheless, in these solutions the gas temperature may be too high within $`10`$ kpc of the center, i.e. inside the central cD galaxy. If so, it may be necessary to consider flows with a radially variable conductive suppression factor $`f`$.

Flows on galactic scales also benefit from the combined effects of AGN heating and conduction, but the agreement is less satisfactory than for cluster flows. The parameters that provide marginally acceptable, but not ideal, flow solutions on galaxy/group scales – ($`\epsilon `$,$`f`$) = (0,0.5),(0,1) and ($`5\times 10^3`$,1) – require that $`f`$ be much closer to the Spitzer value than for cluster flows. Such large values $`f0.51`$ – at the higher temperatures of clusters – would be inconsistent with the “cold fronts” observed in clusters (e.g. Markevitch et al. 2000) and may promote cluster isothermality as described by Loeb (2002). Values of $`f\begin{array}{c}>\hfill \\ \hfill \end{array}0.5`$ also exceed the predictions of Narayan & Medvedev (2001). Perhaps the field geometry is more favorable to thermal conduction in undisturbed parts of the flows, which we emphasize here, than in the cold fronts where fields transverse to the thermal gradient are expected (Ettori & Fabian 2000). Nevertheless, higher values of the conductivity coefficient $`f`$ are required in galaxies/groups to achieve the same benefits as in hotter cluster gas and it can be doubted if such large conductivities are physically acceptable.

Another difficulty with our galaxy/group heated flows are the high hot gas iron abundances, several times higher than observed. The high iron abundance occurs in the computed flows when the gas density is lowered (by heating) to fit the observed profile. The Type Ia supernova rate that we assume, SNu$`(t_n)=0.06`$ SNu, is already much lower than the combined rate observed in E and S0 galaxies, $`(0.16\pm 0.05)h_{70}^2`$ (Cappellaro et al. 1999). Therefore, to fully accept the beneficial effects of heating plus conduction in galactic flows, it is also necessary to hypothesize some means of reducing the computed iron abundance: some (but not all of) the iron produced in Type Ia may cool before entering the hot gas, some distributed cooling is present to remove the excess iron, etc. Flows on the galaxy/group scale are much more constrained by X-ray observations than cluster scale flows, although the latter have received by far the most attention from observers.

Studies of the evolution of hot gas in elliptical galaxies at UC Santa Cruz are supported by NASA grant NAG 5-8049 and NSF grants AST-9802994 and AST-0098351 for which we are very grateful. FB is supported in part by grant MURST-Cofin 00.


format: text

text:
# Decoherence and Records for the Case of a Scattering Environment

## I Introduction

The notion of decoherence plays a central role in studies of emergent classicality and the foundations of quantum theory generally Zur1 . While it is usually regarded as signifying the loss of quantum coherence for a system of interest, it may also usefully be regarded as an indication of the degree to which information about the system is stored somewhere in the system or in its immediate environment GeH1 ; Hal1 . It is in this way that decoherence is related to “generalized measurements”. An important application of these ideas is in quantum cosmology Har1 . There, in applying quantum theory to the very early universe, there are no actual measuring device to measure what was happening. The process of decoherence, however, guarantees that measurements we make in the present are correlated with alternatives in the past.

These ideas are perhaps most transparent when formulated in terms of the decoherent histories approach to quantum theory GeH1 ; GeH2 ; Gri ; Omn ; Har1 ; Hal2 . Other approaches to decoherence, such as Zurek’s “einselection” approach Zur1 ; Zur2 ; Zur4 , related density matrix approaches JoZ  or quantum state diffusion GiP ; HaZ , may be equally useful for analyzing these issues, but will not be explored here. It is the aim of this paper, continuing in the spirit of Ref.Hal1 , to investigate the connection between decoherence and information storage. To fix ideas, we briefly review the decoherent histories approach (although the general results of this paper are by no means specific that approach).

In the decoherent histories approach GeH1 ; GeH2 ; Gri ; Omn ; Har1 ; Hal2 , probabilities are assigned to histories via the formula,

$$p(\alpha _1,\alpha _2,\mathrm{})=\mathrm{Tr}\left(C_{\underset{¯}{\alpha }}\rho C_{\underset{¯}{\alpha }}^{}\right)$$
(1)
where $`C_{\underset{¯}{\alpha }}`$ denotes a time-ordered string of projectors at times $`t_1\mathrm{}t_n`$,

$$C_{\underset{¯}{\alpha }}=P_{\alpha _n}e^{\frac{i}{\mathrm{}}H(t_nt_{n1})}P_{\alpha _{n1}}\mathrm{}e^{\frac{i}{\mathrm{}}H(t_2t_1)}P_{\alpha _1}$$
(2)
and $`\underset{¯}{\alpha }`$ denotes the string $`\alpha _1,\alpha _2,\mathrm{}\alpha _n`$. We are interested in sets of histories which satisfy the condition of decoherence, which is that decoherence functional

$$D(\underset{¯}{\alpha },\underset{¯}{\alpha }^{})=\mathrm{Tr}\left(C_{\underset{¯}{\alpha }}\rho C_{\underset{¯}{\alpha }^{}}^{}\right)$$
(3)
is zero when $`\underset{¯}{\alpha }\underset{¯}{\alpha }^{}`$. Decoherence implies the weaker condition that $`\mathrm{Re}D(\underset{¯}{\alpha },\underset{¯}{\alpha }^{})=0`$ for $`\underset{¯}{\alpha }\underset{¯}{\alpha }^{}`$, and this is equivalent to the requirement that the above probabilities satisfy the probability sum rules.

The stronger condition of decoherence is the more interesting one since it is related to the existence of records – information storage about the histories somewhere in the system. More precisely, if the initial state is pure, decoherence means that there exist a set of alternatives at the final time $`t_n`$ which are perfectly correlated with the alternatives $`\alpha _1\mathrm{}\alpha _n`$ at times $`t_1\mathrm{}t_n`$ GeH2 ; GeH3 . This follows because, with a pure initial state $`|\mathrm{\Psi }`$, the decoherence condition implies that the states $`C_{\underset{¯}{\alpha }}|\mathrm{\Psi }`$ are an orthogonal set. It is therefore possible to introduce a projection operator $`R_{\underset{¯}{\beta }}`$ (which is generally not unique) such that

$$R_{\underset{¯}{\beta }}C_{\underset{¯}{\alpha }}|\mathrm{\Psi }=\delta _{\underset{¯}{\alpha }\underset{¯}{\beta }}C_{\underset{¯}{\alpha }}|\mathrm{\Psi }$$
(4)
It follows that the extended histories characterized by the chain $`R_{\underset{¯}{\beta }}C_{\underset{¯}{\alpha }}|\mathrm{\Psi }`$ are decoherent, and one can assign a probability to the histories $`\underset{¯}{\alpha }`$ and the records $`\underset{¯}{\beta }`$, given by

$$p(\alpha _1,\alpha _2,\mathrm{}\alpha _n;\beta _1,\beta _2\mathrm{}\beta _n)=\mathrm{Tr}\left(R_{\beta _1\beta _2\mathrm{}\beta _n}C_{\underset{¯}{\alpha }}\rho C_{\underset{¯}{\alpha }}^{}\right)$$
(5)
This probability is then zero unless $`\alpha _k=\beta _k`$ for all $`k`$, in which case it is equal to the original probability $`p(\alpha _1,\mathrm{}\alpha _n)`$. Hence either the $`\alpha `$’s or the $`\beta `$’s can be completely summed out of Eq.(5) without changing the probability, so the probability for the histories can be entirely replaced by the probability for the records at a fixed moment of time at the end of the history:

$$p(\underset{¯}{\alpha })=\mathrm{Tr}\left(R_{\underset{¯}{\alpha }}\rho (t_n)\right)=\mathrm{Tr}\left(C_{\underset{¯}{\alpha }}\rho C_{\underset{¯}{\alpha }}^{}\right)$$
(6)
Conversely, the existence of records $`\beta _1,\mathrm{}\beta _n`$ at some final time perfectly correlated with earlier alternatives $`\alpha _1,\mathrm{}\alpha _n`$ at $`t_1,\mathrm{}t_n`$ implies decoherence of the histories.

These issues are most usefully investigated in the context of particular models, and it then becomes possible to ask some more precise questions: Which dynamical variables in the environment store the information about the decoherent histories? Or what is essentially the same thing, how are the “pointer basis” variables stored in the environment? How is the amount of decoherence related to the amount of information stored?

Ref.Hal1  investigated these questions in the context of the quantum Brownian model (QBM), which consists of a particle of large mass $`M`$ moving in a potential $`V(x)`$ and linearly coupled to a bath of harmonic oscillators. The total system is therefore described by the action,

$`S_T[x(t),q_n(t)]`$ $`=`$ $`{\displaystyle 𝑑t\left[\frac{1}{2}M\dot{x}^2V(x)\right]}`$ (7)
$`+`$ $`{\displaystyle \underset{n}{}}{\displaystyle 𝑑t\left[\frac{1}{2}m_n\dot{q}_n^2\frac{1}{2}m_n\omega _n^2q_n^2c_nq_nx\right]}`$
In the traditional discussion of decoherence in this model, it is shown that for a continuum of oscillators in a thermal state, the influence functional or density matrix become approximately diagonalized in position. This may be seen, for example, through the master equation for the reduced density matrix $`\rho (x,y)`$ of the distinguished system CaL , which in the high temperature limit is

$`{\displaystyle \frac{\rho }{t}}`$ $`=`$ $`{\displaystyle \frac{i\mathrm{}}{2m}}\left({\displaystyle \frac{^2\rho }{x^2}}{\displaystyle \frac{^2\rho }{y^2}}\right){\displaystyle \frac{i}{\mathrm{}}}\left(V(x)V(y)\right)\rho `$ (8)
$``$ $`\gamma (xy)\left({\displaystyle \frac{\rho }{x}}{\displaystyle \frac{\rho }{y}}\right){\displaystyle \frac{2m\gamma kT}{\mathrm{}^2}}(xy)^2\rho `$
In Ref.Hal1 , the issue of how the information about position is stored in the environment was addressed. The system is linear in the oscillators, so the classical and quantum dynamics coincide for the environment. Classically, the equations of motion of the environment of oscillators are,

$$m_n\ddot{q}_n+m_n\omega _n^2q_n=c_nx(t)$$
(9)
The solution to this equation, with fixed $`p_n(0)`$, $`q_n(0)`$ is

$`q_n(\tau )`$ $`=`$ $`q_n(0)\mathrm{cos}\omega _n\tau +{\displaystyle \frac{p_n(0)}{m_n\omega _n}}\mathrm{sin}\omega _n\tau `$
$``$ $`{\displaystyle \frac{c_n}{m_n\omega _n}}{\displaystyle _0^\tau }𝑑tx(t)\mathrm{sin}\omega _n(\tau t)`$ (10)
$`p_n(\tau )`$ $`=`$ $`p_n(0)\mathrm{cos}\omega _n\tau m_n\omega _nq_n(0)\mathrm{sin}\omega _n\tau `$
$``$ $`c_n{\displaystyle _0^\tau }𝑑tx(t)\mathrm{cos}\omega _n(\tau t)`$ (11)
where $`p_n=m\dot{q}_n`$. From this solution, one can see that at the final time $`\tau `$, the positions and momenta of the environment of oscillators depend on the particle’s trajectory $`x(t)`$ via the temporally non-local quantities

$`X_n^s`$ $`=`$ $`{\displaystyle _0^\tau }𝑑tx(t)\mathrm{sin}\omega _n(\tau t)`$ (12)
$`X_n^c`$ $`=`$ $`{\displaystyle _0^\tau }𝑑tx(t)\mathrm{cos}\omega _n(\tau t)`$ (13)
These are essentially the Fourier modes of the particle’s trajectory $`x(t)`$. Hence, each oscillator stores a single Fourier mode of the trajectory, and therefore by using a large number of oscillators, information about many Fourier modes is stored from which the approximate trajectory may be recovered. Furthermore, since it is the Fourier modes that are naturally registered in the environment, rather than positions at each moment of time, decoherence is in fact most clearly seen in terms of the variables (12), (13), rather than position, as shown in Ref.Hal1 . The variables are non-local in time so it can only be seen at the level of an influence functional expressed in path integral language, rather than a master equation. Hence, in this model, it was possible to see exactly how the environment stored the information about the system’s trajectory in configuration space. Furthermore, a detailed quantitative estimate of the information content was also carried out in Ref.Hal1 .

Although a very illustrative model, the quantum Brownian motion model is not the most relevant model for decoherence in physically interesting situations. Far more physically significant is the case in which the environment is a set of light particles which interact with the distinguished particle by a scattering process. The resulting master equation, first derived by Joos and Zeh, is very similar in form to the QBM case, Eq.(8) JoZ ; Dio ; GaF . But the dynamics of the environment, and therefore the means of information storage, are rather different.

The aim of this paper is to investigate the connection between decoherence and records in the case of decoherence by a scattering environment. In some ways it is simpler, since, in the usual assumption of widely separated timescales for system and environment dynamics, each environmental particle scatters briefly off the distinguished particles, and moves freely thereafter, carrying some information about the distinguished particles. This process can therefore be described by a Markovian master equation, and the process of information storage and decoherence may be described in a moment by moment manner (unlike the quantum Brownian motion case, where the environment oscillators store information about the entire history). This means in fact that we do not need to make use of the full machinery of the decoherence functional – it is sufficient in fact to look at the evolution of the reduced density operator.

In the quantum Brownian motion model case, the system variables the environment measures were actually identified quite simply, from the classical equations of motion. In the scattering case, the system variables measured most directly by the environment are also determined quite easily, by examining simple scattering processes. In particular, suppose we consider the scattering of some light particles off a dilute gas of a set of more massive particles with coordinates $`𝐪_j`$. Then it follows quite straightforwardly from simple scattering theory (and we will in fact demonstrate this) that the scattering amplitude is proportional to the Fourier transform of the number density of the massive particles,

$$N_𝐤=\underset{j}{}e^{i𝐤𝐪_j}$$
(14)
This means that, loosely speaking, for a known interaction potential, measurements of the initial and final momenta of the scattering environment determine the number density of the distinguished system.

Of course, the number density is closely related to position, which is normally held to be the preferred basis in these calculations. But, following the lead of the oscillator model, we expect decoherence to look simplest in terms of the dynamical variables which are most simply and directly stored in the environment. Our aim is therefore to give a derivation of the master equation which emphasize the central role played by the number density. We have found that the derivation is in fact most transparent in terms of non-relativistic many body quantum field theory, where the number density appears very naturally. We will give an alternative and more general derivation of the master equation, using many body theory, which brings out the role of local number density more clearly, hence showing the connection with records.

It is pertinent at this stage to mention the Lindblad form of the master equation Lin , which is the most general possible form a master equation can take under the assumption that the evolution is Markovian (a condition well-satisfied in a wide variety of interesting models). The Lindblad master equation is

$$\frac{d\rho }{dt}=i[H,\rho ]\frac{1}{2}\underset{j=1}{\overset{n}{}}\left(\{L_j^{}L_j,\rho \}2L_j\rho L_j^{}\right)$$
(15)
Here, $`H`$ is the Hamiltonian of the distinguished subsystem (sometimes modified by terms depending on the $`L_j`$) and the $`n`$ operators $`L_j`$ model the effects of the environment. The master equation of quantum Brownian motion, for example, is of this form with

$$L=\left(\frac{4m\gamma kT}{\mathrm{}^2}\right)^{\frac{1}{2}}x+i\left(\frac{\gamma }{2mkT}\right)^{\frac{1}{2}}p$$
(16)
as described in Refs.Dio2 ; HaZ . (Actually, the master equation (8) is not strictly of the Lindblad form, and as a consequence can suffer from a violation of positivity Amb . However, the difference between Eq.(8) and the Lindblad form with $`L`$ given by Eq.(16) is of order $`1/T`$ which does not matter for high temperatures).

The Lindblad operators $`L_j`$ determine the basis in which the density operator tends to become approximately diagonal, or what is essentially the same, the sets of variables describing an approximately decoherent set of histories. This may be seen from the formal solution to the Lindblad equation DGHP . Consider the case of a single Lindblad operator $`L=L_R+iL_I`$, where $`L_R`$, $`L_I`$ are hermitian. Divide the finite time interval $`[0,t]`$ into $`K`$ subintervals, so that $`t=K\delta t`$, and let $`\delta t0`$, $`K\mathrm{}`$, holding $`t`$ constant. Then, the formal solution to the Lindblad equation is obtained by taking the limit $`\delta t0`$, $`K\mathrm{}`$ (with $`t`$ fixed) of the expression,

$`\rho (t)`$ $`=`$ $`\left({\displaystyle \frac{\delta t}{\pi }}\right)^K{\displaystyle d^2ł_1\mathrm{}d^2ł_K}`$ (17)
$`\times `$ $`{\displaystyle \underset{m=1}{\overset{K}{}}}\mathrm{exp}\left({\displaystyle \frac{\delta t}{2}}(\mathrm{}_m^{}L\mathrm{}_mL^{})\right)\mathrm{exp}\left({\displaystyle \frac{\delta t}{2}}|L\mathrm{}_m|^2\right)\mathrm{exp}\left({\displaystyle \frac{i}{\mathrm{}}}H^{}\delta t\right)\rho (0)`$
$`\times `$ $`{\displaystyle \underset{m=1}{\overset{K}{}}}\mathrm{exp}\left({\displaystyle \frac{i}{\mathrm{}}}H^{}\delta t\right)\mathrm{exp}\left({\displaystyle \frac{\delta t}{2}}|L\mathrm{}_m|^2\right)\mathrm{exp}\left({\displaystyle \frac{\delta t}{2}}(\mathrm{}_m^{}L\mathrm{}_mL^{})\right)`$
Here, $`H^{}=H+\frac{i\mathrm{}}{4}[L,L^{}]`$, and the $`\mathrm{}_m`$ are complex numbers at the discrete moments of time labelled by $`m`$. We use the notation

$$|L\mathrm{}_m|^2\left(L_R\mathrm{Re}\mathrm{}_m\right)^2+\left(L_I\mathrm{Im}\mathrm{}_m\right)^2$$
(18)
The ordering of the operators at each moment of time is irrelevant in the limit $`\delta t0`$ (although the operators at different times are time-ordered, according to increasing $`m`$). That this is the solution is readily verified by explicit computation. The solution has the form of a “measurement process” of the $`L`$’s, continuous in time, with “feedback” via the terms $`(\mathrm{}_m^{}L\mathrm{}_mL^{})`$ Caves . That is, one can see that the effect of the environment is to produce a tendency towards diagonality in $`L`$.

We shall show that a many-body theory derivation of the master equation in the case of a scattering environment leads to a master equation of the Lindblad form (under the assumption that the environment dynamics are much faster than the system dynamics), and that the Lindblad operators are essentially the local number density. The previous forms of the master equation are recovered in the one-particle sector for the system of massive particles.

This work grew out of a more ambitious programme, in the context of the decoherent histories approach to quantum theory, which aims to give a very general account of emergent classicality. In particular, it is asserted that at sufficiently coarse-grained scales, the local densities (number, momentum, energy) define a a set of habitually decohering variables, even in the absence of an environment GeH1 ; Hal3 . This is because they are locally conserved, and therefore slowly varying when coarse-grained over sufficiently large volumes, and thus are expected to be approximately decoherent (because exactly conserved quantities are exactly decoherent in the histories approach). This is therefore a different mechanism for decoherence than the usual mechanism of decoherence through an environment. Hence, in order to close the gap between the familiar system–environment models and the less familiar hydrodynamic models without an obvious environment, it is useful to rewrite the system–environment models in terms of local densities as we do here.

In Section 2, we briefly review many body field theory, and carry out the simple scattering calculation leading to the result the the scattering particles effectively measure the local number density, Eq.(14).

In Section 3, we use many body field theory to derive the master equation for the system, using a slow motion limit for the gas of massive particles. As anticipated it has the Lindblad form in with the Lindblad operators proportional to the Fourier-transformed number density $`N_𝐤`$.

In Section 4 we show that our master equation reduces to an earlier result of Gallis and Fleming GaF  in the one-particle sector for the gas of massive particles. (This is essentially the same as the master equation of Joos and Zeh JoZ  but the comparison with Gallis and Fleming is more direct).

The master equation of Sections 3 and 4 does not have any dissipation and is analogous to the Lindblad equation of quantum Brownian motion with $`L`$ proportional to $`x`$. In Section 5, we go beyond the slow-motion limit to derive a master equation with dissipative terms.

We summarize and conclude in Section 6.

## II Many Body Field Theory

The dynamics of a many body system is very conveniently handled using many body quantum field theory. We now set up the formalism of many body field theory Zub ; FeWa  which we will use to derive the master equation. We consider a set of non-relativistic system particles described by a field $`\psi (𝐱)`$ interacting through a potential $`\varphi (𝐱)`$ with an environment described by a field $`\chi (𝐱)`$. The total system is described by the Hamiltonian

$`H`$ $`=`$ $`{\displaystyle d^3x\left(\frac{1}{2M}\psi ^{}(𝐱)\psi (𝐱)+\frac{1}{2m}\chi ^{}(𝐱)\chi (𝐱)\right)}`$ (19)
$`+`$ $`{\displaystyle \frac{1}{2}}{\displaystyle d^3xd^3x^{}\psi ^{}(𝐱)\psi (𝐱^{})\varphi (𝐱𝐱^{})\chi ^{}(𝐱^{})\chi (𝐱)}`$
(For simplicity we set $`\mathrm{}=1`$ hereafter). In this language, the number densities $`N(𝐱)`$ and $`n(𝐱)`$ of the system and environment fields are

$`N(𝐱)`$ $`=`$ $`\psi ^{}(𝐱)\psi (𝐱)`$ (20)
$`n(𝐱)`$ $`=`$ $`\chi ^{}(𝐱)\chi (𝐱)`$ (21)
(This is the field-theoretic version of Eq.(14)).

The above relations are also more conveniently written in terms of $`a_𝐤`$ and $`b_𝐤`$, the annihilation operators for the system and environment, respectively, and the Hamiltonian then is

$`H`$ $`=`$ $`{\displaystyle \underset{𝐪}{}}\left(E_𝐪a_𝐪^{}a_𝐪+\omega _𝐪b_𝐪^{}b_𝐪\right)`$
$`+`$ $`{\displaystyle \frac{1}{2V}}{\displaystyle \underset{𝐤_1^{}+𝐤_2^{}=𝐤_1+𝐤_2}{}}\nu (𝐤_2^{}𝐤_2)a_{𝐤_1}^{}b_{𝐤_2}^{}a_{𝐤_1^{}}b_{𝐤_2^{}}`$
where $`E_𝐪=𝐪^2/2M`$, $`\omega _𝐪=𝐪^2/2m`$, $`V`$ is the spatial volume of the system (which we assume is in a box) and

$$\nu (𝐤)=d^3xe^{i𝐤𝐱}\varphi (𝐱)$$
(23)
The Fourier transformed number densities are

$`N_𝐤`$ $`=`$ $`{\displaystyle \underset{𝐪}{}}a_𝐪^{}a_{𝐪+𝐤}`$ (24)
$`n_𝐤`$ $`=`$ $`{\displaystyle \underset{𝐪}{}}b_𝐪^{}b_{𝐪+𝐤}`$ (25)
and one may see that the Hamiltonian has the more concise form

$`H`$ $`=`$ $`{\displaystyle \underset{𝐪}{}}\left(E_𝐪a_𝐪^{}a_𝐪+\omega _𝐪b_𝐪^{}b_𝐪\right)+{\displaystyle \frac{1}{2V}}{\displaystyle \underset{𝐤}{}}\nu (𝐤)N_𝐤n_𝐤`$ (26)
$`=`$ $`H_0+H_{int}`$ (27)
From these relations we see that the environment couples to the number density of the system. It is this feature of many body field theory that makes it the appropriate medium for the derivation of the master equation emphasizing the role of number density.

The $`S`$-matrix is

$$S=T\mathrm{exp}\left(i_{\mathrm{}}^{\mathrm{}}𝑑tH_{int}(t)\right)$$
(28)
where

$$H_{int}(t)=\frac{1}{2V}\underset{𝐤}{}\nu (𝐤)N_𝐤(t)n_𝐤(t)$$
(29)
and here

$`N_𝐤(t)`$ $`=`$ $`{\displaystyle \underset{𝐪}{}}a_𝐪^{}a_{𝐪+𝐤}e^{i(E_𝐪E_{𝐪+𝐤})t}`$ (30)
$`n_𝐤(t)`$ $`=`$ $`{\displaystyle \underset{𝐪}{}}b_𝐪^{}b_{𝐪+𝐤}e^{i(\omega _𝐪\omega _{𝐪+𝐤})t}`$ (31)
We may now use this formalism to look at a simple scattering situation to determine how the environment stores information about the system. In the quantum Brownian motion case, the nature of information storage was determined in essence by solving the classical equations of motion. A similar strategy works here too. Let us suppose the distinguished system is classical, and consider what happens when the environment scatters off it. Suppose the environment starts in an initial momentum state $`|𝐤_0`$ and scatters into a final state $`|𝐤_f`$. The scattering amplitude for this process, to first order, is

$`𝐤_f|S|𝐤_0`$ $`=`$ $`{\displaystyle \frac{i}{2V}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑t{\displaystyle \underset{𝐤}{}}\nu (𝐤)N_𝐤(t)𝐤_f|n_𝐤(t)|𝐤_0`$ (32)
$`=`$ $`{\displaystyle \frac{i}{2V}}\nu (𝐤){\displaystyle 𝑑t}N_𝐤(t)e^{i(\omega _{𝐤_f}\omega _{𝐤_0})t}`$
where $`𝐤=𝐤_f𝐤_0`$. This simple result shows that a single scattering event by the environment stores information about the Fourier transform (in space and time) of the number density. It is in this sense that the number density has a preferred status – this is the variable that is measured most directly by the environment and is the exact analogy of the relations Eqs.(12), (13) in the quantum Brownian case.

The measured variables above are of course non-local in time, involving a temporal Fourier transform of the number density, so cannot in fact be compatible with a Markovian master equation. Under a reasonable slow motion assumption, the system timescale is much slower than the environment timescale, and we may ignore the time-dependence in $`N_𝐤(t)`$, yielding

$$𝐤_f|S|𝐤_0=\frac{i}{2V}\nu (𝐤)N_𝐤\delta (\omega _{𝐤_f}\omega _{𝐤_0})$$
(33)
This corresponds more directly to a Markovian master equation, as we shall see.

It remains to briefly sketch the connection between these results and the discussion in Section 1 of records in the decoherent histories approach. We imagine that the environment consists of a very large number of particles which scatter off the system particles. Each scattering event consists of an incoming environment particle with momentum $`𝐤_0`$ scattered into a final state of momentum $`𝐤_f`$, as outlined above. After the scattering event, which is essentially instantaneous (on the timescale of system dynamics), we may imagine that the scattered particle propagates freely and may be measured at any time in the future. Therefore, the records in the decoherent histories approach consist of projections at the end of this histories onto the momenta of all the scattered environment particles, from which the number densities $`N_𝐤`$ of the system at a series of times may be retrodicted.

## III Derivation of the Master Equation

Following the method first used by Joos and Zeh JoZ , we may derive the master equation for the reduced density operator $`\rho `$ of the system by considering the scattering of the environment off the system, to second order in interactions. We assume that the system and environment are initally uncorrelated, so the total density operator is

$$\rho _T=\rho _0\rho _{}$$
(34)
We also assume that each scattering event takes place on a timescale which is extremely short compared to the timescale of system dynamics. This means that in an interval of time $`\mathrm{\Delta }t`$ which is long for the environment but short for the system, we may write,

$$\rho _T(t+\mathrm{\Delta }t)=S\rho _T(t)S^{}$$
(35)
(where we are using the interaction picture). Expanding (35) to second order, the $`S`$-matrix may be written,

$$S=1+iU_1U_2$$
(36)
where

$$U_1=_{\mathrm{}}^{\mathrm{}}𝑑tH_{int}(t)$$
(37)
and

$$U_2=\frac{1}{2}𝑑t_1𝑑t_2\mathrm{T}\left(H_{int}(t_1)H_{int}(t_2)\right)$$
(38)
The requirement of unitarity, $`S^1=S^{}`$, implies that $`U_1=U_1^{}`$ and

$$U_2+U_2^{}=U_1^2$$
(39)
We will therefore write

$$U_2=\frac{1}{2}U_1^2+iB$$
(40)
where $`B=B^{}`$, so we now have

$$S=1+i(U_1B)\frac{1}{2}U_1^2$$
(41)
Inserting this in (35), we obtain

$$\frac{d\rho _T}{dt}\mathrm{\Delta }t=i[U_1B,\rho _T]+U_1\rho _TU_1\frac{1}{2}U_1^2\rho _T\frac{1}{2}\rho _TU_1^2$$
(42)
We now trace Eq.(42) over the environment to obtain the master equation for the system density operator $`\rho `$. As is usual in this sort of model, we assume that the environment is so large that its state is essentially unaffected by the interaction with the system. Since the total density operator starts out in the factored state (34), this then means that, to a good approximation, $`\rho _T`$ persists in the approximately factored form $`\rho \rho _{}`$, and we may insert this in the right-hand side of Eq.(42) Gar . We thus obtain the preliminary form for the master equation

$$\frac{d\rho }{dt}\mathrm{\Delta }t=i[\mathrm{Tr}_{}(U_1\rho _{})\mathrm{Tr}_{}(B\rho _{}),\rho ]+\mathrm{Tr}_{}\left(U_1\rho _TU_1\frac{1}{2}U_1^2\rho _T\frac{1}{2}\rho _TU_1^2\right)$$
(43)
We now work out these terms in more detail. We first consider the simple but useful slow motion approximation, in which we ignore the time-dependence of $`N_𝐤(t)`$. This implies that

$$U_1\frac{1}{2V}\underset{𝐤}{}\nu (𝐤)N_𝐤\underset{𝐪}{}b_𝐪^{}b_{𝐪𝐤}2\pi \delta (\omega _𝐪\omega _{𝐪𝐤})$$
(44)
The important terms for decoherence are the final three terms on the right-hand side of (44). When traced, these give,

$$\mathrm{Tr}_{}\left(U_1\rho _TU_1\frac{1}{2}U_1^2\rho _T\frac{1}{2}\rho _TU_1^2\right)=\underset{\mathrm{𝐤𝐤}^{}}{}c(𝐤,𝐤^{})\left(N_𝐤^{}\rho N_𝐤\frac{1}{2}N_𝐤N_𝐤^{}\rho \frac{1}{2}\rho N_𝐤N_𝐤^{}\right)$$
(45)
where

$$c(𝐤,𝐤^{})=\nu (𝐤)\nu (𝐤^{})\underset{\mathrm{𝐪𝐪}^{}}{}\delta (\omega _𝐪\omega _{𝐪𝐤})\delta (\omega _𝐪^{}\omega _{𝐪^{}𝐤^{}})b_𝐪^{}b_{𝐪𝐤}b_𝐪^{}^{}b_{𝐪^{}𝐤^{}}_{}$$
(46)
We will take the environment to be a thermal state, which is diagonal in the momentum states. It follows that

$$b_𝐪^{}b_{𝐪𝐤}b_𝐪^{}^{}b_{𝐪^{}𝐤^{}}_{}\delta _{𝐪,𝐪^{}𝐤^{}}\delta _{𝐪^{},𝐪𝐤}$$
(47)
This implies $`𝐤=𝐤^{}`$, and also that the two delta-functions are the same in Eq.(46). We then interpret the square of the delta-function in the usual way,

$`\left[\delta (\omega _𝐪\omega _{𝐪𝐤})\right]^2`$ $`=`$ $`\delta (0)\delta (\omega _𝐪\omega _{𝐪𝐤})`$ (48)
$`=`$ $`{\displaystyle \frac{\mathrm{\Delta }t}{2\pi }}\delta (\omega _𝐪\omega _{𝐪𝐤})`$
We now have

$$c(𝐤,𝐤^{})=\delta _{𝐤,𝐤^{}}c(𝐤)\frac{\mathrm{\Delta }t}{2\pi }$$
(49)
where

$`c(𝐤)`$ $`=`$ $`{\displaystyle \frac{1}{2\pi }}\left|\nu (𝐤)\right|^2{\displaystyle \underset{𝐪}{}}\delta (\omega _𝐪\omega _{𝐪𝐤})b_𝐪^{}b_{𝐪𝐤}b_{𝐪𝐤}^{}b_𝐪_{}`$ (50)
$`=`$ $`{\displaystyle \frac{1}{2\pi }}\left|\nu (𝐤)\right|^2{\displaystyle \underset{𝐪}{}}\delta (\omega _𝐪\omega _{𝐪𝐤})b_𝐪^{}b_𝐪_{}\left(b_{𝐪𝐤}^{}b_{𝐪𝐤}_{}+1\right)`$
The terms involving environment averages have the usual thermal form (for a bosonic environment),

$$b_𝐪^{}b_𝐪_{}=\frac{1}{e^{\beta (\omega _𝐪\mu )}1}$$
(51)
where $`\beta =1/kT`$ with $`T`$ temperature, and $`\mu `$ is the chemical potential.

The form Eq.(49) means that the important terms in the master equation are of the Lindblad form,

$$\mathrm{Tr}_{}\left(U_1\rho _TU_1\frac{1}{2}U_1^2\rho _T\frac{1}{2}\rho _TU_1^2\right)=\mathrm{\Delta }t\underset{𝐤}{}c(𝐤)\left(N_𝐤\rho N_𝐤^{}\frac{1}{2}N_𝐤^{}N_𝐤\rho \frac{1}{2}\rho N_𝐤^{}N_𝐤\right)$$
(52)
where we have used the fact that $`N_𝐤^{}=N_𝐤`$. The remaining two terms in Eq.(43) clearly just modify the unitary dynamics of the system. First we have

$$\mathrm{Tr}_{}\left(U_1\rho _{}\right)=\frac{1}{2V}\underset{𝐤}{}\nu (𝐤)N_𝐤\underset{𝐪}{}b_𝐪^{}b_{𝐪𝐤}_{}2\pi \delta (\omega _𝐪\omega _{𝐪𝐤})$$
(53)
Clearly from the term $`b_𝐪^{}b_{𝐪𝐤}_{}`$ this expression will be zero unless $`𝐤=0`$, and therefore it is proportional to $`N`$, the total particle number operator (although the overall coefficient will need to be regularized). This therefore contributes a term to the master equation of the form $`[N,\rho ]`$. We assume that there is a fixed number of system particles so it is reasonable to take this term to be zero.

The other remaining term in Eq.(43) involves the time ordering terms in $`U_2`$ and is a bit more complicated to evaluate. Fortunately, the detailed form of this expression is not needed here, and it can in fact be easily shown that this term has the form

$$\mathrm{Tr}_{}(B\rho _{})=\mathrm{\Delta }t\underset{𝐤}{}d(𝐤)N_𝐤N_𝐤^{}$$
(54)
for some coefficient $`d(𝐤)`$ which we will not need. Inserting all these results in Eq.(43), the factors of $`\mathrm{\Delta }t`$ all drop out, and we obtain, in the Schrödinger picture,

$$\frac{d\rho }{dt}=i[H_0\underset{𝐤}{}d(𝐤)N_𝐤N_𝐤^{},\rho ]+\underset{𝐤}{}c(𝐤)\left(N_𝐤\rho N_𝐤^{}\frac{1}{2}N_𝐤^{}N_𝐤\rho \frac{1}{2}\rho N_𝐤^{}N_𝐤\right)$$
(55)
As desired, this is the Lindblad form with the Lindblad operators given by

$$L_𝐤=c^{\frac{1}{2}}(𝐤)N_𝐤$$
(56)
We have therefore produced a derivation of the master equation for a scattering environment which shows very clearly the connection between the preferred basis (diagonalization in the Lindblad operators) and the information storage about the system, as indicated by the simple scattering calculation, Eq.(33).

It is interesting to note that the decoherence effect is second order in interactions, but we were able to anticipate it from the simple first order calculation, Eq.(33). The reason for this is the relationship Eq.(40), which shows that the important part of the second order terms is the square of the first order terms, and this is a consequence of unitarity.

## IV Comparison with Previous Works

It is useful to check that the master equation we have derived reproduces known results when we restrict to the one-particle sector for the system. We will compare the results of this to the derivation of Gallis and Fleming GaF  (which is essentially the same as Joos and Zeh JoZ  and Diosi Dio ).

In the one-particle sector we may work with a density matrix $`\rho (𝐤,𝐤^{})=𝐤|\rho |𝐤^{}`$, or equivalently $`\rho (𝐱,𝐲)`$ in the position representation. We use the relations

$`[N_𝐪,a_𝐤^{}]`$ $`=`$ $`a_{𝐤𝐪}^{}`$ (57)
$`[N_𝐪,a_𝐤]`$ $`=`$ $`a_{𝐤+𝐪}`$ (58)
These relations imply that

$`N_𝐪\rho (𝐤,𝐤^{})N_𝐪`$ $`=`$ $`\rho (𝐤𝐪,𝐤^{}𝐪)`$ (59)
$`N_𝐪N_𝐪\rho (𝐤,𝐤^{})`$ $`=`$ $`\rho (𝐤,𝐤^{})`$ (60)
$`\rho (𝐤,𝐤^{})N_𝐪N_𝐪`$ $`=`$ $`\rho (𝐤,𝐤^{})`$ (61)
In the position representation, this means

$$N_𝐤\rho (𝐱,𝐲)N_𝐤=e^{i𝐤.(𝐱𝐲)}\rho (𝐱,𝐲)$$
(62)
The master equation for the one-particle density operator $`\rho (𝐱,𝐲)`$ is then

$$\frac{\rho (𝐱,𝐲)}{t}=i𝐱|[H_0,\rho ]|𝐲F(𝐱𝐲)\rho (𝐱,𝐲)$$
(63)
where

$$F(𝐱𝐲)=\frac{\pi }{2V^2}d^3qd^3k\left|\nu (𝐤)\right|^2n_𝐪(n_{𝐪𝐤}+1)\delta (\omega _𝐪\omega _{𝐪𝐤})(1e^{i𝐤(𝐱𝐲)})$$
(64)
Note that the term involving the coefficient $`d(𝐤)`$ in Eq.(55) drops out because $`[N_𝐤N_𝐤^{},\rho ]=0`$ in the one-particle sector.

To compare this with the Gallis-Fleming result GaF , we first introduce the quantity

$$f(𝐤,𝐤^{})=\frac{m}{2\pi }\nu (𝐤𝐤^{})$$
(65)
(which appears in the usual Born approximation to first order scattering). Then, letting $`𝐤𝐤+𝐪`$ in (64), we get

$$F(𝐫)=\frac{4\pi ^2}{m^2}d^3qd^3k\left|f(𝐪,𝐤)\right|^2n_𝐪(n_𝐤+1)\delta (\omega _𝐪\omega _𝐤)(1e^{i(𝐪𝐤)𝐫})$$
(66)
The delta-function implies that $`𝐪^2=𝐤^2`$, and we find that

$$F(𝐫)=\frac{4\pi ^2}{m^2}𝑑qq^3n_q(n_q+1)𝑑\mathrm{\Omega }𝑑\mathrm{\Omega }^{}|f(𝐪,𝐤)|^2(1e^{i(𝐪𝐤)𝐫})$$
(67)
This in fact agrees with Gallis and Fleming if we identify $`q/m`$ as their $`v(q)`$, the speed of the incoming particles, and $`4\pi ^3q^2n_q(n_q+1)`$ as the density of particles with speed $`q`$. In the one particle sector there is therefore agreement with earlier work. (At least up to an overall numerical factor which we could not rectify. However, we have also spotted some small and probably insignificant numerical errors in Ref.GaF ).

Mention should also be made of the master equations derived by Unruh and Zurek, which used a field as an environment for a particle UnZu , and Anastopoulos and Zoupas AnZo , which used a photon field and as an environment for a spinor field. Also of relevance is the general account of the derivation of master equations given by Omnès Omn1 . These works are rather different to the present paper.

## V Beyond the Slow Motion Approximation

The derivation above assumed, in essence, that the system dynamics are infinitely slow. Not surprisingly, the resulting master equation does not involve dissipation, since, in the approximation used, the system is essentially at rest for the timescale of a single scattering event. It is analogous to the master equation of quantum Brownian motion with the Lindblad operator $`L`$ proportional to $`x`$, Eq.(16). To get a more realistic equation with dissipation we therefore need to go beyond the slow motion approximation.

Because the local number density is a locally conserved quantity, it obeys a continuity equation of the form,

$$\dot{N}_𝐤=i𝐤𝐏_𝐤$$
(68)
where $`𝐏_𝐤`$ is the local momentum density,

$$𝐏_𝐤=\underset{𝐪}{}\left(𝐪+\frac{1}{2}𝐤\right)a_𝐪^{}a_{𝐪+𝐤}$$
(69)
It is reasonable to expect that the master equation will involve this operator when we go beyond the infinitely slow limit. We now briefly repeat the derivation of the master equation, this time allowing a slow time-dependence in $`N_𝐤(t)`$.

We have

$$U_1=\frac{1}{2V}\underset{𝐤}{}\nu (𝐤)\underset{𝐪}{}b_𝐪^{}b_{𝐪𝐤}𝑑tN_𝐤(t)e^{i(\omega _𝐪\omega _{𝐪𝐤})t}$$
(70)
To take into account the time-dependence of $`N_𝐤(t)`$, we write,

$$N_𝐤(t)=N_𝐤+t\dot{N}_𝐤+\mathrm{}$$
(71)
where $`\dot{N}_𝐤`$ is given in terms of the momentum density, Eq.(69). Inserting this in Eq.(70), the factor of $`t`$ may be rewritten in terms of a delta-function derivative, yielding,

$$U_1=\frac{1}{2V}\underset{𝐤}{}\nu (𝐤)\underset{𝐪}{}b_𝐪^{}b_{𝐪𝐤}\left(N_𝐤\delta (\omega _𝐪\omega _{𝐪𝐤})i\dot{N}_𝐤\delta ^{}(\omega _𝐪\omega _{𝐪𝐤})+\mathrm{}\right)$$
(72)
We now use this expression for $`U_1`$ in the derivation of the master equation. So for example, we get, in place of Eq.(45),

$$\mathrm{Tr}_{}\left(U_1^2\rho _{}\right)=\mathrm{\Delta }t\underset{𝐤}{}c(𝐤)N_𝐤^{}N_𝐤\rho +i\underset{𝐤}{}b(𝐤)\left(N_𝐤\dot{N}_𝐤^{}N_𝐤^{}\dot{N}_𝐤\right)\rho +\mathrm{}$$
(73)
where $`c(𝐤)`$ is given by Eq.(50) and

$$b(𝐤)=\left|\nu (𝐤)\right|^2\underset{𝐪}{}\delta (\omega _𝐪\omega _{𝐪𝐤})\delta ^{}(\omega _𝐪\omega _{𝐪𝐤})b_𝐪^{}b_𝐪_{}\left(b_{𝐪𝐤}^{}b_{𝐪𝐤}_{}+1\right)$$
(74)
This coefficient may in fact be shown to be simply related to $`c(𝐤)`$. The delta-function derivative may be dealt with by noting the formal relation

$$\delta (x)\delta ^{}(x)=\frac{1}{2}\frac{}{x}\left[\delta (x)\right]^2$$
(75)
Now note that

$$\omega _𝐪\omega _{𝐪𝐤}=\frac{1}{2m}\left(2𝐤𝐪𝐤^2\right)$$
(76)
It follows that the delta function derivatives may be expressed in terms of derivatives with respect to $`q_i`$ as

$$\delta ^{}(\omega _𝐪\omega _{𝐪𝐤})=\frac{2m}{𝐤^2}k_i\frac{}{q_i}\delta (\omega _𝐪\omega _{𝐪𝐤})$$
(77)
Inserting these relations in Eq.(74) and integrating by parts yields,

$$b(𝐤)=\frac{1}{2}\left|\nu (𝐤)\right|^2\underset{𝐪}{}\left[\delta (\omega _𝐪\omega _{𝐪𝐤})\right]^2\frac{2m}{𝐤^2}k_i\frac{}{q_i}\left(b_𝐪^{}b_𝐪_{}\left(b_{𝐪𝐤}^{}b_{𝐪𝐤}_{}+1\right)\right)$$
(78)
where we will interpret the square of the delta-function as in Eq.(48).

Now for simplicity work in the high temperature limit, so

$$b_{𝐪𝐤}^{}b_{𝐪𝐤}_{}1$$
(79)
and

$$b_𝐪^{}b_𝐪_{}e^{\mu \beta }e^{\beta \omega _𝐪}$$
(80)
It follows that

$`k_i{\displaystyle \frac{}{q_i}}b_𝐪^{}b_𝐪_{}`$ $`=`$ $`\beta k_i{\displaystyle \frac{\omega _𝐪}{q_i}}b_𝐪^{}b_𝐪_{}`$ (81)
$`=`$ $`\beta {\displaystyle \frac{𝐤𝐪}{m}}b_𝐪^{}b_𝐪_{}`$
and we arrive at the very simple result,

$$b(𝐤)=\frac{\beta }{2}\mathrm{\Delta }tc(𝐤)$$
(82)
It is not difficult to see that we then arrive at a master equation which is once again of the Lindblad form, but this time with Lindblad operators of the form

$`L_𝐤`$ $`=`$ $`c^{\frac{1}{2}}(𝐤)\left(N_𝐤i{\displaystyle \frac{\beta }{2}}\dot{N}_𝐤\right)`$ (83)
$`=`$ $`c^{\frac{1}{2}}(𝐤)\left(N_𝐤{\displaystyle \frac{\beta }{2}}𝐤𝐏_𝐤\right)`$
(up to terms of order $`\beta ^2`$, which can be dropped in the approximation we are using). This is clearly closely analogous to the QBM result, Eq.(16). (A closely analogous formula appears in Diosi’s paper Dio ).

## VI Summary and Discussion

We have given a derivation of the master equation describing a many-body system interacting with a reasonably general class of environments. The form of the master equation emphasizes the central role of the local number density, which is the system variable measured most directly by the environment in a scattering situation.

We did not in fact give a specific form for the interaction between the system and environment since it was not necessary to illustrate the general points we are making. Some specific forms for this interaction are discussed elsewhere GaF ; JoZ .

The derivation reduces to familiar results of Gallis and Fleming GaF , Diosi Dio , and Joos and Zeh JoZ , when we restrict to the one-particle sector of the many-body field theory. The many-body derivation confers some advantages of the usual derivations (which consider scattering theory in quantum mechanics) in that they avoid essentially classical assumptions about fluxes of scattering particles. Our derivation also has the possibility of being extended to a low temperature regime (and to Bose-Einstein condensation, for example) and to fermionic environments, although we do not discuss this here.

## VII Acknowledgements

We are grateful to Lajos Diosi for useful conversations. P.D. was supported by PPARC.


text:
# R&D ON THE GEM READOUT OF THE TESLA TPC

## 1 Introduction

A Time Projection Chamber (TPC, ) is proposed to be the main tracking device for the detector at the TESLA collider .
Our purpose is to study the TPC readout using Gas Electron Multipliers (GEMs, ) instead of conventionally used wire chambers to produce the gas amplification. GEMs offer several advantages: $`\stackrel{}{E}\times \stackrel{}{B}`$ effects degrading the spatial resolution are strongly reduced, and the amount of required material in the endcap is low; ion feedback into the drift volume leading to field distortion is naturally suppressed to a level of a few percent. Due to the absence of ion tails in the pulse shape the use of GEMs improves the intrinsic two track resolution in the drift direction .
In the following our present TPC test setups and the measurements are presented.

## 2 Setups and Measurements

### 2.1 Large prototype

The setup is shown in figure 1, left side. It consists of the field cage, the endplate with the GEM module and two scintillators forming the trigger system . The readout is done using an $`11.4\mathrm{M}Hz`$ Flash ADC. The chamber volume is filled with a gas mixture of $`93\%`$ argon, $`2\%`$ carbon dioxide and $`5\%`$ methane (TESLA TDR ). Cosmic muons are used as ionising particles. No magnetic field is applied ($`B=0\mathrm{T}`$).

The endplate contains the GEM module consisting of two ”standard” GEMs : A copper coated capton foil with holes of a double conical shape (inner(outer) diameter: $`55(70)\mu \mathrm{m}`$) and a pitch of $`140\mu \mathrm{m}`$. Each foil has a size of $`10\times \mathrm{\hspace{0.17em}10}\mathrm{c}m^2`$. They are mounted in cascade with a gap of $`1.5\mathrm{m}m`$ between the two GEMs (transfer gap) and between the GEM and the array of readout pads, respectively (induction gap). The corresponding fields are $`E_{\mathrm{induction}}=E_{\mathrm{transfer}}1.3\mathrm{kV}/\mathrm{cm}`$.
Our standard readout pads are rectangles of a size of $`2\times 6\mathrm{mm}^2`$, as proposed in the TESLA TDR. According to simulations  more sophisticated pad geometries may lead to a better point resolution due to larger charge sharing. Therefore we also performed measurements using chevron shaped pads. Since the studies are made without magnetic field, the pad sizes are scaled in $`x`$ according to the width of the charge cloud, governed by the diffusion coefficient $`𝒟(B)`$. In the drift volume, $`𝒟(B=0\mathrm{T})6𝒟(B_{\mathrm{T}ESLA}=4\mathrm{T})`$ holds true for our gas. Thus, we use rectangles with an extension in $`x`$ of $`14\mathrm{mm}`$ and the corresponding chevron size (fig. 1, right).
Measurements: Some already intensively studied GEM properties  are verified for a TPC: The gas amplification depends exponentially on the applied voltages across the GEMs and the stability of the gain is satisfactory: The gain variation equals $`1.9\%`$ over a period of 70 hours. Taking into account the correction for atmospheric pressure, we expect to reduce the variation to a value $`<1\%`$, which is needed to pin down the relative error $`\frac{\sigma _{\mathrm{d}E/\mathrm{d}x}}{\mathrm{d}E/\mathrm{d}x}`$ below $`5\%`$. The point resolution $`\sigma _x`$ (respectively $`\sigma _z`$) is determined calculating the residuals $`x_{\mathrm{m}eas}x_{\mathrm{f}it}`$, $`x_{\mathrm{m}eas}`$ being the centre of charge of the 3 d electron cluster assuming a linear charge distribution and $`x_{\mathrm{f}it}`$ the coordinate calculated from the linear track fit, see figure 2. The four left plots contain results of the small rectangular pads: $`x`$ and $`z`$ resolution as a function of the drift length (upper plots): Apart from very short drift distances, the charge cloud is larger than the pad size, and the resolution is dominated by diffusion. In $`x`$, this leads to an overall worse result than expected for TESLA ($`\sigma _{x,\mathrm{T}ESLA}140\mu m`$). The lower two plots show the degradation of the resolution with increasing track angle, due to the prolongation of the projection of the charge cloud on the corresponding coordinate ($`x`$ and $`z`$). The right plot in figure 2 contains the results using large pads. The resolution $`\sigma _x`$ is dominated by the pad size rather than by diffusion. For small drift distances the point resolution with chevrons is better than with rectangles. This effect of the pad geometry gets smaller with increasing drift distance. The corresponding resolution at TESLA can be estimated by scaling down the measured resolution by the factor corresponding to the ratio of our pad sizes.

### 2.2 Small prototype (”Mini TPC”)

A side view sketch of the ”mini TPC”  can be seen in figure 3. It is built to measure the ion feedback and features a short drift distance ($`16\mathrm{m}m`$) making it possible to run at a moderate cathode voltage ($`U_{\mathrm{c}ath}1.5\mathrm{k}V`$) which allows the measurement of the cathode current $`I_{\mathrm{c}ath}`$ (in the order of few $`\mathrm{n}A`$). Dividing $`I_{\mathrm{c}ath}`$ by the sum of all positive currents gives the ion feedback . The used gas mixture is the same as in the previous case (TDR gas), the ionisation is produced by a radioactive iron source. In the plot of figure 3 the ion feedback is shown as a function of the drift field. Both induction and transfer field are adjusted to $`1\mathrm{k}V/cm`$, the GEM voltages are $`U_{\mathrm{G}EM2}=390\mathrm{V}`$ and $`U_{\mathrm{G}EM1}=370\mathrm{V}`$. For a drift field of $`250\mathrm{V}/m`$, as foreseen for the TESLA detector, the feedback is $`4\%`$.

## 3 Future Plans

The future measurement programme includes GEM operation in a magnetic field up to $`5\mathrm{T}`$ and TPC studies in the DESY electron test beam ($`6\mathrm{G}eV`$). The determination of the point resolution under more realistic conditions will be possible and it is planned to study the two track resolution. Resolution studies will also be carried out with different GEM tower geometries. Simulations will accompany the measurement programme.

## 4 Summary and Conclusion

Studies on the GEM readout for the TESLA TPC are performed. The gain variation without correction for atmospheric pressure is lower than $`2\%`$ which is promising for $`\mathrm{d}E/\mathrm{d}x`$ measurements. Resolution studies in $`x`$ with $`B=0\mathrm{T}`$ show different behaviour of rectangular and chevron pads favouring chevrons especially for short drift distances. The $`z`$ resolution is $`0.51.0\mathrm{m}m`$. It has to be taken into account that resolution measurements depend strongly on the diffusion — thus on the gas mixture and on the geometry of the GEM tower. Finally, for a double GEM structure, the ion feedback is determined to be in the order of a few percent.

## Acknowledgement

I would like to thank T. Behnke, N. Ghodbane, T. Kuhl, T. Lux and F. Sefkow for their support and contributions.


text:
# Stellar yields in CNO from rotating stellar models

## 1. Rotation and the CO–core masses

In Fig. 1, the variations of the CO–core masses as a function of the initial mass are indicated for various initial metallicities and for rotating and non–rotating stellar models. In stars with initial masses below 30–40 M, both at low and high metallicity, the CO cores are enhanced by rotation. As an example, at solar metallicity a rotating 20 M stellar model with $`v_{\mathrm{ini}}=300`$ km s-1 has a CO–core mass equal to that of a non–rotating 26 M stellar model. Thus one expects that rotating models will inject more carbon and oxygen than their non–rotating counterparts.

The increase of the CO–cores due to rotation is less marked at low metallicity. Paradoxically this is a consequence of the more efficient rotational mixing at low $`Z`$ (Meynet & Maeder 2002b). Indeed the more efficient mixing continuously replenishes in hydrogen the H–burning shell. As a consequence the H–burning shell progresses less rapidly outwards and thus prevents the He–core mass and therefore the CO–core mass to increase as much as at higher metallicity.

For stars above 30–40 M, at solar metallicity, the CO–cores in rotating models are smaller than in non–rotating ones. Due to rotational mixing, these models enter into the Wolf–Rayet phase at an earlier stage of their evolution. During the Wolf–Rayet phase, which lasts longer in rotating models, the stars loose mass at a higher rate, producing smaller final masses. As we shall see below, the CNO yields for these high mass stars are little affected by rotation.

## 2. Rotation and mass loss by stellar winds

In the framework of the radiative driven wind theory, it is possible to obtain an expression for the ratio of the mass loss rates from a rotating and a non–rotating star lying at the same position in the HR diagram (Maeder & Meynet 2000b). Typically for a 20 M star model at the end of the MS phase, the enhancement factor due to rotation is at most 1.67 at break–up while for a 60 M stellar model, in the same conditions, the factor amounts to 3.78. Rotation does not only increase the quantities of mass lost by stellar winds but also induces anisotropies in the winds. Typically, polar winds are expected for fast rotating hot stars, which, as a consequence, will loose small amounts of angular momentum and thus reach more easily the break–up limit (Maeder 2002). Other effects induced by rotation may change the quantity of mass lost by stellar winds: the increase of the MS lifetimes, the bluer evolutionary tracks, the easier Wolf–Rayet star formation.

Fig. 2 shows the final masses obtained for different initial masses, velocities and metallicities. At low metallicity, for stars with initial masses below 40 M and with $`\upsilon _{\mathrm{ini}}`$= 300 km s-1, rotating and non–rotating stellar models end with final masses close to their initial values. Only in the case of the 60 M there is a significant difference. One expects that the difference will increase when the mass increases and/or when the initial velocity increases. Concerning this last point it is worthwhile to recall here that there are some indications that the relative number of fast rotators increases when the metallicity decreases (Maeder et al. 1999). This fact, together with the effect of rotation on both mixing and mass loss, may have interesting consequences for the nucleosynthesis expected from the first stellar generations in the Universe.

At solar metallicity, mass loss has an important impact on the values of the final masses. All stars with initial masses above 40 M end with final masses between 10–15 M. In general, as expected, rotating models produce smaller final masses. One notes also that the new mass loss rates used in the recent computations are smaller than those used in the stellar grids of Schaller et al. (1992). This has some impact on the carbon yields from the most massive stars as we shall see in Sect. 4.

## 3. Effects of rotation on the carbon and oxygen yields

Stellar yields at the metallicities $`Z`$ = 0.00001 and 0.004, for 4He, 12C, 14N and 16O from rotating and non–rotating models, in the mass range between 2 and 60 M, are discussed in Meynet & Maeder (2002b). Let us recall here the main results. Figure 3 shows the yields in 12C for non–rotating and rotating stellar models ($`v_{\mathrm{ini}}`$ = 300 km s-1) for two metallicities and for mass ranges between 2 and 120 M. For the metallicities considered here, the effects of the stellar winds are small. We see that the yields in carbon are generally enhanced in the rotating models. The enhancement factors are equal to 1.3 and 2.4 for the 20 M models at $`Z`$ = 10-5 and 0.004 respectively. Larger carbon yields in rotating models are a consequence of the greater CO–core masses in rotating models (see Sect. 1 above). Interestingly, the yields from the rotating models are not very different from those computed by Maeder (1992) for $`Z`$ = 0.001 computed without rotation but with a moderate overshooting, while the rotating models were computed without overshooting. Rotation, by enlarging the CO core masses, acts in this respect as an overshoot.

Figure 4 shows the yields in 16O. As for carbon, the yields in oxygen are enhanced by rotation. For the 20 M, the enhancement factor is 1.8 at $`Z`$ = 10-5 and 1.6 at $`Z`$ = 0.004. Again the yields from rotating models (without overshooting) are very similar to those obtained by Maeder (1992).

At solar metallicity, the main difference with respect to the situation found at low metallicity essentially concerns the high mass star range (stars with initial masses above about 40 M). For this mass range, mass loss by stellar winds becomes the dominant effect. As a consequence the yields in carbon and oxygen in these stars are much less affected by rotation. As a numerical example, the yields in carbon and oxygen of a 60 M star model differ by less than 11% between the rotating and the non–rotating cases.

It was proposed by Maeder (1992) that massive stars (going through a WC phase) may be important contributors of carbon at high metallicity, an idea further sustained by various chemical evolution models (Prantzos et al. 1994; Carigi 2000; Gustafsson et al. 1999). Since the mass loss rates have been reduced since the work of Maeder (1992), the present yields of carbon from WC stars are smaller, and those of oxygen greater. However these stars remain important sources of carbon at solar and higher metallicities.

## 4. Effects of rotation on the nitrogen yields

Nitrogen is produced by transformation of carbon and oxygen through CNO cycle in H–burning zones. Nitrogen is said to have a primary origin when the carbon (oxygen) used for its synthesis is produced by the star itself. The secondary channel corresponds to the case when the carbon transformed into nitrogen is the one initially present in the star. When nitrogen is primary, its abundance in the interstellar medium evolves in lockstep with other primary elements like, for instance, oxygen. In that case, the N/O ratio remains constant when the metallicity increases. When nitrogen is mainly produced by the secondary channel, the N/O ratio is expected to increase steeply with the increasing metallicity.

If, at low metallicity, observations clearly require some source of primary nitrogen, the identification of these sources (massive or intermediate mass stars) remains controversial. Interestingly, there are some observational features which may give some hints: let us suppose that we observe starbursts in galaxies. If primary nitrogen is produced by massive stars, one expects to observe a very small scatter of the N/O ratios in starbursts of different ages, since both nitrogen and oxygen are released at the same time 111However the observation of a small scatter of the N/O ratios at low metallicity does not necessarily imply that primary nitrogen is produced by massive stars. Indeed the systems may be sufficiently old for having allowed intermediate mass stars to have released their primary nitrogen. In that case, low star formation must be invoked in order to explain the observed low metallicity.. In case primary nitrogen is released by intermediate mass stars, there will be some time delay between the oxygen and the nitrogen release. Thus there is some chance to observe systems which have already been enriched in oxygen by the massive stars but which are in the process of being enriched in primary nitrogen by intermediate mass stars. One expects in that case some scatter in the observed values of the N/O ratios. Recently, numerous Damped Lyman Alpha (DLA) systems have been observed with values of the N/O ratios well below the plateau level. This strongly suggests that intermediate mass stars are responsible for the primary nitrogen production.

What are the predictions of the theoretical stellar models ? In the classical scenario primary nitrogen is produced during the thermal pulse phase of AGB stars undergoing HBB (see van den Hoek and Groenewegen 1997; Marigo 2001). The quantities of primary nitrogen produced and expelled depends on various parameters as the mass loss along the AGB and the strength of the third dredge–up. Recently, we have proposed a new scenario involving rotational diffusion: carbon and oxygen produced in the He–core migrate by rotational diffusion in the H–burning shell where they are transformed into primary nitrogen. In contrast to the classical scenario, this process occurs in the whole mass range and is thus not restricted to the intermediate mass range.

The left panel of Fig. 5 shows the chemical abundances as a function of the lagrangian mass in a rotating 3 M star model at $`Z=0.00001`$ at the beginning of the Thermal Pulse–AGB phase. The most striking feature is the very high 14N abundance obtained at the border of the He–core and in the outer convective zone. Its abundance is more than 70 times the initial metal content of the star ! The right panel shows the 14N stellar yields normalised to the initial metal content obtained for different initial mass stars at various metallicities. At very low metallicity, all rotating models produce quantities of new 14N which are multiples of the initial metal content. However, for the initial velocities considered here ($`v_{\mathrm{ini}}`$ = 300 km s-1), the intermediate mass stars are the main producers. For higher metallicities, the secondary channel for nitrogen production is the dominant one.

Why the process of primary nitrogen production only works at low metallicity ? At low metallicity rotational mixing of chemical species is more efficient due to the steeper gradients of angular velocity inside the stars. These steeper gradients arise as a result of less efficient transport by meridional circulation of the angular momentum in the outer shells (Meynet & Maeder 2002ab). In metal poor stars, the H–burning shell is also nearer from the He–burning core. This shortens the timescale for diffusion between these two zones.

Left panel of Fig. 6 shows the evolution of the N/O ratio as a function of the O/H ratio in the framework of a very simple model for the chemical evolution of galaxies (closed box model with instantaneous recycling). The hatched zones in Fig. 6 represent regions where many observed points are found (see caption). At low metallicity, only the rotating models involving the contribution of the intermediate mass stars can reproduce the plateau level (at about Log (N/O) = -1.7). Massive stars alone may contribute significantly only if they rotate much faster than the typical rotational velocities found at solar metallicity. While the evolution in the N/O versus O/H plane mainly depends on rotation, the evolution in the C/O versus O/H plane is mainly sensitive to the mass range of the stars contributing to the enrichments of the ISM. Thus, the combination of the diagrams C/O vs O/H, more sensitive to the mass interval, and of the diagram N/O vs. O/H, more sensitive to rotation, may be particularly powerful to disentangle the two effects of rotation and mass interval, and to specify the properties of the star populations responsible for the early chemical evolution of galaxies.

The DLAs present also new and interesting constraints on the source of primary nitrogen at low metallicity. Recent discussions of this point may be found in the paper of Pettini et al. (2002) as well as in the contributions by Henry, Molaro in the present volume. As explained above, the numerous DLAs observed with N/O ratios below the plateau level favour a primary nitrogen production by intermediate mass stars. Secondly the time delay required to explain the frequency of the systems with low N/O ratios points towards a relative long time delay between the release of oxygen and that of nitrogen (of about 700 Myr according to Pettini et al. 2002). This is well in agreement with the predictions of the rotating models (see Fig. 7).

Finally some authors suggest that the N/O values of DLAs are distributed mainly along two plateaux. One is at Log(N/O) $``$ -1.7, as shown on the left panel of Fig. 6. A second plateau seems to be present at about 0.7 dex below the first one (see the papers by Henry and Molaro in this volume). A possible interpretation of the lower plateau is that it is produced by stars with initial masses above 8 M (see Molaro in this volume). Interestingly, the present massive fast rotating stellar models would be in agreement with this viewpoint. As shown on the left panel of Fig. 6, the contribution of rotating stars with masses superior to 8 M and $`v_{\mathrm{ini}}=400`$ km s-1 predict a N/O value at $`Z`$ = 10-5 about 0.7 dex below the high plateau value.

## 5. Conclusion

We focused our discussion here on the main isotopes of carbon, oxygen and nitrogen. However rotation may also affect the other less abundant CNO isotopes. In particular, 13C and to a lesser extent 17O may also be abundantly produced in rotating models and have a primary production channel at low metallicity. As an example, at $`Z`$ = 0.00001, the amounts of 13C synthesized in rotating intermediate mass star models may be nearly 10 000 times the amount produced in non–rotating models. Interestingly some primary 13C production seems to be required by the observations (Prantzos et al. 1996). Important productions of primary 22Ne are also predicted by the rotating models.

## References

Angulo, C., Arnould, M., Rayet, M., et al. (the NACRE Collaboration) 1999, Nucl. Phys., A656, 3

Carigi, L. 2000, Revista Mexicana de Astronomia y Astrofisica, 36, 171

Gustafsson, B., Karlsson, T., Olsson, E., Edvardsson, B., Ryde, N. 1999, A&A 342, 426

Heger, A., Langer, N. 2000, ApJ 544, 1016

Henry, R.B.C., Edmunds, M.G., Koeppen, J. 2000, ApJ 541, 660

van den Hoek L.B., Groenewegen M.A.T. 1997, A&ASS 123, 305

Maeder, A. 1992, A&A 264, 105

Maeder, A. 2002, A&A 392, 575

Maeder, A., Grebel, E.K., Mermilliod, J.–C. 1999, A&A 346, 459

Maeder, A. & Meynet, G. 2000a, ARA&A, 38, 143

Maeder, A. & Meynet, G. 2000b, A&A 361, 159

Maeder, A. & Meynet, G. 2001, A&A 373, 555

Maeder, A. & Zahn, J.–P. 1998, A&A 334, 1000

Marigo P. 2001, A&A 370, 194

Meynet, G., Maeder, A. 2002a, A&A 381L, 25

Meynet, G., Maeder, A. 2002b, A&A 390, 561

Pettini, M., Ellison, S.L., Bergeron, J., Petitjean, P. 2002, A&A 391, 21

Prantzos, N., Aubert, O., Audouze, J. 1996, A&A 309, 760

Prantzos, N., Vangioni–Flam, E., Chauveau, S. 1994, A&A 285, 132

Schaller, G., Schaerer, D., Meynet, G., Maeder, A. 1992, A&AS 96, 269

Zahn, J.-P. 1992, A&A 265, 106.