Skip to content
Snippets Groups Projects
Commit f6647fd7 authored by Soeren Peters's avatar Soeren Peters
Browse files
parent a1b8c8da
No related branches found
No related tags found
1 merge request!282Feature/rename collision k17
......@@ -84,64 +84,64 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
{
auto& distribution = parameter.distribution;
real f_000 = distribution[d000];
real f_P00 = distribution[dP00];
real f_M00 = distribution[dM00];
real f_0P0 = distribution[d0P0];
real f_0M0 = distribution[d0M0];
real f_00P = distribution[d00P];
real f_00M = distribution[d00M];
real f_PP0 = distribution[dPP0];
real f_MM0 = distribution[dMM0];
real f_PM0 = distribution[dPM0];
real f_MP0 = distribution[dMP0];
real f_P0P = distribution[dP0P];
real f_M0M = distribution[dM0M];
real f_P0M = distribution[dP0M];
real f_M0P = distribution[dM0P];
real f_0PP = distribution[d0PP];
real f_0MM = distribution[d0MM];
real f_0PM = distribution[d0PM];
real f_0MP = distribution[d0MP];
real f_PPP = distribution[dPPP];
real f_MPP = distribution[dMPP];
real f_PMP = distribution[dPMP];
real f_MMP = distribution[dMMP];
real f_PPM = distribution[dPPM];
real f_MPM = distribution[dMPM];
real f_PMM = distribution[dPMM];
real f_MMM = distribution[dMMM];
real f000 = distribution[d000];
real fP00 = distribution[dP00];
real fM00 = distribution[dM00];
real f0P0 = distribution[d0P0];
real f0M0 = distribution[d0M0];
real f00P = distribution[d00P];
real f00M = distribution[d00M];
real fPP0 = distribution[dPP0];
real fMM0 = distribution[dMM0];
real fPM0 = distribution[dPM0];
real fMP0 = distribution[dMP0];
real fP0P = distribution[dP0P];
real fM0M = distribution[dM0M];
real fP0M = distribution[dP0M];
real fM0P = distribution[dM0P];
real f0PP = distribution[d0PP];
real f0MM = distribution[d0MM];
real f0PM = distribution[d0PM];
real f0MP = distribution[d0MP];
real fPPP = distribution[dPPP];
real fMPP = distribution[dMPP];
real fPMP = distribution[dPMP];
real fMMP = distribution[dMMP];
real fPPM = distribution[dPPM];
real fMPM = distribution[dMPM];
real fPMM = distribution[dPMM];
real fMMM = distribution[dMMM];
////////////////////////////////////////////////////////////////////////////////////
//! - Define aliases to use the same variable for the moments (m's):
//!
real& m_111 = f_000;
real& m_211 = f_P00;
real& m_011 = f_M00;
real& m_121 = f_0P0;
real& m_101 = f_0M0;
real& m_112 = f_00P;
real& m_110 = f_00M;
real& m_221 = f_PP0;
real& m_001 = f_MM0;
real& m_201 = f_PM0;
real& m_021 = f_MP0;
real& m_212 = f_P0P;
real& m_010 = f_M0M;
real& m_210 = f_P0M;
real& m_012 = f_M0P;
real& m_122 = f_0PP;
real& m_100 = f_0MM;
real& m_120 = f_0PM;
real& m_102 = f_0MP;
real& m_222 = f_PPP;
real& m_022 = f_MPP;
real& m_202 = f_PMP;
real& m_002 = f_MMP;
real& m_220 = f_PPM;
real& m_020 = f_MPM;
real& m_200 = f_PMM;
real& m_000 = f_MMM;
real& m111 = f000;
real& m211 = fP00;
real& m011 = fM00;
real& m121 = f0P0;
real& m101 = f0M0;
real& m112 = f00P;
real& m110 = f00M;
real& m221 = fPP0;
real& m001 = fMM0;
real& m201 = fPM0;
real& m021 = fMP0;
real& m212 = fP0P;
real& m010 = fM0M;
real& m210 = fP0M;
real& m012 = fM0P;
real& m122 = f0PP;
real& m100 = f0MM;
real& m120 = f0PM;
real& m102 = f0MP;
real& m222 = fPPP;
real& m022 = fMPP;
real& m202 = fPMP;
real& m002 = fMMP;
real& m220 = fPPM;
real& m020 = fMPM;
real& m200 = fPMM;
real& m000 = fMMM;
//////////////////////////////////////////////////////(unsigned long)//////////////////////////////
//! - Calculate density and velocity using pyramid summation for low round-off errors as in Eq. (J1)-(J3) \ref
......@@ -182,57 +182,57 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
//!
////////////////////////////////////////////////////////////////////////////////////
// Z - Dir
forwardInverseChimeraWithK(f_MMM, f_MM0, f_MMP, vvz, vz2, c36o1, c1o36);
forwardInverseChimeraWithK(f_M0M, f_M00, f_M0P, vvz, vz2, c9o1, c1o9);
forwardInverseChimeraWithK(f_MPM, f_MP0, f_MPP, vvz, vz2, c36o1, c1o36);
forwardInverseChimeraWithK(f_0MM, f_0M0, f_0MP, vvz, vz2, c9o1, c1o9);
forwardInverseChimeraWithK(f_00M, f_000, f_00P, vvz, vz2, c9o4, c4o9);
forwardInverseChimeraWithK(f_0PM, f_0P0, f_0PP, vvz, vz2, c9o1, c1o9);
forwardInverseChimeraWithK(f_PMM, f_PM0, f_PMP, vvz, vz2, c36o1, c1o36);
forwardInverseChimeraWithK(f_P0M, f_P00, f_P0P, vvz, vz2, c9o1, c1o9);
forwardInverseChimeraWithK(f_PPM, f_PP0, f_PPP, vvz, vz2, c36o1, c1o36);
forwardInverseChimeraWithK(fMMM, fMM0, fMMP, vvz, vz2, c36o1, c1o36);
forwardInverseChimeraWithK(fM0M, fM00, fM0P, vvz, vz2, c9o1, c1o9);
forwardInverseChimeraWithK(fMPM, fMP0, fMPP, vvz, vz2, c36o1, c1o36);
forwardInverseChimeraWithK(f0MM, f0M0, f0MP, vvz, vz2, c9o1, c1o9);
forwardInverseChimeraWithK(f00M, f000, f00P, vvz, vz2, c9o4, c4o9);
forwardInverseChimeraWithK(f0PM, f0P0, f0PP, vvz, vz2, c9o1, c1o9);
forwardInverseChimeraWithK(fPMM, fPM0, fPMP, vvz, vz2, c36o1, c1o36);
forwardInverseChimeraWithK(fP0M, fP00, fP0P, vvz, vz2, c9o1, c1o9);
forwardInverseChimeraWithK(fPPM, fPP0, fPPP, vvz, vz2, c36o1, c1o36);
////////////////////////////////////////////////////////////////////////////////////
// Y - Dir
forwardInverseChimeraWithK(f_MMM, f_M0M, f_MPM, vvy, vy2, c6o1, c1o6);
forwardChimera( f_MM0, f_M00, f_MP0, vvy, vy2);
forwardInverseChimeraWithK(f_MMP, f_M0P, f_MPP, vvy, vy2, c18o1, c1o18);
forwardInverseChimeraWithK(f_0MM, f_00M, f_0PM, vvy, vy2, c3o2, c2o3);
forwardChimera( f_0M0, f_000, f_0P0, vvy, vy2);
forwardInverseChimeraWithK(f_0MP, f_00P, f_0PP, vvy, vy2, c9o2, c2o9);
forwardInverseChimeraWithK(f_PMM, f_P0M, f_PPM, vvy, vy2, c6o1, c1o6);
forwardChimera( f_PM0, f_P00, f_PP0, vvy, vy2);
forwardInverseChimeraWithK(f_PMP, f_P0P, f_PPP, vvy, vy2, c18o1, c1o18);
forwardInverseChimeraWithK(fMMM, fM0M, fMPM, vvy, vy2, c6o1, c1o6);
forwardChimera( fMM0, fM00, fMP0, vvy, vy2);
forwardInverseChimeraWithK(fMMP, fM0P, fMPP, vvy, vy2, c18o1, c1o18);
forwardInverseChimeraWithK(f0MM, f00M, f0PM, vvy, vy2, c3o2, c2o3);
forwardChimera( f0M0, f000, f0P0, vvy, vy2);
forwardInverseChimeraWithK(f0MP, f00P, f0PP, vvy, vy2, c9o2, c2o9);
forwardInverseChimeraWithK(fPMM, fP0M, fPPM, vvy, vy2, c6o1, c1o6);
forwardChimera( fPM0, fP00, fPP0, vvy, vy2);
forwardInverseChimeraWithK(fPMP, fP0P, fPPP, vvy, vy2, c18o1, c1o18);
////////////////////////////////////////////////////////////////////////////////////
// X - Dir
forwardInverseChimeraWithK(f_MMM, f_0MM, f_PMM, vvx, vx2, c1o1, c1o1);
forwardChimera( f_M0M, f_00M, f_P0M, vvx, vx2);
forwardInverseChimeraWithK(f_MPM, f_0PM, f_PPM, vvx, vx2, c3o1, c1o3);
forwardChimera( f_MM0, f_0M0, f_PM0, vvx, vx2);
forwardChimera( f_M00, f_000, f_P00, vvx, vx2);
forwardChimera( f_MP0, f_0P0, f_PP0, vvx, vx2);
forwardInverseChimeraWithK(f_MMP, f_0MP, f_PMP, vvx, vx2, c3o1, c1o3);
forwardChimera( f_M0P, f_00P, f_P0P, vvx, vx2);
forwardInverseChimeraWithK(f_MPP, f_0PP, f_PPP, vvx, vx2, c3o1, c1o9);
forwardInverseChimeraWithK(fMMM, f0MM, fPMM, vvx, vx2, c1o1, c1o1);
forwardChimera( fM0M, f00M, fP0M, vvx, vx2);
forwardInverseChimeraWithK(fMPM, f0PM, fPPM, vvx, vx2, c3o1, c1o3);
forwardChimera( fMM0, f0M0, fPM0, vvx, vx2);
forwardChimera( fM00, f000, fP00, vvx, vx2);
forwardChimera( fMP0, f0P0, fPP0, vvx, vx2);
forwardInverseChimeraWithK(fMMP, f0MP, fPMP, vvx, vx2, c3o1, c1o3);
forwardChimera( fM0P, f00P, fP0P, vvx, vx2);
forwardInverseChimeraWithK(fMPP, f0PP, fPPP, vvx, vx2, c3o1, c1o9);
////////////////////////////////////////////////////////////////////////////////////
//! - Setting relaxation rates for non-hydrodynamic cumulants (default values). Variable names and equations
//! according to <a href="https://doi.org/10.1016/j.jcp.2017.05.040"><b>[ M. Geier et al. (2017),
//! DOI:10.1016/j.jcp.2017.05.040 ]</b></a>
//! => [NAME IN PAPER]=[NAME IN CODE]=[DEFAULT VALUE].
//! - Trace of second order cumulants \f$ C_{200}+C_{020}+C_{002} \f$ used to adjust bulk
//! viscosity:\f$\omega_2=OxxPyyPzz=1.0 \f$.
//! - Third order cumulants \f$ C_{120}+C_{102}, C_{210}+C_{012}, C_{201}+C_{021} \f$: \f$ \omega_3=OxyyPxzz
//! \f$ set according to Eq. (111) with simplifications assuming \f$ \omega_2=1.0\f$.
//! - Third order cumulants \f$ C_{120}-C_{102}, C_{210}-C_{012}, C_{201}-C_{021} \f$: \f$ \omega_4 = OxyyMxzz
//! \f$ set according to Eq. (112) with simplifications assuming \f$ \omega_2 = 1.0\f$.
//! - Third order cumulants \f$ C_{111} \f$: \f$ \omega_5 = Oxyz \f$ set according to Eq. (113) with
//! simplifications assuming \f$ \omega_2 = 1.0\f$ (modify for different bulk viscosity).
//! - Fourth order cumulants \f$ C_{220}, C_{202}, C_{022}, C_{211}, C_{121}, C_{112} \f$: for simplification
//! all set to the same default value \f$ \omega_6=\omega_7=\omega_8=O4=1.0 \f$.
//! - Fifth order cumulants \f$ C_{221}, C_{212}, C_{122}\f$: \f$\omega_9=O5=1.0\f$.
//! - Sixth order cumulant \f$ C_{222}\f$: \f$\omega_{10}=O6=1.0\f$.
//! - Trace of second order cumulants \f$ C{200}+C{020}+C{002} \f$ used to adjust bulk
//! viscosity:\f$\omega2=OxxPyyPzz=1.0 \f$.
//! - Third order cumulants \f$ C{120}+C{102}, C{210}+C{012}, C{201}+C{021} \f$: \f$ \omega3=OxyyPxzz
//! \f$ set according to Eq. (111) with simplifications assuming \f$ \omega2=1.0\f$.
//! - Third order cumulants \f$ C{120}-C{102}, C{210}-C{012}, C{201}-C{021} \f$: \f$ \omega4 = OxyyMxzz
//! \f$ set according to Eq. (112) with simplifications assuming \f$ \omega2 = 1.0\f$.
//! - Third order cumulants \f$ C{111} \f$: \f$ \omega5 = Oxyz \f$ set according to Eq. (113) with
//! simplifications assuming \f$ \omega2 = 1.0\f$ (modify for different bulk viscosity).
//! - Fourth order cumulants \f$ C{220}, C{202}, C{022}, C{211}, C{121}, C{112} \f$: for simplification
//! all set to the same default value \f$ \omega6=\omega7=\omega8=O4=1.0 \f$.
//! - Fifth order cumulants \f$ C{221}, C{212}, C{122}\f$: \f$\omega9=O5=1.0\f$.
//! - Sixth order cumulant \f$ C{222}\f$: \f$\omega{10}=O6=1.0\f$.
//!
////////////////////////////////////////////////////////////////////////////////////
//! - Calculate modified omega with turbulent viscosity
......@@ -274,39 +274,39 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
//!
////////////////////////////////////////////////////////////
// 4.
real c_211 = m_211 - ((m_200 + c1o3) * m_011 + c2o1 * m_110 * m_101) * oneOverRho;
real c_121 = m_121 - ((m_020 + c1o3) * m_101 + c2o1 * m_110 * m_011) * oneOverRho;
real c_112 = m_112 - ((m_002 + c1o3) * m_110 + c2o1 * m_101 * m_011) * oneOverRho;
real cm211 = m211 - ((m200 + c1o3) * m011 + c2o1 * m110 * m101) * oneOverRho;
real cm121 = m121 - ((m020 + c1o3) * m101 + c2o1 * m110 * m011) * oneOverRho;
real cm112 = m112 - ((m002 + c1o3) * m110 + c2o1 * m101 * m011) * oneOverRho;
real c_220 = m_220 - (((m_200 * m_020 + c2o1 * m_110 * m_110) + c1o3 * (m_200 + m_020)) * oneOverRho - c1o9 * (drho * oneOverRho));
real c_202 = m_202 - (((m_200 * m_002 + c2o1 * m_101 * m_101) + c1o3 * (m_200 + m_002)) * oneOverRho - c1o9 * (drho * oneOverRho));
real c_022 = m_022 - (((m_002 * m_020 + c2o1 * m_011 * m_011) + c1o3 * (m_002 + m_020)) * oneOverRho - c1o9 * (drho * oneOverRho));
real cm220 = m220 - (((m200 * m020 + c2o1 * m110 * m110) + c1o3 * (m200 + m020)) * oneOverRho - c1o9 * (drho * oneOverRho));
real cm202 = m202 - (((m200 * m002 + c2o1 * m101 * m101) + c1o3 * (m200 + m002)) * oneOverRho - c1o9 * (drho * oneOverRho));
real cm022 = m022 - (((m002 * m020 + c2o1 * m011 * m011) + c1o3 * (m002 + m020)) * oneOverRho - c1o9 * (drho * oneOverRho));
////////////////////////////////////////////////////////////
// 5.
real c_122 =
m_122 - ((m_002 * m_120 + m_020 * m_102 + c4o1 * m_011 * m_111 + c2o1 * (m_101 * m_021 + m_110 * m_012)) +
c1o3 * (m_120 + m_102)) *
real cm122 =
m122 - ((m002 * m120 + m020 * m102 + c4o1 * m011 * m111 + c2o1 * (m101 * m021 + m110 * m012)) +
c1o3 * (m120 + m102)) *
oneOverRho;
real c_212 =
m_212 - ((m_002 * m_210 + m_200 * m_012 + c4o1 * m_101 * m_111 + c2o1 * (m_011 * m_201 + m_110 * m_102)) +
c1o3 * (m_210 + m_012)) *
real cm212 =
m212 - ((m002 * m210 + m200 * m012 + c4o1 * m101 * m111 + c2o1 * (m011 * m201 + m110 * m102)) +
c1o3 * (m210 + m012)) *
oneOverRho;
real c_221 =
m_221 - ((m_200 * m_021 + m_020 * m_201 + c4o1 * m_110 * m_111 + c2o1 * (m_101 * m_120 + m_011 * m_210)) +
c1o3 * (m_021 + m_201)) *
real cm221 =
m221 - ((m200 * m021 + m020 * m201 + c4o1 * m110 * m111 + c2o1 * (m101 * m120 + m011 * m210)) +
c1o3 * (m021 + m201)) *
oneOverRho;
////////////////////////////////////////////////////////////
// 6.
real c_222 = m_222 + ((-c4o1 * m_111 * m_111 - (m_200 * m_022 + m_020 * m_202 + m_002 * m_220) -
c4o1 * (m_011 * m_211 + m_101 * m_121 + m_110 * m_112) -
c2o1 * (m_120 * m_102 + m_210 * m_012 + m_201 * m_021)) *
real cm222 = m222 + ((-c4o1 * m111 * m111 - (m200 * m022 + m020 * m202 + m002 * m220) -
c4o1 * (m011 * m211 + m101 * m121 + m110 * m112) -
c2o1 * (m120 * m102 + m210 * m012 + m201 * m021)) *
oneOverRho +
(c4o1 * (m_101 * m_101 * m_020 + m_011 * m_011 * m_200 + m_110 * m_110 * m_002) +
c2o1 * (m_200 * m_020 * m_002) + c16o1 * m_110 * m_101 * m_011) *
(c4o1 * (m101 * m101 * m020 + m011 * m011 * m200 + m110 * m110 * m002) +
c2o1 * (m200 * m020 * m002) + c16o1 * m110 * m101 * m011) *
oneOverRho * oneOverRho -
c1o3 * (m_022 + m_202 + m_220) * oneOverRho - c1o9 * (m_200 + m_020 + m_002) * oneOverRho +
(c2o1 * (m_101 * m_101 + m_011 * m_011 + m_110 * m_110) +
(m_002 * m_020 + m_002 * m_200 + m_020 * m_200) + c1o3 * (m_002 + m_020 + m_200)) *
c1o3 * (m022 + m202 + m220) * oneOverRho - c1o9 * (m200 + m020 + m002) * oneOverRho +
(c2o1 * (m101 * m101 + m011 * m011 + m110 * m110) +
(m002 * m020 + m002 * m200 + m020 * m200) + c1o3 * (m002 + m020 + m200)) *
oneOverRho * oneOverRho * c2o3 +
c1o27 * ((drho * drho - drho) * oneOverRho * oneOverRho));
......@@ -315,19 +315,19 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
//!
////////////////////////////////////////////////////////////
// 2.
real mxxPyyPzz = m_200 + m_020 + m_002;
real mxxMyy = m_200 - m_020;
real mxxMzz = m_200 - m_002;
real mxxPyyPzz = m200 + m020 + m002;
real mxxMyy = m200 - m020;
real mxxMzz = m200 - m002;
////////////////////////////////////////////////////////////
// 3.
real mxxyPyzz = m_210 + m_012;
real mxxyMyzz = m_210 - m_012;
real mxxyPyzz = m210 + m012;
real mxxyMyzz = m210 - m012;
real mxxzPyyz = m_201 + m_021;
real mxxzMyyz = m_201 - m_021;
real mxxzPyyz = m201 + m021;
real mxxzMyyz = m201 - m021;
real mxyyPxzz = m_120 + m_102;
real mxyyMxzz = m_120 - m_102;
real mxyyPxzz = m120 + m102;
real mxyyMxzz = m120 - m102;
////////////////////////////////////////////////////////////////////////////////////
// incl. correction
......@@ -339,34 +339,35 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
//! DOI:10.1016/j.camwa.2015.05.001 ]</b></a> Note that the division by rho is omitted here as we need rho times
//! the gradients later.
//!
real Dxy = -c3o1 * omega * m_110;
real Dxz = -c3o1 * omega * m_101;
real Dyz = -c3o1 * omega * m_011;
real dxux = c1o2 * (-omega) * (mxxMyy + mxxMzz) + c1o2 * OxxPyyPzz * (m_000 - mxxPyyPzz);
real Dxy = -c3o1 * omega * m110;
real Dxz = -c3o1 * omega * m101;
real Dyz = -c3o1 * omega * m011;
real dxux = c1o2 * (-omega) * (mxxMyy + mxxMzz) + c1o2 * OxxPyyPzz * (m000 - mxxPyyPzz);
real dyuy = dxux + omega * c3o2 * mxxMyy;
real dzuz = dxux + omega * c3o2 * mxxMzz;
////////////////////////////////////////////////////////////////////////////////////
switch (turbulenceModel)
{
case TurbulenceModel::None:
case TurbulenceModel::AMD: //AMD is computed in separate kernel
break;
case TurbulenceModel::Smagorinsky:
turbulentViscosity.value = calcTurbulentViscositySmagorinsky(turbulentViscosity.SGSconstant, dxux, dyuy, dzuz, Dxy, Dxz , Dyz);
break;
case TurbulenceModel::QR:
turbulentViscosity.value = calcTurbulentViscosityQR(turbulentViscosity.SGSconstant, dxux, dyuy, dzuz, Dxy, Dxz , Dyz);
break;
default:
break;
switch (turbulenceModel) {
case TurbulenceModel::None:
case TurbulenceModel::AMD: // AMD is computed in separate kernel
break;
case TurbulenceModel::Smagorinsky:
turbulentViscosity.value =
calcTurbulentViscositySmagorinsky(turbulentViscosity.SGSconstant, dxux, dyuy, dzuz, Dxy, Dxz, Dyz);
break;
case TurbulenceModel::QR:
turbulentViscosity.value =
calcTurbulentViscosityQR(turbulentViscosity.SGSconstant, dxux, dyuy, dzuz, Dxy, Dxz, Dyz);
break;
default:
break;
}
////////////////////////////////////////////////////////////
//! - Relaxation of second order cumulants with correction terms according to Eq. (33)-(35) in
//! <a href="https://doi.org/10.1016/j.jcp.2017.05.040"><b>[ M. Geier et al. (2017),
//! DOI:10.1016/j.jcp.2017.05.040 ]</b></a>
//!
mxxPyyPzz += OxxPyyPzz * (m_000 - mxxPyyPzz) - c3o1 * (c1o1 - c1o2 * OxxPyyPzz) * (vx2 * dxux + vy2 * dyuy + vz2 * dzuz);
mxxPyyPzz += OxxPyyPzz * (m000 - mxxPyyPzz) - c3o1 * (c1o1 - c1o2 * OxxPyyPzz) * (vx2 * dxux + vy2 * dyuy + vz2 * dzuz);
mxxMyy += omega * (-mxxMyy) - c3o1 * (c1o1 + c1o2 * (-omega)) * (vx2 * dxux - vy2 * dyuy);
mxxMzz += omega * (-mxxMzz) - c3o1 * (c1o1 + c1o2 * (-omega)) * (vx2 * dxux - vz2 * dzuz);
......@@ -376,9 +377,9 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
// mxxMyy += -(-omega) * (-mxxMyy);
// mxxMzz += -(-omega) * (-mxxMzz);
//////////////////////////////////////////////////////////////////////////
m_011 += omega * (-m_011);
m_101 += omega * (-m_101);
m_110 += omega * (-m_110);
m011 += omega * (-m011);
m101 += omega * (-m101);
m110 += omega * (-m110);
////////////////////////////////////////////////////////////////////////////////////
// relax
......@@ -388,8 +389,8 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
//! <a href="https://doi.org/10.1016/j.jcp.2017.05.040"><b>[ M. Geier et al. (2017),
//! DOI:10.1016/j.jcp.2017.05.040 ]</b></a>
//!
real wadjust = Oxyz + (c1o1 - Oxyz) * KERNEL_ABS(m_111) / (KERNEL_ABS(m_111) + quadricLimitD);
m_111 += wadjust * (-m_111);
real wadjust = Oxyz + (c1o1 - Oxyz) * KERNEL_ABS(m111) / (KERNEL_ABS(m111) + quadricLimitD);
m111 += wadjust * (-m111);
wadjust = OxyyPxzz + (c1o1 - OxyyPxzz) * KERNEL_ABS(mxxyPyzz) / (KERNEL_ABS(mxxyPyzz) + quadricLimitP);
mxxyPyzz += wadjust * (-mxxyPyzz);
wadjust = OxyyMxzz + (c1o1 - OxyyMxzz) * KERNEL_ABS(mxxyMyzz) / (KERNEL_ABS(mxxyMyzz) + quadricLimitM);
......@@ -415,16 +416,16 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
////////////////////////////////////////////////////////////////////////////////////
//! - Compute inverse linear combinations of second and third order cumulants
//!
m_200 = c1o3 * (mxxMyy + mxxMzz + mxxPyyPzz);
m_020 = c1o3 * (-c2o1 * mxxMyy + mxxMzz + mxxPyyPzz);
m_002 = c1o3 * (mxxMyy - c2o1 * mxxMzz + mxxPyyPzz);
m_210 = ( mxxyMyzz + mxxyPyzz) * c1o2;
m_012 = (-mxxyMyzz + mxxyPyzz) * c1o2;
m_201 = ( mxxzMyyz + mxxzPyyz) * c1o2;
m_021 = (-mxxzMyyz + mxxzPyyz) * c1o2;
m_120 = ( mxyyMxzz + mxyyPxzz) * c1o2;
m_102 = (-mxyyMxzz + mxyyPxzz) * c1o2;
m200 = c1o3 * (mxxMyy + mxxMzz + mxxPyyPzz);
m020 = c1o3 * (-c2o1 * mxxMyy + mxxMzz + mxxPyyPzz);
m002 = c1o3 * (mxxMyy - c2o1 * mxxMzz + mxxPyyPzz);
m210 = ( mxxyMyzz + mxxyPyzz) * c1o2;
m012 = (-mxxyMyzz + mxxyPyzz) * c1o2;
m201 = ( mxxzMyyz + mxxzPyyz) * c1o2;
m021 = (-mxxzMyyz + mxxzPyyz) * c1o2;
m120 = ( mxyyMxzz + mxyyPxzz) * c1o2;
m102 = (-mxyyMxzz + mxyyPxzz) * c1o2;
//////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////
......@@ -434,23 +435,23 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
//! to Eq. (43)-(48) <a href="https://doi.org/10.1016/j.jcp.2017.05.040"><b>[ M. Geier et al. (2017),
//! DOI:10.1016/j.jcp.2017.05.040 ]</b></a>
//!
c_022 = -O4 * (c1o1 / omega - c1o2) * (dyuy + dzuz) * c2o3 * factorA + (c1o1 - O4) * (c_022);
c_202 = -O4 * (c1o1 / omega - c1o2) * (dxux + dzuz) * c2o3 * factorA + (c1o1 - O4) * (c_202);
c_220 = -O4 * (c1o1 / omega - c1o2) * (dyuy + dxux) * c2o3 * factorA + (c1o1 - O4) * (c_220);
c_112 = -O4 * (c1o1 / omega - c1o2) * Dxy * c1o3 * factorB + (c1o1 - O4) * (c_112);
c_121 = -O4 * (c1o1 / omega - c1o2) * Dxz * c1o3 * factorB + (c1o1 - O4) * (c_121);
c_211 = -O4 * (c1o1 / omega - c1o2) * Dyz * c1o3 * factorB + (c1o1 - O4) * (c_211);
cm022 = -O4 * (c1o1 / omega - c1o2) * (dyuy + dzuz) * c2o3 * factorA + (c1o1 - O4) * (cm022);
cm202 = -O4 * (c1o1 / omega - c1o2) * (dxux + dzuz) * c2o3 * factorA + (c1o1 - O4) * (cm202);
cm220 = -O4 * (c1o1 / omega - c1o2) * (dyuy + dxux) * c2o3 * factorA + (c1o1 - O4) * (cm220);
cm112 = -O4 * (c1o1 / omega - c1o2) * Dxy * c1o3 * factorB + (c1o1 - O4) * (cm112);
cm121 = -O4 * (c1o1 / omega - c1o2) * Dxz * c1o3 * factorB + (c1o1 - O4) * (cm121);
cm211 = -O4 * (c1o1 / omega - c1o2) * Dyz * c1o3 * factorB + (c1o1 - O4) * (cm211);
//////////////////////////////////////////////////////////////////////////
// 5.
c_122 += O5 * (-c_122);
c_212 += O5 * (-c_212);
c_221 += O5 * (-c_221);
cm122 += O5 * (-cm122);
cm212 += O5 * (-cm212);
cm221 += O5 * (-cm221);
//////////////////////////////////////////////////////////////////////////
// 6.
c_222 += O6 * (-c_222);
cm222 += O6 * (-cm222);
////////////////////////////////////////////////////////////////////////////////////
//! - Compute central moments from post collision cumulants according to Eq. (53)-(56) in
......@@ -460,41 +461,41 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
//////////////////////////////////////////////////////////////////////////
// 4.
m_211 = c_211 + c1o3 * ((c3o1 * m_200 + c1o1) * m_011 + c6o1 * m_110 * m_101) * oneOverRho;
m_121 = c_121 + c1o3 * ((c3o1 * m_020 + c1o1) * m_101 + c6o1 * m_110 * m_011) * oneOverRho;
m_112 = c_112 + c1o3 * ((c3o1 * m_002 + c1o1) * m_110 + c6o1 * m_101 * m_011) * oneOverRho;
m211 = cm211 + c1o3 * ((c3o1 * m200 + c1o1) * m011 + c6o1 * m110 * m101) * oneOverRho;
m121 = cm121 + c1o3 * ((c3o1 * m020 + c1o1) * m101 + c6o1 * m110 * m011) * oneOverRho;
m112 = cm112 + c1o3 * ((c3o1 * m002 + c1o1) * m110 + c6o1 * m101 * m011) * oneOverRho;
m_220 =
c_220 + (((m_200 * m_020 + c2o1 * m_110 * m_110) * c9o1 + c3o1 * (m_200 + m_020)) * oneOverRho - (drho * oneOverRho)) * c1o9;
m_202 =
c_202 + (((m_200 * m_002 + c2o1 * m_101 * m_101) * c9o1 + c3o1 * (m_200 + m_002)) * oneOverRho - (drho * oneOverRho)) * c1o9;
m_022 =
c_022 + (((m_002 * m_020 + c2o1 * m_011 * m_011) * c9o1 + c3o1 * (m_002 + m_020)) * oneOverRho - (drho * oneOverRho)) * c1o9;
m220 =
cm220 + (((m200 * m020 + c2o1 * m110 * m110) * c9o1 + c3o1 * (m200 + m020)) * oneOverRho - (drho * oneOverRho)) * c1o9;
m202 =
cm202 + (((m200 * m002 + c2o1 * m101 * m101) * c9o1 + c3o1 * (m200 + m002)) * oneOverRho - (drho * oneOverRho)) * c1o9;
m022 =
cm022 + (((m002 * m020 + c2o1 * m011 * m011) * c9o1 + c3o1 * (m002 + m020)) * oneOverRho - (drho * oneOverRho)) * c1o9;
//////////////////////////////////////////////////////////////////////////
// 5.
m_122 = c_122 + c1o3 *
(c3o1 * (m_002 * m_120 + m_020 * m_102 + c4o1 * m_011 * m_111 + c2o1 * (m_101 * m_021 + m_110 * m_012)) +
(m_120 + m_102)) * oneOverRho;
m_212 = c_212 + c1o3 *
(c3o1 * (m_002 * m_210 + m_200 * m_012 + c4o1 * m_101 * m_111 + c2o1 * (m_011 * m_201 + m_110 * m_102)) +
(m_210 + m_012)) * oneOverRho;
m_221 = c_221 + c1o3 *
(c3o1 * (m_200 * m_021 + m_020 * m_201 + c4o1 * m_110 * m_111 + c2o1 * (m_101 * m_120 + m_011 * m_210)) +
(m_021 + m_201)) * oneOverRho;
m122 = cm122 + c1o3 *
(c3o1 * (m002 * m120 + m020 * m102 + c4o1 * m011 * m111 + c2o1 * (m101 * m021 + m110 * m012)) +
(m120 + m102)) * oneOverRho;
m212 = cm212 + c1o3 *
(c3o1 * (m002 * m210 + m200 * m012 + c4o1 * m101 * m111 + c2o1 * (m011 * m201 + m110 * m102)) +
(m210 + m012)) * oneOverRho;
m221 = cm221 + c1o3 *
(c3o1 * (m200 * m021 + m020 * m201 + c4o1 * m110 * m111 + c2o1 * (m101 * m120 + m011 * m210)) +
(m021 + m201)) * oneOverRho;
//////////////////////////////////////////////////////////////////////////
// 6.
m_222 = c_222 - ((-c4o1 * m_111 * m_111 - (m_200 * m_022 + m_020 * m_202 + m_002 * m_220) -
c4o1 * (m_011 * m_211 + m_101 * m_121 + m_110 * m_112) -
c2o1 * (m_120 * m_102 + m_210 * m_012 + m_201 * m_021)) *
m222 = cm222 - ((-c4o1 * m111 * m111 - (m200 * m022 + m020 * m202 + m002 * m220) -
c4o1 * (m011 * m211 + m101 * m121 + m110 * m112) -
c2o1 * (m120 * m102 + m210 * m012 + m201 * m021)) *
oneOverRho +
(c4o1 * (m_101 * m_101 * m_020 + m_011 * m_011 * m_200 + m_110 * m_110 * m_002) +
c2o1 * (m_200 * m_020 * m_002) + c16o1 * m_110 * m_101 * m_011) *
(c4o1 * (m101 * m101 * m020 + m011 * m011 * m200 + m110 * m110 * m002) +
c2o1 * (m200 * m020 * m002) + c16o1 * m110 * m101 * m011) *
oneOverRho * oneOverRho -
c1o3 * (m_022 + m_202 + m_220) * oneOverRho - c1o9 * (m_200 + m_020 + m_002) * oneOverRho +
(c2o1 * (m_101 * m_101 + m_011 * m_011 + m_110 * m_110) +
(m_002 * m_020 + m_002 * m_200 + m_020 * m_200) + c1o3 * (m_002 + m_020 + m_200)) *
c1o3 * (m022 + m202 + m220) * oneOverRho - c1o9 * (m200 + m020 + m002) * oneOverRho +
(c2o1 * (m101 * m101 + m011 * m011 + m110 * m110) +
(m002 * m020 + m002 * m200 + m020 * m200) + c1o3 * (m002 + m020 + m200)) *
oneOverRho * oneOverRho * c2o3 +
c1o27 * ((drho * drho - drho) * oneOverRho * oneOverRho));
......@@ -503,9 +504,9 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
//! <a href="https://doi.org/10.1016/j.camwa.2015.05.001"><b>[ M. Geier et al. (2015),
//! DOI:10.1016/j.camwa.2015.05.001 ]</b></a>
//!
m_100 = -m_100;
m_010 = -m_010;
m_001 = -m_001;
m100 = -m100;
m010 = -m010;
m001 = -m001;
////////////////////////////////////////////////////////////////////////////////////
//! - Chimera transform from central moments to well conditioned distributions as defined in Appendix J in
......@@ -516,67 +517,67 @@ __host__ __device__ void runK17CompressibleNavierStokes(CollisionParameter& para
//!
////////////////////////////////////////////////////////////////////////////////////
// X - Dir
backwardInverseChimeraWithK(m_000, m_100, m_200, vvx, vx2, c1o1, c1o1);
backwardChimera( m_010, m_110, m_210, vvx, vx2);
backwardInverseChimeraWithK(m_020, m_120, m_220, vvx, vx2, c3o1, c1o3);
backwardChimera( m_001, m_101, m_201, vvx, vx2);
backwardChimera( m_011, m_111, m_211, vvx, vx2);
backwardChimera( m_021, m_121, m_221, vvx, vx2);
backwardInverseChimeraWithK(m_002, m_102, m_202, vvx, vx2, c3o1, c1o3);
backwardChimera( m_012, m_112, m_212, vvx, vx2);
backwardInverseChimeraWithK(m_022, m_122, m_222, vvx, vx2, c9o1, c1o9);
backwardInverseChimeraWithK(m000, m100, m200, vvx, vx2, c1o1, c1o1);
backwardChimera( m010, m110, m210, vvx, vx2);
backwardInverseChimeraWithK(m020, m120, m220, vvx, vx2, c3o1, c1o3);
backwardChimera( m001, m101, m201, vvx, vx2);
backwardChimera( m011, m111, m211, vvx, vx2);
backwardChimera( m021, m121, m221, vvx, vx2);
backwardInverseChimeraWithK(m002, m102, m202, vvx, vx2, c3o1, c1o3);
backwardChimera( m012, m112, m212, vvx, vx2);
backwardInverseChimeraWithK(m022, m122, m222, vvx, vx2, c9o1, c1o9);
////////////////////////////////////////////////////////////////////////////////////
// Y - Dir
backwardInverseChimeraWithK(m_000, m_010, m_020, vvy, vy2, c6o1, c1o6);
backwardChimera( m_001, m_011, m_021, vvy, vy2);
backwardInverseChimeraWithK(m_002, m_012, m_022, vvy, vy2, c18o1, c1o18);
backwardInverseChimeraWithK(m_100, m_110, m_120, vvy, vy2, c3o2, c2o3);
backwardChimera( m_101, m_111, m_121, vvy, vy2);
backwardInverseChimeraWithK(m_102, m_112, m_122, vvy, vy2, c9o2, c2o9);
backwardInverseChimeraWithK(m_200, m_210, m_220, vvy, vy2, c6o1, c1o6);
backwardChimera( m_201, m_211, m_221, vvy, vy2);
backwardInverseChimeraWithK(m_202, m_212, m_222, vvy, vy2, c18o1, c1o18);
backwardInverseChimeraWithK(m000, m010, m020, vvy, vy2, c6o1, c1o6);
backwardChimera( m001, m011, m021, vvy, vy2);
backwardInverseChimeraWithK(m002, m012, m022, vvy, vy2, c18o1, c1o18);
backwardInverseChimeraWithK(m100, m110, m120, vvy, vy2, c3o2, c2o3);
backwardChimera( m101, m111, m121, vvy, vy2);
backwardInverseChimeraWithK(m102, m112, m122, vvy, vy2, c9o2, c2o9);
backwardInverseChimeraWithK(m200, m210, m220, vvy, vy2, c6o1, c1o6);
backwardChimera( m201, m211, m221, vvy, vy2);
backwardInverseChimeraWithK(m202, m212, m222, vvy, vy2, c18o1, c1o18);
////////////////////////////////////////////////////////////////////////////////////
// Z - Dir
backwardInverseChimeraWithK(m_000, m_001, m_002, vvz, vz2, c36o1, c1o36);
backwardInverseChimeraWithK(m_010, m_011, m_012, vvz, vz2, c9o1, c1o9);
backwardInverseChimeraWithK(m_020, m_021, m_022, vvz, vz2, c36o1, c1o36);
backwardInverseChimeraWithK(m_100, m_101, m_102, vvz, vz2, c9o1, c1o9);
backwardInverseChimeraWithK(m_110, m_111, m_112, vvz, vz2, c9o4, c4o9);
backwardInverseChimeraWithK(m_120, m_121, m_122, vvz, vz2, c9o1, c1o9);
backwardInverseChimeraWithK(m_200, m_201, m_202, vvz, vz2, c36o1, c1o36);
backwardInverseChimeraWithK(m_210, m_211, m_212, vvz, vz2, c9o1, c1o9);
backwardInverseChimeraWithK(m_220, m_221, m_222, vvz, vz2, c36o1, c1o36);
distribution[dP00] = f_P00;
distribution[dM00] = f_M00;
distribution[d0P0] = f_0P0;
distribution[d0M0] = f_0M0;
distribution[d00P] = f_00P;
distribution[d00M] = f_00M;
distribution[dPP0] = f_PP0;
distribution[dMM0] = f_MM0;
distribution[dPM0] = f_PM0;
distribution[dMP0] = f_MP0;
distribution[dP0P] = f_P0P;
distribution[dM0M] = f_M0M;
distribution[dP0M] = f_P0M;
distribution[dM0P] = f_M0P;
distribution[d0PP] = f_0PP;
distribution[d0MM] = f_0MM;
distribution[d0PM] = f_0PM;
distribution[d0MP] = f_0MP;
distribution[d000] = f_000;
distribution[dPPP] = f_PPP;
distribution[dPMP] = f_PMP;
distribution[dPPM] = f_PPM;
distribution[dPMM] = f_PMM;
distribution[dMPP] = f_MPP;
distribution[dMMP] = f_MMP;
distribution[dMPM] = f_MPM;
distribution[dMMM] = f_MMM;
backwardInverseChimeraWithK(m000, m001, m002, vvz, vz2, c36o1, c1o36);
backwardInverseChimeraWithK(m010, m011, m012, vvz, vz2, c9o1, c1o9);
backwardInverseChimeraWithK(m020, m021, m022, vvz, vz2, c36o1, c1o36);
backwardInverseChimeraWithK(m100, m101, m102, vvz, vz2, c9o1, c1o9);
backwardInverseChimeraWithK(m110, m111, m112, vvz, vz2, c9o4, c4o9);
backwardInverseChimeraWithK(m120, m121, m122, vvz, vz2, c9o1, c1o9);
backwardInverseChimeraWithK(m200, m201, m202, vvz, vz2, c36o1, c1o36);
backwardInverseChimeraWithK(m210, m211, m212, vvz, vz2, c9o1, c1o9);
backwardInverseChimeraWithK(m220, m221, m222, vvz, vz2, c36o1, c1o36);
distribution[dP00] = fP00;
distribution[dM00] = fM00;
distribution[d0P0] = f0P0;
distribution[d0M0] = f0M0;
distribution[d00P] = f00P;
distribution[d00M] = f00M;
distribution[dPP0] = fPP0;
distribution[dMM0] = fMM0;
distribution[dPM0] = fPM0;
distribution[dMP0] = fMP0;
distribution[dP0P] = fP0P;
distribution[dM0M] = fM0M;
distribution[dP0M] = fP0M;
distribution[dM0P] = fM0P;
distribution[d0PP] = f0PP;
distribution[d0MM] = f0MM;
distribution[d0PM] = f0PM;
distribution[d0MP] = f0MP;
distribution[d000] = f000;
distribution[dPPP] = fPPP;
distribution[dPMP] = fPMP;
distribution[dPPM] = fPPM;
distribution[dPMM] = fPMM;
distribution[dMPP] = fMPP;
distribution[dMMP] = fMMP;
distribution[dMPM] = fMPM;
distribution[dMMM] = fMMM;
macroscopicValues.rho = drho;
macroscopicValues.vx = vvx;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment