# Position dependent radiation fields near accretion disks

Kara Smith Department of Physics & Astronomy
University of Nevada, Las Vegas
4505 S. Maryland Pkwy
Las Vegas, NV, 89154-4002, USA  Nevada Center for Astrophysics
University of Nevada, Las Vegas
4505 S. Maryland Pkwy.
Las Vegas, NV 89154, USA  Daniel Proga Department of Physics & Astronomy
University of Nevada, Las Vegas
4505 S. Maryland Pkwy
Las Vegas, NV, 89154-4002, USA  Nevada Center for Astrophysics
University of Nevada, Las Vegas
4505 S. Maryland Pkwy.
Las Vegas, NV 89154, USA  Randall Dannen Department of Physics & Astronomy
University of Nevada, Las Vegas
4505 S. Maryland Pkwy
Las Vegas, NV, 89154-4002, USA  Nevada Center for Astrophysics
University of Nevada, Las Vegas
4505 S. Maryland Pkwy.
Las Vegas, NV 89154, USA  Sergei Dyda Department of Astronomy
University of Virginia
530 McCormick Rd
Charlottesville, VA 22904, USA  Tim Waters Theoretical Division, Los Alamos
National Laboratory,
NM 87545, USA

###### Abstract

In disk wind models for active galactic nuclei (AGN) outflows, high-energy
radiation poses a significant problem wherein the gas can become overionized,
effectively disabling what is often inferred to be the largest force acting on
the gas: the radiation force due to spectral line opacity. Calculations of
this radiation force depend on the magnitude of ionizing radiation, which can
strongly depend on the position above a disk where the radiation is
anisotropic. As our first step to quantify the position and direction
dependence of the radiation field, we assumed free streaming of photons and
computed energy distributions of the mean intensity and components of flux as
well as energy-integrated quantities such as mean photon energy. We find a
significant dependence of radiation field properties on position, but this
dependence is not necessarily the same for different field quantities. A key
example is that the mean intensity is much softer than the radial flux at many
points near the disk. Because the mean intensity largely controls ionization,
this softening decreases the severity of the overionization problem. The
position dependence of mean intensity implies the position dependence of gas
opacity, which we illustrate by computing the radiation force a fluid element
feels in an accelerating wind. We find that in a vertical accelerating flow,
the force due to radiation is not parallel to the radiation flux. This
misalignment is due to the force’s geometric weighting by both the velocity
field’s directionality and the position dependence of the mean intensity.

galaxies: active - methods: numerical - hydrodynamics - radiation: dynamics

††journal: ApJ

##  1 Introduction

The key properties of active galactic nuclei (AGN), such as compactness, high
luminosity, time variability, spectral energy distribution (SED), and spectral
features can be well-explained as consequences of disk accretion onto a
supermassive black hole (SMBH). Spectral features such as broad emission lines
(BELs) have for decades been used as a basis for classifying AGN, while broad
absorption lines (BALs) in the ultraviolet (UV) are a spectacular
manifestation for gas not only accreting onto an SMBH but also flowing
outward. Even though it has been shown the gas responsible for line
emission/absorption is irradiated by the SMBH accretion disk (e.g., see
Krolik, 1999, for a textbook review), the structure and dynamics of the region
where these broad lines are formed (the broad line region, BLR) remains
elusive.

Spectral features indicate that these accretion disks commonly produce
outflows of matter in the form of jets and winds. The latter is one of the
fundamental mechanisms by which the central SMBH interacts with its host
galaxy (e.g., Giustini & Proga, 2019; Laha et al., 2020, and references
therein). Therefore, the gas outflows are important not only for the inner
workings of AGN but also for AGN feedback (Fabian, 2012). The mechanisms
launching and driving these accretion disk winds include thermal driving,
magnetic driving, and radiation pressure due to photon scattering by spectral
lines (‘line driving’; e.g., see for a review Proga, 2007, and references
therein). While each of these three mechanisms requires its own methods of
study and computation, a commonality between them is that each depends on the
radiation fields of the accretion disk as the source of ionizing radiation.

One major problem for any AGN wind model is overionization of the gas, which
can both suppress or inhibit outflows and fail to reproduce observed line
profiles, namely lines produced from C iv, Si v, and N v. Even a clumpy
magnetically driven outflow not shielded from the powerful AGN radiation will
be overionized and not reproduce observed BALs (de Kool & Begelman, 1995). In
a line driving scenario, the problem is much more severe; an overionized gas
cannot be accelerated by line driving at all (e.g., Stevens & Kallman, 1990;
Proga et al., 2000; Dannen et al., 2019).

Estimates of the ionizing radiation based on constraints of the SED and
luminosity obtained by a distant observer show that for the so-called
unobscured AGN, without shielding, the wind would be overionized (Murray et
al., 1995). However, a fluid element in the wind is not a distant observer; it
sees the radiation source as a nearby observer. Regardless of which driving
mechanism is favored or dominant, understanding the radiation fields seen by
the fluid element in the wind is fundamental to determining if the wind is as
overionized as the distant observer implies.

To understand the radiation fields as seen by the outflow, we compute mean
intensity, flux, mean photon energy, and the line force as a function of
position. We find that a simple multi-temperature blackbody treatment, which a
distant observer sees, is woefully inadequate for describing the radiation
very near to the disk. We have developed a code capable of accurately
computing radiation fields near the disk and can do so inexpensively in terms
of computational resources. In this paper, we use our code to explore
variables that contribute to line driving.

The paper is organized as follows. In §2, we give a basic overview of the
theory behind radiation fields and line-driven winds from accretion disks. In
§3, we discuss the parametrization of the disk and the computation of
radiation fields. In §4, we discuss the results of the study; including
comparisons between SEDs at different points relative to the disk (§4.1), the
sensitivity of mean photon energy, $\left\langle h\nu\right\rangle$ and mean
ionizing intensity, $J_{\rm{ion}}$, to changes in the coronal radiation
(§4.3), and finally, the differences in the direction between the radiation
flux and the radiation force for a vertical velocity field (§4.2). The paper
closes with a summary and discussion of the results in §5.

##  2 Preliminaries

###  2.1 From Stars to Disks: Necessary Changes in Perspective.

Models of outflows in AGN are informed by the relatively well-understood winds
from O stars (for a textbook review see Lamers & Cassinelli, 1999, and
references therein). However, one must remember that stars are simpler than
disks when discussing the radiation field. For example, stellar photospheres
are relatively uniform in temperature; therefore, under free-streaming
conditions, the specific intensity, $I_{\nu}$, above the star is direction
independent inside the solid angle subtended by the star, $\Omega_{*}$.
Further, when free-streaming is applicable, $I_{\nu}$ equals the intensity of
the stellar photosphere, $I_{\nu,0}$. Therefore, calculations that require
integration over $\Omega_{*}$ can be simplified, e.g., when calculating the
mean intensity, $J_{\nu}=(4\pi)^{-1}\oint_{\Omega_{*}}I_{\nu}\ d\Omega$, one
can remove the term for specific intensity from the geometric integral,
yielding the well-known relation $J_{\nu}=I_{\nu,0}\mathcal{W}(r)$, where
$\mathcal{W}(r)=\frac{1}{2}(1-\sqrt{1-(r_{*}/r)^{2}})$ is the geometric
dilution factor and $r_{*}$ is the stellar radius. For $r\approx r_{*}$,
$\mathcal{W}(r)\approx 1/2$ and for large distances,
$\mathcal{W}(r)\approx(r_{*}/2r)^{2}$ (see Lamers & Cassinelli, 1999).

