Two-dimensional-one-dimensional splitting scheme for the numerical solution of problems of transport of multicomponent suspensions using θ coordinates

. The work considers the spatial-three-dimensional problem of suspension transport, which takes into account many parameters and processes (multicomponent fractional composition of suspension, particle sedimentation rate, suspension distribution, intensity of sources of suspension distribution, etc.). For the basic equation of this problem, a symmetric representation of the convective terms is used, which makes it possible to ensure the unconditional skew symmetry of the convective transport operator. The methodology for constructing additive circuits (splitting circuits) is outlined, which makes it possible to reduce the solution of the original problem to a sequential (or parallel for multiprocessor computers) solution of two-dimensional and one-dimensional analogues. The feasibility of this approach is determined, first of all, by the significant differences in space-time scales for the difference operators of diffusion transfer in the horizontal and vertical directions, as well as their spectra. The constructed splitting scheme is monotonic, and each of the difference equations is non-degenerate if the grid Peclet number is less than one.


Introduction
The article discusses a 3D model of diffusion-convection-sedimentation of suspended matter in an aquatic environment [1] and its numerical implementation.This model has a coordinate in the vertical direction, which is a set of topologically identical layers of constant or variable thickness.Due to the dependence of the -coordinate on depth, the model makes it possible to study the process of resuspension in the shallow coastal zone.
The numerical implementation of the model is based on additive splitting schemes along geometric directions.Real three-dimensional turbulence in coastal systems can be considered locally homogeneous and isotropic only on scales significantly smaller than the depth of the water area, since the stratification of water in a well-mixed shallow reservoir is limited by the maximum vertical size of turbulent eddies [2][3][4][5][6].It is believed that for the latter reason, the process of turbulent mixing can be divided into horizontal turbulent exchange and vertical turbulent diffusion.Operators of convective and diffusion transport in horizontal and vertical directions for problems of suspended matter transport have significantly different physical and spectral properties.In connection with the above, to construct economical algorithms for the operational forecast of suspended matter distribution, it is advisable to use two-dimensional-one-dimensional splitting schemes [7].It was previously shown for spatial-three-dimensional initial-boundary value problems for equations of parabolic type that the use of these schemes makes it possible to obtain algorithms for solving parabolic mesh equations that are low-cost in terms of the number of arithmetic operations, which are economical in terms of the number of operations and the time required to carry out exchanges between processors [8][9][10].

Initial-boundary value problem of suspension transport
We will use a rectangular Cartesian coordinate system Oxyz , where is the axis Ох passes along the surface of the undisturbed water surface and is directed towards the sea.Let hH  is the total depth of the water area; H is the depth when the surface of the reservoir is undisturbed;  is the elevation of the free surface relative to the geoid.
In the model, we will carry out the transition from the z -coordinate system to the  - coordinate system, for which we use a cartesian coordinate system in the horizontal plane, and as a vertical variable we use a dimensionless variable   , 0;1.

 
In the  -coordinate system, the water column is divided into the same number of layers at each point, regardless of depth, so using the "new" coordinate system solves some of the problems associated with adding and subtracting layers [11].
When moving to the  -coordinate system, we use the formula: where 0 a   on the free surface of the reservoir (upper boundary), In what follows, instead of expressions (1), we will use

