Stiff Source Research Articles

Staggered Semi-Implicit Hybrid Finite Volume/Finite Element Schemes for Turbulent and Non-Newtonian Flows

This paper presents a new family of semi-implicit hybrid finite volume/finite element schemes on edge-based staggered meshes for the numerical solution of the incompressible Reynolds-Averaged Navier–Stokes (RANS) equations in combination with the k−ε turbulence model. The rheology for calculating the laminar viscosity coefficient under consideration in this work is the one of a non-Newtonian Herschel–Bulkley (power-law) fluid with yield stress, which includes the Bingham fluid and classical Newtonian fluids as special cases. For the spatial discretization, we use edge-based staggered unstructured simplex meshes, as well as staggered non-uniform Cartesian grids. In order to get a simple and computationally efficient algorithm, we apply an operator splitting technique, where the hyperbolic convective terms of the RANS equations are discretized explicitly at the aid of a Godunov-type finite volume scheme, while the viscous parabolic terms, the elliptic pressure terms and the stiff algebraic source terms of the k−ε model are discretized implicitly. For the discretization of the elliptic pressure Poisson equation, we use classical conforming P1 and Q1 finite elements on triangles and rectangles, respectively. The implicit discretization of the viscous terms is mandatory for non-Newtonian fluids, since the apparent viscosity can tend to infinity for fluids with yield stress and certain power-law fluids. It is carried out with P1 finite elements on triangular simplex meshes and with finite volumes on rectangles. For Cartesian grids and more general orthogonal unstructured meshes, we can prove that our new scheme can preserve the positivity of k and ε. This is achieved via a special implicit discretization of the stiff algebraic relaxation source terms, using a suitable combination of the discrete evolution equations for the logarithms of k and ε. The method is applied to some classical academic benchmark problems for non-Newtonian and turbulent flows in two space dimensions, comparing the obtained numerical results with available exact or numerical reference solutions. In all cases, an excellent agreement is observed.

Cerebrospinal fluid dynamics coupled to the global circulation in holistic setting: Mathematical models, numerical methods and applications.

This paper presents a mathematical model of the global, arterio‐venous circulation in the entire human body, coupled to a refined description of the cerebrospinal fluid (CSF) dynamics in the craniospinal cavity. The present model represents a substantially revised version of the original Müller‐Toro mathematical model. It includes one‐dimensional (1D), non‐linear systems of partial differential equations for 323 major blood vessels and 85 zero‐dimensional, differential‐algebraic systems for the remaining components. Highlights include the myogenic mechanism of cerebral blood regulation; refined vasculature for the inner ear, the brainstem and the cerebellum; and viscoelastic, rather than purely elastic, models for all blood vessels, arterial and venous. The derived 1D parabolic systems of partial differential equations for all major vessels are approximated by hyperbolic systems with stiff source terms following a relaxation approach. A major novelty of this paper is the coupling of the circulation, as described, to a refined description of the CSF dynamics in the craniospinal cavity, following Linninger et al. The numerical solution methodology employed to approximate the hyperbolic non‐linear systems of partial differential equations with stiff source terms is based on the Arbitrary DERivative Riemann problem finite volume framework, supplemented with a well‐balanced formulation, and a local time stepping procedure. The full model is validated through comparison of computational results against published data and bespoke MRI measurements. Then we present two medical applications: (i) transverse sinus stenoses and their relation to Idiopathic Intracranial Hypertension; and (ii) extra‐cranial venous strictures and their impact in the inner ear circulation, and its implications for Ménière's disease.

ChemNODE: A neural ordinary differential equations framework for efficient chemical kinetic solvers

