Kinematic and dynamic analysis of dam break flow impact on vertical walls using weakly compressible SPH

This article presents the kinematic and dynamic analysis of a dam break flow based on data obtained from numerical solutions by the smoothed particle hydrodynamics (SPH) method. The method and original algorithms necessary for correct pressure evaluation are thoroughly described. The pressure evaluation method consists of data reading using virtual sensors and filtration in the time domain using the weight function. A simple convergence study showing the independence of the evaluated parameters of spatial resolution is presented together with the validation of the introduced methods and algorithms using a simple hydrostatic problem and experimental data available in literature. We focus on two parameters that describe the problem: the distance of the downstream vertical wall from the edge of the liquid column and the column’s height to width ratio. We found that the impact can be divided into three consecutive phases characterized by specific kinematic (flow patterns) and dynamic (exerted pressure and forces) behavior and different roles of the investigated parameters during these phases. During the early stages of an impact, the column’s distance from the vertical wall plays a major role. The dependence between the column distance and the force peak in this stage was identified in the form of a power function. In the second stage, when a rolling wave emerges, the vertical wall position influences the shape of the wave and the pressure distribution on the wall. The total force is greater in this phase for lower column height to width ratios due to the higher total momentum of the liquid. In the third stage, when the rolling wave impacts the liquid surface, the employed methodology with the two-dimensional solution and free-surface approach seems to reach its limits of applicability. More complex modeling would be necessary to capture this phase of the impact properly. c © 2021 University of West Bohemia.


