3D simulations of rising magnetic flux tubes in a compressible rotating interior: The effect of magnetic tension
LeibnizInstitut für Astrophysik Potsdam (AIP), An der Sternwarte 16, 14482 Potsdam, Germany
email: y.fournier@aip.de
Received: 1 November 2016
Accepted: 4 July 2017
Context. Longterm variability in solar cycles represents a challenging constraint for theoretical models. Meanfield BabcockLeighton dynamos that consider noninstantaneous rising flux tubes have been shown to exhibit longterm variability in their magnetic cycle. However a relation that parameterizes the risetime of nonaxisymmetric magnetic flux tubes in terms of stellar parameters is still missing.
Aims. We aim to find a general parameterization of the risetime of magnetic flux tubes for solarlike stars.
Methods. By considering the influence of magnetic tension on the rise of nonaxisymmetric flux tubes, we predict the existence of a control parameter referred as Γ_{α1}^{α2}. This parameter is a measure of the balance between rotational effects and magnetic effects (buoyancy and tension) acting on the magnetic flux tube. We carry out two series of numerical experiments (one for axisymmetric rise and one for nonaxisymmetric rise) and demonstrate that Γ_{α1}^{α2} indeed controls the risetime of magnetic flux tubes.
Results. We find that the risetime follows a power law of Γ_{α1}^{α2} with an exponent that depends on the azimuthal wavenumber of the magnetic flux loop.
Conclusions. Compressibility does not impact the rise of magnetic flux tubes, while nonaxisymmetry does. In the case of nonaxisymmetric rise, the tension force modifies the force balance acting on the magnetic flux tube. We identified the three independent parameters required to predict the risetime of magnetic flux tubes, that is, the stellar rotation rate, the magnetic flux density of the flux tube, and its azimuthal wavenumber. We combined these into one single relation that is valid for any solarlike star. We suggest using this generalized relation to constrain the risetime of magnetic flux tubes in BabcockLeighton dynamo models.
Key words: stars: magnetic field / magnetohydrodynamics (MHD) / magnetic fields / Sun: magnetic fields
© ESO, 2017
1. Introduction
Solar and stellar dynamo theory is a complex field that has seen a lot of progress in recent decades (e.g., review by Charbonneau 2014). It has also benefited from comparisons with diskresolved stellar spot data (e.g. review by Strassmeier 2009), but still lacks a universal solution for all stars possessing convective envelopes. In many stars, differential rotation is likely responsible for the conversion of poloidal into toroidal fields. However the regeneration of poloidal magnetic field is a more subtle issue (Rädler et al. 2003; Miesch 2005; Charbonneau 2010).
At least two types of dynamo model offer an explanation. The turbulent dynamo theory suggests that there is a net effect from the turbulent electromotive force on large scales of which one part – the αeffect – provides the necessary regeneration (Brun et al. 2004; Käpylä et al. 2010, 2012; Brown et al. 2010; Racine et al. 2011). The BabcockLeighton dynamo model (BLdynamo) on the contrary, describes a regeneration mechanism taking place near the surface where magnetic flux emergence plays a major role (Babcock 1962; Leighton 1969). In many of those dynamo solutions, strong toroidal magnetic fields reside near the tachocline. We suppose that these fields form magnetic flux tubes. The concentration of magnetic flux is a special issue and is beyond the scope of this paper. We focus on the properties of the rise of such magnetic structures from the tachocline to the surface. Hence, the present work addresses the concept of a BLdynamo.
BLdynamos require active regions to regenerate the poloidal field. The formation of active regions is also a debated issue (Cheung & Isobe 2014). We recognize two main concepts here: the magnetic flux tube model and the local formation model. The latter received an interesting incentive through the Negative Effective Magnetic Pressure Instability (NEMPI) which could be responsible for the formation of active regions directly at the surface (Warnecke et al. 2013, and references therein). Within the same concept, we should also underline the local convective model (Rempel & Cheung 2014) which suggests the formation of active regions due to local convective motion and granulation. In contrast to the local formation, the magnetic flux tube models all require coherent magnetic structures preceding the emergence. The present document discusses this scenario.
Some mechanisms have been suggested to form magnetic flux tubes. They could be the manifestation of the concentration of magnetic flux due to turbulence in the bulk of the convection zone (Nelson et al. 2014). Alternatively, they could form at the tachocline due to an entropic instability in a magnetic layer sitting at the bottom of the convection zone (Cattaneo & Hughes 1987; Cattaneo et al. 1989; Matthews et al. 1995; Hughes et al. 1997; Schüssler & Rempel 2002). We suppose that flux tubes form as a result of the destabilization of a magnetic layer in an unstable equilibrium against some kind of RayleighTaylor instability (Wissink et al. 2000; Fan 2001).
In the present paper, we follow and discuss the idea of Schüssler (1980) and assume that stellar cycles as well as the dynamo can be explained by rising flux tubes formed in the tachocline. We study rising flux tubes in compressible rotating stellar interiors with highly resolved, threedimensional numerical simulations. The issue of rising flux tubes was first addressed by Parker (1955) who described the buoyant instability and the fact that flux tubes could rise as coherent structures through the turbulent convection zone. This idea became really attractive, however, after Spruit (1981) derived the equations for the thinfluxtube approximation. The approximation turned out to be a very useful description for numerical experiments with a high level of precision. This idea became enriched by an increasing amount of physics and a complex model was presented recently by Weber & Fan (2015) which incorporates most of the relevant physics into the thinfluxtube approximation. This latter publication concluded a long series of work which addressed various aspects of the rise of magnetic flux tubes; the influence of initial conditions (Spruit & van Ballegooijen 1982; MorenoInsertis 1983; Yoshimura 1985; MorenoInsertis et al. 1992; Fan et al. 1994), and the need of rotation to reproduce the observations (Choudhuri & Gilman 1987; van Ballegooijen 1983; Schüssler & Solanki 1992; Caligari et al. 1994, 1996; DeLuca et al. 1997; Granzer et al. 2000; Granzer 2004).
In particular, Schüssler & Solanki (1992) showed that in order to reproduce large polar spots on shortperiod stars, the model needed to take rotation into account. They also showed that several properties of the rise of a flux tube scale with the magnetic Rossby number, (1)where is the Alfvén speed of a magnetic field B in a medium with density ρ, while μ_{0} is the magnetic permeability, H_{P} is the pressure scale height, and Ω is the angular velocity. This is an important result if one wants to learn from other stars, as it defines the regime of the rise. One of the aims of the present work is to discuss this result and to try to generalize it to the nonaxisymmetric case.
Finally, it has been kept in mind that the thinfluxtube approximation has limits (Cheung et al. 2006). This is why Boussinesq as well as anelastic rising flux tubes (“thick flux tubes”) were simulated in axisymmetry (twodimensional (2D); MorenoInsertis 1983, 1997; Choudhuri & Gilman 1987; Chou & Fisher 1989; Fan et al. 1998). Getting rid of the thin aspect of the flux tube introduced the need for twist to maintain the coherence of the magnetic structures along their rise (Browning & Priest 1983; MorenoInsertis & Emonet 1996; Longcope & Klapper 1997; MorenoInsertis 1997).
Later, Fan (2008) showed that the azimuthal asymmetry of the rise of a magnetic flux tube influences its dynamics. A variety of papers have been published on nonaxisymmetric (threedimensional (3D)) studies (Jouve & Brun 2009; Weber et al. 2011; Fan et al. 2013; Pinto & Brun 2013; Weber & Fan 2015).
We present the first fully compressible nonaxisymmetric model of rising flux tubes in a rotating stellar interior. We conducted parameter studies which enable us to extend the control parameter (1) for the axisymmetric case to the 3D case. The original idea consists in challenging the flux tube theory, and incorporating our results into a meanfield dynamo model, which could be validated by observations and could help to constrain the flux tube dynamo theory. In Sect. 2 we first describe the numerical experiments. Then in Sect. 3 we derive the scaling parameter from first principles and make a prediction for the nonaxisymmetric rise. Before studying the nonaxisymmetric case we first validate the 2D setup in Sect. 4 by comparing our results with the literature. In Sect. 5, we show that our theoretical predictions for nonaxisymmetric rises are valid. Finally, we discuss the limits of our results and possible improvements in Sect. 6. In the Appendix, we enclose the derivation of the stratified interior and the full list of our simulations.
2. Equations and setup
Our study relies on 2D and 3D numerical simulations. Since most of the presented simulations differ by just one parameter, we often refer to some standard cases. For each simulation we specify the modified parameter presupposing that all others are taken from the standard model. From the end of this section we refer to a fiducial 2D simulation (STD2D), and a standard 3D simulation (STD3D).
2.1. Compressible MHD set of equations
We are interested in rising magnetic flux tubes in stellar interiors. While in the present paper, we focus on solarlike stars, our model is constructed such that we could cover more than this specific case and extend our study to other dwarfs and red giants in the future. We solve the fully compressible, resistive MHD equations. One could argue that in solarlike stars the anelastic approximation is sufficient and saves a lot of computational time. But in case of red giants this assumption may not hold anymore.
The MHD equations as we solve them can be written in the following compact form: (2)(3)(4)(5)(6)Here, u, P_{tot}, B, g, f, T, F_{rad}, k_{B}, m and μ are the velocity field, the total pressure being the sum of the thermal (P) and the magnetic pressure (P_{m} = B^{2}/ 2μ_{0}), magnetic flux density, gravitational acceleration, external acceleration (namely Coriolis), the temperature, the radiative flux, the Boltzmann constant, the atomic mass unit and the mean molecular weight, respectively. Apart from the usual symbols, ∂_{t}(.) and I are the partial time derivative and the identity matrix, respectively.
We solve the equations in a fraction of a spherical shell with varying azimuthal extent depending on the needs of various nonaxisymmetries (Fig. 1). This system is solved on a spherical grid with the parallelized NIRVANA code, described in detail by Ziegler (2011)^{1}. The spherical coordinate system is (r,θ,φ).
Fig. 1 Sketch of the numerical setup. The yellow contour shows the buoyant part of the flux tube with a lack of entropy. This section of the flux tube will buoyantly rise toward the surface and emerge as an active region. The orange contour surfaces show the neutrally buoyant parts of the flux tube, with an excess of entropy. The latter sections will remain rooted at the bottom of the simulated domain shown by the wedgelike shape of solid lines. The azimuthal morphology shown here is valid only in 3D, for an azimuthal wavenumber of m = 8. In 2D, the flux tube is buoyant everywhere, and the computational domain is a meridional plane. The direction of rotation is also indicated. 

