text:
# Deep inside the core of Abell 1795: the Chandra view

## 1 INTRODUCTION

Abell 1795 (Abell, Corwin & Olowin 1989) is a nearby rich cD galaxy cluster well studied at optical, radio and X-ray wavelengths. In X-rays, it shows a relaxed structure (Buote & Tsai 1996) and a surface brightness strongly peaked in the centre, where cooler gas is present (Fabian 1994, Briel & Henry 1996, Allen & Fabian 1997, Ikebe et al. 1999, Allen et al. 2001, Tamura et al. 2001). Strong optical emission lines around the cD galaxy (Cowie et al. 1983), an excess in the blue optical light (Johnstone, Fabian & Nulsen 1987; Cardiel, Gorgas & Aragon-Salamanca 1997) and emission in ultra-violet (Mittaz et al. 2001) due to the formation of massive stars, perhaps triggered by the central radio source 4C 26.42 (McNamara et al. 1996), support a scenario in which the central peak of cooler gas is due to a massive cooling flow (Edge et al. 1992, Fabian 1994).

In a previous paper (Fabian et al. 2001b), we have presented a spatial analysis in three X-ray colors of the Chandra observation of the core of A1795. In the present work, we analyze the spectral characteristics of the intracluster medium, studying the metal enrichment, deprojecting the physical properties and determining the form of the gravitational potential within the cluster core.

