A&A 377, 175-191 (2001)
DOI: 10.1051/0004-6361:20011075

Multicomponent radiatively driven stellar winds

II. Gayley-Owocki heating in multitemperature winds of OB stars

J. Krticka1,2 - J. Kubát2


1 - Ústav teoretické fyziky a astrofyziky PrF MU, Kotlárská 2, 611 37 Brno, Czech Republic
2 - Astronomický ústav, Akademie ved Ceské republiky, 251 65 Ondrejov, Czech Republic

Received 29 May 2001 / Accepted 24 July 2001

Abstract
We show that the so-called Gayley-Owocki (Doppler) heating is important for the temperature structure of the wind of main sequence stars cooler than the spectral type O6. The formula for Gayley-Owocki heating is derived directly from the Boltzmann equation as a direct consequence of the dependence of the driving force on the velocity gradient. Since Gayley-Owocki heating deposits heat directly on the absorbing ions, we also investigated the possibility that individual components of the radiatively driven stellar wind have different temperatures. This effect is negligible in the wind of O stars, whereas a significant temperature difference takes place in the winds of main sequence B stars for stars cooler than B2. Typical temperature differences between absorbing ions and other flow components for such stars is of the order $10^3\,{\rm K}$. However, in the case when the passive component falls back onto the star, the absorbing component reaches temperatures of order $10^6\,{\rm K}$, which allows for emission of X-rays. Moreover, we compare our computed terminal velocities with the observed ones. We found quite good agreement between predicted and observed terminal velocities. The systematic difference coming from the using of the so called "cooking formula'' has been removed.

Key words: hydrodynamics - stars: mass-loss - stars: early-type - stars: winds, outflows


1 Introduction

Since the founding of the theory of radiatively-driven stellar winds by Lucy & Solomon (1970) and Castor et al. (1975, hereafter CAK) many of the initial assumptions introduced by these authors were examined. To the most important ones belong the radial streaming approximation (Friend & Abbott 1986; Pauldrach et al. 1986), the wind stability (Abbott 1980; Owocki & Rybicki 1984), the limitations of the Sobolev approach (Poe et al. 1990; Owocki & Puls 1999), the thermal structure of the wind (Drew 1989) and many others.

Another important assumption, studied already at the dawn of the radiatively-driven stellar wind theory by Castor et al. (1976) is the condition of the one-component flow. They discussed encounters which transfer momentum received by absorbing ions (typically C, N, O, etc.) to passive, nonabsorbing ions, mainly hydrogen and helium. They showed that for the high-density winds, such encounters are not important for the overall dynamics of the wind and that high-density winds can be considered as one-component. However, for the low-density winds Springmann & Pauldrach (1992, hereafter SP) showed that momentum transfer between absorbing and nonabsorbing plasma influences the wind thermal balance and even the wind dynamics. Thus, for the low-density winds the flow is essentially multicomponent. They proposed that the so-called "ion-runaway" may occur. Based on the simplified theory of the multicomponent flow many interesting results were seen. Porter & Drew (1995) re-examined the model of a wind-compressed disk in the presence of dynamical decoupling of absorbing ions and passive plasma, Porter & Skouza (1999) showed the possibility of formation of pulsating shells around stars with low-density radiatively driven wind, and Hunger & Groote (1999) explained the H/He abundance anomalies in Bp stars on the basis of helium decoupling.

The first detailed numerical models of multicomponent radiatively driven stellar winds were presented by Babel (1995, 1996). However, Krticka & Kubát (2000, hereafter KK0) showed that due to the functional dependence of the radiative force decoupling does not occur. Moreover, Krticka & Kubát (2001, hereafter KKI) using nonisothermal multicomponent models, concluded that winds of B stars are frictionally heated such that the possibility of decoupling of absorbing ions from the passive plasma is excluded.

The solar wind is well-known to possess large temperature differences between electrons and protons. Such differences were obtained also by Bürgi (1992), who used the three-component models of the solar wind. So the natural question arises, whether similar temperature differences exist in the radiatively-driven stellar wind or, in other words, whether the assumption of equal temperatures of all wind components is acceptable. In this paper we intend to answer this question.

Any effect which deposits heat separately on to individual component of the flow may influence our results. Thus, we shall include the effect of Doppler heating, introduced in the stellar wind domain by Gayley & Owocki (1994, hereafter GO). Because it arises from the dependence of the radiative force on the velocity via the Doppler effect it deposits heat directly to the absorbing ion component and thus, it can trigger the temperature difference between absorbing and passive ions.

Proper treatment of ionization balance may be important for the correct description of decoupling of individual components of the flow. Thus, we decided to compute electrical charges of individual components using adequate ionization balance formulas.

2 Boltzmann equation with velocity-dependent force

The procedure of the derivation of the hydrodynamic equations from the Boltzmann equation for particle distribution function Fs of the particle s is thoroughly described in a number of textbooks. However, it is commonly assumed that the Boltzmann equation can be written in the form (we use the Einstein summation law)

$\displaystyle \frac{\partial F_s}{\partial t} + \xi_{sh}\frac{\partial F_s}{\pa... ..._s}{\partial \xi_{sh}} = \left(\frac{{\rm d}F_s}{{\rm d}t}\right)_{{\rm coll}},$     (1)

i.e. it is assumed that the force fs is independent of velocity, which is acceptable when we consider gravitational or electrical forces. However, the latter assumption is not valid for the case of a wind driven by radiative force coming from line absorption, which is strongly dependent on a velocity gradient. Therefore, we rederive here hydrodynamic equations without an assumption of the force independent of velocity. In this case, the Boltzmann equation for the one-particle distribution function Fs of particles of type s is written as
 
$\displaystyle \frac{\partial F_s}{\partial t} +\xi_{sh} \frac{\partial F_s}{\pa... ...h}}{m_s}F_s\right) = \left(\frac{{\rm d}F_s}{{\rm d}t}\right)_{{\rm coll}}\cdot$     (2)

Here $\xi_{sh}$ (h=1,2,3) are the velocity components of individual particles of a type s with mass ms, and fsh are the components of an external force acting on them. The right-hand side term expresses the effect of collisions. By a definition, the integral of the distribution function over the velocity space is a number density ns of s-particles,
 
$\displaystyle n_s=\int {\rm d}\vec{\xi}_s F_s.$     (3)

Now, the usual way to obtain hydrodynamic equations is the following. One multiplies the Boltzmann Eq. (2) by multipliers ms, $m_s\xi_{sh}$ and $m_s\xi_{sh}\xi_{sk}$and integrates it over the velocity space. For the discussion of the role of velocity dependent force in the Boltzmann equation we confine us to the left-hand side of the Boltzmann equation. The right-hand side (i.e. the collisional term) remains unaffected by the presence of such forces, so we assume in this section that the gas is collisionless, i.e. that the right hand side of the Boltzmann equation is zero.

2.1 Continuity equation

Multiplying the Boltzmann Eq. (2) by ms and integrating over the velocity space we obtain the continuity equation

 
$\displaystyle \frac{\partial }{\partial t}\left(n_s m_s\right)+\frac{\partial }{\partial x_k}\left(n_s m_s v_{sk}\right)=0,$     (4)

where vsk are the components of the mean velocity of particles s, and
$\displaystyle n_s m_s v_{sk}= m_s \int {\rm d}\vec{\xi}_s \xi_{sk} F_s.$     (5)

Here, the force term disappeared through integration by parts and the collisional term is also zero. Apparently, velocity-dependent external force does not change the continuity equation.

2.2 Momentum equation

