diff --git a/applications/solvers/combustion/fireFoam/createPyrolysisModel.H b/applications/solvers/combustion/fireFoam/createPyrolysisModel.H index e5b8d4e457..4c93ad2c6a 100644 --- a/applications/solvers/combustion/fireFoam/createPyrolysisModel.H +++ b/applications/solvers/combustion/fireFoam/createPyrolysisModel.H @@ -1,6 +1,3 @@ Info<< "Creating pyrolysis model" << endl; -autoPtr pyrolysis -( - regionModels::pyrolysisModels::pyrolysisModel::New(mesh) -); +regionModels::pyrolysisModels::pyrolysisModelCollection pyrolysis(mesh); diff --git a/applications/solvers/combustion/fireFoam/fireFoam.C b/applications/solvers/combustion/fireFoam/fireFoam.C index 37ae08ec86..e95a4318c7 100644 --- a/applications/solvers/combustion/fireFoam/fireFoam.C +++ b/applications/solvers/combustion/fireFoam/fireFoam.C @@ -35,7 +35,7 @@ Description #include "turbulenceModel.H" #include "basicReactingCloud.H" #include "surfaceFilmModel.H" -#include "pyrolysisModel.H" +#include "pyrolysisModelCollection.H" #include "radiationModel.H" #include "SLGThermo.H" #include "solidChemistryModel.H" @@ -84,7 +84,7 @@ int main(int argc, char *argv[]) surfaceFilm.evolve(); - pyrolysis->evolve(); + pyrolysis.evolve(); if (solvePrimaryRegion) { diff --git a/applications/solvers/combustion/fireFoam/readPyrolysisTimeControls.H b/applications/solvers/combustion/fireFoam/readPyrolysisTimeControls.H index 95dc78ee5d..c8ba13c886 100644 --- a/applications/solvers/combustion/fireFoam/readPyrolysisTimeControls.H +++ b/applications/solvers/combustion/fireFoam/readPyrolysisTimeControls.H @@ -29,6 +29,6 @@ Description \*---------------------------------------------------------------------------*/ -scalar maxDi = pyrolysis->maxDiff(); +scalar maxDi = pyrolysis.maxDiff(); // ************************************************************************* // diff --git a/applications/solvers/combustion/fireFoam/solidRegionDiffusionNo.H b/applications/solvers/combustion/fireFoam/solidRegionDiffusionNo.H index 2d443c50bf..4aac1f390e 100644 --- a/applications/solvers/combustion/fireFoam/solidRegionDiffusionNo.H +++ b/applications/solvers/combustion/fireFoam/solidRegionDiffusionNo.H @@ -1 +1 @@ -scalar DiNum = pyrolysis->solidRegionDiffNo(); +scalar DiNum = pyrolysis.solidRegionDiffNo(); diff --git a/src/regionModels/pyrolysisModels/Make/files b/src/regionModels/pyrolysisModels/Make/files index 6badede5bc..0a9f942f78 100644 --- a/src/regionModels/pyrolysisModels/Make/files +++ b/src/regionModels/pyrolysisModels/Make/files @@ -3,5 +3,6 @@ pyrolysisModel/pyrolysisModel.C pyrolysisModel/pyrolysisModelNew.C reactingOneDim/reactingOneDim.C noPyrolysis/noPyrolysis.C +pyrolysisModel/pyrolysisModelCollection.C LIB = $(FOAM_LIBBIN)/libpyrolysisModels diff --git a/src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.C b/src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.C index 607b757d02..fe55a6bc27 100644 --- a/src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.C +++ b/src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.C @@ -40,6 +40,8 @@ namespace pyrolysisModels defineTypeNameAndDebug(noPyrolysis, 0); addToRunTimeSelectionTable(pyrolysisModel, noPyrolysis, mesh); +addToRunTimeSelectionTable(pyrolysisModel, noPyrolysis, dictionary); + // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // @@ -57,11 +59,39 @@ bool noPyrolysis::read() } +bool noPyrolysis::read(const dictionary& dict) +{ + if (pyrolysisModel::read(dict)) + { + // no additional info to read + return true; + } + else + { + return false; + } +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // noPyrolysis::noPyrolysis(const word& modelType, const fvMesh& mesh) : - pyrolysisModel(mesh) + pyrolysisModel(modelType, mesh), + solidChemistry_(solidChemistryModel::New(regionMesh())), + solidThermo_(solidChemistry_->solidThermo()) +{} + + +noPyrolysis::noPyrolysis +( + const word& modelType, + const fvMesh& mesh, + const dictionary& dict +): + pyrolysisModel(modelType, mesh, dict), + solidChemistry_(solidChemistryModel::New(regionMesh())), + solidThermo_(solidChemistry_->solidThermo()) {} @@ -73,60 +103,44 @@ noPyrolysis::~noPyrolysis() // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // +void noPyrolysis::preEvolveRegion() +{ + //Do nothing +} + + +void noPyrolysis::evolveRegion() +{ + //Do nothing +} + const volScalarField& noPyrolysis::rho() const { - FatalErrorIn("const volScalarField& noPyrolysis::rho() const") - << "rho field not available for " << type() << abort(FatalError); - return volScalarField::null(); + return (solidThermo_.rho()); } const volScalarField& noPyrolysis::T() const { - FatalErrorIn("const volScalarField& noPyrolysis::T() const") - << "T field not available for " << type() << abort(FatalError); - return volScalarField::null(); + return (solidThermo_.T()); } const tmp noPyrolysis::Cp() const { - FatalErrorIn("const tmp& noPyrolysis::Cp() const") - << "Cp field not available for " << type() << abort(FatalError); - - return tmp - ( - new volScalarField - ( - IOobject - ( - "noPyrolysis::Cp", - time().timeName(), - primaryMesh(), - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - primaryMesh(), - dimensionedScalar("zero", dimEnergy/dimVolume/dimTime, 0.0) - ) - ); + return (solidThermo_.Cp()); } const volScalarField& noPyrolysis::kappa() const { - FatalErrorIn("const volScalarField& noPyrolysis::kappa() const") - << "kappa field not available for " << type() << abort(FatalError); - return volScalarField::null(); + return (solidThermo_.kappa()); } const volScalarField& noPyrolysis::K() const { - FatalErrorIn("const volScalarField& noPyrolysis::K() const") - << "K field not available for " << type() << abort(FatalError); - return volScalarField::null(); + return (solidThermo_.K()); } @@ -140,7 +154,7 @@ const surfaceScalarField& noPyrolysis::phiGas() const // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -} // End namespace surfaceFilmModels +} // End namespace pyrolysisFilmModels } // End namespace regionModels } // End namespace Foam diff --git a/src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.H b/src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.H index b4f34299ee..9179ad1977 100644 --- a/src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.H +++ b/src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.H @@ -73,6 +73,15 @@ protected: //- Read control parameters from dictionary virtual bool read(); + //- Read control parameters from dictionary + virtual bool read(const dictionary& dict); + + //- Reference to the solid chemistry model + autoPtr solidChemistry_; + + //- Reference to solid thermo + basicSolidThermo& solidThermo_; + public: @@ -85,6 +94,15 @@ public: //- Construct from type name and mesh noPyrolysis(const word& modelType, const fvMesh& mesh); + //- Construct from type name and mesh and dict + noPyrolysis + ( + const word& modelType, + const fvMesh& mesh, + const dictionary& dict + ); + + //- Destructor virtual ~noPyrolysis(); @@ -111,6 +129,16 @@ public: //- Return the total gas mass flux to primary region [kg/m2/s] virtual const surfaceScalarField& phiGas() const; + + + // Evolution + + //- Pre-evolve region + virtual void preEvolveRegion(); + + //- Evolve the pyrolysis equations + virtual void evolveRegion(); + }; diff --git a/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModel.C b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModel.C index 05c684c043..cccaf1edf0 100644 --- a/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModel.C +++ b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModel.C @@ -40,6 +40,7 @@ namespace pyrolysisModels defineTypeNameAndDebug(pyrolysisModel, 0); defineRunTimeSelectionTable(pyrolysisModel, mesh); +defineRunTimeSelectionTable(pyrolysisModel, dictionary); // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * // @@ -91,14 +92,33 @@ void pyrolysisModel::constructMeshObjects() // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // +void pyrolysisModel::readPyrolysisControls() +{ + filmCoupled_ = readBool(coeffs_.lookup("filmCoupled")); + reactionDeltaMin_ = + coeffs_.lookupOrDefault("reactionDeltaMin", 0.0); +} + + bool pyrolysisModel::read() { if (regionModel1D::read()) { - filmCoupled_ = readBool(coeffs_.lookup("filmCoupled")); - reactionDeltaMin_ = - coeffs_.lookupOrDefault("reactionDeltaMin", 0.0); + readPyrolysisControls(); + return true; + } + else + { + return false; + } +} + +bool pyrolysisModel::read(const dictionary& dict) +{ + if (regionModel1D::read(dict)) + { + readPyrolysisControls(); return true; } else @@ -134,6 +154,26 @@ pyrolysisModel::pyrolysisModel(const word& modelType, const fvMesh& mesh) } +pyrolysisModel::pyrolysisModel +( + const word& modelType, + const fvMesh& mesh, + const dictionary& dict +) +: + regionModel1D(mesh, "pyrolysis", modelType, dict), + filmCoupled_(false), + filmDeltaPtr_(NULL), + reactionDeltaMin_(0.0) +{ + if (active_) + { + read(dict); + constructMeshObjects(); + } +} + + // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // pyrolysisModel::~pyrolysisModel() @@ -175,7 +215,7 @@ scalar pyrolysisModel::maxDiff() const // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -} // End namespace surfaceFilmModels +} // End namespace pyrolysisModels } // End namespace regionModels } // End namespace Foam diff --git a/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModel.H b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModel.H index 11a7c4ef67..90ca858763 100644 --- a/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModel.H +++ b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModel.H @@ -70,6 +70,9 @@ private: //- Construct fields void constructMeshObjects(); + //- Read pyrolysis controls + void readPyrolysisControls(); + //- Disallow default bitwise copy construct pyrolysisModel(const pyrolysisModel&); @@ -93,9 +96,12 @@ protected: // Protected Member Functions - //- Read control parameters from dictionary + //- Read control parameters virtual bool read(); + //- Read control parameters from dictionary + virtual bool read(const dictionary& dict); + public: @@ -105,17 +111,31 @@ public: // Declare runtime constructor selection table - declareRunTimeSelectionTable - ( - autoPtr, - pyrolysisModel, - mesh, - ( - const word& modelType, - const fvMesh& mesh - ), - (modelType, mesh) - ); + declareRunTimeSelectionTable + ( + autoPtr, + pyrolysisModel, + mesh, + ( + const word& modelType, + const fvMesh& mesh + ), + (modelType, mesh) + ); + + declareRunTimeSelectionTable + ( + autoPtr, + pyrolysisModel, + dictionary, + ( + const word& modelType, + const fvMesh& mesh, + const dictionary& dict + ), + (modelType, mesh, dict) + ); + // Constructors @@ -125,12 +145,58 @@ public: //- Construct from type name and mesh pyrolysisModel(const word& modelType, const fvMesh& mesh); + //- Construct from type name and mesh and dictionary + pyrolysisModel + ( + const word& modelType, + const fvMesh& mesh, + const dictionary& dict + ); + + //- Return clone + autoPtr clone() const + { + notImplemented("autoPtr clone() const"); + return autoPtr(NULL); + } + // Selectors - //- Return a reference to the selected surface film model + //- Return a reference to the selected pyrolysis film model static autoPtr New(const fvMesh& mesh); + //- Return a reference to a named selected pyrolysis model + static autoPtr New + ( + const fvMesh& mesh, + const dictionary& dict + ); + + //- Return pointer to new pyrolysis created on freestore from Istream + class iNew + { + const fvMesh& mesh_; + + public: + + iNew(const fvMesh& mesh) + : + mesh_(mesh) + {} + + autoPtr operator()(Istream& is) const + { + keyType key(is); + dictionary dict(is); + + return autoPtr + ( + pyrolysisModel::New(mesh_, dict) + ); + } + }; + //- Destructor virtual ~pyrolysisModel(); diff --git a/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelCollection.C b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelCollection.C new file mode 100644 index 0000000000..615aa3d84b --- /dev/null +++ b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelCollection.C @@ -0,0 +1,180 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "pyrolysisModelCollection.H" +#include "volFields.H" + +namespace Foam +{ +namespace regionModels +{ +namespace pyrolysisModels +{ + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +defineTemplateTypeNameAndDebug(IOPtrList, 0); + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +pyrolysisModelCollection::pyrolysisModelCollection +( + const fvMesh& mesh +) +: + IOPtrList + ( + IOobject + ( + "pyrolysisZones", + mesh.time().constant(), + mesh, + IOobject::MUST_READ, + IOobject::NO_WRITE + ), + pyrolysisModel::iNew(mesh) + ), + mesh_(mesh) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + + +void pyrolysisModelCollection::preEvolveRegion() +{ + forAll(*this, i) + { + this->operator[](i).preEvolveRegion(); + } +} + + +void pyrolysisModelCollection::evolveRegion() +{ + forAll(*this, i) + { + this->operator[](i).evolveRegion(); + } +} + + +void pyrolysisModelCollection::evolve() +{ + forAll(*this, i) + { + pyrolysisModel& pyrolysis = this->operator[](i); + + if (pyrolysis.active()) + { + if (pyrolysis.primaryMesh().changing()) + { + FatalErrorIn("pyrolysisModelCollection::evolve()") + << "Currently not possible to apply " + << pyrolysis.modelName() + << " model to moving mesh cases" << nl<< abort(FatalError); + } + + // Pre-evolve + pyrolysis.preEvolveRegion(); + + // Increment the region equations up to the new time level + pyrolysis.evolveRegion(); + + // Provide some feedback + if (pyrolysis.infoOutput()) + { + Info<< incrIndent; + pyrolysis.info(); + Info<< endl << decrIndent; + } + } + } +} + + +void pyrolysisModelCollection::info() const +{ + forAll(*this, i) + { + this->operator[](i).info(); + } +} + + +scalar pyrolysisModelCollection::maxDiff() const +{ + scalar maxDiff = 0.0; + forAll(*this, i) + { + if (maxDiff < this->operator[](i).maxDiff()) + { + maxDiff = this->operator[](i).maxDiff(); + } + + } + return maxDiff; +} + + +scalar pyrolysisModelCollection::solidRegionDiffNo() const +{ + scalar regionDiNum = 0.0; + scalar totalDiNum = 0.0; + + forAll(*this, i) + { + const pyrolysisModel& pyrolysis = this->operator[](i); + + if (pyrolysis.regionMesh().nInternalFaces() > 0) + { + surfaceScalarField KrhoCpbyDelta + ( + pyrolysis.regionMesh().surfaceInterpolation::deltaCoeffs() + * fvc::interpolate(pyrolysis.K()) + / fvc::interpolate(pyrolysis.Cp()*pyrolysis.rho()) + ); + + regionDiNum = + max(KrhoCpbyDelta.internalField()) + *pyrolysis.time().deltaTValue(); + } + + if (regionDiNum > totalDiNum) + { + totalDiNum = regionDiNum; + } + } + + return totalDiNum; +} + + +} // End namespace pyrolysisModels +} // End namespace regionModels +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelCollection.H b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelCollection.H new file mode 100644 index 0000000000..84e28ef473 --- /dev/null +++ b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelCollection.H @@ -0,0 +1,122 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Class + Foam::pyrolysisModelCollection + +Description + A centralized pyrolysis collection. + + Container class for a set of pyrolysis with functions implemented + to loop over the functions for each type. + +SourceFiles + pyrolysisModelCollection.C + +\*---------------------------------------------------------------------------*/ + +#ifndef pyrolysisModelCollection_H +#define pyrolysisModelCollection_H + +#include "IOPtrList.H" +#include "pyrolysisModel.H" +#include "fvc.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// Forward declaration of friend functions and operators +class fvMesh; + +namespace regionModels +{ +namespace pyrolysisModels +{ + +/*---------------------------------------------------------------------------*\ + Class pyrolysisModelCollection Declaration +\*---------------------------------------------------------------------------*/ + +class pyrolysisModelCollection +: + public IOPtrList +{ + // Private data + + //- Reference to the finite volume mesh this zone is part of + const fvMesh& mesh_; + + + // Private Member Functions + + //- Disallow default bitwise copy construct + pyrolysisModelCollection(const pyrolysisModelCollection&); + + //- Disallow default bitwise assignment + void operator=(const pyrolysisModelCollection&); + + +public: + + + // Constructors + + //- Construct from mesh + pyrolysisModelCollection(const fvMesh&); + + + // Member Functions + + //- Pre-evolve regions + void preEvolveRegion(); + + //- Evolve the pyrolysis equation regions + void evolveRegion(); + + //- Evolve regions + void evolve(); + + //- Provide some feedback from pyrolysis regions + void info() const; + + //- Return max diffusivity allowed in the solid + scalar maxDiff() const; + + //- Mean diffusion number of the solid regions + scalar solidRegionDiffNo() const; + +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace pyrolysisModels +} // End namespace regionModels +} // End namespace Foam + + +#endif + +// ************************************************************************* // diff --git a/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelNew.C b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelNew.C index 8daa194082..98861a0f5d 100644 --- a/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelNew.C +++ b/src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelNew.C @@ -75,6 +75,33 @@ autoPtr pyrolysisModel::New(const fvMesh& mesh) } +autoPtr pyrolysisModel::New +( + const fvMesh& mesh, + const dictionary& dict +) +{ + + const word modelType = dict.lookup("pyrolysisModel"); + + Info<< "Selecting pyrolysisModel " << modelType << endl; + + dictionaryConstructorTable::iterator cstrIter = + dictionaryConstructorTablePtr_->find(modelType); + + if (cstrIter == dictionaryConstructorTablePtr_->end()) + { + FatalErrorIn("pyrolysisModel::New(const fvMesh&, const dictionary&)") + << "Unknown pyrolysisModel type " << modelType + << nl << nl << "Valid pyrolisisModel types are:" << nl + << dictionaryConstructorTablePtr_->sortedToc() + << exit(FatalError); + } + + return autoPtr(cstrIter()(modelType, mesh, dict)); +} + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace surfaceFilmModels diff --git a/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.C b/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.C index 942094d838..1ce85f10b1 100644 --- a/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.C +++ b/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.C @@ -48,19 +48,40 @@ namespace pyrolysisModels defineTypeNameAndDebug(reactingOneDim, 0); addToRunTimeSelectionTable(pyrolysisModel, reactingOneDim, mesh); +addToRunTimeSelectionTable(pyrolysisModel, reactingOneDim, dictionary); // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // +void reactingOneDim::readReactingOneDimControls() +{ + const dictionary& solution = this->solution().subDict("SIMPLE"); + solution.lookup("nNonOrthCorr") >> nNonOrthCorr_; + time_.controlDict().lookup("maxDi") >> maxDiff_; + + coeffs().lookup("radFluxName") >> primaryRadFluxName_; + coeffs().lookup("minimumDelta") >> minimumDelta_; +} + + bool reactingOneDim::read() { if (pyrolysisModel::read()) { - const dictionary& solution = this->solution().subDict("SIMPLE"); - solution.lookup("nNonOrthCorr") >> nNonOrthCorr_; - time_.controlDict().lookup("maxDi") >> maxDiff_; + readReactingOneDimControls(); + return true; + } + else + { + return false; + } +} - coeffs().lookup("radFluxName") >> primaryRadFluxName_; - coeffs().lookup("minimumDelta") >> minimumDelta_; + +bool reactingOneDim::read(const dictionary& dict) +{ + if (pyrolysisModel::read(dict)) + { + readReactingOneDimControls(); return true; } else @@ -420,6 +441,108 @@ reactingOneDim::reactingOneDim(const word& modelType, const fvMesh& mesh) } +reactingOneDim::reactingOneDim +( + const word& modelType, + const fvMesh& mesh, + const dictionary& dict +) +: + pyrolysisModel(modelType, mesh, dict), + solidChemistry_(solidChemistryModel::New(regionMesh())), + solidThermo_(solidChemistry_->solidThermo()), + kappa_(solidThermo_.kappa()), + K_(solidThermo_.K()), + rho_(solidThermo_.rho()), + Ys_(solidThermo_.composition().Y()), + T_(solidThermo_.T()), + primaryRadFluxName_(dict.lookupOrDefault("radFluxName", "Qr")), + nNonOrthCorr_(-1), + maxDiff_(10), + minimumDelta_(1e-4), + + phiGas_ + ( + IOobject + ( + "phiGas", + time().timeName(), + regionMesh(), + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + regionMesh(), + dimensionedScalar("zero", dimMass/dimTime, 0.0) + ), + + phiHsGas_ + ( + IOobject + ( + "phiHsGas", + time().timeName(), + regionMesh(), + IOobject::NO_READ, + IOobject::NO_WRITE + ), + regionMesh(), + dimensionedScalar("zero", dimEnergy/dimTime, 0.0) + ), + + chemistrySh_ + ( + IOobject + ( + "chemistrySh", + time().timeName(), + regionMesh(), + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + regionMesh(), + dimensionedScalar("zero", dimEnergy/dimTime/dimVolume, 0.0) + ), + + QrCoupled_ + ( + IOobject + ( + primaryRadFluxName_, + time().timeName(), + regionMesh(), + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + regionMesh() + ), + + Qr_ + ( + IOobject + ( + "QrPyr", + time().timeName(), + regionMesh(), + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + regionMesh(), + dimensionedScalar("zero", dimEnergy/dimArea/dimTime, 0.0), + zeroGradientFvPatchVectorField::typeName + ), + + lostSolidMass_(dimensionedScalar("zero", dimMass, 0.0)), + addedGasMass_(dimensionedScalar("zero", dimMass, 0.0)), + totalGasMassFlux_(0.0), + totalHeatRR_(dimensionedScalar("zero", dimEnergy/dimTime, 0.0)) +{ + if (active_) + { + read(dict); + } +} + + // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // reactingOneDim::~reactingOneDim() @@ -560,6 +683,8 @@ void reactingOneDim::preEvolveRegion() void reactingOneDim::evolveRegion() { + Info<< "\nEvolving pyrolysis in region: " << regionMesh().name() << endl; + const scalarField mass0 = rho_*regionMesh().V(); solidChemistry_->solve @@ -591,7 +716,7 @@ void reactingOneDim::evolveRegion() void reactingOneDim::info() const { - Info<< "\nPyrolysis: " << type() << endl; + Info<< "\nPyrolysis in region: " << regionMesh().name() << endl; Info<< indent << "Total gas mass produced [kg] = " << addedGasMass_.value() << nl diff --git a/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.H b/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.H index 1b47cd15c4..ed1bd9fa0d 100644 --- a/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.H +++ b/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.H @@ -66,6 +66,9 @@ private: //- Disallow default bitwise assignment void operator=(const reactingOneDim&); + //- Read model controls + void readReactingOneDimControls(); + protected: @@ -154,6 +157,9 @@ protected: //- Read control parameters from dictionary bool read(); + //- Read control parameters from dict + bool read(const dictionary& dict); + //- Update submodels void updateFields(); @@ -193,6 +199,16 @@ public: //- Construct from type name and mesh reactingOneDim(const word& modelType, const fvMesh& mesh); + //- Construct from type name, mesh and dictionary + reactingOneDim + ( + const word& modelType, + const fvMesh& mesh, + const dictionary& dict + ); + + + //- Destructor virtual ~reactingOneDim(); diff --git a/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/pyrolysisProperties b/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/pyrolysisZones similarity index 66% rename from tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/pyrolysisProperties rename to tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/pyrolysisZones index ed2bd1c315..6253de1741 100644 --- a/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/pyrolysisProperties +++ b/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/pyrolysisZones @@ -11,27 +11,35 @@ FoamFile format binary; class dictionary; location "constant"; - object pyrolysisProperties; + object pyrolysisZones; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -active true; +1 +( + pyrolysis + { + active true; -pyrolysisModel reactingOneDim; + pyrolysisModel reactingOneDim; -regionName panelRegion; + regionName panelRegion; -reactingOneDimCoeffs -{ - filmCoupled false; + reactingOneDimCoeffs + { + filmCoupled false; - radFluxName Qr; + radFluxName Qr; - minimumDelta 1e-8; + minimumDelta 1e-8; - moveMesh false; -} + reactionDeltaMin 1e-8; -infoOutput true; + moveMesh false; + } + + infoOutput true; + } +) // ************************************************************************* //