# Warped discs and the directional stability of jets in Active Galactic Nuclei

## 1 Introduction

The angular momentum of gas accreting onto massive black holes in Active Galactic Nuclei (AGN) is likely to change with time, either as a result of mergers funneling fresh gas towards the nucleus (Hernquist & Mihos 1995), or because of radiation or disc wind driven warping instabilities in an existing accretion flow (Pringle 1996, 1997; Maloney, Begelman & Pringle 1996; Maloney, Begelman & Nowak 1998). In either case, the angular momentum of gas at large radius in the disc will be misaligned with the rotation axis of the black hole, while at small radius the combined action of viscosity and differential precession induced by the Lense-Thirring effect (Lense & Thirring 1918) leads to alignment of the disc and hole angular momenta (Bardeen & Petterson 1975; Kumar & Pringle 1985). Achieving this alignment requires that the black hole exert a torque on the disc gas, and implies that an equal and opposite torque act on the black hole itself. Over time, this causes a change in the spin axis of the hole towards alignment with the large angular momentum reservoir provided by the disc at large radius.

Whether the timescale for alignment is long or short compared to the lifetime of an AGN has implications for our understanding of the accretion history of these objects. If the black hole is spinning, then jets, regardless of whether they derive power from the aligned inner disc or directly from the hole itself (Blandford & Znajek 1977; see also Ghosh & Abramowicz 1997; Livio, Ogilvie & Pringle 1999), trace the spin of the black hole. If the alignment timescale is long, then we expect that the jet direction reflects the initial spin of the black hole established during (or prior) to the formation of the host galaxy (Rees 1978). Moreover the jet direction is expected to be stable over time, irrespective of variations in the angular momentum of accreting gas. Conversely, if the timescale is short then jets define the angular momentum of inflowing material, and the observation of systems where they appear to be stable over long time periods of $`10^7`$ \- $`10^8`$ yr (Alexander & Leahy 1987; Liu, Pooley & Riley 1992; Scheuer 1995) implies a constancy in the average angular momentum of accreting gas.

The rate of realignment was computed by Scheuer & Feiler (1996), who found an approximate analytic solution to the equations governing the evolution of a warp. The timescale depends on two poorly known viscosities, $`\nu _1`$ and $`\nu _2`$, which correspond to the $`(R,\varphi )`$ and $`(R,z)`$ components of the shear as defined by Papaloizou & Pringle (1983). If $`\nu _1\nu _2`$ (i.e. the timescale for warp evolution is comparable to that for the evolution of the surface density) then the formula of Scheuer & Feiler (1996) is similar to that previously assumed by Rees (1978). Subsequently, Natarajan & Pringle (1998) revisited the problem, and pointed out that for the parameters of AGN discs the assumption that $`\nu _1\nu _2`$ is at variance with the results of detailed analyses of the hydrodynamics of viscous discs, which suggest instead that $`\nu _2\nu _1`$. Recalculating the alignment timescale with this modification they found that it was short — much less than any reasonable estimate of AGN lifetime — and used this to argue that the spin of black holes in AGN, and the direction of jets produced in the inner regions of the accretion flow, rapidly adjust to trace the angular momentum of gas at large radius.

The analysis in Natarajan & Pringle (1998) employed the same torque formula as Scheuer & Feiler (1996), which was derived under the restrictive conditions of small amplitude warps and a disc viscosity that is constant with radius (or, equivalently, that in a steady-state the disc surface density is independent of radius). Since AGN may well harbour strongly warped discs, and almost certainly have more complex surface density profiles, our goal in this paper is to generalise their results. Section 2 below sets out the equations for the evolution of a warped disc, which we solve numerically in Section 3 for comparison with the approximate analytic solution. The resulting torque is calculated in Section 4, and the alignment timescale derived for a simple AGN disc model in Section 5. Section 6 discusses the implications of our results.

## 2 Warped disc equations

### 2.1 Governing equations

The evolution of a thin, warped accretion disc in a Keplerian potential was studied by Papaloizou & Pringle (1983). For a disc with a surface density profile $`\mathrm{\Sigma }(R,t)`$, with angular momentum at radius $`R`$ parallel to the unit vector $`\widehat{l}`$, the evolution can be expressed most compactly in terms of an equation for the angular momentum density $`𝐋=(GMR)^{1/2}\mathrm{\Sigma }\widehat{l}`$. Adopting the simplest description (Pringle 1992; see also Papaloizou & Pringle 1983), the time evolution is given by,

