Chinese
Adv Search
Home | Accepted | Article In Press | Current Issue | Archive | Special Issues | Collections | Featured Articles | Statistics

2021,  3 (2):   105 - 117   Published Date:2021-4-20

DOI: 10.1016/j.vrih.2020.12.003
1 Introduction2 Related work3 Methods3.1 NB-MultiFLIP3.2 Velocities transfer from particles to grid3.2.1 PIC method for two-phase flow3.2.2 FLIP method for two-phase flow3.2.3 APIC method for two-phase flow3.3 Velocities transfer from grid to particles3.3.1 PIC method for two-phase flow3.3.2 FLIP method for two-phase flow3.3.3 APIC method for two-phase flow3.4 Escaped particle handling4 Results5 Conclusion and future studies

Abstract

Background
The interaction of gas and liquid can produce many interesting phenomena, such as bubbles rising from the bottom of the liquid. The simulation of two-phase fluids is a challenging topic in computer graphics. To animate the interaction of a gas and liquid, MultiFLIP samples the two types of particles, and a Euler grid is used to track the interface of the liquid and gas. However, MultiFLIP uses the fluid implicit particle (FLIP) method to interpolate the velocities of particles into the Euler grid, which suffer from additional noise and instability.
Methods
To solve the problem caused by fluid implicit particles (FLIP), we present a novel velocity transport technique for two individual particles based on the affine particle-in-cell (APIC) method. First, we design a weighed coupling method for interpolating the velocities of liquid and gas particles to the Euler grid such that we can apply the APIC method to the simulation of a two-phase fluid. Second, we introduce a narrowband method to our system because MultiFLIP is a time-consuming approach owing to the large number of particles.
Results
Experiments show that our method is well integrated with the APIC method and provides a visually credible two-phase fluid animation.
Conclusions
The proposed method can successfully handle the simulation of a two-phase fluid.

Content

