Exploring the effects of Delta Baryons in magnetars

The AMM of a particle is a deviation from the magnetic moment of that particle, as predicted by the “classical” tree-level prediction. Historically, the term anomalous was used to describe the deviation from the Dirac equation prediction for a system of fermions under the influence of a magnetic field, the magnetic moment, and thus, refers to fundamental particles. Dipole moment, on the other hand, is used for composite particles, such as baryons, since their value depends on the configuration of quarks and gluons inside it, and thus, are not anomalous in the strict sense. As commonly used in the literature and for simplicity, in this work we use the term AMM in all cases.

The energy spectrum of baryons with an AMM can be empirically determined, but a theoretical derivation of their values from first principles is yet an unaccomplished task. The AMMs of nucleons are measured to a very high precision, with errors of the order of $10^{-9}$ (Greene et al., 1979; Schneider et al., 2017), but the same does not apply to heavier baryons. Measurements of the hyperon AMMs are precise to an order of $10^{-2}$ (Zyla et al., 2020), while $\Delta$s are experimentally determined only for the positively charged $\Delta^{++}$ and $\Delta^{+}$. For the $\Delta^{+}$, there is a single measurement of $\mu_{\Delta^{+}}$/$\mu_{N}=2.7^{+1.0}_{-1.3}\pm 1.5$ that comes from the $\gamma p\rightarrow p\pi\gamma^{\prime}$ reaction (Kotulla et al., 2002), while for the $\Delta^{++}$ there are several measurements coming from the $\pi^{+}p\rightarrow\pi^{+}p\gamma$ bremsstrahlung cross section, with values in the range $\mu_{\Delta^{++}}$/$\mu_{N}=3.7-7.5$  (Zyla et al., 2020). These measurements include systematic uncertainties, but additional theoretical uncertainties lead to errors $\sim\pm 3$. Complementary to experimental results, lattice quantum chromodynamics (LQCD) has been able to extract AMM values for $\Delta$ baryons. The values utilized in this paper are based on the predictions from LQCD provided in (Cloet et al., 2003) that lie within the experimental uncertainties of the experimentally measured AMMs.

Different properties of baryons considered in this study are shown in Tab. 1. The AMM strength coefficients $\kappa_{b}$ are related to the magnetic moments $\mu_{b}$ through the relation

which depends on the baryon charge $q_{b}$, the nuclear magneton $\mu_{N}=e/2M_{p}$, with $e$ being the electron charge, and the ratio of the proton mass $M_{p}$ to the baryon mass $M_{b}$. Although the expression (1) is derived for spin-1/2 fermions in the non-relativistic regime, it is still commonly employed to the description of the spin-3/2 particles (see López Castro and Mariano, 2002, 2001). This subject is controversial, as the Rarita-Schwinger equation with minimal coupling predicts a gyromagnetic ratio of 2/3, while low energy/optical theorems predict a value of 2. For a more in-depth discussion we refer to Deser et al. (2000), which studies a generic non-minimal electromagnetic coupling in the Rarita-Schwinger formalism. In this work, we also account for the leptons (electron and muon) AMMs, given by $\kappa_{e}/\mu_{B_{e}}=1.15965\times 10^{-3}$ and $\kappa_{\mu}/\mu_{B_{\mu}}=1.16592\times 10^{-3}$, respectively, with $\mu_{B_{l}}=e/2M_{l}$, for $l=\{e,\mu\}$.

To describe baryon-dense matter subject to a strong magnetic field, we must start from a Lagrangian density describing how the particles interact with each other and with the external electromagnetic field. The photons are simply described by the massless Proca Lagrangian density, followed by a term containing the electromagnetic interaction for charged baryons and leptons, and a term describing the AMMs of baryons $b$ and leptons $l$

where $F_{\mu\nu}=\partial_{\mu}A_{\nu}-\partial_{\nu}A_{\mu}$, $\sigma^{\mu\nu}=\frac{i}{2}\left[\gamma^{\mu},\gamma^{\nu}\right]$, and the vector potential $A_{\mu}=(0,0,Bx,0)$ is chosen such that the magnetic field is parallel to the $z$-axis. Leptons are described additionally by the (free with respect to the strong force) Dirac Lagrangian density

