Realizable high-order finite-volume schemes for quadrature-based moment methods

Dilute gas–particle flows can be described by a kinetic equation containing terms for spatial transport, gravity, fluid drag and particle–particle collisions. However, direct numerical solution of kinetic equations is often infeasible because of the large number of independent variables. An alternative is to reformulate the problem in terms of the moments of the velocity distribution. Recently, a quadrature-based moment method was derived for approximating solutions to kinetic equations. The success of the new method is based on a moment-inversion algorithm that is used to calculate non-negative weights and abscissas from the moments. The moment-inversion algorithm does not work if the moments are non-realizable, which might lead to negative weights. It has been recently shown [14] that realizability is guaranteed only with the 1st-order finite-volume scheme that has an inherent problem of excessive numerical diffusion. The use of high-order finite-volume schemes may lead to non-realizable moments. In the present work, realizability of the finite-volume schemes in both space and time is discussed for the 1st time. A generalized idea for developing realizable high-order finite-volume schemes for quadrature-based moment methods is presented. These finite-volume schemes give remarkable improvement in the solutions for a certain class of problems. It is also shown that the standard Runge–Kutta time-integration schemes do not guarantee realizability. However, realizability can be guaranteed if strong stability-preserving (SSP) Runge–Kutta schemes are used. Numerical results are presented on both Cartesian and triangular meshes.     

A high order kinetic flux-vector splitting method for the reduced five-equation model of compressible two-fluid flows

We present a high order kinetic flux-vector splitting (KFVS) scheme for the numerical solution of a conservative interface-capturing five-equation model of compressible two-fluid flows. This model was initially introduced by Wackers and Koren (2004) [21]. The flow equations are the bulk equations, combined with mass and energy equations for one of the two fluids. The latter equation contains a source term in order to account for the energy exchange. We numerically investigate both one- and two-dimensional flow models. The proposed numerical scheme is based on the direct splitting of macroscopic flux functions of the system of equations. In two space dimensions the scheme is derived in a usual dimensionally split manner. The second order accuracy of the scheme is achieved by using MUSCL-type initial reconstruction and Runge–Kutta time stepping method. For validation, the results of our scheme are compared with those from the high resolution central scheme of Nessyahu and Tadmor [14]. The accuracy, efficiency and simplicity of the KFVS scheme demonstrate its potential for modeling two-phase flows.     

Nonuniform time-step Runge–Kutta discontinuous Galerkin method for Computational Aeroacoustics

With many superior features, Runge–Kutta discontinuous Galerkin method (RKDG), which adopts Discontinuous Galerkin method (DG) for space discretization and Runge–Kutta method (RK) for time integration, has been an attractive alternative to the finite difference based high-order Computational Aeroacoustics (CAA) approaches. However, when it comes to complex physical problems, especially the ones involving irregular geometries, the time step size of an explicit RK scheme is limited by the smallest grid size in the computational domain, demanding a high computational cost for obtaining time accurate numerical solutions in CAA. For computational efficiency, high-order RK method with nonuniform time step sizes on nonuniform meshes is developed in this paper. In order to ensure correct communication of solutions on the interfaces of grids with different time step sizes, the values at intermediate-stages of the Runge–Kutta time integration on the elements neighboring such interfaces are coupled with minimal dissipation and dispersion errors. Based upon the general form of an explicit p-stage RK scheme, a linear coupling procedure is proposed, with details on the coefficient matrices and execution steps at common time-levels and intermediate time-levels. Applications of the coupling procedures to Runge–Kutta schemes frequently used in simulation of fluid flow and acoustics are given, including the third-order TVD scheme, and low-storage low dissipation and low dispersion (LDDRK) schemes. In addition, an analysis on the stability of coupling procedures on a nonuniform grid is carried out. For validation, numerical experiments on one-dimensional and two-dimensional problems are presented to illustrate the stability and accuracy of proposed nonuniform time-step RKDG scheme, as well as the computational benefits it brings. Application to a one-dimensional nonlinear problem is also investigated.     

A high-order gas-kinetic Navier–Stokes flow solver

