A&A 440, 185-198 (2005)
DOI: 10.1051/0004-6361:20052706

Accretion in dipole magnetic fields: flow structure and X-ray emission of accreting white dwarfs

J. B. G. Canalle1,2 - C. J. Saxton1,3 - K. Wu1 - M. Cropper1 - G. Ramsay1


1 - Mullard Space Science Laboratory, University College London, Holmbury St. Mary, Dorking, Surrey RH5 6NT, UK
2 - State University of Rio de Janeiro, Rua São Francisco Xavier, 524/3023-D, CEP 20559-900, Rio de Janeiro, RJ, Brazil
3 - Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany

Received 16 January 2005 / Accepted 12 March 2005

Abstract
Field-channelled accretion flows occur in a variety of astrophysical objects, including T Tauri stars, magnetic cataclysmic variables and X-ray pulsars. We consider a curvilinear coordinate system and derive a general hydrodynamic formulation for accretion onto stellar objects confined by a stellar dipole magnetic field. The hydrodynamic equations are solved to determine the velocity, density and temperature profiles of the flow. We use accreting magnetic white-dwarf stars as an illustrative example of astrophysical applications. Our calculations show that the compressional heating due to the field geometry is as important as radiative cooling and gravity in determining the structure of the post-shock flow in accreting white-dwarf stars. The generalisation of the formulation to accretion flows channelled by higher-order fields and the applications to other astrophysical systems are discussed.

Key words: accretion, accretion disks - hydrodynamics - shock waves - stars: magnetic fields - stars: novae, cataclysmic variables - stars: pre-main-sequence

1 Introduction

Accretion is a common phenomenon in astrophysical systems ranging from young stellar objects, interacting binaries, galaxies to galaxy clusters. When the magnetic-field stress is larger than the ram pressure of the accreting material, the flow is confined to follow the magnetic-field lines. The accretion hydrodynamics in these systems are therefore dependent on the magnetic-field geometry.

Field-channelled accretion flow can occur in young stellar objects (e.g. Gullbring et al. 2000; Lamzin et al. 2001; Königl 1991; Hartman et al. 1994; Li 1996; Calvet & Gullbring 1998; Koldoba et al. 2002; Stelzer & Schmitt 2004; Lamzin 1998; Romanova et al. 2003), neutron star accretion from inter-stellar medium (e.g. Toropina et al. 2003), and interacting binaries containing a white dwarf or a neutron star (e.g. Koldoba et al. 2002; Li et al. 1996; Arons 1993; Kryukov et al. 2000; Lovelace et al. 1995; Elsner & Lamb 1977; Ghosh & Lamb 1978). In white-dwarf and neutron-star binaries the magnetic field of the compact star can also be strong enough to affect the cooling processes. Despite the fact that these stars may have complex magnetic-field structures, the dipole field component is important, as it has a longer range than the higher-order field components and dominates in the regions sufficiently far from the star.

Here we investigate the accretion flow onto stellar objects in the regime $B^2 \gg 8\pi \rho v^2$ (where B is the magnetic field, and $\rho$ and v are the density and velocity of the accreting material respectively), so that the flow is strictly confined by the magnetic field. We apply a curvilinear coordinate system that is natural to the dipole-field geometry for the hydrodynamic formation. We solve the hydrodynamic equations for an accretion flow in which the cooling function has a power-law dependence on the temperature and density and thereby obtain the flow-velocity, density and temperature profiles.

We apply the model to the post-shock flow in magnetic cataclysmic variables (mCVs), which are close binary systems containing a magnetic white dwarf accreting from a red-dwarf companion star (see Warner 1995). We calculate the temperature, velocity and density structure of the post-shock emission region and model the X-ray line and continuum emission. We compare our results to those obtained by the plane-parallel model (Chevalier & Imamura 1982; Cropper et al. 1999; Wu 1994; Wu et al. 1994) which is generally used in spectral analysis of X-rays emitted from mCVs.

The paper is organised as follows: in Sect. 2 we derive the hydrodynamic equations in curvilinear coordinates; in Sect. 3 we present the treatment of the boundary conditions and discuss briefly the numerical scheme that we use; and in Sects. 4 to 6 we show an illustrative example - accretion onto a magnetic white dwarf stars - and present hydrodynamic structure and spectral calculations. In Sect. 7 we consider the use of this formulation in some other astrophysical applications. Our summary and conclusions comprise Sect. 8.

2 Formulation

2.1 Coordinate systems

For the analytic study of field-channelled accretion flow, we choose a coordinate system with one component along the magnetic-field lines. In this representation, the flow is at most 2-dimensional (2D); and it can often be reduced to 1-dimensional (1D) form when the system has a special geometric symmetry. In the present study we consider accretion channelled by a dipole field, which is axi-symmetric and, it can be shown, essentially 1D. We define a curvilinear coordinate system (u, w, $\varphi $) for the dipolar field as follows. The first coordinate u is defined by the magnetic-field lines generated by a point magnetic dipole. The dipole is oriented in the z direction and located at the center of the accreting star (the origin of the coordinate systems), where (x,y,x) are the usual Cartesian coordinates. As a static magnetic field is the gradient of a potential, we use the equipotential curves of the dipole field as second coordinate w. The third coordinate is the azimuthal coordinate, $\varphi $, the same as that in the conventional spherical coordinate system (r, $\theta$, $\varphi $). The unit vectors $\hat{\vec{u}}$, $\hat{\vec{w}}$ and $\hat{\vec{\varphi}}$ are orthogonal to each other, as shown in Fig. 1.

In spherical coordinates, the field lines and equipotential surfaces of a dipole field satisfy $r = a_1 \sin^2 \theta$ and $r^2 = a_2 \cos \theta$ respectively, where a1 is the radius at which the field line intercepts the mid-plane of the dipole and a2 is the radius at which the equipotential surface meets the polar axis. Thus, we have

 \begin{displaymath}u = \frac{\sin^2 \theta}{r}~{\rm and}
\end{displaymath} (1)


 \begin{displaymath}w = \frac{\cos \theta}{r^2} \cdot
\end{displaymath} (2)

It follows that $\theta = \cos^{-1} (wr^2)$ and the radial location r is a function of u and v, as given in Eq. (A.7) in the appendices. Also,

 \begin{displaymath}\left[{\begin{array}{ccc}
x\ & y
\end{array}}\right]
=
\...
...ray}{ccc}
\cos~\varphi\ & \sin~\varphi
\end{array}}\right] ,
\end{displaymath} (3)


 
z = w r3 . (4)

(See Appendix A for the coordinate transformations and their derivations.)

2.2 Hydrodynamics equations


  \begin{figure}
\par\includegraphics[angle=270,width=8.4cm,clip]{2706f01.ps}
\end{figure} Figure 1: The projection of the ( u, w, $\varphi $) cordinates on the xz-plane. The stellar magnetic dipole is located at the origin and the symmetry axis of the dipole is the z axis. The central circle respresents the surface of the star. The solid lines are the constant-u curves and the dashed lines are the constant-w curves. The unit vectors  $\hat{\vec{u}}$ and  $\hat{\vec{w}}$ are the unit normals to the constant-u and constant-w curves respectively. The $\varphi $ coordinate, whose normal is perpendicular to the xz-plane, is not shown.
Open with DEXTER

The hydrodynamics of the flow is described by the conservation equations

 \begin{displaymath}\frac{{\rm d}}{{\rm d}t}\rho
+(\nabla\cdot\rho\vec{ v})=0 ,
\end{displaymath} (5)


 \begin{displaymath}\frac{{\rm d}\vec{ v}}{{\rm d}t}+
\frac{1}{\rho}\nabla(P_{\rm e}+P_{\rm i})
=\vec{ g}+\vec{ f}_{\rm rad} ,
\end{displaymath} (6)


 \begin{displaymath}\frac{{\rm d}P_{\rm i}}{{\rm d}t}-\gamma\frac{P_{\rm i}}{\rho}
\frac{{\rm d}}{{\rm d}t}\rho
=-(\gamma - 1)\Gamma_{\rm ei} ,
\end{displaymath} (7)


 \begin{displaymath}\frac{{\rm d}}{{\rm d}t}P_{\rm e}
-\gamma\frac{P_{\rm e}}{\r...
...c{ q}+\vec{ F}_{\rm rad})
-\Lambda_{\rm h}-\Gamma_{\rm ei}] ,
\end{displaymath} (8)