Disks are flattened structures that are not uniform in temperature. Not only
the distance, $r$ but also the viewing angle, $\theta$, measured from the
normal of the disk are necessary variables to consider in describing any
observed properties of the disk. A nearby observer is taken to be a fluid
element with a radius less than the disk’s outer radius, $r_{\rm out}$. For
fluid elements very near to the disk, the mean intensity does approach
$J_{\nu}\rightarrow\frac{1}{2}I_{\nu,0}$, as it does for stars. However, the
specific intensity above the disk is not direction-independent within the
solid angle subtended by the disk, $\Omega_{D}$. Instead, the disk
photospheric temperature and specific intensity are typically a decreasing
function of the radius along the disk, $r_{D}$. In a standard accretion disk
theory (Shakura & Sunyaev, 1973), the disk photospheric temperature,
$T(r_{D})$, follows a $r_{D}^{-3/4}$ scaling for large $r_{D}$ (a similar
scaling can also be derived for an irradiated disk) Assuming that the disk
radiates as a blackbody at every point, $I_{\nu,0}(r_{D})$ follows the famous
$r_{D}^{-3}$ scaling for large $r_{D}$, so near the disk, $J_{\nu}$ scales as
$1/r^{3}$ rather than $1/r^{2}$.

If the disk were of constant specific intensity, then the radial radiation
flux for a large $r$ and given $\theta$ would be
$F_{\nu,r}=I_{\nu,0}A_{D}\cos\theta/r^{2}$, where $A_{D}$ is the area of the
disk. The factor of $\cos\theta$ accounts for the geometrical foreshortening
effect. Because disks have non-uniform specific intensity profiles, we must
integrate them over the entire disk surface. A simple way to do this is to
assume that the angular distance between points on the disk is negligible and
only to integrate over radii (see, e.g., Frank et al., 2002). This is only a
good approximation for the radial flux as seen by a very distant observer,
though, so we denote it $F_{\nu,dist}$:

$F_{\nu,dist}=\frac{2\pi\cos\theta}{r^{2}}\int_{r_{in}}^{r_{out}}I_{\nu,0}(r_{D})\cdot
r_{D}\ dr_{D},$ (1)

where $r_{in}$ is the radius of the disk’s inner edge. This estimate of flux
is not adequate for computing properties from the point of view of a nearby
fluid element, as non-radial flux components comprise a significant portion of
the total flux near the disk. In fact, as we approach arbitrarily close to the
disk, the radiation flux is dominated by local radiation. Each point on the
disk radiates a single-temperature blackbody, not the multi-temperature
blackbody produced by Equation (1). Even if a fluid element receives
substantial radiation from the entire disk, assuming negligible angular
differences between points on the disk can lead to inaccurate estimates.
Finally, the scaling in Equation (1) is $1/r^{2}$, but the flux should be
proportional to $1/r^{3}$ near the disk, as noted above. Therefore, the
distant observer approximation (Equation 1) cannot be applied to points close
to the disk. However, it is still useful for benchmarking numerical
calculations, as demonstrated by previous studies (e.g., see (Proga et al.,
1998)).

###  2.2 Line Driven Disk Winds

One of the primary quantities derived from the radiation field is the
radiation force due to lines, $\mathbfcal{F}^{rad}$. Winds can be launched
when continuum photons interact with gas through scattering in spectral lines
(Castor et al., 1975, hereafter CAK). CAK expresses the radial component of
the radiation force, $\mathcal{F}^{rad}_{r}$, in terms of the radial radiation
flux, $F_{r}$, as

$\mathcal{F}^{rad}_{r}=M(t)\frac{\sigma_{e}}{c}F_{r},$ (2)

where $M(t)$ is the force multiplier, which is a function of optical depth
parameter $t$, $\sigma_{e}$ is the electron scattering coefficient, and $c$ is
the speed of light. CAK defined $t$ using the Sobolev approximation (Rybicki &
Hummer, 1978), $t\equiv\sigma_{e}v_{th}\rho|dv_{r}/dr|^{-1}$, where $v_{th}$
is the thermal velocity of the gas, and $\rho$ is the density. Notice that
this definition of t assumes a radial flow; additionally, CAK’s theory uses
radial force and flux as analogs for their total values because the theory is
developed for optically thin gas around a spherically symmetric star. CAK’s
theory has been successfully extended to a disk (Proga et al., 1998, 1999,
hereafter PSD98 and PSD99). However, the gas velocity is no longer radial, so
instead.

$t\equiv\rho\sigma_{e}v_{th}\left|\frac{dv_{\ell}}{d\ell}\right|^{-1},$ (3)

where $\ell$ is the distance from the emitting point on the disk along a given
line of sight defined by the unit vector $\hat{\boldsymbol{n}}$ and $v_{\ell}$
is the velocity along the line of sight:
$v_{\ell}=\mathbf{v}\cdot\hat{\boldsymbol{n}}$. To compute the radiation
force, we must integrate the intensity over all such lines of sight, and
because $dv_{\ell}/d\ell$ is direction dependent, $\mathcal{M}(t)$ cannot be
removed from this integral. Thus, we have an expression for the radiation
force,

$\mathbfcal{F}^{rad}=\oint_{\Omega_{D}}\mathcal{M}(t)\cdot\frac{\sigma_{e}I}{c}\hat{\boldsymbol{n}}d\Omega,$
(4)

where $I$ is the frequency integrated specific intensity of the disk. In a
sense then, Equation (2) is the “distant observer” analog of Equation (4),
where $I$ is non-zero for only one direction: purely radial. CAK developed a
correction to account for the changes in velocity along a non-radial line of
sight. This correction, elaborated upon by Friend & Abbott (1986) and
Pauldrach et al. (1986), uses the symmetry of stars to integrate over all
lines of sight analytically. The line force remains purely radial, though,
also due to symmetry. For a disk, we cannot use such symmetry; instead, we
must use fast numerical methods to account for sufficiently many lines of
sight to the disk to sample $\Omega_{D}$. One might assume that, generally,
the line force has the same direction as the radiation flux – as is done to
derive Equation (2) – even when that direction is not purely radial, and that
the biggest term of $dv_{\ell}/d\ell$ suffices to compute the underlying
velocity shear tensor (Proga et al., 1998; Nomura et al., 2020). In section
§4.2, however, we show that the force and flux are not aligned near a disk of
non-uniform brightness.

The CAK line-driven wind theory informed earlier studies of line-driven winds
from disks around cataclysmic variables (Proga et al., 1998, 1999, CVs;) and
in AGN (Murray et al., 1995; Proga & Kallman, 2004). In CVs, launching and
driving winds are robust, as the winds are optically thin to continuum photons
and the radiation is dominated by UV, similar to the O stars modeled by CAK.
In AGN, however, the acceleration of a wind is less robust. The wind can be
accelerated to high velocities only if both the disk continuum radiation is
strong in the UV and the wind can be shielded from the ionizing radiation from
the inner disk and the X-ray corona (Dyda et al., 2023). Because AGN winds are
much denser than those of O stars and CVs, gas launched from large radii can
be accelerated when shielded by a dense flow launched from smaller radii that
fail to escape the disk due to overionization (Proga et al., 2000, hereafter
PSK).

