Simulation of Incompressible Free Surface Flow Using the Volume Preserving Level Set Method

This study aims to develop a numerical scheme in collocated Cartesian grids to solve the level set equation together with the incompressible two-phase flow equations. A seventh-order accurate upwinding combined compact difference (UCCD7) scheme has been developed for the approximation of the first-order spatial derivative terms shown in the level set equation. Developed scheme has a higher accuracy with a three-point grid stencil to minimize phase error. To preserve the mass of each phase all the time, the temporal derivative term in the level set equation is approximated by the sixth-order accurate symplectic Runge-Kutta (SRK6) scheme. All the simulated results for the dam-break, Rayleigh-Taylor instability, bubble rising, two-bubble merging, and milkcrown problems in two and three dimensions agree well with the available numerical or experimental results.     

Development of a High-Resolution Scheme for Solving the PNP-NS Equations in Curved Channels

A high-order finite difference scheme has been developed to approximate the spatial derivative terms present in the unsteady Poisson-Nernst-Planck (PNP) equations and incompressible Navier-Stokes (NS) equations. Near the wall the sharp solution profiles are resolved by using the combined compact difference (CCD) scheme developed in five-point stencil. This CCD scheme has a sixth-order accuracy for the second-order derivative terms while a seventh-order accuracy for the first-order derivative terms. PNP-NS equations have been also transformed to the curvilinear coordinate system to study the effects of channel shapes on the development of electroosmotic flow. In this study, the developed scheme has been analyzed rigorously through the modified equation analysis. In addition, the developed method has been computationally verified through four problems which are amenable to their own exact solutions. The electroosmotic flow details in planar and wavy channels have been explored with the emphasis on the formation of Coulomb force. Significance of different forces resulting from the pressure gradient, diffusion and Coulomb origins on the convective electroosmotic flow motion is also investigated in detail.     

An Interface-Capturing Method for Resolving Compressible Two-Fluid Flows with General Equation of State

In this study, a stable and robust interface-capturing method is developed to resolve inviscid, compressible two-fluid flows with general equation of state (EOS). The governing equations consist of mass conservation equation for each fluid, momentum and energy equations for mixture and an advection equation for volume fraction of one fluid component. Assumption of pressure equilibrium across an interface is used to close the model system. MUSCL-Hancock scheme is extended to construct input states for Riemann problems, whose solutions are calculated using generalized HLLC approximate Riemann solver. Adaptive mesh refinement (AMR) capability is built into hydrodynamic code. The resulting method has some advantages. First, it is very stable and robust, as the advection equation is handled properly. Second, general equation of state can model more materials than simple EOSs such as ideal and stiffened gas EOSs for example. In addition, AMR enables us to properly resolve flow features at disparate scales. Finally, this method is quite simple, time-efficient and easy to implement.     

Development of a Volume of Fluid Method for Computing Interfacial Incompressible Fluid Flows

This study is aimed to develop a volume of fluid (VOF) method to capture the free surface flow. The incompressible two-phase flow is computed by second-order Adams-Bashforth algorithm with a uniform staggered Cartesian grid arrangement. The tangent of hyperbola for interface capturing (THINC) scheme and weighted linear interface calculation (WLIC) based geometrical reconstruction procedure have been implemented in the operator-splitting method for the VOF method. The proposed VOF method preserves mass well, and the interface normal vector can be easily estimated from the level set (LS) function. The LS function, which is a continuous signed distance function around the interface, is represented by solving the re-initialization equation. Numerical results using the present scheme are compared with experimental data and other numerical results in the Rayleigh-Taylor instability, dam-break flow, travelling solitary wave, Kelvin-Helmholtz instability, rising bubble and merging bubble problems. We also present numerical results in detail between computations made with the proposed VOF method and computations made with the conventional LS method.     

An Implementation of MAC Grid-Based IIM-Stokes Solver for Incompressible Two-Phase Flows