Introduction
A simple dam break is a problem of permanent interest since not all of its essential features have been already described suitably. It has many modifications, but the general idea remains the same: liquid is initially kept in a tank, surrounded by walls that keep it at rest; then suddenly, some of these walls or their parts disappear. This naturally leads to a discharge of the liquid from its original place in the form of a surge. Obstacles of various shapes can stand in its way, breaking it or changing the flow direction. Flow structures formed by these interactions may appear quite impressive and the assessment of hydrodynamic forces acting on the impacted experimental objects is crucial because these effects must be considered when designing engineering structures that have to withstand similar dynamic loads regularly or constructions that might encounter similar incidents during their life, for example, coastal and off-shore structures or buildings in potentially flooded areas. These are the main reasons why a dam break problem remains so popular among researchers for decades, and they keep investigating it from different perspectives.
Many works on a dam break problem deal mostly with the kinematics of the flow. Height of the liquid surface or position of the surge front is investigated in most of these contributions, and they differ in the employed methods and used problem geometry. Martin et al. [30] published their investigation of the collapse of rectangular and cylindrical water columns of different sizes, including measured data of the surge propagation and height of the collapsing column. A set of non-dimensional variables was proposed for the problem to compare the results of different scales. The kinematic data provided by Martin et al. in [30] were used numerous times for the validation of methods of free-surface flow problems, for example, the volume of fluid for the finite volume method by Hirt and Nichols [18] and the moving particle semiimplicit by Koshizuka and Oka [24]. A dam break experiment was performed in the latter work, as well. Stansby et al. [42] identified mushroom-like jets shortly after the water release in their experiment with wetted bed and these jets appeared in the numerical simulation, as well. Nsom [36] investigated the tip of the surge for high viscosity liquid and derived a theoretical relation for the surge front position. Jánosi et al. [23] studied the friction and drag caused by the bottom of the tank. Cruchaga et al. [12] carried out a work dealing with the kinematics of the flow after a dam break of water and shampoo, where the experimental results served for the validation of the finite element method. Hu and Sueyoshi [20] combined the experiment with the numerical solutions using two different methods (the constrained interpolation profile and the moving particle semi-implicit). Ozmen-Cagatay and Kocaman [37] measured the surface evolution for dry and wetted beds and compared the shallow water model with a commercially available RANS solver. Feizi Khankandi et al. [14] experimentally investigated different shapes of the canal and their influence on the surface height. Bellos and Hrissanthou [5] solved a dam break for the case with a triangular hump on a flat bottom, both experimentally and computationally. A very similar arrangement was used by Ozmen-Cagatay et al. [38] in their experimental and numerical investigation. Recently, Hu et al. [21] tested various numerical approaches and used data previously published by Ozmen-Cagatay and Kocaman [37] and Feizi Khankandi et al. [14]. Another recent work was published by Yang et al. [45], who numerically investigated the turbulence and trajectories of concrete fluid particles. Finally, Wang et al. [43] experimented with rectangular and triangular channels.
Several authors focused primarily on the dynamic effects of the flow in their investigations. Zhou et al. [47] investigated the dynamics of the flow after a dam break in their experiment by measuring pressure values in time at different heights of a vertical wall impacted by a surge. They used these values for validation of their computational method for the prediction of water flow on a deck of a ship and its dynamic effects. Colicchio et al. [11] used these data as a reference for their comparative study, in which they solved a dam break by three different numerical methods (the boundary element method, the smoothed particle hydrodynamics, and the level-set method). Colagrossi and Landrini [10] investigated the influence of the surrounding gaseous phase on the solution using SPH. Bukreev [6] investigated the force acting on a vertical wall after the impact of a wave generated by a dam break and compared the results with theory. Another work combining an experimental investigation with a numerical solution was published by Wemmenhove et al. [44]. The experimental setup was very similar to that used by Zhou et al. [47], but the numerical investigation employing the finite element method included a flexible impacted vertical wall. Marrone et al. [29] simulated several different dam break cases using SPH in both 2D and 3D and compared them to each other. Abadie and Mokrani [1] pointed out the influence of mesh for pressure evaluation of liquid impacts. Lobovský et al. [28] published a very thorough work primarily dealing with the experimental investigation of the dynamics of the flow by measuring the pressure at different heights of the vertical wall. Laigle and Labbe [25] used SPH for the investigation of mudflows and tested their model on a dam break. Mirauda et al. [32] solved a dam break problem using SPH with two different ways of modelling the wall boundary conditions and compared their results with the experiments of Lobovský et al. [28] and with the numerical results found using three modifications of SPH method investigated by Cercos-Pita et al. [9]. The resulting flow after a solitary wave impacts on a vertical wall and accompanied dynamic effects were investigated by Paprota et al. [39] using SPH and a potential-flow-based semi-analytical method. A similar work, but with an impacting liquid body with a free surface, was published by Pugliese Carratelli et al. [40], who used SPH and FVM.
Another category of papers investigates the effect on some objects obstructing the dam break flow. Objects of different shapes are used and the total force exerted by the flow is evaluated. Bukreev and Zykov [7] placed a flat plate into the flow caused by a dam break and measured the total force. Cummins et al. [13] simulated a dam break in a rectangular tank with a rigid column placed inside of it using the SPH method. The total force exerted on the column was evaluated and compared with experimental data. Fondelli et al. [16] solved a dam break in a tank containing a rigid prismatic block using FVM with VOF. The total force exerted on the block was evaluated and compared with experimental results. A similar case was investigated by Aureli et al. [4], who used three different numerical approaches and compared the results with the experiment.
The articles published so far mainly deal with the kinematics of the flow field for various geometries and employed different experimental or computational methods. Fewer contributions investigated the dynamic effects induced by the flow. Pressure values obtained from the numerical solutions were compared with experimental data only for validation at particular points given by the placement of the sensors in the experiment. A similar description holds for the evaluation of total forces exerted on objects in the flow. This work exploits the strengths of the numerical solution, which allows to easily capture more information than most experiments, and to analyze the collected data from both kinematic and dynamic viewpoints because these two perspectives are naturally connected. We analyzed a relatively wide range of similar yet different cases in this work. This approach revealed some patterns in the flow depending on the chosen parameters: distance of the downstream vertical wall from the edge of the liquid column and the column's width to height ratio. Another objective of this paper is to describe the used weakly compressible SPH method in detail, especially the algorithms that are necessary for correct pressure evaluation because the preceding authors usually omit this essential part.
The article expands on the paper [22] published by the authors previously. There is only one variable parameter (vertical wall position) and only three computed and analyzed variants. The current paper works with more variants, and therefore it can bring conclusions that are generally applicable.

