Simulation Of Compressible Flows Research Articles

A parallel modular computing environment for three-dimensional multiresolution simulations of compressible flows

Numerical investigation of compressible flows faces two main challenges. In order to accurately describe the flow characteristics, high-resolution nonlinear numerical schemes are needed to capture discontinuities and resolve wide convective, acoustic and interfacial scale ranges. The simulation of realistic three-dimensional (3D) problems with state-of-the-art finite-volume method (FVM) based on approximate Riemann solvers with weighted nonlinear reconstruction schemes requires the usage of high-performance computing (HPC) architectures. Efficient compression algorithms reduce computational and memory load. Fully adaptive multiresolution (MR) algorithms with LTS have proven their potential for such applications. While modern central processing units (CPUs) requires multiple levels of parallelism to achieve peak performance, the fine-grained MR mesh adaptivity results in challenging compute/communication patterns. Moreover, local time stepping (LTS) incurs strong data dependencies which challenge a parallelization strategy.We address these challenges with a block-based MR algorithm, where arbitrary cuts in the underlying octree are possible. This allows for a parallelization on distributed-memory machines via the Message Passing Interface (MPI). We obtain neighbor relations by a simple bit-logic in a modified Morton Order. The block-based concept allows for a modular setup of the source code framework in which the building blocks of the algorithm, such as the choice of the Riemann solver or the reconstruction stencil, are interchangeable without loss of parallel performance. We present the capabilities of the modular framework with a range of test cases and scaling analysis with effective resolutions beyond one billion cells using O(104) cores.

Temporally adaptive conservative scheme for unsteady compressible flow

An adaptive conservative time integration scheme (ACTI) for compressible flow simulations is devised. A finite volume method is used, in which every cell employs a time step size which is some power of two times smaller than a global time step. By first advancing those cells with smaller time steps, this allows to use the already computed accumulated fluxes when the adjacent cells with larger time steps catch up. Conservation and periodic synchronization are thus guaranteed by construction; both important advantages of the new scheme (especially in the presence of shock waves) compared to previously proposed adaptive time stepping methods for compressible flow. Accuracy and computational speedup are demonstrated for challenging 1D and 2D test cases. Comparisons of ACTI results with an analytical solution and with results obtained without ACTI show excellent agreement. Although maximum wave speeds and/or cell sizes vary by several orders of magnitude throughout the computational domains, the CFL numbers employed by ACTI in each cell vary not more than by approximately a factor of two. Without ACTI, on the other hand, the CFL numbers are tiny in most of the cells. Since this results in a dramatic reduction of required flux computations, high speedup gains are achieved due to ACTI; in particular, if only few cells experience a very small time step, while the time step size limitations are less severe for the majority of cells.

An implicit lattice Boltzmann flux solver for simulation of compressible flows

The lattice Boltzmann flux solver (LBFS) is a novel Boltzmann-typed solver under the finite volume framework, which has many advantages over traditional lattice Boltzmann equation (LBE) solvers. However, the existing versions of LBFS are explicit in temporal discretization. Therefore, they are of low computational efficiency and are mostly applied to the simulation of flows around simple geometries. In this paper, an implicit lattice Boltzmann flux solver for compressible flows is proposed, which can efficiently march in time with a large CFL number, and therefore facilitates the simulation of flows around complex geometries in practical engineering applications. Aiming to obtain the time-accurate results of unsteady flows, the implicit LBFS is further combined with the dual time-stepping technique to efficiently simulate unsteady flows, in which a sub-iteration is applied on the pseudo time level. In the sub-iteration, the Jacobian matrix is preconditioned by the lower-upper symmetric Gauss-Seidel (LUSGS) method and solved by the generalized minimal residual method (GMRES). The accuracy and efficiency of the proposed solver are demonstrated by various simulations in a wide Mach number regime. To further demonstrate its performance in practical engineering applications, the proposed solver is also used to simulate the flow over complex aircraft models.

Global and Local Conservation of Linear and Quadratic Invariants in the Simulation of Compressible Flow

Global and Local Conservation of Linear and Quadratic Invariants in the Simulation of Compressible Flow

Modification alpha formalism of Shakura–Sunyaev for the coefficient of turbulent viscosity in an astrophysical disk of finite thickness

