Influence of mesh size on accuracy results in seismically excited water tanks

<!-- metadata commented in wiki content


<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
'''Influence of mesh size on accuracy results in seismically excited water tanks'''</div>

<div id="_GoBack" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Alejandro Hernández–Martínez<sup> a</sup>, Hugo Hernández–Barrios<sup> b</sup>, Adrián David García–Soto<sup> a</sup> and Jesús G. Valdés–Vázquez<sup>a</sup></div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<sup>a </sup>Departamento de Ingeniería Civil, División de Ingenierías, Campus Guanajuato, Universidad de Guanajuato. Av. Juárez 77, Zona Centro, Guanajuato, Gto., México.</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
alejandro.[mailto:hernandez@ugto.mx hernandez@ugto.mx], [mailto:adgarcia@ugto.mx adgarcia@ugto.mx]</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<sup>b </sup>Facultad de Ingeniería Civil, Universidad Michoacana de San Nicolás de Hidalgo</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Edificio de posgrado de ingeniería civil, Ciudad Universitaria, Morelia, Mich., México.</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
[mailto:hugohernandezbarrios@yahoo.com.mx hugohernandezbarrios@yahoo.com.mx]</div>
-->

==Abstract==

The aim of this article is to present the influence of mesh spacing in the error parameters computed for rectangular tanks seismically excited. The potential flow theory is considered using the finite difference method to set a linear numerical formulation due to its simplicity and relatively low computation time, compared with those methods that try to solve directly the Navier–Stokes equations. Wave height on the free surface, as well as pressure, shear force and overturning moment at the tanks wall base are used as parameters for the error evaluation.

'''Keywords''': tank, sloshing, error, mesh spacing

==Introduction==