In PSK, the radiation transfer of the ionizing radiation was treated in a
simplified way, namely by accounting only for the attenuation of direct
radiation. Such a treatment results in a lower limit for the mean intensity of
ionizing radiation $J_{ion}\equiv\int_{\nu_{H}}^{\infty}J_{\nu}\ d\nu$, where
$\nu_{H}$ corresponds to 13.6 eV and $J_{\nu}$ is the mean intensity as a
function of frequency. It does not account for scattered radiation found by
Sim et al. (2010) and Higginbottom et al. (2014). The latter studies used
Monte Carlo photoionization calculations and showed that radiation scattering
can allow ionizing photons to propagate ‘around’ the shielding inner flow and
thereby increase the ionization of the outer outflow. The techniques used by
Sim et al. (2010) and Higginbottom et al. (2014) are computationally
expensive, and as such, it is not yet feasible to use them in the context of
radiation-hydrodynamical simulations, wherein we would need to compute line
distributions for each grid point at each time step. It is necessary, then, to
find computationally inexpensive and sufficiently accurate estimators of the
ionization state of the gas.

The main three quantities that can be used to estimate ionization are the
above mentioned $J_{\rm{ion}}$, gas number density $n$, and mean photon
energy, $\left\langle h\nu\right\rangle\equiv\int_{0}^{\infty}h\nu J_{\nu}/J\
d\nu$, where $h$ is Planck’s constant. In optically thin plasmas, the
ionization state is determined by the ratio of the first two quantities, as
shown by Tarter et al. (1969), referred to as the ionization parameter,
$\xi\equiv(4\pi)^{2}J_{\rm{ion}}/n$. The mean photon energy is vital for high
values of $\xi$, where the main heating is Compton heating, which heats the
gas to the Compton temperature, $T_{C}\equiv\left\langle
h\nu\right\rangle/(4k_{B})$. If $T_{C}>10^{6}$ K, the gas will undergo runaway
heating and could even become thermally unstable (e.g., Kallman & McCray,
1982; Dyda et al., 2017; Proga et al., 2022). We find that, without
attenuation of the central source, when the central source has just 5% the
luminosity of the disk, all positions above the disk will suffer from runaway
heating effects, and, even at only 1%, a significant portion of points will
overheat as well.

##  3 Calculating Radiation and Radiation Force

To evaluate the properties of the radiation field, it is helpful to start with
the fundamental equation governing the coupling between radiation and matter
as a function of position along a line of sight, $\ell$, and time, $t$, which
is the radiation transfer (RT) equation:

$\frac{1}{c}\frac{\partial I_{\nu}}{\partial
t}+\hat{\boldsymbol{n}}\cdot\nabla I_{\nu}=j_{\nu}-\chi_{\nu}I_{\nu},$ (5)

where $j_{\nu}$ is the emissivity and $\chi_{\nu}$ the total opacity. Here, we
consider a wind with conditions similar to those studied by CAK. Specifically,
we focus on the time-independent solutions and assume that the wind emission
is negligible compared to the photospheric emission so that we may set
$j_{\nu}$ to zero. In addition, we consider that the opacity due to
scattering, $\chi^{sc}_{\nu}$, and absorption, $\chi^{abs}_{\nu}$, (i.e.,
$\chi_{\nu}=\chi^{sc}_{\nu}+\chi^{abs}_{\nu}$ with units $\rm cm^{-1}$) for
continuum processes. This corresponds to the situation where continuum photons
emitted at the photosphere can stream freely through the wind. Again,
following CAK, we will account for RT in spectral lines by using the Sobolev
approximation. This approximation takes advantage of the fact that when a
photon is emitted at a certain frequency, $\nu$, in the atmosphere, it
interacts with the outflowing gas through a specific line only when the
frequency is very close to the line center frequency, $\nu_{0}$. This
resonance point, which is determined by the Doppler shift, occurs when the
photon frequency matches the line center frequency. Under these assumptions,
$I_{\nu}$ can be treated as constant along any line of sight, and we need only
calculate the intensity of the photosphere. So, for both the corona and the
disk, we let $I_{\nu}=I_{\nu,0}$ for continuum photons that do not interact
via lines, and $I_{\nu}=I_{\nu,0}~{}{\rm exp}(-\tau_{\nu})$ for the photons
that do interact, where $\tau_{\nu}$ is the optical depth due to the
absorption line in the wind.

Because $I_{\nu}=I_{\nu,0}$ for continuum photons, calculations of moments of
intensity can be simplified and, therefore, sped up. To do this, we use the
methods described in PSD99 and adopted by Dyda et al. (2017). We use a
spherical-polar coordinate system $(r_{1},\theta_{1},\phi_{1})$ centered on a
point in the wind, $W$ Let the point $W^{\prime}$ be the point on the disk
vertically below $W$. The coordinate system’s zenith is directed along
$\overline{WC}$ – where $C$ is the center of the disk – and $\phi_{1}$ = 0 the
half-plane opposite $\overline{WC}$ from $W^{\prime}$. We also introduce a
shorthand $\mu=\cos\theta_{1}$. We can write the expressions for
$\hat{\boldsymbol{n}}$ and $\Omega$ in $(r_{1},\theta_{1},\phi_{1})$-space as

$\hat{\boldsymbol{n}}=(\mu,\ -\sqrt{1-\mu^{2}}\cos\phi_{1},\
\sqrt{1-\mu^{2}}\sin\phi_{1})$ (6)

and

$d\Omega=\text{sin}\theta_{1}\ d\phi_{1}d\theta_{1}=-\ d\phi_{1}d\mu.$ (7)

Equations for the zeroth moment, the mean intensity, $J_{\nu}$, and the first
moment, or the flux, $F_{\nu}$, are quite simple in this coordinate system:

$J_{\nu}(r,\theta)=\frac{1}{4\pi}\oint_{\Omega(r,\theta)}I_{\nu,0}(r,\theta,\mu,\phi_{1})\
d\Omega$ (8)

and

$\mathbf{F}_{\nu}(r,\theta)=\oint_{\Omega(r,\theta)}I_{\nu,0}(r,\theta,\mu,\phi_{1})\hat{\boldsymbol{n}}(\mu,\phi_{1})d\Omega.$
(9)

In Equations (8) and (9), we explicitly show $J_{\nu}$, $\mathbf{F}_{\nu}$,
$\Omega$, $I_{\nu,0}$, and $\hat{\boldsymbol{n}}$ as functions of $r$,
$\theta$, $\mu$ and $\phi_{1}$. For the sake of brevity, we omit these
explicit dependencies in the rest of this paper.

To compute integrals over d$\Omega$ numerically, we use the Gauss-Legendre
quadrature, in which we divide the angular $d\mu d\phi_{1}$-space into a
discrete $\Delta\mu\Delta\phi_{1}$-sum:

$J_{\nu}=\frac{1}{4\pi}\sum_{i}\sum_{j}I_{\nu,0}(\mu_{i},{\phi_{1}}_{j})\omega_{i}\omega_{j}$
(10)

and

$\mathbf{F}_{\nu}=\sum_{i}\sum_{j}I_{\nu,0}(\mu_{i},{\phi_{1}}_{j})\hat{\boldsymbol{n}}(\mu_{i},{\phi_{1}}_{j})\omega_{i}\omega_{j},$
(11)

where $\mu_{i}$ and ${\phi_{1}}_{j}$ are the roots of Legendre polynomials
while $\omega_{i}$ and $\omega_{j}$ are the associated weights (for a review
of Gauss-Legendre quadrature as well as other methods of numerical integration
see Press et al., 2007). In this paper, we used a 15x15 grid of quadrature
points, which is a relatively low number of grid points, especially compared
to older methods such as the one used in PSD98. With such a grid, on a single
12th Gen Intel i7-12700H Core (i.e., a modern laptop), computing SEDs for
$J_{\nu}$ and $\mathbf{F_{\nu}}$ sampled at 300 frequencies takes 0.1 seconds
per location, ($r$, $\theta$). Using Gauss-Laguerre Quadrature (again, see
Press et al., 2007) with 50 frequency quadrature points, we can accurately
compute $J$ and $\mathbf{F}$ in about 0.02 seconds per location. This
computational speed makes feasible the goal of incorporating these methods in
a time-dependent radiation-hydrodynamical simulation.