The foundation for the development of modern compressible flow solver is based on the Riemann solution of the inviscid Euler equations. The high-order schemes are basically related to high-order spatial interpolation or reconstruction. In order to overcome the low-order wave interaction mechanism due to the Riemann solution, the temporal accuracy of the scheme can be improved through the Runge–Kutta method, where the dynamic deficiencies in the first-order Riemann solution is alleviated through the sub-step spatial reconstruction in the Runge–Kutta process. The close coupling between the spatial and temporal evolution in the original nonlinear governing equations seems weakened due to its spatial and temporal decoupling. Many recently developed high-order methods require a Navier–Stokes flux function under piece-wise discontinuous high-order initial reconstruction. However, the piece-wise discontinuous initial data and the hyperbolic-parabolic nature of the Navier–Stokes equations seem inconsistent mathematically, such as the divergence of the viscous and heat conducting terms due to initial discontinuity. In this paper, based on the Boltzmann equation, we are going to present a time-dependent flux function from a high-order discontinuous reconstruction. The theoretical basis for such an approach is due to the fact that the Boltzmann equation has no specific requirement on the smoothness of the initial data and the kinetic equation has the mechanism to construct a dissipative wave structure starting from an initially discontinuous flow condition on a time scale being larger than the particle collision time. The current high-order flux evaluation method is an extension of the second-order gas-kinetic BGK scheme for the Navier–Stokes equations (BGK-NS). The novelty for the easy extension from a second-order to a higher order is due to the simple particle transport and collision mechanism on the microscopic level. This paper will present a hierarchy to construct such a high-order method. The necessity to couple spatial and temporal evolution nonlinearly in the flux evaluation can be clearly observed through the numerical performance of the scheme for the viscous flow computations.     

Adaptive Runge–Kutta discontinuous Galerkin methods using different indicators: One-dimensional case

In this paper, we systematically investigate adaptive Runge–Kutta discontinuous Galerkin (RKDG) methods for hyperbolic conservation laws with different indicators which were based on the troubled cell indicators studied by Qiu and Shu [J. Qiu, C.-W. Shu, A comparison of troubled-cell indicators for Runge–Kutta discontinuous Galerkin mehtods using weighted essentially non-osillatory limiters, SIAM J. Sci. Comput. 27 (2005) 995–1013]. The emphasis is on comparison of the performance of adaptive RKDG method using different indicators, with an objective of obtaining efficient and reliable indicators to obtain better performance for adaptive computation to save computational cost. Both h-version and r-version adaptive methods are considered in the paper. The idea is to first identify “troubled cells” by different troubled-cell indicators, namely those cells where limiting might be needed and discontinuities might appear, then adopt an adaptive approach in these cells. A detailed numerical study in one-dimensional case is performed, addressing the issues of efficiency (less CPU cost and more accurate), non-oscillatory property, and resolution of discontinuities.     

Optimal Runge–Kutta smoothers for the p-multigrid discontinuous Galerkin solution of the 1D Euler equations

This work presents a family of original Runge–Kutta methods specifically designed to be effective relaxation schemes in the numerical solution of the steady state solution of purely advective problems with a high-order accurate discontinuous Galerkin space discretization and a p-multigrid solution algorithm. The design criterion for the construction of the Runge–Kutta methods here developed is different form the one traditionally used to derive optimal Runge–Kutta smoothers for the h-multigrid algorithm, which are designed in order to provide a uniform damping of the error modes in the high-frequency range only. The method here proposed is instead designed in order to provide a variable amount of damping of the error modes over the entire frequency spectrum. The performance of the proposed schemes is assessed in the solution of the steady state quasi one-dimensional Euler equations for two test cases of increasing difficulty. Some preliminary results showing the performance for multidimensional applications are also presented.     

High order conservative Lagrangian schemes with Lax–Wendroff type time discretization for the compressible Euler equations

In this paper, we explore the Lax–Wendroff (LW) type time discretization as an alternative procedure to the high order Runge–Kutta time discretization adopted for the high order essentially non-oscillatory (ENO) Lagrangian schemes developed in 3 and 5. The LW time discretization is based on a Taylor expansion in time, coupled with a local Cauchy–Kowalewski procedure to utilize the partial differential equation (PDE) repeatedly to convert all time derivatives to spatial derivatives, and then to discretize these spatial derivatives based on high order ENO reconstruction. Extensive numerical examples are presented, for both the second-order spatial discretization using quadrilateral meshes [3] and third-order spatial discretization using curvilinear meshes [5]. Comparing with the Runge–Kutta time discretization procedure, an advantage of the LW time discretization is the apparent saving in computational cost and memory requirement, at least for the two-dimensional Euler equations that we have used in the numerical tests.     