1 Introduction
The interaction of gases and liquids can produce exciting phenomena, such as bubbles rising from the bottom of the liquid and crown-shaped splashes being created after a droplet drops into a pool. Reproducing such eye-attracting scenarios on a computer is an essential topic in computer graphics. However, this is an arduous task because two fluid phases are introduced into the system. We propose a stable and less noisy algorithm for the simulation of a two-phase fluid.
Many researchers have made significant efforts in this area. Some methods, such as FLIP, introduce only one type of fluid into the system, which cannot animate the complicated interaction of gases and liquids. To address this problem, Boyd and Bridson proposed the MultiFLIP method, which introduces two types of particles, i.e., gas and liquid, into the system[1]. Although MultiFLIP can handle a two-phase fluid simulation, MultiFLIP adopted the FLIP method to transport the velocities between particles and a Euler grid. FLIP is designed to remove the dissipation of the particle-in-cell transport method, but it suffers from noise and instability when simulating a fluid. Jiang et al. presented the affine particle-in-cell (APIC) method to not only remove dissipation but also achieve a less noisy and stable simulation. MultiFLIP is a time-consuming and memory-wasted algorithm because of the large number of particles involved in the system[2]. Based on the assumption that interesting phenomena usually occur at the interface of the gases and liquids, Lyu introduced the narrowband method to a MultiFLIP simulation to reduce the number of particles and achieve a higher performance.
To achieve a less noisy and more stable two-phase fluid simulation, we apply APIC to the MultiFLIP algorithm and obtain the MultiAPIC method which uses the APIC method to transport the velocities between particles and the Euler grid. Unlike one-phase fluid simulation, MultiAPIC maintains two respective velocities on the Euler grid for both gases and liquids. These two respective velocities are merged before the projection step, which is used to obtain a divergence-free velocities field. The merging method handles the influences between gas and liquid. Therefore, the application of the APIC method to the MultiFLIP algorithm depends on the velocities coupling method of the gases and liquids.
Furthermore, our algorithm is implemented using the narrowband method to enhance the performance of MultiAPIC. As the tenet of the narrowband method, all interesting phenomena are presented around the interface of the gases and liquids. We modified the narrowband method such that it is suitable for MultiAPIC because the narrowband method introduced by Ferstl et al. was originally used for the FLIP method[3].
In summary, our contributions are as follows:
(1) We designed a weighed coupling method for interpolating the velocities of liquid and gas particles to the Euler grid so that we can apply the APIC method to the simulation of a two-phase fluid.
(2) We introduced a modified narrowband method to our simulation system because too many particles involved in the system consume both time and memory.
We combine these contributions to obtain a narrowband MultiAPIC (NB-MultiAPIC).
2 Related work
We divided the reviews into three categories. First, we present several different models for a two-phase fluid simulation. We then show the development of transport methods between particles and the Euler grid. Finally, we review the techniques used by narrowband methods.
Two-phase flow. Over the last few decades, researchers have developed methods to simulate two-phase fluids. Popularized by Foster et al., numerical methods for two-phase fluids have received significant attention[4]. The method proposed by Liu et al., which presents a boundary-capturing approach to construct the variable coefficient of a Poisson equation, is a well-known method[5]. This method considers the coefficients and solutions to be discontinuous at the gas-liquid interface and uses the ghost fluid method to capture the boundary conditions. Based on this method, Hong and Kim proposed the use of discontinuous fluids to animate small-scale details of incompressible two-phase fluids by modeling the surface tension of both free and bubble surfaces and considering the discontinuities in the velocity gradient field[6]. Later, the approach was extended by Losasso et al. to simulate a multiphase flow[7]. Kim et al. added a volume control method to preserve the volume of the bubbles[8]. Ando et al. proposed a stream function solver for a two-phase simulation at the price of solving a vector Poisson system instead of the variable Poisson equation, which can obtain the divergence-free velocities directly[9]. Apart from the developments of the above Euler models, some Lagrangian methods have been enhanced to solve the two-phase fluid simulation problem. Smoothed particle hydrodynamics (SPH) has been used to model multiphase flows, such as those developed by Ren et al. and Yang et al.[10,11]. In recent years, several new models based on learning algorithms have been proposed. Ma et al. used neural networks and a direct numerical simulation method to find the closure terms for a simple model of the average flow[12].
Velocity transport. Foster and Metaxas first introduced PIC techniques to computer graphics using liquid simulations[4]. To remove the dissipation, Zhu and Bridson presented a blending method for FLIP and PIC. Subsequently, some researchers proposed several extensions of FLIP to enhance the accuracy of the velocity transport method[13]. Jens et al. coupled high-resolution FLIP with low-resolution implicit SPH[14]. FLIP suffers from more noise and instability when transporting the velocities of particles and the Euler grid. To achieve a stable velocity transport method and prevent this loss of information during the velocity transport step, Jiang et al. introduced the APIC method by allowing for the exact conservation of angular momentum across the transfers between particles and grid[2]. To significantly improve the energy and vorticity conservation, Fu et al. developed an extension of the original APIC, which is a polynomial particle-in-cell method[15]. Boyd and Bridson developed MultiFLIP by coupling the velocities of the gases and liquids to simulate two-phase fluids[1].
Narrowband Method. Ferstl et al. introduced a narrowband method by removing most of the particles to enhance the performance of the fluid simulation[3]. Later, Sato et al. extended the narrowband method and obtained an extended narrowband method, which only maintains the particles in the area where active small-scale details occur[16]. Lyu et al. introduced the narrowband method and extended a narrowband method to the MultiFLIP algorithm to enhance the simulation performance because both gas and liquid particles are involved in the system[17].

Algorithm 1 NB-MultiFLIP vs. NB-MultiAPIC

Require: Steps are colored differently for NB-MultiAPIC and NB-MultiFLIP. Steps in red are unique to NB-MultiAPIC, while steps in blue are unique to NB-MultiFLIP. Other steps, which are not colored, are shared by both methods.

1: Advect the position and velocities of gas and liquid particles

2: Track and construct the gas–liquid interface based on the particles and previous interface

