diff --git a/applications/solvers/combustion/SLFMFoam/rhoEqn.H b/applications/solvers/combustion/SLFMFoam/AMC/AMC.C similarity index 71% rename from applications/solvers/combustion/SLFMFoam/rhoEqn.H rename to applications/solvers/combustion/SLFMFoam/AMC/AMC.C index b69701e74..f76b9563d 100644 --- a/applications/solvers/combustion/SLFMFoam/rhoEqn.H +++ b/applications/solvers/combustion/SLFMFoam/AMC/AMC.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2016 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -21,30 +21,29 @@ License You should have received a copy of the GNU General Public License along with OpenFOAM. If not, see . -Global - rhoEqn - -Description - Solve the continuity for density. - \*---------------------------------------------------------------------------*/ +#include "AMC.H" + +#include "error.H" + +// * * * * * * * * * * * * * Functions * * * * * * * * * * * * // + +//Amplitude Mapping Closure (from KIVA) +Foam::tmp Foam::CMC::AMC(const scalarField& eta) { - fvScalarMatrix rhoEqn - ( - fvm::ddt(rho) - + fvc::div(phi) - == - parcels.Srho(rho) - + fvOptions(rho) - ); + tmp tRes(new scalarField(eta.size(), 0.0)); - rhoEqn.solve(); + scalarField& Res = tRes.ref(); - fvOptions.correct(rho); + for(label i=0 ; i. + +Function + Foam::CMC::AMC + +Description + +SourceFiles + AMC.C + +\*---------------------------------------------------------------------------*/ + +#ifndef AMC_H +#define AMC_H + +#include "tmp.H" +#include "scalarField.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace CMC +{ + +/*---------------------------------------------------------------------------*\ + Function AMC Declaration +\*---------------------------------------------------------------------------*/ + + const scalar pi = 3.141592; + + const scalar spi = Foam::sqrt(pi); + + //Amplitude Mapping Closure (from KIVA) + tmp AMC(const scalarField& eta); + + + //Amplitude Mapping Closure (from KIVA) + //Define exp(-2*(erf^-1(2*eta - 1))^2) + inline scalar AMC(const scalar eta); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace CMC +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "AMCI.H" + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/AMC/AMCI.H b/applications/solvers/combustion/SLFMFoam/AMC/AMCI.H new file mode 100644 index 000000000..ad663725b --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/AMC/AMCI.H @@ -0,0 +1,54 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + + +//Amplitude Mapping Closure (from KIVA) +//Define exp(-2*(erf^-1(2*eta - 1))^2) +inline Foam::scalar Foam::CMC::AMC(const scalar eta) +{ + const scalar a0 = (2.0*eta - 1.0); + scalar a1 = 0.5; + scalar slope; + scalar da = GREAT; + + while(mag(da) > SMALL) + { + slope = (2.0/spi)*Foam::exp(-1.0*Foam::pow(a1,2.0)); + da = (a0 - Foam::erf(a1))/slope; + a1 = a1 + da; + } + + return Foam::exp(-2.0*Foam::pow(a1,2.0)); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunction.C b/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunction.C new file mode 100644 index 000000000..b6a13aef6 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunction.C @@ -0,0 +1,120 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 "BetaFunction.H" + +#include "error.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +// const scalar Foam::CMC::BetaFunction::staticData(); + + +// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void Foam::CMC::BetaFunction::_setParameters() +{ + scalar gamma = mf_*(1.0 - mf_)/(mfVar_ + SMALL) - 1.0; + + if (gamma < SMALL || (mfVar_)/(mf_*(1 - mf_) + SMALL) < 0.001) + { + fdelta_ = true; + } + else + { + a_ = max(mf_*gamma, 0.0); + b_ = max((1.0 - mf_)*gamma, 0.0); + + if(a_ < 1.0) + { + delta0_ = true; + } + if(b_ < 1.0) + { + delta1_ = true; + } + + _limitAB(); + } +} + + +void Foam::CMC::BetaFunction::_limitAB() +{ + scalar fmax = 1.0/(1.0 + (b_ - 1.0)/(a_ - 1.0)); + + if(a_ > 500.0) + { + a_ = 500.0; + b_ = (a_ - 1.0 - fmax*(a_ - 2.0))/fmax; + } + else if(b_ > 500.0) + { + b_ = 500.0; + a_ = (1.0 + fmax*(b_ - 2.0))/(1.0 - fmax); + } +} + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::CMC::BetaFunction::BetaFunction(const scalar mf, const scalar mfVar) +: + mf_(mf), + mfVar_(mfVar), + a_(0.0), + b_(0.0), + fdelta_(false), + delta0_(false), + delta1_(false) +{ + _setParameters(); +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::CMC::BetaFunction::~BetaFunction() +{} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * // + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunction.H b/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunction.H new file mode 100644 index 000000000..81c2d04f5 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunction.H @@ -0,0 +1,199 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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::BetaFunction + +Description + +SourceFiles + BetaFunctionI.H + BetaFunction.C + BetaFunctionIO.C + +\*---------------------------------------------------------------------------*/ + +#ifndef BetaFunction_H +#define BetaFunction_H + +#include "scalar.H" +#include "Switch.H" +#include "scalarField.H" +#include "tmp.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// Forward declaration of classes +class Istream; +class Ostream; + +// Forward declaration of friend functions and operators +namespace CMC +{ + class BetaFunction; +} +Istream& operator>>(Istream&, CMC::BetaFunction&); +Ostream& operator<<(Ostream&, const CMC::BetaFunction&); + + +/*---------------------------------------------------------------------------*\ + Class BetaFunction Declaration +\*---------------------------------------------------------------------------*/ + +namespace CMC +{ + +class BetaFunction +{ + // Private data + + //- Mean mixture fraction + scalar mf_; + + //- Mixture fraction variance + scalar mfVar_; + + //- Beta pdf parameter alpha + scalar a_; + + //- Beta pdf parameter beta + scalar b_; + + //- Flag beta pdf special case: forced-delta function + Switch fdelta_; + + //- delta function at pure oxidizer stream + Switch delta0_; + + //- delta function at pure fuel stream + Switch delta1_; + + + // Private Member Functions + + //- Calculate alpha and beta and set special case flags + void _setParameters(); + + //- Limit too large value of beta pdf parameters + void _limitAB(); + + //- Disallow default bitwise copy construct + // BetaFunction(const BetaFunction&); + + //- Disallow default bitwise assignment + // void operator=(const BetaFunction&); + + +public: + + // Static data members + + //- Static data staticData + // static const scalar staticData; + + + // Constructors + + //- Construct from components + BetaFunction(const scalar mf, const scalar mfVar); + + //- Construct from Istream + BetaFunction(Istream&); + + //- Construct as copy + // BetaFunction(const BetaFunction&); + + + //- Destructor + ~BetaFunction(); + + + // Member Functions + + tmp etaFunc(const UList& eta) const + { + return pow(eta, a_-1.0)*pow(1.0-eta, b_-1.0); + } + + scalar etaFunc(const scalar eta) const + { + return Foam::pow(eta, a_-1.0)*Foam::pow(1.0-eta, b_-1.0); + } + + + // Access + + Switch fdelta() const {return fdelta_;} + + Switch delta0() const {return delta0_;} + + Switch delta1() const {return delta1_;} + + scalar mf() const {return mf_;} + + scalar mfVar() const {return mfVar_;} + + scalar alpha() const {return a_;} + + scalar beta() const {return b_;} + + // Check + + // Edit + + // Write + + + // Member Operators + + // void operator=(const BetaFunction&); + + + // Friend Functions + + // Friend Operators + + // IOstream Operators + + friend Istream& ::Foam::operator>>(Istream&, BetaFunction&); + friend Ostream& ::Foam::operator<<(Ostream&, const BetaFunction&); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace CMC +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "BetaFunctionI.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunctionI.H b/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunctionI.H new file mode 100644 index 000000000..eca8d9ea3 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunctionI.H @@ -0,0 +1,58 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunctionIO.C b/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunctionIO.C new file mode 100644 index 000000000..9d568fa18 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaFunction/BetaFunctionIO.C @@ -0,0 +1,65 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 "BetaFunction.H" +#include "IOstreams.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::CMC::BetaFunction::BetaFunction(Istream& is) +{ + // Check state of Istream + is.check("Foam::CMC::BetaFunction::BetaFunction(Foam::Istream&)"); +} + + +// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * // + +Foam::Istream& Foam::operator>>(Istream& is, CMC::BetaFunction&) +{ + // Check state of Istream + is.check + ( + "Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::CMC::BetaFunction&)" + ); + + return is; +} + + +Foam::Ostream& Foam::operator<<(Ostream& os, const CMC::BetaFunction&) +{ + // Check state of Ostream + os.check + ( + "Foam::Ostream& Foam::operator<<(Foam::Ostream&, " + "const Foam::CMC::BetaFunction&)" + ); + + return os; +} + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGrid.C b/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGrid.C new file mode 100644 index 000000000..71b421187 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGrid.C @@ -0,0 +1,126 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 "BetaGrid.H" + +#include "error.H" +#include "dictionary.H" +#include "SubList.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +// const scalar Foam::CMC::BetaGrid::staticData(); + + +// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::CMC::BetaGrid::BetaGrid(const dictionary& dict) +: + etaCut_(dict.lookup("detailedEta")), + nSpacings_(dict.lookup("detailedN")), + etaSpace_(sum(nSpacings_)+3) +{ + // Check for input list size + if (etaCut_.size() != nSpacings_.size() + 1) + { + FatalErrorInFunction + << "Number of grid interval points does not match number of grid sizes" + << abort(FatalError); + } + + // Check for evenness of number of spacings in each interval + forAll (nSpacings_, i) + { + if ((nSpacings_[i] % 2) == 1) + { + FatalErrorInFunction + << "Number of spacings in intervals should be even" + << abort(FatalError); + } + } + + forAll(nSpacings_, i) + { + const label ns = nSpacings_[i]; + const labelList::subList prev(nSpacings_, i); + const label baseI = sum(prev) + 1; + + const scalar low = etaCut_[i]; + const scalar upp = etaCut_[i+1]; + const scalar delta = (upp - low) / scalar(ns); + + for (label j = 0; j < ns; j++) + { + etaSpace_[baseI+j] = low + delta * j; + } + } + etaSpace_[0] = 0.0; + etaSpace_[etaSpace_.size()-2] = etaCut_.last(); + etaSpace_[etaSpace_.size()-1] = 1.0; +} + + +Foam::CMC::BetaGrid::BetaGrid(const BetaGrid&) +{} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::CMC::BetaGrid::~BetaGrid() +{} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * // + +void Foam::CMC::BetaGrid::operator=(const BetaGrid& rhs) +{ + // Check for assignment to self + if (this == &rhs) + { + FatalErrorInFunction + << "Attempted assignment to self" + << abort(FatalError); + } +} + +// * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * // + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGrid.H b/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGrid.H new file mode 100644 index 000000000..fc17e0ab2 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGrid.H @@ -0,0 +1,161 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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::BetaGrid + +Description + +SourceFiles + BetaGridI.H + BetaGrid.C + BetaGridIO.C + +\*---------------------------------------------------------------------------*/ + +#ifndef BetaGrid_H +#define BetaGrid_H + +#include "scalarField.H" +#include "labelList.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// Forward declaration of classes +class Istream; +class Ostream; + +// Forward declaration of friend functions and operators +namespace CMC +{ + class BetaGrid; +} +class dictionary; + +Istream& operator>>(Istream&, CMC::BetaGrid&); +Ostream& operator<<(Ostream&, const CMC::BetaGrid&); + + +/*---------------------------------------------------------------------------*\ + Class BetaGrid Declaration +\*---------------------------------------------------------------------------*/ + +namespace CMC +{ + +class BetaGrid +{ + // Private data + + //- Piecewise uniform interval points + scalarField etaCut_; + + //- Sizes of piecewise uniform intervals + labelList nSpacings_; + + //- Mixture fraction eta grid points + scalarField etaSpace_; + + + // Private Member Functions + + //- Disallow default bitwise copy construct + BetaGrid(const BetaGrid&); + + //- Disallow default bitwise assignment + void operator=(const BetaGrid&); + + +public: + + // Static data members + + //- Static data staticData + // static const scalar staticData; + + + // Constructors + + //- Construct from components + BetaGrid(const dictionary& dict); + + //- Construct from Istream + BetaGrid(Istream&); + + //- Construct as copy + // BetaGrid(const BetaGrid&); + + + //- Destructor + ~BetaGrid(); + + + // Member Functions + + // Access + const labelList& N() const {return nSpacings_;} + + const scalarField& etaCut() const {return etaCut_;} + + const scalarField& etaSpace() const {return etaSpace_;} + + // Check + + // Edit + + // Write + + + // Member Operators + + // void operator=(const BetaGrid&); + + + // Friend Functions + + // Friend Operators + + // IOstream Operators + + friend Istream& ::Foam::operator>>(Istream&, BetaGrid&); + friend Ostream& ::Foam::operator<<(Ostream&, const BetaGrid&); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace CMC +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "BetaGridI.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGridI.H b/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGridI.H new file mode 100644 index 000000000..eca8d9ea3 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGridI.H @@ -0,0 +1,58 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGridIO.C b/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGridIO.C new file mode 100644 index 000000000..d567d459b --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaGrid/BetaGridIO.C @@ -0,0 +1,65 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 "BetaGrid.H" +#include "IOstreams.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::CMC::BetaGrid::BetaGrid(Istream& is) +{ + // Check state of Istream + is.check("Foam::CMC::BetaGrid::BetaGrid(Foam::Istream&)"); +} + + +// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * // + +Foam::Istream& Foam::operator>>(Istream& is, CMC::BetaGrid&) +{ + // Check state of Istream + is.check + ( + "Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::CMC::BetaGrid&)" + ); + + return is; +} + + +Foam::Ostream& Foam::operator<<(Ostream& os, const CMC::BetaGrid&) +{ + // Check state of Ostream + os.check + ( + "Foam::Ostream& Foam::operator<<(Foam::Ostream&, " + "const Foam::CMC::BetaGrid&)" + ); + + return os; +} + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegrator.C b/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegrator.C new file mode 100644 index 000000000..fc16c5ffb --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegrator.C @@ -0,0 +1,180 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 "BetaIntegrator.H" + + +#include "error.H" +#include "SubList.H" +#include "SubField.H" +#include "interpolateXY.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +// const scalar Foam::CMC::BetaIntegrator::staticData(); + + +// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +//extended Simpson's rule (Numerical recipes, 2nd Ed. p.128) +//for equally spaced and even N intervals (or odd N+1 points) +Foam::scalar Foam::CMC::BetaIntegrator::simps +(const scalar xl, const scalar xh, const label n, const UList& fx) + const +{ + scalar evensum(0.0), oddsum(0.0), sum(0.0); + + scalar h = (xh - xl)/scalar(n); + + for(label i=0 ; i& f) const +{ + scalar total = 0.0; + + const labelList &N(bg_.N()); + const scalarField &etaCut(bg_.etaCut()); + + forAll(N, i) + { + const labelList::subList prev(N, i); + const label baseI = sum(prev) + 1; + + const scalarField::subField subInterval(f, N[i]+1, baseI); + + total += simps(etaCut[i], etaCut[i+1], N[i], subInterval); + } + + return total; +} + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::CMC::BetaIntegrator::BetaIntegrator(const BetaFunction& bf, const BetaGrid &bg) +: + etaSpace_(bg.etaSpace()), + pdfNum_(bg.etaSpace().size(), 0.0), + pdfDen_(0.0), + bf_(bf), + bg_(bg) +{ + const scalarField::subField domain(etaSpace_, etaSpace_.size()-2, 1); + scalarField::subField range(pdfNum_, etaSpace_.size()-2, 1); + + if (bf_.fdelta()) + { + pdfDen_ = 1.0; + } + else + { + range = bf_.etaFunc(domain); + + if (!bf_.delta0()) + { + pdfNum_.first() = bf_.etaFunc(etaSpace_.first()); + } + + if (!bf_.delta1()) + { + pdfNum_.last() = bf_.etaFunc(etaSpace_.last()); + } + + pdfDen_ = sumSimps(pdfNum_) + + Foam::pow(etaSpace_[1], bf_.alpha())/(bf_.alpha() + SMALL) + + Foam::pow(etaSpace_[1], bf_.beta())/(bf_.beta() + SMALL); + } +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::CMC::BetaIntegrator::~BetaIntegrator() +{} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + +// beta-pdf weighted integration for given mf, mfVar and f +Foam::scalar Foam::CMC::BetaIntegrator::betaIntegrate(const scalarField& f) const +{ + scalar result = 0.0; + + if(bf_.fdelta()) + { + result = interpolateXY(bf_.mf(), etaSpace_, f); + } + else + { + scalar num = sumSimps (f * pdfNum_) + + f.first()*Foam::pow(etaSpace_[1], bf_.alpha())/(bf_.alpha() + SMALL) + + f.last()*Foam::pow(etaSpace_[1], bf_.beta())/(bf_.beta() + SMALL); + result = num/pdfDen_; + } + + return result; +} + +// * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * // + +void Foam::CMC::BetaIntegrator::operator=(const BetaIntegrator& rhs) +{ + // Check for assignment to self + if (this == &rhs) + { + FatalErrorInFunction + << "Attempted assignment to self" + << abort(FatalError); + } +} + +// * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * // + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegrator.H b/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegrator.H new file mode 100644 index 000000000..692bff582 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegrator.H @@ -0,0 +1,168 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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::CMC::BetaIntegrator + +Description + +SourceFiles + BetaIntegratorI.H + BetaIntegrator.C + BetaIntegratorIO.C + +\*---------------------------------------------------------------------------*/ + +#ifndef BetaIntegrator_H +#define BetaIntegrator_H + +#include "scalarField.H" +#include "BetaFunction.H" +#include "BetaGrid.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// Forward declaration of classes +class Istream; +class Ostream; + +// Forward declaration of friend functions and operators +namespace CMC +{ + class BetaIntegrator; +} +Istream& operator>>(Istream&, CMC::BetaIntegrator&); +Ostream& operator<<(Ostream&, const CMC::BetaIntegrator&); + + +/*---------------------------------------------------------------------------*\ + Class BetaIntegrator Declaration +\*---------------------------------------------------------------------------*/ + +namespace CMC +{ + +class BetaIntegrator +{ + // Private data + + //- Quadrature points + const scalarField &etaSpace_; + + //- Beta pdf numerators + scalarField pdfNum_; + + //- Beta pdf denominator + scalar pdfDen_; + + //- Description of data_ + BetaFunction bf_; + + //- Description of data_ + const BetaGrid &bg_; + + + // Private Member Functions + + //- Calculate integral using Simpson's rule + scalar simps(const scalar xl, const scalar xh, const label n, const UList& fx) const; + + //- Calculate integral using Simpson's rule + scalar sumSimps(const UList& f) const; + + //- Disallow default bitwise copy construct + BetaIntegrator(const BetaIntegrator&); + + //- Disallow default bitwise assignment + void operator=(const BetaIntegrator&); + + +public: + + // Static data members + + //- Static data staticData + // static const scalar staticData; + + + // Constructors + + //- Construct from components + BetaIntegrator(const BetaFunction& bf, const BetaGrid &bg); + + //- Construct as copy + // BetaIntegrator(const BetaIntegrator&); + + + //- Destructor + ~BetaIntegrator(); + + + // Member Functions + + scalar betaIntegrate(const scalarField& f) const; + + // Access + + scalar pdfDen () const {return pdfDen_;} + + // Check + + // Edit + + // Write + + + // Member Operators + + // void operator=(const BetaIntegrator&); + + + // Friend Functions + + // Friend Operators + + // IOstream Operators + + friend Istream& ::Foam::operator>>(Istream&, BetaIntegrator&); + friend Ostream& ::Foam::operator<<(Ostream&, const BetaIntegrator&); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace CMC +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "BetaIntegratorI.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegratorI.H b/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegratorI.H new file mode 100644 index 000000000..eca8d9ea3 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegratorI.H @@ -0,0 +1,58 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegratorIO.C b/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegratorIO.C new file mode 100644 index 000000000..881819730 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/BetaIntegrator/BetaIntegratorIO.C @@ -0,0 +1,56 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 "BetaIntegrator.H" +#include "IOstreams.H" + +// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * // + +Foam::Istream& Foam::operator>>(Istream& is, CMC::BetaIntegrator&) +{ + // Check state of Istream + is.check + ( + "Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::CMC::BetaIntegrator&)" + ); + + return is; +} + + +Foam::Ostream& Foam::operator<<(Ostream& os, const CMC::BetaIntegrator&) +{ + // Check state of Ostream + os.check + ( + "Foam::Ostream& Foam::operator<<(Foam::Ostream&, " + "const Foam::CMC::BetaIntegrator&)" + ); + + return os; +} + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/EEqn.H b/applications/solvers/combustion/SLFMFoam/EEqn.H deleted file mode 100644 index 8239d4462..000000000 --- a/applications/solvers/combustion/SLFMFoam/EEqn.H +++ /dev/null @@ -1,33 +0,0 @@ -{ - volScalarField& he = thermo.he(); - - fvScalarMatrix EEqn - ( - mvConvection->fvmDiv(phi, he) - + ( - he.name() == "e" - ? fvc::div(phi, volScalarField("Ekp", 0.5*magSqr(U) + p/rho)) - : fvc::div(phi, volScalarField("K", 0.5*magSqr(U))) - ) - - fvm::laplacian(turbulence->alphaEff(), he) - == - rho*(U&g) - + parcels.Sh(he) - + radiation->Sh(thermo) - + combustion->Sh() - + fvOptions(rho, he) - ); - - EEqn.relax(); - - fvOptions.constrain(EEqn); - - EEqn.solve(); - - fvOptions.correct(he); - thermo.correct(); - radiation->correct(); - - Info<< "T gas min/max = " << min(T).value() << ", " - << max(T).value() << endl; -} diff --git a/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructure.C b/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructure.C new file mode 100644 index 000000000..3d37bca4c --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructure.C @@ -0,0 +1,206 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 "FlameStructure.H" + +#include "fileName.H" +#include "IFstream.H" +#include "hashedWordList.H" +#include "interpolateXY.H" + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::SLFM::FlameStructure::FlameStructure(const string& NstDir, scalarList NstList, const scalarField &etaGrid, const hashedWordList &spTable) +: + NstDir_(NstDir), + NstList_(NstList), + eta_(etaGrid), + species_(spTable), + Y_(spTable.size(), scalarFieldArray1d(NstList_.size(), scalarField(eta_.size(), 0.0))), + W_(spTable.size(), scalarFieldArray1d(NstList_.size(), scalarField(eta_.size(), 0.0))), + T_(NstList_.size(), scalarField(eta_.size(), 0.0)), + Q_(NstList_.size(), scalarField(eta_.size(), 0.0)), + h_(NstList_.size(), scalarField(eta_.size(), 0.0)), + rho_(NstList_.size(), scalarField(eta_.size(), 0.0)), + Rgas_(NstList_.size(), scalarField(eta_.size(), 0.0)) +{ + for (label n = 0; n < NstList_.size(); n++) + { + scalar N = NstList_[n]; + + fileName fname = "sdr"+Foam::name(N)+".inp"; + IFstream fin(NstDir_/fname); + + string gbg; + label etamax_SLFM, NoSpecies; + + //INPUT FILE... + fin.getLine(gbg); + + //No.GRID... + fin.getLine(gbg); + + //Number of eta-point and species + fin >> etamax_SLFM >> NoSpecies; + + //blank line + fin.getLine(gbg); + + //PRESSURE... + fin.getLine(gbg); + + //1.00... + fin.getLine(gbg); + + //MIXTURE... + fin.getLine(gbg); + + //read eta space (from SLFM library) + scalarField etaValue_SLFM(etamax_SLFM, 0.0); + for(label j = 0 ; j < etamax_SLFM ; j++) + { + fin>>etaValue_SLFM[j]; + } + + //blank line + fin.getLine(gbg); + + //INITIAL... + fin.getLine(gbg); + + //species loop + scalarField temp(etamax_SLFM, 0.0); + for(label i = 0 ; i < NoSpecies ; i++) + { + word spName; + //H2... (species name) + fin >> spName; + const label spI (species_[spName]); + + // Info << spI << spName << endl; + + //read species mass fraction (from SLFM library) + for(label j=0 ; j> Yj; + temp[j] = max(0.0, Yj); + } + + Y_[spI][n] = interpolateXY(etaGrid, etaValue_SLFM, temp); + + //read reaction rate (from SLFM library) + for(label j = 0 ; j < etamax_SLFM ; j++) + { + fin>>temp[j]; + } + + W_[spI][n] = interpolateXY(etaGrid, etaValue_SLFM, temp); + + //blank line + fin.getLine(gbg); + } + + //read temperature + { + //H2... (species name) + fin.getLine(gbg); + + //read temperature + for(label j = 0 ; j < etamax_SLFM ; j++) + { + scalar Yj = 0.0; + fin >> Yj; + temp[j] = max(0.0, Yj); + } + + T_[n] = interpolateXY(etaGrid, etaValue_SLFM, temp); + + //read heat source + for(label j = 0 ; j < etamax_SLFM ; j++) + { + fin>>temp[j]; + } + + Q_[n] = interpolateXY(etaGrid, etaValue_SLFM, temp); + + //blank line + fin.getLine(gbg); + } + + //DENSITY... + fin.getLine(gbg); + + //read density (from SLFM library) + for(label j = 0 ; j < etamax_SLFM ; j++) + { + fin>>temp[j]; + } + //Info<>temp[j]; + } + + h_[n] = interpolateXY(etaGrid, etaValue_SLFM, temp); + } +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::SLFM::FlameStructure::~FlameStructure() +{} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * // + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructure.H b/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructure.H new file mode 100644 index 000000000..96a5ba8d8 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructure.H @@ -0,0 +1,215 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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::SLFM::FlameStructure + +Description + +SourceFiles + FlameStructureI.H + FlameStructure.C + FlameStructureIO.C + +\*---------------------------------------------------------------------------*/ + +#ifndef FlameStructure_H +#define FlameStructure_H + +#include "error.H" +#include "scalarField.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// Forward declaration of classes +class Istream; +class Ostream; + +class hashedWordList; + +namespace SLFM +{ + class FlameStructure; +} + +// Forward declaration of friend functions and operators +Istream& operator>>(Istream&, SLFM::FlameStructure&); +Ostream& operator<<(Ostream&, const SLFM::FlameStructure&); + +namespace SLFM +{ + +/*---------------------------------------------------------------------------*\ + Class FlameStructure Declaration +\*---------------------------------------------------------------------------*/ + +class FlameStructure +{ + // Private data + typedef List scalarFieldArray1d; + + typedef List scalarFieldArray2d; + + typedef List scalarFieldArray3d; + + //- Species mass fraction + word NstDir_; + + //- Species mass fraction + const scalarList NstList_; + + //- Species mass fraction + const scalarField eta_; + + //- Species mass fraction + const hashedWordList &species_; + + //- Species mass fraction + scalarFieldArray2d Y_; + + //- Species production rate + scalarFieldArray2d W_; + + //- Temperature + scalarFieldArray1d T_; + + //- Heat source + scalarFieldArray1d Q_; + + //- Enthalpy + scalarFieldArray1d h_; + + //- Density + scalarFieldArray1d rho_; + + //- Specific gas constant for mean W + scalarFieldArray1d Rgas_; + + + // Private Member Functions + + //- Disallow Construct null + FlameStructure(); + + //- Disallow default bitwise copy construct + // FlameStructure(const FlameStructure&); + + //- Disallow default bitwise assignment + // void operator=(const FlameStructure&); + + +public: + + // Static data members + + //- Static data staticData + // static const dataType staticData; + + + // Constructors + + //- Construct from components + FlameStructure(const string& NstDir, scalarList NstList, const scalarField &etaGrid, const hashedWordList& spTable); + + //- Construct from Istream + // FlameStructure(Istream&); + + //- Construct as copy + // FlameStructure(const FlameStructure&); + + + // Selectors + + //- Select null constructed + // static autoPtr New(); + + + //- Destructor + ~FlameStructure(); + + + // Member Functions + + // Access + + //- Select null constructed + scalarFieldArray2d& Y() {return Y_;} + + //- Select null constructed + scalarFieldArray2d& W() {return W_;} + + //- Select null constructed + scalarFieldArray1d& T() {return T_;} + + //- Select null constructed + scalarFieldArray1d& Q() {return Q_;} + + //- Select null constructed + scalarFieldArray1d& rho() {return rho_;} + + //- Select null constructed + scalarFieldArray1d& h() {return h_;} + + //- Select null constructed + scalarFieldArray1d& Rgas() {return Rgas_;} + + // Check + + // Edit + + // Write + + + // Member Operators + + // void operator=(const FlameStructure&); + + + // Friend Functions + + // Friend Operators + + // IOstream Operators + + friend Istream& operator>>(Istream&, FlameStructure&); + friend Ostream& operator<<(Ostream&, const FlameStructure&); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace SLFM +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "FlameStructureI.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructureI.H b/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructureI.H new file mode 100644 index 000000000..eca8d9ea3 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructureI.H @@ -0,0 +1,58 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * // + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructureIO.C b/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructureIO.C new file mode 100644 index 000000000..75b6748f5 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/FlameStructure/FlameStructureIO.C @@ -0,0 +1,72 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 "FlameStructure.H" +#include "IOstreams.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + + /* +Foam::SLFM::FlameStructure::FlameStructure(Istream& is) +: + base1(is), + base2(is), + member1(is), + member2(is) +{ + // Check state of Istream + is.check("Foam::SLFM::FlameStructure::FlameStructure(Foam::Istream&)"); +} + */ + + +// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * // + +Foam::Istream& Foam::operator>>(Istream& is, SLFM::FlameStructure&) +{ + // Check state of Istream + is.check + ( + "Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::SLFM::FlameStructure&)" + ); + + return is; +} + + +Foam::Ostream& Foam::operator<<(Ostream& os, const SLFM::FlameStructure&) +{ + // Check state of Ostream + os.check + ( + "Foam::Ostream& Foam::operator<<(Foam::Ostream&, " + "const Foam::SLFM::FlameStructure&)" + ); + + return os; +} + + +// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/Make/files b/applications/solvers/combustion/SLFMFoam/Make/files index 0df9202d9..3349e8db0 100644 --- a/applications/solvers/combustion/SLFMFoam/Make/files +++ b/applications/solvers/combustion/SLFMFoam/Make/files @@ -1,3 +1,17 @@ +FlameStructure/FlameStructure.C +FlameStructure/FlameStructureIO.C + +./AMC/AMC.C + +./BetaFunction/BetaFunction.C +./BetaFunction/BetaFunctionIO.C + +./BetaGrid/BetaGrid.C +./BetaGrid/BetaGridIO.C + +./BetaIntegrator/BetaIntegrator.C +./BetaIntegrator/BetaIntegratorIO.C + SLFMFoam.C EXE = $(FOAM_APPBIN)/SLFMFoam diff --git a/applications/solvers/combustion/SLFMFoam/Make/options b/applications/solvers/combustion/SLFMFoam/Make/options index f7f22123e..3e21a2bcd 100644 --- a/applications/solvers/combustion/SLFMFoam/Make/options +++ b/applications/solvers/combustion/SLFMFoam/Make/options @@ -1,28 +1,19 @@ EXE_INC = \ -I. \ + -I./AMC \ + -I./BetaFunction \ + -I./BetaGrid \ + -I./BetaIntegrator \ -I$(LIB_SRC)/finiteVolume/lnInclude \ - -I${LIB_SRC}/meshTools/lnInclude \ + -I$(LIB_SRC)/meshTools/lnInclude \ -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \ -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \ - -I$(LIB_SRC)/lagrangian/basic/lnInclude \ - -I$(LIB_SRC)/lagrangian/intermediate/lnInclude \ - -I$(LIB_SRC)/lagrangian/coalCombustion/lnInclude \ - -I$(LIB_SRC)/lagrangian/distributionModels/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \ -I$(LIB_SRC)/transportModels/compressible/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/properties/liquidProperties/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/properties/liquidMixtureProperties/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/properties/solidProperties/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/properties/solidMixtureProperties/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/thermophysicalFunctions/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \ -I$(LIB_SRC)/ODE/lnInclude \ - -I$(LIB_SRC)/regionModels/regionModel/lnInclude \ - -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \ -I$(LIB_SRC)/combustionModels/lnInclude \ -I$(LIB_SRC)/sampling/lnInclude @@ -32,24 +23,12 @@ EXE_LIBS = \ -lmeshTools \ -lturbulenceModels \ -lcompressibleTurbulenceModels \ - -llagrangian \ - -llagrangianIntermediate \ - -llagrangianTurbulence \ -lspecie \ -lcompressibleTransportModels \ -lfluidThermophysicalModels \ - -lliquidProperties \ - -lliquidMixtureProperties \ - -lsolidProperties \ - -lsolidMixtureProperties \ - -lthermophysicalFunctions \ -lreactionThermophysicalModels \ - -lSLGThermo \ -lchemistryModel \ - -lradiationModels \ -lODE \ - -lregionModels \ - -lsurfaceFilmModels \ -lcombustionModels \ -lfvOptions \ -lsampling diff --git a/applications/solvers/combustion/SLFMFoam/Mixturefraction.H b/applications/solvers/combustion/SLFMFoam/Mixturefraction.H new file mode 100644 index 000000000..f4c7a205d --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/Mixturefraction.H @@ -0,0 +1,12 @@ +// mixture fraction equation +fvScalarMatrix mfEqn +( + + fvm::div(phi, mf) + - fvm::laplacian((1/Sc)*turbulence->mut(), mf) //let 1/Sc = 1.47 +); + + + mfEqn.relax(); + mfEqn.solve(); + diff --git a/applications/solvers/combustion/SLFMFoam/MixturefractionVar.H b/applications/solvers/combustion/SLFMFoam/MixturefractionVar.H new file mode 100644 index 000000000..bcd78f1df --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/MixturefractionVar.H @@ -0,0 +1,17 @@ +//mixture fraction variance equation +//cf) mfVar should not be ZERO +SDR = turbulence->epsilon() * mfVar / turbulence->k(); +volVectorField Gradmf = fvc::grad(mf); + +fvScalarMatrix mfVarEqn +( + fvm::div(phi, mfVar) + - fvm::laplacian((1/Sc) * turbulence->mut(), mfVar) + == + 2 * ((1/Sc) * turbulence->mut()) * (Gradmf & Gradmf) + - 2 * rho * SDR +); + +mfVarEqn.relax(); + +mfVarEqn.solve(); diff --git a/applications/solvers/combustion/SLFMFoam/SLFMFoam.C b/applications/solvers/combustion/SLFMFoam/SLFMFoam.C index ba96da7ff..0dea1ccf0 100644 --- a/applications/solvers/combustion/SLFMFoam/SLFMFoam.C +++ b/applications/solvers/combustion/SLFMFoam/SLFMFoam.C @@ -22,24 +22,36 @@ License along with OpenFOAM. If not, see . Application - simpleReactingParcelFoam + SLFMFoam Description - Steady state solver for compressible, turbulent flow with reacting, - multiphase particle clouds and optional sources/constraints. + Steady Laminar Flamelet Model(SLFM) solver for turbulent combustion. SLFM + assumes a turbulent flame composed of thin stretched laminar flamelets. + Each flamelet library is generated by external program and imported by + the solver. + +References + A.Y. Klimenko, R.W. Bilger, Progress in Energy and Combustion Science 25 (1999) 595-687 + N. Peters, Turbulent Combustion, Cambridge University Press (2000) + +Contact + POSTECH combustion lab. + huh@postech.ac.kr \*---------------------------------------------------------------------------*/ #include "fvCFD.H" #include "turbulentFluidThermoModel.H" -#include "basicReactingMultiphaseCloud.H" #include "rhoCombustionModel.H" -#include "radiationModel.H" -#include "IOporosityModelList.H" #include "fvOptions.H" -#include "SLGThermo.H" #include "simpleControl.H" +#include "FlameStructure/FlameStructure.H" +#include "interpolateXY.H" +#include "AMC.H" +#include "BetaFunction.H" +#include "BetaGrid.H" +#include "BetaIntegrator.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // int main(int argc, char *argv[]) @@ -57,6 +69,16 @@ int main(int argc, char *argv[]) turbulence->validate(); + //SLFM + Info<<"Read SLFM library"<correct(); - runTime.write(); + if(runTime.write() == true) + { + #include "updateYi.H" //update species + rho.write(); + } Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" << " ClockTime = " << runTime.elapsedClockTime() << " s" diff --git a/applications/solvers/combustion/SLFMFoam/UEqn.H b/applications/solvers/combustion/SLFMFoam/UEqn.H index a29d4b8c5..733ff9277 100644 --- a/applications/solvers/combustion/SLFMFoam/UEqn.H +++ b/applications/solvers/combustion/SLFMFoam/UEqn.H @@ -7,7 +7,6 @@ + turbulence->divDevRhoReff(U) == rho()*g - + parcels.SU(U) + fvOptions(rho, U) ); fvVectorMatrix& UEqn = tUEqn.ref(); diff --git a/applications/solvers/combustion/SLFMFoam/YEqn.H b/applications/solvers/combustion/SLFMFoam/YEqn.H deleted file mode 100644 index 20d55a8d5..000000000 --- a/applications/solvers/combustion/SLFMFoam/YEqn.H +++ /dev/null @@ -1,53 +0,0 @@ -tmp> mvConvection -( - fv::convectionScheme::New - ( - mesh, - fields, - phi, - mesh.divScheme("div(phi,Yi_h)") - ) -); - -{ - combustion->correct(); - dQ = combustion->dQ(); - label inertIndex = -1; - volScalarField Yt(0.0*Y[0]); - - forAll(Y, i) - { - if (Y[i].name() != inertSpecie) - { - volScalarField& Yi = Y[i]; - - fvScalarMatrix YEqn - ( - mvConvection->fvmDiv(phi, Yi) - - fvm::laplacian(turbulence->muEff(), Yi) - == - parcels.SYi(i, Yi) - + combustion->R(Yi) - + fvOptions(rho, Yi) - ); - - YEqn.relax(); - - fvOptions.constrain(YEqn); - - YEqn.solve(mesh.solver("Yi")); - - fvOptions.correct(Yi); - - Yi.max(0.0); - Yt += Yi; - } - else - { - inertIndex = i; - } - } - - Y[inertIndex] = scalar(1) - Yt; - Y[inertIndex].max(0.0); -} diff --git a/applications/solvers/combustion/SLFMFoam/createClouds.H b/applications/solvers/combustion/SLFMFoam/createClouds.H deleted file mode 100644 index 954b74e06..000000000 --- a/applications/solvers/combustion/SLFMFoam/createClouds.H +++ /dev/null @@ -1,9 +0,0 @@ -Info<< "\nConstructing reacting cloud" << endl; -basicReactingMultiphaseCloud parcels -( - "reactingCloud1", - rho, - U, - g, - slgThermo -); diff --git a/applications/solvers/combustion/SLFMFoam/createFieldRefs.H b/applications/solvers/combustion/SLFMFoam/createFieldRefs.H index 5842906a6..6600f4d79 100644 --- a/applications/solvers/combustion/SLFMFoam/createFieldRefs.H +++ b/applications/solvers/combustion/SLFMFoam/createFieldRefs.H @@ -1,2 +1,2 @@ -const volScalarField& T = thermo.T(); const volScalarField& psi = thermo.psi(); +//const volScalarField& T = thermo.T(); diff --git a/applications/solvers/combustion/SLFMFoam/createFields.H b/applications/solvers/combustion/SLFMFoam/createFields.H index 03413049b..a5383e293 100644 --- a/applications/solvers/combustion/SLFMFoam/createFields.H +++ b/applications/solvers/combustion/SLFMFoam/createFields.H @@ -10,8 +10,6 @@ autoPtr combustion rhoReactionThermo& thermo = combustion->thermo(); thermo.validate(args.executable(), "h", "e"); -SLGThermo slgThermo(mesh, thermo); - basicSpecieMixture& composition = thermo.composition(); PtrList& Y = composition.Y(); @@ -26,6 +24,7 @@ if (!composition.contains(inertSpecie)) } volScalarField& p = thermo.p(); +volScalarField& T = thermo.T(); volScalarField rho ( @@ -33,9 +32,7 @@ volScalarField rho ( "rho", runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE + mesh ), thermo.rho() ); @@ -95,6 +92,11 @@ autoPtr turbulence // Set the turbulence into the combustion model combustion->setTurbulence(turbulence()); + +Info<< "Creating field kinetic energy K\n" << endl; +volScalarField K("K", 0.5*magSqr(U)); + + Info<< "Creating multi-variate interpolation scheme\n" << endl; multivariateSurfaceInterpolationScheme::fieldTable fields; @@ -104,20 +106,257 @@ forAll(Y, i) } fields.add(thermo.he()); -volScalarField dQ +/***************************************************************/ +//SLFM-related Fields +Info<<"Creating field mf, mfVar, SDR\n"< Neta(etamax); + +for(label j=0 ; j scalarFieldArray1d; +typedef List scalarFieldArray2d; +typedef List scalarFieldArray3d; + +SLFM::FlameStructure slfmLibrary(runTime.constant()/NstFolder, NstList, etaGrid, composition.species()); + +scalarFieldArray2d& Y_SLFM (slfmLibrary.Y()); +scalarFieldArray2d& W_SLFM (slfmLibrary.W()); +scalarFieldArray1d& T_SLFM (slfmLibrary.T()); +scalarFieldArray1d& Q_SLFM (slfmLibrary.Q()); +scalarFieldArray1d& rho_SLFM (slfmLibrary.rho()); +scalarFieldArray1d& h_SLFM (slfmLibrary.h()); +scalarFieldArray1d& Rgas_SLFM (slfmLibrary.Rgas()); + +for(label n=0 ; n stoi) +upperSum = Foam::pow(stretchFac,upperN)-1.0; +upperSum = upperSum/(stretchFac-1.0); + +//summation of geometric progression (eta < stoi) +lowerSum = Foam::pow(stretchFac,lowerN)-1.0; +lowerSum = lowerSum/(stretchFac-1.0); + +//set reference delta eta +upperDelta = (1.0-stoiMF)/upperSum; +lowerDelta = stoiMF/lowerSum; + +//set boundary and stoichiometric value of eta +etaValue[0] = 0.0; +etaValue[lowerN] = stoiMF; +etaValue.last() = 1.0; + +//set internal value of eta (eta > stoi) +for(label j=0 ; j. - -\*---------------------------------------------------------------------------*/ - -{ - volScalarField& rDeltaT = trDeltaT.ref(); - - const dictionary& pimpleDict = pimple.dict(); - - // Maximum flow Courant number - scalar maxCo(readScalar(pimpleDict.lookup("maxCo"))); - - // Maximum time scale - scalar maxDeltaT(pimpleDict.lookupOrDefault("maxDeltaT", GREAT)); - - // Smoothing parameter (0-1) when smoothing iterations > 0 - scalar rDeltaTSmoothingCoeff - ( - pimpleDict.lookupOrDefault("rDeltaTSmoothingCoeff", 0.1) - ); - - // Damping coefficient (1-0) - scalar rDeltaTDampingCoeff - ( - pimpleDict.lookupOrDefault("rDeltaTDampingCoeff", 0.2) - ); - - // Maximum change in cell temperature per iteration - // (relative to previous value) - scalar alphaTemp(pimpleDict.lookupOrDefault("alphaTemp", 0.05)); - - - Info<< "Time scales min/max:" << endl; - - // Cache old reciprocal time scale field - volScalarField rDeltaT0("rDeltaT0", rDeltaT); - - // Flow time scale - { - rDeltaT.ref() = - ( - fvc::surfaceSum(mag(phi))()() - /((2*maxCo)*mesh.V()*rho()) - ); - - // Limit the largest time scale - rDeltaT.max(1/maxDeltaT); - - Info<< " Flow = " - << gMin(1/rDeltaT.primitiveField()) << ", " - << gMax(1/rDeltaT.primitiveField()) << endl; - } - - // Reaction source time scale - { - volScalarField::Internal rDeltaTT - ( - mag - ( - parcels.hsTrans()/(mesh.V()*runTime.deltaT()) - + combustion->Sh()() - ) - /( - alphaTemp - *rho() - *thermo.Cp()()() - *T() - ) - ); - - Info<< " Temperature = " - << gMin(1/(rDeltaTT.field() + VSMALL)) << ", " - << gMax(1/(rDeltaTT.field() + VSMALL)) << endl; - - rDeltaT.ref() = max - ( - rDeltaT(), - rDeltaTT - ); - } - - // Update tho boundary values of the reciprocal time-step - rDeltaT.correctBoundaryConditions(); - - // Spatially smooth the time scale field - if (rDeltaTSmoothingCoeff < 1.0) - { - fvc::smooth(rDeltaT, rDeltaTSmoothingCoeff); - } - - // Limit rate of change of time scale - // - reduce as much as required - // - only increase at a fraction of old time scale - if - ( - rDeltaTDampingCoeff < 1.0 - && runTime.timeIndex() > runTime.startTimeIndex() + 1 - ) - { - rDeltaT = max - ( - rDeltaT, - (scalar(1.0) - rDeltaTDampingCoeff)*rDeltaT0 - ); - } - - Info<< " Overall = " - << gMin(1/rDeltaT.primitiveField()) - << ", " << gMax(1/rDeltaT.primitiveField()) << endl; -} - - -// ************************************************************************* // diff --git a/applications/solvers/combustion/SLFMFoam/setSDR.H b/applications/solvers/combustion/SLFMFoam/setSDR.H new file mode 100644 index 000000000..c2e416861 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/setSDR.H @@ -0,0 +1,79 @@ +scalarField AMCcoeff(CMC::AMC(etaValue)); + +//cell loop + +forAll(mf, cellI) +{ + + + scalar jl(0), jh(0), vl(0), vh(0); + scalar jfac(0), vfac(0); + scalar C1coeff(0); + + //find eta-index, factor + if(mf[cellI] < etaValue[1]) + { + jl = 0; + jh = 1; + } + else if(mf[cellI] > etaValue[etamax-2]) + { + jl = etamax-2; + jh = etamax-1; + } + else + { + jl = label( interpolateXY(mf[cellI], etaValue, etaIndex) ); + jh = jl+1; + } + jfac = max(0.0, (mf[cellI]-etaValue[jl])/(etaValue[jh]-etaValue[jl])); + + //find var-index, factor + scalar scaledVar = min(0.99999, mfVar[cellI]/(mf[cellI]*(1.0-mf[cellI])+SMALL)); + if(scaledVar < varValue[1]) + { + vl = 0; + vh = 1; + } + else if(scaledVar > varValue[NVar-1]) + { + vl = NVar-1; + vh = NVar; + } + else + { + vl = label( interpolateXY(scaledVar, varValue, varIndex) ); + vh = vl+1; + } + vfac = max(0.0, (scaledVar-varValue[vl])/(varValue[vh]-varValue[vl])); + + + //Bi-linear interpolation on j and v + //Numerical recipes, 2nd Ed. p.117 + C1coeff + = C1table[jl][vl]*(1-jfac)*(1-vfac) + + C1table[jh][vl]*(jfac)*(1-vfac) + + C1table[jh][vh]*(jfac)*(vfac) + + C1table[jl][vh]*(1-jfac)*(vfac); + + jlc[cellI] = jl; + jhc[cellI] = jh; + jfc[cellI] = jfac; + vlc[cellI] = vl; + vhc[cellI] = vh; + vfc[cellI] = vfac; + + + scalar Coeff = AMCcoeff[lowerN]; + + if (Equilibrium == true) + { + Neta[lowerN][cellI] = 0; //Equilibrium + } + else + { + Neta[lowerN][cellI] = SDR[cellI]*Coeff*C1coeff; //CSDR at stoichiometic m.f. + } +} + + diff --git a/applications/solvers/combustion/SLFMFoam/setVarSpace.H b/applications/solvers/combustion/SLFMFoam/setVarSpace.H new file mode 100644 index 000000000..d0413aabc --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/setVarSpace.H @@ -0,0 +1,23 @@ +//set normalized variance space +{ + scalar sum(0), del(0), fac(1.075); + + sum = Foam::pow(fac,NVar) - 1.0; + sum = sum/(fac-1.0); + + del = 1.0/sum; + + varValue[0] = 0.0; + varValue[NVar] = 1.0; + + for(label v=0 ; v NstList[NstList.size()-2]) + { + nhc[cellI] = NstList.size()-2; //just before extinction + nfc[cellI] = 1; + } + else + { + nlc[cellI] = label( interpolateXY(Nst, NstList, NstIndex) ); + nhc[cellI] = nlc[cellI]+1; + nfc[cellI] = (Nst-NstList[nlc[cellI]]) + /(NstList[nhc[cellI]]-NstList[nlc[cellI]]); + } + + //update Favre mean temp. and rho(=1/RT) + scalar jl = jlc[cellI], jh = jhc[cellI]; + scalar vl = vlc[cellI], vh = vhc[cellI]; + scalar nl = nlc[cellI], nh = nhc[cellI]; + scalar jfac = jfc[cellI], vfac = vfc[cellI], nfac = nfc[cellI]; + + scalar T00 = Ttable[jl][vl][nl]*(1-jfac)+Ttable[jh][vl][nl]*jfac; + scalar T01 = Ttable[jl][vl][nh]*(1-jfac)+Ttable[jh][vl][nh]*jfac; + scalar T10 = Ttable[jl][vh][nl]*(1-jfac)+Ttable[jh][vh][nl]*jfac; + scalar T11 = Ttable[jl][vh][nh]*(1-jfac)+Ttable[jh][vh][nh]*jfac; + + scalar T0 = T00*(1-vfac)+T10*vfac; + scalar T1 = T01*(1-vfac)+T11*vfac; + + T[cellI] = T0*(1-nfac)+T1*nfac; + + scalar R00 = Rtable[jl][vl][nl]*(1-jfac)+Rtable[jh][vl][nl]*jfac; + scalar R01 = Rtable[jl][vl][nh]*(1-jfac)+Rtable[jh][vl][nh]*jfac; + scalar R10 = Rtable[jl][vh][nl]*(1-jfac)+Rtable[jh][vh][nl]*jfac; + scalar R11 = Rtable[jl][vh][nh]*(1-jfac)+Rtable[jh][vh][nh]*jfac; + + scalar R0 = R00*(1-vfac)+R10*vfac; + scalar R1 = R01*(1-vfac)+R11*vfac; + + scalar Rgas = R0*(1-nfac)+R1*nfac; + rho[cellI] = p[cellI]/(Rgas*T[cellI]); +} + +//correct processor boundary value of T and rho +T.correctBoundaryConditions(); +rho.correctBoundaryConditions(); + +//update outlet boundary value of rho +const fvPatchList& patches = mesh.boundary(); +forAll(patches, patchI) +{ + const fvPatch& curPatch = patches[patchI]; + if(curPatch.name() == outletName) + { + forAll(curPatch, faceI) + { + label cellI = curPatch.faceCells()[faceI]; + rho.boundaryFieldRef()[patchI][faceI] = rho[cellI]; + } + } +} + diff --git a/applications/solvers/combustion/SLFMFoam/updateYi.H b/applications/solvers/combustion/SLFMFoam/updateYi.H new file mode 100644 index 000000000..09a38c766 --- /dev/null +++ b/applications/solvers/combustion/SLFMFoam/updateYi.H @@ -0,0 +1,56 @@ +//force the saving of the old-time values +//important!! + +forAll(rho, cellI) +{ + scalar& Nst = Neta[lowerN][cellI]; + + //find Nst index + if(Nst < NstList.first()) + { + nlc[cellI] = 0; + nfc[cellI] = 0; + } + else if(Nst > NstList[NstList.size()-2]) + { + nhc[cellI] = NstList.size()-2; //just before extinction + nfc[cellI] = 1; + } + else + { + nlc[cellI] = label( interpolateXY(Nst, NstList, NstIndex) ); + nhc[cellI] = nlc[cellI]+1; + nfc[cellI] = (Nst-NstList[nlc[cellI]]) + /(NstList[nhc[cellI]]-NstList[nlc[cellI]]); + } + + //update Favre mean temp. and rho(=1/RT) + scalar jl = jlc[cellI], jh = jhc[cellI]; + scalar vl = vlc[cellI], vh = vhc[cellI]; + scalar nl = nlc[cellI], nh = nhc[cellI]; + scalar jfac = jfc[cellI], vfac = vfc[cellI], nfac = nfc[cellI]; + + + forAll(Y, yi) + { + scalar y00 = Ytable[yi][jl][vl][nl]*(1-jfac)+Ytable[yi][jh][vl][nl]*jfac; + scalar y01 = Ytable[yi][jl][vl][nh]*(1-jfac)+Ytable[yi][jh][vl][nh]*jfac; + scalar y10 = Ytable[yi][jl][vh][nl]*(1-jfac)+Ytable[yi][jh][vh][nl]*jfac; + scalar y11 = Ytable[yi][jl][vh][nh]*(1-jfac)+Ytable[yi][jh][vh][nh]*jfac; + + scalar y0 = y00*(1-vfac)+y10*vfac; + scalar y1 = y01*(1-vfac)+y11*vfac; + + scalar ygas = y0*(1-nfac)+y1*nfac; + + Y[yi][cellI] = ygas; + } +} + +//correct processor boundary value of T and rho +forAll(Y, yi) +{ + Y[yi].correctBoundaryConditions(); +} + +