Instability of Magnetic Equilibria in Barotropic Stars

Instability of Magnetic Equilibria in Barotropic Stars

J.P. Mitchell , J. Braithwaite , A. Reisenegger, H. Spruit, J.A. Valdivia, and N. Langer
Argelander-Institut, University of Bonn, Auf dem Hügel 71, 53121 Bonn, Germany
Instituto de Astrofísica, Facultad de Física, Pontificia Universidad Católica de Chile, Av. Vicuña Mackenna 4860, 7820436 Macul, Santiago, Chile
Max-Planck-Institut für Astrophysik, Karl-Schwarzschild-Str. 1, D-85748 Garching, Germany
Departamento de Física, Facultad de Ciencias, Universidad de Chile, Casilla 653, Santiago, Chile
Accepted . Received 1; in original form

In stably stratified stars, numerical magneto-hydrodynamics simulations have shown that arbitrary initial magnetic fields evolve into stable equilibrium configurations, usually containing nearly axisymmetric, linked poloidal and toroidal fields that stabilize each other. In this work, we test the hypothesis that stable stratification is a requirement for the existence of such stable equilibria. For this purpose, we follow numerically the evolution of magnetic fields in barotropic (and thus neutrally stable) stars, starting from two different types of initial conditions, namely random disordered magnetic fields, as well as linked poloidal-toroidal configurations resembling the previously found equilibria. With many trials, we always find a decay of the magnetic field over a few Alfvén times, never a stable equilibrium. This strongly suggests that there are no stable equilibria in barotropic stars, thus clearly invalidating the assumption of barotropic equations of state often imposed on the search of magnetic equilibria. It also supports the hypothesis that, as dissipative processes erode the stable stratification, they might destabilize previously stable magnetic field configurations, leading to their decay.

(magnetohydrodynamics)MHD-stars; stars: magnetic fields; stars: neutron; stars: white dwarfs; stars
pagerange: Instability of Magnetic Equilibria in Barotropic StarsLABEL:lastpagepubyear: 2014

1 Introduction

There are several kinds of stars, namely upper main sequence stars, white dwarfs, and neutron stars, that, contrary to the Sun, have magnetic fields that are organized on large scales and persist unchanged over long time scales. Much of their interior is stably stratified, due to (a) gradients of entropy in white dwarfs and radiative envelopes of main sequence stars; and (b) gradients of composition (relative abundances of particle species) in the case of neutron stars. Clearly, no dynamo action is taking place in these stars, so their magnetic fields must be in stable hydromagnetic equilibrium states in which the Lorentz force is balanced by small perturbations to the otherwise spherical pressure and gravitational forces (with solid shear stresses in neutron star crusts possibly also playing a role).

It has long been known that purely toroidal (azimuthal) and purely poloidal (meridional) magnetic field configurations are unstable (Tayler, 1973; Markey & Tayler, 1973; Wright, 1973). However, 3-dimensional magneto-hydrodynamics (MHD) simulations of self-gravitating balls of highly conducting gas (Braithwaite & Spruit, 2004; Braithwaite & Nordlund, 2006) have shown their magnetic field to evolve naturally from an initially random configuration to a nearly axisymmetric structure in which poloidal and toroidal components of comparable amplitudes stabilize each other, persisting for many Alfvén times. Further work has suggested that the stable stratification of the fluid plays an important role in stabilizing these structures (Braithwaite, 2009; Akgün et al., 2013), and Reisenegger (2009) has conjectured that there are no stable magnetic configurations in stars that are barotropic, since they are not stably stratified.

On the other hand, attempting to construct plausible axisymmetric magnetic equilibria for these stars, several authors (Tomimura & Eriguchi, 2005; Yoshida et al., 2006; Haskell et al., 2008; Akgün & Wasserman, 2008; Kiuchi & Kotake, 2008; Ciolfi et al., 2009; Lander & Jones, 2009; Duez & Mathis, 2010; Pili et al., 2014) have imposed a barotropic equation of state, forcing the magnetic field structure to satisfy the non-linear and thus highly non-trivial Grad-Shafranov equation (Grad & Rubin, 1958; Shafranov, 1966). This assumption is not only very restrictive and unjustified for any of the stars likely to contain hydromagnetic equilibria (only very massive, radiation-dominated stars or extremely cold white dwarfs might come close to being barotropic); it might even be incompatible with having stable equilibria if the conjecture mentioned above is correct.

On the other hand, stable stratification will in the long term be eroded by non-ideal MHD processes such as heat diffusion (in main sequence and white dwarf stars), beta decays, and ambipolar diffusion (the latter two in neutron stars; see Reisenegger 2009). If the above conjecture is true, these processes would destabilize the magnetic field on these long time scales, making it decay (unless stabilized, e.g., by the solid neutron star crust).

Thus, it is important to clarify whether stable magnetic equilibria can exist in barotropic stars. A previous study by Lander & Jones (2012), based on a perturbative approach, gives a negative answer. The present work explores the same question more extensively, using full MHD simulations, in which we have evolved initially disordered magnetic fields as well as the ordered, axially symmetric, twisted-torus magnetic equilibria written analytically by Akgün et al. (2013), in both barotropic and (for comparison) stably stratified stars, in order to see whether stable equilibria can be reached.