We make use of two different relativistic models to describe the still widely unknown strong interaction between baryons. The first model is a non-linear version of the Walecka model, where the baryon interactions are mediated by the $\sigma$, $\omega$, $\rho$ and $\phi$ mesons, within the mean field approximation. We choose the recently proposed L3$\omega\rho$ parametrization  (Lopes, 2022), which includes an additional $\omega\rho$ interaction that allows the correct prediction of slope of the symmetry energy, neutron-star radii and tidal deformabilities. The $\phi$ meson (with hidden strangeness) couples only to the hyperons, allowing a higher maximum mass to be reproduced for neutron stars, thus circumventing the well-known hyperon puzzle (Chamel et al., 2013), with an effect similar to the higher-order $\omega^{4}$ self interaction, also included. In addition to satisfying standard astrophysical constraints from LIGO/VIRGO and NICER, the model satisfies nuclear ground-state properties of finite nuclei and bulk properties of infinite nuclear matter. It is also consistent with the PREX-2 results for the symmetry energy of $L=106\pm 37$ (Reed et al., 2021), though at the lower end of the error band.

The model Lagrangian density is written as $\mathcal{L}=\mathcal{L}_{b}+\mathcal{L}_{\rm m}$, where the first term is the (interacting) Dirac Lagrangian density for nucleons, hyperons, and $\Delta$s, and the second term accounts for the self interaction among scalar and vector mesons¹¹1Spin-3/2 baryons are in fact described by the Rarita-Schwinger Lagrangian density $$\mathcal{L}_{RS}=-\sum_{b=\Delta}\frac{1}{2}\bar{\psi}_{\mu\,b}\left(\epsilon^{\mu\eta\lambda\nu}\gamma_{5}\gamma_{\eta}\partial_{\lambda}-im_{b}\sigma^{\mu\nu}\right)\psi_{\nu\,b}~,$$ (4) where $\epsilon^{\mu\eta\lambda\nu}$ is the Levi-Civita symbol, $\sigma^{\mu\nu}=\frac{i}{2}\left[\gamma^{\mu},\gamma^{\nu}\right]$, and $\psi_{\mu\,b}$ is a vector-valued spinor with additional components (when compared to the four component spinor in the Dirac equation). Nonetheless, its equation of motion can be written compactly as $\left(i\gamma^{\nu}\partial_{\nu}-m\right)\psi_{\mu}=0$, see de Paoli et al. (2013).

and

where $I_{3\,b}$ is the baryon isospin 3rd component, given in Table  1. The mass of the baryons is modified by the medium, giving rise to an effective mass $M_{b}^{*}=M_{b}-g_{\sigma b}\sigma$. The fittings of the self couplings $\lambda$ and $\kappa$, and the couplings between the mesons $i=\{\sigma,\omega,\rho,\phi\}$ and baryons $b$, defined in terms of the ratios $x_{ib}=g_{ib}/g_{iN}$, are discussed in detail in Lopes (2022).

Relevant for this work, the scalar meson couplings are fitted to reproduce the hyperon potential depth $U_{\Lambda}=-28$ MeV (for isospin-symmetric matter at saturation) and the remaining relative strength of the coupling constants are determined by SU(3) symmetry group arguments, as proposed by Lopes and Menezes (2014), determining the complete hyperon-meson coupling scheme from a single free parameter, $\alpha_{v}$. Despite the value of $\alpha_{v}$, hyperons are always present in the neutron-star matter and the sequence of hyperon thresholds are always the same, with an inversely proportional relationship between $\alpha_{v}$ and the stiffness of the EoS. In this work, we choose to use $\alpha_{v}=0.5$, which results in values for the additional potentials $U_{\Sigma}=+21.8$ MeV and $U_{\Xi}=+35.3$ MeV, and a stiffer EoS with respect to the values $\alpha_{v}=0.75$ or 1.0 that are considered in Lopes and Menezes (2014). Though the potential for the $\Xi^{-}$mesons is repulsive in the parametrization used, recent observational constraints predict it to be attractive (Fabbietti et al., 2021; Acharya et al., 2019; Friedman and Gal, 2021), but to reproduce such an attractive potential we would need an extra free parameter in the meson couplings Lopes (2022).