The global version of the gyrokinetic turbulence code GENE

The understanding and prediction of transport due to plasma microturbulence is a key open problem in modern plasma physics, and a grand challenge for fusion energy research. Ab initio simulations of such small-scale, low-frequency turbulence are to be based on the gyrokinetic equations, a set of nonlinear integro-differential equations in reduced (five-dimensional) phase space. In the present paper, the extension of the well-established and widely used gyrokinetic code GENE [F. Jenko, W. Dorland, M. Kotschenreuther, B.N. Rogers, Electron temperature gradient driven turbulence, Phys. Plasmas 7 (2000) 1904–1910] from a radially local to a radially global (nonlocal) version is described. The necessary modifications of both the basic equations and the employed numerical methods are detailed, including, e.g., the change from spectral methods to finite difference and interpolation techniques in the radial direction and the implementation of sources and sinks. In addition, code verification studies and benchmarks are presented.     

Dynamic analysis of double-row self-aligning ball bearings due to applied loads,internal clearance,surface waviness and number of balls

In this paper, a three degrees of freedom (dof) model was established for a double-row self-aligning ball bearing (SABB) system, and was applied to study the dynamic behavior of the system during starting process and constant speed rotating process. A mathematical model was developed concerning stiffness and damping characteristics of the bearing, as well as three-dimensional applied load, rotor centrifugal force, etc. Balls and races were all considered as nonlinear springs, and the contact force between ball and race was calculated based on classic Hertzian elastic contact deformation theory and deformation compatibility theory. The changes of each ball?s contact force and loaded angle of each row were taken into account. In order to solve the nonlinear dynamical equilibrium equations of the system, these equations were rewritten as differential equations and the fourth order Runge–Kutta method was used to solve the equations iteratively. In order to verify accuracy of the dynamical model and correctness of the numerical solution method, a kind of SABB-BRF30 was chosen for case studies. The effects of several important governing parameters, such as radial and axial applied loads, normal internal, inner and outer races waviness, and number of balls were investigated. These parametric studies led to a complete characterization of the shaft-bearing system vibration transmission. The research provided a theoretical reference for new type bearing design, shaft-bearing system kinetic analysis, optimal design, etc.     

High-order parabolic beam approximation for aero-optics

The parabolic beam equations are solved using high-order compact differences for the Laplacians and Runge–Kutta integration along the beam path. The solution method is verified by comparison to analytical solutions for apertured beams and both constant and complex index of refraction. An adaptive 4th-order Runge–Kutta using an embedded 2nd-order method is presented that has demonstrated itself to be very robust. For apertured beams, the results show that the method fails to capture near aperture effects due to a violation of the paraxial approximation in that region. Initial results indicate that the problem appears to be correctable by successive approximations. A preliminary assessment of the effect of turbulent scales is undertaken using high-order Lagrangian interpolation. The results show that while high fidelity methods are necessary to accurately capture the large scale flow structure, the method may not require the same level of fidelity in sampling the density for the index of refraction. The solution is used to calculate a phase difference that is directly compared with that commonly calculated via the optical path difference. Propagation through a supersonic boundary layer shows that for longer wavelengths, the traditional method to calculate the optical path is less accurate than for shorter wavelengths. While unlikely to supplant more traditional methods for most aero-optics applications, the current method can be used to give a quantitative assessment of the other methods as well as being amenable to the addition of more physics.     

Splitting multisymplectic integrators for Maxwell’s equations

In the paper, we describe a novel kind of multisymplectic method for three-dimensional (3-D) Maxwell’s equations. Splitting the 3-D Maxwell’s equations into three local one-dimensional (LOD) equations, then applying a pair of symplectic Runge–Kutta methods to discretize each resulting LOD equation, it leads to splitting multisymplectic integrators. We say this kind of schemes to be LOD multisymplectic scheme (LOD-MS). The discrete conservation laws, convergence, dispersive relation, dissipation and stability are investigated for the schemes. Theoretical analysis shows that the schemes are unconditionally stable, non-dissipative, and of first order accuracy in time and second order accuracy in space. As a reduction, we also consider the application of LOD-MS to 2-D Maxwell’s equations. Numerical experiments match the theoretical results well. They illustrate that LOD-MS is not only efficient and simple in coding, but also has almost all the nature of multisymplectic integrators.     