Solving for detailed chemical kinetics remains one of the major bottlenecks for computational fluid dynamics simulations of reacting flows using a finite-rate-chemistry approach. This has motivated the use of neural networks to predict stiff chemical source terms as functions of the thermochemical state of the combustion system. However, due to the nonlinearities and multi-scale nature of combustion, the predicted solution often diverges from the true solution when these machine learning models are coupled with a computational fluid dynamics solver. This is because these approaches minimize the error during training without guaranteeing successful integration with ordinary differential equation solvers. In the present work, a novel neural ordinary differential equations approach to modeling chemical kinetics, termed as ChemNODE, is developed. In this machine learning framework, the chemical source terms predicted by the neural networks are integrated during training, and by computing the required derivatives, the neural network weights are adjusted accordingly to minimize the difference between the predicted and ground-truth solution. A proof-of-concept study is performed with ChemNODE for homogeneous autoignition of hydrogen-air mixture over a range of composition and thermodynamic conditions. It is shown that ChemNODE accurately captures the chemical kinetic behavior and reproduces the results obtained using the detailed kinetic mechanism at a fraction of the computational cost.

Existence and uniqueness of generalized solutions to hyperbolic systems with linear fluxes and stiff sources

Motivated by the modeling of boiling two-phase flows, we study systems of balance laws with a source term defined as a discontinuous function of the unknown. Due to this discontinuous source term, the classical theory of partial differential equations (PDEs) is not sufficient here. Restricting to a simpler system with linear fluxes, a notion of generalized solution is developed. An important point in the construction of a solution is that the curve along which the source jumps, which we call the boiling curve, must never be tangent to the characteristics. This leads to exhibit sufficient conditions which ensure the existence and uniqueness of a solution in two different situations: first when the initial data is smooth and such that the boiling curve is either overcharacteristic or subcharacteristic; then with discontinuous initial data in the case of Riemann problems. A numerical illustration is given in this last case.

A universal centred high-order method based on implicit Taylor series expansion with fast second order evolution of spatial derivatives

In this paper, a centred universal high-order finite volume method for solving hyperbolic balance laws is presented. The scheme belongs to the family of ADER methods where the Generalized Riemann Problems (GRP) is a building block. The solution to these problems is carried through an implicit Taylor series expansion, which allows the scheme to work very well for stiff source terms. The series expansion is high order for the state and requires the evolution in time of spatial derivatives. A Taylor expansion of second order based on a linearization around the resolved state, is proposed for evolving spatial derivatives. A von Neumann stability analysis is carried out to investigate the range of CFL values for which stability and accuracy are balanced. The scheme implements a centred, low dissipation approach for dealing with the advective part of the system which profits from small CFL values. Numerical tests demonstrate that the present scheme can solve, efficiently, hyperbolic balance laws in both conservative and non-conservative form. The scheme is proven to be well-balanced, the C-property, a classical assessment of well-balancing of non-conservative schemes, is numerically demonstrated. An empirical convergence rate assessment shows that the expected theoretical orders of accuracy are achieved up to the fifth order.

Asymptotic preserving schemes on conical unstructured 2D meshes

AbstractIn this article, we consider the hyperbolic heat equation. This system is linear hyperbolic with stiff source terms and satisfies a diffusion limit. Some finite volume numerical schemes have been proposed which reproduce this diffusion limit. Here, we extend such schemes, originally defined on polygonal meshes, to conical meshes (using rational quadratic Bezier curves). We obtain really new schemes that do not reduce to the polygonal version when the conical edges tend to straight lines. Moreover, these schemes can handle curved unstructured meshes so that geometric error on initial data representation is reduced and geometry of the domain is improved. Extra flux coming from conical edge (through his midedge point) has a deep impact on the stabilization when compared with the original polygonal scheme. Cross‐stencil phenomenon of polygonal scheme has disappeared, and issue of positivity for the diffusion problem (although unresolved on distorted mesh and/or with varying cross‐section) has been in some sense improved.

Maximum-Principle-Preserving Local Discontinuous Galerkin Methods for Allen-Cahn Equations

