# Imaging thick accretion disks and jets surrounding black holes

Zhenyu Zhang1 Yehui Hou1 Minyong Guo2,3∗ Bin Chen1,4,5

###### Abstract

Based on the horizon-scale magnetofluid model developed in [1], we investigate
the millimeter-wave images of a geometrically thick accretion disk or a funnel
wall jet around a Kerr black hole. By employing the numerical method to solve
the null geodesic and radiative transfer equations, we obtain the optical
appearances at various observational angles and frequencies, generated by the
thermal synchrotron radiation within the magnetofluid. For the thick disk, we
specifically examine the impact of emission anisotropy on images, concluding
that anisotropic synchrotron radiation could play an important role in the
observability of the photon ring. For the funnel wall jet, we find that both
the outflow and inflow funnel walls exhibit annular structures on the imaging
plane. The outflow jet yields a brighter primary image than the photon ring,
whereas the inflow jet does not. Based on our investigation, the inflow funnel
wall model can not be ruled out by current observations of M87*.

1Department of Physics, Peking University, No.5 Yiheyuan Rd, Beijing 100871,
P.R. China

2 Department of Physics, Beijing Normal University, Beijing 100875, P. R.
China

3 Key Laboratory of Multiscale Spin Physics, Ministry of Education, Beijing
100875, P. R. China

4Center for High Energy Physics, Peking University, No.5 Yiheyuan Rd, Beijing
100871, P. R. China

5 Collaborative Innovation Center of Quantum Matter, No.5 Yiheyuan Rd, Beijing
100871, P. R. China

$\ast$ Corresponding author<EMAIL_ADDRESS>
##  1 Introduction

In recent years, the horizon-scale photographs of supermassive black holes
within the millimeter wave band have been viewed with great anticipation for
exploring black hole physics [2, 3, 4]. It is believed that the observed
bright rings are due to the so-called photon ring formed by multi-level images
caused by the strong gravitational lensing effect [5, 6, 7, 8, 9]. However,
since the observations are sensitive to the astronomical environments around
black holes, the extent to which the currently photographed bright rings are
correlated with the photon rings is still a controversial topic [10, 11, 12,
13]. Therefore, it is essential to investigate more deeply the physical
properties of the astronomical environments and their contributions to the
black hole images.

It is widely recognized that a supermassive black hole accretes hot,
magnetized plasma, thereby forming a luminous accretion disk. The light
emitted by the thermal synchrotron electrons within this disk is believed to
be the most probable source in the frequency band of black hole image
observations. Moreover, for a high-spin black hole, the electromagnetic energy
can potentially power a relativistic jet [14]. The base of the jet is
enveloped by what is known as the “funnel wall”, capable of generating
substantial thermal synchrotron radiation. Consequently, both the disk and the
jet can serve as light sources, producing horizon-scale images of black holes.

In the field of theoretical research, many previous studies mainly focus on
the gravitational lensing effect stemming from spacetime geometry, where the
light sources are modeled to be of simple structures and often as
geometrically thin disks [15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27,
28, 29, 30, 31, 32]. However, for a low-luminosity active galactic nuclei
(LLAGN) like the M87*, the horizon-scale accretion flow is mostly
geometrically thick [33, 2]. In this case, it is essential to have precise
knowledge of various factors such as the streamlines, particle number density,
electron temperature, and the magnetic field structure. This crucial physical
information can be obtained by solving the general relativistic
magnetohydrodynamic (GRMHD) equations. Given the analytical challenges posed
by the GRMHD equations, researchers typically resort to numerical techniques
to ascertain the physical properties of the magnetofluid [34, 35, 36, 37, 38,
39, 40, 41]. However, there are shortcomings in numerical simulation, for
example, the relationship between the initial and final states of simulation
is somewhat ambiguous, and the numerical approach is computationally
demanding.

To balance computationally expensive numerical simulations and over-simplified
emission models, semi-analytic models have been utilized in some literature.
The basic idea behind these models is to circumvent the GRMHD simulation and
give the distributions of number density, temperature, and magnetic field in
other ways. For example, most of these works considered radiatively
inefficient accretion flow (RIAF) models [42, 43, 44, 45, 46, 47, 48, 49, 50,
51], where the vertically-averaged electron density and temperature have a
nearly power-law dependence on radius [52]. Another model that was often
utilized is the magnetized torus solution [53]. This model, characterized by
its purely toroidal configuration, was frequently adopted as the initial
condition for GRMHD simulations. Interestingly, some studies have also
considered the images of the torus [54, 55]. After establishing the emission
profiles for the accreted plasma, black hole images can be obtained by using
general relativistic radiative transfer (GRRT) methods.

Recently, we introduced a new analytical model of magnetofluids surrounding a
rotating black hole within the framework of GRMHD[1]. Assuming that the fluid
acceleration at the event horizon scale is dominated by gravitation, we
derived a comprehensive set of explicit expressions for thermodynamics and the
associated magnetic field. We showed that our model can effectively describe
the morphology of both thick disks and jets at the horizon scale and provides
a direct pathway for the study of black hole imaging. In the present work,
given that our model provides the essential physical parameters for the
emission profile of the magnetofluid within the GRMHD framework, we would like
to utilize the model to study the optical appearances of accretion disks and
funnel wall jets near black holes, by using a GRRT method to generate the
images numerically. Our investigation for the images of the disk will not be
limited to analyzing the characteristics of intensity maps under various
observational angles and frequencies, but will specifically probe the
influence of emission anisotropy on the images.

The remaining parts of the paper are organized as follows. In sec. 2, we
concisely review the analytical model for thick accretion disks and jets,
developed in [1]. In sec. 3 we introduce thermal synchrotron radiation and
radiative transport, which are fundamental in magnetohydrodynamic imaging. In
Sec. 4, we present the imaging results for accretion disks and jets,
accompanied by comprehensive discussions. Sec. 5 provides a summary and
discussion of our work. In addition, we will work in the unit with
$GM=c=4\pi\epsilon_{0}=1$, where $M$ is the mass of the black hole.

##  2 Analytic model for thick accretion disks and jets

Recently, we have proposed a new analytical model for the magnetohydrodynamic
fluid in [1], trying to capture key features of thick disks and FWJs.
Specifically, we considered the fluid to be electrically neutral and postulate
complete ionization of the plasma into electrons and protons. We focused on
high-temperature systems with ultra-relativistic electrons in the adiabatic
limit. To simplify the analysis, we introduced a linear relationship between
the temperatures of electrons and ions, represented by a constant factor,
$z=T_{\text{ion}}/T_{\text{e}}$, and considered both the non-relativistic and
ultra-relativistic limits of ions. With the relativistic Euler equations for
ideal, adiabatic magnetofluids, we expressed the temperature and pressure as
power functions of the particle number density, with the exponents determined
by $z$. Then, we demonstrated that in the scenarios where both the sound speed
and Alfv́en velocity are sub-relativistic, the streamlines follow geodesics in
the leading-order approximation, namely, the ballistic approximation.
Furthermore, for a stationary, axisymmetric, ballistic fluid configuration in
Kerr spacetime, we derived explicit expressions for the fluid thermodynamics
under specific conditions regarding the streamlines. Especially, one of the
conditions leads to the conical solution, which is of particular importance in
the study of thick disks and FWJs. For the detailed discussions on the model,
please refer to [1].

In this study, we would like to use this model to investigate the images of
thick accretion disks and funnel wall jets encircling Kerr black holes. In a
Kerr spacetime, the metric reads

$\displaystyle\mathrm{d}s^{2}=-\mathopen{}\mathclose{{}\left(1-\frac{2r}{\Sigma}}\right)\mathrm{d}t^{2}+\frac{\Sigma}{\Delta}\mathrm{d}r^{2}+\Sigma\mathrm{d}\theta^{2}+\mathopen{}\mathclose{{}\left(r^{2}+a^{2}+\frac{2ra^{2}}{\Sigma}\mathrm{sin}^{2}\theta}\right)\mathrm{sin}^{2}\theta\mathrm{d}\phi^{2}-\frac{4ra}{\Sigma}\mathrm{sin}^{2}\theta\mathrm{d}t\mathrm{d}\phi\,,$
(2.1)

in Boyer-Lindquist (BL) coordinates. Note that we have set $GM=1$ where $M$ is
the mass of the black hole. In Eq. (2.1) we have defined

$\displaystyle\Delta=r^{2}-2r+a^{2}\,,\quad\Sigma=r^{2}+a^{2}\cos^{2}\theta\,,$
(2.2)

where $a$ is the spin parameter of the black hole. From $\Delta=0$, the outer
and inner radii of the horizons of a Kerr black hole can be obtained as

$\displaystyle r_{\pm}=1\pm\sqrt{1-a^{2}}\,.$ (2.3)