where ${\rm d}/{\rm d}t\equiv \partial/\partial t+ \vec{ v}\cdot\nabla$, and $P_{\rm e}$, $P_{\rm i}$, $\vec{v}$ and $\rho$ are the electron pressure, ion pressure, velocity and density. $\gamma$ is the adiabatic index of the gas, $\vec{g}$ is the gravitational acceleration, $\vec{f}_{\rm rad}$ is the acceleration due to the radiative force, $\vec{F}_{\rm rad}$ is the radiative flux, $\Gamma_{\rm ei}$ is the electron-ion exchange rate, and $\Lambda_{\rm h}$ is the heating function.

For stationary flows, the time derivatives are zero. For one temperature flows, the electron and the ion energy equations can be combined by eliminating $\Gamma_{\rm ei}$, yielding

\begin{displaymath}(\vec{ v}.\nabla)P - \frac{\gamma P}{\rho}(\vec{ v}.\nabla)\r...
...nabla\cdot(\vec{ q}+\vec{ F}_{\rm rad})
- \Lambda_{\rm h}] ,
\end{displaymath} (9)

where $P = P_{\rm e} + P_{\rm i}$ is the total pressure. We ignore the details of radiative transfer for simplicity and assume that the total effect of heating and energy loss can be described by a cooling function $\Lambda$ that depends only on the local hydrodynamical variables. We also ignore the radiative force and set $\vec{f}_{\rm rad} = 0$. Under these assumptions, we have a set of three conservation equations

 \begin{displaymath}\nabla \cdot (\rho \vec{ v}) = 0 ,
\end{displaymath} (10)


 \begin{displaymath}(\vec{ v}\cdot \nabla) \vec{ v} + \frac{1}{\rho}\nabla P
= \vec{ g} ,
\end{displaymath} (11)


 \begin{displaymath}(\vec{ v}\cdot \nabla) P
- \frac{\gamma P}{\rho}(\vec{ v}\cdot \nabla )\rho
= -(\gamma -1)\Lambda .
\end{displaymath} (12)

This set of hydrodynamic equations is closed after we specify an equation of state for the gas.

2.3 Dipole-field channelled accretion flow

As we have assummed by default that the flow follows the field lines, the velocity has only the component along the $\hat w$ direction, i.e.

 \begin{displaymath}\vec{ v}= v ~ \hat{\vec w} .
\end{displaymath} (13)

If we ignore the effect of stellar rotation, the acceleration is due to only the gravitational force, which is radial, i.e. $\vec{g}= -g_{*}\hat{\vec r}$, where g*=GM*/R*2, M* and R* and the mass and radius of the accreting star, and G is the gravitational constant. For accretion flow along a dipole stellar magnetic field, the gravitational acceleration has two components, i.e.

 \begin{displaymath}\vec{ g}= g_{u}~ \hat{\vec u} + g_{w}~\hat{\vec w} ,
\end{displaymath} (14)

in which

 \begin{displaymath}g_{u} = -g_{*}~ \frac{Q_{11}}{r(u,w)^{2}} ,
\end{displaymath} (15)


 \begin{displaymath}g_{w} = -g_{*}~ \frac{Q_{12}}{r(u,w)^{2}} ,
\end{displaymath} (16)

where r(u,w) is given by Eq. (A.7). The terms Q11 and Q12 above are the elements of the transformation matrix  ${\mathsf Q}$ (given in Appendix A).

The hydrodynamic equations decomposed into the orthogonal components in the (u, w, $\varphi $)-coordinate system are

 \begin{displaymath}\frac{\partial}{\partial w}(h_{1}h_{3}\rho~v)=0 ,
\end{displaymath} (17)


 \begin{displaymath}\frac{v}{h_{2}}\frac{\partial v}{\partial w} + \frac{1}{h_{2} \rho }
\frac{\partial P}{\partial w} =g_{w} ,
\end{displaymath} (18)


 \begin{displaymath}\frac{v}{h_{2}}\frac{\partial P}{\partial w}
- \frac{\gamma}...
...rho}\frac{\partial \rho}{\partial w}
= - (\gamma -1)\Lambda .
\end{displaymath} (19)

In the equations above h1, h2 and h3 are the metric elements of the (u, w, $\varphi $)-coordinate system, and their explicit expression are derived in Appendix A. As we have assumed that the flow strictly follows the magnetic-field lines (and omitted an explicit consideration of the plasmas and MHD effects, e.g. the tension of the field lines and the magnetic stress), we can ignore the momentum and force components in u and $\varphi $ directions. The three relevant hydrodynamic equations in this study are therefore Eqs. (17)-(19), which are respectively the mass-continuity, momentum-conservation and energy-conservation equations along the magnetic-field lines.

Direct integration of the mass-continuity Eq. (17) yields

 \begin{displaymath}h_{1}~h_{3}~\rho~v = C ,
\end{displaymath} (20)

where C is a constant (related to the local accretion rate and the geometry of the field) to be determined. It follows that

 \begin{displaymath}\frac{\partial}{\partial w} (h_{1} ~ h_{3} ~ \rho ~ v^{2})
=h_{1}~ h_{3} ~ \rho ~ v \frac{\partial v }{\partial w} ,
\end{displaymath} (21)

and the momentum-conservation equation becomes

 \begin{displaymath}\frac{\partial}{\partial w} (h_{1} ~ h_{3} ~ \rho ~ v^{2})
...
...partial P}{\partial w}
=h_{1}~ h_{2} ~ h_{3}~ \rho ~ g_{w} .
\end{displaymath} (22)

To simplify the mathematics we consider a variable $\xi$, as in Cropper et al. (1999), which is defined as

 \begin{displaymath}\xi \equiv v + \frac{P}{\rho ~ v} ,
\end{displaymath} (23)

and use it to replace P as a dependent variable in the hydrodynamic equations. We also consider a geometrical function

\begin{displaymath}{\mathcal H}(u,w) \equiv \frac{\partial}{\partial w}\ \ln~(h_1 h_3) ,
\end{displaymath} (24)

which describes the variation of the cross-section area of the flow-flux tube. (Numerical values of ${\mathcal H}$ for different field lines are shown in Appendix B.) We may then rewrite the momentum-conservation equation as

 \begin{displaymath}\frac{\partial \xi}{\partial w}
= {{g_w h_2}\over{v}} +{\mathcal H}(\xi-v)
.
\end{displaymath} (25)

Using the mass-continuity and momentum-conservation equations in terms of the ${\mathcal H}$ function and the variable $\xi$, we rewrite the energy-conservation equation as

 \begin{displaymath}\frac{\partial v}{\partial w}
={{-h_2}\over{\gamma(\xi-v)-v}...
...da
+{{\gamma{\mathcal H}}\over{h_2}}(\xi-v)v
+g_w
}\right].
\end{displaymath} (26)

We now have a set of coupled differential equations (Eqs. 26 and 25) that determine the accretion flow channelled by a dipole magnetic field. When the explicit form of the cooling function $\Lambda$ is given and the boundary conditions are specified, Eqs. (25) and (26) can be solved along each field line (constant u and $\varphi $) to obtain the velocity, density and temperature profiles.

2.4 Cooling function

In this study, we simply consider an ideal-gas law

 \begin{displaymath}P = \frac{\rho~k_{\rm B}~ T}{ \mu~ m_{\rm H}} ,
\end{displaymath} (27)

for the equation of state of the gas (where T is the thermal temperature of the gas, $\mu$ is the molecular mass of the gas, $m_{\rm H}$ is the Hydrogen atomic mass and $k_{\rm B}$ is the Boltzmann constant) and ignore the microscopic plasma effects due to the magnetic field. In a variety of astrophysical situations, the cooling function of an accretion flow can be approximated by

 \begin{displaymath}\Lambda = A \rho^{\alpha} T^{\beta} ,
\end{displaymath} (28)

where $\alpha$ and $\beta$ are the power-law indices for the density and temperature. When $\alpha = 2$ and $\beta = 1/2$, $\Lambda$ is the cooling function for the optically thin free-free emission, which is an important radiative-loss process in accreting compact stars and in some young stellar objects. From the mass-continuity Eq. (20) and the equation of state for the gas (27), we can express the cooling function in terms of the hydrodynamic variables used in our formulation:

 \begin{displaymath}\Lambda = A ~\left(\frac{C}{h_1~h_3~ v}\right)^{\alpha}
\left[{v(\xi-v)}\right]^\beta .
\end{displaymath} (29)

Expressing the cooling function in terms of the hydrodynamic variables, the energy conservation is
 
                                    $\displaystyle {{{\rm d}v}\over{{\rm d}w}}$ = $\displaystyle {{-h_{2}}\over{\gamma(\xi-v) - v}}