The $\Delta$ couplings are treated more freely, as their behavior is not well known. The scarce information present in the literature, such as transport models (Cozma and Tsang, 2021) and quasi-elastic scattering of electrons off nuclei (Bodek and Cai, 2020), allows us to infer that the nucleon-$\Delta$ potential is slightly more attractive than the nucleon-nucleon one, so that, $U_{N}-30{\rm~MeV}\lesssim U_{\Delta}\lesssim U_{N}$, which implies $x_{\sigma\Delta}$ is greater than $1$. Also, the vector coupling is constrained by LQCD results as respecting the relation $0\lesssim x_{\sigma\Delta}-x_{\omega\Delta}\lesssim 0.2$, and no constraint is put in the $x_{\rho\Delta}$ value (Wehrberger et al., 1989; Ribes et al., 2019; Raduta, 2021). Early investigations on the effect of these parameters were made in de Oliveira et al. (2000, 2007) and their role in the stellar particle composition and maximum-mass was studied considering $x_{\sigma\Delta}=1.0$ and $1.1$, within two classes of relativistic mean-field models in Dexheimer et al. (2021a). Following the previous study, we analyse the scenarios with $x_{\sigma\Delta}=x_{\omega\Delta}=1.0$ and $x_{\sigma\Delta}=x_{\omega\Delta}=1.2$, keeping $x_{\rho\Delta}=1.0$, that generates, respectively, potentials $U_{\Delta}=-66.25$ MeV (equal to the nucleon potential) and $-79.50$ MeV.

The second model we use in this work is the chiral mean-field (CMF) model, which is based on a nonlinear realization of the chiral sigma model. As in all chiral models, the masses of the baryons are generated (and not just modified) by the medium. In this way, at large temperatures and/or densities they decrease allowing chiral symmetry to be restored, in agreement with LQCD results (Aarts et al., 2017). In this work, we restrict ourselves to the hadronic version of the model (with leptons) developed by Dexheimer and Schramm (2008), and disregard the possibility of phase transitions to quark matter. We add an additional $\omega\rho$ interaction to be in better agreement with data for the slope of the symmetry energy, neutron-star radii, and tidal deformabilities (Dexheimer et al., 2015, 2019). We also add a higher-order $\omega^{4}$ interaction to reproduce more massive neutron stars (Dexheimer et al., 2021b).

The mean-field model Lagrangian density has the terms

and

where the effective mass of baryons is $M_{b}^{*}=g_{\sigma b}\sigma+g_{\delta b}I_{3\,b}\delta+g_{\zeta b}\zeta+M_{0_{b}}$, including additional corrections given by the scalar-isovector $\delta$ and scalar-isoscalar (with hidden strangeness) $\zeta$ mesons, and a small bare mass correction $M_{0}$. The couplings were fitted to reproduce hadronic vacuum masses, decay constants, nuclear saturation properties, and to reach more than $2.1$ M_⊙ stars. See Roark and Dexheimer (2018) for a complete list of coupling constants. We follow the SU(3) and SU(6) coupling schemes for the scalar and vector couplings of mesons and baryons. In this way, there are only two free parameters left: one fitted to reproduce for symmetric matter at saturation the potential $U_{\Lambda}=-27$ MeV and another one (${x_{\omega\Delta}}=g_{\omega\Delta}/g_{\omega N}=1.25$) fitted to reproduce under the same conditions $U_{\Delta}=-64$ MeV $\sim U_{N}$. They result additionally in $U_{\Sigma}=6$ MeV and $U_{\Xi}=-17$ MeV. A much larger $x_{\omega\Delta}$ would suppress all $\Delta$s, while a much lower value would suppress all hyperons.

