text:
# 1 Introduction

## 1 Introduction

In spite of the great success in the quantum field theory, there is still no consistent scenario in which the quantum gravity is included. The most promising scenario in this direction would be the string theory, where the consistency is ensured by the requirement of additional dimensions. Initially the sizes of extra dimensions had been assumed to be as small as $`M_p^1`$, but later it has been observed that there is no reason to believe such a tiny compactification radius. The idea of the large extra dimension may solve or weaken the hierarchy problem. Denoting the volume of the $`n`$-dimensional compact space by $`V_n`$, the observed Planck mass is obtained by the relation $`M_p^2=M_{}^{n+2}V_n`$, where $`M_{}`$ denotes the fundamental scale of gravity. If one assumes more than two extra dimensions, $`M_{}`$ may be close to the TeV scale without conflicting any observable bounds.

Although this new idea inspired creativity in many physicists to lead them to a new paradigm of phenomenology, a drastic modification is required for the conventional cosmological scenarios. Models of inflation and baryogenesis are especially sensitive to this low fundamental scale. To avoid extreme fine tuning, we should reconstruct the conventional scenarios of the standard cosmology. This requires inclusion of novel ideas that are quite different from the conventional one that was used for the models with large fundamental scale $`M_{}M_p`$, where $`M_p`$ denotes the Planck scale. This makes it difficult to construct a specific model for the early evolution of the Universe. For example, if one puts the inflaton fields on the brane, their masses are required to be unnaturally small. On the other hand, in generic cases the mass of the inflaton is bounded from below to achieve successful reheating. Thus it seems quite difficult to construct a model for successful inflation driven by a field on the brane.

A way to avoid these difficulties is put forward by Arkani-Hamed et al., where inflation is assumed to occur before the stabilization of the internal dimensions. In this case, however, the late oscillation of the radion is a problem, which will be solved by the second weak inflation.

One can find other ways to solve problems of inflation in models with large extra dimensions. Due to some dynamical mechanisms, the extra dimension may be stabilized before the Universe exited from inflation. If the stabilization of the internal dimensions occurred before the end of inflation, inflation cannot be induced by the fields on the brane, since their energy densities are highly suppressed. In this case one may use the bulk field rather than a field on the brane.

One may find another interesting possibility, “brane inflation”, which uses the interbrane distance as the inflaton.

In this paper we focus our attention to the second possibility, where the bulk field drives hybrid inflation after the radion stabilization.

We first make a brief review of the hybrid bulk field inflation, then we will examine the effect of the Kaluza-Klein modes.222 In ref., interesting aspects of bulk inflation after the stabilization of the radion field is studied by using the KK modes in the 4D effective theory. See also ref.. In ref., it is claimed that the phase transition becomes so slow that it is impossible to produce the required density perturbations after inflation with normal parameter values. However, we find that the phase transition is fast because of the huge number of the destabilized Kaluza-Klein modes at the end of inflation. The number of the destabilized Kaluza-Klein modes becomes up to about $`O(M_p^2/M_{}^2)`$. Although the problem of the slow phase transition seems to be solved by the Kaluza-Klein modes, another problem is induced by the excited Kaluza-Klein modes. The overproduction of the excited Kaluza-Klein modes is the problem, because they efficiently emit Kaluza-Klein gravitons when they decay into lower excited modes.

## 2 Hybrid inflation due to bulk scalar fields

In ref., it is argued that a single field inflationary model cannot provide adequate density perturbations, either with the inflaton on the brane or with the bulk field inflation. Thus it is natural to invoke a hybrid inflationary model in higher dimensions with the potential for the bulk field,

$$V(\varphi _5,\sigma _5)_{5D}=\lambda ^2M_{}\left(\sigma _{5o}^2\frac{1}{M_{}}\sigma _5^2\right)^2+\frac{m_\varphi ^2}{2}\varphi _5^2+\frac{g^2}{M_{}}\sigma _5^2\varphi _5^2,$$
(2.1)
where $`\varphi _5`$ is the inflaton field with a chaotic initial condition. The coupling constants $`\lambda `$ and $`g`$ are assumed to be $`O(1)`$. To keep things as simple as possible, here we have assumed five dimensional theory with only one extra dimension. Extensions to higher dimensional theories are straightforward. As the higher dimensional fields $`\varphi _5`$ and $`\sigma _5`$ has a mass of dimension $`3/2`$, the interaction terms are nonrenormalizable. Higher dimensional fields $`\varphi _5,\sigma _5`$ are related to the effective four dimensional fields $`\varphi ,\sigma `$ by a scaling

$$\varphi =\sqrt{R}\varphi _5,\sigma =\sqrt{R}\sigma _5,$$
(2.2)
which leads to couplings suppressed by the four dimensional Planck mass $`M_p`$, where $`M_p^2=M_{}^3R`$. Here $`R`$ denotes the volume of the extra dimensions. After dimensional reduction, the effective four dimensional potential for the 0-modes reads;

$`V(\varphi ,\sigma )_{4D}`$ $``$ $`\lambda ^2{\displaystyle \frac{M_{}^2}{M_p^2}}\left({\displaystyle \frac{M_p^2\sigma _0^2}{M_{}^2}}\sigma ^2\right)^2`$ (2.3)
$`+{\displaystyle \frac{m_{\varphi ^2}}{2}}\varphi ^2+g^2\left({\displaystyle \frac{M_{}}{M_p}}\right)^2\varphi ^2\sigma ^2.`$
Including the Kaluza-Klein modes, the potential becomes333 We have assumed that the nonrenormalizable terms in the five dimensional potential (2.3) are obtained by integrating out massive modes. These massive modes, which mediate the interactions $`\sigma _5^4`$ or $`\varphi _5^2\sigma _5^2`$, are assumed to have no Kaluza-Klein number.

$`V(\varphi _n,\sigma _n)_{KK}`$ $``$ $`\lambda ^2{\displaystyle \frac{M_{}^2}{M_p^2}}\left({\displaystyle \frac{M_p^2\sigma _0^2}{M_{}^2}}(\sigma ^2+{\displaystyle \underset{n}{}}\sigma _n^2)\right)^2`$ (2.4)
$`+{\displaystyle \frac{m_{\varphi ^2}}{2}}(\varphi ^2+{\displaystyle \underset{n}{}}\varphi _n^2)+g^2\left({\displaystyle \frac{M_{}}{M_p}}\right)^2(\varphi ^2+{\displaystyle \underset{n}{}}\varphi _n^2)(\sigma ^2+{\displaystyle \underset{n}{}}\sigma _n^2)`$
$`+{\displaystyle \underset{n}{}}{\displaystyle \frac{n^2}{R^2}}\varphi _n^2+{\displaystyle \underset{n}{}}{\displaystyle \frac{n^2}{R^2}}\sigma _n.`$
Let us first ignore the Kaluza-Klein modes in the effective potential and discuss the hybrid bulk field inflation induced by the potential eq.(2.3). In ref., the parameter values $`\sigma _oM_{}=10^5`$GeV, $`m_\varphi 10`$GeV were considered. With these values, the density perturbation

$$\frac{\delta \rho }{\rho }\left(\frac{g}{2\lambda ^{3/2}}\right)\times 10^5$$
(2.5)
were obtained. Then the reheating temperature is calculated within the perturbation theory. In ref., the decay rates of the inflaton field are already discussed. For the decay of the inflaton field into two Higgs fields on the brane, the decay width becomes

$$\mathrm{\Gamma }_{HH}\frac{M_{}^4}{32\pi M_p^2m_\varphi }.$$
(2.6)
On the other hand, their decay into excited Kaluza-Klein modes are highly suppressed,

$$\mathrm{\Gamma }_{\varphi \varphi \varphi _n\varphi _n}\lambda ^2\frac{M_{}^2}{M_p^2}.$$
(2.7)
Thus it seems appropriate to calculate the reheating temperature only by the decay into Higgs fields on the brane. Then the reheating temperature becomes $`100`$ MeV, which satisfies the requirement for the successful big bang nucleosynthesis, and at the same time solves the problem of overproduction of Kaluza-Klein gravitons.

In this model, however, a serious problem was reported in ref.. The authors of ref. claimed that the phase transition becomes extremely slow and it is impossible to reproduce the present Universe with the above parameter values. Let us first make a brief review of the argument. The slope of the potential (2.3) in the $`\varphi `$ direction is given by

$$\frac{dV}{d\varphi }=\left[2g^2\left(\frac{M_{}}{M_p}\right)^2\sigma ^2+m_\varphi ^2\right]\varphi ,$$
(2.8)
and the inflaton $`\varphi `$ rolls down the potential with the Hubble parameter

$$H=\sqrt{\frac{8\pi }{3}}\frac{\lambda \sigma _0^2}{M_{}}.$$
(2.9)
The inflaton $`\varphi `$ slowly rolls down the potential until it reaches the critical value

$$\varphi _c=\frac{\lambda }{g}\left(\frac{\sigma _0M_p}{M_{}}\right),$$
(2.10)
where the potential in the $`\sigma `$ direction is destabilized. The evolution of the inflaton $`\varphi `$ is given by

$$\varphi =\varphi _cexp\left[\frac{1}{\sqrt{24\pi }\lambda }\left(\frac{m_\varphi }{\sigma _0}\right)^2M_{}t\right]$$
(2.11)
where we take the time variable $`t`$ as $`t=0`$ at $`\varphi =\varphi _c`$. In the $`\sigma `$ direction, the slope of the potential (2.3) is

$$\frac{dV}{d\sigma }=\left(\frac{M_{}}{M_p}\right)^2\left[\lambda ^2\sigma ^2+2g^2(\varphi ^2\varphi _c^2)\right]\sigma .$$
(2.12)
At $`t=0`$, the initial value of $`\sigma `$ is expected to be about $`\sigma _{ini}HM_{}`$. The phase transition lasts until the first term in eq.(2.8) comes to dominate the evolution. It happens when

$$2g^2\left(\frac{M_{}}{M_p}\right)^2\sigma ^2m_\varphi ^2,$$
(2.13)
where $`\sigma `$ becomes about $`\sigma _{end}10^{14}`$GeV. Solving the equation for $`\sigma `$, one obtains large e-foldings $`N_e10^5`$. During this period, $`\sigma `$ evolves from $`\sigma _{ini}10^5`$GeV to $`\sigma _{end}10^{14}`$GeV.

