In numerical linear stability investigations, the rates of change of the kinetic and thermal energy of the perturbation flow are often used to identify the dominant mechanisms by which kinetic or thermal energy is exchanged between the basic and the perturbation flow. Extending the conventional energy analysis for a single-phase Boussinesq fluid, the energy budgets of arbitrary infinitesimal perturbations to the basic two-phase liquid–gas flow are derived for an axisymmetric thermocapillary bridge when the material parameters in both phases depend on the temperature. This allows identifying individual transport terms and assessing their contributions to the instability if the basic flow and the critical mode are evaluated at criticality. The full closed-form energy budgets of linear modes have been derived for thermocapillary two-phase flow taking into account the temperature dependence of all thermophysical parameters. The influence of different approximations to the temperature dependence on the linear stability boundary of the axisymmetric flow in thermocapillary liquid bridges is tested regarding their accuracy. The general mechanism of symmetry breaking turns out to be very robust.

In lubrication problems, which concern thin film flow, cavitation has been considered as a fundamental element to correctly describe the characteristics of lubricated mechanisms. Here, the well-posedness of a cavitation model that can explain the interaction between viscous effects and micro-bubbles of gas is studied. This cavitation model consists of a coupled problem between the compressible Reynolds partial differential equation (PDE) (that describes the flow) and the Rayleigh–Plesset ordinary differential equation (that describes micro-bubbles evolution). A simplified form without bubbles convection is studied here. This coupled model seems never to be studied before from its mathematical aspects. Local times existence results are proved and stability theorems are obtained based on the continuity of the spectrum for bounded linear operators. Numerical results are presented to illustrate these theoretical results.

We investigate the Cauchy problem of the viscous liquid-gas two-phase flow model in ℝ3. Under the assumption that the initial data is close to the constant equilibrium state in the framework of Sobolev space H2(ℝ3), the Cauchy problem is shown to be globally well-posed by an energy method. If additionally, for 1 ⩽ p < 6/5, Lp-norm of the initial perturbation is bounded, the optimal convergence rates of the solutions in Lq-norm with 2 ⩽ q ⩽ 6 and optimal convergence rates of their spatial derivatives in L2-norm are also obtained by combining spectral analysis with energy methods.

This paper presents topology optimization of capillary, the typical two-phase flow with immiscible fluids, where the level set method and diffuse-interface model are combined to implement the proposed method. The two-phase flow is described by the diffuse-interface model with essential no slip condition imposed on the wall, where the singularity at the contact line is regularized by the molecular diffusion at the interface between two immiscible fluids. The level set method is utilized to express the fluid and solid phases in the flows and the wall energy at the implicit fluid-solid interface. Based on the variational procedure for the total free energy of two-phase flow, the Cahn-Hilliard equations for the diffuse-interface model are modified for the two-phase flow with implicit boundary expressed by the level set method. Then the topology optimization problem for the two-phase flow is constructed for the cost functional with general formulation. The sensitivity analysis is implemented by using the continuous adjoint method. The level set function is evolved by solving the Hamilton-Jacobian equation, and numerical test is carried out for capillary to demonstrate the robustness of the proposed topology optimization method. It is straightforward to extend this proposed method into the other two-phase flows with two immiscible fluids.

With the use of temporal derivative of flux function, a two-stage temporal discretization has been recently proposed in the design of fourth-order schemes based on the generalized Riemann problem (GRP) [21] and gas-kinetic scheme (GKS) [28]. In this paper, the fourth-order gas-kinetic scheme will be extended to solve the compressible multicomponent flow equations, where the two-stage temporal discretization and fifth-order WENO reconstruction will be used in the construction of the scheme. Based on the simplified two-species BGK model [41], the coupled Euler equations for individual species will be solved. Each component has its individual gas distribution function and the equilibrium states for each component are coupled by the physical requirements of total momentum and energy conservation in particle collisions. The second-order flux function is used to achieve the fourth-order temporal accuracy, and the robustness is as good as the second-order schemes. At the same time, the source terms, such as the gravitational force and the chemical reaction, will be explicitly included in the two-stage temporal discretization through their temporal derivatives. Many numerical tests from the shock-bubble interaction to ZND detonative waves are presented to validate the current approach.

