General Relativistic Mean-Field Dynamo Model for Proto-Neutron Stars

Abstract

Neutron stars, and magnetars in particular, are known to host the strongest magnetic fields in the Universe. The origin of these strong fields is a matter of controversy. In this preliminary work, via numerical simulations, we study, for the first time in non-ideal general relativistic magnetohydrodynamic (GRMHD) regime, the growth of the magnetic field due to the action of the mean-field dynamo due to sub-scale, unresolved turbulence. The dynamo process, combined with the differential rotation of the (proto-)star, is able to produce an exponential growth of any initial magnetic seed field up to the values required to explain the observations. By varying the dynamo coefficient we obtain different growth rates. We find a quasi-linear dependence of the growth rates on the intensity of the dynamo. Furthermore, the time interval in which exponential growth occurs and the growth rates also seems to depend on the initial configuration of the magnetic field.

1. Introduction

A neutron star (NS) is the smallest and densest star in the Universe, with radius of the order of 10 km and mass between about 1.1 and 2.17 solar masses (e.g., [1,2,3]). Together with black holes, it is one of the possible compact objects in which the core of a giant star—a star with a mass greater than eight solar masses—can collapse at the end of its life cycle.

NSs show a vast and diverse phenomenology. Their magnetic fields are between ${10}^8$ and ${10}^{15}$ times stronger than Earth’s one. Intense magnetic fields can deform the star, making it more or less oblate [4,5,6,7,8]), and give observable signals of gravitational waves (GWs), especially during merging events like GW170817, the first gravitational wave signal due to the merging of two neutron stars observed by the LIGO and Virgo detectors on 17 August 2017 [2,9,10].

A newly born, fast rotating and ultra-magnetized proto-NS can play as central engine driving the energetic emission of both Long and Short Gamma-Ray Bursts (GRBs) [9,11,12,13,14,15], the luminous electromagnetic events known to occur in the Universe. Once formed, rapidly rotating and strongly magnetized NSs would lose rotational kinetic energy in a short time, on a timescale of seconds or less. A relativistic plasma, driven the combination of the rotation with the magnetic field, flows away from the star, and its emission of X-ray and $\gamma$-ray in the photosphere can reproduce the observational characteristics of a GRB [11].

The origin of these strong magnetic fields of neutron stars—and in particular those of magnetars, which have magnetic fields of the order of ${10}^{14}$–${10}^{15}$ gauss (G)—is a matter of controversy and many evolutionary scenarios have been proposed so far. In the original magnetar model [16] the strength of the magnetic field is the result of the growth of the seed field amplified by the dynamo mechanism—a process through which a rotating, convecting, and electrically conducting fluid can generate a magnetic field by self-inductive action converting kinetic energy into magnetic energy—in the proto-neutron star (PNS), the remnant of the supernova explosion which represents the first phase of life of a NS. The amplification of the magnetic field by a turbulent dynamo can occur either by a convective dynamo (e.g., [17,18]) or the magnetorotational instability (MRI, e.g., [19,20,21,22]).

Most of the studies done so far have focused on the so-called mean-field dynamo theory (e.g., [18,23,24,25]), in which the magnetic field is amplified by an electromotive force generated by the coupling between velocity and magnetic field small scale turbulent fluctuations (which we do not solve and must be 3D to generate dynamo)—almost certainly induced by MRI. However, all the studies conducted so far have been done using classical dynamo equations. A fully covariant generalization of equations in the general relativistic magnetohydrodynamical (GRMHD) regime was first proposed by Bucciantini & Del Zanna in [26] using the so-called 3 + 1 formalism.

In this preliminary work we examine for the first time with numerical simulations in non-ideal GRMHD the amplification of a seed field by the mean-field dynamo in the so-called kinematic approximation, in which the star is assumed to be in stationary hydrodynamic equilibrium. As a result, the velocity field is taken as pre-assigned and time-independent. Our goal is to verify the feasibility of the proposed mechanism and to estimate the average growth rate of the magnetic field, investigating its dependence on the intensity of the dynamo.

