Polarized radiation from an accretion shock in accreting millisecond pulsars using exact Compton scattering formalism

Pulse profiles of accreting millisecond pulsars can be used to determine neutron star (NS) parameters, such as their masses and radii, and therefore provide constraints on the equation of state of cold dense matter. Information obtained by the Imaging X-ray Polarimetry Explorer (IXPE) can be used to decipher pulsar inclination and magnetic obliquity, providing ever tighter constraints on other parameters. In this paper, we develop a new emission model for accretion-powered millisecond pulsars based on thermal Comptonization in an accretion shock above the NS surface. The shock structure was approximated by an isothermal plane-parallel slab and the Stokes parameters of the emergent radiation were computed as a function of the zenith angle and energy for different values of the electron temperature, the Thomson optical depth of the slab, and the temperature of the seed blackbody photons. We show that our Compton scattering model leads to a significantly lower polarization degree of the emitted radiation compared to the previously used Thomson scattering model. We computed a large grid of shock models, which can be combined with pulse profile modeling techniques both with and without polarization included. In this work, we used the relativistic rotating vector model for the oblate NS in order to produce the observed Stokes parameters as a function of the pulsar phase. Furthermore, we simulated the data to be produced by IXPE and obtained constraints on model parameters using nested sampling. The developed methods can also be used in the analysis of the data from future satellites, such as the enhanced X-ray Timing and Polarimetry mission.