Multiplying the Boltzmann Eq. (2) by $m_s\xi_{sh}$ and integrating over the velocity space we obtain momentum equation

 
$\displaystyle \frac{\partial }{\partial t}\left(n_s m_s v_{sh}\right)+ \frac{\p... ...(n_s m_s v_{sh} v_{sk}+p_{s,hk}\right)- \int {\rm d}\vec{\xi}_s f_{sh} F_s = 0,$     (6)

where ps,hk are the components of the momentum transfer tensor, and
$\displaystyle n_s m_s v_{sh} v_{sk}+p_{s,hk} = m_s \int {\rm d}\vec{\xi}_s \xi_{sh} \xi_{sk} F_s.$     (7)

Again, the velocity term was integrated by parts. Note, that due to the velocity dependence of the external force the term containing fshcannot be moved in front of the integral.

2.3 Energy equation

Multiplying the Boltzmann Eq. (2) by $m_s\xi_{sk}\xi_{sk}$and integrating over the velocity space we obtain energy equation

 
$\displaystyle \frac{\partial }{\partial t}\left(\frac{1}{2}n_s m_s v_s^2+\frac{... ...s,ki}+\frac{1}{2}p_{s,kki}\right)- \int {\rm d}\vec{\xi}_s \xi_{sk}f_{sk}F_s=0,$     (8)

where ps=1/3 ps,kk is the scalar hydrostatic pressure and the last term was simplified using integration by parts. The equation for the temperature can be derived by subtracting the momentum Eq. (6) dot-multiplied by vsh from the Eq. (8),
 
$\displaystyle \frac{\partial }{\partial t}\left(\frac{3}{2} p_s\right)+\frac{\p... ...int {\rm d}\vec{\xi}_s f_{sh} F_s-\int {\rm d}\vec{\xi}_s \xi_{sh}f_{sh}F_s =0.$     (9)

Note that the last two terms cancel if the force does not depend on the velocity. For this "standard'' case the external force does not give any contribution to net heating.

2.3.1 Gayley-Owocki (Doppler) heating

Let us explore now the last term on the left hand side of Eq. (9) for the case of the force caused by absorption of radiation in spectral lines. This force depends on particle velocity through the velocity dependence of the line absorption coefficient owing to the Doppler effect. Therefore, Gayley & Owocki (1994, hereafter GO) termed the heating effect by Doppler heating, but terming it Gayley-Owocki heating (or GO heating in the abbreviated form) might be more appropriate. Let us denote the heating term in the comoving fluid-frame as

$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \int {\rm d}{\vec{\xi}}_{{\rm i}} \xi_{{\rm i},h}f_{{\rm i},h}{F}_{{\rm i}}.$     (10)

Here the index ${\rm i}$stands for "absorbing ions''. Although this term is written in the comoving fluid-frame, the same expression holds in the non-relativistic case also in the observer frame. Let us assume complete redistribution and an angle independent opacity and emissivity in the atomic frame. However, the emissivity in the comoving fluid-frame is generally angle dependent due to the Doppler effect. Therefore, the radiative force acting on an atom in the comoving fluid-frame is
 
$\displaystyle f_{{\rm i},h} = \frac{m_{\rm i}}{c} \int_0^\infty {\rm d} \nu \oi... ...ga \left[\kappa(\vec{n},\nu) I(\vec{n},\nu) - \epsilon(\vec{n},\nu)\right] n_h,$     (11)

where $\kappa({\vec n},\nu)$ and $\epsilon({\vec n},\nu)$ are the absorption and emission coefficients per unit mass, which can be expressed as


$\displaystyle \kappa({\vec n},\nu) = \kappa\, \varphi\left(\nu-\frac{\nu_0}{c}{\vec{\xi}}_{{\rm i}} \!\cdot\!{\vec n}\right),$     (12a)

$\displaystyle \epsilon({\vec n},\nu) = \epsilon\, \varphi\left(\nu-\frac{\nu_0}{c}{\vec{\xi}}_{{\rm i}} \!\cdot\!{\vec n}\right),$     (12b)


where $\varphi(\nu)$ is the absorption (emission) profile in the atomic frame. After some rearrangement we obtain

$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \frac{{m}_{{\rm i}}}{c} \int_{0}^{\inf... ...\varphi\left(\nu-\frac{\nu_0}{c}{\vec{\xi}}_{{\rm i}} \!\cdot\!{\vec n}\right).$     (13)

If for the calculation of the last integral we choose one of the velocity axes parallel to the direction of n, then the integrals over other two axes vanishes (we integrate odd function) and the Gayley-Owocki heating formula becomes
$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \frac{{m}_{{\rm i}}}{c} \int_{0}^{\inf... ...m i}}({\xi}_{{\rm i}}) \varphi \left(\nu-\frac{\nu_0}{c}{\xi}_{{\rm i}}\right).$     (14)

where ${\xi}_{{\rm i}}$is a velocity component in the direction of ${\vec n}$. If we rewrite photon-line-of-sight velocity component as ${\xi}_{{\rm i}}=wv_{{\rm th,i}}$and use frequency displacement from the line center in Doppler units $x=\left(\nu-\nu_0\right)/\Delta\nu_{\rm D}$, then the GO heating takes the form
$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \frac{{m}_{{\rm i}}}{c} \int_{-\infty}... ...nfty v_{{\rm th,i}}^2\,{\rm d}w\,w\, {F}_{{\rm i}}(v_{{\rm th,i}}w)\,\psi(w-x),$     (15)

where we introduced the intensity $\tilde{I}({\vec n},x)$ as $\tilde{I}({\vec n},x)\, {\rm d}x = I({\vec n},\nu)\, {\rm d}\nu$ and the function
$\displaystyle \psi(x)=\varphi(x\Delta\nu_{\rm D})$     (16)

is normalized according to Castor (1974) as
$\displaystyle \int_{-\infty}^{\infty} \psi(x) {\rm d}x=1.$     (17)

Note that the thermal speed $v_{{\rm th,i}}$ is really ionic because it comes from the velocity distribution of absorbing ions (shall not be interchanged with $v_{{\rm th}}$, which comes from normalization of force multipliers). We neglect absorption in the resonance wings of the profile and approximated $\psi(w-x)\approx\delta(w-x)$. Finally, we assume that the velocity distribution is given by the Maxwellian velocity distribution,
$\displaystyle {F}_{{\rm i}}(v_{{\rm th,i}}w) =\frac{{n}_{{\rm i}}}{v_{{\rm th,i}}\sqrt\pi}\,{\rm e}^{-w^2}.$      

The latter assumption was made purely due to simplicity reasons. Relaxing it could lead to interesting effects especially if the number of collisions is not sufficient to maintain an equilibrium (cf., e.g., Scudder 1994; Cranmer 1998). However, we postpone the analysis of the non-Maxwellian effects to a future paper.

Thus, the GO heating formula takes the form of

$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \frac{{m}_{{\rm i}} {n}_{{\rm i}} v_{{... ...\,x \phi(x) \left[\kappa \tilde I({\vec n},x)-\Delta\nu_{\rm D}\epsilon\right],$     (18)

where
 
$\displaystyle \phi(x)=\frac{1}{\sqrt\pi}{\rm e}^{-x^2}.$     (19)

Due to the symmetry of the absorption profile (19) the product $x \phi(x)$ is an odd function and, thus, after integration over x, vanishes. This means that in the case of complete redistribution the process of emission gives no direct contribution to GO heating, i.e.
$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \frac{{m}_{{\rm i}} {n}_{{\rm i}} v_{{... ...ega\, \int_{-\infty}^{\infty} {\rm d}x \,x \phi(x) \kappa \tilde I({\vec n},x).$     (20)

In the static medium the product $x \phi(x) \tilde I({\vec n},x)$ is an odd function of x, thus, there is no GO heating effect in static stellar atmospheres. In the particular case of a spherically-symmetric stellar wind the expression for the Gayley-Owocki heating takes the form of
 
$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \frac{2\pi\kappa{\rho}_{{\rm i}}v_{{\r... ..._{-1}^1{\rm d}\mu \int_{-\infty}^{\infty} {\rm d}x\, x \phi(x) \tilde I(\mu,x),$     (21)

where $\mu=\cos\theta$, which was actually used by GO.

2.3.2 Formula for Gayley-Owocki heating in the stellar wind domain

In the case of a two-level atom without continuum the solution of the transfer equation in the Sobolev approximation is (Rybicki & Hummer 1978; Owocki & Rybicki 1985, GO)

 
$\displaystyle \tilde I(\mu,x)=\tilde I_{\rm c} \left\{\frac{\beta_{\rm c}}{\beta}+\left[D(\mu)- \frac{\beta_c}{\beta}\right] {\rm e}^{-\tau_\mu\Phi(x)}\right\},$     (22)

where $\tilde I_{\rm c}$ is the core intensity, $D(\mu)$ is unity for $\mu>\mu_*$ and zero otherwise ( $\mu_*=\left(1-R_*^2/r^2\right)^{1/2}$), core penetration and escape probabilities are given by
  
$\displaystyle \beta_{\rm c}=\frac{1}{2}\int_{\mu_*}^{1}{\rm d}\mu\frac{1-{\rm e}^{-\tau_\mu}}{\tau_\mu},$     (23)
$\displaystyle \beta=\frac{1}{2}\int_{-1}^{1}{\rm d}\mu\frac{1-{\rm e}^{-\tau_\mu}}{\tau_\mu},$     (24)

respectively, and
$\displaystyle \Phi(x)=\int_{x}^{\infty}{\rm d}x'\phi(x').$     (25)

The Sobolev optical depth $\tau_\mu$ is given by (Castor 1974; Rybicki & Hummer 1978)
 
$\displaystyle \tau_\mu=\frac{{\rho}_{{\rm i}}\kappa v_{{\rm th}}r}{{\mbox{\eufont Y}}_{{\rm i}}{{v_{\rm r}}}_{{\rm i}} \left(1+\sigma\mu^2\right)},$     (26)

where the variable $\sigma $ was introduced by Castor (1974)
 
$\displaystyle \sigma=\frac{{\rm d}\ln {{v_{\rm r}}}_{{\rm i}} }{{\rm d}\ln r} -1.$     (27)

Inserting the solution of the transfer Eq. (22) into the expression for the GO heating Eq. (21) we obtain
 
$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \frac{2\pi\kappa{\rho}_{{\rm i}}v_{{\r... ...right] \int_{-\infty}^{\infty} {\rm d}x\, x \phi(x) {\rm e}^{-\tau_\mu\Phi(x)}.$     (28)

The effect of line ensemble is usually described using the concept of a line-strength distribution function (CAK, Abbott 1982; Puls et al. 2000)
 
$\displaystyle {\rm d}N(\kappa)=-N_0 \left(\frac{{\rho}_{{\rm e}}/\left(W {m}_{{... ...cm}^{-3}}\right)^{\delta} \kappa^{\alpha-2}{\rm d}\kappa\frac{{\rm d}\nu}{\nu},$     (29)

where normalization constant N0 is taken in the form of
 
$\displaystyle N_0=\frac{c k}{v_{{\rm th}}}\left(1-\alpha\right)\alpha\, \sigma_{\rm e}^{1-\alpha}.$     (30)

The GO heating formula for this line ensemble can be obtained by the integration of the heating term for one line (Eq. (28)) over the CAK distribution function Eq. (29). In this case it takes the form of
$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \frac{{\rho}_{{\rm i}}v_{{\rm th}}v_{{... ...right] \int_{-\infty}^{\infty} {\rm d}x\, x \phi(x) {\rm e}^{-\tau_\mu\Phi(x)},$     (31)

where $ \tilde I_{\rm c} =\Delta\nu_{\rm D}L/4\pi^2 R_*^2$ was used. Finally, applying substitution $y=\kappa{\rho}_{{\rm i}}v_{{\rm th}}r/{\mbox{\eufont Y}}_{{\rm i}}{{v_{\rm r}}}_{{\rm i}} $preceding equation becomes
 
$\displaystyle {Q}_{{\rm i}}^{{\rm GO}} = \frac{{\rho}_{{\rm i}}v_{{\rm th,i}}k ... ...}{\sigma_{\rm e}{\rho}_{{\rm i}}v_{{\rm th}}r}\right)^{\alpha} G(\sigma,\mu_*),$     (32)

where the function $G(\sigma,\mu_*)$ is given by the triple integration
 
$\displaystyle G(\sigma,\mu_*)=\int_0^{\infty}{\rm d}y\,y^{\alpha-1}\int_{-1}^1{... ...\infty} {\rm d}x\, x \phi(x) \exp\left(-\frac{y\Phi(x)}{{1+\sigma\mu^2}}\right)$     (33)

and in the integrals for $\beta_{\rm c}$ and $\beta $ Eqs. (23), (24) the Sobolev depth shall be computed using
$\displaystyle \tau_\mu=\frac{y}{1+\sigma\mu^2}$     (34)

instead of Eq. (26). Contrary to the radiative force formula (Castor 1974) the GO heating formula depends on the absorption profile. For the determination of the GO heating term we selected Gaussian profile (which comes from the Maxwellian velocity distribution).

  
3 Model equations

3.1 Basic assumptions

We assume that the stationary, spherically symmetric stellar wind consists of three components, namely absorbing ions, nonabsorbing hydrogen atoms and ions, and electrons, denoted by subscripts ${\rm i}$, ${\rm p}$, ${\rm e}$, respectively. Each of them is described by a density $\rho_a$, radial velocity ${ {v_{\rm r}}}_a$, temperature Ta, electrical charge qa=e za(where e is an elementary charge and za denotes the ionization degree - may have a non-integer value), and particle mass ma. Subscript a stands for $a={\rm i,p,e}$. Contrary to our previous models (KKI), we allow for different temperature of each component and for radial changes of electrical charge. We assume that chemical composition is given by the factor z*, which is a stellar metallicity relative to the solar value.

3.2 Continuity equations

In the case of a stationary spherically symmetric stellar wind each component is described by the continuity Eq. (4) in the form of



$\displaystyle \frac{1}{r^2}\frac{{\rm d}}{{\rm d} r}\left(r^2\rho_a{ {v_{\rm r}}}_a\right) = S_a, $     (35a)



where the term Sa accounts for radial change of mass-loss rate of individual components due to the ionization. Whereas for all types of ions the mass-loss rate is constant through the wind and thus number of these particles is conserved ( ${S}_{{\rm i}}={S}_{{\rm p}} =0$), we account for the possibility of variation of electron number via ionization and recombination. Because the total electric charge is conserved,

$\displaystyle \sum_a \frac{{\rm d}}{{\rm d} r}\left(r^2q_a\frac{\rho_a}{m_a}{ {v_{\rm r}}}_a\right)=0$      

and continuity equation holds separately for all components except electrons, we obtain electron continuity equation in the form of



$\displaystyle \frac{1}{r^2}\frac{{\rm d}}{{\rm d} r}\left(r^2{\rho}_{{\rm e}}{{... ...{{m}_{{\rm e}}}{m_a} \rho_a{ {v_{\rm r}}}_a \frac{{\rm d} z_a}{{\rm d}r}\cdot $     (35b)



Although inclusion of a term ${S}_{{\rm e}}$ into the electron continuity equation does not significantly alter the model, it is important to obtain well converged model.

3.3 Momentum equations

In the case of stationary spherically symmetric stellar wind the momentum equation Eq. (6) has the form of

 
$\displaystyle { {v_{\rm r}}}_a\frac{{\rm d} { {v_{\rm r}}}_a}{{\rm d}r}= {g}_{a... ... {v_{\rm r}}}_b-{ {v_{\rm r}}}_a}{\vert{ {v_{\rm r}}}_b-{ {v_{\rm r}}}_a\vert},$     (36)

where square of isothermal sound speed is aa2=kTa/ma, E is a charge separation electric field. Gravitational acceleration has the form $g={G{\mbox{\eufont M}}}/{r^2}$, where $\mbox{\eufont M}$ is the stellar mass and G is the gravitational constant. The radiative acceleration acting on free electrons can be expressed as
$\displaystyle {g}_{{\rm e}}^{{\rm rad}}=\frac{{m}_{{\rm p}}}{{m}_{{\rm e}}}\frac{G\Gamma{\mbox{\eufont M}}}{r^2},$     (37)

where the ratio of the radiative force caused by absorption of radiation by free electrons and gravitational force is
$\displaystyle \Gamma=\frac{\sigma_{\rm e} L}{4\pi c G{\mbox{\eufont M}}}$      

L is the stellar luminosity and $\sigma_{\rm e}$ is the mass scattering coefficient of the free electrons (do not confuse it with Thomson scattering cross section).

The radiative acceleration acting on absorbing ions is taken in the form of Castor et al. (1975)

 
$\displaystyle {g}_{{\rm i}}^{{\rm rad}}= \frac{1}{{\mbox{\eufont Y}}_{{\rm i}}}... ...}_{{\rm i}}} \frac{{\rm d}{ {v_{\rm r}}}_{{\rm i}}}{{\rm d}r} \right)^{\alpha},$     (38)

with force multipliers k, $\alpha $, $\delta $ after Abbott (1982). The finite disk correction factor (Friend & Abbott 1986; Pauldrach et al. 1986) is
$\displaystyle f=\frac{\left(1+\sigma\right)^{\alpha+1}- \left(1+\sigma\mu_{\rm ... ...lpha+1\right)\left(1-\mu_{\rm c}^2\right)\sigma\left(1+\sigma\right)^{\alpha}},$     (39)

where W is a dilution factor. The thermal speed $v_{{\rm th}}=\sqrt{2k{T}_{{\rm i}}/{m}_{{\rm p}}}$ is, owing to the Abbott's normalization, hydrogen thermal speed. However, because $v_{{\rm th}}$ physically describes ionic thermal speed, it depends on ionic temperature. Finally, ${\mbox{\eufont Y}}_{{\rm i}}$ is the photospheric ratio of the metallic ion density to the passive plasma density. We selected the same value as in KKI, namely ${\mbox{\eufont Y}}_{{\rm i}}=0.0127$ (this value corresponds to the solar ratio of sum of densities of C, N, O, Fe to the density of bulk plasma).

Constant of friction evaluated using Fokker-Planck approximation (cf. Burgers 1969) has the following form:

$\displaystyle {K}_{ab}={n}_a{n}_b\frac{4\pi {q}_a^2{q}_b^2}{k T_{ab}}\ln\Lambda,$     (40)

where na and nb are number densities of individual components and mean temperature of both components
$\displaystyle T_{ab}=\frac{m_aT_b+m_bT_a}{m_a+m_b}\cdot$     (41)

The Coulomb logarithm is of the form
$\displaystyle \ln\Lambda= \ln\left[\frac{3k{T}_{{\rm e}}}{e^2}\left(\frac{k{T}_{{\rm e}}}{4\pi n e^2}\right)^{1/2}\right],$     (42)

where n is the particle density ( $n={n}_{{\rm p}}+{n}_{{\rm e}}+{n}_{{\rm i}}$). Finally, the Chandrasekhar function G(x), defined in terms of the error function ${\rm erf}(x)$ (Dreicer 1959) is
$\displaystyle G(x)=\frac{1}{2 x^2}\left({\rm erf}(x)-\frac{2x}{\sqrt\pi}\exp\left(-x^2\right)\right).$     (43)

The argument of the Chandrasekhar function is
$\displaystyle x_{ab}=\frac{\vert{ {v_{\rm r}}}_b-{ {v_{\rm r}}}_a\vert}{\alpha_{ab}},$     (44)

where
$\displaystyle \alpha_{ab}^2=\frac{2k\left(m_aT_b+m_bT_b\right)} {m_am_b}\cdot$     (45)

3.4 Energy equation

Energy Eq. (9) in the case of a stationary, spherically symmetric multicomponent flow has the form of (cf. Burgers 1969)

 
$\displaystyle \frac{3}{2}k{ {v_{\rm r}}}_a\frac{\rho_a}{m_a}\frac{{\rm d}T_a}{{\rm d}r}$ + $\displaystyle a_a^2\rho_a \frac{1}{r^2}\frac{{\rm d}}{{\rm d}r}\left(r^2 { {v_{\rm r}}}_a\right)= Q_a^{{\rm rad}}$  
  + $\displaystyle \frac{1}{\sqrt\pi}\sum_{b\neq a}K_{ab} \frac{2k\left(T_b-T_a\right)}{m_a+m_b} \frac{\exp\left(-x_{ab}^2\right)}{\alpha_{ab}}$  
  + $\displaystyle \sum_{b\neq a}\frac{m_b}{m_a+m_b}K_{ab}G(x_{ab})\vert{ {v_{\rm r}}}_b-{ {v_{\rm r}}}_a\vert.$ (46)

Two terms on the left-hand side stand for advection and adiabatic cooling, respectively. The right-hand side terms describe radiative heating/cooling, heat exchange by encounters of different particles caused by unequal temperatures of the components, and frictional heating.

There are two sources of radiative heating/cooling. First source are bound-free and free-free transitions and the second is Gayley-Owocki heating/cooling.

3.4.1 Radiative energy term of bound-free and free-free transitions

Bound-free and free-free transitions (which will be called "classical'' radiative transitions) deposit energy directly on electrons. Therefore, this classical radiative energy term should be considered in an electron energy equation. We decided to estimate the radiative heating/cooling term $Q^{{\rm rad}}$ using two mechanisms only, hydrogen Lyman bound-free and free-free transitions. The detailed form of heating and cooling in the above mentioned transitions is nearly the same as in KKI and will not be repeated here (see also Kubát et al. 1999). The only difference is that the temperature in these equations is now the electron temperature. $J_{\nu}$ at the base of the wind is taken as an emergent radiation from a spherically symmetric static hydrogen model atmosphere for a corresponding stellar type (Kubát 2001).

3.4.2 Gayley-Owocki heating/cooling

Contrary to bound-free and free-free transitions Gayley-Owocki heating/cooling deposits energy directly to absorbing ions. GO heating/cooling term has the following form (GO, Eq. (32))

 
$\displaystyle {Q}_{{\rm i}}^{{\rm rad}}\equiv{Q}_{{\rm i}}^{{\rm GO}} ={\rho}_{... ...}{\sigma_{\rm e}{\rho}_{{\rm i}}v_{{\rm th}}r}\right)^{\alpha} G(\sigma,\mu_*),$     (47)

where the function $G(\sigma,\mu_*)$ is given by Eq. (33). Here $v_{{\rm th,i}}$ is a thermal speed of driving ions, thus, $v_{{\rm th,i}}=\sqrt{2k{T}_{{\rm i}}/{m}_{{\rm i}}}$. Inspecting the GO heating term (47) we conclude that wind is heated via the Doppler effect by direct radiation whereas the Gayley-Owocki cooling is caused by the diffuse radiation. The sign of GO heating depends on the value of a variable $\sigma $. In the case where $\sigma=0$, when the expansion is locally isotropic, the GO heating and cooling vanishes. When the $\sigma $ is positive, the cooling term dominates and the wind is cooled by the GO heating. In the opposite case, when the $\sigma $ is negative, the heating term dominates and thus, the wind is heated.

3.5 Charge separation electric field

The equation for charge separation electric field can be obtained directly from the third Maxwell equation, which in the case of spherical symmetry can be written as

 
$\displaystyle \frac{1}{r^2}\frac{{\rm d}}{{\rm d}r}\left(r^2 E\right) =4\pi\sum_aq_an_a.$     (48)

3.6 Ionic charge

The ionization structure of stellar wind should be derived using time consuming NLTE calculations (e.g. Pauldrach et al. 1994). Because we want to determine only a mean charge of selected elements, we use a simpler approximate method. As described by Mihalas (1978, Eq. (5.46) therein), the ionization equilibrium in stellar winds can be approximated by

 
$\displaystyle \frac{n_{a,j}}{n_{a,j+1}}\approx \frac{1}{2} \frac{{n}_{{\rm e}}}... ...R}}{{T}_{{\rm e}}}\right)^{1/2}\exp\left(\frac{\chi_{a,j}}{k T_{\rm R}}\right),$     (49)

where $\chi_{a,j}$ is the ionization potential, Ua,j is the partition function and $T_{\rm R}$ is the radiation temperature (we set $T_{\rm R}=\frac{3}{4}T_{\rm eff}$). Partition function approximations are taken from Smith & Dworetsky (1988). Electrical charge of absorbing and passive ions are then derived using a formula
$\displaystyle z_a=\frac{\sum_j j n_{a,j}}{\sum_jn_{a,j}}\cdot$     (50)

Clearly, electron electrical charge is ${z}_{{\rm e}}=-1$ everywhere.

3.7 Critical points

Critical points are points where derivatives of variables cannot be determined directly from differential equations. The derivation of critical point conditions for our set of equations is simpler than that of KKI because we use different temperatures of each component and correct momentum equation for electrons. Due to these generalizations, the set of critical point equations is not as coupled as it was in KKI and thus, the obtained critical point conditions are simpler.

We write model equations in a simplified form, where we explicitly write only terms containing derivatives of individual variables and other terms are included into the terms $F_{\rm i}$. Thus, the continuity Eqs. (35a), (35b) are



$\displaystyle \rho_a\frac{{\rm d}{ {v_{\rm r}}}_a}{{\rm d} r}+{ {v_{\rm r}}}_a\frac{{\rm d}\rho_a}{{\rm d} r}= F_{a,1}(r,\rho_a,{ {v_{\rm r}}}_a). $     (51a)



In the electron continuity Eq. (35b) we neglected the derivatives of ionic charge because their contribution to electron continuity equation is only marginal. However, inclusion of such term influences critical point and regularity conditions for electrons only, which will not be used (see bellow).

Similarly we can rewrite momentum Eq. (36). In the momentum equations of absorbing ions we shall linearize a term containing the velocity gradient. Note that because model equations are not quasi-linear (i.e. linear with respect to the derivatives of the independent variables), the mathematically more correct method would employ some form of transformation to the quasi-linear form (cf. Courant & Hilbert 1962). However, because the results are essentially the same in this case, we present analysis of critical points in a simplified form. Thus, momentum equations are



$\displaystyle { {v_{\rm r}}}_a\frac{{\rm d} { {v_{\rm r}}}_a}{{\rm d}r}- \frac{... ...}{T_a}\frac{{\rm d}T_a}{{\rm d}r}= F_{a,2}(r,E,\rho_b, {v_{\rm r}}_b,T_b,z_b). $     (51b)



Similarly, due to the dependence of the Doppler term on the velocity gradient (in the Sobolev approximation) we shall write energy equations in the form of



$\displaystyle \frac{3}{2}{ {v_{\rm r}}}_a\rho_a\frac{{a}_a^2}{T_a}\frac{{\rm d}... ... { {v_{\rm r}}}_a}{{\rm d}r}= F_{a,3}(r,\rho_b, {v_{\rm r}}_b,T_b,z_{\rm b}). $     (51c)



The system of equations is closed by the equation for charge separation electric field, which has a simple form,



$\displaystyle \frac{{\rm d}E}{{\rm d}r}=F_{4}(r,\rho_b,z_b). $(51d)



The system of Eq. (51) can be simplified by inserting the derivatives of density from the Eq. (51a) and derivatives of temperature (51c) into the momentum Eq. (51b). We obtain modified linearized momentum equations

 
$\displaystyle \left[{ {v_{\rm r}}}_a- \frac{\partial {g}_{a}^{{\rm rad}}}{\part... ...{{\rm d} { {v_{\rm r}}}_a}{{\rm d}r}= F_{a}(r,E,\rho_b, {v_{\rm r}}_b,T_b,z_b),$     (52)

where Fa is a combination of Fa,i from Eqs. (51a), (51b), and (51c). Because the charge separation electric field Eq. (51d) does not introduce any physically interesting critical point, the system of Eq. (52) consists of only three independent critical point equations. Each of them will be discussed separately in the following subsections.

3.7.1 Critical point of passive ions

For the passive plasma the critical point condition (52) has a simple form

 
$\displaystyle {{v_{\rm r}}^2}_{{\rm p}}=\frac{5}{3}{a}_{{\rm p}}^2.$     (53)

In order to obtain continuous solution of model equations, at this point should the quantity $F_{a}(r,E,\rho_b, {v_{\rm r}}_b,T_b,z_b)$ vanish (regularity condition). This requirement is fulfilled if
 
    $\displaystyle \frac{10{a}_{{\rm p}}^2}{3r}-g+\frac{{q}_{{\rm p}}}{{m}_{{\rm p}}... ...b}-{ {v_{\rm r}}}_{\rm p}}{\vert{ {v_{\rm r}}}_{b}-{ {v_{\rm r}}}_{\rm p}\vert}$  
    $\displaystyle ~~~~~-\frac{2}{3}\frac{1}{{\rho}_{{\rm p}}{{v_{\rm r}}}_{{\rm p}}... ...m_{{\rm p}}+m_b} \frac{\exp\left(-x_{{{\rm p}}b}^2\right)}{\alpha_{{{\rm p}}b}}$  
    $\displaystyle ~~~~~- \frac{2}{3}\frac{1}{{\rho}_{{\rm p}}{{v_{\rm r}}}_{{\rm p}... ...}b}G(x_{{{\rm p}}b})\vert{ {v_{\rm r}}}_b-{ {v_{\rm r}}}_{{\rm p}}\vert=0. \,\,$ (54)

This equation is a generalization of well-known regularity condition for the coronal wind.

3.7.2 Critical point of absorbing ions

Critical point condition Eq. (52) for absorbing ions has form

$\displaystyle { {v_{\rm r}}}_{{\rm i}}- \frac{\partial {g}_{{\rm i}}^{{\rm rad}... ...^{{\rm rad}}}{\partial\left({\rm d}{{v_{\rm r}}}_{{\rm i}} /{\rm d}r\right)}=0.$     (55)

However, as was discussed by KKI, this critical point condition is not met anywhere in the wind. Therefore, to obtain CAK type solution of the stellar wind we use another condition, which is connected to generalized Abbott waves (see KKI). Such a condition can be derived by multiplying the critical point conditions Eq. (52) by density of corresponding component and summing them
$\displaystyle \sum_a \rho_a \frac{{\rm d} { {v_{\rm r}}}_a}{{\rm d}r} \left[{ {... ... i}} /{\rm d}r\right)}\right]= \sum_a F_{a}(r,E,\rho_b, {v_{\rm r}}_b,T_b,z_b).$     (56)

To assure that common equation of motion (obtained by summing of individual momentum Eq. (36)) does not depend on the derivatives of any variables, a condition
 
$\displaystyle \sum_a \rho_a \frac{{\rm d} { {v_{\rm r}}}_a}{{\rm d}r} \left[{ {... ... rad}}}{\partial\left({\rm d}{{v_{\rm r}}}_{{\rm i}} /{\rm d}r\right)}\right]=0$     (57)

shall be fulfilled. This condition we use to fix the mass-loss rate.

3.7.3 Critical point of electrons

The last critical point condition for electrons has again a very simple form

 
$\displaystyle { {v_{\rm r}}}_{{\rm e}}^2=\frac{5}{3}{a}_{{\rm e}}^2.$     (58)

This critical point condition has a similar form to the condition of nonabsorbing ions Eq. (53). Thus, regularity condition for electrons resembles the regularity condition for passive ions Eq. (54),
 
$\displaystyle \frac{10{a}_{{\rm e}}^2}{3r}$ + $\displaystyle {g}_{{\rm e}}^{{\rm rad}}-g+\frac{{q}_{{\rm e}}}{{m}_{{\rm e}}}E+... ...}}_b-{ {v_{\rm r}}}_{\rm e}}{\vert{ {v_{\rm r}}}_b-{ {v_{\rm r}}}_{\rm e}\vert}$  
  - $\displaystyle \frac{2}{3}\frac{1}{{\rho}_{{\rm e}}{{v_{\rm r}}}_{{\rm e}} } \fr... ...m_{{\rm e}}+m_b} \frac{\exp\left(-x_{{{\rm e}}b}^2\right)}{\alpha_{{{\rm e}}b}}$  
  - $\displaystyle \frac{2}{3}\frac{1}{{\rho}_{{\rm e}}{{v_{\rm r}}}_{{\rm e}} } \su... ..._{{{\rm e}}b}G(x_{{{\rm e}}b})\vert{ {v_{\rm r}}}_b-{ {v_{\rm r}}}_{\rm e}\vert$  
  - $\displaystyle \frac{2}{3}\frac{1}{{\rho}_{{\rm e}}{{v_{\rm r}}}_{{\rm e}} }{Q}_{{\rm e}}^{{\rm rad}}=0.$ (59)

However, numerical tests showed that electron regularity condition (59) is approximately fulfilled at the electron critical point if the zero current condition is used as a boundary condition. Thus, this condition should not necessarily be included into the set of model equations.

3.8 Boundary conditions

3.8.1 Boundary conditions for temperatures

We assume that the flow at the inner boundary is in radiative equilibrium and that the boundary temperature of all components is the same, thus, we write boundary condition for temperatures in the form of


$\displaystyle T_{a}={T}_{{\rm e}},\quad a={\rm i, p},$     (60a)


$\displaystyle Q^{{\rm rad}}(R_{*})=0. $     (60b)


Boundary values of ionic charges can be directly obtained from the condition of ionization equilibrium (49).

3.8.2 Boundary condition for velocity

Conditions (54), (57) can be generally used to fix the boundary values of model quantities. However, inclusion of two inner conditions directly into model equations sometimes leads to numerical problems. Therefore, we use a more secure method, which gives essentially the same results.

We start to calculate our models at the passive plasma critical point. Consequently, the boundary condition for the passive plasma velocity is the critical point condition Eq. (53). Boundary condition for the velocity of absorbing ions may be obtained from the passive plasma regularity condition Eq. (54). Because we suppose equal boundary temperatures of each component Eq. (60a), the regularity condition may be simplified

 
$\displaystyle \frac{10{a}_{{\rm p}}^2}{3r}-g+\frac{{q}_{{\rm p}}}{{m}_{{\rm p}}... ...b}-{ {v_{\rm r}}}_{\rm p}}{\vert{ {v_{\rm r}}}_{b}-{ {v_{\rm r}}}_{\rm p}\vert}$$\displaystyle ~~~~- \frac{2}{3}\frac{1}{{\rho}_{{\rm p}}{{v_{\rm r}}}_{{\rm p}}... ...{\rm p}}b}G(x_{{{\rm p}}b})\vert{ {v_{\rm r}}}_b-{ {v_{\rm r}}}_{\rm p}\vert=0.$ (61)

The boundary value of electron velocity is chosen to fulfil the electron regularity condition Eq. (59) at the electron critical point Eq. (58). As was already mentioned, this condition is approximately satisfied if the zero current condition

 
$\displaystyle {n}_{{\rm e}}{ {v_{\rm r}}}_{{\rm e}}={z}_{{\rm p}}{n}_{{\rm p}}{ {v_{\rm r}}}_{{\rm p}} + {z}_{{\rm i}} {n}_{{\rm i}}{ {v_{\rm r}}}_{{\rm i}}$     (62)

is used as boundary condition for electron velocity. The latter condition (62) was applied in our models.

3.8.3 Boundary conditions for density

We write the boundary condition for the passive plasma density in the same form as in KKI,

 
$\displaystyle z_*{\rho}_{{\rm p}}(R_{*})= \frac{1}{{\mbox{\eufont Y}}_{{\rm i}}... ...}) \frac{{{v_{\rm r}}}_{{\rm i}} (R_{*})}{{{v_{\rm r}}}_{{\rm p}} (R_{*})}\cdot$     (63)

Here we only newly introduced the relative abundance z* which accounts for different chemical composition in stellar atmospheres.

The boundary value of ionic density is determined numerically to obtain CAK-type solution (see Sect. 3.9). Boundary electron density is calculated from the condition of quasi-neutrality

 
$\displaystyle {n}_{{\rm e}}=\sum_{a\neq{\rm e}}z_a n_a.$     (64)

3.8.4 Boundary conditions for electric field

Because we have no critical point condition to determine the intensity of the electric field at the stellar surface, we used the condition of neutrality, which simply sets the gradient of the electric field at the stellar surface to zero (cf. Eq. (48)).

   
3.9 Numerical method

We apply the Henyey method (Henyey et al. 1964), which is a modification of the well-known Newton-Raphson method to solve equations described here together with the appropriate boundary conditions. We use essentially the same method as KKI, except that the vector of variables at each grid point d has the form of

$\displaystyle \vec{\psi}_d =\left(\rho_{a,d},v_{{\rm r}a,d},T_{a,d},z_{a,d}, E_{\rm d},\Delta v_{{\rm r},d}\right), \quad a={\rm e},{\rm i},{\rm p}$     (65)

where the velocity difference
$\displaystyle \Delta {v_{\rm r}}={{v_{\rm r}}}_{{\rm i}} -{{v_{\rm r}}}_{{\rm p}}$     (66)

may be added to the set of variables to assure better convergence of the models.

First of all we search for the boundary density $\rho_{0}={\rho}_{{\rm i}}(R_{*})$. We compute several wind models (each of them is a result of several Newton-Raphson iterative steps) for the region near the star for different values of $\rho_{0}$ (for more details see KK0, KKI). We select such value of $\rho_{0}$ which allows wind model to pass smoothly through the point defined by the Eq. (57) and to obtain CAK-type solution. After the appropriate value of $\rho_{0}$ is chosen, we compute the wind model downstream of the point defined by the Eq. (57) again using several Newton-Raphson iterative steps.

The detailed method of calculation of Gayley-Owocki heating/cooling term is given in Appendix A.

4 Results of calculations


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524f1.eps}\\ [4mm] \includegraphics[width=8.4cm,clip]{MS1524f2.eps} \end{figure} Figure 1: Upper panel: comparison of the radial wind velocity of one-component (dotted line) and three-component radiatively driven stellar wind models of an O6 star. Radiatively accelerated ions are denoted using dashed-dotted line, passive plasma with full line and electrons with dashed line. Notice that all curves are very nearly the same. Lower panel: comparison of temperature stratification of one-component and three-component models. Assignment of all curves is the same as for velocity. Notice that the curves describing temperature of individual components of a three-component model coincide, which is not true if we compare one and three component models.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524f3.eps}\\ [4mm] \includegraphics[width=8.4cm,clip]{MS1524f4.eps} \end{figure} Figure 2: The same as Fig. 1 for an O8 star. There is no significant difference in the temperatures of each component. The wind temperature of the one-component and three-component models differ due to the Gayley-Owocki heating.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524f5.eps}\\ [4mm] \includegraphics[width=8.4cm,clip]{MS1524f6.eps} \end{figure} Figure 3: The same as Fig. 1 for a B0 star. For this star Gayley-Owocki heating (in the outer parts of the wind) and cooling (in the inner parts of the wind) effects are important for temperature structure. Note that the frictional heating is negligible in this case and that the temperatures of particular components are nearly the same.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524f7.eps}\\ [4mm] \includegraphics[width=8.4cm,clip]{MS1524f8.eps} \end{figure} Figure 4: The same as Fig. 1 for a B2 star. Gayley-Owocki heating and cooling is greater than for a B0 star - cf. Fig. 3. Frictional heating is negligible. The temperatures of particular components are nearly the same.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524f9.eps}\\ [4mm] \includegraphics[width=8.4cm,clip]{MS1524fa.eps} \end{figure} Figure 5: The same as Fig. 1 for a B3 star. The wind is heated by both frictional and Gayley-Owocki heating in the outer parts of the wind. The temperatures of absorbing ions and electrons are nearly the same whereas the ionic temperature slightly differs mainly in the outer parts of the wind.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524fb.eps}\\ [4mm] \includegraphics[width=8.4cm,clip]{MS1524fc.eps} \end{figure} Figure 6: The same as Fig. 1 for a B4 star. The wind is heated by both frictional and Gayley-Owocki heating in the outer parts of the wind. Notice that the effect of heating is more pronounced than for the case of a B3 star (Fig. 5) and the temperature at its maximum is larger than for a cooler star. The temperatures of absorbing ions and electrons are nearly the same.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524fd.eps} \end{figure} Figure 7: The dependence of the ionic charge on radius for different wind models. Accelerated ions are denoted using dashed-dotted line.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.5cm,clip]{MS1524fe.eps}\\ [4mm] \includegraphics[width=8.5cm,clip]{MS1524fg.eps} \end{figure} Figure 8: The same as Fig. 1 for a B5 star. The passive component decouples just above the stellar surface and subsequently falls back onto the star. Note that the temperature of the ionic component reaches the value of the order of $10^6~{\rm K}$. The temperatures of absorbing ions and electrons are nearly the same.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524fh.eps}\\ [4mm] \includegraphics[width=8.4cm,clip]{MS1524fi.eps} \end{figure} Figure 9: The same as Fig. 1 for a main sequence star $\tau $ Sco. The wind is heated both by frictional and Gayley-Owocki heating.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524fj.eps}\\ [4mm] \includegraphics[width=8.4cm,clip]{MS1524fk.eps} \end{figure} Figure 10: The same as Fig. 1 for a giant star $\epsilon $ CMa. The wind is heated both by frictional and Gayley-Owocki heating. Note that the wind temperature at its maximum nearly reaches the effective temperature of the star.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.4cm,clip]{MS1524fl.eps}\\ [4mm] \includegraphics[width=8.4cm,clip]{MS1524fm.eps} \end{figure} Figure 11: The same as Fig. 1 for giant star $\beta $ CMa. The wind is heated both by frictional and Gayley-Owocki heating.
Open with DEXTER

