Software package for predicting possible scenarios for changing the geometry of the bottom of shallow water reservoirs using high-performance computing Текст научной статьи по специальности «Математика»

UDC 519.6 10.23947/2587-8999-2019-2-2-118-135

Software package for predicting possible scenarios for changing the geometry of the bottom of shallow water reservoirs using high-performance computing*

A. I. Sukhinov1, A. E. Chistyakov1, V. V. Sidoryakina2, E. A. Protsenko2, S. V. Protsenko1

1DonStateTechnicalUniversity,Rostov-on-Don,Russia

2Chekhov Taganrog Institute Taganrog branch of Rostov State University of Economics, Taganrog, Russia

The article is devoted to the study of the model of transport and sedimentation of suspended solids in the coastal zone. The model takes into account the following processes: advection transport due to the movement of the aqueous medium, microturbulent diffusion and gravitational sedimentation of particles of the suspension, as well as a change in the geometry of the bottom caused by the sedimentation of particles of the suspension or the rise of particles of bottom sediments. The article presents the results of a study of the correctness of the initial-boundary-value problem corresponding to the constructed model. Software package has been developed for predicting possible scenarios for changing the geometry of the bottom of reservoirs in shallow water using high-performance computing.

Keywords: distributed computing, high-performance computing, parallel programming, mathematical model, dynamics of sea sediments, bottom relief.

Introduction. The study of dynamic effects in the coastal zones of water bodies is very relevant in connection with the study of sedimentation and sedimentation processes, the generation and evolution of accumulative forms, the determination of the morpholithodynamic regime of the adjacent coastal zone [1]-[4]. The most widely accumulative forms are presented in the conditions of a gentle slope of the underwater and surface parts of the coast, i.e. on shallow shores, where conditions of shallow water and deformed waves dominate at a considerable distance from the coast [5]-[6]. Under the action of waves and currents, the material begins to move, undergoing sorting in accordance with the shape and mass of individual fractions, as well as mechanical and chemical changes. As a result, a flow of suspended solids is formed in the coastal zone, contributing to the formation of sediment.

Non-stationary spatially 3D models of suspension transport have been introduced into the practice of mathematical modeling relatively recently and have been verified by numerically solving a number of model and some real problems [7]-[10]. The authors (Sukhinov A.I., Sukhinov A.A., Chistyakov A.E., Protsenko E.A., Degtyareva E.E., Sidoryakina V.V.) previously proposed spatially 3D model of suspension transport, taking into account the following physical parameters and processes: advection transport due to the movement of the aqueous medium, microturbulent diffusion and gravitational sedimentation of particles of the suspension, as well as a

*The research is done with the financial support from Russian Foundation for Basic Research, Project No. 19-01-00701.

change in the geometry of the bottom caused by the sedimentation of particles of the suspension or the rise of particles of bottom sediments. In this work, an analytical study of the correctness of the initial-boundary-value problem corresponding to the constructed model is carried out. Based on the developed parallel algorithms implemented as a complex of programs, numerical experiments were performed for model problems of bottom sediment transport and bottom topography transformation, the results of which are consistent with real physical experiments.

Physical statement of the problem. Comprehensive studies of aquatic ecosystems are an integral part of environmental management. An important place in these studies is the study of the movement of suspensions, which are the starting material for the formation of bottom sediments.

In oceanology, particles of various origin, passively suspended in sea water and having sizes from 0.5 p,m to 1 mm, are considered to be suspensions. Particles may consist of organic and inorganic substances. Inorganic minerals consist mainly of clay minerals (silica, alumina, montmorillonite, illite, etc.) and non-clay minerals (quartz, mica, etc.). Organic materials can exist in the form of plants and bacteria.

The main factors of weighing, redistribution and transport of bottom material is the combined effect of waves and currents. A particle involved in the flow moves in the direction of the water flow and, under the influence of variables in magnitude and direction of pulsating velocities, simultaneously continuously makes vertical movements (rises - falls) (Fig. 1).

Fig. 1. Scheme of movement of suspended and entrained sediments

In the process of vertical mixing, the particle can sink to the bottom and remain on it until the moment when the lifting force is sufficient to detach it.

Suspension particles differ in size, density, area and, therefore, physico-chemical activity, different residence times in water and sedimentation rate. The mass of suspended solids in the volume of water is characterized by their concentration. The concentration of suspended particles is associated with seasonal factors and the hydrochemical regime of the reservoir, and also depends on anthropogenic factors.

Since it is very difficult to obtain data on the concentration of suspended matter in a pond from field data and in each case the available data may be incomplete, this problem often requires solving using modern methods of mathematical modeling.

Continuous 3Dmodelofdiffusion-convection-suspensionandthecorrespondinginitial-

boundaryvalueproblem. Let us consider a continuous mathematical model of suspended matter propagation in an aqueous medium, taking into account diffusion and convection of suspended

matter, the effect of gravity on suspension, the presence of a bottom and a free surface. We will use a rectangular Cartesian coordinate system Oxyz, where the axis Ox passes on undisturbed water surface and is directed towards the sea, the axis Oz directed vertically down. Let be h=H+q- total water depth, [m]; H - depth with undisturbed surface of the reservoir, [m]; q- elevation of the free surface relative to the geoid (sea level), [m].