Introduction
Rapidly rotating accreting neutron stars (NSs) can be used to study the properties of the extremely dense matter inside their cores.Some of these are accretion-powered millisecond pulsars (AMPs), which show X-ray pulsations at the stellar spin frequency.Their emission profiles can be modeled, accounting for light bending and other relativistic effects, and from these profiles, the information about the mass and radius of the NS can be extracted(see e.g., Pechenick et al. 1983;Miller & Lamb 1998;Poutanen & Gierliński 2003;Poutanen & Beloborodov 2006;Morsink et al. 2007;Lo et al. 2013;Watts et al. 2016;Salmi et al. 2018;Bogdanov et al. 2019).The mass and radius constraints can be translated to the equation of state (EOS) constraints for matter inside the NS core (e.g., see Lindblom 1992;Lattimer 2012;Hebeler et al. 2013;Baym et al. 2018).
This technique was recently applied to the case of rotationpowered millisecond pulsars using data from the Neutron Star Interior Composition Explorer (NICER, Miller et al. 2019;Riley et al. 2019;Miller et al. 2021;Riley et al. 2021;Salmi et al. 2022).In these sources, the hot regions on the NS surface are heated by a magnetospheric return current and the emission escaping the NS surface can be computed using selfconsistent atmosphere models (see e.g., Zavlin et al. 1996;Ho & Lai 2001;Heinke et al. 2006;Ho & Heinke 2009;Haakonsen et al. 2012;Salmi et al. 2020).For AMPs, these models are not adequate because the emitted radiation is Comptonized in an accretion "shock" formed above the hot region.In reality, a standard gasdynamical shock might not exist at all, but the kinetic energy of incoming particles is released in the surface layer through Coulomb collisions and plasma instabilities (Zel'dovich & Shakura 1969;Suleimanov et al. 2018;Salmi et al. 2020).Comptonization of soft photons generated within the layer as well as by the underlying NS surface is then responsible for the power-law spectra extending up to 100 keV as observed from AMPs (Poutanen & Gierliński 2003;Gierliński & Poutanen 2005;Falanga et al. 2005aFalanga et al. ,b, 2007Falanga et al. , 2011Falanga et al. , 2012)).Electron scattering also causes the radiation to be significantly polarized.Variations of polarization with pulsar phase can be used to determine geometrical parameters such as inclination angle and magnetic obliquity (Viironen & Poutanen 2004), which in their turn allow us to improve constraints on the NS mass and radius.
Previously, polarized emission models for AMPs used the Thomson scattering approximation in an optically thin NS atmosphere (see Sunyaev & Titarchuk 1985;Viironen & Poutanen 2004;Salmi et al. 2021).However, the polarization degree due to Compton scattering on hot electrons strongly depends on the electron temperature (Nagirner & Poutanen 1994;Poutanen 1994); for example, for 100 keV electrons, polarization is less than 50% of that for Thomson scattering.On the other hand, the models for nonpolarized radiation -in the context of pulse profile modeling-typically employ approximate formulas for the anisotropy of the radiation and empirical models for the Comptonized spectra (see e.g., Poutanen & Gierliński 2003;Leahy et al. 2008;Steiner et al. 2009;Salmi et al. 2018).Self-consistent accretion-heated atmosphere models have also been developed (see e.g., Zampieri et al. 1995;Deufel et al. 2001;Suleimanov et al. 2018), although these are relatively computationally expensive and have a large number of free parameters.
In this work, we applied the formalism for Compton scattering in a hot slab (Nagirner & Poutanen 1993;Poutanen & Svensson 1996) to compute the Stokes parameters of the emergent radiation as a function of energy and emission angle using only three model parameters: the electron temperature and optical depth of the hot slab on top of the NS surface, and the temperature of the seed blackbody photons coming from the NS.By employing pulse profile modeling for polarized radiation from rapidly rotating oblate NSs (Poutanen 2020;Loktev et al. 2020), we also simulated data for the Imaging X-ray Polarimeter Explorer (IXPE, Weisskopf et al. 2022) and updated the NS geometry parameter constraints previously predicted in Salmi et al. (2021) using the Thomson scattering model.These methods can also be applied when analyzing the observations from IXPE1 or from future X-ray polarimetric missions such as the enhanced Xray Timing and Polarimetry mission (eXTP; Zhang et al. 2019;Watts et al. 2019).The atmosphere look-up tables produced here can also be combined with any AMP pulse profile modeling codes. 2 .
The remainder of this paper is structured as follows.In Sect.2, we present the theory that describes the formation of polarized radiation in a Comptonizing slab above the NS surface.We then apply this theory to obtain the intensity and polarization of the radiation escaping from the slab.We then describe the method to obtain pulse profiles using our new radiation model.We compare our resulting spectra and pulse profiles to those obtained using previous models.We also discuss pulse profiles obtained for different NS parameters.In Sect.3, we generate synthetic data and apply our fitting routine to determine the NS parameters from the data.We discuss the applications of our model and predictions for the upcoming observations in Sect. 4. In Sect.5, we conclude and summarize our findings.

Radiative transfer equation
We consider a simple NS atmosphere model in which the atmosphere is a plane parallel slab consisting of electrons and lying above an optically thick source that radiates as a blackbody.The hot slab has the Thomson optical depth τ T = σ T n e H, where H is the vertical height, σ T is the Thomson cross section, and n e is the electron concentration.The electron gas in the atmosphere is considered to be isotropic and isothermal with the electron temperature T e .The momentum distribution of electrons is given by relativistic Maxwellian distribution characterized by the dimensionless temperature Θ e = kT e /m e c 2 : where γ is the electron Lorentz factor and K 2 is the modified Bessel function of the second kind.This distribution is normalized to unity: where p = γ 2 − 1 is the electron dimensionless momentum.The radiation field can be described by the Stokes vector I(τ, x, µ), where µ is the cosine of the zenith angle, an angle between the normal to the slab and the direction of photon propagation, and x = E/m e c 2 is the photon energy measured in the units of the electron rest mass.
The Stokes vector usually contains four components I, Q, U, and V, but because of the azimuthal symmetry and the absence of sources of circular polarization, the latter two are equal to zero.We therefore use just two Stokes parameters I(τ, x, µ) = (I, Q) T , where superscript T denotes the transposed vector.
The photon distribution at the bottom of the slab is considered to be Planckian of the temperature T bb .The incident radiation (for zeroth scattering order, n = 0) at the bottom of the slab is given by the Stokes vector where Θ bb = kT bb /m e c 2 .In order to describe the propagation of polarized radiation through the hot electron slab, we solve the radiative transfer equation (RTE) (Nagirner & Poutanen 1994;Poutanen & Svensson 1996): where dτ = σ T n e dz is the Thomson optical depth, S is the source function (also a Stokes vector), and σ(x) is the dimensionless Compton scattering cross-section (in units of the Thomson cross-section σ T ).We solve the RTE using the iterative scattering method of Poutanen & Svensson (1996), where the intensity is represented as a series expansion in scattering orders.The Stokes vector of unscattered radiation at all optical depths is From the known Stokes vector at the n-th scattering order, we find the source function for the following scattering order as where R is the 2 × 2 azimuth averaged redistribution matrix describing Compton scattering by isotropic electrons (see Appendix A1 in Poutanen & Svensson 1996).If the source function is known, the Stokes vector is found via the formal solution of the RTE: Iterations proceed until the required accuracy of the total Stokes vector I = n I n is achieved so that the maximal contribution of the next scattering is less than 1% of the total spectrum in all energies and angles.Finally, we can compute the polarization degree (PD) of the emergent radiation as  The positive values of P imply that the polarization vector lies in the meridional plane defined by the slab normal and the photon momentum, while for the negative values, the polarization vector is perpendicular to that plane, as in the case of optically thick electron scattering atmosphere (Chandrasekhar 1960;Sobolev 1963).

Spectro-polarimetric properties of accretion shocks
We can now apply the derived formalism to predict the beaming patterns, energy spectra, and polarization of the accretion shock at the NS surface.We take as a fiducial set of parameters:  2021), we note that the PD differs significantly while the spectrum remains mainly unchanged.This effect is most obviously seen in the angular dependence of the PD.The absolute values of the PD are reduced by a factor of two for the photon energies of 2, 5, and 8 keV, and increased at higher energies.For Compton scattering, the PD changes the sign already at the photon energy of 8 keV and is positive at higher values, while for Thomson scattering the sign changes at photon energies above 18 keV.Last but not least, the maximum in the absolute value of the PD for the photon energies higher than 8 keV has shifted from µ = 0 to higher values of the cosine of the zenith angle.Keeping in mind the energy range of the IXPE satellite (2-8 keV), we can say that compared to the previously developed model for AMPs that used Thomson scattering, the Compton scattering model predicts a lower PD.This result is in agreement with our expectations, as for hot electrons the relativistic aberration becomes important (Poutanen 1994).In the case of cold electrons, the PD of the Thomson singly scattered radiation is highly dependent on the scattering angle, reaching 100% at 90 • .For relativistic electrons (but still in the Thomson scattering regime), the same happens in the electron rest frame; however, because of the relativistic motion of the electron, the photon in the external frame is preferentially scattered in the direction of the electron motion, resulting in zero final PD independent of the scattering angle (Bonometto et al. 1970;Nagirner & Poutanen 1993).In our case of 50 keV electrons, the electrons are not yet relativistic and the maximum PD is reduced to 60% from the Thomson case.We can predict here that with the higher electron temperature, the PD will decrease further.This will be examined in the following section.