We computed several wind models for different stellar spectral types. Parameters of individual model stars are listed in Table 1. Main sequence stellar parameters are taken from Harmanec (1988). For $\tau $ Sco, $\epsilon $ CMa, and $\beta $ CMa we used the same parameters as SP and Cassinelli et al. (1995, 1996), respectively. Note, that here we used slightly different parameters of $\tau $ Sco than KKI. Force multipliers were adopted from Abbott (1982). The parameters of absorbing ions were selected in the following way. Because main sequence models were computed mainly for a demonstration of particular effects, it was sufficient to choose an ion which is simple enough and which describes the basic line driving. So, we selected a carbon atom with ${m}_{{\rm i}}=12{m}_{{\rm p}}$ as a driving ion for them. On the other hand, a driving ion of individual giant models was selected more carefully with respect to the stellar type. For $\epsilon $ and $\beta $ CMa we selected iron as a driving ion. The effective temperature of $\tau $ Sco is higher, thus, we again selected carbon as a driving ion for this star. We stress that the selection of driving ions does not influence the amount of radiative force. However, it affects the thermal balance of the wind (via the frictional and Gayley-Owocki heating). For comparison purposes we also computed nonisothermal one-component models (see KKI) of these stars' winds with the same stellar and wind parameters (however, without GO heating).


 

 
Table 1: Adopted parameters of model stars. $\mbox{\eufont M}$ is the stellar mass in units of a solar mass, R* is the stellar radius in units of solar radius, $T_{\rm eff}$ is the star's effective temperature, k, $\alpha $, and $\delta $ are radiative force multipliers, and z* is the metallicity.
Stellar Stellar parameters Wind parameters  
type $\mbox{\eufont M}$ R* $T_{\rm eff}$ z* k $\alpha $ $\delta $ $\frac{{m}_{{\rm i}}}{{m}_{{\rm p}}}$
(star) $[{\mbox{\eufont M}}_{\odot}]$ $[R_{\odot}]$ $[{\rm K}] $          
O6 31.65 9.85 41700 1.0 0.174 0.606 0.120 12.0
O8 21.66 7.51 35600 1.0 0.166 0.607 0.120 12.0
B0 14.57 5.80 29900 1.0 0.156 0.609 0.120 12.0
B2 8.62 4.28 23100 1.0 0.377 0.537 0.091 12.0
B3 6.07 3.56 19100 1.0 0.477 0.506 0.089 12.0
B4 5.12 3.26 17200 1.0 0.365 0.509 0.105 12.0
B5 4.36 3.01 15500 1.0 0.235 0.511 0.12 12.0
$\tau $ Sco 19.60 5.50 33000 0.3 0.113 0.604 0.095 12.0
$\epsilon $ CMa 15.2 16.2 21000 0.18 0.135 0.561 0.092 55.8
$\beta $ CMa 15.5 11.6 23250 0.39 0.125 0.564 0.099 55.8