On nonlinear cascades and resonances in the outer magnetosphere

The paper addresses nonlinear phenomena that control the interaction between plasma flow (solar wind) and magnetic barrier (magnetosphere). For the first time we demonstrate that the dominant solar wind kinetic energy: (i) excites boundary resonances and their harmonics which modulate plasma jets under the bow shock; (ii) produces discrete three-wave cascades, which could merge into a turbulent-like one; (iii) jet produced cascades provide the effective anomalous plasma transport inside and out of the magnetosphere; (iv) intermittency and multifractality characteristics for the statistic properties of jets result in a super-ballistic turbulent transport regime. Our results could be considered as suggestive for the space weather predictions, for turbulent cascades in different media and for the laboratory plasma confinement (e.g., for fusion devices).     

Effect of one-way clutch on the nonlinear vibration of belt-drive systems with a continuous belt model

This study focuses on the nonlinear steady-state response of a belt-drive system with a one-way clutch. A dynamic model is established to describe the rotations of the driving pulley, the driven pulley, and the accessory shaft. Moreover, the model considers the transverse vibration of the translating belt spans for the first time in belt-drive systems coupled with a one-way clutch. The excitation of the belt-drive system is derived from periodic fluctuation of the driving pulley. In automotive systems, this kind of fluctuation is induced by the engine firing harmonic pulsations. The derived coupled discrete–continuous nonlinear equations consist of integro-partial-differential equations and piece-wise ordinary differential equations. Using the Galerkin truncation, a set of nonlinear ordinary differential equations is obtained from the integro-partial-differential equations. Applying the Runge–Kutta time discretization, the time histories of the dynamic response are numerically solved for the driven pulley and the accessory shaft and the translating belt spans. The resonance areas of the coupled belt-drive system are determined using the frequency sweep. The effects of the one-way clutch on the belt-drive system are studied by comparing the frequency–response curves of the translating belt with and without one-way clutch device. Furthermore, the results of 2-term and 4-term Galerkin truncation are compared to determine the numerical convergence. Moreover, parametric studies are conducted to understand the effects of the system parameters on the nonlinear steady-state response. It is concluded that one-way clutch not only decreases the resonance amplitude of the driven pulley and shaft's rotational vibration, but also reduces the resonance region of the belt's transverse vibration.     

Local DG method using WENO type limiters for convection–diffusion problems

The local discontinuous Galerkin (LDG) method is a spatial discretization procedure for convection–diffusion equations, which employs useful features from high resolution finite volume schemes, such as the exact or approximate Riemann solvers serving as numerical fluxes and limiters, which is termed as Runge–Kutta LDG (RKLDG) when TVD Runge–Kutta method is applied for time discretization. It has the advantage of flexibility in handling complicated geometry, h-p adaptivity, and efficiency of parallel implementation and has been used successfully in many applications. However, the limiters used to control spurious oscillations in the presence of strong shocks are less robust than the strategies of essentially non-oscillatory (ENO) and weighted ENO (WENO) finite volume and finite difference methods. In this paper, we investigated RKLDG methods with WENO and Hermite WENO (HWENO) limiters for solving convection–diffusion equations on unstructured meshes, with the goal of obtaining a robust and high order limiting procedure to simultaneously obtain uniform high order accuracy and sharp, non-oscillatory shock transition. Numerical results are provided to illustrate the behavior of these procedures.     

A general strategy for the optimization of Runge–Kutta schemes for wave propagation phenomena

We analyze optimized explicit Runge–Kutta schemes (RK) for computational aeroacoustics, and wave propagation phenomena in general. Exploiting the analysis developed in [S. Pirozzoli, Performance analysis and optimization of finite-difference schemes for wave propagation problems, J. Comput. Phys. 222 (2007) 809–831], we rigorously evaluate the performance of several time integration schemes in terms of appropriate error and cost metrics, and provide a general strategy to design Runge–Kutta methods tailored for specific applications. We present families of optimized second- and third-order Runge–Kutta schemes with up to seven stages, and describe their implementation in the framework of Williamson’s 2N$2N$-storage formulation [J.H. Williamson, Low-storage Runge–Kutta schemes, J. Comput. Phys. 35 (1980) 48–56]. Numerical simulations of the 1D linear advection equation and of the 2D linearized Euler equations are performed to demonstrate the validity of the theory and to quantify the improvement provided by optimized schemes.     

