Application of Smoothed Particle Hydrodynamics Method for Tsunami Force Modeling on Building with Openings

Smoothed Particle Hydrodynamics (SPH) serves as a numerical technique extensively employed for simulating free surface flow. The computational intricacy of the SPH method arises from the numerous computations of a particle's properties, derived from interactions with surrounding particles. To address this complexity, experts developed DualSPHysics. This study employs the SPH method, specifically the DualSPHysics application, for tsunami modeling. To accurately represent tsunami characteristics, precise numerical parameters are essential for numerical modeling. This research provides valuable insights into optimizing numerical parameters for accurate SPH simulations. Therefore, the research aims to identify the exact values of crucial parameters in DualSPHysics model. Validation of numerical calculations involves comparing the tsunami forces, as simulated by DualSPHysics, with secondary data obtained from physical experiments results. The findings highlight the significance of particle size (dp) as a crucial numerical parameter in DualSPHysics modeling. A smaller particle size contributes to model’s accuracy. The determination of the particle size must account for model’s shortest characteristic (s). According to simulations those have been carried out, it is recommended to set the maximum limit value of dp/s at 1/3.67 to achieve precise calculation. Furthermore, the DualSPHysics simulation demonstrates a reduction in force due to the opening configuration (n). `


Introduction
Tsunami waves can cause damage to coastal infrastructure.The destructive power of tsunami waves on an infrastructure is influenced by many factors, including; characteristics of tsunami waves (force, wave speed and height), building characteristics (shape, pores, and building dimensions), and environmental conditions around the infrastructure.around infrastructure.The attributes of a structural system have a significant impact on the ability of the structure to withstand the effects of tsunamis, earthquakes, and wind loading [1].An open system that allows water to flow with minimal resistance is one of the structural attributes of buildings that have shown good performance in overcoming the impact of past tsunamis [2] which is proven by Lukkunaprasit et al. [3].The interaction between tsunami waves and porous buildings was physically modeled.The tsunami wave forces acting on buildings with various openings variations were compared.The results showed that the pores significantly reduce the tsunami force acting on the building.
Tsunami is a natural phenomenon that is difficult to predict.One of the efforts in studying tsunami waves is through numerical modeling.In terms of cost and time efficiency, numerical modeling is often more profitable than physical models.With advances in technology and recent applications it is possible to model wave phenomena by numerical methods.The Smoothed Particle Hydrodynamics (SPH) method developed by Gingold and Monaghan (1977) is one of them [4].The SPH method has been widely used to model fluid dynamics, especially free surface flow.Besides its wide application, this method has good accuracy, stability and adaptability in numerical modeling.Nevertheless, SPH method also has its limitations.This method requires relatively long computation time due to its calculation properties.To address this computing problems DualSPHysics was developed.

164
DualSPHysics was developed based on method of Smoothed Particle Hydrodynamics (SPH), a numerical approach for modeling the movement of fluid surfaces.SPH employs a meshless particle system based on the Lagrangian model, treating the fluid as a collection of interconnected particles.These particles interact and integrate with one another in accordance with the principles of fluid conservation.[5].SPH implementation with DualSPHysics makes it possible to simulate real life problem such as tsunami.This study employs the SPH method, specifically the DualSPHysics application, for tsunami modeling.To accurately represent tsunami characteristics, precise numerical parameters are essential for numerical modeling.This research provides valuable insights into optimizing numerical parameters for accurate SPH simulations.Therefore, the research aims to identify the exact values of crucial parameters in DualSPHysics model.In this paper, the tsunami force acting on building with opening which calculated using the SPH method were studied.

Smoothed Particle Hydrodynamics (SPH) Method
Smoothed Particle Hydrodynamics (SPH) is a particlebased method without the need for a mesh structure, following the Lagrangian approach.Its origins can be traced back to the late 1970s when Lucy [6], Gingold, and Monaghan [4] initially proposed it to model astrophysical phenomena.This approach utilizes a group of particles positioned at irregular intervals to depict the movement of a fluid.Each of these particles possesses specific characteristics like energy, velocity, pressure, density, and mass.In this method, the properties of each particle for a specific time increment are approximated based on the properties of nearby particles within a defined influence area.The SPH formulation, being unaffected by arbitrary particle arrangement, enables this method to handle significant deformations without requiring surface treatments.This capability is widely regarded as the most appealing aspect of this approach [7].This method offers several advantages in fluid simulations, which can be summarized as follows: 1. SPH allows for the simulation of free surface flow, enabling the modeling of fluid behavior at boundaries and interfaces.. 2. Through ongoing advancements, SPH method has achieved notable improvements in adaptability, stability, and accuracy, making it a reliable choice for fluid simulations.
3. The versatility of SPH is evident in its broad applicability across different scales, ranging from small systems to large scale and even astronomical scale scenarios.It can effectively handle both discrete and continuous fluid systems 4. SPH applications cover various physical aspects of fluids such as viscosity, external force, internal force, density, and so on.