4.1 Winds with Gayley-Owocki heating

As was shown, e.g., by SP, high density winds are well coupled. For such winds the effect of frictional heating is negligible. Similarly, due to a high wind density the heat exchange between individual components is capable to maintain the same temperature for all components. However, for spectral types cooler than O6 a subtle effect of GO heating/cooling influences the temperature structure. This behaviour is displayed in Figs. 1 and 2 for wind models of O6 and O8 stars. There is a large velocity gradient near the star, the variable $\sigma $(Eq. (27)) is positive and thus, the wind is very slightly cooled by the GO cooling compared to the model without GO heating/cooling effects. On the other hand, in the outer parts of the wind the velocity gradient is lower, the variable $\sigma $ is negative and thus GO heating dominates (cf. Eqs. (32), (33)). At the outermost parts of the wind the temperatures of models with and without GO heating are again nearly the same mainly due to the lowering of a stellar angular diameter.

As was discussed by GO, for stars with a lower density wind GO heating and GO cooling are much more evident. This can be seen in Figs. 3, 4 for the case of B0 and B2 stars. As was shown by KKI, for such stars the frictional heating is negligible and thus, changes in temperature stratification are caused only by GO heating/cooling. Due to the changed temperature stratification another effects becomes important. Because the radiative force in the CAK parametrization (Eq. (38)) depends explicitly on the thermal velocity, higher wind temperature causes its lowering. Lower radiative force in the outer parts of the wind (above the critical point) leads to lowering the outflow velocity (cf. Vink et al. 1999). However, the description of the dependence of the radiative force on the temperature via the Eq. (38) is only approximate. Further calculations needed for better quantitative understanding the dependence of the radiative force on the temperature are currently under way and will be reported in future paper(s).