In the approximation of one-fluid hydrodynamics, a closed system of Favre-averaged magneto-hydrodynamic equations is formulated, intended for the numerical simulation of compressible turbulent flows of electrically conductive media in the presence of a magnetic field. Special emphasis is paid to the method of obtaining, within the framework of irreversible thermodynamics, the constitutive relations for the turbulent flux heat and the total (kinetic plus magnetic) tensor of turbulent stresses. A new approach to modeling the coefficient of turbulent kinematic viscosity for an astrophysical disk is proposed, which takes into account the influence of an external and generated magnetic field, as well as the processes of convective heat transfer on turbulence in a stratified layer of finite thickness, and thereby modifies the Shakura–Sunyaev alpha formalism developed by for a thin disk and widely used in the astrophysical literature.

A non-oscillatory face-centred finite volume method for compressible flows

This work presents a face-centred finite volume (FCFV) paradigm for the simulation of compressible flows. The FCFV method defines the unknowns at the face barycentre and uses a hybridisation procedure to eliminate all the degrees of freedom inside the cells. In addition, Riemann solvers are defined implicitly within the expressions of the numerical fluxes. The resulting methodology provides first-order accurate approximations of the conservative quantities, i.e. density, momentum and energy, as well as of the viscous stress tensor and of the heat flux, without the need of any gradient reconstruction procedure. Hence, the FCFV solver preserves the accuracy of the approximation in presence of distorted and highly stretched cells, providing a solver insensitive to mesh quality. In addition, FCFV is capable of constructing non-oscillatory approximations of sharp discontinuities without resorting to shock capturing or limiting techniques. For flows at low Mach number, the method is robust and is capable of computing accurate solutions in the incompressible limit without the need of introducing specific pressure correction strategies. A set of 2D and 3D benchmarks of external flows is presented to validate the methodology in different flow regimes, from inviscid to viscous laminar flows, from transonic to subsonic incompressible flows, demonstrating its potential to handle compressible flows in realistic scenarios.

Oscillation-free nodal discontinuous spectral element method for the simulation of compressible multicomponent flows

The nodal discontinuous spectral element method (DSEM) is improved to avoid pressure and velocity oscillations from material interfaces in the solution of multi-component compressible flows. A novel double flux method compatible with the DSEM scheme is developed to account for variations in the specific heats ratio at material interfaces. The double flux method consists of a modified Harten-Lax-van Leer-Contact (HLLC) Riemann solver and an averaging technique to freeze the temperature-dependent specific heats ratio during the solution update to remove oscillations. The double flux model is combined with a conservative method of calculating the convective flux function to avoid unnecessary energy conservation error at shock fronts. Numerical diffusion terms are also introduced for capturing a shock wave, contact discontinuity, and material interface, without inducing oscillations at the material interface. A series of one- and two-dimensional tests are conducted to evaluate the efficacy of the numerical scheme in capturing sharp discontinuities, as well as vortices, without introducing significant conservation error. According to the benchmark results, the new approach maintains the pressure, velocity, and temperature equilibriums across material interfaces. Artificial diffusion is added to the equation in a controlled manner that enables a higher Courant-Friedrichs-Lewy (CFL) number for high polynomial orders. The present scheme's performance is validated by comparing the spatiotemporal evolution of characteristic points in the shock-bubble interaction and Richtmyer-Meshokov's instability problems.

A Semi-Lagrangian Godunov-Type Method without Numerical Viscosity for Shocks