Let in the closure area G={0<x<L ,0<y<L ,0<z<H(x,y)} there are suspended particles that are

at (x,y,z) and at the time t have a concentration c=c(x,y,z,t), [mg / l]; t— temporary variable, [s]. We will also use the notation L = max H(x,y).

z 0<x<L,0<y<Ly v v '

The behavior of suspended particles will be described by the following system of equations:

Qc+d(uc) | d(vc) | d((w+wg)c) ^ dt dx dy dz \dx2 dy

^d2c . d2c d f dc\ p

dz\ vdz) ' (1)

QH s

-=—w c,

dt p g

where u,v,w- vector components U fluid velocity, [m / s]; wg- hydraulic particle size or sedimentation rate, [m / s]; juk, /uv— coefficients of horizontal and vertical turbulent diffusion of

particles, respectively, [m2 / s]; F — power sources of particles; s— bottom porosity.

The terms on the left side (except for the time derivative) of the first equation of system (1) describe the advective transport of particles due to the inertial motion of the aquatic environment and sedimentation under the action of gravity. The terms on the right-hand side describe the diffusion of the suspension. The vertical diffusion coefficient is chosen different from the horizontal diffusion coefficient due to the fact that the effect of the difference of these coefficients is often observed in different environments and can be caused by various factors.

Add to system (1) the initial and boundary conditions (assuming that the deposition of particles on the bottom is irreversible).

As initial conditions at time t=0 accept

c(x,y,z,0)=c0(x,y,z); (2)

H(x,y,0)=HQ(x,y). (3)

We set boundary conditions on the edges. ABCDAlOClDl (set streams of suspensions both towards the coast, and along the coast)

-on the edges S1 =AAOB(x=0,0<y<L ,0<z<L) , S2 =AADp{y=Ly,0<x<L,0<z<L) and

S =BOCf (y=0,0<x<L ,0<z<L )

c=c*, where c* =c*(x,y,z,t),te[0,T]; (4)

-on the edges S4=DDff (x=L,0<y<L ,0<z<L) and S5=AOCfi, (z=0,0<x<L,0<y<L)

c=0; (5)

-on a surface S6=ABCD (z=H (x,y,t),0<x<L ,0<y<L )

de K do K _

(6)

dn u dz u

i v 1 v

The boundary condition (5) takes place with relatively small bottom slopes:

fdH , max I — I +

* Mi dx

2

dH

□ 1.

The following condition of nondegeneracy of the solution domain is set for all (x,y,t) under which the initial boundary problem is posed:

H(x,y,t)>h0 =const>0, 0<t<T. (7)

When studying combined sediment and sediment transport models, it is possible to increase the concentration of suspended particles in the bottom layer due to the rise of sediment particles when the shear stress exceeds a certain critical value [11]-[16]. Then instead of the boundary condition (6) we will consider the boundary condition of the form

3c

—=ac,a=const >0. (8)

3z

Linearization of the initial-boundary value problem of transport and sedimentation of suspensions.In order to create a linearized model on a time interval 0<t<T build a uniform grid wz

in steps t , ie many points coT={tn =m,n=0,1,...,N,Nt=T] .

Functions Cn)(x,y,ztn^) and H{n)(x,ytn_J we define at each step of the time grid wz. If a n=1, then as c®(x,y,z,t0), H(1)(x,y,t0) it is enough to take the functions of the initial condition, that is c(1)(x,y,z,0)=c0 (x, y,z), H(1)(x,y,t0)=H0(x,y) respectively. If n=2,...,N, then functions Cn)(x,y,z,tn_ 1)=Cn_1)(x,y,z,^_1) assumed to be known, since problem (1) - (6) for the previous time interval is assumed to be solved t <t<t ,.

n - 2 n-1

System (1) in the gap tn <t<tn we write in the form:

dCn) d{uQ">) d(ven) d((w+wg)Cnj_ fden deA df den

I+F

dt dx dy dz ^ dx2 dy2 ) dz\ v dz

dH(n) s ,, ■=—w e(n)

(9)

dt p

and supplement it with the initial conditions:

d1{x,y,z,t0)=e0(x,y,z),On{x,y,z,tJ^^^^^^x,y,zn n=2,. ..,N. (10)

m{x,y,t0)=Hlx,y)HnyO=Hn{x,yn n=2,...,N. (11)

The boundary conditions (4) - (6) are assumed to be satisfied for all time intervals t <t<t ,n=1,2,...,N.

n—1 n7 ' ' '

Defining the function e(n){x,y,z,tn ^)=On—1){x,y,z,tn ^ in the time interval tn_<t<tn, can find function H(n)(x,y,tn^). To this end, we integrate both sides of the second equation of system (9) over the variable t,t ,<t<t . Will get

7 n— 1 n '

'r dH" s \

¡z^-dt=—w [cm. (i2)

dt p gC

From equality (12) it is not difficult to get

H(n)=H("-1)-pw£ I cf")dt. (13)

p g "=k

We introduce at each time step tnl<t<t,n-l,2,...,Nregion

G— ={0<x<Lx,0<y<Ly,0<z<Hn-1( x,ytn — .

We have a chain of linear initial-boundary value problems for each time layer, where for the interval tn <i<tn, n-1,2,.. ,,/V view system is considered

dcn , d(udn)), d(vcn), d((w+wg)c<n>)_ ^ ^ f dO | dOl | d f | p ^ ^

dt dx dy dz ^ dx2 dy2 J dz\ v dz

(x,y,z)GG-i , Gn-1={0<x<L,0<y<Ly,0<z<Hn-1(x,y,tn_^

P N %

H(n)=H(n-1)-pw £ I cn)dt, n=1,2,... ,N. (15)

p g"=\-1

with initial conditions:

d"{x,y,ztn—-X){x,y,ztJ, ( 16)

H("{x,ytn-Hn)(x,y,tn- (17)

Note that at each time step, the boundary surfaces will change (except the face S5 ). Considering the time span tn_ <t<tn, we carry out the task of boundary conditions on the edges of the

region Gn-1:

-on the edges Sl^-l{x=0,0<y<Ly,0<z<H"-\0,ytn_)) , S^-Lfi^Lfi^H^L/- and SUy=0,0<x<Lc,0<z<H(n-1)(x,0,tn-1))

cn)=c', где c=c*(x^y,z,t),te[tn-1,tn] ; (18)

-on the edges S^(x=Lx,0<y<Ly,0<z<Hin-1)(Lx,y,tnA)) and Sn1(z=0,0<x<L,0<y<L)=ApcPi

c(n)=0 ; (19)

-on a surface S^(z=H(n-1)(x,ytn_),0<x<Lx,0<y<L„)

dc(n) w dé") w

dc =—gc" or —=—gc{n). (20)

dn dz ^

The boundary condition (8) will be replaced by the following

dc(n)

-=adn),a=const >0. (21)

dz y J

Thus, it is assumed that the bottom relief within a given time step when calculating the

distribution of concentrations of suspended matter does not change and is taken from the previous

time layer. Initially at this time step tn_<2<Xn,«=1,2,..^N the initial boundary value problem is solved

for the convection-diffusion equation (14) with a fixed bottom relief function H(n-1), and then an

update (recalculation) of the relief function is performed H(n) in accordance with equality (15). The determination of the conditions of existence, uniqueness and continuous dependence of the solution on the input data of the problem is carried out on a fixed time layer in these assumptions and subject to condition (7).

We will not investigate in this paper the existence of solutions of the initial boundary value problems (14) - (20) and (14) - (19), (21). Questions of the existence of solutions of initial-boundary value problems for equations of parabolic type with lower derivatives (diffusion-convection equations) are considered, for example, in monographs [17]-[18].

Investigation of the uniqueness of the solution of the initial-boundary problem of suspension transport.Consider the initial boundary value problem (14) - (20), formulated for an arbitrary time layer tn<t<tn.

Multiply the left and right side of equation (14) by the function Cn) and get:

ddn)

Cn)__i-c(n)

dt

f

d(uC")) djvd")) c((w+wg)c(n))

dx

Cy

dz

(d2c(n) d2Cn)) ,,d( dC") -+- '■ " 1

dx2 dy2

+Cn)—I ft— \+dn)F . (22) dz\ v dz

The left side of equality (22) can be transformed as follows:

dx dy ^ n ^ ^ '

dCn)

c(n)--i-c(n)

dt

dz \ 2 dt

2 dt 2 VV ' '

where U=||«,v,mh-vvJ .

In view of (23), equation (22) will be written as

(23)

2 dt

dl\({C'"jll)

>=AC<">

'd2Cn) d2c(n)^ -+-

dx2 dy2

d ( dCn)

+d,n"—I ft

dz\ v dz

|+c(n)F.

(24)

Then we integrate both sides of equation (24) over time on the interval tn_<t<tn and, after this, by spatial variables in the region G^. In the first term, the order of integration is changed by the Fubini theorem [19]. Will get

= i iif.fc(n)M

(d2c(n) d2c(n)N

dG

d ( dc(n)

dx2 dy

The first term on the left side of (25) is obviously equal to

dt +fj IG dt +fj fF)FdGL dt.

(25)

fff1i f ^Xg- 1=fjf^(cn))Xx,y,zl)-(cn))Xx,y,zl-1))dGn-1.

