text:
# Iron Abundance Profiles of 12 Clusters of Galaxies Observed With BeppoSAX

## 1. Introduction

The presence of metals in the hot gas contained within clusters of galaxies provides important clues as to how metals are deposited into the intracluster medium (ICM) from the constituent galaxies of the cluster. The two most likely mechanisms for injecting metals into the ICM are through ram pressure stripping of the gas from galaxies by the ICM (e.g., Gunn & Gott 1972; Gaetz, Salpeter, & Shaviv 1987) and galactic winds (e.g., Ostriker & Yahil 1973; Metzler & Evrard 1994). Recent work has suggested that negative radial abundance gradients are a common feature in clusters. Analysis of ASCA data have indicated that in some clusters the metal abundance declines from $``$40%–60% of solar in the center to $``$20% of solar at a distance of 0.5–1 Mpc (Ikebe et al. 1997; Ezawa et al. 1997; Sarazin, Wise, & Markevitch 1998; Finoguenov, David, & Ponman 1999; Dupke & White 2000a,b)

Most of the above studies, however, dealt with groups or cooler clusters. The large, energy-dependent Point Spread Function (PSF) of ASCA preferentially scatters hard photons, making spatially-resolved spectroscopy problematic for hot clusters. Whereas the problem is less severe for cooler ($`<5`$ keV) clusters, substantial errors in the radial temperature and abundance profiles can occur if the PSF is not dealt with properly (see Takahashi et al. 1995). As such, little observational work has been done on the radial abundance profiles of hotter clusters. There has also been correspondingly little theoretical work performed on the predicted abundance profiles of clusters from hydrodynamical simulations. White (1999) has derived PSF-corrected abundance profiles for a large sample of clusters observed with ASCA, but the deconvolution technique used to correct for the PSF unfortunately led to rather large uncertainties in the derived abundance profiles. As a result, only weak constraints could be placed on the presence of abundance gradients from these data.

BeppoSAX is better-suited to perform spatially-resolved spectroscopy on hot clusters. The PSF of BeppoSAX is one-half that of the ASCA GIS, and more importantly, is only weakly dependent on energy. The abundance profiles of five higher temperature clusters have already been derived from BeppoSAX, finding at least marginal evidence for abundance gradients in A2029, A2256, A3266, and PKS0745-191, with no evidence for a gradient in A2319 (Molendi & De Grandi 1999; Molendi, De Grandi, & Fusco-Femiano 2000; De Grandi & Molendi 1999a,b; Molendi et al. 1999). In this Paper, we increase this sample size by analyzing BeppoSAX archival data for 12 clusters of galaxies, and derive radial abundance profiles for each cluster, nine of which have temperatures greater than 5 keV. Throughout this paper, we assume $`H_0=50`$ km s-1 Mpc-1 and $`q_0=0.5`$.

## 2. Sample and Data Reduction

