Calculation of the Total Flow Components in the Models of Wind Fluid Motion Текст научной статьи по специальности «Математика»

Original Russian Text © V. S. Kochergm, S. V. Kochergm, S. N. Sklyar, 2023, published in MORSKOY GIDROFIZICHESKIY ZHURNAL, Vol. 39, Iss. 3 (2023)

Original article

Calculation of the Total Flow Components in the Models of Wind Fluid Motion

V. S. Kochergin 1 S. V. Kochergin S. N. Sklyar 2

1 Marine Hydrophysical Institute of RAS, Sevastopol, Russian Federation 2 American University in Central Asia, Bishkek, Kyrgyz Republic H vskocher@gmail.com

Abstract

Purpose. The object of the work is to construct an effective numerical method to solve the problem for a stream function and to determine subsequently the total flow components in the models of fluid wind motion in a reservoir. Its efficiency was analyzed, and the work of the proposed difference approximations was illustrated using a class of test problems with the known analytical solutions. Methods and Results. The difference scheme and the corresponding computational algorithm were constructed based on the projection variant of the integro-interpolation method, which permitted (within a single approach) to solve the problem for the stream function, to calculate its derivatives, and to determine subsequently the total flow horizontal components.

Conclusions. The used discretization method permits to preserve automatically the most important features of the initial differential model at switching to its discrete analog. In particular, its application makes it possible to reproduce correctly the behavior of the stream function derivatives, and hence, the behavior of the total flow horizontal components in the areas of its highest intensity.

Keywords: wind currents, stream function, total flow component, singularly perturbed problem, analytical solution, difference scheme

Acknowledgements: The study was carried out within the framework of the state task on theme FNNN-2021-0005 "Complex interdisciplinary studies of oceanological processes which determine functioning and evolution of the ecosystems of coastal zones of the Black and Azov seas".

For citation: Kochergin, V.S., Kochergin, S.V. and Sklyar, S.N., 2023. Calculation of the Total Flow Components in the Models of Wind Fluid Motion. Physical Oceanography, 30(3), pp. 274-287. doi:10.29039/1573-160X-2023-3-274-287

DOI: 10.29039/1573-160X-2023-3-274-287

© V. S. Kochergin, S. V. Kochergin, S. N. Sklyar, 2023

© Physical Oceanography, 2023

Introduction

When solving many problems of hydrothermodynamics, special integral functions are used. In particular, when modeling the fluid dynamics in a reservoir 1, the integral stream function or the level surface is usually used. The corresponding problems belong to the class of singularly perturbed ones [1] and can form areas with large solution gradients, the so-called boundary 2 or internal transition layers. It is known that in these problems there are increased requirements for difference

1 Marchuk, G.I. and Sarkisyan, A.S., 1988. Mathematical Modeling of the Ocean Circulation. Moscow: Nauka, 302 p. (in Russian).

2 Vasilyeva, A.B. and Butuzov, V.F., 1973. [Asymptotic Expansions of Solutions of Singularly Perturbed Equations]. Moscow: Nauka, 272 p. (in Russian).

ISSN 1573-160X PHYSICAL OCEANOGRAPHY VOL. 30 ISS. 3 (2023)

The content is available under Creative Commons Attribution-NonCommercial 4.0 International (CC BY-NC 4.0) License

schemes that implement them numerically [2]. Additional difficulties appear when it is required to calculate the derivatives of the solution of a singularly perturbed boundary value problem, namely, it needs to be carried out when determining the horizontal components of the total flow.

All of the above explains the relevance of the problem being solved. These problems can be solved if special discretization methods are applied, which allow auto-saving of the most important properties during the transition from the original differential model to its discrete counterpart. The present study is devoted to the implementation of one of these methods. Comparison of the numerical solution with the available exact analytical analogue makes it possible to judge the accuracy of the algorithm used. The discretization method used in the paper permitted to correctly reproduce the behavior of the stream function derivatives, and hence the horizontal components of the total flow.