In this paper, a novel implementation of immersed interface method combined with Stokes solver on a MAC staggered grid for solving the steady two-fluid Stokes equations with interfaces. The velocity components along the interface are introduced as two augmented variables and the resulting augmented equation is then solved by the GMRES method. The augmented variables and/or the forces are related to the jumps in pressure and the jumps in the derivatives of both pressure and velocity, and are interpolated using cubic splines and are then applied to the fluid through the jump conditions. The Stokes equations are discretized on a staggered Cartesian grid via a second order finite difference method and solved by the conjugate gradient Uzawa-type method. The numerical results show that the overall scheme is second order accuracy. The major advantages of the present IIM-Stokes solver are the efficiency and flexibility in terms of types of fluid flow and different boundary conditions. The proposed method avoids solution of the pressure Poisson equation, and comparisons are made to show the advantages of time savings by the present method. The generalized two-phase Stokes solver with correction terms has also been applied to incompressible two-phase Navier-Stokes flow.     

Conservative Semi-Lagrangian Finite Difference WENO Formulations with Applications to the Vlasov Equation

In this paper, we propose a new conservative semi-Lagrangian (SL) finite difference (FD) WENO scheme for linear advection equations, which can serve as a base scheme for the Vlasov equation by Strang splitting [4]. The reconstruction procedure in the proposed SL FD scheme is the same as the one used in the SL finite volume (FV) WENO scheme [3]. However, instead of inputting cell averages and approximate the integral form of the equation in a FV scheme, we input point values and approximate the differential form of equation in a FD spirit, yet retaining very high order (fifth order in our experiment) spatial accuracy. The advantage of using point values, rather than cell averages, is to avoid the second order spatial error, due to the shearing in velocity (v) and electrical field (E) over a cell when performing the Strang splitting to the Vlasov equation. As a result, the proposed scheme has very high spatial accuracy, compared with second order spatial accuracy for Strang split SL FV scheme for solving the Vlasov-Poisson (VP) system. We perform numerical experiments on linear advection, rigid body rotation problem; and on the Landau damping and two-stream instabilities by solving the VP system. For comparison, we also apply (1) the conservative SL FD WENO scheme, proposed in [22] for incompressible advection problem, (2) the conservative SL FD WENO scheme proposed in [21] and (3) the non-conservative version of the SL FD WENO scheme in [3] to the same test problems. The performances of different schemes are compared by the error table, solution resolution of sharp interface, and by tracking the conservation of physical norms, energies and entropies, which should be physically preserved.     

A Phase-Field Model Coupled with Lattice Kinetics Solver for Modeling Crystal Growth in Furnaces

In this study, we present a new numerical model for crystal growth in a vertical solidification system. This model takes into account the buoyancy induced convective flow and its effect on the crystal growth process. The evolution of the crystal growth interface is simulated using the phase-field method. A semi-implicit lattice kinetics solver based on the Boltzmann equation is employed to model the unsteady incompressible flow. This model is used to investigate the effect of furnace operational conditions on crystal growth interface profiles and growth velocities. For a simple case of macroscopic radial growth, the phase-field model is validated against an analytical solution. The numerical simulations reveal that for a certain set of temperature boundary conditions, the heat transport in the melt near the phase interface is diffusion dominant and advection is suppressed.     

Development of a Combined Compact Difference Scheme to Simulate Soliton Collision in a Shallow Water Equation

In this paper a three-step scheme is applied to solve the Camassa-Holm (CH) shallow water equation. The differential order of the CH equation has been reduced in order to facilitate development of numerical scheme in a comparatively smaller grid stencil. Here a three-point seventh-order spatially accurate upwinding combined compact difference (CCD) scheme is proposed to approximate the first-order derivative term. We conduct modified equation analysis on the CCD scheme and eliminate the leading discretization error terms for accurately predicting unidirectional wave propagation. The Fourier analysis is carried out as well on the proposed numerical scheme to minimize the dispersive error. For preserving Hamiltonians in Camassa-Holm equation, a symplecticity conserving time integrator has been employed. The other main emphasis of the present study is the use of u−P−α formulation to get nondissipative CH solution for peakon-antipeakon and soliton-anticuspon head-on wave collision problems.     

Construction,Analysis and Application of Coupled Compact Difference Scheme in Computational Acoustics and Fluid Flow Problems