In this article a comparison study of the numerical methods for compressible two-phase flows is presented. Although many numerical methods have been developed in recent years to deal with the jump conditions at the fluid-fluid interfaces in compressible multiphase flows, there is a lack of a detailed comparison of these methods. With this regard, the transport five equation model, the modified ghost fluid method and the cut-cell method are investigated here as the typical methods in this field. A variety of numerical experiments are conducted to examine their performance in simulating inviscid compressible two-phase flows. Numerical experiments include Richtmyer-Meshkov instability, interaction between a shock and a rectangle SF6 bubble, Rayleigh collapse of a cylindrical gas bubble in water and shock-induced bubble collapse, involving fluids with small or large density difference. Based on the numerical results, the performance of the method is assessed by the convergence order of the method with respect to interface position, mass conservation, interface resolution and computational efficiency.

We propose a fully conservative and less oscillatory multi-moment scheme for the approximation of hyperbolic conservation laws. The proposed scheme (CIP-CSL3ENO) is based on two CIP-CSL3 schemes and the essentially non-oscillatory (ENO) scheme. In this paper, we also propose an ENO indicator for the multimoment framework, which intentionally selects non-smooth stencil but can efficiently minimize numerical oscillations. The proposed scheme is validated through various benchmark problems and a comparison with an experiment of two droplets collision/separation. The CIP-CSL3ENO scheme shows approximately fourth-order accuracy for smooth solution, and captures discontinuities and smooth solutions simultaneously without numerical oscillations for solutions which include discontinuities. The numerical results of two droplets collision/separation (3D free surface flow simulation) show that the CIP-CSL3ENO scheme can be applied to various types of fluid problems not only compressible flow problems but also incompressible and 3D free surface flow problems.

The Shan-Chen multiphase lattice Boltzmann model (LBM) coupled with Carnahan-Starling real-gas equation of state (C-S EOS)was proposed to simulate three-dimensional (3D) cavitation bubbles. Firstly, phase separation processes were predicted and the inter-phase large density ratio over 2×104 was captured successfully. The liquid-vapor density ratio at lower temperature is larger. Secondly, bubble surface tensions were computed and decreased with temperature increasing. Thirdly, the evolution of creation and condensation of cavitation bubbles were obtained. The effectiveness and reliability of present method were verified by energy barrier theory. The influences of temperature, pressure difference and critical bubble radius on cavitation bubbles were analyzed systematically. Only when the bubble radius is larger than the critical value will the cavitation occur, otherwise, cavitation bubbles will dissipate due to condensation. According to the analyses of radius change against time and the variation ratio of bubble radius, critical radius is larger under lower temperature and smaller pressure difference condition, thus bigger seed bubbles are needed to invoke cavitation. Under higher temperature and larger pressure difference, smaller seed bubbles can invoke cavitation and the expanding velocity of cavitation bubbles are faster. The cavitation bubble evolution including formation, developing and collapse was captured successfully under various pressure conditions.

The level set method is one of the most successful methods for the simulation of multi-phase flows. To keep the level set function close the signed distance function, the level set function is constantly reinitialized by solving a Hamilton-Jacobi type of equation during the simulation. When the fluid interface intersects with a solid wall, a moving contact line forms and the reinitialization of the level set function requires a boundary condition in certain regions on the wall. In this work, we propose to use the dynamic contact angle, which is extended from the contact line, as the boundary condition for the reinitialization of the level set function. The reinitialization equation and the equation for the normal extension of the dynamic contact angle form a coupled system and are solved simultaneously. The extension equation is solved on the wall and it provides the boundary condition for the reinitialization equation; the level set function provides the directions along which the contact angle is extended from the contact line. The coupled system is solved using the 3rd order TVD Runge-Kutta method and the Godunov scheme. The Godunov scheme automatically identifies the regions where the angle condition needs to be imposed. The numerical method is illustrated by examples in three dimensions.

Semi-Lagrangian (S-L) methods have no CFL stability constraint, and are more stable than the Eulerian methods. In the literature, the S-L method for the level-set re-initialization equation was complicated, which may be unnecessary. Since the re-initialization procedure is auxiliary, we propose to use the first-order S-L scheme coupled with a projection technique to improve the accuracy at the grid points just adjacent to the interface. Standard second-order S-L method is used for evolving the level-set convection equation. The implementation is simple, including on the block-structured adaptive mesh. The efficiency of the S-L method is demonstrated by extensive numerical examples including passive convection of interfaces with corners/kinks/large deformation under given velocity fields, a geometrical flow with topological changes, simulations of bubble/ droplet dynamics in incompressible two-phase flows. In terms of accuracy it is comparable to the other existing methods.