For both models, the equations of motion for the mesonic fields are obtained from the Lagrangian densities via the Euler-Lagrange equations. Under charge neutrality and $\beta$-equilibrium conditions, we can write the chemical potential of a baryon as a relation between the chemical potential of the neutron and the electron, $\mu_{n}$ and $\mu_{e}$, respectively, and its electric charge, i.e.,

At low (effectively zero) temperature, the Fermi energy spectra of baryons is

while for leptons it is simply $E_{F\,l}^{\ast}=\mu_{e}$.

In the presence of a magnetic field, the Fermi momentum (squared) can be calculated from the difference between the Fermi energy (squared) and

1. 1.

   the square of the effective mass modified by the AMM for particles that are not electrically charged $(q_{b}=0)$,

   $$k_{F,b}^{2}(s)={E_{F\,b}^{\ast}}^{2}-\left(M_{b}^{\ast}-s\kappa_{b}B\right)^{2};$$

   (11)

2. 2.

   the square of the effective mass modified by Landau quantization and AMM for particles that are electrically charged $(q_{b}\neq 0)$,

   $$k_{F,b}^{2}(\nu,s)={E_{F\,b}^{\ast}}^{2}-\left(\sqrt{M_{b}^{\ast 2}+2\nu|q_{b}|B}-s\kappa_{b}B\right)^{2}.$$

   (12)

For the momentum of leptons, $M^{*}$ is simply $M$. In the latter case, the Fermi momentum refers to the local direction of the magnetic field, hereafter defined as the $z$-axis. In the transverse direction to the local magnetic field, the Fermi momentum is restricted to discrete values $2\nu|q_{b}|B$, where the Landau levels $\nu$ relate to the orbital angular momentum $n$ via the relation

where $n=0,1,2...$ . For particles with spin $1/2$, the first Landau level ($\nu=0$) is occupied by a single spin projection: $s=+1$ for $q_{b}>0$ and $s=-1$ for $q_{b}<0$. The second level ($\nu=1$) and above are occupied by both spin projections $s=\{\pm 1\}$. For the spin-3/2 positively charged $\Delta$s, the first level ($\nu=0$) is occupied by the spin projections s$=\{+3,+1\}$, the second level ($\nu=1$) by $s=\{+3,\pm 1\}$, and hereafter all spin states are occupied. For the negatively charged $\Delta^{-}$ spin projection, signs are reversed for the lowest levels.

At zero temperature, there is a maximum Landau level allowed, beyond which the momentum in Eq. (12) becomes imaginary given by

The number density for each baryon is also defined separately for electrically charged and uncharged particles, respectively,

as well as the scalar densities,

The expressions for energy density and pressure are different for each model and can be obtained from the energy-momentum tensor for matter (discussed in the following).

For spherically symmetric neutron stars, given an EoS $P(\varepsilon)$, the global structure can be obtained by solving the Tolman-Oppenheimer-Volkoff (TOV) equations of hydrostatic equilibrium

where $M(r)$ is the stellar mass contained within the radius $r$ and $\nu(r)$ is a gravitational potential for the line element in spherical coordinates

The TOV equations cannot be applied to describe the structure of the magnetars we study in this work because the spherical symmetries assumed in Eq. (22) will not hold. This is due to the strong magnetic fields we infer for such objects, which produce highly deformed stellar shapes. Instead, the stellar structure must be determined by solving equations in General Relativity describing the stationary configuration for the fluid, coupled with Einstein field equations. The energy-momentum tensor, which contains the information on the matter properties of stars, enters the stellar structure equations as the source of the Einstein equations. Neglecting the coupling to the electromagnetic field, one generally assumes a perfect fluid and the energy-momentum tensor takes the form

where $\varepsilon$ denotes the (matter) energy density, $P$ the pressure, and $u^{\mu}$ the fluid four-velocity.

The EoS then relates pressure and energy density to the relevant thermodynamic quantities. In Chatterjee et al. (2015), the general expression for the energy-momentum tensor in the presence of an electromagnetic field was derived, starting from a microscopic Lagrangian including interactions between matter and the electromagnetic field

