Dispersed multiphase flows are dynamic systems of two or more phases, at least one of which consists of fine particles in any state of matter dispersed in a fluid. They are ubiquitous in nature and technology. A few examples of their occurence in nature are – from large to small scales – weather patterns, ocean dynamics and blood flow. On the other side, propulsion systems, power generation, heat exchange and manufacturing technologies highlight the importance of multiphase flows in technical and industrial applications. Numerical modeling of such systems is particularly challenging in terms of the mathematical description of the dispersed phase and the associated physical phenomena such as turbulence and various mechanisms of mass, momentum and energy exchange.

QBMMs are based on the macroscopic description in terms of a set of moments associated with the NDF. The governing moment transport equations can be derived from the underlying PBE. Presuming a single internal coordinate \(\xi\) in the presence of advection and diffusion with drift function \(\mathcal{A}(\xi)\) and diffusivity \(\mathcal{D}(\xi)\), the evolution of the NDF \(n\) over time \(t\) is governed by the PBE

            \(\frac{\partial n(\xi, t)}{\partial t} = - \frac{\partial }{\partial \xi} \left[ \mathcal{A}(\xi) \, n(\xi, t) \right] \; + \; \frac{\partial^2 }{\partial \xi^2} \left[ \mathcal{D}(\xi) \, n(\xi, t) \right]\).

Applying the definition of the \(k\)th moment \(m_k = \int_\Omega \xi^k n(\xi) \mathrm{d} \xi\) to the PBE yields the \(k\)th moment transport equation

            \(\frac{\mathrm{d} m_k(t)}{\mathrm{d} t} = - k \int_\Omega \xi^{k-1} \, \mathcal{A}(\xi) \, n(\xi, t) \mathrm{d} \xi \; + \; k(k-1) \int_\Omega \xi^{k-2} \, \mathcal{D}(\xi) \, n(\xi, t) \mathrm{d} \xi\).

Thus, without knowledge of n , the moment equations are unclosed. This closure problem is solved by QBMMs employing quadrature rules, i.e.

            \(\frac{\mathrm{d} m_k(t)}{\mathrm{d} t} \approx - k \sum\limits_j w_j \, \xi_j^{k-1} \, \mathcal{A}(\xi_j) \; + \; k(k-1) \sum\limits_j w_j \, \xi_j^{k-2} \, \mathcal{D}(\xi_j)\).

If the quadrature nodes \(\xi_j\) and weigths \(w_j\) are those of the Gaussian quadrature associated with the known moment sequence this method is referred to as the quadrature method of moments (QMOM) [1], which is the basis of various derived methods, see [2]. All have in common that the quadrature formula can be computed solely from moment information. We extended the QMOM making use of the anti-Gaussian quadrature formulas due to [3].

In cooperation with the Geophysical Fluid Dynamics group of the FU Berlin, we are investigating the theoretical framework of Shockless Explosion Combustion (SEC), which represents a completely new operating concept for gas turbines. Core of the SEC process is the controlled quasi-homogeneous auto-ignition of an air-fuel mixture, which should lead to a dynamic approximation of an energetically more favorable constant volume combustion (CVC) [2]. The numerical investigation of the flows occuring here involves very small time and length scales, which also requires significant effort from the computational fluid dynamics (CFD) solver. Therefore, we are developing our own CFD solver and pay special attention to efficient numerical methods and parallelization of the simulations. In the following, we give a rough overview of the numerical methods we use and focus especially on methods for adaptive mesh refinement (AMR).

The general basis for the numerical solution of problems we discuss here are the multi-species reactive Euler Equations. These are given in their three-dimensional conservative Form as
    \begin{gather}
            \frac{\partial \mathbf{U}}{\partial t}
            + \frac{\partial \mathbf{F}(\mathbf{U})}{\partial x}
            + \frac{\partial \mathbf{G}(\mathbf{U})}{\partial y}
            + \frac{\partial \mathbf{H}(\mathbf{U})}{\partial z}
            = \mathbf{S}(\mathbf{U}) \,, \label{eq:eulerEquation} \\
             \mathbf{U} = \begin{bmatrix}
                \rho \\ \rho u \\ \rho v \\ \rho w \\ \rho E \\ \rho Y_s
                \end{bmatrix} ,\,
             \mathbf{F} = \begin{bmatrix}
                \rho u \\ \rho u^2+p \\ \rho uv \\ \rho uw \\ u(\rho E+p) \\ \rho u Y_s
            \end{bmatrix} ,\,
             \mathbf{G} = \begin{bmatrix}
                \rho v \\ \rho uv \\ \rho v^2+p \\ \rho vw \\ v(\rho E+p) \\ \rho v Y_s
            \end{bmatrix} ,\,
             \mathbf{H} = \begin{bmatrix}
                \rho w\\ \rho uw \\ \rho vw \\ \rho w^2+p \\ w(\rho E+p) \\ \rho w Y_s
            \end{bmatrix} ,\, \label{eq:eulerComponents}
    \end{gather}
where \(\rho\) is the density, \(\mathbf{u}=(u,v,w)\) the velocity, \(E=e+||\mathbf{u}||^2/2\) the total energy, \(e\) the internal energy, \(p\) the pressure and \(Y_s,\,s=1,\ldots,N_s-1\) are the mass fractions of the reactive species that are normalized to sum to unity. The equations above are closed by the caloric perfect gas approximation \(p=\rho(\gamma-1)e\), where \(\gamma=c_p/c_v\) is the ratio of specific heat. The term \(\mathbf{S}(\mathbf{U})\) is a placeholder for any possible source term that accounts, e.g. external forces, chemical reactions or any injection of mass, momentum, energy or species. To compute the numerical solution to the above equations, we utilize fully conservative finite volume schemes discretised on block-structured Cartesian grids. Due to the usage of the external C++ library AMReX [5] in our CFD solver, we can locally refine the grid in regions of interest through adaptive mesh refinement (AMR) techniques [1]. To ensure stability of any small or irregular cells that occur due embedding complex geometries on the structured grids, we use a conservative cut-cell method [3].