3: Advect the velocities of gas and liquid on the Euler grid

4: Map the velocities of gas and liquid particles to the grid

 [NB-MultiAPIC] Transfer by mass weighted APIC

 [NB-MultiFLIP] Transfer by FLIP

5: Update the velocities on the grid

6: Bump the particles near the interface and resample the particles

7: Handle the escaped particles

8: Add external forces

9: Solve the variable Poisson equation to gain the divergence-free velocities

10: Transport the velocities of the Euler grid to the gas and liquid particles

 [NB-MultiAPIC] Transfer by APIC

 [NB-MultiFLIP] Transfer by FLIP

3 Methods
Our method is an extension of MultiFLIP and inherits its main framework. To better understand the idea of our method, we briefly introduce MultiFLIP before we describe the modified details of MultiAPIC.
3.1 NB-MultiFLIP
In this section, we briefly describe the details of NB-MultiFLIP. NB-MultiFLIP is a two-phase fluid simulation algorithm that introduces FLIP and the narrowband method into the system. To animate a more creative two-phase fluid scene, NB-MultiFLIP adopts an inviscid, incompressible Navier-Stokes equation,
t u + u u = - 1 ρ p + 1 ρ g ,     u = 0 ,
where
u
is the velocity on the Euler grid,
p
is the pressure,
ρ
is the density, and
g
is an external force such as gravity.
The procedure for NB-MultiFLIP is presented in Algorithm 1. The key components of the NB-MultiFLIP algorithm steps are interface tracking, particle bumping, escaped particle handling, and velocity coupling.
(1) Because two types of particles are involved in the system, it is important to construct a distinct interface between the gases and liquids in the interface tracking steps.
(2) Particles near the interface need to be bumped back to their own side to prevent the system from becoming a mixed soup of gas and liquid particles.
(3) Escaped particles are particles that move into the other side (such as gas particles moving into the body of liquids and liquid particles moving into the body of gases) and run far away from the gas-liquid interface. Because these particles represent some small-scale features of fluids, NB-MultiFLIP treats them specially and constructs a small interface around these particles.
(4) NB-MultiFLIP samples two respective velocities on the Euler grid for both a gas and liquid. NB-MultiFLIP couples these two velocities before solving the variable Poisson equation. The coupling method has a significant influence on the results.
Furthermore, NB-MultiFLIP adopts a narrowband method to enhance the performance of the simulation. The basic idea of the narrowband method is to remove the particles far away from the interface. For the domain that is not covered by particles, NB-MultiFLIP uses the Euler grid to sample the properties of the gases and liquids. The key technique of NB-MultiFLIP is coupling the velocities on the Euler grid and the velocities of particles before the divergence-free projection step. More details regarding NB-MultiFLIP can be found in[17].
3.2 Velocities transfer from particles to grid
Our method samples two respective velocities for both a gas and liquid. Because we aim to simulate a two-phase fluid, we need to combine the influence of the gas and liquid. Therefore, our method couples these two respective velocities of a gas and liquid before solving the variable Poisson equation to achieve a divergence-free velocity:
1 ρ p = u * ,
where
u *
is the combined grid velocity field of the gas and liquid velocity fields. The coupling method plays an important role in merging the influence of the gas and liquid. In real life, liquids are more likely to push the gases away when liquids and gases combine together because liquids are heavier than gases. However, the coupling method provided by MultiFLIP treats the velocities of the gas and liquid equally, resulting in the motion of the liquid being restrained by the gas. Therefore, to reduce the influence of the gas on the liquid, we introduce three different types of velocity transfer methods for a two-phase fluid simulation. We adopt different mass weights for the gas and liquid when transferring velocities from the particles to the grid. Furthermore, to enhance the performance of our simulation, we introduced the narrowband method into our simulation system. When coupling velocities of the gas and liquid, we use the velocities on the grid to represent the velocities of the fluid for those grid cells not covered by the particles. In the following, we provide detailed information on discretizing the velocities of the gas and liquid when we use the PIC, FLIP, and APIC methods for velocity transfer.
3.2.1 PIC method for two-phase flow
The particle-in-cell method is a simple and stable velocity transfer method that interpolates the velocities of the particles onto a grid with a kernel function:
m a i * = p ω a i p m p ,
m a i * u ̂ a i = p ω a i p m p u a p ,
where
ω a i p = N ( x p - x a i )
is the trilinear interpolation weight,
m p
is the mass of a particle,
x p
is the position of a particle,
x a i
is the location of the center of the grid cell, and
u a p
denotes the velocity of a particle in direction
a
. In our experiment, we used 1 for the mass of the gas particles and 1000 for the mass of the liquid particles. Equations (3) and (4) are used to compute the mass and momentum of the grid faces, respectively. Using these two quantities, we can obtain the interpolated velocities through
u ̂ a i = m a i * u ̂ a i / m a i *
.
Unlike the original narrowband method, which combines the velocities interpolated from particles and velocities on the grid based on a distance function
ϕ
, our method combines these two velocities based on the number of particles using:
m a i = m a x 0 ,    ρ f r a c a i V -   m a i * + m a i * ,
m a i u a i * = m a x 0 ,    ρ f r a c a i V -   m a i * u g a i + m a i * u ̂ a i ,
where
V = x 3
is the volume of a grid cell,
f r a c a i
is the fraction of the grid cell face covered by the liquid, and
u g a i
is the component of advected velocity in direction
a
. We used 1 for the density
ρ
of gas and 1000 for the density of the liquid. Finally, we obtain the velocity for the divergence-free solver through
u a i * = m a i u a i * / m a i
.
3.2.2 FLIP method for two-phase flow
The FLIP method was designed to remove the dissipation of the PIC method. Unlike PICs, it does not store the external properties of the particles. Therefore, we used the same method as the PIC to transfer velocities from the particles to the grid.
3.2.3 APIC method for two-phase flow
The PIC loses information during a velocity transfer because it dissipates a significant amount of angular momentum. Based on the velocity transfer method introduced by MultiFLIP, the velocities are transferred from the particles to the Euler grid through the following:
N a i = p ω a i p ,
N a i u ̂ a i * = p ω a i p ( u a p + ( c p a ) T ( x a i - x p ) )
u ˙ a i * = u g a i        ϕ > k × Δ x u ̂ a i *         ϕ k × Δ x ,
where
c p a
is a vector stored per particle in direction
a
, which represents the angular momentum of a particle. Before solving the projection equation, the velocities of the gas and liquid are coupled using
u a i * = f r a c a i × u ˙ a i _ l * + ( 1 - f r a c a i ) × u ˙ a i _ a * ,
where
u ˙ a i _ l *
and
u ˙ a i _ a *
are the transferred velocities on the Euler grid.
However, the liquid may flow in an abnormal path without considering the weights of the gas and liquid. To solve this problem, APIC also interpolates the angular momentum from the particles to the grid when calculating the momentum of the grid faces as follows:
m a i = m a x 0 ,    ρ f r a c a i V -   m a i * + p ω a i p m p ,
m a i u a i * = m a x 0 ,    ρ f r a c a i V -   m a i * u g a i + p ω a i p m p ( u a p + ( c p a ) T ( x a i - x p ) ) .
Finally, we obtain the velocity for the divergence-free solver through
u a i * = m a i u a i * / m a i
. With the angular momentum, our method can not only remove the dissipation of the fluid but also provides a stable and less noisy two-phase fluid simulation algorithm.
3.3 Velocities transfer from grid to particles
We used particles to track the details of the gas-liquid interface. After solving the variable Poisson equation, we need to transfer the divergence-free velocities from the grid to the particles. Similar to the transfer from particles to the grid, we will provide details on how to discretize the velocity transfer using PIC, FLIP, and APIC methods.
3.3.1 PIC method for two-phase flow
The PIC method simply interpolates the velocities on the grid into particles using:
u a p = i ω a i p u a i .
This simple handling of velocities leads to fluid dissipation. The PIC method restrains the motion of both the gas and liquid. Thus, few researchers have used it to transfer velocities when simulating a fluid.
3.3.2 FLIP method for two-phase flow
The FLIP method for a two-phase flow is based on the original FLIP method introduced by Zhu and Bridson. The velocities of the particles are updated using:
u a p = u a p _ i n f + u i n t e r p   - u i n t e r p * ,
where
u a p _ i n f
is the final velocity of a particle in the previous step;
u i n t e r p
is the interpolated velocity of the divergence-free velocity on the grid, which is the result of the variable Poisson equation solver; and
u i n t e r p *
is the interpolated velocity after particle bumping.
3.3.3 APIC method for two-phase flow
Apart from the velocity, the APIC method also stores the external information of the particles to conserve the angular momentum through the following:
u a p = i ω a i p u a i ,
c p a = i ω a i p u a i ,
where
c p a
is a vector in direction
a
.
After solving the divergence-free projection step, a coupled divergence-free velocity field is obtained. The MultiFLIP method suggests that the velocities are updated on the grid faces whose liquid (gas) fraction is not equal to zero. However, the handling method of MultiFLIP smears the influence of the gas into the liquid, which causes the liquid to splash too many particles, which move along an abnormal path, as shown in Figure 7. To obtain a better result, our method updates the gas velocities on the grid faces that are covered by gas (
f r a c a i < 0.5 )
and the liquid velocities on the grid faces that are covered by liquid (
f r a c a i 0.5 )
.
After updating the velocities on the grid, our method extrapolates the liquid (gas) velocities to the grid cells that are not covered by liquid (gas). Finally, we transfer the velocities from the Euler grid to the particles.
3.4 Escaped particle handling
Escaped particles run far away from the body of the fluid. These particles represent small-scale details of the fluid. MultiFLIP constructs a small gas-liquid interface around these particles so that no external operators are adopted to handle the escaped particles. However, the escaped particles of the liquid are significantly influenced by the gas. Therefore, we treat the escaped particles of the liquid as ballistic particles. When we use the APIC method to transfer the velocities between the particles and the grid, we only store the velocities of the liquid-escaped particles. Although the angular momentum of the liquid-escaped particles is not recorded, the results are acceptable because the number of escaped particles is insufficient to have much influence.
4 Results
We ran our experiments on a desktop PC equipped with 3.2 GHz i5-3470 processors and 8 GB of memory. The experiments were divided into three groups. The first, second, and third groups are simulated using NB-MultiFLIP, NB-MultiPIC, and NB-MultiAPIC, respectively. All the groups were simulated using a narrowband method. The results of all the experiments, including the number of particles and runtime, can be found in Table 1.
Experimental results including number of particles and runtime
Method Number of Particles (Thousand) Avg. Time/Timestep (s)
G a s
L i q u i d
T o t a l
Tension NB-MultiPIC 133.45 101.97 235.42 2.22
NB-MultiFLIP 133.31 102.01 235.33 2.26
NB-MultiAPIC 126.98 126.03 253.01 2.23
Drop NB-MultiPIC 539.50 509.24 1048.73 11.04
NB-MultiFLIP 562.15 550.19 1112.34 13.20
NB-MultiAPIC 528.33 516.71 1045.04 12.54
Box glugging NB-MultiPIC 98.96 83.61 182.57 1.65
NB-MultiFLIP 112.88 101.66 214.54 2.10
NB-MultiAPIC 121.27 104.69 225.96 2.20
Tube NB-MultiPIC 232.36 230.20 462.39 4.75
NB-MultiFLIP 271.23 255.65 526.88 5.13
NB-MultiAPIC 280.55 273.57 554.12 6.05
Dam break NB-MultiPIC 250.15 234.90 485.05 4.41
NB-MultiFLIP 272.57 253.65 526.21 5.87
NB-MultiAPIC 261.21 241.19 502.39 5.00
Bubble NB-MultiPIC 135.50 145.99 281.49 3.89
NB-MultiFLIP 140.66 159.68 300.34 4.28
NB-MultiAPIC 139.01 155.05 294.06 4.02
Figure 1 presents a cube-shaped liquid oscillating in a zero-gravity zone. The results show that NB-MultiAPIC and NB-MultiPIC can simulate two-phase fluid scenes similar to NB-MultiFLIP. The scenes simulated using these three methods are nearly the same. NB-MultiPIC can also preserve the energy of the liquid as in other methods when there is a slight rotation of the scene. NB-MultiFLIP does not generate more noise than other methods in this case.
Figure 2 shows a droplet dropping into a water pool and creating a crown on the gas-liquid interface. As shown in Figure 2, NB-MultiPIC restrains the motion of water and cannot form an intact crown after the droplet drops into the water pool. Although NB-MultiFLIP does not restrain the motion of the liquid, NB-MultiFLIP has more noise than the other transfer methods. Figure 3 shows the number of particles escaped during the simulation of a water droplet. NB-MultiFLIP splashes more escaped particles than our method. As time passes, the number of escaped particles decreases. However, after frame 250, the amplitude of the blue line does not decay, which means that MultiFLIP may produce more noise than our method. NB-MultiAPIC does not splash extra-escaped particles after frame 250, resulting in a normal liquid simulation similar to that of a liquid found in the real world.
Figure 4 and Figure 5 show a dam break with obstacles and a bubble rising in water with obstacles, respectively. The results show that our method can handle the interaction between the fluid and obstacles. NB-MultiFLIP and NB-MultiAPIC splash more particles than NB-MultiPIC. This result meets our expectation that the NB-MultiPIC restrains the motion of the liquid.
Figure 6 and Figure 7 demonstrate the liquid flow from a box container to another through tubes in the horizontal and vertical directions. Although the movement of the liquid is intended when the liquid is simulated using NB-MultiAPIC, more escaped particles run away from the body of the liquid when simulating a fluid using the FLIP transfer method. This is because NB-MultiFLIP is noisier than the other approaches.
Figure 8 shows liquid flowing from one box container to another through tubes placed in the vertical direction with different weights when coupling the velocities of the gas and liquid. The demonstration shows that the particles of the liquid are carried by particles of gas when the weights of the liquids are small. The smaller the weight of the liquid, the more chaotic the system.
Figure 9 shows the velocity transfer method from the Euler grid to the particles of MultiFLIP and our method. The first two subfigures are the results simulated using MultiFLIP, which indicates that the velocity transfer method suggested by MultiFLIP smears the influence of gas into the liquid and makes the particles of the liquid flow along the path of the gas. In addition, NB-MultiAPIC can solve this problem. NB-MulitAPIC does not generate too many abnormal escaped particles and results in a more accurate scene of the liquid and gas interaction.
In summary, we compared the performances of three different velocity transfer methods. As shown in Table 1, there are no obvious differences between the runtime of these three transfer methods. NB-MultiPIC restrains the motion of the liquids and suffers from dissipation. Although NB-MultFLIP can preserve the energy of the liquid, NB-MultiFLIP is the noisiest of these transfer methods and may splash too many escaped particles carried by the gas and flow along abnormal paths. NB-MultiAPIC not only maintains the kinetic energy of the liquid it also remains stable. NB-MultiAPIC can generate more accurate results than the other methods without reducing the efficiency of the simulation system.
5 Conclusion and future studies
We introduced a technique to apply the APIC method to MultiFLIP. The key components of the technique are the coupling method between the gas and liquid, as well as the coupling of particles and the grid based on a narrowband method. In addition, we used ballistic particles to represent the small-scale details of the liquid. The kinetic energy of the escaped particles is absorbed by the main body of the liquid after the escaped particles return to the main body of the liquid. The experiments show that NB-MultiAPIC can not only remove dissipation but also has less noise. NB-MultiPIC suffers from loss of kinetic energy when simulating a two-phase fluid, and NB-MultiFLIP is an unstable method and has more noise than the other approaches.
Although our method can handle the simulation of a two-phase fluid, there may be more than two phases in daily life, and our method cannot reproduce scenes that involve more than two phases of flow. Therefore, in a future study, we would like to extend the MultiAPIC/MultiFLIP method to simulate a multiphase flow.