<span id='OLE_LINK2'></span><span id='OLE_LINK3'></span>Damage in water storage tanks caused by earthquakes can be a serious problem, especially for the ones that are used to supply drinking water to people. Design standards and guidelines such as [1–4], employ expressions based on models proposed by Housner <span id='cite-_Ref472620840'></span>[[#_Ref472620840|[5]]] or Veletsos & Yang <span id='cite-_Ref473453008'></span>[[#_Ref473453008|[6]]] to estimate the response of liquid sloshing and to design tanks. These models are based on simple assumptions and it is not uncommon that they lead to unrepresentative results of real phenomena.

In recent years, a wide variety of computational fluid dynamics (CFD) techniques have been proposed to solve the sloshing on tanks. Several approaches such as, finite difference method (FDM), finite volume method (FVM), finite difference method (FEM), and boundary element method (BEM), among many variations and modifications have been used to address the problem. All works mainly use two types of approaches to represent the fluid flow in the tank: (a) methods that try to solve directly the Navier–Stokes equations, and (b) methods that employ a potential flow theory.

Navier–Stokes is second order, nonlinear and unstable differential equations system; this is the main reason why the methods based on this approach are complex and tend to require a lot of computational effort. FEM and/or BEM are techniques commonly used for this approach; also, hybrid methodologies have been proposed, such as Lin ''et al''. <span id='cite-_Ref473454009'></span>[[#_Ref473454009|[7]]]. On the other hand, the use of potential flow theory facilitates the establishment of a solution to the fluid motion in the tank. To use this approach, various simplifying assumptions to the Navier–Stokes equations are required, such as: the fluid is homogeneous, incompressible, inviscid and irrotational. Thus, the use of potential flow theory is usually less complicated and results can be obtained in shorter time.

One of the CFD methods used to propose a solution for sloshing in tanks is the FDM, due to its simplicity and relatively easy way to create the equations systems to be solved. Several solutions have been proposed using the FDM, such as those established by Chen ''et al''. <span id='cite-_Ref473454402'></span>[[#_Ref473454402|[8]]] and Hernández–Barrios ''et al''. [9, 10]. However, a systematic and extensive study to investigate the most advantageous mesh spacing using this method with adequate precision is lacking. Consequently, the main objective of this paper is to determinate the required mesh spacing to expect certain level of error in the parameters of interest using the FDM to solve rectangular tanks subjected to earthquake ground motions.

==Numerical Model==

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_Hernandez-Martinez_571879942-image1.jpeg|282px]] </div>

<div id="_Ref473892242" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 1 2D physical domain tank</div>

The equations of motion for the liquid–tank system depicted in <span id='cite-_Ref473892242'></span>[[#_Ref473892242|Figure 1]], which uses an inertial reference system (that moves with the tank) and considering a velocity potential ''&#x03a6;'', are described below. Domain–wide continuity of the liquid is met using the Laplace equation defined by Eq. <span id='cite-ZEqnNum755399'></span>[[#ZEqnNum755399|(1)]].

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum755399'></span>  <math>\frac{{\partial }^2\phi }{\partial x^2}+\frac{{\partial }^2\phi }{\partial y^2}=</math><math>0</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (1)
|}


On the walls, Eq. <span id='cite-ZEqnNum897311'></span>[[#ZEqnNum897311|(2)]] below must be satisfied.

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum897311'></span>  <math>\frac{\partial \phi }{\partial x}=0</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (2)
|}


And, on the bottom of the tanks the Eq. <span id='cite-ZEqnNum964230'></span>[[#ZEqnNum964230|(3)]] must be met.

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum964230'></span>  <math>\frac{\partial \phi }{\partial y}=0</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (3)
|}


On the surface ( <math display="inline">y=\delta </math> ) the movement of the fluid is governed by dynamic and kinematic conditions defined by Eqs. <span id='cite-ZEqnNum541660'></span>[[#ZEqnNum541660|(4)]] and <span id='cite-ZEqnNum503657'></span>[[#ZEqnNum503657|(5)]], respectively.

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum541660'></span>  <math>\frac{\partial \phi }{\partial t}+\frac{1}{2}\left[{\left(\frac{\partial \phi }{\partial x}\right)}^2+\right. </math><math>\left. {\left(\frac{\partial \phi }{\partial y}\right)}^2\right]+</math><math>\left\{\begin{array}{cc}
x & y
\end{array}\right\}\left\{\begin{array}{c}
a_{cx}\left(t\right)\\
a_{cy}\left(t\right)
\end{array}\right\}+g\delta =0</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (4)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum503657'></span>  <math>\frac{\partial \delta }{\partial t}+\frac{\partial \phi }{\partial x}\frac{\partial \delta }{\partial x}-</math><math>\frac{\partial \phi }{\partial y}=0</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (5)
|}


Pressure is estimated from Bernoulli equation, expressed by Eq. <span id='cite-ZEqnNum380733'></span>[[#ZEqnNum380733|(6)]], where ''&#x03c1;'' represents the fluid density;  <math display="inline">a_{cx}(t)</math> and  <math display="inline">a_{cy}(t)</math> represent respectively, the horizontal and vertical acceleration induced at tank base.

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum380733'></span>  <math>p=-\rho \left[\frac{\partial \phi }{\partial t}+\left\{\begin{array}{cc}
x & y
\end{array}\right\}\left\{\begin{array}{c}
a_{cx}\left(t\right)\\
a_{cy}\left(t\right)
\end{array}\right\}+\frac{1}{2}\left[{\left(\frac{\partial \phi }{\partial x}\right)}^2+\right. \right. </math><math>\left. \left. {\left(\frac{\partial \phi }{\partial y}\right)}^2\right]+\right. </math><math>\left. gy\right]</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (6)
|}


Discretizing in terms of second–order finite differences, and in correspondence with the scheme shown in <span id='cite-_Ref473892222'></span>[[#_Ref473892222|Figure 2]], on the walls the Eq. <span id='cite-ZEqnNum897311'></span>[[#ZEqnNum897311|(2)]] is expressed according to Eqs. <span id='cite-ZEqnNum176844'></span>[[#ZEqnNum176844|(7)]] and <span id='cite-ZEqnNum688163'></span>[[#ZEqnNum688163|(8)]].

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum176844'></span>  <math>-3{\phi }_{i,j}+4{\phi }_{i+1,j}-{\phi }_{i+2,j}=0\mbox{ }\mbox{ }i=</math><math>1,\mbox{ }j\in \left[1,j_{max}\right]</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (7)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum688163'></span>  <math>{\phi }_{i-2,j}-4{\phi }_{i-1,j}+3{\phi }_{i,j}=0\mbox{ }\mbox{ }i=</math><math>i_{max},\mbox{ }j\in \left[1,j_{max}\right]</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (8)
|}


<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_Hernandez-Martinez_571879942-image13.jpeg|306px]] </div>

<div id="_Ref473892222" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 2 Computational domain</div>

On the bottom of the tank, the condition established by Eq. <span id='cite-ZEqnNum964230'></span>[[#ZEqnNum964230|(3)]] is expressed according Eq. <span id='cite-ZEqnNum248952'></span>[[#ZEqnNum248952|(9)]].

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum248952'></span>  <math>-3{\phi }_{i,j}+4{\phi }_{i,j+1}-{\phi }_{i,j+2}=0\mbox{ }\mbox{ }i\in \left[2,i_{max-1}\right],\mbox{ }j=</math><math>1</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (9)
|}


For interior mesh points, the Laplace equation defined by Eq. <span id='cite-ZEqnNum755399'></span>[[#ZEqnNum755399|(1)]] in its discretized form is expressed by Eq. <span id='cite-ZEqnNum718595'></span>[[#ZEqnNum718595|(10)]].

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum718595'></span>  <math>C_1{\phi }_{i,j-1}+{\phi }_{i-1,j}-2\left(1+C_1\right){\phi }_{i,j}+</math><math>{\phi }_{i+1,j}+C_1{\phi }_{i,j+1}=0\mbox{ }\mbox{ }i\in \left[2,i_{max-1}\right],\mbox{ }j\in \left[2,i_{max-1}\right]</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (10)
|}


where  <math display="inline">C_1={\left(\Delta x/\Delta y\right)}^2</math>

On the surface, combining Eqs. <span id='cite-ZEqnNum541660'></span>[[#ZEqnNum541660|(4)]] and <span id='cite-ZEqnNum503657'></span>[[#ZEqnNum503657|(5)]], neglecting the higher order terms and using an implicit Crank–Nicholson scheme in time, the discretized expressions is shown in Eq. <span id='cite-ZEqnNum758334'></span>[[#ZEqnNum758334|(11)]]. Semi–viscous ''&#x03bc;'' parameter is stated in 0.5% to simulate water viscosity, as recommended by Chen ''et al.'' <span id='cite-_Ref473544347'></span>[[#_Ref473544347|[8]]].

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum758334'></span>  <math>\mbox{ }\mbox{ }C_3{\phi }_{i,j-2}^{\left(n+1\right)}-</math><math>4C_3{\phi }_{i,j-1}^{\left(n+1\right)}+\left(1+3C_3+\right. </math><math>\left. \frac{1}{2}\mu \mbox{ }\Delta t\right){\phi }_{i,j}^{\left(n+1\right)}=</math><math>-C_3{\phi }_{i,j-2}^{\left(n\right)}+\mbox{ }4C_3{\phi }_{i,j-1}^{\left(n\right)}\mbox{ }+</math><math>\mbox{ }\left(1-3C_3-\frac{1}{2}\mu \mbox{ }\Delta t\right){\phi }_{i,j}^{\left(n\right)}\mbox{ }\mbox{ }-</math><math>\frac{1}{2}\Delta t\left(a_{cx}^{\left(n+1\right)}+\right. </math><math>\left. a_{cx}^{\left(n\right)}\right)-\Delta t\mbox{ }g\mbox{ }{\delta }_{i,j}^{\left(n\right)}i\in \left[2,i_{max-1}\right],\mbox{ }j=</math><math>j_{max}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (11)
|}


where  <math display="inline">C_3=\frac{g{\left(\Delta t\right)}^2}{8\mbox{ }\Delta y}</math>

From the kinematic condition expressed by Eq. <span id='cite-ZEqnNum503657'></span>[[#ZEqnNum503657|(5)]], neglecting the higher order terms and discretizing in finite difference terms, the wave height can be obtained from Eq. <span id='cite-ZEqnNum179600'></span>[[#ZEqnNum179600|(12)]].

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum179600'></span>  <math>\begin{array}{c}
{\delta }_{i,j}^{\left(n+1\right)}={\delta }_{i,j}^{\left(n\right)}+\frac{\Delta t}{2}\left(\frac{{\phi }_{i,j-2}^{\left(n\right)}-4{\phi }_{i,j-1}^{\left(n\right)}+3{\phi }_{i,j}^{\left(n\right)}}{2\Delta y}\right)+\frac{\Delta t}{2}\left(\frac{{\phi }_{i,j-2}^{\left(n+1\right)}-4{\phi }_{i,j-1}^{\left(n+1\right)}+3{\phi }_{i,j}^{\left(n+1\right)}}{2\Delta y}\right)\\
i\in \left[1,i_{max}\right],\mbox{ }j=j_{max}
\end{array}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (12)
|}


From Eq. <span id='cite-ZEqnNum380733'></span>[[#ZEqnNum380733|(6)]], pressure can be estimated by Eq. <span id='cite-ZEqnNum608540'></span>[[#ZEqnNum608540|(13)]].

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <span id='ZEqnNum608540'></span>  <math>p_{i,j}^{\left(n+1\right)}=-\rho \left[\frac{{\phi }_{i,j}^{\left(n+1\right)}-{\phi }_{i,j}^{\left(n\right)}}{\Delta t}+\right. </math><math>\left. \frac{1}{2}\left\{\begin{array}{cc}
x & y
\end{array}\right\}\left\{\begin{array}{c}
a_{cx}^{\left(n+1\right)}+a_{cx}^{\left(n\right)}\\
a_{cy}^{\left(n+1\right)}+a_{cy}^{\left(n\right)}
\end{array}\right\}+gy\right]\mbox{ }\mbox{ }i\in \left[1,i_{max}\right],\mbox{ }j\in \left[1,j_{max}\right]</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (13)
|}


<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_Hernandez-Martinez_571879942-image21-c.png|504px]] </div>

<div id="_Ref473892490" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 3 Volcán de Cerro Prieto accelerogram (February 7, 1987) – Mexicali BC, Mexico</div>

To assess the numerical model accuracy previously described, two very different seismic records are selected. The first one is the record Volcán de Cerro Prieto (VCP) showed in <span id='cite-_Ref473892490'></span>[[#_Ref473892490|Figure 3]], which is a seismic accelerogram with a great peak ground motion, close to 1.5 times gravity in its horizontal component. This record is selected to evaluate the response to a strong impulsive seismic acceleration. Additionally, this record has a strong vertical peak ground motion near to 0.7 times gravity. The second one is the record of Central de Abastos (CDA) which is shown in <span id='cite-_Ref473892509'></span>[[#_Ref473892509|Figure 4]]. This record is selected because tends to cause high sloshing on tanks because it has a long dominant period, characteristic of soft soils.

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_Hernandez-Martinez_571879942-image22-c.png|510px]] </div>

<div id="_Ref473892509" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 4 Central de Abastos accelerogram (September 19, 1985) – Mexico City</div>

==Results==

The numerical model previously described is used to evaluate wave height on the free surface, pressure, shear force and overturning moment at tank wall base. Wave height and pressure are obtained directly from Eqs. <span id='cite-ZEqnNum179600'></span>[[#ZEqnNum179600|(12)]] and <span id='cite-ZEqnNum608540'></span>[[#ZEqnNum608540|(13)]] respectively, while shear force and overturning moment are obtained using integration Simpson’s rule at the mesh point on the walls.

Initially, the maximum values of aforementioned parameters are computed with a gross mesh. Subsequently the number on mesh points in ''y'' direction is incremented by one. Seeking a regular mesh in ''x'' direction, the number of mesh divisions is determined by searching that ''C<sub>1</sub>'' in Eq. <span id='cite-ZEqnNum718595'></span>[[#ZEqnNum718595|(10)]] is the closest to the unity, being always greater than one. The maximum values of each parameter are compared to those obtained from the previous mesh, where the estimated errors are the difference. <span id='cite-_Ref473893901'></span>[[#_Ref473893901|Figure 5]] shows the error evolution for a 100 m³ tank (2''a'' = 6.25m,  <math display="inline">H=</math> 3.1 m) subjected to VCP accelerogram. Similar analyses are performed for tank capacities from 100 m³ to 1500 m³ with typical dimensions showed in <span id='cite-_Ref473893250'></span>[[#_Ref473893250|Table 1]].

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_Hernandez-Martinez_571879942-image24-c.jpeg|492px]] </div>

<div id="_Ref473893901" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 5 Error evolution for a 100 m³ tank</div>

Many previous works report results just considering horizontal acceleration. Therefore, two sets of analyses for each accelerogram are computed, one set just with horizontal acceleration and another one that considers horizontal and vertical acceleration at the same time. <span id='cite-_Ref473893250'></span>[[#_Ref473893250|Table 1]] shows the maximum values obtained for VCP accelerogram for both sets of analyses, where ''a<sub>cx</sub>'' implies that results are obtained with horizontal acceleration only, and ''a<sub>cx</sub> ''+ ''a<sub>cy</sub>'' means that horizontal and vertical acceleration are considered. The presented values in <span id='cite-_Ref473893250'></span>[[#_Ref473893250|Table 1]] correspond to the densest mesh computed for each case because they are considered the most accurate results. As can be noted, wave heights are practically the same for both sets of analysis. However, for those results computed with horizontal and vertical acceleration, pressure, shear force, and overturning moment at the wall base are on average 17.8%, 14.4% and 3.3% greater respect to those obtained with horizontal acceleration only. On the other hand, if the same procedure is computed for CDA accelerogram, the obtained results are practically the same with or without vertical acceleration, having differences of 1.3% for pressure at most. This occurs because pressure, shear, and moment on VCP record are dominated by the impulsive induced motion, while for CDA accelerogram the motion is dominated by the sloshing effect, where vertical acceleration has a little influence. Additionally, vertical acceleration is considerably small for CDA accelerogram as can be shown in <span id='cite-_Ref473892509'></span>[[#_Ref473892509|Figure 4]].

Note from <span id='cite-_Ref473893901'></span>[[#_Ref473893901|Figure 5]] that for the same grid spacing, several error values are obtained for each parameter. Considering that error evolution showed in <span id='cite-_Ref473893901'></span>[[#_Ref473893901|Figure 5]] is different for each tank geometry, a set of grid spacing values is obtained to achieve the same error. <span id='cite-_Ref474411948'></span>[[#_Ref474411948|Figure 6]] shows the distribution of required grid spacing to achieve several error levels for all tank dimensions and VPC accelerogram. Thin line joins the minimum and maximum grid spacing values, central mark indicates the average value, while the thick line joins the ± standard deviation around average. On <span id='cite-_Ref474411948'></span>[[#_Ref474411948|Figure 6]] only wave height and pressure are considered, because these parameters are obtained directly from the numerical model. Shear force and overturning moment are obtained from numerical integration, adding a numerical error, which leads to a smaller mesh spacing to achieve the same error. This can be observed in <span id='cite-_Ref474411970'></span>[[#_Ref474411970|Figure 7]], where required mesh spacing for all parameters is showed.

<div id="_Ref473893250" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Table 1 Result values for Volcán de Cerro Prieto accelerogram</div>

{| style="width: 100%;margin: 1em auto 0.1em auto;border-collapse: collapse;" 
|-
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|Capacity

(m³)
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|Depth
{| class='formulaSCP' style='width: 100%;'
|-
|
{| style='margin:auto;width: 100%; text-align:center;'
|-
| <math display="inline">H</math>
|}
| style='width: 5px;text-align: right;white-space: nowrap;' | (m) 
|}
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|Width
{| class='formulaSCP' style='width: 100%;'
|-
|
{| style='margin:auto;width: 100%; text-align:center;'
|-
| <math display="inline">2a</math>
|}
| style='width: 5px;text-align: right;white-space: nowrap;' | (m) 
|}
|  colspan='2'  style="border-top: 1pt solid black;text-align: center;"|Wave Height

(m)
|  colspan='2'  style="border-top: 1pt solid black;text-align: center;"|Pressure

(kPa)
|  colspan='2'  style="border-top: 1pt solid black;text-align: center;"|Shear

(kN)
|  colspan='2'  style="border-top: 1pt solid black;text-align: center;"|Moment

(kN–m)
|-
|  style="border-bottom: 1pt solid black;text-align: center;width: 7%;"|<math display="inline">a_{cx}</math> 
|  style="border-bottom: 1pt solid black;text-align: center;width: 11%;"|<math display="inline">a_{cx}+a_{cy}</math> 
|  style="border-bottom: 1pt solid black;text-align: center;width: 6%;"|<math display="inline">a_{cx}</math> 
|  style="border-bottom: 1pt solid black;text-align: center;width: 10%;"|<math display="inline">a_{cx}+a_{cy}</math> 
|  style="border-bottom: 1pt solid black;text-align: center;width: 7%;"|<math display="inline">a_{cx}</math> 
|  style="border-bottom: 1pt solid black;text-align: center;width: 11%;"|<math display="inline">a_{cx}+a_{cy}</math> 
|  style="border-bottom: 1pt solid black;text-align: center;width: 7%;"|<math display="inline">a_{cx}</math> 
|  style="border-bottom: 1pt solid black;text-align: center;width: 11%;"|<math display="inline">a_{cx}+a_{cy}</math> 
|-
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|100
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|3.1
|  style="border-top: 1pt solid black;text-align: center;"|5.25
|  style="border-top: 1pt solid black;text-align: center;"|0.370
|  style="border-top: 1pt solid black;text-align: center;"|0.370
|  style="border-top: 1pt solid black;text-align: center;"|56.59
|  style="border-top: 1pt solid black;text-align: center;"|67.44
|  style="border-top: 1pt solid black;text-align: center;"|106.55
|  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|122.69
|  style="border-top: 1pt solid black;text-align: center;"|186.94
|  style="border-top: 1pt solid black;text-align: center;"|193.27
|-
|  style="border-bottom: 1pt solid black;text-align: center;"|6.25
|  style="border-bottom: 1pt solid black;text-align: center;"|0.400
|  style="border-bottom: 1pt solid black;text-align: center;"|0.400
|  style="border-bottom: 1pt solid black;text-align: center;"|58.35
|  style="border-bottom: 1pt solid black;text-align: center;"|69.13
|  style="border-bottom: 1pt solid black;text-align: center;"|110.27
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|126.41
|  style="border-bottom: 1pt solid black;text-align: center;"|191.74
|  style="border-bottom: 1pt solid black;text-align: center;"|198.07
|-
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|300
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|3.6
|  style="border-top: 1pt solid black;text-align: center;"|8.25
|  style="border-top: 1pt solid black;text-align: center;"|0.496
|  style="border-top: 1pt solid black;text-align: center;"|0.496
|  style="border-top: 1pt solid black;text-align: center;"|69.05
|  style="border-top: 1pt solid black;text-align: center;"|81.53
|  style="border-top: 1pt solid black;text-align: center;"|152.26
|  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|174.02
|  style="border-top: 1pt solid black;text-align: center;"|301.16
|  style="border-top: 1pt solid black;text-align: center;"|311.07
|-
|  style="border-bottom: 1pt solid black;text-align: center;"|10.25
|  style="border-bottom: 1pt solid black;text-align: center;"|0.495
|  style="border-bottom: 1pt solid black;text-align: center;"|0.495
|  style="border-bottom: 1pt solid black;text-align: center;"|70.40
|  style="border-bottom: 1pt solid black;text-align: center;"|82.83
|  style="border-bottom: 1pt solid black;text-align: center;"|155.55
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|177.31
|  style="border-bottom: 1pt solid black;text-align: center;"|305.00
|  style="border-bottom: 1pt solid black;text-align: center;"|314.91
|-
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|500
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|3.8
|  style="border-top: 1pt solid black;text-align: center;"|9.25
|  style="border-top: 1pt solid black;text-align: center;"|0.545
|  style="border-top: 1pt solid black;text-align: center;"|0.545
|  style="border-top: 1pt solid black;text-align: center;"|73.39
|  style="border-top: 1pt solid black;text-align: center;"|86.55
|  style="border-top: 1pt solid black;text-align: center;"|171.05
|  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|195.30
|  style="border-top: 1pt solid black;text-align: center;"|354.95
|  style="border-top: 1pt solid black;text-align: center;"|366.61
|-
|  style="border-bottom: 1pt solid black;text-align: center;"|14.25
|  style="border-bottom: 1pt solid black;text-align: center;"|0.440
|  style="border-bottom: 1pt solid black;text-align: center;"|0.440
|  style="border-bottom: 1pt solid black;text-align: center;"|75.16
|  style="border-bottom: 1pt solid black;text-align: center;"|88.25
|  style="border-bottom: 1pt solid black;text-align: center;"|175.56
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|199.80
|  style="border-bottom: 1pt solid black;text-align: center;"|360.25
|  style="border-bottom: 1pt solid black;text-align: center;"|371.90
|-
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|800
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|4.35
|  style="border-top: 1pt solid black;text-align: center;"|11.30
|  style="border-top: 1pt solid black;text-align: center;"|0.483
|  style="border-top: 1pt solid black;text-align: center;"|0.407
|  style="border-top: 1pt solid black;text-align: center;"|84.55
|  style="border-top: 1pt solid black;text-align: center;"|99.60
|  style="border-top: 1pt solid black;text-align: center;"|226.04
|  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|257.82
|  style="border-top: 1pt solid black;text-align: center;"|532.22
|  style="border-top: 1pt solid black;text-align: center;"|549.70
|-
|  style="border-bottom: 1pt solid black;text-align: center;"|16.30
|  style="border-bottom: 1pt solid black;text-align: center;"|0.377
|  style="border-bottom: 1pt solid black;text-align: center;"|0.341
|  style="border-bottom: 1pt solid black;text-align: center;"|86.03
|  style="border-bottom: 1pt solid black;text-align: center;"|101.02
|  style="border-bottom: 1pt solid black;text-align: center;"|230.29
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|262.06
|  style="border-bottom: 1pt solid black;text-align: center;"|538.80
|  style="border-bottom: 1pt solid black;text-align: center;"|556.29
|-
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|1000
|  rowspan='2' style="border-top: 1pt solid black;text-align: center;"|4.5
|  style="border-top: 1pt solid black;text-align: center;"|12.30
|  style="border-top: 1pt solid black;text-align: center;"|0.370
|  style="border-top: 1pt solid black;text-align: center;"|0.370
|  style="border-top: 1pt solid black;text-align: center;"|87.81
|  style="border-top: 1pt solid black;text-align: center;"|103.37
|  style="border-top: 1pt solid black;text-align: center;"|243.00
|  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|277.00
|  style="border-top: 1pt solid black;text-align: center;"|590.41
|  style="border-top: 1pt solid black;text-align: center;"|609.77
|-
|  style="border-bottom: 1pt solid black;text-align: center;"|18.30
|  style="border-bottom: 1pt solid black;text-align: center;"|0.332
|  style="border-bottom: 1pt solid black;text-align: center;"|0.332
|  style="border-bottom: 1pt solid black;text-align: center;"|89.12
|  style="border-bottom: 1pt solid black;text-align: center;"|104.62
|  style="border-bottom: 1pt solid black;text-align: center;"|246.86
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|280.86
|  style="border-bottom: 1pt solid black;text-align: center;"|596.63
|  style="border-bottom: 1pt solid black;text-align: center;"|615.98
|-
|  rowspan='2' style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|1500
|  rowspan='2' style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|4.6
|  style="border-top: 1pt solid black;text-align: center;"|16.30
|  style="border-top: 1pt solid black;text-align: center;"|0.373
|  style="border-top: 1pt solid black;text-align: center;"|0.373
|  style="border-top: 1pt solid black;text-align: center;"|90.85
|  style="border-top: 1pt solid black;text-align: center;"|106.71
|  style="border-top: 1pt solid black;text-align: center;"|257.25
|  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|292.78
|  style="border-top: 1pt solid black;text-align: center;"|635.93
|  style="border-top: 1pt solid black;text-align: center;"|656.60
|-
|  style="border-bottom: 1pt solid black;text-align: center;"|20.30
|  style="border-bottom: 1pt solid black;text-align: center;"|0.305
|  style="border-bottom: 1pt solid black;text-align: center;"|0.305
|  style="border-bottom: 1pt solid black;text-align: center;"|91.17
|  style="border-bottom: 1pt solid black;text-align: center;"|107.03
|  style="border-bottom: 1pt solid black;text-align: center;"|258.22
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|293.75
|  style="border-bottom: 1pt solid black;text-align: center;"|637.39
|  style="border-bottom: 1pt solid black;text-align: center;"|658.07
|}


<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_Hernandez-Martinez_571879942-image29-c.jpeg|450px]] </div>

<div id="_Ref474411948" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 6 Grid spacing distribution for wave height and pressure for VCP accelerogram</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_Hernandez-Martinez_571879942-image30-c.jpeg|438px]] </div>

<div id="_Ref474411970" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 7 Grid spacing distribution for all parameters and VCP accelerogram</div>

<span id='cite-_Ref474411987'></span>[[#_Ref474411987|Table 2]] shows the recommended values of mesh separation &#x0394;''y'' to get certain error level for both records considering horizontal and vertical acceleration at the same time. Suggestions are based on the average minus the standard deviation values showed in Figures 6 and 7 for VCP record. A same criterion is used for CDA record. As it can be observed from <span id='cite-_Ref474411987'></span>[[#_Ref474411987|Table 2]], for small errors, similar recommended values are obtained for both records, while mesh spacing is a little more conservative for VPC accelerogram for larger error levels. Thus, values shown in <span id='cite-_Ref474411987'></span>[[#_Ref474411987|Table 2]] can be used as reference to create the mesh for the desired error level when perform the described numerical model.

<div id="_Ref474411987" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Table 2 Recommended mesh spacing to achieve the expected error</div>

{| style="width: 74%;margin: 1em auto 0.1em auto;border-collapse: collapse;" 
|-
|  rowspan='3' style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|Error
|  colspan='4'  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|Recommended mesh spacing (m)
|-
|  colspan='2'  style="border-top: 1pt solid black;border-bottom: 1pt solid black;border-right: 1pt solid black;text-align: center;"|VCP record (impulsive motion)
|  colspan='2'  style="border-top: 1pt solid black;border-left: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|CDA record (sloshing effect)
|-
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|Wave height & pressure
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;border-right: 1pt solid black;text-align: center;"|All

parameters
|  style="border-top: 1pt solid black;border-left: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|Wave height &

pressure
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|All parameters
|-
|  style="border-top: 1pt solid black;text-align: center;"|10–1
|  style="border-top: 1pt solid black;text-align: center;"|0.45
|  style="border-top: 1pt solid black;text-align: center;"|0.13
|  style="border-top: 1pt solid black;text-align: center;"|0.50
|  style="border-top: 1pt solid black;text-align: center;"|0.20
|-
|  style="text-align: center;"|10–2
|  style="text-align: center;"|0.20
|  style="text-align: center;"|0.06
|  style="text-align: center;"|0.25
|  style="text-align: center;"|0.08
|-
|  style="text-align: center;"|10–3
|  style="text-align: center;"|0.08
|  style="text-align: center;"|0.04
|  style="text-align: center;"|0.08
|  style="text-align: center;"|0.03
|-
|  style="text-align: center;"|10–4
|  style="text-align: center;"|0.03
|  style="text-align: center;"|0.02
|  style="text-align: center;"|0.03
|  style="text-align: center;"|0.02
|-
|  style="border-bottom: 1pt solid black;text-align: center;"|10–5
|  style="border-bottom: 1pt solid black;text-align: center;"|0.02
|  style="border-bottom: 1pt solid black;text-align: center;"|0.01
|  style="border-bottom: 1pt solid black;text-align: center;"|0.02
|  style="border-bottom: 1pt solid black;text-align: center;"|0.01
|}


Since reasonably similar recommended values are obtained even for two records that induce very different behavior in the liquid–tank system, it is expected that comparable results can be obtained for accelerograms which leads to an intermediate behavior respect to those here presented.

==Conclusions==

The suggested values of mesh separation &#x0394;''y'' presented in Table 2 can be used as a guide to perform analysis of seismically excited water tanks. This can be done even of thanks are not rectangular, because normally an irregular mesh is transformed into a rectangular system to apply the FDM.

The forces magnitude on the walls appears to be directly related to the peak ground acceleration. Similarly, the influence of vertical acceleration on the results tends to be proportional to the magnitude of peak vertical acceleration, i.e.: if peak vertical acceleration is small, results tend not to be affected. Moreover, a large sloshing effect not necessarily induces large forces on tank boundaries.

Each analysis was performed in a computer with an Intel® Xeon® X5650 [mailto:@ @] 2.67 GHz processor that never required more than two hours to solve the system for the entire accelerogram and meshes with smaller spacing, which have near 400 thousand points. The numerical model was coded in Fortran 95 using LU factorization algorithm with compressed sparse row storage scheme to solve the equations system. Whereas some commercial programs that solve the Navier–Stokes equations may require one of two days of execution to solve a similar problem, one or two hours of computational time (depending on the required precision) for the numerical model may only be needed. Therefore, numerical model presented can be and excellent alternative to simplified model used in design standards and guidelines to analyze and design water storage tanks.

==Acknowledgements==

The authors acknowledge the support of the Universidad de Guanajuato and the Universidad Michoacana de San Nicolás de Hidalgo to develop the presented paper. Also, the authors thank the reviewers for their comments and suggestions, which helped to improve this paper.

==References==

<span id='_Ref472620798'></span>
:[[#cite-_Ref472620798|[1]]] American Waters Works Association (AWWA). ''Steel water storage tanks – Design, construction, maintenance and repair''. McGraw–Hill, 2010.

<span id='_Ref472620808'></span>
:[[#cite-_Ref472620808|[2]]] Priestley MJN, Davidson BJ, Honey GD, Hopkins DC, Martin RJ, Ramsey et al. editors. ''Seismic design of storage tanks – Recommendations from the New Zealand National Society of Earthquake Engineering''. Wellington, New Zealand. 1986.

<span id='_Ref472620815'></span>
:[[#cite-_Ref472620815|[3]]] Comisión Federal de Electricidad (CFE). ''Manual de Diseño de Obras Civiles – Diseño por Sismo''. Instituto de Investigaciones Eléctricas, México. 2015.

:[4] American Concrete Institute (ACI). ''Seismic Design of Liquid–Containing Concrete Structures and Commentary ACI 350.3–06'', USA, 2006.

<span id='_Ref472620840'></span>
:[[#cite-_Ref472620840|[5]]] Housner GW. The dynamic behavior of water storage tanks. ''Bulletin of Seismological Society of America'' 1963; '''53''':381–387.

<span id='_Ref473453008'></span>
:[[#cite-_Ref473453008|[6]]] Veletsos AS, Yang JY. Earthquake response of liquid storage tanks. Proceeding of the 2<sup>nd</sup> Annual Engineering Mechanics Division ASCE. 1977. 1–24

<span id='_Ref473454009'></span>
:[[#cite-_Ref473454009|[7]]] Lin G, Liu J, Li J, Hu Z. A scaled boundary finite element approach for sloshing analysis of liquid storage tanks. ''Engineering Analysis with Boundary Elements'' 2015; '''56''':70–80.

<span id='_Ref473454402'></span><span id='_Ref473544347'></span>
:[[#cite-_Ref473544347|[8]]] Chen W, Haroun MA, Liu F. Large amplitude liquid sloshing in seismically excited tanks. ''Earthquake Engineering & Structural Dynamics'' 1996. '''25''':653–669.

:[9] Hernández–Barrios H, Heredia–Zavoni E, Aldama–Rodríguez AA. Nonlinear sloshing response of cylindrical tanks subjected to earthquake ground motion. ''Engineering Structures'' 2007. '''29''':3364–3376.

:[10] Hernández–Barrios H, Hernández–Martínez A, Valdés–Vázquez JG. Non–linear sloshing effect on storage tanks subjected to high earthquake ground motion. ''Revista Internacional de Métodos Numéricos en Ingeniería'' 2015; '''31''':198–206.