4.2 Different temperatures of individual components

For stars with lower wind density, individual wind components have different temperatures. This is shown in Figs. 5 and 6 for B3 and B4 stars.

Near the stellar surface the wind is relatively dense, the heat exchange between individual components is effective and temperatures of individual components are nearly equal. Compared to the model without GO cooling the wind temperature is slightly lower.

However, this is not the case in the outer parts of the wind. The wind is heated by frictional and GO heating there. Because the wind is more tenuous, the heat exchange between individual components is not so effective as it is near the star and the temperatures of individual components differ. The temperature of absorbing ions is the highest and electron and passive component temperatures are nearly equal.

There are three possible mechanisms which heat (or cool) individual components of the flow selectively and thus, which could make temperatures of particular components different. First, radiative heating/cooling (in our case made by bound-free and free-free transitions) deposits (or picks up) energy from the electron thermal pool. This causes electrons to incline to temperatures given by former models for whole wind - see Drew (1989) for the one-component case and KKI when frictional heating is included. Second, the GO heating cools absorbing ions below the point where $\sigma=0$ (see GO), and heats them above this point. Finally, frictional heating itself deposits thermal energy unevenly. From the functional behaviour of the last right-hand side term of Eq. (46) we can infer that temperature increase is proportional to $\rho_b/\left(m_a+m_b\right)$. Thus, those frictionally heated are mainly low-density components, i.e. electrons and absorbing ions, both via their collisions with passive ions.