In this paper, we study the classical Allen-Cahn equations and investigate the maximum-principle-preserving (MPP) techniques. The Allen-Cahn equation has been widely used in mathematical models for problems in materials science and fluid dynamics. It enjoys the energy stability and the maximum-principle. Moreover, it is well known that the Allen-Cahn equation may yield thin interface layer, and nonuniform meshes might be useful in the numerical solutions. Therefore, we apply the local discontinuous Galerkin (LDG) method due to its flexibility on h-p adaptivity and complex geometry. However, the MPP LDG methods require slope limiters, then the energy stability may not be easy to obtain. In this paper, we only discuss the MPP technique and use numerical experiments to demonstrate the energy decay property. Moreover, due to the stiff source given in the equation, we use the conservative modified exponential Runge-Kutta methods and thus can use relatively large time step sizes. Thanks to the conservative time integration, the bounds of the unknown function will not decay. Numerical experiments will be given to demonstrate the good performance of the MPP LDG scheme.

Simulation of non-Newtonian viscoplastic flows with a unified first order hyperbolic model and a structure-preserving semi-implicit scheme

We discuss the applicability of a unified hyperbolic model for continuum fluid and solid mechanics to modeling non-Newtonian flows and in particular to modeling the stress-driven solid-fluid transformations in flows of viscoplastic fluids, also called yield-stress fluids. In contrast to the conventional approaches relying on the non-linear viscosity concept of the Navier-Stokes theory and representation of the solid state as an infinitely rigid non-deformable solid, the solid state in our theory is deformable and the fluid state is considered rather as a “melted” solid via a certain procedure of relaxation of tangential stresses similar to Maxwell’s visco-elasticity theory. The model is formulated as a system of first-order hyperbolic partial differential equations with possibly stiff non-linear relaxation source terms. The computational strategy is based on a staggered semi-implicit scheme which can be applied in particular to low-Mach number flows as usually required for flows of non-Newtonian fluids. The applicability of the model and numerical scheme is demonstrated on a few standard benchmark test cases such as Couette, Hagen-Poiseuille, and lid-driven cavity flows. The numerical solution is compared with analytical or numerical solutions of the Navier-Stokes theory with the Herschel-Bulkley constitutive model for nonlinear viscosity.

Very high-order asymptotic-preserving schemes for hyperbolic systems of conservation laws with parabolic degeneracy on unstructured meshes

In this paper, we consider the numerical approximation of hyperbolic systems of conservation laws with stiff source terms and parabolic degeneracy in the asymptotic limit. We are more precisely interested in the design of high-order asymptotic-preserving schemes on unstructured meshes. Our approach is based on a very simple modification of the numerical flux associated with the usual HLL scheme and boils down to a sharp control of the underlying numerical diffusion. The strategy allows to capture the correct asymptotic parabolic behaviour and to preserve the high-order accuracy also in the asymptotic limit. Numerical experiments are proposed to illustrate these properties.

Modified shallow water model for viscous fluids and positivity preserving numerical approximation

Shallow water equations are widely used in the simulation of those geophysical flows for which the flow horizontal length scale is much greater than the vertical one. Inspired by the example of lava flows, we consider here a modified model with an additional transport equation for a scalar quantity (e.g., temperature), and the derivation of the shallow water equations from depth-averaging the Navier-Stokes equations is presented. The assumption of constant vertical profiles for some of the model variables is relaxed allowing the presence of vertical profiles, and it follows that the non-linearity of the flux terms results in the introduction of appropriate shape coefficients.The space discretization of the resulting system of hyperbolic partial differential equations is obtained with a modified version of the finite volume central-upwind scheme introduced by Kurganov and Petrova in 2007. The time discretization is based on an implicit-explicit Runge-Kutta method which couples properly the hyperbolic part and the stiff source terms, avoiding the use of a very small time step; the use of complex arithmetic increases accuracy in the implicit treatment of stiff terms. The whole scheme is proved to preserve the positivity of flow thickness and the stationary steady-states.Some numerical experiments are performed to validate the proposed method and to show the incidence on the numerical solutions of shape coefficients introduced in the model.

Stability and bifurcations in low inertia PV rich power networks