Open with DEXTER 
2.2. Rotating adiabaticallystratified stellar interior
We study the rise of magnetic flux tubes in a hydrostatic, adiabatically stratified layer. The design and the analysis of a convective zone is by itself a complex problem and should be addressed separately. We choose to first tackle the issue of a nonconvective environment. In this situation, the radiative flux transports the entire luminosity radially along the unit vector e_{r}: (7)where κ is the thermal conductivity. This layer is also hydrostatic, the gradient of pressure balances gravity (here approximated by a point mass), (8)where M_{⋆} is the stellar mass. We defined the pressure scale height at the top of the domain as follows: (9)where P_{0} and R_{0} are the pressure and the radius at the top of the domain, respectively. The logarithmic temperature gradient is (10)where the special case of ∇ = 0.4, referring to an adiabatic stratification, is denoted by ∇_{ad}. Making use of these four equations one can obtain the analytical expression for the three thermodynamic quantities T, ρ, and P.
Fig. 2 Radial profiles of ρ, T, and P, normalized by their respective values at the top of the domain. The radial profile of the normalized gradient of entropy is a constant (adiabatic) in this configuration. 

Open with DEXTER 
The index 0 refers to the quantities at the top boundary of the domain. We illustrate their profiles in Fig. 2 where we also show the constant entropy gradient profile characterizing adiabatic layers. For more details, please refer to the Appendix. Furthermore, the entropy gradient is defined as (14)where c_{v} and c_{p} are the volumetric heat capacity and the pressure heat capacity, respectively.
A proper treatment of the top boundary condition requires at least a few tens of grid points to resolve the pressure scale height at the top of the domain (0.964 R_{⋆}). In case the pressure scale height is not sufficiently resolved, the hydrostatic equilibrium cannot be maintained. Therefore stratification imposes a minimum radial resolution necessary for the simulation.
2.3. Boundary conditions
Both models, STD2D and STD3D, have similar boundary conditions. As illustrated in Fig. 1, boundary conditions for the thermodynamical variables are a constant thermal flux at the bottom of the convection zone, and a constant temperature T_{0} at the top of the domain. In the latitudinal direction, we choose zero thermal flux condition at the equator and at high latitude. For the momentum equation, all boundaries are stressfree. For the STD3D setup we apply periodic boundary conditions in the azimuthal direction. For the magnetic field, we use a pseudoaxis boundary condition for the highlatitude boundary, a reverse boundary at the equator, and a pseudovacuum at the inner and outer radial boundaries, that is, a verticalfield condition.
2.4. Useful parameters
Throughout this study we make use of several parameters. The plasmaβ is the ratio of gas pressure to magnetic pressure, (15)The rotational Mach number is used to describe rotation, defined as the ratio of the rotation speed to the sound speed c_{s}, (16)where ϖ_{ini} is the axis distance of the initial tube location, (r_{ini},θ_{ini}), that is, ϖ_{ini} = r_{ini}sinθ_{ini}. We also introduce the Lorenz number (Lo) according to (17)We further introduce the relative risetime defined as (18)where τ_{rise} and P_{rot} are the time needed by the flux tube to reach the emergence line (0.95 R_{⋆}) of the computational domain, and the rotation period, respectively.
Numerical experiments can take advantage of dimensionless variables. We present the dimensionless system used in our setup in Table 1.
Dimensionless definition of main quantities.
The dimensionless pressure scale height at the top of the domain is defined by (19)
2.5. Criterion for adaptive mesh refinement (AMR)
Activeregion latitudes span from 0° to 40°. Catching such a wide latitudinal band requires a considerable fraction of the spherical domain. Therefore we need to investigate rising magnetic flux tubes in global simulations. However, flux tubes are small coherent structures in this large domain. Because of the strong stratification of the solar interior and supposing that magnetic flux tubes conserve their internal flux all along their rise, a 60 Mm active region can only result from the emergence of a magnetic flux tube that had a radius less than 1% of the convective zone’s radial extent, when it formed at the bottom of the convective zone. In addition, the numerical conservation of the magnetic flux in the tube forces us to resolve the magnetic flux tube with at least 50 points in diameter.
Adaptive mesh refinement (AMR) allows us to locally add resolution, and therefore meets the requirements of our situation by resolving small structures in large domains (Ziegler 2012). The AMR procedure checks, at each time step, the refinement (derefinement) condition: if a cell fulfills a given criterion, a resolution level is automatically added (removed). In 2D, a cell will be divided into four smaller cells; in 3D into eight cells. We chose the presence of a magnetic field as a refinement condition: if the strength of the magnetic field exceeds 10% of the maximum value in the domain, the grid will be refined. We limit the code to two refinement levels on top of the base level to avoid overrefinement and save computation time.
2.6. Flux tube definition
At the bottom of the prescribed static layer, we introduce a flux tube in nonequilibrium. The flux tube is a twisted torus of constant magnetic field along the azimuthal direction. The magnetic field strength decreases with the distance from the tube’s center at (r_{ini},θ_{ini}). The strength of the azimuthal magnetic field is assumed to be (20)where r_{ft}, R_{ft}, and B_{0} are the distance from (r_{ini},θ_{ini}), the initial radius of the flux tube, and the maximum strength of the magnetic flux tube, respectively. The radius of the flux tube (R_{ft}) defines the torus, which contains about 98% of the magnetic flux initially. B_{0} is the magnetic field strength in the middle of the flux tube, and corresponds to the value used in thinfluxtube approximation models.
While a purely azimuthal magnetic field is by construction divergence free, ensuring the solenoidality of a twisted flux tube is less trivial. We therefore derive the field twist from a vector potential. This component is circular around the tube center, with strength (21)where λ is the twist parameter and A_{φ} is the azimuthal component of the magnetic vector potential. Figure 3 shows the resulting radial profiles of B_{p} and B_{φ}. In the present paper, we fixed all characteristics of the flux tube except B_{0}. Its value is computed from the input parameters.
Fig. 3 Profile of the dimensionless magnetic quantities over the normalized distance from the flux tube center. B_{φ}, B_{p}, and P_{m} are the toroidal component of the twisted magnetic field, the poloidal component of the twisted magnetic field, and the resulting magnetic pressure, respectively. 

