Computational singular perturbation with non-parametric tabulation of slow manifolds for time integration of stiff chemical kinetics

This paper presents a novel tabulation strategy for the adaptive numerical integration of chemical kinetics using the computational singular perturbation (CSP) method. The strategy stores and reuses CSP quantities required to filter out fast dissipative processes, resulting in a non-stiff chemical source term. In particular, non-parametric regression on low-dimensional slow invariant manifolds (SIMs) in the chemical state space is used to approximate the CSP vectors spanning the fast chemical subspace and the associated fast chemical time-scales. The relevant manifold and its dimension varies depending on the local number of exhausted modes at every location in the chemical state space. Multiple manifolds are therefore tabulated, corresponding to different numbers of exhausted modes (dimensions) and associated radical species. Non-parametric representations are inherently adaptive, and rely on efficient approximate-nearest-neighbor queries. As the CSP information is only a function of the non-radical species in the system and has relatively small gradients in the chemical state space, tabulation occurs in a lower-dimensional state space and at a relatively coarse level, thereby improving scalability to larger chemical mechanisms. The approach is demonstrated on the simulation of homogeneous constant pressure H₂–air and CH₄–air ignition, over a range of initial conditions. For CH₄–air, results are shown that outperform direct implicit integration of the stiff chemical kinetics while maintaining good accuracy.     

Explicit time-reversible orbit integration in Particle In Cell codes with static homogeneous magnetic field

A new explicit time-reversible orbit integrator for the equations of motion in a static homogeneous magnetic field – called Cyclotronic integrator – is presented. Like Spreiter and Walter’s Taylor expansion algorithm, for sufficiently weak electric field gradients this second order method does not require a fine resolution of the Larmor motion; it has however the essential advantage of being symplectic, hence time-reversible. The Cyclotronic integrator is only subject to a linear stability constraint (ΩΔt < π, Ω being the Larmor angular frequency), and is therefore particularly suitable to electrostatic Particle In Cell codes with uniform magnetic field where Ω is larger than any other characteristic frequency, yet a resolution of the particles’ gyromotion is required. Application examples and a detailed comparison with the well-known (time-reversible) Boris algorithm are presented; it is in particular shown that implementation of the Cyclotronic integrator in the kinetic codes SCEPTIC and Democritus can reduce the cost of orbit integration by up to a factor of ten.     

Integration of large chemical kinetic mechanisms via exponential methods with Krylov approximations to Jacobian matrix functions

Recent trends in hydrocarbon fuel research indicate that the number of species and reactions in chemical kinetic mechanisms is rapidly increasing in an effort to provide predictive capabilities for fuels of practical interest. In order to cope with the computational cost associated with the time integration of stiff, large chemical systems, a novel approach is proposed. The approach combines an exponential integrator and Krylov subspace approximations to the exponential function of the Jacobian matrix. The components of the approach are described in detail and applied to the ignition of stoichiometric methane–air and iso-octane–air mixtures, here described by two widely adopted chemical kinetic mechanisms. The approach is found to be robust even at relatively large time steps and the global error displays a nominal third-order convergence. The performance of the approach is improved by utilising an adaptive algorithm for the selection of the Krylov subspace size, which guarantees an approximation to the matrix exponential within user-defined error tolerance. The Krylov projection of the Jacobian matrix onto a low-dimensional space is interpreted as a local model reduction with a well-defined error control strategy. Finally, the performance of the approach is discussed with regard to the optimal selection of the parameters governing the accuracy of its individual components.     

Static properties of coupled atomic-molecular Bose-Einstein condensates in modified Gross-Pitaevskii approach

We analyze the stationary state properties of an atomic Bose-Einstein condensate coupled to a molecular condensate via a Raman photoassociation process using Gross-Pitaevskii and modified Gross-Pitaevskii model and compare the results. We find that the static properties depend crucially on the Raman detuning parameter, which can be experimentally varied. We also show how the analytical expressions for densities of atoms and molecules of the hybrid system can be found in the Thomas-Fermi approximation where kinetic energies are not included. We have explored the feasibility of maximum molecular BEC formation by varying the Raman detuning.     

A discrete action principle for electrodynamics and the construction of explicit symplectic integrators for linear,non-dispersive media

In this work, we derive a discrete action principle for electrodynamics that can be used to construct explicit symplectic integrators for Maxwell’s equations. Different integrators are constructed depending on the choice of discrete Lagrangian used to approximate the action. By combining discrete Lagrangians in an explicit symplectic partitioned Runge–Kutta method, an integrator capable of achieving any order of accuracy is obtained. Using the von Neumann stability analysis, we show that the integrators greatly increase the numerical stability and reduce the numerical dispersion compared to other methods. For practical purposes, we demonstrate how to implement the integrators using many features of the finite-difference time-domain method. However, our approach is also applicable to other spatial discretizations, such as those used in finite element methods. Using this implementation, numerical examples are presented that demonstrate the ability of the integrators to efficiently reduce and maintain a minimal amount of numerical dispersion, particularly when the time-step is less than the stability limit. The integrators are therefore advantageous for modeling large, inhomogeneous computational domains.     