However, the Kaluza-Klein modes are important in this case. Seeing the effective potential eq.(2.4), one can easily find that not only the 0-mode $`\sigma `$ but also the Kaluza-Klein modes $`\sigma _n`$ are destabilized to contribute the phase transition. The number of the destabilized modes at the end of inflation is estimated by the potential (2.4),

$$N_{KK}(R\lambda \sigma _0)\left(\frac{M_p}{M_{}}\right)^2.$$
(2.14)
In the presence of a large number of destabilized channels near $`\sigma \sigma _n0`$, the phase transition becomes inevitably fast.

In the original calculation of ref., the initial value of $`\sigma `$ is $`\sigma _{ini}HM_{}10^5`$GeV, which is much smaller than the critical value $`\sigma _{end}`$ in eq.(2.13). On the other hand, in the case when $`\sigma ^2`$ in eq.(2.8) is replaced by the sum of the huge number of the Kaluza-Klein modes, the condition (2.13) is already satisfied as soon as the phase transition starts at $`\varphi \varphi _c`$. In this case, each Kaluza-Klein mode $`\sigma _n`$ may have the initial value of $`\sigma _nH`$, if the mass is smaller than the Hubble parameter.

Although the problem of the slow phase transition is avoided by the huge numbers of the Kaluza-Klein modes, there arises another serious problem. As is discussed in ref., the reheating after the bulk inflation is not a problem when only the 0-mode oscillates and decays after inflation. In this model, however, the decay products of the excited Kaluza-Klein modes will dominate the Universe after inflation, because the inflation ends with the oscillation and the production of the excited Kaluza-Klein modes. Then the Kaluza-Klein gravitons are produced when the excited Kaluza-Klein modes decay into lower modes. The decay rate is enhanced by the large number of accessible modes in the final state, thus the excited modes are very short lived.

Let us examine the condition for the Kaluza-Klein gravitons to decay safe before nucleosynthesis. The decay width of the Kaluza-Klein gravitons into fields on the brane is estimated in ref.,

$$\mathrm{\Gamma }\frac{E^3}{M_p^2}.$$
(2.15)
Here $`E`$ denotes the energy of the graviton propagating in the bulk. In the most optimistic case, when $`EM_{}`$, the Kaluza-Klein gravitons may decay before nucleosynthesis if the fundamental scale is larger than $`10^6`$GeV.

## 3 Conclusions and Discussions

In this paper, we have discussed the possibility of successful hybrid inflation due to the bulk scalar field, in models with large extra dimensions. Bulk field inflation is already discussed in papers, where the problem of Kaluza-Klein gravitons and the problem of the slow phase transition are discussed. In ref., it is argued that the production rate of the Kaluza-Klein graviton is so small that the energy of the inflaton is safely drained into the standard model fields on the brane. However, including the Kaluza-Klein interactions, the production of the Kaluza-Klein gravitons becomes efficient. The problem of the slow phase transition is discussed in ref.. However, because of the huge number of excited Kaluza-Klein states that becomes unstable at the end of inflation, the phase transition becomes fast. The remaining problem is the overproduction of the Kaluza-Klein gravitons, which puts a lower bound to the fundamental scale. Even in the most optimistic case, the bound becomes $`M_{}>10^6`$GeV.

The significant effect of the Kaluza-Klein excited states, which we have discussed for the bulk inflation, is generically important in cosmological models that utilizes the phase transition of the bulk field.

## 4 Acknowledgment

We wish to thank K.Shima for encouragement, and our colleagues in Tokyo University for their kind hospitality.


format: text

text:
# Hydrodynamic Models of Line-Driven Accretion Disk Winds III: Local Ionization Equilibrium

## 1 Introduction

Ultraviolet observations of Cataclysmic Variables (CVs) show that virtually all nonmagnetic CVs with high mass accretion rates ($`4\times 10^{10}M_{\mathrm{}}\mathrm{yr}^1`$) exhibit P Cygni profiles in the resonance lines (e.g. Córdova & Mason, 1982; Greenstein & Oke, 1982; Guinan & Sion, 1982; Klare et al., 1982; Vitello & Shlosman, 1993; Prinja & Rosen, 1995; Friedjung, Selvelli, & Cassatella, 1997; Gänsicke, Beuermann, & Thomas, 1997; Knigge et al., 1997; Prinja et al., 2000). These profiles consist of blue-shifted absorption and redshifted emission, with the emission most prominent in high inclination systems (e.g. Hutchings, 1980; Krautter et al., 1981; Holm, Panek, & Shiffer III, 1982; Córdova & Mason, 1985; Mason et al., 1995; Cheng et al., 1997), and the absorption appearing at low and intermediate inclinations. The most straightforward explanation for these observations is in terms of a wind originating from the accretion disk. The wind is driven by radiation pressure from the disk UV continuum absorbed by line transitions, analogous to the radiation driven winds in early type stars. Support for this scenario came from Drew (1987), who calculated theoretical C IV 1549 Å line-profiles from one-dimensional kinematic CV wind models, for a wide variety of wind and disk parameters. Drew (1987) found that the dependence on the form of the UV resonance lines with inclination angle could be accounted for by a bipolar wind emerging from an accretion disk.

One-dimensional dynamical models for line-driven winds (LDW) in CVs were developed independently by Vitello & Shlosman (1988) and Kallman (1988). Further efforts came through two-dimensional kinematic modeling (Shlosman & Vitello, 1993; Knigge, Woods, & Drew, 1995), which succeeded in showing consistency between the assumed polar geometry of a disk wind and observed profiles; but such models were unable to self-consistently calculate the wind dynamics including the effects of rotation and the anisotropic disk radiation field.

Two-dimensional disk wind dynamical models were developed by Icke (1980), but these did not take into account the radiation pressure due to line absorption. For a typical white dwarf mass of $`0.6M_{\mathrm{}}`$ (Leibert, 1980; Silvestri et al., 2001) the disk luminosity required to produce such a wind without line radiation pressure (assuming radiation pressure due to continuum scattering only) would require a luminosity of $`10^5L_{\mathrm{}}`$ (Icke, 1980), several magnitudes above the observed luminosity of CVs which vary between $`0.01L_{\mathrm{}}`$ and $`10L_{\mathrm{}}`$ (Patterson, 1984). Icke (1980) did find that, for sufficiently high disk luminosities, a biconical disk wind would form. Icke (1981) suggested that biconical winds were a general property of accretion disk winds independent of the wind driving mechanism.

Two-dimensional dynamical models of isothermal line-driven disk winds were presented by us Pereyra, Kallman, & Blondin (1997, hereafter Paper 1). Results from Paper 1 showed, in analogy with line-driven winds from early type stars, that terminal velocities are approximately independent of the luminosity of the disk, although increments in luminosity produce increments in mass-loss rate, and that rotational forces are important in the study of winds from accretion disks. They cause the velocity streamlines to collide, which reduce the speed and increase the density of the wind producing an enhanced density region. The highest absorption occurs in the enhanced density region where density is increased relative to a spherically diverging wind with the same mass loss rate and the velocity is roughly half the terminal velocity. This density increase is necessary in order to produce at least marginally optically thick lines.

We developed (Pereyra, 1997; Pereyra, Kallman, & Blondin, 2000, hereafter Paper 2) an adiabatic wind model (rather than isothermal), driven by a standard accretion disk (Shakura & Sunyaev, 1973; Lynden-Bell & Pringle, 1974) rather than an isothermal disk. The hydrodynamic models developed of line-driven accretion disk winds (LDADWs) in CVs included the radial structure of an optically thick accretion disk with the corresponding radiation fields and surface temperature distributions. The corresponding energy conservation equation was also implemented, including the adiabatic heating and cooling effects due to compression and expansion. From the computational models, assuming single scattering and constant ionization, we calculated theoretical line profiles for the C IV 1550 Å line, and found that the line profiles obtained were consistent with observations in their general form and strong dependence with inclination angle.

Proga, Stone, & Drew (1998) also developed two-dimensional models for these systems, obtaining similar results to the models presented here. They assumed an isothermal wind and constant values for the line radiation parameters throughout the wind. They found similar wind velocities and wind mass loss rates as we do here, although they assumed disk luminosities about an order of magnitude greater. A difference we find with the results of the models presented by Proga, Stone, & Drew (1998) is that they obtain unsteady flows characterized by large amplitude fluctuations in velocity and density, while for our previous models (Paper 1; Pereyra 1997; Paper 2) and in the models presented here we find a steady flow. We note that these large fluctuations (rather than the steady flow we find) could account for the difference in luminosity required by the Proga models, with respect to our models, to obtain similar wind mass loss rates. One possibility suggested by Proga, Stone, & Drew (1998) was that the spatial resolution of the models we presented in Paper 1 was too low to adequately resolve the wind. In paper 2 we increased significantly the spatial resolution of our models and obtained similar steady results when compared to our earlier spatial resolutions used in Paper 1; thus showing that differences in spatial resolution of the models was not generating the difference between our steady disk wind models and the “intrinsically unsteady” disk winds reported by Proga, Stone, & Drew (1998). Another possibility is that the difference in the unsteady flows of the Proga models with respect to the steady flows we find is due to the boundary condition treatment in the Proga models, which appears to lead to instabilities at the base of the wind and in turn to the strong fluctuations reported.

It is well known that for early type stars, when a Sobolev treatment is used for the line radiation pressure (as we have done in our models and as has been applied by Proga, Stone, & Drew (1998) for accretion disk winds in CVs), a steady wind solution is obtained (e.g. Owocki & Puls, 1999). The treatment of boundary conditions in numerical models for line driven winds in early type stars was discussed by Owocki, Cranmer, & Blondin (1994) and Cranmer & Owocki (1996). They find that the treatment of the lower (base of wind) boundary condition is “somewhat problematic”. Further they find that if they used an arbitrarily high density, similar to the boundary conditions used by Proga, Stone, & Drew (1998), they would find oscillations at the base of the wind in both density and velocity. In both papers, in addition to other details in the boundary condition treatment (Owocki, Cranmer, & Blondin, 1994; Cranmer & Owocki, 1996), they lowered the density at the inner boundary until they found steady conditions at the base and were simultaneously able to resolve acceleration in the subsonic region. Therefore a conclusion in the area of numerical modeling that we derive is that by varying the treatment of boundary conditions at the base of a line driven wind, one can obtain strong fluctuations in density and velocity which are not necessarily physically intrinsic but are rather a consequence of the numerical details of the model.