This paper is organized as follows. In Section 2, we briefly present the clasical dynamo theory and its fully covariant generalization proposed by Bucciantini & Del Zanna. In Section 3, we describe our model and its initialization. In Section 4 we show and discuss the results of our numerical simulations. Finally, we conclude in Section 5.

2. The Mean-Field Dynamo

The dynamo theory describes the process through which a rotating, convecting, and electrically conducting fluid can generate a magnetic field by self-inductive action converting kinetic energy into magnetic energy and maintain it over astronomical time scales. In particular, motion through a magnetic field induces an electric field that generates a current according to Ohm’s law. The result is a highly nonlinear system, in which this current in turn produces a magnetic field which on the one hand generates a new electric field through Faraday’s law and on the other can oppose motion by means of the Lorentz force. However, the onset of dynamo action can be studied in the linear approximation, i.e., the velocity field is assumed to be given (kinematic problem).

We will concentrate on the process known as mean-field dynamo in the magnetohydrodynamical (MHD) regime, which has been applied to a large number of astrophysical contexts, as the Sun [27], neutron stars (e.g., [23]), and accretion disks (e.g., [28,29]).

In the following we assume a signature $(-,+,+,+)$ for the space-time metric and we use Greek letters (running from 0 to 3) for 4D space-time components, and Latin letters (running from 1 to 3) for 3D spatial components. Moreover, we set $G=c=M_{\odot }=1$ (where $M_{\odot }$ is the solar mass). Finally, for the electromagnetic quantities we make use of the Lorentz-Heaviside units, where ${\epsilon }_0={\mu }_0=1$.

2.1. Classical Theory

In the classical MHD, the time evolution of the electric field can be neglected in the Ampère-Maxwell equation, which can thus replaced by the Ampère’s law $\mathbf{\nabla }\times \mathbf{B}=\mathbf{J}$. Substituting Ohm’s law

$$\mathbf{J}=\sigma {\mathbf{E}}^{\prime }{\mathbf{E}}^{\prime }=\mathbf{E}+\mathbf{U}\times \mathbf{B}$$

(1)

where $\sigma$ is the electric conductivity, U is the fluid velocity and ${\mathbf{E}}^{\prime }$ is the electric field in the frame comoving with the fluid—into the Faraday’s law of induction, one can write a single evolution equation for B, which is called the induction equation:

$${\partial }_t\mathbf{B}=\mathbf{\nabla }\times (\mathbf{U}\times \mathbf{B}-{\eta }_d\mathbf{J})$$

(2)

where ${\partial }_t=\partial /\partial t$ and ${\eta }_d={\sigma }^{-1}$ is the magnetic diffusivity. As mentioned in Section 1, in the mean-field dynamo theory, if small scale turbulent fluctuations of the magnetic field are well developed, the coupling between velocity and magnetic field fluctuations generates an electromotive force that can be written as

$$\langle \delta \mathbf{v}\times \delta \mathbf{B}\rangle ={\alpha }_{turb}\mathbf{B}-{\beta }_{turb}\mathbf{J}$$

(3)

The first term is the so-called $\alpha$-term, describing the creation of poloidal fields starting from toroidal fields, while the second term describes the strongly enhanced dissipation of the field by the turbulence. In general both ${\alpha }_{turb}$ and ${\beta }_{turb}$ will be tensors, however we will focus here on the isotropic case where they can be dealt with as scalars. Thanks to the combined action of these two terms, this electromotive force amplifies the initial magnetic field according to the mechanism known as $\alpha \Omega$-dynamo, where the $\Omega$ effect refers to the amplification of the toroidal field by shear (i.e., differential rotation). It is possible to prove (e.g., [26]) that the comoving electric field can be rewritten as

$${\mathbf{E}}^{\prime }=\xi \mathbf{B}+\eta \mathbf{J}$$

(4)