DualSPHsyics
The computational intricacy of the SPH method arises from the numerous computations of a particle's properties, derived from interactions with surrounding particles, so the computation time required is relatively long.From the root of the problem, the SPH community has conducted various researches to speed up the simulation of this method.According to Dominguez et al., in the last 10 years, the use of hardware such as graphics processing units (GPUs) for scientific computing has been common.SPH in particular, this method is compatible with GPUaccelerated simulations thanks to vector lists of particles and their interactions.Multi-GPU based SPH codes are now available [8].
In collaborative research by Universidade de Vigo, University of Manchester, University of Lisbon, Università di Parma, Flanders Hydraulics Research, Universitat Politécnica de Catalunya, and New Jersey Institute of Technology, the open source code DualSPHysics was developed to address SPH computing problems.DualSPHysics is a set of C++, CUDA and Java code developed to solve problems in real engineering practice.This code was developed to study surface free flow phenomena where the Eulerian method is difficult to apply, such as waves or the impact of dambreaks on offshore structures.More details on the DualSPHysics program can be seen in [8].DualSPHysics is a free, opensource SPH code released online.DualSPHysics is an openly accessible, SPH software package that has been made available on the internet as an open-source code.

Tsunami force on building with openings
The effect of the openings of the building on the acting tsunami force were studied by Triatmadja and Nurhasanah [9].Wave force on cubes with several variations of opening were observed Figure 1.The opening configuration (n) in the building is defined as Equation 1.