Feldmeier & Shlosman (1999) and Feldmeier, Shlosman, & Vitello (1999) developed two-dimensional stationary models for CV disk winds where they solve the dynamics under the assumption that the wind streamlines lie on straight cones. The values of wind mass loss rates obtained through the models presented by Feldmeier, Shlosman, & Vitello (1999) for an isothermal disk (equation  of their paper) are consistent with those found in Paper 1. In addition Feldmeier, Shlosman, & Vitello (1999) find that in an isothermal disk (equations - of their paper) the wind mass loss rate scales as $`\dot{M}_{wind}L_{disk}^{1/\alpha }`$, where $`L_{disk}`$ is the disk luminosity and $`\alpha `$ is line force multiplier parameter. This is the same scaling law reported in Paper 1. Further, the biconical geometry and wind tilt angle in the models of Feldmeier, Shlosman, & Vitello (1999) are roughly consistent with the velocity structures we find through full 2.5D hydrodynamic modeling of our earlier models (Paper 1; Pereyra 1997; Paper 2) and in the models presented here. We further note that Feldmeier et al. suggest that detailed non-LTE calculations of the line-force multiplier may be important in the modeling of CV disk winds. This work is an additional step in that direction.

In spite of these encouraging results, a more realistic model which would include a detailed calculation of the UV spectrum from dynamically calculated density and velocity fields is desirable in order to make specific comparisons with observations. For instance, a maximum flux greater than $`1.7`$ times the continuum has been reported for the emission component of P Cygni line profiles of C IV 1550Å from some cataclysmic variables (e.g. Córdova & Mason, 1982), while from our previous models (Pereyra 1997; Paper 2) we find that the emission maximum of the C IV lines stays below 1.3 times the continuum flux for those angles where the absorption component is observed. This discrepancy may be due to some of the simplifications in our previous models or it may indicate a C IV line-emitting region other than the disk wind itself. Furthermore, in recent observational results, Froning et al. (2001) find relatively narrow absorption lines ($`700\mathrm{km}\mathrm{s}^1`$) in the FUV spectra of the cataclysmic variable U Gem. As Froning et al. (2001) have indicated, such lines cannot arise from the disk photosphere. In their work they suggest the possibility of an outer disk chromosphere to account for these lines, but it is also possible that some of these narrow FUV lines could be accounted for by an accretion disk wind.

An additional motivation for this work is that the models developed here can also be extended to the study of LDADWs in Quasars (QSOs), where the local ionization equilibrium will play an important in the overall dynamics (Murray et al., 1995; Hillier et al., 2002). Evidence for winds in QSOs is found in the broad absorption lines (BALs) observed in approximately 10% of the QSOs (e.g. Weymann, 1997). If the wind in QSOs is generated from an accretion disk, then the presence of BALs in only a fraction of QSOs can be accounted for as a viewing-angle effect (Turnshek, 1984); similar to the case of CVs where the accretion disk wind produces P Cygni profiles when observed “face-on” and emission lines (without absorption components) when seen “edge-on” (Pereyra, Kallman, & Blondin, 2000). The x-ray luminosity in QSOs are generally comparable to the UV/optical luminosities (e.g. Tananbaum et al., 1979; Grindlay et al., 1980; Mushotzsky, Done, & Pounds, 1993; Laor, 1997; George et al., 2000). Drew & Boksenberg (1984) found that such high x-ray luminosities would ionize the wind to a point where the populations of ions responsible for the BAL observed in QSOs would be too low to produce observable absorption lines. Murray et al. (1995) found that, with an appropriate x-ray shielding mechanism, LDADWs can produce wind densities and velocities consistent with BALs observed in QSOs. Strong changes in the ionization balance can produce observable effects in line formation in QSOs; as well as significant changes in the line radiation force parameters which will play an important role in the overall wind dynamics (Hillier et al., 2002).

In our previous models (Pereyra 1997; Paper 2) we had assumed constant ionization equilibrium throughout the wind and single scattering in the line profile calculations. In this work we extend our previous models to include radiative heating and cooling and ionization balance in the wind, and a calculation of the line radiation force parameters at each point in space at each time step. Our results are similar to those of Paper 2, and therefore act to confirm the general results found in our previous simpler models. They are a further step towards a model capable of self-consistently producing theoretical spectrum which may be compared in detail with observations; they can also be extended to the study of LDADWs in QSOs/AGN, where the local ionization equilibrium will play a crucial role in the overall dynamics (Murray et al., 1995; Hillier et al., 2002).

In §2 we discuss the radial structure of the accretion disk and the radiation field as implemented in our models. In §3 we present the ionization balance calculations applied in this work. We derive the expression used for the treatment of the line radiation pressure in §4. The treatment of radiative heating and cooling in this work is discussed in §5. In §6 we present and discuss the hydrodynamic calculations. The results of the models are presented and discussed in §7. We present a summary and the conclusions of this work in §8.

## 2 Accretion Disk

A key ingredient in the models presented in this work is the radiation field emerging from the accretion disk. The boundary layer where the disk intercepts the white dwarf and the white dwarf photosphere may also contribute to the radiation field. The white dwarf by itself (without considering the effects of the accreting mass) will have a luminosity of the order of $`0.01L_{\mathrm{}}`$ (Leibert, 1980), while the accretion disk will have typical luminosities of the order of $`L_{\mathrm{}}`$ (Warner, 1987). As we indicated in Paper 2, there is still uncertainty concerning the existence and spectrum of the boundary layer (Vrtilek et al., 1994; Mauche, Raymond, & Mattei, 1995; Long et al., 1996; Szkody et al., 1999; Mauche & Raymond, 2000; Froning et al., 2001). In the models presented here we have chosen not to explore the boundary layer issue in depth, but instead to consider the non-boundary layer scenario. We assume in this work that the radiation field is generated by the accretion disk alone, neglecting any boundary layer or white dwarf radiation. We assume here the standard accretion disk (Shakura & Sunyaev, 1973; Lynden-Bell & Pringle, 1974). In paper 2 we discussed some of the observational evidence for the existence of an accretion disk in CVs. We note that in this work we study the winds originating in nonmagnetic CVs in which an accretion disk is formed.

In the standard accretion disk model, the disk is assumed to be in a steady state. Shear stresses transport angular momentum outwards as the material of the gas spirals inwards. Conservation of angular momentum leads to the following expression:

$$W2\pi r^2\dot{M}_a\omega r^2=C(\mathrm{constant})$$
(1)
where $`r`$ is the radius, $`\dot{M}_a`$ is the mass accretion rate, $`\omega `$ is the corresponding angular velocity. and $`W`$ is defined by:

$$W\underset{z_0}{\overset{z_0}{}}w_{r\varphi }𝑑z$$
(2)
where $`z_0`$ is the half thickness of the disk and $`w_{r\varphi }`$ is shear stress between adjacent layers. With the additional assumption that the shear stresses are negligible in the inner disk radius:

$$W=\frac{\dot{M}_a}{2\pi r^2}(\omega r^2\omega _0r_0^2)$$
(3)
where $`\omega _0`$ is the angular velocity at the inner disk radius $`r_0`$. In CVs the inner disk radius corresponds to the white dwarf radius $`R_{wd}`$, i.e. the accretion disk extends from the inner Lagrangian point of the binary system down to the white dwarf surface.

Taking into account the work done by shear stresses, and assuming that, as the mass accretes inwards, the gravitational energy lost is emitted locally:

$$Q=\frac{1}{4\pi r}\frac{d}{dr}\left[\dot{M}_a\left(\frac{\omega ^2r^2}{2}\frac{GM}{r}\right)W2\pi \omega r^2\right]$$
(4)
where $`Q`$ is the radiated energy per area of the disk surface. Assuming that the disk material follows Keplerian orbits, i.e.

$$\omega =\left(\frac{GM_{wd}}{r^3}\right)^{1/2}$$
(5)
where $`G`$ is the gravitational constant and $`M_{wd}`$ is the white dwarf mass; from equations (3)-(4) one finds:

$$Q=\frac{3\dot{M}_{accr}GM_{wd}}{8\pi r^3}\left[1\left(\frac{r_o}{r}\right)^{1/2}\right]$$
(6)
The function $`Q(r)`$ (equation ) was originally derived by Shakura & Sunyaev (1973) for binary systems with an accretion disk about a black hole. Popham & Narayan (1995) developed detailed models of boundary layers and accretion disks in cataclysmic variables and found that their models supported the same function $`Q(r)`$ for the accretion disk in CVs.

Further, assuming that the disk is emitting locally as a blackbody, the radial temperature distribution of the disk will be given by:

$$T=\left\{\frac{3\dot{M}_{accr}GM_{wd}}{8\pi r^3\sigma _s}\left[1\left(\frac{r_o}{r}\right)^{1/2}\right]\right\}^{1/4}$$
(7)
where $`\sigma _s`$ is the Stefan-Boltzmann constant. Equation (6) is implemented in our model in the calculation of radiation flux throughout the wind, and consequently in the calculation of the continuum and line radiation pressure force. We assume local blackbody radiation. In our model we also assume that the temperature at the base of the wind is equal to the disk surface temperature for each radius. We wish to note here that the blackbody assumption of equation 7 will not strictly hold. Radiative diffusion of the energy from the inner parts of the disk to the disk photosphere will produce deviations from the blackbody assumption, limb darkening effects will produce deviations from the projected area in the angular dependence of the emitted flux, line scattering and electron scattering will produce deviations from the blackbody emission, the relatively high velocity gradients in the rotational shears as well as the wind velocity gradients will affect the line scattering; thus a detailed calculation of the spectrum emitted from the accretion disk would requires a 2.5D radiative transfer calculation of the disk atmosphere coupled with the corresponding velocity gradients.