Grid of models
The model we consider for the accretion shock has three parameters: the optical depth of the Comptonizing plasma τ T , the dimensionless electron temperature Θ e , and the dimensionless characteristic photon temperature Θ bb (associated with T e and T bb , respectively).To determine these parameters from the upcoming IXPE data, we need a dense grid of models.For the optical depth of the shock τ T , we chose the lower and higher value of 0.5 and 3.5, respectively, and varied the parameter in steps of 0.1.For the electron temperature, we varied Θ e from 0.04 to 0.2 (i.e., T e from ∼20 to 100 keV) with a step of 0.004 (i.e., ∼2 keV).For the seed photon temperature, the dimensionless Θ bb was varied in the interval (1-3)×10 −3 (i.e., T bb =0.5-1.5 keV) in steps of 2 × 10 −4 (i.e., ∼0.1 keV).Altogether, we computed 13981 models.
The impact of the NS atmosphere parameters on the PD can be seen in Fig. 3.For the electron temperature, we see the expected decrease in the PD with T e as electrons become more relativistic.For the two other parameters, the main feature is that with increasing the parameter value, the change of sign of the PD happens at higher energies.Within the mentioned intervals of parameter values, the impact is most significant with increasing the optical depth of the atmosphere and least significant with changes to the seed photon energy.For the IXPE observational range of 2-8 keV, this means that, in general, PD increases with the increase in these parameters.
In order to use the calculated Stokes parameters more effectively, we wrote a 3D interpolation routine over this grid of parameters, so that for any arbitrary combination of these parameters within the intervals mentioned, we can obtain the intensity and PD of the emission as functions of energy and zenith angle.We tested that the results of this interpolation match the direct calculation results for any combination of the parameters within the grid to within 8% (within 1% in the IXPE energy range).

