A&A 370, 635-648 (2001)
A. Bardou1 - B. von Rekowski1 - W. Dobler1 - A. Brandenburg1,2 - A. Shukurov1
1 - Department of Mathematics, University of Newcastle, Newcastle upon Tyne NE1 7RU, UK
2 - Nordita, Blegdamsvej 17, 2100 Copenhagen, Denmark
Received 23 November 2000 / Accepted 15 February 2001
We consider the effect of vertical outflows on the mean-field dynamo in a thin disk. These outflows could be due to winds or magnetic buoyancy. We analyse both two-dimensional finite-difference numerical solutions of the axisymmetric dynamo equations and a free-decay mode expansion using the thin-disk approximation. Contrary to expectations, a vertical velocity can enhance dynamo action, provided the velocity is not too strong. In the nonlinear regime this can lead to super-exponential growth of the magnetic field.
Key words: accretion, accretion disks - magnetic fields - galaxies: magnetic fields
Magnetic fields play an important role in accretion disks. Inside the disk, for example, they allow the magneto-rotational instability to develop (Balbus & Hawley 1998). This instability can lead to the development of turbulence which in turn transports angular momentum radially outwards, thus allowing the matter to accrete. On the other hand, large scale magnetic fields are often considered necessary for launching winds or jets from accretion disks. The origin of such large scale magnetic fields is still unclear: they may be advected from the surrounding medium or be generated by a dynamo inside the disk. Advection appears unlikely as turbulence leads to enhanced viscosity and magnetic diffusivity, so that the two are of the same order (Prandtl number of order unity; Pouquet et al. 1976). Thus, turbulent magnetic diffusivity can compensate the dragging of the field by viscously induced accretion flow (van Ballegooijen 1989; Lubow et al. 1994; Heyvaerts et al. 1996). Dynamo action is a plausible mechanism for producing magnetic fields in accretion disks (Pudritz 1981; Stepinski & Levy 1988). However, dynamo magnetic fields can affect winds from accretion disks (Campbell et al. 1998; Campbell 1999; Brandenburg et al. 2000). On the other hand, the wind can also affect the dynamo. In particular, an enhancement by wind was suggested by Brandenburg et al. (1993) in the context of galactic dynamos. In the same context, Elstner et al. (1995) also consider the effects of shear in the wind.
In this paper, we consider the generation of magnetic field by a dynamo acting in a thin accretion disk. We provide a study of the effects of vertical velocities. We show how vertical velocities can enhance the dynamo, allowing a larger growth rate and leading to super-exponential growth of the magnetic field in the nonlinear regime.
The vertical velocities considered can have several origins. They can, for example, be due to a wind emanating from the disk or to magnetic buoyancy inside the disk. Note that we are interested here in magnetic buoyancy as generating vertical outflows. We neglect the associated dynamo effect (Moss et al. 1999).
The paper is organized as follows. In Sect. 2, we describe the two complementary approaches used and specify our choice of parameters. In Sect. 3, we discuss the linear results without vertical outflows. Results with vertical outflows are discussed in Sect. 4 for the linear regime and in Sect. 5 for the nonlinear regime. Section 6 presents our conclusions.
We use two approaches to solve this equation. Firstly, a 2D numerical simulation based on a finite-difference scheme is applied (Sect. 2.1). Secondly, we develop a 1D free-decay mode expansion (Sect. 2.2). Details of these two complementary approaches are given below. Each of these two methods has its own advantages. The finite-difference simulation allows investigation of the nonlinear behaviour and applies to 2D geometry. On the other hand, the free-decay mode expansion gives the full eigenstructure and provides a tool for qualitative analysis. Furthermore, boundary conditions corresponding to an ideal vacuum can be implemented easily in the 1D model, whereas our 2D model where the disk is embedded in a halo can provide only an approximation to a vacuum outside the disk. Our two approaches are therefore complementary.
The problem is governed by three dimensionless numbers which appear
in the dimensionless form of Eq. (1).
Let h, the half-thickness of the disk, be a characteristic length-scale,
a characteristic value of ,
Vz0 a characteristic vertical velocity and
value of .
We also define the rotational shear
The dimensionless numbers are the
azimuthal magnetic Reynolds
Our computational domain contains the disk embedded in a halo: , , where the disk occupies the region , , R=1.5 being the outer radius of the disk and h=0.15 its half-thickness. This results in a disk aspect ratio of h/R = 0.1 and the upper boundary of the domain at . Using a resolution of grid points and a uniform grid, the resulting mesh spacing is . Equation (5) is solved using a sixth order finite-difference scheme in space and a third order Runge-Kutta time advance scheme.
We control the symmetry of our solutions (dipolar or quadrupolar) by using
appropriate initial conditions.
On the boundaries of the computational domain (i.e. far away from the disk)
we impose the conditions
that the normal component of ,
vanishes together with the normal
derivative of the tangential components,
is antisymmetric about
the midplane and vanishes outside the disk. We adopt
As appropriate for accretion disks, we adopt a softened Keplerian angular
velocity profile in r in the whole
computational domain (including the halo as well as the disk), assuming
to be z-independent:
The turbulent magnetic diffusivity is given by an interpolation between
in the disk midplane and
in the halo,
To obtain the magnetic Reynolds numbers of Eqs. (2)-(4) we choose .
In accretion disks, the dynamo number is approximately constant
with radius. In this model this is achieved by setting
We consider solutions with exponential time dependence
with respect to the
orthonormal base formed
by the free-decay modes (Galerkin expansion),
After making this expansion, the differential equations become an algebraic set of
equations forming an eigenvalue problem of 2(N+1) equations
and Ynm are defined by
The rotation of the thin disk is assumed to be close to Keplerian and independent of z for the finite-difference model (see Sect. 2.1). Only the assumption of the independence of on z is needed in the free-decay mode expansion.
We use two configurations: one where the magnetic diffusivity is the same in the disk and in the halo around the disk ( ) and one with a low conductivity in the halo or a vacuum around the disk. Note that assuming a vacuum or low conductivity in the halo above the disk is a rather good approximation for galaxies (Sokoloff & Shukurov 1990; Poezd et al. 1993). It also allows for a connection with established theory. Note finally that the low-conductivity approximation ( ) is a good approximation for potential fields ( ) and hence for a particular case of force-free fields ( ).
|Figure 1: Real part of the growth rate of the magnetic field as a function of , with , obtained with the free-decay-mode expansion with . Asterisks denote non-oscillatory solutions whereas diamonds denote oscillatory solutions. Solid lines are for quadrupolar modes and dotted lines for dipolar modes. Signs above the curves indicate the sign of|
|Open with DEXTER|
The boundary conditions of the free-decay mode expansion imply a vacuum above the disk. Technically, to achieve in this model, we restrict the dynamo action to an inner layer of the disk. The outer part of the disk is thus playing the role of the halo. The results simply need to be rescaled for a proper comparison.
In the finite-difference model, we study the case and the case of a low conductivity halo with . Note that, with this model, to approximate a vacuum in the halo we have to choose as high as possible. But, this choice is restricted because the higher this ratio, the smaller the time step.
In this paper, the model with is referred to as the "disk+halo'' model whereas the model with a vacuum or low conductivity above the disk is referred to as the "disk+vacuum'' model.
In both models, the profile for
is taken to be
In this work, the two control parameters are the dynamo number and the vertical magnetic Reynolds number . For reference, we give here possible values for these two parameters. While dynamo numbers are rather low in galactic disks, (Ruzmaikin et al. 1988), in accretion disks, they can be as high as , depending on the Mach number of the turbulence, see e.g. Pudritz (1981). The value of obviously depends on the origin of the vertical outflow. With magnetic buoyancy, a rough estimate can be obtained assuming that the vertical motions occur at the Alfvén speed. Estimating the Alfvén speed from magnetic equipartition, one gets a vertical magnetic Reynolds number of order unity.
In the following, the figures obtained with the finite-difference model are labelled with "FD'' and the ones obtained with the free-decay mode (Galerkin) expansion with "GE''.
We first compare the results of the free-decay mode expansion (with ) with the results published by Stepinski & Levy (1991). We thus take a linear dependence of on z, . The calculated growth rates are plotted in Fig. 1. The agreement is good but a discrepancy is observed at large where the increase of the growth rate with is faster according to Stepinski & Levy than in our free-decay mode expansion for both quadrupolar and dipolar symmetry ( ). This discrepancy seems to be due to an insufficient number of grid points in Stepinski & Levy's model. Indeed, if we increase the truncation level of the free-decay mode expansion, no change is observed in the results whereas if we decrease the truncation level, we recover the results of Stepinski & Levy. Note that the results discussed in this paper are for where the agreement is excellent.
|Figure 2: a-d) Real part of the growth rate of the magnetic field as a function of for . Asterisks denote non-oscillatory solutions and diamonds denote oscillatory solutions. Solid lines indicate quadrupolar modes and dashed lines indicate dipolar modes. a) and c) are for the finite-difference model (FD) and b) and d) are for the free-decay mode expansion model (GE) with a) and b) for the "disk+halo'' model and c) and d) for the "disk+vacuum'' model. Note the approximate symmetry of graphs a) and b) with respect to the vertical axis. Note also that the dipolar mode in c) is the "forgotten mode'', not seen in d)|
|Open with DEXTER|
|Figure 3: a-d) Real part of the growth rate of the magnetic field as a function of for the dominant quadrupolar mode and for different dynamo numbers. All modes are non-oscillatory. a) and c) are for the finite-difference model (FD) and b) and d) are for the free-decay mode expansion model (GE) with a) and b) for the "disk+halo'' model and c) and d) for the "disk+vacuum'' model. The models of panels a) to d) correspond to the ones in Fig. 2|
|Open with DEXTER|
|Figure 4: a-d) Vertical structure of the non-oscillatory quadrupolar magnetic field for . Panels a) and c) are for the finite-difference model (FD) while panels b) and d) are for the free-decay mode expansion (GE). Panels a) and b) give the radial structure and c) and d) give the azimuthal structure of the magnetic field. In the finite-difference model, eigenfunctions are taken at r=1 and are normalized such that for each . Note that in the free-decay mode expansion, we plot and , i.e. the time-independent part of the magnetic field|
|Open with DEXTER|
|Figure 5: a-d) Magnetic field configuration obtained using the finite-difference model with the linear "disk+halo'' configuration and for the dominant non-oscillatory quadrupolar mode. Contours represent poloidal magnetic field lines with solid lines being clockwise and dotted lines being counter-clockwise. The greyscales are the intensity and direction of the azimuthal magnetic field; light grey (dark grey) signifies ( ). The parameters are a) , , b) , , c) , , d) , (where is maximum). Comparing a) and b) shows the effect of increasing the dynamo number. b- d) show the effect of increasing|
|Open with DEXTER|
We now want to compare the finite-difference simulation with the free-decay mode expansion. In particular, we want to test our method to mimic, in the case of the free-decay mode expansion, a model for starting with the "disk+vacuum'' model (see Sect. 2). For an infinitely extended halo with , the growth rates of dipolar and quadrupolar modes should be exchanged when reverses sign, i.e. the graph of as a function of should be symmetric with respect to the vertical axis (Proctor 1977). The finite-difference simulation shows this symmetry with a good precision (Fig. 2a). The results of the free-decay mode expansion (Fig. 2b) are less symmetric. This is not surprising since the ratio between the height of the computational domain and the height of the disk is about 7 in the finite-difference simulation whereas it is only 2 in the free-decay mode expansion.
We now discuss our results for the "disk+halo'' and "disk+vacuum'' configuration. Let us first consider the "disk+halo'' configuration (Figs. 2a and b). For moderate dynamo number , the leading mode is quadrupolar and non-oscillatory if ( ) and it is dipolar and non-oscillatory if ( ). In the following, we thus study the effect of vertical motions for in quadrupolar symmetry and for in dipolar symmetry. Above , the symmetry of the leading mode changes and it becomes oscillatory. The critical dynamo number is .
Let us now consider the "disk+vacuum'' configuration (Figs. 2c and d). The results obtained with the two approaches are rather different. With the free-decay mode expansion, the dipolar mode is always oscillatory and subcritical for . On the other hand, with the finite-difference simulation, the dipolar mode becomes supercritical for and is non-oscillatory for all . This is because the dominant dipolar mode in the finite-difference simulation is the "forgotten mode'' i.e. a dipolar mode corresponding to n=0 in Eqs. (25) and (26). The "forgotten mode'' cannot be obtained with our 1D calculation (free-decay mode expansion). However, we reproduced the corresponding profile with a modified 1D model which takes into account the finite aspect ratio of the disk. The model is based on the assumption of separability of r- and z-dependence and yields equations similar to the thin-disk Eqs. (12), (13), but with additional correction terms of order (h/r)2. Thus, the dipolar modes of Figs. 2c and 2d are of a different nature. As for the quadrupolar mode, the results obtained with the two approaches are not so different. The branch for is rather similar in both cases with a critical dynamo number . The branch for is all subcritical with the finite-difference simulation whereas it becomes supercritical for with the free-decay mode expansion. However, we have checked that, in the finite-difference model, if we take a steeper transition between and , the branch for moves upwards, resembling the one for the free-decay mode expansion. The qualitative behaviour for the case of low conductivity in the halo is in many respects similar to that of the "disk+halo'' configuration (Figs. 2a and c). But there are quantitative differences. There is no intersection point where the symmetry of the leading mode would change for . The diagram does not look symmetric anymore.
The results of the finite-difference model presented here are in agreement with the accretion disk models of Rüdiger et al. (1995) and von Rekowski et al. (2000). Note, however, that these models have an inner radius and that at this radius a perfectly conducting boundary condition is used, which leads to a vertical current layer on the boundary and potentially lowers the excitation threshold for magnetic fields.
We now discuss the magnetic field structure. Figure 4 shows the eigenfunctions for Br and for and for different in the "disk+halo'' configuration. With this choice of the dynamo number, a maximum in the curve occurs around . An obvious feature is that the vertical velocity pushes the zeros of Br and outwards. For the dynamo to work, it is essential that Br changes sign in the disk since only then can magnetic flux of the opposite sign to that in the midplane leave through the disk surface (Ruzmaikin et al. 1988). When takes large values, the radial component of the magnetic field Br has no zero inside the disk and thus the dynamo cannot work ( ). This is the technical reason for the decrease of with . The physical reason for it is investigated in Sect. 4.3. The increase of the growth rate at small will be explained in Sect. 4.3. We note however that when we decrease the truncation level of the free-decay mode expansion down to N=0 (one doublet), the maximum disappears.
The magnetic field configuration in 2D geometry is shown in Fig. 5. It shows both the effect of increasing the dynamo number and the effect of increasing the vertical magnetic Reynolds number. When one decreases the dynamo number from to , in the absence of any vertical velocity, the magnetic field becomes more concentrated at the outer edge of the disk (Figs. 5a and 5b). If the vertical magnetic Reynolds number is then increased from , the magnetic field is radially redistributed over the whole disk (Fig. 5c) and then advected towards the halo (Fig. 5d). From panel (b) to (d), the growth rate increases as shown by the curve for in Fig. 3a. Thus, the magnetic field diffusion becomes less efficient as the vertical scale of the magnetic field increases with . This leads to an increase of the growth rate up towards its maximum at (Fig. 5d). Note that in the halo the poloidal magnetic field lines become more and more aligned with lines of constant angular velocity as increases. Note also that the mode structure of Fig. 5d (where is maximum for ) is very similar to the mode structure for , where that mode has attained its maximum growth rate.
We now study the linear behaviour with vertical velocity for positive where the dipolar non-oscillatory mode is dominant. Figure 6 shows that there is no maximum in the growth rates; at of order 1 or even less the dynamo is switched off.
The magnetic field structure in 2D geometry is shown in Fig. 7. The magnetic field structures for and and 0.75 show that the field is advected vertically outwards and the poloidal field lines are aligned with the contours lines very quickly.
|Figure 6: Real part of the growth rate of the magnetic field as a function of for the dominant dipolar mode and for different dynamo numbers in the "disk+halo'' model. Asterisks denote non-oscillatory modes whereas diamonds denote oscillatory ones|
|Open with DEXTER|
|Figure 7: a) and b) Magnetic field configuration obtained using the finite-difference model with the linear "disk+halo'' configuration and for the dominant non-oscillatory dipolar mode. The parameters are a) , , b) , . The notation is the same as in Fig. 5|
|Open with DEXTER|
The general picture is as follows. The positive contribution to comes from and whereas dissipation has a negative contribution. In the presence of a vertical velocity, an additional negative contribution comes from the stretching of the magnetic field . But, at the same time, the stretching reduces dissipation by increasing the scale of the magnetic field (see Eq. (40)). Thus, depending on the relative values of the reduction of dissipation and of the negative stretching term, the growth rate is enhanced or not (a maximum occurs or not). Obviously, when the magnetic field is completely stretched, dissipation cannot be reduced anymore. However, this does not necessarily correspond to the maximum of the growth rate as the terms and evolve themselves with (which quantifies the vertical velocity). These different quantities are plotted in Fig. 8 in a case where a maximum occurs and in a case with no maximum.
We note that dissipation increases again after being reduced by the stretching of the magnetic field in the vertical direction. This increased dissipation is due to the formation of boundary layers near z=h. We also note that the term decreases whereas increases. The evolution of these terms is linked to the disappearance of the positive area of Br (with the convention of sign of Fig. 4). Finally, note that at large , the growth rate always decreases as the stretching terms become dominant.
|Figure 8: a-d) Contributions of the different terms in Eqs. (36) and (37) to the growth rate in the "disk+vacuum'' model. Lines are labelled with the subscripts used for the different terms in Eqs. (36) and (37). a) and c) represent the radial and azimuthal terms for , for which a maximum occurs, b) and d) represent the radial and azimuthal terms for , for which no maximum occurs|
|Open with DEXTER|
In short, the growth rate can increase with , thus forming a maximum, because dissipation is reduced by the stretching of the magnetic field.
We define an effective vertical magnetic Reynolds number
|Figure 9: Magnetic energy and effective vertical magnetic Reynolds number as a function of time for the "disk+halo'' model with magnetic buoyancy (cf. Eq. (42)). and . Super-exponential growth occurs before saturation sets in. The mode is quadrupolar and non-oscillatory. The initial magnetic field strength is . The steep increase for is due to the generation of azimuthal field from the initial large scale poloidal one at the chosen high|
|Open with DEXTER|
In this subsection we assume that the only nonlinearity is due to magnetic buoyancy
and thus appears in the vertical velocity term only. In particular, we do not consider
here any dependence of
In this case the dimensionless form of the induction equation in the steady state reads
The magnetic field is extended
into the whole halo and almost completely stretched.
Therefore, adding a wind term such that
|Figure 10: Magnetic energy as a function of for the "disk+halo'' model with -quenching. All modes are quadrupolar and non-oscillatory|
|Open with DEXTER|
|Figure 11: a-c) Magnetic field configuration for the "disk+halo'' model with -quenching: the effect of vertical velocity. Poloidal field (lines) and azimuthal field (grey scales). a) , ; b) , ; c) , . a-c) following the lower curve in Fig. 10. All modes are quadrupolar and non-oscillatory|
|Open with DEXTER|
Figure 10 shows the magnetic energy as a function of for the three dynamo numbers , and . From this we can see the following features:
Perfect agreement between the two curves cannot be expected, because the averaging in Eq. (43) was somewhat arbitrary; the agreement could perhaps be improved by choosing a different definition of . The fact that with magnetic buoyancy the magnetic energies are larger for a given value of can therefore not be interpreted as a relative enhancement of dynamo action.
|Figure 12: Magnetic energy as a function of defined in Eq. (43) for the models with -quenching with wind and -quenching combined with magnetic buoyancy. and . All modes are quadrupolar and non-oscillatory. With -quenching alone,|
|Open with DEXTER|
Disk dynamos have been studied intensively over the past 20 years, but only recently the effects of outflows and winds have been included (e.g. Campbell 1999; Dobler et al. 1999; Brandenburg et al. 2000). It was immediately clear that these changes in the boundary conditions affect the structure of the solutions. The present work is an attempt to provide a more thorough survey of solutions with wind for different values and sign of the dynamo number. One of the effects that we have highlighted is the enhancement of the kinematic dynamo growth rate for intermediate strengths of the wind. In the nonlinear case this corresponds to a phase shortly before saturation where the growth becomes super-exponential.
In the present work we have deliberately isolated the dynamo problem from the wind problem. In reality the two aspects together constitute a more complicated self-consistent dynamo problem that needs to be looked at in more detail in future. One step in that direction was taken by Brandenburg et al. (2000) where the feedback on the flow from the Lorentz force was taken fully into account. A more formidable task would be to obtain self-consistent dynamo action in three dimensions without assuming an -effect, thus obtaining dynamo action and outflows in a self-consistent manner from first principles. In the context of accretion disks, this would require simulating turbulent dynamo action, where the turbulence itself is driven by dynamo-generated magnetic fields (Brandenburg et al. 1995). This would be of great interest not only from the astrophysical point of view, but also from that of dynamo theory in general, because, as was recently realized, the catastrophic quenching of the standard -dynamo could be alleviated in the presence of open boundary conditions which may substantially improve the working of the dynamo (Blackman & Field 2000). This is simply because the -effect generates helical fields, but at the same time magnetic helicity is conserved in a closed domain. Thus, being able to get rid of helical fields through open boundaries may be vital for the dynamo (Brandenburg & Dobler 2001). Furthermore, the production of fields at large scales could be facilitated by the production - and subsequent loss through boundaries - of small scale fields with opposite sign of magnetic helicity. This latter aspect could not be addressed with mean-field models and would require a direct simulation of accretion disk turbulence in a global model.
We acknowledge financial support from PPARC (Grant PPA/G/S/1997/00284) and from the Leverhulme Trust (Grant F/125/AL). We thank D. Sokoloff for helpful discussions. The use of the PPARC supported parallel computer GRAND in Leicester and the PPARC funded Compaq MHD Cluster in St Andrews is acknowledged.