Section 2 gives a brief introduction of the effect that stratification or its absence is expected to have on the stability of magnetic equilibria. Section 3 explains the numerical scheme used for our simulations. In Section 4, we explain how the disordered fields are created, and show the results of their evolution. In Section 5, we present the axially symmetric equilibria we use, as well as their evolution. The conclusions from our results can be found in Section 6.

2 Stable stratification vs. barotropy.

The physical effect of stable stratification can be illustrated by the thought experiment in which a fluid element is taken from a given position and moved vertically against gravity to a new position. Then, it is allowed to achieve pressure equilibrium with its new surroundings, without allowing it to change its entropy or chemical composition (which generally occur on much longer time scales). If the entropy or chemical composition vary with depth in the star, they will be different inside and outside the displaced fluid element, causing its density to be different, and thus creating a restoring force that is quantified by the squared Brunt-Väisälä (or buoyancy) frequency (force per unit mass per unit displacement)


Here, is the local acceleration of gravity, is the fluid mass density, is its pressure, and the subscripts “eq” and “ad” refer to the equilibrium profile existing in the star and to the changes produced in the adiabatic displacement of the fluid element, respectively. If , the restoring force tends to move the fluid element back to its initial position, causing the fluid to be stably stratified, whereas causes a runaway fluid element, resulting in a convective instability. If, on the other hand, the fluid is barotropic, there is a one-to-one relation between pressure and density, (i. e. because the specific entropy and/or composition are uniform throughout the fluid or adjust to an equilibrium as fast as the fluid element can move), the two partial derivatives will be equal, and the fluid will be neutrally stable.

In a hydromagnetic equilibrium state, the forces on a given fluid element must be balanced,


where is the effective gravitational potential (including a centrifugal contribution in a rotating star), is the electric current density, and is the magnetic field. From this, it follows that


In the barotropic case, and must be parallel everywhere, so the left-hand side vanishes, forcing the right-hand side to vanish as well, and thus imposing a strong constraint (three scalar equations) on the magnetic field configuration. This constraint is much weaker for the stably stratified case. In order to understand the latter, we make the astrophysically well justified assumption that the last term in Eq. (2) is much smaller than the other two, so the thermodynamic variables can be written as and , where and are their values in an unmagnetized star, and and are small perturbations caused by the Lorentz force, which can be regarded as independent because of the third, relevant thermodynamic variable entering their relation (small perturbations of specific entropy or chemical composition; see also Reisenegger 2009; Mastrano et al. 2011). In this case, , so the left-hand side of Eq. (3) is generally non-zero, only being constrained to be a horizontal vector field (perpendicular to ). Thus, only the vertical component of the curl on the right-hand side is required to vanish (a single scalar equation), a much weaker constraint. The particular, axially symmetric case is discussed in Section 5.

In realistic stars, the characteristic Alfvén frequency associated with magnetically induced motions (e. g., the growth rate of magnetic instabilities) is much smaller than , thus the buoyancy will produce a powerful restoring force, largely impeding magnetically induced vertical (radial) motions, and thus strongly constraining possible magnetic instabilities. If, on the other hand, we had (or ), there would be no such constraints on the motion, and the magnetic field could decay much more easily. This general idea, as well as more detailed arguments based on the same physics (Braithwaite, 2009; Akgün et al., 2013) motivate the conjecture that there will be stable magnetic equilibria only in stably stratified stars (Reisenegger, 2009).

3 The Models

In this section we describe the setup of the model, starting with a description of the numerical code.

3.1 The numerical scheme.

We use a three-dimensional Cartesian grid-based MHD code developed by Nordlund & Galsgaard (1995) (see also Gudiksen & Nordlund 2005), which has been used in many astrophysical contexts such as star formation, stellar convection, sunspots, and the intergalactic medium (e.g. Padoan et al. (2007); Braithwaite (2012); Collet et al. (2011); Padoan & Nordlund (1999)). It has a staggered mesh, so that different variables are defined at different locations in each grid box, improving the conservation properties. The third-order predictor-corrector time-stepping procedure of Hyman (1979) is used, and interpolations and spatial derivatives are calculated to fifth and sixth order respectively. The high order of the discretization is a bit more expensive per grid point and time step, but the code can be run with fewer grid points than lower order schemes, for the same accuracy. Because of the steep dependence of computing cost on grid spacing (4th power for explicit 3D) this results in greater computing economy. We model the star as a “star in a box”, by modeling a self-gravitating ball of gas inside of the cubic computational domain. The simulations described here, unless stated otherwise, are run at a resolution of , and the star has a radius of grid spacings (see below).

For stability, high-order “hyper-diffusive” terms are employed for thermal, kinetic, and magnetic diffusion. These are an effective way of smoothing structure on small, badly-resolved scales close to the spatial Nyquist frequency whilst preserving structure on larger length scales. This results in a lower “effective” diffusivity compared to standard diffusion, and is more computationally efficient than achieving the same by increasing the resolution. In this study we are interested in physical processes taking place on a dynamic timescale, and so the diffusion coefficients are simply set to the lowest value needed to reliably prevent unpleasant numerical effects (‘zig-zags’), so that the diffusion timescale is as long as possible. All three diffusivities are equal, i.e. the Prandtl and magnetic Prandtl numbers are both set to . The code contains a Poisson solver to calculate the self-gravity.

