An Equation of State for Magnetized Neutron Star Matter and Tidal Deformation in Neutron Star Mergers

The extreme properties of neutron stars (NSs) not only open up many possibilities related to the composition, structure, and dynamics of stable cold matter in the observable universe but also to the matter interaction at the fundamental level (Rezzolla et al. 2018). Almost every aspect of them, be it mass, radius, rotational fRequency, or magnetic field, represents matter at extreme conditions. The structure of the NS depends on the nuclear equation of state (EOS), which was poorly known so far (Baiotti 2019; Perot et al. 2019; Malik et al. 2019; Hebeler et al. 2013). Although the observation of high-mass pulsars like PSR  J1614–2230 (M = 1.908 ± 0.016M_⊙; Arzoumanian et al. 2018) and PSR J0348–0432 (M = 2.01 ± 0.04 M_⊙; Antoniadis et al. 2013) severely constrains the EOS as well as their interaction, it also questions the possible presence of exotic matter in them. Apart from the observational constraints from high-mass stars, the observed gravitational wave event GW 170817 of a binary NS merger in August 2017 disfavors some of the stiff EOSs (Abbott et al. 2017; Malik et al. 2018). The precise knowledge of the NS radius can also constrain the behavior of EOSs. Of late, NICER (Neutron star Interior Composition Explorer mission) has come up with a measurement of the radius, ${12.71}_{-1.19}^{+1.14}\,\mathrm{km}$, for the NS with mass ${1.34}_{-0.16}^{+0.15}\,{{\rm{M}}}_{\odot }$ (Riley et al. 2019), and other independent analyses show that the radius is ${13.02}_{-1.06}^{+1.21}\,\mathrm{km}$ for NSs with mass $M={1.44}_{-0.14}^{+0.15}\,{M}_{\odot }$ (Miller et al. 2019). However, the empirical estimates of the radius of a canonical NS ($M=1.4{M}_{\odot }$) should be ${R}_{1.4}=(11.9\pm 1.22)\,$ km (Lattimer & Lim 2013). Recently, Fattoyev et al. (2018) and Malik et al. (2018), using the extracted bounds on the NS tidal deformability of the GW 170817 event, suggested that R_1.4 ≲ 13.76 km.

