E-ISSN: 2528-388X P-ISSN: 0213-762X INERSIA Vol. No. December 2023 Application of Smoothed Particle Hydrodynamics Method for Tsunami Force Modeling on Building with Openings Putra Anggita*. Radianta Triatmadja, and Nur Yuwono Department of Civil and Environmental Engineering. Universitas Gadjah Mada. Yogyakarta 55281. Indonesia ABSTRACT Keywords: DualSPHsyics SPH Tsunami 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. 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 . as a crucial numerical parameter in DualSPHysics modeling. A smaller particle size contributes to modelAos accuracy. The determination of the particle size must account for modelAos shortest characteristic . 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 . This is an open access article under the CCAeBY license. Introduction Tsunami waves can cause damage to coastal The destructive power of tsunami waves on an infrastructure is influenced by many factors, including. characteristics of tsunami waves . orce, wave speed and heigh. , building characteristics . hape, pores, and building dimension. , 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 . 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 . which is proven by Lukkunaprasit et al. 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 Tsunami is a natural phenomenon that is difficult to 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 . is one of them . 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 Nevertheless. SPH method also has its This method requires relatively long computation time due to its calculation properties. address this computing problems DualSPHysics was *Corresponding author E-mail: putraanggita@mail. https://dx. org/10. 21831/inersia. Received November 1st, 2022. Revised November 20th, 2023. Accepted November 20th, 2023 Available online December 31st, 2023 Putra Anggita, et al. INERSIA. Vol. No. December 2023 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. 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 In this paper, the tsunami force acting on building with opening which calculated using the SPH method were studied. The versatility of SPH is evident in its broad applicability across different scales, ranging from small systems to large scale and even astronomical scale It can effectively handle both discrete and continuous fluid systems SPH applications cover various physical aspects of fluids such as viscosity, external force, internal force, density, and so on. 2 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 According to Dominguez et al. , in the last 10 years, the use of hardware such as graphics processing units (GPU. 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 . Literature Review 1 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 . Gingold, and Monaghan . initially proposed it to model astrophysical 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 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 The SPH formulation, being unaffected by arbitrary particle arrangement, enables this method to handle significant deformations without requiring surface This capability is widely regarded as the most appealing aspect of this approach . In collaborative research by Universidade de Vigo. University of Manchester. University of Lisbon. University di Parma. Flanders Hydraulics Research. Universitat Politycnica de Catalunya, and New Jersey Institute of Technology, the open source code DualSPHysics was developed to address SPH computing DualSPHysics is a set of C . CUDA and Java code developed to solve problems in real engineering 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 . 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. This method offers several advantages in simulations, which can be summarized as follows: 3 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 . Wave force on cubes with several variations of opening were observed Figure 1. The opening configuration . in the building is defined as Equation 1. SPH allows for the simulation of free surface flow, enabling the modeling of fluid behavior at boundaries and Through ongoing advancements. SPH method has achieved notable improvements in adaptability, stability, and accuracy, making it a reliable choice for fluid n=Ao/Ab . Where Ao is the area of the openings and Ab the front area of the box. INERSIA. Vol. No. December 2023 Putra Anggita, et al. In this context, where r represents the position vector. W denotes the weighting function . , and h represents the smoothing length, which governs the extent of the influence domain E. The smoothing length . is a variable that relies on both space and time and is utilized to determine the domain of influence for neighboring The value of h must be greater than the initial distance between the particles . In addition, as a general guideline, dp can be set to approximately 1/10 of the modelAos shortest characteristic length . Figure 1. Tsunami models on building with various openings studied in Triatmadja and Nurhasanah . In this study, tsunami waves were modeled through a 24 y 45 y 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. Equation 2 can be transformed into discrete as Equation 3, this formulation is an approximation of a function on a particle a. ya ya. = Oc ycoyca yca ycOycayca yuUyca yca where mb is the particle mass. Ab the particle density and Wab the weighting function . 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. Figure 2. Force reduction due to openings presented in Triatmadja and Nurhasanah . SPH Theory ONya. = Oc ycoyca 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 . uch as position, pressure, velocity, and densit. for a particle by interpolating the values from neighboring particles. represent the transition from a continuous medium . he flui. to a discrete medium . he particle. , 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 . The key characteristics of SPH method, which relies on integral interpolant, can be summarized as follows. yca yayca ONycOycayca yuUyca 2 Kernel function The selection of the weghting function . 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 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 5 to Equation 9. Positivity: c Oe yc A . Ou 0 inside the domain E 1 Integral interpolant Compact support: SPH method performs calculations using integral interpolation known as kernel approximation. The underlying principle of approximating any function A. is as follows ya. = OE ya. cA)ycO. c Oe ycA. yccycA . A E. = ycO. c Oe yc . Normalization: O ya. c A )ycO. c Oe yc A . yccyc A = 1 0 inside the domain E Putra Anggita, et al. INERSIA. Vol. No. December 2023 Delta function behavior: where superscripts and are coordinate directions. g is gravitational acceleration. E is particle stress. v is particle e is internal energy per unit mass. A is shear strain rate (A OI 0. and ij is Monaghan's artificial viscosity. 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 . lim ycO. c Oe yc A . yccyc A = yu. c Oe yc A ) hIe0 Monotonically decreasing behavior, . A E. c Oe yc . 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 . is taken into account for particle a. The description of the particle and the domain of influence is shown in Figure 3. From the set of Equation . the particle force equation can be derived as Equation 11 and Equation 12 . ycyycyceycycycycyce yaycn ycyycn ycyyc = yuN Oc ycoyc ycycn ycyc 2 yu ycO. yc ycycnycycaycuycycnycyc yaycn = Oe Oc ycoyc yc Where rij = xi - xj, is the fluid viscosity coefficient. The pressure . is calculated by the constitutive Equation 13. ycyycn = ya . uUycn Oe yuU0 ) Where K is the fluid stiffness value and A0 is the initial density value. Based on equations . , the acceleration of particle i can be derived as Equation 14. Figure 3. Particle and domain influence presented in Dominguez et al. Kernel smoothing and particle approximation can be used to discretise partial differential equations. The SPH formulation is derived by spatially discretising the NavierStokes equations, so that they become a set of Ordinary Differential Equations that can be solved through time 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 . ycu yuiycOycnyc yccyuUycn = Oc ycoyc . cuycnyu Oe ycuycyu ) yu , yccyc yc=1 yuiycu ycaycn = ycyycyceycycycycyce yuUycn . aycn ycycnycycaycuycycnycyc yaycn yaycnyceycuycyceycycuycayco ) . where Fiexternal represents external forces such as body force and force due to contact. 3 Governing Equation The governing equations are described by the NavierStokes equations for a compresible fluid. The continuity and momentum equations in Lagrangian form can be written as Equation 15 and Equation 16. ycn dyuU = OeyuUyu. dyc yuyu yuyu ycu yuiycOycnyc yccycycnyu yua yua = Oc ycoyc ( ycn 2 Oe yc 2 ) yu , yccyc yuUycn yuUyc yuiycuycn yc=1 . yuyu yuyu ycu yccycycnyu yua yuiycO yuiycO yua = Oc ycoyc ( ycn . yuycnyc Oe yc . yuycnyc ) , yccyc yuUycn yuUyc yuiycu yuUycn yuUyc yuiycu yc=1 ycn ycy yuiycOycnyc , yuUycn = Oe ycn 2 ycnyc Ocycuyc=1 ycoyc . cycn Oe ycyc ) = Oe yuycE e f where d is the total or material derivative, v the vector velocity. A the density. P the pressure, e the dissipation terms and f the accelerations due to external forces, such as the gravity. ycu yco yccycuycnyu yc . cycn Oe ycyc )ycOycnyc , = ycycn yuA Oc yccyc yc=1 yuUyc INERSIA. Vol. No. December 2023 Putra Anggita, et al. 4 Tsunami Surge Force on Building with Openings can be modified to account for the effect of the opening . , the modified equation is written as Equation 19. ya = yayce yuU. Oe y. yuUyaAEaycO 2 The surge celerity caused by a dam break that represents a tsunami front can be expressed in the equation as Equation ycO = yaOoyciEa0 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. Where C is the celerity coefficient, g the gravitational acceleration and h0 the depth of the basin. 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. 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. ya = yaycn yuUyaAEaycO 2 where Ci is the impact coefficient that takes momentum change after the impact into consideration. Equation . Figure 4. Research flowchart. Putra Anggita, et al. INERSIA. Vol. No. December 2023 Figure 5. Dam break by Triatmadja and Nurhasanah . The DualSPHysics model in this study as Figure 5 was modeled based on physical model conducted by Triatmadja and Nurhasanah . 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 . 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. Table 1. DualSPHysics constants and parameters Parameter and constant Gamma value. Type and value Rho. A . g/m. 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 . and artificial viscosity (). The DualSPHysics constant and parameter executions are presented in Table 1. CFL number Coefficient value Coef. of sound Kernel type Wendland Step Algorithm Verlet Viscosity Artificial Particle size, dp . 01 & 0. Simulating cases with very small particles diameter . 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 SPH numerical simulation results with variation of dp on various opening configuration . are shown in Figure 6 to Figure 11. Figure 6. Comparison of force acting on the building without opening. Figure 7. Comparison of force acting on the building with n = 0. Results INERSIA. Vol. No. December 2023 Putra Anggita, et al. Figure 8. Comparison of force acting on the building with n = 0. Figure 9. Comparison of force acting on the building with n = 0. Figure 10. Comparison of force acting on the building with n = 0. Figure 11. Comparison of force acting on the building with n = 0. Maximum surge force due to tsunami wave hitting the building with a particle size of 0. 02 m were higher than that 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. 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. addition, simulations with smaller particle showed better Based on Figure 12, the experimental results showed a significant reduction in force due to openings configuration . 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 Figure 13 showed the relationship between n value and , where n is the opening configuration, is the force reduction due to opening configuration . F0 is the maximum force acting on building without opening, and F is the maximum force acting on building with opening. Figure 12. Maximum force acting on the building with various opening configuration . Putra Anggita, et al. INERSIA. Vol. No. December 2023 Figure 13. Relationship between force reduction and n 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. Figure 14. Building model with n = 0. Selection of the size of particle . 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 1/10 of the modelAos shortest characteristic length . For instance, in the building model with n = 0. the s value is 0. 037 m, which can serve as a reference for determining an appropriate dp value (Figure . Thus the maximum value of the comparison between dp and s can be determined as a reference for the DualSPHsyics 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 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 presented the results of DualSPHysics simulations with variaton of artificial viscosity value. Figure 15. Surge force profile of wave with = 1, 0. 1, and 0. 01 values INERSIA. Vol. No. December 2023 Putra Anggita, et al. Figure 16. Surge force profile of wave with = 0. 01, 0. 001, and 0. 0001 values Similar surge force profiles were obtained on simulations with = 0. 0001, 0. 001 and 0. 01 (Figure . 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 The more viscous the fluid, the slower the wave Further adjustments and calibration of DualSPHysics model against physical model are needed due to limitations in research. Conclusion The size of the particles . is a critical numerical parameter in DualSPHysics calculations. A smaller particle diameter leads to improved accuracy in the SPH When determining the value of dp, it is essential to take into account modelAos shortest characteristic . 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. References