$`{\displaystyle \frac{𝐋}{t}}`$ $`=`$ $`{\displaystyle \frac{3}{R}}{\displaystyle \frac{}{R}}\left[{\displaystyle \frac{R^{1/2}}{\mathrm{\Sigma }}}{\displaystyle \frac{}{R}}(\nu _1\mathrm{\Sigma }R^{1/2})𝐋\right]`$ (1)   
$`+`$ $`{\displaystyle \frac{1}{R}}{\displaystyle \frac{}{R}}\left[(\nu _2R^2\left|{\displaystyle \frac{\widehat{l}}{R}}\right|^2{\displaystyle \frac{3}{2}}\nu _1)𝐋\right]`$   
$`+`$ $`{\displaystyle \frac{1}{R}}{\displaystyle \frac{}{R}}\left({\displaystyle \frac{1}{2}}\nu _2R|𝐋|{\displaystyle \frac{\widehat{l}}{R}}\right)+{\displaystyle \frac{\stackrel{}{\omega }_p\times 𝐋}{R^3}}.`$   
Here $`\nu _1`$ and $`\nu _2`$ are viscosities acting on the $`(R,\varphi )`$ and $`(R,z)`$ components of the shear, as defined by Papaloizou & Pringle (1983). Neither is equal to the shear viscosity ‘$`\nu `$’ that enters into the Navier-Stokes equations, and the general relation between these three quantities is extremely complex (Ogilvie 1999). For the parameters of AGN discs, the simplest assumptions suggest that $`\nu _2\nu _1`$ (Natarajan & Pringle 1998), though for most of this paper we will allow both viscosities to be general power laws with radius,

$`\nu _1`$ $`=`$ $`\nu _{10}R^\beta `$   
$`\nu _2`$ $`=`$ $`\nu _{20}R^\beta .`$ (2)   
For thin, vertically averaged accretion disc models, the shear viscosity $`\nu `$ is a power-law in radius in a given opacity regime (e.g. Frank, King & Raine 1992). Later, the index $`\beta `$ will be chosen to correspond to the scaling expected at radii in the disc where the bulk of the realignment torque acts. However, for a general warped disc, the radial run of both $`\nu _1`$ and $`\nu _2`$ will also depend on the shape of the warp (Ogilvie 1999). This effect is not taken into account in this paper. Our assumption that $`\nu _1`$ and $`\nu _2`$ are simple power-laws in radius, independent of warp shape, is therefore strictly valid in the limit of small warps, and is an approximation for larger amplitude warps.

The range of behaviour allowed by equation (1) is considerably more diverse than simple diffusion of $`𝐋`$. Nonetheless, we will make frequent use of the characteristic timescales for diffusion of surface density and warp, which we define conventionally as,

$$t_{\nu _1}=\frac{R^2}{\nu _1},t_{\nu _2}=\frac{R^2}{\nu _2}.$$
(3)   
The final term on the right hand side of equation (1) represents the effect of Lense-Thirring precession on the disc. For a black hole with mass $`M`$, we take $`\stackrel{}{\omega }_p=\omega _𝐩(0,0,1)`$, where

$$\omega _𝐩=2ac\left(\frac{GM}{c^2}\right)^2$$
(4)   
where $`a`$ is the dimensionless angular momentum of the black hole.

Several assumptions and simplifications are necessary to derive this equation, and these are set out in detail in Papaloizou & Pringle (1983) (see also Pringle 1999; Ogilvie 1999). Particularly important is the boundary between the diffusive behaviour described by equation (1) and wave-like evolution, which occurs when the usual Shakura & Sunyaev (1973) viscosity parameter $`\alpha H/R`$, where $`H`$ is the disc scale height (Papaloizou & Lin 1995). The observational constraints on $`\alpha `$ in AGN discs are weak (Siemiginowska & Czerny (1989) and Siemiginowska, Czerny & Kostyunin (1996) suggest $`\alpha =10^110^2`$), but in the ionized parts of the disc it is plausible to assume that $`\alpha 0.1`$, as in better studied disc systems (eg. Cannizzo 1993). In any case, theoretical estimates for $`H/R`$ at the relevant radii (Collin-Souffrin & Dumont 1990) suggest that in AGN $`\alpha H/R`$, so that we are safely in the diffusive regime. In this respect AGN discs are very different from protostellar discs, for example, where wave-like behaviour is likely to be important (e.g. Larwood et al. 1996). We also note that although Papaloizou & Pringle (1983) considered only small amplitude warps describable using linear perturbation theory, Ogilvie (1999) has shown that strongly warped discs can also be described using a similar approach.

### 2.2 The Scheuer & Feiler analytic solution

Scheuer & Feiler (1996) obtained an approximate steady state analytic solution to equation (1) for the case where $`\beta =0`$, in the limit where the angle of misalignment of the disc and the hole spin was small. In this case, the surface density obeys the usual relation for a planar disc with accretion rate $`\dot{M}`$ and zero torque inner boundary condition imposed at $`R_{\mathrm{in}}`$,

$$\nu _1\mathrm{\Sigma }=\frac{\dot{M}}{3\pi }\left(1\sqrt{\frac{R_{\mathrm{in}}}{R}}\right).$$
(5)   
The shape of the disc is given by $`\widehat{l}=(l_x,l_y,l_z)`$, where

$`l_x=Ke^\varphi \mathrm{cos}\varphi `$   
$`l_y=Ke^\varphi \mathrm{sin}\varphi `$ (6)   
and $`\varphi =2(\omega _p/\nu _2R)^{1/2}`$. The constant $`K`$ depends on the disc inclination $`i`$ at large radius, at large radius $`K=\mathrm{sin}i`$.

### 2.3 Numerical solutions

The above analytic solution applies only for $`\beta =0`$, and as derived assumes that the disc warp is of small amplitude. Relaxing these assumptions requires a steady state numerical solution of equation (1).

For $`\beta <1/2`$, steady state solutions can be calculated efficiently by solving directly the ordinary differential equations (ODEs) obtained by setting $`\mathrm{\Sigma }/t=0`$ in equation (1). We give a full description of our methods in the Appendix, but briefly we have found it simplest to use equations expressed in terms of the surface density $`\mathrm{\Sigma }(R)`$ and unit tilt vector $`\widehat{l}=𝐋/|𝐋|`$, which can readily be solved iteratively using a finite-difference technique described by Pereyra (1979). We use the Numerical Algorithms Group’s implementation, and impose boundary conditions of $`\mathrm{\Sigma }(R_{\mathrm{in}})=0`$, $`\mathrm{\Sigma }(R_{\mathrm{out}})=\mathrm{\Sigma }_{\mathrm{out}}`$, $`\widehat{l}(R=R_{\mathrm{in}})=(0,0,1)`$, $`\widehat{l}(R=R_{\mathrm{out}})=\widehat{l}_{\mathrm{out}}`$. These boundary conditions assume that the inner disc is aligned with the spin axis of the hole, we therefore have to check post facto that this is indeed a valid assumption for each model calculated.

For $`\beta 1/2`$, we have been unable to obtain a converged numerical solution to the steady state ODEs using these methods. In this regime, we instead evolve the time dependent equation (1) from an arbitrary initial condition (of a flat, uniformly tilted disc) until a steady state is obtained, using the numerical method described by Pringle (1992). This is straightforward but computationally expensive, since many viscous times at the outer edge of the disc are required to obtain a steady state solution. Defining the viscous time of the disc by $`t_{\nu _1}=R_{\mathrm{out}}^2/\nu _1`$, we find that more than 10 viscous times of evolution are necessary. Moreover, for either method capturing the torque accurately requires a large range between the inner and outer disc radii, of at least $`R_{\mathrm{out}}/R_{\mathrm{in}}\stackrel{>}{}10^4`$. For these reasons, we are only able to obtain solutions using the time dependent code for high values of $`\beta `$ ($`\beta =1.5`$), when the viscous time is a weak function of increasing radius.

## 3 Steady-state disc shape

### 3.1 Comparison with the analytic solution

Figure 1 shows the numerical steady-state solution for the $`\beta =0`$ case, and compares it to the approximate analytic solution. We adopt units in which $`R_{\mathrm{in}}=1`$, and take $`\nu _{10}=\nu _{20}=1`$. We choose $`\omega _p=200`$, which gives an inner aligned region extending out to around $`10^2R_{\mathrm{in}}`$. This is roughly appropriate for plausible AGN parameters, though we defer until later consideration of detailed disc models. Convergence to the limiting inclination at infinity is slow, so it is necessary to impose the outer boundary condition of $`i=0.1`$ at a large radius, $`R=10^5R_{\mathrm{in}}`$. At this radius $`l_x`$ and $`l_y`$ are chosen to facilitate comparison with the analytic solution of equation (6).

Evidently, the analytic solution provides an excellent description of the disc shape. The inner disc is aligned to the rotation axis of the hole out to close to $`10^2R_{\mathrm{in}}`$, and warps to around half the limiting inclination by $`10^3R_{\mathrm{in}}`$. Moving inward from infinity, the warp is twisted by an angle of $`\pi `$ before being flattened effectively into the aligned plane.

### 3.2 Aligned radius