The present trend of inserting increasingly more solar photovoltaic (PV) sources into the electricity grid leads to a significant reduction in mechanical inertia. Inertia represents energy reserve in the grid, that inherently and instantaneously supports frequency stability. Therefore, recent studies have re-examined the frequency stability of the grid. Most of these consider scenarios where the participation from grid-following inverters remains relatively low. More importantly, these studies include an infinite bus in their analyses. The infinite bus acts as a very stiff voltage source which has a significantly stabilising influence, the more so as all control refers back to the infinite bus. Here, the frequency stability of a grid with significant generation from grid-tied PV inverters is considered without reference to an infinite bus. In the proposed simplified grid model, all the synchronous generators (SGs) are lumped as one large SG complete with classically operating controls. Similarly, all the PV generators are lumped as one large, quasi-instantaneous, non-linear power source. In this simplified network, the feasible operating region is identified using bifurcation techniques. It transpires that the stable operating region is bounded by a locus of Hopf bifurcations linked (inter-alia) to the fraction of power generated by grid-tied PV inverters.

Toward the implementation of a multi-component framework in a density-based flow solver for handling chemically reacting flows

PurposeThe purpose of this study is to present and demonstrate a numerical method for solving chemically reacting flows. These are important for energy conversion devices, which rely on chemical reactions as their operational mechanism, with heat generated from the combustion of the fuel, often gases, being converted to work.Design/methodology/approachThe numerical study of such flows requires the set of Navier-Stokes equations to be extended to include multiple species and the chemical reactions between them. The numerical method implemented in this study also accounts for changes in the material properties because of temperature variations and the process to handle steep spatial fronts and stiff source terms without incurring any numerical instabilities. An all-speed numerical framework is used through simple low-dissipation advection upwind splitting (SLAU) convective scheme, and it has been extended in a multi-component species framework on the in-house density-based flow solver. The capability of solving turbulent combustion is also implemented using the Eddy Dissipation Concept (EDC) framework and the recent k-kl turbulence model.FindingsThe numerical implementation has been demonstrated for several stiff problems in laminar and turbulent combustion. The laminar combustion results are compared from the corresponding results from the Cantera library, and the turbulent combustion computations are found to be consistent with the experimental results.Originality/valueThis paper has extended the single gas density-based framework to handle multi-component gaseous mixtures. This paper has demonstrated the capability of the numerical framework for solving non-reacting/reacting laminar and turbulent flow problems. The all-speed SLAU convective scheme has been extended in the multi-component species framework, and the turbulent model k-kl is used for turbulent combustion, which has not been done previously. While the former method provides the capability of solving for low-speed flows using the density-based method, the later is a length-scale-based method that includes scale-adaptive simulation characteristics in the turbulence modeling. The SLAU scheme has proven to work well for unsteady flows while the k-kL model works well in non-stationary turbulent flows. As both these flow features are commonly found in industrially important reacting flows, the convection scheme and the turbulence model together will enhance the numerical predictions of such flows.

Deep Residual Networks for Flamelet/progress Variable Tabulation with Application to a Piloted Flame with Inhomogeneous Inlet

ABSTRACT In this work, a deep neural network is presented which is trained on flamelet/progress variable (FPV) tables and validated in a combustion large eddy simulation (LES) of the Sydney/Sandia flame with inhomogeneous inlets. Using data scaling and transformation techniques, as well as novel network architectures as the residual skip connection we are able to store all combustion relevant quantities in one single network, thus reducing effectively the memory footprint compared to the FPV tables, while keeping data retrieval times similar to table interpolation. The accuracy of the deep neural networks (DNN) is compared to its FPV counterpart and achieves excellent agreement. The DNN is also able to accurately predict the stiff progress variable source term in the thin reaction zone. Finally, the DNN thermochemistry representation is validated in combustion simulations of a 2D laminar premixed flame and a 3D LES of the Sydney/Sandia piloted jet flame with inhomogeneous inlet. The DNN results are in very good agreement with the conventional tabulated FPV simulations and show a promising way to efficiently reduce the storage size of high dimensional pretabulated thermochemical state space in reactive flow simulations. Abbreviations: ANN: Artificial neural network; API: Application programming interface; CFD: Computational fluid dynamics; CFL: Courant-Friedrichs-Lewy; CPU: Central processing unit; DL: Deep learning; DNN: Deep neural network; FGM: Flamelet generated manifold; FLOPS: Floating point operation per second; FPV: Flamelet/progress variable; GPU: Graphics processing unit; HPC: High-performance computing; ILDM: Intrinsic low-dimensional manifold; LES: Large eddy simulation; MAE: Mean absolute error; MPI: Message passing interface; MSE: Mean squared error; PDF: Probability density function; SGS: Subgrid scale; SOM: Self-organizing map; TCI: Turbulence–chemistry interaction.