3.2 Timescales, numerical acceleration scheme

The presence of several very different timescales in stellar MHD problems presents computational difficulties. There is the sound crossing time (governing deviations from pressure equilibrium), the Alfvén crossing time (on which a magnetic field evolves towards equilibrium), the time scale on which the magnetic field evolves under magnetic diffusion, and the rotation period of the star. We simplify the problem by assuming a non-rotating star, and by making sure the magnetic diffusion inherent in the numerical method is small enough that it does not significantly affect the evolution of the field, while still allowing for magnetic reconnection. The timescales to be included are then the sound travel time and the Alfvén time. For realistic Ap star field strengths of about a few kG, is of order a few years, which is 5 orders of magnitude longer than the sound crossing time. To cover a range like this, we make use of the fact that in a star close to hydrostatic equilibrium the evolution of the magnetic field of a given configuration depends on the field amplitude only through a factor in time scale. That is to say, if a field has the initial state , a field with initial amplitude, , where is a constant, evolves approximately as


In other words, the time axis scales in proportion to the Alfvén crossing time. This approximation is valid as long as the Alfvén crossing time is sufficiently long compared with the sound crossing time that the evolution takes place close to hydrostatic equilibrium, and at the same time sufficiently short compared to the magnetic diffusion time scale.

We make use of this scaling to speed up the computation. The field strength of the configuration is multiplied by a time dependent factor , chosen so as to keep the Alfvén crossing time approximately constant between time steps. With Eq. (4) the resulting (unphysical) field is related to the physical field by , where is related to by . For the numerical implementation see the Appendix. This acceleration scheme (Braithwaite & Nordlund 2006) makes it possible to follow the decay of an unstable configuration to very low field strengths, by maintaining an artificially high amplitude magnetic field, and consequently shorter Alfvén time, considerably decreasing computational needs and making such studies feasible.

Since the Alfvén speed is not uniform within the star, it is necessary to decide on some suitable average for the definition of ; we use


where , , , , and are the radius, mass, average density, average magnetic field, and total magnetic energy of the star. Similarly, the sound crossing time is defined as


where is the total thermal energy content of the star and is the adiabatic index, further discussed in Section 3.3.

The initial field strength of the configuration , or equivalently the value of the initial Alfvén time , is an adjustable parameter. Its choice involves a compromise. A lower value of the initial field strength increases the separation between and . The approximation made is then better, but computationally more expensive. The value used here corresponds to .

3.3 Implementing stably stratified vs. barotropic stellar models

The code assumes a chemically uniform, monatomic, classical ideal gas, whose specific entropy is , and


This partial derivative as given corresponds to an adiabatic change as discussed in § 2.

Simulations so far (e. g., Braithwaite & Spruit 2004; Braithwaite & Nordlund 2006; Braithwaite 2009) have taken an initial density profile inside the star corresponding to a polytrope , with


where (here ) is the usual “polytropic index”, and the value was chosen to roughly match the radiative zones of Ap stars. For these values, since , we have and , so the star is stably stratified. We use this model as a reference against which to compare the new barotropic model for which we choose an initial profile with (), so now , , and .



Figure 1: (a)Initial dimensionless specific entropy, with zero point at its central value, as a function of . The red dashed line is for the stably stratified model (), and the solid black line is for the barotropic model (). (b) The magnetic diffusivity versus radius for the stably stratified model.

The profiles of the specific entropies of both models can be seen in Figure 1(a). The star in either model is surrounded by a low-density, poorly conducting atmosphere, which causes the magnetic field outside of the star to relax to a potential field. Figure 1(b) shows the magnetic diffusivity profile for the stably stratified model. The atmosphere has a high temperature, so that its density does not become too small towards the edges of the computational domain.

Note, however, that the numerical code contains heat diffusion. Since the chosen profiles correspond to a radially decreasing temperature, the diffusion will reduce the temperature gradient, increasing the specific entropy gradient and thus making the fluid increasingly stably stratified. In order to counteract this effect and keep the fluid barotropic, we introduce an ad-hoc term in the evolution equation for the internal energy per unit volume , namely,


where is the local temperature, is the initial value of the specific entropy inside the star, and is a timescale (to be chosen) at which this “entropy term” forces the star back to its initially isentropic structure.

4 Disordered initial fields

The setup of the disordered initial field is done in the same way as in Braithwaite & Spruit (2004) and Braithwaite & Nordlund (2006). This consists in assigning random phases to locations in three-dimensional wavenumber space to wavenumbers from that which corresponds to the size of the star up to wavenumbers corresponding to about four grid spacings. The amplitudes are scaled as a function of wavenumber as :


where is the amplitude for wavenumber , and and are two randomly created phases. A reverse Fourier transformation is performed to obtain a scalar field. Two more scalar fields are produced in the same way and these become the three components of a vector potential. The vector potential is then multiplied by , where was set equal to roughly a quarter of the radius of the star. This was done to concentrate the field in the inner quarter of the star, so that the field dies off sufficiently quickly outside the star.

