- Email: [email protected]
Sloshing waves and resonance modes of fluid in a 3D tank by a time-independent finite difference method
ARTICLE IN PRESS Ocean Engineering 36 (2009) 500–510
Contents lists available at ScienceDirect
Ocean Engineering journal homepage: www.elsevier.com/locate/oceaneng
Sloshing waves and resonance modes of fluid in a 3D tank by a time-independent finite difference method Chih-Hua Wu, Bang-Fuh Chen Department of Marine Environment and Engineering, National Sun Yat-sen University, Kaohsiung 804, Taiwan
a r t i c l e in f o
a b s t r a c t
Article history: Received 25 June 2008 Accepted 31 January 2009 Available online 13 February 2009
A 3D time-independent finite difference method is developed to solve for wave sloshing in a threedimensional tank excited by coupled surge and sway motions. The 3D equations of fluid motion are derived in a moving coordinate system. The three-dimensional tank, with an arbitrary depth and a square base, is subjected to a range of excitation frequencies with motions that exhibit multiple degrees of freedom. For demonstration purposes the numerical scheme is validated by a benchmark study. Five types of sloshing waves were observed when the tank is excited by various excitation frequencies. A spectral analysis identified the resonant frequencies of each type of wave and the results show a strong correlation between resonant modes and the occurrence of the sloshing wave types. The method can be used to simulate fluid sloshing in a 3D tank with six-degrees of freedom. & 2009 Elsevier Ltd. All rights reserved.
Keywords: Fully nonlinear free surface wave Sloshing fluids 3D tank Finite difference method
1. Introduction Free surface sloshing in a moving container is associated with various engineering problems, such as tank trucks on highways, liquid oscillations in large storage tanks caused by earthquakes, sloshing of liquid cargo in ocean-going vessels and the motion of liquid fuel in aircraft and spacecraft. It is known that partially filled tanks are prone to violent sloshing under certain motions, especially when near-resonant excitation occurs. The large liquid movement creates highly localized impact pressures on tank walls, which in turn cause structural damage and may even create moments that affect the stability of the vehicle, which carries the container. Sloshing waves in moving tanks have been studied numerically, theoretically and experimentally in the past several decades and many significant phenomena have been considered in those studies, especially the linear and nonlinear effects of sloshing for both inviscid and viscous liquids. Most reported studies involved tanks excited by limited excitation directions and with a fixed excitation frequency throughout the excitation. In reality, as the tank is excited by accelerations due to an earthquake or waves, the excitation directions include multi-degrees of freedom (surge/ sway/heave/pitch/roll/yaw) and the excitation frequency also varies with time. The potential formulation of the problem is often used in studying sloshing, for example Ockendon et al. (1996) among others. The most distinguished analytical works are Faltinsen’s series of studies (Faltinsen, 1978; Faltinsen and
Corresponding author. Tel./fax: +886 7 5255065.
E-mail address: [email protected] (B.-F. Chen). 0029-8018/$ - see front matter & 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.oceaneng.2009.01.020
Timokha, 2001, 2002) for sloshing fluids in 2D tanks and Faltinsen et al. (2005) where their asymptotic modal system is extended to model nonlinear sloshing in a 3D rectangular tank. Besides the potential flow approaches, many numerical studies (computational fluid dynamics) of the problems with primitive variables have been made. These have focused particularly on the fully nonlinear free surface effects. Many papers give successful examples for two-dimensional sloshing (see, for instance, Chen and Chiang, 1999a; Celebi and Akyildiz, 2002; Turnbull et al., 2003; Aliabadi et al., 2003; Frandsen, 2004; Chen and Nokes, 2005; Wu, 2007) and for three-dimensional sloshing (Wu et al., ¨ nal, 2005, 2006; Kim et al., 2007; 1998; Kim, 2001; Akyildiz and U Lee et al., 2007a, b). The reported techniques for handling a wavyfree surface include VOF, SOLA, SURF and also the s-transformation technique to stretch the grid from the bed to surface. However, the 3D numerical simulation of fluid sloshing in a tank is still very limited. In the present study, a 3D tank with different ratios of depth/excitation amplitude, multiple degrees of freedom of excitation and excitation frequencies are considered. In this three-dimensional model, we develop a 3D time-independent finite difference method to incorporate the incompressible and inviscid Navier–Stokes equations, fully nonlinear kinematic and dynamic free surface conditions in the analysis of the seismic response of sloshing fluid in a rectangular tank with a square base. The time-varying moving boundary is mapped onto a timeindependent domain through proper transformation functions and a special finite difference approximation is made in order to overcome the difficulty of maintaining the accuracy of the finite difference expression for the second derivative when the difference mesh is stretched near the boundary. The main focus of this paper is the simulation of a 3D tank undergoing different
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
combinations of motions with varying excitation directions, and not only transient but nearly steady-state phenomena are recorded and studied. Hill (2003) compared the transient waves with the steady-state waves and showed that the maximum transient response can far exceed the steady-state response of the basin. The transient effect on moving vehicles is more dangerous than that of steady-state condition. Section 2 introduces the equations of motion that are written in a moving coordinate system attached to the accelerating tank. The proper coordinate transformation functions are adopted to map the time-dependent domain into a fixed unit cubic. The proposed finite difference method is developed in Section 3. Section 4 presents the detailed results and provides comprehensive discussion of all the phenomena found in this study.
2. Mathematical formulation
and the kinematic and dynamic free surface conditions are @h @h @h þu þw ¼v @t @x @z
2 ð2a_ b_ b_ g_ 2 Þx 2ðb_ w g_ vÞ
(1)
@v @v @v @v 1 @p € þu þv þw ¼ gy Y C ðg€ x a€ zÞ @t @x @y @z r @y ð2b_ g_ a_ 2 g_ 2 Þy 2ðg_ u a_ wÞ
(2)
p¼0
(6)
The Poisson equation can be derived by taking partial derivatives of (1), (2) and (3) with respect to x, y and z, respectively, and summing the results. @2 p @2 p @2 p @ @u @u @u u þv þw þ 2 þ 2 ¼ r 2 @x @x @y @z @x @y @z @ @v @v @v r u þv þw @y @x @y @z @ @w @w @w r u þv þw @z @x @y @z 2 2r½a_ b_ þ b_ g_ þ a_ g_ Þ ða_ 2 þ b_ þ g_ 2 Þ @v @w 2r a_ @z @y @w @u @u @v þb_ þ g_ @x @z @y @x
2
ð2a_ g_ a_ 2 b_ Þz 2ða_ v b_ uÞ
(3)
where u, v and w are velocity components in x-, y- and z-directions, X€ C , Y€ C and Z€ C are the acceleration components of the tank in x, y and z directions; (a_ , b_ , g_ ) and (a€ , b€ , g€ ) are the corresponding angular velocities and accelerations with respect to the x, y and z axes, p is the pressure, r is the fluid density, and gx,gy and gz are the components of the acceleration due to gravity. The continuity equation for incompressible flow is @u @v @w þ þ ¼0 @x @y @z
(4)
x b1 ðy; zÞ b2 ðy; zÞ b1 ðy; zÞ z b3 ðx; yÞ zn ¼ b4 ðx; yÞ b3 ðx; yÞ
yn ¼ 1
y þ d0 hðx; z; tÞ (8)
where the instantaneous water depth, h(x,z,t), is a single-valued function measured from tank bottom, d0 represents the vertical distance between the still water surface and the tank bottom, b1(y,z) and b2(y,z) are horizontal distances from the z-axis to the west (W1) and east (W2) walls, respectively, and b3(x,y) and b4(x,y) are horizontal distance from the x-axis to the south (W3) and north (W4) walls, respectively (see Fig. 1). Based on the assumption of h(x,z,t), the phenomena of wave breaking, run-up and tank roof impacts are not considered in the present study. With the above equations, one can map the west wall to x* ¼ 0
θ
Y
X
Yaw
Z Surge
X Y
D
Pitch
Z A Sway
E
h (x, z, t) b1 (y, z)
N
b4 (x, y) b2(y, z) b3 (x, y)
W3
F
X
W4
Z
W1
C
G
I
Roll
(7)
Many finite difference (FDM) and finite volume methods have been reported to solve the free surface displacement of sloshing fluid in tanks. A brief review of those methods can be found in Chen and Nokes (2005). In the present study, the irregular boundary, such as the time-varying fluid surface, non-vertical walls and non-horizontal bottom, can be mapped onto a square by the proper coordinate transformations (Hung and Wang, 1987; Chen et al., 1999b) as follows: xn ¼
@w @w @w @w 1 @p € þu þv þw ¼ gz Z C ða€ y b€ xÞ @t @x @y @z r @z
Heavy
(5)
and
The momentum equations are written based on a moving coordinate system (the coordinate system is chosen to move with the tank, including surge, sway, heave, yaw, roll and pitch motions, see Fig. 1). @u @u @u @u 1 @p € þu þv þw ¼ gx X C ðb€ z g€ yÞ @t @x @y @z r @x
501
E L
Fig. 1. Definition sketch of the tank motion.
B
W2
B
ARTICLE IN PRESS 502
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
Fig. 2. The concepts of transformation coordinate systems and stretch of the grid. The stretching constants presented in (x*, y*, z*) coordinate are l1 ¼ l2 ¼ l3 ¼ 0.5, k1 ¼ k2 ¼ k3 ¼ 2.
and east wall to x* ¼ 1, the north wall to z* ¼ 0 and south wall to z* ¼ 1, the free surface to y* ¼ 0 and the tank bottom to y* ¼ 1. In this way the computational domain is transformed to a fixed unit cubic domain. The coordinates (x*, y*, z*) can be further transformed such that the layers near the boundary are stretched to capture the sharp local velocity gradients. The following exponential functions provide these stretching transformations: n ðxn 1Þ
X ¼ l1 þ ðxn l1 Þek1 x Z ¼ l3 þ ðzn l3 Þek3 z
n ðyn 1Þ
Y ¼ l2 þ ðyn l2 Þek2 y
n ðzn 1Þ
(9)
The constants ki and li, i ¼ 1,2,3, control the mesh size and stretching in X-, Y- and Z-directions, respectively. The concepts of transformation coordinate system and grid stretching are presented in Fig. 2. The geometry of the flow field and the meshes in the computational domain (X–Y–Z system) become timeindependent throughout the computational analysis. The dimensional parameters are normalized in the following way: u v w p U ¼ pffiffiffiffiffiffiffiffi V ¼ pffiffiffiffiffiffiffiffi W ¼ pffiffiffiffiffiffiffiffi P ¼ rgd0 gd0 gd0 gd0 rffiffiffiffiffi g h xc y zc H¼ T¼t Xc ¼ Y c ¼ c Zc ¼ d0 d0 d0 d0 d0 sffiffiffiffiffi a b g a_ d0 Ya ¼ Yb ¼ Yg ¼ YaT ¼ 2p 2p 2p 2p g sffiffiffiffiffi sffiffiffiffiffi b_ d0 g_ d0 a€ d0 YbT ¼ YgT ¼ YaTT ¼ 2p g 2p g 2p g
YbTT ¼
b€ d0 2p g
YgTT ¼
g€ d0 2p g
With the aforementioned transformations and dimensionless variables, Eqs. (1)–(7) can be written in dimensionless form. For example, the momentum equation in x-direction is: U T þ ðC 10 C 13 U X þ C 11 C 14 U Y þ C 12 C 15 U Z Þ þ C 1 C 13 UU X þ C 2 C 14 UU Y þ C 3 C 15 UU Z þ C 4 C 13 VU X þ C 5 C 14 VU Y þ C 6 C 15 VU Z þ C 7 C 13 WU X þ C 8 C 14 WU Y þ C 9 C 15 WU Z ¼ Gx ðC 1 C 13 P X þ C 2 C 14 P Y þ C 3 C 15 PZ Þ 2p 4p2 ðzYbTT yYgTT Þ ð2YaT YbT Y2bT Y2gT Þx d0 d0 (10) 4pðYbT W YgT VÞ X CTT
The dimensionless forms of the other equations are omitted in the text. In Eq. (10), the coefficients C1–C15 are due to coordinate transformations and can be found in Chen et al. (1999b). PX denotes a partial derivative of P with respect to X; UT is the partial derivative of U with respect to dimensionless time T; XCTT, YCTT and ZCTT are dimensionless ground accelerations in x, y and z directions; the other terms have similar meanings.
3. Computational algorithm In this three-dimensional analysis, the fluid flow is solved in a unit cubic mesh in the transformed flow domain. All the computations use the dimensionless equations in X–Y–Z coordinate system. The difference equations for space derivatives use the central difference approximation, except at the boundaries were forward or backward differences are employed. A staggered grid system is used in the analysis. That is, the pressure P is defined at the center of a cell, whereas the velocity components U, V and W are calculated 0.5DX, 0.5DY and 0.5DZ behind, below or forward of the cell center. The Crank–Nicholson second-order finite difference scheme and the Gauss–Seidel Point successive over-relaxation iterative procedure are used to calculate the velocity and pressure, respectively. However, it is sensitive to numerical solution for the convective term and the second-order central difference upwind scheme is applied in the present numerical scheme. The numerical scheme is similar to that described in Chen and Nokes (2005) and the detailed iteration procedures are omitted here.
4. Results and discussions The natural modes of a 3D tank with square base can be obtained by solving a linearized natural sloshing standing wave problem where the circular frequency can be expressed as qffiffiffiffiffiffiffiffiffiffiffiffiffiffi o2i;j ¼ g li;j tanhðli;j d0 Þ; li;j ¼ p i2 þ j2 (11) and i, j are the natural mode’s components in the x and z directions. The FFT analysis of ChiChi earthquake is presented in Fig. 3 and a wide range of excitation frequencies between 0.4o1 and 2.5o1 is considered in this study. The accuracy of the numerical results significantly depends on the spatial grid resolution and the selected time step. Thus, it is important to verify the appropriate time step for a particular spatial grid before large numerical problems are simulated. In the present study, the time step is restricted by Dto minfDxmin =jui;j;k j; Dymin =jvi;j;k j; Dzmin =jwi;j;k jg, which implies that a fluid particle cannot move more than one cell in a single time step. The selection of the spatial grid was also based upon a parametric study. Fig. 4 shows the sloshing displacement at point A with various mesh sizes selections. As depicted in the figure, the results of mesh size (20 10 20) agree well with those of Kim (2001), whose mesh size was (30 30 30). The results of finer meshes and smaller time-step selection present minor difference from that of mesh size ¼ (20 10 20). The influence of stretching constants is illustrated in Fig. 4(c). The mass balance is shown in Fig. 4(d) and good numerical accuracy is noted. Therefore, the mesh size of (20 10 20), DT ¼ 0.006, l1 ¼ l2 ¼ l3 ¼ 0.5 and k1 ¼ k2 ¼ k3 ¼ 2 are used in this study. The simulation of the cases
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
with strong shallow-fluid sloshing is beyond the limits of the present numerical model. The proposed numerical scheme can be applied in a non-resonant oscillation with a0/d0 ¼ 2.2 and the d0/L is limited to 0.1. In a near-resonant oscillation, however, the limit of a0/d0 becomes 0.05 in the present study.
503
benchmark tests. The comparison between the present numerical results and the experimental and theoretical results reported by Faltinsen et al. (2005) are shown in Fig. 5 and the agreement is very good. Fig. 6 further compares the sloshing displacements of fluid in a tank excited by vertical motion with those reported by Frandsen (2004). The agreement is also very good.
4.1. Benchmark tests In order to validate the accuracy of numerical simulation in the present study, the results obtained from the present numerical model are compared with those reported in the literature for
Fig. 3. The spectral analysis of 921, 1999 ChiChi earthquake in Taiwan o1: The fundamental natural frequency of a tank with d0/L ¼ d0/B ¼ 0.25.
4.2. Horizontal ground motion In a real marine environment, the excitation forcing of the tank might involve multiple directions (surge, sway, heave, pitch, yaw and roll). In the present study, we focus on horizontal excitation only. The tank, in general, is not excited exactly in surge (x-) or sway (z-) directions. Rather it is excited in a coupled surge and sway motion. We define y as the excitation angle of the horizontal ground motion (see Fig. 1) that represents the various coupled surge/sway directions. For longitudinal excitation (surge motion, y ¼ 01), the work of Faltinsen et al. (2005) indicated four different kinds of waves, namely planar wave, swirling wave, irregular wave and square-like wave, could be generated through this type of forcing. As the tank excitation angle y is equal to 451 (diagonal motion), most of the four kinds of waves can be observed while the stable square-like waves are seldom seen. In the present study, horizontal excitations with various excitation angles and
Fig. 4. The wave history at corner A of the tank under diagonal motion, the ratio d0/L ¼ d0/B ¼ 0.25, the ratio of excitation displacement a0/L ¼ 0.0093, ox ¼ oz ¼ 0.99o1 the effect of (a) mesh size; (b) time step; (c) stretching constants; (d) mass balance. Vnew: instant water volume; Vinit: initial water volume.
Fig. 5. The comparison of present result with experimental and theoretical results reported by Faltinsen et al. (2005), d0/L ¼ d0/B ¼ 0.5. Excitation displacement a0/ L ¼ 0.0078, (a) surge motion, ox ¼ 1.037o1; (b) diagonal motion, ox ¼ oz ¼ 1.115, o1.
ARTICLE IN PRESS 504
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
excitation frequencies are simulated. The tank is a rigid with a square base (base-side length ¼ L) and d0/L ¼ 0.25. The amplitude of the ground displacement is 0.005L. The wave trains of different
sloshing waves are presented below and the detailed sloshing wave patterns and the associated resonance modes are also discussed separately in the following subsections.
1000 800 600 h/ ah
4.2.1. Diagonal waves Diagonal wave sloshing in a tank was firstly investigated by Miles (1994). Under diagonal forcing (excitation angle y ¼ 451) with the excitation frequency ox ¼ oz ¼ 0.9o1, (o1 ¼ o1,0), Fig. 7(a) plots the time histories of surface displacement at points A, C, D and E, and the free surface displacement at point E versus that at point F is depicted in Fig. 7(b). The coordinate of point A is (0.5,0,0.5) and others can be easy to be located in similar ways. The surface profile of wall W1 and W3 are shown in Fig. 7(c) and (d) and both of these figures present identical results due to diagonal forcing. The profiles of W3 and W4 are nearly antisymmetric, those of W1 and W2 have similar behaviors and the results of W2 and W4 are, therefore, omitted in the figure. As shown in Fig. 7(b), the relationship between free surface elevation at points E and F shows a perfect linear relationship. The consequence is due to almost linear and anti-symmetric
present Frandsen (2004)
400 200 0 -200 -400 40
50
60
70
80
90
t*ω1 Fig. 6. The wave history on tank’s corner (point A) under heave motion, the ratio d0/L ¼ d0/B ¼ 0.5, a0/L ¼ 0.026, oy ¼ 2o1, t: real time, ah: the initial perturbation ¼ 0.0001215 m (Frandsen, 2004).
0.15 0.3 0.2
Elevation at F
Point A Point C Point E Point D
H 0.1 0.0
0.05 0.00 -0.05 -0.10
0
40
20
60 T
80
0. 5 0. 0. 4 0. 3 0. 2 1 -0 0 -0 .1 -0 .2 . -0 3 . 4 -0 .5
120
60 50 40 30 T 20 10
H0
W1
100
0
-0.15 -0.15 -0.10 -0.05 0.00 0.05 0.10 0.15 Elevation at E
60 50 40 30 T 20 10
H0 -0 -0 .5 . -0 4 -0 .3 -0 .2 .1 0 0. 0. 1 2 0. 3 0. 0. 4 5
-0.1 -0.2
0.10
W3
0
Fig. 7. The wave history (7a, 7b) and wave profiles (7c, 7d) of a diagonal wave. ox ¼ oz ¼ 0.9o1, y ¼ 451.
Fig. 8. The wave history of a single-directional wave. ox ¼ oz ¼ 0.4o1, y ¼ 51.
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
sloshing displacement at points E and F. The displacements at points A and C are nearly anti-symmetric, while the displacements at points B and D are almost negligibly small and they are omitted from the figure. 4.2.2. Single-directional waves In this section, the tank is also under coupled surge–sway motion, the excitation angle is 51, and the excitation frequency is 0.4o1, far from the first fundamental frequency. The occurrence of single-directional waves can be obtained in a tank excited with excitation frequencies from 0.1o1 to 0.7o1. The histories of free surface displacement at points A, C, E and D are plotted in Fig. 8(a) and the displacement at point E versus that at point F is shown in Fig. 8(b). As shown in Fig. 8(b), the relationship between the displacements at points E and F is perfectly linear. Fig. 8(b) implies that the wave in the tank is at an angle of line 51 to the x-axis. We call this kind wave a single-directional wave. The sloshing displacements of point A and D show double troughs and peaks and the reason for the occurrence of double peaks (troughs) will be discussed in the section of spectral analysis. 4.2.3. Square-like waves Under the coupled surge–sway motion, when we increase the excitation frequency to 1.5o1, the waves become square-like. The
sloshing histories of surface displacements at points A, C, E and D are plotted in Fig. 9(a), the free surface displacement at point E versus that at point F is shown in Fig. 9(b), and the surface wave profile of W1 and W3 with different excitation angle versus time are depicted in Fig. 9(c)–(f). The upper two surface profiles are for an excitation angle of 51 and the lower two are for an excitation angle of 451. The relationship of sloshing displacement at points E and F is no longer linear but the phenomenon remains linear-like. Besides, for an excitation angle of 451, the surface profiles of W3 show uniform peaks between AE during T ¼ 63–67 that implies a terraced plane occurs.
4.2.4. Swirling waves When the sloshing waves move along the tank walls in a clockwise or counterclockwise direction, these kinds of waves are called ‘‘swirling’’ waves. In general, the swirling waves were found when the excitation frequencies of the tank were close to the first natural frequency. In the present study, we simulated a tank excited by coupled surge–sway motion with excitation angle y ¼ 51 and an excitation frequency ¼ 0.97o1. The corresponding results are shown in Fig. 10. The wave history of points A (up to T ¼ 3700), E and D are depicted in Fig. 10(a) and the relationship between sloshing displacement at points E and F is plotted in Fig. 10(c). As shown in Fig. 10(a), a ‘‘beating’’ phenomenon is found and a switching of rotation direction for swirling waves can be
Point A Point C Point E Point D
0.10 0.05 H 0.00 -0.05 0
H
H
20
40 T
60
505
80
H
H
Fig. 9. The wave history and wave profiles of a square-like wave. ox ¼ oz ¼ 1.5o1. (a)–(d): y ¼ 51; (e) and (f): y ¼ 451.
ARTICLE IN PRESS 506
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
clearly identified in Fig. 10(c). In general, a beat wave is the product of two waves with slightly different frequencies added together. When a tank excites near the first natural frequency, the beating phenomenon occurs. The beating does not disappear over a long period of time since the damping is very small of the fluid motion inside a smooth tank with no internal members to block the flow, and no heavy impact on the tank ceiling. Under the assumption of invisid fluid, no artificial damping is added and the numerical damping is very small in the present numerical scheme and the beating phenomenon still presents up to T ¼ 3700 (see Fig. 10(b)). The pattern depicted in
Fig. 10(c) is not as circular as theoretical predictions reported by Faltinsen et al. (2003). This might be due to the transient beating waves cause different sloshing displacements passing through points E and F. The snapshots of free surface for counterclockwise swirling waves are presented in Fig. 11. We increase the water depth in the tank to d0/L ¼ 1 and the corresponding swirling waves are presented in Fig. 12 that shows a nearly circular relationship between the elevation of points E and F even in a very early transient state. The depth might be key effect of the circular relationship between the sloshing displacement of points E and F.
Fig. 10. The wave history of a swirling wave. ox ¼ oz ¼ 0.97o1, y ¼ 51.
Fig. 11. The free surface profile of a swirling wave. ox ¼ oz ¼ 0.97o1, y ¼ 51.
Fig. 12. The wave history of a swirling wave. The ratio d0/L ¼ d0/B ¼ 1, a0/L ¼ 0.005 ox ¼ oz ¼ 1.03o1, y ¼ 51.
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
4.2.5. Irregular waves Waves that slosh irregularly inside the tank are named ‘‘irregular’’ waves or ‘‘chaotic’’ waves. The irregular waves can be found as the tank is under an excitation frequency close to the first natural frequency. In the present study, irregular waves also can be found when the tank’s excitation frequency is far away from the first resonant frequency. Fig. 13 presents the results of the tank excited at y ¼ 51 with an excitation frequency of 2.3o1. The elevation of irregular waves at points A, D and E, are shown in Fig. 13(a). The relationship between free surface displacements at points E and F is shown in Fig. 13(b), and a chaotic relation is observed. The surface profiles for an excitation angle of 51 are shown in Fig. 13(c) and (d), and the lower two are surface profiles (Fig. 13(e) and (f)) with an excitation angle of 451. The surface displacements along four walls all show an irregular pattern for excitation angles of 51 and 451, while the corresponding profiles are omitted in the text. 4.3. Spectral analysis In this section, spectral analyses are made to obtain the dominant resonant frequencies of each type of sloshing wave mentioned in the previous sections. The results of transient and steady-state conditions are presented. According to the assumption of invicid fluid in the present study, the damping is contributed by the numerical damping due to using the upwind scheme in the computational domain. However, it is hard to obtain the steady-state result because of the small numerical
507
damping. It is well known that resonant sloshing results from the first fundamental frequency, called ‘‘the primary resonance’’. If the tank is excited with a frequency away from the first fundamental mode of the tank, the secondary resonance also can be triggered due to the effects of the other natural modes of the system. Fig. 14 presents the spectral analyses for each type. For a excitation angle of 51, and an excitation frequency of 0.4o1 (Fig. 14(a)), single-directional waves occur. Two resonant peaks are identified in their spectra, one corresponding to the first fundamental frequency and the other to the excitation frequency. The latter is dominant in both transient and steady-state conditions. As we know, the beating phenomenon occurs as the excitation frequency, oe, of the tank close to the first natural frequency and the beating period is equal to j2p=ðoe o1 Þj. For a non-resonant excitation, the excitation frequency in addition to the first natural frequency will also affect the resonant frequency of the tank and the sloshing period of the fluid in the tank is also equal to j2p=ðoe o1 Þj. In the present study, the dimensionless period is equal to 14.6. Fig. 15 illustrates sloshing displacement at point A, the tank displacement and the harmonic sine history based on the first fundamental mode o1. As shown in the figure, the first double troughs (peaks) occur around T ¼ 15 and 29, respectively. The trend of tank displacement and that of the harmonic history clearly demonstrate the combined effects of the first fundamental mode and excitation mode (0.4o1). The sloshing displacement at point A, therefore, consists of primary resonance (due to the first fundamental mode) and the excitation mode.
Fig. 13. The wave history and wave profiles of a irregular wave. ox ¼ oz ¼ 2.3o1. (a)–(d): y ¼ 51; (e) and (f): y ¼ 451.
ARTICLE IN PRESS 508
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
Fig. 14. The spectral analysis at point A of (a) single-direction waves (0.4o1), (b) square-like waves (1.5o1), (c) swirling waves (0.97o1) and (d) irregular waves (2.3o1), y ¼ 51.
considering the nature of the mode shape of the resonant modes, as discussed in Section 4.4. 4.4. Wave patterns and resonance modes
Fig. 15. The influence of resonant modes for single-directional waves.
For excitation angles of 51, and excitation frequency increasing to 1.5o1 (Fig. 14(b)), square-like waves appear and the dominant resonant frequency is the first fundamental frequency during transient period with another resonant peak at the excitation frequency. However, when it reaches steady state, the dominant one becomes the excitation frequency. For the same excitation angle but the excitation frequency is much nearer the first fundamental frequency (Fig. 14(c)), a few little peaks appear that indicates the secondary resonance corresponds to 2ox and o2,2. The secondary resonance is most likely related to the occurrence of the swirling waves since the spectra of the other cases of swirling waves all present resonant peaks corresponding to 2ox and o2,2. As the excitation frequency increases further to 2.3o1 (o3,2), irregular waves occur. Although the dominant resonant frequency is the excitation frequency, several secondary resonances are present and they correspond to modes o1,0, o3,0(2.13o1) and o5,0 (2.75o1)during the transient period. These are also likely to be related to the occurrence of irregular waves. As it comes to stable condition (T ¼ 1800–2500), all the odd modes of natural frequencies disappear and only the excitation frequency exists. The said phenomena may be clarified by
Faltinsen et al. (2000) presented two subclasses of wave patterns exist. The first one consists of a two-dimensional Stokes wave pattern, and the corresponding waves are called planar waves. They are, 1 ; mX0 fð1Þ m ðxÞ ¼ cos mp x þ 2 1 ; nX0 (12) fð2Þ n ðyÞ ¼ cos npr y þ 2r The second subclass is a three-dimensional wave pattern given by the product of two Stokes waves ð2Þ fð1Þ m ðxÞfn ðyÞ
(13)
Faltinsen also suggested the mixed modes and the combination of the two Stokes modes to represent three-dimensional wave patterns. As stated in the previous sub-section, the major resonant frequency can be identified from the spectral analysis of each type of sloshing wave. The spectral analysis of a square-like wave presents two peaks, and they correspond to the first fundamental frequency and excitation frequency. The combination of o1,0 with 1.5 cos(301) weight and o0,1 with 1.5 sin(301) weight is depicted in Fig. 16(a) and clear the terraced planes can be found in Fig. 16(b). The spectral analysis of irregular waves indicates the resonant peaks of o1,0, o3,0 and o3,2, and the mode shape from the addition of those three resonant modes and the corresponding surface wave contours are depicted in Fig. 16(c) and (d). The odd modes (o1,0, o3,0) vanish when the sloshing comes to the steady state and only o3,2 remains. The wave pattern o3,2 is shown on the top right-hand corner of Fig. 16(c). The effect of odd modes o1,0, o3,0 on an irregular wave during transient state and that of only o3,2 during steady state can be noted from the surface wave contours.
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
509
Fig. 16. Wave patterns and resonance node (a) square-like waves, y ¼ 301; (b) irregular waves, y ¼ 51.
Fig. 16 provides clear evidence of the correlation between sloshing surface patterns and resonant modes.
5. Conclusion A three-dimensional time-independent finite difference method is developed to study the sloshing waves and resonant modes of fluid in a 3D tank. Only coupled surge–sway motions with various excitation angles under a wide range of excitation frequencies are presented. The following conclusions are reached. 1. The method is used to solve a fully nonlinear sloshing fluid in a three-dimensional tank with square base. The numerical scheme is validated by rigorous benchmark tests. 2. Five types of sloshing waves are observed and they are closely related to the excitation frequency. When the sloshing frequency is close to the first fundamental mode, swirling wave patterns occur. The square-like and irregular wave patterns were also observed in the literature (Faltinsen et al., 2003) when the tank is under near-resonant excitation. In the present study, the square-like wave and irregular wave also can be triggered when the excitation frequency is far away from the first fundamental frequency with an excitation angle of 51. 3. The resonant peaks of spectral analyses were shown to correspond to the primary mode and excitation frequency for the diagonal, single-directional and square-like waves. For the swirling wave, the secondary resonance of mode o2,2 occur in addition to the primary resonance. For irregular waves, the secondary resonance of odd modes o3,0 and o5,0 also appear along with the resonant peaks of the primary mode and the excitation frequency.
4. The surface patterns are generated by the combining of fundamental modes that reflect resonance. The plots of mode shape of the resonant modes and the surface contour of the wave pattern demonstrate a clear correlation between the occurrence of each type of sloshing wave and corresponding resonant modes. 5. The method developed can be used to simulate 3D fluid sloshing in a tank excited by motions including all six-degrees of freedom. Acknowledgement This study is supported by National Science Council under a Grant NSC 96-2221-E-110-107-MY3. Reference Aliabadi, S., Johnson, A., Abedi, J., 2003. Comparison of finite element and pendulum models for simulation of sloshing. Computers and Fluids 23, 535–545. ¨ nal, E., 2005. Experimental investigation of pressure distribution on a Akyildiz, H., U rectangular tank due to the liquid sloshing. Ocean Engineering 32, 1503–1516. ¨ nal, N.E., 2006. Sloshing in a three-dimensional rectangular tank: Akyildiz, H., U numerical simulation and experimental validation. Ocean Engineering 33, 2135–2149. Celebi, M.S., Akyildiz, H., 2002. Nonlinear modeling of liquid sloshing in moving rectangular tank. Ocean Engineering 29, 1527–1553. Chen, B.F., Chiang, S.W., 1999a. Complete 2D and fully nonlinear analysis of ideal fluid in tanks. Journal of Engineering Mechanics 125 (1), 70–78. Chen, B.F., Yuan, Y.S., Lee, J.F., 1999b. Three dimensional nonlinear hydrodynamic pressures by earthquakes on dam faces with arbitrary reservoir shapes. Journal of Hydraulic Research 37, 163–187. Chen, B.F., Nokes, R., 2005. Time-independent finite difference analysis of fully non-linear and viscous fluid sloshing in a rectangular tank. Journal of Computational Physics 209, 47–81.
ARTICLE IN PRESS 510
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
Faltinsen, O.M., 1978. A numerical non-linear method of sloshing in tanks with two dimensional flow. Journal of Ship Research 22 (3), 193–202. Faltinsen, O.M., Rognebakke, O.F., Lukovsky, I.A., Timokha, A.N., 2000. Multidimensional modal analysis of nonlinear sloshing in a rectangular tank with finite water depth. Journal of Fluid Mechanics 407, 201–234. Faltinsen, O.M., Timokha, A.N., 2001. Adaptive multimodal approach to nonlinear sloshing in a rectangular tank. Journal of Fluid Mechanics 432, 167–200. Faltinsen, O.M., Timokha, A.N., 2002. Asymptotic modal approximation of nonlinear resonant sloshing in a rectangular tank with small fluid depth. Journal of Fluid Mechanics 470, 319–357. Faltinsen, O.M., Rognebakke, O.F., Timokha, A.N., 2003. Resonance three-dimensional nonlinear sloshing in a square-base basin. Journal of Fluid Mechanics 487, 1–42. Faltinsen, O.M., Rognebakke, O.F., Timokha, A.N., 2005. Classification of threedimensional nonlinear sloshing in a square-base tank with finite depth. Journal of Fluids and Structures 20, 81–103. Frandsen, J.B., 2004. Sloshing motions in excited tanks. Journal of Computational Physics 196, 53–87. Hill, D.F., 2003. Transient and steady-state amplitudes of forced waves in rectangular basins. Physics of Fluids 15 (6), 1576–1587. Hung, T.K., Wang, M.H., 1987. Nonlinear hydrodynamic pressure on rigid dam motion. Journal of Engineering Mechanics, ASCE 113 (4), 482–499.
Kim, Y., 2001. Numerical simulation of sloshing flows with impact load. Applied Ocean Research 23, 53–62. Kim, Y., Nam, B.W., Kim, D.W., Kim, Y.S., 2007. Study on coupling effects of ship motion and sloshing. Ocean Engineering 34, 2176–2187. Lee, S.J., Kim, M.H., Lee, D.H., Kim, J.W., Kim, Y.H., 2007a. The effects of LNG-tank sloshing on the global motions of LNG carriers. Ocean Engineering 34, 10–20. Lee, D.H., Kim, M.H., Kwon, S.H., Kim, J.W., Lee, Y.B., 2007b. A parametric sensitivity study on LNG tank sloshing loads by numerical simulations. Ocean Engineering 34, 3–9. Miles, J.W., 1994. Faraday waves: rolls versus squares. Journal of Fluid Mechanics 269, 353–371. Ockendon, J.R., Ockendon, H., Waterhouse, D.D., 1996. Multi-mode resonances in fluids. Journal of Fluid Mechanics 315, 317–344. Turnbull, M.S., Borthwick, A.G.L., Eatock Taylor, R., 2003. Numerical wave tank based on a s-transformation finite element inviscid flow solver. International Journal for Numerical Methods in Fluids 42, 641–663. Wu, G.X., Ma, Q.W., Taylor, R.E., 1998. Numerical simulation of sloshing waves in a 3D tank based on a finite element methods. Applied Ocean Research 20, 337–355. Wu, G.X., 2007. Second-order resonance of sloshing in a tank. Ocean Engineering 34, 2345–2349.
Contents lists available at ScienceDirect
Ocean Engineering journal homepage: www.elsevier.com/locate/oceaneng
Sloshing waves and resonance modes of fluid in a 3D tank by a time-independent finite difference method Chih-Hua Wu, Bang-Fuh Chen Department of Marine Environment and Engineering, National Sun Yat-sen University, Kaohsiung 804, Taiwan
a r t i c l e in f o
a b s t r a c t
Article history: Received 25 June 2008 Accepted 31 January 2009 Available online 13 February 2009
A 3D time-independent finite difference method is developed to solve for wave sloshing in a threedimensional tank excited by coupled surge and sway motions. The 3D equations of fluid motion are derived in a moving coordinate system. The three-dimensional tank, with an arbitrary depth and a square base, is subjected to a range of excitation frequencies with motions that exhibit multiple degrees of freedom. For demonstration purposes the numerical scheme is validated by a benchmark study. Five types of sloshing waves were observed when the tank is excited by various excitation frequencies. A spectral analysis identified the resonant frequencies of each type of wave and the results show a strong correlation between resonant modes and the occurrence of the sloshing wave types. The method can be used to simulate fluid sloshing in a 3D tank with six-degrees of freedom. & 2009 Elsevier Ltd. All rights reserved.
Keywords: Fully nonlinear free surface wave Sloshing fluids 3D tank Finite difference method
1. Introduction Free surface sloshing in a moving container is associated with various engineering problems, such as tank trucks on highways, liquid oscillations in large storage tanks caused by earthquakes, sloshing of liquid cargo in ocean-going vessels and the motion of liquid fuel in aircraft and spacecraft. It is known that partially filled tanks are prone to violent sloshing under certain motions, especially when near-resonant excitation occurs. The large liquid movement creates highly localized impact pressures on tank walls, which in turn cause structural damage and may even create moments that affect the stability of the vehicle, which carries the container. Sloshing waves in moving tanks have been studied numerically, theoretically and experimentally in the past several decades and many significant phenomena have been considered in those studies, especially the linear and nonlinear effects of sloshing for both inviscid and viscous liquids. Most reported studies involved tanks excited by limited excitation directions and with a fixed excitation frequency throughout the excitation. In reality, as the tank is excited by accelerations due to an earthquake or waves, the excitation directions include multi-degrees of freedom (surge/ sway/heave/pitch/roll/yaw) and the excitation frequency also varies with time. The potential formulation of the problem is often used in studying sloshing, for example Ockendon et al. (1996) among others. The most distinguished analytical works are Faltinsen’s series of studies (Faltinsen, 1978; Faltinsen and
Corresponding author. Tel./fax: +886 7 5255065.
E-mail address: [email protected] (B.-F. Chen). 0029-8018/$ - see front matter & 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.oceaneng.2009.01.020
Timokha, 2001, 2002) for sloshing fluids in 2D tanks and Faltinsen et al. (2005) where their asymptotic modal system is extended to model nonlinear sloshing in a 3D rectangular tank. Besides the potential flow approaches, many numerical studies (computational fluid dynamics) of the problems with primitive variables have been made. These have focused particularly on the fully nonlinear free surface effects. Many papers give successful examples for two-dimensional sloshing (see, for instance, Chen and Chiang, 1999a; Celebi and Akyildiz, 2002; Turnbull et al., 2003; Aliabadi et al., 2003; Frandsen, 2004; Chen and Nokes, 2005; Wu, 2007) and for three-dimensional sloshing (Wu et al., ¨ nal, 2005, 2006; Kim et al., 2007; 1998; Kim, 2001; Akyildiz and U Lee et al., 2007a, b). The reported techniques for handling a wavyfree surface include VOF, SOLA, SURF and also the s-transformation technique to stretch the grid from the bed to surface. However, the 3D numerical simulation of fluid sloshing in a tank is still very limited. In the present study, a 3D tank with different ratios of depth/excitation amplitude, multiple degrees of freedom of excitation and excitation frequencies are considered. In this three-dimensional model, we develop a 3D time-independent finite difference method to incorporate the incompressible and inviscid Navier–Stokes equations, fully nonlinear kinematic and dynamic free surface conditions in the analysis of the seismic response of sloshing fluid in a rectangular tank with a square base. The time-varying moving boundary is mapped onto a timeindependent domain through proper transformation functions and a special finite difference approximation is made in order to overcome the difficulty of maintaining the accuracy of the finite difference expression for the second derivative when the difference mesh is stretched near the boundary. The main focus of this paper is the simulation of a 3D tank undergoing different
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
combinations of motions with varying excitation directions, and not only transient but nearly steady-state phenomena are recorded and studied. Hill (2003) compared the transient waves with the steady-state waves and showed that the maximum transient response can far exceed the steady-state response of the basin. The transient effect on moving vehicles is more dangerous than that of steady-state condition. Section 2 introduces the equations of motion that are written in a moving coordinate system attached to the accelerating tank. The proper coordinate transformation functions are adopted to map the time-dependent domain into a fixed unit cubic. The proposed finite difference method is developed in Section 3. Section 4 presents the detailed results and provides comprehensive discussion of all the phenomena found in this study.
2. Mathematical formulation
and the kinematic and dynamic free surface conditions are @h @h @h þu þw ¼v @t @x @z
2 ð2a_ b_ b_ g_ 2 Þx 2ðb_ w g_ vÞ
(1)
@v @v @v @v 1 @p € þu þv þw ¼ gy Y C ðg€ x a€ zÞ @t @x @y @z r @y ð2b_ g_ a_ 2 g_ 2 Þy 2ðg_ u a_ wÞ
(2)
p¼0
(6)
The Poisson equation can be derived by taking partial derivatives of (1), (2) and (3) with respect to x, y and z, respectively, and summing the results. @2 p @2 p @2 p @ @u @u @u u þv þw þ 2 þ 2 ¼ r 2 @x @x @y @z @x @y @z @ @v @v @v r u þv þw @y @x @y @z @ @w @w @w r u þv þw @z @x @y @z 2 2r½a_ b_ þ b_ g_ þ a_ g_ Þ ða_ 2 þ b_ þ g_ 2 Þ @v @w 2r a_ @z @y @w @u @u @v þb_ þ g_ @x @z @y @x
2
ð2a_ g_ a_ 2 b_ Þz 2ða_ v b_ uÞ
(3)
where u, v and w are velocity components in x-, y- and z-directions, X€ C , Y€ C and Z€ C are the acceleration components of the tank in x, y and z directions; (a_ , b_ , g_ ) and (a€ , b€ , g€ ) are the corresponding angular velocities and accelerations with respect to the x, y and z axes, p is the pressure, r is the fluid density, and gx,gy and gz are the components of the acceleration due to gravity. The continuity equation for incompressible flow is @u @v @w þ þ ¼0 @x @y @z
(4)
x b1 ðy; zÞ b2 ðy; zÞ b1 ðy; zÞ z b3 ðx; yÞ zn ¼ b4 ðx; yÞ b3 ðx; yÞ
yn ¼ 1
y þ d0 hðx; z; tÞ (8)
where the instantaneous water depth, h(x,z,t), is a single-valued function measured from tank bottom, d0 represents the vertical distance between the still water surface and the tank bottom, b1(y,z) and b2(y,z) are horizontal distances from the z-axis to the west (W1) and east (W2) walls, respectively, and b3(x,y) and b4(x,y) are horizontal distance from the x-axis to the south (W3) and north (W4) walls, respectively (see Fig. 1). Based on the assumption of h(x,z,t), the phenomena of wave breaking, run-up and tank roof impacts are not considered in the present study. With the above equations, one can map the west wall to x* ¼ 0
θ
Y
X
Yaw
Z Surge
X Y
D
Pitch
Z A Sway
E
h (x, z, t) b1 (y, z)
N
b4 (x, y) b2(y, z) b3 (x, y)
W3
F
X
W4
Z
W1
C
G
I
Roll
(7)
Many finite difference (FDM) and finite volume methods have been reported to solve the free surface displacement of sloshing fluid in tanks. A brief review of those methods can be found in Chen and Nokes (2005). In the present study, the irregular boundary, such as the time-varying fluid surface, non-vertical walls and non-horizontal bottom, can be mapped onto a square by the proper coordinate transformations (Hung and Wang, 1987; Chen et al., 1999b) as follows: xn ¼
@w @w @w @w 1 @p € þu þv þw ¼ gz Z C ða€ y b€ xÞ @t @x @y @z r @z
Heavy
(5)
and
The momentum equations are written based on a moving coordinate system (the coordinate system is chosen to move with the tank, including surge, sway, heave, yaw, roll and pitch motions, see Fig. 1). @u @u @u @u 1 @p € þu þv þw ¼ gx X C ðb€ z g€ yÞ @t @x @y @z r @x
501
E L
Fig. 1. Definition sketch of the tank motion.
B
W2
B
ARTICLE IN PRESS 502
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
Fig. 2. The concepts of transformation coordinate systems and stretch of the grid. The stretching constants presented in (x*, y*, z*) coordinate are l1 ¼ l2 ¼ l3 ¼ 0.5, k1 ¼ k2 ¼ k3 ¼ 2.
and east wall to x* ¼ 1, the north wall to z* ¼ 0 and south wall to z* ¼ 1, the free surface to y* ¼ 0 and the tank bottom to y* ¼ 1. In this way the computational domain is transformed to a fixed unit cubic domain. The coordinates (x*, y*, z*) can be further transformed such that the layers near the boundary are stretched to capture the sharp local velocity gradients. The following exponential functions provide these stretching transformations: n ðxn 1Þ
X ¼ l1 þ ðxn l1 Þek1 x Z ¼ l3 þ ðzn l3 Þek3 z
n ðyn 1Þ
Y ¼ l2 þ ðyn l2 Þek2 y
n ðzn 1Þ
(9)
The constants ki and li, i ¼ 1,2,3, control the mesh size and stretching in X-, Y- and Z-directions, respectively. The concepts of transformation coordinate system and grid stretching are presented in Fig. 2. The geometry of the flow field and the meshes in the computational domain (X–Y–Z system) become timeindependent throughout the computational analysis. The dimensional parameters are normalized in the following way: u v w p U ¼ pffiffiffiffiffiffiffiffi V ¼ pffiffiffiffiffiffiffiffi W ¼ pffiffiffiffiffiffiffiffi P ¼ rgd0 gd0 gd0 gd0 rffiffiffiffiffi g h xc y zc H¼ T¼t Xc ¼ Y c ¼ c Zc ¼ d0 d0 d0 d0 d0 sffiffiffiffiffi a b g a_ d0 Ya ¼ Yb ¼ Yg ¼ YaT ¼ 2p 2p 2p 2p g sffiffiffiffiffi sffiffiffiffiffi b_ d0 g_ d0 a€ d0 YbT ¼ YgT ¼ YaTT ¼ 2p g 2p g 2p g
YbTT ¼
b€ d0 2p g
YgTT ¼
g€ d0 2p g
With the aforementioned transformations and dimensionless variables, Eqs. (1)–(7) can be written in dimensionless form. For example, the momentum equation in x-direction is: U T þ ðC 10 C 13 U X þ C 11 C 14 U Y þ C 12 C 15 U Z Þ þ C 1 C 13 UU X þ C 2 C 14 UU Y þ C 3 C 15 UU Z þ C 4 C 13 VU X þ C 5 C 14 VU Y þ C 6 C 15 VU Z þ C 7 C 13 WU X þ C 8 C 14 WU Y þ C 9 C 15 WU Z ¼ Gx ðC 1 C 13 P X þ C 2 C 14 P Y þ C 3 C 15 PZ Þ 2p 4p2 ðzYbTT yYgTT Þ ð2YaT YbT Y2bT Y2gT Þx d0 d0 (10) 4pðYbT W YgT VÞ X CTT
The dimensionless forms of the other equations are omitted in the text. In Eq. (10), the coefficients C1–C15 are due to coordinate transformations and can be found in Chen et al. (1999b). PX denotes a partial derivative of P with respect to X; UT is the partial derivative of U with respect to dimensionless time T; XCTT, YCTT and ZCTT are dimensionless ground accelerations in x, y and z directions; the other terms have similar meanings.
3. Computational algorithm In this three-dimensional analysis, the fluid flow is solved in a unit cubic mesh in the transformed flow domain. All the computations use the dimensionless equations in X–Y–Z coordinate system. The difference equations for space derivatives use the central difference approximation, except at the boundaries were forward or backward differences are employed. A staggered grid system is used in the analysis. That is, the pressure P is defined at the center of a cell, whereas the velocity components U, V and W are calculated 0.5DX, 0.5DY and 0.5DZ behind, below or forward of the cell center. The Crank–Nicholson second-order finite difference scheme and the Gauss–Seidel Point successive over-relaxation iterative procedure are used to calculate the velocity and pressure, respectively. However, it is sensitive to numerical solution for the convective term and the second-order central difference upwind scheme is applied in the present numerical scheme. The numerical scheme is similar to that described in Chen and Nokes (2005) and the detailed iteration procedures are omitted here.
4. Results and discussions The natural modes of a 3D tank with square base can be obtained by solving a linearized natural sloshing standing wave problem where the circular frequency can be expressed as qffiffiffiffiffiffiffiffiffiffiffiffiffiffi o2i;j ¼ g li;j tanhðli;j d0 Þ; li;j ¼ p i2 þ j2 (11) and i, j are the natural mode’s components in the x and z directions. The FFT analysis of ChiChi earthquake is presented in Fig. 3 and a wide range of excitation frequencies between 0.4o1 and 2.5o1 is considered in this study. The accuracy of the numerical results significantly depends on the spatial grid resolution and the selected time step. Thus, it is important to verify the appropriate time step for a particular spatial grid before large numerical problems are simulated. In the present study, the time step is restricted by Dto minfDxmin =jui;j;k j; Dymin =jvi;j;k j; Dzmin =jwi;j;k jg, which implies that a fluid particle cannot move more than one cell in a single time step. The selection of the spatial grid was also based upon a parametric study. Fig. 4 shows the sloshing displacement at point A with various mesh sizes selections. As depicted in the figure, the results of mesh size (20 10 20) agree well with those of Kim (2001), whose mesh size was (30 30 30). The results of finer meshes and smaller time-step selection present minor difference from that of mesh size ¼ (20 10 20). The influence of stretching constants is illustrated in Fig. 4(c). The mass balance is shown in Fig. 4(d) and good numerical accuracy is noted. Therefore, the mesh size of (20 10 20), DT ¼ 0.006, l1 ¼ l2 ¼ l3 ¼ 0.5 and k1 ¼ k2 ¼ k3 ¼ 2 are used in this study. The simulation of the cases
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
with strong shallow-fluid sloshing is beyond the limits of the present numerical model. The proposed numerical scheme can be applied in a non-resonant oscillation with a0/d0 ¼ 2.2 and the d0/L is limited to 0.1. In a near-resonant oscillation, however, the limit of a0/d0 becomes 0.05 in the present study.
503
benchmark tests. The comparison between the present numerical results and the experimental and theoretical results reported by Faltinsen et al. (2005) are shown in Fig. 5 and the agreement is very good. Fig. 6 further compares the sloshing displacements of fluid in a tank excited by vertical motion with those reported by Frandsen (2004). The agreement is also very good.
4.1. Benchmark tests In order to validate the accuracy of numerical simulation in the present study, the results obtained from the present numerical model are compared with those reported in the literature for
Fig. 3. The spectral analysis of 921, 1999 ChiChi earthquake in Taiwan o1: The fundamental natural frequency of a tank with d0/L ¼ d0/B ¼ 0.25.
4.2. Horizontal ground motion In a real marine environment, the excitation forcing of the tank might involve multiple directions (surge, sway, heave, pitch, yaw and roll). In the present study, we focus on horizontal excitation only. The tank, in general, is not excited exactly in surge (x-) or sway (z-) directions. Rather it is excited in a coupled surge and sway motion. We define y as the excitation angle of the horizontal ground motion (see Fig. 1) that represents the various coupled surge/sway directions. For longitudinal excitation (surge motion, y ¼ 01), the work of Faltinsen et al. (2005) indicated four different kinds of waves, namely planar wave, swirling wave, irregular wave and square-like wave, could be generated through this type of forcing. As the tank excitation angle y is equal to 451 (diagonal motion), most of the four kinds of waves can be observed while the stable square-like waves are seldom seen. In the present study, horizontal excitations with various excitation angles and
Fig. 4. The wave history at corner A of the tank under diagonal motion, the ratio d0/L ¼ d0/B ¼ 0.25, the ratio of excitation displacement a0/L ¼ 0.0093, ox ¼ oz ¼ 0.99o1 the effect of (a) mesh size; (b) time step; (c) stretching constants; (d) mass balance. Vnew: instant water volume; Vinit: initial water volume.
Fig. 5. The comparison of present result with experimental and theoretical results reported by Faltinsen et al. (2005), d0/L ¼ d0/B ¼ 0.5. Excitation displacement a0/ L ¼ 0.0078, (a) surge motion, ox ¼ 1.037o1; (b) diagonal motion, ox ¼ oz ¼ 1.115, o1.
ARTICLE IN PRESS 504
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
excitation frequencies are simulated. The tank is a rigid with a square base (base-side length ¼ L) and d0/L ¼ 0.25. The amplitude of the ground displacement is 0.005L. The wave trains of different
sloshing waves are presented below and the detailed sloshing wave patterns and the associated resonance modes are also discussed separately in the following subsections.
1000 800 600 h/ ah
4.2.1. Diagonal waves Diagonal wave sloshing in a tank was firstly investigated by Miles (1994). Under diagonal forcing (excitation angle y ¼ 451) with the excitation frequency ox ¼ oz ¼ 0.9o1, (o1 ¼ o1,0), Fig. 7(a) plots the time histories of surface displacement at points A, C, D and E, and the free surface displacement at point E versus that at point F is depicted in Fig. 7(b). The coordinate of point A is (0.5,0,0.5) and others can be easy to be located in similar ways. The surface profile of wall W1 and W3 are shown in Fig. 7(c) and (d) and both of these figures present identical results due to diagonal forcing. The profiles of W3 and W4 are nearly antisymmetric, those of W1 and W2 have similar behaviors and the results of W2 and W4 are, therefore, omitted in the figure. As shown in Fig. 7(b), the relationship between free surface elevation at points E and F shows a perfect linear relationship. The consequence is due to almost linear and anti-symmetric
present Frandsen (2004)
400 200 0 -200 -400 40
50
60
70
80
90
t*ω1 Fig. 6. The wave history on tank’s corner (point A) under heave motion, the ratio d0/L ¼ d0/B ¼ 0.5, a0/L ¼ 0.026, oy ¼ 2o1, t: real time, ah: the initial perturbation ¼ 0.0001215 m (Frandsen, 2004).
0.15 0.3 0.2
Elevation at F
Point A Point C Point E Point D
H 0.1 0.0
0.05 0.00 -0.05 -0.10
0
40
20
60 T
80
0. 5 0. 0. 4 0. 3 0. 2 1 -0 0 -0 .1 -0 .2 . -0 3 . 4 -0 .5
120
60 50 40 30 T 20 10
H0
W1
100
0
-0.15 -0.15 -0.10 -0.05 0.00 0.05 0.10 0.15 Elevation at E
60 50 40 30 T 20 10
H0 -0 -0 .5 . -0 4 -0 .3 -0 .2 .1 0 0. 0. 1 2 0. 3 0. 0. 4 5
-0.1 -0.2
0.10
W3
0
Fig. 7. The wave history (7a, 7b) and wave profiles (7c, 7d) of a diagonal wave. ox ¼ oz ¼ 0.9o1, y ¼ 451.
Fig. 8. The wave history of a single-directional wave. ox ¼ oz ¼ 0.4o1, y ¼ 51.
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
sloshing displacement at points E and F. The displacements at points A and C are nearly anti-symmetric, while the displacements at points B and D are almost negligibly small and they are omitted from the figure. 4.2.2. Single-directional waves In this section, the tank is also under coupled surge–sway motion, the excitation angle is 51, and the excitation frequency is 0.4o1, far from the first fundamental frequency. The occurrence of single-directional waves can be obtained in a tank excited with excitation frequencies from 0.1o1 to 0.7o1. The histories of free surface displacement at points A, C, E and D are plotted in Fig. 8(a) and the displacement at point E versus that at point F is shown in Fig. 8(b). As shown in Fig. 8(b), the relationship between the displacements at points E and F is perfectly linear. Fig. 8(b) implies that the wave in the tank is at an angle of line 51 to the x-axis. We call this kind wave a single-directional wave. The sloshing displacements of point A and D show double troughs and peaks and the reason for the occurrence of double peaks (troughs) will be discussed in the section of spectral analysis. 4.2.3. Square-like waves Under the coupled surge–sway motion, when we increase the excitation frequency to 1.5o1, the waves become square-like. The
sloshing histories of surface displacements at points A, C, E and D are plotted in Fig. 9(a), the free surface displacement at point E versus that at point F is shown in Fig. 9(b), and the surface wave profile of W1 and W3 with different excitation angle versus time are depicted in Fig. 9(c)–(f). The upper two surface profiles are for an excitation angle of 51 and the lower two are for an excitation angle of 451. The relationship of sloshing displacement at points E and F is no longer linear but the phenomenon remains linear-like. Besides, for an excitation angle of 451, the surface profiles of W3 show uniform peaks between AE during T ¼ 63–67 that implies a terraced plane occurs.
4.2.4. Swirling waves When the sloshing waves move along the tank walls in a clockwise or counterclockwise direction, these kinds of waves are called ‘‘swirling’’ waves. In general, the swirling waves were found when the excitation frequencies of the tank were close to the first natural frequency. In the present study, we simulated a tank excited by coupled surge–sway motion with excitation angle y ¼ 51 and an excitation frequency ¼ 0.97o1. The corresponding results are shown in Fig. 10. The wave history of points A (up to T ¼ 3700), E and D are depicted in Fig. 10(a) and the relationship between sloshing displacement at points E and F is plotted in Fig. 10(c). As shown in Fig. 10(a), a ‘‘beating’’ phenomenon is found and a switching of rotation direction for swirling waves can be
Point A Point C Point E Point D
0.10 0.05 H 0.00 -0.05 0
H
H
20
40 T
60
505
80
H
H
Fig. 9. The wave history and wave profiles of a square-like wave. ox ¼ oz ¼ 1.5o1. (a)–(d): y ¼ 51; (e) and (f): y ¼ 451.
ARTICLE IN PRESS 506
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
clearly identified in Fig. 10(c). In general, a beat wave is the product of two waves with slightly different frequencies added together. When a tank excites near the first natural frequency, the beating phenomenon occurs. The beating does not disappear over a long period of time since the damping is very small of the fluid motion inside a smooth tank with no internal members to block the flow, and no heavy impact on the tank ceiling. Under the assumption of invisid fluid, no artificial damping is added and the numerical damping is very small in the present numerical scheme and the beating phenomenon still presents up to T ¼ 3700 (see Fig. 10(b)). The pattern depicted in
Fig. 10(c) is not as circular as theoretical predictions reported by Faltinsen et al. (2003). This might be due to the transient beating waves cause different sloshing displacements passing through points E and F. The snapshots of free surface for counterclockwise swirling waves are presented in Fig. 11. We increase the water depth in the tank to d0/L ¼ 1 and the corresponding swirling waves are presented in Fig. 12 that shows a nearly circular relationship between the elevation of points E and F even in a very early transient state. The depth might be key effect of the circular relationship between the sloshing displacement of points E and F.
Fig. 10. The wave history of a swirling wave. ox ¼ oz ¼ 0.97o1, y ¼ 51.
Fig. 11. The free surface profile of a swirling wave. ox ¼ oz ¼ 0.97o1, y ¼ 51.
Fig. 12. The wave history of a swirling wave. The ratio d0/L ¼ d0/B ¼ 1, a0/L ¼ 0.005 ox ¼ oz ¼ 1.03o1, y ¼ 51.
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
4.2.5. Irregular waves Waves that slosh irregularly inside the tank are named ‘‘irregular’’ waves or ‘‘chaotic’’ waves. The irregular waves can be found as the tank is under an excitation frequency close to the first natural frequency. In the present study, irregular waves also can be found when the tank’s excitation frequency is far away from the first resonant frequency. Fig. 13 presents the results of the tank excited at y ¼ 51 with an excitation frequency of 2.3o1. The elevation of irregular waves at points A, D and E, are shown in Fig. 13(a). The relationship between free surface displacements at points E and F is shown in Fig. 13(b), and a chaotic relation is observed. The surface profiles for an excitation angle of 51 are shown in Fig. 13(c) and (d), and the lower two are surface profiles (Fig. 13(e) and (f)) with an excitation angle of 451. The surface displacements along four walls all show an irregular pattern for excitation angles of 51 and 451, while the corresponding profiles are omitted in the text. 4.3. Spectral analysis In this section, spectral analyses are made to obtain the dominant resonant frequencies of each type of sloshing wave mentioned in the previous sections. The results of transient and steady-state conditions are presented. According to the assumption of invicid fluid in the present study, the damping is contributed by the numerical damping due to using the upwind scheme in the computational domain. However, it is hard to obtain the steady-state result because of the small numerical
507
damping. It is well known that resonant sloshing results from the first fundamental frequency, called ‘‘the primary resonance’’. If the tank is excited with a frequency away from the first fundamental mode of the tank, the secondary resonance also can be triggered due to the effects of the other natural modes of the system. Fig. 14 presents the spectral analyses for each type. For a excitation angle of 51, and an excitation frequency of 0.4o1 (Fig. 14(a)), single-directional waves occur. Two resonant peaks are identified in their spectra, one corresponding to the first fundamental frequency and the other to the excitation frequency. The latter is dominant in both transient and steady-state conditions. As we know, the beating phenomenon occurs as the excitation frequency, oe, of the tank close to the first natural frequency and the beating period is equal to j2p=ðoe o1 Þj. For a non-resonant excitation, the excitation frequency in addition to the first natural frequency will also affect the resonant frequency of the tank and the sloshing period of the fluid in the tank is also equal to j2p=ðoe o1 Þj. In the present study, the dimensionless period is equal to 14.6. Fig. 15 illustrates sloshing displacement at point A, the tank displacement and the harmonic sine history based on the first fundamental mode o1. As shown in the figure, the first double troughs (peaks) occur around T ¼ 15 and 29, respectively. The trend of tank displacement and that of the harmonic history clearly demonstrate the combined effects of the first fundamental mode and excitation mode (0.4o1). The sloshing displacement at point A, therefore, consists of primary resonance (due to the first fundamental mode) and the excitation mode.
Fig. 13. The wave history and wave profiles of a irregular wave. ox ¼ oz ¼ 2.3o1. (a)–(d): y ¼ 51; (e) and (f): y ¼ 451.
ARTICLE IN PRESS 508
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
Fig. 14. The spectral analysis at point A of (a) single-direction waves (0.4o1), (b) square-like waves (1.5o1), (c) swirling waves (0.97o1) and (d) irregular waves (2.3o1), y ¼ 51.
considering the nature of the mode shape of the resonant modes, as discussed in Section 4.4. 4.4. Wave patterns and resonance modes
Fig. 15. The influence of resonant modes for single-directional waves.
For excitation angles of 51, and excitation frequency increasing to 1.5o1 (Fig. 14(b)), square-like waves appear and the dominant resonant frequency is the first fundamental frequency during transient period with another resonant peak at the excitation frequency. However, when it reaches steady state, the dominant one becomes the excitation frequency. For the same excitation angle but the excitation frequency is much nearer the first fundamental frequency (Fig. 14(c)), a few little peaks appear that indicates the secondary resonance corresponds to 2ox and o2,2. The secondary resonance is most likely related to the occurrence of the swirling waves since the spectra of the other cases of swirling waves all present resonant peaks corresponding to 2ox and o2,2. As the excitation frequency increases further to 2.3o1 (o3,2), irregular waves occur. Although the dominant resonant frequency is the excitation frequency, several secondary resonances are present and they correspond to modes o1,0, o3,0(2.13o1) and o5,0 (2.75o1)during the transient period. These are also likely to be related to the occurrence of irregular waves. As it comes to stable condition (T ¼ 1800–2500), all the odd modes of natural frequencies disappear and only the excitation frequency exists. The said phenomena may be clarified by
Faltinsen et al. (2000) presented two subclasses of wave patterns exist. The first one consists of a two-dimensional Stokes wave pattern, and the corresponding waves are called planar waves. They are, 1 ; mX0 fð1Þ m ðxÞ ¼ cos mp x þ 2 1 ; nX0 (12) fð2Þ n ðyÞ ¼ cos npr y þ 2r The second subclass is a three-dimensional wave pattern given by the product of two Stokes waves ð2Þ fð1Þ m ðxÞfn ðyÞ
(13)
Faltinsen also suggested the mixed modes and the combination of the two Stokes modes to represent three-dimensional wave patterns. As stated in the previous sub-section, the major resonant frequency can be identified from the spectral analysis of each type of sloshing wave. The spectral analysis of a square-like wave presents two peaks, and they correspond to the first fundamental frequency and excitation frequency. The combination of o1,0 with 1.5 cos(301) weight and o0,1 with 1.5 sin(301) weight is depicted in Fig. 16(a) and clear the terraced planes can be found in Fig. 16(b). The spectral analysis of irregular waves indicates the resonant peaks of o1,0, o3,0 and o3,2, and the mode shape from the addition of those three resonant modes and the corresponding surface wave contours are depicted in Fig. 16(c) and (d). The odd modes (o1,0, o3,0) vanish when the sloshing comes to the steady state and only o3,2 remains. The wave pattern o3,2 is shown on the top right-hand corner of Fig. 16(c). The effect of odd modes o1,0, o3,0 on an irregular wave during transient state and that of only o3,2 during steady state can be noted from the surface wave contours.
ARTICLE IN PRESS C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
509
Fig. 16. Wave patterns and resonance node (a) square-like waves, y ¼ 301; (b) irregular waves, y ¼ 51.
Fig. 16 provides clear evidence of the correlation between sloshing surface patterns and resonant modes.
5. Conclusion A three-dimensional time-independent finite difference method is developed to study the sloshing waves and resonant modes of fluid in a 3D tank. Only coupled surge–sway motions with various excitation angles under a wide range of excitation frequencies are presented. The following conclusions are reached. 1. The method is used to solve a fully nonlinear sloshing fluid in a three-dimensional tank with square base. The numerical scheme is validated by rigorous benchmark tests. 2. Five types of sloshing waves are observed and they are closely related to the excitation frequency. When the sloshing frequency is close to the first fundamental mode, swirling wave patterns occur. The square-like and irregular wave patterns were also observed in the literature (Faltinsen et al., 2003) when the tank is under near-resonant excitation. In the present study, the square-like wave and irregular wave also can be triggered when the excitation frequency is far away from the first fundamental frequency with an excitation angle of 51. 3. The resonant peaks of spectral analyses were shown to correspond to the primary mode and excitation frequency for the diagonal, single-directional and square-like waves. For the swirling wave, the secondary resonance of mode o2,2 occur in addition to the primary resonance. For irregular waves, the secondary resonance of odd modes o3,0 and o5,0 also appear along with the resonant peaks of the primary mode and the excitation frequency.
4. The surface patterns are generated by the combining of fundamental modes that reflect resonance. The plots of mode shape of the resonant modes and the surface contour of the wave pattern demonstrate a clear correlation between the occurrence of each type of sloshing wave and corresponding resonant modes. 5. The method developed can be used to simulate 3D fluid sloshing in a tank excited by motions including all six-degrees of freedom. Acknowledgement This study is supported by National Science Council under a Grant NSC 96-2221-E-110-107-MY3. Reference Aliabadi, S., Johnson, A., Abedi, J., 2003. Comparison of finite element and pendulum models for simulation of sloshing. Computers and Fluids 23, 535–545. ¨ nal, E., 2005. Experimental investigation of pressure distribution on a Akyildiz, H., U rectangular tank due to the liquid sloshing. Ocean Engineering 32, 1503–1516. ¨ nal, N.E., 2006. Sloshing in a three-dimensional rectangular tank: Akyildiz, H., U numerical simulation and experimental validation. Ocean Engineering 33, 2135–2149. Celebi, M.S., Akyildiz, H., 2002. Nonlinear modeling of liquid sloshing in moving rectangular tank. Ocean Engineering 29, 1527–1553. Chen, B.F., Chiang, S.W., 1999a. Complete 2D and fully nonlinear analysis of ideal fluid in tanks. Journal of Engineering Mechanics 125 (1), 70–78. Chen, B.F., Yuan, Y.S., Lee, J.F., 1999b. Three dimensional nonlinear hydrodynamic pressures by earthquakes on dam faces with arbitrary reservoir shapes. Journal of Hydraulic Research 37, 163–187. Chen, B.F., Nokes, R., 2005. Time-independent finite difference analysis of fully non-linear and viscous fluid sloshing in a rectangular tank. Journal of Computational Physics 209, 47–81.
ARTICLE IN PRESS 510
C.-H. Wu, B.-F. Chen / Ocean Engineering 36 (2009) 500–510
Faltinsen, O.M., 1978. A numerical non-linear method of sloshing in tanks with two dimensional flow. Journal of Ship Research 22 (3), 193–202. Faltinsen, O.M., Rognebakke, O.F., Lukovsky, I.A., Timokha, A.N., 2000. Multidimensional modal analysis of nonlinear sloshing in a rectangular tank with finite water depth. Journal of Fluid Mechanics 407, 201–234. Faltinsen, O.M., Timokha, A.N., 2001. Adaptive multimodal approach to nonlinear sloshing in a rectangular tank. Journal of Fluid Mechanics 432, 167–200. Faltinsen, O.M., Timokha, A.N., 2002. Asymptotic modal approximation of nonlinear resonant sloshing in a rectangular tank with small fluid depth. Journal of Fluid Mechanics 470, 319–357. Faltinsen, O.M., Rognebakke, O.F., Timokha, A.N., 2003. Resonance three-dimensional nonlinear sloshing in a square-base basin. Journal of Fluid Mechanics 487, 1–42. Faltinsen, O.M., Rognebakke, O.F., Timokha, A.N., 2005. Classification of threedimensional nonlinear sloshing in a square-base tank with finite depth. Journal of Fluids and Structures 20, 81–103. Frandsen, J.B., 2004. Sloshing motions in excited tanks. Journal of Computational Physics 196, 53–87. Hill, D.F., 2003. Transient and steady-state amplitudes of forced waves in rectangular basins. Physics of Fluids 15 (6), 1576–1587. Hung, T.K., Wang, M.H., 1987. Nonlinear hydrodynamic pressure on rigid dam motion. Journal of Engineering Mechanics, ASCE 113 (4), 482–499.
Kim, Y., 2001. Numerical simulation of sloshing flows with impact load. Applied Ocean Research 23, 53–62. Kim, Y., Nam, B.W., Kim, D.W., Kim, Y.S., 2007. Study on coupling effects of ship motion and sloshing. Ocean Engineering 34, 2176–2187. Lee, S.J., Kim, M.H., Lee, D.H., Kim, J.W., Kim, Y.H., 2007a. The effects of LNG-tank sloshing on the global motions of LNG carriers. Ocean Engineering 34, 10–20. Lee, D.H., Kim, M.H., Kwon, S.H., Kim, J.W., Lee, Y.B., 2007b. A parametric sensitivity study on LNG tank sloshing loads by numerical simulations. Ocean Engineering 34, 3–9. Miles, J.W., 1994. Faraday waves: rolls versus squares. Journal of Fluid Mechanics 269, 353–371. Ockendon, J.R., Ockendon, H., Waterhouse, D.D., 1996. Multi-mode resonances in fluids. Journal of Fluid Mechanics 315, 317–344. Turnbull, M.S., Borthwick, A.G.L., Eatock Taylor, R., 2003. Numerical wave tank based on a s-transformation finite element inviscid flow solver. International Journal for Numerical Methods in Fluids 42, 641–663. Wu, G.X., Ma, Q.W., Taylor, R.E., 1998. Numerical simulation of sloshing waves in a 3D tank based on a finite element methods. Applied Ocean Research 20, 337–355. Wu, G.X., 2007. Second-order resonance of sloshing in a tank. Ocean Engineering 34, 2345–2349.