Another effect which influences the temperature balance is the heat exchange between components, described by the second right-hand side term of Eq. (46). Clearly, the heat exchange depends mainly on the product of number densities of components. Thus, similarly to the differences in velocity, temperature differentiation takes place mainly for the low density wind. As discussed above, in such winds absorbing ions can be heated more than other two components. On the other hand, due to their large number densities, electrons and passive ions will share nearly the same temperature. All these effects influence models in Figs. 5 and 6.

We determine the charge separation field directly from the Maxwell Eq. (48). However, this has only marginal effect on the wind models, because our models tend to fulfil quasi-neutrality,

$\displaystyle {n}_{{\rm e}}\approx{n}_{{\rm p}} + {z}_{{\rm i}} {n}_{{\rm i}}.$      

The dependence of ionic charge for selected wind models is given in Fig. 7. Clearly, because the ionization equilibrium depends mainly on radiative temperature, the ionic charge is nearly constant through the wind.

Many B stars exhibit UV-excess (e.g. Cassinelli et al. 1995, 1996; Morales et al. 2001). Note, that frictional heating and Gayley-Owocki heating of the stellar wind could be one of the possible explanations (Babel 1995).

Frictional heating in one-component models:
The frictional heating term in the energy equation of a one-component wind can be roughly estimated without using multicomponent models. In the ionic equation of motion (36) all terms without radiative force and frictional term (corresponding to encounters between passive component and ions) are small. We shall neglect them in this paragraph to derive an approximate expression for frictional heating. The radiative acceleration term can be expresses as
 
$\displaystyle {g}_{s}^{{\rm rad}}\approx \frac{1}{{\rho}_s} K_{s{\rm p}}G(x_{s{\rm p}}).$     (67)

Here the subscript s denotes absorbing ions (note, that in these equations we allow for more than one type of absorbing ions). The frictional term in the energy equation of one-component wind can be approximated as
 
$\displaystyle Q^{{\rm fric}}$ $\textstyle \approx$ $\displaystyle \displaystyle\sum_{{{\rm absorbing} \,\, {\rm ions}}} K_{s{\rm p}}G(x_{s{\rm p}})\vert{ {v_{\rm r}}}_{\rm p}-{ {v_{\rm r}}}_s\vert$  
  $\textstyle \approx$ $\displaystyle \displaystyle\sum_{{{\rm absorbing} \,\, {\rm ions}}} {\rho}_s{g}_{s}^{{\rm rad}}\left({ {v_{\rm r}}}_s-{ {v_{\rm r}}}_{\rm p}\right),$ (68)