We obtain an explicit expression for $I_{\nu,0}$ following PSK’s model of
radiation from the disk and corona. We separate the intensity into that from
each surface: the intensity of the disk, $I_{\nu,D}$, and that of the corona,
$I_{\nu,*}$. For the disk, we assume a geometrically thin, optically thick
disk with a temperature profile given by Shakura & Sunyaev (1973). The disk
extends from some inner radius, $r_{\rm in}$, to some outer radius, $r_{\rm
out}$. We take $r_{in}$ to be $r_{ISCO}=6GM/c^{2}$, where $M$ is the mass of
the central black hole. Because the exact shape of the central source is not
well known, we assume a spherical extended corona with a radius also equal to
$r_{ISCO}$. The total luminosity of the disk, $L_{D}$, is given by

$L_{D}=\frac{GM\dot{M}}{2r_{in}},$ (12)

where $\dot{M}$ is the mass accretion rate. We parametrize the accretion rate
in terms of the Eddington fraction, $\Gamma\equiv\dot{M}/\dot{M}_{Edd}$, where
$\dot{M}_{Edd}=8\pi cr_{in}/\sigma_{e}$ is the Eddington accretion rate.

We also assume the corona radiates with a power-law spectrum – i.e. the
luminosity $L_{\nu,*}\propto\nu^{\beta}$ for some constant $\beta$ – when
$0.1\ keV<h\nu<511\ keV$ (Sim et al., 2010). We normalize the corona
luminosity by setting, as an input parameter, the ratio $x\equiv L_{*}/L_{D}$.
We include further details about the computation of the intensity of the
corona in §3.1. The accretion disk system is parameterized by $M$, $\Gamma$,
$r_{\rm out}$, $x$, and $\beta$. We list the values used in our calculations
in Table 1.

Table 1: Values of the accretion disk and corona parameters used in this paper. Parameter | Value(s)
---|---
$M$ | $10^{8}M_{\odot}$
$\Gamma$ | 0.7
$r_{out}$ | 100 $r_{\rm in}$
$x$ | 0.001-0.25
$\beta$ | -1.1

###  3.1 Evaluation of Intensity

To develop an expression for the intensity of the disk, we begin with the
surface brightness of the disk, $Q$. The total surface brightness is the sum
of two effects:

$Q=Q_{int}+Q_{irr},$ (13)

where $Q_{int}$ is the intrinsic brightness of the disk due to internal
processes and $Q_{irr}$ is the brightness due to irradiation by an outside
source (in our case, the black hole corona). The intrinsic brightness is given
by the Shakura-Sunyaev model and the brightness due to irradiation is given by
PSD98:

$Q_{int}=Q_{0}\cdot\frac{r_{\rm in}^{3}}{r_{D}^{3}}\left(1-\sqrt{\frac{r_{\rm
in}}{r_{D}}}\right)$ (14)
$Q_{irr}=\frac{xQ_{0}}{3\pi}\left[\text{arcsin}\left(\frac{r_{\rm
in}}{r_{D}}\right)-\frac{r_{\rm in}}{r_{D}}\sqrt{1-\left(\frac{r_{\rm
in}}{r_{D}}\right)^{2}}\right],$ (15)

where $Q_{0}=3GM\dot{M}/8\pi r_{\rm in}^{3}$. Assuming that the disk radiates
as a blackbody at every point, we can use the Stefan-Boltzmann law to write
temperature as a function of brightness and, thereby, as a function of radius
along the disk.

$T(r_{D})=\left(\frac{Q_{int}(r_{D})+Q_{irr}(r_{D})}{\sigma}\right)^{1/4},$
(16)

where $\sigma$ is the Stefan-Boltzmann constant.

As per our assumption that the disk radiates as a blackbody, we take the
intensity from a surface element at radius $r_{D}$ to be the blackbody
intensity of that surface element.

$I_{\nu,D}(r_{D})\equiv B_{\nu}(T(r_{D})),$ (17)

where $r_{D}$ is given by

$\displaystyle r_{D}^{2}=r_{1}^{2}+r^{2}-2rr_{1}\mu$ (18) $\displaystyle
r_{1}=\frac{r\cos\theta}{\mu\cos\theta-\sqrt{1-\mu^{2}}\cos\phi_{1}\sin\theta},$

where $r_{1}$ is the distance from $W$ to the disk along the line of sight.
Connecting back to Equations (8) and (9), $I_{\nu,D}$ is a function of $r$,
$\theta$, $\mu$ and $\phi_{1}$ through the fact that $r_{D}$ is a function of
those four variables.

We then model the corona’s specific luminosity according to a power-law
spectrum with cutoffs at $h\nu_{X}=0.1\ keV$ and $h\nu_{\gamma}=511\ keV$.

$L_{\nu,*}=\begin{cases}xL_{0}\,\nu^{\beta}&\nu_{X}<\nu<\nu_{\gamma}\\\
0&otherwise\end{cases}$ (19)

where $L_{0}$ is a normalization constant such that $x=L_{*}/L_{D}$. The
intensity of the corona is then given in terms of the luminosity,

$I_{\nu,*}=\frac{L_{\nu,*}}{4\pi^{2}r_{*}^{2}},$ (20)

where $r_{*}$ is the radius of the corona. For our purposes $r_{*}=r_{\rm
in}$.

###  3.2 Radiation Force due to Lines

Computing the radiation force due to a line requires the local line absorption
coefficient, $\kappa_{L}$ – which is related to the opacity in Equation 5 by
$\chi_{l,\nu}=\rho\kappa_{L}\phi(\nu)$, where $\rho$ is the density of the gas
and $\phi(\nu)$ is the assumed line absorption profile – and the local
intensity, $I_{\nu}$.

$d\mathbfcal{F}_{L}^{rad}=\frac{\Delta\nu_{D}\kappa_{L}}{c}I_{\nu}\hat{\boldsymbol{n}}d\Omega,$
(21)

where $\Delta\nu_{D}\equiv\nu_{0}v_{th}/c$ is the Doppler width of the line.
To integrate Equation (21) over d$\Omega$, we cannot simply use
$I_{\nu}=I_{\nu,0}$, as this is not a continuum process; we must include an
escape probability to account for the interaction of photons with spectral
lines along a line of sight (Owocki & Puls, 1996). For a derivation of this
escape probability from $\phi(\nu)$, see Lamers & Cassinelli (1999). With the
escape probability included, we have an integral for the radiation force due
to lines in terms of the intensity of the disk photosphere.

$\mathbfcal{F}^{rad}_{L}=\oint_{\Omega}\frac{\sigma_{e}\Delta\nu_{D}}{c}\frac{1-e^{-\tau_{L}}}{\tau_{L}}I_{\nu,0}\hat{\boldsymbol{n}}d\Omega,$
(22)

where $\tau_{\nu}$ is the local optical depth, which we use the Sobolev
approximation in three dimensions (see Rybicki & Hummer, 1978) to compute.

$\tau_{L}=\int_{r}^{\infty}\rho\kappa_{L}\phi(\nu)\
dr\mathrel{\stackrel{{\scriptstyle\makebox[0.0pt]{\mbox{\tiny
Sob}}}}{{=}}}\rho\kappa_{L}v_{th}\left|\frac{dv_{\ell}}{d\ell}\right|^{-1},$
(23)

where $\phi(\nu)$ is the aforementioned assumed line absorption profile.