From this potential field, the magnetic field is then computed by taking the curl, and then its magnitude is scaled so as to obtain the desired total magnetic energy, which was times the thermal energy.

We used four different barotropic models, each including different initially disordered field configurations and an entropy timescale value roughly equal to the sound crossing time. We also include one model where a disordered magnetic field was put into a stably stratified model as a comparison. The works of Braithwaite & Spruit (2004) and Braithwaite & Nordlund (2006) have shown that disordered fields in a stably stratified model can reach stable equilibria. The evolution of the total magnetic energy of all of these models can be found in Figure 2. The conclusion is that the stably stratified model reaches a stable equilibrium in a few Alfvén timescales after which, the field decays only on the long diffusive timescale, which in a main-sequence star is of the order yr. On the other hand, none of the barotropic models reach such a stable equilibrium. Figure 3 contains 3-D renderings of the magnetic field for both a stably stratified model and a barotropic model after a time of 70. In the stably stratified model, there is a visible twisted torus structure, while the barotropic model still has a disordered magnetic field configuration. This suggests that stable stratification is a crucial ingredient for reaching a stable MHD equilibrium from a disordered field.

Figure 2: Total magnetic energy relative to initial total magnetic energy versus time, given in initial Alfvén times, for models that have an initially disordered magnetic field configuration. The solid violet curve is the evolution of a stably stratified model that reaches a stable equilibrium, and is shown for comparison. The dashed-dotted gray line shows the evolution of the stably stratified model where the model is kept static, to show how the system will evolve under just diffusive processes. All other curves are the evolution of barotropic models which began with initially different disordered fields, none of which reach stable equilibria.

(a) (b)

Figure 3: Magnetic field configurations after a time of roughly 70 for models with initially disordered fields in a stably stratified star (panel (a)) and a barotropic star (panel (b)). Lines shown in red are outside of the star, those that are blue and green are inside the star. Panel (a) shows that the field has evolved into a large twisted torus inside the star. It was shown that this configuration decays on a diffusive timescale and is stable. Panel (b) shows that the barotropic model still has a disordered magnetic field, and much of the magnetic flux is at large radii or even outside of the star.

We have also investigated the effect of using different values of on the evolution of a disordered initial field. Simulations were carried out with four different values of , as well as a simulation where the entropy term was not used, which is akin to being equal to infinity. Figure 4 shows the evolution of the total magnetic energy for these models, as well as the stably stratified one. It is obvious that increasing the value of has little effect on the evolution, although after a few , models with a larger entropy timescale decay slightly more slowly. The reason for this is that the longer timescale allows for a slight positive entropy gradient to evolve, thus a small buoyant restoring force will act to slow down the rise of the magnetic field. It should be noted, however, that even in the case where the entropy term is not utilized in the code, the formation of a stable stratification does not occur quickly enough for a stable equilibrium to be created.

Figure 4: Evolution of the total magnetic energy relative to its initial value, for models all starting from the same initially disordered magnetic field configuration, and evolved with values of 1, 5, 13, and 40 , where is the entropy timescale defined in Eq. 9, depicted as a short-dashed magenta, short-dashed dotted blue, solid red, and long-dashed dotted green curves, as well as a model in which the entropy term was not used, which is plotted as a long-dashed gold curve. To compare, the evolution of the random magnetic field in a stably stratified model is also shown with the solid violet curve. As the timescale is increased in the barotropic models, the magnetic energy decays more slowly at late times. However, none of the simulations of barotropic models reaches a stable equilibrium.

5 Axisymmetric equilibria

Akgün et al. (2013) wrote down analytic expressions of linked poloidal-toroidal magnetic field configurations that correspond to hydromagnetic equilibria in stably stratified stars. (A particular form of these was also used by Mastrano et al. (2011)) They have not been shown to be stable, and they are not equilibria in the barotropic case, but they might be stable equilibria (or close to these) in the stably stratified case and a good approximation to the best candidates for stability in the barotropic case.

In order to motivate their functional form, consider a general, axisymmetric field written as a sum of poloidal and toroidal components,


where and are characteristic values of the two components, and and are dimensionless scalar functions that depend on , the radial coordinate, and , the polar angle in the star. Since there cannot be an azimuthal component of the Lorentz force (which, in axial symmetry, cannot be balanced by pressure gradients and gravity), and must be parallel, thus one can write .

The latter condition (a special case of the condition of the vanishing vertical component of discussed in Section 2) is the only condition required to be satisfied by and in order to yield a hydromagnetic equilibrium in a stably stratified star. Besides this, these functions are arbitrary and need not satisfy, e. g., the commonly assumed Grad-Shafranov equation, which is physically justified only for barotropic fluids and otherwise represents an arbitrary, additional constraint (e. g., Reisenegger 2009). Of course, in neither case there is a guarantee for the stability, and thus for the astrophysical relevance, of the constructed equilibria.

For the poloidal part, Akgün et al. (2013) choose


where x is the dimensionless radial coordinate, , for the stellar radius, and