Computational method
The weakly compressible smoothed particle hydrodynamics (SPH) method was employed for the numerical solution of the dam break problem in this work. It is a mesh-free Lagrangian particle method that is particularly suitable for transient multiphase problems. These properties make it a good choice for the dam break problem. The free surface approach is used, which means that the gaseous phase surrounding the liquid is omitted. This simplification is justified by the fact that its influence is negligible and there is also a positive effect on the side of computational performance [10].
This method was first used for free-surface flow by Monaghan [34] and since then it has been used many times for similar problems, both in 2D and 3D [4,10,11,13,29]. For a more detailed description of the method, see for example the monograph by Liu and Liu [26].
In the SPH method, the spatial derivatives of the field variables are approximated using the so-called weight function W . Liu et al. [27] name eight properties which a function has to possess to be a good weight function (it is normalized, compactly supported, non-negative, monotonically decreasing, degenerating to the Dirac delta function, symmetric, smooth, and applicable to more dimensions). However, not any function having these properties is suitable for being employed as a weight function. From the practical point of view, the crucial parameters are accuracy (the ability of a function to approximate the field variables and their derivatives correctly), and stability (the ability of a function to provide a numerically stable solution).
Many different weight functions have been proposed and tested. Hongbin and Xin [19] formulated criteria for weight function evaluation and tested ten different functions. They found the Gaussian kernel and the Q-spline function as the best choice. Yang et al. [46] proposed a new cosine weight function and compared it with the Gaussian kernel and the cubic spline function. No substantial differences were found in the performance of the three functions.
In this work, the truncated Gaussian kernel function was employed. It can be written as where d is the number of spatial dimensions, h is the so-called smoothing length, and R is the particle distance normalized by the smoothing length. Smoothing length is taken the same as the initial particle spacing. Truncation of the kernel does not affect the solution, but it improves the computational performance significantly by reducing the number of interacting particles. The continuity equation and the momentum equation are partial differential equations, which are transformed into ordinary differential equations with respect to time using a standard SPH procedure. The fluid is considered compressible and with no physical viscosity. In particular where t is the time, is the density, m is the mass, v is the velocity vector, p is the pressure, and f is the external force intensity per unit mass. Indices i and j denote particles, and Π is the artificial viscosity that serves for numerical stabilization.
Since the fluid is modelled as compressible, an equation of state is necessary to close the system of equations (2) and (3). The adopted equation of state was first used in SPH in [34]. It is usually written in the form where 0 is the reference density, c is the speed of sound, and γ the exponent of the equation. It applies for all fluid particles, and the particle identification indices are not written. The exponent is usually chosen as 7 for liquids. The value of the speed of sound is not the physical one, its choice is based on the specific problem. It should be chosen about ten times higher than the highest velocity occurring in the problem. That makes the maximal Mach number of the problem so low that the fluid can be regarded as incompressible.
Monaghan's artificial viscosity is used in the computation and it can be written as where α is the artificial viscosity coefficient, and x is the position vector. Symbols with both i and j indices are taken as average values from the mentioned particles. The pressure field is very noisy in weakly compressible SPH and it is unsuitable for pressure evaluation, as stated by Gomez-Gesteira et al. [17]. Several ways how to deal with this noise were proposed. The first one is to add artificial mass diffusion into the governing equations, which was done by Antuono et al. [3], Ferrari et al. [15] or Molteni and Colagrossi [33]. Another way is to re-initialize the density field after a certain number of integration steps (20-50). There are two variants: the zero order filter proposed by Shepard [41] and the first order moving least squares [10]. Shepard filtering is used in this work since it is easier for implementation and less computationally expensive. Moving least squares conserve the total energy slightly better [10], but the main purpose of the filtering is to eliminate the noise, which the Shepard filtering does well enough for the pressure evaluation. Re-initialized density is found for every particle using the formula A free surface boundary is formed naturally by the particles. The pressure value is zero on the free surface. Since the liquid is considered inviscid, the wall is a free-slip boundary condition. The model proposed by Adami et al. [2] with a slight modification is employed in this work. Walls consist of three layers of so-called dummy particles. These particles remain fixed in space and their velocity is zero. Pressure values of these particles are extrapolated from the fluid using formula where the indices w and f denote the wall particle and fluid particles, respectively. The density of the dummy particles is evaluated from their pressure values using the equation of state (4). The dummy particles take part in the governing equations (2) and (3), just like liquid particles. The slight modification lies in the application of the artificial viscosity. The coefficient α in equation (5) is taken ten times smaller for fluid-wall interaction than for fluid-fluid interaction, which brings the model closer to the true free-slip condition. The original paper by Adami et al. [2] does not mention artificial viscosity regarding the proposed wall boundary condition method. Various explicit schemes can be employed for the integration of the governing equations (2) and (3). A modification of the leapfrog algorithm was employed in [26]. The standard leap frog cannot be used because the acceleration of the particles depends not only on the particle positions, but on the velocity and density, as well. One time step of the scheme can be written as v where the superscript denotes the time step number and Δt is the time step length.
Since an explicit scheme is used, the solution is conditionally stable. The Courant-Friedrichs-Lewy condition has to be satisfied. Maximal value of a time step for SPH can be written There are other stability criteria, but the CFL condition is the most restrictive from the conditions presented by Morris et al. [35] when used in the solved cases. Some details of pressure evaluation are rarely mentioned in the papers describing SPH. Laigle and Labbe [25] used rectangular virtual pressure sensors. The average pressure value of the particles present in the sensor area is recorded every N time steps. Cummins et al. [13] mention the necessity of signal filtering in the time domain using the triangular-weighted doublesided moving average. Recently, Meringolo et al. [31] proposed a promising method of pressure signal filtering based on wavelet transform.
In this work, the pressure is evaluated by sensors fixed in space. The value at a certain time step is a weighted average of the surrounding fluid particles. The dummy particles are not considered in the pressure evaluation because they carry no new information since their parameters are extrapolated from the fluid particles. The truncated Gaussian kernel function defined by (1) is taken as the weight function. The used formula stands where the index S denotes the sensor location. The reading of the sensor depends on its influence domain. The smoothing length of the weight function is chosen the same as for the simulation itself. This choice keeps the number of influencing fluid particles in 2D simulations between 10 and 20, which appears to be a reasonable amount. Larger influence domain leads to incorrect reading since particles from too far out are considered while a smaller one gives too noisy readings.
Raw signal from the sensors is still very noisy due to high-frequency oscillations in the pressure field. The weight function in the time domain is used here for the suppression of high frequencies. This filtering requires the definition of the smoothing period τ that replaces the smoothing length h. The choice of this smoothing period depends on the solved problem since it should reduce the numerical noise but not the physical phenomena occurred. The time filtering formula stands where the subscripts m and n denote time samples. There is no phase delay due to the symmetry of the Gaussian kernel function, but filtering must be carried out in post-processing. The chosen sampling period of five integration steps proved to be sufficiently fine to avoid signal aliasing.
For the total force affecting the wall, the pressure value of each sensor is multiplied by the corresponding area and these elementary forces are simply summed up. In the mathematic notation, this can be written as The SPH program, based on the Fortran code provided in [26], was used for the computations presented in this work. Post-processing is carried out in MATLAB.