Reference

1.

Boyd L, Bridson R. MultiFLIP for energetic two-phase fluid simulation. ACM Transactions on Graphics, 2012, 31(2): 1–12 DOI:10.1145/2159516.2159522

2.

Jiang C, Schroeder C, Selle A, Teran J, Stomakhin A. The affine particle-in-cell method. ACM Transactions on Graphics, 2015, 34(4): 1–10 DOI:10.1145/2766996

3.

Ferstl F, Ando R, Wojtan C, Westermann R, Thuerey N. Narrow band FLIP for liquid simulations. Computer Graphics Forum, 2016, 35(2): 225–232 DOI:10.1111/cgf.12825

4.

Foster N, Metaxas D. Realistic animation of liquids. Graphical Models and Image Processing, 1996, 58(5): 471–483 DOI:10.1006/gmip.1996.0039

5.

Liu X D, Fedkiw R P, Kang M. A boundary condition capturing method for poisson's equation on irregular domains. Journal of Computational Physics, 2000, 160(1): 151–178 DOI:10.1006/jcph.2000.6444

6.

Hong J M, Kim C H. Discontinuous fluids. ACM Transactions on Graphics, 2005, 24(3): 915–920 DOI:10.1145/1073204.1073283

7.

Losasso F, Shinar T, Selle A, Fedkiw R. Multiple interacting liquids. ACM Transactions on Graphics, 2006, 25(3): 812–819 DOI:10.1145/1141911.1141960