where the are constants. Outside the star, this gives an exact dipole. For the internal field, three of the four constants can be obtained from the condition that all components of the magnetic field are continuous across the stellar boundary (implying and , where the primes denote derivatives with respect to ) and that the current density vanishes at the surface (implying ). The fourth constant, here chosen as , is a free parameter, which can be tuned to change the size of the closed field line region of the equilibria, as seen in Figure 5.

Figure 5: Cross sectional view of the poloidal field in the equilibria from Akgün et al. (2013), with the parameter =-100 and -10 (left and right, respectively), controlling the volume of the closed field lines within the star.

The closed field line region is important because the toroidal component of the field, described by the function , can only exist in the volume where the poloidal field lines are closed within the star. This motivates the choice:


where the exponent , so that the current due to the toroidal field decreases smoothly to zero at the stellar surface.

5.1 Stability tests in stably stratified models.

We do not know a priori whether these equilibria will be stable in any type of star. As a first step, we verified their stability in a stably stratified model. If the magnetic field configurations are not stable in the stably stratified model, it seems likely that they will not be stable in a barotropic model.

We utilize the aforementioned equilibria of Akgün et al. (2013) with the toroidal exponent of Eq. 14 set equal to 2, as was done in their analytic work, and from Eq. 13 equal to -100. The reason for choosing this value of was two-fold. The equilibria found to evolve from random field configurations in the stably stratified models were found to have larger tori than the equilibria of Akgün et al. (2013) in which the value of was set to zero. Secondly, we chose to be a large amplitude, so as to be spread out the toroidal flux across a larger area of the star. This would have the effect of decreasing the local magnetic energy density in the torus, as compared to a configuration with a smaller closed field line region with the same fraction of the toroidal to total magnetic energy. In spreading out the toroidal energy, the local value of , defined as , in the torus would be larger, and it would thus be less buoyant.

We ran a number of simulations with various initial poloidal field energies relative to total magnetic energy , ranging between 0.0015 and 0.93. Such a large range was used, because at low values of this parameter, the Tayler instability (Tayler, 1973) is expected to set in, and at high values, instability along the “neutral line” may occur (Markey & Tayler, 1973; Wright, 1973). Models with initial between 0.008 and 0.8 were found to be stable, and to decay on a diffusive timescale, while all others decayed more quickly due to instabilities. Figure 6 shows the evolution of the magnetic energy for some of the models. To try to diagnose the source of instabilities, we analyzed the evolution of the =0-3 azimuthal modes, which were determined by performing a root-mean-square integration of the component of the velocity in the meridional plane as:


where is the amplitude of the theta component of the velocity for a particular -mode, is the area of the star in the meridional plane, and is the component of the velocity. We then calculated the component of the kinetic energy in each of these modes, and plotted these values normalized to the initial total magnetic energy versus time in Figure 7. Here, it is evident that the model with of 0.5 is stable to all of the modes, as the curves for each mode evolve to a steady value. The model with initial of 0.87 was unstable to m=2 which can be seen by the sharp peak in the m=2 curve at a time of roughly 2-5 after which this model calms to a new equilibrium as the m=2 mode flattens out. The model with initial of 0.0025 was unstable to the m=1 mode which begins to set in at a time of roughly 4, as can be seen by the peak in the figure. It was then found that all configurations with initial above 0.8 were unstable to the m=2 mode, while those with initial below 0.008 were unstable to the m=1 mode. In the m=1 instability, the toroidal component of the field can be thought of as a stack of rings on top of one another; as the instability occurs, the rings slip with respect to each other. In the m=2 instability, the region along the neutral line, which is the space where the poloidal field goes to zero, experiences tension, which results in kinking of the neutral line, stretching any toroidal field that may be present along the neutral line. This kinking of the neutral line bends the toroidal field lines into a shape that is similar to the seam of a tennis ball. Both the =1 and =2 unstable configurations end up reaching another stable equilibrium after some time. In the case of the =1 unstable model, the equilibrium is reached at a much later time. However, it can be seen that the =2 unstable model reaches a new “tennis ball” shaped equilibria after roughly 8 when it begins to decay on a diffusive timescale, and the kinetic energy amplitude in the =2 mode drops off and becomes stable in Figure 7. These results are similar to those of Braithwaite (2009), who studied the stability of twisted torus configurations in stably stratified stars.

Figure 6: Evolution of the total magnetic energy given in initial Alfvén times, for an initially ordered magnetic field with equal to 0.0025, 0.50, and 0.87 as the dashed red curve, solid black curve, and dotted blue curve respectively, in a stably stratified star.
Figure 7: component of the kinetic energy in each of the azimuthal m=0-3 modes relative to the total initial magnetic energy versus time for stably stratified models with initially ordered magnetic fields and values of: 0.0025 (top left), 0.5 (top right), and 0.87 (bottom). The dashed black curve represents the m=0 mode, the solid gray curve the m=1 mode, the dotted blue curve is the m=2 mode, and the dashed-dotted red curve is the m=3 mode.

5.2 Stability tests in barotropic models

With the knowledge that these equilibria are stable in a stably stratified star for initial values between about 0.008 and 0.8, we now investigate whether such stable equilibria can exist in barotropic stars.

5.2.1 Existence of stable equilibria?