From the BeppoSAX Science Data Center (SDC) archive (available at http://www.sdc.asi.it/sax\_main.html) we have obtained data for 12 clusters of galaxies. Within this sample, eight of them possess cooling flows (A85, A496, A1795, A2029, A2142, A2199, A3562, and 2A0335+096) and four do not (A2163, A2256, A2319, and A3266). We define a cluster as having a cooling flow if the cooling rate is greater than $`20M_{}`$ yr-1; see Peres et al. 1998 and White, Jones, & Forman 1997. All of the clusters are at low redshift ($`0.03z0.09`$) except for A2163, which has a redshift of $`z=0.203`$. We analyze data taken with the Medium Energy Concentrator Spectrometer (MECS) onboard BeppoSAX. The MECS (see Boella et al. 1999 for details) is composed of two identical gas scintillation proportional counters (three detectors before 1997 May 9) that are sensitive in the 1.3–10.5 keV energy range. The event files for all 12 clusters were subjected to the standard screening criteria of the BeppoSAX SDC.

Since we are interested in deriving the radial abundance profiles for these clusters, it is important to account for scattering owing to the (PSF) of the MECS instrument. Fortunately, the detector + telescope PSF of BeppoSAX is nearly independent of energy. This is because the Gaussian PSF of the MECS detector improves with increasing energy, while the PSF of the grazing incidence Mirror Unit degrades with increasing energy (D’Acri, De Grandi, & Molendi 1998), leading to a partial cancellation when these two effects are combined.

To correct for the PSF we have used the routine effarea, available as part of the SAXDAS 2.0 suite of BeppoSAX data reduction programs. The program effarea convolves the surface brightness profile of the cluster (which has been determined using ROSAT PSPC data; see Mohr, Mathiesen, & Evrard 1999 and Ettori & Fabian 1999) with the PSF of the MECS to determine the amount of contamination in the spatial region in question at each energy from other regions of the cluster. Thus, an energy-dependent correction vector is formed. When this correction vector is multiplied by the observed spectrum, the PSF-corrected energy spectrum is obtained for the region. In practice, this information is incorporated into the auxiliary response file (the .arf file), which is subsequently used in the spectral fitting. This task also corrects for vignetting. A more complete description of the task is given in Molendi (1998) and D’Acri et al. (1998). This correction does not appear to affect the spectrum greatly; D’Acri et al. (1998) and Kaastra, Bleeker, & Mewe (1998) found only small changes between the uncorrected and corrected temperature profiles for Virgo and A2199, respectively. In addition, D’Acri et al. (1998) found that the correction vectors amounted to 5% or less for energies above 3 keV outside of the innermost bin (the innermost bin loses some flux via scattering but gains very little from photons scattered in from greater off-axis radii). Our inner three spatial bins are identical in angular extent to the bins of D’Acri et al. (1998), so we assume that our correction vectors are as small as theirs. In addition, we find only a small difference between our corrected and uncorrected abundance profiles for our sample.

We extracted spectra in four annular regions for each cluster, with inner and outer radii of $`0^{}2^{}`$, $`2^{}4^{}`$, $`4^{}6^{}`$, and $`6^{}9^{}`$. At $`9^{}`$ the telescope entrance window support structure (the strongback) strongly absorbs X-rays, making spectroscopy difficult in this area of the detector. In addition, for off-axis angles greater than $`10^{}`$, the departure of the PSF from radial symmetry becomes noticeable (Boella et al. 1997). With this in mind, we ended the radial profiles at $`9^{}`$. At this radius, the profiles extended to 17%–33% of the virial radius, where $`r_{virial}=3.9(T/10\mathrm{keV})^{1/2}`$ Mpc (see Evrard, Metzler, & Navarro 1996), for all clusters except A2163, for which the profile extended to 55% of $`r_{virial}`$. We also extracted one global spectrum ($`0^{}9^{}`$) and also two spectra covering the two regions $`00.075r_{virial}`$ and $`0.0750.173r_{virial}`$. Background was obtained from the deep blank sky data provided by the SDC. We used the same region filter to extract the background as we did the data, so that both background and data were affected by the detector response in the same manner. The energy channels were regrouped to contain at least 25 counts.

For each cluster, XSPEC Version 11.0 was used to fit the spectrum. The MECS2 and MECS3 (and MECS1 when available) data were fit separately, but with the same temperature, metallicity, and normalization. We assumed a MEKAL model with an absorption component fixed at the Galactic value, allowing the temperature and metallicity to vary. For the iron-to-hydrogen ratio we use a value of $`4.68\times 10^5`$ (Anders & Grevesse 1989) since this is the value assumed for most previous abundance determinations, although it should be noted that more recent determinations of the Fe/H ratio point to a value of $`3.24\times 10^5`$ (Ishimaru & Arimoto 1997). Because of a systematic shift of 45–50 eV at 6.6 keV in the MECS channel–to–energy conversion, we have also allowed the redshift to vary (F. Fiore, private communication). This systematic shift was evident in our sample; when the redshift was allowed to vary, the measured redshift was less than the optically-determined redshift in all 12 clusters, and inconsistent with the optically-determined redshift at the 90% confidence level for eight of them. A modest decrease in the reduced $`\chi ^2`$ also occurred for most of the clusters when the redshift was allowed to vary. However, freeing the redshift did not affect the values obtained for the temperature and metallicity by more than 5% (and in most instances much less), so the results do not depend significantly on this shift. We also fit only the 3.0–10.5 keV part of the spectrum because of excess flux found below 3.0 keV that might be the result of uncertainties in the calibration of the MECS instruments; see Irwin & Bregman (2000) for a more complete description. All quoted errors are 90% confidence levels for one interesting parameter ($`\mathrm{\Delta }\chi ^2=2.71`$) unless otherwise noted. Below we present the results from the abundance analysis. The results of the temperature analysis are given in Irwin & Bregman (2000). One additional cluster (A3562) not analyzed in Irwin & Bregman (2000) is presented here.

## 3. Abundance of the Hot Gas

### 3.1. Global Abundances

For hot clusters it is the Fe-K$`\alpha `$ line that primarily determines the metal abundance measurement. The fits to the global spectra were adequate, with $`\chi _\nu ^2<1.20`$ in all cases. For the eight clusters harboring cooling flows, the weighted average emission-weighted global abundance was $`0.37\pm 0.04`$ (the error is given by the weighted standard deviation for the subsample). The four clusters lacking cooling flows had an average emission-weighted global abundance of $`0.27\pm 0.02`$. The ratio of cooling flow to non-cooling flow global abundance is $`1.4\pm 0.2`$, in good agreement with the value of 1.6 obtained by Allen & Fabian (1998) with a sample of 30 clusters observed with ASCA when they used a single component thermal model.

### 3.2. Radial Abundance Profiles

The radial metal abundance profile for each cluster is shown in Figure 1. Four of the clusters have had abundance gradients determined with ASCA, and are generally consistent within the errors with the profiles determined here: A496 (Dupke & White 2000b; Finoguenov et al. 1999), A2199 (Finoguenov et al. 1999), A2029 (Sarazin et al. 1998; Finoguenov et al. 1999) and 2A0335+096 (Kikuchi et al. 1999). In addition, analysis of this same BeppoSAX data for A2029, A2256, A2319, and A3266 by other authors are in agreement with our results (Molendi & De Grandi 1999; Molendi et al. 1999; De Grandi & Molendi 1999), with the exception of A2256 for which Molendi et al. (2000) found a significantly decreasing profile whereas ours is flat. However, our inability to find an abundance gradient is the result of the short observing time of the data we used; the total observing time of the Molendi et al. (2000) observation was 2.6 times longer than ours, and they were able to place much tighter constraints on the abundance than we were. In both studies, the 1$`\sigma `$ error bars overlapped one another.

In general, the abundance profiles decline with radius. Such a trend was also found in clusters observed by ASCA (Finoguenov et al. 1999; Kikuchi et al. 1999; Dupke & White 2000a), although there are notable exceptions to this trend such as A1060 (Tamura et al. 1996). These previous studies focused primarily on cooler clusters ($`kT<5`$ keV). This study confirms the ubiquity of negative abundance gradients in hotter clusters (nine of our 12 clusters have $`kT>5`$ keV), as illustrated in Figure 2, where all 12 abundance profiles, normalized to the global abundance for each cluster, are plotted versus radius in units of the virial radius. A strong negative gradient in the abundance profile is evident out to 25% of the virial radius. However, there are a few points quite discordant with this trend at large radii. The clusters responsible for these discordant points are A1795, A2142, and A2163.

We investigated the possibility that the use of single component temperature models in the central regions of cooling flow clusters has biased the determination of the abundance in these regions. The addition of a cooling flow component did not lead to a significant improvement in the fits, and the cooling rate was only weakly constrained. This is a result of excluding data below 3.0 keV, where the effect of cooling flows on the spectra is most pronounced. When we fixed the cooling rate at values determined by other satellites, the best-fit abundance value was actually slightly higher than in the single component case by a few percent, indicating that the central enhancement in the abundance is real.

We have fit a linear model of the form $`(Z/Z_{mean})=ab(r/r_{virial})`$ to the data. While including all values, the fit is not very good, with $`\chi ^2=72.2`$ for 47 degrees of freedom, and can be excluded at $`>98\%`$ confidence level. Removal of the four discordant data points yielded an excellent fit, with $`\chi ^2=46.1`$ for 43 degrees of freedom for values of $`a=1.212`$ and $`b=2.705`$ (2.149–3.250 at 90% confidence). Although it is not clear why four of the clusters show a $``$2$`\sigma `$ departure from the observed trend outside $``$25% of the virial radius, the derived slope is nonetheless an accurate representation of the data inside of 25% of the virial radius. Excluding non-cooling flow clusters in addition to the discordant points led to the same best-fit slope.

### 3.3. Comparison To Theoretical Models

Surprisingly little theoretical work has been performed to date to quantify the behavior of abundance gradients in clusters. One such study is that of Metzler & Evrard (1997), who used gas dynamical simulations which included the effects of galactic winds to generate an ensemble of 18 realizations of clusters spanning a wide range of temperatures. Assuming that the iron present in the gas resulted from the ejection of metals from galaxies, a composite iron abundance profile for the 18 realizations was determined. It was found to decrease significantly with increasing radius. This resulted from the spatial distribution of the galaxies being more centrally condensed than the original primordial X-ray gas. Since the ejected gas traces the galaxies, an abundance gradient is created relative to the primordial gas.

The composite iron abundance profile from the 18 three-dimensional realizations of Metzler & Evrard (1997) is shown superposed on the BeppoSAX data in Figure 2 as a solid line. The composite profile was normalized assuming an average global abundance of 40% of the solar value. Metzler & Evrard (1997) scaled each cluster in radial units of $`r_{170}`$ before merging the profiles, where $`r_{170}`$ is the radius within which the mean density is 170 times the critical density. This is very close to our adopted value of $`r_{virial}`$. The normalized abundance profile of Metzler & Evrard (1997) matches the data very well, ranging from a value of 1.3 in the center to about 0.6 at 25% of the virial radius.

We note that recent work by Finoguenov et al. (1999) and Dupke & White (2000b) has shown that significant radial gradients exist in the elemental abundance ratios in clusters. This would argue against the idea that a difference in the spatial distribution between the galaxies and the gas is responsible for the observed abundance gradient if it is assumed that both Type Ia and Type II elements were released into the ICM under the same conditions. The observational evidence implies that Type Ia and Type II elemental contamination is not the same in the inner and outer parts, with Type Ia ejecta being dominant in the center. Dupke & White (2000b) propose that the observed abundance ratios can be explained by a phase of early, vigorous protogalactic wind driven by Type II supernovae that spread metals uniformly throughout the cluster, followed by a phase of less vigorous Type Ia supernovae wind that was suppressed by the central dominant galaxy in the center of the cluster. Thus, the measured iron abundances originate from a combination of Type II supernovae (having no radial gradient) and Type Ia supernovae (having a radial gradient). The gradient we measure is the linear combination of these two effects.

### 3.4. Cooling Versus Non-cooling Flow Clusters

Inspection of the abundance profiles of the eight cooling flow clusters reveals a radially decreasing trend in all of them. The gradient is quite significant in A85, A496, A2199, A3562, and 2A0335+096, and less significant in A1795, A2029, and A2142. This seems to be a universal trend in cooling flow clusters since all nine cooling flow clusters in the sample of Finoguenov et al. (1999), as well as all three clusters in the samples of Kikuchi et al. (1999) and Dupke & White (2000a) show the same trend, at least marginally. No clear trend exists among the non-cooling flow clusters. The profile of A2256 is flat, although a decline with radius cannot be ruled out owing to the large errors (a decline was seen with a longer BeppoSAX observation by Molendi et al. 2000; see Figure 1). A2319, A3266, and A2163 show some evidence for a negative gradient although the profiles are formally consistent with a flat profile. Results with ASCA for clusters with little or no cooling flows are also mixed. Whereas Coma, A401, A1060, and A2670 show a flat abundance profile, A399 has an abundance gradient (Watanabe et al. 1997; Fujita et al. 1996; Tamura et al. 1996; Finoguonov et al. 1999). It is possible that mergers that could be responsible for the disruption of the cooling flow in clusters can also destroy abundance gradients. The degree to which the abundance gradient remains intact may be an indicator of how recent and/or violent the merger was. All four non-cooling flow clusters are merger candidates. Of the five cooling flow clusters with significant gradients, only A85 is a merger candidate. Of the three cooling flow clusters with only weak evidence for an abundance gradients two are recent merger candidates. Oegerle & Hill (1994) conclude that the velocity distribution of galaxies within A1795 suggest that the cluster underwent a merger with a relatively small subcluster such that the cooling flow was not disrupted. Oegerle, Hill, & Fitchett (1995) also find evidence for subclustering from the galaxy velocity distribution in A2142. ROSAT PSPC and HRI images of A2142 also suggest that a merger is in process (Henry & Briel 1996). It is possible that mergers have mixed the gas of these two clusters enough to decrease the significance of an abundance gradient without destroying the cooling flow. A2029 does not show evidence for a merger, although the lack of a significant abundance gradient may be the result of poor photon statistics rather than a lack of a gradient. Further studies of A2029 will be needed to confirm this.

As suggested by Fabian et al. (1994) and Allen & Fabian (1998), the presence of abundance gradients in cooling flow clusters provides a reasonable explanation as to why cooling flow clusters have higher emission-weighted metal abundances than their non-cooling flow counterparts. If the cooling gas at the center of clusters has a higher abundance than the outer regions, the strongly peaked X-ray surface brightness distribution characteristic of cooling flow clusters will preferentially weight the global emission-weighted abundance towards the central value, with the low-abundance outer regions contributing relatively little to the global value. Such an effect is not present in non-cooling flow clusters owing to the broader, less-peaked X-ray surface brightness distribution of these clusters. But is this the only reason for the higher measured metallicity in cooling flow clusters, or is there really an abundance difference between cooling flow and non-cooling flow clusters? To test this, we determined abundances in two regions for each cluster, a circle with a radius of $`0.075r_{virial}`$ and an annulus with inner and outer radii of $`0.075r_{virial}`$ and 0.173$`r_{virial}`$. The outer radius was set by the need to extend the annulus no further than $`9^{}`$ for any cluster, and the inner radius was set to be larger than the cooling radius for any of the cooling flow clusters.

The results are shown in Table 1. Inside of $`0.075r_{virial}`$ the cooling flow subsample had an average abundance of $`0.42\pm 0.06`$, while the non-cooling flow subsample had an average of $`0.33\pm 0.04`$. In the outer region, the cooling flow subsample had an average abundance of $`0.30\pm 0.02`$, while the non-cooling flow subsample had an average of $`0.24\pm 0.03`$. Thus, there appears to be an abundance gradient in non-cooling flow clusters as well as cooling flow clusters. The presence of abundance gradients in non-cooling flow clusters (especially in hot clusters) would seem to rule out claims that previously observed abundance gradients were artifacts of incorrect spectral modeling due to either the presence of cooling gas at the center of clusters or uncertainties in modeling the Fe-L line complex for cooler clusters.

Although the average abundances of the cooling flow and non-cooling flow clusters are marginally consistent within the errors within the cooling region, the data suggests that there is a discrepancy between the abundance of cooling and non-cooling flow clusters outside of the cooling radius. Of course, any conclusions drawn from a subsample of four clusters should be viewed with caution. A larger sample of clusters will be necessary to confirm if the abundance differences between cooling flow and non-cooling flow clusters are real. This result is in the opposite sense expected if higher metallicity gas from the cooling flow region has been mixed in with the lower metallicity gas of the outer regions as the result of a merger in a present-day non-cooling flow cluster. This trend also appeared in a sample of clusters observed by ASCA. Fukazawa et al. (1998) analyzed 26 clusters with temperatures of at least 3.0 keV, and found global abundances after having excised the inner cooling regions. The average abundance for the 16 cooling flow clusters was $`0.29\pm 0.07`$, while the average for the ten non-cooling flow clusters was $`0.23\pm 0.05`$. Although at a lower significance than our BeppoSAX sample, the trend is in the same direction. Fukazawa et al. (1998) did not claim to extract their spectra in units of the virial radius, which is the unit of distance that allows the most direct comparison of radial properties for clusters of varying sizes and temperatures. This might wash out some of the information on the abundance gradient when the profiles are combined, and lead to a result of lower statistical significance.

Allen & Fabian (1998) present the interesting scenario in which metallicity gradients can be explained (at least in part) by a significant fraction of the metals in clusters residing on grains in cluster cores. Since the lifetime of a grain to sputtering in the hot ICM is inversely proportional to the electron density of the hot gas (e.g., Draine & Salpeter 1979), more metals would be released into the gas phase in the high density cooling flow regions of clusters than in lower density regions, providing the dust grains are $`10`$ microns. By extending this argument outside the core region, it might be possible to explain the difference in metallicity at large radii between cooling flow and non-cooling flow clusters. If the density of gas outside the core is greater in cooling flow clusters than in non-cooling flow clusters, more metals will be released into the gas phase for the former. A difference in density among cooling and non-cooling flow clusters outside of 200 kpc was observed by White et al. (1997), who deprojected the surface brightness profiles of 207 clusters observed with Einstein Observatory and found that the electron densities of clusters with cooling rates greater than $`10M_{}`$ yr-1 were larger than for clusters with cooling rates less than $`10M_{}`$ yr-1 out to a radius of 1 Mpc. For a 7 keV cluster, our $`0.0750.173`$$`r_{virial}`$ bin corresponds to 240–560 kpc. At a radius of 450 kpc, White et al. (1997) found that the cooling flow clusters had an electron density of $``$$`7\times 10^4`$ cm-3, whereas the non-cooling flow clusters had an electron density of $``$$`3\times 10^4`$ cm-3. For a grain lifetime of $`2\times 10^6`$ $`a/n_e`$ yrs (where $`a`$ is the size of the grain in microns and $`n_e`$ is the electron density of the hot gas in cm-3), grains with a size of around 2 microns would survive in non-cooling flow clusters for a Hubble time while being sputtered in cooling flow clusters. While rather large, dust grains of this size are not unreasonable. In addition, this argument assumes that the dust was deposited into the ICM at the time the cluster was created. In reality, dust is deposited into the ICM throughout the history of the cluster, allowing smaller dust grains to survive in low density non-cooling flow clusters while still being sputtered in high density coling flow clusters. This in turn can lead to a higher abundance in the outer regions of cooling flow clusters since more of the grains have been sputtered away owing to the higher electron density.

However, more recent determinations of the density profiles of clusters determined with ROSAT PSPC data indicate that the non-cooling flow clusters in our sample do not have a lower electron density in the 0.075–0.173 $`r_{virial}`$ region than the cooling flow clusters (J. Mohr, private communication). At 10% of the virial radius, the average electron density for the non-cooling flow clusters A2256, A2319, and A3266 was $`1.80\pm 0.20\times 10^3`$ cm-3, while for the cooling flow clusters A85, A496, A1795, A2029, A2142, A2199, and A3562 the average electron was $`1.92\pm 0.46\times 10^3`$ cm-3 Thus, it does not appear likely that the difference in abundances outside the cooling radius can be the result of an increased sputtering of grains in cooling flow clusters.

Since non-cooling flow clusters are usually merger candidates, one possible explanation for the difference in abundances between cooling and non-cooling flow clusters is the mixing of gas via mergers from larger ($`>0.175r_{virial}`$) radii rather than from gas interior to this region. If the iron abundance of the gas continues to decline with increasing radius beyond $`0.2r_{virial}`$ as the simulations of Metzler & Evrard (1997) imply, the mean metallicity in our 0.075–0.173 $`r_{virial}`$ bin will be lowered upon mixing with this low metallicity gas during a merger. A2142, a cooling flow cluster that is thought to be undergoing a merger and shows only weak evidence for an abundance gradient, also has the lowest metallicity in the 0.075–0.173 $`r_{virial}`$ bin among cooling flow clusters in our sample, consistent with the merger-mixing theory. However, A1795, which is also believed to have recently experienced a merger has a high metallicity value in the 0.075–0.173 $`r_{virial}`$ bin, although the uncertainty is rather large. A larger sample of clusters with abundances determined out to larger radii will be needed to explore this topic further.

## 4. Conclusions

We have analyzed 12 clusters with the MECS instrument onboard BeppoSAX and find evidence for a negative abundance gradient in most of them. This gradient is present in all eight cooling flow clusters, and to a lesser significance in the non-cooling flow clusters. Since $`kT>5`$ keV for nine of the 12 clusters, this work extends the ubiquity of abundance gradients to hot clusters. The slope of the abundance gradient is in good agreement with cluster simulations performed by Metzler & Evrard (1997). In our sample cooling flow clusters had a higher metallicity than non-cooling flow clusters both in the inner regions of the cluster and the outer regions out to several hundred kpc; further work with a larger sample will be needed to confirm this result for clusters in general. Outside the cooling region, it seems unlikely that this difference can be explained by the more efficient sputtering of dust grains in cooling flow clusters, since the densities of cooling flow and non-cooling flow clusters are similar in this region. The mixing of low metallicity gas during a merger may be responsible for this difference.

We thank the anonymous referee for many useful suggestions and comments. We thank Chris Metzler for making his theoretical abundance profiles available to us. We also thank Joe Mohr for kindly providing us his ROSAT PSPC density profiles. JAI thanks Renato Dupke and Gus Evrard for many useful comments and conversations. This research has made use of data obtained through the BeppoSAX Science Data Center and the High Energy Astrophysics Science Archive Research Center Online Service, provided by the NASA/Goddard Space Flight Center. This work has been supported by Chandra Fellowship grant PF9-10009, awarded through the Chandra Science Center. The Chandra Science Center is operated by the Smithsonian Astrophysical Observatory for NASA under contract NAS8-39073.


format: text

text:
# Strong-Coupling Theory for Counter-Ion Distributions

## Abstract

The Poisson-Boltzmann approach gives asymptotically exact counter-ion density profiles around charged objects in the weak-coupling limit of low valency and high temperature. In this paper we derive, using field-theoretic methods, a theory which becomes exact in the opposite limit of strong coupling. Formally, it corresponds to a standard virial expansion. Long-range divergences, which render the virial expansion intractable for homogeneous bulk systems, are shown to be renormalizable for the case of inhomogeneous distribution functions by a systematic expansion in inverse powers of the coupling parameter. For a planar charged wall, our analytical results compare quantitatively with extensive Monte-Carlo simulations.

Recent years have witnessed a revival of the interest in classical charged systems. Specific attention has been paid to the failure of the Poisson-Boltzmann (PB) approach, which is known to give reliable results only in the limit of low-valency ions or high temperatures. Corrections to PB have been attributed to correlations between ions, or, more precisely, correlated ion-density fluctuations, and, if present, additional non-electrostatic interactions between ions. These corrections are particularly important for the interaction between macroscopic similarly charged objects, where they can lead to attractions. In as much as the PB approach is accurate for weakly charged systems, no systematic theory is available for the distribution of counter-ions around charged objects in the opposite limit of high-valency ions; moreover, it was not clear whether such a limit exists and whether it is physically meaningful. In this paper we show, using field-theoretic methods, that while PB corresponds to the asymptotically exact theory in the weak-coupling limit, (corresponding to low-valency ions or high temperatures), our novel strong-coupling theory becomes asymptotically exact in the opposite limit of high-valency ions or low temperatures and constitutes a physically sound limit. For the case of a planar charged wall, we give explicit results for the asymptotic density profile in the strong-coupling limit. We also have performed extensive Monte-Carlo (MC) simulations of this system. The resulting density profiles agree for weak and for strong coupling with predictions from PB theory and our strong-coupling theory, respectively. The strong-coupling limit is experimentally easily reached at room temperatures with highly charged walls and/or multivalent counter ions and thus relevant from the application point of view.

To proceed, consider the Hamiltonian of a system of $`N`$ ions of valency $`q`$ at an impenetrable and oppositely charged wall of number density of surface charges $`\sigma _s`$,

$$\frac{}{k_BT}=\underset{j<k}{}\frac{\mathrm{}_Bq^2}{|𝐫_j𝐫_k|}+2\pi q\mathrm{}_B\sigma _s\underset{j=1}{\overset{N}{}}z_j$$
(1)
where $`\mathrm{}_B=e^2/4\pi \epsilon k_BT`$ is the Bjerrum length which measures the distance at which two unit charges interact with thermal energy; in water $`\mathrm{}_B0.7nm`$. The dielectric constant is assumed to be homogeneous throughout the system. The Gouy-Chapman length, $`\mu =1/2\pi q\sigma _s\mathrm{}_B`$, measures the distance from the wall at which the potential energy of an isolated ion reaches thermal energy. Rescaling all lengths by $`\mu `$ according to $`𝐫=\mu \stackrel{~}{𝐫}`$, the Hamiltonian reads

$$\frac{}{k_BT}=\underset{j<k}{\overset{N}{}}\frac{\mathrm{\Xi }}{|\stackrel{~}{𝐫}_j\stackrel{~}{𝐫}_k|}+\underset{j=1}{\overset{N}{}}\stackrel{~}{z}_j$$
(2)
and thus only depends on the coupling parameter $`\mathrm{\Xi }=2\pi q^3\mathrm{}_B^2\sigma _s`$. Using the fact that the typical distance of an ion from the wall is in reduced units $`\stackrel{~}{z}1`$ (which holds both in the weak and strong-coupling limits, as we will show below) the confinement energy, the second term in Eq.(2), is of order unity per particle. The typical distance between ions scales as $`\stackrel{~}{r}\mathrm{\Xi }^{1/2}`$ (for $`\mathrm{\Xi }>1`$) or like $`\stackrel{~}{r}\mathrm{\Xi }^{1/3}`$ (for $`\mathrm{\Xi }<1`$ and assuming a liquid-like structure) and thus the repulsive energy between two ions at this distance, as determined by the first term in Eq.(2), is of order $`\mathrm{\Xi }^{1/2}`$ for $`\mathrm{\Xi }>1`$ or $`\mathrm{\Xi }^{2/3}`$ for $`\mathrm{\Xi }<1`$. It thus follows on rather general grounds that for a large coupling parameter, $`\mathrm{\Xi }>1`$, the ionic structure is dominated by mutual repulsions, suggesting crystallization, while for weak coupling $`\mathrm{\Xi }<1`$ these repulsions should be rather unimportant. In Fig.1 we show ion-distribution snapshots from Monte-Carlo simulations for $`\mathrm{\Xi }=0.1,10,10^4`$. For small $`\mathrm{\Xi }`$, Fig.1a, repulsive ion-ion interactions are indeed not playing a dominant role, the ion distribution is rather disordered and mean-field theory should work, since each ion moves in a diffuse cloud of neighboring ions. For large $`\mathrm{\Xi }`$, on the other hand, ion-ion repulsions are strong and ion-ion distances are large compared to the distance from the wall. In Fig.1c the ions form a flat layer on the charged wall. A two-dimensional one-component plasma is known to crystallize for plasma parameter $`\mathrm{\Gamma }125`$. From the definition of the two-dimensional plasma parameter we obtain the relation $`\mathrm{\Xi }=2\mathrm{\Gamma }^2`$. This leads to a crystallization threshold (in units of our coupling parameter) of $`\mathrm{\Xi }31000`$ meaning that all three systems in Fig.1 are below the crystallization threshold. However, mean-field theory is expected to break down, at least for the system shown in Fig.1c, because each ion moves, though confined by its immediate neighbors in the lateral directions, almost independently from the other ions along the vertical direction (which constitutes the soft mode). Let us now see how these notions can be made concrete in a field-theoretic formulation.

The partition function of a system of $`N`$ counter-ions interacting via Coulomb interactions $`v(𝐫)=1/r`$ among themselves and with a distribution $`\sigma (𝐫)`$ is given by

$$𝒵_N=\frac{1}{N!}\underset{j=1}{\overset{N}{}}d𝐫_j\mathrm{exp}\left\{q^2\mathrm{}_B\underset{j<k}{}v(𝐫_j𝐫_k)+q\mathrm{}_Bd𝐫\sigma (𝐫)\underset{j}{}v(𝐫𝐫_j)+\underset{j}{}h(𝐫_j)\right\}$$

where the field $`h`$ has been added to calculate density distributions later on. We implicitly assume the configurational integral to span the upper half space ($`z>0`$) only. At this point we employ a Hubbard-Stratonovitch transformation, similar to previous implementations of a field theory for charged system, which leads to

$$𝒵_N=\frac{1}{N!}\frac{𝒟\varphi }{𝒵_v}\mathrm{exp}\left\{\frac{1}{2\mathrm{}_B}d𝐫d𝐫^{}\varphi (𝐫)v^1(𝐫𝐫^{})\varphi (𝐫^{})+ıd𝐫\sigma (𝐫)\varphi (𝐫)\right\}\left[d𝐫\mathrm{e}^{h(𝐫)ıq\varphi (𝐫)}\right]^N$$

where we introduced the notation $`𝒵_v=\sqrt{detv}`$. For simplicity, we perform a Legendre transformation to the grand-canonical ensemble, $`𝒬=_N\lambda ^N𝒵_N`$, and thereby introduce the fugacity $`\lambda `$. The inverse Coulomb operator follows from Poisson’s law as $`v^1(𝐫)=^2\delta (𝐫)/4\pi `$. The grand-canonical partition function can therefore be written as

$$𝒬=\frac{𝒟\varphi }{𝒵_v}\mathrm{exp}\left\{d𝐫\left[\frac{[\varphi (𝐫)]^2}{8\pi \mathrm{}_B}ı\sigma (𝐫)\varphi (𝐫)\lambda \theta (z)\mathrm{e}^{h(𝐫)ıq\varphi (𝐫)}\right]\right\}$$
(3)
where $`\theta (z)=1`$ for $`z>0`$ and zero otherwise. Next we rescale the action, similarly to our rescaling analysis of the Hamiltonian in the beginning. All lengths are rescaled by the Gouy-Chapman length, $`𝐫=\mu \stackrel{~}{𝐫}`$, the fluctuating field is rescaled by the valency, $`\varphi =\stackrel{~}{\varphi }/q`$. We also use the explicit form $`\sigma (𝐫)=\sigma _s\delta (z)`$ for the fixed charge distribution. The result is

$$𝒬=\frac{𝒟\stackrel{~}{\varphi }}{𝒵_v}\mathrm{exp}\left\{\frac{1}{8\pi \mathrm{\Xi }}d\stackrel{~}{𝐫}\left[[\stackrel{~}{\varphi }(\stackrel{~}{𝐫})]^24ı\delta (\stackrel{~}{z})\stackrel{~}{\varphi }(\stackrel{~}{𝐫})4\mathrm{\Lambda }\theta (\stackrel{~}{z})\mathrm{e}^{h(\stackrel{~}{𝐫})ı\stackrel{~}{\varphi }(\stackrel{~}{𝐫})}\right]\right\}$$
(4)
where the rescaled fugacity $`\mathrm{\Lambda }`$ is defined by

$$\mathrm{\Lambda }=2\pi \lambda \mu ^3\mathrm{\Xi }=\frac{\lambda }{2\pi \sigma _s^2\mathrm{}_B}.$$
(5)
The expectation value of the counter-ion density, $`\rho (\stackrel{~}{𝐫})`$, follows by taking a functional derivative with respect to the generating field $`h`$, $`\rho (\stackrel{~}{𝐫})=\delta \mathrm{ln}𝒬/\delta h(\stackrel{~}{𝐫})\mu ^3`$, giving rise to

$$\frac{\rho (\stackrel{~}{𝐫})}{2\pi \mathrm{}_B\sigma _s^2}=\mathrm{\Lambda }\mathrm{e}^{ı\stackrel{~}{\varphi }(\stackrel{~}{z})}.$$
(6)
The normalization condition for the counter-ion distribution, $`\mu d\stackrel{~}{z}\rho (\stackrel{~}{z})=\sigma _s/q`$, which follows directly from the definition of the grand-canonical partition function, leads to

$$\mathrm{\Lambda }_0^{\mathrm{}}d\stackrel{~}{z}\mathrm{e}^{ı\stackrel{~}{\varphi }(\stackrel{~}{z})}=1.$$
(7)
This is an important equation since it shows that the expectation values of the fugacity term in Eq.(4) is bounded and of the order of unity per unit area. Let us first repeat the saddle-point analysis, which, because of the structure of the action in Eq.(4), should be valid for $`\mathrm{\Xi }1`$. The saddle-point equation reads

$$\frac{\mathrm{d}^2\stackrel{~}{\varphi }(\stackrel{~}{z})}{\mathrm{d}\stackrel{~}{z}^2}=2ı\mathrm{\Lambda }\mathrm{e}^{ı\stackrel{~}{\varphi }(\stackrel{~}{z})}$$
(8)
with the boundary condition $`\mathrm{d}\stackrel{~}{\varphi }(\stackrel{~}{z})/\mathrm{d}\stackrel{~}{z}=2ı`$ at $`\stackrel{~}{z}=0`$. The solution of this differential equation is

$$ı\stackrel{~}{\varphi }(\stackrel{~}{z})=2\mathrm{ln}\left(1+\mathrm{\Lambda }^{1/2}\stackrel{~}{z}\right)$$
(9)
while the boundary condition leads to $`\mathrm{\Lambda }=1`$, which shows that the saddle-point approximation is indeed valid in the limit $`\mathrm{\Xi }1`$. Combining Eqs.(6) and (9), the density distribution of counter ions is given by the well-known result

$$\frac{\rho (\stackrel{~}{𝐫})}{2\pi \mathrm{}_B\sigma _s^2}=\frac{1}{(1+\stackrel{~}{z})^2}.$$
(10)
Let us now consider the opposite limit, when the coupling constant $`\mathrm{\Xi }`$ is large. In this case, the saddle-point approximation breaks down, since the prefactor in front of the action in Eq.(4) becomes small. This has been seen in a systematic loop-wise expansion around the saddle point, where the corrections to the saddle-point solution are proportional to powers of $`\mathrm{\Xi }`$ and thus become large. However, from the field-theoretic partition function Eq.(4), it is self-evident what has to be done in this limit. Since the fugacity term is bounded, as evidenced by Eq.(7), one can expand the partition function (and also all expectation values) in powers of $`\mathrm{\Lambda }/\mathrm{\Xi }`$. Upon Legendre transformation to the canonical ensemble, this gives the standard virial expansion. The normalization condition Eq.(7) can be solved by an expansion of the fugacity as $`\mathrm{\Lambda }=\mathrm{\Lambda }_0+\mathrm{\Lambda }_1/\mathrm{\Xi }+\mathrm{}`$, which leads to an expansion of the density profile with the small parameter $`1/\mathrm{\Xi }`$. While the standard virial expansion fails for homogeneous bulk systems because of infra-red divergences, these divergences are renormalized for the present case of inhomogeneous distribution functions via the normalization condition Eq.(7) as we will now demonstrate. To leading order in this expansion, the rescaled density is

$$\frac{\rho (\stackrel{~}{𝐫})}{2\pi \mathrm{}_B\sigma _s^2}=\mathrm{\Lambda }\mathrm{exp}\left\{\frac{\mathrm{\Xi }}{2}v(0)+\frac{1}{2\pi }d𝐫v(𝐫\stackrel{~}{𝐫})\delta (z)\right\}.$$
(11)
From the normalization condition Eq.(7) we obtain

$$\mathrm{\Lambda }_0=\mathrm{exp}\left\{\frac{\mathrm{\Xi }}{2}v(0)\frac{1}{2\pi }d𝐫v(𝐫)\delta (z)\right\}$$
(12)
and thus the density distribution is to leading order given by

$$\frac{\rho (\stackrel{~}{𝐫})}{2\pi \mathrm{}_B\sigma _s^2}=\mathrm{e}^{\stackrel{~}{z}}.$$
(13)
The density profile in the strong coupling limit $`\mathrm{\Xi }\mathrm{}`$ is thus given by a simple exponential, while the profile in the weak-coupling (Poisson-Boltzmann) limit $`\mathrm{\Xi }0`$, Eq.(10), is given by a power law. Both distributions give the same density at the wall, which is a trivial consequence of the contact-value theorem since the net force on the wall is zero. Just to avoid confusion at this stage, we stress that the exponential density profile Eq.(13) has nothing to do with the Debye-Hückel approximation. It is true that an exponential density profile (though with a different density contact value and thus violating the contact-value theorem) also follows from linearizing the differential equation Eq.(8). However, the linearized solution can never be a more faithful representation of the true density profile than the full non-linear solution, and secondly, nothing in the saddle-point equation indicates that linearization should be valid for large values of $`\mathrm{\Xi }`$, because Eq.(8) does not depend on $`\mathrm{\Xi }`$. It is the saddle-point approach itself (the framework within which the Debye-Hückel approximation can be formulated) which becomes invalid when $`\mathrm{\Xi }`$ becomes large.

In Fig.2 we show counter-ion density profiles obtained using Monte Carlo simulations for various values of the coupling parameter $`\mathrm{\Xi }`$. As can be seen, the PB density profile Eq.(10) is only realized for $`\mathrm{\Xi }<1`$, while the strong-coupling profile Eq.(13) is indeed the asymptotic solution and agrees with simulation results for $`\mathrm{\Xi }>10^4`$. Experimentally, a coupling parameter $`\mathrm{\Xi }=100`$, which is already quite close to the strong-coupling limit (see Fig.2 inset), is reached with divalent ions for a surface charged density $`\sigma _s3.6nm^2`$, which is feasible with compressed charged monolayers, and with trivalent counter ions for $`\sigma _s1nm^2`$, which is a typical value. The strong-coupling limit is therefore experimentally accessible.

In our simulations we employ periodic boundary conditions in the lateral directions, while the ions are unconfined vertically. Forces and energies are calculated following the approach by Lekner and Sperb where the sum over the periodic images is transformed into a rapidly converging expansion in terms of Bessel functions. Since the number of particles in the simulation, $`N`$, is related to the linear box size (in units of the Gouy-Chapman length $`\mu `$) by $`\stackrel{~}{L}=L/\mu =\sqrt{2\pi N\mathrm{\Xi }}`$, effects due to the finite number of particles and effects due to the finite lateral box size are connected. In Fig.3 we show MC density profiles for a fixed coupling constant $`\mathrm{\Xi }=0.1`$ and for various number of particles. For all but the smallest systems with $`5`$ and $`15`$ particles finite-size effects are negligible.

In summary, we have derived within a field-theoretic framework the strong-coupling theory for counter-ion distributions, valid in the limit of large valencies and/or low temperatures. It corresponds to a standard virial expansion of the field-theoretic action, which gives meaningful results because long-range divergences, which spoil the free-energy virial expansion for bulk systems and lead to non-analytic terms, are subtracted by a fugacity renormalization. This theory thereby complements the Poisson-Boltzmann approach, which works for the oppposite limit of weak coupling. We note that an exponential density profile (albeit with a different prefactor) has been derived in the strong-coupling limit by Shklovskii using a heuristic model where ions bound to the wall are in chemical equilibrium with free ions.

Next-leading order results for the density profile within our strong-coupling theory, and comparison of these corrections with Monte-Carlo simulations, will be presented in a future publication. Also, generalization to other geometries, including the pressure between two charged walls, the presence of salt, dielectric boundaries, and hard-core interactions are straightforward within our theory and will be published soon. For the case of a cylinder with a radius $`R`$ and linear charge density $`\tau `$ we obtain that the strong-coupling theory is valid for $`\mathrm{\Xi }_c1`$ with the cylindrical coupling parameter defined as $`\mathrm{\Xi }_c=q^3\mathrm{}_B^2\tau /R`$. The Manning condensation threshold is reflected by the fact that the density distribution is not normalizable for $`\tau \mathrm{}_B<1`$. For $`\tau \mathrm{}_B>1`$ only part of the counter ions are condensed and the normalization condition Eq.(7) is modified. The result for the density profile in the limit $`\mathrm{\Xi }_c1`$ is

$$\rho (r)=\frac{(\tau \mathrm{}_Bq1)^2}{\pi \mathrm{}_Bq^2R^2}\left(\frac{r}{R}\right)^{2\tau \mathrm{}_Bq}$$
(14)
and thus exhibits purely algebraic decay, in contrast to the PB prediction but in agreement with Shklovskii’s results. Finally, we note that the strong-coupling limit not necessarily entails unrealistically small Gouy-Chapman lengths. Writing the Gouy Chapman length as $`\mu =\mathrm{}_Bq^2/\mathrm{\Xi }`$, it is clear that $`\mu `$ can be fixed at a moderate value for highly valent ions even when $`\mathrm{\Xi }1`$. In this case, however, one would have to worry about the effect of discrete surface charges.

\***

We thank Felix Csajka, Henri Orland, Christian Seidel, and Julian Shillcock for useful discussions. AGM acknowledges financial support from FCT through the grant Praxis XXI/BD/13347/97.


text:
# Screening of Hydrodynamic Interactions in Semidilute Polymer Solutions: A Computer Simulation Study

## Abstract

We study single–chain motion in semidilute solutions of polymers of length $`N=1000`$ with excluded–volume and hydrodynamic interactions by a novel algorithm. The crossover length of the transition from Zimm (short lengths and times) to Rouse dynamics (larger scales) is proportional to the static screening length. The crossover time is the corresponding Zimm time. Our data indicate Zimm behavior at large lengths but short times. There is no hydrodynamic screening until the chains feel constraints, after which they resist the flow: “Incomplete screening” occurs in the time domain.

The dynamics of polymer chains in solution  has been the subject of long–standing theoretical investigations, even for the simple case of flexible uncharged chains in good solvent. While in the dilute limit the validity of Zimm scaling predictions  is generally accepted, as confirmed by experiments  and computer simulation , the theoretical understanding becomes much more involved as soon as one considers finite concentrations . This is so due to the complicated interplay between excluded–volume interactions, hydrodynamic interactions, and entanglements. Only for the opposite limit of dense melts, where the first two interactions are fully screened, there exists a well–controlled description in terms of the Rouse or reptation model . However, the details of the crossover, the underlying mechanism of the screening of hydrodynamic interactions, and the concentration dependence of the screening length have been a subject of considerable debate.

In this Letter, we present the first computer simulation study which is able to contribute to the resolution of these questions. Experiments, such as light scattering  or non–equilibrium methods  usually focus on collective concentration fluctuations, while single–chain motions are only accessible by labeling techniques (neutron  or light  scattering). Computer simulations can in principle analyze both types of motion; however, for reasons of statistical accuracy we had to confine ourselves to single–chain motion.

We study the equilibrium fluctuations of a three–dimensional semidilute system of flexible bead–spring polymer chains with full excluded volume interactions, coupled to a hydrodynamic background to fully take into account hydrodynamic interactions, using an efficient method which we have recently developed and tested . The polymer system is simulated by Langevin stochastic dynamics, the solvent by a stochastic D3Q18 lattice Boltzmann model , and a point–particle coupling is introduced via a monomeric friction coefficient. The present work uses the same model with the same parameters but in the semidilute regime. One particular advantage, without which the study would have been unfeasible, is the fact that the lattice Boltzmann solvent does not alter the good solvent statistics of the chain conformations. We therefore first equilibrated the multi–chain system without the computationally expensive solvent, using a combination of stochastic dynamics and slithering–snake Monte Carlo moves comprising several momomers (roughly one blob, see below). This run produced a set of configurations, which were afterwards coupled to the solvent. For further details of the model, we refer the reader to Ref. .

Semidilute systems are characterized by a very low monomer concentration $`c`$, which is nevertheless large enough to induce strong overlap of the coils. The static conformations are well understood  in terms of the “blob size” $`\xi _S`$, i. e. the typical mesh size of the temporary network. On scales below $`\xi _S`$, the chains are self–avoiding walks (SAWs) characterized by the scaling law $`RaN^\nu `$, where $`a`$ is the monomer size, $`R`$ the chain extension, $`N`$ the degree of polymerization and $`\nu 0.59`$. The concentration dependence of $`\xi _S`$ results from $`c\xi _S^3(\xi _S/a)^{1/\nu }`$. On scales above $`\xi _S`$ the excluded–volume interaction is screened, and the chains are random walks (RWs, $`RaN^{1/2}`$). The overall chain is thus a RW of blobs with $`R^2\xi _S^2N/(\xi _S/a)^{1/\nu }`$. This picture implies that rather long chains are necessary in order to clearly observe both regimes; guided by the idea of having roughly $`30`$ blobs of $`30`$ monomers each available, we chose $`N=1000`$, and varied $`\xi _S`$ by studying the concentration values $`c=0.00837,0.0397,0.0734,0.134,0.231`$ for the statics, and the latter three values for the dynamics (data are always given in the Lennard–Jones unit system of Ref. ). The number of chains $`M=50`$ was kept fixed; this is large enough to ensure that even the most concentrated system does not exhibit self–overlap due to the periodic boundary conditions. Our data for the static chain conformations are in perfect agreement with the blob scaling picture , as have been those of previous extensive Monte Carlo simulations .

Dynamic scaling for a single chain which exhibits no special length scale except $`a`$ and $`R`$ implies $`\tau R^z`$, where $`\tau `$ is the conformational relaxation time, and $`z=3`$ for the Zimm model (applicable to dilute solutions without chain overlap where hydrodynamic interactions are fully developed), while $`z=4`$ for the RW Rouse model (applicable to dense melts where hydrodynamic interactions are fully screened; we do not consider reptation–like slowing down, which occurs only for longer chains, and does not play any role for our simulation data). In terms of the chain diffusion constant $`D`$ this corresponds, via $`D\tau R^2`$, to $`DR^1`$ for Zimm (chain as a Stokes sphere) and $`DN^1`$ for Rouse (monomers as independent friction centers). Furthermore, the scaling of lengths with the corresponding times implies a subdiffusive behavior for the single–monomer mean square displacement $`\mathrm{\Delta }r^2t^{2/z}`$, and a $`k^zt`$ behavior for the single–chain dynamic structure factor $`S(k,t)=N^1_{ij}\mathrm{exp}\left(i\stackrel{}{k}(\stackrel{}{r}_i(t)\stackrel{}{r}_j(0))\right)`$ in the scaling regime of intermediate length scales (between $`a`$ and $`R`$) and time scales (between $`\tau _0`$, the microscopic time for monomer relaxation, and $`\tau `$).

For a semidilute system, one expects a crossover between these cases. Indeed, our data for $`\mathrm{\Delta }r^2`$ do exhibit a crossover from a Zimm–like $`t^{2/3}`$ behavior at short times to $`t^{1/2}`$ at longer times. The behavior at short length and time scales is thus Zimm–like. The pure Zimm model  predicts that the decay rate, i. e., in the given context, the prefactor $`A`$ of the law $`\mathrm{\Delta }r^2=At^{2/3}`$ should only depend on solvent viscosity and not on concentration. We nevertheless found a weak concentration dependence of $`A`$ (roughly $`20\%`$ within the given concentration range, see below). Figure 1 thus studies $`\mathrm{\Delta }r^2/(At^{2/3})=f(t/t_c)`$, where $`t_c`$ is the concentration–dependent crossover time, which again is the Zimm relaxation time of a dynamic crossover length (or hydrodynamic screening length), $`t_c\xi _H^3`$, and $`f(x)x^{1/6}`$ for large $`x`$. We find a very good data collapse assuming that $`\xi _H\xi _S`$ or $`t_cc^{3\nu /(3\nu 1)}=c^{2.3}`$, as done in Fig. 1. The assumptions $`\xi _Hc^1`$  and $`\xi _Hc^{1/2}`$  produced significantly poorer collapses. The simulation thus confirms de Gennes’ prediction $`\xi _H\xi _S`$ , as also observed in experiments (except for weak corrections ).

$`S(k,t)`$ was studied for the most dilute system $`c=0.0734`$. We prefer scaling plots of the raw data using asymptotic exponents over fits to functional forms derived from approximate theories. From $`S(k,0)`$ we estimated the crossover wavenumber as $`k_c0.45`$ and the scaling regime as $`0.2<k<1.5`$. From $`\mathrm{\Delta }r^2`$ we read off $`t_c10^3`$; the nonuniversal regime $`t<30`$ was discarded. Figures 2 and 3 show that for short times $`tt_c`$ the decay can be described quite well by Zimm scaling, regardless of wavenumber, while for $`tt_c`$ there is a simultaneous smooth crossover to Rouse dynamics for those wavenumbers which have not yet fully decayed, i. e. for $`k<k_c`$. Note that the initial Zimm regime of these wavenumbers can be easily overlooked in the representation of the inset of Fig. 2.

The physical picture which results from this observation is thus free Zimm motion up to the crossover time (on all length scales), after which screening sets in, leading to Rouse–like motion. Hence the most important finding of our simulation is that hydrodynamic screening must necessarily be viewed as a dynamic time–dependent phenomenon. We consider this to be the logical consequence of the (correct) original treatment by de Gennes  (see below). Nevertheless, this has so far been overlooked in the literature, the main reason being that single–chain motion on length scales beyond $`\xi `$ is not accessible to standard scattering experiments  which are sensitive to collective concentration fluctuations: On scales $`k\xi 1`$ the overall solution is homogeneous, and one observes a simple diffusive decay $`\mathrm{exp}(D_{coop}k^2t)`$ with $`D_{coop}\xi ^1`$. Accordingly, single–chain motion on scales beyond $`\xi `$ was not treated explicitly in Ref. . The experiments on labeled chains  produced data which are fully consistent with our view, but were interpreted incorrectly (see below). The rest of the paper will be devoted to further discussion, and elucidation of the underlying mechanism.

Hydrodynamic interaction is the presence of long–ranged correlations in the stochastic displacements of a system of Brownian particles, caused by fast diffusive momentum transport through the surrounding fluid. They can be calculated by solving the stationary Stokes equation around a system of spheres , resulting in a complicated multipole expansion which contains many–body terms representing the multiple scattering of the flow. In the dilute limit, however, it is sufficient to just consider the leading–order pair interaction, which decays like $`r^1`$, where $`r`$ is the interparticle distance (Oseen tensor). Conversely, screened hydrodynamic interactions are described by a Yukawa–like decay $`r^1\mathrm{exp}(r/\xi _H)`$ defining the hydrodynamic screening length $`\xi _H`$. Such an interaction occurs for Darcy flow through a porous medium, where fixed frictional obstacles with friction constant $`\zeta `$ exert a force $`\zeta \stackrel{}{u}`$ on the flow with velocity $`\stackrel{}{u}`$. Denoting the obstacle concentration with $`c`$, the flow is described, on scales beyond the typical interparticle distance, by a modified Stokes equation $`\rho \stackrel{}{u}/t=\eta ^2\stackrel{}{u}\zeta c\stackrel{}{u}`$, where $`\rho `$ is the fluid density and $`\eta `$ the viscosity. This implies $`\eta \xi _H^2=\zeta c`$.

The simplest approach to hydrodynamic screening in polymer solutions just replaces the $`r^1`$ Oseen interaction by a screened Yukawa–like interaction, leading to uncorrelated displacements of monomers whose distance exceeds $`\xi _H`$. The resulting motion of the chain is Zimm–like on short length and time scales and Rouse–like on length scales beyond $`\xi _H`$ for all times . The Darcy flow thus produces the desired crossover.

Unfortunately, this picture generates as many questions as it answers. In particular, the obstacles must be the mobile polymer chains themselves, whereas strict Darcy flow requires fixed obstacles. Moreover, momentum is present and being transported infinitely far in polymeric as well as in simple fluids. This fundamental conceptual difficulty was recognized by Richter et al. . In their “incomplete screening” model they proposed that the hydrodynamic interaction should cross over to a second $`r^1`$ regime on very large scales, but with the overall viscosity as a prefactor. For the single–chain short–time behavior this model also predicts Rouse–like motion. However, this regime is now restricted to length scales $`\xi _Hk^1\xi _H\eta _{solution}/\eta _{solvent}`$. On larger scales there is an additional Zimm regime. Richter et al.  used this to interpret the mixture of Rouse– and Zimm–like signals in their scattering data. Similar arguments were used by Martin , who observed Zimm scaling on all length scales in the initial decay rate of dynamic light scattering of labeled chains.

The simple model and the more refined version by Richter et al.  are at variance with both our data and our theoretical arguments for the short-time behavior, see below. It should be noted that the “incomplete screening” model must have fundamental conceptual difficulties, since $`\eta _{solution}`$ appears in the short–time dynamics, although $`\eta _{solution}`$ is established only on time scales beyond the overall chain relaxation.

It is therefore clear that a consistent theoretical description has to study the dynamics of the coupled polymer–solvent system. The first attempt by Freed and Edwards  considered a multiple scattering series of the flow around the monomers, which is in spirit somewhat similar to the multipole expansion of Ref. . After some approximations an effective Darcy equation results, with $`\xi _Hc^1`$. Later this scheme was shown to be inadequate; in Ref.  it was proven that a system of phantom chains (which do not interact with each other, but to which the original approach  should apply as well) does not exhibit any hydrodynamic screening whatsoever. This is consistent with new results for colloidal suspensions, where possible hydrodynamic screening has been discussed recently; the result of large–scale computer simulations  was the absence of screening.

With respect to hydrodynamic screening in polymer solutions we can thus draw the following conclusions: (i) neither does the presence of higher–order terms of the multipole expansion  at finite concentrations of scattering centers lead to screening; (ii) nor can such terms be of any importance in the semidilute limit, where one can reach arbitrarily small monomer concentrations $`c`$, while still keeping the chains at strong overlap; (iii) as this “colloidal” mechanism of screening does not apply, the underlying physics must rather be polymer–specific; and (iv) the mechanism must lead to a time–delayed screening.

Concerning the short–time behavior, we note that the semidilute system is governed by a Kirkwood diffusion equation , with a pure Oseen–type $`r^1`$ diffusion tensor, which describes the short–time diffusive behavior, and a force term due to connectivity, excluded volume, and entanglements. Within this formalism, it can be shown rigorously  that the initial decay rates of correlation functions are only governed by the diffusion tensor and the statistics of the chain conformations. In particular, considering the initial decay rate of the single–chain structure factor , one obtains the same formula as for an isolated chain in solvent — the effect of the other chains is merely the modification of the conformations. Zimm chains, however, have always $`z=3`$ independently of chain statistics; the initial decay rate is given by $`\mathrm{\Gamma }(k)(k_BT/\eta )k^3`$, while the fractal dimension only enters the prefactor . For systems in the semidilute limit we thus conclude, in accordance with our simulation results and the experimental data by Martin , that for short times the single–chain dynamics is Zimm–like, independently of length scales.

In his pioneering 1976 paper  de Gennes noticed that the decisive mechanism for screening is the connectivity and the strong coupling to the temporary matrix (expressed in terms of an elastic gel, which is physically more appropriate than a rigid porous medium): After the time needed for a blob to move its own size, which is the blob’s Zimm time $`t_c`$, the chain will, on average, feel the constraints by the temporary matrix. From then on it is unable to follow the flow, but rather lags behind, and starts to produce Darcy–type frictional resistance. As the blob can be envisioned as a Stokes sphere with radius $`\xi `$ and friction coefficient $`\eta \xi `$, the obstacles which produce the Darcy flow are not the monomers but rather the blobs. Hence, the hydrodynamic screening length is given by $`\eta \xi _H^2\eta \xi _Sc_{blob}\eta \xi _S\xi _S^3`$, i. e. the hydrodynamic screening length is, apart from prefactors, identical to the static screening length, $`\xi _H\xi _S`$. This argument  makes the picture fully self–consistent. On length scales beyond $`\xi `$, and time scales beyond $`t_c`$, the semidilute solution is just a Rouse melt of blobs, while the conformations within the blobs are already fully relaxed. In the Rouse regime, momentum transport is no longer described by a simple Navier–Stokes equation. It rather occurs mainly along the chain backbones, due to the connectivity forces. This results in a very efficient randomization of a locally applied “kick”.

It should be noted that our simulated system deviates somewhat from that ideal scenario. The most dilute system has a density of $`9\%`$ of a typical dense melt, and thus one must expect that higher–order terms in the multipole expansion  do play a role. We believe that these are the main source of the $`c`$ dependence of the prefactor of the initial $`t^{2/3}`$ law in $`\mathrm{\Delta }r^2`$. We expect finite size effects  with respect to the linear box size $`L`$ to be rather small, since for our data $`(kL)^1<0.06`$.

To summarize, we have presented the first computer simulation study which was able to study the dynamic crossover from Zimm to Rouse behavior in semidilute polymer solutions. This was made possible by a novel algorithm, whose essential feature is the replacement of the solvent by a Navier–Stokes background which is coupled dissipatively to the monomers. Our results are fully consistent with de Gennes’ scaling picture , and emphasize the fact that hydrodynamic screening is a dynamic effect which becomes relevant only after the crossover time. Incomplete screening thus indeed occurs, however not on large length scales , but on short time scales. Any theoretical description which builds upon a screened hydrodynamic interaction depending only on distance, but disregards the time dependence, cannot describe the phenomena correctly.

We gratefully acknowledge fruitful discussions with K. Kremer, and thank the Max Planck Society for generous allocation of Cray T3E computer time at RZG.


text:
# References

Recent investigations have brought evidence of a deep connection between two-dimensional (2D) gravitating systems and conformal mechanics . The most natural context to test this conjecture is the Anti-de Sitter/conformal field theory (AdS/CFT) correspondence  in two spacetime dimensions . In fact, for $`D=2`$ this correspondence essentially states that gravity on AdS2 should be described by a conformally invariant quantum mechanics.

Most of the progress in this direction has been achieved in the context of 2D dilaton gravity, mainly because of the simplicity of the model . For these models the conformal symmetry (generated by a Virasoro algebra) has a natural interpretation in terms of the asymptotic symmetries of the gravitational system. Moreover, the central charge of the algebra, whose value is crucial for calculating the statistical entropy of 2D black holes, can be calculated using the deformation algebra of the boundary of AdS2 .

Despite the simplicity of the model, the various attempts to identify the conformal quantum mechanics that should be dual to gravity on AdS2, and to calculate the entropy of 2D black holes by counting states of the CFT, met only partial success . The conformal mechanics involved could not be identified and a mismatch of a $`\sqrt{2}`$ factor between thermodynamical and statistical entropy was found.

The puzzle has become even more intricate in view of the results of Ref. . There it was shown that in $`D=2`$ the AdS/CFT correspondence has a realization in terms of a two-dimensional CFT, which is essentially a theory of open strings with Dirichlet boundary conditions. Moreover, it was shown that the degeneracy of states of the 2D CFT explains correctly the thermodynamical entropy of 2D black holes.

In this letter we show that the boundary dynamics induced by AdS2 dilaton gravity can be described by a De Alfaro-Fubini-Furlan (DFF) model  coupled to an external source with conformal dimension two. The external source encodes the information about the gauge symmetries of the 2D gravity system. Alternatively, the dynamics of the boundary fields admits an equivalent description in terms of a mechanical system with anholonomic constraints. In both cases, the mechanical system is invariant under the action of the full one-dimensional conformal group generated by a Virasoro algebra, which also appears as asymptotic symmetry algebra of AdS2. We compute the central charge of this algebra and find perfect agreement between statistical and thermodynamical entropy of AdS2 black holes.

Our starting point is the Jackiw-Teitelboim (JT) model , with action

$$I=\frac{1}{2}d^2x\sqrt{g}\eta [R+2\lambda ^2],$$
(1)
where $`\eta `$ represents the dilaton. Two-dimensional anti-de Sitter space, or more generally black holes in AdS2, are solutions of this model ,

$$ds^2=(\lambda ^2r^2a^2)dt^2+(\lambda ^2r^2a^2)^1dx^2,\eta =\eta _0\lambda x,$$
(2)
where $`a`$ is given in terms of the black hole mass $`M`$, $`a^2=2M/\lambda `$. The thermodynamical black hole entropy is given by

$$S=4\pi \sqrt{\frac{\eta _0M}{2\lambda }}.$$
(3)
The asymptotic symmetries of AdS2 are defined as the transformations which leave the asymptotic form of the metric invariant, i. e. they preserve the large $`r`$ behavior

$`g_{tt}`$ $`=`$ $`\lambda ^2r^2+\gamma _{tt}(t)+𝒪\left({\displaystyle \frac{1}{r^2}}\right),`$
$`g_{tr}`$ $`=`$ $`{\displaystyle \frac{\gamma _{tr}(t)}{\lambda ^3r^3}}+𝒪\left({\displaystyle \frac{1}{r^5}}\right),`$
$`g_{rr}`$ $`=`$ $`{\displaystyle \frac{1}{\lambda ^2r^2}}+{\displaystyle \frac{\gamma _{rr}(t)}{\lambda ^4r^4}}+𝒪\left({\displaystyle \frac{1}{r^6}}\right),`$ (4)
where the fields $`\gamma _{\mu \nu }`$ parametrize the first sub-leading terms in the expansion and can be interpreted as deformations of the boundary.

The asymptotic form (4) is preserved by infinitesimal diffeomorphisms $`\chi ^\mu (x,t)`$ of the form

$`\chi ^t`$ $`=`$ $`ϵ(t)+{\displaystyle \frac{\ddot{ϵ}(t)}{2\lambda ^4r^2}}+{\displaystyle \frac{\alpha ^t(t)}{r^4}}+𝒪\left({\displaystyle \frac{1}{r^5}}\right),`$
$`\chi ^r`$ $`=`$ $`r\dot{ϵ}(t)+{\displaystyle \frac{\alpha ^r(t)}{r}}+𝒪\left({\displaystyle \frac{1}{r^2}}\right),`$ (5)
where $`ϵ(t)`$ and $`\alpha ^\mu (t)`$ are arbitrary functions, the $`\alpha ^\mu `$ describing pure gauge diffeomorphisms. In Ref.  it was shown that the symmetries (5) generate a Virasoro algebra.

The asymptotic behavior of the scalar field $`\eta `$, compatible with the transformations (5), must take the form

$$\eta =\eta _0\left(\lambda \rho (t)r+\frac{\gamma _\eta (t)}{2\lambda r}\right)+𝒪\left(\frac{1}{r^3}\right),$$
(6)
where $`\rho `$ and $`\gamma _\eta `$ play a role analogous to that of the $`\gamma _{\mu \nu }`$. Introducing the new fields, invariant under the pure gauge diffeomorphisms parametrized by $`\alpha ^\mu `$,

$`\beta ={\displaystyle \frac{1}{2}}\rho \gamma _{rr}+\gamma _\eta ,`$
$`\gamma =\gamma _{tt}{\displaystyle \frac{1}{2}}\gamma _{rr},`$ (7)
the equations of motion following from the action (1) yield in the limit $`r\mathrm{}`$111At first sight, it seems not necessary to require (9), since it comes from the leading term in the stress tensor component $`T_{rt}`$, which is of order $`1/x^2`$ . However, $`T_{\mu \nu }`$, transforming as a tensor, is clearly not invariant under coordinate transformations. In fact, in the light-cone coordinates used below, (9) originates from an order $`1`$ term, so requiring (9) is really necessary for consistency.

$$\lambda ^2\ddot{\rho }=\rho \gamma \beta ,$$
(8)
$$\dot{\rho }\gamma +\dot{\beta }=0.$$
(9)
Eqs. (8) and (9) determine a mechanical system with anholonomic constraint, since the one-form

$$\omega \gamma d\rho +d\beta $$
(10)
is not exact. The Lagrange equations of the first kind for the fields $`\phi _i=\{\rho ,\beta ,\gamma \}`$ read

$$F_im_i\ddot{\phi }_i+\mathrm{\Lambda }\omega _i=0,$$
(11)
where $`F_i`$ is the force that can be derived from a potential $`U`$, $`F_i=_iU`$, $`m_i`$ denote the masses of the fields, $`\mathrm{\Lambda }`$ is a Lagrange multiplier, and the $`\omega _i`$ are the components of the one-form $`\omega `$. If we choose

$$m_\rho =\lambda ^1,m_\beta =m_\gamma =0,$$
(12)
and

$$U=\lambda \beta \rho ,$$
(13)
the Lagrange equations (11) yield (8), together with the Lagrange multiplier $`\mathrm{\Lambda }=\lambda \rho `$. Before we proceed, we note that from Eqs. (8) and (9), one gets the conservation law

$$T+U=\frac{1}{2}\lambda ^1\dot{\rho }^2+\lambda \beta \rho =\mathrm{const}.$$
(14)
Notice that $`T+U`$ is essentially the mass of the black holes considered in .
The boundary fields $`\phi _i`$ span a representation of the full infinite dimensional group generated by the Killing vectors (5). In fact, under the asymptotic symmetries (5), they transform as

$`\delta \rho `$ $`=`$ $`ϵ\dot{\rho }\dot{ϵ}\rho ,`$
$`\delta \beta `$ $`=`$ $`ϵ\dot{\beta }+\dot{ϵ}\beta +{\displaystyle \frac{\ddot{ϵ}\dot{\rho }}{\lambda ^2}},`$
$`\delta \gamma `$ $`=`$ $`ϵ\dot{\gamma }+2\dot{ϵ}\gamma {\displaystyle \frac{\stackrel{\mathrm{}}{ϵ}}{\lambda ^2}}.`$ (15)
The above transformations are easily recognized as (anomalous) transformation laws for conformal fields of weights $`1,1,2`$ respectively. We are interested in the transformation laws of the equation of motion (8), the constraint (9), and the conserved charge (14). Using (15), we get

$`\delta [\gamma \dot{\rho }+\dot{\beta }]=`$ (16)
$`ϵ{\displaystyle \frac{d}{dt}}[\gamma \dot{\rho }+\dot{\beta }]+2\dot{ϵ}[\gamma \dot{\rho }+\dot{\beta }]+\ddot{ϵ}[\gamma \rho +\beta +{\displaystyle \frac{\ddot{\rho }}{\lambda ^2}}],`$
$`\delta [\gamma \rho +\beta +{\displaystyle \frac{\ddot{\rho }}{\lambda ^2}}]=`$ (17)
$`ϵ{\displaystyle \frac{d}{dt}}[\gamma \rho +\beta +{\displaystyle \frac{\ddot{\rho }}{\lambda ^2}}]+\dot{ϵ}[\gamma \rho +\beta +{\displaystyle \frac{\ddot{\rho }}{\lambda ^2}}],`$
$$\delta [T+U]=ϵ\frac{d}{dt}[T+U].$$
(18)
We see that the constraint transforms like a conformal field of weight 2 with anomaly term. The conformal weights of the equation of motion and the conserved charge are 1 and 0 respectively, and anomalies are absent. The above equations imply that on-shell the constraint and the equation of motion are invariant under the transformations (15).

Alternatively, we may describe the dynamical system (8) in terms of the DFF model  of conformal mechanics, coupled to an external source. To this aim, we start from the conservation law (14), i. e. $`T+U=c`$. Introducing the new field $`q=\sqrt{\rho /\lambda }`$, which has conformal dimension $`1/2`$, and eliminating $`\beta `$ from (14) by means of (8), we arrive at the equation

$$\ddot{q}\frac{g}{q^3}=\frac{\lambda ^2}{2}\gamma q,$$
(19)
with $`g=c/(2\lambda )`$, whereas from Eq. (14) follows

$$\frac{\dot{q}^2}{2}+\frac{g}{2q^2}=\frac{\lambda ^2}{4}\beta .$$
(20)
One can easily check the equivalence of the system (8), (9) with (19), (20). Moreover, (19) is easily recognized as the equation of motion for the DFF model coupled to an external source $`\gamma `$.

At this point it is straightforward to write down an action for the dynamical system (19). It is given by

$$I=𝑑t\left[\frac{1}{2}\dot{q}^2\frac{g}{2q^2}+\frac{1}{4}\lambda ^2\gamma q^2\right].$$
(21)
(21) resembles very much the IR-regularized action proposed by DFF , the only difference consists in the fact that the external source $`\gamma `$, which couples to the field $`q`$, is not constant, but represents an operator of conformal dimension two. Note that in the calculation of $`\delta I`$, $`\gamma `$, being an external source, is not varied. One can easily show that the action is (up to a total derivative) invariant under the conformal transformations (15).

The dynamical system described by Eqns. (8), (9) (or equivalently by (19)), defines a one-dimensional conformal field theory (CFT1). The invariance group of the model coincides with the group of asymptotic symmetries of AdS2, and can be realized as the diff1 group describing time-reparametrizations $`\delta t=ϵ(t)`$. In analogy with 2D CFT, one would like to identify the stress-energy tensor $`T_{tt}`$ associated with the CFT1. This analogy suggests that $`T_{tt}`$ is proportional to the constraints (9),

$$T_{tt}=\lambda (\dot{\rho }\gamma +\dot{\beta }).$$
(22)
(The constant of proportionality has been chosen to ensure that $`T_{tt}`$ has the dimensions of a mass squared). In fact, from Eqns. (8), (9) and (16) it is evident that $`T_{tt}`$ plays the same role as the holomorphic $`T_{++}`$ (and antiholomorphic $`T_{}`$) stress-energy tensor play in 2D CFT. A more compelling argument leading to Eq. (22) relies on the identification of $`T_{tt}`$ as the boundary value for the stress-energy tensor of a 2D CFT.

To do this, we choose the conformal gauge

$$ds^2=e^{2\omega }dx^+dx^{}.$$
(23)
Then the action (1) takes the Liouville-like form

$$I=d^2x\left(_+\omega _{}\eta +_{}\omega _+\eta \frac{1}{2}\lambda ^2\eta e^{2\omega }\right).$$
(24)
The action must be complemented by the constraints (the equations of motion for the missing components of the metric)

$$T_{\pm \pm }=_\pm ^2\eta 2_\pm \eta _\pm \omega =0,$$
(25)
where $`T_{\pm \pm }`$ denote the components of the stress-energy tensor.

For $`\eta \mathrm{}`$ the potential term in the action (24) goes to zero and the model becomes an exact 2D CFT . In fact, defining the new fields $`X,Y`$,

$$\omega =XY,\eta =X+Y,$$
(26)
the action and the constraints become

$`I`$ $`=`$ $`2{\displaystyle d^2x\left(_+X_{}X_+Y_{}Y\right)},`$ (27)
$`T_{\pm \pm }`$ $`=`$ $`_\pm ^2X+_\pm ^2Y2_\pm X_\pm X+2_\pm Y_\pm Y=0.`$ (28)
Taking $`\eta \mathrm{}`$ we reach the boundary of AdS2, which in light-cone coordinates is located at $`x^+=x^{}`$. One can show that $`T_{++}|_{boundary}=T_{}|_{boundary}=T_{tt}`$, with $`T_{tt}`$ given by Eq. (22).

Let us now show that $`T_{tt}`$ generates the diff1 group. Introducing the charges

$$\widehat{J}=ϵT_{tt},$$
(29)
and using the transformation law (16), one has

$$\delta _ϵT_{tt}=[\widehat{J},T_{tt}]=ϵ\dot{T_{tt}}+2\dot{ϵ}T_{tt}.$$
(30)
Expanding in Fourier modes,

$$T_{tt}=L_me^{im\lambda t},ϵ(t)=a_me^{im\lambda t},$$
(31)
and using Eq. (30), one finds that $`L_m`$ generates a Virasoro algebra,

$$[L_m,L_m]=(mn)L_{m+n}.$$
(32)
Eq. (22) does not give the most general form of the CFT1 stress-energy tensor. We have the freedom to add to it a constant term (which we choose proportional to the black hole mass $`M`$) and an improvement term,

$$T_{tt}|_{impr}=T_{tt}+\lambda M+b\ddot{\rho },$$
(33)
where $`b`$ is an arbitrary constant. The improvement term is a total derivative and does not affect the status of $`T_{tt}`$ as generator of the diff1 group. Though the constant $`b`$ is at the CFT level undetermined, it can be computed using the underlying gravitational dynamics. Comparing (33) with the expression for $`T_{tt}`$ found in Ref. , one gets $`b=2\eta _0`$.

We can now calculate the central charge $`C`$ associated with the Virasoro algebra (32). A naive computation will give the value $`C=24\eta _0`$ found in Ref. . This can be shown explicitly using the general transformation law of a CFT1 stress-energy tensor,

$$\delta T_{tt}=ϵ\dot{T}_{tt}+2\dot{ϵ}T_{tt}+\frac{C}{12}\stackrel{\mathrm{}}{ϵ}.$$
(34)
Using Eqns. (15) and (16), one can easily show that $`T_{tt}`$ given by Eq. (33) follows the transformation law (34), with $`C=24\eta _0`$. This value of the central charge, once inserted in the Cardy formula, $`S=2\pi \sqrt{\frac{CL_0}{6}}`$ , produces a black hole entropy, which differs by a $`\sqrt{2}`$ factor from the thermodynamical value (3).

However in this computation (and in those of Refs. ,  as well), one only considers the contribution of the $`r\mathrm{}`$ boundary of AdS2. The black hole solution (2) has also an inner boundary  (located at $`r=0`$ for the ground state or at the horizon for the generic black hole), which can give a contribution to the central charge. That this inner boundary can be crucial for understanding the black hole entropy has been shown in .

There is a simple way to compute this contribution. We need to change the coordinates, from the Schwarzschild $`(r,t)`$ frame used in Eq. (2) to the conformal frame, where the vacuum and black hole solutions have, respectively, the form ,

$$ds^2=\frac{1}{\lambda ^2x^2}(dt^2+dx^2),$$
(35)
$$ds^2=\frac{a^2}{\mathrm{sinh}^2(a\lambda \sigma )}(d\tau ^2+d\sigma ^2).$$
(36)
The key point is that in the conformal frame the inner boundary (and the horizon) is pushed to $`x=\mathrm{}`$, whereas the the timelike $`r=\mathrm{}`$ boundary of AdS2 is now located at $`x=0`$. The information about the existence of the inner boundary is now encoded in the coordinate transformation

$$t=\frac{1}{a\lambda }e^{a\lambda \tau }\mathrm{cosh}(a\lambda \sigma ),x=\frac{1}{a\lambda }e^{a\lambda \tau }\mathrm{sinh}(a\lambda \sigma ),$$
(37)
that maps the vacuum (35) into the black hole solution (36) .

Evaluating these transformations on the $`x=\sigma =0`$ boundary, where our one- dimensional conformal field theory lives, we find

$$t=\frac{1}{a\lambda }e^{a\lambda \tau }.$$
(38)
Thus the vacuum and black hole solutions correspond to different time variables on the boundary. Moreover, for $`\mathrm{}<\tau <\mathrm{}`$, $`0<t<\mathrm{}`$, part of the ”history” seen by the vacuum observer cannot be seen by the black hole observer. It follows that there will be a term in the entropy describing the entanglement of states, which has the form of a contribution $`C_{ent}`$ to the central charge. This contribution can be calculated using a method similar to that employed in Ref.  ($`C_{ent}`$ is interpreted as a Casimir energy).

The transformation law of the stress-energy tensor under a general change of coordinates $`t=t(\tau )`$ is given by the transformation (34) in its finite form,

$$T_{\tau \tau }=\left(\frac{dt}{d\tau }\right)^2T_{tt}\frac{C_{ent}}{12}\left(\frac{dt}{d\tau }\right)^2\{\tau ,t\},$$
(39)
where $`\{\tau ,t\}`$ is the Schwarzian derivative. Applied to the transformation (38), Eq. (39) gives

$$T_{\tau \tau }=(a\lambda t)^2T_{tt}\frac{C_{ent}}{24}a^2\lambda ^2.$$
(40)
From Eq. (33) it follows that by fixing the diffeomorphisms invariance, on-shell we can always have $`T_{\tau \tau }=\lambda M`$. In fact, on-shell the term proportional to the constraints is zero and we can always choose $`\rho =const`$. Because $`T_{tt}`$ in Eq. (40) refers to the vacuum, we have $`T_{tt}=0`$, and (40) becomes

$$\lambda M=\frac{C_{ent}}{24}a^2\lambda ^2.$$
(41)
The coordinate transformation (38) maps the ground state into the black hole with mass $`M=a^2\eta _0\lambda /2`$ , which inserted into Eq. (41) gives

$$C_{ent}=12\eta _0.$$
(42)
Notice that the entanglement contribution is negative, yielding a total central charge

$$C_{tot}=C+C_{ent}=12\eta _0.$$
(43)
Using this value of the central charge in Cardy’s formula, one finds perfect agreement with the thermodynamical entropy of the 2D black hole (3).

Let us now discuss the physical interpretation of the dynamical system (19) and its relationship with 2D dilaton gravity. The equation of motion (19) describes a mechanical system coupled to an external source. Alternatively, one can think of $`\gamma `$ as a time-dependent coupling constant, appearing in the harmonic oscillator potential.

Because $`\gamma `$ is arbitrary (the only constraint on it is the transformation law (15)), the dynamical system is essentially non-deterministic. Moreover, because the external source is time-dependent, the energy is not conserved and we have a time-dependent Hamiltonian, whose evolution is not completely fixed by the dynamics of the system.

Strictly speaking we have to deal with an ensemble of Hamiltonians. One has thus a strong analogy with disordered systems in statistical mechanics. The main difference is that, whereas in the case of disordered systems we have a probability distribution for the couplings, in our case they are arbitrary functions, for which we only give the transformation law under the conformal group. Moreover, in our case $`\gamma `$ is a smooth function of $`t`$.

From the point of view of the 2D gravitational theory the meaning of the source $`\gamma `$ is clear. The fields $`\beta `$ and $`\gamma `$ describe deformations of the boundary of AdS2, generated by 2D bulk diffeomorphisms, whereas the field $`\rho `$ describes deformations of the dilaton. Thus, the function $`\gamma `$ encodes the information about the gauge symmetry of the 2D gravity theory. The non-deterministic nature of the dynamical system (19) is a consequence of the gauge freedom of the gravitational dynamics. This indicates an interesting relationship between gauge symmetries and non-deterministic dynamical systems.

In Ref. , it was pointed out that the non constant value of the dilaton ($`\rho 0`$ in terms of boundary fields) breaks the $`SL(2,R)`$ isometry group of AdS2, and that the origin of the central charge (and hence of the black hole entropy) can be traced back to this breaking. Different values of $`\rho `$ represent different vacua of the 2D gravity theory which break $`SL(2,R)`$. Indeed the conformal transformations (15), (which are the boundary counterpart of the 2D diffeomorphisms) map all these vacua one into the other. Moreover, because the energy (14) is invariant (on-shell) under conformal transformations, all these vacua are degenerate in energy.

Now, the crucial point is that there is a one-to-one correspondence between these vacua and the solutions of the dynamical system (19). Note that also the (IR-regularized) DFF model (21) breaks the $`SL(2,R)`$ symmetry of the original model if we take the source $`\gamma `$ constant. Introducing a time-dependent external source, transforming as conformal field of weight 2, we reinstate the full conformal symmetry. From the point of view of 2D gravity we are considering different $`\rho `$-dependent vacua.

The above considerations indicate a natural way to explain statistically the entropy of 2D black holes. This entropy can be interpreted in terms of the degeneracy of the $`\rho `$-vacua. The energy (14) is invariant under conformal transformations, so that we can calculate the entropy by counting the independent excitations in the configuration space of vacua.

From the point of view of the dynamical system this degeneracy is encoded in the external source $`\gamma `$. The quadratic mass-temperature and mass-entropy dependence , which is typical of a 2D CFT and which, in principle, could be used to rule out the duality of AdS2 gravity with a conformal quantum mechanics, is presumably related to the fact that the conformally invariant mechanical system in question is not a usual mechanical system but a DFF model coupled to an external source.