where $\mu_{0}$ is the vacuum permeability, $g^{\mu\nu}$ the metric tensor, and $x$ is the magnetisation. The electromagnetic field tensor has been expressed as $F_{\mu\nu}=\epsilon_{\alpha\beta\mu\nu}u^{\beta}B^{\alpha}~,$ with $\epsilon_{\alpha\beta\mu\nu}$ being the four-dimensional Levi-Civita symbol (Gourgoulhon, 2012). Assuming an isotropic medium and a magnetisation parallel to the magnetic field, the magnetisation tensor $M_{\mu\nu}$ can be written as

with the magnetization four-vector defined as $a_{\mu}=\frac{x}{\mu_{0}}B_{\mu}$. In the absence of magnetisation, i.e. for $x=0$, this expression reduces to the standard magnetohydrodynamics form for the energy-momentum tensor (c.f. Gourgoulhon, 2012).

Strong magnetic fields result in an anisotropy of the energy momentum tensor and break spherical symmetry, such that with increasing strength of the magnetic field, the shape of a magnetar departs more and more from a spherical shape. Interpreting the spatial elements of the fluid rest frame energy-momentum tensor as pressures, then there is a difference induced by the orientation of the magnetic field, commonly referred to as “parallel” and “perpendicular” pressures. Several earlier works tried to compute the mass-radius relations of strongly magnetised neutron stars through a first approach using isotropic TOV equations (Rabhi et al., 2009; Ferrer et al., 2010; Strickland et al., 2012; Lopes and Menezes, 2012; Dexheimer et al., 2013; Casali et al., 2014). In these works, the components of the macroscopic energy-momentum tensor in the fluid rest frame are used to obtain the energy density $\varepsilon$, parallel ($P_{\parallel}$) and perpendicular ($P_{\perp}$) pressures. In Heaviside-Lorentz natural units, the pure electromagnetic contribution to the energy-momentum tensor, which is anisotropic, has values of $B^{2}/2$ and $-B^{2}/2$ in the perpendicular and parallel directions to the local magnetic field, respectively. However, this approach can drastically overestimate the mass of neutron stars (as shown in Fig. 3 of Gomes et al., 2017).

Several works obtained the global structure models of magnetars by solving coupled Einstein-Maxwell equations, taking into account the anisotropy of the stress-energy tensor (Bocquet et al., 1995; Cardall et al., 2001; Oron, 2002; Ioka and Sasaki, 2004; Kiuchi and Yoshida, 2008; Yasutake et al., 2010; Frieben and Rezzolla, 2012; Yoshida et al., 2012; Pili et al., 2014; Thapa et al., 2020; Tsokaros et al., 2022). In these studies either a perfect fluid, a polytropic EoS, or a realistic EoS was assumed, but do not take into account the magnetic field modifications due to its quantisation.

Ideally, to explore magnetic field effects such as Landau quantisation and AMM on the global properties of the star, one must solve the coupled Einstein-Maxwell equations, along with a magnetic field dependent EoS. In  Chatterjee et al. (2015) and Franzon et al. (2016), global numerical models for magnetars were obtained by consistently solving Einstein-Maxwell equations with magnetic field dependent quark EoS. It was however explicitly demonstrated by Chatterjee et al. (2015, 2021) that the maximum mass of a neutron star is minimally modified due to the magnetic field dependence of the microscopic EoS, even for the highest magnetic fields. Therefore in this work, we assume a non-magnetic ($B=0$) matter contribution to the EoS to compute global neutron-star models and the magnetic field enters structure calculations only through the dominant pure electromagnetic field contribution. Although it remains to be checked explicitly in future work, the effects of Landau quantisation and AMMs are not expected to sensibly affect the results of this study. Note however, that this is not the case for microscopic properties of matter, as discussed in the following.

We start our discussion with the single-particle interaction potential for the $\Delta$ baryons in dense nuclear matter, which is a measure of the energy cost of adding one particle $b$ in a $b$-less matter with the given condition. For the L3$\omega\rho$ model, it can be written as

and, for the CMF model,