Hereafter, we assume the values of the cosmological parameters to be $`H_0=50`$ km s-1 Mpc-1, $`\mathrm{\Omega }_\mathrm{m}=1`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0`$. At the nominal redshift of A1795 ($`z=0.0632`$, Girardi et al. 1998), an angular scale of 1 arcsec ($`2`$ ACIS pixels) corresponds to a physical length of 1.65 kpc. For the currently popular cosmology of $`(H_0,\mathrm{\Omega }_\mathrm{m},\mathrm{\Omega }_\mathrm{\Lambda })=(70,0.3,0.7)`$, distances change by a factor of 0.74, luminosity by 0.742 = 0.54 and the age of the Universe at the cluster redshift by 1.06.

All the errors quoted are $`1\sigma `$ (68.3 per cent level of confidence) unless stated otherwise.

## 2 Chandra datasets

In Table 1, we present a summary of the Chandra (Weisskopf et al. 2000) observations which were both made with the back-illuminated CCD S3 of the Advanced Camera for Imaging and Spectroscopy (ACIS) instrument. We have used the CIAO software (version 2.1) to clean further the revision-1 level-2 events files provided by the standard pipeline processing that filters for ASCA grade classifications of 0,2,3,4,6 and defines Good Time Intervals for a nominal total exposure given in parenthesis in Table 1. The FP–120°dataset has been analyzed with the August 2001 release of of ACIS response products (Fits Encoded Function and gain files; see at http://asc.harvard.edu/cal/Links/Acis/acis/ Cal\_prods/matrix/matrix.html), whereas the FP–110°events file is calibrated according to the previous release of December 1999. Therefore the original FP–120°dataset has been reprocessed with acis\_process\_events routine in CIAO 2.1.3 with the gain files available in CALDB 2.7. We have also corrected for time intervals with (i) bad attitude solutions, and (ii) high background with flares in the lightcurve with relative value larger than 20 per cent the mean count rate. This process reduces the effective exposure by about 5 per cent in the two datasets.

## 3 Imaging analysis

The optical position of the cD (which hosts the radio source 4C+26.42) is (RA, Dec) 13h 48m 52$`\stackrel{s}{.}`$43, 26°35′34$`\stackrel{}{.}`$0 (J2000). This is within the $`5\times 5`$ arcsec2 bright envelope to the north of the X-ray emitting wake, also known as an H$`\alpha `$ emitter, discussed in our previous work (Fabian et al. 2001b). The X-ray centre for our analysis is chosen at coordinates (RA, Dec) 13h 48m 52$`\stackrel{s}{.}`$54, 26°35′25$`\stackrel{}{.}`$3 (J2000). This is consistent with the centroid of the X-ray emission determined in three different X-ray colour images. The behaviour of the centroid and the chosen centre are plotted in Fig. 1.

The choice of the centre is also determined to make circular annuli a good representation of the cluster surface brightness as shown in Fig. 2. This is relevant in the deprojection analysis applied in the present work that assumes a spherical symmetric emission projected on the sky in circular regions.

The azimuthally averaged brightness profile is shown in Fig. 3. Two main breaks are evident at about 25 and 300 kpc, respectively, when the profile is compared with a smoothed distribution provided, e.g., from a $`\beta `$model (e.g. Cavaliere & Fusco-Femiano 1976) fitted over the radial range 0–400 kpc (with best-fit parameters $`r_\mathrm{c}=49\pm 1`$ kpc, $`\beta =0.510\pm 0.002`$, $`\chi ^2`$ = 971, 243 degrees of freedom). The presence of breaks in the brightness distribution is now routinely observed in Chandra data of galaxy clusters (e.g. Allen, Ettori, Fabian 2001, Markevitch et al. 2001, Vikhlinin et al. 2001, Mazzotta et al. 2001) and is mainly associated with discontinuities in the gas density profile probably due either to the aggregation of subclumps or to the “sloshing” of the gas decoupled from dark matter into the changing gravitational potential (as suggested for the break in the surface brightness about 70″southward the central galaxy in A1795 from Markevitch, Vikhlinin & Mazzotta 2001).

It is worth noting that a $`\beta `$model is generally not a good representation of the brightness profile of cooling flow clusters, like A1795, over large radial ranges due to the peaked central emission (but see improvements in modelling this kind of profiles obtained from Rosat PSPC observations in Mohr et al. 1999 and Ettori 2000). When the central 200 kpc radius region is not included in the fit, Ettori & Fabian (1999) measure from a Rosat PSPC observation that extends up to 1.5 Mpc a core radius of 250 kpc and a $`\beta `$ of 0.75, consistent with results from fits that include a model for the central excess (Mohr et al. 1999, Ettori 2000). Due to the strong positive correlation between core radius and $`\beta `$ parameter, the evidence of a smaller core in the brightness profile as observed from Chandra implies a lower value of $`\beta `$ as estimated. Moreover, it is interesting to note the dependence of these parameters from the radial range considered. For example, we find that, from a statistical point of view, the best modelling with a $`\beta `$model in the central 400 kpc of A1795 is obtained over the range 0–60 kpc ($`\chi ^2`$/dof = 45/34, probability P=0.09), with best-fit parameters $`r_\mathrm{c}=23\pm 2`$ kpc and $`\beta =0.353\pm 0.011`$ that show a significant lower value of $`r_\mathrm{c}`$ and a corresponding lower value of $`\beta `$. We remind, however, that the $`\beta `$model is a convenient fitting function that assumes both an isothermal gas and the hydrostatic equilibrium between the galaxy population, the intracluster plasma and the underlying gravitational potential. In the cooling flow regions, the assumption of isothermality is definitely not correct and also the hydrostatic condition might be violated.

In Fig. 3, we show a map of the differences, $`(IM)/\sqrt{M}`$, between the unsmoothed exposure corrected image, $`I`$, and the two-dimensional distribution of a smoothed model, $`M`$, that we choose to be a $`\beta `$model with parameters constrained over the radial range 0–400 kpc. Note that positive excesses are present corresponding to (i) the central 20 arcsec–radius structure (i.e. cD galaxy plus filament) and (ii) extended emission 50 arcsec from the centre towards North-East. Negative values are predominant in the South-East direction.

The evidence for regions with an excess/deficit in emission along the North-South direction, but not aligned, becomes clear when the brightness profiles extracted in these sectors are compared to the emission from westward of the centre (Fig. 4). Starting at about 50 kpc, the excess in the North-East and the deficit in the South-East become evident with deviations respective to the West region of about $`10\sigma `$. This suggests that the core is in a non-relaxed status with indication of merging activity. In Sect. 5.4, we discuss more on the dynamical status of the cluster core.

### 3.1 X-ray colour profile

To study the spatial distribution of the photons with respect to their energy, we extract images in two X-ray ‘colours’ (0.5–1.5 keV, 1.5–7 keV), correct them by the respective exposure maps and subtract the background estimated from blank fields (see further comments in Sect. 5.2). The surface brightness profiles are then obtained in bins of 20 physical pixels (equal to 9.8 arcsec) and the color ratios determined (see Fig. 5)

The colour ratio rapidly decreases and flattens moving outwards, indicating a hardening in the photon counts. We have quantified this in a non-parametric way, calculating a weighted average (and relative error) of the colour ratio for all the points above a given radius. Then, a distribution of the number of points enclosed within a given colour ratio $`\pm 3\sigma `$ is obtained. The values corresponding to (16, 50, 84) per cent of the cumulative function of this distribution is plotted in Fig. 5. The break between the “most–probable” value and the rapid increase of the profile in the centre is located at $`r=72_{32}^{+40}`$ arcsec ($`119_{53}^{+66}`$ kpc).

This excess of the relative amount of soft ($`<`$1.5 keV) photons with respect to the harder counts population traces the region with a lower gas temperature. Whether this gas is cooling or just following the behaviour of the gravitational potential is discussed in the following sections.

## 4 Spectral analysis

The spectra were extracted from the processed events-2 files, selecting in PI space for energies between 0.5 and 7 keV to minimize the effects both of the uncertainties in the calibration of the observed soft X-rays and of the cosmic rays that dominate the background components at energies larger than 5 keV (see at the page http://asc.harvard.edu/ cal/Links/ Acis/acis/ Cal\_prods/bkgrnd/ current/index.html maintained by M. Markevitch). The background spectrum has been extracted from the same region of the CCD in exam from blank field exposures (namely aciss\_B\_s3\_bg\_evt\_060600.fits and aciss\_C\_s3\_bg\_evt\_191000.fits for FP–110°and FP–120°, respectively) adapted to our observations through the routine make\_acisbg provided to the Chandra community by M. Markevitch. As presented in Table 2, the total counts rate estimated in each ring in the adopted 0.5–7 keV energy range is contaminated from the background by a percentage between less than 1 per cent (in the inner annulus) and 4 per cent (in the outer ring), making this correction lower than statistical uncertainties.

Redistribution Matrix Files (RMF) and Auxiliary Response Files (ARF) were made by averaging over a $`32\times 32`$ pixel-grid of calibration files covering the back–illuminated chip 7, using the CIAO tools mkrmf and mkarf and using weighting factors equal to the number of counts in the source in the region covered by the calibration.

In the present analysis, we use both single–phase and multi–phase models, considering that we are mapping the very central part of the cluster core. The single–phase model just assumes an emission from an optically–thin plasma (Mekal –Kaastra 1992, Liedhal et al. 1995– in XSPEC v. 11.1.0 –Arnaud 1996) absorbed by a column density that we have left free to vary (hereafter, model abs1T; the reference column density is the Galactic value of $`1.2\times 10^{20}`$ cm-2 from radio HI maps in Dickey & Lockman 1990).

The multi–phase model combines thermal emission from the ambient, outer gas with a continuous distribution of gas states represented by a cooling–flow model (Johnstone et al. 1992). An absorbing column density fixed to the Galactic value is put in front of the combined model. The cooling–flow component is also intrinsically absorbed by uniformly distributed amount of hydrogen atoms per cm2 at the cluster redshift that is left free to vary during the search for a minimum $`\chi ^2`$ (absCF model). The temperature and metallicity of the cooling–flow model are fixed to be equal to the ones of the thermal component. All the absorption of X-rays due to the interstellar medium has been parametrized using the Tübingen-Boulder model (tbabs in XSPEC v. 11.1.0; Wilms, Allen & McCray 2000).

Two representative spectra are shown in Figure 6 with the respective abs1T model. In Table 2, we quote the best-fit results for the gas temperature and metallicity obtained with a single-phase plasma model (abs1T) and with a multi-phase component (absCF). The minimun $`\chi ^2`$ and the number of degrees of freedom for these models are quoted in Table 3. In Figure 7, we plot the best-fit values obtained with the model abs1T for some interesting quantities and, where available, the region of uncertainty at the 90 per cent confidence level from the analysis of ASCA data in Allen (2000). Steep temperature and metallicity profiles are present within the inner 300 kpc. The emission-weighted gas temperature increases from 3.4 to 6.5 keV, whereas the metal abundance decreases with radius from about 0.7 times the solar photospheric values in Anders & Grevesse (1989) to $``$0.2. The column density is consistent with the Galactic value above 50 kpc, but it is higher by a factor of 2 in the inner two annuli. A remarkable good agreement is present between the redshift measured independently in each ring and the optical estimate from Girardi et al. (1998).

In Fig. 8 and Table 3, we show that (i) the null-hypothesis probability on the goodness of the fits is larger than 0.1 on 6 out of nine annuli for both the abs1T and absCF models with the worst reduced $`\chi ^2`$ of about 1.27 in the second and sixth annulus, and (ii) that the cooling flow model, absCF, is equivalent to the single-phase model, abs1T, from a statistical point of view when the F (Fisher) test is considered (e.g. Bevington & Robinson 1992; on how a single-phase gas can mimic a multi-phase medium through projection effects see, e.g., Ettori 2001).

In the same Fig. 8, we compare also the spectral results for the FP–110°and FP–120°datasets. Gas temperatures, metallicities and column densities show differences generally between 1 and 3 $`\sigma `$ and a clear trend with temperature and metal abundance estimates higher and column density values lower in the FP–120°dataset.

## 5 Deprojection analysis of the spectral and imaging data

In this Section, we discuss the constraints on the physical quantities obtained from the deprojection analysis applied to both the best-fit results of the spectral analysis presented in the previous section (hereafter spectral deprojection) and the surface brightness profile (spatial deprojection).

The deprojection analysis makes the assumption that the emission is spherically symmetric and the X-ray emitting plasma is in hydrostatic equilibrium with the underlying gravitational potential, which is likely to be reasonable on intermediate to larger scales considering that the sound crossing time in the inner 100 kpc is less than $`10^8`$ years.

Considering that we are not able to distinguish between single phase and cooling flow models in the fitting analysis of our spectra as discussed in Section 4, we adopt hereafter for simplicity the single phase description to infer the properties of the intracluster medium.

### 5.1 Spectral deprojection

The physical quantities constrained from the projected spectra need to be translated to their 3-dimensional values to allow us to investigate the proper characteristics of the X-ray emitting plasma. Fitting a thermal model to a spectrum obtained collecting X-ray counts in rings provides, for each annulus, (i) an estimate for the Emission Integral, $`EI=n_\mathrm{e}n_\mathrm{p}𝑑V=0.82n_\mathrm{e}^2𝑑V`$, through the normalisation $`K`$ of the model, $`K=\frac{10^{14}}{4\pi d_{\mathrm{ang}}^2(1+z)^2}EI`$ (see Mekal model in XSPEC); (ii) a direct measurement of the emission-weighted gas temperature, $`T_{\mathrm{ring}}`$, metal abundance, $`Z_{\mathrm{ring}}`$, and luminosity, $`L_{\mathrm{ring}}`$. The purpose of the deprojection is, for example, to recover the value of the gas temperature in shells, $`T_{\mathrm{shell}}T_i`$, that is defined as

$$T_{\mathrm{ring}}T_j=\frac{_{i,\mathrm{outer}\mathrm{shell}}^{i=j}T_iw_{ij}}{_{i,\mathrm{outer}\mathrm{shell}}^{i=j}w_{ij}}$$
(1)
where $`w_{ij}=L_i\times Vol(i,j)/Vol(i)=ϵ_iVol(i,j)`$ provides the luminosity for a given shell $`i`$ with volume $`Vol(i)`$ weighted by the part of this volume projected on the ring $`j`$, $`Vol(i,j)`$. Using this notation, it is simple to note that $`L_{\mathrm{ring}}L_j=_{i,\mathrm{outer}\mathrm{shell}}^{i=j}ϵ_iVol(i,j)=_{i,\mathrm{outer}\mathrm{shell}}^{i=j}w_{ij}`$.

From Kriss, Cioffi & Canizares (1983; see also McLaughlin 1999), the volume of each shell observed through each ring adopted in the spectral analysis can be evaluated and a matrix, $`\mathrm{𝐕𝐨𝐥}`$, can be built with components equal to the parts the volume of the shells (rows $`i`$) seen at each ring (or annuli; column $`j`$).

The deprojected physical quantities can be then obtained through the following matrix products (shown by the symbol $`\mathrm{\#}`$):

$$\begin{array}{c}n_\mathrm{e}=\left[(\mathrm{𝐕𝐨𝐥}^T)^1\mathrm{\#}(EI/0.82)\right]^{1/2}\hfill \\ ϵ=(\mathrm{𝐕𝐨𝐥}^T)^1\mathrm{\#}L_{\mathrm{ring}}\hfill \\ ϵT_{\mathrm{shell}}=(\mathrm{𝐕𝐨𝐥}^T)^1\mathrm{\#}(L_{\mathrm{ring}}T_{\mathrm{ring}})\hfill \\ ϵZ_{\mathrm{shell}}=(\mathrm{𝐕𝐨𝐥}^T)^1\mathrm{\#}(L_{\mathrm{ring}}Z_{\mathrm{ring}}),\hfill \end{array}$$
(2)
where $`(\mathrm{𝐕𝐨𝐥}^T)^1`$ indicates that the matrix is firstly transposed and then inverted. The emission due to the shells projected along the line of sight but with the corresponding annuli outside the field-of-view is taken into account with an edge correction factor estimated assuming a power law distribution of the emission proportional to $`r^4`$ (cfr. equation A8 in the appendix of McLaughlin 1999).

We have applied this technique to the single-phase results of our spectral analysis (model abs1T in Sect. 4). In rings where part of the flux was masked for the presence of a point-source, we correct the normalisation $`K`$ by the relative amount of area not considered. (We note that, if we fit simultaneously the nine spectra with a combination of single phase models properly weighted with the correspondent volume of the shell observed at each ring, we measure a $`\chi ^2`$ of 2604.3 with 2278 degrees of freedom that is statistically undistinguishable from the sums of the results in Table 2 of $`\chi ^2=`$2576.4 and 2261 degrees of freedom).

The deprojected temperature gradient can then be compared with the output of the spatial deprojection code to constrain the gravitational potential and the mass deposition rate. The error bars come from 100 Monte-Carlo simulations obtained from scattering the original projected input with respect to their Gaussian error.

The deprojected gas temperature and metallicity are slightly steeper than the observed, projected profiles (cf. Fig. 7). The gas temperature rises from $`2.6_{0.5}^{+0.3}`$ to $`6.5_{0.3}^{+0.4}`$ keV at 270 kpc. The profile can be fitted ($`\chi ^2`$=10.5, 7 d.o.f., P=0.16) with a power law: $`T_{\mathrm{gas}}=2.74_{0.23}^{+0.41}\times (r/10\mathrm{k}\mathrm{p}\mathrm{c})^{0.27(0.06,+0.04)}`$ keV. The metal (iron) abundance decreases from a central value of about 0.8 times the solar abundance to $`0.18(\pm 0.07)`$.

### 5.2 Spatial deprojection

To reach more firm conclusions on the state of the intracluster plasma in the inner region, we have combined the results of the spectral deprojection with the spatial deprojection of the X-ray surface brightness profile, which has been done applying the code developed from Fabian et al. (1981; see details on the code in White, Jones & Forman 1997). The surface brightness profile has been obtained in the 0.5–7 keV energy range with a bin-size of 9.8 arcsec (i.e. 20 times the original ACIS bin-size). Considering that the cluster emission covers the whole CCD, we have estimated a background count-rate of about 1.5 $`\times 10^6`$ cts s-1 arcsec-2 from the blank-field used in the spectral analysis (Sect. 4). We detect signal up to 4.4 arcmin ($``$ 440 kpc) from the assigned centre with $`5\sigma `$ confidence.

The results from the spectral deprojection analysis are used to constrain the process of the spatial deprojection in the way described below (and similar to the procedure adopted in Allen, Ettori & Fabian 2001, Schmidt, Allen & Fabian 2001). In the spatial deprojection, there are 5 free parameters: (i) the column density, $`N_\mathrm{H}`$, fixed to the averaged value observed in the spectral analysis of about $`1.6\times 10^{20}`$ cm-2; (ii) the metallicity fixed to 0.3 times the solar value in accordance with the mean behaviour of the deprojected abundance profile; (iii) the outer pressure, fixed to the value of $`27.6\times 10^4`$ K cm-3, extrapolated from the pressure profile estimated from the deprojected spectral temperature and gas density profiles; (iv-v) the scale (core) radius and the velocity dispersion describing the assumed potential law \[in the present case either an isothermal sphere (Binney & Tremaine 1987, eqn. 4-124b) or a Navarro, Frenk & White profile (NFW, 1997)\]. The last two parameters are obtained from a $`\chi ^2`$ minimization between the gas temperature profile obtained from spectral deprojection of the single-phase analysis results (cf. Sect. 4) and the output of the spatial deprojection code properly rebinned (cf. eqn.1) to match the number of bins of the spectral profile. (In the spatial deprojection, the errors are obtained from the distribution of the results of the deprojection of 100 Monte-Carlo simulations of the original surface brightness profile.)

The best-fit between the two temperature profiles are obtained for the following values of the two free parameters of the potential laws: $`r_\mathrm{c}=0.10_{0.01}^{+0.01}`$ Mpc, $`\sigma _{\mathrm{IS}}=\sqrt{(4/9)\pi G\rho _0r_\mathrm{c}^2}=658_{11}^{+14}`$ km s-1 for the isothermal sphere; $`r_\mathrm{s}=0.49_{0.13}^{+0.11}`$ Mpc, $`\sigma _{\mathrm{NFW}}=\sqrt{(GM_{200})/(2r_{200})}=\sqrt{50}H_0r_{200}=837_{85}^{+69}`$ km s-1 when a NFW potential law is assumed. Note that $`\rho _0`$ and $`r_\mathrm{c}`$ define the central matter density and the core radius, respectively, in the isothermal sphere, whereas $`M_{200}`$ indicates the total mass within the radius $`r_{200}`$ where an overdensity of 200 with respect to the background is reached. Note also that the two velocity dispersions are related by the equation $`2\sigma _{\mathrm{NFW}}^2=(d\mathrm{ln}\rho /d\mathrm{ln}r)\sigma _{\mathrm{IS}}^2`$ (cfr. eqn.4-127b in Binney & Tremaine 1987) and tend to be equal when regions extending at $`r10r_\mathrm{c}`$ ca be considered (see Figure 4-8 in Binney & Tremaine 1987; in our case, we are mapping regions within $`4r_\mathrm{c}`$). However, we observe a lower $`\chi ^2`$ for a potential with a flat profile in the core, namely $`\chi ^2=`$ 7.1 (with 7 degrees of freedom, P=0.41) for the isothermal sphere against 16.5 (P=0.02) for a NFW profile.

In Fig. 9, we show the gas temperature and density profiles from the spectral deprojection analysis compared to those obtained from the best-fit results in the spatial deprojection. It is worth noting that the good agreement between the two gas density profiles is due to the dominance of the density in the X-ray emissivity that allows an almost temperature-independent deprojection, whereas we have required that the two temperature profiles match, searching for a minimun $`\chi ^2`$ in comparing these profiles with the parameters of the gravitational potential left free to vary. Moreover, this comparison shows the self-consistency of our results and provides us with a cooling time and a mass deposition rate profiles with resolution of the imaging data (see discussion in Section 5.5).

### 5.3 Cluster mass profile in the core

In our spectral analysis, we are able to resolve in 5 radial bins the gas temperature profile in the inner 100 kpc radius with a relative uncertainty less than 5 (10) per cent for the projected (deprojected) temperature. This enables us to investigate the mass profile in the cluster core assuming that the intracluster medium is there in hydrostatic equilibrium with the dark matter potential.

We make use of the gas temperature and electron density values estimated from the deprojected best-fit spectral results obtained with a single-phase modeli (see Section 5.1). Then, we apply the equation of hydrostatic equilibrium between the gravitational potential and the intracluster plasma to estimate the mass profile, $`M_\mathrm{X}`$,

$$M_\mathrm{X}=\frac{Tr^2}{G\mu m_\mathrm{p}}\left(\frac{d\mathrm{ln}T+d\mathrm{ln}n}{dr}\right).$$
(3)
In doing this, we interpolate the gas temperature profile with a polynomial of 2nd order to smooth over any fluctuations for the results from spectral deprojection. Moreover, we estimate the mass profile using a smoothed temperature and density profiles obtained from the best-fit results of the spatial deprojection (see Section 5.2). We assume a relative error propagated from the results of the spectral deprojection of about 15 per cent consistent with 1.5 times the uncertainty on the gas temperature measurements. In Fig. 10, we compare these dark matter profiles with the best-fit results obtained fitting the outer region of the cluster emission observed with Rosat PSPC (for details see Ettori & Fabian 1999) both with a $`\beta `$model (Cavaliere & Fusco-Femiano 1976; $`r_\mathrm{c}`$ =0.25 Mpc, $`\beta `$ = 0.75) and a gas model obtained from the NFW potential law ($`r_\mathrm{s}`$ =0.77 Mpc, $`\eta `$ = 10.70).

To assess the shape of the underlying gravitational potential in the central region, we describe the dark matter profile with a power law expression, $`\rho _{\mathrm{grav}}=\rho _0(r/r_0)^\alpha `$, integrate it over the volume

$$M_{\mathrm{grav}}=_0^R4\pi \rho _{\mathrm{grav}}r^2𝑑r=\frac{4\pi \rho _0r_0^3}{3\alpha }\left(\frac{R}{r_0}\right)^{3\alpha }$$
(4)
and fit this mass profile to the one from the spectral deprojection results (see Fig. 10) between 10 and 100 kpc, where the discrepancy among the different forms of the potential is more significant. We measure $`\alpha =0.59_{0.17}^{+0.12}`$ (in the range 0.27–0.81 at the 90 per cent confidence level; Fig. 11), that indicates a remarkable flat profile suggesting the presence of a core. In particular, a value of 1 (typical for cuspy dark matter profile like NFW) is excluded at more than 3 $`\sigma `$. In other words, we expect four times more than the observed mass in the inner bin, if we assume a power law index of 1 and a fixed mass at 100 kpc. On the other hand, the deviation from a flat core (i.e. $`\alpha =`$ 0) makes the mass profile lie apart from that expected given a King profile. The presence of this “shoulder” has been detected also by Xu et al. (1998) from ASCA and Rosat PSPC data and is now resolved with Chandra . Finally, from the presence of a shear in the $`H_\alpha `$ filament southward the cD galaxy, Hu, Cowie and Wang (1985) conclude that the central gravitational potential is not too deep. Following their arguments on the tidal shearing of the filament in the cluster potential, we obtain an upper limit on the central mass (assuming constant density) consistent with our results and that excludes any dark matter density profile sharper than the NFW form below 10 kpc (see Fig. 10).

It is worth noticing that this result in not function of the X-ray center adopted in our analysis (see discussion in Section 3). As shown in Fig. 10, moving the center of the circular annuli used to collect the spectra to the cD galaxy does not change the mass profile. Moreover, the presence of a multi-phase gas could affect our estimate of the gravitational potential. Gunn & Thomas (1996) show that, given the same emissivity, the single phase assumption underestimates the total mass by 20–40 per cent. Considering that the evidence for a truly multi-phase gas is weak and the gas is not certainly described from a steady-state cooling flow model (see, e.g., the discrepancy between the deposition rate obtained from spatial and spectral deprojection in Section 5.5), the correction factors for A1795 would be smaller than these. Adopting these upper limits and increasing by 40 per cent the mass values in the two inner bins do not change the results significantly (cf. Fig. 11).

### 5.4 On the dynamical state of the core

The analysis of the surface brightness distribution in Section 3 suggests that the intracluster medium has an unrelaxed nature that might introduce a kinetic pressure component in the hydrostatic equilibrium equation so raising the total effective mass in the central core. This could possibly explain the flattening of the gravitating mass profile. Assuming that underlying potential is described from NFW and $`\mathrm{\Delta }M`$ is the difference with respect to the observed mass estimates, a bulk motion with velocity (assumed independent from radial position) of $``$ 300 km s-1 is required from the equation $`d(v^2\rho _{\mathrm{gas}})/dr=(G\rho _{\mathrm{gas}}\mathrm{\Delta }M)/r^2`$. This value is lower than the gas sound speed of about 800 km s-1 and definitely subsonic.

It is worth noting that the shape of the gravitational potential is consistent with the detected motion of the cD galaxy. Oegerle & Hill (1994) measure a cD velocity redshifted relative to the cluster of $``$ 150 km s-1. The emission-line velocity map in Hu et al. (1985) shows that the H$`\alpha `$/X-ray filament has the same velocity of the cluster and is blueshifted from the cD velocity by still 150 km s-1 (from recent H$`\alpha `$ observation, this can be considered as un upper limit; C. Crawford, priv. comm.). In Fig. 12, the cD velocity is estimated as $`v_{\mathrm{cD}}=(2GM_{\mathrm{grav}}/r)^{0.5}`$, where $`G`$ is the gravitational constant and $`M_{\mathrm{grav}}`$ is the total gravitating mass within the radius $`r`$ calculated from eq.4. Fixing the gravitational mass at $`r`$= 100 kpc, we can change the power law index $`\alpha `$ and investigate the behaviour of the cD velocity. As shown in Fig. 12, a cD velocity larger than 400 km s-1 is required for potential well described by a dark matter density profile steeper than $`r^1`$. Given the observed velocity offset along the line of sight, the presence of a density profile flatter than a NFW profile is then required.

If we model the motion of the dominant galaxy as an harmonic oscillation at the bottom of the potential well (e.g. Lazzati & Chincarini 1998), we assume that the cD galaxy passed through the X-ray cluster centroid at the time zero with maximal velocity and, then, underwent an oscillatory motion described from the equations, $`x_{\mathrm{cD}}=A\mathrm{cos}(\omega t+\varphi ),v_{\mathrm{cD}}=A\omega \mathrm{sin}(\omega t+\varphi )`$, where the phase $`\varphi =\pi /2`$ and $`\omega =(GM_{\mathrm{grav}}/r^3)^{0.5}5\times 10^{16}`$. Solving the above two equations, we find that the cD galaxy crossed the centroid about 1.3 $`\times 10^8`$ years ago with velocity of $`400`$ km s-1. Considering our interpretation of the X-ray filament pointing southward as a cooling wake (see Fabian et al. 2001b for a discussion), we can use this age as an upper limit on the cooling time of the densest, coolest phases in the central plasma.

### 5.5 On the cooling flow model

As discussed in Section 4, we cannot conclude that the gas in the core of A1795 is multi-phase. However, this gas would cool radiatively in about $`10^9`$ years (3.8 $`\times 10^8`$ yrs in the central 10 kpc radius, with 10th and 90th percentile of 3.2 and 4.4 $`\times 10^8`$ yrs, respectively), approaching the estimated age of the Universe at the cluster redshift of $`1.2\times 10^{10}`$ years at about 200 kpc (Fig. 13). This can be considered as an upper limit on the region where cooling is taking place. In particular, the physical extension of the central region, where a X-ray soft component dominates the emission, can be obtained from the X-ray colour profile shown in Fig. 5, where the profile flattens at about 100 kpc. This is consistent with the break observed in the $`\dot{M}(r)`$ profile discussed below and presented in Fig. 13.

To maintain the cluster core in pressure equilibrium as energy is lost to radiation, and if there is not balancing from any heating sources, a flow of material occurs inward and deposition of cool gas mass takes place in the central region of the cluster. In this cooling flow scenario (e.g. Fabian 1994), the amount of gas deposited within a given radius, $`\dot{M}(<r)`$, can be estimated and used to parametrize the luminosity originating from the cooling materials. We measure $`\dot{M}(<r)`$ in two ways (Fig. 13): in the spatial deprojection, it is provided from the luminosity associated to the gas that either (i) cools completely in a shell or (ii) passes the shell changing temperature under the action of the cluster potential (cf. White, Jones & Forman 1997); in the spectral analysis, $`\dot{M}(<r)`$ is estimated as the amount of gas that cools from the ambient temperature to zero through the modelling of the spectra accumulated in concentric circular region with a multi-phase model (Johnstone et al. 1992) absorbed by a column density intrinsic to the galaxy cluster that is required to reproduce properly the soft part of the spectra.

The cumulative values of $`\dot{M}(<r)`$ obtained in the present analysis (Fig. 13) enclose the constraints available from the ASCA analysis in Allen et al. (2001), where no spatial resolution of $`\dot{M}(<r)`$ was possible. The deviation between the spatial and spectral results on the $`\dot{M}`$ profile, apart from the central two bins where the assumed spherical geometry is probably inappropriate given the observed cooling wake, is due to limitations on the validity of the steady–state cooling flow model. The dynamical scenario in which cooling flows establish and evolve considers merging with infalling substructures that interrupt the subsonic flow of material. These mergers are been shown to be responsible for several features in the surface brightness and temperature distribution (see, e.g., Markevitch et al. 2001, Vikhlinin et al. 2001, Mazzotta et al. 2001) and can be energetically relevant, thus throwing into question the assumption of steady–state mass accretion underlying the standard cooling flow model (e.g. Nulsen 1986). Therefore, we adopt hereafter the spectral estimates of the mass deposition rate as more reliable and use these as reference values.

Recent analysis of XMM-Newton data of A1795 (Tamura et al. 2001, Molendi & Pizzolato 2001) does not show detectable emission from gas cooling below 1–2 keV. Our lower limit from the deprojection of the best-fit gas temperature in the central 20 kpc radius is 1.8 keV at 90 per cent confidence level. When an isobaric cooling flow component is considered to model the XMM-Newton Reflection Grating Spectrometers (RGS) spectra of the central region, an upper limit (90 per cent level of confidence) of 150 $`M_{}`$ yr-1 is obtained (Tamura et al. 2001; note that these authors estimate that about 80 per cent of the emission observed in the RGS spectra comes from regions enclosed within $``$ 60″). This value is larger than our value of 74 $`M_{}`$ yr-1 (90 per cent confidence limit of 89 $`M_{}`$ yr-1) measured within 75 arcsec from the X-ray centre. The Chandra -determined cooling flow therefore appears to be consistent with the present XMM-Newton RGS constraint. Within the central 200 kpc, we estimate a deposition rate of about 100 $`M_{}`$ yr-1 ($`<`$ 121 $`M_{}`$ yr-1 at 90 per cent c.l.), consistent with the early ASCA result of Fabian et al. (1994b) and lower than the estimate in Allen et al. (2001), where the assumption of an isothermal gas in the core (due to the limited spatial resolution of ASCA ) causes an overestimate of the integrated deposition rate (up to a factor of $``$ 3 in the case of A2390 discussed in Allen, Ettori & Fabian 2001). When a cooling flow model with a low temperature cut-off is considered (e.g. Tamura et al. 2001), good fits can be obtained for the inner 3 annuli with lower estimates of $`T_{\mathrm{cut}\mathrm{off}}`$ 1.4 keV. With respect to absCF, this model preserves the same number of degrees-of-freedom and provides a better $`\chi ^2`$ just for the second annulus that encloses a region between 10 and 25 arcsec ($`\mathrm{\Delta }\chi ^2`$= –0.7, 3.3, –0.8 for the inner three annuli, respectively, when compared with the results in Table 2). This is qualitatively in agreement with the results from XMM-Newton EPIC analysis in Molendi & Pizzolato (2001).

Note that we model the cluster emission with absCF, in which the cool emission is suppressed with an intrinsic absorption of about $`10^{21}`$ particle cm-2. If we adopt the same model that describes the RGS spectra and consider a region of 75 arcsec radius (80 per cent of the emission observed in the RGS spectra comes from a projected radius of about 60 arcsec), fixing the outer thermal component to 6.4 keV with metallicity of 0.4 solar and the Galactic column density to 3 $`\times 10^{20}`$ cm-2, and not including any intrinsic absorption, we measure a normalization for the cooling flow component of 140 $`M_{}`$ yr-1 (90 per cent confidence limit range of 133–152 $`M_{}`$ yr-1) with a considerably worst $`\chi ^2`$ of 654 (387 d.o.f.) with respect to 535 (384) obtained with absCF. On the other hand, a cooling flow model with a low temperature cut-off provides a better fit ($`\chi ^2=`$ 505, 383 d.o.f.; F-test probability of 3 $`\times 10^3`$) with $`T_{\mathrm{cut}\mathrm{off}}`$ 1.4 keV (in the range 1.3–1.6 keV at 90 per cent c.l.) and normalization $`\dot{M}`$ 324 $`M_{}`$ yr-1 (range at the 90 per cent c.l.: 271–469 $`M_{}`$ yr-1). The intrinsic absorption discussed above can be interpreted as one method to suppress the line emission from gas below $``$ 1.5 keV. Other methods and related issues are discussed by Peterson et al. (2001) and Fabian et al. (2001a).

### 5.6 Distribution of the metals

The X-ray emitting plasma is about 5 times more massive than the stars in galaxies (David et al. 1990, White et al. 1993) and, thus, is the major reserve of both baryons and heavy metals contained in galaxy clusters, since the Iron abundance is typically between 0.3 and 0.5 times the solar value (Mushotzky & Loewenstein 1997). Studies on the correlation between Fe mass and light coming from E and S0 galaxies (Arnaud et al. 1992, Renzini 1997) conclude that larger amount of iron resides in the intracluster medium than inside galaxies and its enrichment originated through releases from early-type galaxies.

The processes that preferably enrich the plasma are (i) (proto)galactic winds (De Young 1978, Metzler & Evrard 1994), that occur at early times and are characterized from supernova (SN) type II ejecta with large abundance of $`\alpha `$ elements, and (ii) ram pressure stripping (Gunn & Gott 1972), that takes place on longer time scales due to the continuous accretion of fields galaxies in the cluster potential well and produces mostly SN Ia ejecta. A negative gradient in the metallicity profiles is observed preferentially in cooling flow clusters (Allen & Fabian 1998, Irwin & Bregman 2001, De Grandi & Molendi 2001) and can be explained as enrichment of the inflowing gas by SNe Ia and stellar mass loss in the outer parts of the central dominant galaxy (Reisenegger, Miralda-Escudé & Waxman 1996).

What discriminates between these two main processes is the different elemental mass yields: SN Ia ejecta tend to be rich in Ni, whereas SN II ejecta present larger ratio between $`\alpha `$ elements (e.g., O, Mg, Ar, Ca, S, Si) and Fe. Recent evidence of iron gradients in cluster cores with decreasing ratio between $`\alpha `$ elements and Fe moving inward suggest that, while the global intracluster metal abundances are consistent with SN II ejecta (Mushotzky & Loewenstein 1997), SN Ia productions are dominant in the central cluster regions (Ishimaru & Arimoto 1997, Fukuzawa et al. 1998, Finoguenov, David & Ponman 2000, Allen et al. 2001, Dupke & Arnaud 2001).

To resolve and quantify the contribution from the single heavy elements, we fit the annular spectra with an absorbed Vmekal model (abs1TZ in Table 3) using the solar elemental abundances from Grevesse & Sauval (1998) that update the commonly-used values in Anders & Grevesse (1989; note that the results on the azimuthally-averaged metallicity in Table 2 and Fig. 7 refer at these values for a direct comparison with previous work) and differ from these mainly in the abundance of Iron, Oxygen and Sulphur by a factor of 0.68, 0.79 and 1.32, respectively. We fix Helium to the cosmic abundance value and consider 5 groups of metals: O, C, N, Na, Mg and Al tied together; Ne; Si; S; Fe tied with Ca, Ar and Ni. In Fig. 14, we show the profiles for Ne, Si, S and Fe. We omit to plot the constraints on the value of the O abundance that is approximately flat around 0.4 times solar. In general, these fits show a null-hypothesis probability higher than the corresponding single abundance fit (right panel in Fig. 14).

To estimate the Fe mass fraction originated from SN Ia and SN II ejecta, we perform a $`\chi ^2`$ minimization at each radial bin of the merit function

$$\chi ^2=\underset{i}{}\frac{(ratio_imod_i)^2}{err_i^2}$$
(5)
where $`ratio_i`$ = \[O/Fe, Ne/Fe, S/Fe, Si/Fe\], $`err_i`$ is the 1 $`\sigma `$ error on $`ratio_i`$ propagated from the uncertainties in the fit and $`mod`$ is equal to $`f\times ratio_{\mathrm{SNIa}}+(1f)\times ratio_{\mathrm{SNII}}`$. The elemental abundance ratios for SN Ia (0.037, 0.006, 0.585, 0.538 for O/Fe, Ne/Fe, S/Fe, Si/Fe, respectively) and SN II (3.82, 2.69, 2.29, 3.53) yields are from Nomoto et al. (1997a, 1997b; for SN Ia yields, we refer to deflagration model, W7) and summarized in Table 4 in Dupke & Arnaud (2001). The radial profile of $`f_{\mathrm{SNIa}}`$ is plotted at the bottom of the left panel in Fig. 14. The fraction in number of the SNe Ia is above 65 per cent everywhere within 200 kpc.

This is consistent with the higher null-hypothesis probability that we observe when the abundance ratios are fixed to the SN Ia production values instead of SN II values (right panel in Fig. 14).

A supernova produces 0.74 $`M_{}`$ of iron when Type Ia (Thielemann, Nomoto & Hashimoto 1996), $``$ 0.12 $`M_{}`$ when Type II (Gibson, Loewenstein & Mushotzky 1997). These values imply a total iron mass of $`M_{\mathrm{Fe}}=0.74f_{\mathrm{SNIa}}N_{\mathrm{SNIa}}+0.12(1f_{\mathrm{SNIa}})N_{\mathrm{SNII}}`$, where $`N_{\mathrm{SNIa}}`$ and $`N_{\mathrm{SNII}}`$ are the number of SN Ia and SN II, respectively, and are estimated comparing this calculation of the iron mass with what measured from the observed iron abundance, $`Z_{\mathrm{Fe}}`$, plotted in Fig. 7. In the latter case, $`M_{\mathrm{Fe}}=A_{\mathrm{Fe}}y_{\mathrm{Fe},}Z_{\mathrm{Fe}}M_\mathrm{H}`$, where $`A_{\mathrm{Fe}}=55.8`$, $`y_{\mathrm{Fe},}=3.16\times 10^5`$ from Grevesse & Sauval (1998) and $`M_\mathrm{H}=M_{\mathrm{gas}}/1.33`$ (we are assuming $`\mu =0.6`$, $`n_\mathrm{e}=1.21n_\mathrm{p}`$), and rises with radius from $`3\times 10^7M_{}`$ in the inner 10 kpc up to few times $`10^9M_{}`$ at $`r>100`$ kpc. From the compilation on the rates of SNe Type Ia and II in Madau, Della Valle & Panagia (1998), a conversion factor between 1 and 5 is obtained from $`N_{\mathrm{SNII}}`$ to $`N_{\mathrm{SNIa}}`$. Using a central value of 3.5 (see also Iwamoto et al. 1999), we require a number of SN II that ranges between $`10^8`$ and $`4\times 10^{10}`$ moving outward.

The amount of the total iron produced in SNIa is then $`M_{\mathrm{Fe},\mathrm{SNIa}}/M_{\mathrm{Fe}}\left[1+3.5\times (0.12/0.74)\times (f_{\mathrm{SNIa}}^11)\right]^1(0.81)`$ in the central 200 kpc.

Assuming that about 80 $`M_{}`$ of star formation are required to generate a supernova Type II (e.g. Thomas & Fabian 1990), we conclude that the cumulative star formation that took place during the cluster history is $`0.72\times 10^{12}M_{}`$ within the central 100-200 kpc, consistent with the expected amount of gas deposited during the existence of the cooling flow. From the $`H_\alpha `$ luminosity of about $`10^{42}`$ erg s-1 coming from the central 12 kpc radius (Cowie et al. 1983), we infer a present star formation rate (SFR) of about 9 $`M_{}`$ yr-1 \[Kennicutt 1983; Crawford et al. (1999) indicate a visible SFR of about 2 $`M_{}`$ yr-1 within 5 kpc-radius from optical observation of the central dominant galaxy\], that has to be compared with an inferred cumulative value of about $`9\times 10^9`$ $`M_{}`$ on the same region in exam. These values can agree with constant star formation rate over the last Gyr.

The deposition of material, star formation activity and residuals in form of absorbing dust (as suggested for central regions of cooling flow clusters; e.g. Voit & Donahue 1995, Fabian et al. 1994a, Allen 2000) can (i) explain the observed excess in the column density measurements in the inner 50 kpc of about 2 with respect to the Galactic value interpolated from radio maps, and (ii) be responsible for emerging reprocessed emission in the far infrared. From IRAS scans over the central 4 arcmin of A1795, Allen et al. (2001) infer a far-infrared luminosity of $`<3.3\times 10^{44}`$ erg s-1, that puts an upper limit on the present SFR of 15 $`M_{}`$ yr-1 (from relation in Kennicutt 1998), in agreement with the previous estimates.

Finally, this star formation activity releases energy during supernova explosions providing an amount of thermal energy per gas particle, $`E`$

$$E=\frac{\eta E_{\mathrm{SN}}N_{\mathrm{SNII}}}{M_{\mathrm{gas}}/(\mu m_\mathrm{p})}0.4\left(\frac{\eta }{0.1}\right)\left(\frac{N_{\mathrm{SNII}}}{10^{10}}\right)\left(\frac{10^{12}M_{}}{M_{\mathrm{gas}}}\right)\mathrm{keV},$$
(6)
where a kinetic energy release by one SN, $`E_{\mathrm{SN}}`$, of about 1.3 $`\times 10^{51}`$ erg (e.g. Iwamoto et al. 1999), a factor $`\eta `$ that represents the efficiency of this kinetic energy in heating the ICM through galactic winds (and assumed here equal to 10 per cent) and typical values for the inner part of A1795 are adopted. When we compare $`E`$ to the thermal energy per particle measured (cf. Fig. 7), we conclude that processes from SNe explosions are responsible for about 8 per cent of that energy.

## 6 Summary and Conclusions

The entropy of the gas increases monotonically moving outwards almost proportionally to the radius, with a best-fit $`\mathrm{log}(K)=0.09+0.97\mathrm{log}(r)`$ (see Fig. 15). This behaviour traces the shock heating (more significant at larger radius) at which the intracluster gas has been subjected during its history of accretion in the gravitational potential well, whereas the efficiency of the cooling has removed the central gas at lower entropy (cf. Tozzi & Norman 2001).

Together with the estimated cooling time due to energy loss by radiation of less than $`10^9`$ years and a decrease by a factor of 3 in the gas temperature, this seems to indicate that the central 200 kpc are undergoing significant cooling with a gas mass deposition rate of about 100 $`M_{}`$ yr-1 in the absence of any heating process. This would require an intrinsic absorption of about $`10^{21}`$ cm-2. However, a simple single phase model with absorption in excess of the Galactic value by a factor of 2, within the inner 50 kpc, is consistent with the observed spectra.

The gravitational potential follows a NFW profile between 60 and 300 kpc, becoming flatter within 60 kpc with a power law index of about $`0.6`$ once a power law description is adopted for the total matter density. The shape of the potential is in agreement with the motion of the central dominant galaxy and suggests that the central cluster region is not relaxed.

This dark matter distribution, combined with the integrated gas mass within 300 kpc, provides a gas mass fraction ranging between 10 and 20 per cent (Fig. 16), in a good agreement with the estimate of $`0.16\pm 0.01`$ at 500 kpc from Rosat PSPC brightness profile and ASCA temperature in Ettori & Fabian (1999) and completely consistent with that generally observed in galaxy clusters (Evrard 1997, Ettori & Fabian 1999, Mohr et al. 1999).

Finally, the present Chandra observation of the core of A1795 allows us to resolve radially the contribution from Oxygen, Neon, Sulphur, Silicon and Iron to the plasma metallicity and, consequently, to estimate between 80 and 100 per cent the amount of Iron produced in supernova Type Ia events within the central 200 kpc.

## ACKNOWLEDGEMENTS

We thank all the staff members involved in the Chandra project for such an excellent performance of the instruments. ACF is grateful to NASA for the opportunity to participate to the project as an InterDisciplinary Scientist. ACF and SWA thank the Royal Society for support. We thank the referee for a careful reading of the manuscript and for several suggestions that improved the presentation of this work.


format: text

text:
# On the duration of the subsonic propeller state of neutron stars in wind-fed mass-exchange close binary systems

## 1 Introduction

The sequence of states which a magnetized neutron star in a wind-fed mass-exchange binary system follows as it spins down from the initially very short periods can be expressed in the form of the following chain: ejector $``$ propeller $``$ accretor. This classification, first suggested by Shvartsman (sh70 (1970)), reflects three different evolutionary stages and three different mechanisms of energy release responsible for the neutron star emission.

The spindown of a neutron star in the state of ejector is governed by the canonical spin-powered pulsar mechanism. The spindown power dominates the star energy budget and is spent predominantly to the generation of magneto-dipole waves and particle acceleration. The pulsar-like spindown ceases when the pressure of the material being ejected by the neutron star can no longer balance the pressure of the surrounding gas. The stellar wind plasma penetrates into the neutron star light cylinder and interacts with the star magnetosphere. This corresponds to the neutron star state transition: ejector $``$ propeller.

Neutron star in the state of propeller is spinning down due to the interaction between its fast rotating magnetosphere and the surrounding material. Davies et al. (dfp79 (1979)) and Davies & Pringle (dp81 (1981)) have shown that during this state the star magnetosphere is surrounded by a spherical quasi-static envelope in which the plasma temperature is of the order of the free-fall temperature,

$$T(R)T_{\mathrm{ff}}(R)=\frac{GM_{\mathrm{ns}}m_\mathrm{p}}{kR}.$$
(1)
Here $`M_{\mathrm{ns}}`$ is the mass of the neutron star, $`m_\mathrm{p}`$ is the proton mass and $`k`$ is the Boltzmann constant. The rotational energy loss by the neutron star is convected up through the envelope by the turbulent motions and lost through its outer boundary. The structure of the envelope and the spindown rate of the neutron star depend on the value of the ratio:

$$\kappa =\frac{\omega R_\mathrm{m}}{V_\mathrm{s}(R_\mathrm{m})},$$

where $`R_\mathrm{m}`$ is the magnetospheric radius, $`\omega =2\pi /P_\mathrm{s}`$ is the neutron star angular velocity and $`V_\mathrm{s}(R_\mathrm{m})`$ is the sound speed at the magnetospheric radius, which according to Eq. (1) is of the order of the free-fall velocity, $`V_{\mathrm{ff}}`$:

$$V_\mathrm{s}(R_\mathrm{m})V_{\mathrm{ff}}(R_\mathrm{m})=\sqrt{\frac{2GM_{\mathrm{ns}}}{R_\mathrm{m}}}.$$
(2)
On this basis Davies et al. (dfp79 (1979)) distinguished three sub-states of the propeller state111According to their classification case ‘a’ corresponds to the pulsar-like spindown, i.e. the ejector state.: (b) very rapid rotator ($`\kappa 1`$); (c) supersonic propeller ($`\kappa 1`$) and (d) subsonic propeller ($`\kappa <1`$).

In cases ‘b’ and ‘c’ the magnetospheric radius of the neutron star exceeds its corotational radius,

$$R_{\mathrm{cor}}=\left(\frac{GM_{\mathrm{ns}}}{\omega ^2}\right)^{1/3}.$$

This means that the neutron star in these cases is in the centrifugal inhibition regime and hence, the effective accretion onto its surface is impossible.

In the case ‘d’ the magnetospheric radius is smaller than the corotational radius. In this situation the plasma being penetrated from the base of the envelope into the star magnetic field is able to flow along the magnetic field lines and to accrete onto the star surface. However the effective plasma penetration into the magnetosphere does not occur, if the magnetospheric boundary is interchange stable. According to Arons & Lea (al76 (1976)) and Elsner & Lamb (el76 (1976)) the onset condition for the interchange instability of the magnetospheric boundary reads

$$T<T_{\mathrm{cr}}0.3T_{\mathrm{ff}}.$$

This means that the neutron star can change its state from the subsonic propeller to accretor only when the cooling of the envelope plasma (due to the radiation and convective motions) dominates the energy input to the envelope due to the propeller action by the neutron star.

Investigating this particular situation Davies & Pringle (dp81 (1981), hereafter DP81) have shown that the energy input to the envelope dominates the energy losses until the spin period of the star reaches the break period, $`P_{\mathrm{br}}`$. Assuming the following values of the neutron star parameters: the magnetic moment $`\mu =10^{30}\mu _{30}\mathrm{G}\mathrm{cm}^3`$ and the mass $`m=1(M_{\mathrm{ns}}/M_{\mathrm{}})`$, and putting the strength of the stellar wind (in terms of the maximum possible accretion rate) $`\dot{M}_0=10^{15}M_{15}\mathrm{g}\mathrm{s}^1`$ they estimated the value of the break period as 60 s.

However, putting the same values of parameters and following the same method of calculations I found the value of $`P_{\mathrm{br}}`$ to be of the order of 450 s, i.e. by a factor of 7.5 larger than that previously estimated in DP81. In this letter I present the calculations and show that this result forced us to change some basic conclusions about the origin of the long periods X-ray pulsars.

## 2 Break period

According to the picture presented by Davies & Pringle the magnetosphere of the neutron star in the state of subsonic propeller is surrounded by the adiabatic ($`pR^{5/2}`$) spherically symmetrical plasma envelope. Until the energy input to the envelope dominates the energy losses the temperature of the envelope plasma is of the order of the free-fall temperature, $`TT_{\mathrm{ff}}`$ and, correspondingly, the sound speed is of the order of the free-fall velocity, $`V_\mathrm{s}V_{\mathrm{ff}}`$. Under this condition the height of the homogeneous atmosphere through out the envelope is comparable to $`R`$ and thus, the envelope is extended from the magnetospheric radius,

$$R_\mathrm{m}r_\mathrm{m}\left(\frac{\mu ^2}{\dot{M}_0\sqrt{2GM_{\mathrm{ns}}}}\right)^{2/7}=$$
(3)
$$=1.2\times 10^9\mathrm{cm}\mu _{30}^{4/7}\dot{M}_{15}^{2/7}m^{1/7},$$

up to the accretion radius of the neutron star,

$$R_\alpha \frac{2GM_{\mathrm{ns}}}{V_{\mathrm{rel}}^2}=2.7\times 10^{10}\mathrm{cm}mV_8^2.$$

Here, $`V_8`$ is the relative velocity between the neutron star and the stellar wind plasma, $`V_{\mathrm{rel}}`$, expressed in units of $`10^8\mathrm{cm}\mathrm{s}^1`$ and $`\dot{M}_0`$ is strength of the stellar wind which is expressed following DP81 in terms of the maximum possible accretion rate:

$$\dot{M}_0=\pi R_\alpha ^2\rho _{\mathrm{}}V_{\mathrm{rel}},$$

where $`\rho _{\mathrm{}}`$ is the density of the stellar wind plasma. Within the considered picture the envelope is quasi-static, i.e. the mass flux through the envelope is almost zero. In this situation the physical meaning of the parameter $`\dot{M}_0`$ is the rate by which the stellar wind plasma overflows the outer edge of the envelope compressing the envelope plasma.

The interaction between the fast rotating magnetosphere and the base of the envelope leads to the turbulization of the envelope plasma. The velocity of the convection motions at the magnetospheric boundary is obviously limited as

$$V_\mathrm{t}\omega R_\mathrm{m}.$$
(4)
Under the condition $`R_\mathrm{m}<R_{\mathrm{cor}}`$ the maximum velocity of the convective motions is smaller than $`V_\mathrm{s}`$ and hence, the Mach number at the base of the envelope is $`M_{\mathrm{Mach}}=V_\mathrm{t}/V_\mathrm{s}<1`$.

The rate of energy loss by the neutron star and, correspondingly, the energy input to the envelope in this case can be expressed as (see Eqs. 3.2.1 and 3.2.2 in DP81)

$$L_\mathrm{d}=2\pi R_\mathrm{m}^2\rho V_\mathrm{t}^3.$$

As it has been argued by Davies & Pringle the cooling of the envelope plasma is determined by the combination of convective motions and bremsstralung radiation. In order to evaluate the energy balance between the heating and cooling processes they introduced the convective efficiency parameter (see also Cox & Guili cc68 (1968)):

$$\mathrm{\Gamma }=\frac{\mathrm{Excess}\mathrm{heat}\mathrm{content}\mathrm{of}\mathrm{convective}\mathrm{blob}}{\mathrm{Energy}\mathrm{radiated}\mathrm{in}\mathrm{the}\mathrm{lifetime}\mathrm{of}\mathrm{a}\mathrm{blob}}.$$

Under the conditions of interest this parameter can be expressed as (for discussion see DP81, page 221)

$$\mathrm{\Gamma }=M_{\mathrm{Mach}}^2\left[\frac{V_\mathrm{t}t_{\mathrm{br}}}{R}\right]=\frac{V_\mathrm{t}^3t_{\mathrm{br}}}{V_\mathrm{s}^2R},$$
(5)
where $`t_{\mathrm{br}}`$ is the bremsstralung cooling time:

$$t_{\mathrm{br}}=6.3\times 10^4\left[\frac{T}{10^9\mathrm{K}}\right]^{1/2}\left[\frac{n}{10^{11}\mathrm{cm}^3}\right]^1\mathrm{s}.$$
(6)
Here, $`n`$ is the number density of the envelope plasma which at the base of the envelope can be evaluated as

$$n(R_\mathrm{m})=\frac{\mu ^2}{4\pi kT_{\mathrm{ff}}(R_\mathrm{m})R_\mathrm{m}^6}.$$
(7)
As it has been shown in DP81 the cooling of the envelope during the subsonic propeller state occurs first at its inner radius. Thus, the energy input to the envelope due to the neutron star propeller action dominates the radiative losses if $`\mathrm{\Gamma }(R_\mathrm{m})1`$.

Combining Eqs. (17) I find the condition $`\mathrm{\Gamma }(R_\mathrm{m})1`$ to be satisfied if the spin period of the neutron star is $`P_\mathrm{s}P_{\mathrm{br}}`$, where the break period is

$$P_{\mathrm{br}}450\mu _{30}^{16/21}\dot{M}_{15}^{5/7}m^{4/21}\mathrm{s}.$$
(8)
This value of the break period exceeds the value of $`P_{\mathrm{br}}`$ presented in DP81 by a factor of 7.5 (see Eq. 4.8 in DP81). Hence the natural question about the reason of this inconsistency arises. One of the most possible reasons is that Davies & Pringle have mistakenly used the value of the magnetospheric radius: $`4.4\times 10^8\mu _{30}^{4/7}\dot{M}_{15}^{2/7}m^{1/7}`$ cm (see Eq. 3.2.3 of their paper) instead of the correct value which is expressed in their paper by Eq. (2.5), i.e. $`10^9\mu _{30}^{4/7}\dot{M}_{15}^{2/7}m^{1/7}`$ cm. Taking into account that $`P_{\mathrm{br}}(R_\mathrm{m})R_\mathrm{m}^{5/2}`$ one finds that the correct value of $`P_{\mathrm{br}}`$ should be larger than that derived in DP81 by a factor of 7.7, i.e. very close to the value of the break period obtained in this letter.

## 3 Discussion

One of the main astrophysical reasons for the investigation of the spindown of neutron stars is the existence of X-ray sources which display pulses with long periods (in excess of 100 s). On the basis of their calculations Davies & Pringle suggested that the periods of neutron stars spinning down due to propeller mechanism can be as long as 100 s only if the stars are situated in the weak stellar wind, i.e. $`\dot{M}_0<4\times 10^{14}\mathrm{g}\mathrm{s}^1`$. They also pointed out that in this case however it is difficult to account for a substantial population of long period pulsators.

In the light of the recalculated value of the break period obtained in this paper (Eq. 8) I find that the propeller mechanism can be responsible for the long spin period of a neutron star even if it is situated in the essentially stronger stellar wind:

$$\dot{M}_08\times 10^{15}\mu _{30}^{16/15}m^{4/15}P_{100}^{7/5}\mathrm{g}\mathrm{s}^1,$$
(9)
where $`P_{100}`$ is the observed spin period of the neutron star expressed in units of 100 s. The corresponding spindown time-scale of the neutron star in the state of subsonic propeller is

$$\tau _\mathrm{d}10^5\mu _{30}^2mI_{45}P_{100}\mathrm{yr},$$
(10)
i.e. smaller that the characteristic evolutionary time-scale of the early spectral type supergiants. Here $`I_{45}`$ is the moment of inertia of the neutron star expressed in units of $`10^{45}\mathrm{g}\mathrm{cm}^2`$.

## 4 Conclusion

The value of the break period at which the spinning down neutron star changes its state from the subsonic propeller to accretor obtained by Davies & Pringle (dp81 (1981)) is underestimated by a factor of 7.5. The incorporation of the re-estimated value of the break period into the spindown scenario suggested by Davies & Pringle shows that the propeller mechanism can be responsible for the origin of the long period X-ray pulsators even if the strength of the stellar wind, in which a neutron star is situated, is in excess of $`10^{15}\mathrm{g}\mathrm{s}^1`$. The analysis of the spin evolution of a neutron star situated in the strong stellar wind will be presented in a forthcoming paper.

###### Acknowledgements.

I acknowledge the support of the Follow-up program of the Alexander von Humboldt Foundation. The work was partly supported by the Federal program “INTEGRATION” under the grant KO 232.


text:
# Vortex glass transition in a random pinning model

\[

## Abstract

We study the vortex glass transition in disordered high temperature superconductors using Monte Carlo simulations. We use a random pinning model with strong point-correlated quenched disorder, a net applied magnetic field, longrange vortex interactions, and periodic boundary conditions. From a finite size scaling study of the helicity modulus, the RMS current, and the resistivity, we obtain critical exponents at the phase transition. The new exponents differ substantially from those of the gauge glass model, but are close to those of the pure three-dimensional XY model.

\]

The magnetic field-temperature phase diagram for vortices in disordered high temperature superconductors has been the focus of considerable interest in the past few years, centered around e.g. the suggestion of M. P. A. Fisher of a possible vortex glass phase with vanishing linear resistance . Numerical simulation is often an important tool in investigations of phase transitions in vortex systems. For example, recent simulations have provided valuable information about the phase diagram for the case of weak disorder or weak magnetic fields . Here evidence for a first order transition separating a Bragg glass phase , i.e. a dislocation-free solid with algebraically decaying translational order, from a vortex liquid phase was obtained. However, in the opposite limit of strong disorder or strong magnetic fields, the existence of a vortex glass phase, i.e. a vortex solid which is topologically disordered in terms of frozen in dislocations, and the critical properties of the vortex glass transition, have not been settled. In this paper we study these issues by Monte Carlo simulations and finite size scaling.

Simulations of disordered superconductors have often used a so-called gauge glass model , where the disorder is modeled as a random vector potential added to the phase difference of the superconducting order parameter. This disorder corresponds to spatially random magnetic flux present in the system. Interestingly, recent simulation results for the gauge glass  gives similar values for the critical exponents as obtained in certain recent experiments where a vortex glass has been reported . However, the gauge glass model has two features that are not particularly realistic, namely, it does not assume any pinning directly affecting the vortex core energy, and, secondly, it is completely isotropic and does not contain any net magnetic field. It is at present unknown if these details are important for the universality class of the glass transition, i.e. if they modify the critical exponents. In particular, there is a possibility of anisotropic scaling produced by the net field, such that the correlation length diverges with different exponents along and perpendicular to the field. Some of these and other issues have been addressed in the literature.

An XY model with a random coupling constant and a net magnetic field has been simulated using open boundary conditions . However, periodic boundary conditions are preferable when bulk properties are studied, to avoid any influence of the sample surface. Also the effect of screening of the vortex interaction has been considered. Gauge field fluctuations lead to screening of the vortex interaction. This screening is usually rather weak in the high temperature superconductors. In models for the strong screening limit, simulations give evidence for the absence of a vortex glass phase at finite temperatures . A vortex glass scenario without any thermodynamic phase transition has also been suggested .

In this paper we consider a random pinning model that contains all the pieces necessary to describe the static universal critical properties of the vortex glass critical point: long-range vortex interactions, strong vortex pinning, a net applied magnetic field, and periodic boundary conditions. The vortex-vortex interaction is a full 3D longrange interaction without screening, that becomes applicable when the bare screening length is much longer than the vortex spacing, i.e. in the strong field limit. The vortex pinning corresponds to uncorrelated quenched point disorder, implemented as a position dependent core energy. This is equivalent to random-$`T_c`$ disorder. The vortex-vortex interaction and the disorder are isotropic in the model, but the applied net magnetic field breaks the spatial symmetry. We include the possibility of anisotropic scaling by allowing for different correlation length exponents in the directions parallel and perpendicular to the magnetic field. The dynamic universality class assumed here is that of relaxation dynamics of the vortex lines, and we stress that there are more dynamic universality classes possible .

The random pinning model is defined by the Hamiltonian

$$H=\frac{1}{2}\underset{𝐫,𝐫^{}}{}V(𝐫𝐫^{})𝐪(𝐫)𝐪(𝐫^{})+\frac{1}{2}\underset{𝐫,\mu }{}D_\mu (𝐫)q_\mu (𝐫)^2$$
(1)
The model is defined on a simple cubic lattice with $`\mathrm{\Omega }=L\times L\times L_z`$ sites, using periodic boundary conditions in all three directions. The vortex line variables are specified by an integer vector field $`𝐪(𝐫)`$, whose $`\mu =x,y,z`$ component is the vorticity on the link from $`𝐫`$ to $`𝐫+𝐞_\mu `$. The partition function is $`Z=\mathrm{Tr}\mathrm{exp}(H/T)`$, where $`T`$ is the temperature, and $`\mathrm{Tr}`$ denotes the sum over all possible integers $`q_\mu `$, subject to the constraint $`𝐪=0`$ on all sites, i.e., the vortex lines have no open ends. An applied net magnetic field is included as a fixed number of vortex lines penetrating the system in the $`z`$-direction. Three different fillings of the applied magnetic field are considered here, i.e., number of vortex lines per link in the $`z`$-direction: $`f=1/2,1/4,1/\sqrt{10}`$. For the irrational filling we used the integer number of vortex lines closest to $`fL^2`$. The figures below are for $`f=1/2`$. The long range vortex-vortex interaction is given by

$$V(𝐫)=\frac{K}{\mathrm{\Omega }}\underset{𝐤}{}\frac{e^{i𝐤𝐫}}{_\mu (22\mathrm{cos}k_\mu )}$$
(2)
where $`K=4\pi ^2J`$ (we set $`J=1`$). On each link in the system is a short-range point-correlated random pinning energy with a uniform distribution in the interval $`0D_\mu (𝐫)K`$. Hence the lattice constant of the discretization lattice in our model corresponds physically to a characteristic length scale for variations in the disorder energy landscape.

The Monte Carlo (MC) trial moves are attempts to insert vortex loops with random orientation on randomly selected plaquettes of the lattice. A MC sweep consists of one attempt on average to insert a loop on every plaquette. The attempts are accepted with probability $`1/(1+\mathrm{exp}\mathrm{\Delta }E/T)`$, where $`\mathrm{\Delta }E`$ is the energy change for inserting the loop. For half and quarter fillings the initial vortex configurations consist of straight lines along the $`z`$-direction, arranged in a regular lattice in the $`xy`$-plane. For filling $`f=1/\sqrt{10}`$ the initial configuration has straight lines placed at random. For equilibration about $`10^5`$ MC sweeps are used, followed by equally many sweeps for collecting data. The results were averaged over up to 2000 samples of the disorder. Thermal averages are denoted by $`\mathrm{}`$ and disorder averages by $`[\mathrm{}]`$. To avoid systematic errors in the calculation of squares of expectation values, two replicas of the system with the same disorder are used.

Superconducting coherence in the vortex line system can be detected by calculating the helicity modulus. One way of doing this  is to add a term $`H_Q=\frac{K}{2\mathrm{\Omega }}𝐐^\mathrm{𝟐}`$ to the Hamiltonian, where $`Q_\mu `$ is the total projected area of vortex loops added during the simulation. The helicity modulus in the direction $`\mu `$ is then given by

$$Y_\mu =1\frac{K}{\mathrm{\Omega }T}[Q_\mu ^2Q_\mu ^2]$$
(3)
and the RMS current density is defined as $`J_\mu =\frac{K}{\mathrm{\Omega }}[Q_\mu ^2]^{1/2}`$. The linear resistivity, $`\rho `$, is obtained from the Kubo formula for the resistance:

$$R_\mu =\frac{1}{2T}\underset{t=t_0}{\overset{t_0}{}}[V_\mu (t)V_\mu (0)]$$
(4)
where $`t`$ is MC time, and the voltage is $`V_\mu \mathrm{\Delta }Q_\mu `$ is the net change in the projected vortex loop area during a sweep. In the calculation of $`\rho `$ the $`H_Q`$ term is not included in $`H`$. The summation time $`t_0`$ is chosen large enough that the resistivity is independent of $`t_0`$, but much shorter than the length of the simulation.

We use a generalization of the Fisher-Fisher-Huse  scaling ansatz to analyze our MC data. At the glass transition temperature $`T_c`$ the correlation length in the $`xy`$-planes, $`\xi `$, and in the $`z`$-direction, $`\xi _z`$, and the correlation time, $`\tau `$, are assumed to diverge as $`\xi |TT_c|^\nu ,\xi _z\xi ^\zeta `$, and $`\tau \xi ^z`$. The anisotropic finite size scaling ansatz  for the helicity modulus is

$`Y_x`$ $`=`$ $`L^{3d\zeta }f_x(L^{1/\nu }(TT_c),L_z/L^\zeta )`$ (5)
$`Y_z`$ $`=`$ $`L^{1d+\zeta }f_z(L^{1/\nu }(TT_c),L_z/L^\zeta ),`$ (6)
where $`d=3`$ is the spatial dimensionality and $`f_\mu `$ are scaling functions (scaling functions will from now on be suppressed in the equations). The current density scales as $`J_xL^{2d\zeta },J_zL^{1d}`$. The linear resistivity $`\rho =E/J`$, where $`E`$ is the electric field, obeys

$$\rho _xL^{d3+\zeta z},\rho _zL^{d1\zeta z},$$
(7)
We did a number of tests of equilibration of our simulations. We followed the standard procedure of calculating the “Hamming distance” between two replicas with identical disorders , for some selected system sizes $`L,L_z`$. We also increased the number of sweeps for equilibration with a factor of 10 for a few selected parameter values, and obtained no deviations from the data in the figures below. The linear resistivity (shown below in Fig. 4) gives an estimate of the relaxation time $`t_0`$, which is the time where the curves saturate. The values for $`t_0`$ gives a rough estimate of the required equilibration time. The equilibration times used in our simulations are $`8t_0`$.

The first step is to verify that the model has a vortex glass phase, instead of a Bragg glass phase  that is expected for low fillings and weak disorder. Figure 1 shows a typical snapshot of a sample configuration for $`T=0.5(T_c)`$, i.e. deep in the glass phase. This calculation was done using an exchange MC algorithm  with 9 uniformly spaced temperatures in $`[0.5,4.5]`$. We also computed the structure function $`S(𝐪)`$, and obtain no essential difference between $`S`$ in the vortex liquid phase and the glass phase, with no indication of a Bragg glass phase . This demonstrates that the low temperature phase of our model, for the filling and disorder strength considered, is indeed a glass where the vortex lines are frozen in random positions.

To locate the critical temperature of the vortex glass to liquid transition, and determine the critical exponents, MC data for the helicity modulus in the $`x`$ and $`z`$ directions is analyzed by the finite size scaling form in Eqs. (5),(6). MC data for different $`T,L,L_z`$ for $`f=1/2`$ is plotted in Fig. 2. To determine the anisotropy exponent $`\zeta `$, a sequence of system sizes $`L_z`$ are examined for each $`L`$. In the relation $`L_z=cL^\zeta `$, both $`c`$ and $`\zeta `$ are varied until data curves for $`Y_x`$ and $`Y_z`$ for different system sizes $`L=6,8,10,12`$ become approximately independent of $`L`$ at a common temperature, which is our estimate for $`T_c`$. Fig. 2 shows the best crossing obtained in this procedure, which gives $`T_c=4.5\pm 0.1,\zeta =1\pm 0.1`$. The same results are obtained from the RMS currents (inset). Equally good fits are obtained for $`c`$ in the interval $`1.5<c<2`$. The errors are estimated as the interval outside which considerably poorer scaling is obtained.

To determine the correlation length exponent $`\nu `$ we use fits to Eqs. (5),(6) to obtain a data collapse for system sizes $`L=612`$ over an entire temperature interval around $`T_c`$. As a measure of the quality of the collapse we define the RMS fit error $`\mathrm{\Delta }=\left[_{L,T,\mu }(LY_\mu (L,T)f_\mu (x))^2\right]^{1/2}`$, where $`x=L^{1/\nu }(TT_c)`$, and $`f_\mu `$ is estimated by e.g. a cubic polynomial fit to the MC data. We did several different types of fits, all giving similar results. Figure 3 shows the best data collapse in a two parameter fit where both $`T_c`$ and $`\nu `$ are varied independently in the $`x`$ and $`z`$-directions. The best fit is obtained for $`T_c=4.5\pm 0.1,\nu =0.7\pm 0.1`$. A good data collapse is obtained, except for the smallest system size, where deviations are obtained below $`T_c`$. However for larger sizes scaling gets better, suggesting that the deviation for small system size is a finite size effect. The inset shows the RMS fit error as a function of $`\nu `$ for fixed $`T_c=4.5`$. Nearly identical results are obtained by instead analyzing data for the RMS currents. The same result for $`\nu `$ is also obtained from an analysis of the derivative of the helicity moduli wrt $`T`$. For the lower fillings, $`f=1/4,1/\sqrt{10}`$, we find similar exponents but with larger error bars.

So far we have presented results for static quantities, and we now consider dynamics. The dynamic critical exponent is obtained from the linear resistivity in Eq. (4), for $`f=1/2`$, $`T=T_c=4.5,L_z=2L`$. Fig. 4 shows finite size scaling data collapses according to Eq. (7) of MC data for $`\rho _x,\rho _z`$ vs. the total MC integration time $`t_0`$ in the Kubo formula. The inset shows the RMS error in a power law fit to the MC data curves, and the best fit is obtained for $`z=1.5\pm 0.2`$. This gives a resistivity exponent in $`\rho t^s`$ of $`s=\nu (z1)0.3`$.

We will now discuss the values obtained for the critical exponents of the random pinning model: $`\zeta =1\pm 0.1,\nu =0.7\pm 0.1,z=1.5\pm 0.2`$. Our correlation length exponent $`\nu `$ is close to the limiting value at a disordered fixed point allowed by the inequality : $`\nu 2/d`$. Within the precision obtained, we can not distinguish this exponent from that of the pure zero field 3D XY model, i.e. $`\nu 0.67`$. Also our dynamic exponent $`z1.5`$ is consistent with the zero field 3D XY model with MC dynamics for vortex loops . If confirmed, these results suggest a scenario where the vortex glass transition in the random pinning model for high fillings and strong disorder belongs to the zero field 3D XY universality class. In particular, this also suggests that the glass transition is driven by similar effectively isotropic, closed vortex loop fluctuations as in the zero field case, on top of a glassy groundstate.

Finally we compare our results with other models for the vortex glass transition and with experiments. The critical exponents obtained here for the random pinning model differ considerably from those of the gauge glass model : $`\nu 1.39,z4.2,s=\nu (z1)4.5`$, and also from a random coupling 3D XY model with an applied field and open boundary conditions : $`\nu 2.2,z3.3,s5.3`$. Hence, these models appear to belong to different universality classes than the random pinning model. Experiments on (K,Ba)BiO3 give $`s3.9`$, and experiments on untwinned proton irradiated YBa2Cu3O7-δ give $`s5.3`$ . These experiments show consistency with expected vortex glass behavior for tilting the magnetic field away from the $`c`$-axis. Unexpectedly, the exponents disagree considerably with our results. Further work is needed in order to clarify the reasons for this discrepancy.

In summary, we have observed a finite temperature continuous vortex glass transition by simulations and finite size scaling analysis of a three-dimensional random pinning model. The critical exponents, $`\zeta 1,\nu 0.7,z1.5`$, are surprisingly close to those of the zero-field pure 3D XY model, but disagree with the gauge glass model and with some experiments. These results motivate further theoretical and experimental work.

We acknowledge very useful discussions with Pierre Le Doussal, Peter Olsson, Stephen Teitel, and Peter Young. This work was supported by the Swedish Natural Science Research Council, STINT, and PDC.


text:
# 1 Introduction

## 1 Introduction

It is generally accepted that under extreme conditions (high density and /or temperature) hadronic matter undergoes two phase transitions, characterized respectively by deconfinement and restoration of chiral symmetry. The QCD vacuum would then be described by a weakly interacting gas of massless quarks and gluons. Major theoretical and experimental efforts have been dedicated to heavy-ion physics looking for signatures of the quark gluon plasma . Special attention has also been given to neutron stars, which provide a natural laboratory to study hadronic matter at high densities.

Theoretical concepts on the behavior of matter at high densities have been mainly developed from model calculations. The Nambu-Jona-Lasinio \[NJL\] type models have been extensively used over the past years to describe low energy features of hadrons and also to investigate restoration of chiral symmetry with temperature or density .

Studies of flavor asymmetric matter at high densities are also motivated by the interest in investigating the behavior of mesons. As a matter of fact  the charge multiplets of mesons, that are degenerated in vacuum or in symmetric matter ($`\rho _u=\rho _d=\rho _s`$) are expected to have a splitting in flavor asymmetric matter. In particular, the masses of kaons (antikaons) would increase (decrease) with density and a similar effect would occur for $`\pi ^{}`$ and $`\pi ^+`$ in neutron matter.

A slight raising of the $`K^+`$ mass and a lowering of the $`K^{}`$ mass  seems to be compatible with the analysis of data on kaonic atoms  and with the results of KaoS and FOPI collaborations at GSI . From the theoretical point of view, the driving mechanism for the mass splitting is attributed mainly to the selective effects of the Pauli principle, although, in the case of $`K^{}`$, the interaction with the $`\mathrm{\Lambda }(1405)`$ resonance plays a significant role as well, at least as low densities are concerned .

A model aiming at describing hadronic behavior in the medium should account for the great variety particle-hole excitations that the medium can exhibit, some of them with the same quantum numbers as the hadrons under study . It has been established, within the framework of NJL models, the presence, in flavor asymmetric media, of low energy pseudoscalar modes, which are excitations of the Fermi sea . Such studies were carried out in quark matter simulating nuclear matter $`(\rho _u=\rho _d,\rho _s=0)`$, for charged kaons, and neutron matter without beta equilibrium, for charged pions. The role of the ’t Hooft interaction was not taken into account in the last case. The combined effect of density and temperature, as well as the effect of vector interaction, was discussed for the case of charged kaons. . Only densities below $`3\rho _0`$ were considered.

The aim of this work is to explore the high density region ($`\rho >\rho _c=2.25\rho _0`$) of quark matter simulating neutron matter in weak equilibrium, having in mind the study of the behavior of kaonic (charged and neutral) and pionic (charged) excitations. The work is carried out in the framework of the SU(3) NJL model. Since for high densities the quarks are supposed to be deconfined the lack of confinement of NJL model is not a drawback in this region.

## 2 Formalism

We work in a flavor SU(3) NJL type model with scalar-pseudoscalar pieces and a determinantal term, the ’t Hooft interaction, which breaks the $`U_A(1)`$ symmetry. We use the following Lagrangian:

$$\begin{array}{ccc}\hfill & =& \overline{q}(i\gamma ^\mu _\mu \widehat{m})q+\frac{1}{2}g_S_{a=0}^8[(\overline{q}\lambda ^aq)^2+(\overline{q}i\gamma _5\lambda ^aq)^2]\hfill \\ & +& g_D\{\text{det}[\overline{q}(\mathrm{\hspace{0.17em}1}+\gamma _5)q]+\text{det}[\overline{q}(\mathrm{\hspace{0.17em}1}\gamma _5)q]\}\hfill \end{array}$$
(1)
The model parameters, the bare quark masses $`m_d=m_u,m_s`$, the coupling constants and the cutoff in three-momentum space, $`\mathrm{\Lambda }`$, are fitted to the experimental values of the masses of pseudoscalar mesons and $`f_\pi `$. Here we use the same parameterization as , $`\mathrm{\Lambda }=602.3`$ MeV, $`g_S\mathrm{\Lambda }^2=3.67`$, $`g_D\mathrm{\Lambda }^5=12.39`$, $`m_u=m_d=5.5`$ MeV and $`m_s=140.7`$ MeV.

The six quark interaction can be put in a form suitable to use the bosonization procedure (see ):

$$_𝒟=\frac{1}{6}g_DD_{abc}(\overline{q}\lambda ^cq)[(\overline{q}\lambda ^aq)(\overline{q}\lambda ^bq)3(\overline{q}i\gamma _5\lambda ^aq)(\overline{q}i\gamma _5\lambda ^bq)]$$
(2)
with: $`D_{abc}=d_{abc},a,b,cϵ\{1,2,\mathrm{..8}\},\text{(structure constants of SU(3))},D_{000}=\sqrt{\frac{2}{3}},D_{0ab}=\sqrt{\frac{1}{6}}\delta _{ab}`$.

The usual procedure to obtain a four quark effective interaction from this six quark interaction is to contract one bilinear $`(\overline{q}\lambda _aq)`$. Then, from the two previous equations, an effective Lagrangian is obtained:

$`L_{eff}`$ $`=`$ $`\overline{q}(i\gamma ^\mu _\mu \widehat{m})q`$ (3)
$`+`$ $`S_{ab}[(\overline{q}\lambda ^aq)(\overline{q}\lambda ^bq)]+P_{ab}[(\overline{q}i\gamma _5\lambda ^aq)(\overline{q}i\gamma _5\lambda ^bq)],`$
where:

$`S_{ab}=g_S\delta _{ab}+g_DD_{abc}<\overline{q}\lambda ^cq>,`$
$`P_{ab}=g_S\delta _{ab}g_DD_{abc}<\overline{q}\lambda ^cq>.`$ (4)
By using the usual methods of bosonization one gets the following effective action:

$`I_{eff}=`$ $`i`$ $`Tr\mathrm{ln}(i_\mu \gamma _\mu \widehat{m}+\sigma _a\lambda ^a+i\gamma _5\varphi _a\lambda ^a)`$ (5)
$``$ $`{\displaystyle \frac{1}{2}}(\sigma _aS_{ab}^{}{}_{}{}^{1}\sigma _b+\varphi _aP_{ab}^{}{}_{}{}^{1}\varphi _b).`$
In order to introduce medium effects, we use the thermal Green functions in the standard way  and restrict ourselves to $`T=0`$. The following gap equations are obtained:

$$M_i=m_i2g_S<\overline{q}_iq_i>\mathrm{\hspace{0.17em}2}g_D<\overline{q}_jq_j><\overline{q}_kq_k>$$
(6)
with $`i,j,k`$ cyclic and $`<\overline{q}_iq_i>=4M_iI_{1,i}`$ are the quark condensates.

The poles in the kaon propagator are obtained from the condition $`(\mathrm{\hspace{0.17em}1}K_PJ_{PP})=\mathrm{\hspace{0.17em}0}`$ with:

$$J_{PP}=4(I_{1,u}+I_{1,s}+[\omega ^2(M_uM_s)^2]I_2(\omega )),$$
(7)
where $`K_P=g_S+g_D<\overline{d}d>`$, $`\omega =q_0\pm (\mu _s\mu _u)`$ is the excitation energy, which, for kaons at rest is equal to $`\pm m_{K^\pm }`$ and the integrals are given by:

$$I_{1,i}=\frac{N_c}{4\pi ^2}_{\lambda _i}^\mathrm{\Lambda }\frac{k^2}{E_i}𝑑k,$$
(8)
where $`\lambda _i`$ is the Fermi momentum, and

$`I_2(\omega )`$ $`=`$ $`{\displaystyle \frac{N_c}{4\pi ^2}}\left\{{\displaystyle _{\lambda _u}^\mathrm{\Lambda }}{\displaystyle \frac{k^2}{[(E_u+\omega )^2E_s^2]E_u}}𝑑k+{\displaystyle _{\lambda _s}^\mathrm{\Lambda }}{\displaystyle \frac{k^2}{[(E_s\omega )^2E_u^2]E_s}}𝑑k\right\}.`$ (9)
Similar expressions are obtained for $`K^0,\overline{K}^0`$, by replacing $`ud`$ and for $`\pi ^+,\pi ^{}`$, by replacing $`sd`$.

## 3 Discussion and conclusions

At variance with , we consider here the case of asymmetric quark matter with strange quarks in chemical equilibrium maintained by weak interactions and with charge neutrality, by imposing the following constraints on the chemical potentials of quarks and electrons and on its densities:

$$\mu _d=\mu _s=\mu _u+\mu _e,\text{ and }\frac{2}{3}\rho _u\frac{1}{3}(\rho _d+\rho _s)\rho _e=0,$$
(10)
with

$$\rho _i=\frac{1}{\pi ^2}(\mu _i^2M_i^2)^{3/2}\theta (\mu _i^2M_i^2),\text{ and }\rho _e=\frac{\mu _e^3}{3\pi ^2}.$$
(11)
As discussed by several authors, this version of the NJL model exhibits a first order phase transition . As shown in , by using a convenient parameterization the model may be interpreted as having a mixing phase — droplets of light $`u,d`$ quarks with a density $`\rho _c=2.25\rho _0`$ surrounded by a non-trivial vacuum — and, above the critical density, a quark phase with partially restored chiral symmetry. As a matter of fact, the energy per particle of the quark system has two minimum, corresponding to the zeros of the pressure (see Fig. 1-a)), the minimum at $`\rho _c=2.25\rho _0`$ being an absolute minimum.

The dynamical quark masses and chemical potentials are plotted in Fig. 1-b) as functions of $`\rho _B`$. As shown in , chiral symmetry is partially restored in the SU(2) sector, the $`u\text{and}d`$ quarks masses decreasing sharply. An important point to be noticed is that, due to the ’t Hooft contribution in the gap equations, the mass of the strange quark decreases smoothly and becomes lower than the chemical potential at densities above 3.8$`\rho _0`$, which we will denote from now on as $`\rho _s^{}`$. A more pronounced decrease is then observed, which is due to the presence of $`s`$ quarks in this regime. As it will be discussed below, this fact is related to a change in the behavior of different observables, as compared to the region of lower densities. At $`\rho =\rho _s^{}`$ (Fig. 1-b)) the difference between the chemical potentials $`\mu _d=\mu _s`$ has a maximum and, as the density increases there is a tendency to the restoration of flavor symmetry that can be seen as well in quark densities ( Fig. 1-c)) and mesonic masses (Fig. 2).