Problem definition
This article presents solutions for simple two-dimensional dam-break problems with solid vertical walls downstream. We investigate three different aspect ratios of the liquid column and three different positions of the vertical wall. The dimensions are chosen to match the experiment of Lobovský et al. [28] for easy validation of the method and for mutual comparison of the results. A scheme with relevant dimensions and validation pressure probes is shown in Fig. 1. The results are presented in non-dimensional variables for easier comparison with outputs from other works. Non-dimensional time, horizontal coordinate, vertical coordinate, pressure, force, and velocity together with the liquid column aspect ratio are defined by where x and y are the coordinates defined in Fig. 1, y 0 is the initial liquid column height, x 0 is the initial liquid column width, z 0 is the transversal dimension of the problem (unity for two dimensions), g is the gravity acceleration vector, and F is the total force exerted on the vertical wall.

Convergence study
To find a suitable spatial resolution for subsequent simulations, the results obtained with a different number of particles were compared. Kinematics of the flow proved to be resolution independent even for relatively coarse resolutions, while pressure evaluation appeared to be much more sensitive to this parameter. Therefore, the dam break with the initial column height y 0 = 0.3 m was simulated using three different spatial resolutions (initial configurations 100×50, 150 × 75, and 200 × 100 particles), and the pressure readings at sensors P1-P4 placed as shown in Fig. 1 were evaluated. The results are shown in Fig. 2. All physical and numerical parameters were kept constant for all three resolutions, except for the time step and domains of dependency of the sensors, which were scaled in accordance with the resolution. It is seen that increasing resolution suppresses oscillations in the signals. Overall, the signal magnitudes are very similar. The biggest differences were identified for sensors P1 and P2 immediately after the impact. It is because the tip of the surge is formed by a relatively small number of particles, so its shape is more dependent on the resolution than the rest of the flow. The effect of artificial viscosity is found to be very small because the surge front reaches the wall practically at the same time even though the effects of artificial viscosity is scaled with spatial resolution. The finest of the compared resolutions was used for the subsequent simulations since it was found sufficient for the desired purpose and further refinement would make the solutions too computationally expensive.