n=Ao/Ab
(1 Where Ao is the area of the openings and Ab the front area of the box. 165 Figure 1.Tsunami models on building with various openings studied in Triatmadja and Nurhasanah [9].
In this study, tsunami waves were modeled through a 24 × 1.45 × 1.5 m flume which is equipped with a dambreak based wave generator.From this study it can be concluded that the reduction of tsunami force acting on the building is not linearly related to the value of n.This is shown in Figure 2.

SPH Theory
The fluid is considered as a collection of particles or nodes referred to as particles.The equations of motion in fluid dynamics are integrated individually for each point using the Lagrangian approach.This process involves calculating the values of essential physical variables (such as position, pressure, velocity, and density) for a particle by interpolating the values from neighboring particles.To represent the transition from a continuous medium (the fluid) to a discrete medium (the particles), a kernel function is employed.This function has a limited range, determined by the distance required to resolve particle interactions, and it effectively describes the transformation from the fluid's continuous nature to the discrete nature of the particles [10].The key characteristics of SPH method, which relies on integral interpolant, can be summarized as follows.

Integral interpolant
SPH method performs calculations using integral interpolation known as kernel approximation.The underlying principle of approximating any function A(r) is as follows In this context, where r represents the position vector, W denotes the weighting function (kernel), and h represents the smoothing length, which governs the extent of the influence domain Ω.The smoothing length (h) is a variable that relies on both space and time and is utilized to determine the domain of influence for neighboring particles.The value of h must be greater than the initial distance between the particles (dp).In addition, as a general guideline, dp can be set to approximately 1/10 of the model's shortest characteristic length [11].
Equation 2 can be transformed into discrete as Equation 3, this formulation is an approximation of a function on a particle a.
where mb is the particle mass, ρb the particle density and Wab the weighting function (kernel).
In SPH method, the derivative of a function is calculated analytically, this is one of the advantages of this method.
In contrast to other methods, the derivative of a function is calculated using the distance between a point particle and the particles around it.The derivative of this interpolation can be obtained by ordinary differentiation shown in Equation 4.

Kernel function
The selection of the weghting function (kernel) is very important in the SPH model because it affects the performance of the SPH model.Several conditions that must be met such as positivity, compact support and normalization.In addition, the value of Wab must fulfill monotonically decreasing behavior with increasing distance between particle a and the surrounding particles and act as a behavior of the delta function.Above conditions can be summarized in the form of Equation 5to Equation 9.
Kernel function relies on two factors: the smoothing length, h, and the non-dimensional distance between particles, q = r/h, where r represents the distance between particles a and b.The parameter h determines the scope of the region surrounding particle a, within which the contributions of all neighboring particles are considered.
It essentially controls the size of the area in which the influence of other particles (domain) is taken into account for particle a.The description of the particle and the domain of influence is shown in Figure 3. Kernel smoothing and particle approximation can be used to discretise partial differential equations.The SPH formulation is derived by spatially discretising the Navier-Stokes equations, so that they become a set of Ordinary Differential Equations that can be solved through time integration.By substituting the SPH approximation for a function and its derivatives into the partial differential equations governing the physics of fluid flow, the discretisation of the flow-forming equations can be written as Equation 10 where superscripts α and β are coordinate directions; g is gravitational acceleration; σ is particle stress; v is particle velocity; e is internal energy per unit mass; ε is shear strain rate (ε ≅ 0.5) and Πij is Monaghan's artificial viscosity.In the analysis of the interaction between waves and structures, Πij can prevent the non-physical shock of the solution results in the collision region and effectively prevent the non-physical penetration of particles when they are close to each other.The role of artificial viscosity is to smoothen the shocks over several particles and to allow the simulation of viscous dissipation, the transformation of kinetic energy into heat.Therefore, to consider artificial viscosity, an artificial viscous pressure term Πij is added [13].
From the set of Equation ( 10) the particle force equation can be derived as Equation 11and Equation 12 [14].
Where rij = xi -xj, μ is the fluid viscosity coefficient.The pressure (pi) is calculated by the constitutive Equation 13.
Where K is the fluid stiffness value and ρ0 is the initial density value.Based on equations ( 11) and ( 12), the acceleration of particle i can be derived as Equation 14.
where Fi external represents external forces such as body force and force due to contact.

Governing Equation
The governing equations are described by the Navier-Stokes equations for a compresible fluid.The continuity and momentum equations in Lagrangian form can be written as Equation 15 and Equation 16. and where d is the total or material derivative, v the vector velocity, ρ the density, P the pressure, Γ the dissipation terms and f the accelerations due to external forces, such as the gravity.167

Tsunami Surge Force on Building with Openings
The surge celerity caused by a dam break that represents a tsunami front can be expressed in the equation as Equation 17.
Where C is the celerity coefficient, g the gravitational acceleration and h0 the depth of the basin.
Based on Newton's second law and assuming a constant surge celerity just before the wave hitting the wall, the force equation due to the impact on the wall of a building can be written as Equation 18.
where Ci is the impact coefficient that takes momentum change after the impact into consideration.Equation ( 18) can be modified to account for the effect of the opening [8], the modified equation is written as Equation 19.
Where Cf is the bulk coefficient of the force that takes into account the drag and impact forces on the building and δ is the force reduction coefficient which is a function of the opening of the building.

Method
Tsunami wave modeling with the DualSPHsyics application conducted in this study aims to determine the correct numerical parameter values in the implementation of the Smoothed Particle Hydrodynamics (SPH) method in accordance with real conditions in the field.For the implementation of the research to be more focused, the research flow chart in Figure 4.The DualSPHysics model in this study as Figure 5 was modeled based on physical model conducted by Triatmadja and Nurhasanah [9 ].The dam break system is used as the basis for tsunami wave generation in the experiment.The flume is divided into two parts, namely upstream and downstream of the channel.The upstream part of the channel functions as a reservoir (basin) with a certain depth and was divided by a sluice gate.Meanwhile, the downstream part of the channel is conditioned as a surface that has not been crossed by a tsunami with a lower water level than the upstream part.In this study, debris flow was not considered in the tsunami model.

Results
Numerical parameter calibration test was carried out in this study to understand the effect of the selection of numerical parameter on the results of the DualSPHsyics simulation.Thus the experimental data can be approximated by numerical simulation.Numerical parameters which were tested are the particle size (dp) and artificial viscosity (α).
The DualSPHysics constant and parameter executions are presented in Table 1.requires GPUs with very high specifications.In addition, the simulation will take a relatively long time.In this case, it took eight hours of DualSPHysics calculation time to simulate a tsunami wave with a duration of four seconds and a particle size of 0.01 m.Thus, simulations with particle size smaller than 0.01 m were not carried out in this study.SPH numerical simulation results with variation of dp on various opening configuration (n) are shown in Figure 6 to Figure 11.Maximum surge force due to tsunami wave hitting the building with a particle size of 0.02 m were higher than that of 0.01 m as shown in Figure 5 to Figure 6.The comparison of the tsunami maximum force between DualSPHsyics simulations and experimental results summarized in Figure 5 to Figure 7.
Based on Figure 12, the experimental results showed a significant reduction in force due to openings configuration (n).Meanwhile, in the DualSPHsyics simulation, the force reductions were not as significant.In general, the results of the DualSPHsyics simulations had the same pattern as the experimental results.However, there were fairly large deviations from the experimental results, both in the DualSPHsyics simulation with dp 0.01 and 0.02 m.In addition, simulations with smaller particle showed better accuracy.
Figure 13 showed the relationship between n value and δ, where n is the opening configuration, δ is the force reduction due to opening configuration (n), F0 is the maximum force acting on building without opening, and F is the maximum force acting on building with opening.The relationship between n value and δ showed a curve that tends to be linear compared to experimental results and value deviations can be found in several points.The simulation presented a good agreement with experimental results in the range between n = 0.4 and n = 0 where the reductions in force due to opening were fairly similar to the experimental values.Selection of the size of particle (dp) holds significant importance in DualSPHysics modeling.It has a notable impact on the overall simulation, including accuracy, wave behavior, and time of simulation.dp can be set to approximately of the model's shortest characteristic length (s).For instance, in the building model with n = 0.4, the s value is 0.037 m, which can serve as a reference for determining an appropriate dp value (Figure 14).Thus the maximum value of the comparison between dp and s can be determined as a reference for the DualSPHsyics simulation.According to DualSPHysics simulations those have been carried out, it is recommended to set the maximum limit value of dp/s at 1/3.67 to achieve precise calculation.
Further simulations were carried out by varying the value of artificial viscosity (α), to see the wave force patterns acting on the building from various values of α.In this case, particle size of 0.01 m was used and solid building acted as the obstacle.Figure 15   Similar surge force profiles were obtained on simulations with α = 0.0001, 0.001 and 0.01 (Figure 16).Based on this results, it can be seen that the simulation results showed good stability with α values less than or equal to 0.01 and showed good agreements with experimental results.The maximum values of the surge force are closest to the experimental results.Whereas, between α = 0.1 and 1, the surge patterns indicated significant differences.This shows that the difference in α value affects the viscosity of the fluid, the higher the α value, the more viscous the fluid will be.The more viscous the fluid, the slower the wave propagation.Further adjustments and calibration of DualSPHysics model against physical model are needed due to limitations in research.

Conclusion
The size of the particles (dp) is a critical numerical parameter in DualSPHysics calculations.A smaller particle diameter leads to improved accuracy in the SPH simulation.When determining the value of dp, it is essential to take into account model's shortest characteristic (s).Setting the maximum limit value of dp/s at 1/3.67 is recommended to achieve precise calculation.Furthermore, the value of the artificial viscosity coefficient (α) affects the viscosity of the fluid.In this study, α values in the range of 0.01 to 0.0001 showed a good agreement with the physical model and can be used as a reference on simulation using DualSPHsyics.

Figure 3 .
Figure 3. Particle and domain influence presented in Dominguez et al. [8]

Figure 6 .
Figure 6.Comparison of force acting on the building without opening.

Figure 7 .Figure 8 .
Figure 7.Comparison of force acting on the building with n = 0.075.

Figure 9 .
Figure 9.Comparison of force acting on the building with n = 0.4.

Figure 10 .
Figure 10.Comparison of force acting on the building with n = 0.6.

Figure 11 .
Figure 11.Comparison of force acting on the building with n = 0.81.

Figure 12 .
Figure 12.Maximum force acting on the building with various opening configuration (n).

Figure 13 .
Figure 13.Relationship between force reduction and n
presented the results of DualSPHysics simulations with variaton of artificial viscosity value.

Table 1 .
DualSPHysics constants and parameters