High fidelity modeling and simulation of moderately dense sprays at relatively low cost is still a major challenge for many applications. For that purpose, we introduce a new multi-fluid model based on a two-size moment formalism in sections, which are size intervals of discretization. It is derived from a Boltzmann type equation taking into account drag, evaporation and coalescence, which are representative of the complex terms that arise in multi-physics environments. The closure of the model comes from a reconstruction of the distribution. A piecewise affine reconstruction in size is thoroughly analyzed in terms of stability and accuracy, a key point for a high-fidelity and reliable description of the spray. Robust and accurate numerical methods are then developed, ensuring the realizability of the moments. The model and method are proven to describe the spray with a high accuracy in size and size-conditioned variables, resorting to a lower number of sections compared to one size moment methods. Moreover, robustness is ensured with efficient and tractable algorithms despite the numerous couplings and various algebra thanks to a tailored overall strategy. This strategy is successfully tested on a difficult 2D unsteady case, which proves the efficiency of the modeling and numerical choices.

The modified ghost fluid method (MGFM), due to its reasonable treatment for ghost fluid state, has been shown to be robust and efficient when applied to compressible multi-medium flows. Other feasible definitions of the ghost fluid state, however, have yet to be systematically presented. By analyzing all possible wave structures and relations for a multi-medium Riemann problem, we derive all the conditions to define the ghost fluid state. Under these conditions, the solution in the real fluid region can be obtained exactly, regardless of the wave pattern in the ghost fluid region. According to the analysis herein, a practical ghost fluid method (PGFM) is proposed to simulate compressible multi-medium flows. In contrast with the MGFM where three degrees of freedomat the interface are required to define the ghost fluid state, only one degree of freedomis required in this treatment. However, when these methods proved correct in theory are used in computations for the multi-medium Riemann problem, numerical errors at the material interface may be inevitable. We show that these errors are mainly induced by the single-medium numerical scheme in essence, rather than the ghost fluid method itself. Equipped with some density-correction techniques, the PGFM is found to be able to suppress these unphysical solutions dramatically.

A mathematical multi-zone ice accretion model used in the numerical simulation of icing on airfoil surface based on three water states, namely, continuous film, rivulets and beads is studied in this paper. An improved multi-zone roughness model is proposed. According to the flow state of liquid water and film flow, rivulets flow governing equations are established to calculate film mass distribution, film velocity, rivulet wetness factor and rivulet mass distribution. Force equilibrium equations of droplet are used to establish the critical conditions of water film broken into rivulets and rivulets broken into beads. The temperature conduction inside the water layer and ice layer is considered. Using the proposed model ice accretion on a NACA0012 airfoil profile with a 4° angle of attack under different icing conditions is simulated. Different ice shapes like glaze ice, mixed ice and rime ice are obtained, and the results agree well with icing wind tunnel experiment data. It can be seen that, water films are formed on the surface, and heights of the films vary with icing time and locations. This results in spatially-temporally varying surface roughness and heat transfer process, ultimately affects the ice prediction. Model simulations indicate that the process of water film formation and evolution cannot be ignored, especially under glaze ice condition.

In the first part, we study the convergence of discrete solutions to splitting schemes for two-phase flow with different mass densities suggested in [Guillen-Gonzalez, Tierra, J.Comput.Math. (6)2014]. They have been formulated for the diffuse interface model in [Abels, Garcke, Grün, M3AS, 2012, DOI:10.1142/S0218202511500138] which is consistent with thermodynamics. Our technique covers various discretization methods for phase-field energies, ranging from convex-concave splitting to difference quotient approaches for the double-well potential. In the second part of the paper, numerical experiments are presented in two space dimensions to identify discretizations of Cahn-Hilliard energies which are ϕ-stable and which do not reduce the acceleration of falling droplets. Finally, 3d simulations in axial symmetric geometries are shown to underline even more the full practicality of the approach.

A front tracking method combined with the real ghost fluid method (RGFM) is proposed for simulations of fluid interfaces in two-dimensional compressible flows. In this paper the Riemann problem is constructed along the normal direction of interface and the corresponding Riemann solutions are used to track fluid interfaces. The interface boundary conditions are defined by the RGFM, and the fluid interfaces are explicitly tracked by several connected marker points. The Riemann solutions are also used directly to update the flow states on both sides of the interface in the RGFM. In order to validate the accuracy and capacity of the new method, extensive numerical tests including the bubble advection, the Sod tube, the shock-bubble interaction, the Richtmyer-Meshkov instability and the gas-water interface, are simulated by using the Euler equations. The computational results are also compared with earlier computational studies and it shows good agreements including the compressible gas-water system with large density differences.