A series of models have been run with the same equilibria used in barotropic models. For these models, was set to be equal to the sound crossing time of the star. The initial values used were between 0.03 and 0.96. The motivation for using values of above 0.9 comes from the works by Tomimura & Eriguchi (2005); Ciolfi et al. (2009); Lander & Jones (2009); Armaza et al. (2014), where magnetic equilibria have been calculated in axial symmetry for barotropic stars, and none of which obtained 0.9. The rest of the spectrum of values was chosen based on the results of Section 5.1, where for a stably stratified star, the equilibria were stable for initial values of 0.008 to roughly 0.8. In addition, Ciolfi & Rezzolla (2013) have calculated magnetic equilibria in barotropic stars and found equilibria with values as low as roughly 0.11.

Figure 8 shows the evolution of the total magnetic energy for barotropic models, as well as a stably stratified model with 0.5 for comparison. It is evident from the figure that none of the equilibria are stable in the barotropic models. An interesting point, however, is that fields with a initial values between 0.33 and 0.2 decay much more slowly at earlier times than all other barotropic models. The reason for the slower decay of the magnetic field for these models can be seen in Figure 9, where the effective magnetic radius, , defined as:


is plotted. It is evident that for magnetic configurations with an initial 0.33, increases, caused by the torus rising in the star. For values below 0.33, the magnetic tension in the torus results in its contraction, causing a decrease of . However, these models that initially experience a contraction of the torus still show that the value will increase after roughly 6 , as the torus rises out of the star.

Figure 8: Evolution of the total magnetic energy relative to its initial value for a stably stratified star with an initially ordered field with 0.5 as the solid violet curve, and a series of curves for barotropic stars with initial values of: 0.03 (red dashed), 0.14 (black dashed-dotted), 0.26 (blue short-dashed), 0.5 (green short-dashed dotted ), 0.74 (gold dashed two-dotted), and 0.96 (gray solid). The stably stratified star’s magnetic field decays on a diffusive timescale (see Section 5.1), while all of the barotropic stars’ magnetic fields decay faster.
Figure 9: Evolution of , the effective magnetic radius (Eq. 16), for the barotropic models discussed in Figure 8. All models with an initial value less than about 0.33 experience an initial decrease in and then subsequent increase as the torus initially contracts due to tension and then rises outwards. All models with an initial above 0.33 simply experience an increase in as the tori simply rise out of the star. In all models the final state is the same, the magnetic field rises out of the star.

5.2.2 Effect of

We also investigated the effect of varying for a magnetic configuration with an initial value of 0.5 in Figure 10. It is evident that a longer allows the magnetic field to decay more slowly. As the entropy timescale becomes larger, the star, due to numerical diffusion, will evolve an entropy gradient, which will impede the rise of the torus. However, even without the entropy adjustment term, the equilibria are still not stable.

Figure 10: Evolution of the total magnetic energy relative to its initial value for models of barotropic stars with an initial value of 0.5, each with a different value. 1, 5, 13, and 40 are plotted as the dashed red, long-dashed dotted black, dotted blue, and short-dashed dotted green curves respectively. In addition, one model where the entropy forcing term was not included is plotted as a gold dashed double-dotted curve, and the comparison model of a stably stratified star is plotted as the solid violet curve. Regardless of the value of used in the barotropic models, the configuration is always unstable.

5.2.3 Variations of the free parameters in the axially symmetric equilibria.

It was shown in Section 5.2.1 that, regardless of the initial fraction used, all the magnetic configurations were unstable. However, the initial fraction of the total magnetic energy in the poloidal component is just one of the free parameters that can be investigated in the axisymmetric configurations. The effect of varying other free parameters, namely the value of in Eq. 14, which affects the distribution of the toroidal component, the value of in Eq. 13, which controls the size of the closed field line region, and the radius of the “neutral line” were also studied. We found that variations to these parameters had no effect on the final state of the simulations, as all configurations were found to be unstable to the same effect seen in Section 5.2.1, where the torus flowed radially out of the star, subsequently decaying in the atmosphere.

5.3 Instabilities

We have analyzed the instability of the -modes in the barotropic models, with the same method as the stably stratified models (See Section 5.1). Figure 11 shows the evolution of the kinetic energies for each of the m=0-3 modes for models with four different initial values. In all cases, the m=0 mode is the dominant instability, quickly rising to a value of a few hundredths of a percent of the total initial magnetic energy, at which point it remains at this range. Even in the cases where and , which in the stably stratified star were found to be unstable to the m=2 and m=1 modes respectively, the m=0 mode is the dominant instability in the barotropic models. To help visualize what the m=0 mode looks like, Figure 12 shows that the cause of the m=0 mode is the rise of the torus, driven by a combination of its own buoyancy and the pressure of the enclosed poloidal field. Thus, the stable stratification seems to play a very strong role in suppressing the m=0 mode, as it always dominates in the absence of stable stratification.