To compute the total radiation force due to lines, we follow CAK and define a
variable $\eta\equiv\kappa_{L}/\sigma_{e}$ so that $\tau_{L}=\eta t$.
Substituting in $t$ and $\eta$ into Equation (21), we get an equation of
radiation force due to a single line:

$\mathbfcal{F}^{rad}_{L}=\oint_{\Omega}\frac{\sigma_{e}\Delta\nu_{D}}{c}\cdot\frac{1-e^{-\eta
t}}{t}\cdot I_{\nu,0}\hat{\boldsymbol{n}}d\Omega.$ (24)

Then, by summing $\mathbfcal{F}^{rad}_{L}$ over all lines, we obtain the full
radiation force due to all lines:

$\mathbfcal{F}^{rad}=\sum_{\text{lines}}\oint_{\Omega}\frac{\sigma_{e}\Delta\nu_{D}}{c}\frac{1-e^{-\eta
t}}{t}\cdot I_{\nu,0}\hat{\boldsymbol{n}}d\Omega.$ (25)

Figure 1: Comparisons of the energy distribution of blackbody radiation
(black-dashed curves) and single line force multipliers (colored vertical
bars) for four blackbody temperatures: $T_{r}$ = $10^{5}$ K, $2\times 10^{5}$
K, $4\times 10^{5}$ K, $6\times 10^{5}$ K. The vertical bars represent the
1,000 strongest individual line force multipliers – normalized to the maximum
force multiplier, $M_{l,max}$ – for fixed values of $\log\xi=2$ and
$t=10^{-10}$ (Eq. 3) as a function of energy. These 1,000 strongest lines have
been binned into quartiles, red being the upper quartile, yellow the second,
green the third, and blue the weakest set of lines. The black dashed line in
each panel represents the assumed blackbody SED (temperature indicated in each
panel) that we used for both determining the ionization balance and computing
$M$. The minimum force multiplier, $M_{l,min}$ is also listed in each panel.

###  3.3 The Distribution of Lines

We can use photoionization codes such as XSTAR  (Kallman & Bautista, 2001) to
compute the abundances of ions and opacities of lines for any given SED, which
allows us to compute the sum in Equation (25). But, as mentioned before, doing
so in the context of radiation-hydrodynamical simulations is computationally
expensive. So, instead, we again follow CAK and describe a statistical
distribution of lines, where the number of lines of a given opacity and
frequency is a function of that opacity and frequency. To recover the power
law given in CAK, $M(t)=kt^{-\alpha}$ for some constant $k$, we have

$dN=N_{0}\eta^{\alpha-2}\frac{1}{\nu}\ d\eta d\nu,$ (26)

where $\alpha$ is the CAK parameter. Owocki et al. (1988, hereafter OCR), add
an exponential drop-off to limit the effect of very strong driving lines. With
this modification, we have

$dN=N_{0}\eta^{\alpha-2}e^{-\eta/\eta_{\rm max}}\frac{1}{\nu}d\eta d\nu,$ (27)

where $\eta_{max}$ is the maximum opacity for which lines are physically
relevant.

We are interested in finding the distribution of lines relative to energy. CAK
introduced a frequency distribution $1/\nu$, which allows for the frequency
dependency of the force multiplier to cancel out entirely. However, based on
photoionization calculations using the methods outlined by Dannen et al.
(2019), we see that the line distribution is not a simple power-law. For
example, in Figure 1, we show that the distribution of lines actually traces
the energy distribution of the mean radiation intensity. This informs our
second modification to the line distribution, which is to use $J_{\nu}/J$ to
weight our lines by frequency rather than the classic $1/\nu$. Combining the
modifications of OCR with our own, we arrive at

$dN=N_{0}\eta^{\alpha-2}e^{-\eta/\eta_{max}}\frac{J_{\nu}}{J}\ d\eta d\nu,$
(28)

We can replace the sum in Equation (25) with an integral over $dN$ and
analytically integrate over $d\eta$ to find our, now $\nu$ dependent, force
multiplier.

$\mathcal{M}(t,\nu)=\frac{\nu J_{\nu}}{J}\cdot
kt^{-\alpha}\cdot\left[\frac{(\tau_{max}+1)^{1-\alpha}-1}{\tau_{max}^{1-\alpha}}\right],$
(29)

where $\tau_{max}=\eta_{max}t$. Accounting for the SED dependence of the line
distribution, not only must the force multiplier be integrated over the solid
angle of the disk but also over all frequencies. We therefore have

$\mathbfcal{F}^{rad}=\int_{0}^{\infty}\oint_{\Omega}\mathcal{M}(t,\nu)\cdot\frac{\sigma_{e}I_{\nu,0}}{c}\hat{\boldsymbol{n}}d\Omega
d\nu$ (30)

as the full expression for the radiation force due to lines.

Figure 2: Comparison of the energy distributions of the disk radiation field
properties at three locations: (100 $r_{in}$, 40∘), (12 $r_{in}$, 85∘), (40
$r_{in}$, 80∘), top to bottom. Each panel shows $4\pi J_{\nu}$ (blue),
$F_{\nu,r}$ (orange), $F_{\nu,\theta}$ (green), $F_{\nu,dist}$ (red dashed),
and the blackbody spectrum of the point projected vertically onto the disk
(black dashed). The area below $F_{\nu,dist}$ is shaded to make it more
apparent at what energies a spectral property is greater or less than that
seen by a distant observer.  Figure 3: Illustrating the misalignment of the
line force and the radiation flux. In a vertical velocity field above a disk,
the line force (black) is not parallel to the radiation flux (red), as it
would be in a spherically symmetric wind case.

##  4 Results

###  4.1 Local vs Distant Radiation Fields

We compute energy distributions $J_{\nu}(r,\theta)$, $F_{\nu,r}(r,\theta)$,
and $F_{\nu,\theta}(r,\theta)$ for $r$ ranging from $2r_{in}$ to $100r_{in}$
and $\theta$ ranging from $0^{\circ}$ to $90^{\circ}$. In Figure 2 we show
example SEDs for three locations, $(r,\theta)$, above the disk: (100, 40),
(40, 80), (12, 85). We selected these three locations as illustrative
representations of the three distinct classes of spectra we find above the
disk.

Specifically, the top panel of Figure 2 shows the spectra at large radius and
intermediate theta, ($r=100\ r_{in},\ \theta=40^{\circ}$). In this regime,
both the spectra of $4\pi J_{\nu}$ (blue solid line) and $F_{\nu,r}$ (orange
solid line) overlap with each other. Moreover, they agree well with the
approximation for the flux seen by a distant observer given in Equation (1),
$F_{\nu,dist}(r,\theta)$ (red dashed line). Agreement between these three SEDs
is consistent with the expectations discussed in §2. First, the $\theta$
component of the flux (green solid line) will disagree with $4\pi J_{\nu}$,
$F_{\nu,r}$, and $F_{\nu,dist}$ because it is weighted by $n_{\theta}$. And
further, $F_{\nu,\theta}$ will be much weaker than $F_{\nu,r}$, which we also
see, $F_{\nu,\theta}$ peaks an order of magnitude lower than $F_{\nu,r}$. For
reference, we plot a curve for a blackbody spectrum with the temperature of
the point projected down onto the disk. We see that $F_{\nu,\theta}$ peaks at
similar energies to the peak of the blackbody spectrum, but the spectrum is
broader, again due to the entire disk contributing to the flux.