Pulse profiles and polarization
We calculated the pulse profiles from the known beaming pattern and PD exactly as in Salmi et al. (2021).After obtaining the beaming pattern and PD of the emission on top of the NS atmosphere, we transport them from the hot emitting regions of the star (called hot spots) to the observer following the formalism described in Morsink et al. (2007), AlGendy & Morsink (2014), Poutanen (2020), andLoktev et al. (2020).We consider an oblate, rapidly rotating star in the Schwarzschild metric (hence, the effects of rotation on space-time are neglected).Gravitational redshift and the light-bending effect, as well as the effect of the fast rotation of the star, are accounted for in the polarization angle (PA) transfer.The only difference in approach is the reference radius used in time-delay calculation.While Loktev et al. ( 2020) use the exact radius of the star at the center of the spot, we follow the path described in Poutanen & Beloborodov (2006) and applied in Salmi et al. (2021) and use the equatorial radius as a reference radius.
We again return to the results from Salmi et al. (2021) and compare the pulse profiles coming from the model with Compton scattering in the atmosphere with previously obtained ones using Thomson scattering.Figure 4 presents this comparison for the case of two antipodal spots, both combined and individually from each spot, for the NS with atmospheric parameters from the fiducial set (inclination i = 60 • and magnetic co-latitude θ = 20 • as in Salmi et al. 2021).We see that the flux does not change significantly, which is in agreement with the beaming patterns not changing either.The PD is now decreasing from 2 to 5 keV and then slowly increasing towards 8 keV.From inspection of the corresponding lines in Fig. 1, and accounting for the gravitational redshift, it is clear that the PD is following the same pattern for the corresponding energies of the emitted photons.The 'jump' by 90 • in the PA for the 8 keV photons comes from the PD changing sign.
The main conclusion here is that we expect a significant drop in the observed PD compared to that predicted by Salmi et al. (2021).The impact of this reduction in PD on the possibility of observing the AMPs is further discussed in Sects.3 and 4.

Generating the synthetic data
The developed model can be further used to analyze the data.Polarimetric measurements of the AMPs are yet to be carried out.In the absence of observations, we used the ixpeobssim framework (Baldini et al. 2022) (version 16) with the pcube algorithm to generate several samples of synthetic data imitating the data to be obtained from the IXPE observations of an AMP.We substituted the source with the model described above and obtained the event list for a specified observing time of 600 ks for a source with a flux of 100 mCrab.In this part of our study, we do not investigate the energy dependency of PD, and collect all the data in one energy bin.These data sets are in exactly the same form as publicly available IXPE data products, so we will be able to apply the developed routines to the real data as soon as they become available.
As mentioned above, applying the Compton scattering formalism reduces the expected PD as compared to the previous studies where Thomson scattering in the NS atmosphere was considered.We synthesized the data for the same parameter set as used in Fig. 4 (see Table 1, Model 1) and compared the PD to the so-called minimal detectable polarization (MDP) at 99% significance (see Weisskopf et al. 2010 for a definition of the MDP and Kislat et al. 2015 for further details) values for the corresponding phase bins.Figure 5a illustrates that the obtained values of the PD are below the MDP.According to Weisskopf et al. (2010), this means that we cannot have a proper measurement of the PD in this case.From this point, there are three main ways to improve the quality of the potential measurement: (i) we can increase the number of counts to reduce the MDP values, (ii) study different energy bands to see if the PD values are higher there, and/or (iii) choose a more promising source, which in terms of simulated data means selecting different source parameters.
Increasing the number of counts by prolonging the exposure time is impossible with AMPs, as they are transient sources observable only during the outbursts, which typically last for a few weeks (see e.g., Patruno & Watts 2021).Also, observations lasting longer than 600 ks are unrealistic for scheduling reasons.The two other options are explored below.
Studying different energy bands is relevant as, according to our predictions from Sect.2.1, the PD changes sign within the 2-8 energy band.This means that averaging over all energies decreases the total PD value.Figure 5b illustrates the results of analyzing the 2-5 keV energy band instead of the 2-8 keV represented in Fig. 5a.The blue line corresponding to the model prediction is now significantly higher and closer to the MDP values.However, the improvement is slight, and we note that the PD values are still below the MDP values in this case.
In our search for a more optimistic (from the point of view of potential observations) yet still valid set of both geometrical Fig. 3. Angular dependency of the PD for different atmosphere parameters.For each plot, we vary only one parameter (marked on the plot) while other parameters are from the fiducial set (τ T = 1, T e = 50 keV, T bb = 1 keV).Cyan, blue, green, orange, and red solid lines show the PD on photon energies of 2, 5, 8, 12, and 18 keV, respectively.
Table 1.Parameter sets for all the computed models and the most probable values obtained in the fitting.