Open with DEXTER 
We assume that flux tubes form at the tachocline. Therefore, we set the initial depth r_{ini} = 0.71 R_{⋆}. Furthermore, we chose θ_{ini} = 20° because the thinfluxtube theory demonstrated that such simulated axisymmetric flux tubes emerge at the observed active latitudes.
We set the initial size of the flux tube such that it meets two requirements: the magnetic flux should correspond to the flux observed in large emergence regions and the diameter of the flux tube should contain more than 50 grid points. Such a resolution is sufficient to ensure the numerical diffusion to be negligible along the flux tube’s rise; and the flux tube will conserve about 90% of its initial flux. The resolution we have chosen represents a compromise between realistic physics and realistic computation time (see Table 2).
In the simulations presented here, magnetic flux tubes rise sufficiently slowly to consider them in pressure equilibrium; the presence of magnetic pressure lowers the thermal pressure inside the radius of the flux tube compared to the thermal pressure of its surroundings.
The thermal state of the flux tube controls how the pressure lack determines the other thermodynamical quantities, ρ and T. The initial thermal state of a flux tube has been shown to have great impact on the dynamics of the rise (MorenoInsertis 1983). We illustrate two extreme situations here: (a) an isothermal flux tube, and (b) a neutrally buoyant flux tube.
In order to discuss these two different thermal states and for more convenience, from now on, we use the indices ()_{i} and ()_{e} for quantities inside the flux tube at r_{ft} = 0 and outside the flux tube at r_{ft}> 2R_{ft}, respectively. As shown in Fig. 4, the isothermal case is accompanied by a lack of density. In such a case, the flux tube is buoyant. The neutrally buoyant case, however, consists of a cool flux tube. The lack of thermal pressure exclusively applies to temperature, allowing the densities inside and outside the flux tube to be equal. We note that even if such a flux tube is buoyantly neutral, it will still rise due to conduction of heat inside the tube. However such a flux tube rises on a much longer time scale.
Fig. 4 Profiles of the magnetic pressure and the thermodynamical quantities for a magnetic flux tube in two extreme hydrostatic states, a in an isothermal state and b in a neutrally buoyant state. 

Open with DEXTER 
2.6.1. The thermal state of the axisymmetric case: STD2D
In the axisymmetric case (azimuthal wavenumber m = 0) the flux tube remains a torus during the totality of its ascent. The flux tube is assumed initially isothermal; hence, it is buoyant everywhere. We denote the internal density by , where ()^{(m)} indicates the wavenumber of the most unstable mode. The internal density can be simply written as (22)We note that the internal density varies with the location of the flux tube, since ρ_{e} comes from the stellar model and is a function of radius r.
Fixed parameters for the numerical setup and the initial conditions of the STD2D case.
2.6.2. The thermal state of the nonaxisymmetric case: STD3D
In the nonaxisymmetric case, the flux tube is initially an axisymmetric torus. While rising, it evolves into an asymmetric Ωshaped loop exhibiting some writhe. Since we envisage that flux tubes rise in the form of loops, we introduce a periodic perturbation along the azimuthal coordinate with an amplitude of ρ_{e} − ρ_{iso}, so that an aziuthal section of the tube is buoyant (see the yellow contour in Fig. 1). The apex of the loop rises, while the feet of the loop remain almost neutrally buoyant and stay around r_{ini}. It is important to note that we do not impose any artificial anchoring on the feet of the loop. The resulting rising Ωshaped loop may writhe and tilt.
In contrast to the axisymmetric flux tube, we break the symmetry by an entropic wave similar to the one described by Jouve & Brun (2009). The dependence of the internal density on the azimuthal direction is denoted by and is defined by (23)For m = 0, the definition is equivalent to the STD2D thermal state.
We are simulating a spherical wedge with an azimuthal extent of π/ 4 which is oneeighth of the full azimuthal range. In order to obtain a single apex in the computed domain, we therefore choose to give a wave number m = 8 to our entropic wave.
The resolution, stratification, and flux tube radius are set to account for the same constraints as in the 2D case. However since 3D simulations are much more time consuming (320 times longer at the same resolution), we had to reduce the resolution by four, leading to lower stratification (≈11) and a ten times larger flux tube radius. We think that even at this large increase in radius and at this large decrease in resolution, the drag force is kept sufficiently low so that we are still in a kinematic regime, the viscous force does not dominate. Furthermore, AMR enables us to keep reasonable computing times, retaining again about 90% of the initial magnetic flux.
We summarize the standard parameter values of the STD3D model in Table 3.
Fixed parameters for the numerical setup and the initial conditions of the STD3D case.
3. The effect of local magnetic tension in nonaxisymmetric rise
We aim to find a parameter that controls the risetime of magnetic flux tubes in rotating stellar interiors from compressible numerical experiments. In contrast to anelastic or thinfluxtube simulations, compressible experiments suffer from an upper limit for the β parameter. Compressible codes are not making use of a background state to solve the MHD equations, hence Δρ/ρ has to be larger than the discretization error, whereas anelastic simulations ensure this by construction. For the setup’s resolution, β is limited to β_{max} ≈ 200 for STD2D and β_{max} ≈ 50 for STD3D. In the solar case, however, β is expected to be of the order of 10^{5}. Hence, it is crucial to be able to define the regime of the buoyant rise in order to scale the results to higher β and compare our results with other simulations and observations.
Choudhuri & Gilman (1987) underlined that the ratio of the buoyant force to the Coriolis effect controls the regime of an axisymmetric buoyant rise, (24)where Δρ and v_{rise} are the lack of density inside the flux tube and the rise velocity of the tube, respectively. When this ratio exceeds unity, the regime of the rise is buoyancy dominated; when the ratio becomes less than unity the regime is rotation dominated.
As Schüssler & Solanki (1992) pointed out for axisymmetric simulations, the rise velocity, v_{rise}, corresponds to the buoyant velocity: (25)with Δρ/ρ = 1 /β according to Eq. (22) and l = H_{P} being the local pressure scale height. Replacing v_{rise} by v_{buoy} in Eq. (24), we can rewrite the ratio as (26)where is the free fall velocity and ℳ_{ff} = v_{ff}/c_{s} is the free fall Mach number. In the present series the latter remains constant and is approximately unity. The variables in Eq. (26) are v_{A} and ϖΩ. Therefore, the ratio is proportional to the Lorentz number. Two simulations with different parameters but the same Lorentz number will deliver the same solution; this conclusion being valid exclusively for axisymmetric simulations.
In a nonaxisymmetric rise, magnetic tension will alter the ratio of forces. Due to the high β in stellar interiors, magnetic tension never dominates, but because it is directed inwards it acts against the buoyant rise. The magnetic tension reduces the rise velocity which alters the Coriolis effect.
Magnetic tension can be approximated by (27)where ℛ is the curvature radius of the rising flux tube.
In the axisymmetric case the curvature radius, Sℛ, is the distance of the flux tube from the rotation axis (ϖ) and it is sufficiently large to make the tension force negligible. The regime depends on two independent variables, β and ℳ_{rot}, which control the buoyant force and the Coriolis effect, respectively. In the nonaxisymmetric case, flux tubes rise in the form of Ωloops. The curvature radius, ℛ, of the Ωloops, naturally connects to the azimuthal wavenumber, m. Nonaxisymmetry introduces a new degree of freedom, that requires an additional parameter: m.
As shown in Fig. 5, ℛ does not only depend on the azimuthal wavenumber m (4 in that case), but also on the nature of the regime of the rise. Flux tubes that rise in a rotationdominated regime (red line) have a smaller curvature radius than flux tubes rising in a buoyancy dominated regime (blue line).
It becomes clear that ℛ depends on β, ℳ_{rot}, and m, and will here be modeled with the ansatz (28)where f_{1} and f_{2} are two functions of the azimuthal wavenumber m. We can already constrain these functions. For instance, increasing β or ℳ_{rot} pushes the system toward the rotationdominated regime (toward the red line in Fig. 5), and decreases ℛ. Hence f_{1} and f_{2} both have to be negative. A further constraint concerns the axisymmetric case, where the curvature radius is a constant and equals ϖ. In that specific case, ℛ depends neither on β nor on ℳ_{rot}, so f_{1} and f_{2} are both zero.
Fig. 5 Sketch representing the approximate shape of rising flux tubes depending on the regime of their rise. The light blue flux tube is in a buoyancy dominated regime. The dark blue flux tube is in a transitional regime, where F_{buoy} ≈ F_{corio}. The red flux tube is in a rotation dominated regime. 