where for the determination of the drift velocity from Eq. (67) an approximation of the Chandrasekhar function
$\displaystyle G(x)\approx \frac{2}{3\sqrt{\pi}}x$     (69)

can be used. We computed a one-component model for a B3 star where the frictional term (68) was inserted into the energy equation of one-component wind model (KKI, Eq. (40c)). Although the temperature was overestimated by $10\%$ (compared to the correct three-component model) in the region where frictional heating is important, the difference of the terminal velocity was only marginal. Note, that the frictional heating approximation Eq. (68) can be used even when there is more than one absorbing ion component.

4.3 Backfalling of hydrogen to the stellar surface

In the case of the wind with the lowest densities the absorbing ions are not able to accelerate sufficiently the passive component of the wind. Thus, the passive component is not dragged out of the atmosphere and falls back onto the stellar surface (see Fig. 8 for a model of a B5 star wind). Such reaccretion should be studied using hydrodynamical calculations (Porter & Skouza 1999). Moreover, Babel (1996) showed that the hydrostatic solution for passive plasma and the wind solution for absorbing ions exists. Probably, this type of solution is common for low-density winds, because Dworetsky & Budaj (2000) whilst studying Ne abundances in peculiar HgMn stars showed that in these stars the radiatively driven stellar wind with hydrogen mass loss rate larger than $10^{-14}\,{\mbox{\eufont M}}_{\odot}/{\rm yr}$ is not present.

Decoupling of velocities of absorbing and passive components is accompanied by decoupling of temperatures of these components (see Fig. 8). This effect is caused by the dependence of the amount of heat transferred between individual components on the velocity difference (see Eq. (46)). Note that the absorbing component attains temperatures sufficient to produce X-rays. This effect can help to explain enhanced X-ray emission observed in mid- and late-B stars (Cohen et al. 1997) which cannot be regarded as a consequence of the standard radiation driven wind-shock mechanism. Another model for X-ray emission based on the shock decoupling was given by Porter & Drew (1995).

4.4 Wind models of particular stars

As was shown by KKI, the heating effect is pronounced in the wind of $\tau $ Sco. Thus, we decided to recompute the wind model with the inclusion of GO heating. This model is shown in the Fig. 9. Similar effects as in wind models of main sequence stars occur in wind models of giants. This can be seen in Figs. 10 and 11 for the wind models of $\beta $ CMa and $\epsilon $ CMa, respectively. For all these stars both frictional and GO heating are important for the temperature structure of the outer parts of the wind.

For the star $\tau $ Sco we used lower than observed value of metallicity. Contrary to Kilian (1994) who determined the value z*=0.6 we reduced the metallicity to z*=0.3to enable larger frictional heating. This change reflects mainly the uncertainties of our model, because, e.g., our model with z*=0.5 and the metallic component described by iron ions instead of carbon ions yields nearly the same velocity and temperature stratification.

Unfortunately, existing measurements of the terminal velocity for this star do not allow us to verify our models precisely. Abbott (1978) and Lamers & Rogerson (1978) determined $v_{\infty}=2000~{\rm km}\,{\rm s}^{-1}$, whereas Lamers et al. (1995) measured $v_{\infty}=1000~{\rm km}\,{\rm s}^{-1}$. However, all of them claim that their values are uncertain. Larger values of $v_{\infty }$ are supported also by a detailed UV-fit of Hamann (1981).

For $\epsilon $ CMa we used metallicity z*=0.18, a value estimated by Gies & Lambert (1992). Similarly to $\tau $ Sco, available determinations of terminal velocity have lower quality. Abbott (1978) determined $v_{\infty}=700~{\rm km}\,{\rm s}^{-1}$and Lamers et al. (1995) measured $v_{\infty}=600~{\rm km}\,{\rm s}^{-1}$. However, it is not clear whether the apparent discrepancy of theoretical and observational terminal velocities is caused by the models or is due to the inaccurate measurements.

According to Gies & Lambert (1992) we reduced the metallicity of $\beta $ CMa to the value z*=0.39. To our knowledge, there is no measured terminal velocity for this star available in the literature. Note that for both $\epsilon $ CMa and $\beta $ CMa enhanced wind temperature can help to explain observed UV-excess (Cassinelli et al. 1995, 1996).

5 Comparison of terminal velocities


 

 
Table 2: Stellar and wind parameters of O6-B5 stars selected from LSL. $\dot{\mbox{\eufont M}}$ is computed mass-loss rate. Terminal velocities measured by LSL (Col. $v_{\infty }$ (LSL)) are compared with theoretical values obtained by LSL using a "cooking formula'' of KPPA (Col. $v_{\infty }$ (KPPA)) and with predicted ones computed with an assumption of a nonisothermal wind model (Col. $v_{\infty }$ (predicted)).
HD Stellar parameters Wind parameters   Terminal velocities
number $\mbox{\eufont M}$ R* $T_{\rm eff}$ z* k $\alpha $ $\delta $ $\dot{\mbox{\eufont M}}$ $v_{\infty }$ (LSL) $v_{\infty }$ (KPPA) $v_{\infty }$ (predicted)
  $[{\mbox{\eufont M}}_{\odot}]$ $[R_{\odot}]$ $[{\rm K}] $         $[{\mbox{\eufont M}}_{\odot}\,{\rm yr}^{-1}]$ $[{\rm km}\,{\rm s}^{-1}]$ $[{\rm km}\,{\rm s}^{-1}]$ $[{\rm km}\,{\rm s}^{-1}]$