Early attempts to calculate deviations of the actual disk spectrum from the local-blackbody assumption, applied locally spectra from stellar atmosphere codes (Herter et al., 1979; Kiplinger, 1979, 1980; Mayo, Wickramasinghe, & Whelan, 1980; Kolykhalov & Sunyaev, 1984). More recently, local plane parallel radiative transfer models implementing disk vertical structure have been developed (Diaz, Wade, & Hubeny, 1996; Hubeny & Hubeny, 1997, 1998; Wade & Hubeny, 1998; Hubeny et al., 2000). However, to our knowledge a self-consistent 2.5D radiative transfer model coupled with velocity shears of an accretion disk atmosphere has not yet been developed. Also, there still remains uncertainties regarding the variation of viscosity in the disk which in turn effects the vertical structure. Thus, due to its simplicity, the assumption of local blackbody spectrum has often been used to model accretion disk and their spectra (e.g. Shlosman & Vitello, 1993; Knigge, Woods, & Drew, 1995; Proga, Stone, & Drew, 1998; Pereyra, Kallman, & Blondin, 2000). In future work we shall consider the radiative transfer effects that deviate the disk spectrum from an ideal local-blackbody, and attempt to couple them with hydrodynamic calculations.

It is generally accepted that the mechanism for producing the shear stress in accretion disks is magnetic turbulence generated by the Balbus-Hawley magnetorotational instability (MRI) (Balbus & Hawley, 1991). Three-dimensional magnetohydrodynamic simulations (e.g. Armitage, 1998; Machida, Hayashi, & Matsumoto, 2000; Hawley, Balbus, & Stone, 2001; Hawley & Krolick, 2001) suggest that the standard assumptions (e.g. Keplerian orbits; zero stress at inner radii) may not hold and that therefore corrections to equation (6) may arise. We believe that the standard Shakura-Sunyaev disk includes much of the physics involved in accretion disks, such as: gravity from a large central mass determining the average velocity fields in the disk, angular momentum transport, conversion of gravitational energy to radiative emission, and monotonic increase of temperature with decreasing radius (except of radius close to the inner disk radii) and local blackbody emission which accounts for the relatively “flat” continuum spectrum observed in systems where accretion disk are inferred. We also believe that it is important to explore “simple” models where the number of free parameters is small. More realistic disk models and radiation fields will be introduced at a later date.

## 3 Ionization Equilibrium and Line Opacity

The main wind driving mechanism in our model is the line radiation force generated as the disk emission encounters ionized gas above the disk. The wind material is thus driven from the disk photosphere achieving velocities of the order of a few thousand of $`\mathrm{km}\mathrm{s}^1`$. The line radiation force will depend on the available line opacity distribution of the wind material, which in turn is a direct function of the ion populations. Therefore the ion populations become an important ingredient in our models. In our previous work (Paper 1, Pereyra 1997, Paper 2), in an approach similar to that presented by Castor, Abbott, & Klein (1975) for OB star winds, we had assumed constant ionization equilibrium throughout the wind. Further we had assumed that the ion populations were similar to those of early type stars. These two assumptions allowed us to apply typical OB star line radiation force parameters (Abbott, 1982) throughout the wind, significantly simplifying the numerical treatment. In this work we include the local ion population calculations simultaneously with the hydrodynamics as we describe below.

The ion populations of the wind material is initially at local thermodynamic equilibrium (LTE), since the wind starts at the disk photosphere. As the wind material rises, the ion populations are no longer in LTE. The material is irradiated by the disk, which is a non-Plank source (equations -). Due to the relatively low wind density, the ion populations above the disk are largely determined by photoionization and radiative recombination processes, which in turn are dependent on the incoming disk radiation. It is difficult to model the collisional and radiative processes which determine the ion populations, varying in space and time, while simultaneously solving the hydrodynamics of the system. In this work we assume that the ion populations in the disk wind are determined by photoionization and radiative recombination processes, similar to a photoionized nebula irradiated by an accretion disk. To illustrate the validity of this approximation we shall estimate the ratio between the radiative photoionization rate $`R_{ij}`$ and the collisional ionization rate $`C_{ij}`$. Under the condition that $`h\nu _okT_R`$ and $`h\nu _okT_e`$, where $`h`$ is Plank’s constant, $`h\nu _o`$ is the minimum energy required to ionize an atom from state $`i`$ to state $`j`$, $`k`$ is Boltzmann’s constant, $`T_R`$ is the temperature of the ionizing radiation, and $`T_e`$ the electron temperature, the ratio can be approximated by (Mihalas, 1978):

$$\frac{R_{ij}}{C_{ij}}\frac{4(2\pi ^3k)^{\frac{1}{2}}h\nu _o^3}{3m_e^{\frac{1}{2}}e^2c^3}\left(\frac{WT_R}{n_eT_e^{\frac{1}{2}}}\right)\mathrm{exp}\left[h\nu _o\left(\frac{1}{kT_e}\frac{1}{kT_R}\right)\right]$$
(8)
where $`m_e`$ is the electron mass, $`e`$ the electron charge, $`c`$ the speed of light, $`W`$ the dilution factor of the radiation from its source, and $`n_e`$ is the electron density. Taking for the CV case: $`n_e10^8`$ (this density is a typical value for the CV disk wind models we presented in Paper 2), $`T_RT_e30,000\mathrm{K}`$, $`W1/2`$, and $`h\nu _o13.6\mathrm{ev}`$, from equation (8) we find:

$$\frac{R_{ij}}{C_{ij}}1.4\times 10^8.$$
(9)
As we have indicated above, in order to calculate the force due to line scattering, one in principle must consider the contribution to all lines which in turn will depend on the ion populations at a particular point in space in a particular time. In the calculation of ion population and line opacity we follow a similar scheme as that applied by Abbott (1982) for OB stars. In this work we assume that photoionization and radiative recombination processes determine the ion populations, therefore the ionization balance in the wind is determined through Saha’s equation for a gas photoionized by an external source:

$$\frac{N_{j+1}n_e}{N_j}=2W\frac{U_{j+1}}{U_j}\left(\frac{2\pi m_ekT_e}{h^2}\right)^{3/2}\frac{T_R}{T_e}\mathrm{exp}\left(\frac{X_{Ij}}{kT_R}\right)$$
(10)
where $`N_j`$ and $`N_{j+1}`$ are the atom densities at stage $`j`$ and at the following ionization stage $`j+1`$, $`U_j`$ and $`U_{j+1}`$ are the partition functions for the corresponding ionization stages, and $`X_{Ij}`$ is the ionization potential for the atom at stage $`j`$. Equation (10) has been applied to several astrophysical systems, e.g. stellar winds (Abbott & Lucy, 1985; Schmutz et al., 1990; Lucy & Abbott, 1993).

In the work presented here, we have assumed that the radiation is primarily coming from the accretion disk which presents a temperature dependent on radius (equation ). Equation (10) assumes that the photoionization rates are much greater than collisional ionization rates, thus:

$$N_jI_jN_{j+1}R_{j+1}$$
(11)
where $`I_j`$ is photoionization rate for an atom at ionization stage $`j`$, and $`R_{j+1}`$ is the radiative recombination rate for the atom at ionization state $`j+1`$. Therefore, from equation (10):

$$I_j=R_{j+1}\frac{2W}{n_e}\frac{U_{j+1}}{U_j}\left(\frac{2\pi m_ekT_e}{h^2}\right)^{3/2}\frac{T_R}{T_e}\mathrm{exp}\left(\frac{X_{Ij}}{kT_R}\right).$$
(12)
In the accretion disk, the radiation is emitted at different temperatures (equation ), thus the total photoionization rate $`I_{Tj}`$ due to the disk radiation must equal the sum of all the contributions of photoionization rates due to the radiation emitted at each section of the disk surface:

$$I_{Tj}=R_{j+1}\frac{2}{n_e}\frac{U_{j+1}}{U_j}\left(\frac{2\pi mkT_e}{h^2}\right)^{3/2}\frac{1}{T_e}_{disk}W(r,\varphi )T_R(r)\mathrm{exp}\left(\frac{X_{Ij}}{kT_R(r)}\right)r𝑑\varphi 𝑑r$$
(13)
where $`W(r,\varphi )`$ is the dilution factor per source area, given by:

$$W(r,\varphi )=\frac{1}{4\pi }\frac{\widehat{z}(\stackrel{}{r}_p\stackrel{}{r})}{|\stackrel{}{r}_p\stackrel{}{r}|}\frac{1}{|\stackrel{}{r}_p\stackrel{}{r}|^2}$$
(14)
i.e.:

$$W(r,\varphi )=\frac{1}{4\pi }\frac{z}{[R^22Rr\mathrm{cos}\varphi +r^2+z^2]^{3/2}}$$
(15)
where $`\stackrel{}{r}_p(R,0,z)`$ is the position where the dilution factor is evaluated and $`\stackrel{}{r}(r\mathrm{cos}\varphi ,r\mathrm{sin}\varphi ,0)`$ is the radiation source position on the disk.

Therefore, under the assumption that the photoionization and radiative recombination processes dominate over the collisional processes, we have:

$$\frac{N_{j+1}n_e}{N_j}=2\frac{U_{j+1}}{U_j}\left(\frac{2\pi mkT_e}{h^2}\right)^{3/2}\frac{1}{T_e}_{disk}W(r,\varphi )T_R(r)\mathrm{exp}\left(\frac{X_{Ij}}{kT_R(r)}\right)r𝑑\varphi 𝑑r$$
(16)
which is the equation that we implement in this work for the calculation of ionization balance for the line radiation force.

Once the ion populations at different stages is calculated for each atom through the resolution of equation (16), the population for each energy level for a given ion stage of a given atom is obtained through Boltzmann’s equations:

$$\frac{N_{ij}}{N_{0j}}=\frac{g_{ij}}{g_{0j}}\mathrm{exp}\left(\frac{X_{ij}}{kT_e}\right)$$
(17)
where $`N_{ij}`$ and $`N_{0j}`$ are the densities of a given atom at a given ionization stage $`j`$ at energy level $`i`$ and at the lowest energy level respectively, $`g_{ij}`$ and $`g_{0j}`$ are the degeneracies of energy level $`i`$ and at the lowest energy level respectively, and $`X_{ij}`$ is the difference in energy between level $`i`$ and the lowest energy level. We wish to note here that the LTE excited population assumption in equation (17), while allowing the local ionization equilibrium calculations to stay numerically manageable, will tend to be more accurate for the lower energy ion states which are the most crucial for the line radiation force. For early type stars Abbott (1982) found that roughly 90% of the line radiation force was due to lines arising from ground or metastable states.

After the ion populations at different stages at different energy levels for a given atom have been obtained, the line opacity of a given line for a given atom is calculated by the expression:

$$\kappa _L=\frac{\pi e^2}{m_ec}g_{lj}f_{luj}\frac{N_{lj}/g_{lj}N_{uj}/g_{uj}}{\rho }\frac{1}{\mathrm{\Delta }\nu _D}$$
(18)
where $`g_{lj}`$ and $`g_{uj}`$ are the degeneracies of the lower energy level $`l`$ and the upper energy level $`u`$ of ionization stage $`j`$ for a given line of a given atom, $`f_{luj}`$ is the absorption oscillator strength corresponding to the transition between lower energy level $`l`$ and upper energy level $`u`$ of ionization stage $`j`$ for a given line of a given atom, $`N_{lj}`$ and $`N_{uj}`$ are the ion densities of ionization stage $`j`$, lower energy level $`l`$, and upper energy level $`u`$ respectively, $`\rho `$ is the total mass density of the gas, and $`\mathrm{\Delta }\nu _D`$ is the Doppler width of the given line. The Doppler width of the line is given by:

$$\mathrm{\Delta }\nu _D=\frac{v_{th}}{c}\nu _L$$
(19)
where $`v_{th}`$ is the thermal velocity and the line frequency $`\nu _L`$ is given by:

$$\nu _L=\frac{X_{uj}X_{lj}}{h}.$$
(20)
## 4 Radiation Force

Once the line opacity $`\kappa _L`$ has been calculated through the above expressions for each atomic line at each ionization stage, the line radiation force is given by:

$$\stackrel{}{f}_{rad}=\underset{lines}{}\frac{\kappa _L\mathrm{\Delta }\nu _D}{c}\left[\frac{1}{\tau _L(\widehat{n})}_0^{\tau _L(\widehat{n})}e^\tau ^{}𝑑\tau ^{}\right]I(\widehat{n})\widehat{n}𝑑\mathrm{\Omega }$$
(21)
where $`I(\widehat{n})`$ is the radiation intensity, and $`\tau _L(\widehat{n})`$ is the line optical depth.

The line optical depth $`\tau _L(\widehat{n})`$ for a given line at a given direction $`\widehat{n}`$, under the Sobolev approximation, is given by (Rybicki & Hummer, 1978):

$$\tau _L(\widehat{n})=\frac{\rho \kappa _Lv_{th}}{|\widehat{n}\stackrel{}{}(\widehat{n}\stackrel{}{v})|}.$$
(22)
As in Paper 1 and Paper 2, we assume in this work that the velocity gradient is primarily in the “$`z`$” direction and approximate the line optical depth by:

$$\tau _L(\widehat{n})\frac{\rho \kappa _Lv_{th}}{\left|\widehat{n}\widehat{z}dv_z/dz\right|}.$$
(23)
Integrating equation (21) in optical depth:

$$\stackrel{}{f}_{rad}=\underset{lines}{}\frac{\kappa _L\mathrm{\Delta }\nu _D}{c}\left[\frac{1\mathrm{exp}\left(\tau _L(\widehat{n})\right)}{\tau _L(\widehat{n})}\right]I(\widehat{n})\widehat{n}𝑑\mathrm{\Omega }$$
(24)
therefore:

$$\stackrel{}{f}_{rad}=\left[\frac{\sigma _0^{\mathrm{}}I(\widehat{n})𝑑\nu }{c}\times \underset{lines}{}\frac{I(\widehat{n})\mathrm{\Delta }\nu _D}{_0^{\mathrm{}}I(\widehat{n})𝑑\nu }\left(\frac{1\mathrm{exp}\left(\frac{\kappa _L}{\sigma }\left[\frac{\sigma }{\kappa _L}\tau _L(\widehat{n})\right]\right)}{\frac{\sigma }{\kappa _L}\tau _L(\widehat{n})}\right)\right]\widehat{n}𝑑\mathrm{\Omega }$$
(25)
where $`\sigma `$ is Thomson cross section per unit mass.

In their work on radiation driven winds in OB stars, Castor, Abbott, & Klein (1975) and Abbott (1982), found that the sum over lines in equation (25) could be approximated by the following expression:

$$\underset{lines}{}\frac{I(\widehat{n})\mathrm{\Delta }\nu _D}{_0^{\mathrm{}}I(\widehat{n})𝑑\nu }\left(\frac{1\mathrm{exp}\left(\frac{\kappa _L}{\sigma }\left[\frac{\sigma }{\kappa _L}\tau _L(\widehat{n})\right]\right)}{\left[\frac{\sigma }{\kappa _L}\tau _L(\widehat{n})\right]}\right)k^{}(\widehat{n})\left(\frac{\sigma }{\kappa _L}\tau _L(\widehat{n})\right)^{\alpha ^{}(\widehat{n})}$$
(26)
where the $`k^{}(\widehat{n})`$ and $`\alpha ^{}(\widehat{n})`$ are constants for a given direction $`\widehat{n}`$. Defining $`t`$:

$$t=\frac{\sigma }{\kappa _L}\tau _L(\widehat{n})$$
(27)
we have:

$$\underset{lines}{}\frac{I(\widehat{n})\mathrm{\Delta }\nu _D}{_0^{\mathrm{}}I(\widehat{n})𝑑\nu }\left(\frac{1\mathrm{exp}\left(\frac{\kappa _L}{\sigma }t\right)}{t}\right)k^{}(\widehat{n})t^{\alpha ^{}(\widehat{n})}.$$
(28)
From equation (28) we obtain the $`k^{}(\widehat{n})`$ and the $`\alpha ^{}(\widehat{n})`$ in a form similar to the Abbott (1982) calculations for OB stars. We calculate the left hand side of equation (28) varying the $`t`$ parameter from $`10^1`$ to $`10^7`$, using the line opacity $`\kappa _L`$ calculated in the form described in the previous section for each line, and using the disk emission described by equations (6)-(7) to calculate the intensities. From these values, the parameters $`k^{}(\widehat{n})`$ and $`\alpha ^{}(\widehat{n})`$ are obtained through linear regression by varying the $`t`$ parameter over two order of magnitudes at a time within the fore mentioned interval, and taking the highest values found for each parameter. We note here that the selection of the highest values within intervals of two orders magnitudes is somewhat arbitrary. To our knowledge this is the first attempt to include the spatial and time dependence of the line radiation force parameters coupled consistently with the hydrodynamics of accretion disk winds, and we believe that the process includes much of the physics involved in line driving from a distributed radiation source: line-scattering, photon momentum generating force on gas, dilution factor, line optical depth in a moving media, line-opacity distribution, energy distribution in spectrum, and ion populations. This approach allows us to explore the effect of variable ion populations on the disk wind rather that assuming a constant ionization equilibrium as we did in Paper 2. In future work we shall develop a more accurate treatment for the line radiation force parameters, including a systematic determination of the appropriate intervals in $`t`$ in which to derive the line force parameters, as well as a study of the directional dependence of these parameters, rather than assume them isotropic as we do in this work (see equation ).

From equation (25) and equation (28), we have:

$$\stackrel{}{f}_{rad}\left[\frac{\sigma _0^{\mathrm{}}I(\widehat{n})𝑑\nu }{c}k^{}(\widehat{n})t^{\alpha ^{}(\widehat{n})}\right]\widehat{n}𝑑\mathrm{\Omega }.$$
(29)
As mentioned above, in this work we use values of the force multiplier parameters $`k`$ and $`\alpha `$ which vary in time and space but are independent of direction. As we see below this will allow us to implement the numerical integrations on intensity $`I(\widehat{n})`$ at each spatial point of the computational grid before implementing the hydrodynamic calculations rather than during the implementation of the hydrodynamic calculations, thus maintaining the computations manageable.

The values of $`k`$ and $`\alpha `$ used in this work are calculated through the expression:

$$kt^\alpha \frac{\left|\left[_0^{\mathrm{}}I(\widehat{n})𝑑\nu k^{}(\widehat{n})t^{\alpha ^{}(\widehat{n})}\right]\widehat{n}𝑑\mathrm{\Omega }\right|}{\left|\left[_0^{\mathrm{}}I(\widehat{n})𝑑\nu \right]\widehat{n}𝑑\mathrm{\Omega }\right|}.$$
(30)
We calculate the right hand side of equation (30) varying the $`t`$ parameter from $`10^1`$ to $`10^7`$, using the values of $`k^{}(\widehat{n})`$ and $`\alpha ^{}(\widehat{n})`$ calculated through equation (28) and using the disk emission described by equations (6)-(7) to calculate the intensities. From these values, the parameters $`k`$ and $`\alpha `$ are obtained through linear regression. We note here that the angle averaging is done in an ad hoc manner in order to reduce computational complexity. As we indicated above, we believe that this approach includes much of the physics of line driving within an accretion disk wind. In future work, we shall develop a more systematic treatment of the line radiation force parameters that will include a study of the directional dependence of these parameters, rather than assume averaged angle values as we do here.

With these values of $`k`$ and $`\alpha `$, we approximate the line radiation force by:

$$\stackrel{}{f}_{rad}\left[\frac{\sigma _0^{\mathrm{}}I(\widehat{n})𝑑\nu }{c}kt^\alpha \right]\widehat{n}𝑑\mathrm{\Omega }.$$
(31)
Therefore, from equation (27):

$$\stackrel{}{f}_{rad}=\left[\frac{\sigma _0^{\mathrm{}}I(\widehat{n})𝑑\nu }{c}k\left(\frac{\sigma }{\kappa _L}\tau _L(\widehat{n})\right)^\alpha \right]\widehat{n}𝑑\mathrm{\Omega }.$$
(32)
Substituting equation (23):

$$\stackrel{}{f}_{rad}=\left[\frac{\sigma _0^{\mathrm{}}I(\widehat{n})𝑑\nu }{c}k\left(\frac{1}{\rho \sigma v_{th}}\left|\widehat{n}\widehat{z}\frac{dv_z}{dz}\right|\right)^\alpha \right]\widehat{n}𝑑\mathrm{\Omega }.$$
(33)
Thus:

$$\stackrel{}{f}_{rad}=\frac{\sigma }{c}k\left(\frac{1}{\rho \sigma v_{th}}\left|\frac{dv_z}{dz}\right|\right)^\alpha \left[\left|\widehat{n}\widehat{z}\right|^\alpha _0^{\mathrm{}}I(\widehat{n})𝑑\nu \right]\widehat{n}𝑑\mathrm{\Omega }.$$
(34)
We wish to note that the $`k`$ and $`\alpha `$ parameters in the above equation depend on the line opacity distribution (equations -), which in turn depends on density and temperature (equations -). Since the density and temperature will vary in space and time, the $`k`$ and $`\alpha `$ parameters will also depend on space and time. Also, we wish to note here that in this work we assume azimuthal symmetry and thus have only two independent spatial coordinates “$`r`$” and “$`z`$”.

