Effect of magnetized phonons on electrical and thermal conductivity of neutron star crust
We study electrical and thermal conductivities of degenerate electrons emitting and absorbing phonons in a strongly magnetized crystalline neutron star crust. We take into account modification of the phonon spectrum of a Coulomb solid of ions caused by a strong magnetic field. Boltzmann transport equation is solved using a generalized variational method. The ensuing three-dimensional integrals over the transferred momenta are evaluated by two different numerical techniques, the Monte-Carlo method and a regular integration over the first Brillouin zone. The results of the two numerical approaches are shown to be in a good agreement. An appreciable growth of electrical and thermal resistivities is reported at quantum and intermediate temperatures ( is the ion plasma temperature) in a wide range of chemical compositions and mass densities of matter even for moderately magnetized crystals ( and are the ion cyclotron and plasma frequencies). This effect is due to an appearance of a soft () phonon mode in the magnetized ion Coulomb crystal, which turns out to be easier to excite than acoustic phonons characteristic of the field-free case. These results are important for modelling magneto-thermal evolution of neutron stars.
keywords:dense matter – stars: neutron.
Magnificent seven is a group of isolated neutron stars emitting quasi-thermal X-ray and optical radiation and situated at distances below 500 pc. They have rotation periods in the 3–11 s range and magnetic fields and characteristic ages, estimated from the magneto-dipole braking formula, of (1–3) G and (1–4) yrs, respectively. So, these are fairly old nearby neutron stars, which, presumably, are seen solely due to their residual heat, but which all have magnetic fields higher than those of typical rotation-powered pulsars. Quite naturally, these stars were viewed as prime examples of objects, which are heated additionally by decay of their magnetic fields (e.g., Kaspi, 2010). Further evidence in favor of this neutron star heating mechanism was presented by Pons et al. (2007), who have shown that there was a correlation between surface temperature and dipole magnetic field for such diverse objects as soft gamma repeaters, anomalous X-ray pulsars, isolated neutron stars, and rotation-powered pulsars.
Vigano et al. (2013) have performed detailed modelling of magneto-thermal evolution of neutron stars. In these simulations, magnetic field decay was able to provide enough heat to explain even the hottest sources such as magnetars. Clearly, in order to obtain effective decay of the magnetic field the electrical conductivity of matter must be relatively low. To achieve that Vigano et al. (2013) have assumed that the nuclear pasta layer at the boundary of the inner crust with the core possessed a high level of disorder. Accordingly, they used impurity scattering with a very high impurity parameter ( in the rest of the crust) to model this situation.
Later on such a high value of the impurity parameter was justified in a molecular dynamics simulation of the nuclear pasta layer by Horowitz et al. (2015). For instance, in the lasagna phase, these authors have discovered a very peculiar defect, which looks like a spiral staircase leading from one lasagna “floor” to the next and which may be responsible for .
While this conjecture is certainly exciting, it brings about a number of questions. Firstly, one may wonder, how come there is so little disorder at lower and higher densities and so much disorder in the narrow density range from to g cm. Secondly, the high disorder layer then must be present also in weakly magnetized neutron stars, where it is not required by observations. Thirdly, as shown by Kaminker et al. (2007), it is more economical to place a heat source of magnetars in more shallow layers in order to avoid excessive energy loss via neutrinos, though these considerations did not take into account proportional worsening of the electron thermal conductivity in the high disorder layer. Ultimately, the question becomes whether it is possible to organize a fast field decay by more conventional means?
In a crystalline neutron star crust charge and heat are transported by degenerate electrons. The main mechanism impeding their transport is emission and absorption of lattice vibrations (phonons) by electrons. Phonons in a magnetized Coulomb crystal of ions (or “magnetized phonons”) are fundamentally different from the field-free case (Usov, Grebenschikov & Ulinich, 1980; Nagai & Fukuyama, 1982, 1983; Baiko, 2009). In the absence of the field, the phonons are acoustic, i.e. their frequency depends linearly on the wave vector. In the presence of the field, a soft mode with a quadratic dependence of the frequency on the wave vector appears.
The appearance of the soft mode can be understood by referring to a simple problem of a charged oscillator with frequency (mass and charge ) brought into a uniform magnetic field (e.g., Landau & Lifshitz, 1994). Its frequency becomes
where . Assuming that , i.e. it represents an ion frequency in an acoustic phonon with sound speed , the square root can be expanded at small wave vectors , which results in the quadratic dependence . Clearly, the mode is softer for stronger magnetic fields. This is precisely what one obtains if a lattice dynamics problem is solved for the Coulomb crystal in the magnetic field.
The soft mode has a dramatic effect on the thermodynamic properties of the crystal at low temperatures. Due to its low frequency it is much easier to excite, which produces, for instance, an increase of the crystal heat capacity per particle by a factor of at temperature and phonon magnetization (cf. Fig. 2 of Baiko, 2009). In this case, is the ion plasma frequency, is the ion plasma temperature (in what follows ), is the ion cyclotron frequency, while , , and are ion charge, mass, and number density, respectively. It appears that the soft mode will affect the effective rate of phonon emission/absorption by electrons in a comparable way. In fact, this rate [cf. Eq. (17) of Baiko & Yakovlev (1995), also see Yakovlev & Urpin (1980); Raikh & Yakovlev (1982)] is expressed via functions and , which are given by similar averages over the phonon spectrum as the thermodynamic quantities and which are expected to grow strongly upon inclusion of the magnetic field. These simple considerations indicate that the problem of electron interaction with magnetized phonons in the Coulomb crystal deserves a serious study.
2 Kinetic equation
Let us write the standard stationary Boltzmann kinetic equation (see Sect. 5 for a discussion of applicability limits):
In this case, is the electron momentum distribution function, which depends also on position but is independent of the electron spin orientation, is the electron velocity for momentum , and is the Lorentz force:
, , and being electron charge, electric and magnetic fields, respectively. On the right hand side of Eq. (2), is the transition probability per unit time of electron with momentum to a state with momentum with an absorption of a phonon with momentum summed over primed and averaged over non-primed electron spin states. Other refer to probabilities of the other three possible processes of the same kind, while is the phonon momentum distribution function. The summation is over all , all and over all phonon modes at given . For brevity, the phonon mode index is suppressed here.
We linearize the Boltzmann equation assuming weak deviation of the electron distribution from the local equilibrium
and also assuming equilibrium phonon distribution
where and are local electron chemical potential and temperature, respectively, while and are electron energy and phonon frequency.
For a strongly degenerate system it is customary (e.g., Ziman, 1960) to assume that the electron distribution deviates noticeably from the local equilibrium one only near the Fermi surface:
where is a new unknown function.
On the right-hand side of Eq. (2), the transition probabilities are determined by relativistic electron scattering probability off a potential
The potential is a sum of screened Coulomb potentials of all ions (labeled by index ) minus the potential of the static lattice, in which all ions are fixed at their lattice nodes . Furthermore, is the ion charge number, is the inverse Thomas-Fermi screening length, and is the operator of ion displacement, which, upon quantization of the ion motion, is given by (Usov et al., 1980; Baiko, 2009):
In this case, is the total number of ions, and are phonon creation and annihilation operators, the sum is over all phonon modes, while vectors are analogous to phonon polarization vectors in a non-magnetized crystal, but have different orthogonality and normalization properties (see Baiko, 2009, for details).
We restrict ourselves to the one-phonon approximation, in which only the first order term in is kept in Eq. (7). Multi-phonon processes become important closer to the crystal melting temperature (Baiko et al., 1998), whereas we are mostly focused on lower temperatures. Then various on the right-hand side of Eq. (2) differ only by the energy-conserving delta-functions and the linearized kinetic equation reads:
In this case, and are electron Fermi velocity and momentum, , and , where is a reciprocal lattice vector chosen in such a way that is in the first Brillouin zone. Let us note, that the electron distribution responds to and to the electric field in exactly the same way, which allows us to set without any loss of generality.
3 Solution by variational principle
Typically, in a magnetic field, one adopts the relaxation time approximation for the collision integral on the right-hand side of Eq. (2) or (9), where the relaxation time is taken from the respective non-magnetic problem (e.g., Ziman, 1960; Askerov, 1970; Urpin & Yakovlev, 1980). This approach is well justified if we do not expect the collision probability to depend on and also the collisions are nearly elastic (i.e. the energy difference between the initial and the final electron states is much lower than ). In our case, both of these assumptions are invalid, since the scatterer (i.e. phonons’) properties depend on , and the most interesting effect is anticipated at low temperatures, where the details of the phonon spectrum are important and the transferred energy is of the order of .
In the field-free case the variational principle has been used successfully (e.g., Ziman, 1960; Flowers & Itoh, 1976; Raikh & Yakovlev, 1982) to study the transport properties at low temperatures. There also exists a generalization of the variational principle to the case of non-zero magnetic field (Ziman, 1960). The kinetic equation (9) can be written in a symbolic form as
where is the unknown function, denotes the terms on the left-hand side which drive the system out of the equilibrium (i.e. terms with an electric field or a temperature gradient), is the collision operator, and is the magnetic operator, which contains the momentum derivative of and which is moved from the left-hand side to the right.
The variational solution is looked for in the form of a linear combination of some basis functions :
where are unknown constants. These constants are solution of a system of linear equations
in which angle brackets denote a scalar product
for the charge transport problem and
for the heat transport problem. Then the scalar products involving and become trivial, for instance
By contrast, the scalar products involving the collision operator are not easy to evaluate. They contain 6D-integrals over and , of which the integrals over and can be taken by the standard methods. The remaining 4D-integrals over the solid angles and can be reduced to 3D-integrals over the transferred momentum over a ball of radius as follows
In the process, we have integrated over the azimuthal angle of vectors and with respect to . The presence of this integration allows one to replace without any loss of accuracy
for the electrical conductivity, and make similar replacements for the thermal conductivity ( is the electron Fermi energy).
Combining everything together we obtain for the electrical and thermal conductivities ( and ) the following systems of linear equations, respectively:
with and being the electron gyro-frequency
The quantities can be called collision frequency tensors. They read
and differs by the replacement
In principle, Eqs. (19) can be solved easily for arbitrary orientations of the magnetic field, electric field, and temperature gradient with respect to the crystal axes. However, we shall assume that the crystal always forms in such a way that the magnetic field is directed along the symmetry axes resulting in the minimum zero-point energy (for bcc crystal this would be a direction towards one of the nearest neighbors, Baiko, 2009). At the same time, we would like to average over the azimuthal angle of the electric field or temperature gradient with respect to the magnetic field. This can be done in various ways. For instance, one can average the effective collision frequency or the effective relaxation time, or the kinetic coefficients themselves. It is not immediately clear, which kind of averaging should be preferred. We note, that the averaging procedure may affect the final practical results in a non-trivial manner (cf. Kobyakov & Pethick, 2015).
In this paper we adopt the simplest approach and average the system of linear equations (19) itself. Then for the electrical conductivity we obtain (upper index is understood for all and )
Once the deviation of the distribution function from the local equilibrium is known
where and refer to components of the electric field parallel and perpendicular to the magnetic field, one can calculate the electric current
and determine the components of the electrical conductivity tensor
In this case it is assumed that the -axis is directed along . Thus the perpendicular time determines the conductivity across the magnetic field, the parallel time determines the conductivity along the magnetic field, and determines the Hall conductivity. In the relaxation time approximation one obtains exactly the same formulae with and replaced by the nonmagnetic relaxation time (Urpin & Yakovlev, 1980).
For the heat flux we get:
and the components of the thermal conductivity tensor read
where it is understood that must be used. To obtain from the following replacements have to be made in the expressions (24) for parallel and perpendicular frequencies
4 Numerical calculations
The most difficult part of the problem is to evaluate 3D-integrals Eq. (17) over the transferred momentum over the ball of radius . In Fig. 1 a sphere of radius (“double Fermi sphere”) is shown by the (red) circle, while rhombi represent the first Brillouin zone shifted by all possible reciprocal lattice vectors, which results in a complete filling of the momentum space. The actual (three-dimensional) first Brillouin zone of the bcc lattice is depicted in the inset (a rhombododecahedron). A reciprocal lattice vector connects the center of the sphere with the center of each rhombus. In order to find the phonon frequencies and the -coefficients for a given , we need to subtract the respective reciprocal lattice vector to obtain a phonon wave vector in the first Brillouin zone, and then solve a lattice dynamics problem at this .
Thus, essentially, we need to integrate over the first Brillouin zone. Methods of such integration in a Coulomb solid are well-developed (e.g., Holas method, Holas, 1977; Albers & Gubernatis, 1981; Baiko, Potekhin & Yakovlev, 2001), however, this is only true if we need to integrate over the whole Brillouin zone. Unfortunately, as is evident from Fig. 1, there are numerous incomplete pieces of the Brillouin zones inside the ball. They appear due to the intersections of the double Fermi sphere with the zones. It seems that the problem of enumerating all these pieces and integrating over them is too cumbersome to deal with. Thus, the first approach that comes to mind (see also Raikh & Yakovlev, 1982; Baiko & Yakovlev, 1995) is to use the Monte-Carlo (MC) integration method, i.e. to find the average of the integrand over the ball by randomly selecting points inside the ball.
At low temperatures, though, the integrand develops very sharp and narrow maxima in the vicinity of points . This means that one needs many more MC coin tosses to get an adequate representation of the integrand (in order to reduce temperature 10 times, times more MC steps is required). Another approach to this integration, which is suited better to lower temperatures, is thus desired and it is illustrated in Fig. 2. In this picture, the dash-dotted (red) line is the volume of the double Fermi ball measured in volumes of the Brillouin zone. The dashed (blue) line is the number of Brillouin zones, which lie wholly inside the ball, and the dotted (black) line is the number of Brillouin zones, which intersect the ball (i.e. it is the number of incomplete zones inside the ball). Finally, the solid (green) line is the number of Brillouin zones, which have their centers inside the ball. It is clear that the total volume of the latter zones is a reasonable approximation of the ball volume. We may thus try to integrate over these zones, taking them wholly, and expect this to be a very accurate result at low temperatures as all points will be accounted for. At higher temperatures, the MC method must be more precise because of the proper integration domain used in this case.
Let us compare integration results by the two methods. We define dimensionless functions (cf. Baiko & Yakovlev, 1995) and plot in Fig. 3 the ratio of calculated by the MC method to the same quantity calculated via the Brillouin zone integration (BZI). The ratios are plotted as functions of for several ion charge numbers (dotted), 6 (short-dashed), 10 (dash-dotted), 26 (long-dashed), 40 (solid) and phonon magnetizations (thin, blue), 1 (intermediate thickness, green), 10 (thick, red). In all cases electrons are assumed to be ultrarelativistic. For and 10 the data are based on MC steps, while for , 26, and 40 this number is doubled.
We note that the ratios look rather insensitive to the charge number with the exception of helium (dots). The lowest magnetization curves display a sharp relative decrease of MC integrals at . In this temperature range BZI results are robust. They are insensitive to the number of integration points in the Holas method, provided it is not too small. By contrast, MC results at these low temperatures are very sensitive to the number of MC steps due to the effect described above. At higher magnetizations though, the MC results converge much better, because wider regions around the Brillouin zone center continue to contribute. At , MC calculations, presumably, are more reliable. The deviation of the BZI values can be explained by different geometry of the overall integration domains (ball vs. a set of whole Brillouin zones). At some intermediate temperatures MC and BZI coincide and the ratio becomes 1. In principle, one can use BZI below this point and MC above it. However, the inaccuracy of the BZI results at high temperatures is so insignificant ( for all elements except helium) that we have decided to drop MC calculations altogether and use BZI in the entire temperature range. Only for low- elements such as helium the accuracy of BZI at becomes somewhat insufficient () and, in principle, the MC approach should be preferred.
The same graphs could be constructed for as well as for but they would reveal nothing principally new.
Our final results are presented in Fig. 4. These plots show quantities (directly proportional to the collision frequencies) as functions of for several representative of the outer neutron star crust chemical elements ranging from He to Zr and mass densities ranging from to g cm. Thin (blue) and thick (red) curves show and , respectively. Vertical dotted lines show the melting temperature if it belongs to the displayed range of temperatures (calculations are artificially extended beyond melting for illustrative purposes). The ion magnetization parameter is set to 0.01, 1, and 10 for dashed, dash-dotted, and solid lines, respectively, and is marked near the curves. We note that at given , , and composition, the mass density determines the electron degree of relativity and enters Eq. (21) only through the Thomas-Fermi screening length and the back-scattering suppression factor in big parentheses.
The lowest magnetization curves essentially reproduce earlier field-free calculations (except at extremely low temperatures , where even such a weak magnetic field starts making a difference). Accordingly, for parallel and perpendicular -functions merge. At high temperatures , all curves merge, which means that the scattering is quasi-elastic and the exact phonon spectrum is not important. At intermediate and low temperatures and one observes a divergence of parallel and perpendicular frequencies ( in all cases) as well as their significant growth as compared to the low magnetization case. This produces a proportional decrease of the electrical and thermal conductivities and represents the main result of our work.
The increase of the collision frequencies is not as strong as one would expect from the comparison with the specific heat in the Introduction. This can be explained by a suppression of the amplitude Eq. (8) of ion deviation from its lattice node in a strong magnetic field and a corresponding reduction of the phonon potential Eq. (7), whereas the argument in the Introduction took into account only statistical weight of phonons.
5 Limits of applicability
While solving the transport problem, we have introduced several important simplifications. First of all, we have defined momentum distribution function for electrons , Eq. (2), as opposed to considering distribution of electrons over quantum numbers appropriate in the magnetic field (i.e., longitudinal momentum, Landau level number etc., see, e.g., Yakovlev, 1984). Secondly, we have neglected a magnetic field effect on the electron screening of the ion potential, Eq. (7). Thirdly, we have neglected an effect (not studied to this day) of electron screening (with or without field) on the magnetized Coulomb crystal phonons.
It is well-known that actual kinetic coefficients in a magnetic field oscillate, e.g. as functions of density at given temperature (e.g., Yakovlev, 1984; Potekhin, 1996, 1999). The lower is the temperature, the stronger are the oscillations. However, if more than one Landau level is populated, one may expect that our consideration based on the momentum distribution yields an average over the oscillations value of the kinetic coefficients. The actual coefficients will trace this average more closely as temperature gets higher. The same argument allows us to neglect the effect of electron screening modification by the magnetic field for electrons populating more than one Landau level. The effect of electron screening on magnetized phonons has not been analysed yet. However, in the absence of the field the Coulomb crystal phonons are modified by electron screening in a well-studied way Pollock & Hansen (1973); Baiko (2002) and this makes a barely visible imprint on kinetics Baiko & Yakovlev (1995). While caution is required, we can expect the same conclusion to hold in the presence of the magnetic field as well.
If electrons populate only the ground Landau level, our approach becomes unreliable, and the problem of electron transport must be reconsidered with account of magnetized phonons, inelasticity of scattering, and electron screening. This problem is further complicated by the fact that the Fermi temperature of electrons populating only the ground Landau level drops very rapidly with decrease of density. Consequently, such electrons become non-degenerate at relevant temperatures and densities. Especially problematic also is a quick onset of the strong screening regime, which accompanies the drop of the Fermi temperature and electron kinetic energy.
In Fig. 5 we show the phonon magnetization parameter as a function of density for , , and G by solid (blue), dash-dotted (green), and dashed (red) curves, respectively. This parameter is independent of ion charge and mass numbers. Dotted vertical lines of the same colours show the densities at which the ground Landau level is completely filled for and for the same magnetic fields (these fractions along with respective magnetic fields are written near the lines). For realistic materials, and therefore, our consideration is valid quantitatively at densities higher than the left-most dotted line for each field. This corresponds to , i.e. moderate and weak magnetization of crystal phonons.
As already mentioned, the case of electrons on the ground Landau level, and thus , requires a separate consideration.
We have calculated electrical and thermal conductivities of degenerate electrons emitting and absorbing phonons in a strongly magnetized crystalline neutron star crust. The novel features of our study include (i) an account of modification of the phonon properties in the magnetic field and (ii) an application of a generalized variational method to solve the Boltzmann kinetic equation in the magnetic field with due consideration of the electron-phonon process inelasticity.
Our results apply at weak and intermediate phonon magnetization and indicate a significant growth of electrical and thermal resistivities at quantum and intermediate temperatures in a broad range of chemical compositions and mass densities of matter. This effect is due to an appearance of a soft phonon mode in the magnetized ion Coulomb crystal, which is easier to excite in an interaction with an electron than an ordinary acoustic phonon in the field-free case.
At our results hint at an even stronger increase of the electrical and thermal resistivities, however, a detailed analysis of the problem is needed under the assumption that electrons populate only the ground Landau level.
These results are important for quantitative modelling of cooling of neutron stars and evolution of their magnetic fields. In particular, they imply an accelerated decay of the magnetic field in the outer neutron star crust accompanied by an additional heat release.
The author is grateful to Prof. D.G. Yakovlev for reading the manuscript and making useful remarks. This work was supported by RSF, grant No. 14-12-00316.
- Albers & Gubernatis (1981) Albers R.C., Gubernatis R.C., 1981, LASL report No. LA-8674-MS
- Askerov (1970) Askerov B.M., 1970, Kinetic Effects in Semiconductors. Nauka, Moscow
- Baiko (2002) Baiko D.A., 2002, Phys. Rev. E, 66, 056405
- Baiko (2009) Baiko D.A., 2009, Phys. Rev. E, 80, 046405
- Baiko et al. (1998) Baiko D.A., Kaminker A.D., Potekhin A.Y., Yakovlev D.G., 1998, Phys. Rev. Lett., 81, 5556
- Baiko et al. (2001) Baiko D.A., Potekhin A.Y., Yakovlev D. G., 2001, Phys. Rev. E, 64, 057402
- Baiko & Yakovlev (1995) Baiko D.A., Yakovlev D.G., 1995, Astron. Lett., 21, 702
- Chugunov (2012) Chugunov A.I., 2012, Astron. Lett., 38, 25
- Flowers & Itoh (1976) Flowers E., Itoh N., 1976, Astrophys. J., 206, 218
- Holas (1977) Holas A., 1977, J. Comput. Phys., 23, 150
- Horowitz et al. (2015) Horowitz C.J., Berry D.K., Briggs C.M., Caplan M.E., Cumming A., Schneider A.S., 2015, Phys. Rev. Lett., 114, 031102
- Kaminker et al. (2007) Kaminker A.D., Yakovlev D.G., Potekhin A.Y., Shibazaki N., Shternin P.S., Gnedin O.Y., 2007, Astrophys. Space Sci., 308, 423
- Kaspi (2010) Kaspi V., 2010, PNAS, 107, 7147
- Kobyakov & Pethick (2015) Kobyakov D., Pethick C.J., MNRAS, 449, L110
- Landau & Lifshitz (1994) Landau L.D., Lifshitz E.M., 1994, The Classical Theory of Fields. Butterworth-Heinemann, Boston
- Nagai & Fukuyama ( 1982) Nagai T., Fukuyama H., 1982, J. Phys. Soc. Jap., 51, 3431
- Nagai & Fukuyama ( 1983) Nagai T., Fukuyama H., 1983, J. Phys. Soc. Jap., 52, 44
- Pollock & Hansen (1973) Pollock E.L., Hansen J.P., 1973, Phys. Rev. A, 8, 3110
- Pons et al. (2007) Pons J.A., Link B., Miralles J.A., Geppert U., 2007, Phys. Rev. Lett., 98, 071101
- Potekhin (1996) Potekhin A.Y., 1996, Astron. Astrophys., 306, 999
- Potekhin (1999) Potekhin A.Y., 1999, Astron. Astrophys., 351, 787
- Raikh & Yakovlev (1982) Raikh M.E., Yakovlev D.G., 1982, Astrophys. Space Sci., 87, 193
- Urpin & Yakovlev (1980) Urpin V.A., Yakovlev D.G., 1980, Sov. Astron., 24, 425
- Usov et al. (1980) Usov N.A., Grebenschikov Yu.B., Ulinich F.R., 1980, ZhETF, 78, 296
- Vigano et al. (2013) Vigano D., Rea N., Pons J.A., Perna R., Aguilera D.N., Miralles J.A., 2013, MNRAS, 434, 123
- Yakovlev (1984) Yakovlev D.G., 1984, Astrophys. Space Sci., 98, 37
- Yakovlev & Urpin (1980) Yakovlev D.G., Urpin V.A., 1980, Sov. Astron., 24, 303
- Ziman (1960) Ziman J.M., 1960, Electrons and Phonons. Clarendon, Oxford