ZHANG Yu-xin (张雨新), WAN De-cheng (万德成)
State Key laboratory of Ocean Engineering, School of Naval Architecture, Ocean and Civil Engineering,
Shanghai Jiao Tong University, Shanghai, China, E-mail: iceyuxin@163.com
HINO Takanori
Graduate school of Engineering, Yokohama National University, Yokohama, Japan
Comparative study of MPS method and level-set method for sloshing flows*
ZHANG Yu-xin (张雨新), WAN De-cheng (万德成)
State Key laboratory of Ocean Engineering, School of Naval Architecture, Ocean and Civil Engineering,
Shanghai Jiao Tong University, Shanghai, China, E-mail: iceyuxin@163.com
HINO Takanori
Graduate school of Engineering, Yokohama National University, Yokohama, Japan
(Received January 6, 2014, Revised March 12, 2014)
This paper presents a comparative study of a meshless moving particle semi-implicit (MPS) method and a grid based level-set method in the simulation of sloshing flows. The numerical schemes of the MPS and level-set methods are outlined and two violent sloshing cases are considered. The computed results are compared with the corresponding experimental data for validation. The impact pressure and the deformations of free surface induced by sloshing are comparatively analyzed, and are in good agreement with experimental ones. Results show that both the MPS and level-set methods are good tools for simulation of violent sloshing flows. However, the second pressure peaks as well as breaking and splashing of free surface by the MPS method are captured better than by the level-set method.
moving particle semi-implicit (MPS) method, level-set method, liquid sloshing, impact pressure, free surface
Liquid sloshing is a kind of nonlinear free surface flows, and usually involves some complicated phenomena, such as breaking wave, overturning of free surfaces and splashing liquid. The impact pressure induced by sloshing may damage the structure of the tank, and affect the safety of ship with liquid cargo interacting with waves. Therefore, accurate prediction of impact pressure induced by liquid sloshing is of significant importance in ocean engineering.
CFD has been proved to be an effective tool to compute sloshing flows in the last decades[1-6]. The methods adopted in CFD can be generally divided into two categories: the grid-based method and the meshless method. In the framework of grid-based method, the discretizations of governing equations are implemented based on grid system, and the deformation of free surface is usually described by an extra function which is updated in every time step, such as volume of fluid (VOF)[7,8]and level-set[9,10]methods. Nevertheless, there is no need for grid in using the meshless method, in which fluid is represented by a set of interacting particles. These particles have physical properties, such as mass, momentum, and energy, etc. Motions of particles are described by the Lagrangian approach. Due to the meshless character and Lagrangian nature, the meshless method has two great advantages: (1) fragment and coalescence of fluid can be computed straightforwardly since there is no constant topology relationship between particles, (2) numerical diffusion in the discretization of convection term is eliminated by use of substantial derivative in governing equations. However, the meshless method is of high computation cost and usually suffers from unphysical pressure oscillation. There are both the pros and cons between the grid-based method and the meshless method.
The aim of the present paper is to make comparative investigation on two commonly used methods: a meshless method, the moving particle semi-implicit (MPS) method[11,12], and a grid-based method, the level-set method, in a sloshing context. The present MPS solver is developed at Shanghai Jiao Tong University, which employs four improved numerical schemes to overcome the unphysical pressure oscillation commonly observed in traditional MPS method, such as: (1) momentum conservative gradient model[13], (2) kernel function without singularity[12], (3) mixed source term for pressure Poisson equation (PPE)[13,14], (4) accurate free surface detection method[15]. The effect of these improved schemes on the prediction of impact pressure in sloshing will be analyzed by comparing the results by the present MPS method with that by traditional MPS method. On the other hand, the present level-set method is incorporated in a gridbased solver, i.e. SURF (Solution algorithm for Unstructured RANS with FVM)[16], which is an unstructured finite volume method (FVM) solver dedicated to simulate free surface flow in marine hydrodynamics. The code has been developed since 1994 at the CFD Research Center, National Maritime Research Institute (NMRI) in Japan. By use of unstructured grid, SURF can deal with complex geometry. A single phase level-set is adopted in SURF to make it more efficiently solve the flow field around ship[17]. This code has been commercialized by NMRI themselves and used in many shipyards in Japan.
In the following sections, the numerical schemes of the MPS and level-set methods are first outlined, then comparison analysis are carried out against liquid sloshing problems. The calculated pressures by traditional MPS method and present MPS method are compared to show the improvement on pressure computation. In addition, both MPS and level-set results are compared with experiment to validate both solvers. Discrepancy of the numerical results by the MPS and level-set methods is discussed.
1.1 Governing equations
For incompressible and viscous fluid, governing equations including mass and momentum conservation equations are represented as
whereD/Dt denotes substantial derivative,t is the time,V the velocity vector,ρthe density,pthe pressure,νthe kinematic viscosity, andgthe gravity acceleration. Thanks to Lagrangian nature, convection terms do not turn up in the momentum equations of MPS calculation, thus avoiding numerical diffusion in the discretization of convection terms.
1.2 Particle models
Governing equations are transformed to the equations of particle interactions based on so called particle models, namely the gradient model and Laplace model.
In original MPS method, the gradient operator is modeled as a local weighted average of the gradient vectors between particlei and its neighboring particle j , given as[11]
whereφis an arbitrary scalar function,Dthe number of space dimensions,rthe coordinate vector, W(|r-r |)the kernel function, and n0the initial
ji particle number density.
Equation 3 is not conservative from momentum point of view since two isolated neighbor particles with different pressures will be accelerated to infinity. A conservative form can be obtained by adding2φiinto its right side[13]
The kernel function in Eq.(6) has similar shape with the one in Eq.(5), but without singularity at r=0.
In the MPS method, the density of fluid is described by the so-called particle number densityn, defined as[11]
The Laplacian operator is derived by Koshizuka[11]from the physical concept of diffusion, expressed as:
where λis a parameter, introduced to keep the variance increase equal to that of the analytical solution. Both viscous force∇2Vin Eq.(2) and ∇2Pon the right hand side of the Poisson pressure equation (i.e., Eq.(10) and Eq.(11)) are discretized by using Eq.(8).
1.3 Model of incompressibility
The incompressible condition in traditional MPS method is represented by keeping the particle number density constant. Thus the Poisson equation of pressure (PPE) is expressed as[11]
wheren∗is the particle number density in the temporal field.
The source term in Eq.(10) is solely based on the deviation of the temporal particle number density from the initial value. Pressure obtained from Eq.(10) is prone to oscillate in spatial and temporal domains since the particle number density field is not smooth. To overcome this, Tanaka and Masunaga[13]proposed a mixed source term for PPE, which combines the velocity divergence and the particle number density together. The main part of the mixed source term is the velocity divergence, while the particle number density is used to keep fluid volume constant. This improved PPE is rewritten by Lee[14]as
whereγis a blending parameter with a value between 0 and 1. The range of 0.01≤γ≤0.05is better according to numerical tests conducted by Lee[14], and smallerγimplies smoother pressure field. Therefore, in this paper,γ=0.01 is used for all simulations.
1.4 Free surface detection
Interaction domain is truncated near free surface since there is no particle out of domain for single phase calculation. Thus particles on the free surface have smaller particle number density than those under the free surfaces. Consider this, in traditional MPS method, particles are defined as surface particles when they satisfy[11]
whereβis a parameter with a value between 0.80 and 0.99. For the surface particles, zero pressure condition is imposed. Therefore, the detection of surface particle has significant effect on the pressure field.
Equation (12) has a low accuracy since inner particles with small particle number density may be misjudged as surface particles. Therefore, unreal pressurearound the misjudged particles occurs. This usually causes unphysical pressure oscillation. To avoid this problem, an improved detection method is employed in this paper. This improved detection method defines a vector function as follows[15]
Particles would be judged as surface particles when the magnitudes of theirF satisfy
where αis a parameter, and has a value of 0.9 in this paper,F0is the initial value ofFfor surface particle, also equal to theFof boundary particles without neighboring fluid particle.
The vector Frepresents the asymmetry of neighboring particles. It has a small value at the core of the domain and a large value near the free surface. It should be pointed out that Eq.(14) is not valid for splashed particle which has few neighboring particles, so it is only used for particles with number density between0.8n0and 0.97n0. Particles with number density lower than0.8n0should be surface particles, while those with number density higher than0.97n0should get pressure through the Poisson equation.
1.5 Level-set method in SURF solver
In the level set method, the shape of free surface is defined by a mathematical function φ(x,y,z,t), which is a signed distance from the interface, i.e.,
According to kinematic condition, the fluid particles on a free surface remain on the free surface. Thus the following mathematical form is obtained
1.6 The localized level set method
Thus, the level set function is localized within the bandwidth 2γ2from the interface. The transport equation (Eq.(16)) is modified as
where C(φ)is the cut-off function defined as
In computation the update ofφis performed only in the region of |φ|≤γ2.
1.7 Re-initialization
As the level set function is no longer a distance function after the convection, a re-initialization process is necessary in level set computation. In SURF, the re-initialization process is conducted using the partial differential equation as
where τis the pseudo time,φ0the initial value and S a sign function, defined as
As the solution of Eq.(20) converges,φbecomes the distance function again since∇φ=1.
In a numerical process,S(φ)is approximated as
whereεis a typical grid spacing.
Thus Eq.(20) is rewritten as
2.1 Computational conditions
2.2 Computational results
Before discussing the results of the MPS and level set methods, comparison is first carried out between traditional MPS and the present MPS to show the improvement on smoothing pressure field. Here, traditional MPS employs Eqs.(3), (5), (10) and (12) in computation, while Eqs.(4), (6), (11) and (14) are used in the present MPS computation. Figure 2 shows the pressure time histories at position P2 obtained by experiment, traditional MPS and the present MPS computations. Although the traditional MPS is able to predict the impact phenomena, it cannot predict the value of impact pressure well. Strong oscillations with high frequency are observed in traditional MPS results. However, the pressure curves in the present MPS results are much smoother, and thus seems more reasonable. The improvement on pressure prediction is evident.
Figure 3 compares pressure time histories in Case 1 between experiment, MPS and level set results. It is seen that the liquid impact behavior is well predicted in numerical simulations. The general evolutions of impact pressures obtained numerically by the MPS and level set methods agree with experimental results. Two pressure peaks are observed in each periodicimpact behavior. The first one has a larger value and a shorter duration. It is due to the impact of liquid on the side wall of the tank with a large horizontal velocity. After the first peak, liquid goes up along the wall, and finally falls down due to the gravity. Second peak occurs when the falling liquid hits the underlying liquid. The first pressure peak is overestimated by both the MPS and level-set methods compared with experimental data. The reason for this is that breaking wave occurs when liquid impacts on the side wall (see Fig.4), air is pocketed by liquid and wall. However, the effect of existence of air is not taken into account in single phase simulation, thus higher pressure are caused when the air bubble vanishes. Two phase calculation may simulate such phenomena better. However, this part of work is not included in this paper. There is a pressure oscillation after the first pressure peak in the level set results, which might be due to numerical error involved in treatment of wall boundary. Further improvement on this will be desired, however, it does not affect the general profile of pressure curves.
Figure 5 shows the results of Case 2. Again, the agreements between experiment and numerical simulations are satisfactory. In this case, flow is also violent. Breaking wave is observed. However, it occurs earlier than it did in Case 1, see Fig.6, in which the wave breaks before it reach the side wall. The results in that the values of first pressure peaks at positions P1 and P3 are relatively smaller than the ones in Case 1, and also implies that excitation period has significant effect on the wave in sloshing flow. The results of MPS and level set are compared in Fig.5, and it is seen that the MPS simulation shows a better results since the level set simulation overestimates the first pressure peak at positions P1 and P3, and underestimates the second peak at position P1.
A detailed comparison between MPS and levelset is shown in Fig.7, which illustrates the pressure evolutions in a single impact behavior at position P1 for the two cases. MPS shows a better prediction on the moment and the value of second pressure peak than level-set compared with experimental data. The reason for this might be that the wave breaks after impacting on the side walls of the tank, and MPS is able to compute the breaking free surfaces better than level-set, thus shows a good prediction of the subsequent flow which is related to the second pressure peak.
Figures 8 and 9 show some snapshots in experiment and simulations for Case 1 and Case 2, respectively. It is seen that the flows are violent. Breakingfree surface and splashing liquid are observed when liquid impacts on the side wall of the tank. Compared with experiment, both MPS and level set predict the profiles of free surface well. The general motion of liquid is accurately captured in computation. However, particle method is more flexible to describe the liquid splashing phenomena, and thus the free surfaces in MPS results behave more naturally.
In this paper, comparative study of the MPS method and level-set method in the simulation of sloshing flows was presented. Two violent sloshing flows in 2-D rectangular tank subjected to sway motion were computed. The impact pressures obtained numerically show two peaks during each periodic impact behavior, and are in agreement with experimental obsercations. However, the second peaks by the MPS method are predicted better than by the level-set method. The computed free surfaces agree well with experimental ones, suggesting that both the MPS and the level-set methods are able to simulate the large deformation of free surface in violent sloshing flows. Nevertheless, due to the Lagrangian nature, the MPS method can simulate the breaking wave and the splashing liquid better than the level-set method do.
Acknowledgement
The support of the Center of High Performance Computing (HPC) of Shanghai Jiao Tong University for this work is gratefully acknowledged.
[1] KIM Y. Numerical simulation of sloshing flows with impact load[J]. Applied Ocean Research, 2001, 23(1): 53-62.
[2] LIU D., LIN P. A numerical study of three-dimensional liquid sloshing in tanks[J]. Journal of Computational Physics, 2008, 227(8): 3921-3939.
[3] PAN Xu-jie, ZHANG Huai-xin and LU Yun-tao. Numerical simulation of viscous liquid sloshing by moving-particle semi-implicit method[J]. Journal of Marine Science and Application, 2008, 7(3): 184-189.
[4] HU Chang-hong, YANG Kyung-Kyu and KIM Yonghwan. 3-D numerical simulations of violent sloshing by CIP-based method[J]. Journal of Hydrodynamics, 2010, 22(5Suppl.): 253-258.
[5] LI Yu-long, ZHU Ren-chuan and MIAO Guo-ping et al. Simulation of tank sloshing based on OpenFOAM and coupling with ship motions in time domain[J]. Journal of Hydrodynamics, 2012, 24(3): 450-457.
[6] WU C. H., FALTINSEN O. M. and CHEN B. F. Numerical study of sloshing liquid in tanks with baffles by time-independent finite difference and fictitious cell method[J]. Computers and Fluids, 2012, 63: 9-26.
[7] LIU D., LIN P. Three-dimensional liquid sloshing in a tank with baffles[J]. Ocean Engineering, 2009, 36(2): 202-212.
[8] MING Ping-jian, DUAN Wen-yang. Numerical simulation of sloshing in rectangular tank with VOF based on unstructured grids[J]. Journal of Hydrodynamics, 2010, 22(6): 856-864.
[9] YU K. Level-set RANS method for sloshing and green water simulations[D]. Doctoral Thesis, College Station, Texas, USA: Texas A&M University, 2008.
[10] CHEN Y. G., DJIDJELI K. and PRICE W. G. Numerical simulation of liquid sloshing phenomena in partially filled containers[J]. Computers and Fluids, 2009, 38(4): 830-842.
[11] KOSHIZUKA S., NOBE A. and OKA Y. Numerical analysis of breaking waves using the moving particle semi-implicit method[J]. International Journal for Numerical Methods in Fluids, 1998, 26(7): 751-769.
[12] ZHANG Yu-xin, WAN De-cheng. Numerical Simulation of liquid sloshing in low-filling tank by MPS[J]. Chinese Journal of Hydrodynamics, 2012, 27(1): 100-107(in Chinese).
[13] TANAKA M., MASUNAGA T. Stabilization and smoothing of pressure in MPS method by quasi-compressibility[J]. Journal of Computational Physics, 2010, 229(11): 4279-4290.
[14] LEE B. H., PARK J. C. and KIM M. H. et al. Stepby-step improvement of MPS method in simulating violent free-surface motions and impact-loads[J]. Computer Methods in Applied Mechanics and Engineering, 2011, 200(9-12): 1113-1125.
[15] ZHANG Yu-xin, WAN De-cheng. Application of improved MPS method in sloshing problem[C]. Proceedings of the 23rd National Congress on Hydrodynamics and 10th National Congress on Hydrodynamics. Xi’an, China, 2011, 156-162(in Chinese).
[16] WACKERS J., KOREN B. and RAVEN H. et al. Freesurface viscous flow solution methods for ship hydrodynamics[J]. Archives of Computational Methods in Engineering, 2011, 18(1): 1-41.
[17] HINATSU M., HINO T. Computation of viscous flows around a wigley hull running in incident waves by use of unstructured grid method[C]. Proceedings of the 20th International Offshore and Polar Engineering Conference. Kitakyushu, Japan, 2002.
[18] PENG D., MERRIMAN B. and OSHER S. et al. A PDE-based fast local level set method[J]. Journal of Computational Physics, 1999, 155(2): 410-438.
[19] HINATSU M., TSUKADA Y. and FUKUSAWA R. et al. Experiments of two-phase flows for the joint research[C]. Proceedings of SRI-TUHH Mini-Workshop on Numerical Simulation of Two-Phase Flows. Tokyo, Japan, 2001.
[20] COLAGROSSI A., LUGNI C. and GRECO M. et al. Experimental and numerical investigation of 2D sloshing with slamming[C]. Proceedings of 19th International Workshop on Water Waves and Floating Bodies. Cortona, Italy, 2004.
[21] FALTINSEN O. M., ROGNEBAKKE O. F. and TIMOKHA A. N. Resonant three-dimensional nonlinear sloshing in a square-base basin[J]. Journal of Fluid Mechanics, 2003, 487: 1-42.
10.1016/S1001-6058(14)60065-2
* Project supported by the National Natural Science Foundation of China (Grant Nos. 51379125, 51411130131 and 11272120), the National Key Basic Research Development of China (973 Program, Grant No. 2013CB036103), the High Technology of Marine Research Project of the Ministry of Industry and the Information Technology of China and the Program for Professor of Special Appointment (Eastern Scholar) at Shanghai Institutions of Higher Learning (Grant No. 2013022).
Biography: ZHANG Yu-xin (1985-), Male, Ph. D. Candidate
WAN De-cheng,
E-mail: dcwan@sjtu.edu.cn