Boundary treatment of implicit-explicit Runge-Kutta method for hyperbolic systems with source terms

In this paper, we develop a high order finite difference boundary treatment method for the implicit-explicit (IMEX) Runge-Kutta (RK) schemes solving hyperbolic systems with possibly stiff source terms on a Cartesian mesh. The main challenge is how to obtain the solutions at ghost points resulting from the wide stencil of the interior high order scheme. We address this problem by combining the idea of using the RK schemes at the boundary and an inverse Lax-Wendroff (ILW) procedure in [29,30]. Specifically, we only apply the ILW procedure in the starting stage of the RK method. In the intermediate stages, we solve out the intermediate solutions as well as their first-order spatial derivatives at the boundary by using the RK schemes, which are then used to compute solutions at ghost points by Taylor expansions. Our method is different from the widely used approach for the explicit RK schemes by imposing boundary conditions at intermediate stages [14,29,30] which does not apply to IMEX schemes. In addition, the intermediate boundary conditions are available for explicit RK schemes only up to third order while our method applies to IMEX and explicit RK schemes of arbitrary order. The stability and the desired accuracy of our boundary treatment with IMEX RK method up to fourth order are demonstrated numerically through 1D examples and 2D reactive Euler equations.

Chemistry reduction using machine learning trained from non-premixed micro-mixing modeling: Application to DNS of a syngas turbulent oxy-flame with side-wall effects

A chemistry reduction approach based on machine learning is proposed and applied to direct numerical simulation (DNS) of a turbulent non-premixed syngas oxy-flame interacting with a cooled wall. The training and the subsequent application of artificial neural networks (ANNs) rely on the processing of ‘thermochemical vectors’ composed of species mass fractions and temperature (ANN input), to predict the corresponding chemical sources (ANN output). The training of the ANN is performed aside from any flow simulation, using a turbulent non-adiabatic non-premixed micro-mixing based canonical problem with a reference detailed chemistry. Heat-loss effects are thus included in the ANN training. The performance of the ANN chemistry is then tested a-posteriori in a two-dimensional DNS against the detailed mechanism and a reduced mechanism specifically developed for the operating conditions considered. Then, three-dimensional DNS are performed either with the ANN or the reduced chemistry for additional a-posteriori tests. The ANN reduced chemistry achieves good agreement with the Arrhenius-based detailed and reduced mechanisms, while being in terms of CPU cost 25 times faster than the detailed mechanism and 3 times faster than the reduced mechanism when coupled with DNS. The major potential of the method lies both in its data driven character and in the handling of the stiff chemical sources. The former allows for easy implementation in the context of automated generation of case-specific reduced chemistry. The latter avoids the Arrhenius rates calculation and also the direct integration of stiff chemistry, both leading to a significant CPU time reduction.

MPI-AMRVAC: A parallel, grid-adaptive PDE toolkit

We report on the latest additions to our open-source, block-grid adaptive framework MPI-AMRVAC, which is a general toolkit for especially hyperbolic/parabolic partial differential equations (PDEs). Applications traditionally focused on shock-dominated, magnetized plasma dynamics described by either Newtonian or special relativistic (magneto)hydrodynamics, but its versatile design easily extends to different PDE systems. Here, we demonstrate applications covering any-dimensional scalar to system PDEs, with e.g. Korteweg–de Vries solutions generalizing early findings on soliton behavior, shallow water applications in round or square pools, hydrodynamic convergence tests as well as challenging computational fluid and plasma dynamics applications. The recent addition of a parallel multigrid solver opens up new avenues where also elliptic constraints or stiff source terms play a central role. This is illustrated here by solving several multi-dimensional reaction–diffusion-type equations. We document the minimal requirements for adding a new physics module governed by any nonlinear PDE system, such that it can directly benefit from the code flexibility in combining various temporal and spatial discretization schemes. Distributed through GitHub, MPI-AMRVAC can be used to perform 1D, 1.5D, 2D, 2.5D or 3D simulations in Cartesian, cylindrical or spherical coordinate systems, using parallel domain-decomposition, or exploiting fully dynamic block quadtree-octree grids.