where $\xi =-{\alpha }_{turb}$ and $\eta ={\eta }_d+{\beta }_{turb}$.

2.2. Mean-Field Dynamo in 3 + 1 Resistive GRMHD

As mentioned in Section 1, the first fully covariant generalization of Equation (3) in GRMHD was first proposed by Bucciantini & Del Zanna in [26] (see also [30]). Extending the usual form adopted for resistive GRMHD, the expression for the comoving electric field appears to be

$$e^{\mu }=\xi b^{\mu }+\eta j^{\mu }$$

(5)

where $j^{\mu }$ and $b^{\mu }$ are the comoving conduction current and the comoving magnetic field, respectively. While the covariant form of equations is very important—and quite elegant—from a theoretical point of view, for numerical convenience time and space are better treated if split. Using the so-called $3+1$ formalism—in which the generic line element can be written as

$$\mathrm{d}{s}^2=-{\alpha }^2\mathrm{d}{t}^2+{\gamma }_{ij}(\mathrm{d}{x}^i+{\beta }^i\mathrm{d}t)(\mathrm{d}{x}^j+{\beta }^j\mathrm{d}t)$$

(6)

where $\alpha$ is known as the lapse function, $\beta$ as the shift vector, and ${\gamma }_{ij}$ is the 3-metric, used to raise/lower the indexes of any spatial three-dimensional vector or tensor—Maxwell’s equations become

$$\begin{array}{l}{\partial }_t\left(\sqrt{\gamma }{E}^i\right)-\sqrt{\gamma }{ϵ}^{ijk}{\partial }_j(\alpha B_k-{ϵ}_{klm}{\beta }^l{E}^m)=-\sqrt{\gamma }(\alpha J^i-q{\beta }^i)\\ {\partial }_t\left(\sqrt{\gamma }{B}^i\right)+\sqrt{\gamma }{ϵ}^{ijk}{\partial }_j(\alpha E_k+{ϵ}_{klm}{\beta }^l{B}^m)=0\\ {\partial }_i\left(\sqrt{\gamma }{E}^i\right)=\sqrt{\gamma }q\\ {\partial }_i\left(\sqrt{\gamma }{B}^i\right)=0\end{array}$$

(7)

where $\gamma =det{\gamma }_{ij}$ and ${ϵ}_{ijk}$ are respectively the determinant and the Levi-Civita pseudo-tensor of the 3-metric. We note that we must integrate the electric field over time and the current J appears as a source term. Writing Equation (5) in the 3 + 1 form we find

$$J^i=q{v}^i+{\sigma }_E\Gamma (E^i+{ϵ}^{iJk}{v}_J{B}_k-\left(v_k{E}^k\right)v^i)+{\sigma }_B\Gamma (B^i-{ϵ}^{iJk}{v}_J{E}_k-\left(v_k{B}^k\right)v^i)$$

(8)

where $\Gamma$ is the Lorentz factor, q is the electric charge density, ${\sigma }_E=1/\eta$, and ${\sigma }_B=-\xi /\eta$. Due to the presence of $\eta$ in the denominator, stiff source terms arise in the equation for the evolution of the electric field, requiring some sort of implicit treatment (e.g., [31]).

3. Model and Numerical Set-Up

The fluid configuration of our PNS is in axisymmetric equilibrium and has been modeled via the XNS code [7,32,33] with an uniform grid in both radial (256 points in the range $[0,30]$ (geometrized units), with the star resolved in about half of the points) and $\theta$ (64 points in the range $[0,\pi ]$) directions. A polytropic relationship between pressure and density is, i.e.,

$$p=K{\rho }^{\gamma }$$

(9)

where p is the thermal pressure, $\rho$ is the mass density, and $\gamma$ is the adiabatic index.

Table 1. Parameters of the fluid model in geometrized units ($c=G=M_{\odot }=1$). The corresponding values in cgs units are reported in the text.