\biggl[\ (\gamma-1) A \biggl( \frac{C}{h_1 h_3}\biggr)^{\alpha-1}
\left(\xi-v \right)^\beta v^{\beta-\alpha}$  
    $\displaystyle +{{\gamma{\mathcal H}}\over{h_2}}v(\xi-v) + g_w\biggr] .$ (30)

3 Transonic accretion onto stellar objects

For non-relativistic, spherical (Bondi-Hoyle type) accretion onto stellar objects, there are two classes of physical solutions (see e.g. Frank et al. 1992). In the first class, the flow is subsonic from infinity to the surface of the accreting star. In the second class, the flow is initially subsonic and becomes supersonic at certain radius from the star. As it requires the flow velocity to be zero at the stellar surface, the flow will become subsonic again via a shock, which converts the kinetic energy of the flow to the thermal energy of accreting material.

The situation is the same for flows channelled by a dipole magnetic field, although the dipole field lines, except the two originating from the poles, are all closed at finite distances. (See e.g. Koldoba et al. 2002, for more details.) In this study we consider the second class and focus on determining the density and temperature structures of the region between the shock and the star. We will show how to construct the shock boundary conditions within the hydrodynamic framework that we use and solve the hydrodynamic equations with these boundary conditions.

   
3.1 Boundary conditions

We apply the Rankine-Hugoniot condition to derive the shock boundary condition, assuming that the pre-shock flow is supersonic and cold, such that the Mach number ${\cal M} \rightarrow \infty$, we have

 \begin{displaymath}V_{\rm s}= \frac{1}{4}V_{\rm pre} ,
\end{displaymath} (31)

(for an adiabatic index $\gamma = 5/3$), where $V_{\rm pre}$ and $V_{\rm s}$ are the pre-shock and post-shock velocities respectively. As a reasonable approximation, we may set $V_{\rm pre} = V_{\rm ff_{s}}$, where $V_{\rm ff_s}$ is the free-fall velocity of the gas before encountering the shock. We denote the free-fall velocity at the stellar surface as

 \begin{displaymath}V_*=\sqrt{\frac{2 G~M_{*}}{R_{*}}} ,
\end{displaymath} (32)

and at any general position (u,w) $V_{\rm ff}=V_*/\sqrt{r}$. Hence, we have the boundary condition

 \begin{displaymath}V_{\rm s}
= {{V_*}\over{4\sqrt{r_{\rm s}}}} \cdot
\end{displaymath} (33)

(Here we use the conventions: the shock is located at $r = r_{\rm s}$ and the stellar surface is located at r=1.)

If the collisional energy exchange between electrons and ions is efficient compared to radiative cooling then the electrons and ions have approximately equal temperatures. The shock temperature is given by

 \begin{displaymath}T_{\rm s} =
\frac{3}{8}\frac{\mu~G~M_{*}~m_{\rm H}}{k_{\rm B}(R_{*}+x_{\rm s})}
\end{displaymath} (34)

(see Frank et al. 1992) where the shock height above the stellar surface  $x_{\rm s} \equiv (r_{\rm s}-1)R_{*}$. Recall the definition of the variable $\xi$, Eq. (23), and we have

 \begin{displaymath}\xi_{\rm s}
=V_{\rm s}+\frac{k_{\rm B}T_{\rm s}}{\mu m_{\rm H}V_{\rm s}} ,
\end{displaymath} (35)

which implies that, regardless of the shock location $r_{\rm s}$,

\begin{displaymath}\xi_{\rm s}= V_{\rm {ff}_{\rm s}} .
\end{displaymath} (36)

At the stellar surface, we have the stationary-wall boundary condition for the velocity, that is v=0. The variable $\xi$ specifies the specific momentum flux, and determining its boundary value at the stellar surface is less trivial, as the flow is not necessarily perpendicular to the stellar surface. Nevertheless, knowing only the value of v at the stellar surface together with the values of v and $\xi$ at the shock is sufficient to solve the hydrodynamic equations.

3.2 Modified hydrodynamics equations and numerical technique

The stellar surface corresponds to r=1 in the spherical coordinates, but the functional form of the stellar surface in the $(u, w, \varphi)$ coordinates is not trivial. At the surface the independent variable w has different values for different magnetic field lines. The boundary value of w for a field line is determined (using Eq. (2)) only after we specify the colatitude on the stellar surface at which the field line is anchored (which we shall denote as $\theta _0$).

As the boundary values of w both at the shock and at the stellar surface are not easily specified in terms of simple functions of w and u, we search for an alternative independent variable. The criterion for the variable is that it increases or decreases montonically in the region of interest. As shown above (Sect. 3.1) the velocity v is better defined at the shock boundary and the stellar-surface boundary. It should be smooth in the post-shock region and is monotonic in at least in the region just beneath the shock and in the region just above the stellar surface.

Without losing qualitative generality, hereafter we will use $\alpha = 2$ and $\beta = 1/2$ (Eq. (29), describing Bremsstrahlung cooling) in our calculations. Then the energy-conservation Eq. (30) can be expressed with v the independent variable and w the dependent variable,

 
                              $\displaystyle \frac{{\rm d} w}{{\rm d} v}$ = $\displaystyle - \frac{1}{h_2}
\big[ \gamma (\xi - v) - v \big] ~v^{3/2}
\biggl[ \frac{(\gamma-1)A~C}{h_{1}h_{3}}\sqrt{\xi-v}$  
    $\displaystyle +\frac{\gamma {\mathcal H}}{h_{2}} (\xi-v)~v^{5/2}
+g_{w} v^{3/2} \biggr]^{-1}
.$ (37)

Applying the chain rule for differentiation, we can also derive (from Eqs. (25)) another equation for the variable $\xi$, that is
 
                             $\displaystyle \frac{{\rm d} \xi}{{\rm d} v}$ = $\displaystyle \left[ h_{2} \frac{\partial w}{\partial v}\right]
\left[ \frac{g_{w}}{v} + \frac{{\mathcal H}}{h_{2}}(\xi -v) \right]$  
  = $\displaystyle -\big[ \gamma(\xi-v)-v \big]
\left[ g_{w} v^{1/2} + \frac{{\mathcal H}}{h_{2}} (\xi-v)v^{3/2} \right]$  
    $\displaystyle \times \biggl[ \frac{(\gamma-1)A~C}{h_{1}h_{3}} \sqrt{\xi-v}
\!+\...
...ac{\gamma {\mathcal H}}{h_{2}} (\xi-v)~v^{5/2}
\!+\!g_{w}~v^{3/2} \biggr]^{-1}.$ (38)

The shock-boundary value of v is still not constant among the field lines (i.e. for different choices of u). Consider a new variable $\tau \equiv{v}/{V_{\rm f\!\!f}}$, where $V_{\rm f\!\!f}(r)$ is the local free-fall velocity at the radial distance r from the centre of the star. For the shock-jump conditions that we assume (Eqs. (31) and (33)), the variable $\tau$ has a fixed value at the shock and is constant at 1/4. Thus, $\tau=0$ at the stellar surface if v=0. The free-fall velocities $V_{\rm f\!\!f}=V_*/\sqrt{r}$, and hence

 \begin{displaymath}v=\frac{\tau~ V_*}{\sqrt{r(u,w)}} ,
\end{displaymath} (39)

where r(u,w) is given in Eq. (A.7).

The corresponding hydrodynamic equations, for ${\rm d}w/{\rm d}\tau$ and  ${\rm d}\xi/{\rm d}\tau$, can be obtained using the chain rule of differentiation with