Efficient computation of compressible and incompressible flows

The combination of explicit Runge–Kutta time integration with the solution of an implicit system of equations, which in earlier work demonstrated increased efficiency in computing compressible flow on highly stretched meshes, is extended toward conditions where the free stream Mach number approaches zero. Expressing the inviscid flux Jacobians in terms of Mach number, an artificial speed of sound as in low Mach number preconditioning is introduced into the Jacobians, leading to a consistent formulation of the implicit and explicit parts of the discrete equations. Besides extension to low Mach number flows, the augmented Runge–Kutta/Implicit method allowed the admissible Courant–Friedrichs–Lewy number to be increased from O(1 0 0) to O(1 0 0 0). The implicit step introduced into the Runge–Kutta framework acts as a preconditioner which now addresses both, the stiffness in the discrete equations associated with highly stretched meshes, and the stiffness in the analytical equations associated with the disparity in the eigenvalues of the inviscid flux Jacobians. Integrated into a multigrid algorithm, the method is applied to efficiently compute different cases of inviscid flow around airfoils at various Mach numbers, and viscous turbulent airfoil flow with varying Mach and Reynolds number. Compared to well tuned conventional methods, computation times are reduced by half an order of magnitude.     

A high-order low-dispersion symmetry-preserving finite-volume method for compressible flow on curvilinear grids

A new high-order finite-volume method is presented that preserves the skew symmetry of convection for the compressible flow equations. The method is intended for Large-Eddy Simulations (LES) of compressible turbulent flows, in particular in the context of hybrid RANS–LES computations. The method is fourth-order accurate and has low numerical dissipation and dispersion. Due to the finite-volume approach, mass, momentum, and total energy are locally conserved. Furthermore, the skew-symmetry preservation implies that kinetic energy, sound-velocity, and internal energy are all locally conserved by convection as well. The method is unique in that all these properties hold on non-uniform, curvilinear, structured grids. Due to the conservation of kinetic energy, there is no spurious production or dissipation of kinetic energy stemming from the discretization of convection. This enhances the numerical stability and reduces the possible interference of numerical errors with the subgrid-scale model. By minimizing the numerical dispersion, the numerical errors are reduced by an order of magnitude compared to a standard fourth-order finite-volume method.     

A conservative level-set based method for compressible solid/fluid problems on fixed grids

A three-dimensional Eulerian method is presented for simulating dynamic systems comprising multiple compressible solid and fluid components where internal boundaries are tracked using level-set functions. Aside from the interface interaction calculation within mixed cells, each material is treated independently and the governing constitutive laws solved using a conservative finite volume discretisation based upon the solution of Riemann problems to determine the numerical fluxes. The required reconstruction of mixed cell volume fractions and cut cell geometries is presented in detail using the level-set fields. High-order accuracy is achieved by incorporating the weighted-essentially non-oscillatory (WENO) method and Runge–Kutta time integration. A model for elastoplastic solid dynamics is employed formulated using the tensor of elastic deformation gradients permitting the equations to be written in divergence form. The scheme is demonstrated using selected one-dimensional initial value problems for which exact solutions are derived, a two-dimensional void collapse, and a three-dimensional simulation of a confined explosion.     

The gyrokinetic water-bag modeling in toroidal geometry

The present paper addresses the gyrokinetic water-bag model in toroidal geometry. The previous works were focused on the water-bag concept in magnetized cylindrical plasmas. Here we report on the possibility to improve the water-bag model by taking into account the curvature and gradient drifts. After a presentation of the model, a local linear analysis with some approximations is performed. Interchange and ion temperature gradient instabilities are examined with this new gyro-water-bag model in order to show its ability and its theoretical interest in describing kinetic instabilities in toroidal geometry.     

Adaptive hierarchic transformations for dynamically p-enriched slope-limiting over discontinuous Galerkin systems of generalized equations

We study a family of generalized slope limiters in two dimensions for Runge–Kutta discontinuous Galerkin (RKDG) solutions of advection-diffusion systems. We analyze the numerical behavior of these limiters applied to a pair of model problems, comparing the error of the approximate solutions, and discuss each limiter’s advantages and disadvantages. We then introduce a series of coupled p-enrichment schemes that may be used as standalone dynamic p-enrichment strategies, or may be augmented via any in the family of variable-in-p slope limiters presented.