From ee765506530e8c2513ba5d6f7ea5355e4d2fc94c Mon Sep 17 00:00:00 2001 From: Henry Date: Sun, 29 Sep 2013 22:07:36 +0100 Subject: [PATCH] chemistryModel/laminar/PaSR: Add support for LTS to laminar model and to PaSR by derivation --- src/combustionModels/PaSR/PaSR.C | 108 ++---------------- src/combustionModels/PaSR/PaSR.H | 11 +- src/combustionModels/laminar/laminar.C | 33 +++++- src/combustionModels/laminar/laminar.H | 7 +- .../basicChemistryModel/basicChemistryModel.H | 4 + .../chemistryModel/chemistryModel.C | 38 ++++-- .../chemistryModel/chemistryModel.H | 9 ++ 7 files changed, 93 insertions(+), 117 deletions(-) diff --git a/src/combustionModels/PaSR/PaSR.C b/src/combustionModels/PaSR/PaSR.C index 14ae46a112..dcc3c9e315 100644 --- a/src/combustionModels/PaSR/PaSR.C +++ b/src/combustionModels/PaSR/PaSR.C @@ -35,7 +35,7 @@ Foam::combustionModels::PaSR::PaSR const fvMesh& mesh ) : - Type(modelType, mesh), + laminar(modelType, mesh), Cmix_(readScalar(this->coeffs().lookup("Cmix"))), turbulentReaction_(this->coeffs().lookup("turbulentReaction")), kappa_ @@ -50,21 +50,8 @@ Foam::combustionModels::PaSR::PaSR ), mesh, dimensionedScalar("kappa", dimless, 0.0) - ), - integrateReactionRate_ - ( - this->coeffs().lookupOrDefault("integrateReactionRate", true) ) -{ - if (integrateReactionRate_) - { - Info<< " using integrated reaction rate" << endl; - } - else - { - Info<< " using instantaneous reaction rate" << endl; - } -} +{} // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // @@ -76,28 +63,12 @@ Foam::combustionModels::PaSR::~PaSR() // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // -template -Foam::tmp Foam::combustionModels::PaSR::tc() const -{ - return this->chemistryPtr_->tc(); -} - - template void Foam::combustionModels::PaSR::correct() { if (this->active()) { - const scalar dt = this->mesh().time().deltaTValue(); - - if (integrateReactionRate_) - { - this->chemistryPtr_->solve(dt); - } - else - { - this->chemistryPtr_->calculate(); - } + laminar::correct(); if (turbulentReaction_) { @@ -108,7 +79,7 @@ void Foam::combustionModels::PaSR::correct() tmp tmuEff(this->turbulence().muEff()); const volScalarField& muEff = tmuEff(); - tmp ttc(tc()); + tmp ttc(this->tc()); const volScalarField& tc = ttc(); forAll(epsilon, i) @@ -138,18 +109,7 @@ template Foam::tmp Foam::combustionModels::PaSR::R(volScalarField& Y) const { - tmp tSu(new fvScalarMatrix(Y, dimMass/dimTime)); - - fvScalarMatrix& Su = tSu(); - - if (this->active()) - { - const label specieI = this->thermo().composition().species()[Y.name()]; - - Su += kappa_*this->chemistryPtr_->RR(specieI); - } - - return tSu; + return kappa_*laminar::R(Y); } @@ -157,32 +117,7 @@ template Foam::tmp Foam::combustionModels::PaSR::dQ() const { - tmp tdQ - ( - new volScalarField - ( - IOobject - ( - typeName + ":dQ", - this->mesh().time().timeName(), - this->mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - this->mesh(), - dimensionedScalar("dQ", dimEnergy/dimTime, 0.0), - zeroGradientFvPatchScalarField::typeName - ) - ); - - if (this->active()) - { - volScalarField& dQ = tdQ(); - dQ = kappa_*this->chemistryPtr_->dQ(); - } - - return tdQ; + return kappa_*laminar::dQ(); } @@ -190,44 +125,17 @@ template Foam::tmp Foam::combustionModels::PaSR::Sh() const { - tmp tSh - ( - new volScalarField - ( - IOobject - ( - typeName + ":Sh", - this->mesh().time().timeName(), - this->mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - this->mesh(), - dimensionedScalar("zero", dimEnergy/dimTime/dimVolume, 0.0), - zeroGradientFvPatchScalarField::typeName - ) - ); - - if (this->active()) - { - scalarField& Sh = tSh(); - Sh = kappa_*this->chemistryPtr_->Sh(); - } - - return tSh; + return kappa_*laminar::Sh(); } template bool Foam::combustionModels::PaSR::read() { - if (Type::read()) + if (laminar::read()) { this->coeffs().lookup("Cmix") >> Cmix_; this->coeffs().lookup("turbulentReaction") >> turbulentReaction_; - this->coeffs().lookup("integrateReactionRate") - >> integrateReactionRate_; return true; } else diff --git a/src/combustionModels/PaSR/PaSR.H b/src/combustionModels/PaSR/PaSR.H index 5b4c8d6081..6d5a1ffd92 100644 --- a/src/combustionModels/PaSR/PaSR.H +++ b/src/combustionModels/PaSR/PaSR.H @@ -38,6 +38,8 @@ SourceFiles #ifndef PaSR_H #define PaSR_H +#include "laminar.H" + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // namespace Foam @@ -52,7 +54,7 @@ namespace combustionModels template class PaSR : - public Type + public laminar { // Private data @@ -65,16 +67,9 @@ class PaSR //- Mixing parameter volScalarField kappa_; - //- Integrate reaction rate over the time-step - // using the selected ODE solver - bool integrateReactionRate_; - // Private Member Functions - //- Return the chemical time scale - tmp tc() const; - //- Disallow copy construct PaSR(const PaSR&); diff --git a/src/combustionModels/laminar/laminar.C b/src/combustionModels/laminar/laminar.C index b7a40921e9..a6c712517e 100644 --- a/src/combustionModels/laminar/laminar.C +++ b/src/combustionModels/laminar/laminar.C @@ -25,6 +25,7 @@ License #include "laminar.H" #include "fvmSup.H" +#include "localEulerDdtScheme.H" // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // @@ -76,7 +77,37 @@ void Foam::combustionModels::laminar::correct() { if (integrateReactionRate_) { - this->chemistryPtr_->solve(this->mesh().time().deltaTValue()); + word ddtScheme(this->mesh().ddtScheme("Yi")); + + if (ddtScheme == fv::localEulerDdtScheme::typeName) + { + const scalarField& rDeltaT = + this->mesh().objectRegistry::lookupObject + ( + "rDeltaT" + ); + + if (this->coeffs().found("maxIntegrationTime")) + { + scalar maxIntegrationTime + ( + readScalar(this->coeffs().lookup("maxIntegrationTime")) + ); + + this->chemistryPtr_->solve + ( + min(1.0/rDeltaT, maxIntegrationTime)() + ); + } + else + { + this->chemistryPtr_->solve((1.0/rDeltaT)()); + } + } + else + { + this->chemistryPtr_->solve(this->mesh().time().deltaTValue()); + } } else { diff --git a/src/combustionModels/laminar/laminar.H b/src/combustionModels/laminar/laminar.H index 25e086d81c..b756843394 100644 --- a/src/combustionModels/laminar/laminar.H +++ b/src/combustionModels/laminar/laminar.H @@ -57,12 +57,17 @@ class laminar // using the selected ODE solver bool integrateReactionRate_; +protected: - // Private Member Functions + // Protected Member Functions //- Return the chemical time scale tmp tc() const; +private: + + // Private Member Functions + //- Disallow copy construct laminar(const laminar&); diff --git a/src/thermophysicalModels/chemistryModel/chemistryModel/basicChemistryModel/basicChemistryModel.H b/src/thermophysicalModels/chemistryModel/chemistryModel/basicChemistryModel/basicChemistryModel.H index f7143df63f..55b94ccd5d 100644 --- a/src/thermophysicalModels/chemistryModel/chemistryModel/basicChemistryModel/basicChemistryModel.H +++ b/src/thermophysicalModels/chemistryModel/chemistryModel/basicChemistryModel/basicChemistryModel.H @@ -156,6 +156,10 @@ public: // and return the characteristic time virtual scalar solve(const scalar deltaT) = 0; + //- Solve the reaction system for the given time step + // and return the characteristic time + virtual scalar solve(const scalarField& deltaT) = 0; + //- Return the chemical time scale virtual tmp tc() const = 0; diff --git a/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.C b/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.C index 098083b6fd..1942fb384f 100644 --- a/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.C +++ b/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.C @@ -25,6 +25,7 @@ License #include "chemistryModel.H" #include "reactingMixture.H" +#include "UniformField.H" // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // @@ -742,9 +743,10 @@ void Foam::chemistryModel::calculate() template +template Foam::scalar Foam::chemistryModel::solve ( - const scalar deltaT + const DeltaTType& deltaT ) { CompType::correct(); @@ -795,9 +797,9 @@ Foam::scalar Foam::chemistryModel::solve // initialise timing parameters scalar t = 0; + scalar timeLeft = deltaT[celli]; scalar tauC = this->deltaTChem_[celli]; - scalar dt = min(deltaT, tauC); - scalar timeLeft = deltaT; + scalar dt = min(timeLeft, tauC); // calculate the chemical source terms while (timeLeft > SMALL) @@ -823,17 +825,39 @@ Foam::scalar Foam::chemistryModel::solve dc = c - c0; for (label i=0; i +Foam::scalar Foam::chemistryModel::solve +( + const scalar deltaT +) +{ + // Don't allow the time-step to change more than a factor of 2 + return min + ( + this->solve >(UniformField(deltaT)), + 2*deltaT + ); +} + + +template +Foam::scalar Foam::chemistryModel::solve +( + const scalarField& deltaT +) +{ + return this->solve(deltaT); +} + + template Foam::scalar Foam::chemistryModel::solve ( diff --git a/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.H b/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.H index 232855d620..dcb092c335 100644 --- a/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.H +++ b/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.H @@ -70,6 +70,11 @@ class chemistryModel //- Disallow default bitwise assignment void operator=(const chemistryModel&); + //- Solve the reaction system for the given time step + // of given type and return the characteristic time + template + scalar solve(const DeltaTType& deltaT); + protected: @@ -221,6 +226,10 @@ public: // and return the characteristic time virtual scalar solve(const scalar deltaT); + //- Solve the reaction system for the given time step + // and return the characteristic time + virtual scalar solve(const scalarField& deltaT); + //- Return the chemical time scale virtual tmp tc() const;