ISSN 1573-160X, PHYSICAL OCEANOGRAPHY, Vol. 31, Iss. 2, pp. 149-160 (2024)
Original Russian Text © S. G. Demyshev, 2024, published inMorskoy Gidrofizicheskiy Zhurnal, 2024, Vol. 40, Iss. 2, pp. 165-179
Original article
Finite-Difference Approximation of the Potential Vorticity Equation for a Stratified Incompressible Fluid and an Example of its Application for Modeling the Black Sea Circulation
Part I. Finite-Difference Equation of Potential Vorticity of Ideal Fluid
Purpose. The study is purposed at deriving the discrete equations of absolute and potential vorticity for a three-dimensional stratified incompressible fluid as an exact consequence of the finite-difference equations of sea dynamics in the field of a potential mass force in the adiabatic approximation provided that viscosity and diffusion are absent. The properties of two-dimensional projections of the absolute vorticity equation onto coordinate planes and the three-dimensional potential vorticity equation are analyzed. Methods and results. In order to determine the discrete analogues of absolute and potential vorticity, an additional grid is introduced, where the finite-difference equations for the components both of absolute and potential vorticity are written down. Two-dimensional analogues of the three-dimensional equation of absolute vorticity on the planes (x, y), (y, z) and (x, z) are obtained; they possess the feature of preserving vorticity, energy and enstrophy (square of vorticity). A discrete equation for potential vorticity of a stratified incompressible fluid is derived from the finite-difference system of three-dimensional equations of sea dynamics in the adiabatic approximation at the absence of viscosity and diffusion. Conclusions. In the case of a linear equation of state, the discrete equations of absolute vorticity and potential vorticity which are the exact consequence of finite-difference formulation are obtained. The equation of potential vorticity is of a divergent form, and two-dimensional analogues of the absolute vorticity equation on the planes (x, y), (y, z) and (x, z) have two quadratic invariants that provide preservation of the average wave number.
Keywords: discrete equation, dynamics of sea, kinetic energy, vortex, potential vorticity, Ertel invariant
Acknowledgments: The study was carried out with financial support of the Russian Science Foundation grant 23-27-00141.
For citation: Demyshev, S.G., 2024. Finite-Difference Approximation of the Potential Vorticity Equation for a Stratified Incompressible Fluid and an Example of its Application for Modeling the Black Sea Circulation. Part I. Finite-Difference Equation of Potential Vorticity of Ideal Fluid. Physical Oceanography, 31(2), pp. 149-160.
© 2024, S. G. Demyshev
© 2024, Physical Oceanography
Currently, increasing attention is being paid to the study of vortex flows in seas and oceans based on the analysis of potential vorticity of a stratified fluid, which began with the works [1, 2]. A special role is played by the well-known Ertel's theorem [2], which establishes a connection between the displacements of isopycnal surfaces and the vorticity flow of an ideal fluid. Under the condition of mass force potentiality, incompressibility of fluid and absence of viscosity and diffusion, the water motion is isopycnic in nature, so that a particle of liquid located on the surface of constant density remains on it during movement. A detailed ISSN 1573-160X PHYSICAL OCEANOGRAPHY VOL. 31 ISS. 2 2024 149
S. G. Demyshev
Marine Hydrophysical Institute of RAS, Sevastopol, Russian Federation H demyshev@gmail.com
Abstract
Introduction
The content is available under Creative Commons Attribution-NonCommercial 4.0 International (CC BY-NC 4.0) License
presentation of historical and scientific chronology of potential vorticity (PV) research is given in [3, 4].
Transition to a discrete system of fluid dynamics equations raises the question of correspondence of the resulting difference solution to its continuous analogue and, in particular, satisfaction of PV conservation law in a finite-difference problem. The importance of this approach is that the Ertel potential vorticity equation describes the balance between nonlinear forces in the equations of motion and the density advection equation. Therefore, small errors associated with the inconsistency of discrete equation for potential vorticity with the original finite-difference formulation can lead to incorrect conclusions. In addition, the derivation of equation for PV in a relatively general form (without a quasi-static approximation) provides subsequent studies of any of its simplified versions.
Schemes that ensure satisfaction of a number of conservation laws for one-dimensional shallow water equation are considered in [5]. A new finite-difference scheme has been obtained, where the laws of energy, mass, center of mass, and momentum conservation are observed. Difference analogues, which have a set of linear and quadratic invariants, were studied for two-dimensional equations of flow dynamics. In [6], invariant discretization schemes are constructed for one-dimensional and two-dimensional shallow water equations with periodic boundary conditions. It is indicated that for maintaining invariance a redistribution of grid points is required, i.e., during the process of integrating model equations the grid cannot remain fixed. The resulting schemes conserve mass and momentum, but are not energy efficient.
The works 1 2 are devoted to the development of a new numerical method for spatio-temporal solution for two-dimensional Navier-Stokes equations in the incompressibility approximation. It consists in constructing difference schemes that satisfy conservation laws in an integral formulation, not in a differential form. Therefore, discretization of model equations is carried out on the basis of the balance between energy and mass flows at the boundaries of discrete regions for local and integral variables 1
The work 2 considers an improved version of flow approximation at the boundaries of regions, which provides a systematic and rigorous derivation of conditions for modeling differential equations while maintaining mass and momentum. For two- and three-dimensional Navier-Stokes equations in the incompressibility approximation, their discretization is studied in the formulation of energy - momentum - angular momentum conservation [7]. It is demonstrated that energy is conserved at skew-symmetric linearization of the Navier-Stokes equations, and momentum and angular momentum are conserved at Newtonian linearization. Based on numerical calculations, it is concluded that linearization with two Newton steps at each time step ensures preservation of all three parameters over large time intervals.
1 Scott, J.R. and Chang, S.-C., 1993. A New Flux Conserving Newton's Method Scheme for the Two-Dimensional, Steady Navier-Stokes Equations. NASA Technical Memorandum; 106160. NASA, 50 p. Available at: https://ntrs.nasa.gov/api/citations/19930019437/downloads/19930019437.pdf [Accessed: 29 March 2024].
2 Scott, J.R., 1994. A New Flux-Conserving Numerical Scheme for the Steady, Incompressible Navier-Stokes Equations. NASA Technical Memorandum; 106520. NASA, 49 p. Available at: https://ntrs.nasa.gov/api/citations/19940024968/downloads/19940024968.pdf [Accessed: 29 March 2024].
In [8], a generalization of the Arakawa-Lamb scheme [9] with high-order discretization, but accurate to time approximation, was obtained. It was revealed that the properties of symmetry and the laws of energy and enstrophy conservation are satisfied for periodic problems using high-order summation-by-parts operators.
Significant progress in this direction is due to the results of work [10]. It provides a discrete construction of the Nambu and Poisson brackets, which preserves their antisymmetry property. In turn, this makes it possible to construct numerical models of incompressible fluid flows in two and three dimensions satisfying basic conservation laws. The schemes conserve energy and potential vorticity, as well as potential enstrophy [11], which is considered particularly important since it prevents spurious energy flow into high-wavenumber motions. Based on time-discrete analogues of Nambu brackets, a difference scheme that has two discrete quadratic invariants (energy and potential enstrophy) in time and space is obtained. This demonstrates the possibility of generalizing the Arakawa-Lamb scheme to a time-discrete model without loss of invariants [12].
The present work lies in the field of the abovementioned studies and focuses on derivation of potential vorticity discrete equation for an incompressible stratified ideal fluid as an exact consequence of model difference equations without a quasi-static approximation.
Potential vorticity equation of a stratified incompressible fluid in differential formulation
We are to consider a system of differential equations for a stratified incompressible fluid in a field of potential forces under the conditions of adiabaticity and the absence of diffusion, viscosity and external sources. Then, in the Boussinesq approximation in the Cartesian coordinate system, the fluid motion in the region Q with boundary dQ in the Gromeka-Lamb form is described by the following system of equations:
— + |tU = -— V(P + E) + g (1)
dt p0 p 0
VU = 0, (2)
dT -
— + div(TU) = 0, (3)
dt
^ + div(SU) = 0, (4)
p = f (T, S). (5)
The following notations were introduced: U = (u, v, w) are the components of flow velocity vector along the axes (x, y, z) directed to the east, north, and vertically downward, respectively; g = (0, 0, g) is gravitational acceleration; (T, S, P, p) are seawater temperature, salinity, pressure, and density; p0 = 1 g/cm3 (further we assume pressure and density normalized to P0); f = (0,0, fz) is the Coriolis parameter, where fz = 2a sin^; a is angular velocity of the Earth's rotation; 9 is latitude. Relationship (5) is the equation of state, which in this work is assumed to be linear. In terms of satisfying the conservation laws in a discrete problem, such an assumption is a fundamental simplification, which will be discussed in more detail later.
An absolute vorticity and kinetic energy of motion are introduced into equation (1):
1= (£1 ,£y ,£z),
^ = dw_ dv_ = du__ dw ^ = dv__ ft (6)
dy dz dz dx dx dy
u2 + v2 + w2 (7)
E = Po-2-■ (7)
At z = 0 w = -ctt, at z = H(x, y) w = 0. (8)
At the side walls the impermeability conditions are as follows: for meridional u = 0, for zonal boundary sections v = 0. (9)
For temperature and salinity, the conditions of the absence of flows are set. Initial conditions:
at t = to (T, S) = (T0, S0), u = u0, v = v0, w = w0. (10)
Applying the rot operation to equation (1), for £ we obtain
fUvxG xü) = Vx ( gp ). (11)
dt
Equation (11) explicitly describes the transformation of an important hydrodynamic characteristic - an absolute vorticity. Correct approximation of the properties of this equation is of fundamental importance. For example, if the right side of this equation (barotropic motion) is equal to zero, the Helmholtz theorem conditions are satisfied and vortex filaments and vortex tubes 3 are preserved in the liquid. To fulfill these properties, the finite-difference approximation of the system of equations (1)-(6) in case of constant density must result in a discrete equation of the absolute vorticity with a zero right-hand side. Otherwise, computational features arise in the difference problem. Due to absolute vorticity (6) definition, we have
Vl = 0. (12)
Then equation (11) can be rewritten as follows:
+ Ü(V) - 1 (VÜ ) = Vx( gp). (13)
dt
If two-dimensional motion in (x, y), (x, z), (y, z) planes is considered, then for
each projection the following conservation laws can be specified. In (x, z) plane, for example, equations (2) and (13) are rewritten as
du dw d£y d£y d£y
— + — = 0, —+ u—+ v—^ = 0. dx dz dt dx dz
Then the conservation laws (yy is stream function):
f J£Q=| Jy)2 ,Q=f J(Vyy)2 „Q=0
and the antisymmetry property are satisfied: J (yy, £y) = - J (£y, yy).
If two quadratic conservation laws are satisfied, it can be demonstrated 4 [12] that the average wavenumber does not depend on time. This property means that
3 Kochin, N.E., Kibel, I.A. and Roze, N.V., 1964. Theoretical Hydromechanics. New York: Interscience Publishers, 577 p.
4 Dymnikov, V.P., 1984. [Computational Methods in Geophysical Fluid Dynamics]. Moscow: Department of Computational Mathematics of the USSR Academy of Sciences, 148 p. (in Russian). 152 PHYSICAL OCEANOGRAPHY VOL. 31 ISS. 2 2024
a systematic energy transfer to high-wavenumber motions cannot exist. If this property is not met in the numerical model, a false (computational) energy cascade to minimally resolved (two-step) waves on the grid may occur. The group velocity of two-step waves is zero, so energy accumulates on these scales.
The vorticity equation in the form (13) allows to obtain the potential vorticity equation for a stratified incompressible fluid.
The equations (3)-(5) imply the density advection equation
dp + div(pU) = 0. (14)
Due to equality (12), we have
V[U(Vi)] - V[i (VU)] -V [Vx (g p )] = 0. (15)
From the density advection equation (14) we obtain three equations for density gradients. To do this, we differentiate equation (14) sequentially with respect to x, y, z. We multiply the resulting equations by ix ,iy ,iz, equation (13) - by the corresponding derivatives of p. After simple transformations, considering relations (14) and (15), we obtain the equation for PV [2]:
— + (UV) © = 0, © = | V p, (16)
dt
where © is potential vorticity of a stratified incompressible fluid. Equations (2) and (16) imply an important integral conservation law for a three-dimensional motion of a baroclinic fluid
-f © d Q = 0. (17)
dt Q
Equality (17) means that potential vorticity a is conserved in each particle as it moves. Let us note that equation (16) characterizes the main feature of a stratified incompressible fluid under the following conditions: potential mass force (gravity force) and incompressibility of seawater, absence of viscosity and density diffusion in motion equations [4]. The result obtained by G. Ertel has two remarkable features. If we consider any Lagrangian invariant (for example, -9) instead of density, then
a = \ V9 will also be an invariant [4]. The second feature is that for an incompressible fluid, due to the form of equation (16), an invariant is any ¥ operator
dT
differentiable with respect to a. Let us multiply (16) by-= ¥'. Then we get
d to
dT (to)_n d
= 0, — J ¥(©) dQ = 0, J ¥(©) dQ = const.
dt dt Q
Discrete equation of absolute vorticity and conservation laws
We approximate a basin with an uneven bottom by boxes; their centers correspond to integer values of indices i, j, k (i = i1,...,iN, j = j1,..., jM , k = 1, ..., Ki, j), edges - to i+1/2, j+1/2, k+1/2. Horizontal dimensions of the boxes
(hx, hy) are constant; uneven approximation (hk2 = zk+1/2 - zk-1/2, hk2+112 = zk+1 - zk)
is applied vertically [13, 14].
Difference operators have the following form (for j, k - similarly):
PHYSICAL OCEANOGRAPHY VOL. 31 ISS. 2 2024 153
♦j = *'+1/2, j, k + f-1/2, j, k, Sx^,,j,k = ♦+1/2, j, f -1/2, j •k, V2A;,k =§x^-,j,k +6 y , j, k,
1 1 K K, (18)
{♦}Qk =7rS+i,j,khxh,, (♦)"= 1 Z jk, F = hkhxhy.
Q k i, j F i, j k=1 i, j k=1
Temperature, salinity, and horizontal velocity components are calculated at zk horizons, vertical velocity is calculated at zk+1/2 horizons. The distribution of variables is given in Fig. 1.
^/-1/2,7+1/2, k-m ^ /,7+1/2^-1/2 ---*-
Ш
/+1 /2,/+1/2, Ar—1/2
i,j,k-m/
7,7-1/2,iA-i/:
-.-...ß^AwU-
ifJ-1/2, A-
'Щ--1/2,7+1/2,А+Ш
Щ+Т/У"
1 /2 J-1 /2, A+1 /2 Jp'/ J-1 /2, A+1 /2
¿+1/2, k+l/2
Э /+1/2J+1/2, £
^/+1/2,7+1/2^+1/2
fm
f+1/2,7-1/2,A+1/2
К
Fig. 1. Distribution of variables in box (z, j, k). PV (co) is determined at the box vertices indicated by an asterisk, and the components of absolute vorticity iX j+1/2 k 'ii^^ j k+1/2 ,iZ+U2 j+1/2 k are determined at its edges
According to notations (6), (7) and (18), we write out the finite-difference equations of model (1)-(5) (differential in time) [9, 14]:
dU' -[v,iz+[w,iy=-Sx(Ei+H2,.j/k + w), C19.1)
+1/2, j, k dt
dv,
i, j+1/2, k dt
+
["' ^ li, j+1/2,k ]i„+1/2,k = Л ( ^i, j+1/2, k + W (19-2)
dw,
i, j, k+1/2
dt
-[ u, ]+/2 +[ V,]+,/2 = "6 ^ ,j, k+1/2 + ^ ,j, k+1/2 ) + «Pi ,j, k+1/2. (^3)
(20)
(21)
5 A,,,, k + § ]vi, j, k + 5zW,,,, k = 0
dT
i, j, k dt
+ 6 * ( U , T , k ) + 6 ] ( v, j, T ,j, k ) + 5 z ( Wi j, kT, j, k ) = 0,
dS,
, j, k
dt
+ § * ( U ,j, A, j, k ) + § ] ( v„ j, kSK j, k ) + 5 z ( w, j, kSi, j, k ) = 0, Pi, j, k = aTt, j, k + ßSi, j, k.
(22) (23)
In accordance with difference operators (18), the vorticity components (Fig. 1) have the following form:
£] j+1/2,k+1/2 = 5] (Wi , j+1/2, k+1/2 ) " 5z (Vi,j+1/2,k+1/2 \
4-z1/2, j,k+1/2 5z (ui +1/2, j,k +1/2 ) 5x (Wi +1/2, j,k +1/2 X 41 +1/2, j+1/2,* — 5x (Vi+1/2, j+1/2,* X — 5y (Mi+1/2, j+1/2,* X + f ■z 1/2 ■
From approximation (24) it follows that at points i + 1/2, j + 1/2, k + 1/2 (analogous to (12)) the following is fulfilled:
5 x 4x + 5 y 4y + 5 z 4z = 0.
(25)
We assume that the terms in square brackets on the left side of equations (19.1)-(19.3) are written down in the form
1
[v, 4z] i+1/2, j ,k ~ o ] z 1/2, j ,k S;z 1/2, j ,k
=3
4z
Si z
Vi+1, j, k 4iz 1/2, j, k + Vi j, k ^siz1/2, j, k
[^ 4yz 1/2,j,k - 1/2,j,k 4i+1/2,j,k
Wi+1, j, k 1/2,j ,k + Wi, j, k 4i-1/2, j, k
[u, 4 z]i
j+1/2,k
2ui,j+1/2,k 4i,j+1/2,k
x u + -
,j+1,k 4i, j+1/2,k + Ui ,j ,k 4 i,j-1/2,k
(26)
[W 4x I, j +1/2,k
2w;
i, j+1/2,k
w,
k 4x
i, j+1,kSi, j+1/2, k
+ Wi, j,k 4
i, j ,k*i, j-1/2,k
j+1/2,k
[u, 4y ]
2u 4 y
i,j ,k+1.2 Si,j ,k+1/2
Ui, j ,k+14i, j ,k +1/2
Ui, j ,k 4yj ,k-1/2
[v. 4x ]
i, j ,k +1/2
2v 4 x
, j ,k +1/2 si, j ,k +1/2
v
i, j, k +1s i+1/2, j, k
+ v.
,kSi-1/2, j ,k
Finite-difference representation (26) ensures, in the case of two-dimensional divergence-free motion in (x, y), (x, z), (y, z) planes, the fulfillment of two quadratic conservation laws: kinetic energy and enstrophy (vorticity square), as well as antisymmetry property [9, 14].
Taking into account expression (24) and relations (26), we write down the equations for the components of the absolute vorticity - for 4x at i, j +1/ 2, k +1/2 point, for 4y at i +1/2, j +1/2, k point, and for 4z at i +1/2, j +1/2, k point (analogue of equation (13)):
dr + 5 y ([v,4x ]) + 5 z ([w,4x ]) - 5 y ([u,4y ]) - 5 z ([u,4z ]) - g5y pz, (27.1)
-4 y -z
-4- + 5x([u,4y]) + 5z([w,4y])-5x([v,4x])-5z([v,4z])--g5xp , (27.2) dt
+ 5 x ([u ,4z ]) + dy ([v,xz ]) - 5 x ([w,4x ]) - 5 y ([w,4y ]) - 0. (27.3)
We are to analyze the system of equations (27.1)-(27.3) in terms of fulfillment of conservation laws.
Let us consider, for example, two-dimensional divergence-free motion in (x, z) plane. In this case, continuity equation (20) is transformed to the form
S x,* + S z",,j, * = 0 (28)
and provides the introduction of the flow function: ui+1/2 * = Szvf+1/2 *, w. *+1/2 = Sxvf*+1/2.
From the continuity equation (28), it follows that
Y y __Y . y
S,'+1/2, k+1/2 = °z Vi+1/2, к+1/2 + 5x Vi+1/2, k+1/2 ■
The vorticity equation is simplified and transformed [9] to the form
dК
'1+1/2, *+1/2 1 V j /
-+ _Z JS (Vi
dt 3 5=1
У
К У 2 > S-
where J = 5X «1/2, *+1/2 )5y (Kiy 1/2, *+1/2 ) 5y (Vy 1/2 , к+1/2 )5x (Kiy 1/2, *+1/2 X
■+1/2, * +1/2' Ы+1/2, *+1/2 - У
) = 0,
J 2 = 5x
V/+1/2, *+1/2 5y (КУ+1/2, *+1/2 )
- 5,
Vf+1/2 , *+1/2 Sx (К
i+1/2, *+1/2
(29)
(30)
j 3 = 5 y
x f-
5x (УУ
1+1/2, i+1/2 ""xVYi+1/2, к+1/2
2 )
- 5.
5 „ (УУ
i+1/2, к+1/2 "y VTi +1/2, i+1/2
)
Expression (30) means that equation (29) has two quadratic invariants:
dt Z Ky+1/2,*+1/2hxhZ=o, dt Z
i ,* ш i,*
(5 x V i+1/2 ,*+1/2 )2 + (5z Vi+1/2,*+1/2 )
2
hX = 0,
(31)
j Z (Ki+1/2,*+1/2 dt i ,*
' )2 kK = 0.
Approximation of nonlinear terms in equation (29) provides the antisymmetry property:
J(Vf+1/2, k+1/2 '£f+1/2, k+1/2 ) _ — J(£f+1/2, k+1/2, ¥'+1/2 , k+1/2 ). (32)
A similar situation occurs for the motion in (y, z) and (x,y) planes:
—Z Iх hhk = — Z
dtZ Кj+1/2,*+1/2Kyhz dtZ
j,* ш j,, d x
7 +1/2,*+1/2 ) hyhz
(5y V j+1/2,*+1/2 )2 + (5zVj+1/2,*+1/2 )
2
hyK =
(33)
= -rZ (Kj+1/2,*+1/2 )2hyhz = 0,
at j,*
J(V;+1/2 , *+1/2 x+1/2, *+1/2 ) = - J(Кjj+1/2, *+1/2' V7+1/2, *+1/2 X
К z hh = f(5xVz+1/2J+1/2 )2 + (5yVz+1/2J+1/2 )2
dt z К+1/2,j+1/2hxhy dt Z
dt Z
hxhy =
(34)
i +1/2, j+1/2 )2hxhy = 0
1/2, j+1/2 ' Si+1/2, j+1/2 ) = -J^l 1/2, j+1/2 +1/2, j+1/2 ).
Along with the feature of total energy conservation [5], relations (31)-(34) ensure the presence of two quadratic discrete invariants in the system of equations (1)-(5) with boundary (8)-(9) and initial (10) conditions.
Discrete equation of potential vorticity of a stratified incompressible fluid
We consider the properties of system (27.1)-(27.3) in the case of three-dimensional motion. We write this system in the following form:
-x
4.x
—x
dKx+1/2, j+1/2,* +1/2 +Tx x =5 (-xz ) (35 1)
^t 1 i+1/2, j+1/2,*+1/2 =g5y (pi+1/2, j+1/2,*+1/2 ),
d Ky+1
, ^ 1 i+1/2, j+1/2,*+1/2 -g5x (P i+1/2, j+1/2,*+1/2 ), (35.2)
-z
'cz
d 4i+1/2, j+1/2,*+1/2 , Yz (35 3)
dt i+1/2, j+1/2,*+1/2 V-^'-V
where
(36)
^x, j+1/2,k+1/2 - 5 y ([V,4x ])i, j+1/2,k + 5 z ([W,4x ])i, j+1/2,k -
-5y ([u,4y])i,j+1/2,k - 5z([u,4z]),,.+1/2,k,
^+1/2,j,k+1/2 - 5x ([u,4y ])i +1/2, j,k +1/2 + 5z ([W,4y ])i+1/2,j,k+1/2 --5x ([V 4x ])i+1/2, j,k +1/2 - 5z ([V 4z ])i+1/2,j,k+1/2 ,
^+1/2, j+1/2,k - 5x ([u'4z ])i+1/2, j+1/2,k + 5y ([V'4z ]).-+1/2,j+1/2,k --5x ([W 4x ])i+1/2, j+1/2,k - 5y ([W 4y ])i+1/2,j+1/2,k ■
It is easy to verify that according to relation (25), equations (27.1)-(27.3) at point i + 1/2, j + 1/2, k + and definition (36) we obtain
Sx (x 1/2, j+1/2,*+1/2 ) + ^y ((-+1/2, j+1/2,*+1/2 ) + ^z (i+1/2, j+1/2,*+1/2 ) = (37)
Note that equality (37) is satisfied because density linearly depends on temperature and salinity. We assume that from equations (21)-(23) follows a discrete equation for density, which has the following form:
-PdtjL + 5 x («,, J, k p ,, J, k ) + 5 y (v, j, k p i, j, k ) + 5 z ( w, , j, k p i, j, k ) = 0. (3 8)
We write equation (38) at the points (i, j +1/2, k +1/ 2), (i +1/2, j, k +1/ 2), (i +1/2, J +1/2, k), respectively:
-y
dpi,j+1/2,*+1/2 + 5 ( yz ) + 5 ( yz ) +
dt x (Ui,j+1/2,*+1/2Pi,j+1/2,* +1/2 ) + 5y (Vi,j+1/2,*+1/2Pi,j+1/2,*+1/2 ) +(39 1)
-yz
+5z (Wi,j+1/2,*+1/2Pi, j +1/2,*+1/2 ) = 0
—xz
d p,
+1/2, j ,*+1/2
dt, + 5x (Ui+1/2,j,*+1/2Pi+1/2,j,*+1/2 ) + 5y (Vi+1/2,j,*+1/2Pi+1/2,j,*+1/2 ) +(39 2)
+5z (W+1/2,j,*+1/2Pi +1/2,j,*+1/2 ) = 0
-xy
dPi+1/2, j +1/2,* + ~ (-^ . + ~ (-^ . +
d^ x \Ui+1/2, j+1/2,*P i+1/2, j+1/2,* / + О y I -+1/2, j+1/2,*P i+1/2, j+1/2,* / + (39 3)
—xy
+ 5z (w,+1/2, j+1/2,*P,+1/2,j+1/2,* ) =
y
We differentiate equation (39.1) in the finite-difference sense with respect to x, (39.2) with respect to y and (39.3) - to z. As a result, we get
d5 x (pi+1/2, j+1/2,k +1/2) s / Dx \ n /ja n -^-+ 5x №+1/2, j+1/2,k +1/2 ) - 0, (40.1)
dSy (pi+1/2, j+1/2,k+1/2 ) -xy
dt
dSz (pi+1/2, j +1/2,k+1/2 )
+ 5y (ry-1/2,j+1/2,k+1/2 ) = 0, (40-2)
+ S, (R+v2j+v2mv2 ) = 0, (403)
dt
where the notations R*j+mMm,R^j,k+1/2,R^J+1/2,k are obvious.
Let us mention the properties of the introduced functions in equations (39.1)-(39.3) and (40.1)-(40.3). They have the following form:
~xz ~x+ -x+z
pi,j+1/2,k+1/2 - p/+1/2,j,k+1/2 - pi'+1/2,j+1/2,k - pi+1/2,j+1/2,k+1/2 , (41)
-x -y -z -xyz
Ri,j+1/2,k+1/2 - Ri+1/2,j,k+1/2 - R,i+1/2,j+1/2,k - 5x (ui+1/2, j+1/2,k +1/2pi+1/2,j+1/2,k+1/2 ) +
-xyz
+5y (Vi+1/2,j+1/2,k+1/2pi+1/2,j+1/2,k+1/2 ) + 5z (Wi+1/2,j+1/2,k+1/2pi+1/2,j+1/2,k+1/2 )
We assume that the difference analogue of PV is defined at point i +1/ 2, j +1/ 2, k +1/2 and is written down as follows: rai+1/2 J+1/2 k+1/2 -
4ix1/2, j+1/2,k+1/2 S x (pi+1/2, j+1/2,k +1/2 ) + 4i+1/2, j+1/2,k +1/2 S y (Pi +1/2, j+1/2,k+1/2 ) + (42) z —xy
i+1/2, j+1/2,k +1/2 S z (Pi +1/2, j+1/2,k+1/2 ).
Determination of potential vorticity of a three-dimensional stratified fluid (42) at the box vertices is determined by the fulfillment of discrete equation (25).
To obtain a difference analogue of Ertel's theorem (15), we multiply equations
(35.1)-(35.3) ЬУ Sx (P/+1/2, j+1/2,k+1/2 X Sy (PX+1/2,j+1/2,k+1/2 X Sz (РхУ1/2, j+1/2,k+1/2 ), and
y
the system (40.1)-(40.3) - by Kx+1/2,j+1/2,*+1/2 ,Ky+1/2,j+1/2,*+1// ,Ki+1/2,j+1/2,*+1/2 , respectively. As a result, at point i + 1/2, j + 1/2, * + 1/2 we get
dfrt -x —yz -y —xz -z —xv -x -y -z
— + Yx 5x(Py) +Yy 5y(P ) + Yz 5z(Py) + Kx 5xRx +Ky 5 Ry +Kz 5+R+ =0. dt
Since Kx ,Ky ,Kz satisfy relation (25), Yx, Yy, Yz - relation (37), with regard to equalities (41), at point i + 1/2, j + 1/2, * + 1/2 we obtain finite-difference equation of potential vorticity (differential in time) of a stratified incompressible fluid in a divergent form:
—+ Sx(Yx pyz +£,xRx ) + Sy(Yy pxz + £,yRy ) + 8г(1г pxy ) = 0. (43) dt
Note that formally there is no difficulty in obtaining equation (43) with temporal discretization. Due to additional indexing, difference equations become cumbersome and therefore difficult to read.
Strictly speaking, the form of nonlinear terms in equation (43) does not
correspond to their differential analogue (15). The integral over domain from ^ is
dt
equal to zero when the following relations are satisfied at the boundaries:
YX,j+1/2,k+1/2 = 0, Ri,j+1/2,k+1/2 = 0 at i = i1 , i = iN ,
Yy 1/2,j,k+1/2 = 0, R^jk+1/2 = 0 at j = j1, j = jM, (44)
*il!2,j+j = 0, Rll!2,j+j = 0 at k =1/2, k = Ki j +1/2, From relations (44) it follows that additional boundary conditions, which are absent in the original formulation, are required.
Conclusion
A finite-difference analogue of absolute vorticity is written out for the system of equations of an ideal fluid without a quasi-static approximation. Projections of this equation onto two-dimensional subspaces (x, y), (y, z), (x, z) preserve energy, vorticity, enstrophy and have the antisymmetry property. To obtain the well-known Arakawa-Lamb scheme for the shallow water equations, it is necessary to write out an original difference system of equations for horizontal velocities that differs from (26).
The original result is the obtained discrete equation for potential vorticity of a stratified incompressible fluid as an exact consequence of the original finite-difference system of equations which have a divergent form. In this case, density satisfied the linear equation of state, the approximation of which in this case ensures both total energy conservation and the divergent form of equation for PV. If a nonlinear density dependence on temperature and salinity is applied, a special density approximation at the box edges is required to preserve total energy. In this case, an additional term occurs in the discrete potential vorticity equation; it has no analogue in the differential problem. Another problem arises from the form of advective terms in the PV equation, which are fundamentally different from their differential counterparts. In order for the potential vorticity to be conserved, we require additional boundary conditions under which the PV discrete analogue is an invariant.
REFERENCES
1. Rossby, C.-G., 1940. Planetary Flow Patterns in the Atmosphere. Quarterly Journal of the Royal Meteorological Society, 66(S1), pp. 68-87. https://doi.org/10.1002/j.1477-870X.1940.tb00130.x
2. Ertel, H., 1942. Ein Neuer Hydrodynamischer Wirbelsatz. Meteorologische Zeitschrift, 59(9), pp. 277-281. https://doi.org/10.1127/0941-2948/2004/0013-0451 (in German).
3. Hoskins, B.J., Mclntyre, M.E. and Robertson, A.W., 1985. On the Use and Significance of Isentropic Potential Vorticity Maps. Quarterly Journal of the Royal Meteorological Society, 111(470), pp. 877-946. https://doi.org/10.1002/qj.49711147002
4. Zhmur, V.V., Novoselova, E.V. and Belonenko, T.V., 2021. Potential Vorticity in the Ocean: Ertel and Rossby Approaches with Estimates for the Lofoten Vortex. Izvestiya, Atmospheric and Oceanic Physics, 57(6), pp. 632-641. https://doi.org/10.1134/S0001433821050157
5. Kapsov, E.I., 2019. Numerical Implementation of an Invariant Scheme for One-Dimensional Shallow Water Equations in Lagrangian Coordinates. Keldysh Institute Preprints, (108), pp. 1-28. https://doi.org/10.20948/prepr-2019-108 (in Russian).
6. Bihlo, A. and Popovych, R.O., 2012. Invariant Discretization Schemes for the Shallow-Water Equations. SIAM Journal on Scientific Computing, 34(6), pp. B810-B839. https://doi.org/10.1137/120861187
7. Charnyi, S., Heister, T., Olshanskii, M.A. and Rebholz, L.G., 2019. Efficient Discretizations for the EMAC Formulation of the Incompressible Navier-Stokes Equations. Applied Numerical Mathematics, 141, pp. 220-233. https://doi.org/10.1016/j.apnum.2018.11.013
8. Sorgentone, C., La Cognata, C. and Nordstrom, J., 2015. A New High Order Energy and Enstrophy Conserving Arakawa-Like Jacobian Differential Operator. Journal of Computational Physics, 301, pp. 167-177. https://doi.org/10.1016/jjcp.2015.08.028
9. Arakawa, A. and Lamb, V.R., 1981. A Potential Enstrophy and Energy Conserving Scheme for the Shallow Water Equations. Monthly Weather Review, 109(1), pp. 18-36. https://d0i.0rg/10.1175/1520-0493(1981)109<0018:APEAEC>2.0.C0;2
10. Salmon, R., 2005. A General Method for Conserving Quantities Related to Potential Vorticity in Numerical Models. Nonlinearity, 18(5), pp. R1-R16. https://doi.org/10.1088/0951-7715/18/5/R01
11. Sugibuchi, Y., Matsuo, T. and Sato, S., 2018. Constructing Invariant-Preserving Numerical Schemes Based on Poisson and Nambu Brackets. JSIAM Letters, 10, pp. 53-56. https://doi.org/10.14495/jsiaml.10.53
12. Arakawa, A., 1966. Computational Design for Long-Term Numerical Integration of the Equations of Fluid Motion: Two-Dimensional Incompressible Flow. Part I. Journal of Computational Physics, 1(1), pp. 119-143. https://doi.org/10.1016/0021-9991(66)90015-5
13. Demyshev, S.G., 2004. Energy of the Black Sea Climatic Circulation. Part I: Discrete Equations of the Rate of Change of Kinetic and Potential Energy. Meteorologiya i Gidrologiya, (9), pp. 65-80 (in Russian).
14. Demyshev, S.G., 2005. Numerical Experiments Aimed at the Comparison of Two Finite-Difference Schemes for the Equations of Motion in a Discrete Model of Hydrodynamics of the Black Sea. Physical Oceanography, 15(5), pp. 299-310. https://doi.org/10.1007/s11110-006-0004-2
15. Demyshev, S.G., 2023. Nonlinear Invariants of a Discrete System of the Sea Dynamics Equations in a Quasi-Static Approximation. Physical Oceanography, 30(5), pp. 523-548.
Submitted 09.06.2023; approved after review 25.07.2023;
accepted for publication 18.01.2024.
About the author:
Sergey G. Demyshev, Head of Wave Theory Department, Chief Research Associate, Marine Hydrophysical Institute of RAS (2 Kapitanskaya Str., Sevastopol, 299011, Russian Federation), DSc (Phys.-Math), Scopus Author ID: 6603919865, ResearcherID: C-1729-2016, ORCID ID: 00000002-5405-2282, demyshev@gmail.com
The author has read and approved the final manuscript. The author declares that he has no conflict of interest.