Computational model validation
The proposed pressure evaluation technique was validated on a pressure distribution along the submerged vertical wall. This distribution is depicted in Fig. 3 and compared with the theoretical hydrostatic distribution. Fifty particles in the vertical direction and the same number of pressure sensors were used. A great majority of the sensors gives virtually the same results as the theoretical hydrostatic solution. A small difference between the simulation and the exact solution is near the surface because fewer particles are influencing the sensor. The same happens near the bottom, but the effect is greater there. The sensors are surrounded by particles that are on average in lower depth than the sensor, so the pressure reading is lower. This affects just three to four sensors nearest the bottom and the relative error diminishes with increasing spatial resolution. In the presented case, the relative error is about 5 %. The accuracy of the kinematics of the flow given by the utilized computational method is validated by the comparison of the surge front position with the experimental data previously published by various previously mentioned authors [20,24,28,30]. The result of this comparison is shown in Fig. 4.
The Ritter solution (see, e.g., in Castro-Orgaz and Chanson [8]) implies that there is a linear relation between the surge front position and time (in non-dimensional variables), and there is no influence by other parameters. However, the experimental datasets differ from each other. One reason can be the different experimental setups or measurement techniques. It can be expected that the surge propagation is faster for bigger y 0 because the effect of bed friction is relatively lower compared to inertia. The experiments performed by Martin et al. [30] suggest that it might be true for later phases (T > 2).
If we compare the data published by Koshizuka and Oka [24] and by Lobovský et al. [28], there is only a small difference in y 0 (0.292 m and 0.3 m), but the surge front propagates considerably further in the second case for the same time instances. The initial column height to width ratio A (2 and 0.5) is the main difference between the two compared cases. A possible explanation for this is that the height of the column, which is the source of momentum for the liquid motion, decreases faster for cases with higher ratio A because the overall amount of liquid is smaller.
We consider the data from [28] as the most relevant since they are the most recent, and the methodology of obtaining the data is thoroughly described in this article. Therefore, we make the simulation with the same parameters (y 0 and A) and compare the results primarily with these data. The comparison is also presented in Table 1. Table 1. Surge front position -comparison between simulation and experiment [28] T (1) X simulation (1) X Lobovský et al. [28] (1) relative value (1)  The simulated surge front lags behind the experimental observation in the initial phase of the liquid column collapse. The relative difference is relatively high due to small absolute values (the simulation reaches only 73% distance compared to the experiment for T = 0.5). Later, the current SPH simulation predicts the surge front to be slightly ahead of the experiments (the simulation reaches 107% and 108% distance compared to the experiment for T = 1.5 and T = 2.0, respectively). This is probably because the boundary layer is not properly modeled. However, this fact appears not to have much influence on the kinematics and dynamics of the flow overall since it is mostly determined by the inertia of the fluid and the viscosity plays just a small role.
The pressure evaluation technique is further validated using experimental data provided in the work by Lobovský et al. [28]. The pressure was measured on the vertical wall at four different heights above the bottom of the tank as marked in Fig. 1. Two cases with different initial liquid column heights were evaluated. The results for y 0 = 0.3 m (X = 3.37, A = 0.5) and for y 0 = 0.6 m (X = 1.68, A = 1.0) are shown in shows Fig. 5 and Fig. 6, repectively. The agreement between the numerical and experimental data is obvious. The initial impacts are captured well, both the time of their occurrence and the pressure values. A relatively wide range of results was captured during the experimental campaign. For example, the maximum impact pressure for sensor P1 (non-dimensional according to (17)) was found between 2.55 and 4.45 (percentiles 2.5 and 97.5, respectively) with median 3.11. The experimental results include a certain degree of uncertainty, so they should not be treated as exact.
In the first case, the simulation overpredicted the pressure value for the lowest sensor P1, while the values for sensors P3 and P4 are slightly lower than in the experiment. The value of P2 is in a good agreement with the experiment. Later after the impact, the first three sensors indicate slightly higher pressure, while P4 remains under the experimental data. From T = 4.5 to T = 6.0 the agreement of numerical and experimental data is very good. After that, a slight increase in the experimental values occurs and a similar rise in the numerical values arrives with a noticeable delay.
In the second case, the impact pressure value of the lowest sensor P1 is higher in the simulation than in the experiment. The other three sensors show a good agreement until about T = 4.5, when the simulation prediction values start growing, while the results from the simulation remain almost constant.
Overall, the comparison between the pressure data from the simulations and from the experiment shows that the simulation can be trusted, especially in the early stages of the impact. The results support the idea that the solved flow depends mostly on fluid inertia and viscosity and other factors play just minor roles. Later after the impact, some discrepancies appear, mostly due to the two-dimensional model simplifications.

