Radiatively driven relativistic spherical winds under relativistic radiative transfer Abstract We numerically investigate radiatively driven relativistic spherical winds from the central luminous object with mass M and luminosity L* under Newtonian gravity, special relativity, and relativistic radiative transfer. We solve both the relativistic radiative transfer equation and the relativistic hydrodynamical equations for spherically symmetric flows under the double-iteration processes, to obtain the intensity and velocity fields simultaneously. We found that the momentum-driven winds with scattering are quickly accelerated near the central object to reach the terminal speed. The results of numerical solutions are roughly fitted by a relation of $$\dot{m}=0.7(\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2.6}$$, where $$\dot{m}$$ is the mass-loss rate normalized by the critical one, Γ* the central luminosity normalized by the critical one, τ* the typical optical depth, β* the initial flow speed at the central core of radius R*, and βout the terminal speed normalized by the speed of light. This relation is close to the non-relativistic analytical solution, $$\dot{m} = 2(\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2}$$, which can be re-expressed as $$\beta _{\rm out}^2/2 = (\Gamma _*-1)GM/c^2 R_*$$. That is, the present solution with small optical depth is similar to that of the radiatively driven free outflow. Furthermore, we found that the normalized luminosity (Eddington parameter) must be larger than unity for the relativistic spherical wind to blow off with intermediate or small optical depth, i.e. $$\Gamma _* \gtrsim \sqrt{(1+\beta _{\rm out})^3/(1-\beta _{\rm out})}.$$ We briefly investigate and discuss an isothermal wind. accretion, accretion discs, radiative transfer, relativistic processes, stars: winds, outflows, galaxies: active 1 INTRODUCTION Active astrophysical objects often exhibit outflow phenomena, which are supposed to be driven by thermal, magnetic, and radiative forces. Accretion discs usually play dominant roles in these active objects (Kato, Fukue & Mineshige 2008). Outflows from relativistic objects, such as neutron stars or black holes, have relativistic speeds up to 0.1 c in neutron star winds and up to ∼c in black hole winds. These relativistic outflows are believed to be driven by radiation pressure, when the central object is sufficiently luminous. Furthermore, except for the very vicinity of the central object, the winds are treated as spherically symmetric outflows, even if they emanate from a black hole accretion disc system, in a simple picture (King & Pounds 2003; King 2010; King & Muldrew 2016; see also Poutanen et al. 2007; Middleton et al. 2015). In order to examine optically thick radiatively driven relativistic spherical winds, researchers usually solved the relativistic hydrodynamical equations and relativistic moment ones with closure relations, instead of the relativistic radiative transfer equation (e.g. Ruggles & Bath 1979; Paczyński & Prószyński 1986; Paczyński 1986, 1990; Turolla & Nobili & Calvani 1986; Nobili, Turolla & Lapidus 1994; Akizuki & Fukue 2008, 2009). As already pointed out in Thorne, Flammang & Żytkow (1981), however, the diffusion approximation has a causality problem, whereas the Eddington approximation is also invalid in the relativistic regime (see e.g. Turolla & Nobili 1988; Nobili, Turolla & Zampieri 1991; Turolla, Zampieri & Nobili 1995; Dullemond 1999; Fukue 2005). Some researchers solved not moment equations, but the radiative transfer equation in spherically symmetric outflows in the flat space–time (e.g. Castor 1972; Mihalas, Kunasz & Hummber 1975, 1976a,b; Schinder & Bludman 1989; Zane et al. 1996; Baron & Hauschildt 2004; Chen et al. 2007; Knop, Hauschildt & Baron 2007), and in the curved one (Schinder & Bludman 1989; Zane et al. 1996; Chen et al. 2007; Knop et al. 2007 ). In these researches, however, the velocity laws were usually assumed, and the hydrodynamical equations were not solved. In order to solve the radiatively driven relativistic winds, we should treat both the relativistic hydrodynamical equations and relativistic radiative transfer equation simultaneously. Recently, under the double-iteration processes, both equations were simultaneously solved for the plane-parallel relativistic flows (Fukue 2015; cf. Fukue 2014), and for the spherical relativistic ones (Fukue 2017b; cf. Fukue 2017a). In these studies, particularly in Fukue (2017b) for a spherical case, however, the gravitational field of the central object was ignored for simplicity. Thus, in this paper, for the spherically symmetric case, we simultaneously solve both the relativistic radiative transfer equation and the relativistic hydrodynamical one for spherically symmetric winds under the gravitational effect of the central object. In the next section, we present the basic equations for radiation and gas. In Section 3, we show and discuss numerical solutions of scattering winds, whereas we briefly examine isothermal winds in Section 4. The final section is devoted to concluding remarks. 2 BASIC EQUATIONS Let us consider a steady spherically symmetric wind, which is driven by radiation force in the radial (r) direction from a central luminous object with mass M, radius R*, and luminosity L* (Fig. 1). The central object emits an isotropic uniform intensity I*. Hence, $$L_* = 4\pi R_*^2 \pi I_*$$, while the Eddington luminosity LE of the central object is LE ≡ 4πcGM/(κ0 + σ0), and the normalized central luminosity (Eddington parameter Γ*) is defined as Γ* ≡ L*/LE, where κ0 and σ0 are the absorption and scattering opacities measured in the comoving frame, respectively. Here and in what follows, the subscript 0 denotes the quantities in the comoving frame. Figure 1. View largeDownload slide Relativistic radiatively driven spherical wind in the radial (r) direction from a central luminous object with mass M, radius R*, luminosity L*. Both the spherical coordinates (r, θ, φ) and the impact space coordinates (p, z) are used. The flow velocity v and density ρ are functions of radius r, whereas the radiative intensity I±(p, z) is a function of (p, z); r2 = p2 + z2 and cos θ = μ = z/r. Figure 1. View largeDownload slide Relativistic radiatively driven spherical wind in the radial (r) direction from a central luminous object with mass M, radius R*, luminosity L*. Both the spherical coordinates (r, θ, φ) and the impact space coordinates (p, z) are used. The flow velocity v and density ρ are functions of radius r, whereas the radiative intensity I±(p, z) is a function of (p, z); r2 = p2 + z2 and cos θ = μ = z/r. We use the spherical coordinates (r, θ, φ) as well as the impact space coordinates (p, z), where the z-axis is in the line-of-sight direction of a distant observer (Hummer & Rybicki 1971). The transformation between these coordinates is as follows: r2 = p2 + z2 and cos θ = μ = z/r, where cos θ is the direction cosine from the z-axis. For the spherically symmetric flow, the flow velocity v and density ρ0 are functions of radius r, while the radiative intensity I±(p, z) in the inertial frame is a function of (p, z). As for the numerical processes, the outer radius of the calculation region is Rout; numerically, Rout = 100R*. The tangent rays are set in the logarithmic scale at pin ≤ p ≤ pout; numerically, pin = 0.01R* and pout = 100R*. Except for the gravity term, the basic equations are almost same as those in Fukue (2017b). For the convenience of readers, we briefly repeat them in the followings. 2.1 Relativistic radiative transfer equation The radiative transfer equations in various forms are derived and described in various works of literature (Chandrasekhar 1960; Mihalas 1970; Rybicki & Lightman 1979; Mihalas & Mihalas 1984; Kato, Fukue & Mineshige 1998; Mihalas & Auer 2001; Peraiah 2002; Castor 2004; Kato et al. 2008). The relativistic radiative transfer equation is also found in several works of literature (Mihalas 1980; Kato et al. 1998, 2008; Mihalas & Auer 2001). In the present treatment, the relativistic transfer equation in the mixed frame is   \begin{eqnarray} \pm \frac{\mathrm{\partial} I^{\pm }(p,z)}{\mathrm{\partial} z} &=& \frac{\rho _0}{\gamma ^3 (1-\beta \cdot l)^3} \left[ \frac{j_0}{4\pi } -\left(\kappa _0+\sigma _0\right) I_0 +\sigma _0 J_0 \right] \nonumber \\ &=& -\frac{(\kappa _0+\sigma _0)\rho _0}{\gamma ^3 (1-\beta \cdot l)^3} \left[ \gamma ^4 (1-\boldsymbol\beta \cdot \boldsymbol l)^4 I^{\pm } - S_0 \right], \end{eqnarray} (1)where I±(p, z) are the intensities in the upward and downward directions in the inertial (fixed) frame, respectively, and I0(p, z) and J0(r) the frequency-integrated specific intensity and mean intensity in the comoving (fluid) frames, respectively. In addition, β(r) (=v/c) is the normalized flow speed, $$\gamma =1/\sqrt{1-\beta ^2}$$ the Lorentz factor, and j0(r) the frequency-integrated mass emissivity. Furthermore, ($$\beta \cdot l$$) (=±β cos θ = ±βz/r) is the scalar product between the flow velocity $$\beta$$ (radial direction) and the ray direction $$l$$ (vertical direction) in the inertial frame. Finally, the source function S0(r) is defined in the comoving frame as   $$S_0(r) \equiv \frac{j_0}{4\pi }\frac{1}{\kappa _0+\sigma _0} + \frac{\sigma _0}{\kappa _0+\sigma _0}J_0.$$ (2)In the local thermodynamic equilibrium (LTE) case, j0/(4π) = κ0B0, where B0 ($$=\sigma T_0^4/\pi$$) is the blackbody intensity, T0 being the temperature in the comoving frame. It should be noted that the isotropic scattering is assumed for simplicity. The transformation rules for the radiative quantities are summarized as follows:   \begin{eqnarray} I_0 = \gamma ^4 (1-\beta \cdot l)^4 I, \end{eqnarray} (3)  \begin{eqnarray} J_0 = \gamma ^2 \left(J - 2\beta H + \beta ^2\,K \right), \end{eqnarray} (4)  \begin{eqnarray} H_0 = \gamma ^2 \left[ \left(1 + \beta ^2 \right)H - \beta \left(J + K \right) \right], \end{eqnarray} (5)  \begin{eqnarray} K_0 = \gamma ^2 \left( \beta ^2 J - 2\beta H + K \right), \end{eqnarray} (6)and the Eddington factor f is defined in the comoving frame as   $$f \equiv \frac{K_0}{J_0}.$$ (7) As was derived in Fukue (2017a), the relativistic formal solutions of equation (1) are   \begin{eqnarray} I^+(p,z) &=& {\rm e}^{ \displaystyle G(p,z_*)-G(p,z) - U(p,z_*)+U(p,z)} I^*(p,z_*) \nonumber \\ && + \int ^{z}_{0} \frac{ {\rm e}^{ \displaystyle G(p,\zeta ) - G(p,z) - U(p,\zeta ) + U(p,z)} }{\displaystyle \gamma ^3 \left(1 - \beta \frac{\zeta }{r} \right)^3} \nonumber \\ && \times\, (\kappa _0+\sigma _0)\rho _0 S_0 {\rm d}\zeta , \end{eqnarray} (8)  \begin{eqnarray} I^-(p,z) &=& - \int ^{z}_{z_{\rm out}} \frac{ {\rm e}^{ \displaystyle [G(p,z)-G(p,\zeta )] + [U(p,z)-U(p,\zeta )] }}{\displaystyle \gamma ^3 \left(1 + \beta \frac{\zeta }{r} \right)^3} \nonumber \\ && \times \, (\kappa _0+\sigma _0)\rho _0 S_0 {\rm d}\zeta , \end{eqnarray} (9)where $$z_* = \sqrt{R_*^2-p^2}$$, $$z_{\rm out} =\sqrt{R_{\rm out}^2-p^2}$$, and   \begin{eqnarray} G(p,z) \equiv \int ^z (\kappa _0+\sigma _0)\rho _0 \gamma {\rm d}\zeta , \end{eqnarray} (10)  \begin{eqnarray} U(p,z) \equiv \int ^z (\kappa _0+\sigma _0)\rho _0 \gamma (\beta \cdot l) {\rm d}\zeta . \end{eqnarray} (11) The boundary condition for the downward intensity is   $$I^-(p,z_{\rm out}) = 0 {\rm \ \ \ at\ \ \, } z=z_{\rm out}=\sqrt{R_{\rm out}^2-p^2},$$ (12)whereas that for the upward intensity is   \begin{eqnarray} I^+(p,z_*) &=& I^* {\rm \ \ \ for\ \ \, } p < R_*, \end{eqnarray} (13)  \begin{eqnarray} I^+(p,0) &=& I^-(p,0) {\rm \ \ \ for\ \ \, } p > R_*. \end{eqnarray} (14) Using the radius R*, opacity value κ*(R*), and density ρ*(R*) there, equations (10) and (11) are normalized as   \begin{eqnarray} G(p,z) \equiv \tau _* \int ^z \frac{\kappa _0 + \sigma _0}{\kappa _*} \frac{\rho _0}{\rho _*} \gamma {\rm d}\tilde{\zeta }, \end{eqnarray} (15)  \begin{eqnarray} U(p,z) \equiv \tau _* \int ^z \frac{\kappa _0 + \sigma _0}{\kappa _*} \frac{\rho _0}{\rho _*} \gamma \beta \cdot l{\rm d}\tilde{\zeta }, \end{eqnarray} (16)where $$\tilde{\zeta }=\zeta /R_*$$ and   $$\tau _* = \kappa _* \rho _* R_*,$$ (17)which is a typical optical depth of the wind. For a given velocity profile v and source function S0, the formal solutions (8) and (9) are numerically integrated under the appropriate boundary conditions and parameters, and the refined radiative quantities are obtained. 2.2 Relativistic hydrodynamical equation The relativistic hydrodynamical equations are also found in several works of literature (e.g. Kato et al. 1998, 2008). For the steady spherical wind under the special relativity, the continuity equation becomes   $$4 \pi r^2 \rho _0 \gamma v = 4 \pi r^2 \rho _0 c \gamma \beta = \dot{M},$$ (18)where $$\dot{M}$$ is the mass-loss rate and assumed to be constant. The relativistic equation of motion is   $$c^2 u \frac{{\rm d}u}{{\rm d}r} = c^2 \gamma ^4 \beta \frac{{\rm d}\beta }{{\rm d}r} = -\frac{GM}{r^2} + \frac{\kappa _0+\sigma _0}{c}\gamma 4\pi H_0,$$ (19)where u (=γβ) is the four velocity, and the gas pressure-gradient force is ignored in the present case. We do not use the energy equation for matter. Using the quantities on the central object, radius R*, density ρ*, velocity β* (corresponding Lorentz factor γ*), and luminosity L*, hydrodynamical equations (18) and (19) are respectively normalized as follows:   \begin{eqnarray} \hat{r}^2 \frac{\rho _0}{\rho _*}\frac{\gamma \beta }{\gamma _* \beta _*} &=& 1, \end{eqnarray} (20)  \begin{eqnarray} \dot{m} \gamma ^3 \frac{{\rm d}\beta }{{\rm d}\hat{r}} = \tau _* \frac{\kappa _0 + \sigma _0}{\kappa _*} \frac{\rho _0}{\rho _*} \left( -1 + \Gamma _* \gamma \frac{4\pi H_0}{\pi I_*} \hat{r}^2 \right), \end{eqnarray} (21)where $$\hat{r}=r/R_*$$, and $$\dot{m}$$ ($$\equiv \dot{M}c^2/L_{\rm E}$$) the normalized mass-loss rate, i.e. the critical mass-loss rate is defined as $$\dot{M}_{\rm crit} \equiv L_{\rm E}/c^2$$ in this work. For a given flux profile H0 in the comoving frame, these hydrodynamical equations are numerically integrated under the appropriate boundary conditions and parameters, and the refined velocity profile is obtained. 2.3 Relativistic radiative moment equations In order to analyse the numerical results, we here present the relativistic radiative moment equations (see e.g. Kato et al. 1998, 2008), although we do not solve them in the present numerical calculations. Using the radiative moment quantities in the comoving and inertial frames, the zeroth and first moment equations for the spherically symmetric case respectively become   \begin{eqnarray} \frac{1}{r^2}\frac{{\rm d}}{{\rm d}r}(r^2H) &=& \rho _0 \gamma \left( \frac{j_0}{4\pi } - \kappa _0 J_0 \right) -\rho _0 \gamma (\kappa _0 + \sigma _0)\beta H_0, \nonumber \\ &=& \rho _0 \gamma \left\lbrace\vphantom{\gamma ^2(1+\beta ^2)} \kappa _0 (B_0 - J + \beta H) \right. \nonumber \\ && \left. -\, \sigma _0 \gamma ^2 \beta [(1+\beta ^2)H-\beta (J+K)] \right\rbrace , \end{eqnarray} (22)  \begin{eqnarray} \frac{{\rm d}K}{{\rm d}r} + \frac{3K-J}{r} &=& -\rho _0 \gamma (\kappa _0 + \sigma _0)H_0 +\rho _0 \gamma \beta \left( \frac{j_0}{4\pi } - \kappa _0 J_0 \right) \nonumber \\ &=& -\rho _0 \gamma \left\lbrace\vphantom{(1+\beta ^2)} \kappa _0 (H - \beta B_0 - \beta K) \right. \nonumber \\ && \left. +\, \sigma _0 \gamma ^2 [(1+\beta ^2)H-\beta (J+K)] \right\rbrace , \end{eqnarray} (23)where we use j0/(4π) = κ0B0 in the second lines of both equations. In addition, for a closure relation K0/J0 = f in equation (7), we have a relation   $$(1-f\beta ^2)K = (f-\beta ^2)J + (1-f)2\beta H,$$ (24)in the inertial frame. It should be stressed that if we use this closure relation (24) in the inertial frame, moment equations in the inertial frame have a pathological singularity at β2 = f = 1/3 under the usual Eddington approximation (cf. Fukue 2005). 3 MOMENTUM-DRIVEN SCATTERING WINDS We first examine the purely scattering case with no emission and absorption (j0 = 0 and κ0 = 0), i.e. S0 = J0. In this case, we do not consider any thermal process, and the wind flow is momentum driven. 3.1 Double-iteration processes Similar to the previous studies (Fukue 2015, 2017b), we need double-iteration processes: one is the iteration for radiation field, and the other is that between radiation and gas. We here outline the procedure. Cycle 0 (Initial trial): In the optically thin non-relativistic case, the radiative moments are roughly expressed as J ∝ r−3 near the centre and H ∝ r−2. By considering non-relativistic moment equations for radiation field, we adopt the following analytical solutions:   \begin{eqnarray} J = \tau _* I_* \frac{1}{4\hat{r}^3}, \end{eqnarray} (25)  \begin{eqnarray} H = I_* \frac{1}{4\hat{r}^2}, \end{eqnarray} (26)as initial trial functions for radiation field. For velocity field, we adopt an approximate solution of equation (21) in the non-relativistic limit (see later),   $$\beta = \beta _{\rm out}\sqrt{1-\frac{1}{\hat{r}}},$$ (27)as an initial trial function, where βout is the terminal speed of the flow.  For this initial velocity field and the resultant density distribution, the functions (15) and (16) are tabulated as initial conditions. Cycle 1A (Radiation field):Using initial trial functions (25) and (27) with tabulated functions (15) and (16), we iteratively integrate the relativistic formal solutions (8) and (9) for given parameters, until the solutions of the specific intensities I± converge. Using the resultant intensities and transformation to the comoving frame, we obtain the comoving flux H0(r). Cycle 1B (Velocity field): Using the Eddington flux H0 obtained by Cycle 1A, we solve equation (21) with equation (20) for a given terminal speed βout. As a result, we obtain the refined velocity distribution, and the mass-loss rate $$\dot{m}$$ as an eigenvalue. For this renewed velocity field and the resultant density distribution, the functions (15) and (16) are re-tabulated for the next step. Cycles 2 to L:We then repeat Cycles 1A and 1B, until both the specific intensity and velocity distributions converge. The parameters for this momentum-driven wind are the Eddington parameter Γ*, the typical optical depth τ*, the initial flow speed β* at the central core, the ‘terminal’ speed βout at the outer radius, and the normalized mass-loss rate $$\dot{m}$$. Of these, we set β* = 0.01, as a sufficiently small value. We further set τ* = 3 for a typical case. In addition, the mass-loss rate is determined as an eigenvalue. Hence, we examine the present problem for several values of Γ* and βout. 3.2 Typical solutions The present results are somewhat similar to the previous study without gravity (Fukue 2017b), where the main parameters are τ* and βout. Hence, we consult the previous study on the optical depth dependence, and in this subsection we mainly show the dependence on Γ* and βout. It should be noted that we can obtain the solutions in the limited range of βout for some given Γ*, which is the important difference from the previous study and discussed later in this section. For the case of τ* = 3, examples of numerical solutions are shown in Figs 2 –4. In Fig. 2, the velocity profiles are shown as a function of the normalized radius $$\hat{r}$$ for several values of βout in the case of Γ* = 5 and τ* = 3, whereas the initial trial velocity profiles are plotted by dashed curves for comparison. As shown in Fig. 2, and similar to the previous study (Fukue 2017b), the momentum-driven wind is quickly accelerated near the central object to reach the terminal speed. Figure 2. View largeDownload slide Velocity profiles as a function of radius for several values of βout. The initial trial velocity profiles are plotted by dashed curves for comparison. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5 from bottom to top. Figure 2. View largeDownload slide Velocity profiles as a function of radius for several values of βout. The initial trial velocity profiles are plotted by dashed curves for comparison. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5 from bottom to top. In order to understand this behaviour, similar to the previous study (Fukue 2017b), we again consider the approximate analytical solution of equation (21) in the non-relativistic limit. That is, using equation (20), equation (21) is approximately expressed as   $$\dot{m}\frac{{\rm d}\beta }{{\rm d}\hat{r}} = \frac{\tau _* \beta _*}{\hat{r}^2 \beta } \left( \Gamma _* \frac{4\pi H_0}{\pi I_*}\hat{r}^2 - 1 \right) \sim (\Gamma _* - 1) \frac{\tau _* \beta _*}{\beta } \frac{1}{\hat{r}^2},$$ (28)where we use the relation for the flux in the non-relativistic limit: $$4\pi H_0 4\pi r^2=\pi I_* 4\pi R_*^2 =L_*$$. This equation (28) is easily integrated to yield the velocity profile in the non-relativistic limit,   \begin{eqnarray} \frac{1}{2}\beta ^2 &=& \frac{1}{2}\beta _{\rm out}^2 - (\Gamma _* - 1) \frac{\tau _* \beta _*}{\dot{m}} \frac{1}{\hat{r}}, \end{eqnarray} (29)  \begin{eqnarray} \qquad &=& \frac{1}{2}\beta _{\rm out}^2 - (\Gamma _* - 1) \frac{GM}{c^2r}, \end{eqnarray} (30)where we impose the boundary condition at the outer radius: β = βout at $$\hat{r} = \infty$$. It should be noted that this analytical solution in the non-relativistic limit is just a simple radiatively driven free outflow. We further impose the boundary condition at the inner radius: β = β* ≪ βout at $$\hat{r}=1$$, and we obtain the relation   $$\dot{m} = \frac{2(\Gamma _* - 1)\tau _* \beta _*}{\beta _{\rm out}^2}.$$ (31)Inserting this relation into the approximate solution, we finally obtain a simple analytical solution in the non-relativistic limit:   $$\beta = \beta _{\rm out}\sqrt{1-\frac{1}{\hat{r}}}.$$ (32) For small βout, this analytical solution (32) well reproduces the numerical solutions in Fig. 2. That is, the present solution is quite similar to the radiatively driven outflow in the non-relativistic limit. For large βout, however, it slightly deviates from the numerical ones. We emphasize that this analytical solution does not include any parameters except for βout. In Fig. 3, the radiative quantities normalized by I* are shown as a function of r normalized by R* for several values of βout (=0.1, 0.3, 0.5) in the case of Γ* = 5 and τ* = 3. Overall behaviour of radiative quantities is similar to the previous studies (Fukue 2017a,b). Figure 3. View largeDownload slide Radiative quantities as a function of radius for several values of βout. (a) Normalized mean intensity, (b) normalized flux, and (c) the Eddington factor. The solid curves represent the quantities in the inertial frame, while the dashed ones mean those in the comoving frame. Radius is normalized by R*, while radiative quantities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5. Figure 3. View largeDownload slide Radiative quantities as a function of radius for several values of βout. (a) Normalized mean intensity, (b) normalized flux, and (c) the Eddington factor. The solid curves represent the quantities in the inertial frame, while the dashed ones mean those in the comoving frame. Radius is normalized by R*, while radiative quantities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5. In Fig. 3(a), the mean intensity J is plotted as a function of radius $$\hat{r}$$. Similar to the previous studies (Fukue 2017a,b), the mean intensity in the inertial frame (solid curves) roughly decreases as J ∝ r−2, whereas that in the comoving frame (dashed curves) reduces in the entire region, due to redshift (Fukue 2014, 2017a). In Fig. 3(b), the Eddington flux H is plotted as a function of radius $$\hat{r}$$. Similar to the previous studies (Fukue 2017a,b), the flux in the inertial frame (solid curves) roughly decreases as H ∝ r−2, whereas that in the comoving frame (dashed curves) reduces in the entire region, as the flow velocity increases (Fukue 2017a). In Fig. 3(c), the Eddington factor f is plotted as a function of radius $$\hat{r}$$. Similar to the previous study without gravity (Fukue 2017b), the Eddington factor is almost 1/3 near the core, and increases towards unity as radius becomes large, due to the geometrical and optical depth effects. In Fig. 4, the emergent intensities I+(p, zout) at the outer boundary Rout normalized by I* are plotted as a function of the impact coordinate p normalized by R*. The normalized luminosity Γ* = 5, the typical optical depth τ* is 3, and the values of βout are 0.1, 0.3, 0.5 from outer to inner. Figure 4. View largeDownload slide Emergent intensities in the inertial frame at the outer radius Rout. Radius is normalized by R*, while the intensities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5 from outer to inner. Figure 4. View largeDownload slide Emergent intensities in the inertial frame at the outer radius Rout. Radius is normalized by R*, while the intensities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5 from outer to inner. Since the optical depth is not so large (τ* = 3), the emergent intensity in the core part of p/R* ≤ 1 is almost unity; this part is roughly the incident intensity I*. In the outer envelope of p/R* > 1, on the other hand, the emergent intensity decreases with p due to the scattering effect, and it is less enhanced as the flow speed increases due to the transverse Doppler effect. Overall behaviour of emergent intensities is quite similar to the previous study (Fukue 2017b). 3.3 Relations on the parameters Since the mass-loss rate $$\dot{m}$$ is obtained as an eigenvalue, similar to the previous study (Fukue 2017b), there exist some relations among the normalized luminosity Γ*, the typical optical depth τ*, the normalized mass-loss rate $$\dot{m}$$, and the normalized flow terminal speed βout. The obtained relations between $$\dot{m}$$ and βout for several values of Γ* in the case of τ* = 3 are plotted in Fig. 5. Figure 5. View largeDownload slide Relations between the normalized mass-loss rate $$\dot{m}$$ and the flow terminal speed βout for several values of the normalized luminosity Γ* in the case of τ* = 3 (solid curves with open circles). The dashed lines denote a fitting relation of $$\dot{m} = 0.7 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2.6}$$ for the case of τ* = 3. The dot–dashed lines represent an analytical relation of $$\dot{m}=2 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2}$$ in the non-relativistic regime. The curve with crosses expresses a rough limit of the existence of solutions. Figure 5. View largeDownload slide Relations between the normalized mass-loss rate $$\dot{m}$$ and the flow terminal speed βout for several values of the normalized luminosity Γ* in the case of τ* = 3 (solid curves with open circles). The dashed lines denote a fitting relation of $$\dot{m} = 0.7 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2.6}$$ for the case of τ* = 3. The dot–dashed lines represent an analytical relation of $$\dot{m}=2 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2}$$ in the non-relativistic regime. The curve with crosses expresses a rough limit of the existence of solutions. As is seen in Fig. 5, similar to the plane-parallel case (Fukue 2015) and to the previous spherical one (Fukue 2017b), the mass-loss rate decreases as the flow terminal speed increases, while it increases with the normalized luminosity or with the typical optical depth. Similar to the previous study (Fukue 2017b), these trends are approximated by a fitting relation to the numerical results, and also by an analytical solution, using the back-of-the-envelope calculation. For example, the dashed lines in Fig. 5 denote a fitting relation of $$\dot{m} = 0.7 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2.6}$$ for the case of τ* = 3. On the other hand, the dot–dashed lines in Fig. 5 represent the analytical solution (31) of $$\dot{m}=2(\Gamma _* - 1)\tau _* \beta _* \beta _{\rm out}^{-2}$$ for the case of τ* = 3, which is obtained during the derivation of the approximate solution in the non-relativistic limit. It should be noted that the normalized mass-loss rate $$\dot{m}$$ apparently depends on parameters under the present normalization, but the actual mass-loss rate is not determined in the present way, where the velocity field is calculated. Now, we shall discuss the limited range of solutions. As was already stated, and as seen in Figs 2–5, for given Γ* and τ*, the numerical solutions can be obtained only for the lower values of βout, and it is hard to obtain solutions when the value of βout becomes high. This nature is remarkably different from the results in the previous study without gravity (Fukue 2017b). This behaviour is mainly due to the physical reason, and partly due to the numerical technique. Physically, for large values of the flow speed, the net force on the right-hand side of equation (21) becomes negative even for Γ* > 1, since the comoving flux decreases as the flow speed increases (Fig. 3b). This will be analytically shown below. Numerically, as the flow speed increases, even in the case that the numerical solution may exist, the net force temporary becomes negative during the iteration process, and the solution search is stopped. This limited range is approximately derived as follows. As is seen in Figs 2 and 3, at far from the centre, the flow becomes constant, and the Eddington factor is almost unity due to the geometrical and optical depth effects. In this region of f = 1, K0 = J0 and further K = J. Moreover, the radiative moments in the inertial frame roughly become   \begin{eqnarray} \frac{4\pi H}{\pi I_*} \hat{r}^2 &=& 1, \end{eqnarray} (33)  \begin{eqnarray} \frac{4\pi J}{\pi I_*} \hat{r}^2 &=& 1. \end{eqnarray} (34)Hence, the Eddington flux in the comoving frame becomes   \begin{eqnarray} \frac{4\pi H_0}{\pi I_*} \hat{r}^2 &=& \frac{4\pi }{\pi I_*} \hat{r}^2 \gamma ^2[(1+\beta ^2)H-\beta (J+K)] \nonumber \\ &=& \frac{1-\beta }{1+\beta } \ \ \ \le 1. \end{eqnarray} (35)Inserting this approximate solution into the net force on the right-hand side of equation (21), we have a condition that the net force vanishes:   $$\Gamma _* = \sqrt{\frac{(1+\beta )^3}{1-\beta }} \ \ \ \ge 1.$$ (36)Thus, the normalized luminosity (Eddington parameter) must be larger than unity for the relativistic spherical wind with intermediate or small optical depth. Finally, inserting this condition (36) into the fitting relation of $$\dot{m}=0.7 (\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2.6}$$, we have the limiting condition on the $$\dot{m}$$–βout plane, which is plotted in Fig. 5 (cross marks). As is seen in Fig. 5, this approximate condition well explains the existence domain of the present numerical solutions. This limitation is interpreted as follows. We do not consider the energy equation for matter, and therefore, the kinetic energy flux of matter is ignored compared with the energy flux in the luminosity. However, in the gaseous wind, if the mass-loss rate becomes large and the wind is heavy, the kinetic energy flux of matter cannot be negligible, and the radiation field with a given luminosity cannot accelerate matter to sufficient speed at infinity. In other words, in order to accelerate matter to sufficient speed at infinity, the sufficient luminosity is required, as in equation (35). 4 ISOTHERMAL WINDS We next consider the isothermal case without scattering (σ0 = 0) under LTE (j0/4π = κ0B0), i.e. S0 = B0 = constant. In this case, the wind flow is also momentum driven, although the momentum source is distributed in the entire space. Since the iteration process for radiation field is unnecessary in this isothermal case, we need a single iteration process between radiation and gas. The parameters for this isothermal wind are Γ*, τ*, β*, βout, $$\dot{m}$$, and the value of B0. Of these, we set β* = 0.01 and τ* = 3 for a typical case. The value of B0 is set as B0 = I*. The mass-loss rate is determined as an eigenvalue. Hence, the rests are Γ* and βout. 4.1 Typical solutions Examples of numerical solutions for isothermal winds are shown in Figs 6 –8. Figure 6. View largeDownload slide Velocity profiles as a function of radius for several values of βout in the isothermal case. The initial trial velocity profiles are plotted by dashed curves for comparison. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from bottom to top. Figure 6. View largeDownload slide Velocity profiles as a function of radius for several values of βout in the isothermal case. The initial trial velocity profiles are plotted by dashed curves for comparison. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from bottom to top. Figure 7. View largeDownload slide Radiative quantities as a function of radius for several values of βout in the isothermal case. (a) Normalized mean intensity, (b) normalized flux, and (c) the Eddington factor. The solid curves represent the quantities in the inertial frame, while the dashed ones mean those in the comoving frame. Radius is normalized by R*, while radiative quantities are normalized by I*. The normalized luminosity Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5. Figure 7. View largeDownload slide Radiative quantities as a function of radius for several values of βout in the isothermal case. (a) Normalized mean intensity, (b) normalized flux, and (c) the Eddington factor. The solid curves represent the quantities in the inertial frame, while the dashed ones mean those in the comoving frame. Radius is normalized by R*, while radiative quantities are normalized by I*. The normalized luminosity Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5. Figure 8. View largeDownload slide Emergent intensities in the inertial frame at the outer radius Rout in the isothermal case. Radius is normalized by R*, while the intensities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from top to bottom. Figure 8. View largeDownload slide Emergent intensities in the inertial frame at the outer radius Rout in the isothermal case. Radius is normalized by R*, while the intensities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from top to bottom. In Fig. 6, the velocity profiles are shown as a function of the normalized radius $$\hat{r}$$ for several values of βout in the case of Γ* = 5, τ* = 3, and B0/I* = 1, whereas the initial trial velocity profiles are plotted by dashed curves for comparison. In contrast to the scattering case in Fig. 2, the isothermal wind is gradually accelerated near the central object to reach the terminal speed. This is because the momentum source is uniformly distributed in the entire space. In Fig. 7, the radiative quantities normalized by I* are shown as a function of r normalized by R* for several values of βout (=0.1, 0.3, 0.5) in the case of Γ* = 5, τ* = 3, and B0/I* = 1. In Fig. 7(a), the mean intensity J is plotted as a function of radius $$\hat{r}$$. In contrast to the scattering case (Fig. 3a), the mean intensity in the inertial frame (solid curves) and that in the comoving one (dashed ones) roughly decrease as J ∝ r−1. This behaviour originates from the difference in source distributions. In the scattering wind (Fig. 3a), the radiation source exists only on the central core, whereas the radiation source is uniformly distributed in the isothermal wind (Fig. 7a). Moreover, the mean intensity decreases as the flow speed increases, since the density and the source intensity (emissivity) decrease as the flow speed increases. In Fig. 7(b), the Eddington flux H is plotted as a function of radius $$\hat{r}$$. Due to the same reason, the flux also decreases as H ∝ r−1. In Fig. 7(c), the Eddington factor f is plotted as a function of radius $$\hat{r}$$. In contrast to the scattering case (Fig. 3c), the Eddington factor in the present isothermal wind does not approach unity, but becomes around 1/2 or so. This means that the radiation field of the isothermal wind is between the sufficiently thick isotropic case (f = 1/3) and the intermediately thin spherical case (f ∼ 1). In order to confirm these numerical solutions for the radiative moment quantities, let us examine the relativistic moment equations. When the flow speed is small and of the order of $${\cal O}(\beta ^1)$$, the zeroth moment equation (22), the first one (23), and the closure relation (24) respectively become in the isothermal case without scattering   \begin{eqnarray} \frac{1}{\hat{r}^2}\frac{{\rm d}}{{\rm d}\hat{r}}(\hat{r}^2H) &=& \tau _* \frac{\rho _0}{\rho _*} (B_0 - J + \beta H), \end{eqnarray} (37)  \begin{eqnarray} \frac{{\rm d}K}{{\rm d}\hat{r}} + \frac{3K-J}{\hat{r}} &=& -\tau _* \frac{\rho _0}{\rho _*} (H - \beta K - \beta B_0), \end{eqnarray} (38)  \begin{eqnarray} K &=& fJ + (1-f)2\beta H. \end{eqnarray} (39) At far from the centre, J and H are smaller than B0. Then, using continuity equation (20), equation (37) becomes   $$\frac{1}{\hat{r}^2}\frac{{\rm d}}{{\rm d}\hat{r}}(\hat{r}^2H) \sim \tau _* \frac{\beta _*}{\hat{r}^2 \beta } B_0.$$ (40)At far from the centre, the flow speed is almost constant, and this equation is solved to yield a solution   $$H \sim \frac{\tau _* \beta _*}{\beta }\frac{B_0}{\hat{r}},$$ (41)and the flux roughly varies as $$H \propto \hat{r}^{-1}$$. Similarly, at far from the centre, where B0 is dominated and the flow speed is constant, equation (38) with equation (39) is also solved to yield an approximate solution:   \begin{eqnarray} J &\sim & \tau _* \beta _* \frac{B_0}{\hat{r}} \ \ \ \ \ \ \ \ \ \ \ {\rm for}\ \ \ f=1, \end{eqnarray} (42)  \begin{eqnarray} J &\sim & 6\tau _* \beta _* \frac{B_0}{\hat{r}}\ln \hat{r} \ \ \ {\rm for}\ \ \ f=1/2, \end{eqnarray} (43)and the mean intensity roughly varies as $$J \propto \hat{r}^{-1}$$ at least on the main dependence. In Fig. 8, the emergent intensities I+(p, zout) at the outer boundary Rout normalized by I* are plotted as a function of the impact coordinate p normalized by R*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from top to bottom. In contrast to the scattering case (Fig. 4), the emergent intensity in the outer envelope is quite large. This is also because the source emissivity is uniformly distributed in the entire space in the isothermal wind. Similar to the relations obtained in the scattering case (Fig. 5), in this isothermal wind we can also the relations on the parameters, although we do not discuss them for simplicity. Finally, we should mention the opacity effect. In this isothermal wind, we assume that the absorption opacity κ0 is constant. In the realistic case, the absorption opacity depends on the physical quantities, e.g. the free–free opacity is proportional to $$\rho _0 T_0^{-7/2}$$. Hence, it is   $$\frac{\kappa _0}{\kappa _*} = \frac{\rho _0}{\rho _*} \frac{T_0}{T_*} = \frac{\rho _0}{\rho _*}$$ (44)in the isothermal case. In this case, the results are quite similar to the scattering case. The reason is because the density is concentrated in the centre and quickly decreases towards the outer envelope, and therefore, the emissivity is also concentrated in the central core. 5 CONCLUDING REMARKS In this paper, we have examined the radiatively driven relativistic spherical winds from the central object with mass M and luminosity L* under Newtonian gravity and special relativity, by numerically solving the relativistic radiative transfer equation and relativistic hydrodynamical ones for spherical geometry, simultaneously and consistently. Namely, using the initial trial velocity profile, we iteratively solved the relativistic formal solution of the transfer equation to obtain specific intensities and moment quantities in the inertial and comoving frames. Using the comoving flux obtained, we then solved the relativistic hydrodynamical equations to obtain the refined velocity profile, and the mass-loss rate as an eigenvalue. We repeat these double-iteration procedures until both the intensity and velocity converge. We found that the momentum-driven winds with scattering are quickly accelerated near the central luminous core to reach the terminal speed (Fig. 2), similar to the radiative outflows without gravity (Fukue 2017b). The behaviour of moment quantities and emergent intensity is quite similar to the previous results (figs 3 and 4 in Fukue 2017b). The results of numerical solutions are roughly fitted by a relation of $$\dot{m}=0.7(\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2.6}$$, which is well reproduced by the analytical relation in the non-relativistic limit $$\dot{m} = 2(\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2}$$ (Fig. 5). It should be noted that the normalized mass-loss rate $$\dot{m}$$ apparently depends on parameters under the present normalization, but the actual mass-loss rate is not determined in the present way, where the velocity field is calculated. Furthermore, we found that the normalized luminosity (Eddington parameter) must be larger than unity for the relativistic spherical wind with intermediate or small optical depth: $$\Gamma _* \gtrsim \sqrt{(1+\beta _{\rm out})^3/(1-\beta _{\rm out})}$$ (Fig. 5). We briefly investigate and discuss an isothermal wind. We found that the isothermal winds are gradually accelerated, since the emissivity is uniformly distributed in the entire space. The advantage of the present study is twofold. The first is the self-consistency; in this work, both the radiation and flow fields are obtained simultaneously. The second is the non-singularity; a pathological singularity under the Eddington approximation (e.g. Dullemond 1999; Fukue 2005) does not appear, since we solve the full transfer equation, instead of the truncated moment equations under the closure relation (24). On the other hand, the present solution has several limitations. For example, we do not consider the energy equation for matter, and therefore, the kinetic energy flux of matter is ignored compared with the energy flux in the luminosity. Hence, the present solution implicitly treats the light wind, and cannot treat the heavy wind. Furthermore, the present solution is limited in the intermediate or small optical depth (τ ∼ 10). If the optical depth becomes large, the diffusion time-scale is longer than the advection time-scale, and the present treatment would break down. In this paper, we have considered the gravitational effect, but the gas pressure has not yet been included. If we take into account the gas pressure-gradient force in the equation of motion, we must seek transonic solutions, which is left as a future work. Acknowledgements The author would like to thank an anonymous referee for valuable comments, including the interpretation of solutions. REFERENCES Akizuki C., Fukue J., 2008, PASJ , 60, 337 https://doi.org/10.1093/pasj/60.2.337 CrossRef Search ADS   Akizuki C., Fukue J., 2009, PASJ , 61, 543 https://doi.org/10.1093/pasj/61.3.543 CrossRef Search ADS   Baron E., Hauschildt P. H., 2004, A&A , 427, 987 https://doi.org/10.1051/0004-6361:20040039 CrossRef Search ADS   Castor J. I., 1972, ApJ , 178, 779 https://doi.org/10.1086/151834 CrossRef Search ADS   Castor J. I., 2004, Radiation Hydrodynamics . Cambridge Univ. Press, Cambridge Google Scholar CrossRef Search ADS   Chandrasekhar S., 1960, Radiative Transfer . Dover Press, New York Chen B., Kantowski R., Baron E., Knop S., Hauschildt P. H., 2007, MNRAS , 380, 104 https://doi.org/10.1111/j.1365-2966.2007.11652.x CrossRef Search ADS   Dullemond C. P., 1999, A&A , 343, 1030 Fukue J., 2005, PASJ , 57, 1023 https://doi.org/10.1093/pasj/57.6.1023 CrossRef Search ADS   Fukue J., 2014, PASJ , 66, 73 https://doi.org/10.1093/pasj/psu048 CrossRef Search ADS   Fukue J., 2015, PASJ , 67, 14 Fukue J., 2017a, PASJ , 69, 8 https://doi.org/10.1093/pasj/psw112 CrossRef Search ADS   Fukue J., 2017b, PASJ , 69, 53 Hummer D. G., Rybicki G. B., 1971, MNRAS , 152, 1 https://doi.org/10.1093/mnras/152.1.1 CrossRef Search ADS   Kato S., Fukue J., Mineshige S., 1998, Black-Hole Accretion Disks . Kyoto Univ. Press, Kyoto Kato S., Fukue J., Mineshige S., 2008, Black-Hole Accretion Disks – Towards a New Paradigm . Kyoto Univ. Press, Kyoto King A. R., 2010, MNRAS , 402, 1516 https://doi.org/10.1111/j.1365-2966.2009.16013.x CrossRef Search ADS   King A. R., Muldrew S. I., 2016, MNRAS , 455, 1211 https://doi.org/10.1093/mnras/stv2347 CrossRef Search ADS   King A. R., Pounds K. A., 2003, MNRAS , 345, 657 https://doi.org/10.1046/j.1365-8711.2003.06980.x CrossRef Search ADS   Knop S., Hauschildt P. H., Baron E., 2007, A&A , 463, 315 https://doi.org/10.1051/0004-6361:20065521 CrossRef Search ADS   Middleton M. J., Heil L., Pintore F., Watson D. J., Roberts P., 2015, MNRAS , 447, 3243 https://doi.org/10.1093/mnras/stu2644 CrossRef Search ADS   Mihalas D., 1970, Stellar Atmospheres . Freeman & Co., San Francisco Mihalas D., 1980, ApJ , 237, 574 https://doi.org/10.1086/157902 CrossRef Search ADS   Mihalas D., Auer L. H., 2001, J. Quant. Spectrosc. Radiat. Transfer , 71, 61 https://doi.org/10.1016/S0022-4073(01)00013-9 CrossRef Search ADS   Mihalas D., Mihalas B. W., 1984, Foundations of Radiation Hydrodynamics . Oxford Univ. Press, Oxford Mihalas D., Kunasz P. B., Hummber D. G., 1975, ApJ , 202, 465 https://doi.org/10.1086/153996 CrossRef Search ADS   Mihalas D., Kunasz P. B., Hummber D. G., 1976a, ApJ , 206, 515 https://doi.org/10.1086/154407 CrossRef Search ADS   Mihalas D., Kunasz P. B., Hummber D. G., 1976b, ApJ , 210, 419 https://doi.org/10.1086/154845 CrossRef Search ADS   Nobili L., Turolla R., Zampieri L., 1991, ApJ , 383, 250 https://doi.org/10.1086/170781 CrossRef Search ADS   Nobili L., Turolla R., Lapidus I., 1994, ApJ , 433, 276 https://doi.org/10.1086/174643 CrossRef Search ADS   Paczyński B., 1986, ApJ , 308, L43 https://doi.org/10.1086/184740 CrossRef Search ADS   Paczyński B., 1990, ApJ , 363, 218 https://doi.org/10.1086/169332 CrossRef Search ADS   Paczyński B., Prószyński M., 1986, ApJ , 302, 519 https://doi.org/10.1086/164012 CrossRef Search ADS   Peraiah A., 2002, An Introduction to Radiative Transfer: Methods and Applications in Astrophysics . Cambridge Univ. Press, Cambridge Poutanen J., Lipunova G., Fabrika S., Butkevich A. G., Abolmasov P., 2007, MNRAS , 377, 1187 https://doi.org/10.1111/j.1365-2966.2007.11668.x CrossRef Search ADS   Ruggles C. L. N., Bath G. T., 1979, A&A , 80, 97 Rybicki G. B., Lightman A. P., 1979, Radiative Processes in Astrophysics . Wiley, New York Schinder P. J., Bludman S. A., 1989, ApJ , 346, 350 https://doi.org/10.1086/168015 CrossRef Search ADS   Thorne K. S., Flammang R. A., Żytkow A. N., 1981, MNRAS , 194, 475 https://doi.org/10.1093/mnras/194.2.475 CrossRef Search ADS   Turolla R., Nobili L., 1988, MNRAS , 235, 1273 https://doi.org/10.1093/mnras/235.4.1273 CrossRef Search ADS   Turolla R., Nobili L., Calvani M., 1986, ApJ , 303, 573 https://doi.org/10.1086/164103 CrossRef Search ADS   Turolla R., Zampieri L., Nobili L., 1995, MNRAS , 272, 625 https://doi.org/10.1093/mnras/272.3.625 CrossRef Search ADS   Zane S., Turolla R., Nobili L., Erna M., 1996, ApJ , 466, 871 https://doi.org/10.1086/177561 CrossRef Search ADS   © 2018 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Monthly Notices of the Royal Astronomical Society Oxford University Press