In the present work, a new type of coupled compact difference scheme has been proposed for the solution of computational acoustics and flow problems. The proposed scheme evaluates the first, the second and the fourth derivative terms simultaneously. Derived compact difference scheme has a significant spectral resolution and a physical dispersion relation preserving (DRP) ability over a considerable wavenumber range when a fourth order four stage Runge-Kutta scheme is used for the time integration. Central stencil has been used for the present numerical scheme to evaluate spatial derivative terms. Derived scheme has the capability of adding numerical diffusion adaptively to attenuate spurious high wavenumber oscillations responsible for numerical instabilities. The DRP nature of the proposed scheme across a wider wavenumber range provides accurate results for the model wave equations as well as computational acoustic problems. In addition to the attractive feature of adaptive diffusion, present scheme also helps to control spurious reflections from the domain boundaries and is projected as an alternative to the perfectly matched layer (PML) technique.     

Consistent Forcing Scheme in the Simplified Lattice Boltzmann Method for Incompressible Flows

Considering the fact that the lattice discrete effects are neglected while introducing a body force into the simplified lattice Boltzmann method (SLBM), we propose a consistent forcing scheme in SLBM for incompressible flows with external forces. The lattice discrete effects are considered at the level of distribution functions in the present forcing scheme. Consequently, it is more accurate compared with the original forcing scheme used in SLBM. Through Taylor series expansion and Chapman-Enskog (CE) expansion analysis, the present forcing scheme can be proven to recover the macroscopic Navier-Stokes (N-S) equations. Then, the macroscopic equations are resolved through a fractional step technique. Furthermore, the material derivative term is discretized by the central difference method. To verify the results of the present scheme, we simulate with multiple forms of external force interactions including the space- and time-dependent body forces. Hence, the present forcing scheme overcomes the disadvantages of the original forcing scheme and the body force can be accurately imposed in the present scheme even when a coarse mesh is applied while the original scheme fails. Excellent agreements between the analytical solutions and our numerical results can be observed.     

A Compact Third-Order Gas-Kinetic Scheme for Compressible Euler and Navier-Stokes Equations

In this paper, a compact third-order gas-kinetic scheme is proposed for the compressible Euler and Navier-Stokes equations. The main reason for the feasibility to develop such a high-order scheme with compact stencil, which involves only neighboring cells, is due to the use of a high-order gas evolution model. Besides the evaluation of the time-dependent flux function across a cell interface, the high-order gas evolution model also provides an accurate time-dependent solution of the flow variables at a cell interface. Therefore, the current scheme not only updates the cell averaged conservative flow variables inside each control volume, but also tracks the flow variables at the cell interface at the next time level. As a result, with both cell averaged and cell interface values, the high-order reconstruction in the current scheme can be done compactly. Different from using a weak formulation for high-order accuracy in the Discontinuous Galerkin method, the current scheme is based on the strong solution, where the flow evolution starting from a piecewise discontinuous high-order initial data is precisely followed. The cell interface time-dependent flow variables can be used for the initial data reconstruction at the beginning of next time step. Even with compact stencil, the current scheme has third-order accuracy in the smooth flow regions, and has favorable shock capturing property in the discontinuous regions. It can be faithfully used from the incompressible limit to the hypersonic flow computations, and many test cases are used to validate the current scheme. In comparison with many other high-order schemes, the current method avoids the use of Gaussian points for the flux evaluation along the cell interface and the multi-stage Runge-Kutta time stepping technique. Due to its multidimensional property of including both derivatives of flow variables in the normal and tangential directions of a cell interface, the viscous flow solution, especially those with vortex structure, can be accurately captured. With the same stencil of a second order scheme, numerical tests demonstrate that the current scheme is as robust as well-developed second-order shock capturing schemes, but provides more accurate numerical solutions than the second order counterparts.     

An Indirect-Forcing Immersed Boundary Method for Incompressible Viscous Flows with Interfaces on Irregular Domains