One of the most important and complex effects in compressible fluid flow simulation is a shock-capturing mechanism. Numerous high-resolution Euler-type methods have been proposed to resolve smooth flow scales accurately and to capture the discontinuities simultaneously. One of the disadvantages of these methods is a numerical viscosity for shocks. In the shock, the flow parameters change abruptly at a distance equal to the mean free path of a gas molecule, which is much smaller than the cell size of the computational grid. Due to the numerical viscosity, the aforementioned Euler-type methods stretch the parameter change in the shock over few grid cells. We introduce a semi-Lagrangian Godunov-type method without numerical viscosity for shocks. Another well-known approach is a method of characteristics that has no numerical viscosity and uses the Riemann invariants or solvers for water hammer phenomenon modeling, but in its formulation the convective terms are typically neglected. We use a similar approach to solve the one-dimensional adiabatic gas dynamics equations, but we split the equations into parts describing convection and acoustic processes separately, with corresponding different time steps. When we are looking for the solution to the one-dimensional problem of the scalar hyperbolic conservation law by the proposed method, we additionally use the iterative Godunov exact solver, because the Riemann invariants are non-conserved for moderate and strong shocks in an ideal gas. The proposed method belongs to a group of particle-in-cell (PIC) methods; to the best of the author’s knowledge, there are no similar PIC numerical schemes using the Riemann invariants or the iterative Godunov exact solver. This article describes the application of the aforementioned method for the inviscid Burgers’ equation, adiabatic gas dynamics equations, and the one-dimensional scalar hyperbolic conservation law. The numerical analysis results for several test cases (e.g., the standard shock-tube problem of Sod, the Riemann problem of Lax, the double expansion wave problem, the Shu–Osher shock-tube problem) are compared with the exact solution and Harten’s data. In the shock for the proposed method, the flow properties change instantaneously (with an accuracy dependent on the grid cell size). The iterative Godunov exact solver determines the accuracy of the proposed method for flow discontinuities. In calculations, we use the iteration termination condition less than 10−5 to find the pressure difference between the current and previous iterations.

Two-stage fourth-order gas kinetic solver-based compact subcell finite volume method for compressible flows on triangular meshes

To meet the demand for complex geometries and high resolutions of small-scale flow structures, a two-stage fourth-order subcell finite volume (SCFV) method combining the gas-kinetic solver (GKS) with subcell techniques for compressible flows on (unstructured) triangular meshes was developed to improve the compactness and efficiency. Compared to the fourth-order GKS-based traditional finite volume (FV) method, the proposed method realizes compactness effectively by subdividing each cell into a set of subcells or control volumes (CVs) and selecting only face-neighboring cells for high-order compact reconstruction. Because a set of CVs in a main cell share the same reconstruction, it is more efficient than traditional FV-GKS, where the solution polynomial on each CV needs to be separately reconstructed. Unlike in the single-stage third-order SCFV-GKS, both accuracy and efficiency are improved significantly by two-stage fourth-order temporal discretization, for which only a second-order gas distribution function is needed to simplify the construction of the flux function and reduce computational costs. For viscous flows, it is not necessary to compute the viscous term with GKS. Compared to the fourth-stage Runge–Kutta method, one half of the stage is saved for achieving fourth-order time accuracy, which also helps to improve the efficiency. Therefore, a new high-order method with compactness, efficiency, and robustness is proposed by combining the SCFV method with the two-stage gas-kinetic flux. Several benchmark cases were tested to demonstrate the performance of the method in compressible flow simulations.

Parallel unstructured finite volume lattice Boltzmann method for high-speed viscid compressible flows

Based on the double distribution function Boltzmann-BGK equations, a cell-centered finite volume lattice Boltzmann method on unstructured grids for high-speed viscid compressible flows is presented. In the equations, the particle distribution function is introduced on the basis of the D2Q17 circular function, and its corresponding total energy distribution function is adopted. In the proposed method, the advective term is evaluated by Roe’s flux-difference splitting scheme, and a limiter is used to prevent the generation of oscillations. The distribution functions on the interface are calculated by piecewise linear reconstruction, in which the gradient is computed by the least-squares approach. In order to do large-scale simulations, a parallel algorithm is illustrated. The present method is validated by a flow around the NACA0012 airfoil and a flow past a circular cylinder at high Mach numbers. The results agree well with the published results, which demonstrate that the present method is an efficient numerical method for high-speed viscid compressible flows. The parallel performance results show that the proposed parallel algorithm achieves 90% parallel efficiency on 4800 cores for a problem with [Formula: see text] unstructured triangle cells, which shows the potential to perform fast and high-fidelity simulations of large-scale high-speed viscid compressible flows in complicated computational domains.

Heterogeneous CPU+GPU parallelization for high-accuracy scale-resolving simulations of compressible turbulent flows on hybrid supercomputers