For an annulus in the disc, the local timescales for precession and for transmission of warp are given by,

$$t_p=\frac{R^3}{\omega _p},t_{\nu _2}=\frac{R^2}{\nu _{20}R^\beta }.$$
(7)   
We expect that the disc will be aligned with the spin axis of the hole at radii where the precession timescale is much shorter than the timescale over which the disc can communicate warp inwards. The characteristic radius of the aligned region $`R_{\mathrm{align}}`$ will then scale with the radius where $`t_p=t_{\nu _2}`$, which is given by,

$$R_{\mathrm{align}}=C_1\left(\frac{\omega _p}{\nu _{20}}\right)^{1/(1+\beta )},$$
(8)   
where $`C_1`$ is a constant expected to be of the order of unity which we will use later to fit the numerical results. Note that this expression assumes both that the Lense-Thirring torque is strong enough to ensure that the inner disc is aligned, and that $`\beta >1`$, so that $`t_p`$ increases with radius more rapidly than $`t_{\nu _2}`$.

Figure 2 shows how the disc shape varies with $`\beta `$ for a set of models with fixed $`\omega _p`$ and $`\nu _1=\nu _2=1`$. From the numerical solutions, we compute $`R_{\mathrm{align}}`$ as the radius where $`i`$ first exceeds some small inclination $`\delta `$. With $`C_1=0.165`$ there is excellent agreement between the computed $`R_{\mathrm{align}}`$ and the scaling given by equation (8), i.e. the simple timescale argument given above suffices to give an accurate idea of the scaling of the radius out to which the hole can enforce inner disc alignment.

## 4 Torque

The torque exerted on the misaligned disc as a result of the Lense-Thirring precession is given by the integral,

$$\frac{\mathrm{d}𝐉}{\mathrm{d}t}=\frac{\stackrel{}{\omega }_p\times 𝐋}{R^3}2\pi RdR.$$
(9)   
Of course an oppositely directed torque of the same magnitude is exerted on the black hole. Computing $`\dot{𝐉}`$ for the disc model shown in Fig. 1 we find that the torque agrees with that obtained from the Scheuer & Feiler solution at the level of a few percent (the small discrepancy being primarily due to the ignored inner boundary condition in the analytic solution, which leads to an error in the surface density in the warp region). The torque scales as $`(\omega _p\nu _2)^{1/2}\dot{M}`$, again as found analytically. We have computed additional models with large inclination angles at infinity (up to $`80^{}`$ angle of misalignment) and for those models there is a larger discrepancy between the numerical and analytic results, but still only at the tens of percent level.

### 4.1 Radial dependence

Figure 3 shows the tilt angle and surface density profile for three fiducial models, all normalised to the same accretion rate. The solid line is computed for the parameters used for Figure 1, namely $`\beta =0`$ and $`\nu _{10}=\nu _{20}`$. The short dashed line maintains $`\beta =0`$, but has $`\nu _{10}=1`$, $`\nu _{20}=8`$. Natarajan & Pringle (1998) find that consideration of the hydrodynamics of a viscous disc suggests that in AGN,

$$\frac{\nu _2}{\nu _1}\frac{1}{2\alpha ^2}$$
(10)   
in which case this ratio of viscosities would be appropriate for a plausible Shakura-Sunyaev $`\alpha `$ of around 0.25. Finally the long dashed line illustrates the effect of a declining surface density profile, plotted is a model with $`\beta =1/4`$ and $`\nu _1=\nu _2`$. All three of these models are chosen to have the same viscous timescale (both for the surface density and for warp) at the (arbitrary) choice of inner edge radius. In conventional thin disc models (Shakura & Sunyaev 1973), fixing the viscous time at a given radius corresponds to a fixed disc thickness $`H/R`$ at that radius. We note that in such models the value of $`\beta `$ is fixed by solving for the vertical disc structure, i.e. it is not a ‘parameter’ that can be varied but rather a property of the disc. In Section 5 we discuss the appropriate value for $`\beta `$ in the region of the disc where the dominant torque arises. Also shown for each of the models is where radially the greatest contribution to the torque arises, we plot the components of the integrand in equation (9) (the extra factor of $`R`$ in the plotted quantity takes account of the logarithmic radial scale to give a true impression of where the strongest torque on the disc is found).