An indirect-forcing immersed boundary method for solving the incompressible Navier-Stokes equations involving the interfaces and irregular domains is developed. The rigid boundaries and interfaces are represented by a number of Lagrangian control points. Stationary rigid boundaries are embedded in the Cartesian grid and singular forces at the rigid boundaries are applied to impose the prescribed velocity conditions. The singular forces at the interfaces and the rigid boundaries are then distributed to the nearby Cartesian grid points using the immersed boundary method. In the present work, the singular forces at the rigid boundaries are computed implicitly by solving a small system of equations at each time step to ensure that the prescribed velocity condition at the rigid boundary is satisfied exactly. For deformable interfaces, the forces that the interface exerts on the fluid are computed from the configuration of the elastic interface and are applied to the fluid. The Navier-Stokes equations are discretized using finite difference method on a staggered uniform Cartesian grid by a second order accurate projection method. The ability of the method to simulate viscous flows with interfaces on irregular domains is demonstrated by applying to the rotational flow problem, the relaxation of an elastic membrane and flow in a constriction with an immersed elastic membrane.     

A Discrete Flux Scheme for Aerodynamic and Hydrodynamic Flows

The objective of this paper is to seek an alternative to the numerical simulation of the Navier-Stokes equations by a method similar to solving the BGK-type modeled lattice Boltzmann equation. The proposed method is valid for both gas and liquid flows. A discrete flux scheme (DFS) is used to derive the governing equations for two distribution functions; one for mass and another for thermal energy. These equations are derived by considering an infinitesimally small control volume with a velocity lattice representation for the distribution functions. The zero-order moment equation of the mass distribution function is used to recover the continuity equation, while the first-order moment equation recovers the linear momentum equation. The recovered equations are correct to the first order of the Knudsen number (Kn); thus, satisfying the continuum assumption. Similarly, the zero-order moment equation of the thermal energy distribution function is used to recover the thermal energy equation. For aerodynamic flows, it is shown that the finite difference solution of the DFS is equivalent to solving the lattice Boltzmann equation (LBE) with a BGK-type model and a specified equation of state. Thus formulated, the DFS can be used to simulate a variety of aerodynamic and hydrodynamic flows. Examples of classical aeroacoustics, compressible flow with shocks, incompressible isothermal and non-isothermal Couette flows, stratified flow in a cavity, and double diffusive flow inside a rectangle are used to demonstrate the validity and extent of the DFS. Very good to excellent agreement with known analytical and/or numerical solutions is obtained; thus lending evidence to the DFS approach as an alternative to solving the Navier-Stokes equations for fluid flow simulations.     

An Immersed Interface Method for the Simulation of Inextensible Interfaces in Viscous Fluids

In this paper, an immersed interface method is presented to simulate the dynamics of inextensible interfaces in an incompressible flow. The tension is introduced as an augmented variable to satisfy the constraint of interface inextensibility, and the resulting augmented system is solved by the GMRES method. In this work, the arclength of the interface is locally and globally conserved as the enclosed region undergoes deformation. The forces at the interface are calculated from the configuration of the interface and the computed augmented variable, and then applied to the fluid through the related jump conditions. The governing equations are discretized on a MAC grid via a second-order finite difference scheme which incorporates jump contributions and solved by the conjugate gradient Uzawa-type method. The proposed method is applied to several examples including the deformation of a liquid capsule with inextensible interfaces in a shear flow. Numerical results reveal that both the area enclosed by interface and arclength of interface are conserved well simultaneously. These provide further evidence on the capability of the present method to simulate incompressible flows involving inextensible interfaces.     

Numerical Simulation for Moving Contact Line with Continuous Finite Element Schemes

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.     

Development of an Explicit Symplectic Scheme that Optimizes the Dispersion-Relation Equation of the Maxwell's Equations

In this paper an explicit finite-difference time-domain scheme for solving the Maxwell's equations in non-staggered grids is presented. The proposed scheme for solving the Faraday's and Ampère's equations in a theoretical manner is aimed to preserve discrete zero-divergence for the electric and magnetic fields. The inherent local conservation laws in Maxwell's equations are also preserved discretely all the time using the explicit second-order accurate symplectic partitioned Runge-Kutta scheme. The remaining spatial derivative terms in the semi-discretized Faraday's and Ampère's equations are then discretized to provide an accurate mathematical dispersion relation equation that governs the numerical angular frequency and the wavenumbers in two space dimensions. To achieve the goal of getting the best dispersive characteristics, we propose a fourth-order accurate space centered scheme which minimizes the difference between the exact and numerical dispersion relation equations. Through the computational exercises, the proposed dual-preserving solver is computationally demonstrated to be efficient for use to predict the long-term accurate Maxwell's solutions.     

