Ingeniería 31(1): 98-111, Enero-junio, 2021. ISSN: 2215-2652. San José, Costa Rica
DOI 10.15517/ri.v31i2.45850
Esta obra está bajo una Licencia de Creative Commons. Reconocimiento - No Comercial - Compartir Igual 4.0 Internacional
Flow under long two-dimensional dam sluice gate using WSPH
Flujo en compuerta inferior de represa bidimensional larga por el
método de WSPH
Juan Gabriel Monge-Gapper
School of Mechanical Engineering
University of Costa Rica, San Jose, Costa Rica
Email: juan.mongegapper@ucr.ac.cr
ORCID: 0000-0002-8093-8239
Alberto Serrano-Pacheco
School of Civil Engineering
University of Costa Rica, San Jose, Costa Rica
Email: alberto.serrano@ucr.ac.cr
ORCID: 0000-0002-4136-0841
Received: february 15, 2021 Accepted: June 1, 2021
Abstract
An application of the Weakly compressible Smoothed Particle Hydrodynamics (WSPH) numerical
method is presented here for the case of two-dimensional ow in a long channel with a partially open sluice
gate. The results are compared with an analytical solution provided by shallow water equations (SWE)
and available experimental data for instantaneous discharge depth at the gate. Of particular interest is the
application of this numerical method to a sluice gate case with a high ratio of channel length to depth, which
tends to amplify the effects of the chosen numerical resolution. General model congruence was within 5 %
of a margin error even for relatively low vertical resolution, and some side effects of the equations used to
describe the boundary conditions were identied.
Keywords:
Channel, computational, dam, hydraulics, model, particles, sluice, vena.
Ingeniería 31(1): 98-111, Enero-junio, 2021. ISSN: 2215-2652. San José, Costa Rica DOI 10.15517/ri.v31i2.45850
99
Resumen
Se presenta la aplicación del método numérico de Hidrodinámica de Partículas Suavizadas con
Compresibilidad Articial (WSPH) al caso de ujo bidimensional en un canal abierto largo en el que hay
una compuerta inferior abierta. Se contrastan los resultados del modelado con la solución analítica según las
ecuaciones de ujo poco profundo (SWE) y datos experimentales disponibles para la profundidad instantánea
de descarga en la compuerta. De particular novedad es la aplicación del método a un caso de compuerta con
una alta relación de longitud del canal a su profundidad, lo que tiende a hacer más notorios los efectos de
la resolución en los resultados numéricos. Se observó congruencia general del orden de error medio del 5
% del modelo incluso para resoluciones verticales relativamente bajas y se identicaron los efectos de las
ecuaciones usadas para describir las condiciones de contorno.
Palabras clave:
Canal, compuerta, computacional, hidráulica, modelo, partículas, presa, vena.
MONGE Y SERRANO: Flow under long two-dimensional dam sluice gate using WSPH
100
1. INTRODUCTION
Many analytical one-dimensional hydraulic models that nd practical application involve
ow lengths that are much greater than its depth or width, such as that in channels, rivers and
piping systems. For cases where vertical uid movement within the ow or variations in depth
are small but not negligible, analytical solutions for the shallow water equations (SWE) have
been developed [1] to provide a good approximation if the boundaries are such that they have
an analytical solution. For more complex boundaries, SWE can be used to generate a set of
differential equations to be solved numerically at a much lower computational cost than applying
a general purpose two-dimensional numerical model of uid ow.
However, precision is just one of many properties that a numerical model may exhibit,
and the particular case of ow with a high aspect ratio can give greater insight as to the impact
of highly non-symmetric resolution when using a general-purpose numerical uid ow model
such as WSPH. This is especially attractive for cases with an open boundary since no dynamic
remeshing algorithms or surface detection procedures are needed, as is the case of the Finite
Volume Method (FVM), among many others. The free movement of particles does mean solid
or elastic boundaries need coherent meshless formulation through equations that have numerical
parameters that must be chosen carefully to avoid decient modeling of ow close to discharge
areas, like those at a sluice gate or pipe outlet.
The transient case of a sluice gate opening at the end of a long channel was chosen to set up
the models for this study, because it serves as an excellent showcase for potentially insightful
phenomena. The results were then compared directly with previously published data, including
an analytical solution [1], and the vena contracta locus and depth from physical experiments
[2], which were for a dammed channel with a 50:1 length to depth ratio and constant rectangular
cross section.
Previous studies using SPH to model sluice gate opening [3], [4], successfully have
concentrated on the formation of the hydraulic jump and the effects of the shape and roughness
of the discharge area. Other studies [5], [6], deal with the addition of an obstacle close to the
discharge area and how the shape of the hydraulic jump changes in time. Studies using the
analytical SWE equations to calibrate SPH models [1], [4], do not mention the effect of low
resolution on shallow depths or rely on relatively short lengths to model high aspect ratios.
This study focuses on the effect of using the entire length of the channel according to
boundaries described in previously published studies, and to show the behavior of the shallower
portion of ow during the rst few instants after the opening of the sluice gate over a dry bed, to
illustrate the impact of vertical resolution on the SPH model of this type of ow more specically.
2. A BRIEF OUTLINE OF THE WSPH METHOD
The smoothed particle hydrodynamics (SPH) approach is a meshless discrete method where the
volume of uid is divided into a nite number of particles, usually uniform in size and nature. As
with any other Lagrangian approach to discretization, since the particles are not xed in space, their
101
movement is modelled with the uid dynamics equations for conservation of mass and conservation
of momentum [7]. Any other applicable equations required to model the physical state of the uid
(changes in density, viscosity, temperature, for example) must, of course, be considered as well.
The main dierence of the SPH method in comparison with Newtonian particle methods
1
, such
as the Multiple Particle Method (MPM) and Particle In Cell (PIC) method [8], is that any single
particle interacts both with immediately adjacent neighbors and those within a spatial radius equal
to
h as shown in Fig. 1.
To avoid an unnatural distorted behavior of the particle, the farther away the neighboring par-
ticle is, the less it aects its neighbor. In the SPH method, this is done through a weighting function
called the kernel function W(x,h) [9], which can be any continuous, spatially symmetric function.
However, to guarantee conservation of physical meaning, the value of the scalar or vector eld that
the function is used with must not introduce articial alterations [7]. Also, it should be as compu-
tationally cheap as possible. This function smooths the distribution of the value of any variable
associated to particle dynamics, such as speed, acceleration, density, enthalpy, entropy, or viscosity,
among many others that may be necessary for the particular case. This smoothing helps to stabilize
the computational process and is extensive to the values of any integrals or derivatives that need to
be considered in the conservation of mass and the conservation of momentum equations applicable
to the ow [10]. With a larger spatial radius, more smoothing eect occurs, however, too large of a
value will cancel out any local dynamics and the ow will exhibit numerical clumping. Therefore,
a proper choice for smoothing length h (usually a value between 1.1 and 2.0 times the nominal
particle size) is a matter that requires careful consideration [10], [11].
Figure 1. A visual representation of the smoothing function W(x,h).
Since the central challenge for solving the uid mechanics equations is obtaining the value f(x
i
)
of the dynamic and physical variables pertaining to the state of the ow, the SPH method approxi-
mates them through a numerical scheme known as the particle approximation [7]. The general sta-
tement is that for a particle i, located at position x
i
an approximation of f(x
i
) can be calculated as a
weighted average of the value of f(x) for all N neighboring particles using (1), which can be readily
transformed to computational source code:
1
Often, this term has been used to describe methods where particles are assumed to be rigid, and where forces acting
on any particle occur only when there is direct contact with other particles or eld eects such as gravity. Due to the many
simplications implied by this type of model and practical computational issues such as determining whether a particle
has contact with another, results are seldom acceptable for scientic purposes and numerical stability is dicult to obtain.
MONGE Y SERRANO: Flow under long two-dimensional dam sluice gate using WSPH
102
(1)
According to this denition, the position of any particle must be known, making this scheme
an explicit forward time-step method. The origin and mathematical development of (1) and the
choice for the smoothing function and many other associated parameters can be studied in detail
in monographs on the topic of SPH applied to uid ow [10], [11].
A complete description of all applicable equations of uid ow is beyond of the scope of this
article, so only the main numerical driver for particle dynamics will be mentioned in this section.
This is the application of the particle approximation to the equation of conservation of momentum
[7], shown as (2) for two-dimensional ow in orthogonal coordinates:
(2)
The usual terms for uid ow are: the speed vector u, pressure P, viscosity , plane strain
rate tensor
xy
, and density . However, particle approximation produces additional terms which
include particle mass m, a second appearance of density
and partial derivatives of the smoothing
function. The reason to mention this equation here, aside from its importance, is that the multiplied
density terms
i j
not only affect the value of acceleration or any other variable through (1) but
also drive the variations. In the case of the weakly compressible SPH method (WSPH), this is due
to the form of the equation of state, that gives the uid an articial compressibility to correlate the
variation of pressure as a function of local density
i
. In this study, (3) was used and is known as
the Monaghan equation [7] for local effective pressure affecting particle i:
(3)
The value of
0
is a constant reference density of the uid being modelled; the true speed of
sound in the uid c
0
can be used, but for many subsonic ow phenomena, a much lower value can
be used without considerable loss of accuracy [10], so larger time steps can be used and therefore
lower computational cost is required. However, if this is done, to avoid numerical instability and
particle interpenetration, it is important that c
0
is at least an order of magnitude larger than the highest
expected uid speeds in the case being studied. Empirical constant is conventionally set at 7, but
can be adjusted by the user. Most importantly, the value of
i
evolves according to (1) which, like
all other evolving physical properties in SPH, is a weighted local average of the implied densities
of surrounding particles. Values of
i
far enough from
0
can produce articially high values of
P that can turn into small local numerical instabilities or even articial numerical shock waves all
across the modelled domain.
Therefore, these terms are key to the stable evolution of all variables calculated this way, if an
incompressible uid is to be modeled, one way to prevent numerical instabilities is to add arti-
cial compressibility to the numerical scheme. The addition of such an equation of state for nearly
Ingeniería 31(1): 98-111, Enero-junio, 2021. ISSN: 2215-2652. San José, Costa Rica DOI 10.15517/ri.v31i2.45850
103
incompressible uids is generally known as the weakly compressible SPH method (WSPH) used
for this study [7], [10], [11].
This is the general framework of SPH as a numerical method applicable to solve the equa-
tions of uid dynamics, shown as a starting point. However, it is a complex, well-evolved method
as shown by the considerable number of publications that use SPH to model many different uid
ow phenomena, like those mentioned in comprehensive reviews of the state of the art such as [12]
and up-to-date revisions of future challenges for the application of SPH methods in general [13].
3. METHODOLOGY
As mentioned above, available validation data drove to the selection of the particular dimensions
and boundary conditions for the numerical case setup. References [1], [2], published data for the
case of a long (50:1), horizontal open ow channel with a completely absent gate, a conventional
dam break. These same references also study this setup but with an open gate with a height equal
to 20 % of the initial dam height. The available analytical solutions are for any instant, but vena
contracta size is given for a time point t = 5.0 s after the gate opening for the following three cases:
a dry bed, shallow wet bed and deep wet bed. In this study, the dry bed conguration is of interest
and was used for the case setup. Fig. 2 shows the rectangular arrangement of the two-dimensional
uid that replicates the proportion and size of the physical experiment [2].
Figure 2. Initial particle arrangement for 2D WSPH model of an open dam gate.
In all cases, WSPH two-dimensional x uniformly sized particles were used and arranged in
an evenly spaced rectangular array. The Morris equation for the conservation of momentum [10]
shown in (1) was used to evolve momentum, and the Monaghan equation of state [7] with an arti-
cial speed of sound of c
0
= 80 m/s was used since uid speeds were expected to be well under 10
m/s. Properties of clean water at T = 20 ºC were used (viscosity
= 0.000896 Pa·s and density =
998 kg/m
3
). Cubic form [7], for smoothing function W(x,h) was used in all cases with a smoothing
length of h = 1.1·
x among the simplest typical value, important in terms of computational cost
MONGE Y SERRANO: Flow under long two-dimensional dam sluice gate using WSPH
104
due to the amount of times it needs to be evaluated along with its spatial derivatives, at least once
for each particle and must be recalculated at each time step.
All solid boundaries, including the static sluice gate when applicable, were modelled using a
single line of virtual particles spaced evenly at x/2; their behavior was modeled with the Monaghan
repulsive force equation [14], which affects nearing uid particles only when within a certain range
r0 of the boundary and becomes higher in inverse proportion to the distance between the boundary
particle and the uid particle (this proportionality can be adjusted empirically by the user). No tur-
bulence models or correction factors were applied. Time step was constant and chosen using the
Verlet criteria [11], which is based on expected numerical propagation of data. Iterations were set
to continue up to a physical time of t = 5.0 s. All other signicant parameters varied according to
particle size and boundary dimensions, and are specied in TABLE I for each of the three variants
used for this study.
TABLE I
SUBCASE SETUP FOR CONVENTIONAL DAM BREAK AND 20% SLUICE GATE
Open LR Open HR Gate HR
Fluid 1 m x 35 m 1 m x 35 m 1 m x 50 m
Boundary 2 m x 50 m 2 m x 100 m 2 m x 100 m
Gate None (dam break) None (dam break) Open (0,20 m)
Fluid particles 3 500 14 000 20 000
Virtual particles 4 080 8 160 8 160
Particle size 0.10 m 0.05 m 0.05 m
Time step 1.0·10
-3
s 0.5·10
-3
s 0.5·10
-3
s
Smoothing length 0.110 m 0.055 m 0.055 m
Boundary particle repulsion radius 0.090 m 0.040 m 0.040 m
The rst two setups were intentionally prepared with shorter boundaries to reduce computatio-
nal time and to determine how much of a difference appeared between a 70 % length model (Open
HR) and a 100 % length model (Gate HR) both in terms of computational time and model preci-
sion. The low-resolution model (Open LR) was used to verify evolution to determine whether an
even higher resolution model (>20 000 particles) was necessary.
These models were processed using the authors own source code that implements the WSPH
method, written in C programming language. The rst version of the functional source code, named
Gap-SPH is freely available at a public repository (https://sourceforge.net/projects/gap-sph/) under
GNU GPL v.3 license terms.
4. RESULTS AND DISCUSSION
The rst run, which is with few layers of particles in the vertical direction (only ten rows of
particles, in contrast to 350 columns horizontally), displays a reasonable exact solution except
Ingeniería 31(1): 98-111, Enero-junio, 2021. ISSN: 2215-2652. San José, Costa Rica DOI 10.15517/ri.v31i2.45850
105
for the leading front segment, which has noticeable lag; more than 15 %, in comparison with the
analytical model [1]. This may be due to the fact that this area is made up of a single layer of par-
ticles, as shown in Fig. 3, where a good number of them are not affected by enough neighbors to
actually move at a better pace. Also, the repulsive forces from the boundary virtual particles are
large enough to become dominant over any inertial or viscous terms, forcing an almost perfect, but
articial, alignment of the rst layer of particles. By increasing the vertical resolution of the pro-
blem to twenty rows (a total of 14 000 particles), the quality of the resulting distribution is greatly
improved, as shown in Fig. 4 for the same general area at the chosen reference instant t = 5.0 s.
In this higher resolution model, similar problems to those for the lower resolution model can be
observed at the boundary.
Figure 3. Detail of the water front for dam break with 3500 particles at t = 5.0 s (model OpenLR, 1 m x 35 m).
In these images, it is very important to keep in mind the nature of the model, which is of a high
(1:50) aspect ratio; as such, there is correspondingly high graphical distortion in the vertical scale
in relation to the horizontal scale. Some particles may appear overlapping, but this is not a particle
penetration problem, it is only a consequence of the high horizontal density of the plot, necessary
to show such a long model in a single graph.
Figure 4. Detail of the water front for dam break with 14 000 particles at t = 5.0 s (model OpenHR, 1 m x 35 m).
MONGE Y SERRANO: Flow under long two-dimensional dam sluice gate using WSPH
106
The complete prole for that instant is directly compared in Fig. 5 with the analytical solution
for the SWE developed by [1] and applied to the dam failure with the specied dimensions accor-
ding to [2]. Since the analytical model implies a 1 m x 50 m uid, but the model labeled OpenHR
is for a 1 m x 35 m uid, the data has been purposefully displaced -15 m in the horizontal direction.
Figure 5. Calibration of hydraulic prole for dam break with 14 000 particles at t = 5.0 s.
The apparent rounding at the point of discontinuity at the surface close to x = 20 m cannot be
considered analytical model imprecision or a defect of the numerical model, but rather an irrele-
vant property of the analytical solution; a real linear uid would not display a sharp edge, so it is
not of interest that this is reproduced in the numerical model. The comparatively minor inaccura-
cies related to full channel depth (1 m) that need to be addressed are a small lead (a maximum of
+4,8 % for position x = 27 m) in the upper half of the uid and a comparable lag in the lower half
(a maximum of -5,1 % for position x = 44 m).
However, the differences are not enough to reject the validity of the long channel dam failure
model and additional resolution did not improve precision, so the dam with an opened gate was
modeled with this same resolution but for the full 1 m x 50 m uid. The gate opening measures 0.20
m starting from the bottom of the dam, which amounts to 20 % of the nominal dam height. For this
situation, an analytical solution is also available [1], and the vena contracta depth and position is
checked with experimental data [2]. Figs. 6 to 8 show position elds for the area of interest of the
model for three signicant instants.
Ingeniería 31(1): 98-111, Enero-junio, 2021. ISSN: 2215-2652. San José, Costa Rica DOI 10.15517/ri.v31i2.45850
107
Figure 6. Detail of hydraulic prole for 20 % sluice gate with 20 000 particles at t = 0.5 s.
Figure 7. Detail of hydraulic prole for 20% sluice gate with 20 000 particles at t = 2.5 s.
Figure 8. Detail of hydraulic prole for 20% sluice gate with 20 000 particles at t = 5.0 s.
During the rst few instants (t < 1.2 s), the particles are very regularly aligned in a square array,
almost as if they were static, with the exception of the open discharge section (the open gate) and
the particles immediately adjacent to boundaries, as in Fig. 6. Later, the distribution of particles is
MONGE Y SERRANO: Flow under long two-dimensional dam sluice gate using WSPH
108
more irregular due to the complex dynamics near the gate, and uid accumulates in the discharge
zone, which drops in level at that location as in Figs. 7 and 8.
The effect of vena contracta is denitely present, and its depth and position coincide almost
exactly with the analytical solution as seen for x = 50 m in Fig. 8, as the aperture is 0.2 m, while
the discharge jet is initially 0.12 m deep (analytical solution gives a depth of 0.13 m), then expands
to about 0.16 m deep for a short distance (it is not yet a hydraulic jump as there is not enough uid
downstream to provide such a boundary) and then at an average depth of about 0,11 m before redu-
cing its height in a manner similar to a conventional water head owing over a dry bed.
Figure 9. Calibration for the hydraulic prole for 20% sluice gate with 20 000 particles at t = 5.0 s.
A comparison with the analytical solution for this dataset is shown in Fig. 9 where there is a good
degree of correspondence (2.8 % RMS relative error to full channel depth for the model between 30
m < x < 60 m, with a maximum error of -7.1 % at x = 42 m) except for the area with only one layer
of particles available to model the water front (occurring between 62 m < x < 80 m at t = 5.0 s). As
in the low-resolution dam break case with 3 500 particles, there are not enough numerical energy
transferal processes between particles to properly model the continuum. The error in that area is
large, both in qualitative and quantitative terms, and part of the reason is because of the virtual
particle equations that model the boundary through a strong repulsive force. Since this ends up
being a series of particles in a line, there is not much of a chance of energy transfer thorough vis-
cous shear stress or the increment of local density which would in turn contribute to local pressure
through the equation of state.
As in previous cases, this can be partially improved by increasing the resolution, which obviously
has a high computational cost not only because of the number of particles, but the time step must
be decreased as well. As a reference, the model with 3 500 particles took only two hours of com-
putational time; the model with 14 000 particles, with the same size boundaries and equal physical
time took more than thirty hours of computer time to complete. On the other hand, additional tuning
of the boundary repulsive force could also improve precision somewhat without allowing particle
leaks, because the rst layer of particles has practically no vertical movement but makes up almost
Ingeniería 31(1): 98-111, Enero-junio, 2021. ISSN: 2215-2652. San José, Costa Rica DOI 10.15517/ri.v31i2.45850
109
a quarter of the 0.20 m gate aperture, so the boundary trigger distance r0 may introduce unwanted
distortions to gate ow condition since, to be effective, it must be in the order of magnitude of the
nominal particle size
x.
These results show a good correlation and have readily illustrated many characteristics of
applying a general-purpose numerical method to cases with extreme proportions. This paves the
way for the proposal of future developments for this type of application and gives pointers to tune
applicable numerical model parameters of the WSPH method when modeling long, shallow cases
of open channel ow.
CONCLUSIONS
The application of the WSPH numerical scheme to a classical two-dimensional hydraulics pro-
blem proved to be effective to model ow over a high aspect ratio according to the special cases
explored in this study, where mean RMS error was well under 5 % with a maximum of 7 % in
critical areas in relation to upstream channel depth. Local error was in the order of 10 % (1 % in
relation to upstream channel depth) at the gate discharge where the vena contracta forms. Vertical
resolution and the repulsive force that acts at the boundaries both affect precision, and this means
it can be compensated for in similar models. Other particularly interesting aspects that came from
this study were:
The numerical model exhibited the phenomenon of vena contracta properly, and further
along, a small hydraulic jump in at the gate outlet, both of which are expected according to
the analytical model and the experimental observations. Although precision was acceptable
considering how low vertical resolution was in that area (one to four layers of particles, but
only 10 % local error in predicted depth in relation to analytical and experimental information
for that locus), the repulsive forces from the virtual particles used to model the solid boun-
dary, also contributed to diminish precision.
In the discharge area close to the gate, artifacts similar to physical waves and vortices appea-
red, such as those that appear in open channel ow phenomena. At the entrance to the gate,
since vertical resolution is that of at least twenty particle layers, wavelengths can be mea-
sured to determine whether they are numerical or physical in origin. At the discharge, since
resolution is so low, any periodic irregularities will have a numerical origin and cannot be
considered as a replica of vortexes that may occur in the real uid. This strategy can be very
useful to quickly validate the behavior of a numerical model, since wavelengths in channel
ow cases can be readily calculated with Froude number corrected wavelength equations or
measured in a relatively straightforward manner with physical experiments.
There was excellent correlation between the results of the numerical model, even in the
absence of stabilization algorithms or turbulence models, so it is likely that, for cases with
a high length-to-depth ratio, future studies will not need to assess its impact or consider its
incorporation. A more precise correlation should be obtained through good enough resolution,
MONGE Y SERRANO: Flow under long two-dimensional dam sluice gate using WSPH
110
since the dynamic phenomena related to SWE are so gradual that correctors for vorticity or
impact have little or no effect.
Low resolution in the water front is the main problem in the models studied, and this can
be a serious issue when modeling sediment carry-over. In these cases, interaction with the
boundary layer is key to its accuracy, so a different type of boundary particle will have to be
adopted or developed to properly interact with the uid particles.
Even when sediment transport is not relevant, to improve the accuracy of SPH models of
shallow water problems without requiring computationally costly increments of the number
of particles, future studies may consider creating a new type of boundary particle customized
specially for SWE-type problems.
The SPH scheme is well-adapted to many different open-ow hydraulics problems with
geometric discontinuities such as the general case of ow channels or partially lled piping
systems, and the convenience of not having to necessarily elaborate surface detection algori-
thms or dynamic remeshing, develops great advantages in comparison with mesh-dependent
methods such as FVM or VOF, and even hybrid approaches can be developed to attempt to
meld the positive aspects of each methodology. However, the proposal and modication of
different solid or elastic boundary equations have proven in many instances that it is one of
the most critical aspects of this method that is continually evolving and is far from being a
standard procedure applicable to any case of uid ow.
ACKNOWLEDGEMENTS
The authors would like to thank the anonymous peer reviewers and editors for their insight-
ful suggestions and corrections, which greatly improved readability, style and completeness. We
appreciate the time they invested in reviewing this paper and the consideration of the personal mark
and purpose with which it has been written.
REFERENCES
[1] L. Cozzolino, L. Cimorelli, C. Covelli, R. Della Morte, and D. Pianese, “The analytic solution of the
shallow-water equations with partially open sluice-gates: the dam-break problem”, Advances in Wa-
ter Resources, vol. 80, pp. 90–102, Jun. 2015, doi:10.1016/j.advwatres.2015.03.010.
[2] A. Dena and F. M. Susin, “Hysteretic behavior of the ow under a vertical sluice gate”, Physics of
Fluids, vol. 15, no. 9, pp. 2541–2548, Jul. 2003, doi:10.1063/1.1596193.
[3] M. Khanpour, A. R. Zarrati, and M. Kolahdoozan, “Numerical simulation of the ow under sluice
gates by SPH model”, Scientia Iranica, vol. 21, pp. 1503–1514, Oct. 2014.
[4] C.P. Sun, D.L. Young, L.H. Shen, T.F. Chen, and C.C. Hsian, “Application of localized meshless
methods to 2D shallow water equation problems”, Engineering Analysis with Boundary Elements,
vol. 37, pp. 1339–1350, Nov. 2013, doi:10.1016/j.enganabound.2013.06.010.
Ingeniería 31(1): 98-111, Enero-junio, 2021. ISSN: 2215-2652. San José, Costa Rica DOI 10.15517/ri.v31i2.45850
111
[5] S. Gu, X. Zheng, L. Ren, H. Xie, Y. Huang, J. Wei, and S. Shao, “SWE-SPHysics Simulation of
Dam Break Flows at South-Gate Gorges Reservoir”, Water. vol. 9, art. 387, May 2017, doi:10.3390/
w9060387.
[6] S. Gu, F. Bo, M. Luo, E. Kazemi, Y. Zhang, and J. Wei, “SPH Simulation of Hydraulic Jump on
Corrugated Riverbeds”, Applied Sciences. vol. 9, art. 436, 2019.
[7] J. J. Monaghan, “Smoothed particle hydrodynamics”, Annual Review of Astronomy and
Astrophysics, vol. 30, pp. 543–574, 1992, doi: 10.1146/annurev.aa.30.090192.002551.
[8] G.H. Cottet, “Multi-physics and particle methods”, presented at the Second MIT Conference on
Computational Fluid and Solid Mechanics, Cambridge, United States of America, June 17-20, 2003,
pp. 1296-1298. doi:10.1016/B978-008044046-0.50319-5.
[9] R. A. Gingold, and J. J. Monaghan, “Smoothed particle hydrodynamics: theory and application to
non-spherical stars”. Monthly Notices of the Royal Astronomical Society, vol. 181, pp. 375–389,
1977, doi:10.1093/mnras/181.3.375.
[10] G. R. Liu and M. B. Liu, Smoothed Particle Hydrodynamics: A Meshfree Particle Method. Upper
Saddle River, N.J.: World Scientic, 2003.
[11] D. Violeau, Fluid Mechanics and the SPH Method: Theory and Applications. Oxford: Oxford Uni-
versity Press, Oxford, 2012.
[12] J. J. Monaghan, “Smoothed particle hydrodynamics and its diverse applications”, Annual Review of
Fluid Mechanics, vol. 44, pp. 323–346, January 2012. doi: 10.1146/annurev-uid-120710-101220.
[13] R. Vacondio, et al., “Grand challenges for smoothed particle hydrodynamics numerical schemes”.
Computational Particle Mechanics, 2020, doi:10.1007/s40571-020-00354-1.
[14] J. J. Monaghan, “SPH without a tensile instability”, Journal of Computational Physics, vol. 159, pp.
290–311, April 2000, doi:10.1006/jcph.2000.6439.