In equation (34) we substitute the $`\alpha `$ parameter inside the integral by a fiducial value of $`0.7`$, which is the value we assumed in our previous models presented in Paper 1 and Paper 2. Varying the number density from $`10^6\mathrm{cm}^3`$ to $`10^{10}\mathrm{cm}^3`$ at spatial points throughout our computational grid, and applying the values of disk luminosity $`L_{disk}=L_{\mathrm{}}`$ and white dwarf radii of $`R_{wd}=0.01R_{\mathrm{}}`$ which we use throughout this work, by the way of equation (30) we find values for the $`\alpha `$ parameter which vary from $`0.55`$ to $`0.85`$. We do not expect that the application of the fiducial value for $`\alpha `$ inside the integral of equation (34) will produce a significant difference in the density and velocity wind structures, and further it will allow the evaluation of the integral at all spatial grid points before solving the hydrodynamic equations described in §6, maintaining the computations manageable. Thus,

$$\stackrel{}{f}_{rad}=\frac{\sigma }{c}k(r,z,t)\left(\frac{1}{\rho \sigma v_{th}}\left|\frac{dv_z}{dz}\right|\right)^{\alpha (r,z,t)}\left[\left|\widehat{n}\widehat{z}\right|^{0.7}_0^{\mathrm{}}I(\widehat{n})𝑑\nu \right]\widehat{n}𝑑\mathrm{\Omega }.$$
(35)
We wish to also note that although equation (35) takes into account the line radiation pressure in the “$`r`$” direction, we have included only the derivative in the “$`z`$” direction of the “$`v_z`$” component. As we found in Paper 2, the accretion disk winds in CVs tend to flow in the direction perpendicular to the disk, particularly in the region near the surface of the disk, thus justifying the assumption that the derivatives of the “$`v_r`$” component are negligible. In this work we also assume negligible the radial derivatives of the “$`v_z`$” component and the vertical and radial derivatives of the “$`v_\varphi `$” component in the calculation of the line radiation force. We note here that the disk wind streamlines are “helical” in nature; the wind flow starts out at the base of the wind with a non-zero “$`v_\varphi `$” component due to disk rotation. The “$`v_\varphi `$” component of the velocity fields is calculated throughout our hydrodynamic calculations under the assumption of azimuthal symmetry (see Appendix). We note that even under azimuthal symmetry, line radiation force may have non-zero components along the the “$`\varphi `$” direction; e.g., Gayley & Owocki (2000) for the case of OB stars, assuming azimuthal symmetry, showed that the $`\varphi `$ component of the line radiation force, although presenting modest modifications in the overall wind dynamics, could lead to observable effects. Even though, on one hand, the inclusion of the additional terms in the velocity gradient for the line radiation pressure in the accretion disk wind may possibly produce interesting results, on the other we do not believe that it will significantly affect the overall results currently being implemented. In future work we shall include the effects of the derivatives along both the “$`z`$” and “$`r`$” directions for the different velocity components in the line radiation force, rather than assume “$`dv_z/dz`$” as the dominant term as we have done in this work and in Paper 1 and Paper 2.

Defining $`\stackrel{}{S}(r,z)`$:

$$\stackrel{}{S}(r,z)=\left[\left|\widehat{n}\widehat{z}\right|^{0.7}_0^{\mathrm{}}I(\widehat{n})𝑑\nu \right]\widehat{n}𝑑\mathrm{\Omega }.$$
(36)
Taking the source of radiation to be the accretion disk, we have:

$$S_z(r,z)=_{r_o}^{\mathrm{}}_0^{2\pi }\frac{Q(r^{})}{\pi }\frac{z^{2.7}}{[(r^2+r^2+z^22rr^{}\mathrm{cos}\varphi )^{1/2}]^{4.7}}r^{}𝑑\varphi 𝑑r^{}$$
(37)
and

$$S_r(r,z)=_{r_o}^{\mathrm{}}_0^{2\pi }\frac{Q(r^{})}{\pi }\frac{z^{1.7}(rr^{}\mathrm{cos}\varphi )}{[(r^2+r^2+z^22rr^{}\mathrm{cos}\varphi )^{1/2}]^{4.7}}r^{}𝑑\varphi 𝑑r^{}$$
(38)
where $`S_z`$ and $`S_r`$ are the corresponding components of the vector $`\stackrel{}{S}`$ defined in equation (36) (the $`S_\varphi `$ component is zero due to the axial symmetry of the accretion disk) and $`Q(r^{})`$ is the radiation emission per area of the disk calculated through equation (6).

Thus we have

$$\stackrel{}{f}_{rad}=\frac{\sigma }{c}k(r,z,t)\left(\frac{1}{\rho \sigma v_{th}}\left|\frac{dv_z}{dz}\right|\right)^{\alpha (r,z,t)}[S_z(r,z)\widehat{z}+S_r(r,z)\widehat{r}]$$
(39)
where $`S_z`$ and $`S_r`$ are given by equations (37) and (38).

We also account for the fact that, as discussed by Castor, Abbott, & Klein (1975) for OB star winds, if the velocity gradient increases to sufficiently high values, or the density of the wind decreases to sufficiently low values, the contribution for the radiation pressure due to lines arrives at a maximum value. Similar to the case of OB stars (Abbott, 1982), we found that the maximum value of the line radiation pressure was obtained when $`(1/\rho \sigma v_{th})dv_z/dz10^7`$. In this work we have assumed that the maximum possible value of the line radiation pressure is obtained when $`(1/\rho \sigma v_{th})dv_z/dz=10^7`$ (rather than $`10^8`$ as we did in Paper 2).

Therefore the expression for the total line radiation force per mass for disk winds which we adopt is

$$\stackrel{}{f}_{rad}=\frac{\sigma }{c}k(r,z,t)\left(\mathrm{max}[\frac{1}{\rho \sigma v_{th}}\left|\frac{dv_z}{dz}\right|,\mathrm{\hspace{0.17em}10}^7]\right)^{\alpha (r,z,t)}[S_z(r,z)\widehat{z}+S_r(r,z)\widehat{r}].$$
(40)
In order to calculate the $`k(r,z,t)`$ and $`\alpha (r,z,t)`$ functions self consistently throughout the hydrodynamic calculations described below and simultaneously keep the computation manageable the numerical code was written such that the spatial distribution of the $`k`$ and $`\alpha `$ are calculated for values of number density varying from $`10^6\mathrm{cm}^3`$ to $`10^{10}\mathrm{cm}^3`$ at each computational spatial grid point. A fiducial temperature is taken for these calculations, since as we find in the models presented here, in inner disk region where the line radiation force is most crucial, the temperature does not change significantly due to radiative heating and cooling. During the hydrodynamic calculations, the actual values of $`k`$ and $`\alpha `$ are obtained through geometric interpolation from the above values and the current value of density at each computational grid point.

For the continuum radiation pressure the total radiation flux is calculated throughout the disk wind:

$$\stackrel{}{F}=\left[_0^{\mathrm{}}I(\widehat{n})𝑑\nu \right]\widehat{n}𝑑\mathrm{\Omega }.$$
(41)
Taking the source of radiation to be the accretion disk, we have

$$F_z(r,z)=_{r_o}^{\mathrm{}}_0^{2\pi }\frac{Q(r^{})}{\pi }\frac{z^2}{[(r^2+r^2+z^22rr^{}\mathrm{cos}\varphi )^{1/2}]^4}r^{}𝑑\varphi 𝑑r^{}$$
(42)
and

$$F_r(r,z)=_{r_o}^{\mathrm{}}_0^{2\pi }\frac{Q(r^{})}{\pi }\frac{z(rr^{}\mathrm{cos}\varphi )}{[(r^2+r^2+z^22rr^{}\mathrm{cos}\varphi )^{1/2}]^4}r^{}𝑑\varphi 𝑑r^{}$$
(43)
where $`F_z`$ and $`F_r`$ are the corresponding components of the radiation flux vector $`\stackrel{}{F}`$ originating from the disk (the $`F_\varphi `$ component is zero due to the axial symmetry of the accretion disk) and $`Q(r^{})`$ is the radiation emission per area of the disk calculated through equation (6).

## 5 Radiative Heating and Cooling

Radiative heating and cooling are computed in our model assuming an optically thin disk spectrum. A local balance between heating and cooling is assumed, subject to the constraint of particle number conservation for each element. The ionization processes included in the calculation are photoionization, charge transfer ionization, and collisional ionization. Recombination processes include radiative and dielectronic recombination and charge transfer. The heating terms includes photoionization heating, Compton heating, charge transfer, and collisional de-excitation. Cooling terms include radiative and dielectronic recombination, bremsstrahlung, collisional ionization, collision excitation of bound levels, and (endothermic) charge transfer. Details of the radiative heating and cooling calculations are described in detail by Kallman & McCray (1982).

Disk models have succeeded in accounting for the observed continuum spectrum in many CVs (e.g. Herter et al., 1979; Kiplinger, 1979, 1980; Bath, Pringle, & Whelan, 1980; Mayo, Wickramasinghe, & Whelan, 1980; Krautter et al., 1981; Hassall, 1985; Wade, 1988; la Dous, 1989; Diaz, Wade, & Hubeny, 1996; Wade & Hubeny, 1998). The form of the spectrum applied in this work for the radiative heating and cooling calculations is obtained by integrating contributions over the disk surface calculated from equation (6), assuming local blackbody emission, and viewing the disk at infinity. In Figure 1 we show the disk spectrum used for the radiative heating and cooling calculations. We note that an assumption made in the radiative cooling and heating calculations is that the form of the spectrum is position independent, and thus the rates calculated here will be more accurate higher above the disk. In the future we will introduce more accurate heating and cooling rates into our models.

In Figure 2 and Figure 3 the radiative heating and cooling are shown respectively as a function of temperature for different values of the ionization parameter $`\xi `$ defined as:

$$\xi \frac{4\pi |\stackrel{}{F}|}{n}$$
(44)
where $`F`$ is the radiation flux and $`n`$ the number density.