g-l2[tn-l dt

(26)

Next, we turn to the transformation of the second term of the left-hand side of equality (25). Taking into account the Ostrogradsky-Gauss formula and the boundary conditions (18) - (20), it can be written as [20]

i \t=\ i i 11 (c* %U*Jtyydz\t+\ } f JJ c2wdxdy\t+

/ V \ tf \ t , i \

+

}I ||(c)2(u*,7?)dxdz\it+-1f||(c)\u\njdxdz\lt~11 ||(c'fudydz Irff- (27)

'n-1 VS3,n-1 ) 2tn-1 \S2,n-1 ) 2 tn-1 VS1,n-1 )

A i t ( \ i t i \

1 " in / ,\2 , , | , 1 "

—1II \\(c')2vdxdz dt+111 \\(c)2vdxdz dt+1 |l \\(dn))2wdxdy

2tn-1

'n-1 VS3,n-1

2<n-1

n-1 VS2,n-1

2<n-1

dt

n-1 VS6,n-1

where ¿7* is the known velocity of the aquatic environment on the faces, where the boundary conditions of the first kind are specified; in fact, this is all the side faces, except S^ and top cover S5n_j, on which the suspension concentration is zero, and therefore the flows through them are zero. Let us turn to the transformation of the right side of (25). There is equality

Ill

Gn-=lll

An)

d fddn)\ q v—I — \+Mh—

V

f rs

dd"

hdx1 dx ) ' hdyl dy i

d f ddn)

)

V

dzV v dz

)j

dG =

n-1

d f . Qdn>) d

dc n)

_ dxV dx) dyV dy,

dc n)

-Ill

V

dx

+V

dc n)

V Qy )

+V

d f ,, dd" +—| 6n)v —

dzV v dz

dc n)

dG -

n-1

(28)

dz

dG ,.

n-1

Let be

dcB> , dcn) I

theorem, we have:

lll

dx

d f , ,ddn) <n)-

Vh—I c

dy

d

cv

\+Vu

An),

dz

d

. Then, by virtue of the Ostrogradsky-Gauss

f—I cn)v

dc(n)

hdx\ dx) hdy| dy \ dzV ' v dz = ll Qdxdz+ ll Qdydz + ll Qdxdz+ ll Qdydz +1

dGn_=\\\divQdG= ' Qdxdy+ ll Qdxdy = (29)

S3,n-1 Sl>n-1

= ll Qdxdz+ ll Qdxdz + ll Qdydz+ ll Qdxdy.

S2,n-1 S3,n-1 S1,n-1 S6n-1

Transforming each term from the right-hand side of (29) subject to the conditions on the boundary (18) - (20), we obtain

lll

Vh

df ,,ddn)

\+vh—

d f , dd"'! d

dxV dx) dy| dy \ dz

A").

|+—| cn)v

ddn) ~dz

dG =

= ll c*u,—dxdz+ ll c'v—dxdz + ll c'v—dydz- ll w i&^dxdy.

