diff --git a/applications/solvers/df0DFoam/createFields.H b/applications/solvers/df0DFoam/createFields.H index eefd02e3e..d4b73c372 100644 --- a/applications/solvers/df0DFoam/createFields.H +++ b/applications/solvers/df0DFoam/createFields.H @@ -2,7 +2,8 @@ Info<< "Reading thermophysical properties\n" << endl; -fluidThermo* pThermo = new hePsiThermo(mesh, word::null); +// fluidThermo* pThermo = new hePsiThermo(mesh, word::null); +fluidThermo* pThermo = new heRhoThermo(mesh, word::null); fluidThermo& thermo = *pThermo; thermo.validate(args.executable(), "ha"); diff --git a/applications/solvers/df0DFoam/df0DFoam.C b/applications/solvers/df0DFoam/df0DFoam.C index 8ca8f63f3..7dfd3644b 100644 --- a/applications/solvers/df0DFoam/df0DFoam.C +++ b/applications/solvers/df0DFoam/df0DFoam.C @@ -27,7 +27,8 @@ Description \*---------------------------------------------------------------------------*/ #include "dfChemistryModel.H" #include "CanteraMixture.H" -#include "hePsiThermo.H" +// #include "hePsiThermo.H" +#include "heRhoThermo.H" #ifdef USE_PYTORCH #include diff --git a/applications/solvers/dfHighSpeedFoam/createFields.H b/applications/solvers/dfHighSpeedFoam/createFields.H index b7aa52221..3446fecd6 100644 --- a/applications/solvers/dfHighSpeedFoam/createFields.H +++ b/applications/solvers/dfHighSpeedFoam/createFields.H @@ -2,8 +2,10 @@ Info<< "Reading thermophysical properties\n" << endl; -psiThermo* pThermo = new hePsiThermo(mesh, word::null); -psiThermo& thermo = *pThermo; +// psiThermo* pThermo = new hePsiThermo(mesh, word::null); +// psiThermo& thermo = *pThermo; +rhoThermo* pThermo = new heRhoThermo(mesh, word::null); +rhoThermo& thermo = *pThermo; //move from creatFieldRefs.H to createFields.H //p needed to be created before e diff --git a/applications/solvers/dfHighSpeedFoam/dfHighSpeedFoam.C b/applications/solvers/dfHighSpeedFoam/dfHighSpeedFoam.C index d52239ed5..ee59f7f68 100644 --- a/applications/solvers/dfHighSpeedFoam/dfHighSpeedFoam.C +++ b/applications/solvers/dfHighSpeedFoam/dfHighSpeedFoam.C @@ -32,7 +32,8 @@ Description #include "dfChemistryModel.H" #include "CanteraMixture.H" -#include "hePsiThermo.H" +// #include "hePsiThermo.H" +#include "heRhoThermo.H" #ifdef USE_PYTORCH #include @@ -47,7 +48,8 @@ Description #include "fvCFD.H" #include "dynamicFvMesh.H" -#include "psiThermo.H" +// #include "psiThermo.H" +#include "rhoThermo.H" #include "turbulentFluidThermoModel.H" #include "fixedRhoFvPatchScalarField.H" #include "directionInterpolate.H" diff --git a/applications/solvers/dfLowMachFoam/createFields.H b/applications/solvers/dfLowMachFoam/createFields.H index 911e7b5a3..065d301e9 100644 --- a/applications/solvers/dfLowMachFoam/createFields.H +++ b/applications/solvers/dfLowMachFoam/createFields.H @@ -2,7 +2,8 @@ Info<< "Reading thermophysical properties\n" << endl; -fluidThermo* pThermo = new hePsiThermo(mesh, word::null); +// fluidThermo* pThermo = new hePsiThermo(mesh, word::null); +fluidThermo* pThermo = new heRhoThermo(mesh, word::null); fluidThermo& thermo = *pThermo; thermo.validate(args.executable(), "ha"); diff --git a/applications/solvers/dfLowMachFoam/dfLowMachFoam.C b/applications/solvers/dfLowMachFoam/dfLowMachFoam.C index 437866bbe..72f33f65d 100644 --- a/applications/solvers/dfLowMachFoam/dfLowMachFoam.C +++ b/applications/solvers/dfLowMachFoam/dfLowMachFoam.C @@ -35,7 +35,8 @@ Description #include "dfChemistryModel.H" #include "CanteraMixture.H" -#include "hePsiThermo.H" +// #include "hePsiThermo.H" +#include "heRhoThermo.H" #ifdef USE_PYTORCH #include diff --git a/applications/solvers/dfSprayFoam/createFields.H b/applications/solvers/dfSprayFoam/createFields.H index f9e0393ef..fbcd635f2 100644 --- a/applications/solvers/dfSprayFoam/createFields.H +++ b/applications/solvers/dfSprayFoam/createFields.H @@ -1,7 +1,8 @@ #include "readGravitationalAcceleration.H" Info<< "Reading thermophysical properties\n" << endl; -fluidThermo* pThermo = new hePsiThermo(mesh, word::null); +// fluidThermo* pThermo = new hePsiThermo(mesh, word::null); +fluidThermo* pThermo = new heRhoThermo(mesh, word::null); fluidThermo& thermo = *pThermo; thermo.validate(args.executable(), "ha"); diff --git a/applications/solvers/dfSprayFoam/dfSprayFoam.C b/applications/solvers/dfSprayFoam/dfSprayFoam.C index 213df8938..ff044db6d 100644 --- a/applications/solvers/dfSprayFoam/dfSprayFoam.C +++ b/applications/solvers/dfSprayFoam/dfSprayFoam.C @@ -31,7 +31,8 @@ Description \*---------------------------------------------------------------------------*/ #include "dfChemistryModel.H" #include "CanteraMixture.H" -#include "hePsiThermo.H" +// #include "hePsiThermo.H" +#include "heRhoThermo.H" #include "turbulentFluidThermoModel.H" #ifdef USE_PYTORCH diff --git a/applications/solvers/dfSprayFoam/pEqn.H b/applications/solvers/dfSprayFoam/pEqn.H index e303e0181..19e3a0b92 100644 --- a/applications/solvers/dfSprayFoam/pEqn.H +++ b/applications/solvers/dfSprayFoam/pEqn.H @@ -3,6 +3,10 @@ rho = max(rho, rhoMin); rho = min(rho, rhoMax); rho.relax(); +// Thermodynamic density needs to be updated by psi*d(p) after the +// pressure solution +const volScalarField psip0(psi*p); + volScalarField rAU(1.0/UEqn.A()); surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU)); volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p)); @@ -93,6 +97,7 @@ else p.relax(); // Recalculate density from the relaxed pressure +thermo.correctRho(psi*p - psip0); rho = thermo.rho(); rho = max(rho, rhoMin); rho = min(rho, rhoMax); diff --git a/src/dfCanteraMixture/CanteraMixture.H b/src/dfCanteraMixture/CanteraMixture.H index 3556e4f58..216cdd144 100644 --- a/src/dfCanteraMixture/CanteraMixture.H +++ b/src/dfCanteraMixture/CanteraMixture.H @@ -121,6 +121,15 @@ public: return CanteraGas_->density()/p; } + scalar rho + ( + const scalar& p, + const scalar& T + ) const + { + return CanteraGas_->density(); + } + scalar mu ( const scalar& p, diff --git a/src/dfCanteraMixture/makeThermos.C b/src/dfCanteraMixture/makeThermos.C index f7fffaea4..99ad2c4bc 100644 --- a/src/dfCanteraMixture/makeThermos.C +++ b/src/dfCanteraMixture/makeThermos.C @@ -24,9 +24,11 @@ License \*---------------------------------------------------------------------------*/ #include "basicThermo.H" -#include "psiThermo.H" +// #include "psiThermo.H" +#include "rhoThermo.H" #include "CanteraMixture.H" -#include "hePsiThermo.H" +// #include "hePsiThermo.H" +#include "heRhoThermo.H" #include "addToRunTimeSelectionTable.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -71,7 +73,8 @@ namespace Foam addThermoPhysicsThermo(fluidThermo, CThermo##Mixture##ThermoPhys) -makeThermoPhysicsThermos(hePsiThermo, CanteraMixture, psiThermo); +// makeThermoPhysicsThermos(hePsiThermo, CanteraMixture, psiThermo); +makeThermoPhysicsThermos(heRhoThermo, CanteraMixture, rhoThermo); // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/dfChemistryModel/dfChemistryModel.C b/src/dfChemistryModel/dfChemistryModel.C index 984c9364d..f97c9a4e6 100644 --- a/src/dfChemistryModel/dfChemistryModel.C +++ b/src/dfChemistryModel/dfChemistryModel.C @@ -67,9 +67,12 @@ Foam::dfChemistryModel::dfChemistryModel alpha_(const_cast(thermo.alpha())), T_(thermo.T()), p_(thermo.p()), - rho_(mesh_.objectRegistry::lookupObject("rho")), - mu_(const_cast(dynamic_cast(thermo).mu()())), - psi_(const_cast(dynamic_cast(thermo).psi())), + // rho_(mesh_.objectRegistry::lookupObject("rho")), + rho_(const_cast(dynamic_cast(thermo).rho())), + // mu_(const_cast(dynamic_cast(thermo).mu()())), + // psi_(const_cast(dynamic_cast(thermo).psi())), + mu_(const_cast(dynamic_cast(thermo).mu()())), + psi_(const_cast(dynamic_cast(thermo).psi())), Qdot_ ( IOobject @@ -385,6 +388,8 @@ void Foam::dfChemistryModel::correctThermo() psi_[celli] = mixture_.psi(p_[celli],T_[celli]); + rho_[celli] = mixture_.rho(p_[celli],T_[celli]); + mu_[celli] = mixture_.CanteraTransport()->viscosity(); // Pa-s alpha_[celli] = mixture_.CanteraTransport()->thermalConductivity()/(CanteraGas_->cp_mass()); // kg/(m*s) @@ -418,7 +423,7 @@ void Foam::dfChemistryModel::correctThermo() const volScalarField::Boundary& pBf = p_.boundaryField(); - const volScalarField::Boundary& rhoBf = rho_.boundaryField(); + volScalarField::Boundary& rhoBf = rho_.boundaryFieldRef(); volScalarField::Boundary& TBf = T_.boundaryFieldRef(); @@ -433,7 +438,7 @@ void Foam::dfChemistryModel::correctThermo() forAll(T_.boundaryField(), patchi) { const fvPatchScalarField& pp = pBf[patchi]; - const fvPatchScalarField& prho = rhoBf[patchi]; + fvPatchScalarField& prho = rhoBf[patchi]; fvPatchScalarField& pT = TBf[patchi]; fvPatchScalarField& ppsi = psiBf[patchi]; fvPatchScalarField& ph = hBf[patchi]; @@ -454,6 +459,8 @@ void Foam::dfChemistryModel::correctThermo() ppsi[facei] = mixture_.psi(pp[facei],pT[facei]); + prho[facei] = mixture_.rho(pp[facei],pT[facei]); + pmu[facei] = mixture_.CanteraTransport()->viscosity(); palpha[facei] = mixture_.CanteraTransport()->thermalConductivity()/(CanteraGas_->cp_mass()); @@ -495,6 +502,8 @@ void Foam::dfChemistryModel::correctThermo() ppsi[facei] = mixture_.psi(pp[facei],pT[facei]); + prho[facei] = mixture_.rho(pp[facei],pT[facei]); + pmu[facei] = mixture_.CanteraTransport()->viscosity(); palpha[facei] = mixture_.CanteraTransport()->thermalConductivity()/(CanteraGas_->cp_mass()); diff --git a/src/dfChemistryModel/dfChemistryModel.H b/src/dfChemistryModel/dfChemistryModel.H index 1304eff2b..32bd94f09 100644 --- a/src/dfChemistryModel/dfChemistryModel.H +++ b/src/dfChemistryModel/dfChemistryModel.H @@ -62,7 +62,8 @@ SourceFiles #include "scalarField.H" #include "volFields.H" #include "hashedWordList.H" -#include "psiThermo.H" +// #include "psiThermo.H" +#include "rhoThermo.H" #include "physicoChemicalConstants.H" // for R #include "ChemistryProblem.H" #include "ChemistrySolution.H" @@ -121,7 +122,7 @@ public IOdictionary volScalarField& alpha_; volScalarField& T_; const volScalarField& p_; - const volScalarField& rho_; + volScalarField& rho_; volScalarField& mu_; volScalarField& psi_; // heat release rate, [J/m^3/s] diff --git a/src/dfCombustionModels/FGM/baseFGM/baseFGM.C b/src/dfCombustionModels/FGM/baseFGM/baseFGM.C index 2f2f1879a..1f8cbca14 100644 --- a/src/dfCombustionModels/FGM/baseFGM/baseFGM.C +++ b/src/dfCombustionModels/FGM/baseFGM/baseFGM.C @@ -44,7 +44,8 @@ Foam::combustionModels::baseFGM::baseFGM combustion_(this->coeffs().lookupOrDefault("combustion", false)), solveEnthalpy_(this->coeffs().lookupOrDefault("solveEnthalpy", false)), flameletT_(this->coeffs().lookupOrDefault("flameletT", false)), - psi_(const_cast(dynamic_cast(thermo).psi())), + // psi_(const_cast(dynamic_cast(thermo).psi())), + psi_(const_cast(dynamic_cast(thermo).psi())), Wt_ ( IOobject diff --git a/src/dfCombustionModels/PaSR/PaSR.C b/src/dfCombustionModels/PaSR/PaSR.C index 69e9838ee..30f04e51f 100644 --- a/src/dfCombustionModels/PaSR/PaSR.C +++ b/src/dfCombustionModels/PaSR/PaSR.C @@ -37,7 +37,8 @@ Foam::combustionModels::PaSR::PaSR ) : laminar(modelType, thermo, turb, combustionProperties), - mu_(const_cast(dynamic_cast(thermo).mu()())), + // mu_(const_cast(dynamic_cast(thermo).mu()())), + mu_(const_cast(dynamic_cast(thermo).mu()())), p_(this->thermo().p()), T_(this->thermo().T()), mixingScaleDict_(this->coeffs().subDict("mixingScale")),