From Paper 2 we found that the wind temperature of CV winds range in values of up to $`40,000\mathrm{K}`$. Thus in our work here we have approximated the heating and cooling rates the results presented in Figure 2 and Figure 3 by analytical expressions that agree with the detailed calculation within 20% for the heating rates and within 30% for the cooling rates for temperatures of up to $`50,000\mathrm{K}`$.

For the radiative heating rate we implemented the following expression (in units of $`\mathrm{erg}\mathrm{cm}^3\mathrm{s}^1`$):

$$\mathrm{\Gamma }=1.6\times 10^{21}T^{1/2}\left(1\frac{T}{1.3\times 10^5}\right)$$
(45)
where $`T`$ is the absolute temperature, valid for temperatures of up to $`50,000\mathrm{K}`$.

For the radiative cooling rate we implemented the following expression (in units of $`\mathrm{erg}\mathrm{cm}^3\mathrm{s}^1`$):

$$\mathrm{\Lambda }=1.9\times 10^{29}T^{3/2}$$
(46)
where $`T`$ is the absolute temperature, valid for temperatures of up to $`50,000\mathrm{K}`$.

## 6 Hydrodynamic Model

Our hydrodynamic model uses the Piece-Wise Parabolic Method (PPM) numerical scheme (Colella & Woodward, 1984). We use the 2.5-dimensional hydrodynamic equations of the wind in cylindrical coordinates; i.e., the three-dimensional hydrodynamic equations are reduced to two-dimensional equations by assuming that the derivative of any physical variable with respect to $`\varphi `$ is zero. The origin of the coordinate system is located at the center of the accretion disk and the “z=0” plane is located to coincide with the accretion disk midplane. The hydrodynamic equations are presented in the Appendix.

In this treatment we have implemented the energy equation self consistently. In particular, the adiabatic heating and cooling effects are implemented through the first and second term on the right hand side of equation (A), and the radiative heating and cooling rates are implemented through the last term on the right hand side of equation (A).

As discussed in section 4, our treatment of the line radiation pressure assumes negligible the $`\varphi `$ component of the line radiation force, and therefore the specific angular momentum is conserved within our models. We note that strictly, the azimuthal component of the line radiation force is not required to be zero, and thus the angular momentum is not necessarily conserved, but this is likely a small correction.

As we discovered in Paper 1 and Paper 2, care must be taken with the boundary conditions to avoid numerical instabilities. In this work we treat the boundary conditions in a form equivalent to that discussed in Paper 2. Our spatial grids are also equivalent to those explored and discussed in Paper 2.

In the models presented here we solve the hydrodynamic equations over a range of radii from $`0.02R_{\mathrm{}}`$ ($`2R_{wd}`$ in our models) to $`1.5R_{\mathrm{}}`$. From equation (6) it is found that the energy emission per area $`Q`$ will increase slightly as radius decreases from $`2R_{wd}`$ down to $`1.36R_{wd}`$ and will decrease strongly as radius decreases from $`1.36R_{wd}`$ down to the white dwarf radii ($`R_{wd}`$). In Paper 1 we showed that the wind density tends to decrease with decreasing radiation flux. Thus the mass loss rate for radii less that $`2R_{wd}`$ will not alter significantly the total mass loss rate of the wind. Also, in Paper 2, we presented models with spatial grids with “r” down to $`R_{wd}`$ and found essentially the same results with models with “r” down to $`2R_{wd}`$ with model parameters otherwise identical.

## 7 Results

In Table 1 we present the wind mass loss rates from our previous models and for the models presented here. The second column indicates the dimension of the model, 2.5D are the full 3D hydrodynamic equations under azimuthal symmetry; i.e. the derivatives with respect to the azimuthal angle are assumed zero. 2.5D models allow a consistent treatment of the 3D equations (under azimuthal symmetry) while mathematically having only two independent spatial directions, keeping the calculations more manageable. The third column indicates the radial disk structure assumed. In all models a steady disk is assumed. An “I” indicates an isothermal disk (no radial dependence). An “S” indicates the radial dependence of disk radiation emission and temperature of a standard Shakura-Sunyaev disk. The third and fourth columns indicate the values of the $`k`$ and $`\alpha `$ parameters used in each model. In our previous models (A-F) we had assumed constant ionization equilibrium throughout the wind, and therefore constant values for these two parameters. In our current work, as described above, we implement the spatial dependence (Model G) and the spatial and time dependence (Model H) of these two parameters and find values that vary throughout our spatial grids (and in time for Model H) within the intervals indicated in Table 1.

In implementing local ionization equilibrium in our models, we first run models that include the spatial dependence of the radiation field above the disk but assume a fiducial number density of $`10^8\mathrm{cm}^3`$ (Model G). From this we obtain the values of the “$`k`$” and “$`\alpha `$” line radiation parameters for each point of the computational grid. Results from these runs are presented in Figure 5. The values of $`k`$ we found varied from $`0.03`$ to $`0.5`$ and values of $`\alpha `$ varying from $`0.6`$ to $`0.8`$. The higher values of both $`k`$ and $`\alpha `$ are obtained in the lower inner region of the wind (lower left hand side of Figure 5). The lower values of both $`k`$ and $`\alpha `$ are obtained in the lower outer region of the wind (lower right hand side of Figure 5). At $`z=R_{\mathrm{}}`$, the $`k`$ parameter was found to vary from $`0.16`$ to $`0.26`$ and the $`\alpha `$ parameter was found to vary from $`0.66`$ to $`0.70`$. In this model we find similar terminal velocities with respect to the models presented in Paper 2, but densities values increased by a factor of 2. We find steady wind flows coming from the disk. Wind-mass loss rates calculated in this work are obtained by integrating the mass flux at the start of the wind, at the lower boundary of the spatial computational grid. The total wind mass loss rate we find for this model is $`2\times 10^{11}M_{\mathrm{}}\mathrm{yr}^1`$. This increase in the wind mass loss rate is due to the increase of the $`k`$ in the lower inner region where the radiation flux is stronger.

In order to calculate the $`k(r,z,t)`$ and $`\alpha (r,z,t)`$ functions self consistently throughout the hydrodynamic calculations, as we described in §5, the spatial distributions of the $`k`$ and $`\alpha `$ are calculated in a similar manner as in the model presented in Figure 5, but now for values of number density varying from $`10^6\mathrm{cm}^3`$ to $`10^{10}\mathrm{cm}^3`$ at each computational spatial grid point (Model H). A fiducial electron temperature of $`30,000\mathrm{K}`$ is taken for these calculations. During the hydrodynamic calculations, the actual values of $`k`$ and $`\alpha `$ are obtained through geometric interpolation from the above values and the current value of density at each computational spatial grid point. Results from these models are presented in Figure 6. For the above range of densities at each spatial grid point, we find values of $`k`$ varying from $`0.03`$ to $`0.6`$ and values of $`\alpha `$ varying from $`0.55`$ to $`0.85`$. Again we find steady wind flows coming from the disk. In the lower inner region, due to the increase of density the $`k`$ parameter drops by a factor of $`2`$, producing in turn a decrease of total wind mass rate. The wind mass loss rate for this model is $`4\times 10^{12}M_{\mathrm{}}\mathrm{yr}^1`$. The velocity fields are very similar to the previous models, finding terminal velocities in the order of $`3000\mathrm{km}\mathrm{s}^1`$. We note that the mass loss rate is a approximately a factor of two below the values we found in our earlier models in Paper 2. This difference is due to the detailed calculations of the line radiation force parameters, which when the local density values at the inner region are considered, we find slightly lower values of the $`k`$ line radiation force parameter.

In Figures 7 and 8 we present contour plots of the line radiation force parameter $`k`$ and $`\alpha `$ respectively, once the calculations arrive at a steady solution (Model H). These two figures illustrate the coupling between line force parameters, ionization equilibrium, and wind dynamics within an accretion disk wind. We note that, although the ionization structures are coupled with the hydrodynamics, within typical CV parameters ($`M_{wd}=0.6M_{\mathrm{}}`$; $`R_{wd}=0.01R_{\mathrm{}}`$; $`L_{disk}=L_{\mathrm{}}`$), in the higher wind density areas, the line force parameters do not change significantly; thus justifying the assumption of our earlier models of constant ionization equilibrium (and thus constant line force parameters). We also note that in QSOs/AGN, the presence of strong x-ray emission will generate strong changes in the ionization balance of the wind, and the assumption of constant ionization equilibrium in the wind will no longer hold (Murray et al., 1995; Hillier et al., 2002).

Studies of observed ultraviolet line profiles from CVs have estimated CV wind mass loss rates to be in the order of $`10^{12}`$ to $`10^{11}M_{\mathrm{}}\mathrm{yr}^1`$ and found CV wind terminal speeds in the order of $`3000\mathrm{km}\mathrm{s}^1`$ (e.g. Krautter et al., 1981; Córdova & Mason, 1982; Greenstein & Oke, 1982; Mauche & Raymond, 1987; Prinja & Rosen, 1995; Friedjung, Selvelli, & Cassatella, 1997) thus the models developed in this work finds values for wind mass loss rates and terminal speeds consistent with observations.

In the models presented in this work we also implement radiative heating and cooling. In Figure 10 we present a velocity field graph superimposed with temperature contours for the local ionization equilibrium disk wind model. Comparing with the results from Paper 2 (Figure 9), we find that radiative heating and cooling tends to maintain the wind temperature in the higher density regions near the disk, rather than cooling due to adiabatic expansion as we had found in our earlier models of Paper 2. Therefore the radiative heating and cooling processes are playing an important role in the temperature structure of the wind. We are finding higher overall temperatures than in Paper 2 where we neglected radiative heating and cooling processes and assumed only adiabatic heating cooling processes.

As a consistency check, in Figure 11 we present the results of models with a computational spatial grid that starts at $`r=R_{wd}`$ (rather than $`r=2.0R_{wd}`$ as in the other models discussed previously in this work). In Figure 12 we present results from models with a computational spatial grid starting at $`r=R_{wd}`$ and with an initial velocity at the base of the wind of $`1\mathrm{km}\mathrm{s}^1`$ (rather than $`10\mathrm{km}\mathrm{s}^1`$ as in the other models discussed previously in this work). In both cases we find similar overall results: a steady biconical wind velocity/density structure, with similar velocity and density values. We note that these results although similar, are not identical to those of Figure 6, and thus represent numerical effects since the physical parameters have remained identical. However the deviations from Figure 6 are considerably smaller than the changes of results we have found in earlier work when we physical parameters have been modified as our models have developed (Paper 2), thus we are confident in our overall results. In future work we shall study these numerical effects in greater detail; a precise knowledge of these model limitations may be important to achieve precise comparisons with observations.

