Sloshing in medium sized tanks under earthquake load

This paper presents results obtained by numerical study of sloshing effects in medium sized rectangular tanks subjected to a long-period seismic excitation. The Smoothed Particle Hydrodynamics (SPH) method is used, and the wave displacement and pressure are analysed. The intensity of sloshing effects depends not only on the peak ground acceleration but also on the dominant frequency of ground motion. Excessive sloshing displacements and pressures increase when the natural period of the contained liquid approaches the dominant period of ground excitation.


Introduction
Liquid storage tanks are used to store a wide array of liquids including industrial chemicals, fuel and potable water, all of which is essential to meet the needs of the society as well as those of the industry.Failure of these structures during earthquakes can have significant adverse effects on the environment, in addition to financial losses and disruption of services to the community.The roof of a potable water storage tank collapsed at Christchurch, New Zealand, during the Mw 7.1 Darfield earthquake of 2010, due to uplift forces caused by sloshing, which refers to the free surface displacements of the contained liquid in form of standing waves during earthquakes [1].Extensive sloshing due to ground shaking also caused partial collapse of the reinforced concrete roofs of two reservoirs at the Feng Yuen Water Treatment Plant during the Mw 7.4 Chi-Chi, Taiwan earthquake of 1999 [2].In addition to these recent examples, field reconnaissance reports in the petroleum and gas industry [3][4][5][6][7][8] also suggest that liquid storage tanks are quite susceptible to structural damage during earthquakes.The foundations and walls of liquid storage tanks should be able to resist the base shear, overturning moments, and the hydrodynamic pressures generated by earthquake load, as well as the hydrostatic load.Furthermore, enough freeboard must be provided to avoid spilling and damage to the roof due to impact loads caused by excessive sloshing.Hence, various analytical methods have been developed since 1960s to make sure that these structural requirements are met in the design of liquid storage tanks.Housner's mechanical analogue model [9] is perhaps one of the earliest and most widely used tools for explaining the sloshing dynamics fundamentals in rectangular and cylindrical liquid storage tanks subjected to horizontal ground motion.In Housner's model, the liquid mass is divided into two components, namely the rigid mass which moves in unison with the tank wall, and the convective masses that are attached to the tank wall with springs.The rigid component of the liquid largely accounts for the hydrodynamic pressures acting on the tank wall and foundations, whereas the convective component accounts for the sloshing displacements.Numerical methods capable of solving the contained liquid motion equations, for complex tank geometries, multi-layered fluids, and viscous fluids, have been developed in parallel with the advances in computing technology.Frandsen [10] used the Finite Differences Method (FDM) for sloshing simulations in 2D rectangular tanks.Nakayama and Wachitzsu applied the Finite Element Method (FEM) for nonlinear sloshing [11].However, the use of the FD or FEM in free surface flow modelling requires another numerical treatment such as the Level Set or Volume of Fluid methods to track the free surface profile at each time step [12].Hence, traditional numerical methods require Eulerian meshes to be updated at each time step, which makes the computing expensive.The Smoothed Particle Hydrodynamics (SPH) method is a computationally efficient meshless method initially developed for astrophysical problems [13,14].Later on, its use became more widespread, especially in the fields of the fluid structure interaction modelling, underwater explosions, and ship mechanics.The SPH method is used in this study to investigate the free surface sloshing displacements and the hydrodynamic pressures in a two dimensional rectangular water tank subjected to the long-period earthquake ground motions.The mathematical basis of the SPH model is described in the section entitled "Method and analysis", while the analysis results of sloshing displacements and hydrodynamic pressures in an example tank, subjected to a set of 6 horizontal earthquake acceleration time-history records for 7 different fill levels, are presented in the "Numerical tests" section of this paper.