For all the models, strong torques arise from a broad region that extends from just outside the alignment radius out to one or two orders of magnitude larger radius. The large radial range used in these calculations is thus essential to capture the torque accurately. Increasing the ratio of $`\nu _2`$ relative to $`\nu _1`$ reduces the timescale for diffusion of warp, $`t_{\nu _2}=R^2/\nu _2`$, and allows the warp to push in closer to the black hole. At fixed accretion rate increasing $`\nu _2`$ even by this large factor leaves the surface density profile almost unaltered, and thus the net effect is simply to increase the magnitude of the torque and move it to smaller disc radius. Larger values of $`\beta `$ likewise decrease the warp timescale at large radius and lead to a shrinkage of the aligned region. For the parameters used here, the disc shape is in fact rather similar for the two models in which $`\nu _2>1`$ at large radius. However, for a given accretion rate, models with higher $`\beta `$ have lower surface density at large radii, which tends to decrease $`𝐋`$ in equation (9) and reduce the integrated torque. Both these effects have a significant impact on the calculation of the total torque.

### 4.2 Scaling with $`\beta `$

The numerical results show that the Scheuer & Feiler (1996) solution, equation (6), provides an accurate description of the disc shape and the resultant torque on the black hole. For this solution, the torque integral in equation (9) has a magnitude given by,

$$|\dot{𝐉}|=\frac{\sqrt{2}K\dot{M}}{3\nu _{10}}(GM)^{1/2}(\omega _p\nu _{20})^{1/2}.$$
(11)   
Guided by the numerical results, a straightforward extension of this result to non-zero $`\beta `$ is possible. We first assume that the radius $`R_w`$ out to which we need to integrate equation (9) to find the torque scales as does the aligned radius,

$$R_w=C_2\left(\frac{\omega _p}{\nu _{20}}\right)^{1/(1+\beta )},$$
(12)   
If we additionally assume that the disc shape is primarily a function of $`R_w`$, then the form of the integral in equation (9) suggests the following scaling with $`\beta `$,

$$|\dot{𝐉}|=\frac{K\dot{M}}{3\sqrt{2}\nu _{10}}(GM)^{1/2}\frac{C_2^\beta \omega _p}{\beta +1/2}\left(\frac{\nu _{20}}{\omega _p}\right)^{\frac{2\beta +1}{2\beta +2}},$$
(13)   
which reduces to the formula in equation (11) if $`\beta =0`$. We have allowed ourselves a single adjustable parameter in this expression, to best match the numerical results we adopt $`C_2=0.55`$.

Figure 4 shows how equation (13) compares to the numerical results for choices of parameters testing the scalings implied by the formula. Excellent agreement (at better than the percent level) is obtained for two models with $`\nu _1=\nu _2`$ and varying $`\omega _p`$, while a single model computed with the time-dependent code at $`\beta =1.5`$ also agrees well with the fitting formula. A model with $`\nu _2>\nu _1`$ displays some systematic error, caused by the influence of the increased ratio of $`\nu _2/\nu _1`$ on the surface density profile in the warp region, but equation (13) still provides a reasonable approximation to the numerical values for the torque. The main result is that the torque for a fixed accretion rate (and fixed viscous time at $`R_{\mathrm{in}}`$) drops by roughly an order of magnitude for disc models with $`\beta `$ varying between $`\beta =0`$ and $`\beta =1.5`$, while the dependence on other parameters ($`\omega _p`$, $`\nu _{20}`$) remains close to that found in the $`\beta =0`$ case.

## 5 Alignment timescale

The torque on the black hole due to the interaction with the warped disc causes it to precess and become aligned with the outer disc’s angular momentum vector on a timescale that, to factors of order unity, is just $`t_{\mathrm{align}}=|𝐉|K/|\dot{𝐉}|`$ (e.g. Scheuer & Feiler 1996). $`|\dot{𝐉}|`$ is specified for a given disc model by equation (13), while the angular momentum of a hole of mass $`M`$ and spin parameter $`a`$ is $`|𝐉|=acM(GM/c^2)`$. The resulting expression for $`t_{\mathrm{align}}`$ is similar, but not equivalent, to the estimate given by Natarajan & Pringle (1998), $`t_{\mathrm{align}}=(|𝐉|/|𝐉_{\mathrm{disk}}|)(R_{\mathrm{warp}}^3/\omega _p)`$, where $`𝐉_{\mathrm{disk}}`$ is the angular momentum of the disc within a warp radius $`R_{\mathrm{warp}}`$ defined similarly to equation (12). Our expression takes more complete account of the need to integrate the torque for all radii less than $`R_w`$. For models of the disc where most of the angular momentum within the warp radius lies close to that radius, the final inferred timescale for alignment is similar.

Obtaining an estimate for $`\beta `$ and $`\nu _{10}`$ requires a model for the vertical structure of the disc, which serves the purpose of relating the unknown central disc temperature, that controls the viscosity, to the effective temperature fixed by the requirement that viscous dissipation balance radiative losses. At the radii of the warp, gas pressure and an electron scattering opacity $`\kappa =0.4\mathrm{cm}^2\mathrm{g}^1`$ are dominant. In this regime, Collin-Souffrin & Dumont (1990) obtain a radial profile of column density $`N_{25}`$ (in units of $`10^{25}\mathrm{cm}^2`$) described by,

$`N_{25}=98\alpha ^{4/5}\left({\displaystyle \frac{ϵ}{0.1}}\right)^{3/5}\left({\displaystyle \frac{L/L_{\mathrm{Edd}}}{0.1}}\right)^{2/5}`$   
$`\times \left({\displaystyle \frac{L}{10^{44}\mathrm{erg}/\mathrm{s}}}\right)^{1/5}\left({\displaystyle \frac{R}{10^4R_g}}\right)^{3/5},`$ (14)   
where $`ϵ`$ is the radiative efficiency, $`L=ϵ\dot{M}c^2`$ is the bolometric luminosity, $`L_{\mathrm{Edd}}=4\pi c^3R_g/\kappa `$ is the Eddington luminosity, and $`R_g=2GM/c^2`$. For this model, $`\beta =0.6`$. Although this disc model includes the effects of irradiation of the upper layers of the disc by the central source, we note that the scaling $`\mathrm{\Sigma }\alpha ^{4/5}\dot{M}^{3/5}M^{1/5}R^{3/5}`$ is identical to that of non-irradiated models (e.g. Sincell & Krolik 1998). The normalisation is within a factor of two. Other uncertainties thus dominate over that arising from different treatments of the vertical disc structure.

Figure 5 shows the calculated alignment timescale for a $`10^8M_{}`$ black hole accreting at rates between $`10^2\dot{M}_{\mathrm{Edd}}`$ and $`\dot{M}_{\mathrm{Edd}}`$, where $`\dot{M}_{\mathrm{Edd}}=4\pi cR_g/(ϵ\kappa )`$. We take a radiative efficiency $`ϵ=0.1`$, and assume a maximal Kerr hole, $`a=1`$. For these parameters, $`R_w`$ as defined in equation (12) is $`R_w10^3R_g`$ for reasonable ratios of $`\nu _2/\nu _1`$, so the main warped region of the disc falls self-consistently within the regime where gas pressure and electron scattering opacity are dominant. The Figure shows results calculated assuming that $`\nu _2=\nu _1`$, along with the expectation if $`\nu _2/\nu _1=1/(2\alpha ^2)`$ for a range of plausible $`\alpha `$. The latter models of course predict more rapid realignment, as noted previously (Natarajan & Pringle 1998).

For holes accreting at rates of the order of the Eddington limit, the alignment timescale is found to be $`10^5`$ \- $`10^6`$ yr. This is short compared to lifetimes of AGN, irrespective of what one assumes for the ratio of $`\nu _2/\nu _1`$. Producing the giant structures observed in radio galaxies may require an active epoch of the order of $`10^8`$ yr, while some estimates for the quasar lifetime are similar (e.g. Haehnelt & Rees 1993, though these estimates are subject to considerable uncertainty). In these systems we expect that there was ample time for a rotating black hole, whatever its initial spin axis, to become aligned with the angular momentum of the disc. Thus jets accelerated from the inner disc should be perpendicular to the plane of the outer disc, and constancy in the direction of such jets implies a corresponding stability in the disc plane over time.

Figure 6 plots the alignment timescale as a function of black hole mass, for a range of accretion rates. The alignment timescale is found to be only weakly dependent on the mass of the hole, as expected given the $`M^{1/16}`$ dependence derived by Natarajan & Pringle (1998). For this figure we have assumed that $`\nu _1=\nu _2`$, and that $`\alpha =0.1`$. Shorter alignment timescales are of course predicted if $`\nu _2>\nu _1`$, as shown in Figure 5.

For holes accreting at lower rates (relative to the Eddington limit), the timescale for alignment grows. It is not unambiguously shorter than estimates for the active phase of a hole in a galactic nucleus, and the answer depends crucially on differences in the disc’s response to shear and warp. It remains true that $`t_{\mathrm{align}}t_{\mathrm{grow}}`$, where $`t_{\mathrm{grow}}=M/\dot{M}`$, for any $`\dot{M}`$, so that alignment will always occur if an accretion event lasts long enough to contribute significantly to the hole mass. Perhaps more probably, however, low luminosity AGN may accrete negligible fractions of the hole mass in many brief episodes. In this scenario, it remains possible that the hole spin accumulated over the whole accretion history will be sufficient to control the jet direction during subsequent low luminosity outbursts. The outer disc would then not be expected to be generally perpendicular to the jet direction.

## 6 Discussion

In this paper we have considered the interaction between a misaligned accretion disc and a rotating black hole for the parameters appropriate to Active Galactic Nuclei. A combination of forced precession and disc viscosity allows the black hole to force the inner disc into a plane perpendicular to the hole’s spin axis, but the resulting warped disc exerts a torque that eventually aligns the rotation of the hole with the angular momentum of the outer disc. We have presented steady-state numerical solutions to the full warped disc equations, verified that the analytic solution of Scheuer & Feiler (1996) provides an accurate representation of the disc shape and torque on the hole, and generalised the results for a range of disc models with varying viscosity. The range in models we consider lead to almost an order of magnitude variation in the integrated torque on the hole.

We have computed the timescale for alignment to occur, and find that for holes accreting at rates of the order of the Eddington limit (for a $`10^8M_{}`$ hole, this implies $`L_{\mathrm{bolometric}}=2.5\times 10^{46}\mathrm{erg}/\mathrm{s}`$, so these are luminous AGN) the timescale is short, $`10^5`$ \- $`10^6`$ yr. This is the case for a range of assumptions as to the ratio $`\nu _2/\nu _1`$, which characterizes the disc’s relative response to azimuthal shear and to warp. The short derived timescale implies that the spin of black holes in such AGN rapidly adjusts to match the angular momentum of accreting gas, even if the hole gains only a small ($``$ 10%) of its mass from accretion. For low luminosity AGN the prediction is more model-dependent, but the timescale for alignment in these systems is sufficiently long that the spin axis of the hole may be able to remain stable if individual accretion episodes are brief. Here the assumed value for $`\nu _2/\nu _1`$ is crucial, and efforts to determine this ratio for specific angular momentum transport mechanisms and disc models would be valuable.

Some powerful radio galaxies display jets that appear to have maintained their direction for long periods of time. Our results imply that this cannot be due to any intrinsic stability imparted by the spin of the black hole, but instead must reflect a long term constancy in the angular momentum of the outer regions of the accretion disc. If the gas derives from a reservoir set up by a single accretion event then this constancy would not be surprising. Alternatively, a preferred axis for gas arriving in the nuclear regions might be the consequence of interactions between inflowing gas and the galactic potential.

We have calculated the torque assuming that the disc is able to relax into and maintain a steady-state configuration as realignment occurs. This assumption is justified in the strongly warped inner region, where $`t_{\nu _1}`$ and $`t_{\nu _2}`$ are indeed smaller than the derived $`t_{\mathrm{align}}`$, but fails at large radii of $`10^4`$ \- $`10^5R_g`$ where our solutions show that the disc maintains a small but significant warp. This will not increase the alignment timescale (before the disc reaches a steady-state the torque in our time dependent calculations is larger than the limiting values we solve for), but does suggest that following alignment of the hole the outer disc might retain a modest warp for a long period. For a $`10^8M_{}`$ hole a radius of $`10^4R_g`$ corresponds to $`0.1\mathrm{pc}`$, a scale which is observable with VLBI in nearby AGN such as NGC4258 (eg. Miyoshi et al. 1995). Some of the warps in maser discs could then simply reflect the slow decay of the initial conditions of a past accretion event, and be devoid of a persistent forcing mechanism. Of course NGC4258 itself has an accretion rate much lower than the models considered here (Gammie, Narayan & Blandford 1999), and self-gravity is likely to be important in its disc (Papaloizou, Terquem & Lin 1998), but the main point — that $`t_{\nu _2}`$ at a few tenths of a parsec is likely to be long — remains true.

A related question is whether the disc at $`R10^2R_g`$ is able to attain a steady state, even in the absence of a warp. This region of the disc may be thermally unstable (Lin & Shields 1986; Clarke 1989; Siemiginowska, Czerny & Kostyunin 1996) and prone to a limit cycle involving large changes in the accretion rate. The interplay between such a cycle and the dynamics of a warped disc would be complex and intriguing, though a more promising setting for exploring such a scenario would be in galactic black hole systems, if those possess warped discs.