Kinematic and dynamic analysis
We thoroughly analyzed the kinematics and dynamics of all nine dam break cases and identified three consecutive phases of the flow. For the sake of clarity, the presented figures do not include all analyzed cases but only a chosen sample, which illustrates the similarities and differences between the cases in the given time phase.
The impact force as a function of time in non-dimensional variables was evaluated, and the results are plotted in Fig. 7. This graph shows all nine cases and the whole analyzed period. Therefore, it is referred to multiple times further in the text.

Phase one -initial impact
In all investigated cases, the liquid forms a wedge running towards the downstream vertical wall. The angle between its surface and the horizontal bed evolves from the right angle in the very beginning of the process to a very shallow angle when the surge front reaches further from the initial position. Therefore, the angle at the moment of impact is mostly dependent on the position of the vertical wall, as shown in Fig. 8. There are no distinct gradients in the velocity field near the tip of the surge except for a noticeable deceleration of particles adjacent to the wall, particularly in thin liquid layers. In the cases with the vertical wall further away, the liquid layer consists of a single layer of particles. It affects the precision of the maximal pressure value due to the limited resolution. However, this situation lasts only for a very limited period in a very limited domain, and thus it does not affect the evaluation of the total force.
The pressure instantly rises immediately after the impact. The magnitude of this pressure peak decreases linearly with the vertical wall being further from the starting point, as illustrated in Fig. 9. Since the surge front does not change its velocity significantly, it is mostly the shape Fig. 9. Non-dimensional maximal pressure peak value as a function of non-dimensional wall distance of the liquid tip that affects the maximal pressure value. If the wall is closer, the angle of the impacting liquid wedge is less acute, so the pressure peak is more pronounced.
Local Froude number is a measure, whether there is a jet (high Froude number) or just a refracting wave (low Froude number), as pointed out by Pugliese Caratelli et al. [40]. The local Froude number can be expressed using the previously defined non-dimensional quantities as At the moment of an impact, Y at the tip of the surge is virtually zero, so the Froude number is very high. Therefore, a jet along the vertical wall can be expected. Even later, after the jet emerging, the Froude number remains relatively high (between 3 and 6 for all cases). The jet's behavior is mostly determined by the angle of the impacting surge and thus by the position of the vertical wall. For the cases with a closer wall (X ≤ 3.37), the jet is relatively small in comparison with the amount of the impacting liquid, and it poses significant velocity, see Fig. 10 (left). These cases result in a steep rise in the force that culminates in significant local force maxima, followed by a slight decrease. The closer the wall is, the more distinct the maximum is. If the wall is further away (X ≥ 6.73), the liquid layer on the bottom and the jet on the wall are about the same thickness and the velocity magnitude does not change significantly, see Fig. 10 (right). This scenario also leads to a steep force rise after the impact, but it is followed by an almost constant value of the force for a short period. This force plateau is a result of quasi-steady-state moments when the thickness and velocity of both horizontal and vertical streams are changing very slightly.
We found that the magnitudes of the force maxima are dependent only on the position of the wall X and not on the column aspect ratio A. This finding allowed plotting a graph of these force maxima as a function of the wall position X (Fig. 11). We found an appropriate approximation function for the data Φ = 0.805 1.05X .
It seems to be valid in the investigated range of X and it should be noted that this equation gives Φ = 1 for X = 0, which follows the hydrostatic solution.