Smoothed Particle Hydrodynamics method
The foundation of the SPH method depends on a unity, even, smooth and compact Kernel Function, W(x,y), such that any arbitrary function, j(x,y) is approximated with the given identity: (1) where X = (x,y) is the position variable and S is the influence region, i.e. the smoothed region of the Kernel function.Spiky Kernel [15] is used in this study to overcome numerical challenges, i.e. clustering due to collocation nature of the SPH method.For two dimensional geometries: (2) where α = 10/ph 5 , r is the distance between particles, and h = 4∆x is the smoothing distance.Since the particles adjacent to the boundaries have an initial spacing of ∆x/2 or ∆y/2 with the corresponding boundary, the uniform initial spacing between the particles in x and y directions are determined as ∆x = L/N x , ∆y = H w /N y where ∆x = ∆y.L is the length of the tank in x, and H w is the fill level of the liquid inside the tank.N x and N y are necessarily set to be consistent to give uniformly distributed particles initially, and then the total number of particles is N x N y .j at a given particle is determined after inserting W(r,h) in the Kernel approximation eq. ( 1): (3) where j i = j (r i ), j j = j (r j ), W ij =W(r ij ,h), r ji =r j -r i is the position vector, and r ij =|r ji | is the distance from i to j (Figure 1).m j and ρ j represent the mass and the density of the particle j.Reformulation of the conservation of mass and conservation of momentum by using the Kernel approximation results in the eqs.( 4) and ( 5): where u ij = u i -u j is the velocity vector, P denotes pressure, and μ is the coefficient of dynamic viscosity.f = (ρa G , ρg) denotes the force vector, a G is the horizontal earthquake ground acceleration in the x direction, and g represents the gravitational acceleration with the magnitude of -9.81 m/s 2 in the y direction.The position vector of particles is determined by the kinematics principle, eq. ( 6): The viscosity term includes the second order derivative of the Kernel function, which requires a rather elaborate numerical treatment.As seen in eq. ( 5), the Kernel Function is approximated in this study with the first order derivation term described in [16].Pressure is approximated with the Batchelor's formula [17], suitable for the incompressible flow: (7) where B = 200 • ρ 0 gH w /g and g = 7, keeping the density fluctuations at approximately 1 % [18].The SPH method does not require a specific definition of free surface boundary conditions if no other fluid layers exist in the tank.Solid wall boundaries are assumed to be slipping and impermeable.A number of approaches including the boundary forces (repulsive forces), boundary particles or ghost/image particles, can be used to apply boundary conditions at walls [19].In this study, ghost (image) particles are distributed out of the walls in the range of h, which are regenerated at each time step.The location of the ghost particles is set to be symmetrical with the inner particles, and all information is assigned as being identical with the inner particle, except for the normal velocity component.Since the wall is impermeable, the normal velocity is assigned in the opposite direction.The ghost particle approach by itself is not adequate to properly reflect the impermeability of the wall.Therefore, the collision between a particle and the boundary is also taken into account.The tangential component of velocity remains the same but the normal component of velocity is diminished according to the type of collision.In this study, the collision is considered as plastic, and hence; where the superscript n denotes the direction normal to the boundary.The location of the collided particle is also updated and the normal particle-to-wall-distance is set to be 0.25h.Eqs. 4 and 5 are discretised in time domain with the 4 th order Adams-Bashford Method (AB4).The position of each particle (eq.6) is updated with the Heun's Method.

Numerical tests
The tank used in this study is two dimensional with the base width of 3 m and the height of 3 m.It is partially filled with water up to the fill heights of H = 1 m, 1.5 m, 2 m, 2.25 m, 2.55 m, 2.7 m and 2.85 m, which correspond to the fill levels of 33 %, 50 %, 67 %, 75 %, 85 %, 90 %, and 95 %, respectively.Each particle is subjected to the constant gravity of g = -9.81m/s 2 in the y direction, and to the translational earthquake ground acceleration, a G of Loma Prieta, Whittier, Coyote Lake, Kobe, Imperial Valley and Northridge earthquakes [20], in the x direction, according to the coordinate system given in Figure 1.The time history of acceleration of each earthquake event is given in Figure 2. The mass density of water is taken as ρ 0 = 1000 kg/m 3 and the dynamical viscosity is taken as μ = 0.001 Pa•s.
The water is simulated with (N x x N y ) particles, with an equal and constant mass, where N x = 60 for all fill levels and N y = 20, 30, 40, 45, 48, 51, 54 and 57 for the aforementioned fill levels, respectively.The particles are uniformly distributed at the beginning of each test, at the distance of ∆x.The time stepping is chosen as ∆t = 1.25x10 -4 in order to avoid time-related instabilities.The free surface level location is calculated at the number of nodes uniformly distributed over x by zero-crossing of D 2 ρ.In fact, it vanishes quite above the surface level so that the surface line is determined with the help of a threshold.Therefore D 2 ρ < 350 provided fairly good results for the cases used in this paper.Sloshing is best observed by means of the temporal variation of sloshing displacements at the walls.The temporal variation of sloshing displacements, observed at the left wall, is provided in Figure 3 for the fill level of 67 %.Each case seems similar to others; the oscillating rise-and-fall behaviour is due to the fundamental unsymmetrical standing wave motion.The intensity of the average acceleration of the strong-ground motion is the primary reason for the excessive sloshing free surface displacements.The results provided by the SPH model are quite similar to those presented by the Veletsos' analytical model for rectangular tanks [21].