In the model proposed in [1], the event horizon radius $r_{h}\equiv r_{+}$
serves as the inner boundary of the fluid. When considering the outer boundary
of the fluid, we take it to be at infinity. This means that accretion inflows
can extend from infinity to the event horizon, and the outflows can propagate
indefinitely. Furthermore, the plasma streamlines within the fluid move along
geodesic trajectories. Following the conical solution in [1], let us further
consider the plasma maintaining a constant $\theta$ in its geodesic motion. In
Kerr spacetime, the four-velocity of such a timelike particle can be expressed
as

$\displaystyle
u^{t}=E\bigg{[}1+\frac{2r(r^{2}+a^{2})}{\Delta\Sigma}\bigg{]}-\frac{2arL}{\Delta\Sigma}\,,$
$\displaystyle u^{r}=\sigma_{r}\frac{\sqrt{R}}{\Sigma}\,,\quad
u^{\theta}=0\,,$ $\displaystyle
u^{\phi}=E\frac{2ar}{\Delta\Sigma}+L\frac{\csc^{2}\theta(\Delta-a^{2}\sin^{2}{\theta})}{\Delta\Sigma}\,,$
(2.4)

where

$\displaystyle
L=\pm_{L}a\sqrt{E^{2}-1}\sin^{2}{\theta}\,,\,\,Q=-a^{2}(E^{2}-1)\cos^{4}{\theta}\,,$
(2.5)

and

$\displaystyle
R(r,\theta)=\mathopen{}\mathclose{{}\left(E^{2}-1}\right)r^{4}+2r^{3}+2a^{2}\mathopen{}\mathclose{{}\left(E^{2}-1}\right)\cos^{2}{\theta}r^{2}+2\mathopen{}\mathclose{{}\left[\big{(}aE-L\big{)}^{2}+Q}\right]r-a^{2}Q\,,$
(2.6)

with “$\pm_{L}$” denoting the sign of $L$. Here $E$, $L$ are the conserved
energy and angular momentum per unit mass and $Q$ is the Carter constant per
unit square mass. In addition, $\sigma_{r}=\pm 1$ denotes the direction of
motion, either outward or inward. Note that one can demonstrate that $R\geq 0$
is always true outside the event horizon through simple algebraic
computations. Then, from Eq. (2.5) we can find that $E\geq 1$ must hold if a
particle following a conical geodesic can reach infinity. In such
circumstances, by solving the equation of particle number conservation along
streamlines, the number density of the plasma and electron temperature are
expressed as

$\displaystyle
n_{\text{e}}(r,\theta)=n_{\text{e}}(r_{h},\theta)\sqrt{\frac{R(r_{h},\theta)}{R(r,\theta)}}\,,$
(2.7)

and