\begin{displaymath}\frac{{\rm d}\tau}{{\rm d}v} = \tau \biggl[
{\frac{\sqrt{r(...
...u,w)}{r(u,w)}\biggr)
{\frac{{\rm d}w} {{\rm d}v}} \biggr] ,
\end{displaymath} (40)

and r'w(u,w) given in Appendix A.

The use of $\tau$ as an independent variable provides a convenient way to treat the boundary conditions at both the shock and the stellar surface. However $\tau$ (like v) is not guaranteed to be monotonic in w for all systems and all choices of u. Thus we cannot use it as the independent variable throughout the entire post-shock region. In contrast, the coordinate w is always monotonic along a field line, but it is not practically usable at the stellar surface, where the velocity gradient  $\partial v/\partial w$ approaches infinity.

Given that it is straightforward to use w as the independent variable in the numerical integration in the entire flow except at the stellar surface, we consider a hybrid numerical scheme in our calculations. Therefore we use $\tau$ or v as the variable of integration in the vinicity of the stellar surface, $0\leq v\leq \delta v$ with $\delta \la 10^{-3}$, but we switch to integration in terms of w in the upper parts of the post-shock structure. This algorithm allows us to solve the flow profile along each field line.

We calculate the pressure and velocity profiles of the post-shock structure by numerically integrating a set of differential equations between the upper and lower boundaries. In regions where we treat w as the independent variable, we use Eqs. (25) and (26). Where v or $\tau$ serves as the independent variable then we use equations such as Eqs. (37) and (38) or the equivalent with derivatives in $\tau$.

Neither boundary condition is completely specified from first principles: the position of the shock (expressed in terms of $w_{\rm s}$ or $r_{\rm s}$) is not initially known. Nor do we have foreknowledge of the momentum flux at the stellar surface ($\xi_*$). In practice we choose trial values of  $(\xi_*,w_{\rm s})$and test how the resulting profile conforms to the conditions required at the opposite boundary. In one method, we choose trial values of $\xi_*$, and integrate upwards from the stellar surface until a shock is found. We adjust $\xi_*$ until we match the conditions  $\xi_{\rm s} = V_*/\sqrt{r_{\rm s}}$and  $v_{\rm s} = \frac14 V_*/\sqrt{r_{\rm s}}$simultaneously. Alternatively, we may try values of $r_{\rm s}$ and integrate downwards from the shock to the point where v=0. If $r\neq 1$ when v=0 then we adjust the trial value of $r_{\rm s}$and retry. These tests are applied iteratively in a root-finding routine.

4 Accretion onto a magnetic white dwarf

We now apply our formulation to an astrophysical system as an illustration. We consider a simple case: accretion onto a magnetised white dwarf, in which free-free emission is the dominant cooling mechanism. This kind of system can be found in intermediate polars and polars, which are cataclysmic variables containing a magnetic white dwarf accreting material from a low-mass companion star (for reviews of magnetic cataclysmic variables, see e.g. Cropper 1990; and Warner 1995). The post-shock accretion flow in magnetic cataclysmic variables have been investigated by many workers and analytic and semi-analytic results were obtained (e.g. Chevalier & Imamura 1982; Aizu 1973; Imamura et al. 1996; Wu 1994; Saxton & Wu 2001,1999; Wu et al. 1994; Cropper et al. 1999). Their studies address various issues such as the cooling processes and the effect due to gravity, and the calculations either planar or spherical geometries were considered. (For reviews of the hydrodynamics of post-shock accretion in magnetic cataclysmic variables, see Wu 2000; Beuermann 2004.) Because of the assumed geometry, these studies have not quantified the effects of the curvature of field lines in determining the hydrodynamic structure of shock-heated region.

For the application of our formulation in the accretion onto white dwarfs, we simply set $M_{*} = M_{\rm wd}$ and  $R_{*} = R_{\rm wd}$, where  $M_{\rm wd}$ and $R_{\rm wd}$ and white-dwarf mass and radius, in the hydrodynamic equations and scale the variables accordingly. The white-dwarf mass and radius are not independent - when the mass is specified one can calculate the radius by means of a mass-radius relation. In our calculations, we adopt the Nauenberg (1972) mass-radius relation

 \begin{displaymath}\frac{R_{\rm wd}}{R_{\odot}}=
\frac{0.0225}{\mu_{\rm wd}}\fr...
...qrt{1-(M_{\rm wd}/M_{3})^{4/3}}}
{(M_{\rm wd}/M_{3})^{1/3}} ,
\end{displaymath} (41)

where ${M_{3}} = {5.816}~{M_{\odot}}/{\mu_{\rm wd}^{2}}$ is the Chandrasekhar mass limit. We set the electron mean molecular weight parameter $\mu_{\rm wd} = 2.00$. This is a good approximation for white dwarfs with He, C/O or Ne/Mg composition.

We consider a strong shock for the upper boundary condition and a cool stationary wall for the lower boundary condition. The accreting material is an ideal gas with  $\gamma = 5/3$ and approximately solar abundances (Anders & Ebihara 1982). The mean ionic mass and charge are $\overline{m}=2365.8 m_{\rm e}$ and $\overline{Z}=1.0999$ respectively, and $n_{\rm e}/n_{\rm H}=1.209$, $\overline{Z^2}=1.3912$, $\overline{Z^2/m}=6.007\times10^{23}\ {\rm g}^{-1}$. We omit other radiative transport and microscopic effects at the shock and the stellar-surface boundary (see e.g. Wu & Cropper 2001; Saxton & Wu 2001; Imamura et al. 1996).

4.1 Bremsstrahlung-dominated flows

If the white dwarf has a weak magnetic field or sufficiently high accretion rate then the cooling of the post-shock flow is dominated by free-free emission, and Cyclotron cooling is unimportant. The cooling function is then

 \begin{displaymath}\Lambda = \Lambda_{\rm br} = A \rho^{2} \sqrt{\frac{P}{\rho}}
\end{displaymath} (42)

with A being a constant that depends on the composition of the plasma (see Rybicki & Lightman 1979). For a purely Hydrogen plasma, $A=3.97 \times 10^{16} g_{_{\rm B}}$ in c.g.s. units where $g_{_{\rm B}}\sim 1$ is the Gaunt factor. When we adopt $g_{_{\rm B}}=1.25$ and approximately solar abundances we have $A=5.06 \times 10^{16}$.

4.2 Flows with Bremsstrahlung and Cyclotron cooling

A natural extension of the case for accretion onto weakly magnetic white dwarfs is the accretion onto strong-field systems (i.e. polars, see Cropper 1990) the flow is cooled by emitting both optical Cyclotron radiation and free-free X-rays. We can use a composite cooling function to describe the cooling process, i.e.

\begin{displaymath}\Lambda = A \rho^2 \sqrt{\frac{P}{\rho}}\
\left[{\ 1+\epsilon_{\rm s} f(P,\rho,u,w) \ }\right]
\end{displaymath} (43)


 \begin{displaymath}f\equiv{{\Lambda_{\rm cy}}\over{\Lambda_{\rm br}}}=
\left({\...
...17}{40}}
\left({ {B}\over{B_{\rm s}} }\right)^{\frac{57}{20}}
\end{displaymath} (44)

(see Wu et al. 1994; and Appendix C), where $\epsilon_{\rm s}$ is the parameter that determines the relative efficiency of Cyclotron cooling, and the subscript ${\rm s}$ denotes quantities evaluated at the shock. The penultimate factor in Eq. (44) appears because the Cyclotron cooling power depends on the cross-sectional area of the accretion stream, which is proportional to h1 h3. The last factor describes the local magnetic field strength, given by

\begin{displaymath}B(u,w) = { \sqrt{4-3ur}\over{2r^3} }
B_{\rm p}
,
\end{displaymath} (45)

where $B_{\rm p}$ is its value at the pole.

Then, we replace $\rho$ by C/h1 h3 v and $P/\rho$ by $(\xi-v)v$ throughout the f cooling function, and substitute $\Lambda$ into Eq. (26). The resulting hydrodynamic equations can be solved using the numerical scheme described in Sect. 3.

4.3 Simplification to previously published models

We note that the formulation we derive can reduce to the formulations obtained in the previous studies under certain approximations and restricted conditions. Take, for instance, Eqs. (37) and (38). If we set  h1 = h2 = h3 =1, ${\mathcal H}=0$ and gw =0, and fix $\xi$ to be a constant equal to the free-fall velocity at the white-dwarf surface, then Eq. (38) vanishes and Eq. (37) becomes the same as that in Chevalier & Imamura (1982) and Wu (1994) for the planar flows. If we set u=0, w=1/r2, $h_2=-{\rm d}r/{\rm d}w$, $h_1 h_3 = {\frac12}$, and ${\mathcal H}=0$ then we obtain the same set of two equations of Cropper et al. (1999).

Comparing the formulation of Cropper et al. (1999) with Chevalier & Imamura (1982) and Wu (1994) reveals that the former needs two equations to describe the flow while the later requires only one. It is because in the presence of gravity the quantity $\xi$, which describes the specific momentum flux, is no longer conserved along the field lines. In our formulation, despite the presence of a dipole field, two differential equations are sufficient to determine the hydrodynamics, as in the case of Cropper et al. (1999). We still have two differential equations when we set gw=0 in both Eqs. (37) and (38).