In this paper, we compute a phase field (diffuse interface) model of Cahn-Hilliard type for moving contact line problems governing the motion of isothermal multiphase incompressible fluids. The generalized Navier boundary condition proposed by Qian et al. [1] is adopted here. We discretize model equations using a continuous finite element method in space and a modified midpoint scheme in time. We apply a penalty formulation to the continuity equation which may increase the stability in the pressure variable. Two kinds of immiscible fluids in a pipe and droplet displacement with a moving contact line under the effect of pressure driven shear flow are studied using a relatively coarse grid. We also derive the discrete energy law for the droplet displacement case, which is slightly different due to the boundary conditions. The accuracy and stability of the scheme are validated by examples, results and estimate order.

Newtonian and non-Newtonian dispersed phase droplet formation in non-Newtonian continuous phase in T-junction and cross junction microchannels are investigated by the immiscible lattice BGK model. The effects of the non-Newtonian fluid power-law exponent, viscosity and interfacial tension on the generation of the droplet are studied. The final droplet size, droplet generation frequency, and detachment point of the droplet change with the power-law exponent. The results reveal that it is necessary to take into account the non-Newtonian rheology instead of simple Newtonian fluid assumption in numerical simulations. The present analysis also demonstrates that the lattice Boltzmann method is of potential to investigate the non-Newtonian droplet generation in multiphase flow.

Homogeneous equilibrium model (HEM) has been widely used in cavitating flow simulations. The major feature of this model is that a single equation of state (EOS) is proposed to describe the thermal behavior of bubbly liquid, where both kinematic and thermal equilibrium is assumed between two phases. In this paper, the HEM was coupled with multi-relaxation-time lattice Boltzmann model (MRT-LBM) and the acoustic behavior was simulated. Two approaches were applied alternatively: adjusting speed of sound (Buick, J. Phys. A, 2006, 39:13807-13815) and setting real gas EOS. Both approaches result in high accuracy in acoustic speed predictions for different void (gas) volume of fractions. It is demonstrated that LBM could be successfully applied as a Navier-Stokes equation solver for industrial applications. However, further dissipation and dispersion analysis shows that Shan-Chen type approaches of LBM are deficient, especially in large wave-number region.

The accurate simulation of disperse two-phase flows, where a discrete particulate condensed phase is transported by a carrier gas, is crucial for many applications; Eulerian approaches are well suited for high performance computations of such flows. However when the particles from the disperse phase have a significant inertia compared to the time scales of the flow, particle trajectory crossing (PTC) occurs i.e. the particle velocity distribution at a given location can become multi-valued. To properly account for such a phenomenon many Eulerian moment methods have been recently proposed in the literature. The resulting models hardly comply with a full set of desired criteria involving: 1- ability to reproduce the physics of PTC, at least for a given range of particle inertia, 2- well-posedness of the resulting set of PDEs on the chosen moments as well as guaranteed realizability, 3- capability of the model to be associated with a high order realizable numerical scheme for the accurate resolution of particle segregation in turbulent flows. The purpose of the present contribution is to introduce a multi-variate Anisotropic Gaussian closure for such particulate flows, in the spirit of the closure that has been suggested for out-of-equilibrium gas dynamics and which satisfies the three criteria. The novelty of the contribution is three-fold. First we derive the related moment system of conservation laws with source terms, and justify the use of such a model in the context of high Knudsen numbers, where collision operators play no role. We exhibit the main features and advantages in terms of mathematical structure and realizability. Then a second order accurate and realizable MUSCL/HLL scheme is proposed and validated. Finally the behavior of the method for the description of PTC is thoroughly investigated and its ability to account accurately for inertial particulate flow dynamics in typical configurations is assessed.

A computational fluid dynamics solver based on homogeneous cavitation model is employed to compute the two-phase cavitating flow. The model treats the two-phase regime as the homogeneous mixture of liquid and vapour which are locally assumed to be under both kinetic and thermodynamic equilibrium. As our focus is on pressure wave formation, propagation and its impact on cavitation bubble, the compressibility effects of liquid water have to be accounted for and hence the flow is considered to be compressible. The cavitating flow disturbed by the introduced pressure wave is simulated to investigate the unsteady features of cavitation due to the external perturbations. It is observed that the cavity becomes unstable, locally experiencing deformation or collapse, which depends on the shock wave intensity and freestream flow speed.