text:
# A Bow-Shock Nebula Around a Compact X-ray Source in the Supernova Remnant IC 443

## 1 Introduction

The mixed-morphology Galactic supernova remnant IC 443 (l,b = 189.1°, +3.0°) has been the subject of extensive studies at all wavelength bands (Fesen (1984); Braun & Strom (1986); van Dishoeck, Jansen & Phillips (1993); Asaoka & Aschenbach (1994)). IC 443 is especially well known for its clear interaction with surrounding molecular clouds (Burton et al. (1990); Bocchino & Bykov (2000)). With its large variety of shocked molecular species detected along its SE edge, IC 443 has become a standard laboratory for studying shock chemistry (van Dishoeck et al. 1993 and references therein).

IC 443 is also of interest because it is coincident with the unidentified EGRET source 2EG J0618+2234 (Sturner & Dermer (1995); Esposito et al. (1996)), and thus has stimulated a good deal of theoretical work aiming to explain the production of GeV $`\gamma `$-rays by shell SNRs (Sturner, Dermer & Mattox (1996); Sturner et al. (1997); Gaisser, Protheroe & Stanev (1998); Baring et al. (1999); Bykov et al. (2000)).

In the X-ray band, Petre et al. (1988) found that the bulk of the emission was thermal in origin with a temperature $``$$`10^7`$ K, typical of middle-aged SNRs. However, Wang et al. (1992) found evidence in the Ginga satellite data for a hard X-ray component in the spectrum of IC 443 extending up to 20 keV. Follow-up observations with the ASCA satellite revealed that the bulk of this hard X-ray emission came from a single, unresolved feature located along the edge of the radio shell where the SNR/molecular cloud interaction was strongest (Keohane et al. (1997), hereafter K97). Furthermore, the hard X-ray feature is positionally coincident with a region where the radio spectral index is considerably flatter than the SNR as a whole (i.e $`\alpha <`$ 0.24 vs $`\alpha `$ = 0.42, where $`F_\nu \nu ^\alpha `$) (Green (1986); Kovalenko, Pynzar & Udal’Tsov (1994)).

K97 considered different explanations for the origin of this hard X-ray emission, and favored a model in which synchrotron emission was being produced by TeV electrons in a region of enhanced particle acceleration, resulting from the SNR/molecular cloud interaction. Recent BeppoSAX data analyzed by Bocchino & Bykov (2000) also support the K97 interaction model. An equally viable model, favored by Chevalier (1999, hereafter C99), posits that the hard X-ray feature is synchrotron emission powered by an energetic neutron star. To distinguish between these two competing hypotheses we have undertaken high resolution X-ray and radio observations toward the hard X-ray feature. These new Chandra and VLA observations strongly favor C99’s neutron star hypothesis.

## 2 Observations and Analysis

### 2.1 Chandra X-Ray Observatory

The Chandra X-Ray Observatory performed a short (10 ks) observation of IC 443 on April 10, 2000 during the first cycle of guest observations (A01). The hard X-ray feature was centered on the I3 chip of the Advanced CCD Imaging Spectrometer (ACIS). To produce Chandra event files, we followed the standard CIAO procedures outlined in Version 1.3 of the CIAO Beginner’s Guide for CIAO Release 1.1, using the acisD2000-01gainN0001.fits gain file. From the events file, we extracted high resolution images in both hard ($`E>2.1`$ keV) and soft ($`E<2.1`$ keV) spectral bands. Smoothed versions of these images are presented in Figure 1.

The hard X-ray source is resolved by this high resolution ($``$1″) Chandra observation (Figure 1) as a nebular region of diffuse emission with a cometary tail. Within the nebula there is an unresolved point source that we have designated as CXOU J061705.3+222127. The nebula of hard emission exhibits bow-shock morphology with a width of $``$ 35$`\mathrm{}`$ and a minimum standoff distance of $`r_s`$ $``$ 8.5$`\mathrm{}`$. The X-ray point source lies at the apex of the nebula, located at (epoch 2000) $`\alpha `$ = 06h17m05.31s, $`\delta `$ = 22°21′27″.3, with errors of $`\pm `$2″ in each coordinate. It is especially visible in the soft ($`E<2.1`$ keV) X-ray band (Figure 1).

Likewise, in order to compare the X-ray luminosity of the cometary nebula with statistical studies of pulsar wind nebulae (PWNe) (e.g. Seward & Wang 1988, see §3), we integrated the power-law model of K97 over the Einstein band (0.2-4 keV) – thus deriving L$`{}_{x}{}^{}5\times 10^{33}d^2(1.5\mathrm{kpc})^2ergs^1`$. Inherent in this Lx calculation is the assumption that the cometary nebula is dominated in soft X-rays by the same power-law exibited at hard energies with ASCA. The Chandra flux (over the same range) yields consistent results, albeit with larger uncertainties.

### 2.2 Very Large Array

All radio observations were made on August 26 and December 31, 1997 with the Very Large Array (VLA111The VLA is operated by the National Radio Astronomy Observatory, a facility of the National Science Foundation.) in the C and D arrays, respectively. A log of the observations is summarized in Table 1. The data acquisition and calibration were standard, using J0632+103 as a phase calibrator and 3C 138 as both the flux density and polarization angle calibrator.

The cometary and bow-shock morphologies of the extended source are even more pronounced at radio wavelengths, extending some 2′ to the northeast of the soft X-ray source. Although we display only the $`\lambda `$3.5 cm image in Figure 2, similar structure is visible at $`\lambda `$6 cm and $`\lambda `$20 cm. The peak in the radio emission lies about 6″ to the northeast of the compact X-ray source. There is no radio point source ($`>2`$ mJy) coincident with the X-ray point source. The distribution of linearly polarized emission is also shown in Figure 2. The direction of the E-vectors are generally parallel to the shock normal, as would be expected for evolved SNRs with circumferential magnetic fields. The character of the magnetic field evidently changes substantially in the area dominated by the X-ray emission. Here the magnetic field direction appears to wrap around the “head” of the hard nebula, tracing a bow-shock morphology. The average degree of polarization in this region is 8%. This is likely a lower limit since, owing to the compact size of the nebula, there is likely some beam depolarization occurring. To the northeast, where the field is more ordered, the percentage polarization is uncharacteristically strong for an evolved SNR, exceeding 25% in several locations.

The integrated flux density over the entire length of the cometary nebula is given for three frequencies in Table 1. An approximate value of 230$`\pm `$200 mJy was also obtained from the 327 MHz image presented by Claussen et al. (1997). The increased error in estimating the flux density at low frequencies is due to the uncertainty in subtracting out the background shell emission from the SNR. Within the errors, the radio spectrum is well represented by a flat spectral index $`\alpha _R`$ $``$ 0.0 and a mean flux density of 206 mJy. In contrast, K97 measured an X-ray spectral index of $`\alpha _X`$1.3$`\pm `$0.2, which implies a break in the spectrum near 4$`\times `$1013 Hz. The radio luminosity of the cometary nebula, determined by integrating the spectrum from 10 MHz to 100 GHz, is LR=5.5$`\times `$10$`{}_{}{}^{31}(d/1.5\mathrm{kpc})_{}^{2}`$ erg s-1, or less than 1% of LR for IC 443 as a whole.

### 2.3 Spectra

In addition to the spectral properties of the hard nebula, a spectrum of CXOU J061705.3+222127 was extracted from this short Chandra observation. The surrounding hard nebula was subtracted as background (i.e. an annulus with an inner and outer radius of 4.3″ and 23.6″ respectively). These spectra were binned by a factor of 256 (PHA channels), and the softest three ACIS-I channels were ignored. As discussed above, the point source is significantly softer than the surrounding hard nebula, implying a likely thermal origin for the emission. For this reason, we modeled the spectrum as black body radiation, though this fit was not statistically unique due to large uncertainties in the small data set.