In the case that we are very near the disk, shown in the middle panel of
Figure 2, ($r=12\ r_{in},\ \theta=85^{\circ}$), we see that $F_{\nu,\theta}$
not only peaks at the same energy as the blackbody spectrum, but is very
similar to it at all energies, save for a small, high-energy tail. The
$J_{\nu}$ spectrum is the most similar, in fact, aside from the high energy
tail, $4\pi J_{\nu}=2F_{\nu,\theta}$. It is to be expected that contribution
to $J_{\nu}$ is dominated by $F_{\nu,\theta}$ as radiation from any location
not immediately below the point $(r,\theta)$ will be diminished by geometric
foreshortening. The other two spectra are quite different from the
$F_{\nu,\theta}$, $4\pi J_{\nu}$, and local blackbody spectra. While
$F_{\nu,r}$ spectrum is relatively narrow, it peaks at higher energies than
$J_{\nu}$ and $F_{\nu,\theta}$ because $F_{\nu,r}$ is weighted by $n_{r}$, and
so is dominated by the high energy radiation from the inner disk. On the other
hand, the $F_{\nu,dist}$ spectrum is still a broad multi-temperature
blackbody, as it’s intensity varies uniformly, independent of frequency.

From these two cases, it may seem reasonable to assume that as one moves
further from the disk, $J_{\nu}$ broadens, $F_{\nu,r}$ broadens and increases
in power, and $F_{\nu,\theta}$ broadens slightly but decreases in power.
However, the bottom panel of Figure 2, ($40\ r_{in},\ \theta=80^{\circ}$),
shows that a simple interpolating function between the narrow spectum near the
disk and broad distant observer spectrum does not reproduce the spectra at
intermediate positions. At points with intermediate radii and fairly high
viewing angles, $F_{\nu,r}$ experiences contributions from a larger portion of
the inner disk, broadening its spectrum. The part of the disk closer to the
point in the wind contributes more to the spectrum and is at lower energies.
This leads to a low energy bump in the $F_{\nu,r}$ spectrum. The $J_{\nu}$
spectrum has a more prominent bump that covers both the $F_{\nu,r}$ and
$F_{\nu,\theta}$ peaks.

Figure 4: Position dependence of log${}_{10}(J_{\rm{ion}})$ and
log${}_{10}(\left\langle h\nu\right\rangle)$ for different corona intensities.
The top two rows of panels show contours of $J_{\rm{ion}}$ in units of erg s-1
cm-2, both on its own (top) and scaled by $\mathcal{W}(r)$ (middle). The
bottom panel shows $\left\langle h\nu\right\rangle$ in units of eV. Each
column is computed for a coronal luminosity fraction $x$. The colors are for
visual contrast only. We demonstrate the position dependence of these
quantities as well as the change in that position dependence due to the
luminosity of the corona.  Figure 5: Relations between the energy integrated
properties of the disk radiation field with and without coronal radiation. The
top panels show $\left\langle h\nu\right\rangle$ vs $J$ and the bottom show
$\left\langle h\nu\right\rangle$ vs $J_{\rm{ion}}$. The left column shows each
point colored according to $r$. The right shows each point colored according
to $\theta$. The vertical lines are an analytic estimate for $\left\langle
h\nu\right\rangle$ based on the assumption that $J_{*}/J_{D}\approx x$; dashed
is $x$ = 0.01, dotted is $x$ = 0.05, and dash-dotted is $x$ = 0.25. We do not
show these vertical lines for $x=0$ and $x=0.001$ because the approximation
only holds when $x\gg\left\langle h\nu\right\rangle_{D}/\left\langle
h\nu\right\rangle_{*}\sim 10^{-3}$ for our disk parameters. The maximum value
of $\left\langle h\nu\right\rangle$ in each panel is $\left\langle
h\nu\right\rangle_{*}\approx$ 42 keV.

The position dependences of $J_{\nu}$ are important for our discussion of the
radiation force due to lines, as they determine what lines are available for
line driving, as discussed in §3.3. This means that in some sense
$\mathbfcal{F}^{rad}$ is doubly position dependent. It depends on the position
like SEDs do directly through its integration over $\Omega$. But also, it is
position dependent on what lines are available through $J_{\nu}$. In the
bottom panel of Figure 2, we see that the $J_{\nu}$ spectrum peaks at a very
different position than the $F_{\nu,r}$ spectrum, but relatively close to the
peak of the $F_{\nu,\theta}$ spectrum. This biases the force toward the
$\theta$ direction, as there will be more lines that the radiation comprising
$F_{\nu,\theta}$ can interact with. Such a bias creates a misalignment between
the radiation force and the flux. The velocity of the gas further biases which
component of the flux predominantly contributes to the line force. In the
following subsection, we discuss these two effects and determine their
relative contribution to causing the misalignment.

###  4.2 The Misalignment of $\mathbfcal{F}^{rad}$ and $\mathbf{F}$

For a given velocity field, we can compute $dv_{\ell}/d\ell$ for each line of
sight to the disk, and with gas parameters $k$, $\alpha$, $v_{th}$, and
$\rho$, we are then able to compute the radiation force using Equation (30).
The radiation force is a function of three quantities: the flux from the disk,
the velocity gradient of the accelerating gas, and the lines available for
line driving. The weighting effects of the velocity field cause the force to
be biased toward the direction of greatest acceleration, as it has the
greatest value of $kt^{-\alpha}$. The weighting effects of line availability
cause the force to be more biased by local radiation than the flux alone would
imply. For example, when the point in the wind, $W$, is near enough to the
disk, local radiation is the dominant contributor to $J_{\nu}$ (see Fig. 2).
Because the available lines are determined by $J_{\nu}$, these lines will be
primarily in frequencies produced by the disk directly below $W$, which
creates a vertical bias in the force. The biases due to both the velocity
field and the local SED cause a misalignment between the flux and the force
based on the strength of the directional force multiplier.

Figure 3 shows an example of this misalignment. The red arrows show the
radiation flux purely from the disk ($x=0$). The black arrows show the
radiation force due to this flux in a vertical velocity field, which is the
velocity field we expect to see very near the disk. This velocity field is
defined by a power-law of the vertical distance. In cylindrical coordinates
$(r_{D},z,\phi_{D})$

$v_{z}=v_{0}\left(\frac{z}{r_{D}}\right)^{\delta},$ (31)

where $v_{0}=\sqrt{GM/r_{D}}$, and in spherical coordinates

$v_{z}=\sqrt{\frac{GM}{\sin\theta}}\frac{1}{\tan^{\delta}\theta}\ r^{-1/2}.$
(32)

In a Keplerian disk, $v_{\phi}=v_{0}$ giving us an expression for $\mathbf{v}$
in spherical coordinates:

$\mathbf{v}=v_{z}\left(\cos\theta,-\sin\theta,\tan^{\delta}\theta\right),$
(33)

which we can use to compute $\mathcal{M}(\nu,t)$.

Figure 3 shows the misalignment for $\delta=1/2$. As one moves vertically
upwards from the disk surface, the radiation flux is significantly influenced
by its radial component. But, due to the vertical velocity field, the
radiation force remains close to vertical for much longer than the flux. This
allows the gas to potentially be launched higher than if the force were
parallel to the flux.

As mentioned, this misalignment comes from a mixture of two effects: the
velocity field and the position dependence of the line distribution. To
determine which is the dominant effect, we computed the directions of the
radiation force with the classical $1/\nu$ weighting of frequencies (for which
all frequency dependence of the force multiplier cancels). We found the
difference between the misalignment obtained using Equation (29) and that
obtained using the $1/\nu$ weighting. The points with the largest difference
between the two methods differed by about 10%. We conclude that the
directionality caused by the velocity gradient is the dominant effect causing
this misalignment; but the effect of accounting for the frequency-weighting of
lines is not insignificant.

###  4.3 Position Dependencies of Frequency Integrated Radiation Fields

