A Simple Model of Radiation from a Magnetized Neutron Star: Accreted Matter and Polar Hotspots

A simple and well known model for thermal radiation spectra from a magnetized neutron star is further studied. The model assumes that the star is internally isothermal and possesses dipole magnetic field (B<= 1e14 G) in the outer heat-insulating layer. The heat transport through this layer makes the surface temperature distribution anisotropic; any local surface element is assumed to emit a blackbody (BB) radiation with a local effective temperature. It is shown that this thermal emission is nearly independent of the chemical composition of insulating envelope (at the same taken averaged effective surface temperature). Adding a slight extra heating of magnetic poles allows one to be qualitatively consistent with observations of some isolated neutron stars.


Yury N. Gnedin
This paper is dedicated to the memory of Yury N. Gnedin (1935Gnedin ( -2018. He was my senior colleague at Theoretical Astrophysics Department of Ioffe Insitute from 1971 to 1984, and we often met later, when he moved to Pulkovo Observatory. He was always interested in many scientific fields. For instance, his pet subjects were radiation transfer, neutron stars and magnetic fields (although he never lost scientific interest in a great amount of other things). When neutron stars were discovered in 1967, he was extremely enthusiastic about them and initiated ultrafast spreading of interest to these objects among colleagues. He was the first author of the first publication [1] on neutron stars at the Ioffe Institute after the discovery of these stars. He wrote -with different co-authors -plenty of (now) classical papers on neutron star physics. In particular, they include basic formulation of the radiation transfer problem in a strongly magnetized plasma [2], theoretical prediction of electron cyclotron lines in the spectra of radiation from magnetized neutron stars [3], the first realistic studies of ionization in the magnetized hydrogen atmospheres of neutron stars [4]. Also, he was greatly fond of history.
Although we did not work much together I, as many others, respected Yury Gnedin very much. When he could help, he really did, and he was a Real Gentleman -friendly and open to any one, from a prominent scientist to a hopeless student. I will always remember the wonderful atmosphere, created by Yury Gnedin, and warm hospitality of his family and home. It is amazing that both his sons have grown up into excellent scientists.
I present a small piece of neutron star physics on thermal radiation from magnetized neutron stars. Hopefully, he would like it.

Introduction
I will mostly consider thermal X-ray radiation from the surface of a middle-aged (t 10 6 yr) isolated magnetized neutron star. The radiation is supposed to emerge from warm neutron star interior and emitted from a thin atmosphere or condensed surface (e.g. [5]), where the radiation spectrum is formed. The interiors of the star are thought to be nearly isothermal because of high thermal conductivity of dense matter. Nevertheless, the interiors are thermally insulated from the surface by a thin heat blanketing envelope, where the thermal conduction is poorer and sensitive to the chemical composition and magnetic field. The field redirects the emerging heat flow and results in anisotropic surface temperature distribution. The thermal surface emission can also become anisotropic. The anisotropy of observed surface radiation is used to infer magnetic field strength and geometry, composition of the surface layers, global parameters of the star (such as mass and radius), as well as important parameters of superdense matter in stellar interiors (see, e.g., [5,6] and references therein).
Here I study a simple model of thermal emission from magnetized neutron stars with isothermal interiors, assuming the star emits a blackbody (BB) radiation with local effective surface temperature T s , that varies over stellar surface under the effects of dipole surface magnetic field. This model has been invented long ago by Greenstein and Hartke [7], elaborated in the literature (e.g. [8][9][10][11][12][13]) and reviewed in [6,14]. Here I point out some properties of the surface emission which, to the best of my knowledge, have not been studied in the literature.
Logically, this paper continues the previous one [15] (hereafter Paper I) which shows that simple thermal spectra of magnetized neutron stars can be accurately approximated by two-BB (2BB) models. Sect. 3 outlines theoretical formalism. Sect. 4 is devoted to the effects of chemical composition of the blanketing envelopes and Sect. 5 extends the solution to the case, in which magnetic poles contain additional hot spots. The conclusions are formulated in Sect. 6.
It is important to mention more complicated neutron star emission models. They include sophisticated magnetic field effects on thermal emission, leading to specific spectral, polarization and angular properties of radiation (see, e.g. [6,[16][17][18]), which are not reproduced by the given model.

General formalism
Calculation of radiation spectral flux from a spherical neutron star under formulated assumptions is simple. Here is the summary using the notations of Paper I. Any small surface temperature element emits like a BB, and the flux is obtained by intergating over a vizible part of the surface, taking intro account that emitted quanta propagate in Schwarzschild space-time and demonstrate gravitational redshift of photon energies and light bending. The effects of General Relativity are specified by the compactness parameter x g = r g /R, where r g = 2GM/c 2 is the Schwartschield radius, M is the gravitational star's mass, R is the circumferential radius, G is the gravitational constant and c is the velocity of light. One deals either with local quantities at the stellar surface (e.g., the local surface temperature T s ) or with the quantities detected by a distant observer. The letter ones will be often marked by the symbol '∞'. For instance, T ∞ s = T s 1 − x g is the redshifted surface temperature. As an exclusion, we will denote local (non-redshifted) photon energy by E 0 , and the redshfied energy by Let F ∞ E be a radiative spectral flux density [erg cm −2 s −1 keV −1 ] detected at a distance D R.
It is customary to express F ∞ E as where H ∞ E is the effective flux, that is formally independent of D. It can be calculated as an integral of the radiating surface flux over the visible part of the surface, where σ SB is the Stefan-Boltzmann constant, k B is the Boltzmann constant, dΩ s is a surface solid angle element, and P is the light bending function (e.g. [19], [20], [5] and [21]). It is convenient to integrate over the star's surface and calculate the bolometric luminosity of the star, L s , as well as the averaged non-redshifted effective surface temperature T eff (e.g. [10]), For a uniform surface temperature (T s = T eff ), one immediately gets the standard BB flux, and the bolometric effective flux

Input parameters
Eq.
(2) allows one to compute thermal spectra for any given temperature distribution T s over neutron star surface. We focus on the distribution created by a dipole magnetic field (with the field strength B pole at magnetic poles) due to anisotropic heat transport in a thin (maximum a few hundreds meters) heat blanketing envelope. The input parameters are M, R, chemical composition of the blanketing envelope, and a nonredshifted temperature T b at its bottom (at density ρ b ∼ 10 10 g cm −3 ; see e.g. [14]). The local T s is usually determined by solving local quasistationary radial heat transport within the envelope mediated by an effective radial thermal conductivity. For studying the thermal surface emission it is profitable to use T eff [see Eq. (3)] instead of T b .

2BB representation
According to Paper I, the spectral fluxes H ∞ E , computed from Eq. (2) for iron heat blankets, are accurately fitted by a familiar 2BB model, Here H BB∞ E (T eff ) is given by Eq. (4); 'c' and 'h' refer, respectively, to colder and hotter BB components.
Any fit contains four parameters, which are two effective temperatures T effc and T effh , and two fractions of effective radiating surface areas, s c and s h . Instead of T effc and T effh , it is convenient to introduce two dimensionless parameters p c = T effc /T eff and p h = T effh /T eff , with p c < p h . In Paper I the fits have been done for a number of representative values of M, R, log T eff [K] (from 5.5 to 6.8), log B pole [G] (from 11 to 14), photon energies (0.064 < E 40 keV, removing those E at which the fluxes are negligibly small) and inclination angles i (between line of sight and the magnetic axis). Typical relative fit errors have not exceeded a few per cent, meaning the fits were really good, providing excellent analytic representation of original computed data. Therefore, thermal X-ray spectral fluxes of magnetized neutron stars within the given model are nearly identical to 2BB spectral models.
Moreover, as shown in Paper I, it is sufficient to calculate the fluxes H ∞ E for two observation directions which are (i) pole observation, i = 0, to be denoted as H ∞ E , and (ii) equator observations, i = 90 • , to be denoted as H ⊥∞ E . If these fluxes are known, the radiation flux in any direction i is accurately approximated as We will see (Sect. 5) that that this approximation is more restrictive than Eq. (6). This accurate representation of numerically calculated H i∞ E by a map of four fit parameters (p c , p h , s c , s h ) in space of input parameters has to be taken with the grain of salt. The problem is that the fluxes H i∞ E are close to 1BB fluxes (4), which leads to some degeneracy of fit parameters (with respect to a fit procedure and a choice of grid points). On the other hand, the assumed surface temperature distribution model and the local BB emission model are definitely approximate by themselves, so that too accurate fitting of H i∞ E is actually purely academic. However, these results can be helpful for interpretation of observations, especially because the 2BB model is often used by observers.

Iron and accreted heat blankets
Paper I considered heat blankets made of iron. Here we employ a more general model [10,11,22], in which a heat blanket consists of shells (from top to bottom) of hydrogen, helium, carbon and iron. Light elements (H, He, C) are thought to appear because of accretion of H and/or He and further burning into C. Iron is either formed at neutron star birth or is a final result of carbon burning. The mass of lighter elements ∆M in the heat blanket is treated as a free parameter. The heat blanketing envelope may fully consist of iron (Paper I) or contain some mass ∆M 10 −6 M of lighter elements. Lighter elements affect the surface temperature distribution and thermal surface spectral fluxes.
To simplify the task, we will study fully accreted blankets and conclude on the partly accreted ones in the end of Sect. 4.2.

T s distributions
Our surface temperature distribution T s is axially symmetric and symmetric with respect to the magnetic equator. For illustration, Fig. 2 presents calculated T s for a star with M = 1.4 M and R = 12 km [x g =0.344 and the surface gravity g s = GM/(R 2 1 − x g ) = 1.59 × 10 14 cm s −2 ]. The effective surface temperatures are shown versus colatitude ϑ, which varies from ϑ = 0 at the north pole, to ϑ = 90 • at the magnetic equator, and then to ϑ = 180 • at the south pole. Three panels (a-c) correspond to three values of the magnetic field B pole = 10 12 , 10 13 and 10 14 G at the pole. Each panel displays T s (ϑ) for six values of log T eff [K]= 5.8, 6.0, . . . 6.8. Solid lines are plotted for the heat blankets made of iron, while dashed lines are for the fully accreted blankets. The iron and accretted matter have different thermal conductivities and, hence, different T s (ϑ) distributions at the same T eff . The dotted lines correspond to a non-magnetic star to guide the eye (they give T s = T eff ).
Heat conduction is mainly radiative in the outer non-degenerate layer. In deeper layers, where the electrons become degenerate, heat is mostly transported by electrons. The field affects also thermodynamic properties of the matter, for instance, the pressure (e.g. [23]). The magnetic field effects, which regulate the anisotropy of the surface temperature distribution, are twofold.
The effects of the first type are the classical effects of electron rotation about magnetic field lines. They are especially important near the equator, where the heat, emergent from the stellar interiors, propagates mainly across field lines and is thus suppressed. This may strongly enhance thermal insulation of the equatorial part of the heat blanketing envelope, producing rather narrow equatorial dips of T s (ϑ) (the cold equatorial belt, Fig. 2). The dips for the accreted matter are seen to be much stronger than for the iron matter, especially at lower T eff .
The effects of second type are associated with Landau quantization of electron motion across B-lines. They become more efficient at higher B near magnetic poles, where they make the heat insulating layers more heat transparent and tend to increase T s . These relatively warmer polar 'caps' are wide (Fig. 2); T s increases inside them, if we approach a pole, but not dramatically.
As long as B pole is rather weak (B pole 10 11 G), the classical effects stay more important, producing a pronounced colder equatorial belt, but weakly affecting the polar zones. The surface becomes overall colder, compared with the field-free star with the same internal temperature T b (e.g., see Figs. 22 and 23 in [14]). With increasing B pole , the quantum effects become dominant, creating hot polar zones and making the surface overall hotter (at the same T b ). At B pole 10 12 G, the surface becomes hotter than at B pole = 0; hot polar zones make the equatorial belt less significant.
It is remarkable ( Fig. 2) that outside the equatorial belt, at fixed T eff , the surface temperature profiles T s (ϑ) for iron and accreted heat blankets are mainly close. Moreover, these profiles almost coincide at B pole ≈ 10 13 G. The latter conclusion is independent of specific values of M and R.

Spectral fluxes
Despite of dramatic interplay of the opposite effects and substantial temperature anisotropy, the surface-integrated spectral fluxes H ∞ E of thermal radiation behave smoothly under variations of B pole and/or T eff . This is because of strong effects of General Relativity, which smooth out photon propagation from the surface to the distant observer [8,10,11]. As noted in the cited papers, General Relativity 'hides' magnetic effects on thermal surface emission. Very large magnetic fields, with B pole 10 14 (typical of magnetars), start to affect the surface emission much stronger but such fields are not considered here. They deserve a special study. In particular, their heat-insulating envelopes can become thick and the 1D approach to calculate the heat insulation within them can be questionable.
We have checked that the spectral fluxes, calculated for the accreted envelopes, possess the same properties (Sect. 3.3) as for the iron envelopes. In particular, the flux remains to be close to that calculated in the 1BB approximation at uniform surface temperature T s = T eff (see, e.g., Fig. 1 in Paper I). It is disappointing, meaning that such a flux, taken at T eff as an input parameter, carries little information on the magnetic field. Then one should study deviations from the 1BB approximation. It is important that if one fixes T b and varies B pole , one obtains noticeable T s and flux variations (as demonstrated, for instance, in Figs. 23 and 24 of [14]) but these variations are accompanied by variations of T eff . This would be a good theoretical construction, if one knew T b . However, the observer can measure T eff and wishes to infer T b , which is not easy because at fixed T eff the theoretical spectral flux is slightly sensitive to  All three fluxes for each T eff look close in the logarithmic format of Fig. 3. However, they do not coincide (as clearly seen from Fig. 3 of Paper I). Note that the difference between the fluxes emitted Fit parameters from the star with accreted and iron envelopes is noticeably smaller than the difference between these fluxes and 1BB ones. The difference from 1BB fluxes monotonically increases with growing B pole and E, which is quite understandable. Higher B pole produces stronger anisotropy of T s distribution. Radiation at higher energies is predominantly emitted from hotter places of the surface.
Let us remark that at B pole = 10 13 G the spectral fluxes from the star with iron and accreted heat blankets are almost indistinguishable (but sufficiently different from the 1BB flux). This is because the T s distrbutions for iron and accreted blankets are very close (Fig. 2, Sect. 4.1). They differ only in cold equatorial belts, but the contribution into the fluxes from cold narrow belts appears almost neglegible (also see Paper I). In contrast, the contrbution from hotter surface regions can be important (Sect. 5). At B pole = 10 14 G the situation remains nearly the same as at 10 13 G. Fig. 4 shows similar spectral fluxes but for equator observations (and only for B pole = 10 11 and 10 13 G). In this case the situation is similar to that for pole observations.
The main outcome of Figs. 3 and 4 is that the difference of the fluxes emitted from stars with iron and accreted heat blanketing envelopes is rather small and can be ignored in many applications, especially taking into account approximate nature of heat blanketing models. The results for partly accreted heat blankets would be similar. Accordingly, one can neglect chemical composition of heat blankets for calculating thermal emission from magnetized neutron stars in the adopted model of heat blankets. Therefore, the chemical composition of heat blankets does affect cooling of the neutron stars but almost unaffects thermal emission (under formulated assumptions).

Phase-space maps
As discussed in Sect  Fig. 5 should have been identical as well. However, these curves differ. The difference is especially visible for B pole = 10 11 and 10 12 G, but less visible at higher B pole . The difference is mainly explained by two things. First, weak variations of 2BB fit parameters lead to weaker variations of the fluxes. Second, the 2BB approximation is accurate but not exact by itself, as discussed in the previous section.
It is important that we fit unabsorbed calculated spectral fluxes using unabsorbed 2BB models. Accordingly, as argued in Paper I, one can compare 2BB parameters, inferred from observations (corrected for interstellar absorption and non-thermal radiation component), with the theoretical parameters of unabsorbed 2BB models.

Extra heating of magnetic poles
Sect. 4 describes thermal emission of cooling neutron stars with dipole magnetic field near their surfaces, using the T s model calculated in [10,11]. According to the theory, radiative spectral fluxes should be almost insensitive to the chemical composition of the heat blanketing envelopes; 2BB fits to their thermal X-ray spectral fluxes do not give the difference between the temperatures T effh and T effc of the hotter and colder BB components higher than 20% (for B pole 10 14 G). A brief comparison with observations of isolated middle-aged neutron stars in Paper I shows that there are no reliable candidates for such objects at the moment; the difference between T effh and T effc is actually larger, at least ∼ 50% for RX J1856.5−3754 [24,25], as a promising example.
Here we propose a possible phenomenological extension of the discussed model, which may simplify interpretation of observations of some sources. Specifically, let us assume the presence of additional heating of magnetic poles, which produces polar hotspots and raises the pole temperature. Such an extension has been used, for instance, in Paper I but assuming heating of both poles. Now we consider heating either of both poles or one of them. For simplicity, this heating is assumed to be axially symmetric. It does not violate the axial symmetry of the T s distribution. Note that two equivalent hotspots do not destrpy the initial symmetry of the north and south hemispheres, whereas one hotspot does destroy it. Possible physical justifications for extra heating are mentioned in Sect. 6.
Let T s0 (ϑ) be the basic effective surface temperature used in Sect. 3. We introduce a small phenomenological angle ϑ 0 that deteminies the size of a hotspot. Following Paper I we assume that at ϑ < ϑ 0 (in the nothern hotspot), and T s (ϑ) = T s0 (ϑ) outside the hotspot. The parameter δ specifies an extra temperature enhancement at the magnetic pole; the enhancement smoothly disappears as ϑ → ϑ 0 . This will be our model with one hotspot. Assuming similar temperature enhancement near the second magnetic pole, we will get the model with two hot spots. The presence of spots renormalizes the total effective temperature T eff , Eq. ( 3). Otherwise computation of spectral fluxes is the same as in Sect. 3. These fluxes can also be approximated by 2BB fits (6), and can be analyzed via phase-space maps. In some cases presented below the relative the fit accuracy becomes worse (reaching sometimes ∼ 10 %) but the fit remains robust. We will again take the star with M = 1.4 M and R = 12 km. For certainty, we use the model of iron heat blanket and ϑ 0 = 10 • .

Extra hotspots on both poles
The left panel of Fig. 6 demonstrates calculated spectral fluxes for B pole = 10 14 G and seven values of log T eff [K]=5.6, 5.8,. . . , 6.8 at δ = 0.3. Since these hotspots are identical, observations of the first and second poles (i = 0 and 180 • ) give the same fluxes. Such fluxes are shown by dashed lines. For comparison, the solid lines are calculated at δ = 0, in which case the extra hotspots are absent and the results of Sect. 3 apply. The dotted lines show the fluxes for non-magnetic star. One can see that the presence of hotspots enhances the spectral fluxes at high photon energies. This is expected: enhanced surface temperature of hotspots intensifies generation of high-energy radiation. The higher T eff , the larger photon energies are affected. The solid and dotted lines are the same as in Fig. 3c.
The right panel of Fig. 6 shows the 2BB fit parameters versus δ for the same star with T eff =1 MK. At δ = 0, the results naturally coincide with those in Figs. 5g and 5h. However with increasing δ, the fit parameters become different. The temperature T effh of the hotter BB component and the effective emission surface area s c of the colder component noticeably increase, whereas the emission surface area of the hotter component dramatically falls down. Even with really small hotspot temperature enhancement δ ≤ 0.3 (which gives the fraction of extra luminosity ≤1.4%), one gets absolutely new phase-space portrait with s h s c . Such 2BB fits have been often inferred from observations of cooling isolated neutron stars (see, e.g., Paper I); these sources are usually interpreted as neutron stars, which have small hotspots with noticeably enhanced temperature. Therefore, the theory with hotspots predicts (e.g. Paper I) two types of neutron stars, whose spectra can be approximated by the 2BB models. The first ones are those with smooth surface temperature distributions, created by nonuniform surface magnetic fields (as considered in Sect. 3) to be called spectral 2BB models of smooth magnetic atmospheres. The second sources are those with distinct hotspot BB component to be called 2BB with hotspots. Naturally, there is a smooth transition between them (for instance, by increasing δ in the right panel of Fig. 6). It seems that the observations do not provide good candidates for the sources of the first type (Paper I), but there are some candidates for the sources of the second and intermediate types as we outline later.

Phase resolved spectroscopy in case of two hotspots
Let us present a few calculated lightcurves to illustrate an important problem of pulse fraction. For simplicity, we consider the star as orthogonal rotator with the spin axis perpendicular to the line of sight. Fig. 7 presents several lightcurves for a star with B pole = 10 14 G versus phase angle θ. The lightcurves are normalized by H E (0) (the spectral flux for pole observations at given energy E). The displayed ratio of redshifted or non-redshifted fluxes is the same at a given E, so that we drop the symbol ∞. Note that the normalization flux H E (0) strongly depends on energy by itself. Panel (a) corresponds to δ = 0, in which case the hotspots are actually absent and we have the emission from the smooth magnetic atmosphere. In case (b) the hotspots with δ = 0.3 are available and strongly affect the lightcurves.
Solid lines on each panel in Fig. 7  1 keV are although smooth and resemble those on panel (a). However, at higher energies the lightcurves (b) have shapes with erased dips at nearly equator observations [and then Eq. (7) becomes inaccurate although Eq. (6) works well]. These curves are typical for lightcurves produced by antipodal point sources on the neutron star surface; see, e.g. Fig. 4 in [19] or Fig. 6 in [21]. Clearly, at high energies the star emits from the poles, the hottest places on the surface. Nevertheless, the pulse fraction remains not too high (≤ 30%) in both cases (a) and (b) even at the highest taken energy E ≈ 6.3 keV, at which thermal spectral flux becomes typically very low by itself. Therefore, the enhancement of thermal emission by extra heating of two magnetic poles does not lead to sizable pulse fractions. The dotted lines on Fig. 7 show the same lightcurves as the solid lines but calculated neglecting gravitational bending of light rays. One sees that without the light bending the pulse fraction would be much stronger (as is well known; see e.g. [8]). The difference of solid and dotted curves on Fig. 7 has simple explanation. Without light bending, the equator observations (θ = 90 • or 270 • ) would collect more light emitted from the cold equator. This would lower the observed emission and produce stronger dips of the lightcurves.

Extra hotspot on one pole
Finally, consider the case of a single hotspot, which we put at the north pole. This case is different from the previous one.
The left panel of Fig. 8 presents spectral fluxes calculated for the same values of T eff and B eff as in Fig. 6. The dashed lines are for the north-pole observations (i = 0) at δ = 0.3. These fluxes are the same as in Fig. 6: the observer sees a large fraction of the surface (mainly the northern hemisphere); the southern polar region is unseen.
The solid line in the left panel of Fig. 8 refers to the same case of δ = 0.3 but the star is observed from the south pole (i = 180 • ). Then the observer cannot see the hotspot, so the line exactly coincides with the solid line in Fig. 3c. The dotted line in Fig. 8 is for a non-magnetic star (as in Figs. 3c and 6). Now, with only one extra heated pole, the difference between observations of the north and south poles is pronounced.
The right panel of Fig. 8 presents the maps of 2BB fit parameters versus δ for the same star at T eff = 1 MK for the north-pole, south-pole, and equator observations (to be compared with Fig. 6). The maps the for north-pole and equator observations are nearly the same as those in Fig. 6, but the maps for the south-pole observations are plain. Corresponding fit parameters are independent of δ just because the observer cannot see the emission from the hotspot (which is the only emission which increases with δ).

Phase resolved spectroscopy in case of one hotspot
Phase-resolved lightcurves are plotted in Figs. 9 and 10. They are naturally different from the case of two extra heated magnetic poles. Fig. 9 is analogous to Fig. 7 and shows the light curves at the same T eff , B pole , and extra heating parameter δ = 0.3 (a) and δ = 0.14 (b). With one hotspot, the fractions of extra luminosities are 0.7% (a) and 0.3% (b). The phase-space variations are similar to those in Fig.  7 only at low photon energies, in which case the extra polar heating is almost unnoticeable and the lightcurves have two shallow dips over one rotation period. However with growing E, the extra heating becomes important and increases the depth and width of the dips (such shapes are known to be produced by point sources on neutron star surface, e.g. Fig. 4 in [19]). At E 4 keV in Fig. 9a the pulse fraction becomes close to 100%. In that case the observer would see only one, but very pronounced dip over one rotation, because the extra hotspot would be poorly vizibe in a wide range of θ with the center at θ = 180 • . In Fig. 9b the heating is weaker, the dips are smaller, but nevertheless they are much larger than for two heated spots at δ = 0.3 in Fig. 7. Evidently, one extra-heated magnetic pole is much more profitable than two, for producing large pulse fractions. Fig. 10 is plotted for B pole = 10 11 G to illustrate the effects of B pole on the pulse fraction. Fig. 10a presents the normalized lightcurves for the star with δ = 0.3. Now the pulse fraction is lower than at B pole = 10 14 G. Nevertheless, it remains much higher than it would be if the additional heating were switched on at both poles. Fig. 10b demonstrates that at lower δ = 0.14 at B pole = 10 11 G the pulse fraction is weaker affected by the extra heating but nevertheless reaches about 70% at E ∼ 6 keV.

Discussion and conclusions
Several isolated middle-aged neutron stars which radiation spectra have been (or can be) fitted by 2BB models have been selected in Paper I, based on the recent catalog [26] presented also on the Web: //www. ioffe.ru/astro/NSG/thermal/cooldat.html. These can be cooling neutron stars with dipole magnetic fields and additionally heated polar caps. This selection is illustrative (does not pretend to be complete).
1. XMMU J172054.5−372652 is a neutron star that is probably associated with the SNR G350.1−0.3 [27]; there is no direct evidence of pulsations. Potekhin et al. [26] used archival Chandra data and fitted the X-ray spectrum with a neutron star ( ) atmosphere model assuming M = 1.4M and R = 13 km. They obtained T eff ≈ 2 MK, but did not perform two-component fits which would be interesting.
2. PSR B1055−52 (J1057−5226) is a moderately magnetized middle-aged pulsar. Its characteristic effective magnetic field, reported in the ATNF pulsar catalog [28], is B eff = 1.1 × 10 12 G. Potekhin et al. [26] present the value k B T ∞ eff ≈ 70 eV (T ∞ eff ≈ 0.8 MK) based on the 2BB spectral fit by De Luca et al. [29]; the fit includes also a power-law (PL) non-thermal radiation component. Potekhin et al. [26] have corrected the results of [29] for the distance estimate to the source made by Mignani et al. [30] and reported T effh /T effc ∼ 2.3 with s h s c . It seems worth to try to explain these data with a 2BB model containing hotspots. 3. PSR J1740+1000 possesses the magnetic field B eff = 1.8 × 10 12 G. The 2BB spectral fit was performed in [31]. With the same version of the fit as selected in [26], one has k B T ∞ eff ≈ 70 eV (T ∞ eff ≈ 0.8 MK), T effh /T effc ∼ 2.8 and s h s c , with the same conclusion as for the PSR B1055−52.
4. PSR B1823−13 (J1826−1334), located in the SNR G18.0−00.7, has B eff = 2.8 × 10 12 G. Its X-ray emission is mainly non-thermal [32,33] but has some thermal component. The 1BB+PL fit gives the radius of thermally emitting region R ∞ eff ≈ 5 km, smaller than the expected radius of a neutron star. Adding the second BB component is statistically insignificant with the present data, but might be possible in the future.
5. RX J1605.3+3249 (RBS 1556) is a neutron star studied by many authors (e.g. [34][35][36][37][38]) with contradictory conclusions on its properties (as detailed in [26]). There is no solid evidence of stellar rotation [37]. Its X-ray emission has been analyzed using BB fits and neutron star atmosphere models. Recently Pires et al. [37] analysed the XMM-Newton observations and Malacaria et al. [38] jointly analysed the NICER and XMM-Newton data. These teams improved 2BB fits by adding a broad Gaussian absorption line, in which case they got k B T ∞ effc ∼ 60 eV (T ∞ effc ∼ 0.7 MK), and T effh /T effc ∼ 2, a good opportunity to assume a dipole magnetic field and additionally heated poles.
6. RX J1856.5−3754 is a neutron star with nearly thermal spectrum. It was discovered by Walter et al. [39]. It is rotating with the spin period of ∼7 s; the effective magnetic field is B eff ∼ 1.5 × 10 13 G; magnetic properties are highly debated (e.g. [18,40] and references therein). The spectrum has been measured in a wide energy range, including X-rays, optical and radio, and interpreted with many spectral models, particularly, with the model of thin partially ionized magnetized hydrogen atmosphere on top of solidified iron surface (see, e.g., [5,41]). Note alternative 2BB, 2BB+PL, and 3BB fits constructed by [24] and [25]. The 2BB fits give k B T ∞ effc ∼ 40 eV (T ∞ effc ≈ 0.46 MK), T effh /T effc ∼ 1.6, R effh ∼ 0.5R and R effc ∼ R; they seem closer to the 2BB spectral models with smooth magnetic atmosphere, than 2BB fits for other sources. Extra heating of magnetic poles would simplify this interpretation.
This brief analysis does not reveal any good candidate which would belong to the family of neutron stars with dipole magnetic fields and smooth surface temperature distribution (Sect. 4). This does not mean that such neutron stars do not exist; it can be difficult to identify them because their spectra are close to 1BB spectrum. To interpret the above sources one needs to complicate the model, for instance, by introducing extra heating of magnetic poles. Now let us summarize our results. Following Paper I we have studied simple models (Sect. 3) of thermal spectra emitted from surfaces of isolated neutron stars with dipole surface magnetic fields 10 11 B pole 10 14 G. Such fields make the surface temperature distribution noticeably non-uniform. The model assumes BB emission with a local temperature from any surface element.
In Sect. 4 we have shown that the spectral X-ray fluxes emitted from such neutron stars are almost the same for heat blaketing envelopes composed of iron and fully acctereted matter (at the same average effective surface temperatures T eff ). By measuring the spectral fluxes and T eff , it would be difficult to infer composition of their heat blanketing envelopes, although this composition can noticeably affect neutron star cooling (e.g., [14] and references therein). In the presence of fully accreted matter, the fluxes are accurately fitted by 2BB models (as in Paper I for the iron heat blankets). At M ∼ 1.4 M and R ∼ 12 km, the ratio T effh /T effc of effective temperatures of the hotter to colder BB components cannot be essentially larger than ∼ 1.2. For less compact stars, with smaller M/R, this value can be somewhat higher (see e.g. [8], Paper I). If this ratio is larger from observations, then the model of Sects. 3 and 4 cannot explain the data.
In Sect. 5, following Paper I, we have extended the model by assuming sufficiently weak extra heating of one or both magnetic poles. Increasing extra relative temperature rise δ at the pole(s), does not greatly violate the accuracy of the 2BB spectral approximation, but allows one to obtain larger T effh /T effc . Also, it affects the pulse fraction, which increases with B pole and photon energy E. The pulse fraction can be very high if one magnetic pole is additionally heated, while the other is not.
Let us stress, that there have been many other studies of thermal emission from magnetized neutron stars. The emission produced by anisotropic surface temperature distribution of the star with dipole magnetic fields was analyzed by Page [8]. Quadrupole magnetic field components were added in [9]. A possible presence of toroidal field component in the neutron star crust was studied in [12,13]. Also, there were studies of thermal emission during magnetic field evolution in neutron stars (e.g., [18,40]). These works give a variety of magnetic fields configurations, surface temperature distributions, and phase-resolved spectra of neutron stars to be compared with observations and determine simultaneously both, the surface temperature distribution and magnetic field geometry.
Introducing some extra heating of magnetic poles, we are actually doing the same. Note that hot spots can be produced by pulsar mechanism. Also, the toriodal crustal magnetic fields can significantly widen the cold equatorial belt, creating relatively stronger heated poles [12], which can mimic hotspots. At the first step it might be simpler to introduce phenomenological polar heating, that can be different on the two poles, and to determine the parameters of the extra heating (the T s distribution) by varying the parameters like δ and extra temperature profile [Eq. (8)] at both poles. This would be a simple and physically transparent model of surface temperature distribution. At the next step, one could infer (constrain) the magnetic field geometry. This method can be easily extended to different T s distributions, not necessarily axially symmetric.
Another very important direction would be to go beyond the approximation of local BB emission from any surface patch, and use more realistic emission models (with spectral features and anisotropic radiation, as well as with account of polarization effects, e.g. [42,43]).
Funding: This research was partly (in analyzing the case of accreted matter) supported by the Russian Science Foundation, grant 19-12-00133.

Data Availability Statement:
The data underlying this article will be shared on reasonable request to the corresponding author.