Computational Modeling of Boundary Layer Flashback in a Swirling Stratified Flame Using a LES-Based Non-Adiabatic Tabulated Chemistry Approach

When operating under lean fuel–air conditions, flame flashback is an operational safety issue in stationary gas turbines. In particular, with the increased use of hydrogen, the propagation of the flame through the boundary layers into the mixing section becomes feasible. Typically, these mixing regions are not designed to hold a high-temperature flame and can lead to catastrophic failure of the gas turbine. Flame flashback along the boundary layers is a competition between chemical reactions in a turbulent flow, where fuel and air are incompletely mixed, and heat loss to the wall that promotes flame quenching. The focus of this work is to develop a comprehensive simulation approach to model boundary layer flashback, accounting for fuel–air stratification and wall heat loss. A large eddy simulation (LES) based framework is used, along with a tabulation-based combustion model. Different approaches to tabulation and the effect of wall heat loss are studied. An experimental flashback configuration is used to understand the predictive accuracy of the models. It is shown that diffusion-flame-based tabulation methods are better suited due to the flashback occurring in relatively low-strain and lean fuel–air mixtures. Further, the flashback is promoted by the formation of features such as flame tongues, which induce negative velocity separated boundary layer flow that promotes upstream flame motion. The wall heat loss alters the strength of these separated flows, which in turn affects the flashback propensity. Comparisons with experimental data for both non-reacting cases that quantify fuel–air mixing and reacting flashback cases are used to demonstrate predictive accuracy.     

Push it to the limit: comparing periodic and local approaches to density functional theory for intermolecular interactions

ABSTRACT

With the aim of systematically comparing two popular approaches to density functional theory – all-electron calculations with local basis sets, and periodic calculations employing plane wave basis sets and norm-conserving pseudopotentials – we have computed complete-basis binding energies across the S22 set of intermolecular interactions, a dataset consisting of noncovalent interactions of small- and medium-sized molecules containing first- and second-row atoms, using the Troullier-Martins norm-conserving pseudopotentials with SPW92, a local spin-density approximation; and PBE, a generalised gradient approximation. We have found that it is challenging to reach the basis set limit with these periodic calculations; for the methods and systems examined, a minimum vacuum distance of 30?Å between a system and its nearest images is necessary – unless some form of dipole correction is employed – as is a kinetic energy cutoff of at least 80 Ry. The trends in convergence with respect to vacuum size and kinetic energy cutoff are largely independent of the level of density functional approximation employed. A sense of the impact of each hyperparameter on basis set error provides a foundation for ensuring quality calculations in future studies and allows us to quantify the basis set errors incurred in existing studies on similar systems.     

Domain wall fermions with Majorana couplings

We examine the lattice boundary formulation of chiral fermions with either an explicit Majorana mass or a Higgs-Majorana coupling introduced on one of the boundaries. We demonstrate that the low-lying spectrum of the models with an explicit Majorana mass of the order of an inverse lattice spacing is chiral at tree level. Within a mean-field approximation we show that the systems with a strong Higgs-Majorana coupling have a symmetric phase, in which a Majorana mass of the order of an inverse lattice spacing is generated without spontaneous breaking of the gauge symmetry. We argue, however, that the models within such a phase have a chiral spectrum only in terms of the fermions that are singlets under the gauge group. The application of such systems to non-perturbative formulations of supersymmetric and chiral gauge theories is briefly discussed.     

Turbulent flame visualization using direct numerical simulation

Combustion phenomena are of high scientific and technological interest, in particular for energy generation and transportation systems. Direct Numerical Simulations (DNS) have become an essential and well established research tool to investigate the structure of turbulent flames, since they do not rely on any approximate turbulence models. In this work two complementary DNS codes are employed to investigate different types of fuels and flame configurations. The code is π³ is a 3-dimensional DNS code using a low-Mach number approximation. Chemistry is described through a tabulation, using two coordinates to enter a database constructed for example with 29 species and 141 reactions for methane combustion. It is used here to investigate the growth of a turbulent premixed flame in a methane-air mixture (Case 1). The second code,Sider is an explicit three-dimensional DNS code solving the fully compressible reactive Navier-Stokes equations, where the chemical processes are computed using a complete reaction scheme, taking into account accurate diffusion properties. It is used here to compute a hydrogen/air turbulent diffusion flame (Case 2), considering 9 chemical species and 38 chemical reactions.     