Open with DEXTER 
In the axisymmetric case, ℛ = ϖ and the rise velocity v_{rise} is the buoyant velocity v_{buoy}. In the nonaxisymmetric case, the rise velocity is reduced by the tension force. We make the hypothesis that the rise velocity is a fraction of the buoyant velocity where the fraction is controlled by the curvature radius, (29)where v_{buoy} is the rise velocity of a flux tube buoyantly rising in an axisymmetric manner, and ℛ is the curvature radius of the flux tube at the apex. We introduce a factor k to consider the various effects of the drag force acting on the flux tube, thermal conduction and the twist. The only constraint on k is that it neither depends on β nor on ℳ_{rot}.
The Coriolis effect depends on the rise velocity; the reduction of v_{rise} by the magnetic tension naturally alters the Coriolis effect: (30)where is the altered Coriolis effect (with respect to the axisymmetric case). For the same β and ℳ_{rot}, the Coriolis effect on an axisymmetric flux tube is going to be larger than on a nonaxisymmetric tube, simply because of the reduction of ℛ for higher azimuthal wavenumber.
The regimecontrolling relation for the nonaxisymmetric case becomes: (31)Replacing ℛ, we can rewrite Eq. (31) as a function of β, ℳ_{rot} and m: (32)The regime of the rise is more rotation dominated for lower azimuthal wavenumbers m for given β and ℳ_{rot}.
We can now introduce as (33)with (34) can be seen as a modified Lo that compares the ratio of buoyant force over altered Coriolis effect. Because ℳ_{ff}, and k are constant for a given m and both independent of β and ℳ_{rot}, we can introduce in Eq. (32) and demonstrate that the latter acts as a proxy to the force ratio (35)It can be seen that in the axisymmetric case, where f_{1} and f_{2} are zero, α_{1} and α_{2} both become unity and becomes (36)The definition of recovers the axisymmetric (m = 0) results.
To conclude this section, we identified a dimensionless number , that acts as a proxy to the force ratio. As such controls the regime of the rise for m = 0 and m> 0.
Assuming a given azimuthal wavenumber m for the initial conditions, we predict that two simulations with the same will reveal the same solution. As a result, the relative risetime of a magnetic flux tube in a rotating stellar interior should scale with . In the following Sects. 4 and 5, we will verify this prediction by studying the behavior of two series of numerical experiments.
4. Validation of the setup with a 2D numerical experiment
In agreement with Schüssler & Solanki (1992), we have shown that for the axisymmetric case reduces to the Lorentz number (). We already know that under the thinfluxtube approximation the relative risetime scales with (Choudhuri & Gilman 1987; Schüssler & Solanki 1992). Hence, we are interested in verifying whether this scaling behavior holds in the case of compressible simulations. For this purpose we carried out a series of simulations based on the STD2D setup.
4.1. Parameter study
The main goal of this section is to study the effect of rotation on the axisymmetric rise of a magnetic flux tube. This effect is controlled by . The latter is a composition of ℳ_{rot} and the plasmaβ, with ℳ_{rot} defining the rotational velocity of the star and the plasmaβ the buoyancy of the flux tube. In the case of an isothermal flux tube, β is directly proportional to the lack of density inside the flux tube, that is, the strength of the buoyant force. Hence, the parameter study is carried out in the (β, ℳ_{rot}) parameter space. As seen in Fig. 6, we covered two orders of magnitude for the β parameter and about one order of magnitude for the ℳ_{rot} parameter. This large domain is restricted by physical and numerical limits. Beyond those limits, simulations deliver either unreliable results due to high numerical diffusion (β> 200), or results that are not applicable to stellar interiors, such as β< 1 or rotation velocities being too close to the sound speed (ℳ_{rot} ≈ 1).
Fig. 6 Parameter space of the 2D series, where each point represents one simulation. The filled symbols represent simulations for which is a constant. Regions that are unphysical or numerically inaccessible are delimited by the dashed lines. The lowβ as well as the highℳ_{rot} limits are due to physical constraints. The highβ limit is of numerical nature. 

Open with DEXTER 
4.2. Verification of the scaling behavior
From a derivation of the force balance, we predicted in Sect. 3 that scalable properties of the rise of a magnetic flux tube should scale with . In order to verify our prediction we selected a set of simulations out of the parameter study described above. This subset corresponds to the filled symbols in Fig. 6. In this subseries we fixed and evolved ℳ_{rot} as well as the plasmaβ accordingly. The results are demonstrated in two Figures. First, Fig. 7 shows that the relative rise time is constant for a constant (=1.217), independently of β and ℳ_{rot}. Second, the topleft panel a of Fig. 8 shows the path of two flux tubes with the same but different values of ℳ_{rot} and β. From these plots, we find that over two orders of magnitude in β, neither the risetime nor the path of the flux tube changed significantly. This is a strong evidence that our axisymmetric setup scales with , as predicted by both the thinfluxtube approximation and our theoretical derivation (Sect. 3). The good scaling behavior convinced us that our setup delivers results comparable to hypothetical simulations with β = 10^{5}. It demonstrates the possibility of computing compressible simulations of magnetic flux tubes that rise in the same regime as in the Sun.
Fig. 7 Relative risetime versus plasmaβ also showing the corresponding ℳ_{rot} for a constant . Each point represents a simulation of the STD2D series where we vary β from about 1 to 100. 

Open with DEXTER 
Fig. 8 Crosssections of the magnetic field at various times for five representative simulations, representing the path taken by the flux tube all along its rise. Panel a is a plot of two simulations with the same (same regime): the green contours correspond to a strongly buoyant flux tube in a rapidly rotating interior (labeled a.1 in Table B.1), and the red contours correspond to a weakly buoyant flux tube in a slowly rotating interior (labeled a.2 in Table B.1). Panels b, c, and d show the dependence of the path on , in three different simulations. 

Open with DEXTER 
4.3. Validation by morphological study
In this section we address the question of the influence of on the path taken by the flux tube while rising. This question has been extensively studied (Choudhuri & Gilman 1987; Choudhuri & D’Silva 1990; Schüssler & Solanki 1992; Caligari et al. 1995, 1996; DeLuca et al. 1997; Fan 2008). We now show that our compressible setup reproduces the behavior found in the literature. In the axisymmetric case, the flux tube remains a torus throughout its rise. Furthermore, since our setup does not take into account turbulent convection, the angular momentum of the flux tube remains almost constant. We assume that molecular diffusion is not sufficient to transfer angular momentum from the flux tube to its surroundings. Therefore, an infinitely slowly rising flux tube will follow the path of constant angular momentum on which it lies initially. At finite rise speed, the flux tube follows a more complex path. This is illustrated in the panels b, c and d of Fig. 8. For a high (b) the flux tube rises radially; the flux tube is in a buoyancydominated regime. For a of approximately unity (c) the flux tube rises first radially and then axially; in that case the buoyant force and the Coriolis force acting on the flux tube are of the same order and the flux tube is in a transitional regime. A flux tube in a rotation dominated regime (d) will rise almost axially. These behaviors have been shown to be exclusively axisymmetric (Fan 2008). Since our setup is able to reproduce results from the literature – thinfluxtube results (Schüssler & Solanki 1992) as well as anelastic results (Fan 2008) – for corresponding values of , we are confident that our axisymmetric setup is reliable. It also shows that compressibility does not influence the path of the flux tube within the computational domain (r< 0.964 R_{⋆}).
4.4. General relation for the relative risetime in 2D
By carrying out a set of numerical experiments based on the STD2D setup, varying only , we investigated the influence of the rotationdominated or buoyancydominated character of the regime on the relative risetime. If the dependence is approximated by a power law as shown in Fig. 9, we obtain a relation (37)Unfortunately the proportionality factor relating and F_{buoy}/F_{corio} in Eq. (35) remains unknown and Eq. (37) fails to identify the nature of the regime. Several series would be required to fully identify how the proportionality factor depends on the drag force, the twist and thermal conduction. Instead, making use of a morphologic argument one can estimate the proportionality factor: in panel c of Fig. 8 the magnetic flux tube starts its latitudinal deflection at half of the convective zone; this indicates that for the magnetic flux tube rises in a transitional regime. We assume F_{buoy}/F_{corio} ≈ 1 for this case and estimate the proportionality factor to be about 0.7. In Fig. 9 we indicate the corresponding estimate of F_{buoy}/F_{corio} on the upper axis. Simulations lying on the lefthand side of unity are in a rotationdominated regime, whereas simulations on the righthand side of unity reveal a buoyancydominated regime.
Fig. 9 Relative risetime versus the scaling parameter for the STD2D series. We obtain a selfsimilar function of power − 2. The gray zone illustrates the low scatter of the individual values. The upper axis indicates the corresponding values of which represents the best estimate for F_{buoy}/F_{corio} we could extract from our simulations. 