Along with the radiation force due to lines, other frequency-integrated
properties help determine the state of the gas. The properties of the disk
that we focus on in this study are $\left\langle h\nu\right\rangle$ and
$J_{ion}$. These two quantities largely control the heating and ionization of
the gas surrounding the disk, as discussed in §2. As spectra change with a
position near the disk, so do these frequency-integrated quantities. In Figure
4 we show log${}_{10}(J_{\rm ion})$ (top), log${}_{10}(J_{\rm
ion}/\mathcal{W}(r))$ (middle), and log${}_{10}(\left\langle
h\nu\right\rangle)$ (bottom).

In the case where there is no corona ($x=0$), we see a simple trend: as we
move radially outward, both $J_{ion}$ and $\left\langle h\nu\right\rangle$
decrease. For large enough $r$ and $\theta$, $J_{ion}/\mathcal{W}(r)$ is
radius independent, consistent with a $\mathcal{W}(r)$ scaling. Moving toward
the disk along constant $r<20\ r_{in}$, $J_{ion}$ and $\left\langle
h\nu\right\rangle$ both decrease as $\theta$ increases. At radii larger than
20 $r_{in}$, both increase slightly before decreasing.

When a corona is introduced, the shapes of the contours of constant $J_{ion}$
are mainly unaffected. The values of $J_{ion}$ are increased due to a term
proportional to $\mathcal{W}(r)$ added by the corona. This term never varies
by more than a factor of $\frac{1}{2}\mathcal{W}(r)$ across all $\theta$,
contributing a relatively constant factor at any given radius.

The effect of a corona on $\left\langle h\nu\right\rangle$ is much more
complicated than on $J_{ion}$. Here we introduce a variable $x_{J}\equiv
J_{*}/J_{D}$, where $J_{*}$ and $J_{D}$ are the frequency-integrated mean
intensities of the corona and disk, respectively. We rewrite mean photon
energy as

$\left\langle h\nu\right\rangle=\frac{1}{x_{J}+1}\left\langle
h\nu\right\rangle_{D}+\frac{x_{J}}{x_{J}+1}\left\langle
h\nu\right\rangle_{*},$ (34)

where $\left\langle h\nu\right\rangle_{D}$ and $\left\langle
h\nu\right\rangle_{*}$ are the mean photon energies of the disk and corona,
respectively. This form of the equation is rather useful because $\left\langle
h\nu\right\rangle_{*}$ is a constant; it is approximately 42 keV for our disk
parameters. As $r$ increases, $x_{J}$ approaches $x$ and $\left\langle
h\nu\right\rangle_{D}$ approaches a constant $\left\langle
h\nu\right\rangle_{dist}$ – the mean photon energy of the disk as seen by a
distant observer – which for our disk parameters is $\sim$25 eV.

For constant $\theta\lesssim 40^{\circ}$, as $r$ increases, first
$\left\langle h\nu\right\rangle$ decreases. This is the regime wherein the
radiation from the outer disk becomes more relevant, softening the mean
intensity spectrum. Then, at large $r$, $\left\langle h\nu\right\rangle$
begins increasing with $r$. As $r$ increases, $x_{J}$ approaches $x$, and so
$\left\langle h\nu\right\rangle$ approaches it’s limit: $x/(x+1)\cdot
4.2\times 10^{4}+1/(x+1)\cdot 25$ eV. For our four values of $x$: 0.001, 0.01,
0.05, and 0.25, these limits are 67 eV, 440 eV, 2 keV, and 8.4 keV,
respectively.

The behavior when $\theta\gtrsim 40^{\circ}$ is the opposite of its behavior
in the $x=0$ case. As one moves closer to the disk while maintaining constant
$r$, $\left\langle h\nu\right\rangle$ increases. The cause is that as one
approaches the disk, the mean intensity spectrum experienced loses its high
energy tail due to the foreshortening effect discussed in §4.1, so the total
mean intensity of the disk decreases. The corona does not suffer any
foreshortening and, therefore, becomes the dominant contributor to
$\left\langle h\nu\right\rangle$.

In Figure 5 we show the relation between J, $J_{\rm{ion}}$, and $\left\langle
h\nu\right\rangle$ more directly. We can see that the change in contour shape
present in Figure 4 translates to a reversal in the correlation of $J$ and
$\left\langle h\nu\right\rangle$. Namely, when there is no corona, when one
fixes a value of either $r$ or $\theta$ and allows the other to vary, we find
a positive correlation between $J$ and $\left\langle h\nu\right\rangle$.
However, when a corona is added, most fixed $r$ and $\theta$ yield negative
correlations. Further, from Equation (34), if $x_{J}\gg\left\langle
h\nu\right\rangle_{*}/\left\langle h\nu\right\rangle_{D}\sim 10^{-3}$ we can
ignore the disk contribution to $\left\langle h\nu\right\rangle$ and use
$\left\langle h\nu\right\rangle=x_{J}/(x_{J}+1)\left\langle
h\nu\right\rangle_{*}$. In Figure 5, we show $x/(x+1)\left\langle
h\nu\right\rangle_{*}$, which, if one assumes $x_{J}\approx x$ as it does at
large distances, should be a good estimate of $\left\langle
h\nu\right\rangle$. This benchmark is represented by three vertical lines at
$x/(x+1)\left\langle h\nu\right\rangle_{*}$ for $x$ = 0.01, 0.05, and 0.25
(dashed, dotted, and dash-dotted lines, respectively). In all three cases,
this benchmark overestimates $\left\langle h\nu\right\rangle$ for most points
near the disk, indicating $x_{J}\not\approx x$; rather, $x_{J}<x$, indicating
that near the disk, the coronal contribution to $\left\langle
h\nu\right\rangle$ is less than one would assume from the disk and coronal
luminosities.

##  5 Concluding Remarks

We have investigated how the properties of the radiation field near an
accretion disk change with a position due to purely geometric effects in the
case of a standard, flat Shakura-Sunyaev disk surrounded by optically thin
gas. As explained in §2, quantifying this position dependence is a necessary
prerequisite for accurately modeling important processes in radiation-
hydrodynamical simulations, such as radiative heating and cooling, and the
radiation force due to lines. We described our numerical methods for
computing, in a position-dependent manner, various radiation field
characteristics, including energy distributions of mean intensity and flux,
the line force, the mean ionizing intensity, and the mean photon energy.

Typically, the energy distributions of mean intensity and flux components
differ, and their position dependence is non-trivial. However, certain
distributions are relatively straightforward to estimate in two limiting
cases. Firstly, at the locations near the disk, the $\theta$ component of the
flux is characterized by a single-temperature blackbody distribution. This
flux originates from a small, almost constant temperature disk area directly
below. Secondly, at locations far from the disk, the mean intensity and radial
flux distributions resemble multi-temperature blackbody spectra originating
from the entire disk. Using our methods, we can recover the dependencies in
these two regions very well (e.g., compare the solid green and dashed black
curves in the middle panel of Fig. 2 and the solid blue and orange curves to
the dashed red in the top panel of that figure). In addition, we find that the
geometric dilution factor well captures the radius dependence of the mean
intensity for $r~{}\gtrsim~{}15~{}r_{\rm in}$ and
$\theta~{}\gtrsim~{}40^{\circ}$.

We also illustrated how the direction of the radiation force can be affected
by the position-dependent mean intensity spectrum. In accordance with detailed
photoionization calculations (see Fig. 1), we made the assumption that the
line opacity’s energy distribution scales like $J_{\nu}$ and showed that the
line force can be more perpendicular to the disk than the radiation flux
vector (see Fig. 3).