To further illustrate the importance of the effects due to the field geometry, in the next section we present the numerical results of calculations, in which typical parameters of magnetic cataclysmic variables are used. We will show that the results can differ substantially from those obtained by formulations (e.g. Cropper et al. 1999) without taking account of the geometric effects.

   
5 Hydrodynamic structure


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f02.ps}
\end{figure} Figure 2: Variation of the local mass accretion rate $\rho v$ within the post-shock structure. Different cases of the accretion colatitude, $\theta _0 = 0^\circ , 18^\circ , 36^\circ $, are marked by full, dashed and dotted lines respectively. At the white dwarf surface we have $\rho_* v_* \equiv \dot{m}$ by definition. In this illustrative case $\dot{m}=2.0\ {\rm g}~{{\rm cm}}^{-2}~{\rm s}^{-1}$, B*=0. and $M_*=1.0~M_\odot$. The horizontal axis is the altitude above the white dwarf surface, relative to that of the shock, $r_{\rm s}-1$.
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f03.ps}
\end{figure} Figure 3: Variation of the component of the gravitational acceleration tangential to the magnetic field line. The five curves correspond to cases of hot-spot colatitude $\theta _0=0^\circ , 6^\circ , 12^\circ , 18^\circ , 24^\circ $. The parameters  $(M_*,B_*,\dot{m})$ take the values used in Fig. 2.
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f04.ps}
\end{figure} Figure 4: Post-shock pressure profiles with the parameters ( $M_*,B_*,\dot{m})$ used in Fig. 2. The vertical axis is scaled in terms of the free-fall velocity ( $V_{\rm ff_s}$) and density ( $\rho _{\rm s}$) immediately downstream of the shock. The lower curve (CWRK) is the cylindrical accretion model. The upper curve is the dipole model with $\theta _0=0^\circ $. Varying $\theta _0$ has little effect on the appearance of this curve.
Open with DEXTER

The constant C introduced in (20) is obtained defining the value of $\rho v$ at the bottom of the accretion column, $\rho_* v_*$. Following Cropper et al. (1999) we identify the accretion rate $\dot{m}=\rho_* v_*$, but in the present work the product $\rho v$ is not constant in the accretion column, since the area of the funnel, $\propto h_1h_3$, varies throughout the accretion stream. In our standard illustrative cases, we use $\dot{m}=2.0\ {\rm g}~{{\rm cm}}^{-2}~{\rm s}^{-1}$, a white dwarf mass of  $1.0~M_\odot$, escape velocity $V_*=6.9\times10^{8}~{{\rm cm}}~{\rm s}^{-1}$ and the accretion area is taken to be $3.9\times10^{15}~{{\rm cm}}^2$, exactly 10-3 of the global surface area. Figure 2 shows how the product $\rho v$ decreases with height within the column for several values of accretion colatitude, $\theta _0$.

The gravitational acceleration component in the direction tangential to the magnetic field lines, which we call gw in Eq. (16), is a fraction of the total gravitational acceleration g. Figure 3 shows the relation gw/g. Except for the special case of accretion onto the pole, $\theta _0=0^\circ $, the magnetic field lines are nowhere vertical to the stellar surface. Thus the component gw is not equal to the total acceleration on the surface of the white dwarf.

Scaled relative to the shock height and parameters of the accretion, there are several qualitative differences between the post-shock flows in the present dipolar model and the cylindrical accretion model of Cropper et al. (1999).

The distribution of pressure is shown in Fig. 4. The dipolar accretion model results in proportionally higher pressures throughout the post-shock region.

Corresponding profiles of the density structures are shown in Fig. 5. The dipolar funnel results in greater densities throughout most of the post-shock region (relative to the gas density at the shock), than in the cylindrical accretion case. As in previous models with a power-law expression for Bremsstrahlung cooling (e.g. Chevalier & Imamura 1982; Cropper et al. 1999; Wu 1994), the density necessarily increases asymptotically near the stellar surface.

Figure 6 shows the gas velocity between the white dwarf surface (r=1) and the shock. The velocity profiles show the most significant qualitative difference between dipolar and cylindrical accretion models: in the Bremsstrahlung-dominated case that we have illustrated, the cylindrical model results in a much more constant velocity gradient throughout most of the post-shock region.

Figure 6 also displays a significant quantitative difference between the predicted shock positions in the cylindrical and dipolar models. For the cases shown, with accretion colatitudes $\theta_0\leq 18^\circ$, the shock occurs at a radius around $r_{\rm s}\approx 1.19$. In the equivalent cylindrical accretion model, the shock height is much lower, $r_{\rm s}\approx1.105$, i.e. just over half the height that we would obtain by using the dipole field geometry.

As the magnetic field strength is increased to some tens of ${\rm MG}$, the radiative cooling due to Cyclotron emission becomes large compared to the Bremsstrahlung cooling. As $\epsilon_{\rm s}$ increases, the shock height decreases in both the cylindrical accretion model (Cropper et al. 1999) and our dipole model. For any set of system parameters, $(M_*,\dot{m},\theta_0,\epsilon_{\rm s})$, the dipolar model predicts the shock to occur higher than in the cylindrical model. However the difference is insignificant in high- $\epsilon_{\rm s}$ cases where the shock height is sufficiently low and curvature effects are negligible. For instances in cases with $B_*=10,30,50~{\rm MG}$ we calculate $r_{\rm s}\approx 1.075, 1.016, 1.007$ according to the cylindrical model, but the dipolar model yields $r_{\rm s}\approx 1.100, 1.017, 1.007$ when $\theta _0=0^\circ $.

Figure 7 shows the dependence of shock height upon the colatitude of the accretion spot ($\theta _0$) for several cases of B* appropriate for polars. The maximum shock height generally occurs for accretion along the field line with $\theta _0=0^\circ $. In cases with low B* the variation of $r_{\rm s}$ with $\theta _0$ is less than in the stronger field cases. When the magnetic field is stronger the shock height is smaller compared to the white dwarf radius, and the dipolar funnel effects less significant.

Figure 8 shows equivalent relations between shock height and $\theta _0$, but for a greater accretion rate, $\dot{m}=5.0\ {\rm g}~{{\rm cm}}^{-2}~{\rm s}^{-1}$. Bremsstrahlung cooling is more efficient in this denser flow, and the shock heights are lower than in Fig. 7, as expected from previous studies (e.g. Imamura et al. 1996; Aizu 1973; Wu 1994). Again, as in Fig. 7, $r_{\rm s}$ is more variable in $\theta _0$ in the cases with weaker magnetic field and greater shock height.

Figure 9 shows the temperature as a function of height within the column in our standard illustrative case. For comparison we calculate the temperature profile with the same parameters for the white dwarf and accretion rate, but applied in the cylindrical accretion model. In this specific case the shock temperature in the cylindrical accretion model is  $4~{\rm keV}$ greater than it is in the dipolar accretion model. The dipolar accretion model generally yields a lower shock temperature than cylindrical accretion, because of the greater shock height and thus lower free-fall velocity at the shock. The temperature predictions of cylindrical and dipolar accretion models are much more alike for cases where the shock is closer to the white dwarf surface (and dipolar funnel effects are less significant).

Figure 10 compares the cylindrical and dipolar models' results for the temperature structure of the accretion column, for four different values of the white dwarf mass, in the Bremsstrahlung-dominated limit. The dipolar model generally has a lower shock temperature than the equivalent cylindrical model. The difference between models is greatest for the most massive white dwarfs, as their accretion shocks are higher (and thus more affected by funnel geometry).