K	$\mathit{\gamma }$	${\mathit{\rho }}_{\mathit{c}}$	${\Omega }_{\mathit{c}}$	${\mathit{A}}^2$	${\Omega }_{\mathit{eq}}$	${\mathit{R}}_{\mathit{eq}}$
100	2	$1.28\times {10}^{-3}$	$2.633\times {10}^{-2}$	70	$9.660\times {10}^{-3}$	8.38

shows all parameters defining the fluid configurations. Here ${\rho }_c\approx 7.9\times {10}^{14}$ g/cm${}^3$ is the central density, ${\Omega }_c\approx 5.3$ kHz is the central rotation rate and ${\Omega }_{eq}\approx 2.0$ kHz is the equatorial one, $R_{eq}\approx 12.4$ km is the equatorial radius of the star, and A is a measure of the differential rotation rate and is such that (see [32])

$$A^2({\Omega }_c-\Omega )={\displaystyle \frac{{\mathcal{R}}^2(\Omega +{\beta }^{\varphi })}{{\alpha }^2-{\mathcal{R}}^2{(\Omega +{\beta }^{\varphi })}^2}}$$

(10)

where $\mathcal{R}=\sqrt{{\gamma }_{\varphi \varphi }}$. The corresponding central spin period is $P_c\approx 1.2$ ms $\ll P_{crit}$, where $P_{crit}\sim$ 1–1.5 s is the the critical value corresponding to the upper limit of the rotation period of the PNS which allows the turbulent mean-field dynamo to develop [23]. The choice of this high central rotation rate is due to the desire to maximize the $\Omega$ effect. Note that the equation of state (EoS) and the rotation profile we use are more appropriate for an NS than for a PNS (e.g., [22]). However, the aim of this preliminary work is to verify whether the dynamo model could be a realistic scheme for generating proto-magnetar, taking into account, for the first time, non-ideal effects in general relativity. The results can be easily extended to PNS, for which one should change the EoS and rotation profile.

The study of the time evolution of the system is performed via the ECHO code [34]. Therefore, at the initial time (i.e., at $t=0$), in addition to the configuration of the stellar fluid just described, a small magnetic field (added “by hand”) and the profiles of the $\xi$ and $\eta$ parameters must be inserted (see below).

For the magnetic field we use a purely toroidal initial configuration. As shown in [32], the toroidal field can be written as

$$B^{\varphi }={\displaystyle \frac{K_m{\left({\alpha }^2{\mathcal{R}}^2\rho h\right)}^m}{\alpha {\mathcal{R}}^2}}$$

(11)

where $h=1+\gamma /(\gamma -1)p/\rho$ is the specific enthalpy, $K_m$ is the toroidal magnetization constant, which regulates the strength of the magnetic field (more specifically the magnetic flux through the meridional plane), and $m\ge 1$ is the toroidal magnetization index (related to the distribution of the magnetic field inside the star). We choose $K_m=5\times {10}^{-5}$ and $m=1.5$ so as to have, in physical units, an initial magnetic field of the order of ${10}^{12}$ G, i.e., the typical value inferred for pulsars. The electric field was instead initialized using the relation

$$\mathbf{E}=-\mathbf{v}\times \mathbf{B}$$

(12)

of the ideal MHD. Therefore, the initial electric field is zero with a purely toroidal initial field (as assumed here by us), but not with an initial poloidal field, as we will assume in Section 4.2.

In order to quantify the dynamo action, we can introduce the two characteristic (dimensionless) numbers (e.g., [24])

$$C_{\xi }={\displaystyle \frac{{\xi }_{0}R}{\eta }}{C}_{\Omega }={\displaystyle \frac{\Delta \Omega R^2}{\eta }}$$

(13)

where $\Delta \Omega ={\Omega }_c-{\Omega }_{eq}$, $R=R_{eq}$, and ${\xi }_0$ is the maximum value of $\xi$. These numbers may be regarded as magnetic Reynolds numbers based on the intensity of the $\alpha$-effect and the differential rotation, respectively. The $\xi$ and $\eta$ profiles are chosen so that the diffusion and dynamo processes occur only within the star. While we assume constant $\eta$ within the star, the profile of $\xi$ is similar to that used by Bonanno et al. in [23], namely

