Program Listing for File hydro.cpp
↰ Return to documentation for file (hydro/hydro.cpp
)
// ***********************************************************************************
// Idefix MHD astrophysical code
// Copyright(C) 2020-2022 Geoffroy R. J. Lesur <geoffroy.lesur@univ-grenoble-alpes.fr>
// and other code contributors
// Licensed under CeCILL 2.1 License, see COPYING for more information
// ***********************************************************************************
#include <string>
#include "idefix.hpp"
#include "hydro.hpp"
#include "dataBlock.hpp"
void Hydro::Init(Input &input, Grid &grid, DataBlock *datain) {
idfx::pushRegion("Hydro::Init");
// Save the datablock to which we are attached from now on
this->data = datain;
#if ORDER < 1 || ORDER > 4
IDEFIX_ERROR("Reconstruction at chosen order is not implemented. Check your definitions file");
#endif
#if HAVE_ENERGY
this->gamma = input.GetOrSet<real>("Hydro","gamma",0, 5.0/3.0);
#endif
#ifdef ISOTHERMAL
std::string isoString = input.Get<std::string>("Hydro","csiso",0);
if(isoString.compare("constant") == 0) {
this->haveIsoSoundSpeed = Constant;
this->isoSoundSpeed = input.Get<real>("Hydro","csiso",1);
} else if(isoString.compare("userdef") == 0) {
this->haveIsoSoundSpeed = UserDefFunction;
} else {
IDEFIX_ERROR("csiso admits only constant or userdef entries");
}
#else
// set the isothermal soundspeed, even though it will not be used
this->isoSoundSpeed = -1.0;
#endif
// read Solver from input file
std::string solverString = input.Get<std::string>("Hydro","solver",0);
if (solverString.compare("tvdlf") == 0) {
mySolver = TVDLF;
} else if (solverString.compare("hll") == 0) {
mySolver = HLL;
#if MHD == YES
} else if (solverString.compare("hlld") == 0) {
mySolver = HLLD;
#else
} else if (solverString.compare("hllc") == 0) {
mySolver = HLLC;
#endif
} else if (solverString.compare("roe") == 0) {
mySolver = ROE;
} else {
std::stringstream msg;
#if MHD == YES
msg << "Unknown MHD solver type " << solverString;
#else
msg << "Unknown HD solver type " << solverString;
#endif
IDEFIX_ERROR(msg);
}
// Shock flattening
this->haveShockFlattening = input.CheckEntry("Hydro","shockFlattening")>=0;
// Source terms (always activated when non-cartesian geometry because of curvature source terms)
#if GEOMETRY == CARTESIAN
this->haveSourceTerms = false;
#else
this->haveSourceTerms = true;
#endif
// Check whether we have rotation
int rotation = input.CheckEntry("Hydro","rotation");
if(rotation>=0 ) {
this->haveSourceTerms = true;
this->haveRotation = true;
if(rotation >1 ) IDEFIX_ERROR("Rotation needs a a single rotation velocity (Omega_Z) \
in idefix.ini");
this->OmegaZ = input.Get<real>("Hydro","rotation",0);
}
// Check whether we have shearing box
int shearingbox = input.CheckEntry("Hydro","shearingBox");
if(shearingbox>=0 ) {
this->haveShearingBox = true;
this->haveSourceTerms = true;
if(shearingbox != 1) {
IDEFIX_ERROR("Shearing box needs a scalar value for the shear rate in idefix.ini");
}
this->sbS = input.Get<real>("Hydro","shearingBox",0);
// Get box size
this->sbLx = grid.xend[IDIR] - grid.xbeg[IDIR];
}
// Gravitational potential (deprecated, use [Gravity] block in your input file)
if(input.CheckEntry("Hydro","gravPotential")>=0) {
IDEFIX_DEPRECATED("gravPotential in [Hydro] block is deprecated in the input file. "
"Use a [Gravity] block instead.");
std::string potentialString = input.Get<std::string>("Hydro","gravPotential",0);
if(potentialString.compare("userdef") == 0) {
data->gravity.haveUserDefPotential = true;
data->gravity.havePotential = true;
data->gravity.Init(input,data);
} else {
IDEFIX_ERROR("Unknown type of gravitational potential in idefix.ini. "
"Only userdef is implemented");
}
}
// Body Force (deprecated, use [Gravity] block in your input file)
if(input.CheckEntry("Hydro","bodyForce")>=0) {
IDEFIX_DEPRECATED("bodyForce in [Hydro] block is deprecated in the input file. "
"Use a [Gravity] block instead.");
std::string potentialString = input.Get<std::string>("Hydro","bodyForce",0);
if(potentialString.compare("userdef") == 0) {
data->gravity.haveBodyForce = true;
data->gravity.Init(input,data);
} else {
IDEFIX_ERROR("Unknown type of body force in idefix.ini. "
"Only userdef is implemented");
}
}
// Do we use fargo? (deprecated, use a full fargo block in idefix.ini,
// instead of including it in hydo)
if(input.CheckEntry("Hydro","fargo")>=0) {
IDEFIX_DEPRECATED("Fargo in [Hydro] block is deprecated in the input file. "
"Use a [Fargo] block instead.");
data->haveFargo = true;
data->fargo.Init(input, this->data);
}
// Parabolic terms
// Check whether viscosity is enabled, if so, init a viscosity object
if(input.CheckEntry("Hydro","viscosity")>=0) {
std::string opType = input.Get<std::string>("Hydro","viscosity",0);
if(opType.compare("explicit") == 0 ) {
haveExplicitParabolicTerms = true;
viscosityStatus.isExplicit = true;
} else if(opType.compare("rkl") == 0 ) {
haveRKLParabolicTerms = true;
viscosityStatus.isRKL = true;
} else {
std::stringstream msg;
msg << "Unknown integration type for viscosity: " << opType;
IDEFIX_ERROR(msg);
}
this->viscosity.Init(input, grid, this);
}
// Check whether thermal diffusion is enabled, if so, init a thermal diffusion object
if(input.CheckEntry("Hydro","TDiffusion")>=0) {
std::string opType = input.Get<std::string>("Hydro","TDiffusion",0);
if(opType.compare("explicit") == 0 ) {
haveExplicitParabolicTerms = true;
thermalDiffusionStatus.isExplicit = true;
} else if(opType.compare("rkl") == 0 ) {
haveRKLParabolicTerms = true;
thermalDiffusionStatus.isRKL = true;
} else {
std::stringstream msg;
msg << "Unknown integration type for thermal diffusion: " << opType;
IDEFIX_ERROR(msg);
}
this->thermalDiffusion.Init(input, grid, this);
}
#if MHD == YES
if(input.CheckEntry("Hydro","resistivity")>=0 ||
input.CheckEntry("Hydro","ambipolar")>=0 ||
input.CheckEntry("Hydro","hall")>=0 ) {
//
this->haveCurrent = true;
if(input.CheckEntry("Hydro","resistivity")>=0) {
std::string opType = input.Get<std::string>("Hydro","resistivity",0);
if(opType.compare("explicit") == 0 ) {
haveExplicitParabolicTerms = true;
resistivityStatus.isExplicit = true;
needExplicitCurrent = true;
} else if(opType.compare("rkl") == 0 ) {
haveRKLParabolicTerms = true;
resistivityStatus.isRKL = true;
needRKLCurrent = true;
} else {
std::stringstream msg;
msg << "Unknown integration type for resistivity: " << opType;
IDEFIX_ERROR(msg);
}
if(input.Get<std::string>("Hydro","resistivity",1).compare("constant") == 0) {
this->etaO = input.Get<real>("Hydro","resistivity",2);
resistivityStatus.status = Constant;
} else if(input.Get<std::string>("Hydro","resistivity",1).compare("userdef") == 0) {
resistivityStatus.status = UserDefFunction;
} else {
IDEFIX_ERROR("Unknown resistivity definition in idefix.ini. "
"Can only be constant or userdef.");
}
}
if(input.CheckEntry("Hydro","ambipolar")>=0) {
std::string opType = input.Get<std::string>("Hydro","ambipolar",0);
if(opType.compare("explicit") == 0 ) {
haveExplicitParabolicTerms = true;
ambipolarStatus.isExplicit = true;
needExplicitCurrent = true;
} else if(opType.compare("rkl") == 0 ) {
haveRKLParabolicTerms = true;
ambipolarStatus.isRKL = true;
needRKLCurrent = true;
} else {
std::stringstream msg;
msg << "Unknown integration type for ambipolar: " << opType;
IDEFIX_ERROR(msg);
}
if(input.Get<std::string>("Hydro","ambipolar",1).compare("constant") == 0) {
this->xA = input.Get<real>("Hydro","ambipolar",2);
ambipolarStatus.status = Constant;
} else if(input.Get<std::string>("Hydro","ambipolar",1).compare("userdef") == 0) {
ambipolarStatus.status = UserDefFunction;
} else {
IDEFIX_ERROR("Unknown ambipolar definition in idefix.ini. "
"Can only be constant or userdef.");
}
}
if(input.CheckEntry("Hydro","hall")>=0) {
// Check consistency
if(mySolver != HLL )
IDEFIX_ERROR("Hall effect is only compatible with HLL Riemann solver.");
std::string opType = input.Get<std::string>("Hydro","hall",0);
if(opType.compare("explicit") == 0 ) {
hallStatus.isExplicit = true;
needExplicitCurrent = true;
} else if(opType.compare("rkl") == 0 ) {
IDEFIX_ERROR("RKL inegration is incompatible with Hall");
} else {
std::stringstream msg;
msg << "Unknown integration type for hall: " << opType;
IDEFIX_ERROR(msg);
}
if(input.Get<std::string>("Hydro","hall",1).compare("constant") == 0) {
this->xH = input.Get<real>("Hydro","hall",2);
hallStatus.status = Constant;
} else if(input.Get<std::string>("Hydro","hall",1).compare("userdef") == 0) {
hallStatus.status = UserDefFunction;
} else {
IDEFIX_ERROR("Unknown Hall definition in idefix.ini. Can only be constant or userdef.");
}
}
}
#endif // MHD
// Do we have to take care of the axis?
if(data->haveAxis) {
this->myAxis.Init(grid, this);
this->haveAxis = true;
}
// ALLOCATION SECION ///////////////////
// We now allocate the fields required by the hydro solver
Vc = IdefixArray4D<real>("Hydro_Vc", NVAR,
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
Uc = IdefixArray4D<real>("Hydro_Uc", NVAR,
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
data->states["current"].PushArray(Uc, State::center, "Hydro_Uc");
InvDt = IdefixArray3D<real>("Hydro_InvDt",
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
cMax = IdefixArray3D<real>("Hydro_cMax",
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
dMax = IdefixArray3D<real>("Hydro_dMax",
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
FluxRiemann = IdefixArray4D<real>("Hydro_FluxRiemann", NVAR,
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
#if MHD == YES
Vs = IdefixArray4D<real>("Hydro_Vs", DIMENSIONS,
data->np_tot[KDIR]+KOFFSET, data->np_tot[JDIR]+JOFFSET, data->np_tot[IDIR]+IOFFSET);
#ifdef EVOLVE_VECTOR_POTENTIAL
#if DIMENSIONS == 1
IDEFIX_ERROR("EVOLVE_VECTOR_POTENTIAL is not compatible with 1D MHD");
#else
Ve = IdefixArray4D<real>("Hydro_Ve", AX3e+1,
data->np_tot[KDIR]+KOFFSET, data->np_tot[JDIR]+JOFFSET, data->np_tot[IDIR]+IOFFSET);
data->states["current"].PushArray(Ve, State::center, "Hydro_Ve");
#endif
#else // EVOLVE_VECTOR_POTENTIAL
data->states["current"].PushArray(Vs, State::center, "Hydro_Vs");
#endif // EVOLVE_VECTOR_POTENTIAL
this->emf.Init(input, this);
#endif
// Allocate sound speed array if needed
if(this->haveIsoSoundSpeed == UserDefFunction) {
this->isoSoundSpeedArray = IdefixArray3D<real>("Hydro_csIso",
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
}
if(this->haveCurrent) {
// Allocate current (when hydro needs it)
J = IdefixArray4D<real>("Hydro_J", 3,
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
}
// Allocate nonideal MHD effects array when a user-defined function is used
if(this->resistivityStatus.status == UserDefFunction)
etaOhmic = IdefixArray3D<real>("Hydro_etaOhmic",
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
if(this->ambipolarStatus.status == UserDefFunction)
xAmbipolar = IdefixArray3D<real>("Hydro_xAmbipolar",
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
if(this->hallStatus.status == UserDefFunction)
xHall = IdefixArray3D<real>("Hydro_xHall",
data->np_tot[KDIR], data->np_tot[JDIR], data->np_tot[IDIR]);
// Fill the names of the fields
for(int i = 0 ; i < NVAR ; i++) {
switch(i) {
case RHO:
VcName.push_back("RHO");
break;
case VX1:
VcName.push_back("VX1");
break;
case VX2:
VcName.push_back("VX2");
break;
case VX3:
VcName.push_back("VX3");
break;
case BX1:
VcName.push_back("BX1");
break;
case BX2:
VcName.push_back("BX2");
break;
case BX3:
VcName.push_back("BX3");
break;
#if HAVE_ENERGY
case PRS:
VcName.push_back("PRS");
break;
#endif
default:
VcName.push_back("Vc_"+std::to_string(i));
}
}
for(int i = 0 ; i < DIMENSIONS ; i++) {
switch(i) {
case 0:
VsName.push_back("BX1s");
break;
case 1:
VsName.push_back("BX2s");
break;
case 2:
VsName.push_back("BX3s");
break;
default:
VsName.push_back("Vs_"+std::to_string(i));
}
}
#ifdef EVOLVE_VECTOR_POTENTIAL
#if DIMENSIONS < 3
VeName.push_back("AX3e");
#else
for(int i = 0 ; i < DIMENSIONS ; i++) {
switch(i) {
case 0:
VeName.push_back("AX1e");
break;
case 1:
VeName.push_back("AX2e");
break;
case 2:
VeName.push_back("AX3e");
break;
default:
VeName.push_back("Ve_"+std::to_string(i));
}
}
#endif
#endif
// Initialise boundary conditions
boundary.Init(input, grid, this);
// Init shock flattening
if(haveShockFlattening) {
this->shockFlattening = ShockFlattening(this,input.Get<real>("Hydro","shockFlattening",0));
}
idfx::popRegion();
}
void Hydro::EnrollIsoSoundSpeed(IsoSoundSpeedFunc myFunc) {
if(this->haveIsoSoundSpeed != UserDefFunction) {
IDEFIX_WARNING("Isothermal sound speed enrollment requires Hydro/csiso "
" to be set to userdef in .ini file");
}
#if HAVE_ENERGY
IDEFIX_ERROR("Isothermal sound speed enrollment requires ISOTHERMAL to be defined in"
"definitions.hpp");
#endif
this->isoSoundSpeedFunc = myFunc;
}
void Hydro::EnrollUserDefBoundary(UserDefBoundaryFunc myFunc) {
// This is a proxy for userdef enrollment
boundary.EnrollUserDefBoundary(myFunc);
}
void Hydro::EnrollFluxBoundary(UserDefBoundaryFunc myFunc) {
// This is a proxy for userdef enrollment
boundary.EnrollFluxBoundary(myFunc);
}
void Hydro::EnrollInternalBoundary(InternalBoundaryFunc myFunc) {
// This is a proxy for userdef enrollment
boundary.EnrollInternalBoundary(myFunc);
}
void Hydro::EnrollEmfBoundary(EmfBoundaryFunc myFunc) {
#if MHD == NO
IDEFIX_ERROR("This function can only be used with the MHD solver.");
#endif
this->emfBoundaryFunc = myFunc;
this->haveEmfBoundary = true;
}
// Deprecated function
void Hydro::EnrollGravPotential(GravPotentialFunc myFunc) {
IDEFIX_DEPRECATED("Calling Hydro::EnrollGravPotential is deprecated."
"Use Gravity::EnrollPotential");
data->gravity.EnrollPotential(myFunc);
}
// Deprecated function
void Hydro::EnrollBodyForce(BodyForceFunc myFunc) {
IDEFIX_DEPRECATED("Calling Hydro::EnrollBodyForce is deprecated."
"Use Gravity::EnrollBodyForce");
data->gravity.EnrollBodyForce(myFunc);
}
void Hydro::EnrollUserSourceTerm(SrcTermFunc myFunc) {
this->userSourceTerm = myFunc;
this->haveUserSourceTerm = true;
this->haveSourceTerms = true;
}
void Hydro::EnrollOhmicDiffusivity(DiffusivityFunc myFunc) {
#if MHD == NO
IDEFIX_ERROR("This function can only be used with the MHD solver.");
#endif
if(this->resistivityStatus.status < UserDefFunction) {
IDEFIX_WARNING("Ohmic diffusivity enrollment requires Hydro/Resistivity "
"to be set to userdef in .ini file");
}
this->ohmicDiffusivityFunc = myFunc;
}
void Hydro::EnrollAmbipolarDiffusivity(DiffusivityFunc myFunc) {
#if MHD == NO
IDEFIX_ERROR("This function can only be used with the MHD solver.");
#endif
if(this->ambipolarStatus.status < UserDefFunction) {
IDEFIX_WARNING("Ambipolar diffusivity enrollment requires Hydro/Ambipolar "
"to be set to userdef in .ini file");
}
this->ambipolarDiffusivityFunc = myFunc;
}
void Hydro::EnrollHallDiffusivity(DiffusivityFunc myFunc) {
#if MHD == NO
IDEFIX_ERROR("This function can only be used with the MHD solver.");
#endif
if(this->hallStatus.status < UserDefFunction) {
IDEFIX_WARNING("Hall diffusivity enrollment requires Hydro/Hall "
"to be set to userdef in .ini file");
}
this->hallDiffusivityFunc = myFunc;
}
Hydro::Hydro() {
// Do nothing !
}
real Hydro::GetGamma() {
return(this->gamma);
}
void Hydro::ResetStage() {
// Reset variables required at the beginning of each stage
// (essentially linked to timestep evaluation)
idfx::pushRegion("Hydro::ResetStage");
IdefixArray3D<real> InvDt=this->InvDt;
idefix_for("HydroResetStage",0,data->np_tot[KDIR],0,data->np_tot[JDIR],0,data->np_tot[IDIR],
KOKKOS_LAMBDA (int k, int j, int i) {
InvDt(k,j,i) = ZERO_F;
});
idfx::popRegion();
}
void Hydro::ShowConfig() {
idfx::cout << "Hydro: ";
#if MHD == YES
idfx::cout << "solving MHD equations." << std::endl;
#ifdef EVOLVE_VECTOR_POTENTIAL
idfx::cout << "Hydro: Using EXPERIMENTAL vector potential formulation for MHD." << std::endl;
#endif
#else
idfx::cout << "solving HD equations." << std::endl;
#endif
idfx::cout << "Hydro: Reconstruction: ";
#if ORDER == 1
idfx::cout << "1st order (donor cell)" << std::endl;
#elif ORDER == 2
idfx::cout << "2nd order (PLM Van Leer)" << std::endl;
#elif ORDER == 3
idfx::cout << "3rd order (LimO3)" << std::endl;
#elif ORDER == 4
idfx::cout << "4th order (PPM)" << std::endl;
#endif
#if HAVE_ENERGY
idfx::cout << "Hydro: EOS: ideal with gamma=" << this->gamma << std::endl;
#endif
#ifdef ISOTHERMAL
if(haveIsoSoundSpeed == Constant) {
idfx::cout << "Hydro: EOS: isothermal with cs=" << isoSoundSpeed << "." << std::endl;
} else if(haveIsoSoundSpeed == UserDefFunction) {
idfx::cout << "Hydro: EOS: isothermal with user-defined cs function." << std::endl;
if(!isoSoundSpeedFunc) {
IDEFIX_ERROR("No user-defined isothermal sound speed function has been enrolled.");
}
}
#endif// ISOTHERMAL
idfx::cout << "Hydro: Riemann solver: ";
switch(mySolver) {
case TVDLF:
idfx::cout << "tvdlf." << std::endl;
break;
case HLL:
idfx::cout << "hll." << std::endl;
break;
#if MHD==YES
case HLLD:
idfx::cout << "hlld." << std::endl;
break;
#else
case HLLC:
idfx::cout << "hllc." << std::endl;
break;
#endif
case ROE:
idfx::cout << "roe." << std::endl;
break;
default:
IDEFIX_ERROR("Unknown Riemann solver");
}
if(haveRotation) {
idfx::cout << "Hydro: Rotation ENABLED with Omega=" << this->OmegaZ << std::endl;
}
if(haveShearingBox) {
idfx::cout << "Hydro: ShearingBox ENABLED with shear rate= " << this->sbS
<< " and Lx= " << sbLx << std::endl;
}
if(haveShockFlattening) {
idfx::cout << "Hydro: Shock Flattening ENABLED." << std::endl;
}
if(resistivityStatus.status != Disabled) {
if(resistivityStatus.status == Constant) {
idfx::cout << "Hydro: Ohmic resistivity ENABLED with constant resistivity eta="
<< etaO << std::endl;
} else if(resistivityStatus.status == UserDefFunction) {
idfx::cout << "Hydro: Ohmic resistivity ENABLED with user-defined resistivity function."
<< std::endl;
if(!ohmicDiffusivityFunc) {
IDEFIX_ERROR("No user-defined Ihmic resistivity function has been enrolled.");
}
} else {
IDEFIX_ERROR("Unknown Ohmic resistivity mode");
}
if(resistivityStatus.isExplicit) {
idfx::cout << "Hydro: Ohmic resistivity uses an explicit time integration." << std::endl;
} else if(resistivityStatus.isRKL) {
idfx::cout << "Hydro: Ohmic resistivity uses a Runge-Kutta-Legendre time integration."
<< std::endl;
} else {
IDEFIX_ERROR("Unknown time integrator for Ohmic resistivity");
}
}
if(ambipolarStatus.status != Disabled) {
if(ambipolarStatus.status == Constant) {
idfx::cout << "Hydro: Ambipolar diffusion ENABLED with constant diffusivity xA="
<< xA << std::endl;
} else if(ambipolarStatus.status == UserDefFunction) {
idfx::cout << "Hydro: Ambipolar diffusion ENABLED with user-defined diffusivity function."
<< std::endl;
if(!ambipolarDiffusivityFunc) {
IDEFIX_ERROR("No user-defined ambipolar diffusion function has been enrolled.");
}
} else {
IDEFIX_ERROR("Unknown Ambipolar diffusion mode");
}
if(ambipolarStatus.isExplicit) {
idfx::cout << "Hydro: Ambipolar diffusion uses an explicit time integration." << std::endl;
} else if(ambipolarStatus.isRKL) {
idfx::cout << "Hydro: Ambipolar diffusion uses a Runge-Kutta-Legendre time integration."
<< std::endl;
} else {
IDEFIX_ERROR("Unknown time integrator for ambipolar diffusion");
}
}
if(hallStatus.status != Disabled) {
if(hallStatus.status == Constant) {
idfx::cout << "Hydro: Hall effect ENABLED with constant diffusivity xH="
<< xH << std::endl;
} else if(hallStatus.status == UserDefFunction) {
idfx::cout << "Hydro: Hall effect ENABLED with user-defined diffusivity function."
<< std::endl;
if(!hallDiffusivityFunc) {
IDEFIX_ERROR("No user-defined Hall diffusivity function has been enrolled.");
}
} else {
IDEFIX_ERROR("Unknown Hall effect mode");
}
if(hallStatus.isExplicit) {
idfx::cout << "Hydro: Hall effect uses an explicit time integration." << std::endl;
} else {
IDEFIX_ERROR("Unknown time integrator for Hall effect");
}
}
if(emfBoundaryFunc) {
idfx::cout << "Hydro: user-defined EMF boundaries ENABLED." << std::endl;
}
if(userSourceTerm) {
idfx::cout << "Hydro: user-defined source terms ENABLED." << std::endl;
}
#if MHD == YES
emf.ShowConfig();
#endif
if(viscosityStatus.isExplicit || viscosityStatus.isRKL) {
viscosity.ShowConfig();
}
if(thermalDiffusionStatus.isExplicit || thermalDiffusionStatus.isRKL) {
thermalDiffusion.ShowConfig();
}
if(haveAxis) {
myAxis.ShowConfig();
}
}