In the next section we present calculations of the X-ray spectrum emitted from the post-shock flow, for two represenative cases of M*. We compare the predictions of the dipolar accretion model with those of the planar models, thereby showing the precise effects of curvature.


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f05.ps}
\end{figure} Figure 5: Density profiles of the accretion columns shown in Fig. 4, scaled to the gas density immediately downstream of the shock ( $\rho _{\rm s}$). In the upper curve $\theta _0=0^\circ $, and variation of $\theta _0$ has little effect. The cylindrical accretion model (CWRK) is the lower (dot-dot-dot-dash) curve.
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f06.ps}
\end{figure} Figure 6: Velocity profile of the post-shock structure, normalised to the escape velocity and the stellar surface, V*. The upper curve represents the model of Cropper et al. (1999). The lower curve shows the effects of accretion in a dipole-field accretion funnel, with $\theta _0=0^\circ $. The system parameters are the same as in Fig. 2.
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=8cm,clip]{2706f07.ps}
\end{figure} Figure 7: Effect of the accretion colatitude, $\theta _0$, upon the radial location of the shock, $r_{\rm s}$ (in units of stellar radius), for a white dwarf of mass $M_*=1.0~M_\odot$and $\dot{m}=2.0\ {\rm g}~{{\rm cm}}^{-2}~{\rm s}^{-1}$. Results are shown for cases with magnetic field strength $B_* = 0, 10, 30, 50~{\rm MG}$ at the accretion spot.
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=8cm,clip]{2706f08.ps}
\end{figure} Figure 8: Radial location of the shock as a function of $\theta _0$, with the same conditions as in Fig. 7, except the accretion rate is higher: $\dot{m}=5.0\ {\rm g}~{{\rm cm}}^{-2}~{\rm s}^{-1}$.
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f09.ps}
\end{figure} Figure 9: Variation of the gas temperature in the accretion column as a function of altitude scaled to the shock location. System parameters $(M_*,B_*,\dot{m})$are the same as in Fig. 2. The solid curve is the case of $\theta _0=0^\circ $; results for $\theta _0\le 18^\circ $ appear indistinguishable on this plot. The nearby dotted curve is the extreme case of $\theta =60^\circ $. The upper curve (dot-dot-dot-dashed) portrays results from the cylindrical accretion model (Cropper et al. 1999).
Open with DEXTER


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f10.ps}
\end{figure} Figure 10: Comparison between cylindrical and dipolar geometry for four different values of white dwarf mass. From bottom to top each pair of curves represents: M* = 0.7, 0.8, 0.9 and $1.0~M_\odot$. In each pair the upper (dotted) curve represents the cylindrical model (CWRK) and the lower (full line) curve represents the dipolar model. We take $\theta _0=18^\circ $.
Open with DEXTER

5.1 Comments on curvature effects upon the flow

In Sect. 5 we use an accreting magnetic white-dwarf star as an illustration and demonstrate that the role of the field geometry can be as important as (or even more than) the role of gravity in determining the velocity, density and temperature profiles of the flow. We can elaborate this more clearly by inspecting the meanings of the terms in the hydrdodynamic equations. Each of the terms in the square brackets in the right side of Eq. (26) corresponds to an energy-transport process. The first represents the radiative loss, the second is determined by the geometry of flow, and the third is due to an external force field, which is only gravity here. A characteristic of the second term is the ${\mathcal H}$ function. In the curvilinear coordinates that we use, h1 and h3 are the metric elements corresponding to the coordinates perpendicular to w. Thus ${\mathcal H} (= \partial \ln (h_1 h_3)/\partial w)$ is the measure of the changes in cross-section area of the magnetic-flux tube along the flow, and it plays an important role in determining the efficiency of compressional heating.

It is worth noting that the ${\mathcal H}$ function also equals to  $h_2 \nabla \cdot {\hat b}$, where ${\hat b}$ is the unit vector tangential to the magnetic field, and h2 is the metric for the coordinate component along the field line. We may also use this relation to determine the ${\mathcal H}$ function.

For accretion onto stellar objects, $g_{w} \sim \bar{V}_{\rm ff}^2$ and $\xi \sim \bar{V}_{\rm ff}$, where $\bar{V}_{\rm ff}$ is the mean free-fall velocity. The metric element $h_1 \sim {\cal{O}}(1)$, and hence, the geometry function ${\mathcal H} \sim {\cal{O}}(1)$. From Eq. (38), we can deduce that the effects due to field curvature (geometry) and gravity are often comparable. Generally, the geometry effect is more important for regions close to the stellar surface than for regions further away, because the convergence of the dipole field lines.

We note that Eqs. (25) and (26) are general equations, in principle, applicable for accretion flows channelled by any field geometries. For example, if the magnetic field has planar parallel structure, then we have h1 = h2 = h3 =1. It follows that ${\mathcal H}=0$. In this flow, the accreting material will not be compressed by the magnetic field. If the flow is spherical, then h1=1, h2 = r, $h_3 = r \sin \theta$ (see e.g. Arfken & Weger 2001), and ${\mathcal H} = 2/r$. The flows channelled by dipole field have been shown in the sections above. The generalisation to channelled flows in higher-order fields is possible. What one needs are first to find the functional expressions of the field lines and the equipotential surfaces, and use them to define the coordinates. After the coordinate system is specified, we can derive the metric elements h1, h2 and h3 and determine the ${\mathcal H}$ function. The effects of more complex field geometries depend on the alignment and distribution of the higher-order field components. If the field lines near the accretion spot diverge more rapidly than in the dipolar model (i.e. a more flared accretion funnel) then we expect an accentuation of the effects found in the present work: a higher shock position and lower shock temperature. If the field lines near the accretion spot diverge less rapidly then the accretion flow will be a closer approximation to the older, planar accretion models.

6 X-ray emission

6.1 Spectral calculations


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f11.ps}
\end{figure} Figure 11: Comparison of X-ray spectra calculated using (I) the dipole accretion model with gravity; (II) CWRK accretion model with gravity but without curvature effects; (III) planar accretion model of Wu (1994) without gravity or curvature; (IV) isothermal model using the post-shock temperature and density derived from model (I). The upper panel shows the X-ray spectrum of model (I) directly. The remaining panels present ratios between the dipole model spectrum and the other models. In all cases we set $\dot{m}=2.0\ {\rm{g}}\ {{\rm cm}}^{-2}\ {\rm{s}}^{-1},
B_*=10\ {\rm MG},
M_*=1.0~M_\odot,
\theta_0=0^\circ$.
Open with DEXTER

We use the analytic model to provide electron density and temperature values throughout the post-shock flow. The structure is divided into between $2\times10^3$ and $2\times10^4$ strata, depending on the adaptive integration steps chosen by the routine that solves the set of differential equations. Using these $(n_{\rm e}, T_{\rm e})$ values, an XSPEC subroutine for a MEKAL model of an optically thin plasma (Phillips et al. 1999; Mewe et al. 1995) calculates the X-ray emission spectrum of each layer. The total X-ray spectrum of the accretion column is a sum over all the layers' spectra, weighted by cell volume. We omit the few cells for which $n_{\rm e}>1.0\times10^{18}~{{\rm cm}}^{-3}$, for which the MEKAL model is inapplicable and fails. In reality the gas at these densities is optically thick and merges into the stellar atmosphere.

Figure 11 presents the X-ray spectrum calculated for a high-mass white dwarf ( $1.0~\hbox{$M_\odot$ }$) with a magnetic field strength of  $10~{\rm {MG}}$at the accretion hot-spot. The uppermost panel represents the dipolar accretion model (I). The lower panels are ratio plots, comparing the dipolar model with: (II) the model of Cropper et al. (1999) (planar geometry, including gravity); (III) the model of Wu (1994) (planar, but without gravity); and (IV) an isothermal, homogeneous post-shock model. For photon energies ranging from 0.2 to  $10~{\rm {keV}}$, the spectrum of model (II) is harder than that of the dipolar model (I), and has less line emission. It may be unsurprising that the CWRK model has a harder spectrum, since the shock temperature is higher than in the dipolar accretion model. The planar accretion model without gravity (III) produces a softer spectrum than the dipolar accretion model, and more line emission. The homogeneous shock model (IV) yields a spectrum that is qualitatively very different from any of the inhomogeneous shock models. It shows a much harder spectal slope.

Curvature effects are strongest for those cases where the shock height is significant compared to the stellar radius, and thus the differences between spectra in comparable CWRK and dipolar accretion models is most noticeable for the cases of the most massive white dwarfs. For lower mass white dwarfs the shock height is small compared to the stellar radius, curvature effects are minor, and the dipolar and CWRK models predict similar X-ray spectra. For the case of $M_*=0.5~M_\odot$(Fig. 12) the spectra predicted from these models are almost indistinguishable.


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f12.ps}
\end{figure} Figure 12: Comparison of X-ray spectra calculated according to the four models, labelled as in Fig. 11 but with stellar mass $M_*=0.5~M_\odot$.
Open with DEXTER

6.2 Effects on X-ray spectrum

For white dwarfs of low or medium mass (e.g. $0.5~M_\odot$) the shock temperatures, post-shock structures and X-ray spectra are nearly the same in the dipolar accretion model as in the earlier planar accretion model of Cropper et al. (1999). For white dwarfs of higher mass (e.g. $1.0~M_\odot$), the dipolar accretion model predicts greater shock heights, lower shock temperatures and softer X-ray spectra in the band from $0.2{-}1.0~{\rm {keV}}$. This implies that for a given observed spectrum, the mass estimate M* is necessarily greater when dipole field curvature is taken into account.