Parameter
Injected values Fitting results 1.9 +1.0 −0.8 1.6 +0.9 −0.6 2.0 +0.9 −0.8 2.2 +0.8 −0.9 2.1 +0.8 −0.7 Notes.NS mass, radius, and spin frequency are fixed at M = 1.4 M ⊙ , R = 12 km, ν = 401 Hz.These parameters, as well as Model 1 parameters, are selected based on the SAX J1808.4−3658pulsar and expected values for canonical NSs.The errors represent the 68% credible interval.Additional parameters χ and ∆ϕ are the position angle of the pulsar rotation axis and the phase shift, respectively.As before, i is inclination, θ is co-latitude of the spot, T e and T bb are the temperatures of the electron in the slab and seed photons, respectively, and τ T is the Thomson optical depth of the slab.Values are averaged over two data realizations for each scenario. (a) Columns labeled 1a and 3a correspond to fitting the whole energy band of IXPE, 2-8 keV. (b) Columns labeled 1b and 3b correspond to the data realization generated and fitted in the 2-5 keV energy band.
and atmosphere parameters, we aimed to reproduce the known spectra of the AMPs with the parameters that, coming from our previous investigations, are expected to give a higher PD of the emission.From Fig. 3, we see that the optical depth affects the  2021) (combined from two spots).The NS parameters in both models are those shown in Table 1 of Salmi et al. (2021).
PD strongest from all the atmosphere parameters, and so we aimed to increase it but preserve the overall spectrum behavior presented in Poutanen & Gierliński (2003).We also attempted to increase the electron temperature to match the predictions in Poutanen & Gierliński (2003), but the scenario proposed there gives PD values below the MDP values, and so we chose not to continue with it.We developed two possible adjustments to the fiducial scenario (Model 1): we can choose the geometry such that the angle between the inclination and co-latitude of the spot is almost 90 • (as in Model 2) or use 'more promising' atmosphere parameters (as in Model 3).We can also combine these two (Model 4).The parameter values for these four models are given in the left part of Table 1 in the columns labeled correspondingly 1-4. Figure 5c shows that for Model 4, the PD now exceeds the MDP values.We also checked this for Models 2 and 3, and for Model 2 the PD is slightly higher than the MDP values, while for Model 3 the PD is similar to MDP.In all the models, we only calculated the emission coming from one spot on the surface of the NS.
In addition, we attempted to improve the quality of the results for Model 1 illustrated in Fig. 5a,b by reducing the number of phase bins, as the PD varies less in these cases compared to other models.We find that for Model 1, case a, only reducing the number of bins to one gives us the PD value that barely exceeds the MDP, meaning that we can detect something from such a source, but we would lose the phase variability completely.For Model 1, case b, we can obtain a detection with three phase bins, but only for the bin that includes the peak (in Fig. 5b the corresponding phase shift is 0.5), and again, the PD barely exceeds the MDP.We conclude that reducing the number of bins can support the observations with barely measurable PD, but its impact is not sufficient for us to rely upon it in cases where the source is otherwise unobservable with IXPE.
Nevertheless, the most recent development of the IXPE data analysis shows that for multiple-point observations (i.e., observations split in different phase bins), information about polarization can be extracted even from the data points where PD values do not exceed MDP (González-Caniulef et al. 2023;Suleimanov et al. 2023).As IXPE measures the Stokes parameters q and u rather than PD, we can fit these parameters directly and obtain some constraints on the geometry of the source.This approach is further explored in the following section. .Simulated PD and normalized Stokes q and u profiles for the NS parameter set 1 (panels (a) and (b) are for the 2-8 keV and 2-5 keV energy range, respectively) and set 4 (c) from Table 1.The exposure time is 600 ks in all three cases.The blue lines show the theoretical model, the purple dots are the synthesized observed data, and the red dots are the MDP values for the corresponding phase bins.
Posterior probability distributions for NS parameters for Model 1a (left panels) and Model 4 (right panels) from Table 1 using one data realization.The two-dimensional contours correspond to 68% and 95% credible regions with the corresponding regions colored black and gray.The histograms show the normalized one-dimensional distribution for a given parameter derived from the posterior samples.The mean value and the 68% and 95% credible levels for the parameters are shown by the vertical dashed black, red, and orange lines, respectively.Blue lines indicate the true values of the parameters that were used for the data simulations.

Fitting the data
After synthesizing the data, we applied the multinest (Feroz et al. 2009) multimodal nested sampling algorithm (specifically, the pymultinest package, see Buchner et al. 2014) to obtain new constraints on both the geometrical and atmosphere parameters of the source.As in Salmi et al. (2021), we fit the geometrical parameters, such as inclination, co-latitude, phase shift, and position angle, but unlike there, we also perform a simultaneous fit over the parameters of the atmosphere, such as characteris-tic photon temperature, electron temperature, and optical depth of the thermalization layer.Table 1 presents the most probable values and 68% credible intervals for these parameters.Figure 6 shows the results of the fitting presented for the fiducial set of parameters (Model 1) and for a 'more promising' scenario (Model 4).In both cases, the full energy band of IXPE (2-8 keV) is considered, and the exposure time is 600 ks.
Figure 6 clearly illustrates the improvement of the fitting results if comparing the data synthesized from the Model 1-type of source (left) and the data coming from the Model 4-type of source (right).As we use the uniform prior for all our fittings, in the former case, the constraints originate from our direct fitting of the q and u Stokes parameters.This illustrates the possibility to extract some information even when the MDP values are not reached.In the latter case, we can clearly see the significant constraints coming from the PD values exceeding the MDP by approximately a factor of two.
We can also conclude that the constraints on the atmospheric parameters are notably less restrictive than the ones on the geometrical parameters of the NS.The main reason for obtaining such broad constraints on atmospheric parameters lies in the energy range of the IXPE satellite, as we are unable to detect the whole spectrum and fit it correctly.A solution to this issue could be to use the observational results from several missions simultaneously; for instance, fitting the NuSTAR data in its higher energy range to determine the electron temperature (and consequently the optical depth through the spectral slope).Then, using the values obtained, we can fit the IXPE data with the fixed (or constrained) atmospheric parameters.However, we see that even with large uncertainties on these parameters, we are able to derive constraints on the geometry of the source from the IXPE data only.
Table 1 also shows that the constraints from the Model 3 data are less restrictive than the ones from Model 2, and so we can conclude that the atmospheric parameters are less important in our search for the 'promising' source than the geometric parameters of the NS.We also attempted to improve the constraints by reducing the number of simultaneously fitted parameters (e.g., we fixed the atmospheric parameters and only fitted the geometrical parameters).Decreasing the number of fitted parameters reduces the calculation time, but the credible intervals for the fitted parameters remain the same.
In the fittings presented here, we fixed the mass and radius of the NS.This is done for the simplicity of the calculations, as they are quite time-consuming.Moreover, we only fit the normalized Stokes q and u parameters, and they are not particularly sensitive to the mass or radius.We could also fit the Stokes I parameter (i.e., flux) as well and add mass and radius as free parameters.Adding two more parameters would increase the computational time significantly, but we do not expect any significant impact of this on the credible intervals of other parameters.

Discussion
Our fiducial set of NS parameters is based on the SAX J1808.4−3658, and Fig. 4 shows that, for the parameters selected, our predictions for the PD of the SAX J1808.4−3658-likesource are low.However, with the uncertainties we still have on both the geometry and the atmospheric parameters of the source, there is still a possibility to observe a higher PD.Furthermore, with about 20 AMPs discovered, we can put forward several conclusions about the perspective of observing measurable PD from this class of sources.
It is known that the product of the electron temperature and optical depth for the AMPs is almost invariant (Poutanen 2006).Nevertheless, we have learned that all the atmospheric parameters have their own contribution to the PD values; see Fig. 3. Optical depth is the most important parameter here; we are searching for a source with an atmosphere that is sufficiently optically thick for many Compton scattering events to occur, yet not so thick that the pulsation pattern is excessively influenced by these scatterings.In this sense, objects such as IGR J17591−2342 and IGR J17511−3057 are promising.The former has an estimated optical depth of 1.59-2.3and a pulsed fraction of the emission of 10%-17% (Kuiper et al. 2020;Manca et al. 2023).The latter has a lower optical depth of 1.34, but higher electron and seed photons temperatures (T e = 51 keV, T bb = 1.36 keV) and a pulsed fraction of 14% (Papitto et al. 2010).Constraints on the geometrical parameters of the AMPs are among the goals of the polarimetric observations.So far our main source of information on the possible geometry is the amplitude of the pulsation, which depends approximately on the product of the sines of the inclination and the magnetic co-latitude (Poutanen & Beloborodov 2006).Moreover, while for the inclination we usually have some constraints, little is known about the geometrical configuration of the magnetic field of the AMP.From the emission patterns shown in Fig. 1, the highest PD is received when the line of sight and the normal to the emitting area are orthogonal.One of the most interesting objects to study is Swift J1749.4−2807, for which the inclination is well-constrained, 74.4 • < i < 77.3 • , and the co-latitude is estimated to be θ ≈ 50 • (Altamirano et al. 2010).In this case, we can expect high polarization and we can test whether the inclination obtained from the polarimetry actually matches the orbital one.
We performed all our fittings with 600 ks of exposure time.It is typical for AMPs to go into outburst for a couple of weeks, and so in order to collect enough photons we need either a rather soft spectrum of emission or a very bright source.We can exclude hard sources such as IGR J18245−2452 with a photon index of Γ ≈ 1.4 (Papitto et al. 2013).
Figure 4 illustrates an additional feature of the PD behavior we expect to observe: the observed PD is different between cases where the secondary spot is invisible (e.g., the very beginning of the outburst) and where it is visible (at the later stages of the outburst).In order to get additional information about the geometry of the source, it is crucial to observe the outburst from the earliest to the latest possible stages.
We purposely focused on phase-resolved analysis and did not investigate the possibilities that energy-resolved analysis could open up.IXPE has an energy resolution of 0.57 keV (at 2 keV, Ramsey et al. 2021), and so with sufficient statistics it is possible to study the behavior of PD with energy and obtain more information about the geometry and atmospheric properties of the source.
The model presented in this paper can be improved further.For instance, the isothermal (static) slab approximation of the accretion shock can be replaced with a more physical model where the dynamics of the accreting gas is computed self-consistently with the electron temperature.Also, the model for seed soft photons can be modified to account more accurately for the reprocessing of the hard X-rays at the NS surface, accounting for the polarization properties of this radiation.We intend to work on these developments in the future.

