A direct comparison between volume and surface tracking methods with a boundary-fitted coordinate transformation and third-order upwinding

A Volume Tracking (VT) and a Front Tracking (FT) algorithm are implemented and compared for locating the interface between two immiscible, incompressible, Newtonian fluids in a tube with a periodically varying, circular cross-section. Initially, the fluids are stationary and stratified in an axisymmetric arrangement so that one is around the axis of the tube (core fluid) and the other one surrounds it (annular fluid). A constant pressure gradient sets them in motion. With both VT and FT, a boundary-fitted coordinate transformation is applied and appropriate modifications are made to adopt either method in this geometry. The surface tension force is approximated using the continuous surface force method. All terms appearing in the continuity and momentum equations are approximated using centered finite differences in space and one-sided forward finite differences in time. In each time step, the incompressibility condition is enforced by a transformed Poisson equation, which is linear in pressure. This equation is solved by either direct LU decomposition or a Multigrid iterative solver. When the two fluids have the same density, the former method is about 3.5 times faster, but when they do not, the Multigrid solver is as much as 10 times faster than the LU decomposition. When the interface does not break and the Reynolds number remains small, the accuracy and rates of convergence of VT and FT are comparable. The well-known failure of centered finite differences arises as the Reynolds number increases and leads to non-physical oscillations in the interface and failure of both methods to converge with mesh refinement. These problems are resolved and computations with Reynolds as large as 500 converged by approximating the convective terms in the momentum equations by third-order upwind differences using Lagrangian Polynomials. When the volume of the core fluid or the Weber number decrease, increasing the importance of interfacial tension and leading to breakup of the interface forming a drop of core fluid, the FT method converges faster with mesh refinement than the VT method and upwinding may be required. Finally, examining the generation of spurious currents around a stationary “bubble” in the tube for Ohnesorge numbers between 0.1 and 10 it is found that the maximum velocity remains approximately the same in spite mesh refinements when VT is applied, whereas it is of the same order of magnitude for the coarsest mesh and monotonically decreases with mesh refinement when FT is applied.     

Coupled fluid–structure solver: The case of shock wave impact on monolithic and composite material plates

An unstructured adaptive mesh flow solver, a finite element structure solver and a moving mesh algorithm were implemented in the numerical simulation of the interaction between a shock wave and a structure. In the past, this interaction is mostly considered as one-way in the sense that the shock causes a transient load on the structure while it is reflected uneffected by the impact. A fully coupled approach was implemented in the present work which can account for the effects associated with a mutual interaction. This approach included a compressible flow Eulerian solver of second order accuracy in finite volume formulation for the fluid and a Langargian solver in finite element formulation for the solid structure. A novel implementation of advancing front moving mesh algorithm was made possible with the introduction of a flexible and efficient quad-edge data structure. Adaptive mesh refinement was introduced into the flow solver for improved accuracy as well. Numerical results are further validated by theoretical analysis, experimental data and results from other numerical simulations. Grid dependency study was performed and results showed that the physical phenomena and quantities were independent of the numerical grid chosen in the simulations. The results illuminated complicated flow phenomena and structure vibration patterns, which in order to be detected experimentally require capabilities beyond those of the current experimental techniques. The numerical simulations also successfully modelled the aero-acoustic damping effects on the structure, which do not exist in previous numerical models. Further analysis of the results showed that the mutual interaction is not linear and that the non-linearity arises because the wave propagation in the fluid is not linear and it cascades a non-linear and non-uniform loading on the plate. Non-linearity intensifies when the plate is vibrating at high frequency while the wave propagation speed is low.     

Semi-automatic sparse preconditioners for high-order finite element methods on non-uniform meshes

High-order finite elements often have a higher accuracy per degree of freedom than the classical low-order finite elements. However, in the context of implicit time-stepping methods, high-order finite elements present challenges to the construction of efficient simulations due to the high cost of inverting the denser finite element matrix. There are many cases where simulations are limited by the memory required to store the matrix and/or the algorithmic components of the linear solver. We are particularly interested in preconditioned Krylov methods for linear systems generated by discretization of elliptic partial differential equations with high-order finite elements. Using a preconditioner like Algebraic Multigrid can be costly in terms of memory due to the need to store matrix information at the various levels. We present a novel method for defining a preconditioner for systems generated by high-order finite elements that is based on a much sparser system than the original high-order finite element system. We investigate the performance for non-uniform meshes on a cube and a cubed sphere mesh, showing that the sparser preconditioner is more efficient and uses significantly less memory. Finally, we explore new methods to construct the sparse preconditioner and examine their effectiveness for non-uniform meshes. We compare results to a direct use of Algebraic Multigrid as a preconditioner and to a two-level additive Schwarz method.     

Mesh refinement algorithms in an unstructured solver for multiphase flow simulation using discrete particles

This study developed spray-adaptive mesh refinement algorithms with directional sensitivity in an unstructured solver to improve spray simulation for internal combustion engine application. Inadequate spatial resolution is often found to cause inaccuracies in spray simulation using the Lagrangian–Eulerian approach due to the over-estimated diffusion and inappropriate liquid–gas phase coupling. Dynamic mesh refinement algorithms adaptive to fuel sprays and vapor gradients were developed in order to increase the grid resolution in the spray region to improve simulation accuracy. The local refinement introduced the coarse-fine face interface that requires advanced numerical schemes for flux calculation and grid rezoning with moving boundaries. To resolve the issue in flux calculation, this work implemented the refinement/coarsening algorithms into a collocated solver to avoid tedious interpolations in solving the momentum equations. A pressure correction method was applied to address unphysical pressure oscillations due to the collocation of pressure and velocity. An edge-based algorithm was used to evaluate the edge-centered quantities in order to account for the contributions from all the cells around an edge at the coarse-fine interface. A quasi-second-order upwind scheme with strong monotonicity was also modified to accommodate the coarse-fine interface for convective fluxes. To resolve the issue related to grid rezoning, rezoning was applied to the initial baseline mesh only and the new locations of the refined grids were obtained by interpolating the updated baseline mesh. The time step constraints were also re-evaluated to account for the change resulting from mesh refinement. The present refinement algorithm was used in simulating fuel sprays in an engine combustion chamber. It was found that the present approach could produce the same level of results as those using the uniformly fine mesh with substantially reduced computer time. Results also showed that this approach could alleviate the artifacts related to the Lagrangian discrete modeling of spray drops due to insufficient spatial resolution.     

三维对流扩散方程非等距网格上的四阶紧致格式及其多重网格方法

提出了数值求解三维变系数对流扩散方程非等距网格上的四阶精度19点紧致差分格式,为了提高求解效率,采用多重网格方法求解高精度格式所形成的大型代数方程组。数值实验结果表明本文方法对于不同的网格雷诺数问题,在精确性、稳定性和减少计算工作量方面均明显优于7点中心差分格式。     

Space–time discontinuous Galerkin finite element method for two-fluid flows

A novel numerical method for two-fluid flow computations is presented, which combines the space–time discontinuous Galerkin finite element discretization with the level set method and cut-cell based interface tracking. The space–time discontinuous Galerkin (STDG) finite element method offers high accuracy, an inherent ability to handle discontinuities and a very local stencil, making it relatively easy to combine with local hp-refinement. The front tracking is incorporated via cut-cell mesh refinement to ensure a sharp interface between the fluids. To compute the interface dynamics the level set method (LSM) is used because of its ability to deal with merging and breakup. Also, the LSM is easy to extend to higher dimensions. Small cells arising from the cut-cell refinement are merged to improve the stability and performance. The interface conditions are incorporated in the numerical flux at the interface and the STDG discretization ensures that the scheme is conservative as long as the numerical fluxes are conservative. The numerical method is applied to one and two dimensional two-fluid test problems using the Euler equations.     

A ghost fluid,level set methodology for simulating multiphase electrohydrodynamic flows with application to liquid fuel injection

In this paper, we present the development of a sharp numerical scheme for multiphase electrohydrodynamic (EHD) flows for a high electric Reynolds number regime. The electric potential Poisson equation contains EHD interface boundary conditions, which are implemented using the ghost fluid method (GFM). The GFM is also used to solve the pressure Poisson equation. The methods detailed here are integrated with state-of-the-art interface transport techniques and coupled to a robust, high order fully conservative finite difference Navier–Stokes solver. Test cases with exact or approximate analytic solutions are used to assess the robustness and accuracy of the EHD numerical scheme. The method is then applied to simulate a charged liquid kerosene jet.     

The immersed boundary method for advection–electrodiffusion with implicit timestepping and local mesh refinement

We describe an immersed boundary method for problems of fluid–solute-structure interaction. The numerical scheme employs linearly implicit timestepping, allowing for the stable use of timesteps that are substantially larger than those permitted by an explicit method, and local mesh refinement, making it feasible to resolve the steep gradients associated with the space charge layers as well as the chemical potential, which is used in our formulation to control the permeability of the membrane to the (possibly charged) solute. Low Reynolds number fluid dynamics are described by the time-dependent incompressible Stokes equations, which are solved by a cell-centered approximate projection method. The dynamics of the chemical species are governed by the advection–electrodiffusion equations, and our semi-implicit treatment of these equations results in a linear system which we solve by GMRES preconditioned via a fast adaptive composite-grid (FAC) solver. Numerical examples demonstrate the capabilities of this methodology, as well as its convergence properties.     

A mesh-dependent model for applying dynamic contact angles to VOF simulations

Typical VOF algorithms rely on an implicit slip that scales with mesh refinement, to allow contact lines to move along no-slip boundaries. As a result, solutions of contact line phenomena vary continuously with mesh spacing; this paper presents examples of that variation. A mesh-dependent dynamic contact angle model is then presented, that is based on fundamental hydrodynamics and serves as a more appropriate boundary condition at a moving contact line. This new boundary condition eliminates the stress singularity at the contact line; the resulting problem is thus well-posed and yields solutions that converge with mesh refinement. Numerical results are presented of a solid plate withdrawing from a fluid pool, and of spontaneous droplet spread at small capillary and Reynolds numbers.     

A Method of Lines Based on Immersed Finite Elements for Parabolic Moving Interface Problems

This article extends the finite element method of lines to a parabolic initial boundary value problem whose diffusion coefficient is discontinuous across an interface that changes with respect to time. The method presented here uses immersed finite element (IFE) functions for the discretization in spatial variables that can be carried out over a fixed mesh (such as a Cartesian mesh if desired), and this feature makes it possible to reduce the parabolic equation to a system of ordinary differential equations (ODE) through the usual semi-discretization procedure. Therefore, with a suitable choice of the ODE solver, this method can reliably and efficiently solve a parabolic moving interface problem over a fixed structured (Cartesian) mesh. Numerical examples are presented to demonstrate features of this new method.     

Transient acoustic radiation from impulsively accelerated bodies by the finite element method

Bodies under impulsive motion, immersed in an infinite acoustic fluid, severely put to test any numerical method for the transient exterior acoustic problem. Such problems, in the context of the finite element method (FEM), are not well studied. FE modeling of such problems requires truncation of the infinite fluid domain at a certain distance from the structure. The volume of computation depends upon the extent of this domain as well as the mesh density. The modeling of the fluid truncation boundary is crucial to the economy and accuracy of solution and various boundary dampers have been proposed in the literature for this purpose. The second order damper leads to unsymmetric boundary matrices and this necessitates the use of an unsymmetric equation solver for large problems. The present paper demonstrates the use of FEM with zeroth, first and second order boundary dampers in conjunction with an unsymmetric, out of core, banded equation solver for impulsive motion problems of rigid bodies in an acoustic fluid. The results compare well with those obtained from analytical methods.     

An improved Quiet Direct Simulation method for Eulerian fluids using a second-order scheme

In this paper, a second-order scheme for the Quiet Direct Simulation (QDS) of Eulerian fluids is proposed. The QDS method replaces the random sampling method used in Direct Simulation Monte Carlo (DSMC) methods with a technique whereby particles are moved, have their properties distributed onto a mesh, are destroyed and then are recreated deterministically from the properties stored on the mesh using Gauss–Hermite quadrature weights and abscissas. Particles are permitted to move in physically realistic directions so flux exchange is not limited to cells sharing an adjacent interface as in conventional, direction decoupled finite volume solvers. In this paper the method is extended by calculating the fluxes of mass, momentum and energy between cells assuming a linear variation of density, temperature and velocity in each cell and using these fluxes to update the mass, velocity and internal energy carried by each particle. This Euler solver has several advantages including large dynamic range, no statistical scatter in the results, true direction fluxes to all nearby neighbors and is computationally inexpensive. The second-order method is found to reduce the numerical diffusion of QDS as demonstrated in several verification studies. These include unsteady shock tube flow, a two-dimensional blast wave and of the development of Mach 3 flow over a forward facing step in a wind tunnel, which are compared with previous results from the literature wherever is possible. Finally the implementation of QUIETWAVE, a rapid method of simulating blast events in urban environments, is introduced and the results of a test case are presented.     

Development of Lattice Boltzmann Flux Solver for Simulation of Incompressible Flows

A lattice Boltzmann flux solver (LBFS) is presented in this work for simulation of incompressible viscous and inviscid flows. The new solver is based on Chapman-Enskog expansion analysis, which is the bridge to link Navier-Stokes (N-S) equations and lattice Boltzmann equation (LBE). The macroscopic differential equations are discretized by the finite volume method, where the flux at the cell interface is evaluated by local reconstruction of lattice Boltzmann solution from macroscopic flow variables at cell centers. The new solver removes the drawbacks of conventional lattice Boltzmann method such as limitation to uniform mesh, tie-up of mesh spacing and time interval, limitation to viscous flows. LBFS is validated by its application to simulate the viscous decaying vortex flow, the driven cavity flow, the viscous flow past a circular cylinder, and the inviscid flow past a circular cylinder. The obtained numerical results compare very well with available data in the literature, which show that LBFS has the second order of accuracy in space, and can be well applied to viscous and inviscid flow problems with non-uniform mesh and curved boundary.     

Solid–fluid diffuse interface model in cases of extreme deformations

Diffuse interface methods have been recently proposed and successfully used for accurate compressible multi-fluid computations Abgrall [1]; Kapila et al. [20]; Saurel et al. [30]. These methods deal with extended systems of hyperbolic equations involving a non-conservative volume fraction equation and relaxation terms. Following the same theoretical frame, we derive here an Eulerian diffuse interface model for elastic solid-compressible fluid interactions in situations involving extreme deformations. Elastic effects are included following the Eulerian conservative formulation proposed in Godunov [16], Miller and Colella [23], Godunov and Romenskii [17], Plohr and Plohr [27] and Gavrilyuk et al. [14]. We apply first the Hamilton principle of stationary action to derive the conservative part of the model. The relaxation terms are then added which are compatible with the entropy inequality. In the limit of vanishing volume fractions the Euler equations of compressible fluids and a conservative hyperelastic model are recovered. It is solved by a unique hyperbolic solver valid at each mesh point (pure fluid, pure solid and mixture cell). Capabilities of the model and methods are illustrated on various tests of impacts of solids moving in an ambient compressible fluid.     

Multi-scale Godunov-type method for cell-centered discrete Lagrangian hydrodynamics

This work presents a multi-dimensional cell-centered unstructured finite volume scheme for the solution of multimaterial compressible fluid flows written in the Lagrangian formalism. This formulation is considered in the Arbitrary-Lagrangian–Eulerian (ALE) framework with the constraint that the mesh velocity and the fluid velocity coincide. The link between the vertex velocity and the fluid motion is obtained by a formulation of the momentum conservation on a class of multi-scale encased volumes around mesh vertices. The vertex velocity is derived with a nodal Riemann solver constructed in such a way that the mesh motion and the face fluxes are compatible. Finally, the resulting scheme conserves both momentum and total energy and, it satisfies a semi-discrete entropy inequality. The numerical results obtained for some classical 2D and 3D hydrodynamic test cases show the robustness and the accuracy of the proposed algorithm.     

一种分布式内存环境下的并行二分网格加密算法

给出了一个基于ALBERT(Adaptive multi-Level finite element toolbox using Bisection refinement and Error control by Residual Techniques)的并行二分网格加密算法.其目的是开发基于ALBERT的、适合于分布式内存计算机的并行自适应有限元软件包.首先给出了针对ALBERT的并行化策略,然后重点介绍并行网格加密算法,并证明了并行算法和原有串行算法在加密结果上完全等效.最后,数值实验证明该并行算法能有效地在分布式内存的计算机上执行.     

剪切变稀液滴撞击不同浸润性壁面的数值模拟研究

基于有限元法,采用水平集方法捕捉相界面的移动,构建了液滴撞击固体壁面的数值模型.通过修正的幂律模型描述流体的非牛顿剪切变稀特性,探讨了剪切变稀特性对液滴撞击固体壁面后铺展行为的影响,分析了撞击不同浸润性壁面时剪切变稀特性对液滴撞击壁面行为的影响差异.研究结果表明:随着幂律指数m的减小,液滴撞击过程中的黏性耗散减小,液滴的形貌变化及无量纲参数变化更为显著.接触角为55°的情况下:当m降低至0.85时,液滴铺展过程中开始出现显著区别于牛顿流体液滴的振荡现象;当m降低至0.80时,液滴在回缩过程中会出现中心液膜断裂的情况.接触角为100°时,剪切变稀液滴均会出现振荡行为,振荡幅度随着m的减小而增大.接触角为160°时,牛顿流体液滴与剪切变稀液滴均会在回缩过程中弹起,但剪切变稀液滴的弹起速度更快.此外,基于数值计算结果,本文提出了接触角为55°情况下剪切变稀液滴撞击壁面后的最大无量纲铺展直径预测模型.     

Numerical comparison of Riemann solvers for astrophysical hydrodynamics

The idea of this work is to compare a new positive and entropy stable approximate Riemann solver by Francois Bouchut with a state-of the-art algorithm for astrophysical fluid dynamics. We implemented the new Riemann solver into an astrophysical PPM-code, the Prometheus code, and also made a version with a different, more theoretically grounded higher order algorithm than PPM. We present shock tube tests, two-dimensional instability tests and forced turbulence simulations in three dimensions. We find subtle differences between the codes in the shock tube tests, and in the statistics of the turbulence simulations. The new Riemann solver increases the computational speed without significant loss of accuracy.     

An RKDG finite element method for the one-dimensional inviscid compressible gas dynamics equations in a Lagrangian coordinate

In this paper,Runge-Kutta Discontinuous Galerkin(RKDG) finite element method is presented to solve the onedimensional inviscid compressible gas dynamic equations in a Lagrangian coordinate.The equations are discretized by the DG method in space and the temporal discretization is accomplished by the total variation diminishing Runge-Kutta method.A limiter based on the characteristic field decomposition is applied to maintain stability and non-oscillatory property of the RKDG method.For multi-medium fluid simulation,the two cells adjacent to the interface are treated differently from other cells.At first,a linear Riemann solver is applied to calculate the numerical ?ux at the interface.Numerical examples show that there is some oscillation in the vicinity of the interface.Then a nonlinear Riemann solver based on the characteristic formulation of the equation and the discontinuity relations is adopted to calculate the numerical ?ux at the interface,which suppresses the oscillation successfully.Several single-medium and multi-medium fluid examples are given to demonstrate the reliability and efficiency of the algorithm.     

A RKDG finite element method for the one-dimensional inviscid compressible gas dynamics equations in Lagrangian coordinate

In this paper, Runge-Kutta Discontinuous Galerkin (RKDG) finite element method is presented to solve the one-dimensional inviscid compressible gas dynamic equations in Lagrangian coordinate. The equations are discretized by the DG method in space and the temporal discretization is accomplished by the total variation diminishing Runge-Kutta method. A limiter based on the characteristic field decomposition is applied to maintain stability and non-oscillatory property of the RKDG method. For multi-medium fluid simulation, the two cells adjacent to the interface are treated differently from other cells. At first, a linear Riemann solver is applied to calculate the numerical flux at the interface. Numerical examples show that there is some oscillation in the vicinity of the interface. Then a nonlinear Riemann solver based on the characteristic formulation of the equation and the discontinuity relations is adopted to calculate the numerical flux at the interface, which suppress the oscillation successfully. Several single-medium and multi-medium fluid examples are given to demonstrate the reliability and efficiency of the algorithm.