A heterogeneous parallel algorithm for simulation of compressible turbulent flows and its portable software implementation are presented. The underlying numerical method is based on a family of higher accuracy edge-based reconstruction schemes on unstructured mixed-element meshes. The proposed parallel solution can engage a large number of computing devices of most of the existing computing architectures used in modern supercomputers, including manycore CPUs and GPUs. It is capable of co-execution on both CPUs and accelerators simultaneously. The multilevel parallel algorithm combines: MPI for distributing workload among hybrid cluster nodes and between devices inside nodes; OpenMP for manycore CPUs and other supporting devices, such as Intel Xeon Phi; OpenCL for massively-parallel accelerators, such as GPUs of various vendors, including NVIDIA, AMD, Intel. The main focus is on the adaptation of the numerical method and its computational algorithm to the stream processing parallel paradigm. The very limited device memory inherent in GPU computing is also taken into account. A detailed description of the parallel algorithm is presented, as well as the techniques used for its efficient parallel implementation. Special attention is paid to implicit time integration with its linear solver and calculation of convective fluxes and viscous terms. The use of mixed floating-point precision and overlapping communications and computations is also discussed. Parallel performance is demonstrated in practical applications on different kinds of supercomputers using up to 10 thousand cores and multiple GPUs of comparable overall performance.

Finite element formulation for compressible multiphase flows and its application to pyroclastic gravity currents

The numerical modeling of compressible multiphase flows is of high interest for several engineering applications. In this work, we focus on the study of pyroclastic flows arising from volcanic eruptive events. An accurate evaluation of the effects of this multiphase flow is of crucial importance for the Civil Protection for the preservation of urban settlements in volcanic areas. In this work, we propose a Finite Element formulation for the simulation of pyroclastic flows at conditions of thermal and kinetic equilibrium. This analysis belongs to the class of advection dominated problems, which are known to suffer from numerical instabilities. The required stabilization is provided by using a Variational Multiscale Method, typically used for monophase compressible flows and extended here for the first time to multiphase flows. The stabilized formulation is validated with benchmark problems for compressible flows, which are solved for both monophase and multiphase cases. A throughout comparison of the numerical results of the two different flows is also presented. Moreover, the numerical formulation is applied to the simulation of representative cases of pyroclastic flows considering large-scale computational domains and realistic material properties and initial thermal-kinematic conditions. The numerical analyses presented show the accuracy of the proposed method for the simulation of compressible multiphase flows and its suitability for risk assessment studies of urban settlements prone to be affected by pyroclastic gravity currents.

A novel high-order scheme for numerical simulation of wake flow over helicopter rotors in hover

Accurate prediction of tip vortices is crucial for predicting the hovering performance of a helicopter rotor. A new high-order scheme (we call it WENO-K) proposed by our research group is employed to minimize numerical dissipation and extended to numerical simulation of unsteady compressible viscous flows dominated by tip vortices over hovering rotors. WENO-K is referred to as an adaptively optimized WENO scheme with Gauss-Kriging reconstruction, and its advantage is to reduce dissipation in smooth regions of flow while preserving high-resolution around discontinuities. Here WENO-K scheme is adopted to reconstruct left and right state values within the Roe Riemann solver updating the inviscid fluxes on a structured dynamic overset grid. To minimize the accuracy loss for high-order reconstruction on artificial boundaries of overset grid, a method of multilayer fringes is proposed to carry out interpolation between background grid and blade grid. Massively parallel computing considering automatic load balance on averagely partitioned overset grid is developed to reduce the wall-clock time of an unsteady simulation. Numerical results for Caradonna-Tung (C-T) rotor in hover at the conditions of subsonic and transonic tip Mach numbers show that the thrust coefficient error for the result of WENO-K scheme is no more than 3%. Compared with WENO-JS scheme, WENO-K scheme achieves about 40% improvement on accuracy of predicting rotor thrust with only 4.1% extra computational cost. More importantly, WENO-K scheme can capture more sophisticated unsteady flow structures and resolve tip vortices to a larger wake age with an increment of about 270° compared to WENO-JS scheme.

Development and verification of meshless diffuse approximate method for simulation of compressible flow between parallel plates