Open with DEXTER 
The most relevant information to extract from Fig. 9 is the power of − 2 of the selfsimilar function that relates the relative rise time with . We call this exponent α_{3}. (38)with α_{1} = α_{2} = 1 and α_{3} = − 2, for axisymmetric rises. This result agrees with Jouve et al. (2010): for a given rotation rate the relative risetime is proportional to the inverse square of the magnetic field strength of the flux tube, that is, (39)We emphasize that this result also agrees with Schüssler & Solanki (1992) who underlined the fact that the risetime is a good scaling parameter for the present problem. However, the risetime is an a posteriori quantity. They showed therefore that the ratio of buoyancy over Coriolis force is a reasonable proxy for the risetime, and the risetime can be used as a scaling parameter. In the present work, we choose which, being an a priori quantity, is a more appropriate scaling parameter for the physical problem considered. Nevertheless, their conclusions are compatible with ours. Such an agreement shows that the thinfluxtube approximation is a good approximation in 2D.
5. Extension to the nonaxisymmetric case (3D)
In Sect. 3, we have predicted that in the nonaxisymmetric case, the regime of the rise of the flux tube is not controlled by anymore, but by a , with α_{1} and α_{2} being less than unity. In this section we compute both α_{1} and α_{2} from measurements of the curvature radius ℛ of magnetic flux tubes; we verify that the setup scales with the resulting dimensionless parameter; and confirm the predictions.
5.1. Defining the parameter study
In order to verify the existence of a scaling behavior of the setup, we conducted a parameter study for the nonaxisymmetric case based on the STD3D setup. We carried out this study in the same parameter space (β, ℳ_{rot}) as for the 2D series. In order to visualize this study, Fig. 10 shows the positions of all simulations in the parameter plane; while Fig. 11 illustrates the morphology of four representative flux tubes reaching the surface. We emphasize that, due to the lower resolution of the STD3D model, the numerical limit on β decreased from 200 to 40, as compared to the STD2D case.
Fig. 10 The (β, ℳ_{rot}) parameter space of the STD3D series, for an azimuthal wavenumber m = 8. Each point represents a simulation. Regions that are unphysical or numerically inaccessible are delimited by the dash lines. In a similar manner as for the STD2D series, the lowβ and highℳ_{rot} limits lead to unrealistic physical regimes, while the highβ limit is numerical. 

Open with DEXTER 
Fig. 11 Morphology of the flux tube when it emerges for the four different cases: panels a–d. The simulations are ordered with increasing running from the buoyancydominated regime a through d which is the rotationdominated regime. 

Open with DEXTER 
5.2. Validation of the general scaling parameter
In Sect. 3 we predicted that the regime of the rise is controlled by where α_{1} and α_{2} are functions of f_{1} and f_{2}. In order to validate this idea, we need to compute f_{1} and f_{2} from the simulations. To do so, we first need the curvature radius, ℛ.
Fig. 12 Projection of magnetic flux tubes on a horizontal plane cutting the northern hemisphere of the Sun at z = 0.34 R_{⋆} (with z being the height in cylindrical coordinates). Each point represents the projected position of the maximum entropy inside the flux tube in a given meridional plane. These points are supposed to map the center line of the magnetic flux tube. The radii of the dashed circles are a proxy for the curvature radii ℛ of the flux tubes at the apex. Panels a–d show suitable estimates of ℛ. Panels e and f show cases where the estimation fails. 

Open with DEXTER 
Computing the curvature radius is not a straightforward task. We use a method that consists in projecting the flux tube on a horizontal plane cutting the northern hemisphere of the star at a vertical distance from the equatorial plane of z = 0.34 R_{⋆}. From this projection we construct the circle going through the apex and two points, one of each leg of the loop, being two pressure scale heights deeper than the apex. We consider the radius of such a circle to be a good approximation for ℛ.
In Fig. 12 we plot a few examples of the circles we obtained from this simple method. It is important to note that in the case when the risetime is larger than the growth time of the kink instability (all our flux tubes are kink unstable) the shape of the flux tube becomes more complex than just a simple Ωloop. In such cases the simple method we used to determine ℛ fails; the values for ℛ are too large. Such situations are illustrated in panels e and f of Fig. 12. We exclude these cases from the computation of f_{1} and f_{2} and emphasize them in Table B.2 in italics. The four other panels a−d in Fig. 12 are ordered by decreasing ℛ. This figure shows visually the quality of the method to extract ℛ.
Figure 13 shows the curvature radius versus , with f_{1} and f_{2} being optimized to obtain the smallest residual of a linear fit.
Fig. 13 Dependency of the curvature radius ℛ at the apex, on both β and ℳ_{rot}. We obtain both exponents f_{1} and f_{2} by minimizing the residual of the best fit. Notes: not all simulations are present in this plot. We excluded the simulations for which the method we use to compute the curvature radius gives unreliable results. 

