3 The method for calculating optically thick wind models for WR-stars

We describe the method for solving the equations and for determining the mass-loss rates for optically thick winds of WR-stars. The essence of the analysis is to find under what conditions the optically thick wind models produce the observed mass-loss rates. We will do this for a grid of models of different values of M_*, $L_* = L_{\rm core}$, $R_{\rm s}$ and $R_{\rm hc}$.

The adopted scheme consists of choosing the effective optical depth $\taup_{\rm s}$ at the sonic point. We then make a first estimate of one particular parameter at the sonic point. This parameter is either the radiative luminosity $L_{\rm s}$ (Variant A) or the opacity $\chi_{\rm s}$ (Variant B). We then calculate a self-consistent value of $\dot{M}$ and of v, T, $\chi$, L and their first and second order derivatives at the sonic point, from the above described equations. Knowing these values at the sonic point, we can check if the value of the velocity gradient that is derived from the momentum equation is identical to the one that is needed to satisfy the energy equation. So basically we check that the momentum equation and the energy equation are both satisfied in the transonic region. If this is not the case, we modify our choice of $L_{\rm s}$ or $\chi_{\rm s}$. We developed a numerical method that converges very rapidly and gives the self-consistent models within typically about some tens of iterations.

Up to this point we have only considered the conditions at or near the sonic point. The next step is to find the value of the velocity law exponent $\beta$ that is needed to produce the chosen value of $\taup_{\rm s}$. This is obtained by integrating the optical depth over the whole wind for different values of $\beta$ until the pre-adopted value of $\taup_{\rm s}$ is found. The resulting model then gives, for any particular choice of the stellar parameters M_*, L_*, $R_{\rm s}$, $R_{\rm hc}$ and wind parameter $\mbox{$v_\infty$ }$, and for any pre-chosen value of $\taup_{\rm s}$ (which is basically the same as choosing the gas temperature at the sonic point) the resulting value of $\dot{M}$ as well as the temperature and velocity structure of the supersonic part of the wind. The calculated mass-loss rate can then be compared with the observed values.

In the previous section we described two methods for dealing with the opacity: Variants A and B. We discuss the calculation of the models in these two variants separately. In this section we describe the computations for fixed values of $R_{\rm s}$. These models are named A0 and B0. In Sect. 8 we will describe the models with a variable value of $R_{\rm s}$.

   
3.1 Variant A0

In this variant we assume that L(r) is constant in the wind. The variation of the opacity $\chi(r)$ then follows from the condition that the wind is radiatively driven.

First we adopt the values of $R_{\rm s}$ and $R_{\rm hc}$(the range of these parameters for WR-stars will be justified later in Sect. 6.2) and we adopt a value of the effective optical depth $\taup_{\rm s}$ at this sonic point. We then make a first guess of $L_{\rm s}$, which is only slightly smaller than L_*, that sets the value of $\mbox{$T_{\rm eff}$ }(R_{\rm s})$. This gives the gas temperature at the sonic point $T_{\rm s}$ (Eq. (34)) as well as the sound speed and the flow velocity, $v_{\rm s}$. The mass-loss rate $\mbox{$\dot{M}$ }$ then follows from formula (19).

The initial guess values of other terms at the sonic point are found by the scheme:
(i) $\chi_{\rm s}$ and $({\rm d}T/{\rm d}r)_{\rm s}$ are both found from the combination of Eqs. (16) and (25).
(ii) $({\rm d}v/{\rm d}r)_{\rm s}$ is found from the derivative of the energy equation (10) with the derivative of $\mbox{$L_{\rm adv}$ }$ (Eq. (18)).
(iii) $({\rm d}^2T/{\rm d}r^2)_s$ and $({\rm d} \chi /{\rm d}r)_{\rm s}$ are found from the combination of Eqs. (16) for the temperature gradient and de l'Hopital's rule (26) for the momentum equation.
(iv) $({\rm d}^2v/{\rm d}r^2)_{\rm s}$ is found from the second derivative of the energy Eq. (10).
(v) $({\rm d}^2 \chi/{\rm d}r^2)_{\rm s}$ and $({\rm d}^3T/{\rm d}r^3)_{\rm s}$ are found by combining the regularity condition (Eq. (27)) with the second derivative of the temperature Eq. (16).
(vi) We now have all the required values and gradients at the sonic point and we can test if the momentum and energy equations are both satisfied. To do this we calculate ${\rm d}v/{\rm d}r$ at a point slightly below the sonic point, e.g. at $r=0.99 R_{\rm s}$, from the momentum equation by using Taylor expansions for finding T, $\chi$and v (L is constant). In this calculation we use the values of these physical terms and their first and second derivatives at the sonic point. $({\rm d}^2v/{\rm d}r^2)_{\rm s}$. We then compare this with the value of ${\rm d}v/{\rm d}r$ at that same point, that is required by the energy equation, i.e. by differentiating Eq. (10). If these two values of ${\rm d}v/{\rm d}r$ are not the same, our initial guess of the constant value of L(r) was incorrect. We then choose another value of L(r) and repeat the process untill convergence is reached.
(vii) In the last step we calculate the value of $\beta$ in the supersonic part of the wind that gives the adopted value of $\taup_{\rm s}$ (Eq. (35)). We find $\rho(r)$ from the mass continuity equation with the velocity law of Eq. (38), with w₁ and w₂ derived from $({\rm d}v/{\rm d}r)_{\rm s}$ and $({\rm d}^2v/{\rm d}r^2)_{\rm s}$ and $w_3=v_{\rm s}-w_1-w_2$. The opacity $\chi(r)$ in the supersonic part of the wind, which is needed for the calculation of $\mbox{$\tau^{\prime}$ }(r)$, is derived from the momentum Eq. (36). The temperature structure in the supersonic part is calculated with Eq. (32) with $\mbox{$\tau^{\prime}$ }(r)$ as input. So in practice we calculate both $\mbox{$\tau^{\prime}$ }(r)$ and T(r) by stepwise integration of the Eq. (33), which results in T(r) through Eq. (32), which gives $\chi(r)$ at the next step through Eq. (36).