30614 43.0 27.6 $30\,900$ 1.0 0.158 0.609 0.120 ${8.9\times 10^{-6}}$ $ 1500\,\pm\,200 $ 2241 1450
34656 30.0 9.9 $ 38\,100 $ 1.0 0.171 0.607 0.120 ${1.0\times 10^{-6}}$ $ 2100\,\pm\,100 $ 3598 2590
36861 30.0 12.3 $36\,000$ 1.0 0.167 0.607 0.120 ${1.5\times 10^{-6}}$ $ 2200\,\pm\,300 $ 3084 2230
41117 25.0 43.4 $18\,500$ 1.0 0.410 0.507 0.098 ${2.8\times 10^{-6}}$ $ 500\,\pm\,50 $ 1058 740
43384 19.0 39.8 $16\,300$ 1.0 0.311 0.510 0.112 ${4.0\times 10^{-7}}$ $ 500\,\pm\,100 $ 1015 690
47240 17.0 23.4 $20\,800$ 1.0 0.451 0.514 0.091 ${9.7\times 10^{-7}}$ $ 1000\,\pm\,100 $ 1267 930
51309 11.0 16.3 $16\,700$ 1.0 0.329 0.509 0.109 ${1.7\times 10^{-8}}$ $ 700\,\pm\,100 $ 1309 910
52382 17.0 20.4 $20\,800$ 1.0 0.451 0.514 0.091 ${4.9\times 10^{-7}}$ $ 1200\,\pm\,100 $ 1381 1030
69464 49.0 20.1 $37\,200$ 1.0 0.169 0.607 0.120 ${9.9\times 10^{-6}}$ $ 2100\,\pm\,200 $ 2721 1840
74194 28.0 14.5 $33\,000$ 1.0 0.161 0.608 0.120 ${1.4\times 10^{-6}}$ $ 2000\,\pm\,300 $ 2829 1970
79186 18.0 62.4 $13\,600$ 1.0 0.284 0.519 0.100 ${6.6\times 10^{-7}}$ $ 450\,\pm\,50 $ 772 530
91572 38.0 9.6 $42\,200$ 1.0 0.175 0.606 0.114 ${1.6\times 10^{-6}}$ $ 2400\,\pm\,100 $ 3989 2940
91969 25.0 22.9 $26\,000$ 1.0 0.284 0.568 0.108 ${3.0\times 10^{-6}}$ $ 1500\,\pm\,100 $ 1816 1240
92964 29.0 68.4 $17\,400$ 1.0 0.361 0.509 0.105 ${1.1\times 10^{-5}}$ $ 550\,\pm\,50 $ 807 530
93130 43.0 13.8 $40\,200$ 1.0 0.174 0.606 0.119 ${4.5\times 10^{-6}}$ $ 2500 \,\pm\, 300 $ 3337 2370
96248 25.0 38.9 $20\,800$ 1.0 0.451 0.514 0.091 ${7.3\times 10^{-6}}$ $ 650\,\pm\,50 $ 1116 770
96917 46.0 25.2 $33\,000$ 1.0 0.161 0.608 0.120 ${9.5\times 10^{-6}}$ $ 1800\,\pm\,200 $ 2430 1600
101190 48.0 13.9 $42\,200$ 1.0 0.175 0.606 0.114 ${6.3\times 10^{-6}}$ $ 2800\,\pm\,200 $ 3452 2450
101436 42.0 12.4 $41\,200$ 1.0 0.174 0.606 0.117 ${3.5\times 10^{-6}}$ $ 2700\,\pm\,200 $ 3487 2570
106343 24.0 40.7 $19\,700$ 1.0 0.464 0.506 0.091 ${5.3\times 10^{-6}}$ $ 800\,\pm\,100 $ 1074 730
109867 26.0 38.9 $20\,800$ 1.0 0.451 0.514 0.091 ${6.7\times 10^{-6}}$ $ 1200\,\pm\,200 $ 1144 800
112244 46.0 25.2 $33\,000$ 1.0 0.161 0.608 0.120 ${9.5\times 10^{-6}}$ $ 1600\,\pm\,100 $ 2430 1580
116084 15.0 24.8 $17\,400$ 1.0 0.361 0.509 0.105 ${1.4\times 10^{-7}}$ $ 500\,\pm\,100 $ 1201 830
148379 24.0 40.7 $19\,700$ 1.0 0.464 0.506 0.091 ${5.3\times 10^{-6}}$ $ 500\,\pm\,100 $ 1074 730
151515 41.0 14.9 $ 38\,100 $ 1.0 0.171 0.607 0.120 ${4.1\times 10^{-6}}$ $ 2400\,\pm\,100 $ 3196 2230
151804 70.0 34.0 $34\,000$ 1.0 0.163 0.608 0.120 ${3.6\times 10^{-5}}$ $ 1500\,\pm\,200 $ 2226 1370
152405 25.0 15.3 $30\,500$ 1.0 0.157 0.609 0.120 ${1.0\times 10^{-6}}$ $ 1800\,\pm\,200 $ 2616 1800
152424 52.0 33.4 $30\,500$ 1.0 0.157 0.609 0.120 ${1.5\times 10^{-5}}$ $ 1500\,\pm\,100 $ 2056 1350
154090 26.0 38.9 $20\,800$ 1.0 0.451 0.514 0.091 ${6.7\times 10^{-6}}$ $ 950\,\pm\,50 $ 1144 800
157246 17.0 23.4 $20\,800$ 1.0 0.451 0.514 0.091 ${9.7\times 10^{-7}}$ $ 900\,\pm\,200 $ 1267 930
162978 40.0 16.0 $37\,100$ 1.0 0.169 0.607 0.120 ${4.4\times 10^{-6}}$ $ 2200\,\pm\,200 $ 2939 2090
163758 50.0 20.1 $37\,200$ 1.0 0.169 0.607 0.120 ${9.5\times 10^{-6}}$ $ 2300\,\pm\,200 $ 2765 1890
166596 9.7 9.8 $18\,700$ 1.0 0.419 0.507 0.097 ${9.6\times 10^{-9}}$ $ 700\,\pm\,50 $ 1551 1170
175754 34.0 14.2 $36\,000$ 1.0 0.167 0.607 0.120 ${2.4\times 10^{-6}}$ $ 2000\,\pm\,100 $ 2998 2130
186980 35.0 13.9 $37\,100$ 1.0 0.169 0.607 0.120 ${2.8\times 10^{-6}}$ $ 2100\,\pm\,100 $ 3028 2180
188209 43.0 27.6 $30\,900$ 1.0 0.158 0.609 0.120 ${8.9\times 10^{-6}}$ $ 1600\,\pm\,100 $ 2241 1450
190603 24.0 40.7 $19\,700$ 1.0 0.464 0.506 0.091 ${5.3\times 10^{-6}}$ $ 500\,\pm\,50 $ 1074 730
190864 42.0 14.0 $39\,200$ 1.0 0.173 0.606 0.120 ${3.8\times 10^{-6}}$ $ 2300\,\pm\,200 $ 3338 2340
198478 17.0 36.3 $16\,300$ 1.0 0.311 0.510 0.112 ${3.0\times 10^{-7}}$ $ 550\,\pm\,100 $ 1009 690
204172 23.0 20.0 $26\,000$ 1.0 0.284 0.568 0.108 ${1.8\times 10^{-6}}$ $ 1600\,\pm\,200 $ 1911 1320
206165 19.0 35.8 $18\,500$ 1.0 0.410 0.507 0.098 ${1.6\times 10^{-6}}$ $ 600\,\pm\,50 $ 1011 720
210809 38.0 21.4 $32\,000$ 1.0 0.160 0.608 0.120 ${4.2\times 10^{-6}}$ $ 2000\,\pm\,200 $ 2588 1700
210839 51.0 19.6 $38\,200$ 1.0 0.171 0.607 0.120 ${1.1\times 10^{-5}}$ $ 2200\,\pm\,200 $ 2882 1930
213087 21.0 23.4 $23\,400$ 1.0 0.368 0.541 0.100 ${1.9\times 10^{-6}}$ $ 1400\,\pm\,200 $ 1533 1080
218915 43.0 27.6 $30\,900$ 1.0 0.158 0.609 0.120 ${8.9\times 10^{-6}}$ $ 1800\,\pm\,100 $ 2241 1450



  \begin{figure} \par\includegraphics[width=8.53cm,clip]{MS1524fn.eps} \end{figure} Figure 12: Comparison of predicted and observational (taken from LSL) values of terminal velocities. Vertical lines denote uncertainty of observed values. Straight line is one-to-relation. For comparison, we plotted theoretical terminal velocities computed by LSL using "cooking formula'' (crosses).
Open with DEXTER

In addition, we decided to compare our predicted terminal velocities with that measured by Lamers et al. (1995, hereafter LSL). They found a discrepancy between theoretical values obtained from a "cooking formula'' of Kudritzki et al. (1989, hereafter KPPA) and their experimental values. We computed wind models of O6-B5 stars for which LSL measured the terminal velocity. Parameters of each wind model are given in Table 2. Stellar parameters are taken from LSL, wind parameters are adopted from Abbott (1982). For many stars we found quite a good agreement between observed and predicted terminal velocities (see Fig. 12 for comparison of predicted and observed values). Although some values of predicted terminal velocities miss the measured value significantly (e.g. for the star HD166596), it is evident that the overall agreement between our predicted terminal velocities and the observed ones is much better than that of the "cooking formula'' of KPPA and the systematic difference, which was previously attributed to an overestimation of $\alpha $(by LSL), has been removed.

However, there are still differences between observed and predicted values of $v_{\infty }$. There are three basic reasons for this discrepancy. First, rotation lowers the terminal velocity (cf. Friend & Abbott 1986). However, Petrenz & Puls (2000) using 2D models showed that the influence of the rotation on terminal velocity in many cases is only marginal. Second, our wind models (especially the radiative force) are constrained. Although we included physical processes that have not been included yet (frictional heating, Gayley-Owocki heating, multicomponent nature of the wind), there are still limits. Our treatment of ionization is only approximate, the equilibrium is not determined consistently with radiation field. In addition, our models are not fully consistent with respect to the radiative force, a proper NLTE treatment of the radiative transfer problem would be very useful. This two reasons causes that many of the terminal velocities are not within quoted uncertainties. However, we plan to improve our models in near future.

Another source of differences may come from uncertainties of stellar parameters derived from observations. Note that, e.g., stars HD106343, HD148379, and HD190603 have roughly the same parameters, however different observed terminal velocities.

The "cooking formula'' of KPPA should be consistent with detailed calculations of Pauldrach et al. (1986) with an accuracy about 5%. However, our predicted terminal velocities correspond to those computed by Pauldrach et al. (1986), too. Thus, there is not clear source of discrepancy between terminal velocities observed by LSL and predicted using formula of KPPA. We stress that the effect of frictional or Gayley-Owocki heating on the terminal velocity are negligible for the models described in this section.

6 Conclusions

We computed non-isothermal three-component models of OB star winds with allowing for different temperatures of each component and with inclusion of the Gayley-Owocki (GO) heating/cooling. We showed that temperature differentiation takes place in the winds of B stars starting from spectral type B3. The temperature of absorbing ions is of the order $10^3\,{\rm K}$higher than temperature of other components whereas the temperatures of passive plasma and electrons is nearly equal. The main sources which trigger the temperature differentiation are GO, frictional, and radiative heating.

Another important effect studied in this paper is the GO heating and cooling, which is important mainly for the low density winds. We showed that this effect is a direct consequence of the dependence of the radiative force on the wind velocity. We derived the GO heating formula directly from the Boltzmann equation. More subtle GO cooling operates near the star at the wind base whereas the GO heating affects the flow mainly in outer parts of the wind. These effects become important starting from stellar type O6. Frictional and GO heating provides a possibility for an alternative explanation of UV-excess observed in some B stars.

At the lowest densities either the passive component falls back onto the star or purely metallic wind exists. If the reaccretion takes place then ionic components is frictionally heated to the temperatures of orders millions K creating corona-like region. This effect can explain enhanced X-ray activity in many of B stars.

Finally, we compared our computed terminal velocities with those derived from observation. There is quite good agreement between them. The systematic difference between observed and predicted (by a "cooking formula'' of KPPA) terminal velocities found by LSL was removed. However, we found no effects of frictional or GO heating in our sample.

Acknowledgements

The authors would like to thank Dr. John Porter for pointing their attention to the importance of the effect of Gayley-Owocki (Doppler) heating, Dr. Kenneth Gayley and Prof. Michal Lenc for their comments on the manuscript of this paper. This research has made use of NASA's Astrophysics Data System Abstract Service (Kurtz et al. 2000; Eichhorn et al. 2000; Accomazzi et al. 2000; Grant et al. 2000). This work was supported by a grant GA CR 205/01/0656 and by projects K2043105 and Z1003909.

   
Appendix A: Numerical calculation of G($\sigma $,$\mu _*$)

The function $G(\sigma,\mu_*)$ is computed using numerical quadrature. Firstly, the integral over x can be efficiently computed using a Hermite quadrature formula (cf. Ralston 1965). Quadrature weights and knots were computed using a subroutine IQPACK which is an implementation of method described by Kautsky & Elhay (1982). Satisfactory approximation can be obtained using 20 quadrature points.

For an angle integration we used Legendre quadrature formula with 5 quadrature points. Again, quadrature weights and knots were computed using subroutine IQPACK (Kautsky & Elhay 1982).

Finally, the integration over y was performed using Simpson quadrature rule. The quadrature integral $(0,\infty)$ was approximated by (10-3,105) and the Simpson integration was divided into subintervals of power 10, with 10 quadrature points in each of them.

Numerical tests showed that temperature computed with described approximation of the GO heating has an error less than 1%.

References

 

Copyright ESO 2001