Open with DEXTER 
We find f_{1} = − 0.161 and f_{2} = − 0.145, which leads to α_{1} = 0.793 and α_{2} = 0.855 (see Eq. (34). As predicted in Sect. 3, f_{1} and f_{2} are both negative, which shows that the curvature radius is indeed reduced for less buoyant flux tubes (f_{1}< 0), and in faster rotating environments (f_{2}< 0). It also gives us access to α_{1} and α_{2} which are required to test the validity of Eq. (29). To test this, we could verify that the relative risetime, , scales with . We computed from the simulations the values of α_{1} and α_{2}, and as illustrated in Fig. 14, the relative risetime indeed scales with . We can conclude that the assumption on which our prediction was based is appropriate.
This result has several implications:

The radii of the fitted circles are good proxies for the curvatureradii.

The magnetic tension indeed influences the regime of the rise.

It is possible to simulate the nonaxisymmetric rise of compressible magnetic flux tubes for any solarlike star.
Fig. 14 Relative risetime versus the scaling parameter . As in the 2D case, the function is selfsimilar, but with a different power − 1.39. 

Open with DEXTER 
5.3. Morphology of the flux tube and emergence
The morphology of the nonaxisymmetric rise of magnetic flux tubes has already been discussed for anelastic simulations by Fan (2008) and Jouve et al. (2013), and we show here that our compressible simulations give similar results. In contrast to axisymmetric simulations, where magnetic flux tubes rise radially or parallel to the rotation axis depending on the regime of the rise (see Fig. 8), nonaxisymmetric flux tubes take a more radial path independently of the regime (see Figs. 11 and 15).
Fig. 15 Latitude of emergence versus the scaling parameter . The dashed line represents the initial latitude for all simulations, namely 20°. 

Open with DEXTER 
Alike the axisymmetric case, the angular momentum of the magnetic flux tube has to be conserved. As the latter moves away from the rotation axis, it decelerates. But according to Fan (2008), because of its Ωshape the nonaxisymmetric rising magnetic flux tube builds up a pressure gradient between its apex and its feet. As a result, mass flows along the tube’s center line in the direction opposite to the deceleration. Nonaxisymmetric flux tubes decelerate less than their axisymmetric counterparts (see Fan 2008, their Fig. 5) and the Coriolis force reduces accordingly. As a result, nonaxisymmetric flux tubes always rise radially.
Furthermore, as shown by Jouve et al. (2013), the asymmetry of the loop increases with its azimuthal deflection; Figs. 11 and 16 show similar behavior.
Fig. 16 Azimuthal deflection angle between the initial azimuth of the buoyant part (φ_{ini}) and the azimuth of the emergence region against the scaling parameter . We note that the simulation labeled (f.12) in Table B.2 is not plotted here because we lack emergence data. 

Open with DEXTER 
These morphologic characteristics were already observed in former anelastic studies. Therefore we conclude that our setup is reliable, and we demonstrate that the anelastic approximation is able to reproduce the morphology of compressible magnetic flux tubes.
5.4. General relation for the risetime in 3D
In order to build a BabcockLeighton dynamo that considers noninstantaneous rises (a delayed BabcockLeighton dynamo), we need a relation between the relative risetime () and the initial strength of the magnetic flux tube. In Fig. 14 we plot the relative risetime over . As in the axisymmetric setup the relative risetime follows a power law of , but in the nonaxisymmetric case, with a different α_{3} of, that is: − 1.39. (40)It is usually accepted that the relative risetime is proportional to the inverse square of the magnetic field strength, but as we have just demonstrated, this has to be reviewed. The STD3D series shows that the exponent depends on the azimuthal wavenumber m of the initial perturbation. For a given rotation rate (Ω) we can write in general (41)where α_{u} = α_{1}α_{2}α_{3} is the universal exponent, which depends on the azimuthal wavenumber m. The known values of α_{u} are summarized in Table 4.
Values of the various exponents for both series, with an azimuthal wavenumber of m = 0 and m = 8.
In the particular case of m = 8, the relative risetime follows the initial magnetic strength of the flux tube to the power − 0.942. This might have a significant impact on the solution of delayed BabcockLeighton dynamo.
Again the plot of Fig. 14 fails to identify the nature of the regime. But as in the axisymmetric case we can use a morphologic argument: panel c of Fig. 11 illustrates the less buoyant loop that remains symmetric until it reaches the surface. We suppose that in such a case both Coriolis and buoyant force are comparable. The regime becomes transitional for and therefore the proportionality factor of Eq. (35) is close to 0.4. We indicated the resulting estimate of the force ratio on the upper axis of Fig. 14. Any simulation with is in a rotationdominated regime; simulations with a largerthanunity are in a buoyancydominated regime.
Assuming that the morphological arguments we have used to compute the proportionality factor are acceptable, we can identify the regime of the rise. Figures 9 and 14 suggest that in the rotationdominated domain nonaxisymmetric flux tubes rise faster than axisymmetric ones. But for a given buoyancydominated regime the axisymmetric flux tube will rise faster.
These conclusions disagree with Fan (2008) where the author concludes that in any case nonaxisymmetric flux tubes rise faster than axisymmetric ones. Thanks to the derivation in Sect. 3, we can demonstrate that the author compared flux tubes rising in different regimes. Indeed, we do agree with Fan (2008) that, for a given β and a given rotation rate ℳ_{rot}, a higherm tube will rise faster, not because of nonaxisymmetry, but because its regime has changed and tends towards a buoyancydominated regime.
Additionally, we point out that the nonaxisymmetric α_{3} is closer to zero than the axisymmetric one. This aspect is beyond the scope of the prediction made in Sect. 3 (the prediction was made on α_{1} and α_{2} only). Nevertheless we would like to suggest a physical interpretation. The pressure gradient that emerges between the apex and the feet of a nonaxisymmetric magnetic flux tube is independent of , but depends on stratification and m. We suggest that α_{3} depends on the competition between two mechanisms: the change in path taken by the magnetic flux tube and the importance of the pressure gradient building up between the apex and the feet of the flux tube. In the axisymmetric case the pressure gradient is absent, α_{3} is set exclusively by the variation of paths taken by the flux tube. In the nonaxisymmetric case there is no variation of the paths, so the pressure gradient should be responsible for the value of α_{3}. This leads to the conclusion that in the case of an m = 4 mode the exponent should lie between − 2.0 and − 1.4. This statement agrees with the conclusion drawn in MorenoInsertis (1986), where the author found an exponent of − 1.8. To confirm this interpretation one should run further simulations to extract the dependence of α_{3} on m.
6. Discussion and conclusions
It is widely accepted that solarlike stars maintain their magnetic field by a dynamo process. The choice of the dynamo process is still heavily debated. In BabcockLeighton (or fluxtransport) dynamos, the rise of magnetic flux tubes is an essential ingredient as it links the strong toroidal fields of the tachocline with the nearsurface source term of the dynamo.
In the various implementations of BLdynamos, the risetime of magnetic flux tubes had been assumed to be instantaneous until Jouve et al. (2010), who used a more complete description of the risetime of magnetic flux tubes and discussed the impact of the resulting delayed generation of poloidal field on the dynamo. However, the model describing the risetime was rather simple. We decided to go further and study the risetime of magnetic flux tubes in direct numerical simulations.
A scaling relation for the axisymmetric rise of thin flux tubes was already proposed by Schüssler & Solanki (1992). In the present work, we carried out simulations of nonaxisymmetrically rising magnetic flux tubes in rotating, compressible, and stratified interiors. Nonaxisymmetric perturbations to flux tubes lead to the rise of Ωshaped loops. These loops have a much smaller curvature radius than the corresponding torus in the axisymmetric case. As a result, the tension force reduces notably the rise velocity of the magnetic flux tube. We predicted that this difference changes the scaling relation of the relative risetime depending on the azimuthal wavenumber m of the initial perturbation. This wavenumber controls how strongly the curvature radius is reduced. We present a theoretical approach describing this phenomenon. From this theory, we predict that the parameter controlling the regime of the rise, formerly introduced by Choudhuri & Gilman (1987), needs to be redefined as This dimensionless parameter is proportional to the ratio of the buoyant force over the Coriolis force acting on the flux tube. It defines three regimes: the rotation dominated one when , the transitional regime with , and the buoyancy dominated regime for .
From two series of simulations we carried out with m = 0 and m = 8, we computed α_{1} and α_{2} in both cases and verified that the regime of the rise is indeed controlled by . We then compared our compressible simulations with the previous anelastic and thinfluxtube studies. We found that compressibility neither influences the morphology, nor the risetime of magnetic flux tubes.
Finally, we focused on our main goal, namely extracting the scaling relation of the relative risetime in light of the modified scaling parameter. We found that in contrast to former conclusions, for a given rotation rate Ω, the relative rise time does not necessarily scale with the inverse square of the magnetic field, but that the power law depends on the azimuthal mode of the rising flux tube: which equals (B_{φ}/B_{eq})^{αu} for Ω = const. For m = 0 we confirm α_{3} = − 2, while for m = 8 we found α_{3} = − 1.39. This leads to risetimesthat last less than 0.1 P_{rot} for the buoyancydominated regime, about P_{rot} for the transitional regime, and more than 10 P_{rot} for the rotationdominated regime. This scaling relation can be applied to any main sequence solarlike star, regardless of its rotation period and internal magnetic fields.
Nevertheless, the present model still neglects convection and differential rotation. Our conclusions probably hold for flux tubes being sufficiently magnetic to be unsensitive to convective motion. But convection and differential rotation will influence the above scaling relation. As shown by Weber et al. (2011), weak flux tubes will be advected by convective motions and will rise faster. Nevertheless the present results are a step toward a general description of flux rise that can possibly be tested using solar and stellar cycle observations.
Magnetic tension can explain the variation of α_{1} and α_{2}, but does not address the variation of α_{3}. Furthermore the fact that is a proxy for the force ratio prevents this parameter from identifying the nature of the regime in which the flux tube rises. Explaining why the relative risetime does not scale with the same power law of for various m, and identifying the dependencies of the proportionality factor on viscosity, twist, and thermal conduction, represent a numerical challenge and will require further simulations.
With a complete theory we will be able to construct a universal model for the risetime of magnetic flux tubes in lowmass stars, but the present conclusions are already sufficient to enrich the idea of delayed BabcockLeighton dynamos.
Acknowledgments
We would like to thank Laurène Jouve, Jörn Warnecke, Sydney Barnes and Federico Spada for constructive discussions and comments. This work was supported by the DFG grant Ar 355/91.
References
 Babcock, H. W. 1962, Trans. Int. Astron. Union, Series B, 11, 419 [Google Scholar]
 Brown, B. P., Browning, M. K., Brun, A. S., Miesch, M. S., & Toomre, J. 2010, ApJ, 711, 424 [NASA ADS] [CrossRef] [Google Scholar]
 Browning, P. K., & Priest, E. R. 1983, ApJ, 266, 848 [NASA ADS] [CrossRef] [Google Scholar]
 Brun, A. S., Miesch, M. S., & Toomre, J. 2004, ApJ, 614, 1073 [NASA ADS] [CrossRef] [Google Scholar]
 Caligari, P., Schüssler, M., Stix, M., & Solanki, S. K. 1994, in Cool Stars, Stellar Systems, and the Sun, ed. J.P. Caillault, 64, 387 [Google Scholar]
 Caligari, P., MorenoInsertis, F., & Schüssler, M. M. 1995, ApJ, 441, 886 [NASA ADS] [CrossRef] [Google Scholar]
 Caligari, P., Schüssler, M., Solanki, S. K., Schaerer, D., & Stix, M. 1996, Astrophys. Lett. Comm., 34, 17 [NASA ADS] [Google Scholar]
 Cattaneo, F., & Hughes, D. W. 1987, in NASA Conference Publication 2483, eds. G. Athay, & D. S. Spicer, 101 [Google Scholar]
 Cattaneo, F., Tzihong, C., & Hughes, D. W. 1989, Bull. Amer. Phys. Soc., 34, 1294 [Google Scholar]
 Charbonneau, P. 2010, Liv. Rev. Sol. Phys., 7, 3 [Google Scholar]
 Charbonneau, P. 2014, ARA&A, 52, 251 [NASA ADS] [CrossRef] [Google Scholar]
 Cheung, M. C. M., & Isobe, H. 2014, Liv. Rev. Sol. Phys., 11 [Google Scholar]
 Cheung, M. C. M., MorenoInsertis, F., & Schüssler, M. 2006, A&A, 451, 303 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Chou, D.Y., & Fisher, G. H. 1989, ApJ, 341, 533 [NASA ADS] [CrossRef] [Google Scholar]
 Choudhuri, A. R., & D’Silva, S. 1990, A&A, 239, 326 [NASA ADS] [Google Scholar]
 Choudhuri, A. R., & Gilman, P. A. 1987, ApJ, 316, 788 [NASA ADS] [CrossRef] [Google Scholar]
 DeLuca, E. E., Fan, Y., & Saar, S. H. 1997, ApJ, 481, 369 [NASA ADS] [CrossRef] [Google Scholar]
 Fan, Y. 2001, ApJ, 546, 509 [NASA ADS] [CrossRef] [Google Scholar]
 Fan, Y. 2008, ApJ, 676, 680 [NASA ADS] [CrossRef] [Google Scholar]
 Fan, Y., Fisher, G. H., & McClymont, A. N. 1994, ApJ, 436, 907 [NASA ADS] [CrossRef] [Google Scholar]
 Fan, Y., Zweibel, E. G., & Lantz, S. R. 1998, ApJ, 493, 480 [NASA ADS] [CrossRef] [Google Scholar]
 Fan, Y., Featherstone, N., & Fang, F. 2013, ArXiv eprints [arXiv:1305.6370] [Google Scholar]
 Granzer, T. 2004, Astron. Nachr., 325, 417 [NASA ADS] [CrossRef] [Google Scholar]
 Granzer, T., Schüssler, M., Caligari, P., & Strassmeier, K. G. 2000, A&A, 355, 1087 [NASA ADS] [Google Scholar]
 Hughes, D. W., Wissink, J. G., Matthews, P. C., & Proctor, M. R. E. 1997, in 1st Advances in Solar Physics Euroconference, Advances in Physics of Sunspots, eds. B. Schmieder, J. C. del Toro Iniesta, & M. Vazquez, 118, 66 [Google Scholar]
 Jouve, L., & Brun, A. S. 2009, ApJ, 701, 1300 [NASA ADS] [CrossRef] [Google Scholar]
 Jouve, L., Proctor, M. R. E., & Lesur, G. 2010, A&A, 519, A68 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Jouve, L., Brun, A. S., & Aulanier, G. 2013, ApJ, 762, 4 [NASA ADS] [CrossRef] [Google Scholar]
 Käpylä, P. J., Korpi, M. J., Brandenburg, A., Mitra, D., & Tavakol, R. 2010, Astron. Nachr., 331, 73 [NASA ADS] [CrossRef] [Google Scholar]
 Käpylä, P. J., Mantere, M. J., & Brandenburg, A. 2012, ApJ, 755, L22 [NASA ADS] [CrossRef] [Google Scholar]
 Leighton, R. B. 1969, ApJ, 156, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Longcope, D. W., & Klapper, I. 1997, ApJ, 488, 443 [NASA ADS] [CrossRef] [Google Scholar]
 Matthews, P. C., Hughes, D. W., & Proctor, M. R. E. 1995, ApJ, 448, 938 [NASA ADS] [CrossRef] [Google Scholar]
 Miesch, M. S. 2005, Liv. Rev. Sol. Phys., 2, 1 [Google Scholar]
 MorenoInsertis, F. 1983, A&A, 122, 241 [NASA ADS] [Google Scholar]
 MorenoInsertis, F. 1986, A&A, 166, 291 [NASA ADS] [Google Scholar]
 MorenoInsertis, F. 1997, in 1st Advances in Solar Physics Euroconference. Advances in Physics of Sunspots, eds. B. Schmieder, J. C. del Toro Iniesta, & M. Vazquez, 118, 45 [Google Scholar]
 MorenoInsertis, F., & Emonet, T. 1996, ApJ, 472, L53 [NASA ADS] [CrossRef] [Google Scholar]
 MorenoInsertis, F., Schüssler, M., & FerrizMas, A. 1992, A&A, 264, 686 [NASA ADS] [Google Scholar]
 Nelson, N. J., Brown, B. P., Sacha Brun, A., Miesch, M. S., & Toomre, J. 2014, Sol. Phys., 289, 441 [NASA ADS] [CrossRef] [Google Scholar]
 Parker, E. N. 1955, ApJ, 122, 293 [NASA ADS] [CrossRef] [Google Scholar]
 Pinto, R. F., & Brun, A. S. 2013, ApJ, 772, 55 [NASA ADS] [CrossRef] [Google Scholar]
 Racine, E., Charbonneau, P., Ghizaru, M., Bouchat, A., & Smolarkiewicz, P. K. 2011, ApJ, 735, 46 [NASA ADS] [CrossRef] [Google Scholar]
 Rädler, K.H., Kleeorin, N., & Rogachevskii, I. 2003, Geophys. Astrophys. Fluid Dynamics, 97, 249 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Rempel, M., & Cheung, M. C. M. 2014, ApJ, 785, 90 [NASA ADS] [CrossRef] [Google Scholar]
 Schüssler, M. 1980, Nature, 288, 150 [NASA ADS] [CrossRef] [Google Scholar]
 Schüssler, M., & Rempel, M. 2002, in From Solar Min to Max: Half a Solar Cycle with SOHO, ed. A. Wilson, 508, 499 [Google Scholar]
 Schüssler, M., & Solanki, S. K. 1992, A&A, 264, L13 [NASA ADS] [Google Scholar]
 Spruit, H. C. 1981, A&A, 102, 129 [NASA ADS] [Google Scholar]
 Spruit, H. C., & van Ballegooijen, A. A. 1982, A&A, 106, 58 [NASA ADS] [Google Scholar]
 Strassmeier, K. G. 2009, A&ARv, 17, 251 [NASA ADS] [CrossRef] [Google Scholar]
 van Ballegooijen, A. A. 1983, A&A, 118, 275 [NASA ADS] [Google Scholar]
 Warnecke, J., Losada, I. R., Brandenburg, A., Kleeorin, N., & Rogachevskii, I. 2013, ApJ, 777, L37 [NASA ADS] [CrossRef] [Google Scholar]
 Weber, M. A., & Fan, Y. 2015, Sol. Phys., 290, 1295 [NASA ADS] [CrossRef] [Google Scholar]
 Weber, M. A., Fan, Y., & Miesch, M. S. 2011, ApJ, 741, 11 [NASA ADS] [CrossRef] [Google Scholar]
 Wissink, J. G., Hughes, D. W., Matthews, P. C., & Proctor, M. R. E. 2000, MNRAS, 318, 501 [NASA ADS] [CrossRef] [Google Scholar]
 Yoshimura, H. 1985, PASJ, 37, 171 [NASA ADS] [Google Scholar]
 Ziegler, U. 2011, J. Comp. Phys., 230, 1035 [NASA ADS] [CrossRef] [Google Scholar]
 Ziegler, U. 2012, SIAM J. Sci. Comput., 34, C102 [CrossRef] [Google Scholar]
Appendix A: Adiabatically stratified atmosphere
In this section, we derive the profiles of the three hydrodynamical quantities T, ρ, and P for an adiabatically stratified atmosphere. These profiled are used to mimic a stellar convective zone avoiding convective motions but still reproducing a realistic stellar stratification. In such a case we have to presume that the luminosity of the star is transported exclusively by radiation: Here, Hence, Integrating over radius we obtain a function for T(r), Setting the following boundary condition: T(R_{0}) = T_{0}, In that case the luminosity and κ are the input parameters, but we need to control the stratification. Hence, we rewrite the two latter variables as functions of more meaningful quantities, namely ∇ (the logarithmic temperature gradient) and pressurescale height. Moreover the pressurescale height at the surface is a useful length for our problem since: Moreover the stellar interior is in hydrostatic equilibrium, where which leads to a convenient definition of the pressurescale height: With both ∇ and H_{P0}, we rewrite T(r) in a more convenient way. We first relate ∇ and H_{P0}. which gives, and considering the ideal gas equation of state: Finally, where we immediately identify H_{P0} and rewrite as a function of both: By integration we obtain a convenient temperature profile, where we are able to control both the pressurescale height and the adiabaticity of our layer, Since ∇ is a constant all over the domain it also has the advantage of relating T and P: leading to Reusing the equation of state it becomes clear that:
Appendix B: List of simulations
List of all axisymmetric simulations based on the STD2D setup.
List of all nonaxisymmetric simulations based on the STD3D setup.
All Tables
Fixed parameters for the numerical setup and the initial conditions of the STD2D case.
Fixed parameters for the numerical setup and the initial conditions of the STD3D case.
Values of the various exponents for both series, with an azimuthal wavenumber of m = 0 and m = 8.
All Figures
Fig. 1 Sketch of the numerical setup. The yellow contour shows the buoyant part of the flux tube with a lack of entropy. This section of the flux tube will buoyantly rise toward the surface and emerge as an active region. The orange contour surfaces show the neutrally buoyant parts of the flux tube, with an excess of entropy. The latter sections will remain rooted at the bottom of the simulated domain shown by the wedgelike shape of solid lines. The azimuthal morphology shown here is valid only in 3D, for an azimuthal wavenumber of m = 8. In 2D, the flux tube is buoyant everywhere, and the computational domain is a meridional plane. The direction of rotation is also indicated. 

Open with DEXTER  
In the text 
Fig. 2 Radial profiles of ρ, T, and P, normalized by their respective values at the top of the domain. The radial profile of the normalized gradient of entropy is a constant (adiabatic) in this configuration. 

Open with DEXTER  
In the text 
Fig. 3 Profile of the dimensionless magnetic quantities over the normalized distance from the flux tube center. B_{φ}, B_{p}, and P_{m} are the toroidal component of the twisted magnetic field, the poloidal component of the twisted magnetic field, and the resulting magnetic pressure, respectively. 

Open with DEXTER  
In the text 
Fig. 4 Profiles of the magnetic pressure and the thermodynamical quantities for a magnetic flux tube in two extreme hydrostatic states, a in an isothermal state and b in a neutrally buoyant state. 

Open with DEXTER  
In the text 
Fig. 5 Sketch representing the approximate shape of rising flux tubes depending on the regime of their rise. The light blue flux tube is in a buoyancy dominated regime. The dark blue flux tube is in a transitional regime, where F_{buoy} ≈ F_{corio}. The red flux tube is in a rotation dominated regime. 

Open with DEXTER  
In the text 
Fig. 6 Parameter space of the 2D series, where each point represents one simulation. The filled symbols represent simulations for which is a constant. Regions that are unphysical or numerically inaccessible are delimited by the dashed lines. The lowβ as well as the highℳ_{rot} limits are due to physical constraints. The highβ limit is of numerical nature. 

Open with DEXTER  
In the text 
Fig. 7 Relative risetime versus plasmaβ also showing the corresponding ℳ_{rot} for a constant . Each point represents a simulation of the STD2D series where we vary β from about 1 to 100. 

Open with DEXTER  
In the text 
Fig. 8 Crosssections of the magnetic field at various times for five representative simulations, representing the path taken by the flux tube all along its rise. Panel a is a plot of two simulations with the same (same regime): the green contours correspond to a strongly buoyant flux tube in a rapidly rotating interior (labeled a.1 in Table B.1), and the red contours correspond to a weakly buoyant flux tube in a slowly rotating interior (labeled a.2 in Table B.1). Panels b, c, and d show the dependence of the path on , in three different simulations. 

Open with DEXTER  
In the text 
Fig. 9 Relative risetime versus the scaling parameter for the STD2D series. We obtain a selfsimilar function of power − 2. The gray zone illustrates the low scatter of the individual values. The upper axis indicates the corresponding values of which represents the best estimate for F_{buoy}/F_{corio} we could extract from our simulations. 

Open with DEXTER  
In the text 
Fig. 10 The (β, ℳ_{rot}) parameter space of the STD3D series, for an azimuthal wavenumber m = 8. Each point represents a simulation. Regions that are unphysical or numerically inaccessible are delimited by the dash lines. In a similar manner as for the STD2D series, the lowβ and highℳ_{rot} limits lead to unrealistic physical regimes, while the highβ limit is numerical. 

Open with DEXTER  
In the text 
Fig. 11 Morphology of the flux tube when it emerges for the four different cases: panels a–d. The simulations are ordered with increasing running from the buoyancydominated regime a through d which is the rotationdominated regime. 

Open with DEXTER  
In the text 
Fig. 12 Projection of magnetic flux tubes on a horizontal plane cutting the northern hemisphere of the Sun at z = 0.34 R_{⋆} (with z being the height in cylindrical coordinates). Each point represents the projected position of the maximum entropy inside the flux tube in a given meridional plane. These points are supposed to map the center line of the magnetic flux tube. The radii of the dashed circles are a proxy for the curvature radii ℛ of the flux tubes at the apex. Panels a–d show suitable estimates of ℛ. Panels e and f show cases where the estimation fails. 

Open with DEXTER  
In the text 
Fig. 13 Dependency of the curvature radius ℛ at the apex, on both β and ℳ_{rot}. We obtain both exponents f_{1} and f_{2} by minimizing the residual of the best fit. Notes: not all simulations are present in this plot. We excluded the simulations for which the method we use to compute the curvature radius gives unreliable results. 

Open with DEXTER  
In the text 
Fig. 14 Relative risetime versus the scaling parameter . As in the 2D case, the function is selfsimilar, but with a different power − 1.39. 

Open with DEXTER  
In the text 
Fig. 15 Latitude of emergence versus the scaling parameter . The dashed line represents the initial latitude for all simulations, namely 20°. 

Open with DEXTER  
In the text 
Fig. 16 Azimuthal deflection angle between the initial azimuth of the buoyant part (φ_{ini}) and the azimuth of the emergence region against the scaling parameter . We note that the simulation labeled (f.12) in Table B.2 is not plotted here because we lack emergence data. 

Open with DEXTER  
In the text 