A computational algorithm for a simplified three-dimensional stationary model of fluid wind currents in a reservoir will be constructed. The fluid is assumed to be homogeneous, and the model lacks the advection and horizontal diffusion mechanisms. Such models belong to the Ekman type models [3] and are used in the first approximation to describe the pattern of streams in reservoirs of various nature. In addition, if it is possible to find classes of analytical solutions in such models, then it is convenient to use them for testing numerical methods and corresponding algorithms used in solving general problems of reservoir hydrodynamics.

It is assumed that the problem is considered in a dimensionless setting. Suppose that the reservoir surface in the xy plane is shaped like a rectangle:

its depth H > 0 is constant, the x axis is directed to the east, the y axis - to the north and the z axis - vertically down. In a 3D region

Problem statement

Q0 =[0,r] x [0,q],

(1)

(2)

the following system of equations is considered:

(3)

du dv dw

= 0,

+

+

dx dy dz

supplementing it with boundary conditions

o j du dv

z = o, (x, y) e Qo}: k = _T- ' k dZ = ~Xy , w = °

z = H, (x,y)e Qo}: k^ = , k=, w = 0, { 0 < z <H, (x,y) eôQo} : Unx +Vny = 0.

In (4) there are horizontal components of the total flow

H H

U (x, y) = j u (x, y, z )dz, V (x, y) = j v( x, y, z )dz

(5)

and the following variant of the parametrization of near-bottom friction is accepted:

tx = ^U, xy = ^V, ^ = const > 0.

In accordance with Stommel works [4, 5], it is assumed that l = l0 + Py, k = const > 0.