A meshless numerical model is developed to simulate single-phase, Newtonian, compressible flow in the Cartesian coordinate system. The coupled set of partial differential equations, i.e., mass conservation, momentum conservation, energy conservation, and equation of state is solved by using Diffuse Approximate Method (DAM) and Pressure Implicit with Splitting of Operators (PISO) pressure correction algorithm on an irregular node arrangement. DAM is structured by using the second-order polynomial basis functions and the Gaussian weight function, leading to the weighted least squares approximation on overlapping sub-domains. Implicit time discretization is performed for the predictor step of PISO, while in the corrector steps the equations are discretized explicitly. The numerical model is validated for flow between parallel plates with helium obeying ideal gas law. The solver’s accuracy is assessed by investigating the shape of the Gaussian weight and the number of the nodes in the local subdomains. The calculated velocity, temperature and pressure fields are compared with the Finite Volume Method (FVM) results obtained by OpenFOAM software and show a reasonably good agreement.

Data-based estimation and simulation of compressible pulsating flow with reverse-flow through an orifice

Highly compressible pulsating flows are often encountered in devices where knowledge of the flow rate is required but elimination of pulsations is not an option. The current work is a continuation of a previous investigation that characterized the orifice discharge coefficient Cd as a function of dimensionless groups based on pulsation characteristics. The experimental apparatus has been rebuilt in the current work to mitigate temperature and vibration problems, allowing pressure and ΔP measurements to be made very close to the test section with 159-mm of nylon tubing. Data was acquired for 77 operating conditions spanning a range of pulsation frequencies, mass flow rates and system pressures. They confirm previously reported low Cd's in 0.20 range (calculated from time-averaged pressures) at some high-pressure low-flow operating conditions. Computational Fluid Dynamics (CFD) simulations of 12 of these data points suggest that the low Cd's result from reverse flow. Flow direction changed several times during each pulsation cycle closely tracking the orifice ΔP. A ‘core-and-sheath’ phenomena was observed for reverse-flow operating conditions: a positive core flow surrounded by a sheath of negative flow transitioned to a negative core and positive sheath several times during each pulsation cycle. The simulations also suggested that velocity profiles at the orifice stay stable and similar to steady-state profiles except for periods of rapid transitions. Based on these results a data-based quasi-steady method of estimating pulsating flow has been proposed. A pair of forward and reverse flow Cd's chosen by the data are used to predict instantaneous forward and reverse flows using the steady-state orifice discharge equation for compressible flow. The instantaneous values are then summed up over the pulsation cycle to estimate average mass flow rate. Average prediction errors were within 6%. A previously proposed method where regression was used to model Cd as a function of dimensionless groupings was shown to produce similar results. Both methods are designed to extract information from experimental data in order to overcome theoretical limitations and experimental error. The data is available upon request for further understanding of the flow physics.

A discontinuous Petrov-Galerkin method for compressible Navier-Stokes equations in three dimensions

Application of a Discontinuous Petrov-Galerkin (DPG) method for simulation of compressible viscous flows in three dimensions is presented. The approach enables construction of stable schemes for problems with a small perturbation parameter. The main idea of the method is a weak formulation with a relaxed interelement continuity of the solution. The formulation satisfies the inf-sup condition with the stability constant independent of the small perturbation parameter, which here is the viscosity constant for the compressible Navier-Stokes equations. The DPG discrete formulation uses the specially designed so-called optimal test functions. They do not compromise the inf-sup stability of the continuous formulation. DPG does not use any artificial dissipation for the compressible Navier-Stokes equations. Being a residual minimization method it has got a built-in a posteriori error estimation which allows for mesh adaptivity leading to resolving reliable viscous fluxes, the major difficult task in simulations of viscous flows. We illustrate the method with a few steady state laminar solutions.

An improved third-order HWCNS for compressible flow simulation on curvilinear grids

Due to the very high requirements on the quality of computational grids, stability property and computational efficiency, the application of high-order schemes to complex flow simulation is greatly constrained. In order to solve these problems, the third-order hybrid cell-edge and cell-node weighted compact nonlinear scheme (HWCNS3) is improved by introducing a new nonlinear weighting mechanism. The new scheme uses only the central stencil to reconstruct the cell boundary value, which makes the convergence of the scheme more stable. The application of the scheme to Euler equations on curvilinear grids is also discussed. Numerical results show that the new HWCNS3 achieves the expected order in smooth regions, captures discontinuities sharply without obvious oscillation, has higher resolution than the original one and preserves freestream and vortex on curvilinear grids.