Phase two -bulge and rolling wave
Later, the shape of the jet still depends mostly on the vertical wall position. In general, it tends to be thicker for walls closer downstream and thinner for walls further downstream. The same applies to the horizontal impacting stream as well. At this time, the vertical jet gradually loses its momentum, decelerates, and the fluid eventually starts falling back down, while the horizontal stream that reached the vertical wall at this time is moving up. This results in the formation of liquid bulges on the vertical wall (Fig. 12). The bulge evolves into a rolling wave. The amount of liquid forming the rolling wave depends on the jet thickness, and thus also on the vertical wall position. A thin jet leads to a rolling wave that contains a significant part of the original jet, while thicker jets evolve into rolling waves formed by a relatively smaller amount of liquid and the waves emerge further away from the wall (Fig. 13). The phenomena described in the previous paragraphs cannot be detected in the total forcetime dependency in Fig. 7; the force gradually rises during the relevant period. However, we observed them in the pressure distribution on the vertical wall. At the time when the bulge starts emerging, the pressure rises in the upper part of the wall. It is not too significant for the cases with thicker jets, where rolling waves form further away from the wall, see Fig. 14 (top). However, for the cases with the impacted wall at greater distances, we identified clear spots of pressure rise. These are the spots from which the rolling waves grow, see Fig. 14

(bottom).
It is seen in the force-time graph (Fig. 7) that the overall force is significantly smaller for the cases with the initial column aspect ratio A = 1.0 than for the other two aspect ratios. Shortly after the initial impact, the force on the vertical wall is dependent only on the wall position, and not on the initial column aspect ratio, as pointed out in previous paragraphs. The same applies for the flow field near the wall (Fig. 15).
Later, when the bulge emerges, the shapes of the fluid bodies are still very similar. However, a considerable difference in the horizontal fluid layer thickness can be noted (Fig. 16). A relatively smaller amount of fluid is available for the higher aspect ratio case, and so the liquid level drops more than in the other two scenarios. Less impacting fluid means less momentum and thus lower force.
The situation described in the previous paragraph further develops with the formation of rolling waves. The shapes of the waves themselves are quite similar, but the one with the lower aspect ratio is considerably higher due to the thicker horizontal liquid layer (Fig. 17). Consequently, the overall wall area in contact with the liquid is larger, which leads to a higher total force. In these later phases, the case with A = 0.5 results in a smaller force in comparison with the cases with A = 0.25.