The lower shock temperature does not at first sight provide a route to resolving the high white dwarf masses derived using spectral fits to X-ray data, and which are the subject of some controversy (e.g. Schwope et al. 2002), since even higher mass white dwarfs are required to provide sufficient hard X-ray flux to fit the observed spectrum.

7 Application to other systems

The formulation that we derived is not restricted to accretion in magnetic white dwarfs. It is applicable to a variety of astrophysical systems, provided that the flow is strictly confined by the magnetic field lines, and the radiative heating and cooling processes can be parametrised in terms of the hydrodynamic variables.

An example is accretion in young stellar objects. There are strong observational evidence that accretion flow in T Tauri stars are channelled by the stellar magnetic fields (e.g. Basri et al. 1992). The main differences between the mCV case and the accretion in T Tau stars are: that the cooling process is dominated by line cooling instead of free-free and cyclotron cooling. Moreover, irradiative heating of both pre- and post-shock flow is important in determining the flow hydrodynamics (see e.g. Martin 1996). Under these conditions, we need to relax the assumption of a completely dissociated gas and generalise the hydrodynamical equations to multi-ion and electron flow. In a fully self-consistent treatment the equations of hydrodynamics need to be solved simultaneously with the ionization structure equations (e.g. Lamzin 1998).

Another example of application is the accretion onto a slowly rotating neutron star in a binary system (i.e. X-ray pulsars). The accretion flow is magnetically funnelled, from either the wind of the companion star or an accretion disk (e.g. Koldoba et al. 2002; Lovelace et al. 1995; Elsner & Lamb 1977; Ghosh & Lamb 1978). In cases of disk accretion, the gas is funnelled from the inner edge of the disk, and the radius of this location effectively determines the colatitude of the accretion hot-spot on the stellar surface ($\theta _0$ in our formulation). For details of interaction of the channelled flow and the accretion disk in X-ray pulsars, see e.g. Anzer & Boerner (1983); Li et al. (1996); Spruit & Taam (1990).

8 Summary

We have investigated accretion onto stellar systems in which the flow is strictly channelled by a dipole magnetic field. Such flow occurs in regions close to the stellar surface, where $B \gg 8\pi \rho v^2$. We derive a set of hydrodynamic equations using the curvilinear coordinates natural to the field geometry. The equations are solved to determine the flow velocity, density and temperature profiles. We show that the dipole-field geometry can cause significant compressional heating, and this effect can be comparable to radiative cooling and gravity in determining the structures of the flow near the surface of the accreting star.

The formulation that we have derived is applicable to a variety of astrophysical systems, from white-dwarf stars in magnetic cataclysmic variables to young stellar objects. It would be possible to generalise the formulation for dipole-field channelled flow to the flow channelled by a higher-order field using curvilinear coordinates.

We demonstrate that our analytic model can be efficiently interfaced with spectral software to reproduce fits in analysis of high quality X-ray spectra, obtained by satellites such as Chandra and XMM-Newton. This is a next-generation model for fitting accretion parameters, including white dwarf mass, magnetic field strength, accretion colatitude and mass flux, to specific observed systems.

Acknowledgements

J.B.G.C. is very grateful for the great hospitality of all people of the Mullard Space Science Laboratory (MSSL) and particularly to Gavin Ramsay and Mark Cropper for the invitation to visit MSSL for one year. J.B.G.C. also thanks Roberto Soria, Mat Page and Christian Bridge for discussions. J.B.G.C. acknowledges the support from the Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq) of Brazil and from the State University of Rio de Janeiro (UERJ) for the one-year leave inside the PROCAD program.

Appendix A: Orthogonal basis and metric elements

Let $\hat{\vec{\imath}}$, $\hat{\vec{\jmath}}$ and $\hat{\vec{k}}$ be the standard orthogonal unit vectors in the Cartesian coordinates. The unit vectors of our dipolar curvilinear coordinate system ( $\hat{\vec{u}}, \hat{\vec{w}}, \hat{\vec{\varphi}}$) and ( $\hat{\vec{\imath}}$, $\hat{\vec{\jmath}}$, $\hat{\vec{k}}$) are related to the Cartesian unit vectors via

\begin{displaymath}\left[\begin{array}{c}
\hat{\vec u}\\
\hat{\vec w}\\
\ha...
...\\
\hat{\vec{\jmath}}\\
\hat{\vec{k}}
\end{array}\right],
\end{displaymath} (A.1)

where (see e.g. Arfken & Weger 2001) the transformation matrix is

 \begin{displaymath}{\mathsf U}\!=\!
\left[\begin{array}{ccc}
\frac{\displaystyle...
...tial z}{\displaystyle \partial \varphi}
\end{array}\right] .
\end{displaymath} (A.2)