$$\xi =\left\{\begin{array}{ll}{\xi }_0{f}_{\xi }(r,\theta ){\mathcal{Q}}_{\xi }\left(\mathbf{B}\right)& \mathrm{inside}\mathrm{the}\mathrm{star}\\ 0& \mathrm{in}\mathrm{the}\mathrm{atmosphere}\hfill \end{array}\right.$$

(14)

where

$$f_{\xi }(r,\theta )={\displaystyle \frac{\mathrm{s}}{2}}sin\theta [1+\mathrm{erf}\left({\displaystyle \frac{r-R_c}{\Delta R}}\right)]$$

(15)

with erf the error function,

$$\mathrm{s}=\mathrm{s}\left(\theta \right)=\left\{\begin{array}{ll}+1& \mathrm{for}{\displaystyle \frac{\pi }{2}}-\theta \ge 0\\ -1& \mathrm{for}{\displaystyle \frac{\pi }{2}}-\theta <0\end{array}\right.$$

(16)

and $\Delta R=0.025R$. The function

$${\mathcal{Q}}_{\xi }\left(\mathbf{B}\right)={[1+{\displaystyle \frac{{\mathbf{B}}^2}{B_{eq}^2}}]}^{-1}$$

(17)

is known as quenching function. In the kinematic approximation, quenching is necessary to prevent the infinite growth of the magnetic field due to the absence of the fluid feedback. Here $B_{eq}$ is the equipartition magnetic field with respect to the kinetic energy density of turbulent motions. We assumed

$$B_{eq}^2={\alpha }_{q}p$$

(18)

where ${\alpha }_q$ is a constant and p is the thermal pressure. As mentioned above, relation (14) is similar to the one assumed by Bonanno et al. in [23]. The difference is that we used the sine function with respect to the cosine function to have the maximum dynamo parameter near the equator of the star. The antisymmetry across the equatorial plane is guaranteed by the $\mathrm{s}$ function.

From a physical point of view, relation (14) implies that the star is divided into two regions: an innermost one where the dynamo is negligible and an outermost one where it is maximum. This two regions present two different instabilities: a convective one is active in the inner regions of the star (of radius $R_c$) and a so-called neutron-finger instability, due to the lepton number gradients (see [23,25]), in the outer regions. The two regions are separated by a thin layer—named tachocline—of thickness $\Delta R$ through which the physical properties of the two regions are supposed to vary in a smooth way. In general, $R_c$ varies over time. Taking $R_c$ constant during the entire evolution is equivalent to assume that $R_c$ varies in much longer times than those in which the magnetic field grows exponentially due to the dynamo.

We choose five different values of ${\xi }_0$. The list of runs with the corresponding parameters is reported in

Table 2. Parameters of the fluid model in geometrized units ($c=G=M_{\odot }=1$).

	${\mathit{\xi }}_0$	$\mathit{\eta }$	${\mathit{C}}_{\mathit{\xi }}$	${\mathit{C}}_{\Omega }$	${\mathit{R}}_{\mathit{c}}/\mathit{R}$	${\mathit{\alpha }}_{\mathit{q}}$
Run1	$1.0\times {10}^{-2}$	$1.0\times {10}^{-3}$	$0.084\times {10}^3$	$1.170\times {10}^3$	0.4	$1.0\times {10}^{-5}$
Run2	$5.0\times {10}^{-2}$	$1.0\times {10}^{-3}$	$0.419\times {10}^3$	$1.170\times {10}^3$	0.4	$1.0\times {10}^{-5}$
Run3	$1.0\times {10}^{-1}$	$1.0\times {10}^{-3}$	$0.838\times {10}^3$	$1.170\times {10}^3$	0.4	$1.0\times {10}^{-5}$
Run4	$2.0\times {10}^{-1}$	$1.0\times {10}^{-3}$	$1.676\times {10}^3$	$1.170\times {10}^3$	0.4	$1.0\times {10}^{-5}$
Run5	$3.0\times {10}^{-1}$	$1.0\times {10}^{-3}$	$2.514\times {10}^3$	$1.170\times {10}^3$	0.4	$1.0\times {10}^{-5}$

. The value of ${\alpha }_q$ is chosen so that, in physical units, $B_{eq}\sim {10}^{15}$ G, i.e., the standard value for magnetars. Recall that in this preliminary work we decided to vary only the most important quantity, i.e., the maximum value of $\xi$ (and therefore $C_{\xi }$) for the same star and the same profiles of $\xi$ and $\eta$.

4. Numerical Results

In this section we show the results of our $\alpha \Omega$-dynamo simulations performed via the ECHO code using a third order scheme in time and a fifth order scheme in space. Following the prescription of Tomei et al. in [29], we define the average on the star of any quantity $f=f(r,\theta )$ as

$$\langle f\rangle ={\displaystyle \frac{{\int \int }_{\mathcal{D}}\mathrm{d}r\mathrm{d}\theta \alpha \sqrt{\gamma }f}{{\int \int }_{\mathcal{D}}\mathrm{d}r\mathrm{d}\theta \alpha \sqrt{\gamma }}}$$

(19)

where $\mathcal{D}$ identifies the region inside the star, $\alpha$ and $\gamma =det{\gamma }_{ij}$ are the lapse function and the determinant of the 3-metric defined in Section 2.2, respectively. The presence of the factor $\alpha$ is commonly invoked but still debated. For a better understanding, in this section we report the results in physical units.

Figure 1 shows the time evolution of the average poloidal and toroidal components of the magnetic field in the Run1. We can see that a poloidal field does not arise immediately, but only after a short phase in which the toroidal field remains constant. At this point the dynamo starts and, after a transient phase, supports the exponential amplification of the two components up to ∼3 ms, when the dynamo action is limited by the quenching effect. At the end of the exponential growth phase, the strength of the toroidal component is slightly higher than that of the poloidal one. The time interval in which the dynamo operates is much smaller than that on which the evolution of the neutron-finger instability occurs, which is of the order of seconds [35], in support of the hypothesis of constant $R_c$ during the action of the dynamo.

The time distribution on the $xz$-plane (any meridional cut of the star) of the components of the magnetic field is shown—at four different times—in Figure 2. We can see that the poloidal field arise near the tachocline and in the outer region of the star (where the dynamo parameter is maximum), mainly near the equator (Figure 2a). During the dynamo action, magnetic islands corresponding to the linear dynamo modes migrate towards the poles and the surface of the star while growing in amplitude (Figure 2b,c). This is similar to the migration of sunspots towards the equatorial plane on the surface of the Sun, but in the opposite direction, due to the fact that we have set $\xi =-{\alpha }_{turb}$ (see Section 2). During the post-exponential growth phase, the field is mainly located—along the z-axis—in a region about $0.6$ times the polar diameter (Figure 2d). This evolution is explained as follows. Due to the presence of curls in Maxwell’s equations, the gradients of the various quantities have a strong impact on the time evolution of the system. In fact already at small t the magnetic field only grows where there are the gradients of $\xi$, which are stronger than any other in play. Then the various components of the magnetic field are created, amplified and migrate, but the “imprint” of the gradients of $\xi$ remains.

4.1. Dependence on the $\alpha$-Dynamo Number

We now investigate the dependence of the results on the $\alpha$-dynamo number $C_{\xi }$. The list of runs with the corresponding growth rates for both poloidal (${\omega }_P$) and toroidal (${\omega }_T$) components is reported in

Table 3. Exponential growth rates (in ms${}^{-1}$) for both poloidal (${\omega }_P$) and toroidal (${\omega }_T$) components of the magnetic filed for all runs. The corresponding value of the $\alpha$-dynamo number $C_{\xi }$ is also reported.

	${\mathit{C}}_{\mathit{\xi }}$	${\mathit{\omega }}_{\mathit{P}}$ [ms${}^{-1}$]	${\mathit{\omega }}_{\mathit{T}}$ [ms${}^{-1}$]
Run1	$0.084\times {10}^3$	$2.75\pm 0.02$	$2.67\pm 0.03$
Run2	$0.419\times {10}^3$	$34.6\pm 0.2$	$33.5\pm 0.1$
Run3	$0.838\times {10}^3$	$45.8\pm 0.3$	$45.9\pm 0.2$
Run4	$1.676\times {10}^3$	$152\pm 2$	$146\pm 1$
Run5	$2.514\times {10}^3$	$169\pm 4$	$163\pm 3$

.

Figure 3 shows the dependence of the exponential growth rates with the dynamo number $C_{\xi }$. Blue (red) squares (stars) for the poloidal (toroidal) field. We observe that both rates follow a quasi-linear dependence on the dynamo number, as expected, with a slope of about $(7\pm 1)\times {10}^{-2}$ ms${}^{-1}$. For small values of $C_{\xi }$ there is not much difference between the two rates, while for large values of $C_{\xi }$ a slight difference is visible.

As can be seen from the data shown in Table 3, the most realistic values for the dynamo parameter are those used for the first runs (or lower). The highest values, used in the latest runs, are used for the sole purpose of determining the trend of growth rates on the dynamo parameter $C_{\xi }$ and a possible drastic variation in the slope (or even a non-linear dependence) for these limit values. In the near future we plan to repeat the study for lower values of $C_{\xi }$, in order to confirm or not this quasi-linear trend of growth rates even for less intense dynamo and the weight that the higher values of $C_{\xi }$ had on the estimate of the slope found in this work.

4.2. Dependence on the Initial Configuration of the Magnetic Field

We repeat here Run1 with an initial configuration for the magnetic field different from that used for the previous runs in order to show that the dynamo occurs independently of the initial configuration of the magnetic field and to verify whether or not there is a dependence of the growth rates on the latter. To do this, we consider a purely dipolar instead than toroidal magnetic field. The field is initialized using the analytical solution found by Mastrano et al. in the Newtonian case in [36,37], namely

$${\tilde{\mathbf{B}}}_{dip}=\left\{\begin{array}{ll}{\displaystyle \frac{35}{4}}[(1-{\displaystyle \frac{6}{5}}{d}^2+{\displaystyle \frac{3}{7}}{d}^4)cos\theta {\hat{\mathbf{e}}}_r-(1-{\displaystyle \frac{12}{5}}{d}^2+{\displaystyle \frac{9}{7}}{d}^4)sin\theta {\hat{\mathbf{e}}}_{\theta }]& \mathrm{for}d<1\\ {\displaystyle \frac{2cos\theta {\hat{\mathbf{e}}}_r+sin\theta {\hat{\mathbf{e}}}_{\theta }}{d^3}}& \mathrm{for}d\ge 1\hfill \end{array}\right.$$

(20)

where $d=r/R_p$—with $R_p$ the polar radius of the star—and $\tilde{\mathbf{B}}=\mathbf{B}/B_0$, where $B_0$ sets the strength of the field overall, is continuous at $d=1$. Our definition of d is due to the fact that in our case the star is not completely spherical, but is slightly oblate, as can be seen in Figure 2, and we decided to study the evolution of the system exclusively within the star, where the quantities of interest are significantly appreciable. We choose $B_0={10}^{-9}$ (in geometrized units), corresponding to a maximum value of the magnetic field of the order of ${10}^{12}$ G. Figure 4 shows the time evolution of the average poloidal and toroidal components of the magnetic field. We see that the exponential growth phase begins and ends with a small delay compared to the previous case. Furthermore, just like in the previous case, in the post-exponential growth phase the toroial field has a higher intensity than the poloidal one.

The exponential growth rates are reported in

Table 4. Exponential growth rates (in ms${}^{-1}$) for both poloidal (${\omega }_P$) and toroidal (${\omega }_T$) components of the magnetic filed for Run1. The subscript ${}_t$ (${}_d$) is for the purely toroidal (dipolar) initial configuration of the magnetic field.

	${\mathit{C}}_{\mathit{\xi }}$	${\mathit{\omega }}_{\mathit{P}}$ [ms${}^{-1}$]	${\mathit{\omega }}_{\mathit{T}}$ [ms${}^{-1}$]
Run1${}_t$	$0.084\times {10}^3$	$2.75\pm 0.02$	$2.67\pm 0.03$
Run1${}_d$	$0.084\times {10}^3$	$2.69\pm 0.02$	$2.64\pm 0.01$

, where those for the purely toroidal initial configuration of the magnetic field are also reported. We can see that the growth rates in the case of purely dipolar initial field are lower than those obtained with the purely toroidal initial configuration.

5. Summary and Conclusions

In this work we investigate, for the first time by means of non-ideal axisymmetric GRMHD simulations, the mean-field dynamo process operating in proto-neutron stars. Our simulations are performed in the so-called kinematic approximation, in which the star is assumed to be in hydrostatic equilibrium at all times. Our star is assumed to have an internal structure described by a polytropic equation of state and to be rapidly rotating—with a central spin period of $P_c\sim P_{crit}/1000$, where $P_{crit}$ is the upper limit of the rotation period of the star which allows the dynamo to develop (see [23])—to maximize the $\Omega$-effect. We also assume a purely toroidal initial configuration for the magnetic field with a strength of the order of ${10}^{12}$ G.

Since a proto-neutron star is subject to two substantially different instabilities—a convective one active in the inner regions of the star, where the dynamo is probably negligible, and a neutron-finger instability in the outer regions, where the mean-field dynamo can operate [17,23]—we assume, for the dynamo parameter, the profile defined by Equation (14), similar to that assumed by Bonanno et al. in [23]. With this model we assume that the relevant physical properties of the two regions vary across a thin layer named tachocline. We choose five different values for the maximum value of the dynamo parameter, while the turbulent magnetic diffusivity is assumed constant within the star.

The study of the evolution of electric and magnetic fields under the action of the dynamo is performed in a stationary spacetime in spherical coordinates via the ECHO code [34] in the upgraded version to include non-ideal resistive and dynamo effects in the Ohm’s law [26,30]. We find that the exponential growth rates follow a quasi-linear dependence on the $\alpha$-dynamo number $C_{\xi }$ for both poloidal and toroidal components of the magnetic field, with a slope of about $(7\pm 1)\times {10}^{-2}$ ms${}^{-1}$. Moreover, the dynamo action occurs in a time interval that is much smaller than that on which the evolution of the neutron-finger instability occurs, in support of the hypothesis that the tachocline does not evolve during the action of the dynamo.

To show that the dynamo occurs independently of the initial configuration of the magnetic field and to determine a possible dependence of growth rates on the latter, we perform a further simulation with the lowest value of the dynamo parameter using a purely dipolar (instead than a purely toroidal) initial configuration for the magnetic field, implemented using the solution of Mastrano et al. [36,37]. We find that the initial and final instants of the exponential growth phase are delayed compared to those with the purely tororidal initial configuration. Furthermore, the growth rates are lower than the previous ones.

In the near future we plan to repeat the study for different configurations of the star (lower central rotation rate and/or different equation of state), and/or different profile of the dynamo parameter, and to study the relative importance of the $\alpha \Omega$- and ${\alpha }^2$-dynamo. We also plan to extend the study to the non-kinematic case, thus taking into account the feedback of the fluid. The study performed in this work can also be extended to hypothetical compact objects such as hybrid stars (i.e., neutron stars with a quark core) and quark stars (i.e., neutron star whose interior is totally composed by quark matter [38,39,40].