Let us analyze the results for the masses of $`K^+,K^{},K^0,\overline{K}^0`$ and for $`\pi ^+,\pi ^{}`$ (Fig. 2 a)-c)). In order to understand these results, it is useful to look for the limits of the Dirac and of the Fermi sea continua of $`\overline{q}q`$ excitations with the quantum numbers of the mesons under study. They can be obtained by inspection of the meson dispersion relations (looking for the limits of the regions of poles in the integrals $`I_2(\omega )`$, Eq. 9) and are plotted in Fig. 2 a) and b) (dashed point lines): $`\omega ^{}=\sqrt{M_s^2+\lambda _s^2}+\sqrt{M_{u(d)}^2+\lambda _s^2}`$ is the lower limit of the Dirac continuum, and $`\omega _{up}=\sqrt{M_s^2+\lambda _s^2}\sqrt{M_{}^{2}{}_{u(d)}{}^{}+\lambda _{}^{2}{}_{s}{}^{}}`$, $`\omega _{low}=\sqrt{M_{}^{2}{}_{s}{}^{}+\lambda _{}^{2}{}_{u(d)}{}^{}}\mu _{u(d)}`$ are the upper and lower limits of the Fermi sea continuum. We do not show the corresponding limits for the pions because, in the range of densities studied, the pion modes remain outside the continuum.