The coefficients h1, h2 and h3 are the metric of the curvilinear coordinates system, and they are

 \begin{displaymath}h_{[1,2,3]} = \sqrt{ \left(\frac{\partial x}{\partial [u,w,\v...
...t(\frac{\partial z}{\partial [u,w,\varphi]}\right)^{2} } \cdot
\end{displaymath} (A.3)

The transformation from between dipolar coordinates and spherical coordinates ( $\hat{\vec{r}}, \hat{\vec{\theta}}, \hat{\vec{\varphi}}$) is given by:

\begin{displaymath}\left[\begin{array}{c}
\hat{\vec{r}}\\
\hat{\vec{\theta}}\...
...
\hat{\vec{w}}\\
\hat{\vec{\varphi}}
\end{array} \right] ,
\end{displaymath} (A.4)

where the preceding matrix is

 \begin{displaymath}{\mathsf N}=\left[\begin{array}{ccc}
\sin~\theta ~~\cos~\var...
...theta\\
-\sin~\varphi&\cos~\varphi & 0
\end{array}\right] ,
\end{displaymath} (A.5)

${\mathsf Q} = {\mathsf N} \cdot {\mathsf M}$, and ${\mathsf M}\equiv {\mathsf U}^{t}={\mathsf U}^{-1}$ (the transpose of ${\mathsf U}$). Thus ${\mathsf Q}$ completely prescribes the relationship between ( $\hat{\vec{u}}, \hat{\vec{w}}, \hat{\vec{\varphi}}$) and ( $\hat{\vec{r}}, \hat{\vec{\theta}}, \hat{\vec{\varphi}}$).

With u and w defined according to Eqs. (1) and (2) we have $\sin^2\theta = ur$ and $\cos^2\theta =w^2 r^4$, and as $\cos^{2}~\theta + \sin^{2}~\theta = 1$, we therefore have

 
w2r4 + u r -1 = 0 . (A.6)

This is a quartic equation, which has two real and symmetrical roots and two imaginary roots (e.g. Abramowitz & Stegun 1964). The real root that is of physical relevance is

 \begin{displaymath}r=\sqrt{Y-S_+-S_-}-{u\over{4w^2Y}}
\end{displaymath} (A.7)

where

\begin{displaymath}Y=\sqrt{(S_++S_-)^2+3(S_+-S_-)^2}
\end{displaymath} (A.8)


\begin{displaymath}S_\pm \equiv\left({ -{\textstyle\frac{1}{2}}{W}\pm\sqrt{X} }\right)^{1/3}
,
\end{displaymath} (A.9)


\begin{displaymath}W=-{{u^2}\over{64w^4}}
\ \ \ \ \ {\rm and}
\end{displaymath} (A.10)


\begin{displaymath}X
=\left({1\over{128w^4}}\right)^2
\left({
u^4 + {{256w^2}\over{27}}
}\right)
\cdot
\end{displaymath} (A.11)

The first derivatives of the radius function, $r'_{u}\equiv\partial r/\partial u$and $r'_{w}\equiv\partial r/\partial w$, are obtainable by implicit differentiation of Eq. (A.6):

 \begin{displaymath}r'_{u}=
{{-r}\over{4w^2r^3+u}}
\ {\rm\ \ and}
\end{displaymath} (A.12)


 \begin{displaymath}r'_{w}=
{{-2wr^4}\over{4w^2r^3+u}}
\cdot
\end{displaymath} (A.13)

The relation between (x, y, z) and (u, w, $\varphi $) is given by Eqs. (3) and (4). It follows that

 \begin{displaymath}{\partial\over{\partial u}}
\left[{\begin{array}{ccc}
x\ & ...
...ray}{ccc}
\cos~\varphi\ & \sin~\varphi
\end{array}}\right] ,
\end{displaymath} (A.14)


 \begin{displaymath}\frac{\partial z}{\partial u}
= 3~ w~ r^{2} r'_{u} {\rm ,}
\end{displaymath} (A.15)


 \begin{displaymath}{\partial\over{\partial w}}
\left[{\begin{array}{ccc}
x\ & ...
...ray}{ccc}
\cos~\varphi\ & \sin~\varphi
\end{array}}\right] ,
\end{displaymath} (A.16)


 \begin{displaymath}\frac{\partial z}{\partial w}
= ~ r^{2}(r + 3~w~r'_{w}) {\rm ,}
\end{displaymath} (A.17)


 \begin{displaymath}{\partial\over{\partial\varphi}}
\left[{\begin{array}{ccc}
...
...}{ccc}
-\sin\varphi\ &\cos\varphi\ &0
\end{array}}\right]
.
\end{displaymath} (A.18)

Consider these substitions into Eqs. (A.3): (i) Eqs. (A.14) and (A.15) into h1; (ii) Eqs. (A.16) and (A.17) into h2; and (iii) Eq. (A.18) for h3; then we obtain metric elements explicitly:

 \begin{displaymath}h_{1} = \sqrt{
\frac{r}{4~u}(r + 3~ u~ r'_{u})^{2}
+(3~w~r^{2}~r'_{u})^{2}} {\rm ,}
\end{displaymath} (A.19)


 \begin{displaymath}h_{2} = \sqrt{r^{4}~(r+3~w~r'_{w})^{2}
+{{\textstyle\frac{9}{4}}} u~r r_{w}^{'2}
} {\rm ,}
\end{displaymath} (A.20)


 \begin{displaymath}h_{3} = \sqrt{u~r^{3}} {\rm .}
\end{displaymath} (A.21)

It is also worth noting that $\hat{\vec{u}}, \hat{\vec{w}}$ and $\hat{\vec{\varphi}}$ vary along the corresponding coordinate axes. It can be shown that

 \begin{displaymath}\frac{\partial \hat{\vec u}}{\partial u}
= -\left(\frac{\hat...
...{3}}~
\frac{\partial h_{1}}{\partial \varphi} \right) {\rm ,}
\end{displaymath} (A.22)


 \begin{displaymath}\frac{\partial \hat{\vec w}}{\partial w}
= -\left(\frac{\ha...
...{3}}~
\frac{\partial h_{2}}{\partial \varphi} \right) {\rm ,}
\end{displaymath} (A.23)


 \begin{displaymath}\frac{\partial \hat{\vec{\varphi}}}{\partial \varphi}
= -\le...
...at{\vec w}}{h_{2}}~
\frac{\partial h_{3}}{\partial w} \right)
\end{displaymath} (A.24)

(see e.g. Arfken & Weger 2001).

Appendix B: The ${\mathcal H}$ function

The function ${\mathcal H}(u,w)$ is a purely geometrical function, dependent on the curvature of the coordinate system. In the (u, w, $\varphi $) coordinate systems that we use, it is

 \begin{displaymath}{\mathcal H}= \frac{3~u~r^{5}}{(h_{1}h_{3})^2} \left[
\frac...
..._{7}
+ 3~ \left( w~ r~ r'_{u}\right)^{2} t_{8}\right] {\rm ,}
\end{displaymath} (B.1)

where

 \begin{displaymath}t_{7} = \frac{1}{u} + \frac{2~ r'_{u}}{r}
+ \frac{r''_{uw}}{r'_{w}} {\rm ,}
\end{displaymath} (B.2)


 \begin{displaymath}t_{8 }= \frac{1}{w} + \frac{7}{2}\frac{r'_{w}}{r}
+ \frac{r''_{uw}}{r'_{u}}\cdot
\end{displaymath} (B.3)

In the above expression, $r''_{uw}=\partial r'_{u}/\partial w
= \partial(\partial r/\partial u)/\partial w$, which is obtained by implicit differentiation of Eq. (A.12):
 
r''uw = $\displaystyle -{{\left[{
r'_w + \left({ 8wr^3+12w^2r^2 r'_w }\right) r'_u
}\right]
}\over{ 4w^2r^3 + u }}
\cdot$ (B.4)

We can also relate the change in $\rho~v$ to the ${\mathcal H}(u,w)$ function. The mass-continuity equation implies

\begin{displaymath}h_{1}h_{3}\frac{\partial}{\partial{w}}(\rho~v)+
\rho~v\frac{\partial}{ \partial{w} }( h_{1}h_{3} ) = 0 .
\end{displaymath} (B.5)

Hence,

\begin{displaymath}\frac{\partial}{\partial{w}}{\rm ln}(\rho~v) = -{\mathcal H} .
\end{displaymath} (B.6)

Figure (B.1) shows the ${\mathcal H}(u,w)$ function for some values of $\theta _0$, the angle between the axis of the dipole and the position vector of the foot point of the magnetic-field line.


  \begin{figure}
\par\includegraphics[width=8.4cm,clip]{2706f13.ps}
\end{figure} Figure B.1: From top to bottom, the ${\mathcal H}(u,w)$ function for $\theta _0 = 0^\circ , 6^\circ , 12^\circ , 18^\circ , 24^\circ , 30^\circ $ and $36^\circ $ respectively. The horizontal axis is the radial distance from the centre of the star. In the calculation $\dot{m} = 2.0~ {\rm g}~{{\rm cm}}^{-2}~{\rm s}^{-1}$ and $M_*=1.0~M_\odot$.

   
Appendix C: Cyclotron cooling efficiency $\epsilon_{\mathsf{s}}$

In our illustrative system the Cyclotron emission from the post-shock region is optically thick at frequencies up to some cut-off, $\omega_*$(e.g. Langer et al. 1982; Chanmugam & Wagner 1979; Wada et al. 1980). As derived in Saxton (1999) and the appendices of Cropper et al. (1999), the effective volumetric cooling function is

\begin{displaymath}\Lambda_{\rm cy}
={{\pi D}\over{a}} {{k_{\rm B}T}\over{12\pi^2c^2}} \omega_*^3
,
\end{displaymath} (C.1)

where $\omega_*\approx9.87\omega_{\rm c}
(\Theta/10^7)^{0.05}(T/10^8~{\rm K})^{0.5}$, $\omega_{\rm c}=eB/m_{\rm e}c$ is the Cyclotron frequency, $\Theta=2\pi e n_{\rm e} D / B$, D is the diameter of the accretion flow, $a=\pi D^2/4$ is the cross-sectional area, B is the magnetic field strength, T and $n_{\rm e}$ are the electron number density and temperature respectively. Calculating $\Lambda_{\rm cy}$ in the conditions immediately downstream of the shock, and taking its ratio with the Bremsstrahlung cooling (42) gives the relative efficiency parameter $\epsilon_{\rm s}\equiv \Lambda_{\rm cy,s}/\Lambda_{\rm br,s}$ (Wu 1994). An appropriate series of substitutions leads to Eq. (44).

However there is a practical problem to performing calculations starting from a chosen value of $\epsilon_{\rm s}$, since Eq. (44) depends on functions of the shock position, which is initially unknown. In the earliest cylindrical accretion models (Chevalier & Imamura 1982; Aizu 1973; Wu 1994) it was possible to eliminate the shock height by normalisation, but this is impossible for formulations that include explicitly position-dependent effects such as gravity (Cropper et al. 1999). For computational convenience, we modify the definition by replacing the ${\rm s}$-subscripted quantities in Eq. (44) with reference values at the stellar surface (B*, h1*, h3*), and ideal post-shock pressure and densities  $(P'_{\rm s},\rho'_{\rm s})$ of the zero-gravity, cylindrical model of Wu (1994). In that model the shock was assumed to occur close to the stellar surface. The corresponding efficiency parameter is denoted  $\epsilon_{{\rm s}*}$, and it yields, after solving for the actual shock location, the value of

\begin{displaymath}\epsilon_{\rm s}
=
\epsilon_{{\rm s}*}
\left({
{P_{\rm s}}\ov...
...rm s}}\over{(4-3u)r_{\rm s}^6}}
}\right]^{\frac{57}{40}}\cdot
\end{displaymath} (C.2)

In our main illustrative case, with $M_*=1.0~M_\odot$, $V_*=6.910\times10^{8}~{{\rm cm}}~{\rm s}^{-1}$, $\dot{m}=2.0\ {\rm g}~{{\rm cm}}^{-2}~{\rm s}^{-1}$ and approximately solar abundances, we derive $\epsilon_{\rm s*}=3.2, 72.7, 312$ and $\epsilon_{\rm s }=1.6, 64.6, 297$ for $B_* = 10, 30, 50\ {\rm MG}$ respectively.

References

 

Copyright ESO 2005