Born out of massive interstellar gases over millions of years, the magnetic field present in these compact structures can be very high, although the origin of these high fields is still not well understood. The typical values of the surface magnetic field (B_surf) of NSs range from 10¹² to 10¹⁵ G. It is speculated that the field intensity can be even more at the core. A fraction of the population having the strongest surface magnetic fields of ∼(10¹²–10¹⁵) G are called magnetars and generally belong to soft gamma repeaters (SGRs) or anomalous X-ray pulsars (AXPs; Gomes et al. 2017). Typical examples of such magnetars are 1E 1048.1–5937 and 1E 2259+586 with surface magnetic field B_surf ∼ 10¹⁴ G (Melatos 1999), 4U 0142+61 (B_surf ∼ 10¹⁶ G; Makishima et al. 2014), and SGR 1806−20 (B_surf  ∼ 10¹⁴ G; Kouveliotou et al. 1998), etc. One can refer to the magnetar catalog available online for more examples (Olausen & Kaspi 2014). It is conceivable, though perhaps speculative, that in the interior of magnetars, the magnetic field could be several orders of magnitude larger and therefore is expected to affect the dense matter properties on the scale of QCD (Lai & Shapiro 1991). Therefore, it is desirable to incorporate the magnetic field effects to determine the composition and gross structural properties of these NSs with a high magnetic field. Various authors have incorporated the effects of a magnetic field to account for the properties of NSs: Ciolfi & Rezzolla (2013, 2012), Broderick et al. (2002), Huang et al. (2010), Lopes & Menezes (2012), Casali et al. (2014), Gomes et al. (2014, 2017), Gao et al. (2015), Orsaria et al. (2011), Dexheimer et al. (2014), Franzon & Schramm (2015), Chu et al. (2015), Paulucci et al. (2011), and Strickland et al. (2012). From various investigations, it is well known that strong magnetic fields affect the energy levels of charged particles due to Landau quantization, which may be strong enough to make the pressure of the matter anisotropic, and therefore, spherical symmetry may not be a suitable approximation in order to study the structural properties of NSs (Perez Martinez et al. 2008; Strickland et al. 2012). However, it has been noted that the difference in NS properties (such as mass and radius), calculated separately for the pressure in the parallel and perpendicular directions to the magnetic field, being small, spherical symmetry can still hold (Chu et al. 2015; Huang et al. (2010). In most of these calculations, the divergent vacuum contribution is omitted. Indeed, for vanishing magnetic fields, such a "no sea approximation" leads to a very small difference in the EOS compared to the EOS calculated by taking into account the Dirac sea effect after renormalization (Glendenning 1988, 1989). We shall, however, include here the effects of the magnetic field on the Dirac sea of nucleons. In fact, the inclusion of magnetic field effects for the Dirac vacuum has been the reason for the magnetic catalysis of chiral symmetry breaking in quark matter and has been studied in various effective models like the Nambu–Jona–Lasinio (NJL) models (Menezes et al. 2009; Chatterjee et al. 2011) as well as the quark–meson models (Ferrari et al. 2012; Andersen & Khan 2012; Skokov 2012). Apart from the static structural properties of magnetars, magnetic fields also play an important role in the physics of compact star mergers. The gravitational waves emitted at the late stage of the merger process can possibly be detected directly and are sensitive to the EOS of the dense matter (Read et al. 2013). The magnetic field in a merger process can possibly become extremely large due to magnetorotational instability, and the magnitude can be large enough to affect the EOS of dense matter.

When a long-lived (toroidal) electric current flows in a highly conductive NS matter, the magnetic pressure corresponding to the Lorentz force comes into play. This induces the deformation of stars (Konno et al. 1999). The first discussion of the deformation of magnetized NSs was given by Chandrasekhar and Fermi and by Ferraro. It was discussed in Bocquet et al. (1995) that NSs having magnetic field B > 10¹⁴ G get deformed, and the maximum permitted magnetic field is around 10¹⁸ G for stability. Further, it was observed that the relative deformation is ≈2% when the magnetic field is taken to be approximately 10¹⁵ G. Under these conditions, we consider that it is adequate to assume spherical symmetry for NS structure (Gomes et al. 2019). The anisotropy of magnetic pressure in an NS stimulates a deformation in NSs, which is studied in Schramm & Schramm (2014).

With this motivation, in the present work, we incorporate the effects of a strong magnetic field on the EOS and calculate the NS properties using a model based on a generalization of the sigma model. In such a model, the nucleons are coupled with the σ and pion (${\boldsymbol{\pi }}$) fields along with a potential for the $\sigma ,{\boldsymbol{\pi }}$ fields. Such a model, generalized to include vector mesons coupled to the scalar fields to give masses to vector mesons, was used to study finite temperature aspects of nuclear matter and its applications to NSs (Glendenning 1986). Later on, it was generalized to include the isovector rho meson as well as higher-order nonlinear meson interactions (Sahu & Ohnishi 2000; Jha et al. 2006). The reason was to include the effects of isospin asymmetry as well as to explain the rather high value of nuclear incompressibility that one gets when nonlinear meson interactions are not taken into account. Further, the model has been generalized recently by Malik et al. (2017) to include cross-coupling between isovector and isoscalar mesons. These cross-couplings were found to be instrumental in explaining the density dependence of the nuclear symmetry energy as well as its slope and curvature parameters at the saturation density deduced from a diverse set of experimental data. The EOS was also used to explore the gross structural properties of NSs like mass and radius, which turned out to be consistent with the measurement for the maximum mass while the radius at the canonical mass is within the empirical bounds (Malik et al. 2017). The present investigation of the effects of a strong magnetic field is carried out within this effective chiral model.

The paper is organized as follows. In Section 2, we first recapitulate the essential features of the effective chiral model, and in the next subsection we introduce the magnetic field in the model and calculate in some detail the pressure and energy density in the presence of an external magnetic field. As we shall see, in the presence of a strong magnetic field, the EOS can become anisotropic when the induced magnetization effects become strong. After deriving the EOS, in Section 3, we solve the Tolman–Oppenheimmer–Volkoff (TOV) equations using the EOS thus derived in the presence of a magnetic field. We also consider here the tidal deformation of NSs in the context of NS mergers. The results along with brief discussions are presented in Section 4. Finally, we summarize and present an outlook of the present investigation in Section 5.

We now discuss briefly the salient features of the effective chiral Lagrangian of Malik et al. (2017) to describe dense nuclear matter. The effective Lagrangian of the model interacting through the exchange of the pseudo-scalar meson π, the scalar meson σ, the vector meson ω, and the isovector ρ−meson is given by

$$\begin{eqnarray}&&\begin{array}{l}{ \mathcal L }={\bar{\psi }}_{B}\,\left[\left(i{\gamma }_{\mu }{\partial }^{\mu }-{g}_{\omega }{\gamma }_{\mu }{\omega }^{\mu }-\displaystyle \frac{1}{2}{g}_{\rho }{{\boldsymbol{\rho }}}_{\mu }\cdot {\boldsymbol{\tau }}{\gamma }^{\mu }\right)\right.\\ \left.-{g}_{\sigma }\,\left(\sigma +i{\gamma }_{5}{\boldsymbol{\tau }}\cdot {\boldsymbol{\pi }}\right)\Space{0ex}{2.85ex}{0ex}\right]\,{\psi }_{B}\\ +\displaystyle \frac{1}{2}\left({\partial }_{\mu }{\boldsymbol{\pi }}\cdot {\partial }^{\mu }{\boldsymbol{\pi }}+{\partial }_{\mu }\sigma {\partial }^{\mu }\sigma \right)-\displaystyle \frac{\lambda }{4}{\left({x}^{2}-{x}_{0}^{2}\right)}^{2}\\ -\displaystyle \frac{\lambda b}{6{m}^{2}}{\left({x}^{2}-{x}_{0}^{2}\right)}^{3}-\displaystyle \frac{\lambda c}{8{m}^{4}}{\left({x}^{2}-{x}_{0}^{2}\right)}^{4}\\ -\displaystyle \frac{1}{4}{F}_{\mu \nu }{F}_{\mu \nu }+\displaystyle \frac{1}{2}{{g}_{\omega B}}^{2}{x}^{2}{\omega }_{\mu }{\omega }^{\mu }-\displaystyle \frac{1}{4}{{\boldsymbol{R}}}_{\mu \nu }\cdot {{\boldsymbol{R}}}^{\mu \nu }+\displaystyle \frac{1}{2}{m}_{\rho }^{{\prime} 2}{{\boldsymbol{\rho }}}_{\mu }\cdot {{\boldsymbol{\rho }}}^{\mu }\\ +{\eta }_{1}\,\left(\displaystyle \frac{1}{2}{g}_{\rho }^{2}{x}^{2}{{\boldsymbol{\rho }}}_{\mu }\cdot {{\boldsymbol{\rho }}}^{\mu }\right)+{\eta }_{2}\,\left(\displaystyle \frac{1}{2}{g}_{\rho }^{2}{x}^{2}{{\boldsymbol{\rho }}}_{\mu }\cdot {{\boldsymbol{\rho }}}^{\mu }{\omega }_{\mu }{\omega }^{\mu }\right).\end{array}\end{eqnarray} \tag{ 1 }$$

The first two lines of the above Lagrangian represents the interaction of the nucleon isospin doublet ψ_B with the mesons. In the third and fourth lines, we have the kinetic and nonlinear terms in the pseudo-scalar-isovector pion field "${\boldsymbol{\pi }}$," the scalar field "σ," and higher-order terms of the scalar field in terms of the chiral invariant combination of the two i.e., ${x}^{2}={{\boldsymbol{\pi }}}^{2}+{\sigma }^{2}$. In the fifth line, we have the field strength and the mass term for the vector field "ω" and the isovector field "${\boldsymbol{\rho }}$" meson. The last line contains the cross-coupling terms between ${\boldsymbol{\rho }}$ and ω and also between the ${\boldsymbol{\rho }}$ and σ mesons. g_σ, g_ω, and g_ρ are the usual meson−nucleon coupling strengths of the scalar, vectors, and isovector fields, respectively. Here, we are concerned only with the normal nonpion condensed state of matter, so we take 〈π〉 = 0 and also m_π = 0. The last two terms in the Lagrangian incorporates the effect of cross-couplings between ρ−σ and ρ−ω with coupling strengths ${\eta }_{1}\,{g}_{\rho }^{2}$ and ${\eta }_{2}\,{g}_{\rho }^{2}$, respectively. From our previous work, Malik et al. (2017), where we investigated the role of cross-coupling terms to constrain symmetry energy and NS properties, we concluded that the inclusion of the $\rho \mbox{--}\sigma$ cross-coupling term is sufficient to satisfy the overall properties. Hence, in the present work, we only consider the ρ−σ coupling (Malik et al. 2017) with coupling parameter η₁.

The interaction of the scalar and pseudoscalar mesons with the vector boson generates a dynamical mass for the vector bosons through the spontaneous breaking of the chiral symmetry with the scalar field attaining the vacuum expectation value x₀. Then, the mass of the nucleon (m), the scalar (${m}_{\sigma }$), and the vector meson mass (m_ω) are related to x₀ through

$$\begin{eqnarray}&&m={g}_{\sigma }{x}_{0},\,\,{m}_{\sigma }=\sqrt{2\lambda }{x}_{0},\,\,{m}_{\omega }={g}_{\omega }{x}_{0},\end{eqnarray} \tag{ 2 }$$

where $\lambda =({m}_{\sigma }^{2}-m{\pi }^{2})/2{f}_{\pi }^{2}$ and f_π = x₀ is the pion decay constant reflecting the strength of spontaneous symmetry breaking. Due to the cross-coupling between the ${\boldsymbol{\rho }}$ and σ mesons, there is a contribution to the ${\boldsymbol{\rho }}$-meson mass from the vacuum expectation value of the σ meson: i.e., ${m}_{\rho }^{2}={m}_{\rho }^{{\prime} 2}+{\eta }_{1}{g}_{\rho }^{2}{x}_{0}^{2}$.

To obtain the EOS, we revert to the mean-field procedure, where one assumes the mesonic fields to be classical and uniform mean fields while retaining the quantum nature of the baryonic field, i.e., $\langle \sigma \rangle =\sigma$, $\langle {\omega }_{\mu }\rangle ={\omega }_{0}{\delta }_{\mu 0}$, $\langle {\rho }_{\mu }^{a}\rangle$ = ${\delta }_{\mu 0}{\delta }_{a3}{\rho }_{03}$.

We recall here that this approach has been extensively used to obtain field-theoretical EOSs for high-density matter (Sahu & Ohnishi 2000; Sahu et al. 2004; Jha & Mishra 2008) and gets increasingly valid when the source terms are large Serot & Walecka (1986). The details of the present model and its attributes such as the derivation of the equation of motion of the meson fields and its EOS can be found in Malik et al. (2017). For the sake of completeness, however, we write down the meson field equations in the mean-field ansatz. Moreover, these mean-field equations remain the same even in the presence of a magnetic field, as we will discuss in the next subsection. The mean meson fields, i.e., the vector field (ω), scalar field (σ), and isovector field (${\rho }_{3}^{0}$) are determined by solving the mean field equations, which (with $Y=x/{x}_{0}$) are, respectively, given by

$$\begin{eqnarray}&&\left[{m}_{\omega }^{2}{Y}^{2}+{\eta }_{2}{C}_{\rho }{m}_{\rho }^{2}{\left({\rho }_{3}^{0}\right)}^{2}\right]{\omega }_{0}={g}_{\omega }\rho ,\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&\begin{array}{c}\begin{array}{l}(1-{Y}^{2})-\displaystyle \frac{b}{{m}^{2}{C}_{\omega }}{\left(1-{Y}^{2}\right)}^{2}+\displaystyle \frac{c}{{m}^{4}{C}_{\omega }^{2}}{\left(1-{Y}^{2}\right)}^{3}\\ +\,\displaystyle \frac{2{C}_{\sigma }{m}_{\omega }^{2}{\omega }_{0}^{2}}{{m}^{2}}+\displaystyle \frac{2{\eta }_{1}{C}_{\sigma }{C}_{\rho }{m}_{\rho }^{2}{\left({\rho }_{3}^{0}\right)}^{2}}{{C}_{\omega }{m}^{2}}-\displaystyle \frac{2{C}_{\sigma }{\rho }_{s}}{{mY}}=0,\end{array}\end{array}\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&{m}_{\rho }^{2}\left[1-{\eta }_{1}(1-{Y}^{2}){C}_{\rho }/{C}_{\omega }+{\eta }_{2}{C}_{\rho }{\omega }_{0}^{2}\right]{\rho }_{3}^{0}=\displaystyle \frac{1}{2}{g}_{\rho }({\rho }_{p}-{\rho }_{n}).\end{eqnarray} \tag{ 5 }$$

The quantities ρ and ρ_S are the baryon and the scalar density defined as

$$\begin{eqnarray}&&\rho =\displaystyle \frac{\gamma }{{\left(2\pi \right)}^{3}}\displaystyle \sum _{i=n,p}{\int }_{0}^{{k}_{F}^{i}}d{\boldsymbol{k}},\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&{\rho }_{s}=\displaystyle \frac{\gamma }{{\left(2\pi \right)}^{3}}\displaystyle \sum _{i=n,p}{\int }_{0}^{{k}_{F}^{i}}\displaystyle \frac{{m}^{\star }}{\sqrt{{{m}^{\star }}^{2}+{{\boldsymbol{k}}}^{2}}}d{\boldsymbol{k}},\end{eqnarray} \tag{ 7 }$$

where ${k}_{F}^{i}$ is the Fermi momentum of the nucleon and $\gamma =2$ is the spin degeneracy factor. In the above ${m}^{\star }={g}_{\sigma }x$ is the medium-dependent nucleon mass. ${C}_{\sigma }\equiv {g}_{\sigma }^{2}/{m}_{\sigma }^{2}$ , ${C}_{\omega }\equiv {g}_{\omega }^{2}/{m}_{\omega }^{2}$, and ${C}_{\rho }\equiv {g}_{\rho }^{2}/{m}_{\rho }^{2}$ are the scalar, vector, and isovector coupling parameters that enter into the actual computation. These parameters are given in Table 1.

The energy density and the pressure of the considered model is given by

$$\begin{eqnarray}&&\begin{array}{l}\epsilon =\displaystyle \frac{1}{{\pi }^{2}}\displaystyle \sum _{i=n,p}{\displaystyle \int }_{0}^{{k}_{F}^{i}}{k}^{2}\sqrt{{k}^{2}+{m}^{\star 2}}{dk}+\displaystyle \frac{{m}^{2}}{8{C}_{\sigma }}{\left(1-{Y}^{2}\right)}^{2}\\ -\,\displaystyle \frac{b}{12{C}_{\sigma }{C}_{\omega }}{\left(1-{Y}^{2}\right)}^{3}+\displaystyle \frac{c}{16{m}^{2}{C}_{\sigma }{C}_{\omega }^{2}}{\left(1-{Y}^{2}\right)}^{4}+\displaystyle \frac{1}{2}{m}_{\omega }^{2}{\omega }_{0}^{2}{Y}^{2}\\ +\,\displaystyle \frac{1}{2}{m}_{\rho }^{2}\left[1-{\eta }_{1}(1-{Y}^{2})({C}_{\rho }/{C}_{\omega })+3{\eta }_{2}{C}_{\rho }{\omega }_{0}^{2}\right]{\left({\rho }_{3}^{0}\right)}^{2},\end{array}\end{eqnarray} \tag{ 8 }$$

$$\begin{eqnarray}&&\begin{array}{l}p=\displaystyle \frac{1}{3{\pi }^{2}}\displaystyle \sum _{i\,=\,n,p}{\displaystyle \int }_{0}^{{k}_{F}^{i}}\displaystyle \frac{{k}^{4}}{\sqrt{{k}^{2}+{m}^{\star 2}}}{dk}-\displaystyle \frac{{m}^{2}}{8{C}_{\sigma }}{\left(1-{Y}^{2}\right)}^{2}\\ +\,\displaystyle \frac{b}{12{C}_{\sigma }{C}_{\omega }}{\left(1-{Y}^{2}\right)}^{3}-\displaystyle \frac{c}{16{m}^{2}{C}_{\sigma }{C}_{\omega }^{2}}{\left(1-{Y}^{2}\right)}^{4}+\displaystyle \frac{1}{2}{m}_{\omega }^{2}{\omega }_{0}^{2}{Y}^{2}\\ +\,\displaystyle \frac{1}{2}{m}_{\rho }^{2}\left[1-{\eta }_{1}(1-{Y}^{2})({C}_{\rho }/{C}_{\omega })+{\eta }_{2}{C}_{\rho }{\omega }_{0}^{2}\right]{\left({\rho }_{3}^{0}\right)}^{2}.\end{array}\end{eqnarray} \tag{ 9 }$$

Clearly, in the above, we have ignored the contributions of the Dirac sea of nucleons and have kept the contribution arising from the Fermi sea of nucleons with an effective mass m^⋆ given by the first terms in Equations (8) and (9).

Such a "no-sea approximation" is a reasonable approximation regarding the EOS (Glendenning 1988, 1989). However, in the presence of a magnetic field, the contribution of the Dirac sea can become significant as will be discussed in the next subsection.

2.1. The EOS with Magnetic Field

We shall consider here the effect of a magnetic field on the EOS as given in the previous subsection. We shall consider the magnetic field to be constant and to be in the z direction without loss of generality. Further, we choose here the gauge ${A}_{\mu }={\delta }_{\mu 2}{xB}$, where B is the magnitude of the magnetic field. In the presence of the magnetic field, the nucleon as well as the rho meson kinetic terms will get modified, with the derivative of the fields getting replaced by covariant derivatives. Further, we shall ignore the effects of the anomalous magnetic moments of the nucleons with the magnetic fields. In any case, the effect of such interaction gives negligible corrections to the system EOS as has been shown in Ferrer et al. (2015). In the mean-field approximation, the contribution of the nucleons to the thermodynamic potential Ω_N will depend upon the mass of the nucleons and the external thermodynamic parameters like baryon chemical potential (μ), magnetic field (B), and temperature (T). Because the baryon mass will be determined dynamically by the minimization of the thermodynamic potential, it will be dependent on these parameters implicitly. We can write the thermodynamic potential ${{\rm{\Omega }}}_{N}\,={{\rm{\Omega }}}_{N}({m}^{\star }(\mu ,B,T),\mu ,B,T)$ as

$$\begin{eqnarray}&&{{\rm{\Omega }}}_{N}={{\rm{\Omega }}}_{{N}_{\mathrm{sea}}}+{{\rm{\Omega }}}_{{N}_{\mathrm{med}}},\end{eqnarray} \tag{ 10 }$$

where ${{\rm{\Omega }}}_{{N}_{\mathrm{sea}}}({m}^{\star }(\mu ,B,T),0,B,0)$ is the free energy of the magnetized nucleons from the Dirac sea while ${{\rm{\Omega }}}_{N,\mathrm{med}}$ is the contribution from the Fermi sea. Let us consider charged nucleons, i.e., protons, first. With Landau quantization for the charged nucleons, ${{\rm{\Omega }}}_{{N}_{\mathrm{sea}}}$ is given explicitly as

$$\begin{eqnarray}&&\begin{array}{l}{{\rm{\Omega }}}_{{N}_{\mathrm{sea}}}({m}^{\star }(\mu ,B,T),B,\mu =0,T=0)\\ =\,-\displaystyle \frac{| {qB}| }{{\left(2\pi \right)}^{2}}\displaystyle \sum _{n=0}^{\infty }{\alpha }_{n}\displaystyle \int {{dp}}_{z}{\epsilon }_{n}({p}_{z}),\end{array}\end{eqnarray} \tag{ 11 }$$

where ${\epsilon }_{n}({p}_{z})=\sqrt{{p}_{z}^{2}+2n| {qB}| +{m}^{\star 2}(B,\mu ,T)}$ is the energy of the nucleon with charge q for the nth Landau level and ${\alpha }_{n}=2-{\delta }_{n0}$ is the degeneracy of the Landau level, i.e., all levels except the lowest Landau level is doubly degenerate. Let us note that ${{\rm{\Omega }}}_{{N}_{\mathrm{sea}}}({m}^{\star }(B,\mu ),B,\mu =0,T\,=\,0)$ is not a vacuum term in the strict sense as the nucleon mass still depends on the medium.

The medium contribution to the thermodynamic potential at a given temperature β⁻¹ is given by

$$\begin{eqnarray}&&\begin{array}{l}{{\rm{\Omega }}}_{{N}_{\mathrm{med}}}=-\displaystyle \sum _{n}\displaystyle \frac{| {qB}| }{{\left(2\pi \right)}^{2}\beta }\displaystyle \int {{dp}}_{z}\left[\mathrm{log}(1+{e}^{-\beta ({\epsilon }_{n}-{\mu }^{\star })})\right.\\ \left.+\,\mathrm{log}(1+{e}^{-\beta ({\epsilon }_{n}+{\mu }^{\star })})\right],\end{array}\end{eqnarray} \tag{ 12 }$$

where μ^⋆ is the effective chemical potential of the baryon in the presence of vector mean fields and is given by ${\mu }^{\star }=\mu -{g}_{\omega }{\omega }_{0}-{g}_{\rho }{I}_{3}{\rho }_{0}$. In the zero-temperature limit of the Ω_med, the antibaryonic contribution will vanish and only the particle part will contribute. Using the relation ${\mathrm{lim}}_{\beta \to \infty }\,(1/\beta )\mathrm{log}(1+{e}^{-\beta x})=-x\theta (-x)$, the integrand of Equation (12) becomes $({\epsilon }_{n}-{\mu }^{\star })\theta ({\mu }^{\star }-{\epsilon }_{n})$. The theta function restricts the integration over the variable p_z up to a maximum of ${p}_{F}^{n}=\sqrt{{\mu }^{\star 2}-{m}^{\star 2}-2n| q| B}$ for a given value of Landau level n. Further, the positive value of p_z² restricts the sum over the Landau levels up to a maximum ${n}_{\max }\,=\mathrm{Int}[\tfrac{\sqrt{{\mu }^{\star 2}-{m}^{\star 2}}}{2| q| B}]$. After the integration over p_z, the medium contribution is now given by

$$\begin{eqnarray}&&\begin{array}{l}{{\rm{\Omega }}}_{{\rm{med}}}=\displaystyle -\sum _{n=0}^{{n}_{max}}\displaystyle \frac{{\alpha }_{n}| q| B}{4{\pi }^{2}}\left[\Space{0ex}{3.65ex}{0ex}{\mu }^{\star }{p}_{F}^{n}\right.\\ \left.-\,\left({m}^{\star 2}+2n| q| B\right){\rm{log}}\left(\displaystyle \frac{{p}_{F}^{n}+{\mu }^{\star }}{\sqrt{{m}^{\star 2}+2n| {q}_{i}| B}}\right)\right].\end{array}\end{eqnarray} \tag{ 13 }$$

Now, let us discuss the sea contribution to the free energy given in Equation (11). This integral is divergent. We regularize this with dimensional regularization. This has been used earlier in the context of chiral symmetry breaking in the presence of a magnetic field (Menezes et al. 2009; Chatterjee et al. 2011) as well as in the hadron resonance gas model (Endrödi 2013). Another alternate method often used to regularize such divergent integrals is through a proper time method yielding similar results (Gusynin et al. 1994, 1995a, 1995b; Haber et al. 2014). To regularize ${{\rm{\Omega }}}_{{N}_{\mathrm{sea}}}$, one adds and subtracts a zero magnetic field sea contribution (Menezes et al. 2009). The divergent zero magnetic field part is evaluated in ${\rm{d}}=3-\epsilon ^{\prime}$, while the integral over dp_z, in the presence of a magnetic field, is evaluated in ${\rm{d}}=1-\epsilon ^{\prime}$ with $\epsilon ^{\prime} \to 0$. Such a manipulation results in

$$\begin{eqnarray}&&\begin{array}{l}{{\rm{\Omega }}}_{{N}_{\mathrm{sea}}}(B)-{{\rm{\Omega }}}_{{N}_{\mathrm{sea}}}(B=0)=-\displaystyle \frac{| {qB}{| }^{2}}{2{\pi }^{2}}\left[\Space{0ex}{3.15ex}{0ex}\zeta ^{\prime} (-1,x)\right.\\ -\,\displaystyle \frac{1}{2}({x}^{2}-x)\mathrm{log}(x)+\displaystyle \frac{{x}^{2}}{4}+\displaystyle \frac{1}{12}\mathrm{log}x\\ \left.+\,\displaystyle \frac{1}{12}\left(-\displaystyle \frac{2}{\epsilon ^{\prime} }+{\gamma }_{E}-1+\mathrm{log}\displaystyle \frac{4\pi {\mu }^{2}}{{M}^{2}}\right)\right],\end{array}\end{eqnarray} \tag{ 14 }$$

where ${\gamma }_{E}\simeq 0.577$ is the Euler–Mascheroni constant, μ is the scale related to dimensional regularization, and $\zeta ^{\prime} (-1,x)$ is the derivative of the Riemann–Hurwitz ζ− function ζ(z, x) at z = −1 and is given by Elizalde (1985),

$$\begin{eqnarray}&&\begin{array}{l}\zeta ^{\prime} (-1,x)=-\displaystyle \frac{1}{2}x\mathrm{log}x-\displaystyle \frac{{x}^{2}}{4}+\displaystyle \frac{1}{2}{x}^{2}\mathrm{log}x\\ +\,{x}^{2}{\displaystyle \int }_{0}^{\infty }\displaystyle \frac{2{\tan }^{-1}y+y\mathrm{log}(1+{y}^{2})}{\exp (2\pi {xy})-1}{dy}.\end{array}\end{eqnarray} \tag{ 15 }$$

We have abbreviated here $x\equiv \tfrac{{m}^{{\star }^{2}}}{2| {qB}| }$. Further, ${{\rm{\Omega }}}_{{N}_{\mathrm{sea}}}(B=0)$ is the (divergent) Dirac sea contribution to the thermodynamic potential at a vanishing magnetic field,

$$\begin{eqnarray}&&{{\rm{\Omega }}}_{{N}_{{\rm{sea}}}}(B=0)=\displaystyle -\frac{\gamma }{{\left(2\pi \right)}^{3}}\int d{\boldsymbol{k}}\sqrt{{{\boldsymbol{k}}}^{2}+{m}^{\star 2}}.\end{eqnarray} \tag{ 16 }$$

As noted earlier, because the zero-field vacuum contribution is known to have small effects on the EOS, we shall consider here solely the B-dependent sea contribution. This contribution given in Equation (14) is still divergent, as it has a purely magnetic-field-dependent term $\sim {B}^{2}/\epsilon ^{\prime}$. As explicitly shown in Endrödi (2013), such a divergence is taken care of by adding the pure field contribution B²/2 to the field-dependent contribution of the Dirac vacuum. The contribution is rendered finite by defining the renormalized charge q_r and renormalized magnetic field B_r through Endrödi (2013) and Haber et al. (2014):

$$\begin{eqnarray}&&{B}^{2}={Z}_{q}{B}_{r}^{2},\quad {q}_{r}^{2}={Z}_{q}^{-1}{q}_{r}^{2},\quad {q}_{r}{B}_{r}={qB}.\end{eqnarray} \tag{ 17 }$$

Once this is done, the free energy still depends upon the scale μ of renormalization. However, as explicitly shown in Haber et al. (2014), one can choose the renormalization scale such that the thermodynamic potential can be written only in terms of renormalized quantities so that

$$\begin{eqnarray}&&\begin{array}{l}{{\rm{\Omega }}}_{{N}_{{\rm{sea}}}}+\displaystyle \frac{1}{2}{B}^{2}=\displaystyle \frac{{B}_{r}^{2}}{2}-\displaystyle \frac{| {q}_{r}{B}_{r}| }{{\left(2\pi \right)}^{2}}\\ \,\displaystyle \displaystyle \displaystyle \left[{\zeta }^{{\rm{{\prime} }}}(-1,x)-\displaystyle \frac{1}{2}({x}^{2}-x){\rm{log}}(x)+\displaystyle \frac{{x}^{2}}{4}+\displaystyle \frac{1}{12}{\rm{log}}x\right]\\ \equiv \,\displaystyle \frac{{B}_{r}^{2}}{2}+{{\rm{\Omega }}}_{f}.\end{array}\end{eqnarray} \tag{ 18 }$$

In what follows we shall suppress the subscript "r" from the magnetic field and the charge but it is understood that the field and charges used are the renormalized quantities.

The Dirac vacuum contribution also affects the scalar condensate ${\rho }_{s}$ of Equation (19). The contribution to the scalar density from charged baryons of a given species "i," ${\rho }_{s}^{i}$ is given by

$$\begin{eqnarray}&&\begin{array}{l}{\rho }_{s}^{i}=\langle {\bar{\psi }}_{i}{\psi }_{i}\rangle \\ =\,\displaystyle \sum _{n}{\alpha }_{n}\displaystyle \frac{| {q}_{i}| B}{{\left(2\pi \right)}^{2}}\displaystyle \int {{dp}}_{z}\displaystyle \frac{{m}_{i}^{\star }}{{\epsilon }_{n}}\theta ({\mu }^{\star }-{\epsilon }_{n})-\displaystyle \frac{{m}_{i}^{\star }| {q}_{i}| B}{2{\pi }^{2}}\\ \left[{x}_{i}(1-\mathrm{log}{x}_{i})+\mathrm{log}{\rm{\Gamma }}({x}_{i})+\displaystyle \frac{1}{2}\mathrm{log}\left(\displaystyle \frac{{x}_{i}}{2\pi }\right)\right]\\ \equiv \,{\rho }_{s}^{i,\mathrm{med}}+{\rho }_{s}^{i,\mathrm{field}}.\end{array}\end{eqnarray} \tag{ 19 }$$

The theta function restricts the integration over the variable p_z up to a maximum of ${p}_{F,n}^{i}=\sqrt{{\mu }^{\star 2}-{m}_{i}^{\star 2}-2n| {q}_{i}| B}$ for a given value of n. The positivity of p_z² again restricts the sum over the Landau levels up to a maximum ${n}_{\max }=\mathrm{Int}[\tfrac{\sqrt{{\mu }^{\star 2}-{m}^{\star 2}}}{2| {q}_{i}B}]$. One can perform the integration of p_z analytically to obtain

$$\begin{eqnarray}&&{\rho }_{s}^{i,\mathrm{med}}=\displaystyle \sum _{n=0}^{{n}_{\max }}\displaystyle \frac{| {q}_{i}| B}{2{\pi }^{2}}{\alpha }_{n}{m}_{i}^{\star }\mathrm{log}\left(\displaystyle \frac{{p}_{F,n}^{i}+{\mu }^{\star }}{\sqrt{{m}_{i}^{\star 2}+2n| {q}_{i}| B}}\right).\end{eqnarray} \tag{ 20 }$$

The number density of the charged baryons of a given species similarly is given by

$$\begin{eqnarray}&&{\rho }_{i}=\displaystyle \sum _{n}\displaystyle \frac{2| {q}_{i}| B}{4{\pi }^{2}}\sqrt{{\mu }^{\star 2}-{m}_{i}^{\star 2}-2n| {q}_{i}| B}\end{eqnarray} \tag{ 21 }$$

The meson field equations, in presence of magnetic field are same as given in Equations (3)–(5) except that the scalar density is now given as

$$\begin{eqnarray}&&{\rho }_{s}={\rho }_{s}^{n}+{\rho }_{s}^{p}\end{eqnarray} \tag{ 22 }$$

with the neutron contribution to the scalar density being

$$\begin{eqnarray}&&{\rho }_{s}^{n}=\displaystyle \frac{{m}^{\star }}{2{\pi }^{2}}\left[{k}_{F}^{n}{\mu }_{n}^{\star }-{m}^{\star 2}\mathrm{log}\displaystyle \frac{{\mu }_{n}^{\star }+{k}_{F}^{n}}{{m}^{\star }}\right]\end{eqnarray} \tag{ 23 }$$

while the proton contribution to the scalar density is given by, with ${x}_{p}=\tfrac{{{m}^{\star }}^{2}}{2| {eB}| }$,

$$\begin{eqnarray}&&\begin{array}{l}{\rho }_{s}^{p}=\displaystyle \sum _{n=0}^{{n}_{\max }}\displaystyle \frac{| {eB}| }{2{\pi }^{2}}{\alpha }_{n}{m}^{\star }\mathrm{log}\left(\displaystyle \frac{{k}_{F,n}^{p}+{\mu }_{p}^{\star }}{\sqrt{{m}^{\star 2}+2n| {eB}| }}\right)\\ \ -\ \displaystyle \frac{{m}^{\star }| e| B}{2{\pi }^{2}}\left[{x}_{p}(1-\mathrm{log}{x}_{p})+\mathrm{log}{\rm{\Gamma }}({x}_{p})+\displaystyle \frac{1}{2}\mathrm{log}\left(\displaystyle \frac{{x}_{p}}{2\pi }\right)\right].\end{array}\end{eqnarray} \tag{ 24 }$$

Similarly, the baryon number densities of neutron and protons are given as

$$\begin{eqnarray}&&{\rho }_{n}=\displaystyle \frac{{{k}_{F}^{n}}^{3}}{3{\pi }^{2}}\quad \quad {\rho }_{p}=\displaystyle \sum _{n=0}^{{n}_{\max }}\displaystyle \frac{| {eB}| }{2{\pi }^{2}}{k}_{F,n}^{p},\end{eqnarray} \tag{ 25 }$$

where ${k}_{F,n}^{p}=\sqrt{{\mu }_{p}^{\star 2}-{m}^{\star 2}-2n| e\,B| }$ for a given value of the Landau level n.

Next we write down the EOS, i.e., the energy density and pressure in the present model in the presence of an external magnetic field. The energy density is given by

$$\begin{eqnarray}&&\epsilon ={\epsilon }_{n}+{\epsilon }_{p}+{\epsilon }_{\mathrm{meson}}+\displaystyle \frac{1}{2}{B}^{2}.\end{eqnarray} \tag{ 26 }$$

In the above, the energy density of neutrons _n is given as

$$\begin{eqnarray}&&{\epsilon }_{n}=\displaystyle \frac{1}{8{\pi }^{2}}\left[{k}_{F}^{n}{\mu }_{n}^{\star }(2{m}^{\star 2}+{{k}_{F}^{n}}^{2})-{m}^{\star 4}\mathrm{log}\left(\displaystyle \frac{{\mu }_{n}^{\star }+{k}_{F}^{n}}{{m}^{\star }}\right)\right].\end{eqnarray} \tag{ 27 }$$

The contribution of the protons, on the other hand, arises from the Dirac sea as in Equation (18) as well as the medium, the Fermi sea of protons:

$$\begin{eqnarray}&&\begin{array}{l}{\epsilon }_{p}=-\displaystyle \frac{| {eB}{| }^{2}}{2{\pi }^{2}}\left[{\zeta }^{{\rm{{\prime} }}}(-1,x)-\displaystyle \frac{1}{2}({x}_{p}^{2}-{x}_{p}){\rm{log}}({x}_{p})+\displaystyle \frac{{x}_{p}^{2}}{4}\,+\,\displaystyle \frac{1}{12}{\rm{log}}{x}_{p}\right]\\ \,+\,\displaystyle \frac{| {eB}| }{4{\pi }^{2}}\displaystyle \sum _{n}^{{n}_{max}}\left[{\alpha }_{n}({k}_{F,p}^{n}+{\mu }_{p}^{\star })+{\left({m}_{n}^{p}\right)}^{2}{\rm{log}}\left(\displaystyle \frac{{k}_{F,n}^{p}+{\mu }_{p}^{\star }}{{m}_{n}^{p}}\right)\right],\end{array}\,\end{eqnarray} \tag{ 28 }$$

where we have introduced the mass in the nth Landau label for proton as ${m}_{n}^{p}=\sqrt{{m}^{\star 2}+2n| {eB}| }$. The contribution to the energy density from the mesons arises from the potential terms of the mesons and is given by

$$\begin{eqnarray}&&\begin{array}{l}{\epsilon }_{\mathrm{meson}}=\displaystyle \frac{{m}^{\star 2}}{8{C}_{\sigma }}{\left(1-{Y}^{2}\right)}^{2}-\displaystyle \frac{b}{12{C}_{\sigma }{C}_{\omega }}{\left(1-{Y}^{2}\right)}^{3}\\ +\,\displaystyle \frac{c}{16{m}^{2}{C}_{\sigma }{C}_{\omega }^{2}}{\left(1-{Y}^{2}\right)}^{4}+\displaystyle \frac{1}{2}{m}_{\omega }^{2}{\omega }_{0}^{2}{Y}^{2}\\ +\,\displaystyle \frac{1}{2}{m}_{\rho }^{2}\left[1-{\eta }_{1}(1-{Y}^{2})({C}_{\rho }/{C}_{\omega })+3{\eta }_{2}{C}_{\rho }{\omega }_{0}^{2}\right]{\left({\rho }_{3}^{0}\right)}^{2}.\end{array}\end{eqnarray} \tag{ 29 }$$

Similarly, for the pressure, the negative of the thermodynamic potential can be written as

$$\begin{eqnarray}&&P={P}_{n}+{P}_{p}-{\epsilon }_{\mathrm{meson}}-\displaystyle \frac{1}{2}{B}^{2}\equiv {P}_{0}-\displaystyle \frac{1}{2}{B}^{2}.\end{eqnarray} \tag{ 30 }$$

The contribution of the neutrons to the pressure, using Equation (9) and integrating over the momentum, is given by

$$\begin{eqnarray}&&{P}_{n}=\displaystyle \frac{1}{24{\pi }^{2}}\left[{k}_{F}^{n}{\mu }_{n}^{\star }\left(2{k}_{F}^{n}-3{m}^{\star 2}\right)+3{m}^{\star 4}\mathrm{log}\displaystyle \frac{{k}_{F}^{n}+{\mu }_{n}^{\star }}{{m}^{\star }}\right].\end{eqnarray} \tag{ 31 }$$

The pressure due to the protons, on the other hand, is given by

$$\begin{eqnarray}&&{P}_{p}={P}_{\mathrm{field}}+{P}_{\mathrm{med}}.\end{eqnarray} \tag{ 32 }$$

The magnetic field contribution to the pressure P_field from the Dirac sea is given by, from Equation (18),

$$\begin{eqnarray}&&\begin{array}{l}{P}_{{\rm{field}}}=\displaystyle \frac{| {eB}{| }^{2}}{{\left(2\pi \right)}^{2}}\\ \,\displaystyle \left[{\zeta }^{{\rm{{\prime} }}}(-1,x)-\displaystyle \frac{1}{2}({x}_{p}^{2}-{x}_{p}){\rm{log}}({x}_{p})+\displaystyle \frac{{x}_{p}^{2}}{4}+\displaystyle \frac{1}{12}{\rm{log}}{x}_{p}\right].\end{array}\,\end{eqnarray} \tag{ 33 }$$

The medium contribution to the pressure from the protons, on the other hand, is given by

$$\begin{eqnarray}&&\begin{array}{l}{P}_{{\rm{m}}{\rm{e}}{\rm{d}}}=\displaystyle \sum _{n=0}^{{n}_{max}}\displaystyle \frac{{\alpha }_{n}| eB| }{4{\pi }^{2}}\\ \,\displaystyle \displaystyle \left[{\mu }_{p}^{\star }{k}_{F,n}^{p}-\,\left({m}^{\star 2}+2n| eB| \right){\rm{l}}{\rm{o}}{\rm{g}}\,\left(\displaystyle \frac{{k}_{F,n}^{p}+{\mu }_{p}^{\star }}{{m}^{\star 2}+2n| eB| }\right)\right].\end{array}\end{eqnarray} \tag{ 34 }$$

In order to account for NS matter, one needs to incorporate the charge neutrality and beta equilibrium conditions as well. The charge neutrality conditions are as follows:

$$\begin{eqnarray}&&\displaystyle \sum _{b}{Q}_{b}\,{\rho }_{b}+\displaystyle \sum _{l}{Q}_{l}\,{\rho }_{l}=0,\end{eqnarray} \tag{ 35 }$$

where the index b is summed over nucleons (n, p), while the index l denotes sum over all leptonic states ($e,\mu$). Q_b and Q_l are the electrical charges of baryons and leptons, respectively. ρ_b(n, p) and ρ_l(e, μ) are the total baryon and lepton densities, respectively. Thus, at a given baryon number density ρ_b =(ρ_n + ρ_p), the charge neutrality condition is given by ${\rho }_{p}={\rho }_{e}+{\rho }_{\mu }$. The beta equilibrium condition leads to the chemical potentials of protons, neutrons, electrons, and muons given as μ_n = μ_p + μ_e and μ_e = μ_μ while the number densities of neutrons and protons are given in Equation (25). The lepton number density ρ_l (electron and muons) is given by

$$\begin{eqnarray}&&{\rho }_{l}=\displaystyle \sum _{n}^{{n}_{\max }}{\alpha }_{n}\displaystyle \frac{| {eB}| }{2{\pi }^{2}}\sqrt{{\mu }_{E}^{2}-2n| {eB}| },\end{eqnarray} \tag{ 36 }$$

where we have ignored the mass of electrons, and in the sum above, the maximum number of Landau levels is given as ${n}_{max}={\rm{Int}}[{\mu }_{E}^{2}/(2| {eB}| )]$.

Let us note that the EOS that we have derived given in Equation (30) corresponds to the thermodynamic pressure, i.e., the negative of the thermodynamic potential. However, in the presence of a magnetic field, the hydrodynamic pressure can be highly anisotropic when there is significant magnetization (Chatterjee et al. 2011; Ferrer et al. 2010; Huang et al. 2010; Canuto & Ventura 1977) of the matter. The pressure in the direction of the field P_∥ is the thermodynamic pressure as given in Equation (30). On the other hand, the pressure P_⊥ in the transverse direction of the magnetic field is given by, with P₀ as defined in Equation (30),

$$\begin{eqnarray}&&{P}_{\perp }={P}_{0}-{ \mathcal M }B+\displaystyle \frac{1}{2}{B}^{2},\end{eqnarray} \tag{ 37 }$$

where ${ \mathcal M }=-\partial {\rm{\Omega }}/\partial B$ is the magnetization of the system. Using Equation (33) and Equation (34), the total magnetization can be written as ${ \mathcal M }={{ \mathcal M }}_{\mathrm{med}}+{{ \mathcal M }}_{\mathrm{field}}$ where the magnetization of the medium is given as

$$\begin{eqnarray}&&\begin{array}{l}{{ \mathcal M }}_{{\rm{med}}}=\displaystyle \frac{{\rm{\partial }}{P}_{{\rm{med}}}}{{\rm{\partial }}B}=\displaystyle \sum _{n=0}^{{n}_{max}}\displaystyle \frac{{\alpha }_{n}| e| }{4{\pi }^{2}}\\ \,\displaystyle \displaystyle \left[{\mu }_{p}^{\star }{k}_{F,n}^{p}-\left({m}^{\star 2}+4n| {eB}| \right){\rm{log}}\left(\displaystyle \frac{{k}_{F,n}^{p}+{\mu }_{p}^{\star }}{{m}^{\star 2}+2n| {eB}| }\right)\right].\end{array}\end{eqnarray} \tag{ 38 }$$

On the other hand, the magnetization of the Dirac sea is given as

$$\begin{eqnarray}&&{{ \mathcal M }}_{\mathrm{field}}=\displaystyle \frac{\partial {P}_{\mathrm{field}}}{\partial B}=\displaystyle \frac{{e}^{2}B}{{\pi }^{2}}\,\left[\displaystyle \frac{1}{12}\mathrm{log}{x}_{p}-\displaystyle \frac{1}{24}+\displaystyle \frac{{x}^{3}}{2}{I}_{2}\right],\end{eqnarray} \tag{ 39 }$$

where we have used the expression for $\zeta ^{\prime} (-1,x)$ given in Equation (15) and defined the quantity I₂ as the integral (Chatterjee et al. 2011)

$$\begin{eqnarray}&&{I}_{2}=(2\pi ){\int }_{0}^{\infty }\displaystyle \frac{2{\tan }^{-1}y+y\mathrm{log}(1+{y}^{2})}{(\exp (2\pi {xy})-1)(1-\exp (-2\pi {xy}))}{ydy}.\end{eqnarray} \tag{ 40 }$$

To obtain the gross structural properties of the NS with the EOS as calculated above, we introduce a chemical-potential-dependent magnetic field (Mao et al. 2003; Rabhi et al. 2009; Dexheimer et al. 2012, 2014):

$$\begin{eqnarray}&&B(\mu )={B}_{\mathrm{surf}}+{B}_{0}\left[1-{e}^{\left(\tfrac{-b{\left(\mu -938\right)}^{a}}{938}\right)}\right],\end{eqnarray} \tag{ 41 }$$

where the parameters a and b are chosen to be 2.56 and 4.2426 × 10⁻⁵ MeV^−1.56, respectively, the chemical potential μ is given in MeV, B_surf = 10¹⁵ G is the magnetic field on the surface, and B₀ is the maximum magnetic field at the core. The field increases somewhat mildly with density from its core to the surface. We have taken here a chemical-potential (μ)-dependent magnetic field rather than the density-dependent magnetic field considered earlier, e.g., in Bandyopadhyay et al. (1997, 1998) and Huang et al. (2010). The density-dependent magnetic field can lead to an unphysical discontinuity in the case of a phase transition, e.g., to quark matter as the density changes discontinuously across the phase transition. Although we do not consider a phase transition here, we have taken a μ-dependent magnetic field, which does not give any discontinuity either in density or in the magnetic field.

Although the pressure in the presence of a magnetic field becomes anisotropic, for the range of magnetic fields we consider here, the deformation due to the anisotropic pressure is small and we can assume spherical symmetry for the structure of the NS. We shall give a post facto justification for the same. The equations for the structure of a relativistic spherical and static star composed of a perfect fluid were derived from Einstein's equations by Tolman, Oppenheimer and Volkoff, known as the TOV equations, which are (Weinberg 1972)

$$\begin{eqnarray}&&\displaystyle \frac{{dP}}{{dr}}=-\displaystyle \frac{G}{r}\displaystyle \frac{\left[\varepsilon +P\right]\left[M+4\pi {r}^{3}P\right]}{(r-2{GM})},\end{eqnarray} \tag{ 42 }$$

$$\begin{eqnarray}&&\displaystyle \frac{{dM}}{{dr}}=4\pi {r}^{2}\varepsilon ,\end{eqnarray} \tag{ 43 }$$

with G as the gravitational constant and M(r) as the enclosed gravitational mass. For the specified EOS, these equations can be integrated from the origin as an initial value problem for a given choice of central energy densities (ε_c). The value of r ( = R) where the pressure vanishes defines the surface of the star. We solve these equations to study the structural properties of a static NS using the EOS as derived and given in Equation (26) and Equation (30) for the magnetized, charge-neutral dense nuclear matter.

While simultaneous measurements of the mass and radius of NSs have the potential to constrain the EOS, such measurements are also plagued with uncertainties and model dependence on the radiation mechanisms at the NS surface as well as interstellar absorption. On the other hand, the observation of the inspiralling binary NSs with a gravitational wave detection GW 170817 could provide significant information about the structure of an NS. The tidal distortion of the NSs in a binary system links the EOS to the gravitational wave emission during the inspiral (Malik et al. 2018; Gomes et al. 2019). In the following, we shall estimate this parameter for the EOS for the magnetized nuclear matter.

The tidal deformity parameter λ relates the induced quadrupole moment Q_ij of an NS due to the strong tidal gravitational field ${{ \mathcal E }}_{{ij}}$ of the companion star. This quadrupole deformation in leading order in perturbation is given as (Hinderer 2008)

$$\begin{eqnarray}&&{Q}_{{ij}}=-\lambda {{ \mathcal E }}_{{ij}}.\end{eqnarray} \tag{ 44 }$$

The parameter λ is related to l = 2 and the tidal Love number as ${k}_{2}=\tfrac{3}{2}\lambda {R}^{-5}$, with R being the radius of the NS. One can estimate k₂ perturbatively by estimating the deformation h_αβ of the metric from the spherical metric. We consider here the leading-order static perturbation and axisymmetric perturbation. The deformation of the metric in the Regge–Wheeler gauge can be written as (Hinderer 2008)

$$\begin{eqnarray}&&{h}_{\alpha \beta }={\rm{diag}}\left[\left.-{e}^{\left.{\rm{\Phi }}(r\right)}{H}_{0},{e}^{\left.{\rm{\Lambda }}(r\right)}{H}_{2},{r}^{2}K(r),{r}^{2}{\sin }^{2}\theta K(r\right)\right]{Y}_{20}(\theta ,\phi ),\end{eqnarray} \tag{ 45 }$$

where H₀(r), H₂(r), and K(r) are the perturbed metric functions. It turns out that H₂(r) = −H₀(r) ≡ H(r) using Einstein's equation $\delta {G}_{\alpha }^{\beta }=8\pi \delta {T}_{\alpha }^{\beta }$, while $K^{\prime} (r)=-H^{\prime} (r)-H(r){\rm{\Phi }}^{\prime} (r)$. The logarithmic derivative of the deformation function H(r), i.e., $y(r)=r\tfrac{{H}_{0}^{{\prime} }(r)}{{H}_{0}(r)},$ satisfies the first-order equation (Damour & Nagar 2009)

$$\begin{eqnarray}&&r\displaystyle \frac{{dy}(r)}{{dr}}+y{\left(r\right)}^{2}+y(r)F(r)+{r}^{2}Q(r)=0,\end{eqnarray} \tag{ 46 }$$

with

$$\begin{eqnarray}&&F(r)=\displaystyle \frac{r-4\pi {r}^{3}\left(\epsilon (r)-p(r)\right)}{r-2M(r)},\end{eqnarray} \tag{ 47 }$$

$$\begin{eqnarray}&&\begin{array}{l}Q(r)=\displaystyle \frac{4\pi r\,\left(5\epsilon (r)+9p(r)+\tfrac{\epsilon (r)+p(r)}{\partial p(r)/\partial \epsilon (r)}-\tfrac{6}{4\pi {r}^{2}}\right)}{r-2M(r)}\\ -\,4\,{\left[\displaystyle \frac{M(r)+4\pi {r}^{3}p(r)}{{r}^{2}\left(1-2M(r)/r\right)}\right]}^{2}.\end{array}\end{eqnarray} \tag{ 48 }$$

To calculate the tidal deformation, the equation for the metric perturbation, Equation (46), can be integrated together with the TOV Equations (42) and (43) for a given EOS radially outwards, with the boundary conditions y(0) = 2, $p(0)={p}_{c},$ and $M(0)=0$, where y(0), p_c, and M(0) are the pressure and mass density at the center of the NS, respectively.

The tidal Love number k₂ is related to y_R = y(R) through

$$\begin{eqnarray}&&\begin{array}{l}{k}_{2}=\displaystyle \frac{8{C}^{5}}{5}{\left(1-2C\right)}^{2}\left[2+2C\left({y}_{R}-1\right)-{y}_{R}\right]\\ \ \times \ \left\{2C\left(6-3{y}_{R}+3C(5{y}_{R}-8)\right)\right.\\ +\,4{C}^{3}\left[13-11{y}_{R}+C(3{y}_{R}-2)+2{C}^{2}(1+{y}_{R})\right]\\ {\left.+3{\left(1-2C\right)}^{2}\left[2-{y}_{R}+2C({y}_{R}-1)\right]\mathrm{log}\left(1-2C\right)\right\}}^{-1},\end{array}\end{eqnarray} \tag{ 49 }$$

where C $(\equiv M/R)$ is the compactness parameter of the star of mass M. The dimensionless tidal deformability Λ is defined as (Flanagan & Hinderer 2008; Hinderer 2008; Hinderer et al. 2010; Damour et al. 2012),

$$\begin{eqnarray}&&{\rm{\Lambda }}=\displaystyle \frac{2}{3}{k}_{2}{\left(R/M\right)}^{5}.\end{eqnarray} \tag{ 50 }$$

The observable signature of relativistic tidal deformation will have an effect on the phase evolution of the gravitational wave spectrum from the inspiralling binary NS system. This signal will have cumulative effects of the tidal deformation arising from both stars. Therefore, one can combine the tidal deformabilities and define a dimensionless tidal deformability Λ by taking a weighted average as (Favata 2014)

$$\begin{eqnarray}&&\tilde{{\rm{\Lambda }}}=\displaystyle \frac{16}{13}\left[\displaystyle \frac{({M}_{1}+12{M}_{2}){M}_{1}^{4}{{\rm{\Lambda }}}_{1}+({M}_{2}+12{M}_{1}){M}_{1}^{4}{{\rm{\Lambda }}}_{2}}{{\left({M}_{1}+{M}_{2}\right)}^{5}}\right].\end{eqnarray} \tag{ 51 }$$

In the above, Λ₁ and Λ₂ are the individual tidal deformabilities corresponding to the two components of the NS binary with masses M₁ and M₂, respectively.

As mentioned earlier, we only incorporate the ρ−σ cross-coupling parameterization of the energy density functional (EDF) as given in Equations (8) and (9) in the present effective chiral model given in Equation (1). With this cross-coupling term, the EDF is able to satisfy the nuclear matter properties, specifically the density dependence of the symmetry energy with the available empirical estimates as given in Table 1 as well as also satisfying the present constraints on NS properties such as the NS maximum mass and radius (Malik et al. 2017). These calculations were performed without the effect of a magnetic field. In the following, we show the effects of a strong magnetic field on this EDF. In particular, we investigate the effective mass of nucleons and the relative particle population of charged and uncharged particles in β-equilibrated nuclear matter.

In Figure 1, we show the dependence of the magnetic field on the baryon chemical potential corresponding to a magnetic field at the NS core of ${B}_{0}={10}^{18}\,\mathrm{and}\,{10}^{19}$ G. The formalism adopted here (see Equation (41)) is consistent with the Maxwell equations (Chatterjee et al. 2019; Peres Menezes & Laércio Lopes 2016; Alloy & Menezes 2017). At very high baryon chemical potentials, the value of the magnetic field becomes extreme, comparable to the hadronic scale. Post facto, we have observed that for the central value of a magnetic field in the range ${10}^{15}\mbox{--}{10}^{18}$ G the deformation of the star from spherical symmetry is less than 1%.

Zoom In Zoom Out Reset image size

Let us first discuss the vacuum nucleon mass, i.e., the nucleon mass at zero baryon density as a function of magnetic field. For μ = 0, the Equations (3)–(5) are solved trivially with ${\rho }_{0}^{3}=0$ and ω₀ = 0. The remaining Equation (4) is to be solved for the effective mass ($Y={m}^{\star }/m$). For zero magnetic field, there is a unique solution Y = 1, i.e., ${m}^{\star }=m$. However, because of the Dirac sea response to the magnetic field, the scalar density for protons has a nonvanishing contribution ${\rho }_{s}^{{field}}$ as shown in Equations (19) and (24). The numerical solutions obtained by solving Equation (4) at a nonzero magnetic field for the effective nucleon mass is shown in Figure 2. The present model, being an effective model for nucleon matter, is not expected to be a profound result for a vacuum of strong interactions. Nonetheless, we see a magnetic catalysis for nucleon mass similar to the magnetic catalysis of chiral symmetry breaking for quarks (Chatterjee et al. 2011). For a very large magnetic field (larger than ${eB}\sim 10{m}_{\pi }^{2}$), we do not get any solutions for Equation (4) for the effective nucleon mass. This can be an artifact of the mean field approximations that we use here. Going beyond the mean field with magnetic-field-dependent meson masses can possibly cure this (Haber et al. 2014). In what follows, however, we will ignore such an effect on the meson masses and continue with the mean-field approximation. We note, however, that there is no such limitations at finite baryon density.

Zoom In Zoom Out Reset image size

In Figure 3, in the upper panel, we show the variation of the nucleon mass for constant magnetic field. The effect of magnetic field become significant only when the field strength exceeds about 10¹⁷ G. Due to magnetic catalysis the effective mass in presence of magnetic field is larger than the mass without the magnetic field. On the lower panel, we display the relative populations of the charge particles. The populations of the charged particles are influenced both by magnetic field and the charge neutrality condition. The proton fraction remains sufficiently small until 2.5ρ₀. As the magnetic field is increased protons contribute at still higher density as their masses become heavier with magnetic field.

Zoom In Zoom Out Reset image size

Next, in Figure 4, we show the same variables as in Figure 3 but with a chemical-potential-dependent magnetic field as given in Equation (41). As can be seen in the figure, the proton fraction remains sufficiently small, up to three to four times the nuclear matter density, and the contribution of the neutrons to the pressure remains dominant. As compared to Figure 3, the effect of the magnetic field in this case is milder because the field decreases as the density decreases. The strength of the magnetic field to induce significant changes can be estimated in a straightforward manner. The contributions from the protons become significant when the lowest Landau level (n = 0) is occupied. One can estimate this to happen when $2{eB}\,\gt \sqrt{({\mu }_{p}^{\star 2}-{m}^{\star 2})}$. Equivalently, this corresponds to a magnetic field ${eB}\gt 3.2\times {10}^{19}{({\rho }_{p}/{\rho }_{0})}^{2/3}$ G, where ρ₀ ∼ 0.16 fm⁻³. This is the reason why the effect of a magnetic field is not seen for magnetic fields up to B₀ ∼ 10¹⁸ G. Moreover, with the density dependence as in Equation (41), the effect is seen only at high density. It can be noted that both effective mass and hence the population of different particles do not show any changes for the magnetic field intensity at the NS core within B₀ = 0–10¹⁸ G, whereas there is a noticeable increase in both nucleon effective mass and proton population when B₀ is more than 10¹⁸ G (Broderick et al. 2000). A similar effect can be seen in the relative particle population when ${B}_{0}\gt {10}^{18}$ G. There is an incremental change in the charge particle concentration for the charged species, such as e⁻, p⁺, and μ⁻ at ≈ρ = 3ρ₀, as a result of which the neutron concentration drops at higher densities, in comparison to the no-field case or the NS core field intensity, until ${B}_{0}\leqslant {10}^{18}$ G.

Zoom In Zoom Out Reset image size

The net magnetization in the matter, given in Equations (38) and (39), is shown in Figure 5. It shows that the net magnetization is almost negligible until B = 10¹⁶ G. From B ≥ 10¹⁷ G onwards, it starts increasing and becomes highly oscillatory. The oscillation of the magnetization happens as an outcome of the well-known de Haas–van Alphen effect (Ebert et al. 1994; Shoenberg 1984; Chatterjee et al. 2011) in which the charged particles, due to Landau quantization, exist only in orbitally quantized states in a magnetic field, and as the number of occupied Landau levels changes with the magnetic field, the magnetization becomes oscillatory. The oscillatory behavior is more pronounced with the increasing strength of the magnetic field. The irregularity in the oscillation is due to the medium dependence of the nucleon mass, which itself depends on the magnetic field. This behavior is also consistent with the findings of Broderick et al. (2000) and Chatterjee et al. (2011). We also show the variation of pressure (scaled with the value P₀, the pressure in the absence of a magnetic field) in directions parallel and perpendicular to magnetic field with respect to magnetic field intensity at a baryon density of 4.0ρ₀ in Figure 6. The parallel component of the pressure as given in Equation (30) decreases with magnetic field and even becomes negative as the magnetic field is increased. On the other hand, the perpendicular component as given in Equation (37) monotonically increases with magnetic field. It turns out that at this density (ρ = 4.0ρ₀), the magnetization contribution ${ \mathcal M }B$ is two orders of magnitude smaller compared to the matter contribution. Therefore, the oscillatory behavior of magnetization as seen in Figure 5 is not reflected in the transverse pressure in Figure 6. For ρ = 4.0ρ₀, the difference between P∥ and P_⊥ becomes significant for magnetic field strength beyond B = 10^18.4 G.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Next, in Figure 7, we show the EOS of the magnetized charged neutral matter for different strengths of the core magnetic field B₀. Here we have taken the EOS in Equation (41). Similar to Figure 3 for the effective masses, the pressure does not show any change as the field strength is increased up to B₀ = 10¹⁷ G, the magnetic field at the core of the NS. The effects of the magnetic field become noticeable for B₀ beyond 10¹⁸ G. In the figure, the solid and dashed lines correspond to the pressure in the directions perpendicular and parallel to the magnetic field. Beyond 10¹⁸ G, the difference between these two pressures increases rapidly. As is obvious from the figure, the perpendicular component becomes stiffer while the parallel component becomes softer. For B₀ beyond 10¹⁹ G, the parallel component becomes negative. Later, we shall use TOV Equations (42) and (43) to solve for the mass and radius of a magnetized NS; we keep the magnetic field strength B₀ up to 10¹⁸ G, so that the anisotropy of pressure is not too large and the spherically symmetric TOV equations can be applicable.

Zoom In Zoom Out Reset image size

We then proceed to calculate the mass–radius relationship with the magnetized charge-neutral NS matter, which is shown in Figure 8. The TOV equations are solved by taking the EOS with both parallel and perpendicular pressures. The solid line in the figure corresponds to taking the P_⊥ for pressure in TOV equations. The dashed line corresponds to taking the P_∥ as the pressure. We find that there is no appreciable change in the gravitational mass and radius of the NS when magnetic effects are incorporated for fields straight up to B₀ = 10¹⁷ G compared to the case where the magnetic field is absent. Even with the increase in the central magnetic field up to 10¹⁸ G, there is hardly any change in the gross structural properties of the NS although the particle concentration in the NS matter changes as seen in the lower panel of Figure 4. Increasing the value of B₀ further results in making the anisotropy between the two pressures larger, with even the parallel component of the pressure becoming negative, leading to mechanical instability. For the central value of the magnetic field B₀ up to 10¹⁸ G, the resulting maximum mass of the NS satisfies the maximum mass constraint (M = 2.01 ± 0.04M_⊙; Antoniadis et al. 2013) in both P_⊥ and P_∥ as shown in the figure. We have also shown the results for B₀ = 10¹⁸ G. The maximum mass corresponding to ${P}_{\perp }$ and P_∥ are 1.97 M_⊙ and 1.96 M_⊙ with the asymmetry δ in the masses ($\delta =({M}_{\perp }-{M}_{\parallel })/({M}_{\perp }+{M}_{\parallel })$) to be about 0.5%. The corresponding radii are 11.43 km and 11.47 km, respectively, leading to a deformation in radius less than 1%. This only means that the pressure anisotropy is rather small, leading to a tiny mass asymmetry and hence gives a post facto justification of using the isotropic TOV equations to calculate approximately the gross structural properties of NS. The corresponding radii (R_1.4) of the canonical mass NS (M = 1.4M_⊙) that we obtain in the present model also agrees well with the empirical estimates given by Lattimer & Lim (2013). The overall results obtained from the solutions of TOV equations are tabulated in Table 2. For the central magnetic field of 10¹⁸ G, by taking P_∥ in the TOV equations, we get a lower maximum mass, which is probably expected as the corresponding EOS becomes softer for P_∥.

Zoom In Zoom Out Reset image size

As discussed earlier, the data from the gravitational wave detection from the observation of the inspiralling binary NSs GW 170817 (Abbott et al. 2017) could possibly constrain not only the properties of NSs but also put constraints on the EOS. In Figure 9, we plot the dimensionless tidal deformability Λ (left panel) as given in Equation (50) and the Love number k₂ (right panel) as given in Equation (49) as a function of NS mass for the EOS with different values of the central magnetic field (B₀). Let us note that the Love number k₂ not only depends upon the compactness parameter $C\equiv M/R$ but also on y(R), the value of the logarithmic derivative of the deformation function that depends upon the internal structure of the NS as in Equation (46). The value of k₂ has a peak around 1.0 M_⊙ while it is rather low at higher and lower masses as seen in the right panel of the figure, indicating that the quadrupole deformation is maximum for intermediate-mass ranges for a given EOS. The obtained value of Λ for a 1.4 M_⊙ NS is 520 in the present model without including the magnetic field. Λ and k₂ are also calculated for both the EOSs with pressure perpendicular and parallel to the magnetic field. There is almost no change in their values for ${B}_{0}=0\mbox{--}{10}^{18}$ G with the EOS for both the cases. However, for ${B}_{0}\gt {10}^{18}$, the Λ and k₂ values both increase in the perpendicular direction, and the effect is the opposite for the parallel case. However, the change is very small and for k₂ the effect is only seen for larger NS masses. Thus, the magnetized NS still satisfies the constraint on the tidal deformability parameter Λ bound from the gravitational wave data of GW 170817. The corresponding Love number (k₂) also lies within the acceptable range (Hinderer 2008). The agreement of both these parameters validates the properties of the NS obtained within the model with the inclusion of a magnetic field.

Zoom In Zoom Out Reset image size

Finally, in Figure 10, we plot the tidal deformability parameters Λ₁ and ${{\rm{\Lambda }}}_{2}$, which are linked to the NS binary companion having a high mass M₁ and a low mass M₂ associated with the GW 170817 event. We have plotted them for the EOS for different values of magnetic field considered here. The curves are obtained by varying the high mass (M₁) independently in the range $1.365\lt M/{M}_{\odot }\lt 1.60$ obtained for GW 170817 whereas the low mass (M₂) is determined by keeping the chirp mass (${M}_{\mathrm{chirp}}={\left({M}_{1}{M}_{2}\right)}^{(3/5)}/{\left({M}_{1}+{M}_{2}\right)}^{5}$) fixed at the observed value of 1.188 M_⊙. The long-dashed dark yellow line signifies the 90% probability contour found from this event. The black dotted line, on the other hand, signifies the 50% probability contour. As may be observed from the figure, the EOS obtained from the present chiral model for nucleon matter lies well within the two limits with or without the magnetic field.

Zoom In Zoom Out Reset image size

Let us summarize the salient feature of the present investigation. We have looked into the different effects of a magnetic field on a nuclear matter EOS within the ambit of a chiral model which is a generalized sigma model couple to nucleons. Apart from sigma and pions, the model incorporates ρ and ω mesons with higher-order mesonic fields. The model with cross-coupling between the isovector and scalar field successfully described the symmetry energy parameters (Malik et al. 2017).

In the present calculation incorporating the effects of a magnetic field, we have included the effect of the field on the Dirac sea of protons. This affects the vacuum mass of the nucleon at zero baryon density and temperature. The mass of the nucleons is seen to be increasing with magnetic field, showing the magnetic catalysis effects. For a small field, this increase in mass is seen to be quadratically dependent on the magnetic field, which we take to be homogeneous.

Next, we calculated the effective mass of the nucleons in the medium with nonzero densities in the presence of a constant magnetic field. For a field straight up to 10¹⁷ G, the EOS does not change compared to the case of a vanishing magnetic field. Due to the directional dependence of the magnetic field, the pressure is no longer isotropic. The EOS for the pressure parallel to the magnetic field (P_∥) becomes softer, while that perpendicular to the magnetic field (P_⊥) becomes stiffer. For a small anisotropy in pressure, we have used the TOV equations to obtain the gross structural properties of the NS. Within this approximation, the correction due to the magnetic field remains small for the masses and radii of an NS. The masses of the NS and their radii appear to be within the corresponding acceptance limits.

We next estimated the tidal deformability, the Love number for an NS with and without a magnetic field. For the strength of the magnetic field considered here (B₀ = 10^18.5 G), the tidal deformability parameters do not differ much from the zero magnetic field cases and lie within the acceptable range for the Λ from GW 170817. We have also calculated the tidal deformabilities in the phase space of Λ₁ and Λ₂ associated with the two components of the binary NSs related to GW 170817. The present EOS with or without a magnetic field is consistent with the limits derived from this event.

We have confined our attention to the EOS based on the nucleonic degrees of freedom in this work. The contributions of the exotic degrees of freedom, such as hyperons and kaons, will soften the effect on the EOS, leading to a smaller maximum mass of the stars. The effective chiral model considered for this work is excellent in producing nuclear saturation properties. The nuclear symmetry energy, as well as its slope and curvature parameters at the saturation density, are in harmony with those deduced from a diverse set of experimental data, and the obtained maximum mass of the NS is 1.97 M_⊙ (Malik et al. 2017), which is the lower bound of the observed NS of mass 2 M_⊙. The inclusion of hyperons in this model will further reduce the maximum NS mass irrespective of the choice of the hyperon coupling constant.

N.K.P. would like to acknowledge the Department of Science and Technology (DST), India, for the support DST/INSPIRE Fellowship/2019/IF190058. T.M. would like to thank the kind hospitality provided by the Physical Research Laboratory, Ahmedabad, India, where part of this work was done.