Summary
We present a new model to compute the flux and polarization of the emission coming from the accretion-powered millisecond pulsars.We used the formalism for Compton scattering in a plane-parallel isothermal hot slab to describe the propagation of photons in the accretion shock above the NS surface.Figures 1 and 2 show the result of such calculations for the fiducial set of parameters taken from Salmi et al. (2021).Comparing our results with the previous studies using Thomson scattering formalism shows that, at least for the same NS parameters, our model predicts lower PD values.
We calculated the Stokes parameters of the radiation as a function of two variables, the zenith angle of the emitted photon and its energy, for almost 14,000 combinations of three parameters: electron temperature, the temperature of the seed blackbody photons, and the Thomson optical depth of the slab.The tables of the spectral energy distribution of the Stokes parameters for all the parameters can now be used for further analysis of the accretion-powered millisecond pulsars, both in studying the pulse profiles and the polarimetric properties of the emerged radiation.Figure 3 illustrates one of the ways to use this new data set.We can understand the impact of the NS atmosphere parameters on the angular dependency of the PD from this particular study.
We then used the relativistic rotating vector model for the oblate NS to produce observed Stokes parameters as a function of the pulsar phase.Figure 4 shows the observed fluxes coming from the two antipodal spots of an AMP and the combined total flux as calculated using our model and a total flux calculated using a previously developed Thomson model.A significant drop in the PD between these two models is again present here.
We used the developed model of the AMP emission to generate the synthetic data sets imitating the one that the IXPE satellite could produce after observing an AMP.We note that the NS with the parameters from our fiducial set would emit the light with PD values below the MDP values of IXPE (which still does not exclude the possibility to extract information about the geometry of the source); however, for a different set of NS parameters, the PD is high enough to be detectable by IXPE, as shown in Fig. 5.We then applied the multimodal nested sampling technique to fit the simulated data with our model and obtain constraints on both the geometrical and atmospheric parameters of the NS.The results of the fitting are presented in Fig. 6.The main conclusion of this study is that we can obtain reasonably good constraints on the geometrical parameters from the IXPE data for some of the AMPs; however, the parameters of the atmosphere of the NS need to be constrained from the supporting observations by other X-ray missions such as NuSTAR.
Our findings contribute to the advancement of understanding AMPs and have significant implications for future research in this field.We are waiting for the AMPs to be observed by IXPE so that we may apply the developed formalism to analyze the resulting data.We also hope to see the eXTP mission observing the AMPs, as it would solve several of the issues discussed in this article: the larger effective area will lower the MDP, and so the fainter sources will become observable; the broader energy range will allow us to develop better constraints on the atmospheric parameters; and finally, the wide field instrument of eXTP will be able to catch the transients, such as AMPs, and observe the outbursts from their very beginning.

Fig. 1 .
Fig.1.Intensity of the emergent radiation (left) and the PD (right) as functions of the cosine of the zenith angle.The atmosphere parameters are τ T = 1.0, electron temperature T e =50 keV, and seed photon temperature T bb =1 keV.Black, blue, green, orange, and red solid lines show the model for photon energies of 2, 5, 8, 12, and 18 keV, respectively.The angular dependence of the intensity is normalized so that 1 0 µI(µ)dµ = 1/2.

Fig. 2 .
Fig.2.Intensity of emergent radiation (left) and the PD (right) as functions of photon energy.The same atmosphere parameters are used as in Fig.1.Black, blue, green, orange, red, and magenta solid lines show the model for µ = 0.0, 0.2, 0.4, 0.6, 0.8, and 1.0 respectively.
50 keV, T bb = 1 keV.Figures1 and 2show the spectrum and the PD as functions of the cosine of the zenith angle and energy, respectively, for this set of parameters.The parameters of the source are the same as those used in the Thomson model ofSalmi et al. (2021) (see their Figs.1 and 2), and so these figures can be used to illustrate the impact of using an exact Compton scattering description instead of the Thomson scattering model in the slab of hot electrons.When comparing our results with those ofSalmi et al.  (

Fig. 4 .
Fig. 4. Pulse profiles of the observed flux, PD, and PA for two antipodal spots shown for three different energies (2, 5, 8 keV).The solid black curves correspond to the total flux (PD, PA), the dashed blue lines correspond to the contribution from the primary spot, and the dotted red lines are for the secondary spot.The dash-dotted gray lines correspond to the Thomson model fromSalmi et al. (2021) (combined from two spots).The NS parameters in both models are those shown in Table1ofSalmi et al. (2021).
Fig.5.Simulated PD and normalized Stokes q and u profiles for the NS parameter set 1 (panels (a) and (b) are for the 2-8 keV and 2-5 keV energy range, respectively) and set 4 (c) from Table1.The exposure time is 600 ks in all three cases.The blue lines show the theoretical model, the purple dots are the synthesized observed data, and the red dots are the MDP values for the corresponding phase bins.