Overview

canoP is a computational fluid dynamic (CFD) code leveraging the p4est library for adaptive mesh refinement (AMR) to simulate multiphase flows. It currently implements:

a monofluid model: for a single phase,
a bifluid5eq model: for two phases with the so-called 5 equations homogeneous equilibrium model (see Padioleau 2021),
a bifluid7eq model: for two phases with the so-called 7 equations Baer-Nunziato model (see Chen 2019).

Generic model

canoP is designed to numerically integrate a set of extensive variables \(\boldsymbol{U}\) (referred as “conservative” from now on) with the following generic model:

\[\frac{\partial \boldsymbol{U}}{\partial t}+\vec{\nabla}\cdot\vec{\boldsymbol{F}}(\boldsymbol{U})- \vec{\nabla}\cdot\left(\boldsymbol{D}(\boldsymbol{U}) \vec{\nabla}\boldsymbol{U}\right) = \boldsymbol{S}(\boldsymbol{U})\]

with \(\vec{\boldsymbol{F}}(\boldsymbol{U})\) a hyperbolic flux function, \(-\boldsymbol{D}(\boldsymbol{U})\vec{\nabla}\boldsymbol{U}\) a diffusion flux function, and \(\boldsymbol{S}(\boldsymbol{U})\) a source term.

Note

Part of the source term \(\boldsymbol{S}(\boldsymbol{U})\) can be written as a non-conservative hyperbolic flux by using e.g. in 1D: \(\boldsymbol{F}^-(\boldsymbol{U}) = \int_{x_i}^{x} \boldsymbol{S}(\boldsymbol{U})dx\) and \(\boldsymbol{F}^+(\boldsymbol{U}) = \int_{x}^{x_{i+1}} \boldsymbol{S}(\boldsymbol{U})dx\) on each side of an interface between two volumes centered at \(x_i\) and \(x_{i+1}\).

All the numerical schemes are using a finite volume discretization by approximating the volume averaged of the conservative variables inside controle volumes. The schemes are also second order in space and time and use spatial reconstruction (the so-called MUSCL reconstruction) on a set of primitive variables \(\boldsymbol{P}(\boldsymbol{U})\) and a half time step estimator (the so-called Hancock predictor).

Monofluid model

The conservative variables contain a density, a 2D or 3D momentum, and a total energy, with (we use Einstein convention):

\[ \begin{align}\begin{aligned}\boldsymbol{U}&=(m, m V_i, m E)\\\boldsymbol{P}(\boldsymbol{U})&=(m, V_i, P, T)\\\boldsymbol{\vec{F}}(\boldsymbol{U})&=(m V_i,m V_i V_j+P \delta_{i,j}, (m E+P)V_i)\\-\boldsymbol{D}(\boldsymbol{U})\vec{\nabla}\boldsymbol{U}&=(0, -\tau_{i,j}, -\kappa \partial_i T),\\\boldsymbol{S}(\boldsymbol{U}) &= (0, m g_i, m g_k V_k)\end{aligned}\end{align} \]

with

\[ \begin{align}\begin{aligned}E&=e+V_k V_k/2,\\\tau_{i,j}&=\zeta\partial_k V_k \delta_{i,j} + \mu(\partial_i V_j + \partial_j V_i -2/3 \partial_k V_k \delta_{i,j}),\end{aligned}\end{align} \]

and an equation of state (EOS) giving the pressure/temperature \((P,T)=f(m,e)\), the thermal conductivity \(\kappa\), and the shear and bulk viscosity \(\mu\) and \(\zeta\).

The EOS for the pressure/temperature is a Noble-Abel Stiffened-Gas EOS parametrized with \((\gamma,P_\infty, b, C_v, q, q^\prime)\) following the convention of Le Métayer & Saurel (2016):

\[ \begin{align}\begin{aligned}P&=(\gamma-1)\frac{e-q}{1/m-b}-\gamma P_\infty\\C_v T& = (1/m-b)\frac{e-q}{\gamma-1}\\s(P,T) &= C_v \ln{\frac{T^\gamma}{(P+P_\infty)^{\gamma-1}}} + q^\prime\end{aligned}\end{align} \]

Bifluid5eq model

The conservative variables contain a density, a 2D or 3D momentum, a total energy, a mass fraction, and a volume fraction with

\[ \begin{align}\begin{aligned}\boldsymbol{U}&=(m_g, m_l, m V_i, m E, \alpha_g)\\\boldsymbol{P}(\boldsymbol{U})&=(m_g, m_l, V_i, P, T, \alpha_g)\\\boldsymbol{\vec{F}}(\boldsymbol{U})&=(m_g V_i, m_l V_i, m V_i V_j+P \delta_{i,j}, (m E+P)V_i, 0)\\-\boldsymbol{D}(\boldsymbol{U})\vec{\nabla}\boldsymbol{U}&=(0, 0, -\tau_{i,j}, -\kappa \partial_i T, 0),\\\boldsymbol{S}(\boldsymbol{U}) &= (0, 0, m g_i, m g_k V_k, -V_i \partial_i \alpha_g )\end{aligned}\end{align} \]