8.

Kim B, Liu Y J, Llamas I, Jiao X M, Rossignac J. Simulation of bubbles in foam with the volume control method. ACM Transactions on Graphics, 2007, 26(3): 98 DOI:10.1145/1276377.1276500

9.

Ando R, Thuerey N, Wojtan C. A stream function solver for liquid simulations. ACM Transactions on Graphics, 2015, 34(4): 1–9 DOI:10.1145/2766935

10.

Ren B, Li C F, Yan X, Lin M C, Bonet J, Hu S M. Multiple-fluid SPH simulation using a mixture model. ACM Transactions on Graphics, 2014, 33(5): 1–11 DOI:10.1145/2645703

11.

Yang T, Chang J, Ren B, Lin M C, Zhang J J, Hu S M. Fast multiple-fluid simulation using Helmholtz free energy. ACM Transactions on Graphics, 2015, 34(6): 1–11 DOI:10.1145/2816795.2818117

12.

Ma M, Lu J C, Tryggvason G. Using statistical learning to close two-fluid multiphase flow equations for a simple bubbly system. Physics of Fluids, 2015, 27(9): 092101 DOI:10.1063/1.4930004

13.

Edwards E, Bridson R. A high-order accurate particle-in-cell method. International Journal for Numerical Methods in Engineering, 2012, 90(9): 1073–1088 DOI:10.1002/nme.3356

14.

Cornelis J, Ihmsen M, Peer A, Teschner M. IISPH-FLIP for incompressible fluids. Computer Graphics Forum, 2014, 33(2): 255–262 DOI:10.1111/cgf.12324

15.

Fu C Y, Guo Q, Gast T, Jiang C, Teran J. A polynomial particle-in-cell method. ACM Transactions on Graphics, 2017, 36(6): 1–2 DOI:10.1145/3130800.3130878

16.

Sato T, Wojtan C, Thuerey N, Igarashi T, Ando R. Extended narrow band FLIP for liquid simulations. Computer Graphics Forum, 2018, 37(2): 169–177 DOI:10.1111/cgf.13351

17.

Lyu L, Ren X H, Cao W, Zhu J, Wu E H. Adaptive narrow band MultiFLIP for efficient two-phase liquid simulation. Science China Information Sciences, 2018, 61(11): 114101 DOI:10.1007/s11432-018-9518-3