Phase three -rolling wave impact
The final part of the analysis is dedicated to the impact of the rolling wave on the liquid surface. It is important to remind that the calculations are conducted as free-surface problems, i.e., there is no surrounding gaseous phase. This simplification has a positive effect on the computational performance and there is virtually no difference in results with or without surrounding gas in previously analyzed phases (for gas-liquid density ratio 0.001, which is close to the common air-water density ratio) [10].
As the rolling wave impacts the surface, a gas pocket is trapped under it and affects the flow. This effect cannot be captured with a free-surface approach, and there are differences in both kinematics and dynamics of the flow. For example, different pressure peaks caused by different phenomena were identified when free-surface and surrounding gas approaches were compared [10]. However, even for the surrounding gas approach, there were some discrepancies when it came to comparison with the experimental data. It seems that the flow in this phase is more complex and the two-dimensional modeling cannot capture some phenomena such as droplet formation or surface deformation. Therefore, the results from this phase of the analysis should be approached carefully.
It shows that the total force on the wall peaks sharply at the moment of impact and reaches the global maximum. The only exception is for the case with the vertical wall being the closest (X = 0.84, A = 1.0), where the rolling wave is relatively the smallest, and it impacts far from the wall. In general, with the wall being further away, the relative magnitude of the force peak grows, and it gets more clearly separated from the gradual growth in the previous phase. The force peaks can be identified in the pressure distributions in Fig. 14, as well.
For cases with smaller X, the impacted horizontal liquid layer is relatively thick, it flows slowly, and there is a large amount of liquid between the rolling wave impact point and the vertical wall, see Fig. 18 (top). As a result, the impact influences the total force just slightly. The sudden drop in force after its maximum is caused by a sudden detachment of the liquid from the wall.
For higher values of X, the horizontal layer is thinner, still possesses significant velocity, and there is a relatively small amount of liquid between the wave and the wall. This results in a significant acceleration of the liquid layer under the wave towards the wall, see Fig. 18 (bottom). The high-velocity liquid then causes a substantial pressure rise on the wall. The detachment of the liquid results in the force drop after the peak.
Due to the facts mentioned at the beginning of this subsection, we refrain from a more detailed qualitative analysis of this phase. A more suitable model would be necessary to capture three-dimensional phenomena and air entrapment, which would be necessary for a proper analysis.

Conclusion
This work presented a weakly compressible SPH method for free surface flow together with an original pressure and force evaluation technique suitable for this method. We used this method to solve a dam collapse problem in two dimensions. We compared the obtained results with the available experimental data from the literature and confirmed that the method gives relevant data for both kinematics and dynamics of the flow.
We found out that there are three distinct phases in the impact. In the first one, the dynamic and kinematic variables depend virtually only on the wall distance X. This fact allowed us to formulate a relationship between the force maxima in this first phase and the impacted wall distance. In the second phase of an impact, a higher aspect ratio leads to a lower force, mostly due to a relatively lower amount of liquid and thus overall lower momentum. The third phase begins after the rolling wave impacts the liquid layer lying on the bottom, and the force reaches its global maximum. The total force exerted on the vertical wall depends on both investigated parameters. However, the force peak appears relatively more significant for cases with the wall being further away. In the third phase, the model reached its limits of applicability, so we do not consider the results very reliable.
We believe that the presented kinematic and dynamic analysis performed on numerical data will serve as a reference for future numerical or experimental investigations, or it can aid coastal structure design, especially in its earliest stages because it gives a force load that can be expected under various conditions. Later in the design process, a geometry-specific model should be used.
A relatively wide range of parameters was analyzed in this work, but future works can go beyond the plotted boundaries. Either by a different choice of parameters or by employing a more complex model. A model that is three-dimensional and takes the surrounding gaseous phase into account would increase the credibility of the results, especially in the third phase.