G
is the area in which the process occurs, and which is a parallelepiped Let there be a suspension of a multifractional composition in the region G .To simplify the calculations, we will consider a suspension consisting of three fractions (the problem statement for the general case is presented, for example, in [1,11]).The restriction imposed on the number of fractions does not affect the general idea of obtaining the result presented in the work.
Let us write a system of equations that describes the concentration of individual suspended fractions.We use an equation in which the term describing the advective transport of suspended particles is presented in the so-called symmetric form: .
Here and further in the text where , rr  are the coefficients determining the intensity of transformation of particle of the type r into type    The boundary and initial conditions for equations (1) are formulated as follows: where n u is the projection of the velocity vector onto the outer normal n to the boundary, r c is known concentration values; : 0; :. r b r gr r c wc 3 Construction of a difference scheme To equation ( 8) we add initial conditions of the form: ..., , , , , as well as boundary conditions similar to conditions ( 5) - (7).For all t, : .
Further, the set of internal grid nodes ,, xy    we will denote accordingly as , , , xy    we will denote their cartesian products as On a space-time grid    let us approximate problems ( 8) -( 12) on grids adopted in computational fluid dynamics with setting velocities (convective transfer coefficients) at grid nodes shifted by half a step: for u , grid After introducing the designation nn gr w w w   , let's construct an additive unconditionally stable splitting scheme in geometric directions, focusing on a chain of continuous problems of the formtwo-dimensional and one-dimensional [7,12]: , ( 1), ( 1),   The initial and boundary conditions for equations (13) -(15) logically follow from conditions ( 9) - (12).

w x y h c x y h w x y h hh c x y h c x y ab c x y h x y h h h h c x y c x y h x y h h
, , , , 1,..., .
To the system of difference equations it is necessary to add initial conditions for a chain of initial-boundary value problems of the form (8) for   ,, xy  , as well as approximation of boundary conditions ( 9) - (12).
In addition, we will assume that the values of the velocity vector components of the aquatic environment and the hydraulic size of suspended particles at the grid nodes are known.\

 
with fractional index values: If, with respect to the external normal n the components   ,, u v w of the velocity vector is less than zero, then the boundary conditions (10) are approximated for the surface l  in the following way: (1), ,5 , 0, , , .
Considering the component part of the surface l  , which is formed by points,   0, 0 , 0 1 From ( 25) and ( 26) we obtain Arguing in a similar way, we get For the surface f  , we use the Neumann boundary conditions (11): For surface b  we will use boundary conditions (12) of the third kind when approximating a one-dimensional (vertical) difference problem.We have its difference approximation with an error   On the other hand, we have the formal equality  x y h   from equality (33), we have

C c w x y h c x y h w x y h c x y h h
and, substituting it into equality (33), we obtain Let us move on to approximating the diffusion terms at the boundary of the areas G .
Let 0 x .Then formally on the extended grid   can be written down Taking into account the Neumann condition at the boundary 0 x if the first of conditions ( 24) is satisfied, we can write , , , , .Then from (37) we obtain , (1), 2 0 1 0,5 , ,
The constructed splitting scheme is monotonic, and each of the difference equations is non-degenerate if the grid Peclet number is less than one.

Conclusions
As a result of constructing a scheme for splitting the original problem into two-dimensional and one-dimensional (additive schemes) in coordinates, we come to the need to solve two- dimensional and one-dimensional problems of diffusion-convection that are identical from the point of view of geometry.Each of these problems has a smaller condition number compared to the direct approximation (three-dimensional difference scheme).In addition, with parallel numerical implementation and the use of the Domain Decomposition approach, these schemes lead to savings in time spent on performing exchanges.
The study of the stability and convergence of the constructed difference scheme, as well as its numerical implementation, including a parallel algorithm, is the subject of future work by the authors.

1 r
 and type   1 r  , respectively, 0, 0 rr   ; r  is the power of the external particle source of the type r .

F
will mean that they belong to the class of grid functions.
which are outside the reservoir, the values of the velocity vector components are assumed to be zero.Let us present separately the approximations of the convective terms for the nodes on each of the surfaces ,, l b f    area G .

E3S 2 x
Web of Conferences 458, 03019 (2023) EMMFT-2023 https://doi.org/10.1051/e3sconf/202345803019which can be considered a difference approximation of the convective term at 0 x with an error   Oh along the axis Ox .Along with (25), we can write a difference approximation of the Neumann boundary condition with an error   is no turbulent diffusion outside the reservoir), so one should accept f  condition (12) must be taken into account.Let us formally write the expression for the diffusion term in difference form: