A new volume-of-fluid method with a constructed distance function on general structured grids

A second-order volume-of-fluid method (VOF) is presented for interface tracking and sharp interface treatment on general structured grids. Central to the new method is a second-order distance function construction scheme on a general structured grid based on the reconstructed interface. A novel technique is developed for evaluating the interface normal vector using the distance function. With the normal vector, the interface is reconstructed from the volume fraction function via a piecewise linear interface calculation (PLIC) scheme on the computational domain. Several numerical tests are conducted to demonstrate the accuracy and efficiency of the present method. In general, the new VOF method is more efficient than both the high-order level set and the coupled level set and volume-of-fluid (CLSVOF) methods. The results from the new method are better than those from the benchmark VOF method, particularly in the under-resolved regions, and are comparable to those from the CLSVOF method. Breaking waves over a submerged bump and around a wedge-shaped bow are simulated to demonstrate the application of the new method and sharp interface treatment in a two-phase flow solver on curvilinear grids. The computational results are in good agreement with the available experimental measurements.     

Spectral difference method for compressible flow on unstructured grids with mixed elements

This paper presents the development of a 2D solver for inviscid and viscous compressible flows using the spectral difference (SD) method for unstructured grids with mixed elements. A mixed quadrilateral and triangular grid is first refined using one-level h-refinement to generate a quadrilateral grid while keeping the curvature of boundary edges. The SD method designed for quadrilateral meshes can subsequently be applied for the refined unstructured grid. Results obtained with the SD method for both inviscid and viscous compressible flows compare well with analytical solutions and other published results.     

An all-speed relaxation scheme for interface flows with surface tension

We consider interface flows where compressibility and capillary forces (surface tension) are significant. These flows are described by a non-conservative, unconditionally hyperbolic multiphase model. The numerical approximation is based on finite-volume method for unstructured grids. At the discrete level, the surface tension is approximated by a volume force (CSF formulation). The interface physical properties are recovered by designing an appropriate linearized Riemann solver (Relaxation scheme) that prevents spurious oscillations near material interfaces. For low-speed flows, a preconditioning linearization is proposed and the low Mach asymptotic is formally recovered. Numerical computations, for a bubble equilibrium, converge to the required Laplace law and the dynamic of a drop, falling under gravity, is in agreement with experimental observations.     

A coupled finite volume solver for the solution of incompressible flows on unstructured grids

This paper reports on a newly developed fully coupled pressure-based algorithm for the solution of laminar incompressible flow problems on collocated unstructured grids. The implicit pressure-velocity coupling is accomplished by deriving a pressure equation in a procedure similar to a segregated SIMPLE algorithm using the Rhie–Chow interpolation technique and assembling the coefficients of the momentum and continuity equations into one diagonally dominant matrix. The extended systems of continuity and momentum equations are solved simultaneously and their convergence is accelerated by using an algebraic multigrid solver. The performance of the coupled approach as compared to the segregated approach, exemplified by SIMPLE, is tested by solving five laminar flow problems using both methodologies and comparing their computational costs. Results indicate that the number of iterations needed by the coupled solver for the solution to converge to a desired level on both structured and unstructured meshes is grid independent. For relatively coarse meshes, the CPU time required by the coupled solver on structured grid is lower than the CPU time required on unstructured grid. On dense meshes however, this is no longer true. For low and moderate values of the grid aspect ratio, the number of iterations required by the coupled solver remains unchanged, while the computational cost slightly increases. For structured and unstructured grid systems, the required number of iterations is almost independent of the grid size at any value of the grid expansion ratio. Recorded CPU time values show that the coupled approach substantially reduces the computational cost as compared to the segregated approach with the reduction rate increasing as the grid size increases.     

Hierarchical reconstruction for spectral volume method on unstructured grids

The hierarchical reconstruction (HR) [Y.-J. Liu, C.-W. Shu, E. Tadmor, M.-P. Zhang, Central discontinuous Galerkin methods on overlapping cells with a non-oscillatory hierarchical reconstruction, SIAM J. Numer. Anal. 45 (2007) 2442–2467; Z.-L. Xu, Y.-J. Liu, C.-W. Shu, Hierarchical reconstruction for discontinuous Galerkin methods on unstructured grids with a WENO type linear reconstruction and partial neighboring cells, J. Comput. Phys. 228 (2009) 2194–2212] is applied to a piecewise quadratic spectral volume method on two-dimensional unstructured grids as a limiting procedure to prevent spurious oscillations in numerical solutions. The key features of this HR are that the reconstruction on each control volume only uses adjacent control volumes, which forms a compact stencil set, and there is no truncation of higher degree terms of the polynomial. We explore a WENO-type linear reconstruction on each hierarchical level for the reconstruction of high degree polynomials. Numerical computations for scalar and system of nonlinear hyperbolic equations are performed. We demonstrate that the hierarchical reconstruction can generate essentially non-oscillatory solutions while keeping the resolution and desired order of accuracy for smooth solutions.     

液舱横荡非结构网格高阶格式数值模拟

提出一种液舱横荡数值模拟的方法,将气液两相交界面视为物理间断,通过高阶精度离散格式捕捉间断.根据NVD(Normalized Variable Diagram)实现非结构化网格上高精度离散格式,建立固定网格上自由表面运动模拟方法.在开发的非结构网格有限体积法求解器GTEA(General Transport Equation Analyzer)基础上,实现上述方法.首先对经典的溃坝过程进行模拟,并与文献结果对比验证方法和程序的可信度.对二维矩形液舱在不同激振频率时的横荡进行数值计算,并与实验以及商业软件CFX计算结果进行比较.结果表明方法和软件可以模拟自由面的翻卷、破碎运动现象,对距自由面较深点处流体载荷的计算结果与实验值符合较好,与商业软件CFX相比,在相同计算网格下,算法可以更好的计算次峰值,验证方法正确可行.     

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.     

A spectrally refined interface approach for simulating multiphase flows

This paper presents a novel approach to phase-interface transport based on pseudo-spectral sub-grid refinement of a level set function. In each flow solver grid cell, a set of quadrature points is introduced on which the value of the level set function is known. This methodology allows to define a polynomial reconstruction of the level set function in each cell. The transport is performed using a semi-Lagrangian technique, removing all constraints on the time step size. Such an approach provides sub-cell resolution of the phase-interface and leads to excellent accuracy in the transport, while a reasonable cost is obtained by pre-computing some of the metrics associated with the polynomials. To couple this approach with a flow solver, an converging curvature computation is introduced. First, a second order explicit distance to the sub-grid interface is reconstructed on the flow solver mesh. Then, a least squares approach is employed to extract the curvature from this distance function. This technique is found to combine the high accuracy and good conservation found in the particle level set method with the converging curvature usually obtained with classical high order PDE transport of the level set function. Tests are presented for both transport as well as two-phase flows, that suggest that this technique is capable of retaining the thin liquid structures that are expected in turbulent atomization of liquids.     

Three-dimensional elliptic solvers for interface problems and applications

Second-order accurate elliptic solvers using Cartesian grids are presented for three-dimensional interface problems in which the coefficients, the source term, the solution and its normal flux may be discontinuous across an interface. One of our methods is designed for general interface problems with variable but discontinuous coefficient. The scheme preserves the discrete maximum principle using constrained optimization techniques. An algebraic multigrid solver is applied to solve the discrete system. The second method is designed for interface problems with piecewise constant coefficient. The method is based on the fast immersed interface method and a fast 3D Poisson solver. The second method has been modified to solve Helmholtz/Poisson equations on irregular domains. An application of our method to an inverse interface problem of shape identification is also presented. In this application, the level set method is applied to find the unknown surface iteratively.     

非结构多边形网格滑移线开穴算法

基于非结构任意多边形网格管理体系及束缚滑移技术,针对多介质多次加卸载大变形问题的模拟,发展滑移线开穴算法,并利用金属碰撞分离模型进行验证.给出滑移线上的点和网格的开闭穴状态判别方法及网格开闭情况下滑移线上点的速度计算方法.该滑移计算技术将传统结构网格的开穴算法推广到拉氏非结构多边形网格中,保留拉氏滑移算法无缝连接和较好模拟物质界面的优点,实现了模拟实际问题中物质界面开闭功能.在金属碰撞分离模型中采用多种网格、不同界面形状和多介质等条件进行测试,证明了算法的正确性.     

Parallel re-initialization of level set functions on distributed unstructured tetrahedral grids

Level set functions are employed to track interfaces in various application areas including simulation of two-phase flows and image segmentation. Often, a re-initializing algorithm is incorporated to transform a numerically instable level set function to a signed distance function. In this note, we present a parallel algorithm for re-initializing level set functions on unstructured, three-dimensional tetrahedral grids. The main idea behind this new domain decomposition approach is to combine a parallel brute-force re-initializing algorithm with an efficient way to compute distances between the interface and grid points. Time complexity and error analysis of the algorithm are investigated. Detailed numerical experiments demonstrate the accuracy and scalability on up to 128 processes.     