Ersan Güray, Gökhan Yazıcı
In Table 1, it seems that Loma Prieta has a predominant period of 1.060 s, which is very close to the 2 nd dominant standing wave period of 1.132 s for the 95 % fill level case.This nearresonant mode is directly added to the mode summation in the Veletsos's approach and then the resulting displacements come out to be quite higher.The last column of Table 1, shows the 1 st mode of the Fourier spectrum of the wave histories given in Figure 4. First four cases have similar displacements but the Northridge case results in quite higher sloshing displacements due to the fact that the strong forcing term, caused by the average acceleration, is at least 3 times higher than the other acceleration time-histories.Another important issue that needs to be addressed here is the amplification of the sloshing displacements when the predominant period of the strong ground motion, and the natural vibration period of the contained liquid, are close to each other.The natural period of the contained liquid (T nat ) with the geometrical parameters (L = 3 m and H = 2 m) is approximately 2 seconds [22].Even though the average acceleration of the Imperial Valley excitation has a smaller PGA compared to the Northridge excitation, the peak  Sloshing in medium sized tanks under earthquake load value of the free surface sloshing displacements is quite high since its mean period is closest to the natural period of the tank.A displacement spectrum analysis of strong ground motion also gives an idea on the intensity of sloshing displacements in liquid storage tanks.Figure 4 shows the 0.1 % damped displacement spectrum of each strong ground motion record.Free surface displacements and response displacements behave in a similar manner.Around T=2 s Northridge and Imperial Valley cases refer to higher response displacements when compared to other earthquake acceleration records.Thus free surface displacement amplitudes of dominant waves on the left wall are all higher, as shown in Table 2.Even though the Loma Prieta case has a relatively small spectral response displacement value, its predominant period, given in Table 1, is closer to natural frequencies for all fill levels.The force exerted over the walls is another essential issue in case of tank failure due to sloshing.It is calculated by integration of hydrodynamic pressure over the top wall boundary.The pressure is obtained at uniformly distributed points along the wall using the Kernel approximation: Lateral forces provided by the SPH model are also generally in good agreement with analytical results (Figure 5).In this figure, the upper half of each diagram gives maximum forces; the lower half negative maximum or minimum forces, and the middle region, where the force is close to zero, gives the average value of the lateral force.The lateral pushing effect of sloshing increases according to the fill level.In the cases with high fill levels, sloshed liquid contacts with the roof, leading to a higher lateral push.The fill level is a critical parameter and so the sloshing may result in roof failures due to an intensive sloshing pressure, as pointed out in the introduction.Sloshing waves can reach the top wall Ersan Güray, Gökhan Yazıcı Sloshing in medium sized tanks under earthquake load at any fill level.Maximum impact forces on the roof are given in Table 3 for all fill levels.There is no impact for shallow fill levels since the sloshing wave displacements are under the tank wall height.However, with an increase in fill levels, these waves can contact the top wall.The impact force exerted on the roof clearly depends on the fill level and the spectral displacement of ground excitation at the natural vibration period of the tank.Even at a shallow fill level of 50%, a small impact appears in the Northridge case and, at the fill level of 95 %, the magnitude of the impact force almost reaches the weight of the liquid in the tank.In addition to exerting an impact force over walls, the sloshing also create suction zones.In the case of Northridge earthquake, namely in the case of 95 % fill level case, this effect is quite excessive on the top wall.The uppermost impact occurs suddenly at t = 2.7 s (Figure 6).

Conclusion
The SPH method shows every sign of being a very promising approach to the analysis of sloshing displacements and hydrodynamic pressures in liquid storage tanks.In addition to providing the free surface line, this mesh-free method also allows observation of secondary wave formations, which is quite useful in case of fluid dynamics problems with free surface displacements.Simulations of all these 42 cases were performed on a personal computer equipped with INTEL i3 processor in approximately 4 days.The sloshing wave amplitudes, surface displacements, and lateral force values, as obtained with the Veletsos' analytical model, were found to be in good agreement with the results obtained using the SPH model presented in this study.It was observed that, in addition to peak ground acceleration, the sloshing amplitude is also closely linked to the response spectrum displacement value of the earthquake ground motion.Large sloshing displacements and hydrodynamic pressures are particularly pronounced in simulations using the Northridge earthquake record with the greatest spectral displacement of around T = 2 s, since the tank models used in this study have an approximate natural frequency of approximately 2 s.Another important observation obtained from the simulations was that the tank roofs can be exposed to very large impact forces, sometimes exceeding the weight of the contained liquid, which should be considered carefully by structural engineers.The SPH model presented in this study only considers two dimensional rectangular tank geometries.Therefore, the effects due to interaction with side walls were not investigated.It can also be expanded to handle three-dimensional tank geometries although that may cause a substantial increase in the use of computational resources required for simulations.It should also be pointed out that the SPH formulation presented in this study is derived for tanks with rigid walls.The rigid wall assumption is widely used in the design of rectangular storage tanks, especially the ones with reinforced concrete walls.However, it should be noted that the wall flexibility can alter the distribution of impulsive pressures acting on side walls for tanks with very thin walls [23,24].

Figure 1 .
Figure 1.Rectangular tank geometry and particles, r ji =r j -r i is the position vector from point i to j.

Figure 4 .
Figure 4. Displacement spectra of earthquake acceleration records used in the study.

Figure 5
Figure 5 Comparison of lateral forces exerted on tank due to sloshing liquid