{'fraction_non_alphanumeric': 0.022071307300509338, 'fraction_numerical': 0.002865025466893039, 'mean_word_length': 5.03393085787452, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Hyperbolic decomposition for hydrostatic approximation of free surface flow problems

The article is devoted to the description of a new method of numerical solution for hydrostatic approximation of incompressible hydrodynamic problems with free surfaces and variable density. This technique develops a model of multilayer shallow water, taking into account the mass and momentum exchange between the layers.


Introduction
Multilayer shallow water equations [1] are not always hyperbolic [2], which leads to the fact that initial -boundary value problem becomes incorrect in the process of solution , and the computational algorithms polished for always hyperbolic equations lose stability. The loss of hyperbolicity is believed to occur in situations when conditions for the development of Kelvin -Helmholtz arise instability at the interfaces between layers. In these cases, an intensive mass and momentum exchange between the layers should occur, which is forbidden in the classical multilayer approximation. This disadvantage can be eliminated by including so-called "turbulent" viscosity in the multilayer equations. This viscosity depends in a complex way on the parameters of calculated flow and contains empirical parameters that can be adjusted to different types of flows. The disadvantage of this approach is the critical dependence of the calculation results on the experience of the researcher.
Another approach to solving this problem is to modify the original model based on the derivation of the classical multilayer shallow water equations. We remove the ban on exchange of mass and momentum between the layers [3]. In fact, this means a transition from the Lagrangian description of the motion of boundaries of layers in the vertical direction to the mixed Eulerian -Lagrangian description. In this case, the regularizing mass and momentum flows arise between the layers.
The implementation of this approach uses the method of splitting into physical processes. First, the values of physical variables on the new time step are calculated using the classical equations of multilayer shallow water with variable density. Then, the prescribed vertical coordinates of the interfaces between layers are determined in some way. Finally, these coordinates determine the values of exchange fluxes of mass and momentum between the layers. When calculating the prescribed vertical coordinates of the boundaries of the layers, following considerations are used: the uppermost layer, which is the free boundary of the fluid, remains Lagrangian, the lowermost one is fixed. If the coordinates of intermediate boundaries are found using the condition of equality of the thicknesses of the layers on each vertical line, then such computational grids are called "sigma -grids". Another class of grids often used in computational oceanology is the so-called "Z-grids", when the boundaries of all layers, except the uppermost one, remain fixed. Intermediate options are also possible.
Many different computational algorithms with their own advantages and disadvantages have been developed for the numerical solution of the classical equations of multilayer shallow water with variable density. Usually, computational algorithms are based on the finite volume method, and some variant of downstream transferring of the values is used to calculate convective flows.
It should be noted that high-resolution methods based on solving the Riemann problem and well-proven in the one-layer case are not applicable for the multilayer case, since the Riemann problem in the multilayer case does not have a simple solution.
In the new numerical method, the values of horizontal flows are proposed to be calculated by the method of characteristics, as is done in the CABARET scheme [4,5]. For this, the multilayer model is represented as a sum of one-layer models (hyperbolic decomposition) interacting through reaction forces applied to the interfaces of the layers. The CABARET scheme is used for each of the one-layer models, taking into account the nonlinear correction of the fluxes based on the maximum principle. The condition of incompressibility and the conservation laws of mass and momentum are also satisfied. As a result, conservative and flux variables are found at the new time step without the exchange of mass and momentum between the layers. Since this problem is incorrect, new flux variables are filtered, not violating the conservation laws.
The resulting algorithm has a second order of approximation and maintains a state of rest at a constant density over any bottom topography (well-balanced algorithm). At the second stage of splitting the prescribed vertical coordinates of the layer boundaries are specified according to physical processes. Using new vertical coordinates, the exchange flows between the layers are calculated with second order of approximation.
The properties of the new technique are illustrated on the test problem.
2. The equations of multilayer shallow water with variable density and considering external pressure and bottom relief The equations of multilayer shallow water with variable density without mass and momentum exchange between the layers can be represented as the sum of one-layer equations related to each other by the forces of pressure and reaction acting on their interfaces: here k is layer number counted from the free surface, Z k is coordinates of the upper boundary of k th layer, P 1 = P T is pressure on the free surface Z 1 . The so-called "simple form" of this system appears as follows: where G is matrix of dimension N × N and D -the right side, the specific form of which does not matter to us. It is well known that, even when N = 2, the matrix G can have complex roots  3 and the system is not always hyperbolic. This creates known computational difficulties. For a larger number of layers, the situation is only getting worse. To overcome this difficulty, we can apply the technique we call the "hyperbolic decomposition" of the problem.
If we consider the forces acting on the interface as "external", then this system can be represented as: This leads to a system of independent characteristic equations: Therefore, it is possible to use for the numerical solution of the complete system of multilayer equations algorithms that are well proven in the one-layer case. The incorrectness of the total system does not disappear and arises in the form of high-frequency distortion of the boundaries. They can be regularized by using special filters that do not violate the approximation and lead to the exchange of mass and momentum between the layers, which does not violate the conservation laws.

Features of numerical implementation
The CABARET scheme uses two types of variables: flux-type and conservative-type. The fluxtype variables belong to the middle of the vertical edges of the calculated cells, the conservativetype variables belong to their centers. The computational algorithm consists of the following elements: • the calculation of conservative-type variables at an intermediate time step using conservative difference schemes obtained by the finite volume method (Phase 1) • calculation of flux-type variables at a new time step by extrapolation of local Riemann invariants (Phase 2) • calculation of conservative-type variables at a new time step using finite-volume schemes (Phase 3).
Short-wave perturbations are generated in the areas of loss of hyperbolicity. Computed flux-type variables have to be filtered to regularize the problem.

The exchange of mass and momentum between the layers
Regularization by filtering flux-type variables that does not affect conservative-type values does not always make a full system correct for long enough periods of time. It is necessary to adjust the conservative-type variables, which leads to the exchange of mass and momentum between the layers. In the proposed algorithm, the restructuring of grid is carried out twice at one time step. The first time after the calculation of conservative-type variables on the intermediate layer in phase 1, the second time after phase 3.

Example of calculations. Kurganov's Test
The initial data:    In the work of A. Kurganov [6], the calculation of this problem was made with the same parameters, and similar results were obtained. Further calculation according to the CABARET scheme showed that the solution loses its stability, becomes irregular and the calculation stops (Fig.3). In the calculation on a denser grid with 1601 nodes loss of stability occurs much earlier. On a coarser grid (201 nodes), the instability develops more slowly (Fig. 4), but poor conditionality of the problem inevitably leads to abnormal termination of the calculation (Fig.  5). Such behavior of the solution will most likely be in calculations by A. Kurganov's scheme [6]. The restructuring of the computational grid and, as a consequence, the exchange of mass and momentum between the layers regularizes the problem.

Conclusion
Multilayer hydrostatic model with a free surface, excluding the mass and momentum exchange between the layers, in most cases of interest are incorrect and cannot be used in practical calculations. The restructuring of the computational grid at each time step allows to regularize the problem. Multilayer models should be used to adequately describe the diffusion of mass and momentum associated with the loss of hyperbolicity properties.