Write new closure relationships in Models_NS_AC_Comp.h

The closure relationships are all functions such as source terms of PDEs, interpolation functions, force terms, fluxes, moments of distribution functions, etc. which are needed for the numerical scheme. The example is given on Conservative Allen-Cahn (CAC) equation of kernel NSAC_Comp.

Objective

Write all macroscopic functions needed for setup_collider

Write in file Models_NS_AC_Comp.h all macroscopic quantities or functions which are needed for LBM algorithm.

Examples

Definition of relaxation rate \(\tau_{\phi}\)

The collision rate \(\tau_{\phi}\) is defined with mobility coefficient \(M_{ \phi}\)

(248)\[\underbrace{\tau_{g}}_{\text{tau}}=\frac{1}{2}+3M_{\phi}\frac{\delta t}{\delta x^{2}}\]

Its value tau is defined inside the function tau_PHI by

KOKKOS_INLINE_FUNCTION
real_t tau_PHI(EquationTag1 tag, const LBMState &lbmState) const {
   real_t tau = 0.5 + (e2 * Mphi * dt / SQR(dx));
   return (tau);
}

Moment of order zero

The moment of order zero returned by the function M0_PHI is simply the function \(\phi\). Its value is stored in container lbmState[IPHI]

KOKKOS_INLINE_FUNCTION
real_t M0_PHI(EquationTag1 tag, const LBMState &lbmState) const {
   return lbmState[IPHI];
}

Moment of second order for Cahn-Hilliard equation

The chemical potential \(\mu_{\phi}\) is defined by \(\mu_{\phi}=\underbrace{2H\phi(1-\phi)(1-2\phi)}_{f_{dw}^{\prime}(\phi)}-\zeta\boldsymbol{\nabla}^{2}\phi\). With \(H=12\sigma/W\) and \(\zeta=3/2W\sigma\):

(249)\[\mu_{\phi}=\frac{3}{2}\frac{\sigma}{W}\left[g^{\prime}(\phi)-W^2\boldsymbol{\nabla}^{2}\phi\right]\]

It is computed by the function M2_MU_PHI:

KOKKOS_INLINE_FUNCTION
real_t M2_MU_PHI(EquationTag1 tag, const LBMState &lbmState) const {
        const real_t mu_phi = sigma * 1.5 / W * (g_prime(lbmState[IPHI]) - SQR(W) * lbmState[ILAPLAPHI]);
   return mu_phi;
}

Counter term

(250)\[\boldsymbol{j}_{CT}=M_{\phi}\frac{4}{W}\phi(1-\phi)\frac{\boldsymbol{\nabla}\phi}{\left|\boldsymbol{\nabla}\phi\right|}\]

It is defined in the function RVect<dim> S_ct:

template <int dim>
KOKKOS_INLINE_FUNCTION RVect<dim> S_ct(const LBMState &lbmState) const {
   const real_t norm = sqrt(SQR(lbmState[IDPHIDX]) + SQR(lbmState[IDPHIDY]) + SQR(lbmState[IDPHIDZ])) + NORMAL_EPSILON;
   real_t val_ct = Mphi * counter_term * 4.0 * lbmState[IPHI] * (1.0 - lbmState[IPHI]) / (W * norm);
   real_t force_ct = (1-cahn_hilliard)*val_ct ;

Where counter_term and cahn_hilliard are either 0 or 1

   RVect<dim> term;
   term[IX] = force_ct * lbmState[IDPHIDX];
   term[IY] = force_ct * lbmState[IDPHIDY];
   if (dim == 3) {
      term[IZ] = force_ct * lbmState[IDPHIDZ];
   }
   return term;
}

Source term

(251)\[\mathscr{S}_{\psi}(\psi,\,\overline{\mu})=-\frac{\lambda M_{\psi}}{W^{2}}6\psi(1-\psi)(c_{s}^{co}-c_{l}^{co})(\overline{\mu}(c,\phi)-\overline{\mu}^{eq})\]

The source term is \(\mathscr{S}_{\phi}\) is defined in function S_st

KOKKOS_INLINE_FUNCTION
real_t S_st(EquationTag1 tag, const LBMState &lbmState) const {
   real_t S = p_prime(lbmState[IPHI]) * (c0_co - c1_co) * (mu(lbmState[IC], lbmState[IPHI]) - mu_eq);
   return -lambda * Mphi / SQR(W) * S;
}

where the function \(\mu(c,\phi)\) (mu) is defined by

// =======================================================
//                                     Model for composition
// =======================================================

KOKKOS_INLINE_FUNCTION
real_t mu(real_t c, real_t phi) const { return mu_eq + c - c_co(phi); }

Section author: Alain Cartalade