In isospin-symmetric nuclear matter, all families of baryons experience the same potential, since the meson field $\rho$ (and $\delta$) are null in this situation. In neutron rich matter, particles with positive isospin projections (as the positively charged $\Delta$s) are more bound than their zero- and negative-isospin counterparts, with the largest difference occurring for pure neutron matter, that can be taken as an extrapolation of neutron-star matter in $\beta$-equilibrium. The first two panels of Fig. 1 show how the L3$\omega\rho$ model scalar and vector interactions affect the $\Delta$ potentials. In all cases, the particle potentials eventually become positive as the density increases, corresponding to unbound states, but they stay negative in the relevant interval of densities around nuclear saturation, where their depth determines how much they are bound.

For the L3$\omega\rho$ model, the larger the scalar coupling value (i.e., the parameter $x_{\sigma\Delta}$), the lower the potentials are in the low density regions. Complementary to that, the larger the vector coupling (i.e., the parameter $x_{\omega\Delta}$), the more repulsive the potentials for $\Delta$s are, which reflect in more positive curves in the high density region, where the repulsive channel dominates. Also, it can be seen that the potential depends less on the species of $\Delta$ in the CMF model. Magnetic-field effects are not included, but it was verified that fields up to $B=3\times 10^{18}$ G do not affect these results.

For $B=0$, the matter EoS (namely, $P$ vs. $\varepsilon$) shows a simple monotonically increasing behavior, however its derivatives show interesting features generated, e.g., by changes in particle composition. Next, we discuss the incompressibility modulus (usually referred to simply as compressibility), given by

At saturation density, compressibility values for isospin-symmetric matter can be compared with laboratory data. We find values of 256 MeV and $300$ MeV for the L3$\omega\rho$ and CMF models, respectively. Laboratory values range between $220<K<260$ MeV (Dutra et al., 2014; Lopes, 2022) and $250<K<315$ MeV (Stone et al., 2014).

The top panel of Fig. 2 shows the effect of the inclusion of different particle species in the compressibility for the L3$\omega\rho$ model, in the absence of an external magnetic field. The kinks in the curves are consequence of the onset of new particle species, which are shifted to lower densities by the inclusion of both $\Delta$s and respective stronger scalar interactions. For $x_{\sigma\Delta}=1$, the effective mass of nucleons becomes zero at $n_{B}\sim 0.85$ fm^-3 and, for this reason we lack solutions at higher densities. The bottom panel of Fig. 2 shows that in the CMF model the kinks are much smaller than in the L3$\omega\rho$ model, with the only displacement of the curve occurring at the onset of the first non-nucleon baryon. As a consequence, the different CMF EoSs behave more similarly as the density increases.

The stiffer EoSs are formally the ones with larger values of the speed of sound $v_{s}$, but here we discuss stiffness with respect to $K$, related to $v_{s}$ through $v_{s}^{2}=K/(9\mu)$ ((Dexheimer et al., 2008)). Isospin symmetric matter is softer at low densities, but becomes stiffer at large densities due to the Pauli exclusion principle because, as only nucleonic matter is considered, higher Fermi levels must be occupied (see Fig. 2). The behavior of neutron-star matter (charge neutral and in chemical equilibrium) depends on the composition, but it is always softer than the symmetric matter case after the hyperon or $\Delta$ onsets, as the presence of new Fermi levels turns the EoS softer. Matter with hyperons but no $\Delta$s is stiffer at intermediate densities (than matter with $\Delta$s), however it is softer at large densities, especially in the case of strong scalar interaction (for the L3$\omega\rho$ model). This trend was noticed previously by Dexheimer et al. (2021a), where we showed that the inclusion of $\Delta$s could turn the EoS stiffer (than the cases where they were absent), despite the fact that the new degrees of freedom soften the EoS. This is related to isospin asymmetry, which we discuss in the following.

We define the isospin fraction as the average 3^rd isospin component of a given matter composition, weighted by the relative densities, i.e.,