The application of the QSSA via reaction lumping for the reduction of complex hydrocarbon oxidation mechanisms

The quasi-steady-state approximation (QSSA) has been widely applied for the purposes of chemical kinetic model reduction. Although it is essentially a low-order approximation, it can be shown to lead to significant reductions in the number of fast variables within a mechanism without significant loss of accuracy for model predictions. Due to the couplings between QSSA expressions, the species are commonly solved for using numerical inner iteration techniques. Therefore, although the stiffness of the model system can be reduced, there is a computational overhead in solving the often nonlinear QSSA equations. Greater computational savings can be made where QSS species can be removed from the chemical model via explicit analytical expressions. In many cases these expressions are equivalent to reaction lumping. Where such reaction lumping can be achieved, a reduced mechanism in standard kinetic form can be developed, which contains new lumped reaction rate coefficients, but leads to the removal of QSS species. This paper describes such an approach for mechanisms describing the oxidation of the hydrocarbon fuels n-heptane and cyclohexane, and shows that significant reductions in both species and reactions can be achieved, leading to substantial computational speed-ups. The resulting schemes clearly demonstrate the main atomic flux patterns within the oxidation process. Patterns related to the time-scales of hydrocarbon radical species within alkane oxidation mechanisms are discussed, as well as the potential significance of non-QSS radicals in determining ignition behaviour.     

Analysis of a few numerical integration methods for the Langevin equation

We study the position recurrence relation of several existing numerical integrators for the Langevin equation and use the modified equation approach to analyse their accuracy. We show that for the harmonic oscillator, the BBK integrator converges weakly with order 1 while the vGB82 and Langevin impulse (LÎ)? integrator converge weakly with order 2. We also study a restricted class of velocity definitions—those that lead to explicit starting procedures. We show that some recurrence relations exact for constant force, can achieve the exact virial relation by a proper definition of velocity, extending the result of Pastor et al. on the analysis of BBK integrators in 1988.     

Tabulated chemistry approaches for laminar flames: Evaluation of flame-prolongation of ILDM and flamelet methods

The present study considers the performance of tabulation methods for numerical simulation of complex chemical kinetics in laminar combusting flows and compares their predictions to results obtained by direct calculation. Two tabulation methods are considered: the Flame Prolongation of Intrinsic low-dimensional manifold (FPI) method and Steady Laminar Flamelet Model (SLFM). The FPI method is of current interest as it is a potentially unifying approach capable of dealing with both premixed and non-premixed flames for gaseous fuels. SLFM tabulation methods are popular for non-premixed flames and form a good basis for comparing the performance of the FPI approach. The performance of each method is also evaluated by comparing the results to the direct simulation of the laminar flames using two chemical kinetic schemes: simplified chemistry involving five species and one reaction and detailed chemistry involving 53 species and 325 reaction steps. As part of the evaluation process, the computational cost of each method is also assessed. The laminar flames considered in this study include: freely propagating laminar premixed flames, a two-dimensional axisymmetric methane–air opposed-jet diffusion flame, and a two-dimensional axisymmetric methane–air co-flow diffusion flame. Both tabulation methods are implemented in a parallel adaptive mesh refinement (AMR) framework for solving the complete set of governing partial differential equations. These equations are solved using a fully-coupled finite-volume formulation on body-fitted multi-block quadrilateral mesh. Significant improvements in terms of reduced computational requirements, as measured by both storage and processing time, are demonstrated for the tabulated methods.     

Why many theories of shock waves are necessary: Kinetic functions,equivalent equations,and fourth-order models

We consider several systems of nonlinear hyperbolic conservation laws describing the dynamics of nonlinear waves in presence of phase transition phenomena. These models admit under-compressive shock waves which are not uniquely determined by a standard entropy criterion but must be characterized by a kinetic relation. Building on earlier work by LeFloch and collaborators, we investigate the numerical approximation of these models by high-order finite difference schemes, and uncover several new features of the kinetic function associated with physically motivated second and third-order regularization terms, especially viscosity and capillarity terms.On one hand, the role of the equivalent equation associated with a finite difference scheme is discussed. We conjecture here and demonstrate numerically that the (numerical) kinetic function associated with a scheme approaches the (analytic) kinetic function associated with the given model – especially since its equivalent equation approaches the regularized model at a higher order. On the other hand, we demonstrate numerically that a kinetic function can be associated with the thin liquid film model and the generalized Camassa–Holm model. Finally, we investigate to what extent a kinetic function can be associated with the equations of van der Waals fluids, whose flux-function admits two inflection points.     