We want to point out that only in the final step for the determination of $\beta$ we need a numerical integration. All other steps use straightforward calculations of single values at the critical point. These calculations are exact, insofar as the approximations used to derive the formulae are correct. Therefore the calculations are very fast and the iteration converges rapidly.

   
3.2 Variant B0

In this variant we assume that the opacity can be expressed as $\chi(r) = a \rho / T^{\alpha}$, with a to be determined by the conditions at the sonic point, and $\alpha$ is an adopted constant. For instance, Kramers' opacity law has $\alpha=3.5$. We will adopt several values in the range of $3.0 < \alpha < 5.0$.

The process of calculating self-consistent optically thick radiation driven wind models is very similar to that used in Variant A0. We again start with stellar parameters M_*, L_*, $R_{\rm hc}$ and $R_{\rm s}$ and we adopt $\taup_{\rm s}$. As in Variant A0, we make a first guess of $L_{\rm s} < L_*$. However in Variant B0 L(r) is distance dependent, contrary to the situation in Variant A0. The choice of $L_{\rm s}$ sets the values of $T_{\rm s}$, $v_{\rm s}$ and $\dot{M}$, as described in Variant A0. The steps for solving the equations at the sonic point are:
(i) $\chi_{\rm s}$ and $({\rm d}T/{\rm d}r)_{\rm s}$ are both found from the combination of Eqs. (16) and (25). The value of $\chi_{\rm s}$ then determines the value of a for the adopted opacity dependence.
(ii) $({\rm d}v/{\rm d}r)_{\rm s}$, $({\rm d}L(r)/{\rm d}r)_{\rm s}$ and $({\rm d}^2T/{\rm d}r^2)_{\rm s}$ are found from the application of the de l'Hopital's rule (Eq. (26)) combined with the derivatives of the temperature Eq. (16) and the energy Eq. (19). The value of $({\rm d} \chi /{\rm d}r)_{\rm s}$, which is needed in these solutions, is derived from the definition of $\chi=a \rho / T^{\alpha}$ together with $({\rm d}T/dr)_{\rm s}$ and $({\rm d}v/{\rm d}r)_{\rm s}$.
(iii) $({\rm d}^2v/{\rm d}r^2)_{\rm s}$ is found from the regularity condition (Eq. (27)) with $({\rm d}^3T/{\rm d}r^3)_{\rm s}$ from the second derivative of the temperature Eq. (16) and $({\rm d}^2 \chi/{\rm d}r^2)_{\rm s}$ from the adopted expression for $\chi$.
(iv) $({\rm d}L(r)/{\rm d}r)_{\rm s}$ is found from the energy Eq. (19).
(v) We then check if the energy equation and the momentum equation are both satisfied around the sonic point. To do this we calculate ${\rm d}v/{\rm d}r$ at a point slightly below the sonic point, at $r=0.99 R_{\rm s}$, from the momentum equation by using for the needed terms T, v and L (with $\chi(r)$ being known) the Taylor expansion formula which uses the values of these physical terms and their first and second derivatives at the sonic point. We then compare this with the value of ${\rm d}v/{\rm d}r$ at that same point, that is required by the energy equation, i.e. by differentiating Eq. (10). If this is not the case, our initial guess of L(r) was incorrect. We then choose another value of L(r) and repeat the process until convergence is reached.
(vi) When the model has converged, we derive the parameter $\beta$of the supersonic velocity law that produces the adopted value of $\taup_{\rm s}$. This is done in a similar way as in Variant A0, i.e. with stepwise integration of Eq. (33) for  $\mbox{$\tau^{\prime}$ }(r)$, with the adopted expression for $\chi(r)$.