Figure 11: Evolution of the component of the kinetic energy in each of the azimuthal m=0-3 modes relative to the total initial magnetic energy, for barotropic stars containing initial magnetic field configurations with values equal to: 0.26 (top left), 0.5 (top right), 0.76 (bottom left), and 0.96 (bottom right). In all cases the m=0 mode is the dominant mode of instability.
Figure 12: Meridional projection of the poloidal field lines and toroidal contours (drawn as , where is the cylindrical radius, and the color scale is constant for all snapshots), at times 0.28, 2.07, 3.08, and 9.71 in a barotropic model. At time 0.28 the configuration has not changed much from the initial state. By time 2.07 the torus has drifted towards the radius of the star. After 3.08 the toroidal field has risen into the atmosphere where it decays. By time 9.71 nearly all of the torus has decayed.

In order to see how robustly this mode dominates in de-stabilizing the magnetic equilibria, we have investigated what happens when an initial non-axially symmetric perturbation is added to the system. To do so, we introduced perturbations to the component of the velocity along the equatorial plane within the star. One case included an initial =1 perturbation for a model with initial =0.03, a second model contained an =2 perturbation with initial =0.74, while the third model contained perturbations to =1, 2, and 3 modes for an initial =0.5. The strength of all perturbations was the same, with each mode containing an initial =.01%, roughly of the same order that the =0 mode was seen to reach in the non-perturbed models. We again found that regardless of the initial perturbations added, the =0 mode always became the dominant mode of instability.

In addition to non-axially symmetric perturbations, we have run two simulations with a non-equatorially symmetric perturbation. The perturbation was created by offsetting the magnetic field configuration by 0.016 in the Northern direction, for two initial values of 0.5 and 0.24. Snapshots of the time evolution of the poloidal and toroidal field lines can be seen in Figure 13. The evolution in these northerly perturbed models is similar to their non-perturbed counterparts, with the only difference being that these models have a slight tendency to evolve away from the equator.

Figure 13: Snapshots of the poloidal field (black lines) and toroidal field (drawn as ) (blue color scale), for magnetic field configurations with an initial of 0.5 (left panels), and 0.24 (right panels), which were perturbed axially symmetrically, by moving the configuration .016 above the equator, causing a North-South asymmetry in barotropic models. The snapshots are at time: 0., 1.1, 2.0, 3.3, 5.1, and 7.1 , from top to bottom. By time 1.1 , it is possible to see that the toroidal field in the 0.5 model has risen to slightly larger radii, while the 0.24 model’s torus has contracted to smaller radii. By time 2.0 , the tori in both models begin to move towards the northern pole, whilst continuing their previously seen radial movement. By time 3.3 , the tori in each model have moved noticeably towards higher latitudes, and the magnetic field in the 0.24 model begins to rise to larger radii. By 5.1 , the torus of the 0.24 model has continued its northerly path, while the torus in the 0.5 model has begun to decay in the atmosphere. (Note there is some toroidal component outside the star, this is because the enhanced diffusion utilized in the atmosphere takes a bit of time to make the field relax to a potential field.) By the final time step of 7.1 , the torus of the 0.5 model has decayed significantly, while the torus of the 0.24 model has started to reach the atmosphere where it is now experiencing a fast decay of the torus. In both cases, the magnetic field configurations remain axisymmetric throughout this process.

6 Conclusions

We have conducted two different kinds of numerical experiments.

In the first kind, we started with disordered initial fields. In the previously studied case of stably stratified stellar models, we confirmed that the magnetic field evolved into an ordered, stable, roughly axially symmetric configuration with comparable poloidal and toroidal fields. In the barotropic cases we studied, this never happened, and the magnetic energy decayed much faster than expected from diffusion.

In the second kind, we started with a smooth, axially symmetric magnetic field with poloidal and toroidal components. We confirmed that, in the stably stratified star, there are stable hydromagnetic equilibria for a fairly wide range of values of the initial fraction of poloidal to total magnetic energy, . Outside this range, the field decayed through non-axisymmetric modes, in the toroidally dominated and in the poloidally dominated case. In contrast, in the barotropic case, we found no stable configurations (exploring a fairly large parameter space), and the field always decayed through an axially symmetric () instability in which the toroidal flux rose radially out of the star, being dissipated in the atmosphere.

Although far from a rigorous mathematical proof, this provides strong evidence (added to that previously provided by Lander & Jones 2012) that there are no stable equilibria in barotropic stars. It also strongly supports the intuition that, in a stably stratified star, the buoyancy force strongly constrains the potential instabilities by impeding any substantial radial motions, whereas such motions are not hindered in the barotropic case, and in fact these motions destabilize the magnetic equilibrium.

As argued above, the stabilization provided to the magnetic field by the buoyancy force can be roughly quantified by comparing the Brunt-Väisälä frequency to the Alfvén frequency . It will only be effective if , whereas in the opposite case the star would behave as if it were barotropic, not being able to contain a stable magnetic field. If there are dissipative processes effectively eroding the stable stratification on long time scales, they will lead to a destabilization and eventual decay of the magnetic field, unless it can be stabilized, e. g., by the solid crust of a neutron star. In fact, in the neutron star case, this effect has been argued to act on astrophysically relevant time scales, which does not seem to be the case in white dwarfs and upper main sequence stars (Reisenegger, 2009), with the possible exception of very massive O stars, which are radiation-dominated and thus only weakly stratified.