As usual, two kinds of solutions for mesons masses may be found, corresponding to excitations of the Dirac sea modified by the presence of the medium, and to excitations of the Fermi sea, respectively; however there are significant differences in their behavior at low and high densities. Concerning the Dirac sea excitations, we observe the expected splitting between multiplets: the increase of the masses of $`K^+,K^0\text{and}\pi ^{}`$ with respect to those of $`K^{},\overline{K}^0\text{and}\pi ^+`$ is due to Fermi blocking and is more pronounced for kaons than for pions because there are $`u`$ and $`d`$ quarks in the Fermi sea and therefore there are repulsive effects due to the Pauli principle acting on $`\pi ^+`$. One notices also that at the critical density, $`2.25\rho _0`$ the antikaons enter in the Dirac continuum and became resonances, a behavior already observed in . The new findings are that at $`\rho =\rho _s^{}`$ there is a sharp increase of $`\omega ^{}`$ and decrease of $`\omega _{up}`$ and the antikaons turn again in bound states and, while at low densities the splitting increases with density, at high densities it decreases and the modes show a tendency to become degenerated again.

Below the lower limit of the Fermi sea continuum there are low bound states, denoted by the subscript $`S`$, with quantum numbers of $`K^{},\overline{K}^0,\pi ^+`$, respectively, which are particle-hole excitations of the Fermi sea. In the last case the low energy mode only exists up to $`\rho =\rho _c`$.