Adaptive solution techniques for simulating underwater explosions and implosions

Adaptive solution techniques are presented for simulating underwater explosions and implosions. The liquid is assumed to be an adiabatic fluid and the solution in the gas is assumed to be uniform in space. The solution in water is integrated in time using a semi-implicit time discretization of the adiabatic Euler equations. Results are presented either using a non-conservative semi-implicit algorithm or a conservative semi-implicit algorithm. A semi-implicit algorithm allows one to compute with relatively large time steps compared to an explicit method. The interface solver is based on the coupled level set and volume-of-fluid method (CLSVOF) [M. Sussman, A second order coupled level set and volume-of-fluid method for computing growth and collapse of vapor bubbles, J. Comput. Phys. 187 (2003) 110–136; M. Sussman, E.G. Puckett, A coupled level set and volume-of-fluid method for computing 3D and axisymmetric incompressible two-phase flows, J. Comput. Phys. 162 (2000) 301–337]. Several underwater explosion and implosion test cases are presented to show the performances of our proposed techniques.     

A novel coupled level set and volume of fluid method for sharp interface capturing on 3D tetrahedral grids

We present a new three-dimensional hybrid level set (LS) and volume of fluid (VOF) method for free surface flow simulations on tetrahedral grids. At each time step, we evolve both the level set function and the volume fraction. The level set function is evolved by solving the level set advection equation using a second-order characteristic based finite volume method. The volume fraction advection is performed using a bounded compressive normalized variable diagram (NVD) based scheme. The interface is reconstructed based on both the level set and the volume fraction information. The novelty of the method lies in that we use an analytic method for finding the intercepts on tetrahedral grids, which makes interface reconstruction efficient and conserves volume of fluid exactly. Furthermore, the advection of volume fraction makes use of the NVD concept and switches between different high resolution differencing schemes to yield a bounded scalar field, and to preserve both smoothness and sharp definition of the interface. The method is coupled to a well validated finite volume based Navier–Stokes incompressible flow solver. The code validation shows that our method can be employed to resolve complex interface changes efficiently and accurately. In addition, the centroid and intercept data available as a by-product of the proposed interface reconstruction scheme can be used directly in near-interface sub-grid models in large eddy simulation.     

八叉树网格到非结构混合网格的转换

在数值模拟中, 非结构网格的优势是可以采用相同的数值格式统一处理任意复杂的计算区域, 但在网格生成过程中难度大, 并且不容易控制网格质量。树结构网格可以认为是介于结构网格和非结构网格之间的一种网格, 目前已经有相对成熟的方法快速在复杂区域内生成二维四叉树网格和三维八叉树网格。在实际应用中, 数值方法往往需要在连接协调的非结构网格上做离散, 树结构网格中不同尺寸的网格之间连接不是协调的, 在应用上会受到很多限制。文章实现了树结构网格到非结构混合网格的转换, 这种转换在二维情况下就是将四叉树网格转换为非结构三角形和四边形的混合网格, 三维情况下则将八叉树网格转换为非结构混合网格。这一转换过程的难点在于需要考虑数千种不同的八叉树单元, 并给出能实现连接协调的非结构混合网格划分。可以出现的网格单元包括六面体、三棱柱、金字塔和四面体这4种不同情况。通过特别的分类, 实现了程序的自动生成, 这种程序自动生成技术一方面可以避免人工编写大量程序时的失误, 另一方面也使得对数以千计的不同情况的处理成为可能。通过对几个简单网格的测试, 对网格数据转换方法做了初步的验证。      

基于特征理论的二维可压缩流动的二阶拉氏算法

提出一种求解二维拉氏可压缩流体力学方程的中心型二阶精度有限体积方法.利用特征理论构造网格节点处的局部近似演化算子,算子用来求解网格节点处的速度及压力,利用这些物理量更新节点位置及计算网格界面通量.通过结合一定的重构方案,该方法达到时、空二阶精度,并且形式简单、计算量小,适用于结构网格与非结构网格.典型数值实验表明,本文格式具有良好的收敛性、对称性及鲁棒性,且能自然地求解多物质流动问题.     

基于非结构网格的不可压流动耦合求解器

采用非结构化网格有限容积法求解了不可压N-S方程组,对流项采用GAMMA格式,扩散项采用二阶中心差分格式建立离散方程,用SOAR算法处理压力与速度的耦合关系,得到了一种求解不可压N-S方程的非结构网格耦合求解器。通过方腔顶盖驱动流、后台阶绕流以及方腔自然对流等几个典型的算例,考察了求解器的计算精度及收敛特性,并与SIMPLE算法进行了比较,结果表明该求解器是有效可行的。     