As we had found in our earlier models of Paper 1 and Paper 2, rotational forces are important in the study of winds from accretion disks. They cause the velocity streamlines to collide which results in an enhanced density region. In this region the wind speed is reduced and the wind density is increased. The increase in density caused by the collision of streamlines is important because it permits the appearance of blue-shifted absorption lines as observed in P Cygni profiles of low-inclination CVs.

## 8 Summary and Conclusions

We have developed a 2.5-dimensional hydrodynamic line-driven accretion disk wind model. Our model solves a complete set of adiabatic hydrodynamic partial differential equations, using the PPM numerical scheme and implementing the radial temperature and radiation emission distributions on the surface of an accretion disk. Our models calculate the line radiation force parameters throughout the wind at each computational grid point for each time step through the calculation of line opacities from atomic data, similar to the calculations presented by Abbott (1982) for OB stars. We treat the hydrodynamic equation of energy self consistently including radiative heating and cooling as well as adiabatic expansion and compression.

We find steady wind flows coming from the disk with similar density and velocity structures as our earlier models of Paper 2 in which we had not yet included radiative heating and cooling and in which we had assumed constant ionization equilibrium (therefore constant values for the line radiation force parameters). The steady nature of our solution contrasts with the strong fluctuations in density and velocity found by Proga, Stone, & Drew (1998) for CV winds for physical parameters similar to the ones we used here. A possibility that could account for the difference in results from our models and those of Proga, Stone, & Drew (1998) is the treatment of the lower boundary conditions. In future work we shall explore this possibility in detail.

We note that as our numerical models evolve towards more realistic ones, it may result that a more accurate treatment of the line radiation force (which would include line-overlapping, rather than assuming the Sobolev approximation for all lines) or a more accurate model of the accretion disk (which may result in an unsteady picture of the accretion disk rather than the standard Shakura-Sunyaev disk) may lead to physical instabilities. In our systematic approach to the accretion disk wind problem, we started with an analytic 1D model, which although being our “crudest” one, did allow for an analytic solution which has proven to be invaluable as a test case for our numerical algorithms. From this initial model we have included, on a step by step basis (Paper 1, Pereyra 1997, Paper 2) elements such as a more accurate treatment of the line radiation force (in which we have always applied the Sobolev approximation), radial structure of the accretion disk, adiabatic heating and cooling processes, etc. This approach has allowed us to establish clearly a relationship between the elements we have introduced (both of numerical and physical nature) and our results, as our models have developed into more realistic ones. In this sense we are confident that if intrinsic physical instabilities should appear as we continue our work, we will be able to establish their physical origin.

From our models we calculate wind mass-loss rates and terminal velocities. For typical cataclysmic variable parameter values of disk luminosity $`L_{disk}=L_{\mathrm{}}`$, white dwarf mass $`M_{wd}=0.6M_{\mathrm{}}`$, and white dwarf radii $`R_{wd}=0.01R_{\mathrm{}}`$, we obtain a wind mass-loss rate of $`\dot{M}_{wind}=4\times 10^{12}M_{\mathrm{}}\mathrm{yr}^1`$, and a terminal velocity of $`3000\mathrm{km}\mathrm{s}^1`$. We note that the mass loss rate is a approximately a factor of two below the values we found in our earlier models in Paper 2. This difference is due to the detailed calculations of the line radiation force parameters, which when the local density values at the inner region are considered, we find slightly lower values of the $`k`$ line radiation force parameter.

The results in this work suggest that the details of the ionization equilibrium in CV disk winds, within reasonable parameters, may not have considerable effect over the velocity fields of the wind nor in the wind mass distribution. Although these details have to be taken into account in order to derive an accurate theoretical value of the wind mass loss rate as illustrated by the difference by a factor of two between the models presented here and the earlier ones of Paper 2. The similar wind mass distribution and velocity fields found here with respect to those of Paper 2, would indicate a similar inclination angle dependence of the C IV 1550 Å line. But the decrease in total wind mass loss rate would generate weaker line features. In QSOs, where the ionization balance of the gas material may change significantly, the coupling between the ion populations and the radiative hydrodynamics will have a important effects in the overall wind dynamics; thus a model capable of representing this coupling becomes necessary.

A shortcoming of our models is that, at each spatial point for each time step we assume direction-independent values for the line radiation parameters and we also assume position independent values for the radiative heating and cooling parameters. In future work we plan to include these effects. We also plan to extend our disk wind models to low mass X-ray binaries and QSOs/AGN, where the local ionization equilibrium will play a more crucial role in the overall dynamics.

We thank John Hillier and Dave Turnshek for their useful comments and discussions. We also thank the referee for useful comments and suggestions that have significantly improved the presentation of this work. This work was supported by he National Science Foundation under Grant AST-0071193.

## Appendix A Hydrodynamic Equations

The hydrodynamic equations implemented in the models presented are, the equation of state:

$$P=(\gamma 1)\rho e$$
(A1)
the mass conservation equation:

$$\frac{\rho }{t}+\frac{1}{r}\frac{(r\rho v_r)}{r}+\frac{(\rho v_z)}{z}=0$$
(A2)
the momentum conservation equations:

$`\rho {\displaystyle \frac{v_r}{t}}+\rho v_r{\displaystyle \frac{v_r}{r}}`$ $``$ $`\rho {\displaystyle \frac{v_{\varphi }^{}{}_{}{}^{2}}{r}}+\rho v_z{\displaystyle \frac{v_r}{z}}=`$
$``$ $`\rho {\displaystyle \frac{GM_{wd}}{(r^2+z^2)}}{\displaystyle \frac{r}{(r^2+z^2)^{1/2}}}{\displaystyle \frac{P}{r}}+\rho {\displaystyle \frac{\sigma F_r(r,z)}{c}}`$
$`+`$ $`\rho {\displaystyle \frac{\sigma S_r(r,z)}{c}}k(r,z,t)\left(\mathrm{max}[{\displaystyle \frac{1}{\rho \sigma v_{th}}}\left|{\displaystyle \frac{dv_z}{dz}}\right|,10^7]\right)^{\alpha (r,z,t)}`$
$$\rho \frac{v_\varphi }{t}+\rho v_r\frac{v_\varphi }{r}+\rho \frac{v_\varphi v_r}{r}+\rho v_z\frac{v_\varphi }{z}=0$$
(A4)
$`\rho {\displaystyle \frac{v_z}{t}}+\rho v_r{\displaystyle \frac{v_z}{r}}`$ $`+`$ $`\rho v_z{\displaystyle \frac{v_z}{z}}=`$
$``$ $`\rho {\displaystyle \frac{GM_{wd}}{(r^2+z^2)}}{\displaystyle \frac{z}{(r^2+z^2)^{1/2}}}{\displaystyle \frac{P}{z}}+\rho {\displaystyle \frac{\sigma F_z(r,z)}{c}}`$
$`+`$ $`\rho {\displaystyle \frac{\sigma S_z(r,z)}{c}}k(r,z,t)\left(\mathrm{max}[{\displaystyle \frac{1}{\rho \sigma v_{th}}}\left|{\displaystyle \frac{dv_z}{dz}}\right|,10^7]\right)^{\alpha (r,z,t)}`$
and the energy conservation equation:

$`{\displaystyle \frac{\rho E}{t}}+{\displaystyle \frac{1}{r}}{\displaystyle \frac{r\rho Ev_r}{r}}`$ $`+`$ $`{\displaystyle \frac{\rho Ev_z}{z}}=`$
$``$ $`{\displaystyle \frac{1}{r}}{\displaystyle \frac{rPv_r}{r}}{\displaystyle \frac{Pv_z}{z}}`$
$``$ $`\rho {\displaystyle \frac{GM_{wd}}{(r^2+z^2)}}{\displaystyle \frac{r}{(r^2+z^2)^{1/2}}}v_r+\rho {\displaystyle \frac{\sigma F_r(r,z)}{c}}v_r`$
$`+`$ $`\rho {\displaystyle \frac{\sigma S_r(r,z)}{c}}k(r,z,t)\left(\mathrm{max}[{\displaystyle \frac{1}{\rho \sigma v_{th}}}\left|{\displaystyle \frac{dv_z}{dz}}\right|,10^7]\right)^{\alpha (r,z,t)}v_r`$
$``$ $`\rho {\displaystyle \frac{GM_{wd}}{(r^2+z^2)}}{\displaystyle \frac{z}{(r^2+z^2)^{1/2}}}v_z+\rho {\displaystyle \frac{\sigma F_z(r,z)}{c}}v_z`$
$`+`$ $`\rho {\displaystyle \frac{\sigma S_z(r,z)}{c}}k(r,z,t)\left(\mathrm{max}[{\displaystyle \frac{1}{\rho \sigma v_{th}}}\left|{\displaystyle \frac{dv_z}{dz}}\right|,10^7]\right)^{\alpha (r,z,t)}v_z`$
$`+`$ $`n^2(\mathrm{\Gamma }(T)\mathrm{\Lambda }(T))`$
where $`P`$ is the pressure, $`\gamma `$ is the ratio of specific heats, $`\rho `$ is the density, $`e`$ is the internal energy per mass, $`v_r`$, $`v_\varphi `$, and $`v_z`$ are the corresponding velocity components in cylindrical coordinates, $`G`$ is the gravitational constant, $`M_{wd}`$ is the mass of the white dwarf, $`\sigma `$ is the Thomson cross section per unit mass, $`F_r`$ and $`F_z`$ are the corresponding radiation flux components (see equations  and ), $`c`$ is the speed of light, $`v_{th}`$ is the thermal velocity, and $`k`$ and $`\alpha `$ are the line force multiplier parameters calculated through equation (30), $`S_r`$ and $`S_z`$ are defined through equations (38) and (37) respectively, $`E=v_r^2/2+v_\varphi ^2/2+v_z^2/2+e`$ is the total energy per mass, $`n`$ is the electron number density, $`\mathrm{\Gamma }`$ is the radiative heating rate, and $`\mathrm{\Lambda }`$ is the radiative cooling rate.