¿1 V dy s|j-11 Vh dy & Vh dx y s!,j-1 g ( 7 y

Taking into account (26), (28), (29) and (30), equality (25) takes the form

(30)

n-1

n-1

G

n-1

1JJJW (wt^r J(JJ(^cf u^yydz^t -

'"( il 2 de'} } >-( f\ 2 Qc*\ A -J JJ -(c*)v+c'fth— dxdz dt +| JJ -{c*) v-c*fth— dxdz

J J

+3 JI JJ w^dxdy dt + J

JJJ

ôàn) dx

V "V-v-

-V 2V

v dy j

dc^

JJ

\

dG

dt+ dt=

(31)

- -JJJ( c(n ))2( x^y,z,tn-i)dGn-i+J(JJJc(n)FdGn-i

2G„ ,

dt.

Identity (31) will be fundamental in the study of uniqueness and obtaining an a priori estimate of the norm for the solution of the initial-boundary value problem (14) - (20). In the case of replacing the boundary condition (20) with the boundary condition (21), the quadratic functional (31) changes as follows:

- Jp)2 (x,y,z,0dG-1-J(fJ(1 (cf u+c*fth d-^ydz jdt -

-JIJJI ^o2 v+c! ^xdzU +J | jj[ V)2 v-c*fth ^

A

tn-1 V S3,n-1

!dy,

tn-1VS6,n-1V 2

+ JI JJ^w-a!v Xc{n])2dxdy Idt +J

JJJ

'''dcMV (

dx

'dy j

dc(n)|2 ( dc(n)

dt+

dy

J

2

dG

(32)

dt=

=1JJJ(c(n ))2(x,y,z,tn-1)dGn-1+ J f JJJc(n^FdGn1 Jdt.

2Gn-1 tn-V Gn-1 J

Suppose that equation (14) with the same conditions (16) - (20) satisfies two different solutions of the problem c-c(x,y,z,t),c=cix,}\z,t). For their difference c=c-c2 The following initial boundary problem is valid:

fo | c(i,c) | a(vc) |

dt dx dy dz h

<92c <92c —+—

dx2 dy2

d

dz

f & M —

' v --n

v dzJ

c(x,v,z,0)-0,(x,v,z)eG„-i,

on the edges S ,, S ,, S ,, S. ,,S, ,

© 1,n-1~ 2,n-1~ 3,n-1~ 4,n-1~ 5,n-1

on a surface S

c=c*-c*=0 :

T-=—J\c\~cd=—c ■ dz n n

(33)

(34)

(35)

(36)

For function c equality (33) will take the form in view of equalities (34) - (36)

+J

JJJ

^jfl+M

rdc< vdy J

, ÔC

dG

(37)

dt=0.

Insofar as w >0 and other known quantities under the sign of integrals are positive ^>0, /U >0, then equality (36) is satisfied only under the condition

c(x,y,z,t)=0, (x,y,z)&Gn l, tn l<t<tn, (3 8)

which completes the proof of the uniqueness of the solution of the initial-boundary value problem (14) - (20).

In the case of replacing the boundary condition (20) by the relation (21), instead of the expression (37), we obtain the equality

[\wraV, Y'dxdy

+

sss

dc ,

rde^

\dy J

+ju

dc ~dz

dG

dt+ dt-0.

Demandfulfillmentofinequality

^w -au >0, (x,y,z)eS ,, t <t<t

2 g ' v ' \ ' / 6,n-1 n-1 n

(39)

or

w ,

a<^, (w)^ t„-1<t <„,

2U

(40)

then all the terms in equation (39) are non-negative and equality to zero is possible if and only if c(x, y,z,/)=(), (x, y,z)eG p tn l<t<in, which means the uniqueness of the solution in this case too.

Reasoning is similarly repeated for all layers of the time grid mz. The change in the boundary

conditions associated with the continuous change of the bottom relief depending on the time variable requires additional research and is beyond the scope of this article.

Theorem. Let given a system of equations

ÔÇn) ô(uc(n) I ô (vc(n)) I cj^w+wy^ 11 fd2c(n) d2Çn)\i d(_ ( dc^

dt ôx

ôy

ôz

-=u

ôx2 ôy2 J ôz\ v ôz

+F,

(x,y,z)sGn-i , Gn-i={0<x<Lx,0<y<Ly,0<z^^y/Jl

F N n

H(n)=H(n-i)-Lwg£ j Cn)dt, n=1,2,.. ,N

P

in a simply connected domain Qn =Gn lx(tn l<t<tJ, Gn=(0<x<Lx ,0<y<Ly,0<z<H{n l)(x,y,tn ~)), with a fairly smooth boundary defined by the smoothness of the function z=H(n-1)(x,y), 0<x<Lx, 0<y<L with initial and boundary conditions (16) - (20). Let solution functions <Cn)(x,y,z,tnJ , water velocity

vector |U,v,w+w||r, initial condition Cn-1)(x,y,z,tn_~), right side F(x,y,z,t),boundary condition c*(x,y,z,t), coefficient of vertical turbulent exchange u=u(z),(x,y,z)^Gi satisfy the following smoothness conditions: ön)(x,y,z,tn_)eC2(QnJnC(Q J, gradCn)^C(Ùn ),

Uvw+wJ eO(Qn-1 Mi^) C^xy^,KG ), F(x,yzt)eC(Oj, ^(w^C(Gnl)oC(G ,),

cixyz^CsM^«J, Sn-l=Gn-l\Gn-l, f^C((0<x<Zx,0<y^ z=^{n-1)(x^y))x[tn-l<t<tn]), as wc11 as the conditions of consistency of the boundary and initial conditions, c*(x,y,z,0)=c0(x,y,z),

(x,y,z)eSn-1\(0<x<Z,0<y<L,z-H^xy)), ^=-!*, (0<x<L,,0<y<L ,z=Hn1(x,y)), then the

g

solution to this problem exists and is unique.

Comment. If the boundary condition (20) is replaced by the boundary condition (21), inequality (40) should be added as a sufficient condition for the fulfillment of the previous theorem.

Investigation of the continuous dependence of the solutions of the initial-boundary value problem of suspension transport on the initial, boundary conditions and the function of the right-hand side.The next stage is connected with the study of the continuous dependence of the solution on the functions of the right-hand side, the boundary and initial conditions for the system (14) - (15).

We will assume that

c'>c-—const >0,

0<x<L ,0<y<L ,0<z<H(n)(x,y,t It <t<t. (41)

x7 s y7 V 7*/ 7 n-1/7 n-1 n

For convenience, we introduce the following notation: the union of all parts of the lateral surface the area borders G , denote as S ,, and the lower base of the area G - S ,. By virtue of

n-1 c,n-V n-1 b,n-1 J

the smoothness conditions listed under the conditions of the Theorem, extremes of functions on bounded closed sets are reached. :

M1,n-1-max-fl^ } M2,n-1-

max<!

dcM dcM

dx dy

(42)

i-i 2,n-1 S„_

M3,n-1-mS a_lx{!h}, M4,n-1-Sn ]{lul,H}, M5,n-1-nGin{!h,!v}.

We will focus on equation (31) if the boundary condition (20) is used, and equality (32) in the case of the boundary condition (21). Drawing on the Friedrichs inequality, we have a chain of inequalities:

JJJ

!

v

dc(n) dx

( dc(nN dy y

dc(n)

Yz

A

dG >

<m in{!h,!v}JJJ

G„-( f

((dc(n)I2 (dcwY (dc^Vl

dx

dG >

n-1

>M

5,n-1

nL

V V

1 1

f—+■

1

dy J V dz

JJJ^G1.

(43)

JJ

L2 L2 (H(n-1))2

V x y v ' .

Let us turn to equation (26), from which, by virtue of (42) and (43), we obtain the inequality:

G

n-1

JJJ( c(n 02 ¿GLi+2M5,b_I

V V

\

1 1

+—4

\\

L2 L2 (H(n-1))

JJ

1 IE

dt+

+311 If w (c(n))dxdy dt<JJJc&GnA+M4"-1 |l JJ(cfdSn-1

dt+

(44)

tn-1 VSb,"-1

+2m2,„-^m3,„-1 JI JJ d dS-1 dt+M-1 MWGM

t„-\Sc,"-1 J

From inequality (44) there are two inequalities

M)

t„A G„-1

JJJ(d"))2dG„-1<JJJc02dG„-1+M4,„-1 Jl JJ(c*)2dS"- ldt+

Gn-1 Gn-1 <n-M,"-1 J

+2M2,„-^M3,„-1 Jf JJ ddS"-1V+M JifflFdG„-1^t.

and

JJJ(c("))2dG„-1<M6,„-1fJJJc02dG„-1+ J f JJ {cidSXt+M-M^ Jf JJ WS^t+M-1 'JfJJJWGn-1 dt l

tn-1V Sc,n-1 J 'n-1VSc,n-1 J '"-1V G"-1

(45)

(46)

where M ,-—1

6,n-1 2Mr

w

1 1

■+—+■

\\

v v

l2 l 2 hn-1))'

jj

From the obtained inequalities, it follows the continuous dependence (stability) of the solution of problem (14) - (20) on the functions: the initial condition, the boundary conditions and the right-hand side L for any point in time 0<T<+<», as well as in time integral L2.

Obviously, with the fulfillment of inequality (45) and the conditions of the Theorem, the initial boundary value problem (14) - (19), (20) will also have a solution that depends continuously on the functions: the initial condition, the boundary conditions and the right side in the corresponding norms.

Description of the parallel algorithm. A software package implemented in C++ has been developed for constructing turbulent flows of an incompressible velocity field of the water medium on high-resolution grids for predicting sediment transport and possible scenarios for changing the geometry of the bottom area of shallow water reservoirs. Parallel algorithms implemented in the software package for solving model problems arising in the process of sampling systems of grid equations were developed using MPI technology.

To solve this problem, we used an adaptive modified alternately triangular method of minimal corrections. In parallel implementation, the methods of decomposition of grid areas for computationally time-consuming diffusion-convection problems are used, taking into account the architecture and parameters of a multiprocessor computer system. The decomposition of the calculated two-dimensional region is performed using two spatial variables . The peak performance of a multiprocessor computer system is 18.8 teraflops. As computing nodes are used 128 similar 16-core Blade servers HP ProLiant BL685c, each of which is equipped with four 4-core processors AMD Opteron 8356 2.3 GHz and RAM in the amount of 32GB.

n-1

G

n-1

The software package is used to calculate the geometry of the bottom of shallow water reservoirs and includes the following blocks:

- control unit, it contains a cycle on the time coordinate and calls the functions: calculation of the speed field without taking into account pressure, calculation of the elevation function, speed field refinement, calculation of the depth field and data input / output functions;

- a block for constructing grid equations for the velocity field without taking into account the pressure in accordance with the finite-difference scheme. The coefficients and the right part of the corresponding grid equation presented in canonical form are considered and written to the array;

- block for constructing grid equations for calculating the elevation function;

- block for calculating the velocity field taking into account the elevation of the level (result of this block is the calculation of the values of the velocity vector field on the next time layer);

- block for calculating grid equations by adaptive modified alternately-triangular method of rapid descent;

- block of output values of the velocity field in the file.

Numerical experiments modeling the sediment transport and the dynamics of changes of the bottom topography. After the development of the software package, a series of numerical experiments was performed to simulate the dynamics of changes in the bottom relief of a complex configuration in the coastal zone of the reservoir. The model problems assumed the presence of obstacles on the bottom surface (boulders, underwater breakwaters, breakwaters, dumps, jetties, spurs, etc.) and various irregularities underlying its surface. As an example, the paper presents the results of modeling the dynamics of changes in the bottom for the case when there are obstacles on its surface in the form of pointed structures-intermittent bun. Due to the retention of sediment, boons not only stop the movement of the material carried by the waves along the shore, but also contribute to its deposition. These structures are one of the best means to protect the coast and prevent the invasion of the sea to the mainland.

Fig.2. The geometry of the calculated area at the initial moment of modeling (a-the position of the isolines of the depth function, and-the bottom relief)

The simulation area under consideration has dimensions of 55 m by 55 m horizontally and 2 m vertically (in depth), the peak point rises above sea level up to 1 m.Assume that the liquid is at rest at the initial moment of time. The size of the calculated grid is equal to 110by 110, the step on spatial variables is 0.05 m, the step on time is 0.01 s, the wind speed is 5 m/s and is directed from left to right.

Fig.2 shows the initial position of the isolines of the function of depths and bottom relief, the feature of which is the presence of three bottom bun. Fluctuations in the isoline of the depth function are observed in the Central part of the calculation. These structures have a length of up to 15 m and are located at a distance not exceeding 10 m from each other. The structures are completely submerged in the reservoir, and their maximum height is 1.25 m.

Modeling the process of sediment transport showed that over time there is a «smoothing» of the roughness on the surface, the formation of sediments, a decrease in the depth of the slope of the bottom of the coastal zone and, as a result, a gradual shallowing of the considered area of the reservoir. So, after 5 minutes after the beginning of modeling the isolines of the depth function in the center of the calculated area became undulating, and in the area of the location of the boon-sinuous.

Fig.3.The geometry of the calculated area 5 minutes after the start of the simulation (a - the position of the isolines of the depth function, b - the bottom relief)

S$-

-

4S-

•I -

K

.' - -

: -

19-

o-

0

Fig. 4. The geometry of the calculated area 15 minutes after the start of the simulation (a - the position of the isolines of the depth function, b - the bottom relief)

Sediment deposition occurs in the interbunar compartments, and the pointed peaks of the boon have deformed and taken the form of gentle hills. As a result of these processes, there is a decrease in the depth of the coastal zone and the increase in the beach area (Fig. 3).

More clearly, the results are presented in Fig. 4, when the simulation time was 15 min. the isolines of the depth function acquire a soft wavelike shape over the entire calculation area, including the areas of peak values. There was an active process of long-distance movement of sediments and lowering the depth level. So, during the specified estimated time in the interbunny compartments, the depth decrease was about 0.5 m. the height of the slides is reduced, and the slides themselves become more «smoothed» appearance.

Conclusion.The paper presents a software package for predicting possible scenarios for changing the geometry of the bottom of shallow water reservoirs using high performance computing. Numerical experiments were performed for model problems of bottom sediment transport and bottom relief transformation, the results of which are consistent with real physical experiments. The proposed mathematical model and the developed software package allow us to predict the dynamics of the behavior of the bottom surface, the appearance of sea braids and ridges, their growth and transformation.

References

1. Leontyev, I.O.: Coastal Dynamics: Waves, Moving Streams, Deposits Drifts. GEOS, San Moscow (2001) (in Russian).

2. Liu,X., Qi,S., Huang,Y., ChenY., Du,P.: Predictive modeling in sediment transportation across multiple spatial scales in the Jialing River Basin of China. International Journal of Sediment Research, 30:3, pp. 250-255 (2015).

3. Aksoy, H., Kavvas, M.L.: A review of hillslope and watershed scale erosion andsediment transport models. Catena, 64:2-3, pp. 247-271 (2005).doi:http://dx.doi.org/10.1016/j.catena.2005.08.008.

4. Ouda, M., Toorman, E.A.: Development of a new multiphase sediment transport model for free surface flows. International Journal of Multiphase Flow, 117, pp. 81-102 (2019).

5. Sukhinov, A.I., Sukhinov, A.A.: Reconstruction of 2001 Ecological Disaster in the Azov Sea on the Basis of Precise Hydrophysics Models. Parallel Computational Fluid Dynamics, pp. 231-238, (2005). doi: 10.1016/B978-044452024-1/50030-0.

6. Sukhinov A.I., Sukhinov A.A. Reconstruction of 2001 Ecological Disaster in the Azov Sea on the Basis of Precise Hydrophysics Models. Parallel Computational Fluid Dynamics, Multidisciplinary Applications, Proceedings of Parallel CFD 2004 Conference, Las Palmas de Gran Canaria, Spain, ELSEVIER, AmsterdamBerlin-London-New York-Tokyo, pp. 231-238, (2005). doi: 10.1016/B978-044452024- 1/50030-0.

7. Alekseenko, E., Roux, B., Sukhinov, A., Kotarba, R., Fougere, D.: Coastal hydrodynamics in a windy lagoon. Nonlinear Processes in Geophysics, 20:2, pp. 189-198 (2013). doi: 10.1016/j.compfluid.2013.02.003.

8. Alekseenko, E., Roux, B., Sukhinov, A., Kotarba, R., Fougere, D.: Nonlinear hydrodynamics in a mediterranean lagoon. Computational Mathematics and Mathematical Physics, 57:6, pp. 978-994 (2017). doi: 10.5194/npg-20-189-2013.

9. Francke, T., Lopez-Tarazon J.A., Vericat, D. Bronstert, A., Batalla, R.J.: Flood-based analysis of high-magnitude sediment transport using a non-parametric method. Earth Surface Processes and Landforms, 33:13, pp. 2064-2077 (2008). doi: http://doi.org/10.1002/esp.1654.

10. Karaushev, A.N.: Theory and Methods for River Load Calculation, Leningrad: Gidrometeoizdat, (1977) (in Russian).

11. Alekseevskii, N.I.: Hydrophysics, Moscow: Akademiya, (2006) (in Russian).

12. Goloviznin V.M., Chetverushkin B. N.: New generation algorithms for computational fluid dynamics. Computational Mathematics and Mathematical Physics, 58:8, pp. 1217-1225 (2018). doi: https://doi.org/10.1134/S0965542518080079.

13. Sukhinov, A.I., Chistyakov, A.E., Protsenko, E.A.: Mathematical modeling of sediment transport in the coastal zone of shallow reservoirs. Mathematical Models and Computer Simulations, 6:4, pp.351-363 (2014). doi: 10.1134/S2070048214040097.

14. Sidoryakina, V.V., Sukhinov, A.I.: Well-posedness analysis and numerical implementation of a linearized two-dimensional bottom sediment transport problem. Computational Mathematics and Mathematical Physics, 57:6, pp. 978-994 (2017). doi: 10.7868/ S0044466917060138.

15. Sukhinov, A.I., Sidoryakina, V.V.: Convergence of linearized sequence tasks to the nonlinear sediment transport task solution. Matematicheskoemodelirovanie, 29:11, pp.19-39 (2017).

16. Belotserkovskii, O.M., Gushchin, V.A., Shchennikov, V.V.: Decomposition method applied to the solution of problems of viscous incompressible fluid dynamics. Computational Mathematics and Mathematical Physics, 15, pp. 197-207 (1975).

17. Favorskaya, A.V., Petrov, I.B.: Numerical modeling of dynamic wave effects in rock masses. Doklady Mathematics, 95:3, pp. 287-290 (2017). doi: 10.1134/S1064562417030139.

18. Sukhinov, A.I., Chistyakov, A.E., Shishenya, A.V.: Error estimate for diffusion equations solved by schemes with weights. Mathematical Models and Computer Simulations, 6:3, pp. 324-331(2014).doi: https://doi.org/10.1134/S2070048214030120.

Authors:

Sukhinov Alexander Ivanovich, Don State Technical University (1st Gagarin Square, Rostov-on-Don, Russian Federation), Doctor of Science in Physics and Maths, Professor

Chistyakov Alexander Evgenievich, Don State Technical University (1st Gagarin Square, Rostov-on-Don, Russian Federation), Doctor of Science in Physics and Maths, Associate professor

Sidoryakina Valentina Vladimirovna, Taganrog Institute of A.P. Chekhov (branch) RSUE (Initiative Street, Taganrog, Russian Federation), Candidate of Science in Physics and Maths, Associate professor

ProtsenkoElena Anatolevna, Taganrog Institute of A.P. Chekhov (branch) RSUE (Initiative Street, Taganrog, Russian Federation), Candidate of Science in Physics and Maths, Associate professor

ProtsenkoSofyaVladinirovna, Don State Technical University (1st Gagarin Square, Rostov-on-Don, Russian Federation), postgraduate student

УДК 519.6 10.23947/2587-8999-2019-2-2-118-135

Программный комплекс прогнозирования возможных сценариев изменения геометрии дна мелководных водоемов с использованием высокопроизводительных вычислений* А.И. Сухинов1, А.Е. Чистяков1, В.В. Сидорякина2, Е.А. Проценко2, С.В. Проценко1

1 Донской государственный технический университет, Ростов-на-Дону, РоссийскаяФедерация 2Таганрогский институт имени А.П. Чехова (филиал) Ростовского государственного экономического университета (РИНХ), Таганрог, Российская Федерация

Настоящая работа посвящена изучению модели переноса и осаждения взвешенных веществ в прибрежной зоне. Комплексные исследования водных экосистем являются неотъемлемой частью рационального природопользования. В этих исследованиях важное место занимает изучение перемещения взвесей, являющихся исходным материалом для образования донных осадков.Основными факторами взвешивания, перераспределения и транспорта донного материала является комбинированное воздействие волн и течений. Частица, вовлеченная внутрь потока, движется в направлении водного потока и, находясь под воздействием переменных по величине и направлению пульсационных скоростей, одновременно непрерывно совершает вертикальные движения (поднимается -опускается).Модель учитывает следующие процессы: адвективный перенос, обусловленный движением водной среды, микротурбулентную диффузию и гравитационное осаждение частиц взвеси, а также изменение геометрии дна, вызванное осаждением частиц взвеси или подъемом частиц донных отложений.В статье представлены результаты исследования корректности начально-краевой задачи, соответствующей построенной модели. Разработан программный пакет для прогнозирования возможных сценариев изменения геометрии дна водоемов на мелководье с использованием высокопроизводительных вычислений.

Ключевыеслова:распределенные вычисления, высокопроизводительные вычисления, параллельное программирование, математическая модель, динамика морских наносов, рельеф дна.

*The research is done with the financial support from Russian Foundation for Basic Research, Project No. 19-01-00701.

Авторы:

Сухинов Александр Иванович, Донской государственный технический университет (344000 Ростов-на-Дону, пл. Гагарина, д. 1), доктор физико-математических наук, профессор

Чистяков Александр Евгеньевич, Донской государственный технический университет (344000 Ростов-на-Дону, пл. Гагарина, д. 1), доктор физико-математических наук, доцент

Проценко Елена Анатольевна, Таганрогский институт им. А.П. Чехова (филиал) РГЭУ (РИНЭ) (347936 Таганрог, улица Инициативная, д. 48), кандидат физико-математических наук

Сидорякина Валентина Владимировна, Таганрогский институт им. А.П. Чехова (филиал) РГЭУ (РИНЭ) (347936 Таганрог, улица Инициативная, д. 48), кандидат физико-математических наук, доцент

Проценко Софья Владимировна, Донской государственный технический университет (344000 Ростов-на-Дону, пл. Гагарина, д. 1), аспирант