, Volume 476 (2) – May 1, 2018
9 pages

Publisher
The Royal Astronomical Society
ISSN
0035-8711
eISSN
1365-2966
D.O.I.
10.1093/mnras/sty358
Publisher site
See Article on Publisher Site

### Abstract

Abstract We numerically investigate radiatively driven relativistic spherical winds from the central luminous object with mass M and luminosity L* under Newtonian gravity, special relativity, and relativistic radiative transfer. We solve both the relativistic radiative transfer equation and the relativistic hydrodynamical equations for spherically symmetric flows under the double-iteration processes, to obtain the intensity and velocity fields simultaneously. We found that the momentum-driven winds with scattering are quickly accelerated near the central object to reach the terminal speed. The results of numerical solutions are roughly fitted by a relation of $$\dot{m}=0.7(\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2.6}$$, where $$\dot{m}$$ is the mass-loss rate normalized by the critical one, Γ* the central luminosity normalized by the critical one, τ* the typical optical depth, β* the initial flow speed at the central core of radius R*, and βout the terminal speed normalized by the speed of light. This relation is close to the non-relativistic analytical solution, $$\dot{m} = 2(\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2}$$, which can be re-expressed as $$\beta _{\rm out}^2/2 = (\Gamma _*-1)GM/c^2 R_*$$. That is, the present solution with small optical depth is similar to that of the radiatively driven free outflow. Furthermore, we found that the normalized luminosity (Eddington parameter) must be larger than unity for the relativistic spherical wind to blow off with intermediate or small optical depth, i.e. $$\Gamma _* \gtrsim \sqrt{(1+\beta _{\rm out})^3/(1-\beta _{\rm out})}.$$ We briefly investigate and discuss an isothermal wind. accretion, accretion discs, radiative transfer, relativistic processes, stars: winds, outflows, galaxies: active 1 INTRODUCTION Active astrophysical objects often exhibit outflow phenomena, which are supposed to be driven by thermal, magnetic, and radiative forces. Accretion discs usually play dominant roles in these active objects (Kato, Fukue & Mineshige 2008). Outflows from relativistic objects, such as neutron stars or black holes, have relativistic speeds up to 0.1 c in neutron star winds and up to ∼c in black hole winds. These relativistic outflows are believed to be driven by radiation pressure, when the central object is sufficiently luminous. Furthermore, except for the very vicinity of the central object, the winds are treated as spherically symmetric outflows, even if they emanate from a black hole accretion disc system, in a simple picture (King & Pounds 2003; King 2010; King & Muldrew 2016; see also Poutanen et al. 2007; Middleton et al. 2015). In order to examine optically thick radiatively driven relativistic spherical winds, researchers usually solved the relativistic hydrodynamical equations and relativistic moment ones with closure relations, instead of the relativistic radiative transfer equation (e.g. Ruggles & Bath 1979; Paczyński & Prószyński 1986; Paczyński 1986, 1990; Turolla & Nobili & Calvani 1986; Nobili, Turolla & Lapidus 1994; Akizuki & Fukue 2008, 2009). As already pointed out in Thorne, Flammang & Żytkow (1981), however, the diffusion approximation has a causality problem, whereas the Eddington approximation is also invalid in the relativistic regime (see e.g. Turolla & Nobili 1988; Nobili, Turolla & Zampieri 1991; Turolla, Zampieri & Nobili 1995; Dullemond 1999; Fukue 2005). Some researchers solved not moment equations, but the radiative transfer equation in spherically symmetric outflows in the flat space–time (e.g. Castor 1972; Mihalas, Kunasz & Hummber 1975, 1976a,b; Schinder & Bludman 1989; Zane et al. 1996; Baron & Hauschildt 2004; Chen et al. 2007; Knop, Hauschildt & Baron 2007), and in the curved one (Schinder & Bludman 1989; Zane et al. 1996; Chen et al. 2007; Knop et al. 2007 ). In these researches, however, the velocity laws were usually assumed, and the hydrodynamical equations were not solved. In order to solve the radiatively driven relativistic winds, we should treat both the relativistic hydrodynamical equations and relativistic radiative transfer equation simultaneously. Recently, under the double-iteration processes, both equations were simultaneously solved for the plane-parallel relativistic flows (Fukue 2015; cf. Fukue 2014), and for the spherical relativistic ones (Fukue 2017b; cf. Fukue 2017a). In these studies, particularly in Fukue (2017b) for a spherical case, however, the gravitational field of the central object was ignored for simplicity. Thus, in this paper, for the spherically symmetric case, we simultaneously solve both the relativistic radiative transfer equation and the relativistic hydrodynamical one for spherically symmetric winds under the gravitational effect of the central object. In the next section, we present the basic equations for radiation and gas. In Section 3, we show and discuss numerical solutions of scattering winds, whereas we briefly examine isothermal winds in Section 4. The final section is devoted to concluding remarks. 2 BASIC EQUATIONS Let us consider a steady spherically symmetric wind, which is driven by radiation force in the radial (r) direction from a central luminous object with mass M, radius R*, and luminosity L* (Fig. 1). The central object emits an isotropic uniform intensity I*. Hence, $$L_* = 4\pi R_*^2 \pi I_*$$, while the Eddington luminosity LE of the central object is LE ≡ 4πcGM/(κ0 + σ0), and the normalized central luminosity (Eddington parameter Γ*) is defined as Γ* ≡ L*/LE, where κ0 and σ0 are the absorption and scattering opacities measured in the comoving frame, respectively. Here and in what follows, the subscript 0 denotes the quantities in the comoving frame. Figure 1. View largeDownload slide Relativistic radiatively driven spherical wind in the radial (r) direction from a central luminous object with mass M, radius R*, luminosity L*. Both the spherical coordinates (r, θ, φ) and the impact space coordinates (p, z) are used. The flow velocity v and density ρ are functions of radius r, whereas the radiative intensity I±(p, z) is a function of (p, z); r2 = p2 + z2 and cos θ = μ = z/r. Figure 1. View largeDownload slide Relativistic radiatively driven spherical wind in the radial (r) direction from a central luminous object with mass M, radius R*, luminosity L*. Both the spherical coordinates (r, θ, φ) and the impact space coordinates (p, z) are used. The flow velocity v and density ρ are functions of radius r, whereas the radiative intensity I±(p, z) is a function of (p, z); r2 = p2 + z2 and cos θ = μ = z/r. We use the spherical coordinates (r, θ, φ) as well as the impact space coordinates (p, z), where the z-axis is in the line-of-sight direction of a distant observer (Hummer & Rybicki 1971). The transformation between these coordinates is as follows: r2 = p2 + z2 and cos θ = μ = z/r, where cos θ is the direction cosine from the z-axis. For the spherically symmetric flow, the flow velocity v and density ρ0 are functions of radius r, while the radiative intensity I±(p, z) in the inertial frame is a function of (p, z). As for the numerical processes, the outer radius of the calculation region is Rout; numerically, Rout = 100R*. The tangent rays are set in the logarithmic scale at pin ≤ p ≤ pout; numerically, pin = 0.01R* and pout = 100R*. Except for the gravity term, the basic equations are almost same as those in Fukue (2017b). For the convenience of readers, we briefly repeat them in the followings. 2.1 Relativistic radiative transfer equation The radiative transfer equations in various forms are derived and described in various works of literature (Chandrasekhar 1960; Mihalas 1970; Rybicki & Lightman 1979; Mihalas & Mihalas 1984; Kato, Fukue & Mineshige 1998; Mihalas & Auer 2001; Peraiah 2002; Castor 2004; Kato et al. 2008). The relativistic radiative transfer equation is also found in several works of literature (Mihalas 1980; Kato et al. 1998, 2008; Mihalas & Auer 2001). In the present treatment, the relativistic transfer equation in the mixed frame is   \begin{eqnarray} \pm \frac{\mathrm{\partial} I^{\pm }(p,z)}{\mathrm{\partial} z} &=& \frac{\rho _0}{\gamma ^3 (1-\beta \cdot l)^3} \left[ \frac{j_0}{4\pi } -\left(\kappa _0+\sigma _0\right) I_0 +\sigma _0 J_0 \right] \nonumber \\ &=& -\frac{(\kappa _0+\sigma _0)\rho _0}{\gamma ^3 (1-\beta \cdot l)^3} \left[ \gamma ^4 (1-\boldsymbol\beta \cdot \boldsymbol l)^4 I^{\pm } - S_0 \right], \end{eqnarray} (1)where I±(p, z) are the intensities in the upward and downward directions in the inertial (fixed) frame, respectively, and I0(p, z) and J0(r) the frequency-integrated specific intensity and mean intensity in the comoving (fluid) frames, respectively. In addition, β(r) (=v/c) is the normalized flow speed, $$\gamma =1/\sqrt{1-\beta ^2}$$ the Lorentz factor, and j0(r) the frequency-integrated mass emissivity. Furthermore, ($$\beta \cdot l$$) (=±β cos θ = ±βz/r) is the scalar product between the flow velocity $$\beta$$ (radial direction) and the ray direction $$l$$ (vertical direction) in the inertial frame. Finally, the source function S0(r) is defined in the comoving frame as   $$S_0(r) \equiv \frac{j_0}{4\pi }\frac{1}{\kappa _0+\sigma _0} + \frac{\sigma _0}{\kappa _0+\sigma _0}J_0.$$ (2)In the local thermodynamic equilibrium (LTE) case, j0/(4π) = κ0B0, where B0 ($$=\sigma T_0^4/\pi$$) is the blackbody intensity, T0 being the temperature in the comoving frame. It should be noted that the isotropic scattering is assumed for simplicity. The transformation rules for the radiative quantities are summarized as follows:   \begin{eqnarray} I_0 = \gamma ^4 (1-\beta \cdot l)^4 I, \end{eqnarray} (3)  \begin{eqnarray} J_0 = \gamma ^2 \left(J - 2\beta H + \beta ^2\,K \right), \end{eqnarray} (4)  \begin{eqnarray} H_0 = \gamma ^2 \left[ \left(1 + \beta ^2 \right)H - \beta \left(J + K \right) \right], \end{eqnarray} (5)  \begin{eqnarray} K_0 = \gamma ^2 \left( \beta ^2 J - 2\beta H + K \right), \end{eqnarray} (6)and the Eddington factor f is defined in the comoving frame as   $$f \equiv \frac{K_0}{J_0}.$$ (7) As was derived in Fukue (2017a), the relativistic formal solutions of equation (1) are   \begin{eqnarray} I^+(p,z) &=& {\rm e}^{ \displaystyle G(p,z_*)-G(p,z) - U(p,z_*)+U(p,z)} I^*(p,z_*) \nonumber \\ && + \int ^{z}_{0} \frac{ {\rm e}^{ \displaystyle G(p,\zeta ) - G(p,z) - U(p,\zeta ) + U(p,z)} }{\displaystyle \gamma ^3 \left(1 - \beta \frac{\zeta }{r} \right)^3} \nonumber \\ && \times\, (\kappa _0+\sigma _0)\rho _0 S_0 {\rm d}\zeta , \end{eqnarray} (8)  \begin{eqnarray} I^-(p,z) &=& - \int ^{z}_{z_{\rm out}} \frac{ {\rm e}^{ \displaystyle [G(p,z)-G(p,\zeta )] + [U(p,z)-U(p,\zeta )] }}{\displaystyle \gamma ^3 \left(1 + \beta \frac{\zeta }{r} \right)^3} \nonumber \\ && \times \, (\kappa _0+\sigma _0)\rho _0 S_0 {\rm d}\zeta , \end{eqnarray} (9)where $$z_* = \sqrt{R_*^2-p^2}$$, $$z_{\rm out} =\sqrt{R_{\rm out}^2-p^2}$$, and   \begin{eqnarray} G(p,z) \equiv \int ^z (\kappa _0+\sigma _0)\rho _0 \gamma {\rm d}\zeta , \end{eqnarray} (10)  \begin{eqnarray} U(p,z) \equiv \int ^z (\kappa _0+\sigma _0)\rho _0 \gamma (\beta \cdot l) {\rm d}\zeta . \end{eqnarray} (11) The boundary condition for the downward intensity is   $$I^-(p,z_{\rm out}) = 0 {\rm \ \ \ at\ \ \, } z=z_{\rm out}=\sqrt{R_{\rm out}^2-p^2},$$ (12)whereas that for the upward intensity is   \begin{eqnarray} I^+(p,z_*) &=& I^* {\rm \ \ \ for\ \ \, } p < R_*, \end{eqnarray} (13)  \begin{eqnarray} I^+(p,0) &=& I^-(p,0) {\rm \ \ \ for\ \ \, } p > R_*. \end{eqnarray} (14) Using the radius R*, opacity value κ*(R*), and density ρ*(R*) there, equations (10) and (11) are normalized as   \begin{eqnarray} G(p,z) \equiv \tau _* \int ^z \frac{\kappa _0 + \sigma _0}{\kappa _*} \frac{\rho _0}{\rho _*} \gamma {\rm d}\tilde{\zeta }, \end{eqnarray} (15)  \begin{eqnarray} U(p,z) \equiv \tau _* \int ^z \frac{\kappa _0 + \sigma _0}{\kappa _*} \frac{\rho _0}{\rho _*} \gamma \beta \cdot l{\rm d}\tilde{\zeta }, \end{eqnarray} (16)where $$\tilde{\zeta }=\zeta /R_*$$ and   $$\tau _* = \kappa _* \rho _* R_*,$$ (17)which is a typical optical depth of the wind. For a given velocity profile v and source function S0, the formal solutions (8) and (9) are numerically integrated under the appropriate boundary conditions and parameters, and the refined radiative quantities are obtained. 2.2 Relativistic hydrodynamical equation The relativistic hydrodynamical equations are also found in several works of literature (e.g. Kato et al. 1998, 2008). For the steady spherical wind under the special relativity, the continuity equation becomes   $$4 \pi r^2 \rho _0 \gamma v = 4 \pi r^2 \rho _0 c \gamma \beta = \dot{M},$$ (18)where $$\dot{M}$$ is the mass-loss rate and assumed to be constant. The relativistic equation of motion is   $$c^2 u \frac{{\rm d}u}{{\rm d}r} = c^2 \gamma ^4 \beta \frac{{\rm d}\beta }{{\rm d}r} = -\frac{GM}{r^2} + \frac{\kappa _0+\sigma _0}{c}\gamma 4\pi H_0,$$ (19)where u (=γβ) is the four velocity, and the gas pressure-gradient force is ignored in the present case. We do not use the energy equation for matter. Using the quantities on the central object, radius R*, density ρ*, velocity β* (corresponding Lorentz factor γ*), and luminosity L*, hydrodynamical equations (18) and (19) are respectively normalized as follows:   \begin{eqnarray} \hat{r}^2 \frac{\rho _0}{\rho _*}\frac{\gamma \beta }{\gamma _* \beta _*} &=& 1, \end{eqnarray} (20)  \begin{eqnarray} \dot{m} \gamma ^3 \frac{{\rm d}\beta }{{\rm d}\hat{r}} = \tau _* \frac{\kappa _0 + \sigma _0}{\kappa _*} \frac{\rho _0}{\rho _*} \left( -1 + \Gamma _* \gamma \frac{4\pi H_0}{\pi I_*} \hat{r}^2 \right), \end{eqnarray} (21)where $$\hat{r}=r/R_*$$, and $$\dot{m}$$ ($$\equiv \dot{M}c^2/L_{\rm E}$$) the normalized mass-loss rate, i.e. the critical mass-loss rate is defined as $$\dot{M}_{\rm crit} \equiv L_{\rm E}/c^2$$ in this work. For a given flux profile H0 in the comoving frame, these hydrodynamical equations are numerically integrated under the appropriate boundary conditions and parameters, and the refined velocity profile is obtained. 2.3 Relativistic radiative moment equations In order to analyse the numerical results, we here present the relativistic radiative moment equations (see e.g. Kato et al. 1998, 2008), although we do not solve them in the present numerical calculations. Using the radiative moment quantities in the comoving and inertial frames, the zeroth and first moment equations for the spherically symmetric case respectively become   \begin{eqnarray} \frac{1}{r^2}\frac{{\rm d}}{{\rm d}r}(r^2H) &=& \rho _0 \gamma \left( \frac{j_0}{4\pi } - \kappa _0 J_0 \right) -\rho _0 \gamma (\kappa _0 + \sigma _0)\beta H_0, \nonumber \\ &=& \rho _0 \gamma \left\lbrace\vphantom{\gamma ^2(1+\beta ^2)} \kappa _0 (B_0 - J + \beta H) \right. \nonumber \\ && \left. -\, \sigma _0 \gamma ^2 \beta [(1+\beta ^2)H-\beta (J+K)] \right\rbrace , \end{eqnarray} (22)  \begin{eqnarray} \frac{{\rm d}K}{{\rm d}r} + \frac{3K-J}{r} &=& -\rho _0 \gamma (\kappa _0 + \sigma _0)H_0 +\rho _0 \gamma \beta \left( \frac{j_0}{4\pi } - \kappa _0 J_0 \right) \nonumber \\ &=& -\rho _0 \gamma \left\lbrace\vphantom{(1+\beta ^2)} \kappa _0 (H - \beta B_0 - \beta K) \right. \nonumber \\ && \left. +\, \sigma _0 \gamma ^2 [(1+\beta ^2)H-\beta (J+K)] \right\rbrace , \end{eqnarray} (23)where we use j0/(4π) = κ0B0 in the second lines of both equations. In addition, for a closure relation K0/J0 = f in equation (7), we have a relation   $$(1-f\beta ^2)K = (f-\beta ^2)J + (1-f)2\beta H,$$ (24)in the inertial frame. It should be stressed that if we use this closure relation (24) in the inertial frame, moment equations in the inertial frame have a pathological singularity at β2 = f = 1/3 under the usual Eddington approximation (cf. Fukue 2005). 3 MOMENTUM-DRIVEN SCATTERING WINDS We first examine the purely scattering case with no emission and absorption (j0 = 0 and κ0 = 0), i.e. S0 = J0. In this case, we do not consider any thermal process, and the wind flow is momentum driven. 3.1 Double-iteration processes Similar to the previous studies (Fukue 2015, 2017b), we need double-iteration processes: one is the iteration for radiation field, and the other is that between radiation and gas. We here outline the procedure. Cycle 0 (Initial trial): In the optically thin non-relativistic case, the radiative moments are roughly expressed as J ∝ r−3 near the centre and H ∝ r−2. By considering non-relativistic moment equations for radiation field, we adopt the following analytical solutions:   \begin{eqnarray} J = \tau _* I_* \frac{1}{4\hat{r}^3}, \end{eqnarray} (25)  \begin{eqnarray} H = I_* \frac{1}{4\hat{r}^2}, \end{eqnarray} (26)as initial trial functions for radiation field. For velocity field, we adopt an approximate solution of equation (21) in the non-relativistic limit (see later),   $$\beta = \beta _{\rm out}\sqrt{1-\frac{1}{\hat{r}}},$$ (27)as an initial trial function, where βout is the terminal speed of the flow.  For this initial velocity field and the resultant density distribution, the functions (15) and (16) are tabulated as initial conditions. Cycle 1A (Radiation field):Using initial trial functions (25) and (27) with tabulated functions (15) and (16), we iteratively integrate the relativistic formal solutions (8) and (9) for given parameters, until the solutions of the specific intensities I± converge. Using the resultant intensities and transformation to the comoving frame, we obtain the comoving flux H0(r). Cycle 1B (Velocity field): Using the Eddington flux H0 obtained by Cycle 1A, we solve equation (21) with equation (20) for a given terminal speed βout. As a result, we obtain the refined velocity distribution, and the mass-loss rate $$\dot{m}$$ as an eigenvalue. For this renewed velocity field and the resultant density distribution, the functions (15) and (16) are re-tabulated for the next step. Cycles 2 to L:We then repeat Cycles 1A and 1B, until both the specific intensity and velocity distributions converge. The parameters for this momentum-driven wind are the Eddington parameter Γ*, the typical optical depth τ*, the initial flow speed β* at the central core, the ‘terminal’ speed βout at the outer radius, and the normalized mass-loss rate $$\dot{m}$$. Of these, we set β* = 0.01, as a sufficiently small value. We further set τ* = 3 for a typical case. In addition, the mass-loss rate is determined as an eigenvalue. Hence, we examine the present problem for several values of Γ* and βout. 3.2 Typical solutions The present results are somewhat similar to the previous study without gravity (Fukue 2017b), where the main parameters are τ* and βout. Hence, we consult the previous study on the optical depth dependence, and in this subsection we mainly show the dependence on Γ* and βout. It should be noted that we can obtain the solutions in the limited range of βout for some given Γ*, which is the important difference from the previous study and discussed later in this section. For the case of τ* = 3, examples of numerical solutions are shown in Figs 2 –4. In Fig. 2, the velocity profiles are shown as a function of the normalized radius $$\hat{r}$$ for several values of βout in the case of Γ* = 5 and τ* = 3, whereas the initial trial velocity profiles are plotted by dashed curves for comparison. As shown in Fig. 2, and similar to the previous study (Fukue 2017b), the momentum-driven wind is quickly accelerated near the central object to reach the terminal speed. Figure 2. View largeDownload slide Velocity profiles as a function of radius for several values of βout. The initial trial velocity profiles are plotted by dashed curves for comparison. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5 from bottom to top. Figure 2. View largeDownload slide Velocity profiles as a function of radius for several values of βout. The initial trial velocity profiles are plotted by dashed curves for comparison. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5 from bottom to top. In order to understand this behaviour, similar to the previous study (Fukue 2017b), we again consider the approximate analytical solution of equation (21) in the non-relativistic limit. That is, using equation (20), equation (21) is approximately expressed as   $$\dot{m}\frac{{\rm d}\beta }{{\rm d}\hat{r}} = \frac{\tau _* \beta _*}{\hat{r}^2 \beta } \left( \Gamma _* \frac{4\pi H_0}{\pi I_*}\hat{r}^2 - 1 \right) \sim (\Gamma _* - 1) \frac{\tau _* \beta _*}{\beta } \frac{1}{\hat{r}^2},$$ (28)where we use the relation for the flux in the non-relativistic limit: $$4\pi H_0 4\pi r^2=\pi I_* 4\pi R_*^2 =L_*$$. This equation (28) is easily integrated to yield the velocity profile in the non-relativistic limit,   \begin{eqnarray} \frac{1}{2}\beta ^2 &=& \frac{1}{2}\beta _{\rm out}^2 - (\Gamma _* - 1) \frac{\tau _* \beta _*}{\dot{m}} \frac{1}{\hat{r}}, \end{eqnarray} (29)  \begin{eqnarray} \qquad &=& \frac{1}{2}\beta _{\rm out}^2 - (\Gamma _* - 1) \frac{GM}{c^2r}, \end{eqnarray} (30)where we impose the boundary condition at the outer radius: β = βout at $$\hat{r} = \infty$$. It should be noted that this analytical solution in the non-relativistic limit is just a simple radiatively driven free outflow. We further impose the boundary condition at the inner radius: β = β* ≪ βout at $$\hat{r}=1$$, and we obtain the relation   $$\dot{m} = \frac{2(\Gamma _* - 1)\tau _* \beta _*}{\beta _{\rm out}^2}.$$ (31)Inserting this relation into the approximate solution, we finally obtain a simple analytical solution in the non-relativistic limit:   $$\beta = \beta _{\rm out}\sqrt{1-\frac{1}{\hat{r}}}.$$ (32) For small βout, this analytical solution (32) well reproduces the numerical solutions in Fig. 2. That is, the present solution is quite similar to the radiatively driven outflow in the non-relativistic limit. For large βout, however, it slightly deviates from the numerical ones. We emphasize that this analytical solution does not include any parameters except for βout. In Fig. 3, the radiative quantities normalized by I* are shown as a function of r normalized by R* for several values of βout (=0.1, 0.3, 0.5) in the case of Γ* = 5 and τ* = 3. Overall behaviour of radiative quantities is similar to the previous studies (Fukue 2017a,b). Figure 3. View largeDownload slide Radiative quantities as a function of radius for several values of βout. (a) Normalized mean intensity, (b) normalized flux, and (c) the Eddington factor. The solid curves represent the quantities in the inertial frame, while the dashed ones mean those in the comoving frame. Radius is normalized by R*, while radiative quantities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5. Figure 3. View largeDownload slide Radiative quantities as a function of radius for several values of βout. (a) Normalized mean intensity, (b) normalized flux, and (c) the Eddington factor. The solid curves represent the quantities in the inertial frame, while the dashed ones mean those in the comoving frame. Radius is normalized by R*, while radiative quantities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5. In Fig. 3(a), the mean intensity J is plotted as a function of radius $$\hat{r}$$. Similar to the previous studies (Fukue 2017a,b), the mean intensity in the inertial frame (solid curves) roughly decreases as J ∝ r−2, whereas that in the comoving frame (dashed curves) reduces in the entire region, due to redshift (Fukue 2014, 2017a). In Fig. 3(b), the Eddington flux H is plotted as a function of radius $$\hat{r}$$. Similar to the previous studies (Fukue 2017a,b), the flux in the inertial frame (solid curves) roughly decreases as H ∝ r−2, whereas that in the comoving frame (dashed curves) reduces in the entire region, as the flow velocity increases (Fukue 2017a). In Fig. 3(c), the Eddington factor f is plotted as a function of radius $$\hat{r}$$. Similar to the previous study without gravity (Fukue 2017b), the Eddington factor is almost 1/3 near the core, and increases towards unity as radius becomes large, due to the geometrical and optical depth effects. In Fig. 4, the emergent intensities I+(p, zout) at the outer boundary Rout normalized by I* are plotted as a function of the impact coordinate p normalized by R*. The normalized luminosity Γ* = 5, the typical optical depth τ* is 3, and the values of βout are 0.1, 0.3, 0.5 from outer to inner. Figure 4. View largeDownload slide Emergent intensities in the inertial frame at the outer radius Rout. Radius is normalized by R*, while the intensities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5 from outer to inner. Figure 4. View largeDownload slide Emergent intensities in the inertial frame at the outer radius Rout. Radius is normalized by R*, while the intensities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, and the values of βout are 0.1, 0.3, 0.5 from outer to inner. Since the optical depth is not so large (τ* = 3), the emergent intensity in the core part of p/R* ≤ 1 is almost unity; this part is roughly the incident intensity I*. In the outer envelope of p/R* > 1, on the other hand, the emergent intensity decreases with p due to the scattering effect, and it is less enhanced as the flow speed increases due to the transverse Doppler effect. Overall behaviour of emergent intensities is quite similar to the previous study (Fukue 2017b). 3.3 Relations on the parameters Since the mass-loss rate $$\dot{m}$$ is obtained as an eigenvalue, similar to the previous study (Fukue 2017b), there exist some relations among the normalized luminosity Γ*, the typical optical depth τ*, the normalized mass-loss rate $$\dot{m}$$, and the normalized flow terminal speed βout. The obtained relations between $$\dot{m}$$ and βout for several values of Γ* in the case of τ* = 3 are plotted in Fig. 5. Figure 5. View largeDownload slide Relations between the normalized mass-loss rate $$\dot{m}$$ and the flow terminal speed βout for several values of the normalized luminosity Γ* in the case of τ* = 3 (solid curves with open circles). The dashed lines denote a fitting relation of $$\dot{m} = 0.7 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2.6}$$ for the case of τ* = 3. The dot–dashed lines represent an analytical relation of $$\dot{m}=2 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2}$$ in the non-relativistic regime. The curve with crosses expresses a rough limit of the existence of solutions. Figure 5. View largeDownload slide Relations between the normalized mass-loss rate $$\dot{m}$$ and the flow terminal speed βout for several values of the normalized luminosity Γ* in the case of τ* = 3 (solid curves with open circles). The dashed lines denote a fitting relation of $$\dot{m} = 0.7 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2.6}$$ for the case of τ* = 3. The dot–dashed lines represent an analytical relation of $$\dot{m}=2 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2}$$ in the non-relativistic regime. The curve with crosses expresses a rough limit of the existence of solutions. As is seen in Fig. 5, similar to the plane-parallel case (Fukue 2015) and to the previous spherical one (Fukue 2017b), the mass-loss rate decreases as the flow terminal speed increases, while it increases with the normalized luminosity or with the typical optical depth. Similar to the previous study (Fukue 2017b), these trends are approximated by a fitting relation to the numerical results, and also by an analytical solution, using the back-of-the-envelope calculation. For example, the dashed lines in Fig. 5 denote a fitting relation of $$\dot{m} = 0.7 (\Gamma _*-1) \tau _* \beta _* \beta _{\rm out}^{-2.6}$$ for the case of τ* = 3. On the other hand, the dot–dashed lines in Fig. 5 represent the analytical solution (31) of $$\dot{m}=2(\Gamma _* - 1)\tau _* \beta _* \beta _{\rm out}^{-2}$$ for the case of τ* = 3, which is obtained during the derivation of the approximate solution in the non-relativistic limit. It should be noted that the normalized mass-loss rate $$\dot{m}$$ apparently depends on parameters under the present normalization, but the actual mass-loss rate is not determined in the present way, where the velocity field is calculated. Now, we shall discuss the limited range of solutions. As was already stated, and as seen in Figs 2–5, for given Γ* and τ*, the numerical solutions can be obtained only for the lower values of βout, and it is hard to obtain solutions when the value of βout becomes high. This nature is remarkably different from the results in the previous study without gravity (Fukue 2017b). This behaviour is mainly due to the physical reason, and partly due to the numerical technique. Physically, for large values of the flow speed, the net force on the right-hand side of equation (21) becomes negative even for Γ* > 1, since the comoving flux decreases as the flow speed increases (Fig. 3b). This will be analytically shown below. Numerically, as the flow speed increases, even in the case that the numerical solution may exist, the net force temporary becomes negative during the iteration process, and the solution search is stopped. This limited range is approximately derived as follows. As is seen in Figs 2 and 3, at far from the centre, the flow becomes constant, and the Eddington factor is almost unity due to the geometrical and optical depth effects. In this region of f = 1, K0 = J0 and further K = J. Moreover, the radiative moments in the inertial frame roughly become   \begin{eqnarray} \frac{4\pi H}{\pi I_*} \hat{r}^2 &=& 1, \end{eqnarray} (33)  \begin{eqnarray} \frac{4\pi J}{\pi I_*} \hat{r}^2 &=& 1. \end{eqnarray} (34)Hence, the Eddington flux in the comoving frame becomes   \begin{eqnarray} \frac{4\pi H_0}{\pi I_*} \hat{r}^2 &=& \frac{4\pi }{\pi I_*} \hat{r}^2 \gamma ^2[(1+\beta ^2)H-\beta (J+K)] \nonumber \\ &=& \frac{1-\beta }{1+\beta } \ \ \ \le 1. \end{eqnarray} (35)Inserting this approximate solution into the net force on the right-hand side of equation (21), we have a condition that the net force vanishes:   $$\Gamma _* = \sqrt{\frac{(1+\beta )^3}{1-\beta }} \ \ \ \ge 1.$$ (36)Thus, the normalized luminosity (Eddington parameter) must be larger than unity for the relativistic spherical wind with intermediate or small optical depth. Finally, inserting this condition (36) into the fitting relation of $$\dot{m}=0.7 (\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2.6}$$, we have the limiting condition on the $$\dot{m}$$–βout plane, which is plotted in Fig. 5 (cross marks). As is seen in Fig. 5, this approximate condition well explains the existence domain of the present numerical solutions. This limitation is interpreted as follows. We do not consider the energy equation for matter, and therefore, the kinetic energy flux of matter is ignored compared with the energy flux in the luminosity. However, in the gaseous wind, if the mass-loss rate becomes large and the wind is heavy, the kinetic energy flux of matter cannot be negligible, and the radiation field with a given luminosity cannot accelerate matter to sufficient speed at infinity. In other words, in order to accelerate matter to sufficient speed at infinity, the sufficient luminosity is required, as in equation (35). 4 ISOTHERMAL WINDS We next consider the isothermal case without scattering (σ0 = 0) under LTE (j0/4π = κ0B0), i.e. S0 = B0 = constant. In this case, the wind flow is also momentum driven, although the momentum source is distributed in the entire space. Since the iteration process for radiation field is unnecessary in this isothermal case, we need a single iteration process between radiation and gas. The parameters for this isothermal wind are Γ*, τ*, β*, βout, $$\dot{m}$$, and the value of B0. Of these, we set β* = 0.01 and τ* = 3 for a typical case. The value of B0 is set as B0 = I*. The mass-loss rate is determined as an eigenvalue. Hence, the rests are Γ* and βout. 4.1 Typical solutions Examples of numerical solutions for isothermal winds are shown in Figs 6 –8. Figure 6. View largeDownload slide Velocity profiles as a function of radius for several values of βout in the isothermal case. The initial trial velocity profiles are plotted by dashed curves for comparison. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from bottom to top. Figure 6. View largeDownload slide Velocity profiles as a function of radius for several values of βout in the isothermal case. The initial trial velocity profiles are plotted by dashed curves for comparison. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from bottom to top. Figure 7. View largeDownload slide Radiative quantities as a function of radius for several values of βout in the isothermal case. (a) Normalized mean intensity, (b) normalized flux, and (c) the Eddington factor. The solid curves represent the quantities in the inertial frame, while the dashed ones mean those in the comoving frame. Radius is normalized by R*, while radiative quantities are normalized by I*. The normalized luminosity Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5. Figure 7. View largeDownload slide Radiative quantities as a function of radius for several values of βout in the isothermal case. (a) Normalized mean intensity, (b) normalized flux, and (c) the Eddington factor. The solid curves represent the quantities in the inertial frame, while the dashed ones mean those in the comoving frame. Radius is normalized by R*, while radiative quantities are normalized by I*. The normalized luminosity Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5. Figure 8. View largeDownload slide Emergent intensities in the inertial frame at the outer radius Rout in the isothermal case. Radius is normalized by R*, while the intensities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from top to bottom. Figure 8. View largeDownload slide Emergent intensities in the inertial frame at the outer radius Rout in the isothermal case. Radius is normalized by R*, while the intensities are normalized by I*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from top to bottom. In Fig. 6, the velocity profiles are shown as a function of the normalized radius $$\hat{r}$$ for several values of βout in the case of Γ* = 5, τ* = 3, and B0/I* = 1, whereas the initial trial velocity profiles are plotted by dashed curves for comparison. In contrast to the scattering case in Fig. 2, the isothermal wind is gradually accelerated near the central object to reach the terminal speed. This is because the momentum source is uniformly distributed in the entire space. In Fig. 7, the radiative quantities normalized by I* are shown as a function of r normalized by R* for several values of βout (=0.1, 0.3, 0.5) in the case of Γ* = 5, τ* = 3, and B0/I* = 1. In Fig. 7(a), the mean intensity J is plotted as a function of radius $$\hat{r}$$. In contrast to the scattering case (Fig. 3a), the mean intensity in the inertial frame (solid curves) and that in the comoving one (dashed ones) roughly decrease as J ∝ r−1. This behaviour originates from the difference in source distributions. In the scattering wind (Fig. 3a), the radiation source exists only on the central core, whereas the radiation source is uniformly distributed in the isothermal wind (Fig. 7a). Moreover, the mean intensity decreases as the flow speed increases, since the density and the source intensity (emissivity) decrease as the flow speed increases. In Fig. 7(b), the Eddington flux H is plotted as a function of radius $$\hat{r}$$. Due to the same reason, the flux also decreases as H ∝ r−1. In Fig. 7(c), the Eddington factor f is plotted as a function of radius $$\hat{r}$$. In contrast to the scattering case (Fig. 3c), the Eddington factor in the present isothermal wind does not approach unity, but becomes around 1/2 or so. This means that the radiation field of the isothermal wind is between the sufficiently thick isotropic case (f = 1/3) and the intermediately thin spherical case (f ∼ 1). In order to confirm these numerical solutions for the radiative moment quantities, let us examine the relativistic moment equations. When the flow speed is small and of the order of $${\cal O}(\beta ^1)$$, the zeroth moment equation (22), the first one (23), and the closure relation (24) respectively become in the isothermal case without scattering   \begin{eqnarray} \frac{1}{\hat{r}^2}\frac{{\rm d}}{{\rm d}\hat{r}}(\hat{r}^2H) &=& \tau _* \frac{\rho _0}{\rho _*} (B_0 - J + \beta H), \end{eqnarray} (37)  \begin{eqnarray} \frac{{\rm d}K}{{\rm d}\hat{r}} + \frac{3K-J}{\hat{r}} &=& -\tau _* \frac{\rho _0}{\rho _*} (H - \beta K - \beta B_0), \end{eqnarray} (38)  \begin{eqnarray} K &=& fJ + (1-f)2\beta H. \end{eqnarray} (39) At far from the centre, J and H are smaller than B0. Then, using continuity equation (20), equation (37) becomes   $$\frac{1}{\hat{r}^2}\frac{{\rm d}}{{\rm d}\hat{r}}(\hat{r}^2H) \sim \tau _* \frac{\beta _*}{\hat{r}^2 \beta } B_0.$$ (40)At far from the centre, the flow speed is almost constant, and this equation is solved to yield a solution   $$H \sim \frac{\tau _* \beta _*}{\beta }\frac{B_0}{\hat{r}},$$ (41)and the flux roughly varies as $$H \propto \hat{r}^{-1}$$. Similarly, at far from the centre, where B0 is dominated and the flow speed is constant, equation (38) with equation (39) is also solved to yield an approximate solution:   \begin{eqnarray} J &\sim & \tau _* \beta _* \frac{B_0}{\hat{r}} \ \ \ \ \ \ \ \ \ \ \ {\rm for}\ \ \ f=1, \end{eqnarray} (42)  \begin{eqnarray} J &\sim & 6\tau _* \beta _* \frac{B_0}{\hat{r}}\ln \hat{r} \ \ \ {\rm for}\ \ \ f=1/2, \end{eqnarray} (43)and the mean intensity roughly varies as $$J \propto \hat{r}^{-1}$$ at least on the main dependence. In Fig. 8, the emergent intensities I+(p, zout) at the outer boundary Rout normalized by I* are plotted as a function of the impact coordinate p normalized by R*. The normalized luminosity is Γ* = 5, the typical optical depth τ* = 3, the blackbody intensity B0/I* = 1, and the values of βout are 0.1, 0.3, 0.5 from top to bottom. In contrast to the scattering case (Fig. 4), the emergent intensity in the outer envelope is quite large. This is also because the source emissivity is uniformly distributed in the entire space in the isothermal wind. Similar to the relations obtained in the scattering case (Fig. 5), in this isothermal wind we can also the relations on the parameters, although we do not discuss them for simplicity. Finally, we should mention the opacity effect. In this isothermal wind, we assume that the absorption opacity κ0 is constant. In the realistic case, the absorption opacity depends on the physical quantities, e.g. the free–free opacity is proportional to $$\rho _0 T_0^{-7/2}$$. Hence, it is   $$\frac{\kappa _0}{\kappa _*} = \frac{\rho _0}{\rho _*} \frac{T_0}{T_*} = \frac{\rho _0}{\rho _*}$$ (44)in the isothermal case. In this case, the results are quite similar to the scattering case. The reason is because the density is concentrated in the centre and quickly decreases towards the outer envelope, and therefore, the emissivity is also concentrated in the central core. 5 CONCLUDING REMARKS In this paper, we have examined the radiatively driven relativistic spherical winds from the central object with mass M and luminosity L* under Newtonian gravity and special relativity, by numerically solving the relativistic radiative transfer equation and relativistic hydrodynamical ones for spherical geometry, simultaneously and consistently. Namely, using the initial trial velocity profile, we iteratively solved the relativistic formal solution of the transfer equation to obtain specific intensities and moment quantities in the inertial and comoving frames. Using the comoving flux obtained, we then solved the relativistic hydrodynamical equations to obtain the refined velocity profile, and the mass-loss rate as an eigenvalue. We repeat these double-iteration procedures until both the intensity and velocity converge. We found that the momentum-driven winds with scattering are quickly accelerated near the central luminous core to reach the terminal speed (Fig. 2), similar to the radiative outflows without gravity (Fukue 2017b). The behaviour of moment quantities and emergent intensity is quite similar to the previous results (figs 3 and 4 in Fukue 2017b). The results of numerical solutions are roughly fitted by a relation of $$\dot{m}=0.7(\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2.6}$$, which is well reproduced by the analytical relation in the non-relativistic limit $$\dot{m} = 2(\Gamma _*-1)\tau _* \beta _* \beta _{\rm out}^{-2}$$ (Fig. 5). It should be noted that the normalized mass-loss rate $$\dot{m}$$ apparently depends on parameters under the present normalization, but the actual mass-loss rate is not determined in the present way, where the velocity field is calculated. Furthermore, we found that the normalized luminosity (Eddington parameter) must be larger than unity for the relativistic spherical wind with intermediate or small optical depth: $$\Gamma _* \gtrsim \sqrt{(1+\beta _{\rm out})^3/(1-\beta _{\rm out})}$$ (Fig. 5). We briefly investigate and discuss an isothermal wind. We found that the isothermal winds are gradually accelerated, since the emissivity is uniformly distributed in the entire space. The advantage of the present study is twofold. The first is the self-consistency; in this work, both the radiation and flow fields are obtained simultaneously. The second is the non-singularity; a pathological singularity under the Eddington approximation (e.g. Dullemond 1999; Fukue 2005) does not appear, since we solve the full transfer equation, instead of the truncated moment equations under the closure relation (24). On the other hand, the present solution has several limitations. For example, we do not consider the energy equation for matter, and therefore, the kinetic energy flux of matter is ignored compared with the energy flux in the luminosity. Hence, the present solution implicitly treats the light wind, and cannot treat the heavy wind. Furthermore, the present solution is limited in the intermediate or small optical depth (τ ∼ 10). If the optical depth becomes large, the diffusion time-scale is longer than the advection time-scale, and the present treatment would break down. In this paper, we have considered the gravitational effect, but the gas pressure has not yet been included. If we take into account the gas pressure-gradient force in the equation of motion, we must seek transonic solutions, which is left as a future work. Acknowledgements The author would like to thank an anonymous referee for valuable comments, including the interpretation of solutions. REFERENCES Akizuki C., Fukue J., 2008, PASJ , 60, 337 https://doi.org/10.1093/pasj/60.2.337 CrossRef Search ADS   Akizuki C., Fukue J., 2009, PASJ , 61, 543 https://doi.org/10.1093/pasj/61.3.543 CrossRef Search ADS   Baron E., Hauschildt P. H., 2004, A&A , 427, 987 https://doi.org/10.1051/0004-6361:20040039 CrossRef Search ADS   Castor J. I., 1972, ApJ , 178, 779 https://doi.org/10.1086/151834 CrossRef Search ADS   Castor J. I., 2004, Radiation Hydrodynamics . Cambridge Univ. Press, Cambridge Google Scholar CrossRef Search ADS   Chandrasekhar S., 1960, Radiative Transfer . Dover Press, New York Chen B., Kantowski R., Baron E., Knop S., Hauschildt P. H., 2007, MNRAS , 380, 104 https://doi.org/10.1111/j.1365-2966.2007.11652.x CrossRef Search ADS   Dullemond C. P., 1999, A&A , 343, 1030 Fukue J., 2005, PASJ , 57, 1023 https://doi.org/10.1093/pasj/57.6.1023 CrossRef Search ADS   Fukue J., 2014, PASJ , 66, 73 https://doi.org/10.1093/pasj/psu048 CrossRef Search ADS   Fukue J., 2015, PASJ , 67, 14 Fukue J., 2017a, PASJ , 69, 8 https://doi.org/10.1093/pasj/psw112 CrossRef Search ADS   Fukue J., 2017b, PASJ , 69, 53 Hummer D. G., Rybicki G. B., 1971, MNRAS , 152, 1 https://doi.org/10.1093/mnras/152.1.1 CrossRef Search ADS   Kato S., Fukue J., Mineshige S., 1998, Black-Hole Accretion Disks . Kyoto Univ. Press, Kyoto Kato S., Fukue J., Mineshige S., 2008, Black-Hole Accretion Disks – Towards a New Paradigm . Kyoto Univ. Press, Kyoto King A. R., 2010, MNRAS , 402, 1516 https://doi.org/10.1111/j.1365-2966.2009.16013.x CrossRef Search ADS   King A. R., Muldrew S. I., 2016, MNRAS , 455, 1211 https://doi.org/10.1093/mnras/stv2347 CrossRef Search ADS   King A. R., Pounds K. A., 2003, MNRAS , 345, 657 https://doi.org/10.1046/j.1365-8711.2003.06980.x CrossRef Search ADS   Knop S., Hauschildt P. H., Baron E., 2007, A&A , 463, 315 https://doi.org/10.1051/0004-6361:20065521 CrossRef Search ADS   Middleton M. J., Heil L., Pintore F., Watson D. J., Roberts P., 2015, MNRAS , 447, 3243 https://doi.org/10.1093/mnras/stu2644 CrossRef Search ADS   Mihalas D., 1970, Stellar Atmospheres . Freeman & Co., San Francisco Mihalas D., 1980, ApJ , 237, 574 https://doi.org/10.1086/157902 CrossRef Search ADS   Mihalas D., Auer L. H., 2001, J. Quant. Spectrosc. Radiat. Transfer , 71, 61 https://doi.org/10.1016/S0022-4073(01)00013-9 CrossRef Search ADS   Mihalas D., Mihalas B. W., 1984, Foundations of Radiation Hydrodynamics . Oxford Univ. Press, Oxford Mihalas D., Kunasz P. B., Hummber D. G., 1975, ApJ , 202, 465 https://doi.org/10.1086/153996 CrossRef Search ADS   Mihalas D., Kunasz P. B., Hummber D. G., 1976a, ApJ , 206, 515 https://doi.org/10.1086/154407 CrossRef Search ADS   Mihalas D., Kunasz P. B., Hummber D. G., 1976b, ApJ , 210, 419 https://doi.org/10.1086/154845 CrossRef Search ADS   Nobili L., Turolla R., Zampieri L., 1991, ApJ , 383, 250 https://doi.org/10.1086/170781 CrossRef Search ADS   Nobili L., Turolla R., Lapidus I., 1994, ApJ , 433, 276 https://doi.org/10.1086/174643 CrossRef Search ADS   Paczyński B., 1986, ApJ , 308, L43 https://doi.org/10.1086/184740 CrossRef Search ADS   Paczyński B., 1990, ApJ , 363, 218 https://doi.org/10.1086/169332 CrossRef Search ADS   Paczyński B., Prószyński M., 1986, ApJ , 302, 519 https://doi.org/10.1086/164012 CrossRef Search ADS   Peraiah A., 2002, An Introduction to Radiative Transfer: Methods and Applications in Astrophysics . Cambridge Univ. Press, Cambridge Poutanen J., Lipunova G., Fabrika S., Butkevich A. G., Abolmasov P., 2007, MNRAS , 377, 1187 https://doi.org/10.1111/j.1365-2966.2007.11668.x CrossRef Search ADS   Ruggles C. L. N., Bath G. T., 1979, A&A , 80, 97 Rybicki G. B., Lightman A. P., 1979, Radiative Processes in Astrophysics . Wiley, New York Schinder P. J., Bludman S. A., 1989, ApJ , 346, 350 https://doi.org/10.1086/168015 CrossRef Search ADS   Thorne K. S., Flammang R. A., Żytkow A. N., 1981, MNRAS , 194, 475 https://doi.org/10.1093/mnras/194.2.475 CrossRef Search ADS   Turolla R., Nobili L., 1988, MNRAS , 235, 1273 https://doi.org/10.1093/mnras/235.4.1273 CrossRef Search ADS   Turolla R., Nobili L., Calvani M., 1986, ApJ , 303, 573 https://doi.org/10.1086/164103 CrossRef Search ADS   Turolla R., Zampieri L., Nobili L., 1995, MNRAS , 272, 625 https://doi.org/10.1093/mnras/272.3.625 CrossRef Search ADS   Zane S., Turolla R., Nobili L., Erna M., 1996, ApJ , 466, 871 https://doi.org/10.1086/177561 CrossRef Search ADS   © 2018 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society

### Journal

Monthly Notices of the Royal Astronomical SocietyOxford University Press

Published: May 1, 2018

## You’re reading a free preview. Subscribe to read the entire article.

### DeepDyve is your personal research library

It’s your single place to instantly
that matters to you.

over 12 million articles from more than
10,000 peer-reviewed journals.

All for just $49/month ### Explore the DeepDyve Library ### Unlimited reading Read as many articles as you need. Full articles with original layout, charts and figures. Read online, from anywhere. ### Stay up to date Keep up with your field with Personalized Recommendations and Follow Journals to get automatic updates. ### Organize your research It’s easy to organize your research with our built-in tools. ### Your journals are on DeepDyve Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more. All the latest content is available, no embargo periods. ### Monthly Plan • Read unlimited articles • Personalized recommendations • No expiration • Print 20 pages per month • 20% off on PDF purchases • Organize your research • Get updates on your journals and topic searches$49/month

14-day Free Trial

Best Deal — 39% off

### Annual Plan

• All the features of the Professional Plan, but for 39% off!
• Billed annually
• No expiration
• For the normal price of 10 articles elsewhere, you get one full year of unlimited access to articles.

$588$360/year

billed annually

14-day Free Trial