Numerical Resolution of a Three Temperature Plasma Model

This paper is devoted to the numerical approximation of a three temperature plasma model: one for the ions, one for the electrons and one for the radiation (photons). A reformulation of the model is proposed that allows to build a convex combination-based scheme that unconditionally satisfies a maximum principle, at each sub-iteration of the non-linear iterative process. This yields a very robust scheme that can handle stiff source terms. In addition, the methodology is extended to include the contribution of a radiative flux (Rosseland diffusion approximation) and electronic and ionic conductivities (Spitzer–Harm diffusion approximation). Several numerical results are carried out to demonstrate the interest of the numerical approach.

High Order Asymptotic Preserving Deferred Correction Implicit-Explicit Schemes for Kinetic Models

This work aims to extend the residual distribution (RD) framework to stiff relaxation problems. The RD is a class of schemes which is used to solve hyperbolic system of partial differential equations. Up to our knowledge, it was used only for systems with mild source terms, such as gravitation problems or shallow water equations. What we propose is an IMEX (implicit--explicit) version of the residual distribution schemes, that can resolve stiff source terms, without refining the discretization up to the stiffness scale. This can be particularly useful in various models, where the stiffness is given by topological or physical quantities, e.g. multiphase flows, kinetic models, viscoelasticity problems. Moreover, the provided scheme is able to catch different relaxation scales automatically, without losing accuracy. The scheme is asymptotic preserving and this guarantees that in the relaxation limit, we recast the expected macroscopic behaviour. To get a high order accuracy, we use an IMEX time discretization combined with a Deferred Correction (DeC) procedure, while naturally RD provides high order space discretization. Finally, we show some numerical tests in 1D and 2D for stiff systems of equations.

Transfverter: Imbuing Transformer-Like Properties in an Interlink Converter for Robust Control of a Hybrid AC–DC Microgrid

In a hybrid ac–dc microgrid, stiff voltage sources may appear in either the dc or ac subgrids, which gives rise to multiple operation modes as power dispatch changes. This creates a challenge for designing the interlink converter between the ac and dc subgrids, since the different modes require different interlink controls. To solve this problem, this paper proposes the concept of a transfverter inspired by how transformers link ac grids. Like a transformer, a transfverter can react to the presence of stiff voltage sources on either the dc or ac side and reflect the “stiffness” and voltage stabilizing capability to the other side. A back-to-back converter with droop control is used as the underlying technology to implement this concept. A novel optimization method called model bank synthesis is proposed to find control parameters for the interlink converter that offer the best controller performance across the different microgrid modes without requiring mode changing of the controller. The effectiveness of the proposed solution is validated through both simulation and experiments.

Advanced Control Methods for Power Converters in Distributed Generation Systems and Microgrids

The twenty-two papers in this special section focus on flexible control of power converters which serve as interfaces between the distributed generation (DG) units and the legacy alternating current (ac) grid or the ac or direct current (dc) microgrid (MG), is the key to realization of high penetration of renewable energy in a safe and stable fashion. When connected to the ac legacy grid, these power converters need to provide ancillary services such as frequency and voltage support, harmonic compensation, as well as synthetic inertia emulation. Another emerging solution is to interface the DG units with the ac legacy grid through an intermediate entity called anMG. MG can be based either on ac and dc architecture and can work in both stand-alone and grid-connected modes. Since it is responsible for multiple power converters, an MG has higher operational flexibility than individual units.However, due to a lack of stiff voltage reference source and natural inertia, control of MGs is generally more challenging than control of individual grid-connected power converters.