Unsteady flow simulation of compressible turbulent flow in dual-bell nozzle with movement of extendible section from its initial to working position

The use of extendible nozzles in propulsion systems is one of the ways to increase the geometric degree of expansion of the nozzle and the thrust. Numerical simulation of a supersonic turbulent flow of a viscous compressible gas in an extendible nozzle is considered. The simulations performed take into account the movement of the extendible part of the nozzle from the its initial to the working position. The Reynolds-averaged Navier–Stokes equations and the equations of the SST turbulence model discretized on moving meshes are used for numerical calculations. The unsteady flowfield formed when the nozzle is brought to the working position is studied, and the topological features of the flow (separation, mixing, reattachment, secondary flow), and the spatial and temporal distributions of flow quantities are investigated. The distributions of pressure at different times and the axial force applied to the inner walls of the nozzle during the outflow of a supersonic underexpanded jet from the nozzle are discussed. The distributions of flow quantities computed with unsteady and quasi-state formulations of the problem are compared.

High-order gas-kinetic scheme with parallel computation for direct numerical simulation of turbulent flows

The performance of high-order gas-kinetic scheme (HGKS) has been investigated for the direct numerical simulation (DNS) of isotropic compressible turbulence up to the supersonic regime [9]. Due to the multi-scale nature and coupled temporal-spatial evolution process, HGKS provides a valid tool for the numerical simulation of compressible turbulent flow. Based on the domain decomposition and message passing interface (MPI), a parallel HGKS code is developed for large-scale computation in this paper. The standard tests from the nearly incompressible flow to the supersonic one, including Taylor-Green vortex problem, turbulent channel flow and isotropic compressible turbulence, are presented to validate the parallel scalability, efficiency, accuracy and robustness of parallel implementation. The performance of HGKS for the nearly incompressible turbulence is comparable with the high-order finite difference scheme, including the resolution of flow structure and efficiency of computation. Based on the accuracy of the numerical solution, the numerical dissipation of the scheme in the turbulence simulation is quantitatively evaluated. Meanwhile, based on the kinetic formulation HGKS shows advantage for supersonic turbulent flow simulation with its accuracy and robustness. The current work demonstrates the capability of HGKS as a powerful DNS tool from the low speed to supersonic turbulence study, which is less reported under the framework of finite volume scheme.

Comparison of Boundary Integral and Volume-of-Fluid methods for compressible bubble dynamics

The Boundary Integral Method (BIM) has been widely applied to simulate oscillating bubbles, for its high efficiency and accuracy. A conventional BIM assumes the fluid surrounding the bubble to be inviscid and incompressible. Wang & Blake (J. Fluid Mech., 659, 2010, 191–224) proposed an improved model for bubbles in a weakly compressible flow, which is referred to as CBIM. In this study, an all-Mach method (AMM) implemented in the free software program Basilisk for the simulation of compressible multiphase flows, and using a geometric Volume-of-Fluid (VoF), is employed to study and estimate the accuracy of BIM and CBIM at different Mach numbers. First, for a spherical bubble, an extended Rayleigh-Plesset equation, CBIM and AMM give very close results when Ma≲0.3. However, a deviation between these three schemes gradually becomes evident as Ma increases from 0.3 to 0.6. Second, for the nonspherical deformation of a bubble close to a wall, the results obtained from CBIM and AMM show many similarities, including the evolution of the nonspherical bubble morphology, jet impact velocity, and impact pressure on the wall. Apart from the liquid compressibility, the gas inertia/density is found to be another factor that may affect the applicability of CBIM. In addition, we compare the CBIM and BIM results against an experiment of a spark-generated cavitation bubble, in which the liquid compressibility is found to play a vital role. From the perspective of engineering applications, BIM can reproduce the main features of the bubble dynamics in the first cycle if the initial conditions are set properly. The new findings provide a reference for research of bubble dynamics in both fundamental and applied problems.