as shown in Fig. 3. For nucleonic matter only, $Y_{I_{3}}=0$ means matter with the same amount of protons and neutrons, while $Y_{I_{3}}=-0.5$ means pure neutron matter. The density at which the curves with and without $\Delta$s split marks the appearance of the $\Delta^{-}$s, which increase the isospin asymmetry (turn the isospin fraction more negative). The effect is much larger for the L3$\omega\rho$ model (top panel) than the CMF model (bottom panel), which hints that the amount of $\Delta$s reproduced in each model is different. Both effects generated by the magnetic fields, i.e. Landau quantization and AMM, decrease the isospin asymmetry (less negative $Y_{I_{3}}$) at low and intermediate densities.

A better understanding of the effects of the inclusion of $\Delta$ baryons, magnetic fields, and AMM in neutron-star matter subject to strong magnetic fields can be obtained from Fig. 4. Comparing the top row ($B=0$) with the lower one ($B=3\times 10^{18}$ G), we can see that some of the charged particles are favored when magnetic field effects without AMMs are considered, an effect that is more pronounced for protons, whose onset density is pulled to very low densities for both models. As a consequence, their population becomes more similar to the neutron one in densities below $\sim 0.05$ fm^-3, turning $Y_{I_{3}}$ less negative. The inclusion of AMM enhances this effect. This explains why the isospin asymmetry depends both on the magnetic field and on the AMM in the lower density region, as shown in Fig. 3. The $\Delta^{-}$ threshold (at densities around $0.3$ fm^-3) coincides with the region at intermediate densities beyond which the N+H+$\Delta$ EoSs are softer than the respective N+H EoSs. The $\Lambda$ (and the $\Sigma$ in the CMF model) hyperons appear at larger densities than the $\Delta^{-}$s. The remaining $\Delta$s appear at much larger densities and in amounts that depend on the interactions in the L3$\omega\rho$ model.

To discriminate AMM effects on the particle composition is not trivial, as they depend on the AMM coupling strength and sign, on the particle mass, charge, and density. Additionally, different spin projections are separately enhanced or suppressed, but this cannot be clearly seen in Fig. 4, as it follows the usual convention and shows the sum of all spin projections for each particle. For this reason, we make use of a quantity that reveals the degree of spin polarization, more suited to discuss spin projection asymmetry of fermions.

We define the total spin polarization of a given matter composition, weighted by the relative densities, in analogy to Eq. (29), i.e.,

and shown the results in Fig. 5. For a fixed magnetic field strength, all charged particles are fully spin polarized at low densities: only spin projection 1/2 for protons and spin projection 3/2 for positive $\Delta$s, only spin projection -1/2 for leptons and negative $\Sigma$s, and spin projection -3/2 for negative $\Delta$s. When AMMs are considered, neutral particles obey the same logic, presenting only positive (negative) spin projections according to their positive (negative) sign of $\kappa_{b}$. At intermediate densities, full polarization is broken for more massive particles, but not for leptons and $\Lambda$s. But, regardless, the polarization never goes to zero, meaning that partial spin projection imbalance remains at high densities. Overall, spin polarization fraction is much stronger for the CMF model (bottom panel) than for the L3$\omega\rho$ model (top panel). Full polarization can be understood from Eqs. (11), (12), and (13), which explains why particles with different isospin projections present different momenta and why particles occupying the first Landau level ($\nu=0$) are more abundant when only a few levels are occupied. This happens for strong magnetic fields and low particle densities, or simply less massive particles.

It is a well-established concept that the magnetic field is not constant within neutron stars, but increases towards their centers where the density is larger. But, before we discuss stellar configurations with macroscopic magnetic fields in detail, we study how one more relevant quantity changes as a function of magnetic field strength. The fraction of exotic particles can be defined as the following quantity

for $i=H$ or $\Delta$, shown in Fig. 6. On the left panel for the L3$\omega\rho$ model, the amount of $\Delta$s is slightly reduced at a given density but then increases tremendously at the larger density when the AMM is included, a behaviour quantitatively not reproduced with larger coupling constants, as seen on the middle panel. The amount of hyperons, on the other hand, is not significantly modified by the magnetic field, only slightly decreases in the presence of AMM and is affected by the small fluctuations related to the De Haas-Van Alphen oscillations (De Haas and Van Alphen, 1930). The right panel shows the same qualitative behavior for the CMF model, which has a more clear substitution of hyperons in favor of deltas for higher values of $B$, independently of the density or accounting for the AMM.