$\displaystyle
T_{\text{e}}(r,\theta)=\mathopen{}\mathclose{{}\left\\{\begin{aligned}
&T_{\text{e}}(r_{h},\theta)\bigg{[}\frac{R(r_{h},\theta)}{R(r,\theta)}\bigg{]}^{\frac{(1+z)}{3(2+z)}}\,,\quad&z\ll
z_{c}\\\
&T_{\text{e}}(r_{h},\theta)\bigg{[}\frac{R(r_{h},\theta)}{R(r,\theta)}\bigg{]}^{\frac{1}{6}}\,,\quad&z\gg
z_{c}\\\ \end{aligned}}\right.$ (2.8)

where $n_{\text{e}}(r_{h},\theta)$ and $T_{\text{e}}(r_{h},\theta)$ are the
boundary values of $n_{\text{e}}$ and $T_{\text{e}}$, respectively. Besides,
$z_{c}$ is a characteristic quantity corresponding to
$\Theta_{\text{ion}}\equiv
k_{B}T_{\text{ion}}/m_{\text{ion}}=\Theta_{\text{e}}$ near the horizon. In the
case of $z\ll z_{c}$, $\Theta_{\text{ion}}\ll 1$, indicating that protons are
considered to be non-relativistic particles. Conversely, $z\gg z_{c}$
corresponds to $\Theta_{\text{ion}}\gg 1$, where protons satisfy the
conditions of being ultra-relativistic particles. To facilitate further
analysis, we adopt a Gaussian distribution in the $\theta$ direction for
$n_{\text{e}}(r_{h},\theta)$ and take $T_{\text{e}}(r_{h},\theta)$ to be a
constant. The expressions are as follows:

$\displaystyle
n_{\text{e}}(r_{h},\theta)=n_{i}\exp\bigg{[}-\bigg{(}\frac{\sin{\theta}-\sin{\theta_{J}}}{\sigma}\bigg{)}^{2}\bigg{]}\,,\quad
T_{\text{e}}(r_{h},\theta)=T_{i}\,,$ (2.9)

Here, $\theta_{J}$ represents the mean position in the $\theta$ direction, and
$\sigma$ characterizes the standard deviation of the distribution. From now
on, we set the boundary values of electron temperature and number density as
$T_{i}=10^{11}\,\text{K}$, $n_{i}=10^{12}\text{m}^{-3}$, which is applicable
to the plasma surrounding the M87* [49, 1]. With such a choice of parameters,
we can estimate that $\Theta_{\text{e}}\simeq 15$ and $z_{c}\simeq 122$.
Furthermore, in practical situations, the temperature of protons can hardly
reach the level of being considered highly relativistic particles, and we will
focus solely on the case where $z\ll z_{c}$.

Regarding the electromagnetic field $F_{\mu\nu}$ accompanying the fluid, it
satisfies $F_{\mu\nu}u^{\nu}=0$, which implies the absence of the electric
field in the rest frame of the fluid. As given in [1], the expression for the
magnetic field is as follows:

$\displaystyle
B^{\mu}=\frac{1}{\sqrt{-\text{det}(g_{\mu\nu})}}\frac{F_{\theta\phi}}{u^{r}}\big{(}u_{t}u^{\mu}+\delta_{t}^{\mu}\big{)}\,,\quad\mu=t,r,\theta,\phi\,,$
(2.10)

where $\text{det}(g_{\mu\nu})=-\Sigma^{2}\sin^{2}{\theta}$ is the determinant
of the metric in Kerr spacetime, expressed in BL coordinates, and
$F_{\theta\phi}$ is a function of $\theta$ for conical motions and can be
chosen to be

$\displaystyle
F_{\theta\phi}=\Psi_{0}\,\text{sign}(\cos{\theta})\sin{\theta}\,,$ (2.11)

in subsequent computations, where $\Psi_{0}$ is a constant indicating the
strength of the magnetic field. Furthermore, the other non-zero components of
the corresponding electromagnetic field tensor are

$\displaystyle F_{\phi\theta}=-F_{\theta\phi}\,,\quad F_{r\phi}=-F_{\phi
r}=\frac{u^{\theta}}{u^{r}}F_{\theta\phi},\,\quad F_{r\theta}=-F_{\theta
r}=\frac{u^{\phi}}{u^{r}}F_{\theta\phi}\,.$ (2.12)

##  3 Black hole imaging from thermal synchrotron radiation

In this section, we will briefly revisit the mechanisms of thermal synchrotron
radiation in the fluid medium and introduce the radiative transfer equation as
well as the numerical methods employed in the imaging process.

###  3.1 Thermal synchrotron radiation

As is known to all, electrons in a fluid can emit thermal synchrotron
radiation when they undergo acceleration by electromagnetic fields. It is
helpful for us to start with the general form of the power spectrum of thermal
synchrotron radiations in the fluid rest frames. For a single electron, the
power spectrum is [56, 57]

$P_{\nu}=\sqrt{3}\frac{e^{3}}{m_{\text{e}}}{\mathcal{E}}_{{\text{em}}}F\bigg{(}\frac{\nu}{\nu_{c}}\bigg{)}\,,\quad
F(x)=x\int_{x}^{\infty}K_{5/3}(x^{\prime})dx^{\prime}\,,$ (3.1)

where $e$ is the charge of an electron, $m_{\text{e}}$ represents the mass of
an electron, $\nu$ denotes the radiation frequency, and $\nu_{c}\equiv
3e\gamma^{2}{\mathcal{E}}_{{\text{em}}}/(2\pi m_{\text{e}})$ corresponds to
the characteristic frequency. Here, $\gamma=1/\sqrt{1-v^{2}}$ is the Lorentz
factor and
${\mathcal{E}}_{{\text{em}}}=\sqrt{\mathopen{}\mathclose{{}\left(\vec{E}+\vec{v}\times\vec{B}}\right)^{2}-\mathopen{}\mathclose{{}\left(\vec{v}\cdot\vec{E}}\right)^{2}}$
encodes the curvature of the electron’s orbit caused by the combination of
electric field $\vec{E}$, magnetic field $\vec{B}$, and the velocity vector
$\vec{v}$ of the moving electron. In addition, $K_{5/3}(x^{\prime})$
represents the modified Bessel function, which can be written as

$\displaystyle K_{5/3}(x^{\prime})=\int_{0}^{\infty}e^{-x^{\prime}\cosh
x}\cosh\mathopen{}\mathclose{{}\left(5/3x}\right)dx\,.$ (3.2)

In our study, we focus on relativistic electrons with high temperature,
satisfying

$\displaystyle\Theta_{\text{e}}\equiv\frac{k_{B}T_{\text{e}}}{m_{\text{e}}}\gg
1\,,$ (3.3)

where $k_{B}$ is the Boltzmann constant. In this regime, the radiations
emitted by electrons are predominantly concentrated along the directions of
their motions. Therefore, the emissivity of a single electron can be
approximated as

$\displaystyle J_{\text{e}}\approx
P_{\nu}\delta\mathopen{}\mathclose{{}\left(\Omega-\Omega_{\vec{v}}}\right)\,,$
(3.4)

where $\Omega$ denotes the solid angle along the direction of radiation,
$\Omega_{\vec{v}}$ represents the solid angle along the direction of motion
and $\delta(x)$ refers to the delta function. Taking the thermal distribution
into account, the emissivity of the collection of thermal electrons is given
by

$\displaystyle j_{\nu}=\frac{1}{4\pi}\int d\gamma
d\Omega_{\vec{v}}\frac{dn_{\text{e}}}{d\gamma}J_{e}\approx\frac{\sqrt{3}e^{3}n_{\text{e}}}{8\pi\Theta_{\text{e}}^{3}m_{\text{e}}}{\mathcal{E}}\int_{0}^{\infty}d\gamma\gamma^{2}e^{-\frac{\gamma}{\Theta_{\text{e}}}}F\bigg{(}\frac{4\pi
m_{\text{e}}\nu}{3e\gamma^{2}{\mathcal{E}}}\bigg{)}\,,$ (3.5)

where

$\displaystyle{\mathcal{E}}=\sqrt{\mathopen{}\mathclose{{}\left(\vec{E}+\vec{n}\times\vec{B}}\right)^{2}-\mathopen{}\mathclose{{}\left(\vec{n}\cdot\vec{E}}\right)^{2}}\,,$
(3.6)

with $\vec{n}$ represent the unit vector along the radiation direction, and
$n_{\text{e}}$ denotes the number density of the electrons. Note that we have
used $|\vec{v}|\approx 1$, $\Theta_{\text{e}}\gg 1$ when deriving the relation
(3.5), since the emission is mainly contributed by ultra-relativistic
electrons. With the expansion of $F(x)$, we can obtain two asymptotic
expressions for the integral [58],

$\displaystyle j_{\nu}\approx\begin{cases}\,\frac{2^{4/3}\pi
e^{2}}{3}\frac{n_{\text{e}}\nu}{\Theta_{\text{e}}^{2}}X^{-2/3},&\quad\text{if}\,\,\,\
X\ll 1\,;\\\ \,\\\ \,\frac{2^{-1/2}\pi
e^{2}}{3}\frac{n_{\text{e}}\nu}{\Theta_{\text{e}}^{2}}\exp{\mathopen{}\mathclose{{}\left[-X^{1/3}}\right]},&\quad\text{if}\,\,\,\
X\gg 1\,,\end{cases}$ (3.7)

where

$\displaystyle X\equiv\frac{9\pi
m_{\text{e}}\nu}{e\Theta_{\text{e}}^{2}{\mathcal{E}}}.$ (3.8)

Nevertheless, the expression for $j_{\nu}$ given by Eq. (3.7) is incomplete,
providing only the asymptotic expressions for the cases when $X$ is much
greater than 1 and much less than 1. To obtain the complete expression for
$j_{\nu}$, we employ the approximation method proposed in [58], which involves
connecting the two asymptotic expressions from Eq. (3.7) using a fitting
function. The expression of the approximate fitting function is given by

$\displaystyle\mathrm{j}_{\nu}=\frac{2^{-1/2}\pi
e^{2}}{3}\frac{n_{\text{e}}\nu}{\Theta_{\text{e}}^{2}}\mathopen{}\mathclose{{}\left(1+2^{11/12}X^{-1/3}}\right)^{2}\exp{\mathopen{}\mathclose{{}\left[-X^{1/3}}\right]}\,,$
(3.9)

where we have considered $\Theta_{\text{e}}\gg 1$. It should be noted that the
relative error of $\mathrm{j}_{\nu}$ is within $10\%$ when
$\Theta_{\text{e}}>1.5$, and within $1\%$ when $\Theta_{\text{e}}>5$.

###  3.2 Radiative transfer

Once the emission profile of thermal synchronous radiation at the light source
is determined, in order to obtain the image at the observer’s location, we
need to propagate the radiation information from the light source to the
observer’s screen. Under the geometrical optics approximation, the light moves
along the null geodesic and carries the information of the medium by
accumulating specific intensity. Then, the radiative transport equation of
total intensity $I_{\nu}$ can be expressed as [59]

$\displaystyle\frac{d}{d\lambda}\bigg{(}\frac{I_{\nu}}{\nu^{3}}\bigg{)}=\frac{j_{\nu}-\alpha_{\nu}I_{\nu}}{\nu^{2}}\,,$
(3.10)

where $\lambda$ is the affine parameter of the null geodesic, $\alpha_{\nu}$
is the absorption coefficient when a light ray passes through the fluid.
According to (3.10), when absorption is neglected, the scaling of the
emissivity results in the observed intensity scaling in the same way,
$j_{\nu}\rightarrow\epsilon j_{\nu}$, $I_{\nu_{o}}\rightarrow\epsilon
I_{\nu_{o}}$. However, in the presence of non-zero absorption, such a scale
transformation does not exist. In this work, we take the absorption into
account, and thus, when handling the quantity $j_{\nu}$, it is essential to
retain the constant terms within the expression of $j_{\nu}$. Note that the
absorption coefficient can be determined by Kirchhoff’s law,

$\displaystyle\alpha_{\nu}=\frac{j_{\nu}}{B(T_{\text{e}},\nu)}\,,$ (3.11)

where

$\displaystyle
B(T_{\text{e}},\nu)=\frac{2h\nu^{3}}{\mathrm{e}^{h\nu/k_{B}T_{\text{e}}}-1}\,,$
(3.12)

is the Planck function. We consider the frequency of observed photons, denoted
by $\nu_{o}$, to be less than $400$ GHz, while the temperature of the electron
is estimated to be on the order of $10^{11}$ K. Henceforth, for the photons
away from the horizon by more than $10^{-16}GM/c^{2}$, we find the ratio

$\displaystyle\frac{h\nu}{k_{B}T_{\text{e}}}=\frac{h\nu_{o}}{gk_{B}T_{\text{e}}}\sim\sqrt{-\mathopen{}\mathclose{{}\left(1-\frac{2GM}{rc^{2}}}\right)}\frac{h\nu_{o}}{k_{B}T_{\text{e}}}<0.02\,,$
(3.13)

where $h$ is the Planck’s constant and $g=\nu_{o}/\nu$ is the redshift factor.
Here we employ the static emitter near a Schwarzschild black hole for simple
estimation. This implies that the condition of $h\nu\ll k_{B}T_{\text{e}}$
persistently holds. Thus, the radiative power adheres to the Rayleigh-Jeans
law, which is given by

$\displaystyle\alpha_{\nu}\approx\frac{j_{\nu}}{2\nu^{2}k_{B}T_{\text{e}}}=\frac{j_{\nu}}{k_{B}T_{\text{e}}}\frac{g^{2}}{2\nu_{o}^{2}}\,,$
(3.14)

If the geodesic of light is solved, which we denote as $x^{\mu}(\lambda)$,
then the solution of Eq. (3.10) takes

$I_{\nu_{o}}=\nu_{o}\int^{\lambda_{o}}_{\lambda_{s}}d\lambda\,g^{2}\big{(}x^{\mu}(\lambda),k^{\mu}(\lambda)\big{)}\,j_{\nu}\big{(}x^{\mu}(\lambda),k^{\mu}(\lambda)\big{)}\,\exp\bigg{\\{}-\int^{\lambda_{o}}_{\lambda}d\lambda^{\prime}\nu\alpha_{\nu}(x^{\mu}(\lambda^{\prime}),k^{\mu}(\lambda^{\prime}))\bigg{\\}}\,,$
(3.15)

where $\lambda_{s}$ is the affine parameter of the geodesic integration at the
start point, $\lambda_{o}$ is the affine parameter of the light reaching the
observer and
$k^{\mu}(\lambda)=\dot{x}^{\mu}(\lambda)=dx^{\mu}(\lambda)/d\lambda$ is the
4-momentum of photons. Eq. (3.15) gives the intensity maps in the observer’s
screen. As indicated by Eq. (3.15), it is necessary to obtain the null
geodesic $x^{\mu}(\lambda)$ in order to establish a connection between the
light source and the screen. In this study, we adopt the ray-tracing method
introduced in the Appendix of [60] to numerically solve for the null
geodesics, where a fisheye camera model is employed based on stereographic
projection on the observer’s screen.

Next, for the convenience of numerical calculations, let us rewrite
${\mathcal{E}}$ in a covariant form. Note that in Eqs. (3.6) the vectors
$\vec{n}$ and $\vec{B}$ are evaluated at the tetrad of the local rest frame of
the fluid $\\{e^{\mu}_{(0)}=u^{\mu},e^{\mu}_{(i)},\,i=1,2,3\\}$. We set
${\mathcal{K}}^{\mu}=k^{\mu}|_{\text{fluid}}$ to be the 4-momentum of photons
at the fluid, and the 3-dimensional unit vector $\vec{n}$ in Eq. (3.6) can be
written as

$\displaystyle
n^{\mu}=\frac{{\mathcal{K}}^{\mu}}{u^{\nu}{\mathcal{K}}_{\nu}}+u^{\mu}\,,$
(3.16)

which satisfies $n^{\mu}n_{\mu}=1$,
$n_{(0)}=n_{\mu}e^{\mu}_{(0)}=n_{\mu}u^{\mu}=0$ and
$n_{(i)}=n_{\mu}e^{\mu}_{(i)}={\mathcal{K}}_{(i)}/{\mathcal{K}}^{(0)}$ in the
local rest frame of the fluid. Note that when calculating $n_{(i)}$, we have
used the expression of the frequency,
$\nu=-{\mathcal{K}}_{\mu}e^{\mu}_{(0)}=-{\mathcal{K}}_{\mu}u^{\mu}={\mathcal{K}}^{(0)}$
in the local rest frame of the fluid. Then, it is convenient to define

$\displaystyle f_{\mu}\equiv F_{\mu\nu}{\mathcal{K}}^{\nu}\,,$ (3.17)

and evaluate it at the tetrad

$\displaystyle
f_{(0)}=F_{(0)(a)}{\mathcal{K}}^{(a)}=-E_{(i)}{\mathcal{K}}^{(i)}\,,\quad
f_{(i)}=F_{(i)(a)}{\mathcal{K}}^{(a)}=E_{(i)}{\mathcal{K}}^{(0)}+\tilde{\varepsilon}_{ijk}B_{(k)}{\mathcal{K}}^{(j)}\,,$
(3.18)

where $\tilde{\varepsilon}_{ijk}$ is the volume element of the three-
dimensional spatial surface orthogonal to $u^{\mu}$. We can also obtain

$\displaystyle
f^{\mu}f_{\mu}=f_{(i)}f^{(i)}-f^{2}_{(0)}=\mathopen{}\mathclose{{}\left[\mathopen{}\mathclose{{}\left({\mathcal{K}}^{(0)}\vec{E}+\vec{{\mathcal{K}}}\times\vec{B}}\right)^{2}-\mathopen{}\mathclose{{}\left(\vec{{\mathcal{K}}}\cdot\vec{E}}\right)^{2}}\right]=\nu^{2}{\mathcal{E}}^{2}\,,$
(3.19)

which is independent of the choice of reference frame. Therefore, after
obtaining $x^{\mu}$ and ${\mathcal{K}}^{\mu}$ through numerical ray-tracing
techniques, we can calculate ${\mathcal{E}}$ by considering Eqs. (2.11),
(2.12), (3.17). Subsequently, by substituting Eqs. (2.7), (2.8), (2.9) and
${\mathcal{E}}$ into Eq. (3.9), we can obtain the approximate emissivity
$\mathrm{j}_{\nu}$. Furthermore, the absorption coefficient $\alpha_{\nu}$ can
be obtained easily from Eq. (3.14).

Recall that we have set $T_{i}=10^{11}\,\text{K}$,
$n_{i}=10^{12}\,\text{m}^{-3}$ as the boundary values of the electron
temperature and number density. The electron temperature corresponds to
$\Theta_{\text{e}}\simeq 15$, which provides a high accuracy of approximate
emissivity given by Eq. (3.9). As for the magnetic field, we would like to set
the magnetic field intensity at $(r_{h},\,\theta_{J})$ as the characteristic
magnetic field strength $B_{i}$, that is, $B_{i}\equiv
B(r_{h},\theta_{J})=\sqrt{B^{\mu}B_{\mu}}\,\big{|}_{r=r_{h},\,\theta=\theta_{J}}$.
Furthermore, in the subsequent numerical calculations, we may take $B_{i}=10$
Gauss, which is roughly the order of magnitude of the magnetic field outside
the M87* [61]. Then, for the conical solution of the fluid, the expression of
emissivity can finally be written as

$\displaystyle\frac{\mathrm{j}_{\nu}}{k_{B}T_{i}}=\mathopen{}\mathclose{{}\left(2.94\times
10^{35}}\right)\bigg{[}\frac{n_{i}}{10^{12}\text{m}^{-3}}\bigg{]}\bigg{[}\frac{\nu_{o}}{230\text{GHz}}\bigg{]}\frac{f(X)}{g}\mathrm{e}^{-\big{(}\frac{\sin{\theta}-\sin{\theta_{J}}}{\sigma}\big{)}^{2}}\bigg{[}\frac{R(r_{h},\theta)}{R(r,\theta)}\bigg{]}^{\frac{2-z}{6(2+z)}}\,,$
(3.20)

at $z\ll z_{c}$, where $\nu=\nu_{o}/g$ and
$f(X)\equiv\big{(}1+2^{11/12}X^{-1/3}\big{)}^{2}\mathrm{e}^{-X^{1/3}}$. The
dimensionless quantity $X$ is now given by

$\displaystyle X=\mathopen{}\mathclose{{}\left(1.30\times
10^{2}}\right)\bigg{[}\frac{\nu_{o}}{230\text{GHz}}\bigg{]}\bigg{[}\frac{10\text{Gauss}}{B_{i}}\bigg{]}\frac{1}{g^{2}}\sqrt{\frac{(E^{2}-1)\Sigma(r_{h},\theta_{J})+2r_{h}}{\Sigma(r_{h},\theta_{J})R(r_{h},\theta_{J})\bar{F}_{\mu\alpha}\bar{F}^{\mu\beta}\bar{{\mathcal{K}}}^{\alpha}\bar{{\mathcal{K}}}_{\beta}}}\bigg{[}\frac{R(r,\theta)}{R(r_{h},\theta)}\bigg{]}^{\frac{2(1+z)}{3(2+z)}}\,$
(3.21)

with the renormalized quantities
$\bar{{\mathcal{K}}}^{\mu}={\mathcal{K}}^{\mu}/\nu_{o}$,
$\bar{F}_{\mu\nu}=F_{\mu\nu}/\Psi_{0}$. Lastly, by utilizing Eq. (3.15), we
can acquire the intensity of the image on the observer’s screen.

In addition, it is worth noting that in the unit with $GM=c=1$, we still have
the freedom to fix the thermal dimension. Since we have
$[I_{\nu}]=[j_{\nu}]=[k_{B}T_{i}]$, the observed intensity can be measured in
the unit of $k_{B}T_{i}$. In the following sections, we set $k_{B}T_{i}=1$ for
simplicity.

##  4 Results of images

In this section, we present the images of the black hole illuminated by the
analytical fluid model, based on the methodology outlined in the preceding
sections.

###  4.1 Images of thick disks

This subsection commences with the optical apparance of thick disks
surrounding black holes within our theoretical framework. The disk exhibits
equatorial symmetry, and the boundary condition of the accretion flow is
chosen as

$\displaystyle\theta_{J}=\frac{\pi}{2}\,,\quad\sigma_{r}=-1\,,\quad E=1\,.$
(4.1)

According to Eq. (2.5), when $E$ equals $1$, $L$ becomes $0$, corresponding to
the infalling geodesics. Note that the azimuthal motion does not vanish due to
the frame-dragging effect. This resembles the second scenario discussed in
[49], which also considers an infalling accretion flow. However, in contrast
to the phenomenological power-law functions of $r$ adopted in [49], our model
derives the distributions of particle number density, temperature, and
magnetic field in a self-consistent manner [1].

The coefficient $\sigma$ in Eq. (2.9) signifies the thickness of the disk, and
we set $\sigma=0.2$ in the ensuing computations. Moreover, as previously
mentioned, we focus on the situation where $z\ll z_{c}\simeq 122$. For the
sake of simplicity, we set 111A prevalent formula extensively employed in
GRMHD is
$T_{{\text{ion}}}/T_{\text{e}}=R_{\text{high}}\beta^{2}/(1+\beta^{2})+1/(1+\beta^{2})$,
where $\beta\equiv p_{\text{gas}}/p_{M}$, and $R_{\text{high}}$ is a
phenomenological parameter ranging from 1 to 160 for the M87*. Since
$p_{M}\approx B^{2}_{i}/2\approx 0.72\times 10^{52}$, $p_{\text{gas}}\approx
n_{i}k_{B}T_{i}(1+z)\approx 10^{51}(1+z)$. Thus, for $z\approx 1$ we have
$\beta^{2}\ll 1$, and if $z>6.2$ we may have $\beta^{2}\gg 1$ and $z\approx
R_{\text{high}}$. Therefore, we can choose $z=1$ for jets and $z=20$ for disks
in our model. $z=20$. Additionally, we have set the characteristic magnetic
field strength $B_{i}=10$ Gauss for our subsequent research, a choice informed
by the order of magnitude of the magnetic field outside the M87*.

Figure 1: Upper row: intensity maps depicting a black hole with $a=0.01$,
illuminated by a thick disk, observed at an inclination of $17^{\circ}$. From
left to right, the observational frequencies are $85$ GHz, $230$ GHz, and
$345$ GHz. Lower row: the horizontal (left) and vertical (right) intensity
cuts across various frequencies. Figure 2: Upper row: Intensity maps depicting
a black hole with a spin parameter $a=0.94$. This black hole is illuminated by
a thick disk and is observed at an inclination of $17^{\circ}$. The maps are
arranged from left to right according to their observation frequencies, which
are $85$ GHz, $230$ GHz, and $345$ GHz respectively. Lower row: the horizontal
(left) and vertical (right) intensity cuts for different frequencies. Figure
3: Upper row: intensity maps depicting a black hole with a spin parameter
$a=0.94$. This black hole is illuminated by a thick disk and is observed at an
inclination of $80^{\circ}$. The maps are arranged from left to right
according to their observation frequencies, which are $85$ GHz, $230$ GHz, and
$345$ GHz respectively. Lower row: the horizontal (left) and vertical (right)
intensity cuts for different frequencies.

In Figs. 1 and 2, we present the imaging results for geometrically thick
accretion disks surrounding a black hole with spin parameters $a=0.01$ and
$a=0.94$ respectively. These images are viewed from an observational angle of
$\theta_{o}=17^{\circ}$. In Fig. 3, we also provide the imaging results at the
observational angle $\theta_{o}=80^{\circ}$ with $a=0.94$. The upper row of
each figure illustrates the light intensity distribution of the black hole
images. These images are arranged from left to right according to their
observation frequencies, which are $85$ GHz, $230$ GHz, and $345$ GHz
respectively. The lower row, also ordered from left to right, depicts the
light intensity distribution along the lines where $y=0$ (horizontal cuts) and
$x=0$ (vertical cuts). Specifically, for the light intensity cuts displayed in
the lower row, we have taken the common logarithm of the light intensity in
order to see clearly the change in light intensity.

Within all the intensity cuts diagrams, it can be discerned that each line is
characterized by two prominent peaks. The regions outside these peaks
correspond to the primary image of the accretion disk, with light from these
regions reaching the observer directly from the disk. Conversely, the peaks
themselves represent higher-order images, where the light has circumnavigated
the black hole once or multiple times prior to reaching the observer. These
higher-order images manifest as conspicuous rings within the intensity maps.
These rings are a consequence of the strong gravitational lensing effect,
whereby the gravitational attraction of the black hole warps the path of
light, enabling a portion of it to bypass the black hole and eventually reach
the observer. From the three illustrations, it can be observed that the ring
structure at $85$ GHz is significantly brighter than those at $230$ GHz and
$345$ GHz. The luminosity of the accretion disk image exhibits a pronounced
dependency on the observation frequency; as the observation frequency
increases, the intensity diminishes. This phenomenon is more conspicuous when
the spin is greater, see Fig. 2 for $a=0.94$. However, the position of the
higher-order images remains unchanged across three different frequencies, as
the gravitational lensing effect exerted by the black hole is independent of
the light’s frequency. Moreover, in every image featuring a thick disk, the
ring structure is distinctly visible. The lensing effect enhances the light
intensity near the location of the photon ring by 3-5 orders of magnitude.
This significant enhancement makes the ring structure highly prominent in the
disk model.

Beyond the ring structure, another noteworthy feature is the darkest region
observable from various angles. That region originates from the event horizon
itself. In the case of a geometrically thin disk extending to the horizon, the
primary image of the horizon, also known as the “inner shadow”, has a distinct
outline and might be captured by the EHT [18]. However, for a geometrically
thick disk, the horizon’s outline is obscured by the emission outside the
equatorial plane, making it more challenging to observe. Despite this, an
extremely dark region is still visible to an observer at
$\theta_{o}=17^{\circ}$, as illustrated in Figs. 1 and 2. Furthermore, in Fig.
3 where $\theta_{o}=80^{\circ}$, two distinct black areas emerge, with the
upper one being slightly darker than the lower one, a phenomenon attributed to
the lensing effect.

Figure 4: Upper row: the isotropic intensity maps representing a black hole
with a spin parameter of $a=0.94$. This black hole is illuminated by a thick
disk and observed at an inclination of $17^{\circ}$. The maps are arranged
from left to right according to their observation frequencies, which are $85$
GHz, $230$ GHz, and $345$ GHz respectively. Lower row: the horizontal (left)
and vertical (right) intensity cuts for different frequencies. Figure 5: Upper
row: A comparison of the intensities of $\nu_{o}=85$ GHz generated by both the
averaged and the unaveraged emissivity along the horizontal axis (left) and
the vertical axis (right) at $a=0.94$. Lower row: A comparison of the
intensities of $\nu_{o}=85$ GHz generated by the averaged and the unaveraged
emissivity at $a=0.01$, as well as the anisotropic emissivity as a function of
the pitch angle. The emission position is given by $r=r_{h}+0.5$,
$\theta=\pi/2$. The observation angles for both the upper and lower rows are
$\theta_{o}=17^{\circ}$.

It is particularly important to emphasize that synchrotron radiation is highly
anisotropic, which means that the radiative power varies with the direction of
emission, represented by the pitch angle $\alpha$ which represents the angle
between the wave vector and the magnetic field vector in the fluid rest
frames. Specifically, we can rewrite Eq.(3.6) as
${\mathcal{E}}=\big{|}\vec{n}\times\vec{B}\big{|}=\big{|}\vec{B}\big{|}\sin{\alpha}$
in the ideal MHD fluid. To study the influence of such anisotropy on the black
hole images, it is necessary to consider an isotropic emission profile as a
comparison. Here, we employ an angle-averaged emissivity defined as

$\bar{j}_{\nu}=\frac{1}{4\pi}\int\mathrm{j}_{\nu}d\Omega=\int^{\pi/2}_{0}\mathrm{j}_{\nu}\sin{\alpha}d\alpha\,,$
(4.2)

where $d\Omega=\sin\alpha d\alpha d\zeta$ is the solid angle element and
$\zeta$ denotes the azimuthal angle in the fluid rest frames. We aim to
compare the results generated by Eq. (3.9) and Eq. (4.2). In Fig. 4, we
present the images generated by $\bar{j}_{\nu}$ , where the other parameters
are kept consistent with those in Fig. 2. The distinctions are noticeably
clear. The most conspicuous difference is that at the same frequency, the
brightness of the primary image in Fig. 4 is greater than that in Fig. 2. This
difference can be more clearly observed from the intensity cuts. Regardless of
whether the intensity cuts are along the horizontal or vertical direction, the
changing trend in the primary image of Fig. 2 is obviously steeper than the
one in Fig. 4.

In Fig. 5, we present a comparative chart of intensity cuts from isotropic and
anisotropic emissions observed at $\theta_{o}=17^{\circ}$. To see the details
of the peaks clearly, we restrict to the sections of the intensity cuts with
either $x<0$ or $y<0$. The first row provides a comparison of the intensities
at $\nu_{o}=85$ GHz along the horizontal axis (left) and the vertical axis
(right) at $a=0.94$. The second row presents the case for $a=0.01$. The left
panel displays the comparison of intensities at $\nu_{o}=85$ GHz, generated by
both the averaged and unaveraged emission along the horizontal direction.
Given that the black hole is almost non-rotating at $a=0.01$, the horizontal
and vertical directions are nearly identical, so we only plot the horizontal
intensity cut. Additionally, the right panel illustrates the emissivity
$\mathrm{j}_{\nu}$ as a function of the pitch angle $\alpha$. The emissivity
peaks when $\alpha=\pi/2$, where the radiation direction is perpendicular to
the magnetic field. This fact remains valid even in the high-spin case, as the
emissivity is insensitive to the spin of the black hole. Notably, the
emissivity exhibits a significant increase as the frequency decreases.

According to Fig. 5, it becomes apparent that the primary image generated by
the isotropic emission is brighter within the peaks compared to the
anisotropic scenario. Such a distinction is closely related to the magnetic
field configuration. To see this point more clearly, let us consider a simple
case of a Schwarzschild black hole observed almost face-on, with the disk
residing on the equatorial plane, and a purely radial magnetic field described
by Eq. (2.10). Gravitational lensing imparts a significant radial component to
the three-momentum of the photon emitted from the disk near the horizon and
reaching the observer. Consequently, the emission direction and magnetic field
line maintain a small angle in the emission region. This results in the
suppression of synchrotron emission in the presence of a radial magnetic
field, as shown in the bottom left panel of Fig. 5. This suppression becomes
especially pronounced near the horizon, where the lensing effect intensifies.
In contrast, the photons corresponding to the lensed images form the photon
sphere at a constant radius. Their three-momenta are nearly perpendicular to
the magnetic field in the emission region such that the intensity of lensed
images is less suppressed compared to that of primary images.

Figure 6: The first and third rows represent the intensity maps of a black
hole with $a=0.01$ and $a=0.94$ respectively, illuminated by an outflow jet
and observed at an angle of $\theta_{o}=17^{\circ}$. From left to right, the
observation frequencies are $85$ GHz, $230$ GHz and $345$ GHz. The second and
fourth rows correspond to the horizontal (left) and vertical (right) intensity
cuts at different frequencies for the situations depicted in the first and
third rows, respectively.

###  4.2 Images produced by a funnel wall jet

In this subsection, we explore the features of black hole images illuminated
by the funnel wall jet in our model. The shape of the funnel wall is obtained
by setting

$\displaystyle\theta_{J}=\frac{\pi}{4}\,,\quad\sigma=0.05\,,$ (4.3)

with $\theta_{J}$ denotes the open angle of the funnel, $\sigma$ is small and
represents a geometrically thin wall. In the funnel wall region, we opt for
$z=1$, a choice validated by the GRMHD simulations in [36]. Furthermore, since
both inward and outward flows are observed in the simulations, we consider
both scenarios, represented by $\sigma_{r}=\pm$ respectively. In the case of
the conical model, there remains a free parameter $E$ denoting the velocity of
the flow. In this study, we select $E=1$ to represent a mild flow. As for the
magnetic field, we persist in choosing a characteristic strength of $B_{i}=10$
Gauss.

Figure 7: The first two rows represent the intensity maps of a black hole with
$a=0.01$ and $a=0.94$ respectively, illuminated by an outflow jet and observed
at an angle of $\theta_{o}=90^{\circ}$. From left to right, the observation
frequencies are $85$ GHz, $230$ GHz and $345$ GHz. The left and right images
in the third row respectively represent vertical intensity cuts at different
frequencies for the scenarios illustrated in the first two rows. Figure 8:
Upper row: intensity maps of a Kerr black hole with $a=0.94$ illuminated by an
inflow funnel wall jet observed at $\theta_{o}=17^{\circ}$. From left to
right, the observation frequency are $85$ GHz, $230$ GHz, $345$ GHz. Bottom
row: the horizontal (Left) and vertical (Right) intensity cuts of different
frequencies. Figure 9:  Upper row: intensity maps of a Kerr black hole with
$a=0.94$ illuminated by an inflow funnel wall jet observed at
$\theta_{o}=90^{\circ}$. From left to right, the observation frequency are
$85$ GHz, $230$ GHz, $345$ GHz. Bottom row: the vertical intensity cuts of
different frequencies. We neglect the horizontal intensity cuts as they are
almost invisible.

The imaging results of the outflow are showcased in Figs. 6 and 7,
representing the images at $\theta_{o}=17^{\circ}$ and $\theta_{o}=90^{\circ}$
respectively. From the figures, we notice a prominent feature: the high
brightness of the primary image of the jet base, attributed to the significant
Doppler effect associated with the outflow. In this discussion, we use the
term “ jet base ” to refer to the inner boundary of the jet, that is, the
intersection of the funnel wall and the event horizon. When observed at
$\theta_{o}=17^{\circ}$, the primary image of the jet base appears like a
donut because the observer is situated within the jet cone,
$\theta_{o}<\theta_{J}$. However, at $\theta_{o}=90^{\circ}$, since the
observer is outside the jet cone and positioned on the equatorial plane, two
symmetrically separated images of the jet base can be seen. Both of these
images exhibit increasing brightness and thickness with the increase in
frequency. Besides, at $\theta_{o}=17^{\circ}$, we notice a series of annular
structures. This phenomenon can be attributed to two factors. Firstly, the
funnel wall manifests as a biconical structure. Secondly, the gravitational
lensing effect results in a series of lensed images of the biconical
structure. Similar observations can be made when the observational angle is
$\theta_{o}=90^{\circ}$. At this particular angle, symmetrical lensed images
become visible.

Figure 10: Optical depth as a function of $t$ along a light ray that crosses
the equatorial plane three times. The source is set as a outflow funnel wall,
with $a=0.94$ and $\theta_{o}=17^{\circ}$. The step increase in the optical
depth curves indicates that the light traverses the jet region multiple times.

However, if the emissions near the M87* observed by the EHT originate from the
funnel wall jet, the jet base should not exhibit such high brightness as
depicted in Fig. 6. This suggests that the outflow model might not be
applicable to M87*. Meanwhile, GRMHD simulations indicate that the funnel wall
can manifest an inward flow [39]. In order to closely match the observed
images, we have considered the imaging results of an inflow funnel wall in
Figs. 8 and 9. Both figures present results for $a=0.94$. The former is
observed at $\theta=17^{\circ}$, while the latter is viewed at
$\theta=90^{\circ}$.

In Figure 8, the primary and lensed images of jet base are clearly visible,
along with the distinctive structure of the photon ring, which exhibits the
highest brightness. Interestingly, the brightness first decreases as we move
inward from the photon ring, then increases, before ultimately drops
significantly, a pattern attributable to the primary image of the jet base. In
Figure 9, where the observer is located on the equatorial plane, the full
photon ring is obscured, yet a segment of the ring remains visible. These two
figures suggest that the brightness of the primary image of the jet base is
significantly diminished compared to the photon ring due to the inward flow.
Consequently, the emissions from the inflow funnel wall are not ruled out by
EHT’s observation. The ring-like image in M87* may contain the important
feature of the jet bases in addition to the photon ring [62]. Our imaging
results are also in alignment with the dual cone model presented in [63].

Furthermore, in Figs. 6 and 7, we notice that the intensity cuts of 85 GHz
appear less inclined compared to those at 230 GHz and 345 GHz. As manifested
in the images themselves, the brightness distribution at 85 GHz seems more
smoothed out. This phenomenon may be attributed to the effect of absorption,
prompting us to investigate the optical depth along a light ray. The results
are depicted in Fig. 10, where we concentrate on a specific light ray that
starts at $r=15$, initially heads towards the black hole, then deflects and
escapes to the observer. The data reveals that the influence of optical depth
on radiative transfer significantly varies with frequency. To be specific, the
absorption effect is negligible at $230$ GHz, but at lower frequencies such as
$85$ GHz or $75$ GHz, it begins to play a significant role. Concurrently, the
emissivity is frequency-dependent as well, with a tendency to increase at
lower frequencies. As a result, despite heightened absorption, the observed
intensity at lower frequencies remains brighter. It is essential to recognize
that the above finding hinges on the positioning of the peaks within the
thermal synchrotron radiation spectrum, which in turn depends on the
parameters in the model, i.e., $n_{i}$, $T_{i}$, and $B_{i}$. Therefore, it is
necessary to determine the physical parameters more accurately in our model
based on both the existing and future observational results of M87*.

##  5 Summary and discussion

In this study, we continued our recent work [1] and utilized the analytical
results of the conical solution for the magnetofluid to investigate black hole
imaging from a thick disk and a funnel wall jet. In our study, the spatial
distributions of number density and temperature are explicitly provided by Eq.
(2.7) and Eq. (2.8), respectively, and the magnetic field structure,
exhibiting a highly radial pattern for the conical solution, is described by
Eq. (2.10). We considered the synchrotron radiation of thermal electrons in
the magnetofluid, which was subsequently transferred to a distant observer to
generate the images. In practice, we numerically solved the null geodesic and
radiation transfer equations to obtain the images.

We investigated thick disk images at varying observational angles and
frequencies. Particularly, to examine the impact of emission anisotropy, we
compared the imaging outcomes from normal synchrotron radiation with those
from isotropic (angle-averaged) radiation. We noticed that emission anisotropy
produces more pronounced photon rings at an observational angle of
$\theta_{o}=17^{\circ}$, and in contrast, isotropic emission yields more
distinct primary images. This suggests that for synchrotron radiations within
a radial magnetic field surrounding a black hole, the photon ring structure is
more readily observable for an observer nearly face-on. This insight could be
validated by advancements in the next-generation EHT [64].

Furthermore, we explored the images of a funnel wall jet, considering both
outward and inward flows. Our study shows that the outflow models yield
excessively bright primary images of jet bases, potentially resulting in the
absence of the black hole shadow after blurring, inconsistent with current
observations of M87*. On the contrary, the images illuminated by the emissions
from the inward funnel wall jet are still consistent with the observations.

Building upon our findings, we further affirm the validity of our analytical
model in characterizing the morphology of thick disks and funnel wall jets at
the horizon scale. However, the model still exhibits certain limitations.
First, the model presumes that the magnetofluid is governed by gravity at the
horizon scale, which restricts the permissible parameter space. For instance,
to fully comprehend the dynamics of funnel wall flows, it is crucial to
integrate viscosity into the relativistic Euler equation, a task that presents
considerable challenges. Second, our analytical results are based on a
stationary and axisymmetric fluid configuration, which may not accurately
reflect the myriad of astronomical processes that happen in accretion disks or
funnel wall regions, such as turbulences, magnetic reconnections, and flares.
To procure time-varying results, the application of numerical methods becomes
essential.

## Acknowledgments

The work is partly supported by NSFC Grant No. 12275004, 12205013 and
11873044. MG is also endorsed by ”the Fundamental Research Funds for the
Central Universities” with Grant No. 2021NTST13.

## References

  * [1] Y. Hou, Z. Zhang, M. Guo, and B. Chen, “A new analytical model of magnetofluids surrounding rotating black holes,” arXiv:2309.13304 [gr-qc].
  * [2] Event Horizon Telescope Collaboration, K. Akiyama et al., “First M87 Event Horizon Telescope Results. I. The Shadow of the Supermassive Black Hole,” Astrophys. J. Lett. 875 (2019) L1, arXiv:1906.11238 [astro-ph.GA].
  * [3] Event Horizon Telescope Collaboration, K. Akiyama et al., “First Sagittarius A* Event Horizon Telescope Results. I. The Shadow of the Supermassive Black Hole in the Center of the Milky Way,” Astrophys. J. Lett. 930 no. 2, (2022) L12.
  * [4] R.-S. Lu et al., “A ring-like accretion structure in M87 connecting its black hole and jet,” Nature 616 no. 7958, (2023) 686–690, arXiv:2304.13252 [astro-ph.HE].
  * [5] E. Teo, “Spherical photon orbits around a kerr black hole,” General Relativity and Gravitation 35 (2003) 1909–1926.
  * [6] H. Falcke, F. Melia, and E. Agol, “Viewing the shadow of the black hole at the galactic center,” Astrophys. J. Lett. 528 (2000) L13, arXiv:astro-ph/9912263.
  * [7] S. E. Gralla and A. Lupsasca, “Lensing by Kerr Black Holes,” Phys. Rev. D 101 no. 4, (2020) 044031, arXiv:1910.12873 [gr-qc].
  * [8] A. E. Broderick et al., “The Photon Ring in M87*,” Astrophys. J. 935 (2022) 61, arXiv:2208.09004 [astro-ph.HE].
  * [9] Y. Hou, P. Liu, M. Guo, H. Yan, and B. Chen, “Multi-level images around Kerr–Newman black holes,” Class. Quant. Grav. 39 no. 19, (2022) 194001, arXiv:2203.02755 [gr-qc].
  * [10] S. E. Gralla, “Can the EHT M87 results be used to test general relativity?,” Phys. Rev. D 103 no. 2, (2021) 024023, arXiv:2010.08557 [astro-ph.HE].
  * [11] Z. Younsi, D. Psaltis, and F. Özel, “Black Hole Images as Tests of General Relativity: Effects of Spacetime Geometry,” Astrophys. J. 942 no. 1, (2023) 47, arXiv:2111.01752 [astro-ph.HE].
  * [12] F. Ozel, D. Psaltis, and Z. Younsi, “Black Hole Images as Tests of General Relativity: Effects of Plasma Physics,” Astrophys. J. 941 no. 1, (2022) 88, arXiv:2111.01123 [astro-ph.HE].
  * [13] K. Glampedakis and G. Pappas, “Can supermassive black hole shadows test the Kerr metric?,” Phys. Rev. D 104 no. 8, (2021) L081503, arXiv:2102.13573 [gr-qc].
  * [14] R. D. Blandford and R. L. Znajek, “Electromagnetic extractions of energy from Kerr black holes,” Mon. Not. Roy. Astron. Soc. 179 (1977) 433–456.
  * [15] J.-P. Luminet, “Image of a spherical black hole with thin accretion disk,” Astronomy and Astrophysics, vol. 75, no. 1-2, May 1979, p. 228-235. 75 (1979) 228–235.
  * [16] J. Peng, M. Guo, and X.-H. Feng, “Influence of quantum correction on black hole shadows, photon rings, and lensing rings,” Chin. Phys. C 45 no. 8, (2021) 085103, arXiv:2008.00657 [gr-qc].
  * [17] X.-X. Zeng, G.-P. Li, and K.-J. He, “The shadows and observational appearance of a noncommutative black hole surrounded by various profiles of accretions,” Nucl. Phys. B 974 (2022) 115639, arXiv:2106.14478 [hep-th].
  * [18] A. Chael, M. D. Johnson, and A. Lupsasca, “Observing the Inner Shadow of a Black Hole: A Direct View of the Event Horizon,” Astrophys. J. 918 no. 1, (2021) 6, arXiv:2106.00683 [astro-ph.HE].
  * [19] M. Guerrero, G. J. Olmo, D. Rubiera-Garcia, and D. S.-C. Gómez, “Shadows and optical appearance of black bounces illuminated by a thin accretion disk,” JCAP 08 (2021) 036, arXiv:2105.15073 [gr-qc].
  * [20] Q. Gan, P. Wang, H. Wu, and H. Yang, “Photon ring and observational appearance of a hairy black hole,” Phys. Rev. D 104 no. 4, (2021) 044049, arXiv:2105.11770 [gr-qc].
  * [21] G.-P. Li and K.-J. He, “Observational appearances of a f(R) global monopole black hole illuminated by various accretions,” Eur. Phys. J. C 81 no. 11, (2021) 1018.
  * [22] G. Gyulchev, P. Nedkova, T. Vetsov, and S. Yazadjiev, “Image of the thin accretion disk around compact objects in the Einstein–Gauss–Bonnet gravity,” Eur. Phys. J. C 81 no. 10, (2021) 885, arXiv:2106.14697 [gr-qc].
  * [23] Y. Hou, Z. Zhang, H. Yan, M. Guo, and B. Chen, “Image of a Kerr-Melvin black hole with a thin accretion disk,” Phys. Rev. D 106 no. 6, (2022) 064058, arXiv:2206.13744 [gr-qc].
  * [24] H.-M. Wang, Z.-C. Lin, and S.-W. Wei, “Optical appearance of Einstein-Æther black hole surrounded by thin disk,” Nucl. Phys. B 985 (2022) 116026, arXiv:2205.13174 [gr-qc].
  * [25] Z. Zhang, S. Chen, and J. Jing, “Image of Bonnor black dihole with a thin accretion disk and its polarization information,” Eur. Phys. J. C 82 no. 9, (2022) 835, arXiv:2205.13696 [gr-qc].
  * [26] Y. Meng, X.-M. Kuang, X.-J. Wang, B. Wang, and J.-P. Wu, “Images from disk and spherical accretions of hairy Schwarzschild black holes,” Phys. Rev. D 108 no. 6, (2023) 064013, arXiv:2306.10459 [gr-qc].
  * [27] X.-J. Wang, X.-M. Kuang, Y. Meng, B. Wang, and J.-P. Wu, “Rings and images of Horndeski hairy black hole illuminated by various thin accretions,” Phys. Rev. D 107 no. 12, (2023) 124052, arXiv:2304.10015 [gr-qc].
  * [28] Y.-X. Huang, S. Guo, Y.-H. Cui, Q.-Q. Jiang, and K. Lin, “Influence of accretion disk on the optical appearance of the Kazakov-Solodukhin black hole,” Phys. Rev. D 107 no. 12, (2023) 123009, arXiv:2311.00302 [gr-qc].
  * [29] C. Zhang, Y. Ma, and J. Yang, “Black hole image encoding quantum gravity information,” Phys. Rev. D 108 no. 10, (2023) 104004, arXiv:2302.02800 [gr-qc].
  * [30] I. De Martino, R. Della Monica, and D. Rubiera-Garcia, “Optical appearance of a nonsingular de Sitter core black hole geometry under several thin disk emissions,” Phys. Rev. D 108 no. 12, (2023) 124054, arXiv:2310.11039 [gr-qc].
  * [31] S. Hu, D. Li, C. Deng, X. Wu, and E. Liang, “Influences of tilted thin accretion disks on the optical appearance of hairy black holes in Horndeski gravity,” arXiv:2309.10557 [gr-qc].
  * [32] L. Chakhchi, H. El Moumni, and K. Masmar, “Shadows and optical appearance of a power-Yang-Mills black hole surrounded by different accretion disk profiles,” Phys. Rev. D 105 no. 6, (2022) 064031.
  * [33] L. C. Ho, “The spectral energy distributions of low-luminosity active galactic nuclei,” Astrophys. J. 516 (1999) 672–682, arXiv:astro-ph/9905012.
  * [34] J. C. McKinney, A. Tchekhovskoy, and R. D. Blandford, “General Relativistic Magnetohydrodynamic Simulations of Magnetically Choked Accretion Flows around Black Holes,” Mon. Not. Roy. Astron. Soc. 423 (2012) 3083, arXiv:1201.4163 [astro-ph.HE].
  * [35] M. Moscibrodzka, H. Falcke, and H. Shiokawa, “General relativistic magnetohydrodynamical simulations of the jet in M 87,” Astron. Astrophys. 586 (2016) A38, arXiv:1510.07243 [astro-ph.HE].
  * [36] Event Horizon Telescope Collaboration, K. Akiyama et al., “First M87 Event Horizon Telescope Results. V. Physical Origin of the Asymmetric Ring,” Astrophys. J. Lett. 875 no. 1, (2019) L5, arXiv:1906.11242 [astro-ph.GA].
  * [37] K. Chatterjee, Z. Younsi, M. Liska, A. Tchekhovskoy, S. B. Markoff, D. Yoon, D. van Eijnatten, C. Hesp, A. Ingram, and M. van der Klis, “Observational signatures of disc and jet misalignment in images of accreting black holes,” Mon. Not. Roy. Astron. Soc. 499 no. 1, (2020) 362–378, arXiv:2002.08386 [astro-ph.GA].
  * [38] S. M. Ressler, C. J. White, E. Quataert, and J. M. Stone, “Ab Initio Horizon-Scale Simulations of Magnetically Arrested Accretion in Sagittarius A* Fed by Stellar Winds,” Astrophys. J. Lett. 896 no. 1, (2020) L6, arXiv:2006.00005 [astro-ph.HE].
  * [39] Y. Shen, Y. Hou, Z.-Y. Fan, M. Guo, and B. Chen, “Imaging a Semi-Analytical Jet model Generated by 3D GRMHD Simulation,” arXiv:2311.14954 [astro-ph.HE].
  * [40] R. Fraga-Encinas, M. Moscibrodzka, and H. Falcke, “The Core-shift of Sagittarius A* as a Discriminant between Disk and Jet Emission Models with millimeter-VLBI,” arXiv:2312.12951 [astro-ph.HE].
  * [41] G. N. Wong and L. Arzamasskiy, “Balanced Turbulence and the Helicity Barrier in Black Hole Accretion,” arXiv:2312.16172 [astro-ph.HE].
  * [42] A. E. Broderick and A. Loeb, “Frequency-dependent shift in the image centroid of the black hole at the galactic center as a test of general relativity,” Astrophys. J. Lett. 636 (2006) L109–L112, arXiv:astro-ph/0508386.
  * [43] A. Broderick and A. Loeb, “Imaging the Black Hole Silhouette of M87: Implications for Jet Formation and Black Hole Spin,” Astrophys. J. 697 (2009) 1164–1179, arXiv:0812.0366 [astro-ph].
  * [44] A. E. Broderick, V. L. Fish, S. S. Doeleman, and A. Loeb, “Evidence for Low Black Hole Spin and Physically Motivated Accretion Models from Millimeter VLBI Observations of Sagittarius A*,” Astrophys. J. 735 (2011) 110, arXiv:1011.2770 [astro-ph.HE].
  * [45] H.-Y. Pu, K. Akiyama, and K. Asada, “The Effects of Accretion Flow Dynamics on the Black Hole Shadow of Sagittarius A${}^{\mathclose{}\,\mathopen{}}$,” Astrophys. J. 831 no. 1, (2016) 4, arXiv:1608.03035 [astro-ph.HE].
  * [46] H.-Y. Pu and A. E. Broderick, “Probing the innermost accretion flow geometry of Sgr A* with Event Horizon Telescope,” Astrophys. J. 863 (2018) 148, arXiv:1807.01817 [astro-ph.HE].
  * [47] F. H. Vincent, M. Wielgus, M. A. Abramowicz, E. Gourgoulhon, J. P. Lasota, T. Paumard, and G. Perrin, “Geometric modeling of M87* as a Kerr black hole or a non-Kerr compact object,” Astron. Astrophys. 646 (2021) A37, arXiv:2002.09226 [gr-qc].
  * [48] Y. Chen, Y. Liu, R.-S. Lu, Y. Mizuno, J. Shu, X. Xue, Q. Yuan, and Y. Zhao, “Stringent axion constraints with Event Horizon Telescope polarimetric measurements of M87*,” Nature Astron. 6 no. 5, (2022) 592–598, arXiv:2105.04572 [hep-ph].
  * [49] F. H. Vincent, S. E. Gralla, A. Lupsasca, and M. Wielgus, “Images and photon ring signatures of thick disks around black holes,” Astron. Astrophys. 667 (2022) A170, arXiv:2206.12066 [astro-ph.HE].
  * [50] H.-X. Jiang, C. Liu, I. K. Dihingia, Y. Mizuno, H. Xu, T. Zhu, and Q. Wu, “Shadows of Loop Quantum Black Holes: Semi-analytical Simulations of Loop Quantum Gravity Effects on Sagittarius A* and M 87*,” arXiv:2312.04288 [gr-qc].
  * [51] V. Deliyski, G. Gyulchev, P. Nedkova, and S. Yazadjiev, “Observing naked singularities by the present and next-generation Event Horizon Telescope,” arXiv:2401.14092 [gr-qc].
  * [52] F. Yuan, E. Quataert, and R. Narayan, “Nonthermal electrons in radiatively inefficient accretion flow models of Sagittarius A*,” Astrophys. J. 598 (2003) 301–312, arXiv:astro-ph/0304125.
  * [53] S. S. Komissarov, “Magnetized Tori around Kerr Black Holes: Analytic Solutions with a Toroidal Magnetic Field,” Mon. Not. Roy. Astron. Soc. 368 (2006) 993–1000, arXiv:astro-ph/0601678.
  * [54] F. H. Vincent, W. Yan, M. A. Abramowicz, A. A. Zdziarski, and O. Straub, “A magnetized torus for modeling Sagittarius A* millimeter images and spectra,” Astron. Astrophys. 574 (2015) A48, arXiv:1406.0353 [astro-ph.GA].
  * [55] F. H. Vincent, M. A. Abramowicz, A. A. Zdziarski, M. Wielgus, T. Paumard, G. Perrin, and O. Straub, “Multi-wavelength torus–jet model for Sagittarius A*,” Astron. Astrophys. 624 (2019) A52, arXiv:1902.01175 [astro-ph.HE].
  * [56] C. A. Coulson, L. Landau, E. Lifshitz, and M. Hamermesh, “The classical theory of fields,” The Mathematical Gazette 38 no. 326, (1954) 312\.
  * [57] G. B. Rybicki and A. P. Lightman, Lightman Radiative Processes in Astrophysics.  Lightman Radiative Processes in Astrophysics, 1979.
  * [58] P. K. Leung, C. F. Gammie, and S. C. Noble, “Numerical calculation of magnetobremsstrahlung emission and absorption coefficients,” The Astrophysical Journal 737 no. 1, (2011) 21.
  * [59] R. W. Lindquist, “Relativistic transport theory,” Annals Phys. 37 no. 3, (1966) 487–518.
  * [60] Z. Hu, Z. Zhong, P.-C. Li, M. Guo, and B. Chen, “QED effect on a black hole shadow,” Phys. Rev. D 103 no. 4, (2021) 044057, arXiv:2012.07022 [gr-qc].
  * [61] Event Horizon Telescope Collaboration, K. Akiyama et al., “First M87 Event Horizon Telescope Results. VIII. Magnetic Field Structure near The Event Horizon,” Astrophys. J. Lett. 910 no. 1, (2021) L13, arXiv:2105.01173 [astro-ph.HE].
  * [62] T. Kawashima, K. Toma, M. Kino, K. Akiyama, M. Nakamura, and K. Moriyama, “A Jet-bases Emission Model of the EHT2017 Image of M87*,” Astrophys. J. 909 no. 2, (2021) 168, arXiv:2009.08641 [astro-ph.HE].
  * [63] E. Papoutsis, M. Bauböck, D. Chang, and C. F. Gammie, “Jets and Rings in Images of Spinning Black Holes,” arXiv:2212.06281 [astro-ph.HE].
  * [64] D. Ayzenberg et al., “Fundamental Physics Opportunities with the Next-Generation Event Horizon Telescope,” arXiv:2312.02130 [astro-ph.HE].