On the characteristics-based ACM for incompressible flows

In this paper, the revised characteristics-based (CB) method for incompressible flows recently derived by Neofytou [P. Neofytou, Revision of the characteristic-based scheme for incompressible flows, J. Comput. Phys. 222 (2007) 475–484] has been further investigated. We have derived all the formulas for pressure and velocities from this revised CB method, which is based on the artificial compressibility method (ACM) [A.J. Chorin, A numerical solution for solving incompressible viscous flow problems, J. Comput. Phys. 2 (1967) 12]. Then we analyze the formulations of the original CB method [D. Drikakis, P.A. Govatsos, D.E. Papatonis, A characteristic based method for incompressible flows, Int. J. Numer. Meth. Fluids 19 (1994) 667–685; E. Shapiro, D. Drikakis, Non-conservative and conservative formulations of characteristics numerical reconstructions for incompressible flows, Int. J. Numer. Meth. Eng. 66 (2006) 1466–1482; D. Drikakis, P.K. Smolarkiewicz, On spurious vortical structures, J. Comput. Phys. 172 (2001) 309–325; F. Mallinger, D. Drikakis, Instability in three-dimensional, unsteady stenotic flows, Int. J. Heat Fluid Flow 23 (2002) 657–663; E. Shapiro, D. Drikakis, Artificial compressibility, characteristics-based schemes for variable density, incompressible, multi-species flows. Parts I. Derivation of different formulations and constant density limit, J. Comput. Phys. 210 (2005) 584–607; Y. Zhao, B. Zhang, A high-order characteristics upwind FV method for incompressible flow and heat transfer simulation on unstructured grids, Comput. Meth. Appl. Mech. Eng. 190 (5–7) (2000) 733–756] to investigate their consistency with the governing flow equations after convergence has been achieved. Furthermore we have implemented both formulations in an unstructured-grid finite volume solver [Y. Zhao, B. Zhang, A high-order characteristics upwind FV method for incompressible flow and heat transfer simulation on unstructured grids, Comput. Meth. Appl. Mech. Eng. 190 (5–7) (2000) 733–756]. Detailed numerical experiments show that both methods give almost identical solutions and convergence rates. Both can generate solutions which agree well with published results and experimental measurements. We thus conclude that both methods, being upwind schemes designed for the ACM, have the same performances in terms of accuracy and convergence speed, even though the revised method is more complex with less stringent assumptions made, while the original CB method is simpler due to the use of extra simplifying assumptions.     

Numerical prediction of interfacial instabilities: Sharp interface method (SIM)

We introduce a sharp interface method (SIM) for the direct numerical simulation of unstable fluid–fluid interfaces. The method is based on the level set approach and the structured adaptive mesh refinement technology, endowed with a corridor of irregular, cut-cell grids that resolve the interfacial region to third-order spatial accuracy. Key in that regard are avoidance of numerical mixing, and a least-squares interpolation method that is supported by irregular datasets distinctly on each side of the interface. Results on test problems show our method to be free of the spurious current problem of the continuous surface force method and to converge, on grid refinement, at near-theoretical rates. Simulations of unstable Rayleigh–Taylor and viscous Kelvin–Helmholtz flows are found to converge at near-theoretical rates to the exact results over a wide range of conditions. Further, we show predictions of neutral-stability maps of the viscous Kelvin–Helmholtz flows (Yih instability), as well as self-selection of the most unstable wave-number in multimode simulations of Rayleigh–Taylor instability. All these results were obtained with a simple seeding of random infinitesimal disturbances of interface-shape, as opposed to seeding by a complete eigenmode. For other than elementary flows the latter would normally not be available, and extremely difficult to obtain if at all. Sample comparisons with our code adapted to mimic typical diffuse interface treatments were not satisfactory for shear-dominated flows. On the other hand the sharp dynamics of our method would appear to be compatible and possibly advantageous to any interfacial flow algorithm in which the interface is represented as a discrete Heaviside function.     

Implementation of high order spline interpolations for tracking test particles in discretized fields

A systematic approach for constructing high order spline interpolation methods is proposed for fields known on regular, rectangular grids. These interpolation methods are tested in tracking trajectories of particles submitted to a force that derives from a potential known on a grid. The interplay between the time advancement scheme and the spatial interpolation is studied in detail and it is shown how the order of the trajectory solver is directly affected by the order of the spline interpolation. It is also shown how an interpolation method that preserves topological properties of physical fields can be better exploited with these higher order spline approximations.