In a previous study, Dexheimer et al. (2021a) obtained results on the effects of the inclusion of $\Delta$s in neutron stars for both the $L3\omega\rho$ and CMF models without magnetic fields, using standard TOV equations. As expected from the discussion regarding the hyperon puzzle (Chamel et al., 2013), pure nucleonic stars are always the most massive configuration, and the inclusion of hyperons decreases the maximum mass obtained. It was noticed that the inclusion of $\Delta$s, however, does not modify significantly (and in some cases increases) the maximum stellar mass in relation to the composition with only nucleons and hyperons, and this effect is more obvious for the cases where $\Delta$s are more abundant. It is argued that the addition of $\Delta$s decreases the fraction of nucleons (in fact, neutrons) and hyperons (mostly $\Lambda$s) to create $\Delta$s and some protons, in a way that the overall increase in isospin asymmetry turns the EoS stiffer, even when more species are present.

As described here in section II.3, to compute the effect of the strong magnetic fields on the structure of the magnetars, one must solve the coupled Einstein–Maxwell equations with the equations of state (described in section III.1). For the chosen poloidal field geometry, we solve the Einstein–Maxwell equations within the numerical library LORENE²²2http://www.lorene.obspm.fr using a multi-domain spectral method. In Fig 7, we show the mass radius relations for the $L3\omega\rho$ and the CMF models, with and without $\Delta$s, as a function of equatorial radius for sequences of constant stellar central magnetic field. Despite the fact that, for the choices discussed here for $x_{\sigma\Delta}$ and $x_{\omega\Delta}$ parameters in the $L3\omega\rho$ model, the masses of N+H+$\Delta$ stars never surpass the respective N+H configurations, we still observe that the the maximum mass (shown in Fig. 7) follows the same ordering of a large (and most relevant) portion of Fig. 3 for the compressibility.

In Fig. 7, any differences between the mass-radius curves for the $B=0$ case (solid lines) arise from the differences in the (non-magnetic) EoS, while the differences with magnetic field come from the pure electromagnetic field contribution. We know that the Lorentz force originating from the pure electromagnetic field affects the low density part of the EoS. This is why the maximum mass of very massive stars does not change with increasing magnetic field strength, but the mass and radius of less massive stars increase significantly. For the $L3\omega\rho$ model, the inclusion of $\Delta$s decreases modestly the maximum stellar mass, especially for the larger coupling. However, for the CMF model, we do not see meaningful changes on the mass-radius diagram with the inclusion of $\Delta$s. From Table 2 we see that, keeping the radius of the neutron star fixed (going up vertically in Fig. 7), the increase in the strength of the central magnetic field increases both the central baryon and energy densities, as a larger matter pressure is necessary to balance the Lorentz force. The addition of $\Delta$s decreases both quantities, as these stars are naturally (at B=0) smaller.

At this point, we note that the maximum mass value of the stellar family described by the $L3\omega\rho$ model with $\Delta$s and $x_{\sigma\Delta}=1.0$ (the yellow curve) has not been attained because, as explained earlier, the EoS numerical code stops converging at large densities due to reaching zero nucleon masses. Such behavior indicates that hadronic matter is no longer stable at this point and deconfinement to quark matter must be considered. We leave such analysis to a future work. But, since the trend of the yellow curve is quite obvious, we can conclude that its maximum mass is lower when compared to the other coupling and composition.

Using the full numerical solution, we also study the effect of the EoS on the magnetic field configurations inside a given star. We decompose the magnetic field norm in terms of spherical harmonics

and plot the first four even multipoles ($l=0,2,4,6$) as function of coordinate radius for both the EoS models and coupling strengths in Fig 8. We also plot the profile of the dominant monopolar, spherically symmetric, term ($l=0$) inside the star in Fig 9. For $L3\omega\rho$ model, specially if we include $\Delta$s, the magnetic field norm decreases slightly inside the star but, for CMF model, we do not see any considerable changes.