Analysis of the energy weighted sum rules allows to calculate the strength located in each mode . The sum of the contributions of the three discrete solutions (dashed curves in Fig. 3) deviates from 1 as the density increases, indicating that the modes of the continuum became more relevant with increasing density, a behavior already found in . We can also see from Fig. 3 that above $`\rho _s^{}`$ the lower energy modes are less coupled to the quarks and its strength decrease, which means that they become less relevant and that the opposite occurs for the antikaons that excitations of the Dirac sea. It is a well known feature in many-body physics that the Fermi sea continuum of excitations induces attractive correlations below its lower limit and repulsive above the upper limit . As seen in Fig. 2- a)-b), this continuum becomes very narrow above $`\rho _s^{}`$, which explains the behavior mentioned above.

Finally, we remark that those low energy modes only appear associated to a first order phase transition, as it was shown in . We do not observe kaon condensation, probably due to the simplicity of the model. In fact, the model, in its present form, does not allow to consider color superconductivity, which is though to play a significant role on kaons. On the other side, our treatment is the leading order calculation in the inverse number of colors, $`1/N_C`$, which means that we use the Hartree approximation to calculate the energy of the quark system and RPA to calculate the meson masses. As pointed out by several authors, the quantum fluctuations due to meson-like excitations are expected to be relevant for calculation of the quark masses and certainly will have effects in the EOS. This study can be an useful starting point to investigations along theses lines.

Acknowledgement:

Work supported by FCT, by project PRAXIS/P/FIS/12247/98, by grant SFRH/BD/3296/2000 (Pedro Costa) and GTAE.