with

\[ \begin{align}\begin{aligned}E&=e+V_k V_k/2,\\\tau_{i,j}&=\zeta\partial_k V_k \delta_{i,j} + \mu(\partial_i V_j + \partial_j V_i -2/3 \partial_k V_k \delta_{i,j}),\end{aligned}\end{align} \]

and an equation of state for the mixture giving the pressure/temperature \((P,T)=f(m,e)\), the thermal conductivity \(\kappa\), and the shear and bulk viscosity \(\mu\) and \(\zeta\).

The density and internal specific energy of each phase can be computed using the mass and volume fractions:

\[ \begin{align}\begin{aligned}\rho_g &= m_g/\alpha_g,\quad e_g = m_g e/(m_g+m_l)\\\rho_l &= m_l/(1-\alpha_g), \quad e_l = m_l e/(m_g+m_l)\end{aligned}\end{align} \]

The EOS for the pressure/temperature of the mixture assumes an isobaric closure \(P_g=P_l\) and the EOS of each phase \((P_i,T_i)=f(m_i,e_i)\) is given by a Noble-Abel Stiffened-Gas EOS (see the monofluid model). The thermal conductivity, shear and bulk viscosity of the mixture are assumed to be an average of each phase weighted by the volume fraction.

Bifluid7eq model

The conservative variables contain two mass fractions, two 2D or 3D momenta, two total energies, and a volume fraction with

\[ \begin{align}\begin{aligned}\boldsymbol{U}&=(m_g, m_g V_{g,i}, m_g E_g, m_l, m_l V_{l,i}, m_l E_l, \alpha_g)\\\boldsymbol{P}(\boldsymbol{U})&=(m_g, V_{g,i}, P_g, T_g, m_l, V_{l,i}, P_l, T_l, \alpha_g)\\\boldsymbol{\vec{F}}(\boldsymbol{U})&=(m_g V_{g,i},m_g V_{g,i} V_{g,j}+\alpha_g P_g \delta_{i,j}, (m_g E_g+\alpha_g P_g)V_{g,i},\\&m_l V_{l,i},m_l V_{l,i} V_{l,j}+(1-\alpha_g) P_l \delta_{i,j}, (m_l E_l+(1-\alpha_g)P_l)V_{l,i}, 0)\\-\boldsymbol{D}(\boldsymbol{U})\vec{\nabla}\boldsymbol{U}&=(0, -\tau_{g,i,j}, -\kappa_g \partial_i T_g, 0, -\tau_{l,i,j}, -\kappa_l \partial_i T_l,0),\\\boldsymbol{S}(\boldsymbol{U}) &= (0, m_g g_i + P_I\partial_i\alpha_g, m_g g_k V_{g,k} + P_I U_I\partial_i\alpha_g,\\&0, m_l g_i - P_I\partial_i\alpha_g, m_l g_k V_{l,k}-P_I U_I\partial_i\alpha_g, -U_I \partial_i \alpha_g)\end{aligned}\end{align} \]

with

\[ \begin{align}\begin{aligned}E_g&=e_g+V_{g,k} V_{g,k}/2,\\\tau_{g,i,j}&=\zeta_g\partial_k V_{g,k} \delta_{i,j} + \mu_g(\partial_i V_{g,j} + \partial_j V_{g,i} -2/3 \partial_k V_{g,k} \delta_{i,j}),\\ E_l&=e_l+V_{l,k} V_{l,k}/2,\\\tau_{l,i,j}&=\zeta_l\partial_k V_{l,k} \delta_{i,j} + \mu_l(\partial_i V_{l,j} + \partial_j V_{l,i} -2/3 \partial_k V_{l,k} \delta_{i,j}),\end{aligned}\end{align} \]

The density of each phase can be computed using the mass and volume fractions:

\[ \begin{align}\begin{aligned}\rho_g &= m_g/\alpha_g,\\\rho_l &= m_l/(1-\alpha_g),\end{aligned}\end{align} \]

and an equation of state for each phase giving the pressure/temperature \((P_i,T_i)=f(\rho_i,e_i)\), (given by a Noble-Abel Stiffened-Gas EOS, see the monofluid model), the thermal conductivity \(\kappa_i\), and the shear and bulk viscosity \(\mu_i\) and \(\zeta_i\). The available closure relations for the interface velocity and pressure are (options 1 to 5):

\[ \begin{align}\begin{aligned}U_I &= V_g, \quad P_I = P_l,\\U_I &= V_l, \quad P_I = P_g,\\U_I &= (m_g V_g+m_l V_l)/(m_g+m_l), \quad P_I = \alpha_g P_g + (1-\alpha_g)P_l,\\U_I &= \alpha_g V_g+ (1-\alpha_g) V_l, \quad P_I = \alpha_g P_g + (1-\alpha_g)P_l,\\U_I &= (V_g+ V_l)/2, \quad P_I = (P_g + P_l)/2.\end{aligned}\end{align} \]

Instantaneous relaxation source terms towards \(P_g=P_l\) and \(V_g=V_l\) are available (see Chen 2019).