The tangential wind stress will be given by the formulas |xx = [F cos(rix) + f2 sin(rix)]cos(qmy\ S[xy = [Gj cos(r^x) + G2 sin(rix)]sin(qiy),

in which the following designations are accepted:

(6)

(7)

(8)

ns _ nm _ np

, qm = , qp = ,

nl

r r

l, s = o,1, 2,...; m, p = 1,2,...

(9)

Thus, the wind model contains four real (F1, F2, G1, G2) and four integer ( l, m, s, p ) numerical parameters, the choice of which allows to describe a fairly general wind situation. Moreover, to describe the components of tangential wind stress, a linear combination of expressions (8) can be used, and the problem solution can be found as a linear combination of the corresponding "elementary" solutions. Note that Stommel used a wind model of the form as follows

t„ = -F cos

( ny >

T y = ^

(1o)

which is obtained from expressions (8), (9) when

F = -F, F2 = G1 = G2 = o, l = o, m = 1.

In [6, 7], analytical solutions to problem (1) - (7), (1o) were found. In [8], analytical expressions for flows were obtained from relation (5) in problem (1) -(9), taking into account general wind situation. These analytical solutions are used in the present study to test the proposed computational algorithms. Let us transform

the problem (1) - (9): each equation is integrated in the formula (3) with respect to the variable z in the range from 0 to H, taking into account the boundary conditions (4); then the pressure gradients from the obtained equations are eliminated using the operation of cross differentiation. As a result, we come to the following problem - definition of U (x,y), V (x, y) functions:

- PV = ^, dy dx

( dU dV dy dx

dU dV ( , °

— + —= 0, (^ y)e Qo, dx dy

Unx +Vny = 0, (x,y)e5Q0.

(11)

To solve it, the integral stream function ¥(x, y) is introduced according to the formulas

U =

dY

V = -

dY

dy dx

Then (12) is substituted into (11) and, denoting

' (x, y )-t-dx ■

the following problem to determine Y (x, y) is obtained:

^2

d2xY d2xY

2

dx2 cy2

dY o

"P — = f (x,y), (x,y) e Qo

Y = 0, (x, y) edQo.

(12)

(13)

(14)

Numerical method for solving problem (12) - (14)

To construct it, a projection version of the integro-interpolation method (PVIIM), which was proposed in [9] and studied in [10], is used. This technique, within the unified approach framework, permits to construct a difference scheme for the numerical solution of problem (14) and obtain formulas for approximating the derivatives of this solution. The latter is necessary for U and V calculation according to formulas (12). The works [11-13], where a similar problem was solved, also should be noted.

In Q0 region the computational grid is considered:

mh = {(x, yj) | x =(*-1) Ar, yj =(j -1)

i = 1, n, j = 1, k, Ax = —r—, Ay = ~~~ n -1 k -1

(15)

Let the grid function }, defined on this grid, consist of approximate

values for {* (xt, y.)} quantities - the exact solution to problem (14). The PVIIM

technique is applicable within the framework of successive approximation of the differential equation from expression (14): first with respect to the x variable, then with respect to the y variable. Equation (14) can be rewritten in the following form:

d2* _d¥ ( , d2* , , ^+=/(x,y)-s g (x,y). (16)

In accordance with the PVIIM, on an arbitrary grid cell, the test function [xt, xM ] is considered, equation (16) is multiplied by X (x) and the result is

integrated over the cell [xi, xi+1 ], including by parts (the y variable is perceived as a parameter), as a result, the following integro-difference identity is obtained:

(pX - ^ dX T + ? * f ^ dx? - p fV=i g (^ y)X (x) dx. (17)

^ ' J x x f J x

, £ X+* i PX - „ dX

dx

The test functions X(0)( x) u X (1)( x ) in identity (17) are selected, assuming that they are solutions of the equation

d2 X _5X ( )

but at the same time satisfy various boundary conditions:

X(0)(x,) = 1, X(0)(x,+1 ) = 0,

X«(x,) = 0, X(1) ( x,.+1 ) = 1.

This choice allows to turn into zero the integral on the left side of identity (17). It is obvious that for x e (x, x+1)

BAr/p. _ 0(B(-1

X(0)(x) =' _1 • X(1)(x)= (18)

The integral on the right side of identity (17) will be approximated using the following formulas:

xi+1 xi +1

j g (x, y) X(0)(x) dx xg (x., y ) j X(0)(x) dx,

x * (19)

xi+1 xi +1

j g (x, y) X(1)(x) dx xg (xt+1, y) j X(1)(x) dx.

x. x.

Note that more accurate approximations of the integrals in (17) can lead to a more accurate difference scheme, but we will focus on the variants above

(formula (19)). Let us substitute X = into identity (17) and, using the designations

r=^ «=fi(R)=-^

Y - Y Y - Y

D + Y = +1 D ~ Y = -1

x i . ' x i . '

Ax Ax

as well as necessary calculations, taking into account formulas (18) and (19), we rewrite it in the following form

dY 1+ fi _

-V+ V[1 + R(0 +1)]D+Yi = g(x,,y)Ax—, i = 1,n -1. (20)

After substituting the X = X(1) function into identity (17), we obtain

dYi+1 r, , . 1 -0

- v[1 + R(0 -1)]Dyt = g(x,+1,y)Ax—, i = 1,n -1. (21)

Let us move on to approximation with respect to the y variable. Considering that the g (x, y ) function on the right-hand sides of relations (20) and (21) was determined by formulas (16), these relations are re-written in the following form

V ^ = f ( Xi, y )+ 2 L ^ - v [1 + R (0 +1)1 D+ y, U dy Ax(1 + 0)[ cX [ v !] x i J (22)

- P(y), i = d2Y, , „, 2 f dY

■ i+1 .2

= f ( x+1, y ^dh) f' If+>[1+R ( 0 -1)1 Dx+ Y }-

" cy2 ^ Ax(1 - 0dx v '] x 'J (23) - q (y), i=im-1.

Taking in account the PVIIM for equation (22), it is re-written for arbitrary grid cell [yj, yj+1 ] :

d2Y

V = P (y), y 6 [yj, yj+1 ]. (24)

The integro-difference identity analogue of (17) for equation (24) will look as follows:

dY,. ^ dn ^^T n-Y

dy dy

Test functions n(0), n( ) will be considered as solutions to the problems

yj1 yj+1 d 2 n y>+1

+ j y,v-X^y = j P(y)n(y)dy. (25)

yj yj yj

u ^y? = 0, y e(y,, y,+1),

n(o) ( y j ) = 1, n(o) ( y j+1 ) = °, n(1)(y,) = °, n(1)(yj+1 ) = 1.

These functions are easy to find

n(o )( y ) =

yj+1 - y

Ay

, n(1)(y ) =

y - yj

Ay

(26)

The test functions (26) are substituted into identity (25) and approximations of integrals on the right side of expression (25), similar to approximations from formulas (19), are used:

yj+1

y.,+1

j P (y) n(o)(y) dy xP (yj ) j n(o)(y) dy,

y, yj+1

yj+1

j P(y)n(1)(y)dy xP(yj+1 ) j n(1)(y)dy,

as well as designations

- Y,

Y,,- Y,

D+y = i,j+1 i,j d - Y = - ',j ',j-1

Ay

as a result, finite-difference relations are obtained

dY.. , \ Ay

Ay

"^"âr+D9« =P(yjH, j=1k -1

u- yi,j=p(yj+1 )f, j=û-i

(27)

(28)

Using formula (22), from expressions (27) and (28) the following relations (i = 1,n -1, j = 1,k -1 ) are obtained:

2 Ay

dY

u - uD+Y . +Ay f (x ,y.)

^ dy y hj 2 y ' j'

dY

(29)

-u+uC1+R(e +1)]D+Y,j \,

2Ax [ dx dY

2Ay

u +DY, j +ff ( x, yj+1 )

dY

(3o)

i, j+1

2Ax

u [1 + R ( e +1)] D+ Y,,

dx v -/J

1

1

1

1

Let us carry out similar reasoning for relations (23). As a result, another group of relations (i = 1,n -1; j = 1,k -1 ) is obtained:

1-0

2Ay

= J. 2Ax

m ^yT - D Y'+1, j +ff ( Xi+1 ' ^ )

Ax f ^- m C1+R (0 -« d+y i,j },

(31)

1-

2Ay

<3Y

i+1, j+1

dy + +^f ( x+1, yj+1 )

1 I dY

1 JM —i^ - M[1 + R (0 -1)]D+ Y,,j+1

2Ax

dx

(32)

Relations (29) - (32), written for an arbitrary grid cell [ xi, xM ]x|yj, yj+1 J, i = 1, n -1, j = 1, k -1, provide full information necessary to

construct approximations of both the equation and all necessary derivatives in problem (12) - (14). Indeed, in formula (30) the index j is replaced by j - 1, this can be carried out at the internal nodes of the grid, then the result is added with

3Y..

equation (29), in the end, we eliminate the derivative

1 + (

dy

Ay (DY,j - D- Y,j )+f ( x, yj ) ^+m C1+r (0+1)] d: y „j ]

2

1 Ax

The last formula allows to obtain an approximation of the derivative

dY

_w

dx

which can be used at the internal nodes of the left vertical boundary: dY

= | 1+ R(0 + 1) \D+ Y....+Ax-

dx

= C1+R (0 +1)] Dx+ Y, j + Ax \jy ( D Y, j - D; Y, j ) -1 f ( x, yj )

M

(33)

Similarly, in (31) the index i is substituted by i - 1, and in in formula (32) i is substituted by i - 1 and j by j - 1. The results obtained are added. It provides

dY. .

the opportunity to eliminate the derivative

1 -(

.A. ( d+Y. - D'Y. ) + f ( x. , y. )

AyV y •,j y •,j/ v 1 j )

dy

1 I dY

=Ax f "âf- m ['+R ( 0 -o] DY

2

From the latter expression, we obtain an approximation of the derivative

dY.

—— for the internal nodes of the right vertical boundary:

dx

dY, r, ^ , . 1 -

& =[1+R(0+Ax^ -1 f (x,)—Ay^)

. (34)

Let us consider formulas (33) and (34) at the internal nodes of the region, multiply

1 — 0 1 + 0 the first by ——, the second by ——, add the results and obtain the relation

dY 1 — 0 1 + 0

= [1 + R ( 0 +1)] D+ Y „, + [1 + R (0 — 1)] D-Y,,,, (35)

dY.

which can be used to approximate the derivative —— at the internal nodes of

dx

the grid region. Now, if formula (34) is subtracted from expression (35), then

the derivative —— is eliminated and an approximation of equation (14) at

dx

the internal nodes of the grid area is obtained:

A- {[1 + R (0 +1)] D Y., j —[1 + R (0 —1)] D— Y.,,} +

Ax (36)

(d+Y . — D~Y. .) = f (x ,y.).

Let us move on to constructing an approximation for the derivative -—. In

cy

relation (30), the index j is replaced by j - 1, the result is subtracted from (29), and finally the following approximation is obtained

dY .

h j

1 (D+Y . + D~Y. ), (37)

") \ y •,j y hj J' v '

dy 2

which can be used in the internal nodes of the grid region. Next, in (31) i is replaced by i - 1, the result is added with formula (29), thereby excluding dY.

the derivative —— and obtain

dx

J_

Ay

3Y . Ay / x

u-ij — uD+Y. + f (x ,y.)

* dy 2 J y . j

A- {u [1 + R ( 0 +1)] D+ Y,, j — u [1 + R ( 0 — 1)] D— Yi,.}.

2Ax

From the last relation, the formula for the internal nodes of the lower horizontal boundary is found:

^ = D+Y i, j + {C1+R (0 +1)] D Y i, J-C1+R (0 -1)]D- Y i, J }-

dy 2Ax (38)

- Af ( -, yj ).

Finally, we consider formula (30), replace the index j with j - 1, and in relation (32) replace i with i - 1 and j with j - 1. After adding the results, we obtain

J_

Ay

-M ^Yf+md-y , j +ff ( x, y j+1 ) A-{mC1+R(0 +1)]DYi,j - MC1+R(0 -1)]D- YI,J }.

2Ax

The formula for the internal nodes of the upper horizontal boundary follows from this expression:

dr=DyY',j -2AX+R(0+1)JD+ Y.,j-[1+R(0-1)]D-Y,,j}+A/(x.yj).

(39)

The numerical solution to problem (12) - (14) is found by formulas (36) together with the corresponding boundary conditions. The problem is solved iterative ly by one of the known methods [14]. Further, formulas (33) - (35) are

dY,,

used to calculate the derivatives-—, and formulas (37) - (39) are used to

dx

calculate the derivatives. Formulas (33), (34), (38) and (39) for determining the derivatives at the boundaries of the region (when solving problem (12) - (14) in a rectangle) are significantly simplified due to the boundary condition for the Y function. It should also be noted that the difference scheme (36) was tested in solving problem (14) in [15] and turned out to be more accurate than the known schemes from papers 3,4 and [2, 16].

Results of numerical experiments

We are to illustrate bellow the operation of the proposed numerical methods with the results of the experiments with problem (1), (12) - (14), where the following values of the main parameters are chosen:

r = 11, q = 5, p, = 0,01, p = 1.

The model wind (formulas (8), (9)) is determined by the values

F =1, F2 =0, Gl = -1, G2 =0, l = 0, m = 1,5 = 1, p = 1,

at which a cyclone forms over the water area (Fig. 1).

3 Samarsky, A.A., 1983. The Theory o/Di//erence Schemes. Moscow: Nauka, 616 p. (in Russian).

4 Buleev, N.I., 1983. [Spatial Model o/ Turbulent Exchange]. Moscow: Nauka, 344 p. (in Russian). PHYSICAL OCEANOGRAPHY VOL. 30 NO. 3 2023 283

F i g. 1. Cyclonic wind field

The criterion for the computational algorithm quality was the relative error

lb—oil

Er (b) = ' 1°° (%),

which was calculated in the grid norm

lib = max lb .

Here b is the exact problem solution projected onto the grid (in this case, found using the formulas in [8]); b is an approximate solution to the same problem calculated using the corresponding algorithm. Let Y,^,Vbe the solution to problem (1), (12) - (14) found by the aforementioned method. Fig. 2 shows the analytical field of the stream function, which visually practically does not differ from that found numerically [15].

0 2 4 6 8 10 >'

F i g. 2. Stream function for the given parameters

Fig. 3 shows the analytical field of the total flows, obtained based on the formulas from [8].

I..........

1-;-;-1-1-i-1-i-1-r

0 2 4 6 8 10

F i g. 3. Analytical field of the total flows

The presented field is characterized by an intense boundary layer near the left boundary. The calculated field using expression (36) visually practically does not differ from that shown in Fig. 3, so it is not included in the study. In this case, it is advisable to judge the accuracy of calculations by the relative error values. When calculating total flows U, V using standard formulas

ou=<wu. u,,=+d-y„ j=u«=Dy Yi.,.i=&, (40)

V,J =-D+ Y,J, V,j =-2(D+ Yi,j + D;YU), i = Vn,j = -D-Yn,j, j = IJc (41)

the field shown in Fig. 4 is obtained. 9

4

3 2

1

0

F i g. 4. Field of the total flows calculated by formulas (40) and (41)

After solving problem (14) for the stream function using formula (36), U and V were determined in two ways: by the aforementioned method and using formulas (40), (41). The values of the corresponding errors calculated on grids with different numbers of nodes are given in the table. The values of the relative errors Er (Y) indicate the iterative process convergence due to the comparison of the obtained results with the exact solution. It is easy to check that the U values PHYSICAL OCEANOGRAPHY VOL. 30 NO. 3 2023 285

obtained using formulas (37) - (39) and the U values calculated in accordance with formula (40) are the same, so the Er (u) errors are not given in the table.

Error values at increase in the number of grid nodes

Number of grid nodes Er ( Y ) Er (U ) Er (V ) Er (V )

23x11 0.64259 1.06010 0.05480 98.0685

51x23 0.26472 0.22237 0.02779 95.5197

111x51 0.10534 0.07869 0.01187 90.0591

221x101 0.04036 0.03189 0.00471 80.2305

441x201 0.01247 0.01018 0.00148 63.8351

Conclusion

According to the results of numerical experiments, even with a sufficiently accurate solution to the stream function problem, calculation of the horizontal components of the total flow may turn out to be inefficient if the specifics of the problem are not taken into account. The algorithms used give a fairly good accuracy of reproducing the solution of the equation for the stream function. On a fine grid, the error is hundredths of a percent in the norm used. The technique proposed in this paper allows to solve the stream function problem within the framework of a unified approach and calculate the derivatives of this solution, which guarantees the accuracy of determining the horizontal components of the total flow. The results of the study can be used in simulation of dynamic processes in the sea.

REFERENCES

1. O'Malley Jr., R.E., 1974. Introduction to Singular Perturbations. New York: Academic Press, 206 p.

2. Roache, P.J., 1976. Computational Fluid Dynamics. Albuquerque: Hermosa Publishers, 434 p.

3. Ekman, V.W., 1905. On the Influence of the Earth's Rotation on Ocean-Currents. Arkiv for Matematik, Astronomi och Fysik, 2(11), pp. 1-53.

4. Stommel, H., 1948. The Westward Intensification of Wind-Driven Ocean Currents. EOS, Transactions, American Geophysical Union, 29(2), pp. 202-206. doi: 10.1029/TR029i002p00202

5. Stommel, H., 1965. The Gulf Stream. A Physical and Dynamical Description. Berkeley and Los Angeles, USA: University of California Press, 243 p.

6. Kochergin, V.S., Kochergin, S.V. and Sklyar, S.N., 2020. Analytical Test Problem of Wind Currents. In: T. Chaplina, ed., 2020. Processes in GeoMedia—Volume I. Cham: Springer Geology, pp. 17-25. doi:10.1007/978-3-030-38177-6_3

7. Kochergin, V.S., Kochergin, S.V. and Sklyar, S.N., 2021. Analytical Solution of the Test Three-Dimensional Problem of Wind Flows. In: T. Chaplina, ed., 2021. Processes in GeoMedia - Volume II. Cham: Springer Geology, pp. 65-71. doi:10.1007/978-3-030-53521-6 9

8. Kochergin, V.S., Kochergin, S.V. and Sklyar, S.N., 2022. Analytical Solution of the Equation for the Stream Function in the Model of Flows with Variable Wind Action in Space. Ecological Bulletin o/ Research Centers of the Black Sea Economic Cooperation, 19(1), pp. 16-24. doi:10.31429/vestnik-19-1-16-25 (in Russian).

9. Sklyar, S.N., 1995. A Projective Version of the Integral-Interpolation Method and Its Application for the Discretization of the Singular Perturbation Problems. In: A. S. Alekseev, N. A. Bakhvalov, eds., 1995. Advanced Mathematics, Computations and Applications. Proceedings o/ the International Con/erence AMCA-95. Novosibirsk: NCC Publisher, pp. 380-385.

10. Eremeev, V.N., Kochergin, V.P., Kochergin, S.V. and Sklyar, S.N., 2002. [Mathematical Modeling o/Hydrodynamics o/Deep-SeaBasins]. Sevastopol: ECOSI-Gidrofizika, 238 p. (in Russian).

11. Kochergin, V.P., Sklyar, S.N. and Sultanov, R.K., 1990. On the Problem of Numerical Ocean Hydrodynamics Modelling. Soviet Journal o/ Physical Oceanography, 2(2), pp. 89-98. https://doi.org/10.1007/BF02196859

12. Kochergin, V.P., Sklyar, S.N. and Sultanov, R.K., 1995. On the Calculation of the Barotropic Components of Motion in Models for General Circulation of the Ocean. Physical Oceanography, 6(1), pp. 19-24. https://doi.org/10.1007/BF02197293

13. Kochergin, V.P. and Dunets, T.V., 2001. Computational Algorithm for the Evaluation of Inclinations of the Level in the Problems of the Dynamics of Basins. Physical Oceanography, 11(3), pp. 221-232. https://doi.org/10.1007/BF02508869

14. Il'in, V.P., 2021. Iterative Preconditioned Methods in Krylov Spaces: Trends of the 21st Century. Computational Mathematics and Mathematical Physics, 61(11), pp. 1750-1775. doi:10.1134/S0965542521110099

15. Kochergin, V.S., Kochergin, S.V. and Sklyar, S.N., 2023. Testing of Numerical Methods for Solving the Stream Function Problem in a Model of Stationary Fluid Motion. In: T. Chaplina, ed., 2023. Processes in GeoMedia—Volume VI. Cham: Springer Geology, pp. 513521. doi: 10.1007/978-3-031-16575-7_46

16. Buleev, N.I. and Timukhin, G.I., 1969. Numerical Solution of Hydrodynamics Equations for the Planar Flow of a Viscous, Incompressible Liquid. Izvestiya SO AN SSSR, Series o/ Technical Sciences, 1(3), pp. 14-24 (in Russian).

About the authors:

Vladimir S. Kochergin, Junior Research Associate, Marine Hydrophysical Institute of RAS (2 Kapitanskaya Str., Sevastopol, 299011, Russian Federation), ORCID ID: 0000-0002-6767-1218, ResearcherlD: AAG-4209-2020, vskocher@gmail.com

Sergey V. Kochergin, Senior Research Associate, Marine Hydrophysical Institute of RAS (2 Kapitanskaya Str., Sevastopol, 299011, Russian Federation), Ph.D. (Phys.-Math.), ORCID ID: 0000-0002-3583-8351, ResearcherID: AAAG-4206-2020, svkocher@mail.ru

Sergey N. Sklyar, Head of Department, American University of Central Asia (AUCA) (7/6 Tokombaeva Str., Bishkek, 720060, Kyrgyz Republic), Dr.Sci. (Phys.-Math.), Professor, ORCID ID: 0000-0001-6985-6155, sklyar51@gmail.com

Contribution o/ the co-authors:

Vladimir S. Kochergin - numerical simulation, visualization, analysis of calculation results, text preparation

Sergey V. Kochergin - problem formulation, analysis of calculation results, text preparation

Sergey N. Sklyar - problem formulation, construction of a numerical algorithm, analysis of calculation results

The authors have read and approved the/inal manuscript.

The authors declare that they have no con/lict o/interest.