The best-fit ($`\chi _\nu ^21.1`$) black body temperature and luminosity are kT=0.71$`\pm `$0.08 keV and L=$`(6.5\pm 0.9)\times 10^{31}`$ erg s-1 (i.e. A=$`0.025\pm 0.003`$($`d1.5\mathrm{kpc})^2`$ km2). The 1-5 keV X-ray flux of the point source is $`2\times 10^{13}`$ erg s-1 cm-2. This fits contained two free parameters by assuming K97’s best-fit column density (NH=1.3$`\times 10^{21}`$cm-2).

## 3 Discussion

The soft X-ray point source in IC 443 is one of a growing number of young neutron star candidates associated with supernova remnants (Helfand (1998)). These neutron stars have a diverse range of properties, suggesting that they do not originate from a homogeneous population (Chakrabarty et al. (2000)). The long-period (8-10 s) and claimed high-field (1014-1015 G) of anomalous X-ray pulsars (AXP) and soft gamma-ray repeaters (SGR), differ from the canonical young radio pulsars with periods P$``$15-500 ms and field strengths of B$``$1012-1013 G. The radio quiet neutron stars (RQNs), with their thermal X-ray emission and lack of radio pulsations, may be altogether different from either of these two classes of objects (Pavlov et al. (2000)).

In the absence of detectable pulsations, canonical radio pulsars can be distinguished from AXPs, SGRs, and RQNs by the presence of an extended synchrotron nebula, powered by the energetic relativistic wind (Gaensler, Bock & Stappers (2000)). It is argued that the cometary nebula around the X-ray point source in IC 443 is most likely such a PWN. This nebula has all the expected observational characteristics (C99; Gaensler et al. (2000)): flat spectrum radio emission ($`\alpha _R`$ = 0.1 to 0.3), a steep X-ray spectrum ($`\alpha _X`$ = 1.0 to 1.5) and a high degree of linear polarization ($`>5\%`$). Morphologically it closely resembles PWN detected toward the SNRs W 44 and G 5.4$`1.2`$ (Frail & Kulkarni (1991); Frail et al. (1996)), both of which contain active pulsars. It has been argued (Frail & Kulkarni (1991); Frail et al. (1996)) that these nebulae have developed bow shocks as a result of ram pressure confinement of the pulsar wind, due to the high space velocity of their pulsar through the surrounding medium. We assert that the most tenable explanation for the X-ray point source in IC 443 is that it is a young pulsar which has traveled (ballistically) from its birth place to its present location, near the edge of the decelerating SNR.

Accepting this hypothesis we proceed with inferring the salient physical properties of the neutron star and corresponding PWN. Constraints can be placed on the velocity of the neutron star by a variety of methods. The offset of the object from its origin gives a measure of its transverse velocity VNS. Measuring this offset in the case of IC 443 is made difficult owing to the uncertain location of the explosion center (Asaoka & Aschenbach (1994); Claussen et al. (1997)). We note that the neutron star has traveled a fraction $`\beta `$0.9 of this projected distance (approximated from the distance between the soft point source and the shock), giving V$`{}_{\mathrm{NS}}{}^{}=\beta `$r$`{}_{s}{}^{}/t_s`$. Likewise, the shock velocity of the SNR can be similarly expressed as V$`{}_{s}{}^{}=c_{}`$r$`{}_{s}{}^{}/t_s`$, where $`c_{}`$ is a constant equal to 2/5 for a remnant in the Sedov stage, while C99 argues for a value of $`c_{}=3/10`$ for IC 443 based off of the parameters of his model. Combining both equations for VNS and Vs together, and using the value of V$`{}_{s}{}^{}=100`$ km s-1 favored by C99, gives V$`{}_{\mathrm{NS}}{}^{}=V_s\beta c_{}^{}{}_{}{}^{1}`$=225-300 km s-1, a result that is independent of distance and only weakly dependent on the evolutionary state of the SNR (i.e. $`c_{}`$).

It should be noted that the tail of the hard nebula does not point towards the geometric center of the SNR. However, the blast center and geometric center of an SNR can be quite different in the presence of large-scale density gradients (Dohm-Palmer & Jones (1996); Hnatyk & Petruk (1999)). Such a density gradient, combined with the potential of a cross-wind, could possibly account for the apparent discrepancy.

We are able to make a second approximation of VNS with respect to the local medium by means of the nebula’s bow shock morphology. K97 fit the SNR with an X-ray temperature kT $``$ 1 keV, implying a sound speed c$`{}_{local}{}^{}`$100 km s-1. An approximation of the bow shock angle of 90°-120° gives a minimum Mach number M$``$1.1-1.5 through a standard geometric calculation. A fundamental determination of the speed of sound in this particular region of the SNR gives velocity with respect to the local medium of V$`_{\mathrm{NS}\mathrm{local}}`$=150 M1.5 kT$`{}_{}{}^{0.5}{}_{1\mathrm{keV}}{}^{}`$ km s-1. Accounting for radial expansion of local plasma ($``$100 km s-1, this velocity could be as high as 250 km s-1. Both methods of determining VNS are in good agreement with a model by C99, indicating a remnant age on the order of 30,000 years (i.e. $`\tau _{30}`$1).

The energy source for the synchrotron nebula must ultimately be derived from the particles and field generated by the compact central neutron star and therefore some fraction of the its spindown power Ė is powering this emission. Empirical relations have been derived for known PWN between Ė and the observed X-ray and radio luminosities (Lx and LR) that enable us to estimate the rotation period P and the magnetic field B of the compact object.

For the value of Lx derived in §2 we obtain Ė$``$2$`\times `$$`10^{36}`$ erg s-1 (Seward and Wang 1988). Frail & Scharringhausen (1997) and Gaensler et al. (2000) note that Ė$``$104LR for known radio PWN, which yields a value of Ė$``$6$`\times `$$`10^{35}`$, so in the discussion that follows we adopt Ė=1036 ergs s-1 (i.e. Ė36=1). Equating the pressure of the wind Ė/$`4\pi cr_w^2`$, to the ram pressure due to the neutron star’s motion through the surrounding gas $`\rho `$ V$`{}_{}{}^{2}{}_{\mathrm{PSR}}{}^{}`$, we derive a density of $``$0.1 cm-3 for our adopted values of Ė and VNS. This density is more typical of that expected for the hot X-ray interior of the SNR, and is orders of magnitude below that expected in the dense molecular ring against which the NS and its nebula are projected against (see van Dishoeck et al. 1993).

The spindown luminosity of a pulsar is related to its current period P and dipolar surface field B by Ė $``$ B2P-4, while the characteristic age of the pulsar is expressed as $`\tau _c`$$``$ B-2P2. It is straightforward to show that for a pulsar born spinning rapidly (i.e. P$`<<`$P), which loses energy predominately via magnetic dipole radiation, P=145 ms ($`\tau _{30}`$Ė$`{}_{36}{}^{})^{1/2}`$ and B=3.3$`\times `$1012 G $`\tau _{30}^1`$Ė$`{}_{}{}^{1/2}{}_{36}{}^{}`$. Thus, with our nominal choice of parameters, the period and magnetic field closely match those of other young pulsars (e.g. PSR B1757$`24`$; Frail & Kulkarni (1991)).

In the light of this result, it is worth re-examining the nature of the unidentified EGRET source 2EG J0618+2234 (l,b=189.13°,+3.19°) (Sturner & Dermer (1995); Esposito et al. (1996)) toward IC 443. Nel et al. (1996) note that all known gamma-ray pulsars have a large (greater than 0.5) ratio of $`\dot{E}_{33}/d_{kpc}^{}{}_{}{}^{2}`$. For our derived $`\dot{E}10^{36}`$ erg s-1 and d$``$1.5 kpc, this ratio is over 400, similar to that of the gamma-ray emitting PSR B1706$``$44, strongly suggesting that the unidentified EGRET emission originates from the neutron star or the PWN. The difficulty with this interpretation is that our proposed neutron star and PWN lie several arcminutes outside the formal 95% error radius of the EGRET source. However, since the position of EGRET sources in the Galactic plane are subject to systematic error (e.g. Hunter et al. (1997)), it is entirely possible that this object may be the source of this gamma-ray emission.

## 4 Conclusion

In summary, the hard X-ray feature in IC 443 is best interpreted as a wind nebula powered by a young neutron star. The PWN interpretation is suggested first on morphological grounds by the Chandra and VLA images, which show a comet-shaped nebula with a soft X-ray point source at its apex. The measurement of significant polarization and a flat radio spectrum further strengthens the argument that this is synchrotron emission from a PWN. The X-ray spectrum of the point source is appropriately, but not uniquely fit by a black body model, suggesting that the emission is thermal in origin. The inferred physical properties of the nebula and point source (i.e. $`\tau `$, VNS, Ė, P, B) support our hypothesis that IC 443 and the pulsar were produced approximately 30,000 years ago in a core collapse supernova event.

Despite the evidence in favor of a real physical association between IC 443 and the compact source within the cometary nebula, there are a few puzzles that remain to be explained. As mentioned previously, one might expect the “tail” of the cometary nebula to point back to the origin of the SNR and the pulsar. However, the geometric center of IC 443 lies almost due north of the compact object, and not to the northeast as suggested by the direction of the cometary tail. Given the complex kinematics and distribution of the molecular and neutral gas in the vicinity of IC 443 (Giovanelli & Haynes (1979); van Dishoeck et al. 1993) it may be overly naïve to ascribe the blast center to the geometric center. A final problem concerns the apparent absence of X-ray and radio pulsations from the IC 443 point source. K97 carried out a pulse search using ROSAT and ASCA data, while Kaspi et al. (1996) have searched for radio pulsations. In our view, these null results present no immediate problem for our young neutron star hypothesis since pulsed X-rays from Vela-like pulsars have proven difficult to detect due to their soft spectra and contamination from the surrounding nebular emission (e.g. Becker & Truemper (1997)). Furthermore, it is entirely possible that the radio beam may not intersect our line of sight (Lorimer, Lyne & Camilo (1998)).

We thank D. J. Thompson, B. M. Gaensler, R. Petre, P. Slane, R. A. Chevalier, R.A. Fesen, and J. Kolena for useful discussions, and referee Randall Smith for his thoughtful comments. We also thank G. E. Allen, K. A. Arnaud, J-H. Rho, N. Evans, and M. Johnson for help with data analysis. This research was funded by NASA grant NRA 97-0SS-14.


format: text

text:
# Local Constraints on the Oscillating G Model

## I Introduction

The success of the standard cosmological model in describing the evolution of our Universe, beginning with the era of nucleosynthesis until the present state, has been confronted with serious difficulties resulting from the analysis of cosmological data. At large cosmological scales we find two main problems which are not dealt within the framework of the standard (old) cosmological model. The first one concerns the cosmological dark matter problem according to which the luminous matter (baryonic matter and radiation) content of the Universe represents only a small fraction of the total matter content. In fact, the inflationary models predicted that the total energy density $`\mathrm{\Omega }=1`$, with $`\mathrm{\Omega }`$ given in terms of the critical energy density . This prediction has recently been given further support by observational data resulting from the recent cosmic microwave background (CMB) experiments like Boomerang and Maxima and the high red-shift supernovae (SNIa) measurements , leading to the conclusion that the average energy-density of the Universe is indeed near the critical value. Obviously, these observations have increased the importance of the dark matter problem for the understanding of our Universe.

The second problem is related to the observations that indicate a periodicity of 128$`h^1`$ Mpc (where $`h`$ is the Hubble parameter in units of 100 km s-1 Mpc-1) in the galaxy number distribution, observed in deep pencil beams  in the north and south poles of our galaxy. This shocking discovery would, in its simplest interpretation, indicate that galaxies in the Universe are situated on the surface of concentric spheres with the center situated in our own galaxy. This is in complete contradiction with the basis of modern cosmology: the cosmological principle of homogeneity and isotropy of the Universe. It has been argued , and it seems to be the pervading view among researchers in the field, that this periodicity could be the result of the appearance of an intrinsic length scale in the distribution of matter. However, we have shown  that this explanation is not really satisfactory as there are scenarios of this type that result in a negligible probability for such observation to be obtained in a particular direction.

In a series of works  we have investigated an alternative model based on a massive scalar field which is non-minimally coupled to gravity. The oscillation of the scalar field in cosmic time results in a time-dependent effective gravitational constant. We have shown that this model leads to predictions which are in good agreement with most of the observational data. In fact, although this model was originally proposed  to explain the observed periodicity in the galaxy number distribution, we have shown that it was possible to adjudicate most of the energy density of the Universe to the oscillating massive scalar field which, therefore, could be regarded as candidate for the non-baryonic nature of the cosmological dark energy. That is, this model is able to explain simultaneously both, the problem of the cosmological dark energy and the problem of the periodicity in the galaxy number distribution. We have checked that the model satisfies some of the cosmological constraints. More precisely, we have seen that the model reproduces correctly the primordial nucleosynthesis of 4He, and is consistent with the present value of the energy density of baryonic matter and the age of the Universe. In this work, we will analyze the additional constraints following from local observations, namely, the Viking experiments , which impose bounds on the rate of change in time of the effective gravitational constant and on the effective Brans-Dicke parameter.

In a previous work  we have shown that all but one of the free parameters entering the model are fixed by the cosmological analysis and that with these values it was not possible to satisfy the Brans-Dicke bound. In this work, we analyze the possibility of overcoming this problem by relaxing the single condition freely imposed in our previous cosmological studies. We will show that even with this relaxation it is not possible to satisfy the local constraints and the periodicity observations simultaneously. This result indicates that either the behavior of the scalar field in the presence of local inhomogeneities is different from its behavior at large scales  or that a modified model would be necessary if we want to explain in a unified way the apparent galactic periodicity and the cosmological dark energy.

## II Constraints on the oscillating $`G`$ model

The dynamics of oscillating $`G`$ model is described by the Lagrangian:

$$=\left(\frac{1}{16\pi G_0}+\xi \varphi ^2\right)\sqrt{g}R\sqrt{g}\left[\frac{1}{2}(\varphi )^2+V(\varphi )\right],$$
(1)
where $`G_0`$ is Newton’s gravitational constant, $`\xi `$ stands for the non-minimally coupling constant, $`R`$ is the scalar curvature, $`\varphi `$ is the scalar field and $`V(\varphi )`$ is a scalar potential which in its simplest form is taken as the harmonic potential $`V=m^2\varphi ^2`$, with $`m`$ the mass of the scalar field. If we consider a time-dependent scalar field, the non-minimal coupling results in a time-dependent effective gravitational constant $`G_{\mathrm{eff}}=G_0(1+16\pi G_0\xi \varphi ^2)^1`$. The central feature of the oscillating $`G`$ model is that oscillations in the expectation value of $`\varphi `$ induce oscillations in $`G_{\mathrm{eff}}`$ and this leads to oscillations in the Hubble parameter $`H`$ which manifest themselves in the redshift measurements of distant points of the Universe. In turn, the redshift oscillations give rise to an apparent variation in the density of galaxies. Consequently, a temporal oscillation of the redshift can be mistakenly interpreted as a real spatial periodicity in the galaxy number distribution. This was used in previous works  to explain the observed periodicity of 128$`h^1`$Mpc in the distribution of galaxies in our Universe. To this end, we analyzed the Friedman-Robertson-Walker cosmology with a combination of two non-interacting perfect fluids (radiation and baryonic matter). From the field equations we obtain the following expression for the total effective energy density of the the system (see ):

$$\mathrm{\Omega }_{\mathrm{tot}}=\frac{1}{1+16\pi \xi \varphi _0^2}\left[\mathrm{\Omega }_{\mathrm{mat}}+\frac{4\pi }{3}\dot{\stackrel{~}{\varphi }}_0^2+\frac{4\pi }{3}\dot{\stackrel{~}{\omega }}^2\varphi _0^232\pi \xi \varphi _0\dot{\stackrel{~}{\varphi }}_0\right],$$
(2)
where

$$\dot{\stackrel{~}{\varphi }}_0=\frac{d\varphi }{d\stackrel{~}{t}}|_{\mathrm{today}},\stackrel{~}{t}=tH_0,\stackrel{~}{\omega }=\frac{\omega }{H_0},\stackrel{~}{m}^2=\frac{4\pi }{3}\stackrel{~}{\omega }^2.$$
(3)
In the above equation, a subscript “0” stands for the value of the corresponding quantity at present time $`t=t_0`$. The frequency of oscillation $`\omega =m\sqrt{3/4\pi }`$ is determined by the period of 128$`h^1`$Mpc observed in the pencil beam surveys and turns out to be $`\omega 147H_0`$. Here $`\mathrm{\Omega }_{\mathrm{matt}}=\mathrm{\Omega }_{\mathrm{bar}}+\mathrm{\Omega }_{\mathrm{rad}}`$. Notice that in Eq.(2) and for the present analysis we can neglect the contribution of the photon energy density $`\mathrm{\Omega }_{\mathrm{rad}}`$ because the observations of the cosmic microwave background radiation of 2.725 K implies that $`\mathrm{\Omega }_{\mathrm{rad}}10^3\mathrm{\Omega }_{\mathrm{bar}}`$. Furthermore, the value of $`\mathrm{\Omega }_{\mathrm{bar}}`$ must lie within the range $`[0.01,0.02]h^2`$ determined by the abundance of the light elements other than 4He . Finally, for the total energy density we take the value $`\mathrm{\Omega }_{\mathrm{tot}}=1`$ in accordance with the standard inflationary model and with the recent CMB and SNIa observations . Consequently, Eq.(2) can be interpreted as a constraint relating the initial cosmological values of the scalar field, $`\varphi _0`$ and $`\dot{\stackrel{~}{\varphi }}_0`$, and the coupling parameter $`\xi `$. Another, in some sense, more realistic approach would be to identify $`\mathrm{\Omega }_{\mathrm{matt}}`$ with the total amount of clumped matter in our Universe which would include besides the baryonic component also the so called Cold Dark Matter, leading us to take $`\mathrm{\Omega }_{\mathrm{matt}}0.3`$. However, we will see that even this drastic change of view does not alter our conclusions in a significant way.

A further constraint is imposed by the observed redshift-galaxy-count amplitude $`𝒜_0𝒪(0.5)`$, which for the oscillating $`G`$ model can be approximated by the expression

$$𝒜_0=\frac{16\pi \xi }{\stackrel{~}{\omega }}\left(\stackrel{~}{\omega }^2\varphi _0^2+\dot{\stackrel{~}{\varphi }}_0^2\right).$$
(4)
Here we are considering the additional term $`\dot{\stackrel{~}{\varphi }}_0^2`$ which was set to zero in previous analysis because we want to remove all the arbitrarily imposed conditions on the model in order to examine whether all the constraints can be solved simultaneously. Since the values of $`𝒜_0`$ and $`\stackrel{~}{\omega }`$ are fixed by the pencil beam observations, Eq.(4) represents a constraint between the values of $`\varphi _0`$, $`\dot{\stackrel{~}{\varphi }}_0`$ and $`\xi `$.

We call Eqs.(2) and (4) the global constraints of the oscillating $`G`$ model because the values of $`𝒜_0`$ and $`\stackrel{~}{\omega }`$ are fixed by the large scale observations of the galactic periodicity, and the values of $`\mathrm{\Omega }_{\mathrm{tot}}`$ and $`\mathrm{\Omega }_{\mathrm{mat}}`$ are the result of global cosmological observations.

On the other hand, the Solar System local observations impose an upper bound on the variation of the gravitational constant $`|\dot{G}/(GH)|0.3h^1`$ . For the oscillating $`G`$ model this yields

$$\beta =\frac{\dot{G}_{\mathrm{eff}}}{G_{\mathrm{eff}}H}=\frac{32\pi \xi \varphi _0\dot{\stackrel{~}{\varphi }}_0}{1+16\pi \xi \varphi _0^2},\mathrm{with}|\beta |0.3$$
(5)
It is well known that scalar-tensor models of the kind defined by the Lagrangian (1) can be transformed by means of a conformal transformation into an effective Brans-Dicke theory. Then, such models can be characterized by an effective Brans-Dicke parameter $`\omega _{\mathrm{BD}}^{\mathrm{eff}}`$ which must satisfy the lower bound imposed by the Viking experiments , $`\omega _{\mathrm{BD}}^{\mathrm{eff}}>3000`$. In the case of the oscillating $`G`$ model we obtain

$$\omega _{\mathrm{BD}}^{\mathrm{eff}}=\frac{1+16\pi \xi \varphi _0^2}{128\pi \xi ^2\varphi _0^2},$$
(6)
a constraint that relates $`\varphi _0`$ with $`\xi `$.

Now we proceed to the analysis of the global constraints (2) and (4) and the local constraints (5) and (6). In our previous cosmological studies, we were able to satisfy simultaneously the total energy constraint (2) as well as the nucleosynthesis and age constraints, together with the constraints for the amplitude (4) and the variation of the effective gravitational constant (5) by setting $`\dot{\stackrel{~}{\varphi }}_0=0`$. In fact, in this case the constraint (5) is automatically satisfied ($`\beta =0`$), whereas the constraints (2) and (4), together with the “plateau hypothesis”  that ensures a successful nucleosynthesis, fix the values of the remaining parameters $`\varphi _0`$ ($`10^3`$) and $`\xi `$ ($`6)`$. The evolution of the model with these conditions result in a value for the age of the Universe compatible with the standard bounds .

However, as we have shown in , with these values the oscillating $`G`$ model is unable to satisfy the Brans-Dicke limit (6) with $`\omega _{\mathrm{BD}}^{\mathrm{eff}}>3000`$ (or even the less severe bound $`\omega _{\mathrm{BD}}^{\mathrm{eff}}>500`$). The simplest possibility to overcome this problem is to relax the condition $`\dot{\stackrel{~}{\varphi }}_0=0`$ within the range allowed by the constraints (4) and (5). To this end, we replace the values of $`\varphi _0`$, $`\dot{\stackrel{~}{\varphi }}_0`$ and $`\xi `$ following from the constraints (4), (5) and (6) into the total energy constraint (2). Then we obtain

$$f(\omega _{\mathrm{BD}}^{\mathrm{eff}},\beta ,𝒜_0)=1\mathrm{\Omega }_{\mathrm{mat}}\frac{2b\stackrel{~}{\omega }\omega _{\mathrm{BD}}^{\mathrm{eff}}𝒜_0}{3(a+b)}+\frac{b+2\sqrt{a\stackrel{~}{\omega }𝒜_0b\stackrel{~}{\omega }^2}}{a}=0,$$
(7)
with

$$a=\beta ^2+4\stackrel{~}{\omega }^2,b=\beta ^2+2\stackrel{~}{\omega }𝒜_0+2\sqrt{\stackrel{~}{\omega }^2(𝒜_0^2\beta ^2)\beta ^2\stackrel{~}{\omega }𝒜_0},$$
(8)
a constraint that, for a specific value of $`\mathrm{\Omega }_{\mathrm{mat}}`$, determines the amplitude in terms of the effective Brans-Dicke parameter and the parameter $`\beta `$ (recall that the frequency $`\stackrel{~}{\omega }`$ has been fixed by the period of oscillation). Notice that the initial values $`\varphi _0`$ and $`\dot{\stackrel{~}{\varphi }}_0`$ do not appear at all in Eq.(7). In order to investigate the constraint (7) in a systematic way we have to solve the algebraic equation (7) as $`𝒜_0=𝒜_0(\beta )`$. Actually this is equivalent to solve the differential equation $`df(\beta ,𝒜_0(\beta ))/d\beta =0`$ (i.e., the resulting differential equation $`d𝒜_0/d\beta =F(𝒜_0,\beta )`$) subject to the boundary values $`(\beta ^i,𝒜_0^i)`$ such that $`f(\beta ^i,𝒜_0^i)=0`$, for a fixed $`\omega _{\mathrm{BD}}^{\mathrm{eff}}`$ and $`\mathrm{\Omega }_{\mathrm{mat}}`$. For instance, for $`\dot{\stackrel{~}{\varphi }}_0=0`$, and $`\omega _{\mathrm{BD}}^{\mathrm{eff}}=3000`$, $`\mathrm{\Omega }_{\mathrm{mat}}=0.0236`$ the pair $`(\beta ^i=0,𝒜_0^i0.022)`$ satisfies the constraint (7) as well as the remaining conditions (except of course the order of magnitude in the bound on $`𝒜_0`$). The result of this calculation is plotted in Figure 1 for two different values of $`\mathrm{\Omega }_{\mathrm{mat}}`$ within the range allowed by observations. We see that the range of values $`(\beta ,𝒜_0)`$ that satisfy $`f(\beta ,𝒜_0)=0`$ is extremely narrow and that all of the values for the amplitude within this range are situated well below the lower bound $`𝒜_0𝒪(0.5)`$ imposed by the redshift-galaxy-count observations. The conclusion is that the Brans-Dicke local constraint is not compatible with the observed value for the oscillation amplitude. Further numerical analysis of the constraint (7) show that an increase of the matter density $`\mathrm{\Omega }_{\mathrm{mat}}`$ or of the effective Brans-Dicke parameter leads to even lower values for the amplitude.

We conclude that the relaxation of the original condition $`\dot{\stackrel{~}{\varphi }}_0=0`$ does not allow the oscillating $`G`$ model to satisfy simultaneously both global and local constraints, and that we have to look for further generalizations of this model if we want to consider it as a candidate to explain the apparent galactic periodicity simultaneously with the nature of the non-baryonic dark matter content in the Universe. Needless is to say that had the model succeeded in these tests, then it would be necessary to confront the oscillating $`G`$ model to further tests in light of the recent CMB and SNIa observations.

Finally, it is worthwhile to emphasize that the generalized view on the problem of the galactic periodicity is that perhaps there is no problem at all and that such a “periodicity” is only the result of an excess of power at some characteristic length scales. While this could be the case, the simplest analysis on this matter shows that the existence of a characteristic distance in the large scale distribution is not enough to explain such observations , and therefore serious doubts arise in taking such a comfortable position. Clearly, the observation of galactic periodicity or lack thereoff in directions other than those corresponding to the north and south galactic poles will put an end to the controversy. On the other hand, if the existence of such periodicity in a large number of directions were to be confirmed we would be in the uncomfortable situation of having no model to account for it, and we would need to resort to variations on the oscillating $`G`$ model presented here as the only type of scenario capable of explaining such observations within the context of the cosmological principle.

Acknowledgments

We acknowledge partial support from DGAPA-UNAM Project No. IN121298 and from CONACyT Projects 32551-E and 32272-E.


text:
# 1 Introduction

## 1 Introduction

$`D`$-branes are topological solitons, which charges are described by Grothendieck K-groups . We begin our consideration with instantons  as they ideologically remind $`D`$-branes.

The self-dual equation for the $`SU\left(2\right)`$ instanton has the form $`F_{\mu \nu }=F_{\mu \nu }^{},\left(1.1\right)`$
where Yang-Mills strengths and potentials are defined as follows

$$F_{\mu \nu }=_\mu A_\nu _\nu A_\mu +[A_\mu ,A_\nu ],F_{\mu \nu }^{}=\frac{1}{2}\epsilon _{\mu \nu \lambda \kappa }F_{\lambda \kappa },\epsilon _{0123}=1,$$

$$F_{\mu \nu }=\frac{1}{2i}\sigma ^aF_{\mu \nu }^a\left(x\right),A_\mu =\frac{1}{2i}\sigma ^aA_\mu ^a\left(x\right),$$

$`\sigma ^a/2`$ are generators of the $`SU\left(2\right)`$ group, $`x=(x_0,x_1,x_2,x_3)`$ are coordinates in Euclidean space.

’t Hooft  has found the $`5k`$ parameter solution of the equation (1.1)

$$A_\mu =\frac{1}{2i}\eta _{\mu \nu }_\nu ln\varphi ,$$

where $`\eta _{\mu \nu }=\left(\epsilon _{oa\mu \nu }+\delta _{a\mu }\delta _{o\nu }\delta _{a\nu }\delta _{o\mu }\right)\sigma ^a,`$

$$\varphi =1+\underset{i=1}{\overset{k}{}}\frac{\lambda _i^2}{\left(x_\nu y_{\nu i}\right)^2},k=\frac{1}{16\pi ^2}d^4xTr\left(F_{\mu \nu }F_{\mu \nu }^{}\right),$$

$`\lambda _i`$ and $`y_{\nu i}`$ are parameters, $`k`$ is the topological charge.

It was shown in the works of Atiyah, Hitchin, Drinfeld and Manin  that dynamics of the equation (1.1) is encoded in the principal bundle

| $`SU\left(2\right)`$ | $``$ | $`P`$ |  |
| --- | --- | --- | --- |
|  |  | $``$ | (1.2) |
|  |  | $`S^4`$ |  |

and in the associated vector bundle

| $`\text{}^2`$ | $``$ | $`P`$ |  |
| --- | --- | --- | --- |
|  |  | $``$ | (1.3) |
|  |  | $`S^4`$ |  |

The potential $`A_\mu `$ is a connection of the vector bundle (1.3), therefore it is possible to calculate it by methods of algebraic geometry. For this purpose following the techniques of Horrocks and Barth , it is necessary to calculate the orthonormalized basis of sections $`Y^p`$ of the vector bundle (1.3). Then the potential $`A_\mu `$ will be equal to

$$A_\mu ^{pq}=(Y^p,_\mu Y^q).$$

The topological charge $`k`$ is an element of the homotopy group

$$\pi _3\left(SU\left(2\right)\right)=\text{},$$

characterizing the principal bundle (1.2).

In the next sections we shall apply the algebro-geometric approach to the description of $`D`$-branes.

## 2 Twisted K-theory

Let us consider the following principal bundles describing $`D`$-branes

| $`U\left(n\right)`$ | $``$ | $`P_H`$ |  |
| --- | --- | --- | --- |
|  |  | $``$ | (2.1) |
|  |  | $`X`$ |  |

| $`SU\left(n\right)/\text{}_n`$ | $``$ | $`P_H`$ |  |
| --- | --- | --- | --- |
|  |  | $``$ | (2.2) |
|  |  | $`X`$ |  |

| $`\underset{n\mathrm{}}{lim}SU\left(n\right)/\text{}_n`$ | $``$ | $`P_H`$ |  |
| --- | --- | --- | --- |
|  |  | $``$ | (2.3) |
|  |  | $`X`$ |  |

where $`X`$ is a compact manifold.
This bundles are characterized by the Dixmier-Douady invariant

$$\left[H\right]H^3(X,\text{}),$$

which represents the strength of a Neveu-Schwarz $`B`$-field  interacting with $`D`$-branes

$$H_{\mu \nu \rho }=_\mu B_{\nu \rho }+_\nu B_{\rho \mu }+_\rho B_{\mu \nu }.$$

For the bundle (2.1)
$`\left[H\right]=0,`$ id est $`H_{\mu \nu \lambda }=0,B_{\mu \nu }=0;`$
for the bundle (2.2)
$`n\left[H\right]=0,`$ id est $`H_{\mu \nu \lambda }=0,B_{\mu \nu }0;`$
for the bundle (2.3)
$`\left[H\right]0,`$ id est $`H_{\mu \nu \lambda }0,B_{\mu \nu }0.`$

The vector bundles associated with the principal bundles (2.1), (2.2), (2.3) are given by
$`E_H=P_H\times \text{}^n,`$ where $`Aut\left(\text{}^n\right)=U\left(n\right);`$ (2.4)
$`E_H=P_H\times M_n\left(\text{}\right),`$ where $`Aut\left(M_n\left(\text{}\right)\right)=SU\left(n\right)/\text{}_n;`$ (2.5)
$`E_H=P_H\times 𝒦,`$ where $`Aut\left(𝒦\right)=\underset{n\mathrm{}}{lim}SU\left(n\right)/\text{}_n;`$ (2.6)
$`M_n\left(\text{}\right)`$ is the algebra of $`n\times n`$ matrices, $`𝒦`$ is the algebra of compact operators.
The space of sections of the vector bundle (2.4) is the algebra of continuous functions $`C_0\left(X\right)`$;
the space of sections of the vector bundle (2.5) is the $`C^{}`$-algebra $`C_0(X,E_H)`$ , which is Morita-equivalent to the Azumaya algebra ;
the space of sections of the vector bundle (2.6) is the Rosenberg $`C^{}`$-algebra $`C_0(X,E_H)`$ .

Let us consider the short exact sequences

$$0𝒦C_0(X,E_H)0\left(2.7\right)$$

$$0𝒦^{^{}}C_0(X,E_H)0\left(2.8\right)$$

determining extensions of the algebra $`C_0(X,E_H)`$ by the algebra $`𝒦`$.
If in the commutative diagrams

$$0𝒦C_0(X,E_H)0$$

$$\tau _1$$

$$0𝒦B\left(\right)\stackrel{\pi }{}𝒬\left(\right)0$$

$$0𝒦^{^{}}C_0(X,E_H)0$$

$$\tau _2$$

$$0𝒦B\left(\right)\stackrel{\pi }{}𝒬\left(\right)0$$

(where $`B\left(\right)`$ is the algebra of bounded operators on a Hilbert space $``$, $`𝒬\left(\right)`$ denotes the Calkin algebra) morphisms are related by

$$\tau _2\left(s\right)=\pi \left(U\right)\tau _1\left(s\right)\pi \left(U\right)^{}$$

(s $``$ $`C_0(X,E_H)`$, $`U`$ is a unitary operator on $``$), then two extensions (2.7) and (2.8) are unitarily equivalent.

Let $`Ext(X,\left[H\right])`$ denote the set of unitary equivalence classes of extensions of $`C_0(X,E_H)`$ by $`𝒦`$ modulo splitting extensions. Then twisted K-groups are defined as follows

$$K_j(X,\left[H\right])=Ext(X\times \text{}^j,p_1^{}\left[H\right]),j\text{},$$

where $`p_1:X\times \text{}^jX`$ denotes projection onto the first factor.
The Bott periodicity theorem asserts:

$$K_j(X,\left[H\right])=\{\begin{array}{cc}\hfill K_0(X,\left[H\right])& \text{for even }j,\hfill \\ \hfill K_1(X,\left[H\right])& \text{for odd }j.\hfill \end{array}$$

## 3 Computation of twisted $`K`$-groups

A short exact sequence of $`C^{}`$-algebras

$$0𝒜C_0(X,E_H)0$$

induces a hexagon exact sequence of twisted K-groups

|  |  | $`K_0(X,\left[H\right])`$ | $``$ | $`K_0\left(\right)`$ | $``$ | $`K_0\left(𝒜\right)`$ |
| --- | --- | --- | --- | --- | --- | --- |
| $``$ |  |  |  | $``$ |  | (3.1) |
| $`K_1\left(𝒜\right)`$ | $``$ | $`K_1\left(\right)`$ | $``$ | $`K_1(X,\left[H\right])`$ |  |  |

The hexagon sequence (3.1) enables us to compute twisted K-groups.

Let $`X=S^3`$. We will compute twisted groups $`K_0(S^3,n\left[H\right])`$ and $`K_1(S^3,n\left[H\right])`$ for the case $`H_{\mu \nu \lambda }=0,B_{\mu \nu }0`$ . There is the short exact sequence

$$0C_0\left(S^2\right)𝒦C_0\left(U_1\right)𝒦C_0\left(U_2\right)𝒦$$

$$C_0(S^3,E_{n\left[H\right]})0$$

(where $`\{U_1,U_2\}`$ is the open cover of $`S^3`$ by hemispheres) ,
which induces the hexagon exact sequence of K-groups

| $`K_0(S^3,n\left[H\right])`$ | $``$ | $`K_0\left(U_1\right)K_0\left(U_2\right)`$ | $``$ | $`K_0\left(S^2\right)`$ |  |  |
| --- | --- | --- | --- | --- | --- | --- |
| $``$ |  |  |  | $``$ |  | (3.2) |
| $`K_1\left(S^2\right)`$ | $``$ | $`K_1\left(U_1\right)K_1\left(U_2\right)`$ | $``$ | $`K_1(S^3,n\left[H\right])`$ |  |  |

Since

$$K_0\left(S^2\right)=\text{},K_0\left(U_1\right)K_0\left(U_2\right)=\text{},$$

$$K_1\left(S^2\right)=0,K_1\left(U_1\right)K_1\left(U_2\right)=0,$$

the sequence (3.2) reduces to the exact sequence

$$0K_0(S^3,n\left[H\right])\text{}\text{}K_1(S^3,n\left[H\right])0\left(3.3\right)$$

Reduction of the sequence (3.3) to the exact sequence

$$0\text{}\text{}\text{}_n0$$

yields the following result:

$$K_0(S^3,n\left[H\right])=0,K_1(S^3,n\left[H\right])=\text{}_n.$$

The group $`K_1(S^3,n\left[H\right])=\text{}_n`$ determines topological charges of the $`D6`$-brane for the case $`H_{\mu \nu \lambda }=0,B_{\mu \nu }0.`$
For comparision we present the $`Dp`$-brane spectrum  for the case $`H_{\mu \nu \lambda }=0,B_{\mu \nu }=0`$ :

| $`Dp`$ | $`D9`$ | $`D8`$ | $`D7`$ | $`D6`$ | $`D5`$ | $`D4`$ | $`D3`$ | $`D2`$ | $`D1`$ | $`D0`$ | $`D(1)`$ |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
| $`S^{9p}`$ | $`S^0`$ | $`S^1`$ | $`S^2`$ | $`S^3`$ | $`S^4`$ | $`S^5`$ | $`S^6`$ | $`S^7`$ | $`S^8`$ | $`S^9`$ | $`S^{10}`$ |
| $`K_0(S^{9p})`$ |  | 0 |  | 0 |  | 0 |  | 0 |  | 0 |  |
| $`K_1(S^{9p})`$ | 0 |  | 0 |  | 0 |  | 0 |  | 0 |  | 0 |

Acknowledgement
This material was presented in the lecture given by the author at the Institute of Mathematics (Kiev, Ukraine). The author thanks audience for questions and comments.


text:
# Bulk viscosity in superfluid neutron star cores. II. Modified Urca processes in 𝑛𝑝𝑒𝜇 matter

## 1 Introduction

The bulk viscosity of neutron star matter determines damping of neutron star pulsations which could be excited during neutron star formation or triggered by instabilities developed during neutron star evolution (e.g., Cutler et al. cls90 (1990)). The bulk viscosity can damp gravitational radiation driven instabilities in rotating neutron stars and therefore influence the maximum rotation frequency (e.g., Lindblom l95 (1995), Zdunik 1996, Andersson & Kokkotas Andersson00 (2000)).

In this paper, we continue to study the bulk viscosity in neutron star cores (the layers of density $`\rho 1.5\times 10^{14}`$ g cm-3) assuming that stellar matter consists of neutrons $`n`$, protons $`p`$, electrons $`e`$, and muons $`\mu `$ and taking into account superfluidity of nucleons. Protons and neutrons constitute strongly interacting Fermi liquids while electrons and muons constitute almost ideal, strongly degenerate gases. Electrons are ultrarelativistic; muons may be nonrelativistic or moderately relativistic, depending on density. The bulk viscosity in this matter is mainly produced by neutrino processes of Urca type. These processes are known to be divided into most powerful direct Urca processes and much weaker modified Urca processes.

A direct Urca process consists of two reactions, direct and inverse one,

$$np+l+\overline{\nu }_l,p+ln+\nu _l,$$
(1)
where $`l`$ is either electron or muon and $`\nu _l`$ is an associated neutrino. Direct Urca processes are allowed (Lattimer et al. lpph91 (1991)) for some realistic equations of state at densities higher than certain threshold densities (of several $`\rho _0`$, where $`\rho _0=2.8\times 10^{14}`$ g cm-3 is the standard nuclear matter density). Thus, they may be open in the inner cores of rather massive neutron stars producing large bulk viscosity. The threshold density for the muon process is always higher than for the electron one.

However, for many equations of state, direct Urca processes are forbidden by momentum conservation at any density in a neutron star core. Moreover, at $`\rho 3\rho _0`$ they are prohibited for the majority of equations of state. Thus, they do not operate in the low and medium-mass neutron stars and in the outer cores of almost all neutron star models. If so, the bulk viscosity is determined by the reactions of modified Urca processes

$$n+Np+N+l+\overline{\nu }_l,p+N+ln+N+\nu _l,$$
(2)
where $`N`$ is an additional nucleon required to ensure momentum conservation. For instance, in the $`npe`$ matter one has two modified Urca processes, with $`N=n`$ or $`p`$, which will be referred to as the neutron and proton branches of the modified Urca process (e.g., Friman & Maxwell fm79 (1979), Yakovlev & Levenfish yl95 (1995)). In the $`npe\mu `$ matter, we have four processes corresponding to $`N=n`$ or $`p`$, and $`l=e`$ or $`\mu `$; they will be labeled as $`Nl`$. The rates of the modified Urca processes are 3–5 orders of magnitude lower than those of the direct Urca processes. The modified Urca processes either have no density threshold (as the neutron branch in the $`npe`$ matter) or have much lower thresholds than the direct Urca.

Sawyer (s89 (1989)) and Haensel & Schaeffer (hs92 (1992)) analyzed the bulk viscosity of non-superfluid $`npe`$ matter produced by the neutron branch of the modified Urca process and by the nucleon direct Urca process, respectively.

Recently we (Haensel et al. hly00 (2000), hereafter Paper I) have considered the bulk viscosity due to the direct Urca processes in the $`npe\mu `$ matter with superfluid neutrons and protons. In the present paper we analyze the bulk viscosity of this matter produced by the modified Urca processes.

Let us remind that damping of stellar pulsations is produced also by shear viscosity (calculated for instance by Flowers & Itoh fi76 (1976)). Contrary to the bulk viscosity the shear viscosity grows with decreasing temperature $`T`$ and governs the damping of pulsations of cold neutron stars. Estimates show (e.g., Paper I) that the bulk viscosity determines the damping of pulsations in nonsuperfluid stars at $`T10^8`$ K if direct Urca processes are open in the stellar cores or at $`T10^9`$ K if direct Urca processes are forbidden.

## 2 Bulk viscosity in non-superfluid matter

General expression for the bulk viscosity $`\zeta `$ of nonsuperfluid $`npe\mu `$ matter produced by nonequilibrium (direct or modified) Urca processes was obtained in Paper I. We restrict ourselves to neutrino transparent neutron star cores ($`T10^{10}`$K if direct Urca processes are forbidden). We will focus on the viscosity which describes damping of stellar pulsations with frequency $`\omega 10^3`$$`10^4`$ s-1. The pulsation frequency is typically much higher than the rates of Urca processes. In this high-frequency limit, $`\zeta `$ can be written as a sum of the partial bulk viscosities associated with all Urca processes. For the modified Urca processes ($`Nl`$) in $`npe\mu `$ matter,

$$\zeta =\zeta _{ne}+\zeta _{pe}+\zeta _{n\mu }+\zeta _{p\mu },\zeta _{Nl}=\frac{|\lambda _{Nl}|}{\omega ^2}\left|\frac{P}{X_l}\right|\frac{\eta _l}{n_b},$$
(3)
where $`\zeta _{Nl}`$ is a partial bulk viscosity, $`P`$ is the pressure, $`X_l=n_l/n_b`$ is an electron or muon fraction ($`n_l`$ and $`n_b=n_n+n_p`$ being the lepton and baryon number densities, respectively), $`\eta _l=\mu _n\mu _p\mu _l`$, $`\mu _j`$ is the chemical potential of particle species $`j`$, and $`\lambda _{Nl}`$ determines the difference of the rates of the direct and inverse reactions of an Urca process, Eq. (2): $`\mathrm{\Gamma }_{Nl}\overline{\mathrm{\Gamma }}_{Nl}=\lambda _{Nl}\eta _l`$. As in Paper I, we assume that in the absence of pulsations the matter is in the state of chemical equilibrium \[with respect to the beta and muon reactions (2)\], in which $`\eta _l=0`$. The pulsations are supposed to violate this state only slightly, $`|\eta |T`$. The derivatives in Eq. (3) are taken at chemical equilibrium with $`n_b`$, $`X_e`$ and $`X_\mu `$ as independent variables. In our notations, $`P/X_l`$ and $`\lambda _{Nl}`$ are negative.

As shown in Paper I, a partial bulk viscosity can be rewritten as

$$\zeta _{Nl}=\frac{|\lambda _{Nl}|}{\omega ^2}C_l^2,C_l=n_b\frac{\eta _l}{n_b}.$$
(4)
In Paper I the factors $`C_e`$ and $`C_\mu `$ were calculated and fitted by simple formulae for two model equations of state of the $`npe\mu `$ matter proposed by Prakash et al. (pal88 (1988)). These equations of state correspond to two different forms of the symmetry energy of nucleon matter, $`S(n_b)=\{13[u^{2/3}F(u)]+30F(u)\}`$ MeV, where $`u=n_b/n_0`$ and $`n_0=0.16`$ fm-3 is the saturation nucleon number density. One has $`F(u)=u`$ for model I; and $`F(u)=2u^2/(1+u)`$ for model II.

Thus the problem is reduced to evaluation of the reaction rates $`\mathrm{\Gamma }_{Nl}`$ and the factors $`\lambda _{Nl}`$. For a modified Urca process $`Nl`$ ($`\mathrm{}=c=k_\mathrm{B}=1`$)

$`\mathrm{\Gamma }_{Nl}`$ $`=`$ $`{\displaystyle [\underset{i=1}{\overset{4}{}}\frac{\mathrm{d}p_i}{(2\pi )^3}]\frac{\mathrm{d}p_l}{2\epsilon _l(2\pi )^3}\frac{\mathrm{d}p_\nu }{2\epsilon _\nu (2\pi )^3}}`$ (5)
$`\times (1f_N^{})(1f_p)(1f_l)f_nf_N(2\pi )^4`$
$`\times \delta (E_fE_i)\delta (P_fP_i){\displaystyle \frac{1}{s}}{\displaystyle \underset{\mathrm{spins}}{}}|M|^2,`$
where $`p_i`$ is the nucleon momentum ($`i=n,p,N,N^{}`$), $`p_l`$ and $`\epsilon _l`$ are the electron or muon momentum and energy, while $`p_\nu `$ and $`\epsilon _\nu `$ are the neutrino momentum and energy. The delta functions $`\delta (E_fE_i)`$ and $`\delta (P_fP_i)`$ describe conservation of the energy $`E`$ and the momentum $`P`$ of the particles in initial and final states, $`i`$ and $`f`$; $`|M|^2`$ is the squared matrix element of the reaction, $`s=2`$ is the symmetry factor which excludes double counting of the same collisions of identical particles, and $`f_i=\{1+\mathrm{exp}[(\epsilon _i\mu _i)/T]\}^1`$ is the Fermi-Dirac function.

To proceed further we need the matrix element $`M`$ which depends on the model of nucleon-nucleon interaction. We employ $`M`$ calculated by Friman & Maxwell (fm79 (1979)) using the one-pion-exchange model to describe the long-range part of the $`NN`$ interaction and the Landau Fermi-liquid theory to describe the short-range part. Their result \[their Eq. (39)\] is that the squared matrix element summed over spin states of all particles and averaged over orientations of the emitted neutrinos can be written in the form:

$$\underset{\mathrm{spins}}{}|M|^2=256G_\mathrm{F}^2g_A^2\left(\frac{f^\pi }{m_\pi }\right)^4\frac{\epsilon _\nu }{\epsilon _l}\alpha _N\beta _N.$$
(6)
Here, $`G_\mathrm{F}=1.436\times 10^{49}`$ erg cm3 is the Fermi weak coupling constant, $`g_A=1.26`$ is the axial-vector constant of weak hadron current, $`f^\pi 1`$ is the $`\pi N`$-interaction constant in the $`p`$-state in the one-pion-exchange model potential, and $`m_\pi `$ is the pion mass ($`\pi ^\pm `$). Furthermore, $`\alpha _N`$ describes momentum dependence of the squared matrix element in the Born approximation and $`\beta _N`$ contains various corrections. In principle, these factors are different for the neutron and proton branches of the modified Urca process (since $`n_nn_p`$). Nevertheless, bearing in mind the uncertainties of the employed model, we set $`\alpha _n=\alpha _p`$ and $`\beta _n=\beta _p`$. Friman & Maxwell (fm79 (1979)) adopted the values $`\alpha _n=1.13`$ (calculated at $`\rho =\rho _0`$) and $`\beta _n=0.68`$ (to account for correlation effects). In our approximation the squared matrix element in Eq. (5) can be regarded as constant and can be taken out of the integral.

Further evaluation of $`\mathrm{\Gamma }_{Nl}`$ is standard (e.g., Shapiro & Teukolsky st83 (1983)). It is based on phase-space decomposition, neglecting neutrino momentum in the momentum-conserving delta function and replacing the momenta of other particles $`p_i`$ by their Fermi momenta $`p_{\mathrm{F}i}`$ ($`i=n,p,l`$) whenever possible. Introducing the dimensionless quantities

$$x_i=\frac{\epsilon _i\mu _i}{T},x_\nu =\frac{\epsilon _\nu }{T},\xi =\frac{\eta _l}{T},$$
(7)
we can present the reaction rate in the form

$`\mathrm{\Gamma }_{Nl}`$ $`=`$ $`\mathrm{\Gamma }_{Nl}^{(0)}I;`$ (8)
$`I`$ $`=`$ $`{\displaystyle dx_\nu x_\nu ^2dx_ndx_pdx_Ndx_N^{}dx_l}`$ (9)
$`\times f(x_n)f(x_p)f(x_N)f(x_N^{})f(x_l)`$
$`\times \delta (x_n+x_p+x_N+x_N^{}+x_lx_\nu +\xi ),`$
where $`\mathrm{\Gamma }_{Nl}^{(0)}`$ is some typical reaction rate. We have transformed all the blocking factors $`(1f(x))`$ into the Fermi-Dirac functions $`f(x)`$ by replacing $`xx`$. The dimensionless energy integral $`I`$ for the $`n`$ and $`p`$ branches of the process is the same. For the reactions involving electrons, we have

$`\mathrm{\Gamma }_{ne}^{(0)}`$ $`=`$ $`{\displaystyle \frac{G_\mathrm{F}^2g_A^2m_n^3m_p^{}p_{\mathrm{F}p}}{\pi ^9\mathrm{}^{10}c^8}}\left({\displaystyle \frac{f^\pi }{m_\pi }}\right)^4(k_\mathrm{B}T)^7\alpha _n\beta _n`$ (10)
$``$ $`3.44\times 10^{25}\left({\displaystyle \frac{m_n^{}}{m_n}}\right)^3{\displaystyle \frac{m_p^{}}{m_p}}\left({\displaystyle \frac{n_p}{n_0}}\right)^{1/3}T_9^7{\displaystyle \frac{1}{\mathrm{cm}^3\mathrm{s}}},`$
$`\mathrm{\Gamma }_{pe}^{(0)}`$ $`=`$ $`\mathrm{\Gamma }_{ne}^{(0)}\left({\displaystyle \frac{m_p^{}}{m_n^{}}}\right)^2{\displaystyle \frac{\left(3p_{\mathrm{F}p}+p_{\mathrm{F}e}p_{\mathrm{F}n}\right)^2}{8p_{\mathrm{F}p}p_{\mathrm{F}e}}}\mathrm{\Theta }_{pe}.`$ (11)
Here, $`m_n^{}`$ and $`m_p^{}`$ are the neutron and proton effective masses, respectively, and the step function $`\mathrm{\Theta }_{pe}`$ takes into account the threshold character of the proton reaction branch: it is allowed ($`\mathrm{\Theta }_{pe}=1`$) for $`p_{\mathrm{F}n}<(3p_{\mathrm{F}p}+p_{\mathrm{F}e})`$ and forbidden ($`\mathrm{\Theta }_{pe}=0`$) otherwise (Yakovlev & Levenfish yl95 (1995)). The latter condition holds almost everywhere in the core; it may break at $`\rho \rho _0`$ only for the equations of state with small symmetry energy. Notice that Eqs. (10) and (11) have to be slightly modified if the direct Urca process is allowed (due to modification of the angular integrals; e.g., Yakovlev & Levenfish yl95 (1995)). However these modifications are of no practical interest: if allowed, the direct Urca process fully dominates over the modified Urca processes.

The rates of the neutron and proton branches of the muon modified Urca process are similar. They are obtained from Eqs. (8) – (11) in two steps. First, by introducing the extra factor $`v_{\mathrm{F}\mu }/c=(n_\mu /n_e)^{1/3}`$, where $`v_{\mathrm{F}\mu }`$ is the Fermi velocity of muons. Second, by replacing $`\mathrm{\Theta }_{pe}\mathrm{\Theta }_{p\mu }`$ that corresponds to the replacement $`p_{\mathrm{F}e}p_{\mathrm{F}\mu }`$. The muon step function $`\mathrm{\Theta }_{p\mu }`$ allows the proton branch of the process at $`p_{\mathrm{F}n}<(3p_{\mathrm{F}p}+p_{\mathrm{F}\mu })`$.

The rate $`\overline{\mathrm{\Gamma }}_{Nl}`$ of an inverse modified Urca reaction (lepton capture) differs from the rate of the direct reaction $`\mathrm{\Gamma }_{Nl}`$ only by the argument of delta function in the expression for $`I`$; one should replace $`\xi \xi `$ there. The difference of the direct and inverse reaction rates can be written as

$`\mathrm{\Gamma }_{Nl}\overline{\mathrm{\Gamma }}_{Nl}=\lambda _{Nl}\eta _l=\mathrm{\Gamma }_{Nl}^{(0)}\mathrm{\Delta }I,`$ (12)
$`\mathrm{\Delta }I={\displaystyle _0^{\mathrm{}}}dx_\nu x_\nu ^2\left[J(x_\nu \xi )J(x_\nu +\xi )\right],`$ (13)
$`J(x)={\displaystyle dx_ndx_pdx_Ndx_N^{}dx_lf(x_n)f(x_p)f(x_N)}`$
$`\times f(x_N^{})f(x_l)\delta (x_n+x_p+x_N+x_N^{}+x_lx).`$ (14)
The function $`J(x)`$ is evaluated analytically (e.g., Shapiro & Teukolsky 1983):

$$J(x)=\frac{x^4+10\pi ^2x^2+9\pi ^4}{24(\mathrm{e}^x+1)}.$$
(15)
Thus, the difference $`\mathrm{\Gamma }_{Nl}\overline{\mathrm{\Gamma }}_{Nl}`$ in nonsuperfluid matter is determined by the only parameter $`\xi =\eta _l/T`$. Furthermore, the integral $`\mathrm{\Delta }I`$ is taken analytically:

$`\mathrm{\Delta }I`$ $`=`$ $`{\displaystyle \frac{367\pi ^6}{1512}}\xi (\xi ),`$
$`(\xi )`$ $`=`$ $`1+{\displaystyle \frac{189\xi ^2}{367\pi ^2}}+{\displaystyle \frac{21\xi ^4}{367\pi ^4}}+{\displaystyle \frac{3\xi ^6}{1835\pi ^6}}.`$ (16)
Combined with Eq. (12), this relation gives:

$$|\lambda _{Nl}|=\frac{\mathrm{\Gamma }_{Nl}^{(0)}}{T}\frac{\mathrm{\Delta }I}{\xi }.$$
(17)
As mentioned above, we do not consider large deviations from the chemical equilibrium and restrict ourselves to the case $`|\eta _l|T`$ in which $`1`$. It is instructive to estimate characteristic times of Urca processes. Taking into account that $`\dot{n}_e=\overline{\mathrm{\Gamma }}_{Ne}\mathrm{\Gamma }_{Ne}=\lambda \eta `$, we can introduce, for instance, the time $`\tau _w\delta n_e/(\lambda _eT)`$ required for an electron Urca process to shift the chemical potential difference by $`\eta T`$, within the thermal width of the Fermi levels. In this case the electron number density will be changed by $`\delta n_en_e(T/\mu _e)`$. Adopting typical parameters of neutron star matter we have $`\tau _w20T_9^4`$ s for direct Urca process (Paper I) and $`\tau _wT_9^6`$ months for modified Urca processes. Therefore, in any case $`\tau _w`$ is much longer than typical periods of neutron star pulsations ($`10^3`$ s) which justifies validity of the high-frequency approximation used in our analysis (see above).

Finally, combining Eqs. (4) and (17) we obtain the partial bulk viscosity of the $`npe\mu `$ matter, $`\zeta _{Nl}=\zeta _{Nl0}`$ (the subscript ‘0’ refers to nonsuperfluid case), produced by the modified Urca processes for $`|\eta _l|/T1`$. In standard physical units

$`\zeta _{ne0}`$ $`=`$ $`{\displaystyle \frac{367G_\mathrm{F}^2g_A^2m_n^3m_p^{}p_{\mathrm{F}p}C_e^2}{1512\pi ^3\mathrm{}^{10}c^8\omega ^2}}\left({\displaystyle \frac{f^\pi }{m_\pi }}\right)^4(k_\mathrm{B}T)^6\alpha _n\beta _n`$ (18)
$``$ $`1.49\times 10^{19}\left({\displaystyle \frac{m_n^{}}{m_n}}\right)^3{\displaystyle \frac{m_p^{}}{m_p}}\left({\displaystyle \frac{n_p}{n_0}}\right)^{1/3}\left({\displaystyle \frac{C_e}{100\mathrm{MeV}}}\right)^2`$
$`\times T_9^6\omega _4^2\alpha _n\beta _n\mathrm{g}\mathrm{cm}^1\mathrm{s}^1,`$
$`\zeta _{pe0}`$ $`=`$ $`\zeta _{ne0}\left({\displaystyle \frac{m_p^{}}{m_n^{}}}\right)^2{\displaystyle \frac{\left(3p_{\mathrm{F}p}+p_{\mathrm{F}e}p_{\mathrm{F}n}\right)^2}{8p_{\mathrm{F}p}p_{\mathrm{F}e}}}\mathrm{\Theta }_{pe},`$ (19)
$`\zeta _{n\mu 0}`$ $`=`$ $`\zeta _{ne0}\left({\displaystyle \frac{p_{\mathrm{F}\mu }}{p_{\mathrm{F}e}}}\right)\left({\displaystyle \frac{C_\mu }{C_e}}\right)^2,`$ (20)
$`\zeta _{p\mu 0}`$ $`=`$ $`\zeta _{ne0}({\displaystyle \frac{C_\mu m_p^{}}{C_em_n^{}}})^2{\displaystyle \frac{\left(3p_{\mathrm{F}p}+p_{\mathrm{F}\mu }p_{\mathrm{F}n}\right)^2}{8p_{\mathrm{F}p}p_{\mathrm{F}\mu }}}({\displaystyle \frac{p_{\mathrm{F}\mu }}{p_{\mathrm{F}e}}})\mathrm{\Theta }_{p\mu },`$ (21)
where $`\omega _4=\omega /(10^4\mathrm{s}^1)`$ and $`T_9=T/(10^9\mathrm{K})`$.

For the two equations of state I and II (see above), the density dependence of the partial bulk viscosities is fitted by simple expressions presented in the Appendix. Figure 1 shows the total bulk viscosity of nonsuperfluid matter versus the baryon number density $`n_b`$ for both models. The contributions of the modified and direct Urca processes are included; the latter contribution is calculated using the results of Paper I. The dotted curves show the bulk viscosity for the simplified equations of state based on the same nuclear energy but with the appearance of muons artificially forbidden.

The viscosities for models I and II are similar. The muon thresholds are $`n_b=0.150`$ fm-3 and $`0.152`$ fm-3 for models I and II. At lower densities the muons are absent and the bulk viscosity is produced by the electron modified Urca processes (and the dotted curves are identical to the solid ones). At densities above the thresholds the muon modified Urca processes are switched on and amplify the bulk viscosity. Since $`\zeta _{N\mu }v_{\mathrm{F}\mu }`$, the muon bulk viscosity is switched on without any jump. On the contrary, the bulk viscosity due to the direct Urca processes is switched on in a jump-like manner at their threshold densities. The thresholds for the electron direct Urca process are $`n_b=0.414`$ fm-3 and $`0.302`$ fm-3, for models I and II. As explained in Paper I, the presence of muons lowers the threshold density, mainly due to increasing the number density and Fermi momenta of protons. On the other hand, the muons lower the bulk viscosity produced by the electron direct Urca processes (by decreasing $`n_e`$). At larger $`n_b`$ the total bulk viscosity jumps again, at $`n_b=0.503`$ fm-3 and $`0.358`$ fm-3 for models I and II. This time the jump is associated with switching on the muon direct Urca process. The contribution of the muon direct Urca into the bulk viscosity is even larger than the contribution of the electron direct Urca (Paper I). With this contribution the bulk viscosity becomes larger than in the $`npe`$ matter.

## 3 Bulk viscosity of superfluid matter

### 3.1 Superfluid reduction factors

Now we turn to superfluid reduction of the bulk viscosity produced by the modified Urca processes. As in Paper I, we restrict ourselves to the motions of the $`npe\mu `$ fluid in which all components move with same macroscopic velocity. This assumption reduces equations of stellar pulsations to those of one-fluid hydrodynamics, with a single coefficient of the bulk viscosity; our bulk viscosity $`\zeta `$ is identical to $`\zeta _2`$ of Landau and Lifshitz (1987). The superfluid reduction of the bulk viscosity is similar to the superfluid reduction of the neutrino emissivity. Since the latter problem has been the subject of extensive studies (e.g., Yakovlev et al. yls99 (1999)) we omit technical details.

We adopt the traditional assumption (e.g., Yakovlev et al. yls99 (1999)) that the protons form Cooper pairs due to singlet-state (1S0) pairing (superfluidity of type A, in notations of Paper I). As for the neutrons, they would be thought to undergo either the singlet-state pairing or the triplet-state (3P2) pairing with zero projection of the total angular momentum of the Cooper pair on the quantization axis (superfluidity of type B). Numerous microscopic calculations (see Yakovlev et al. yls99 (1999) for references) show that the singlet-state pairing of neutrons takes place at $`\rho \rho _0`$, while their triplet-state pairing is realized at higher densities. The critical temperatures $`T_{cn}`$ and $`T_{cp}`$ of the neutron and proton superfluids are very model dependent, and we treat them as free parameters.

Microscopically, superfluidity introduces a gap $`\delta `$ into momentum dependence of the nucleon energy, $`\epsilon (p)`$. Near the Fermi surface ($`|pp_\text{F}|p_\text{F}`$), this dependence can be written as

$$\begin{array}{c}\epsilon =\mu \sqrt{\delta ^2+v_\mathrm{F}^2(pp_\mathrm{F})^2}\mathrm{at}p<p_\text{F},\hfill \\ \epsilon =\mu +\sqrt{\delta ^2+v_\mathrm{F}^2(pp_\mathrm{F})^2}\mathrm{at}pp_\text{F}.\hfill \end{array}$$
(22)
In case A the gap is isotropic (independent of orientation of $`p`$), $`\delta _\mathrm{A}=\mathrm{\Delta }_\mathrm{A}(T)`$, while in case B the gap depends on the angle $`\vartheta `$ between $`p`$ and the quantization axis, $`\delta _\mathrm{B}=\mathrm{\Delta }_\mathrm{B}(T)(1+3\mathrm{cos}^2\vartheta )^{1/2}`$. In both cases $`\mathrm{\Delta }(T)`$ is the amplitude which describes temperature dependence of the gap. The amplitude is derived from the standard equation of the Bardeen - Cooper - Schrieffer theory.

For further analysis we introduce the dimensionless quantities ($`k_\mathrm{B}=1`$ as before):

$$v=\frac{\mathrm{\Delta }(T)}{T},\tau =\frac{T}{T_c},y=\frac{\delta }{T}.$$
(23)
The dimensionless gap $`y`$ can be presented in the form: $`y_\mathrm{A}=v_\mathrm{A}`$, $`y_\mathrm{B}=v_\mathrm{B}\sqrt{1+3\mathrm{cos}^2\vartheta }`$. The dimensionless gap amplitude $`v`$ depends only on $`\tau =T/T_c`$. This dependence is accurately fitted as (Levenfish & Yakovlev ly94 (1994)):

$`v_\mathrm{A}`$ $`=`$ $`\sqrt{1\tau }\left(1.4560.157/\sqrt{\tau }+1.764/\tau \right),`$
$`v_\mathrm{B}`$ $`=`$ $`\sqrt{1\tau }\left(0.7893+1.188/\tau \right).`$ (24)
Quite generally, the high-frequency bulk viscosity can be presented in the form

$$\zeta =\underset{Nl}{}\zeta _{Nl0}R^{(N)},$$
(25)
where $`\zeta _{Nl0}`$ is a partial bulk viscosity of nonsuperfluid matter \[given by Eqs. (18)–(21)\] and $`R^{(N)}=R_{Nl}^{(N)}`$ is a factor which describes its superfluid reduction. The factors for the modified Urca processes $`Ne`$ and $`N\mu `$ are equal since the superfluidity affects only nucleons. In non-superfluid matter matter, we have $`R^{(N)}=1`$ and reproduce the results of Sect. 2.

Thus the problem consists in calculating the reduction factors $`R^{(N)}`$. Each factor depends on two parameters, $`v_n`$ and $`v_p`$, dimensionless gap amplitudes of the neutrons and protons. As in Paper I we may assume that the superfluidity affects noticeably only the factor $`\lambda _{Nl}`$ in the expression for the bulk viscosity, Eq. (4). By making the standard angular-energy decomposition in the expressions for the reaction rates $`\mathrm{\Gamma }_{nl}`$ and $`\overline{\mathrm{\Gamma }}_{nl}`$, one can easily show that $`\lambda _{Nl}`$ contains the factor

$$𝒥(v_p,v_n)=4\pi \left[\underset{i=1}{\overset{5}{}}\mathrm{d}\mathrm{\Omega }_i\right]\delta \left(\underset{i=1}{\overset{5}{}}p_i\right)\mathrm{\Delta }I,$$
(26)
where d$`\mathrm{\Omega }_i`$ is the solid angle element in the direction of $`p_i`$ and $`\mathrm{\Delta }I`$ is the energy integral given by Eq. (13). This integral is the main quantity affected by the superfluidity. At $`\xi 1`$ the integrand of Eq. (13) contains the function $`J(x_\nu \xi )J(x_\nu +\xi )2\xi J(x_\nu )/x_\nu `$, where $`J(x_\nu )`$ is given by Eq. (14). Thus, for small deviations from equilibrium one can transform Eq. (13) to:

$`\mathrm{\Delta }I`$ $`=`$ $`4\xi {\displaystyle _0^+\mathrm{}}dx_\nu x_\nu \left[{\displaystyle \underset{i=1}{\overset{5}{}}}{\displaystyle _{\mathrm{}}^+\mathrm{}}dx_if(x_i)\right]`$ (27)
$`\times \delta \left({\displaystyle \underset{i=1}{\overset{5}{}}}x_ix_\nu \right).`$
Here, $`i=n,p,N,N^{},l`$.

Let us now label the quantities in the nonsuperfluid case by the subscript ‘0’. The angular-energy decomposition yields

$$𝒥_0=A_{nl0}\mathrm{\Delta }I_0,\mathrm{\Delta }I_0=\frac{367\pi ^6\xi }{1512}.$$
(28)
In this case $`A_{nl0}`$ is the angular integral containing the momentum conserving delta function, in which the neutrino momentum is neglected and the momenta of all other particles are placed at their Fermi surfaces. For the processes of our study:

$`A_{nl0}`$ $`=`$ $`{\displaystyle \frac{(4\pi )^5}{2p_{\mathrm{F}n}^3}},`$
$`A_{pl0}`$ $`=`$ $`{\displaystyle \frac{(4\pi )^5}{2p_{\mathrm{F}n}^2p_{\mathrm{F}p}}}{\displaystyle \frac{\left(3p_{\mathrm{F}p}+p_{\mathrm{F}l}p_{\mathrm{F}n}\right)^2}{8p_{\mathrm{F}p}p_{\mathrm{F}l}}}\mathrm{\Theta }_{pl}.`$ (29)
Notice that Yakovlev & Levenfish (yl95 (1995)) inaccurately determined $`A_{pl0}`$. Here we present the corrected expression.

Generalization of $`\lambda _{Nl}`$ to the superfluid case can be achieved by introducing the neutron and proton energy gaps into Eq. (26). Then

$`\lambda _{Nl}=\lambda _{Nl0}R^{(N)},R^{(N)}(v_p,v_n)={\displaystyle \frac{𝒥}{A_{Nl0}\mathrm{\Delta }I_0}},`$ (30)
where $`𝒥`$ is given by Eq. (26) with

$`\mathrm{\Delta }I(y_1,y_2)`$ $`=`$ $`4\xi {\displaystyle _0^+\mathrm{}}dx_\nu x_\nu \left[{\displaystyle \underset{i=1}{\overset{5}{}}}{\displaystyle _{\mathrm{}}^+\mathrm{}}dx_if(z_i)\right]`$ (31)
$`\times \delta \left({\displaystyle \underset{i=1}{\overset{5}{}}}z_ix_\nu \right),`$
$`z_i=(\epsilon _i\mu )/T=\mathrm{sign}(x_i)\sqrt{x_i^2+y_i^2}`$ for nucleons, and $`z_l=x_l`$. Equation (30) is the general expression for calculating the reduction factor $`R^{(N)}`$.

### 3.2 Superfluidity of neutrons or protons

First, let the neutrons be normal while the protons undergo Cooper pairing of type A. Since the energy gap is isotropic, the angular integral in Eq. (26) is separated from the energy one, being the same as in the nonsuperfluid case, Eq. (29). Then in Eq. (26) we can put $`y_n=0`$ and $`y_p=v_p`$. Accordingly, the reduction factor depends on the only parameter $`v=v_p`$. For a strong proton superfluidity ($`\tau =T/T_{cp}1`$, $`v1`$), the asymptotes of the reduction factors for the neutron and proton branches of the modified Urca process are:

$`R_{p\mathrm{A}}^{(n)}`$ $`=`$ $`{\displaystyle \frac{126\sqrt{2}}{1835\pi ^{5.5}}}v^{5.5}\mathrm{e}^v={\displaystyle \frac{0.004061}{\tau ^{5.5}}}\mathrm{exp}\left({\displaystyle \frac{1.764}{\tau }}\right),`$ (32)
$`R_{p\mathrm{A}}^{(p)}`$ $`=`$ $`{\displaystyle \frac{6048}{367\pi ^6}}\gamma v^5\mathrm{e}^{2v}={\displaystyle \frac{0.1604}{\tau ^5}}\mathrm{exp}\left({\displaystyle \frac{3.528}{\tau }}\right),`$
$`\gamma `$ $`=`$ $`{\displaystyle \frac{\pi }{16}}\left(21\sqrt{3}{\displaystyle \frac{51}{2}}\mathrm{ln}(\sqrt{3}+2)\right)0.548.`$ (33)
The proton branch is affected by the proton superfluidity stronger than the neutron branch. This is natural (e.g., Yakovlev et al. yls99 (1999)) because in the proton branch we have three protons which belong to the superfluid component of matter, while in the neutron branch there is only one such particle. Under the number of reacting nucleons we mean their total number (in the initial plus final states) in the direct or inverse reaction of an Urca process.

We have calculated $`R_{p\mathrm{A}}^{(n)}`$ and $`R_{p\mathrm{A}}^{(p)}`$ numerically in a wide range of $`v=v_p`$ and proposed the fits which reproduce the numerical results (with the mean errors of $`1`$%) and the asymptotes (32) and (33):

$`R_{p\mathrm{A}}^{(n)}`$ $`=`$ $`{\displaystyle \frac{a^{5.5}+b^{3.5}}{2}}\mathrm{exp}\left(3.245\sqrt{(3.245)^2+v^2}\right),`$
$`a=0.1863+\sqrt{}(0.8137)^2+(0.1310v)^2,`$
$`b=0.1863+\sqrt{}(0.8137)^2+(0.1437v)^2;`$
$`R_{p\mathrm{A}}^{(p)}`$ $`=`$ $`c^5\mathrm{exp}\left(5.033\sqrt{(5.033)^2+(2v)^2}\right),`$
$`c=0.3034+\sqrt{}(0.6966)^2+(0.1437v)^2.`$
Now consider normal protons and superfluid neutrons (superfluidity of type B). Then in Eq. (26) we may put $`z_p=x_p`$ and the reduction factor depends only on $`v=v_n`$. For a strong neutron superfluidity ($`\tau =T/T_{cn}1`$, $`v1`$) the asymptotes of the reduction factors are

$`R_{n\mathrm{B}}^{(n)}`$ $`=`$ $`{\displaystyle \frac{6048}{367\pi ^6}}{\displaystyle \frac{2\gamma }{3\sqrt{3}}}v^4\mathrm{e}^{2v}={\displaystyle \frac{0.00720}{\tau ^4}}\mathrm{exp}\left({\displaystyle \frac{2.376}{\tau }}\right),`$ (36)
$`R_{n\mathrm{B}}^{(p)}`$ $`=`$ $`{\displaystyle \frac{42\sqrt{3}}{1835\pi ^5}}v^5\mathrm{e}^v={\displaystyle \frac{3.066\times 10^4}{\tau ^5}}\mathrm{exp}\left({\displaystyle \frac{1.188}{\tau }}\right).`$ (37)
In this case the neutron branch is affected by the superfluidity stronger than the proton branch.

The calculation of $`R_{n\mathrm{B}}^{(n)}(v_n)`$ for intermediate values $`v_n`$ is difficult since there are three superfluid neutrons with anisotropic gaps. An approximate expression for $`R_{n\mathrm{B}}^{(n)}(v_n)`$ will be proposed in Sect. 3.3. As for the factor $`R_{n\mathrm{B}}^{(p)}`$, we have calculated it numerically for a wide range of $`v`$ and obtained the fit which reproduces the numerical results (with mean error $`1\%`$) and the asymptotes:

$`R_{n\mathrm{B}}^{(p)}`$ $`=`$ $`{\displaystyle \frac{a^5+b^3}{2}}\mathrm{exp}\left(2.110\sqrt{(2.110)^2+v^2}\right),`$ (38)
$`a=0.1973+\sqrt{(0.8027)^2+(0.1257v)^2},`$
$`b=0.1973+\sqrt{(0.8027)^2+(0.1428v)^2}.`$
The factors $`R^{(N)}`$ which describe reduction of the bulk viscosity of $`npe\mu `$ matter by superfluidity of neutrons or protons are plotted in Figs. 2 and 3. For comparison, we present also the reduction factors $`R^{(\mathrm{D})}`$ for the partial bulk viscosity produced by the direct Urca processes (Paper I). Figure 2 displays the reduction factors versus $`T/T_c`$, while Fig. 3 presents them versus $`v`$. As seen from these figures, the reduction factors may be separated into two groups. The first group contains the factors $`R_{p\mathrm{A}}^{(n)},R_{n\mathrm{B}}^{(p)},R_{p\mathrm{A}}^{(\mathrm{D})},R_{n\mathrm{B}}^{(\mathrm{D})}`$ which describe reduction of the reactions with one superfluid nucleon. The other group contains the factors $`R_{p\mathrm{A}}^{(p)}`$ and $`R_{n\mathrm{B}}^{(n)}`$ (asymptote) for the reactions with three superfluid nucleons at once. Notice that the curves inside each group in Fig. 3 are much closer to one another than in Fig. 2. The same situation takes place for the reduction factors of neutrino emissivities (e.g., Yakovlev et al. yls99 (1999)).

### 3.3 Superfluidity of neutrons and protons

If the neutrons and protons are superfluid at once, calculation of the reduction factors is complicated. If however we do not need very accurate results, we can use approximate similarity relations analogous to those suggested by Yakovlev & Levenfish (yl95 (1995)) and described also by Yakovlev et al. (yls99 (1999)) for the reduction of neutrino emissivities in different reactions:

$`R_{\mathrm{AB}}^{(p)}(v_p,v_n)`$ $``$ $`{\displaystyle \frac{R_{\mathrm{AB}}^{(\mathrm{D})}(2v_p,v_n)}{R_{n\mathrm{B}}^{(\mathrm{D})}(v_n)}}R_{n\mathrm{B}}^{(p)}(v_n),`$ (39)
$`R_{\mathrm{AB}}^{(n)}(v_p,v_n)`$ $``$ $`{\displaystyle \frac{R_{\mathrm{AB}}^{(\mathrm{D})}(v_p,2v_n)}{R_{p\mathrm{A}}^{(\mathrm{D})}(v_p)}}R_{p\mathrm{A}}^{(n)}(v_p).`$ (40)
Here, $`R_{p\mathrm{A}}^{(\mathrm{D})},R_{n\mathrm{B}}^{(\mathrm{D})}`$ and $`R_{\mathrm{AB}}^{(\mathrm{D})}`$ are the factors which define reduction of the direct-Urca bulk viscosity by superfluidity of protons, neutrons and protons+neutrons, respectively. These factors were obtained in Paper I. Notice that Eq. (39) becomes exact for normal protons ($`v_p=0`$), while Eq. (40) is exact for normal neutrons. We may assume also that the reduction of the $`N`$-branch of the modified-Urca viscosity by neutron superfluidity can be approximated (for $`v10`$) by the factor: $`R_{n\mathrm{B}}^{(n)}R_{p\mathrm{A}}^{(p)}(v_n)`$.

Finally, Fig. 4 illustrates reduction of the bulk viscosity of $`npe\mu `$ matter (due to the modified Urca processes) with decreasing temperature by superfluidity of neutrons of type B or protons of type A for $`n_b=2n_0`$ and $`\omega =10^4`$ s-1. Thick solid line shows the viscosity of non-superfluid matter (cf. with Fig. 1). Thin dashed lines present the bulk viscosity suppressed by the proton superfluidity at several selected critical temperatures $`T_{cp}`$ indicated near the curves. The dot-and-dashed line shows the effect of neutron superfluidity ($`T_{cn}=10^{10}`$ K) for normal protons. We see that superfluid reduction of the bulk viscosity depends on temperature, superfluidity type, and the critical temperatures $`T_{cn}`$ and $`T_{cp}`$. If $`T_{cn}`$ and $`T_{cp}`$ are not higher than $`10^{10}`$ K for $`n_b2n_0`$, then the superfluid reduction cannot be very large, say, for $`T10^9`$ K. It cannot reach more than two orders of magnitude in the case of superfluidity of one nucleon species or six orders of magnitude if $`n`$ and $`p`$ are superfluid at once for $`T=10^9`$ K at $`T_{cn}=T_{cp}=10^{10}`$ K (Fig. 5 of Paper I). The reduction grows exponentially with further decrease of $`T`$.

## 4 Conclusions

We have calculated the bulk viscosity due to non-equilibrium modified Urca processes in neutron star cores composed of neutrons, protons, electrons, and muons ($`npe\mu `$). We have considered non-superfluid matter (Sect. 2) and described also reduction of the bulk viscosity by superfluidity of neutrons and protons (Sect. 3). In combination with our previous work (Paper I) we have obtained realistic description of the bulk viscosity provided by the modified and direct Urca processes in the $`npe\mu `$ matter of the neutron star cores. The results can be used for studying damping of pulsations in neutron stars and gravitational radiation driven instabilities in rotating neutron stars. Strong superfluidity of neutrons and protons reduces the bulk viscosity and creates favorable conditions for the development of these instabilities.

###### Acknowledgements.

Two of the authors (KPL and DGY) acknowledge hospitality of N. Copernicus Astronomical Center in Warsaw. The authors are grateful to anonymous referee for useful comments. This work was supported in part by the RBRF (grant No. 99-02-18099), INTAS (grant No. 96-0542), and KBN (grant 2 P03D 014 13).

## Appendix

To calculate the bulk viscosity from the equations obtained above one needs to know the factors $`C_e`$ and $`C_\mu `$ and the number densities $`n_n`$, $`n_p`$, $`n_e`$, $`n_\mu `$, and $`n_b`$. Other functions of density such as $`\mathrm{\Gamma }_{Nl}`$ are expressed through the Fermi-momenta of various particles $`p_{\mathrm{F}i}=\mathrm{}(3\pi ^2n_i)^{1/3}`$ and hence through $`n_i`$. The practical expression for calculating $`C_e`$ and $`C_\mu `$ for any given equation of state was obtained in Paper I \[Eq. (19)\]. The number densities $`n_i`$ are usually available for a given equation of state. It is sufficient to specify $`n_p`$, $`n_\mu `$ and $`n_b`$, since $`n_n=n_bn_p`$ and $`n_e=n_pn_\mu `$.

For example, consider two model equations of state I and II used in Paper I and in the present paper. The fits of the quantities $`q_l=(n_l/n_0)^{1/3}(C_l/100\mathrm{MeV})^2`$, $`l=e`$ and $`\mu `$, as functions of the baryon number density $`n_b`$ are given by Eq. (36) of Paper I. They yield practical expressions for $`C_l(n_b)`$ provided the particle number densities are known. Thus, we must fit the number densities $`n_i=n_i(n_b)`$. We present the proton number density as $`n_p=n_{p0}+\mathrm{\Delta }n_p`$, where $`n_{p0}`$ is the proton number density for the case in which creation of muons is artificially forbidden. For the equations of state I and II the muons appear at $`n_b=n_{b\mu }=0.15000`$ and 0.5122831 fm-3, respectively. Therefore, we have $`n_p=n_{p0}=n_e`$ and $`n_\mu =0`$ for $`n_bn_{b\mu }`$.

We have fitted $`n_{p0}`$ in the density range from $`n_b=0.048`$ fm$`{}_{}{}^{3}=0.3n_0`$ to 1.92 fm$`{}_{}{}^{3}=12n_0`$ as

$$n_{p0}=\eta +\alpha n^\beta /(1+\gamma n^\delta ),$$
(A1)
and we have also fitted $`\mathrm{\Delta }n_p`$ and $`n_\mu `$ in the range from $`n_b=n_{b\mu }`$ to 1.92 fm-3 as

$`\mathrm{\Delta }n_p`$ $`=`$ $`a\delta n^{1.5}(1+d\delta n^e)/(1+b\delta n^c),`$
$`n_\mu `$ $`=`$ $`A\delta n^{1.5}(1+D\delta n^E)/(1+B\delta n^C).`$ (A2)
Here, $`nn_b`$, $`\delta n=n_bn_{b\mu }`$, and all number densities are expressed in fm-3. The fit parameters are $`\alpha =0.7295`$, $`\beta =2.503`$, $`\gamma =2.604`$, $`\delta =1.307`$, $`\eta =0`$, $`a=0.1306`$, $`b=1.887`$, $`c=1.509`$, $`A=0.1614`$, $`B=0.1785`$, $`C=1.235`$, $`d=D=0`$ for model I; $`\alpha =3.64`$, $`\beta =3.374`$, $`\gamma =12.07`$, $`\delta =2.096`$, $`\eta =2.673\times 10^5`$, $`a=0.1919`$, $`b=4.946`$, $`c=1.279`$, $`d=0.1558`$, $`e=0.9253`$, $`A=0.2376`$, $`B=3.294`$, $`C=1.468`$, $`D=2.194`$, $`E=1.201`$ for model II. The maximum errors of the fits of $`n_p`$ and $`n_\mu `$ do not exceed 1.3% in the indicated density ranges. The presented fits give all functions versus $`n_b`$. The fit expressions of $`n_b`$ via mass density $`\rho `$ for three versions of each equation of state I and II corresponding to three different values of the compression modulus of saturated nuclear matter, $`K=120`$, 180 and 240 MeV, are given by Eq. (37) of Paper I.