For AGN with rapidly spinning holes, the main warped region occurs at radii of $`10^2R_g`$. During realignment episodes there is a substantial deposition of energy, originating in the spin energy of the hole, into this region of the disc. This is particularly the case if $`\nu _2\nu _1`$, as we have argued is likely in this paper and previously (Natarajan & Pringle 1998). The disc in the vertical structure models we have considered is quite thin at these radii and can easily radiate locally a significantly enhanced luminosity, but substantial changes to the vertical structure might be expected (note that we have assumed that the vertical structure of a warped disc is just like that of a planar disc, which is wrong on several counts). Moreover, the enhanced luminosity would increase the disc flux at optical and infrared wavelengths, which would in any event be raised due solely to the greater covering fraction of the warped disc as seen from the central source. Since the alignment timescale is short, systems of this kind would be rare — most AGN would harbour quite flat central discs — but potentially bright.

## Acknowledgements

We gratefully acknowledge useful discussions with Jim Pringle, Mitch Begelman, Roger Blandford, Julian Krolik and Martin Rees, and thank the referee for a very helpful report. PJA thanks Space Telescope Science Institute, where part of this paper was completed, for support and hospitality.

## Appendix

The ODEs for the steady state shape of the disc are most easily solved by expressing equation (1) in terms of separate equations for the tilt vector $`\widehat{l}`$ and surface density $`\mathrm{\Sigma }`$. For a Keplerian potential these are (Pringle 1992),

$`{\displaystyle \frac{\widehat{𝐥}}{t}}=\left[3\nu _1^{}+{\displaystyle \frac{\mathrm{\Sigma }^{}}{\mathrm{\Sigma }}}\left(3\nu _1+{\displaystyle \frac{1}{2}}\nu _2\right)+{\displaystyle \frac{\nu _2}{R}}\left({\displaystyle \frac{3}{4}}+R^2\left|{\displaystyle \frac{\widehat{𝐥}}{R}}\right|^2\right)\right]{\displaystyle \frac{\widehat{𝐥}}{R}}`$   
$`+{\displaystyle \frac{}{R}}\left({\displaystyle \frac{1}{2}}\nu _2{\displaystyle \frac{\widehat{𝐥}}{R}}\right)+{\displaystyle \frac{1}{2}}\nu _2\left|{\displaystyle \frac{\widehat{𝐥}}{R}}\right|^2\widehat{𝐥}+{\displaystyle \frac{\stackrel{}{\omega }_p\times \widehat{𝐥}}{R^3}},`$ (15)   
where the primes denote derivatives with respect to $`R`$, and

$`{\displaystyle \frac{\mathrm{\Sigma }}{t}}`$ $`=`$ $`{\displaystyle \frac{3}{R}}{\displaystyle \frac{}{R}}\left[R^{1/2}{\displaystyle \frac{}{R}}\left(\nu _1\mathrm{\Sigma }R^{1/2}\right)\right]`$ (16)   
$`+`$ $`{\displaystyle \frac{1}{R}}{\displaystyle \frac{}{R}}\left[\nu _2\mathrm{\Sigma }R^2\left|{\displaystyle \frac{\widehat{𝐥}}{R}}\right|^2\right].`$   
Setting the time derivatives to zero these constitute 3 second order equations for $`\mathrm{\Sigma }`$ and any two components of $`\widehat{l}`$. We solve these with the boundary conditions described in Section 2.3 using the NAG routines D02GAF and D02RAF, which implement the finite difference technique with deferred correction described by Pereyra (1979) (for a general introduction to such methods, see e.g. Press et al. 1992). An iterative approach is employed, in which we first solve equation (16) assuming zero tilt, and then use the resulting surface density profile in the solution of equation (15). The solution for the tilt vector $`\widehat{l}`$ is then recycled for use in equation (16), and we loop until convergence is achieved. Typically only a small number of iterations are required, since the changes to the surface density profile for even quite strong warps vary smoothly with changing disc shape. A finely space finite difference mesh is required to obtain a solution using this scheme, we use between 8000 and 16,000 mesh points evenly spaced in $`\mathrm{log}R`$. Even so, we have found that obtaining solutions with this scheme is still rather difficult, especially when $`\beta >0`$ or $`\nu _1\nu _2`$. For these cases, we start with the easy $`\beta =0`$, $`\nu _1=\nu _2`$ solution, and step towards the desired parameters using the previous solution at each step as the initial guess.

We have been unable to obtain solutions using the above method for large values of $`\beta `$. In this regime, we instead evolve directly equation (1) using the numerical code described by Pringle (1992) until a steady-state is obtained. This code uses an explicit first order finite difference technique, which conserves angular momentum to machine precision. We have modified the boundary conditions to be as close as possible to those described in Section 2.3, so that the results are directly comparable between the two solution methods. The time dependent code is vastly more expensive to run, for the runs described here we are restricted to 100 radial grid points, again spaced evenly in $`\mathrm{log}R`$. Nevertheless this still provides reasonable accuracy for computing the torque.