The wind’s thermodynamics and observational properties, regardless of its
launching mechanism, are determined by its ionization and radiative heating
and cooling rates. The temperature and thermal stability of highly ionized gas
depend on the mean photon energy, $\left\langle h\nu\right\rangle$. Therefore,
we calculated $\left\langle h\nu\right\rangle$ for different cases and found
that the mean intensity is much softer than the radial flux at many points
near the disk. This softening reduces the severity of the overionization
problem as we discussed in the introduction. We also found that near the disk,
solely geometrical factors can cause the value of $\left\langle
h\nu\right\rangle$ to vary by up to a factor of 3 when the central source’s
X-ray luminosity is zero. However, even a small increase in the X-ray
luminosity can increase $\left\langle h\nu\right\rangle$ by several orders of
magnitude. We note that this dependence on the central luminosity of the
X-rays can be reduced by attenuation caused by a bulging disk or an inner disk
wind. We will asses the impact of such intervening gas and the misalignment of
the line force and flux vectors. by carrying out time-dependent radiation-
hydrodynamical simulations of disk winds. To do so, we will incorporate our
radiation field methods in a radiation-MHD code like Athena++ . Our treatment
in this work assumes a geometric optics formalism where radiation propagates
along linear rays. We neglect special and general relativistic effects. For
example, we do not account for Doppler boosting (both Doppler beaming and the
Doppler shift) due to emission from a relativistic disk or compute
gravitational deflection of photons. Using the treatment presented in Equation
(4) of Jiang (2022) we estimated the percent difference in $J$ and $J_{\nu}$
between including and not including Doppler boosting for the three spectra
shown in Figure 2. For the points $(12r_{*},85^{\circ})$,
$(40r_{*},80^{\circ})$, and $(100r_{*},40^{\circ})$, we find the errors in $J$
to be, respectively, 3%, 0.7%, and 0.6%. For $J_{\nu}$, we found the maximum
errors to be 3%, 0.8%, 1%. These errors become progressively worse as one
approaches $r_{isco}$, particularly near $\theta=90^{\circ}$, where they can
grow as large as 50%. However at such positions that these differences are
large the assumption of linear radiation propagation is violated and so a
wholly more detailed relativistic treatment is necessary to accurately
describe radiation fields near the black hole.

Acknowledgments

Support for this work was provided by the National Aeronautics and Space
Administration under TCAN grant 80NSSC21K0496. D.P. and T.W. acknowledge
support from NASA grant HST-GO- 16196.

## References

  * Castor et al. (1975) Castor, J. I., Abbott, D. C., & Klein, R. I. 1975, ApJ, 195, 157, doi: 10.1086/153315
  * Dannen et al. (2019) Dannen, R. C., Proga, D., Kallman, T. R., & Waters, T. 2019, ApJ, 882, 99, doi: 10.3847/1538-4357/ab340b
  * de Kool & Begelman (1995) de Kool, M., & Begelman, M. C. 1995, ApJ, 455, 448, doi: 10.1086/176594
  * Dyda et al. (2017) Dyda, S., Dannen, R., Waters, T., & Proga, D. 2017, MNRAS, 467, 4161, doi: 10.1093/mnras/stx406
  * Dyda et al. (2023) Dyda, S., Davis, S. W., & Proga, D. 2023, arXiv e-prints, arXiv:2310.18557, doi: 10.48550/arXiv.2310.18557
  * Fabian (2012) Fabian, A. C. 2012, ARA&A, 50, 455, doi: 10.1146/annurev-astro-081811-125521
  * Frank et al. (2002) Frank, J., King, A., & Raine, D. J. 2002, Accretion Power in Astrophysics: Third Edition
  * Friend & Abbott (1986) Friend, D. B., & Abbott, D. C. 1986, ApJ, 311, 701, doi: 10.1086/164809
  * Giustini & Proga (2019) Giustini, M., & Proga, D. 2019, A&A, 630, A94, doi: 10.1051/0004-6361/201833810
  * Higginbottom et al. (2014) Higginbottom, N., Proga, D., Knigge, C., et al. 2014, ApJ, 789, 19, doi: 10.1088/0004-637X/789/1/19
  * Jiang (2022) Jiang, Y.-F. 2022, ApJS, 263, 4, doi: 10.3847/1538-4365/ac9231
  * Kallman & Bautista (2001) Kallman, T., & Bautista, M. 2001, ApJS, 133, 221, doi: 10.1086/319184
  * Kallman & McCray (1982) Kallman, T. R., & McCray, R. 1982, ApJS, 50, 263, doi: 10.1086/190828
  * Krolik (1999) Krolik, J. H. 1999, Active galactic nuclei : from the central black hole to the galactic environment
  * Laha et al. (2020) Laha, S., Reynolds, C. S., Reeves, J., et al. 2020, Nature Astronomy, doi: 10.1038/s41550-020-01255-2
  * Lamers & Cassinelli (1999) Lamers, H. J. G. L. M., & Cassinelli, J. P. 1999, Introduction to Stellar Winds
  * Murray et al. (1995) Murray, N., Chiang, J., Grossman, S. A., & Voit, G. M. 1995, ApJ, 451, 498, doi: 10.1086/176238
  * Nomura et al. (2020) Nomura, M., Ohsuga, K., & Done, C. 2020, MNRAS, 494, 3616, doi: 10.1093/mnras/staa948
  * Owocki et al. (1988) Owocki, S. P., Castor, J. I., & Rybicki, G. B. 1988, ApJ, 335, 914, doi: 10.1086/166977
  * Owocki & Puls (1996) Owocki, S. P., & Puls, J. 1996, ApJ, 462, 894, doi: 10.1086/177203
  * Pauldrach et al. (1986) Pauldrach, A., Puls, J., & Kudritzki, R. P. 1986, A&A, 164, 86
  * Press et al. (2007) Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. 2007, Numerical Recipes 3rd Edition: The Art of Scientific Computing, 3rd edn. (USA: Cambridge University Press)
  * Proga (2007) Proga, D. 2007, in Astronomical Society of the Pacific Conference Series, Vol. 373, The Central Engine of Active Galactic Nuclei, ed. L. C. Ho & J. W. Wang, 267
  * Proga & Kallman (2004) Proga, D., & Kallman, T. R. 2004, ApJ, 616, 688, doi: 10.1086/425117
  * Proga et al. (1998) Proga, D., Stone, J. M., & Drew, J. E. 1998, MNRAS, 295, 595, doi: 10.1046/j.1365-8711.1998.01337.x
  * Proga et al. (1999) —. 1999, MNRAS, 310, 476, doi: 10.1046/j.1365-8711.1999.02935.x
  * Proga et al. (2000) Proga, D., Stone, J. M., & Kallman, T. R. 2000, ApJ, 543, 686, doi: 10.1086/317154
  * Proga et al. (2022) Proga, D., Waters, T., Dyda, S., & Zhu, Z. 2022, ApJ, 935, L37, doi: 10.3847/2041-8213/ac87b0
  * Rybicki & Hummer (1978) Rybicki, G. B., & Hummer, D. G. 1978, ApJ, 219, 654, doi: 10.1086/155826
  * Shakura & Sunyaev (1973) Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 500, 33
  * Sim et al. (2010) Sim, S. A., Proga, D., Miller, L., Long, K. S., & Turner, T. J. 2010, MNRAS, 408, 1396, doi: 10.1111/j.1365-2966.2010.17215.x
  * Stevens & Kallman (1990) Stevens, I. R., & Kallman, T. R. 1990, ApJ, 365, 321, doi: 10.1086/169486
  * Tarter et al. (1969) Tarter, C. B., Tucker, W. H., & Salpeter, E. E. 1969, ApJ, 156, 943, doi: 10.1086/150026