Kinetic Monte Carlo simulations of the partial oxidation of methanol on oxygen-covered Cu(1 1 0)

The partial oxidation of methanol to formaldehyde on oxygen-precovered Cu(1 1 0) has been studied using kinetic Monte Carlo simulations. The rates entering the simulation have been derived from density functional theory calculations within the generalized gradient approximation using transition state theory. We demonstrate that kinetic Monte Carlo simulations are a powerful tool to elucidate the microscopic details of the reaction kinetics on surfaces. Furthermore, the comparison of calculated and measured temperature programmed desorption rates allows a genuine assessment of the calculated barrier heights. We find that some of the calculated barriers and adsorption energies have to be adjusted by up to 0.5 eV in order to reproduce the measured desorption spectra. Possible reasons for the discrepancies between experiment and theory are discussed.     

A new explicit multisymplectic integrator for the Kawahara-type equation

We derive a new multisymplectic integrator for the Kawahara-type equation which is a fully explicit scheme and thus needs less computation cost. Multisympecticity of such scheme guarantees the long-time numerical behaviors. Nu- merical experiments are presented to verify the accuracy of this scheme as well as the excellent performance on invariant preservation for three kinds of Kawahara-type equations.     

A patch that imparts unconditional stability to explicit integrators for Langevin-like equations

This paper extends the results in [8] to stochastic differential equations (SDEs) arising in molecular dynamics. It implements a patch to explicit integrators that consists of a Metropolis–Hastings step. The ‘patched integrator’ preserves the SDE’s equilibrium distribution and is accurate on finite time intervals. As a corollary this paper proves the integrator’s accuracy in estimating finite-time dynamics along an infinitely long solution - a first in molecular dynamics. The paper also covers multiple time-steps, holonomic constraints and scalability. Finally, the paper provides numerical tests supporting the theory.     

A time-domain Discontinuous Galerkin method for mechanical resonator quality factor computations

Numerical simulations are becoming increasingly important in the design of micromechanical resonators, in particular for the prediction of complex frequency response in high quality devices where damping is controlled by anchor losses. Frequency based approaches have been shown to predict these accurately, however, they require the solution of eigenvalue problems or the inversion of Helmholtz-type operators which are known to be very difficult for large-scale iterative solvers. We propose using a time-domain approach instead, where a broadband input signal is propagated through the system with a local explicit time-stepper. This is achieved using a new high-order Discontinuous Galerkin (DG) discretization for the linear elasticity equations, in particular a second-order formulation with Compact DG fluxes and a Runge–Kutta time integrator, where the block-diagonal mass matrices allow for efficient, stable, and accurate time stepping. Our solver scales well on distributed parallel computers, even in three spatial dimension and for large problem sizes. The resulting output signal is analyzed using a well-known filter diagonalization method, which is capable of finding accurate frequencies and quality factors for as little as a hundred periods of data. We validate the properties of our scheme on model problems, and demonstrate the feasibility of our proposed analysis process on two high quality factor disk resonators, using an axisymmetric formulation as well as full three dimensional simulations which is shown to scale well.     

Two-spin dynamic autocorrelations and energy diffusion in the Heisenberg paramagnet

We establish the kinetic equation for the two-spin time-dependent correlation function in the Weiss limit. Considering this equation in the limit of small wavenumbers, we obtain an explicit expression for the energy diffusion coefficient. We also discuss some difficulties connected with overcounting problems of the Weiss approximation.     

Kinetic energy driven pairing in cuprate superconductors

Pairing occurs in conventional superconductors through a reduction of the electronic potential energy accompanied by an increase in kinetic energy. In the underdoped cuprates, optical experiments show that pairing is driven by a reduction of the electronic kinetic energy. Using the dynamical cluster approximation we study superconductivity in the two-dimensional Hubbard model. We find that pairing is indeed driven by the kinetic energy and that superconductivity evolves from an unconventional state with partial spin-charge separation, to a superconducting state with quasiparticle excitations.     

Multi-symplectic variational integrators for nonlinear Schrdinger equations with variable coefficients

In this paper, we propose a variational integrator for nonlinear Schrdinger equations with variable coefficients. It is shown that our variational integrator is naturally multi-symplectic. The discrete multi-symplectic structure of the integrator is presented by a multi-symplectic form formula that can be derived from the discrete Lagrangian boundary function. As two examples of nonlinear Schrdinger equations with variable coefficients, cubic nonlinear Schrdinger equations and Gross–Pitaevskii equations are extensively studied by the proposed integrator. Our numerical simulations demonstrate that the integrator is capable of preserving the mass, momentum, and energy conservation during time evolutions. Convergence tests are presented to verify that our integrator has second-order accuracy both in time and space.