2. Model
The model setup is straightforward: a black hole pulls in and crushes the matter from the surrounding space; the particles are then accelerated to near-light velocities and twisted around the black hole, forming a flattened accretion plasma disk in the equatorial plane.
We will use the spacetime metric entirely characterized by the black hole mass parameter and its “spin” (described in our Refs. [
7,
8]; for more details, see also [
9,
10,
11,
12,
13,
14], and bibliographies therein). We will assume that the mass of the accretion disk is negligible compared to the black hole “mass”
M, probably
; the radiative cooling does not strongly affect the dynamics of fluid motion; and the electrons and ions are very weakly coupled by Coulomb interaction and hence ions and electrons plasmas components have different temperatures,
(see Ref. [
15]), and thus it is the electron component that contributes the most to the equation of state of the accretion disk matter. Due to the large difference in the masses of electrons and protons, electrons are highly mobile and provide quasi-neutrality of the plasma. Due to the high conductivity of the plasma, its own magnetic field can be considered as a field “frozen” into the plasma.
Generally speaking, equations of fluid motion in the vicinity of a black hole must be written using the concept of relativistic dynamics. The key points are as follows: We suppose that the space-time near the (non-charged) black hole Sgr A* is described by the Kerr metric—an exact, singular, stationary, and axially symmetric solution of the Einstein–Hilbert equations of the gravitational “field” in vacuum. Next, we introduce the Boyer–Lindquist 4-coordinates,
(it is well known that besides the Boyer–Lindquist coordinate representation, other representations of space-time locations exist). In terms of the Boyer–Lindquist coordinates, the square of interval is written as
with
, i.e., the components of
depend only on the dimensionless combination
and
. Here,
is the Schwarzschild radius,
c is the speed of light,
G is the gravitational constant, and
M is the “mass” of the black hole. The off-diagonal term
in the metric tensor is proportional to the rate of the black hole’s own rotation and to
. For the Minkowski tensor, we use the metric signature
(see, for example, Ref. [
8], and Refs. therein). To satisfy the principle of causality for moving material objects, obviously,
. The four time-space coordinates
give the location of a world-event from the viewpoint of a remote observer. The meaning of space coordinates
is clear once transitioned to the limit
. When the square of the interval becomes
, i.e., at infinity, parameters
may be interpreted as the standard spherical coordinates in flat space-time. As for the parameter
r, strictly speaking, note that it is not the “distance” in the usual meaning from the center of black hole. This is because, for any material object, in the space-time defined by equation
, no central point
exists in the sense of a world-event on a valid world-line.
Next, we consider the motion of the medium far away from the event horizon, i.e., when parameter r is meaningfully greater than . For the flow at , in the expansion of the metric tensor, we may neglect the terms of order and greater. They contribute less than to the components of the metric tensor. Such omission of smaller terms makes our approximation Newtonian or post-Newtonian. The non-diagonal metric-tensor term (describing involvement of the medium in the rotation of space-time in the vicinity of the black hole) gives rise to the “force” analogous to the traditional Coriolis force in the equations for medium flow.
Hence, we write the system of equations of relativistic fluid dynamics in the curved space-time and expand the metric tensor and the fluid energy–momentum tensor into a series with respect to small parameters and ∼. We keep only the leading terms in the equations of fluid motion.
Next, the fluid is presumed to be localized near surface (in a pancake-like accretion disk)—the flows of the disk medium are considered only near this surface. Then, we can transition to cylindrical coordinates and presume that the gas particles orbit near the plane and the vertical component of their velocity .
In this model, when r exceeds the Schwarzschild radius (specifically, ), then the vertical component of gravity acceleration is balanced by the pressure gradient. Here, G is the gravitational constant, c is the speed of light, , and K∼ is the radius depending so-called Kepler parameter: the angular velocity of a test particle in a circular orbit at distance r from a point mass M in Newtonian approximation of gravity.
Gas pressure along the direction perpendicular to the disk plane is determined by the hydrostatic equilibrium, . Generally speaking, pressure may have important contributions from electromagnetic radiation and induced magnetic field. The simplest case, however, is when the pressure is dominated by gas pressure, and the vertical temperature distribution is isothermal, which is roughly appropriate when the disk is optically thick and externally heated. The equation of state of gas/plasma is then , where s is the “isothermal sound speed” (which is not a function of the transversal coordinate z).
If we further assume that , the equation of hydrostatic equilibrium becomes , the solution of which is . Here, the transversal space scale h is introduced, . Parameter has the meaning of density at the disk mid-plane (at ), and parameter h is the characteristic local thickness of the accretion disk.
In order for the accretion disk to be considered thin, it is necessary that . On the other hand, for great distances away from the black hole (), specific relativistic effects may be neglected or parametrized. Thus, we obtain natural bounds: .
Equations of Motion: The equations of motion for an inviscid incompressible fluid express the laws of conservation for mass, momentum, and energy. The flow of gas may be considered incompressible (see, for example, [
16]) when its velocity
v satisfies condition
and the characteristic time scale of the flow change
where
L is the characteristic spatial scale where the flow characteristics change substantially. In this case, the mass conservation law is expressed as
, i.e., the law is not
, but
.
We assume that the domain of the disk where the vortex structure of interest is formed, is characterized by an approximately constant angular velocity
. In a coordinate system rotating with the angular velocity
, where
and
are the axes in the horizontal plane and
z is vertically upwards, the basic equations of motion become
Here,
are the components of the velocity field,
is density,
P is pressure, and
is temperature. The potential
of the force field is
. Furthermore, here,
is the cross product of vectors
and
; and
is the alternating tensor (
, zero for any two indices being equal,
for any even number of permutations from
, and
for any odd number of permutations).
Due to the assumption of an incompressible medium (implying
), the equation of state
turns into the density expression, which only depends on the temperature (and not on the pressure). In fact,
. The coefficient of thermal expansion,
, may be assumed constant and positive. (For gases,
.) Thus, we can set
, where subscript zero denotes the reference values. Due to the fact that
is generally significantly less then one, one may neglect the density variations in all principal terms and hence replace
with the constant value
, except in the “buoyancy” term, which is proportional to
(see, for example, [
16,
17]).
Next, we apply the curl operator
to the linear momentum conservation Equation (
2). This gives
where brackets symbolize cross-product. Since
, term
is small with respect to
at the horizontal scales typical for any hydrodynamical vector
and may be not taken into account.
Consider now temperature
as
, where
is the (constant) baseline temperature of the accretion disk, quantity
is an axially symmetrical part of the temperature distribution that is not time-dependent, and
is the dynamical quantity related to the vortex structures in the fluid. When
, then
, and Equation (
4) can be rewritten in the tensorial form as
where
is the substantial derivative,
is the flow velocity, and
is the gradient operator with components
in the (
) plane.
In the simplest model—in which the vortex structures are realized in a ”thin” flat sheet of an incompressible inviscid fluid (i.e., when
,
,
)—the
z-component of velocity,
w, vanishes and may be dropped in all formulas (see, for example, [
18]). Note also that
because of an axial symmetry of both
and
. Thus, the set of equations for the
z component of vorticity
and for variation of temperature
becomes
where indices
. Equation (
6) permits the introduction of the stream function
:
where tensor
is the antisymmetric unit tensor of the second order,
, and diagonal components are zero. In this case, vorticity
with
.
Temperature Stratification: To describe the effect of temperature stratification—i.e., to show how stationary temperature increases as the distance from the axis of rotation increases—we express (using
and Taylor series expansion) the background distribution of the temperature
in the form
where
is the parameter characterizing the “rapidity” of increase in temperature with distance from the disk rotation axis. (Generally speaking, the question of what shape the temperature profile takes within a black hole’s accretion disk is an open one. Experimental measurements remain challenging, despite significant progress. See, for example, [
19,
20].) When the leading contribution to potential
comes from the centrifugal effects, we can express
. Thus, when
is presumed constant in a band of
r where the vortex hotspots are forming, the set of coupled nonlinear evolution equations becomes
For a general case,
should be replaced:
.
Equations (
11) and (
12) have a transparent physical meaning. Their left sides describe transport of dynamic quantities: vortex
and temperature perturbation
. Their right sides describe “sources” that generate the vortices and temperature perturbations.
In other words, the vortices are generated by the source (the right part of Equation (
11) which is effective (non-zero) only when there exists an inhomogeneous gravity-like force field
with which temperature perturbation
interacts via the equation of state (when
). The quantity
—which characterizes the vortex field in the fluid—is transported by the self-induced flow (the left part of Equation (
11)). On the other hand, in Equation (
12), this temperature perturbation
is transported by the self-generated flow (
); the temperature perturbation
is generated by the source (the right part of Equation (
12)), which is non-zero only when there exists a spatial and time-independent temperature gradient (i.e., when
). The processes are interlinked because the “source” of one dynamic quantity depends on the complex combination involving another dynamic quantity. The set of Equations (
11) and (
12) shows that when temperature stratification is absent (
) and there is no disk rotation (disk
, i.e., centrifugal force is zero), then Equations (
11) and (
12) degenerate into the traditional equations for vortex evolution and transport in a two-dimensional ideal fluid.
Linear Approximation: Assuming that an excess of temperature
above the basic level of temperature is not too large, in view of the link between fields
and
via Equation (
12), we consider a simple linear dependence between the excess of temperature
and the vorticity
that generates this excess:
Obviously, it follows from Equation (
13) and from the meaning of quantities
and
that the dimension of the coefficient of proportionality is
. By substituting Equation (
13) into Equation (
12), we obtain
Hence, we conclude that both the first and second equations in Equations (
14) and (
15) describe the evolution of the same physical quantity. The equations will be consistent when the following condition is imposed on the current function
:
Here, parameter
. The dimension of this quantity is
; i.e.,
R is a space scale factor. Quantity
describes the distribution of generalized vorticity. Its change in time and in space has to satisfy the evolution equation Equation (
16).
The stream function
is found via the Green function approach. In the symbolic integral form, in the boundless space, it is
The subsequent calculation procedure is as follows: (i) The initial vorticity distribution is set; the stream function
is found from Equation (
17); together with it, the non-linear evolution Equation (
16) is numerically solved. (ii) The distribution of vorticity is set in the form of a macrostructure with petals (with constant vorticity inside) whose moving boundaries evolve according to Equation (
16); i.e., we are considering a region bounded by some closed contour (with possibly a rather complex shape) such that quantity
takes a constant value inside and zero outside. For stationary dynamical regimes, this can be accomplished analytically (details of the contour dynamics method and of the operator techniques can be found, for example, in Refs. [
21,
22]), as well as Refs. [
18,
23]); (iii) For a strongly localized vortex, quantity
can be parameterized by the function in form
, i.e., the one with the center at coordinate
, which satisfies the equation of transport
, i.e., when the center of vortex moves according to
. (Here, the symbol “dot” signifies the derivative with respect to time).
Stationary Vortex Structures: Stationary vortex structures—the ones rotating with constant angular velocity
—are simpler to consider in a rotating coordinate system where the structures appear immovable; i.e., when rotation direction is co-aligned with
z axis, the derivative with respect to time becomes
. (The procedure is laid out, for example, in Ref. [
24].) Indeed, when
), then the calculation relying on the properties of Jacobeans produces
. Then, Equation (
12) may be rewritten as:
which is satisfied by the ansatz
The explicit expression of function
F can be found from the obvious fact that temperature-driven flow perturbations must vanish when temperature perturbations vanish themselves. The suitable expressions is function
. Then, temperature fluctuations are expressed via the stream function
Combining this with Equation (
11) where
is taken into account, we obtain the second evolution equation: the one that describes rotation of a stationary vortex structure with angular velocity
caused by the self-induced field of hydrodynamical velocity:
Once Equation (
20) is written out, parameter
becomes
. Here,
is the coefficient of thermal expansion. For almost all physically realizable situations,
(the well-known exception is water in the temperature range between 0 and
C). Parameter
—characterizing the “rapidity” of increase in temperature with distance from the disk rotation axis—can be either positive or negative. When
, the periphery of the accretion disk is heated more than the central zone; when
, the central zone of the disk is hotter than the periphery. Thermal length scale
is determined by the background state of the disk in the framework of the model. Obviously, space scale parameter
R characterizes the “rapidity” of the decrease in temperature (and stream function) with distance from the disk rotation axis. Below, we write
, where
is the Heaviside step-function (equal to unit for positive argument and zero for negative argument), and, leaving the old notations,
.
3. Results
We use the model described above to gain insights into the three-spot structure in the EHT image of Sgr A* accretion disk (
Figure 1A). The bright spots are the zones with higher (relative to the base level) temperature
and (since
) with higher vorticity
. In this framework, thermo-hydrodynamical vorticity-field structures may be conceptually analyzed with two limit approaches: the method of localized vortices and the method of thermo-vorticial spots (blobs).
Localized Vortices: The simplest consideration for a three-spot stationary thermo-vortex structure (symmetric with respect to the rotation axis of a pancake-like thin accretion disk) is to analytically treat the hot zones as narrowly localized formations modeled as 2D delta-functions at the limit case. Doing so would allow us to find the required conditions for the existence of the observed phenomenon and the resulting relationships between key characteristics of the structure. (This approach is not limited to only structures with three spots, but can be generalized to structures with any vortex-number
N.) Hence, we write:
where
characterizes the intensity of one localized vortex, and
is the two-dimensional Dirac function in the
-disk plane. The appearance of the factor
in the right side of Equation (
22) is due to the following reasoning: the delta-function is a dimensional function. Its dimension is
for 2D space. Obviously, the argument of the delta-function must be dimensionless. (Indeed, there is no such thing as
of “one inch” or
of “ten gallons”.) Therefore,
. In view of the meaning of the delta function in the right part of Equation (
22) and of the structure of the left part of Equation (
22), the scale-factor
R is the same
R that was introduced above—the characteristic size of the vortex kernel. Although alternatives may exist, following the Occam’s Razor principle—the problem-solving principle of parsimony that “entities should not be multiplied beyond necessity”—we choose the simplest option. Indeed, among the options to use
R or
R multiplied by dimensionless
, we choose the simple
R.
We further note that and may be expressed as components of the complex quantity in the plane (where l is the radius of a circumscribed circle).
When a vortex is modeled by a delta-function, a weak (logarithmic) singularity appears in the distribution of stream function
(which is the solution of Equation (
22)) when
. If
, or a smoother distribution is adopted for
, any singularities in
function then disappear.
The solution to Equation (
22) for unbounded space is
Here,
,
is the Bessel function of order
,
is the modified Bessel (McDonald) function of order
.
Recall that for small values of the argument, function , function behaves logarithmically: . For large values, function , as known, behaves as , function behaves as . This means that for and ; i.e., when , function exponentially quickly tends to zero, and therefore, the stream function , and consequently the local vortex magnitude and temperature excess , all will also tend to zero. This is the consequence of the initial choice to model the vortices via delta-functions. Numerically, for some characteristic points, and . The derivative of with respect to argument is .
By substituting Equation (
22) into Equation (
21), and then taking into account Equation (
23), and by setting the terms with delta function and their derivatives as equal to zero, we find the set of conditions for the existence of the modeled vortex structure:
In this expression, symbol “prime” in the second summation indicates that the effect of vortex self-action is excluded. When notation
is used, then
and
(where
l is the radius of a circumscribed circle).
To satisfy Equation (
24), the coefficients before every singularity must equal zero. Thus, the following condition for the parameters of the vortex configuration arises:
Recall that
is a dimensional function of four dimensional parameters.
Equation (
25) can be rewritten in a compact dimensionless form, where
:
which says that for any specific
, ratio
is uniquely defined.
Figure 3 plots Equation (
26).
For the presented model, the temperature excess
follows from Equation (
20) and Equation (
23):
The characteristic temperature scale is defined thus as
. When
(disk periphery is cooler), both “anticyclones” (temperatures excesses) and “cyclones” (temperature depressions) are possible.
Figure 1B (on the front page of the article) illustrates the temperature distribution for a system with three localized vortices, calculated using Equation (
27), for a case when
(disk periphery is hotter). In
Figure 1B, the unit for the color scale is the characteristic magnitude of temperature excess
. For the
x and
y axes, the unit is
;
is the angular velocity of rotation of the vortex structure as a whole;
is the characteristic angular velocity of rotation of the accretion disk; and
. In the expression for
, the parameter
(characterizing the “rapidity” of the increase in temperature with distance from the disk rotation axis) may be estimated as
. Here,
is the characteristic external radius of the active part of the accretion disk (obviously not equal to
R) is presumed to be greater than the size of the macrostructure, and
may be estimated as the temperature difference across the disk (from the periphery to the center). In the expression for
, the parameter
is the coefficient of thermal expansion, which (for the model whose equation of state is approximated by the equation for an ideal gas) may be written as
, where
is the averaged temperature of the accretion disk.
Thermo-Vorticial Spots (“Blobs”): Another fruitful approach is to use the concept of thermo-vortical spots, for which the vorticity is constant inside domains bounded by movable boundaries (contours) and is zero outside. Then, the problem of vortex evolution becomes the problem of examination of the movement of the contours and determination of their final macro-configuration. Examples of such macro-configurations are depicted in
Figure 4 and
Figure 5 (discussed below).
Within this framework—called the contour dynamics method (CDM) (see Refs. [
18,
22,
23,
24,
25,
26,
27])—the distribution of vorticity
q for a vortex structure may be written as
. Function
is a two-dimensional Heaviside step-function equal to one inside the domain in consideration and zero outside. Parameter
describes the vorticity inside the spot, which is presumed constant. Obviously, the dimension of this parameter is
.
Using CDM (see the physical foundation and the details on how to perform calculations using CDM in Ref. [
22]), we find the general expression for temperature distribution:
Here, the integral is taken along the contour (whose shape is previously found by solving the problem of the spot configuration);
is the complex variable of the
-position, normalized by the characteristic space-scale
, which is linked to the size of the macro-configuration; parameter
is dimensionless; and domain
D in the complex plane
z is a region filled with the vorticity
. The domain is bounded by a closed contour
; this boundary is expressed in the parametric form
; parameter
s is the contour arc length beginning from some initial point; and derivative
is a unit vector that is obviously tangential to the contour
.
Based on the contour-integral Equation (
28), the temperature distribution
is solved (and plotted in
Figure 4B) for the contour boundary (black line) in
Figure 4A. Shaded gray in
Figure 4A is the constant vorticity distribution
.
The key parameters here are
and
, where
is the averaged temperature of the accretion disk,
is the characteristic radius of the accretion disk, and
is the temperature difference across the disk (from the periphery to the center). Thus, the characteristic level of the temperature excess in this limit case (i.e., in the model of thermo-vorticial spots, not in the model of localized thermo-vortices) can be estimated by the expression
The above-mentioned parameters in Equation (
29) are the very parameters that must be measured experimentally to be able to comprehend the phenomena in images such as
Figure 1A.
Equation (
29) is functionally similar to Equation (
27). The difference is in the index of power dependence of the ratio of the vortex intensity and its angular velocity. This is an insignificant difference considering the fact that the two models describe radically different limit cases.
4. Conclusions
Large-scale long-lived vortices are found in many types of hydrodynamic flows. Large vortices in turbulent flows are called coherent structures. They are observed, for example, in the planetary atmospheres and in the oceans. The horizontal scales of the vortices are much greater than the atmospheric or oceanic thicknesses. In the simplest geophysical context, examples of such structures are the Gulf Stream rings, the vortices shed from coastal currents, the cyclones and anti-cyclones, the Antarctic Polar Vortex, etc. A very well-known example is Jupiter’s Great Red Spot: a huge vortex plunged in the equatorial flow, which has persisted for more than three centuries; the presence of intense small-scale turbulence around it does not destroy it. This prompts the question: under what conditions do large-scale vortex structures form?
In traditional 3D space, the turbulent motion is usually considered to be homogeneous and isotropic. Common sense and the laws of thermodynamics show that it is very difficult to extract energy from a fully chaotic system, and only with some additional specific properties of such systems is this possible to realize. Homogeneous, isotropic turbulence, which does not possess any preferred directions or preferred scales, is extremely symmetric, giving birth to large-scale vortices; self-organization seems to be quasi-improbable in this case. It is evident that the breaking of some structural symmetry is one of the necessary conditions for the possibility of self-organization. It becomes clear that turbulent motion of fluid/gas with broken spherical symmetry can be a candidate for a system where self-organization of 2D turbulence into large-scale 2D vortex structures can take place.
Notably, quasi-2D flows in “thin” pancake-like accretion disks (see, e.g., [
28,
29] on the physics of accretion), where components of hydrodynamical flows that are perpendicular to the disk plane are strongly suppressed, can be a place where large-scale vortices can be self-organized. However, another condition is an insignificant influence of a dissipative process on large-scale motions in the system. As is known, the role of dissipation in a typical hydrodynamical process is characterized by the so-called number of Reynolds, which is determined as the ratio of magnitude of the inertial term to the dissipative term in the equation of fluid motion. For “smooth” flows described by the models of classical hydrodynamics, the introduction of such a number is not a problem; however, for flows of plasma, the determination of the predominant mechanism of dissipation is not apparent. However, the three-spot structure in the EHT image of Sgr A* accretion disk (seen in
Figure 1A) is a clear example of self-organization in plasma.
To examine the observed phenomenon in the Sgr A* disk from the perspective of theoretical hydrodynamics, we first considered a simplified thermo-hydrodynamic model that permits analytical consideration. In the model, the vorticity clusters are approximated via delta-functions. As the result, we established the condition of existence of a regular thermo-vorticial structure (Equation (
26); see also
Figure 2). We also spelled out the relationships that should take place between, on the one hand, the parameters that determine the vortex-structure dynamics—each vortex size (
l), its period of proper revolution (
), and temperature excess
in the vortex—and, on the other hand, the accretion disk characteristics (
and
, i.e.,
) and the angular velocity of the entire accretion disk
.
The necessary conditions for the formation of large quasi-stationary symmetric thermo-vorticial structures in the plasma disk are as follows: (1) the accretion disk has to be pancake-like thin (i.e., , where h is the thickness and is the characteristic radius) and rotate with non-zero angular velocity (); and (2) the disk temperature has to decrease towards the center (i.e., parameter ). A multi-spot thermo-vortex structure forms only when key system parameters fall within the ranges captured by the dimensionless relationship: for a three-spot structure, , where . The temperature of the vortices (i.e., the excess over the base level) is also linked to another parameter of the system, vortex intensity of the hot spot: . Because the bright spots are hot, i.e., , this result also means that the rotation directions of the vortices and the entire accretion disk must co-align.
In the framework of localized vortices, the estimate for temperature excess of the bright hot spots is given by Equation (
27), and in the framework of thermo-vorticial spots (with sharp contour boundaries), by Equation (
29).
For quasi-2D flows in thin layers of ideal fluid, another fruitful approach to analysis also exists (for deeper insight, see, Refs. [
18,
22,
23,
24,
25,
26,
27], and applications and references therein). This approach is called the contour dynamics method (CDM). The gist of the method is that the continuous hydrodynamical velocity distribution may be treated as a set of patches with movable boundaries (contours) and constant vorticity inside. As the underpinnings of the CDM show, such vortex-patch approximation correctly grasps the general tendency in the dynamics and evolution of large-scale flows when the large-scale motions are weakly sensitive to a fine structure in the hydrodynamical velocity field. The contours of the vortex structures within the CDM framework are determined by spatially one-dimensional integro-differential nonlinear equations. Equations of the contour dynamics—which describe the self-induced motion of the vorticity-discontinuity boundaries, or “contours”, in an inviscid, incompressible, two-dimensional fluid with piecewise constant vorticity distribution—may be effectively resolved by either numerical or analytical approaches.
Figure 5 illustrates the variability of shapes that a symmetric stationary macro-structure may take, depending on one of its guiding parameters. The shape may range from a weakly deformed circle to a sharply pronounced three-petal “flower”. (Other parameters define how many petals appear: two, three, etc.)
In view of the visual similarity between the hot-spot arrangement in the Sgr A* accretion disk revealed in
Figure 1A and the thermo-vortex-structure plotted in
Figure 5, we conclude that the observed hot spots in the Sgr A* accretion disk are highly likely to be large-scale quasi-2D quasi-stationary
vortices in their nature.
In conclusion, let us emphasize that further improvement in understanding undoubtedly depends on the progress in numerical simulations, for which solid experimental knowledge of parameters defining the processes is paramount. Specifically, as discussed above, these parameters are the background temperature, its contrast, the size of the macro-structure, the number of petals, the angular velocity of rotation of the structure as a whole, the angular velocity of rotation of the accretion disk, gradients of velocity of the collective flows of plasma in the disk, and the characteristic time of existence of this quasi-stationary structure.