An All Speed Second Order IMEX Relaxation Scheme for the Euler Equations

We present an implicit-explicit finite volume scheme for the Euler equations. We start from the non-dimensionalised Euler equations where we split the pressure in a slow and a fast acoustic part. We use a Suliciu type relaxation model which we split in an explicit part, solved using a Godunov-type scheme based on an approximate Riemann solver, and an implicit part where we solve an elliptic equation for the fast pressure. The relaxation source terms are treated projecting the solution on the equilibrium manifold. The proposed scheme is positivity preserving with respect to the density and internal energy and asymptotic preserving towards the incompressible Euler equations. For this first order scheme we give a second order extension which maintains the positivity property. We perform numerical experiments in 1D and 2D to show the applicability of the proposed splitting and give convergence results for the second order extension.     

Finite Volume Lattice Boltzmann Method for Nearly Incompressible Flows on Arbitrary Unstructured Meshes

A genuine finite volume method based on the lattice Boltzmann equation (LBE) for nearly incompressible flows is developed. The proposed finite volume lattice Boltzmann method (FV-LBM) is grid-transparent, i.e., it requires no knowledge of cell topology, thus it can be implemented on arbitrary unstructured meshes for effective and efficient treatment of complex geometries. Due to the linear advection term in the LBE, it is easy to construct multi-dimensional schemes. In addition, inviscid and viscous fluxes are computed in one step in the LBE, as opposed to in two separate steps for the traditional finite-volume discretization of the Navier-Stokes equations. Because of its conservation constraints, the collision term of the kinetic equation can be treated implicitly without linearization or any other approximation, thus the computational efficiency is enhanced. The collision with multiple-relaxation-time (MRT) model is used in the LBE. The developed FV-LBM is of second-order convergence. The proposed FV-LBM is validated with three test cases in two-dimensions: (a) the Poiseuille flow driven by a constant body force; (b) the Blasius boundary layer; and (c) the steady flow past a cylinder at the Reynolds numbers Re=10, 20, and 40. The results verify the designed accuracy and efficacy of the proposed FV-LBM.     

Uniformity of the fluid flow velocities within hollow fiber membranes of blood oxygenation devices

A finite volume-based computational model was developed to investigate the uniformity of the fluid flow across the hollow fiber membranes in blood oxygenation devices. A two-dimensional annular cross section of a blood oxygenation device including about 3,300 hollow fiber membranes was used in the computation model. The equations governing the steady incompressible laminar flow in the blood oxygenation device were solved numerically and the results were compared with those obtained from the equivalent porous medium approximation. For the porous medium approximation, the Ergun equation was used for evaluating the permeability. The simulation results showed that the fluid molecules spend about six times longer in the fiber bundle region than that in its equivalent porous medium approximation model. The computational model also provides a more detailed fluid flow pattern in the membrane compartment of the blood oxygenator.     

Two-Grid Method for Miscible Displacement Problem by Mixed Finite Element Methods and Mixed Finite Element Method of Characteristics

The miscible displacement of one incompressible fluid by another in a porous medium is governed by a system of two equations. One is elliptic form equation for the pressure and the other is parabolic form equation for the concentration of one of the fluids. Since only the velocity and not the pressure appears explicitly in the concentration equation, we use a mixed finite element method for the approximation of the pressure equation and mixed finite element method with characteristics for the concentration equation. To linearize the mixed-method equations, we use a two-grid algorithm based on the Newton iteration method for this full discrete scheme problems. First, we solve the original nonlinear equations on the coarse grid, then, we solve the linearized problem on the fine grid used Newton iteration once. It is shown that the coarse grid can be much coarser than the fine grid and achieve asymptotically optimal approximation as long as the mesh sizes satisfy $h=H^2$ in this paper. Finally, numerical experiment indicates that two-grid algorithm is very effective.