Note that the numerical models used here, which assume a chemically uniform, classical ideal gas, stably stratified by an entropy gradient, clearly do not directly apply to degenerate stars, such as white dwarfs and neutron stars. However, we have no reason to doubt that the essential physics, in particular the competition between magnetic and buoyancy forces, is the same in these cases, and the (very rough) condition is still required to stabilize a hydromagnetic equilibrium in these cases, even in the neutron star case, where is due to a composition gradient. Additional changes in the neutron star case are the presence of a stabilizing solid crust, as well as superconducting and superfluid regions in the interior, which modify the form of the Lorentz force and the dynamical equations and thus do not allow a direct application of the present results, although some of their features might carry over to this regime.

7 Acknowledgements

This work was supported by the DFG-CONICYT International Collaboration Grant DFG-06, FONDECYT Regular Project 1110213, and Proyecto de Financiamiento Basal PFB-06/2007. Some of the figures were made with VAPOR (

Appendix A Numerical acceleration scheme

The scaling of magnetic fields that evolve from different initial amplitudes (but otherwise identical configuration) with time is described in Section 3.2. The practical implementation in the MHD code is as follows. We need to make a distinction between the field and time in the accelerated code, and the physical field and time reconstructed from it. The MHD induction equation


is first evolved over a time step , where , to yield an intermediate update . The magnetic energy over the numerical volume is measured from the result of this time step:


and the evolution time scale (in the numerical time unit) is calculated:


The intermediate update is modified:


which brings the numerical Alfvén crossing time back to its initial value . The new value of the physical magnetic energy has changed by the amount


and the physical time coordinate by the amount



  • Akgün et al. (2013) Akgün T., Reisenegger A., Mastrano A., Marchant P., 2013, MNRAS, 433, 2445
  • Akgün & Wasserman (2008) Akgün T., Wasserman I., 2008, MNRAS, 383, 1551
  • Armaza et al. (2014) Armaza C., Reisenegger A., Valdivia J. A., Marchant P., 2014, in IAU Symposium Vol. 302. pp 419–422
  • Braithwaite (2009) Braithwaite J., 2009, MNRAS, 397, 763
  • Braithwaite (2012) Braithwaite J., 2012, MNRAS, 422, 619
  • Braithwaite & Nordlund (2006) Braithwaite J., Nordlund Å., 2006, A&A, 450, 1077
  • Braithwaite & Spruit (2004) Braithwaite J., Spruit H. C., 2004, Nature, 431, 819
  • Ciolfi et al. (2009) Ciolfi R., Ferrari V., Gualtieri L., Pons J. A., 2009, MNRAS, 397, 913
  • Ciolfi & Rezzolla (2013) Ciolfi R., Rezzolla L., 2013, MNRAS, 435, L43
  • Collet et al. (2011) Collet R., Hayek W., Asplund M., Nordlund Å., Trampedach R., Gudiksen B., 2011, A&A, 528, A32
  • Duez & Mathis (2010) Duez V., Mathis S., 2010, A&A, 517, A58
  • Grad & Rubin (1958) Grad H., Rubin H., 1958, Int. Atom. Ener. Ag. Conf. Proc., 31, 190
  • Haskell et al. (2008) Haskell B., Samuelsson L., Glampedakis K., Andersson N., 2008, MNRAS, 385, 531
  • Hyman (1979) Hyman J. M., 1979, in Advances in Computer Methods for Partial Differential Equations - III A method of lines approach to the numerical solution of conservation laws. pp 313–321
  • Kiuchi & Kotake (2008) Kiuchi K., Kotake K., 2008, MNRAS, 385, 1327
  • Lander & Jones (2009) Lander S. K., Jones D. I., 2009, MNRAS, 395, 2162
  • Lander & Jones (2012) Lander S. K., Jones D. I., 2012, MNRAS, 424, 482
  • Markey & Tayler (1973) Markey P., Tayler R. J., 1973, MNRAS, 163, 77
  • Mastrano et al. (2011) Mastrano A., Melatos A., Reisenegger A., Akgün T., 2011, MNRAS, 417, 2288
  • Nordlund & Galsgaard (1995) Nordlund Å., Galsgaard K., 1995,
  • Padoan & Nordlund (1999) Padoan P., Nordlund Å., 1999, ApJ, 526, 279
  • Padoan et al. (2007) Padoan P., Nordlund Å., Kritsuk A. G., Norman M. L., Li P. S., 2007, ApJ, 661, 972
  • Pili et al. (2014) Pili A. G., Bucciantini N., Del Zanna L., 2014, MNRAS, 439, 3541
  • Reisenegger (2009) Reisenegger A., 2009, A&A, 499, 557
  • Shafranov (1966) Shafranov V. D., 1966, Reviews of Plasma Physics, 2, 103
  • Tayler (1973) Tayler R. J., 1973, MNRAS, 161, 365
  • Tomimura & Eriguchi (2005) Tomimura Y., Eriguchi Y., 2005, MNRAS, 359, 1117
  • Wright (1973) Wright G. A. E., 1973, MNRAS, 162, 339
  • Yoshida et al. (2006) Yoshida S., Yoshida S., Eriguchi Y., 2006, ApJ, 651, 462
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description