diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files index b982734e8bc4d6ddd8a631a61e3342db949e4973..5e097730fd6259abb845f93137b3abcf73b7b565 100644 --- a/src/finiteVolume/Make/files +++ b/src/finiteVolume/Make/files @@ -342,6 +342,7 @@ $(limitedGradSchemes)/cellMDLimitedGrad/cellMDLimitedGrads.C snGradSchemes = finiteVolume/snGradSchemes $(snGradSchemes)/snGradScheme/snGradSchemes.C $(snGradSchemes)/correctedSnGrad/correctedSnGrads.C +$(snGradSchemes)/faceCorrectedSnGrad/faceCorrectedSnGrads.C $(snGradSchemes)/limitedSnGrad/limitedSnGrads.C $(snGradSchemes)/uncorrectedSnGrad/uncorrectedSnGrads.C $(snGradSchemes)/orthogonalSnGrad/orthogonalSnGrads.C diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/faceCorrectedSnGrad/faceCorrectedSnGrad.C b/src/finiteVolume/finiteVolume/snGradSchemes/faceCorrectedSnGrad/faceCorrectedSnGrad.C new file mode 100644 index 0000000000000000000000000000000000000000..6b04bb2efb4be28927c3b500fd98b2ee4837f615 --- /dev/null +++ b/src/finiteVolume/finiteVolume/snGradSchemes/faceCorrectedSnGrad/faceCorrectedSnGrad.C @@ -0,0 +1,167 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2012 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "faceCorrectedSnGrad.H" +#include "volPointInterpolation.H" +#include "triangle.H" + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +template<class Type> +Foam::fv::faceCorrectedSnGrad<Type>::~faceCorrectedSnGrad() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template<class Type> +Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh> > +Foam::fv::faceCorrectedSnGrad<Type>::fullGradCorrection +( + const GeometricField<Type, fvPatchField, volMesh>& vf +) const +{ + const fvMesh& mesh = this->mesh(); + + GeometricField<Type, pointPatchField, pointMesh> pvf + ( + volPointInterpolation::New(mesh).interpolate(vf) + ); + + // construct GeometricField<Type, fvsPatchField, surfaceMesh> + tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsfCorr + ( + new GeometricField<Type, fvsPatchField, surfaceMesh> + ( + IOobject + ( + "snGradCorr("+vf.name()+')', + vf.instance(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh, + vf.dimensions()*mesh.nonOrthDeltaCoeffs().dimensions() + ) + ); + + Field<Type>& sfCorr = tsfCorr().internalField(); + + const pointField& points = mesh.points(); + const faceList& faces = mesh.faces(); + const vectorField& Sf = mesh.Sf().internalField(); + const vectorField& C = mesh.C().internalField(); + const scalarField& magSf = mesh.magSf().internalField(); + const labelList& owner = mesh.owner(); + const labelList& neighbour = mesh.neighbour(); + + forAll(sfCorr, facei) + { + typename outerProduct<vector, Type>::type fgrad + ( + outerProduct<vector, Type>::type::zero + ); + + const face& fi = faces[facei]; + + vector nf(Sf[facei]/magSf[facei]); + + for (label pi=0; pi<fi.size(); pi++) + { + // Next point index + label pj = (pi+1)%fi.size(); + + // Edge normal in plane of face + vector edgen(nf^(points[fi[pj]] - points[fi[pi]])); + + // Edge centre field value + Type pvfe(0.5*(pvf[fi[pj]] + pvf[fi[pi]])); + + // Integrate face gradient + fgrad += edgen*pvfe; + } + + // Finalize face-gradient by dividing by face area + fgrad /= magSf[facei]; + + // Calculate correction vector + vector dCorr(C[neighbour[facei]] - C[owner[facei]]); + dCorr /= (nf & dCorr); + + // if (mag(dCorr) > 2) dCorr *= 2/mag(dCorr); + + sfCorr[facei] = dCorr&fgrad; + } + + tsfCorr().boundaryField() = pTraits<Type>::zero; + + return tsfCorr; +} + + +template<class Type> +Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh> > +Foam::fv::faceCorrectedSnGrad<Type>::correction +( + const GeometricField<Type, fvPatchField, volMesh>& vf +) const +{ + const fvMesh& mesh = this->mesh(); + + // construct GeometricField<Type, fvsPatchField, surfaceMesh> + tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tssf + ( + new GeometricField<Type, fvsPatchField, surfaceMesh> + ( + IOobject + ( + "snGradCorr("+vf.name()+')', + vf.instance(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh, + vf.dimensions()*mesh.nonOrthDeltaCoeffs().dimensions() + ) + ); + GeometricField<Type, fvsPatchField, surfaceMesh>& ssf = tssf(); + + for (direction cmpt = 0; cmpt < pTraits<Type>::nComponents; cmpt++) + { + ssf.replace + ( + cmpt, + faceCorrectedSnGrad<typename pTraits<Type>::cmptType>(mesh) + .fullGradCorrection(vf.component(cmpt)) + ); + } + + return tssf; +} + + +// ************************************************************************* // diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/faceCorrectedSnGrad/faceCorrectedSnGrad.H b/src/finiteVolume/finiteVolume/snGradSchemes/faceCorrectedSnGrad/faceCorrectedSnGrad.H new file mode 100644 index 0000000000000000000000000000000000000000..a13dd62d18e2c6ccb775d549ea1d98542ec21981 --- /dev/null +++ b/src/finiteVolume/finiteVolume/snGradSchemes/faceCorrectedSnGrad/faceCorrectedSnGrad.H @@ -0,0 +1,157 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2012 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 <http://www.gnu.org/licenses/>. + +Class + Foam::fv::faceCorrectedSnGrad + +Description + Simple central-difference snGrad scheme with non-orthogonal correction. + +SourceFiles + faceCorrectedSnGrad.C + +\*---------------------------------------------------------------------------*/ + +#ifndef faceCorrectedSnGrad_H +#define faceCorrectedSnGrad_H + +#include "snGradScheme.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace fv +{ + +/*---------------------------------------------------------------------------*\ + Class faceCorrectedSnGrad Declaration +\*---------------------------------------------------------------------------*/ + +template<class Type> +class faceCorrectedSnGrad +: + public snGradScheme<Type> +{ + // Private Member Functions + + //- Disallow default bitwise assignment + void operator=(const faceCorrectedSnGrad&); + + +public: + + //- Runtime type information + TypeName("faceCorrected"); + + + // Constructors + + //- Construct from mesh + faceCorrectedSnGrad(const fvMesh& mesh) + : + snGradScheme<Type>(mesh) + {} + + + //- Construct from mesh and data stream + faceCorrectedSnGrad(const fvMesh& mesh, Istream&) + : + snGradScheme<Type>(mesh) + {} + + + //- Destructor + virtual ~faceCorrectedSnGrad(); + + + // Member Functions + + //- Return the interpolation weighting factors for the given field + virtual tmp<surfaceScalarField> deltaCoeffs + ( + const GeometricField<Type, fvPatchField, volMesh>& + ) const + { + return this->mesh().nonOrthDeltaCoeffs(); + } + + //- Return true if this scheme uses an explicit correction + virtual bool corrected() const + { + return true; + } + + //- Return the explicit correction to the faceCorrectedSnGrad + // for the given field using the gradient of the field + tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > + fullGradCorrection + ( + const GeometricField<Type, fvPatchField, volMesh>& + ) const; + + //- Return the explicit correction to the faceCorrectedSnGrad + // for the given field using the gradients of the field components + virtual tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > + correction(const GeometricField<Type, fvPatchField, volMesh>&) const; +}; + + +// * * * * * * * * Template Member Function Specialisations * * * * * * * * // + +template<> +tmp<surfaceScalarField> faceCorrectedSnGrad<scalar>::correction +( + const volScalarField& vsf +) const; + + +template<> +tmp<surfaceVectorField> faceCorrectedSnGrad<vector>::correction +( + const volVectorField& vvf +) const; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace fv + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "faceCorrectedSnGrad.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/faceCorrectedSnGrad/faceCorrectedSnGrads.C b/src/finiteVolume/finiteVolume/snGradSchemes/faceCorrectedSnGrad/faceCorrectedSnGrads.C new file mode 100644 index 0000000000000000000000000000000000000000..ba6d83542eb2cb3c70a152c99d457764382cb0de --- /dev/null +++ b/src/finiteVolume/finiteVolume/snGradSchemes/faceCorrectedSnGrad/faceCorrectedSnGrads.C @@ -0,0 +1,62 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2012 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "faceCorrectedSnGrad.H" +#include "fvMesh.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace fv +{ + makeSnGradScheme(faceCorrectedSnGrad) +} +} + + +template<> +Foam::tmp<Foam::surfaceScalarField> +Foam::fv::faceCorrectedSnGrad<Foam::scalar>::correction +( + const volScalarField& vsf +) const +{ + return fullGradCorrection(vsf); +} + + +template<> +Foam::tmp<Foam::surfaceVectorField> +Foam::fv::faceCorrectedSnGrad<Foam::vector>::correction +( + const volVectorField& vvf +) const +{ + return fullGradCorrection(vvf); +} + + +// ************************************************************************* // diff --git a/src/finiteVolume/fvMesh/fvMesh.H b/src/finiteVolume/fvMesh/fvMesh.H index 97c1b0db589b4305f568d1a34a798dd8c53c820f..7cb9d97a00b0b8105d0f2331a6047dd3034d5c27 100644 --- a/src/finiteVolume/fvMesh/fvMesh.H +++ b/src/finiteVolume/fvMesh/fvMesh.H @@ -312,6 +312,9 @@ public: //- Return face centres as surfaceVectorField const surfaceVectorField& Cf() const; + //- Return face deltas as surfaceVectorField + tmp<surfaceVectorField> delta() const; + // Edit diff --git a/src/finiteVolume/fvMesh/fvMeshGeometry.C b/src/finiteVolume/fvMesh/fvMeshGeometry.C index 5abdb45d118c210ffb10ade57be95fb7b2cbe587..ddd395a722966801a068f41e6306b64bb9148621 100644 --- a/src/finiteVolume/fvMesh/fvMeshGeometry.C +++ b/src/finiteVolume/fvMesh/fvMeshGeometry.C @@ -372,6 +372,53 @@ const surfaceVectorField& fvMesh::Cf() const } +tmp<surfaceVectorField> fvMesh::delta() const +{ + if (debug) + { + Info<< "void fvMesh::delta() : " + << "calculating face deltas" + << endl; + } + + tmp<surfaceVectorField> tdelta + ( + new surfaceVectorField + ( + IOobject + ( + "delta", + pointsInstance(), + meshSubDir, + *this, + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + *this, + dimLength + ) + ); + surfaceVectorField& delta = tdelta(); + + const volVectorField& C = this->C(); + const labelUList& owner = this->owner(); + const labelUList& neighbour = this->neighbour(); + + forAll(owner, facei) + { + delta[facei] = C[neighbour[facei]] - C[owner[facei]]; + } + + forAll(delta.boundaryField(), patchi) + { + delta.boundaryField()[patchi] = boundary()[patchi].delta(); + } + + return tdelta; +} + + const surfaceScalarField& fvMesh::phi() const { if (!phiPtr_) diff --git a/src/thermophysicalModels/basic/basicThermo/basicThermo.C b/src/thermophysicalModels/basic/basicThermo/basicThermo.C index 84e30e31a4767807499e91b55180f838ab3639e3..e2fc874b5ef6ba930a1c85a745cf2032891c6e86 100644 --- a/src/thermophysicalModels/basic/basicThermo/basicThermo.C +++ b/src/thermophysicalModels/basic/basicThermo/basicThermo.C @@ -112,7 +112,7 @@ Foam::basicThermo::basicThermo ( IOobject ( - phasePropertyName("alpha"), + phasePropertyName("thermo:alpha"), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -167,7 +167,7 @@ Foam::basicThermo::basicThermo ( IOobject ( - phasePropertyName("alpha"), + phasePropertyName("thermo:alpha"), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -237,7 +237,7 @@ const Foam::basicThermo& Foam::basicThermo::lookupThermo void Foam::basicThermo::validate ( - const word& app, + const string& app, const word& a ) const { @@ -252,7 +252,7 @@ void Foam::basicThermo::validate void Foam::basicThermo::validate ( - const word& app, + const string& app, const word& a, const word& b ) const @@ -275,7 +275,7 @@ void Foam::basicThermo::validate void Foam::basicThermo::validate ( - const word& app, + const string& app, const word& a, const word& b, const word& c @@ -301,7 +301,7 @@ void Foam::basicThermo::validate void Foam::basicThermo::validate ( - const word& app, + const string& app, const word& a, const word& b, const word& c, diff --git a/src/thermophysicalModels/basic/basicThermo/basicThermo.H b/src/thermophysicalModels/basic/basicThermo/basicThermo.H index 659b6f35ca03512a0d98f5ddef38f37a3e69add7..e652c30ee6cc3216a0227160a9450e99e5223f78 100644 --- a/src/thermophysicalModels/basic/basicThermo/basicThermo.H +++ b/src/thermophysicalModels/basic/basicThermo/basicThermo.H @@ -185,7 +185,7 @@ public: // with energy forms supported by the application void validate ( - const word& app, + const string& app, const word& ) const; @@ -193,7 +193,7 @@ public: // with energy forms supported by the application void validate ( - const word& app, + const string& app, const word&, const word& ) const; @@ -202,7 +202,7 @@ public: // with energy forms supported by the application void validate ( - const word& app, + const string& app, const word&, const word&, const word& @@ -212,7 +212,7 @@ public: // with energy forms supported by the application void validate ( - const word& app, + const string& app, const word&, const word&, const word&, @@ -263,6 +263,14 @@ public: //- Enthalpy/Internal energy [J/kg] virtual const volScalarField& he() const = 0; + //- Enthalpy/Internal energy + // for given pressure and temperature [J/kg] + virtual tmp<volScalarField> he + ( + const volScalarField& p, + const volScalarField& T + ) const = 0; + //- Enthalpy/Internal energy for cell-set [J/kg] virtual tmp<scalarField> he ( diff --git a/src/thermophysicalModels/basic/heThermo/heThermo.C b/src/thermophysicalModels/basic/heThermo/heThermo.C index 8efaf0f5d50ee60a97cb1bbda81fe132b39c3ce1..29634e276936910be17d608bcbe43871a05a4471 100644 --- a/src/thermophysicalModels/basic/heThermo/heThermo.C +++ b/src/thermophysicalModels/basic/heThermo/heThermo.C @@ -174,7 +174,10 @@ Foam::heThermo<BasicThermo, MixtureType>::heThermo ( IOobject ( - BasicThermo::phasePropertyName(MixtureType::thermoType::heName()), + BasicThermo::phasePropertyName + ( + MixtureType::thermoType::heName() + ), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -205,7 +208,10 @@ Foam::heThermo<BasicThermo, MixtureType>::heThermo ( IOobject ( - BasicThermo::phasePropertyName(MixtureType::thermoType::heName()), + BasicThermo::phasePropertyName + ( + MixtureType::thermoType::heName() + ), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -230,6 +236,60 @@ Foam::heThermo<BasicThermo, MixtureType>::~heThermo() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +template<class BasicThermo, class MixtureType> +Foam::tmp<Foam::volScalarField> Foam::heThermo<BasicThermo, MixtureType>::he +( + const volScalarField& p, + const volScalarField& T +) const +{ + const fvMesh& mesh = this->T_.mesh(); + + tmp<volScalarField> the + ( + new volScalarField + ( + IOobject + ( + "he", + mesh.time().timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh, + he_.dimensions() + ) + ); + + volScalarField& he = the(); + scalarField& heCells = he.internalField(); + const scalarField& pCells = p.internalField(); + const scalarField& TCells = T.internalField(); + + forAll(heCells, celli) + { + heCells[celli] = + this->cellMixture(celli).HE(pCells[celli], TCells[celli]); + } + + forAll(he.boundaryField(), patchi) + { + scalarField& hep = he.boundaryField()[patchi]; + const scalarField& pp = p.boundaryField()[patchi]; + const scalarField& Tp = T.boundaryField()[patchi]; + + forAll(hep, facei) + { + hep[facei] = + this->patchFaceMixture(patchi, facei).HE(pp[facei], Tp[facei]); + } + } + + return the; +} + + template<class BasicThermo, class MixtureType> Foam::tmp<Foam::scalarField> Foam::heThermo<BasicThermo, MixtureType>::he ( diff --git a/src/thermophysicalModels/basic/heThermo/heThermo.H b/src/thermophysicalModels/basic/heThermo/heThermo.H index ef9919a99cdb5b3a303780d45045ccfcfe5bd176..f77d017e1396afd3840611449829eb89fbc6653a 100644 --- a/src/thermophysicalModels/basic/heThermo/heThermo.H +++ b/src/thermophysicalModels/basic/heThermo/heThermo.H @@ -161,6 +161,14 @@ public: // Fields derived from thermodynamic state variables + //- Enthalpy/Internal energy + // for given pressure and temperature [J/kg] + virtual tmp<volScalarField> he + ( + const volScalarField& p, + const volScalarField& T + ) const; + //- Enthalpy/Internal energy for cell-set [J/kg] virtual tmp<scalarField> he ( diff --git a/src/thermophysicalModels/basic/psiThermo/psiThermo.C b/src/thermophysicalModels/basic/psiThermo/psiThermo.C index fd4bc962c903bcda157391c428df833b492191d0..dd69890313406b8e9ca182b9d0510aab991766e0 100644 --- a/src/thermophysicalModels/basic/psiThermo/psiThermo.C +++ b/src/thermophysicalModels/basic/psiThermo/psiThermo.C @@ -44,7 +44,7 @@ Foam::psiThermo::psiThermo(const fvMesh& mesh, const word& phaseName) ( IOobject ( - phasePropertyName("psi"), + phasePropertyName("thermo:psi"), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -58,7 +58,7 @@ Foam::psiThermo::psiThermo(const fvMesh& mesh, const word& phaseName) ( IOobject ( - phasePropertyName("mu"), + phasePropertyName("thermo:mu"), mesh.time().timeName(), mesh, IOobject::NO_READ, diff --git a/src/thermophysicalModels/basic/rhoThermo/rhoThermo.C b/src/thermophysicalModels/basic/rhoThermo/rhoThermo.C index 4aa7835caca575346ea7341ee3346b8140a173d5..d34186eea9ee534ceec2376b91acba11c92c8466 100644 --- a/src/thermophysicalModels/basic/rhoThermo/rhoThermo.C +++ b/src/thermophysicalModels/basic/rhoThermo/rhoThermo.C @@ -43,7 +43,7 @@ Foam::rhoThermo::rhoThermo(const fvMesh& mesh, const word& phaseName) ( IOobject ( - phasePropertyName("rhoThermo"), + phasePropertyName("thermo:rho"), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -57,7 +57,7 @@ Foam::rhoThermo::rhoThermo(const fvMesh& mesh, const word& phaseName) ( IOobject ( - phasePropertyName("psi"), + phasePropertyName("thermo:psi"), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -71,7 +71,7 @@ Foam::rhoThermo::rhoThermo(const fvMesh& mesh, const word& phaseName) ( IOobject ( - phasePropertyName("mu"), + phasePropertyName("thermo:mu"), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -95,7 +95,7 @@ Foam::rhoThermo::rhoThermo ( IOobject ( - phasePropertyName("rhoThermo"), + phasePropertyName("thermo:rho"), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -109,7 +109,7 @@ Foam::rhoThermo::rhoThermo ( IOobject ( - phasePropertyName("psi"), + phasePropertyName("thermo:psi"), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -123,7 +123,7 @@ Foam::rhoThermo::rhoThermo ( IOobject ( - phasePropertyName("mu"), + phasePropertyName("thermo:mu"), mesh.time().timeName(), mesh, IOobject::NO_READ, diff --git a/src/thermophysicalModels/basic/rhoThermo/rhoThermos.C b/src/thermophysicalModels/basic/rhoThermo/rhoThermos.C index fcfc3437dd35810afbc88e6f4160879852e23236..beaa43b7350ab95fb0ed3a57535931dbaf2fc7bf 100644 --- a/src/thermophysicalModels/basic/rhoThermo/rhoThermos.C +++ b/src/thermophysicalModels/basic/rhoThermo/rhoThermos.C @@ -30,6 +30,7 @@ License #include "perfectGas.H" #include "incompressiblePerfectGas.H" #include "rhoConst.H" +#include "perfectFluid.H" #include "hConstThermo.H" #include "janafThermo.H" #include "sensibleEnthalpy.H" @@ -101,6 +102,18 @@ makeThermo specie ); +makeThermo +( + rhoThermo, + heRhoThermo, + pureMixture, + constTransport, + sensibleEnthalpy, + hConstThermo, + perfectFluid, + specie +); + makeThermo ( rhoThermo, @@ -200,6 +213,18 @@ makeThermo specie ); +makeThermo +( + rhoThermo, + heRhoThermo, + pureMixture, + constTransport, + sensibleInternalEnergy, + hConstThermo, + perfectFluid, + specie +); + makeThermo ( rhoThermo, diff --git a/src/thermophysicalModels/solidThermo/solidThermo/solidThermo.C b/src/thermophysicalModels/solidThermo/solidThermo/solidThermo.C index 102414e3fbb6ce0618a581b7e83546c07f4860f5..2e12bcf213d0ba2e82b213b55be0891956d95c6c 100644 --- a/src/thermophysicalModels/solidThermo/solidThermo/solidThermo.C +++ b/src/thermophysicalModels/solidThermo/solidThermo/solidThermo.C @@ -50,7 +50,7 @@ Foam::solidThermo::solidThermo ( IOobject ( - "rhoThermo", + phasePropertyName("thermo:rho"), mesh.time().timeName(), mesh, IOobject::NO_READ, @@ -74,7 +74,7 @@ Foam::solidThermo::solidThermo ( IOobject ( - "rhoThermo", + phasePropertyName("thermo:rho"), mesh.time().timeName(), mesh, IOobject::NO_READ, diff --git a/src/thermophysicalModels/specie/equationOfState/perfectFluid/perfectFluid.C b/src/thermophysicalModels/specie/equationOfState/perfectFluid/perfectFluid.C new file mode 100644 index 0000000000000000000000000000000000000000..844f293e2dd580a904cc9fa93ba877dfd8d0df5a --- /dev/null +++ b/src/thermophysicalModels/specie/equationOfState/perfectFluid/perfectFluid.C @@ -0,0 +1,76 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2012 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "perfectFluid.H" +#include "IOstreams.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template<class Specie> +Foam::perfectFluid<Specie>::perfectFluid(Istream& is) +: + Specie(is), + rho0_(readScalar(is)) +{ + is.check("perfectFluid<Specie>::perfectFluid(Istream& is)"); +} + + +template<class Specie> +Foam::perfectFluid<Specie>::perfectFluid(const dictionary& dict) +: + Specie(dict), + rho0_(readScalar(dict.subDict("equationOfState").lookup("rho0"))) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template<class Specie> +void Foam::perfectFluid<Specie>::write(Ostream& os) const +{ + Specie::write(os); + + dictionary dict("equationOfState"); + dict.add("rho0", rho0_); + + os << indent << dict.dictName() << dict; +} + + +// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * // + +template<class Specie> +Foam::Ostream& Foam::operator<<(Ostream& os, const perfectFluid<Specie>& pf) +{ + os << static_cast<const Specie&>(pf) + << token::SPACE << pf.rho0_; + + os.check("Ostream& operator<<(Ostream&, const perfectFluid<Specie>&)"); + return os; +} + + +// ************************************************************************* // diff --git a/src/thermophysicalModels/specie/equationOfState/perfectFluid/perfectFluid.H b/src/thermophysicalModels/specie/equationOfState/perfectFluid/perfectFluid.H new file mode 100644 index 0000000000000000000000000000000000000000..76cf80c4a8fd504ad52bcef4f0f4a55653c20dda --- /dev/null +++ b/src/thermophysicalModels/specie/equationOfState/perfectFluid/perfectFluid.H @@ -0,0 +1,223 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2012 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 <http://www.gnu.org/licenses/>. + +Class + Foam::perfectFluid + +Description + Perfect gas equation of state. + +SourceFiles + perfectFluidI.H + perfectFluid.C + +\*---------------------------------------------------------------------------*/ + +#ifndef perfectFluid_H +#define perfectFluid_H + +#include "autoPtr.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// Forward declaration of friend functions and operators + +template<class Specie> class perfectFluid; + +template<class Specie> +inline perfectFluid<Specie> operator+ +( + const perfectFluid<Specie>&, + const perfectFluid<Specie>& +); + +template<class Specie> +inline perfectFluid<Specie> operator- +( + const perfectFluid<Specie>&, + const perfectFluid<Specie>& +); + +template<class Specie> +inline perfectFluid<Specie> operator* +( + const scalar, + const perfectFluid<Specie>& +); + +template<class Specie> +inline perfectFluid<Specie> operator== +( + const perfectFluid<Specie>&, + const perfectFluid<Specie>& +); + +template<class Specie> +Ostream& operator<< +( + Ostream&, + const perfectFluid<Specie>& +); + + +/*---------------------------------------------------------------------------*\ + Class perfectFluid Declaration +\*---------------------------------------------------------------------------*/ + +template<class Specie> +class perfectFluid +: + public Specie +{ + // Private data + + //- The reference density + scalar rho0_; + +public: + + // Constructors + + //- Construct from components + inline perfectFluid(const Specie& sp, const scalar rho0); + + //- Construct from Istream + perfectFluid(Istream&); + + //- Construct from dictionary + perfectFluid(const dictionary& dict); + + //- Construct as named copy + inline perfectFluid(const word& name, const perfectFluid&); + + //- Construct and return a clone + inline autoPtr<perfectFluid> clone() const; + + // Selector from Istream + inline static autoPtr<perfectFluid> New(Istream& is); + + // Selector from dictionary + inline static autoPtr<perfectFluid> New(const dictionary& dict); + + + // Member functions + + //- Return the instantiated type name + static word typeName() + { + return "perfectFluid<" + word(Specie::typeName_()) + '>'; + } + + + // Fundamental properties + + //- Is the equation of state is incompressible i.e. rho != f(p) + static const bool incompressible = false; + + //- Is the equation of state is isochoric i.e. rho = const + static const bool isochoric = false; + + //- Return density [kg/m^3] + inline scalar rho(scalar p, scalar T) const; + + //- Return compressibility rho/p [s^2/m^2] + inline scalar psi(scalar p, scalar T) const; + + //- Return compression factor [] + inline scalar Z(scalar p, scalar T) const; + + //- Return (cp - cv) [J/(kmol K] + inline scalar cpMcv(scalar p, scalar T) const; + + + // IO + + //- Write to Ostream + void write(Ostream& os) const; + + + // Member operators + + inline void operator+=(const perfectFluid&); + inline void operator-=(const perfectFluid&); + + inline void operator*=(const scalar); + + + // Friend operators + + friend perfectFluid operator+ <Specie> + ( + const perfectFluid&, + const perfectFluid& + ); + + friend perfectFluid operator- <Specie> + ( + const perfectFluid&, + const perfectFluid& + ); + + friend perfectFluid operator* <Specie> + ( + const scalar s, + const perfectFluid& + ); + + friend perfectFluid operator== <Specie> + ( + const perfectFluid&, + const perfectFluid& + ); + + + // Ostream Operator + + friend Ostream& operator<< <Specie> + ( + Ostream&, + const perfectFluid& + ); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "perfectFluidI.H" + +#ifdef NoRepository +# include "perfectFluid.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/thermophysicalModels/specie/equationOfState/perfectFluid/perfectFluidI.H b/src/thermophysicalModels/specie/equationOfState/perfectFluid/perfectFluidI.H new file mode 100644 index 0000000000000000000000000000000000000000..21fc092fe2adca3d3875247001d26cacb61309da --- /dev/null +++ b/src/thermophysicalModels/specie/equationOfState/perfectFluid/perfectFluidI.H @@ -0,0 +1,220 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2012 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "perfectFluid.H" + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +template<class Specie> +inline Foam::perfectFluid<Specie>::perfectFluid +( + const Specie& sp, + const scalar rho0 +) +: + Specie(sp), + rho0_(rho0) +{} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template<class Specie> +inline Foam::perfectFluid<Specie>::perfectFluid +( + const word& name, + const perfectFluid<Specie>& pf +) +: + Specie(name, pf), + rho0_(pf.rho0_) +{} + + +template<class Specie> +inline Foam::autoPtr<Foam::perfectFluid<Specie> > +Foam::perfectFluid<Specie>::clone() const +{ + return autoPtr<perfectFluid<Specie> >(new perfectFluid<Specie>(*this)); +} + + +template<class Specie> +inline Foam::autoPtr<Foam::perfectFluid<Specie> > +Foam::perfectFluid<Specie>::New(Istream& is) +{ + return autoPtr<perfectFluid<Specie> >(new perfectFluid<Specie>(is)); +} + + +template<class Specie> +inline Foam::autoPtr<Foam::perfectFluid<Specie> > +Foam::perfectFluid<Specie>::New +( + const dictionary& dict +) +{ + return autoPtr<perfectFluid<Specie> >(new perfectFluid<Specie>(dict)); +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template<class Specie> +inline Foam::scalar Foam::perfectFluid<Specie>::rho(scalar p, scalar T) const +{ + return rho0_ + p/(this->R()*T); +} + + +template<class Specie> +inline Foam::scalar Foam::perfectFluid<Specie>::psi(scalar, scalar T) const +{ + return 1.0/(this->R()*T); +} + + +template<class Specie> +inline Foam::scalar Foam::perfectFluid<Specie>::Z(scalar, scalar) const +{ + return 1.0; +} + + +template<class Specie> +inline Foam::scalar Foam::perfectFluid<Specie>::cpMcv(scalar, scalar) const +{ + return this->RR; +} + + +// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * // + +template<class Specie> +inline void Foam::perfectFluid<Specie>::operator+= +( + const perfectFluid<Specie>& pf +) +{ + scalar molr1 = this->nMoles(); + + Specie::operator+=(pf); + + molr1 /= this->nMoles(); + scalar molr2 = pf.nMoles()/this->nMoles(); + + rho0_ = molr1*rho0_ + molr2*pf.rho0_; +} + + +template<class Specie> +inline void Foam::perfectFluid<Specie>::operator-= +( + const perfectFluid<Specie>& pf +) +{ + scalar molr1 = this->nMoles(); + + Specie::operator-=(pf); + + molr1 /= this->nMoles(); + scalar molr2 = pf.nMoles()/this->nMoles(); + + rho0_ = molr1*rho0_ - molr2*pf.rho0_; +} + + +template<class Specie> +inline void Foam::perfectFluid<Specie>::operator*=(const scalar s) +{ + Specie::operator*=(s); +} + + +// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * // + +template<class Specie> +inline Foam::perfectFluid<Specie> Foam::operator+ +( + const perfectFluid<Specie>& pf1, + const perfectFluid<Specie>& pf2 +) +{ + scalar nMoles = pf1.nMoles() + pf2.nMoles(); + scalar molr1 = pf1.nMoles()/nMoles; + scalar molr2 = pf2.nMoles()/nMoles; + + return rhoConst<Specie> + ( + static_cast<const Specie&>(pf1) + + static_cast<const Specie&>(pf2), + molr1*pf1.rho0_ + molr2*pf2.rho0_ + ); +} + + +template<class Specie> +inline Foam::perfectFluid<Specie> Foam::operator- +( + const perfectFluid<Specie>& pf1, + const perfectFluid<Specie>& pf2 +) +{ + scalar nMoles = pf1.nMoles() + pf2.nMoles(); + scalar molr1 = pf1.nMoles()/nMoles; + scalar molr2 = pf2.nMoles()/nMoles; + + return rhoConst<Specie> + ( + static_cast<const Specie&>(pf1) + - static_cast<const Specie&>(pf2), + molr1*pf1.rho0_ - molr2*pf2.rho0_ + ); +} + + +template<class Specie> +inline Foam::perfectFluid<Specie> Foam::operator* +( + const scalar s, + const perfectFluid<Specie>& pf +) +{ + return perfectFluid<Specie>(s*static_cast<const Specie&>(pf), pf.rho0_); +} + + +template<class Specie> +inline Foam::perfectFluid<Specie> Foam::operator== +( + const perfectFluid<Specie>& pf1, + const perfectFluid<Specie>& pf2 +) +{ + return pf2 - pf1; +} + + +// ************************************************************************* // diff --git a/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.C b/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.C index ddde795545d66396d6675a73f2acf9fc45f4330b..321c07443b1297f776ac6b0d4580ef8e09b4baa2 100644 --- a/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.C +++ b/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.C @@ -97,6 +97,19 @@ tmp<volScalarField> kOmegaSST::F3() const } +tmp<volScalarField> kOmegaSST::F23() const +{ + tmp<volScalarField> f23(F2()); + + if (F3_) + { + f23() *= F3(); + } + + return f23; +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // kOmegaSST::kOmegaSST @@ -228,6 +241,15 @@ kOmegaSST::kOmegaSST 10.0 ) ), + F3_ + ( + Switch::lookupOrAddToDict + ( + "F3", + coeffDict_, + false + ) + ), y_(mesh_), @@ -289,7 +311,7 @@ kOmegaSST::kOmegaSST / max ( a1_*omega_, - b1_*F2()*F3()*sqrt(2.0)*mag(symm(fvc::grad(U_))) + b1_*F23()*sqrt(2.0)*mag(symm(fvc::grad(U_))) ) ); mut_.correctBoundaryConditions(); @@ -370,6 +392,7 @@ bool kOmegaSST::read() a1_.readIfPresent(coeffDict()); b1_.readIfPresent(coeffDict()); c1_.readIfPresent(coeffDict()); + F3_.readIfPresent("F3", coeffDict()); return true; } @@ -470,7 +493,7 @@ void kOmegaSST::correct() // Re-calculate viscosity - mut_ = a1_*rho_*k_/max(a1_*omega_, b1_*F2()*F3()*sqrt(S2)); + mut_ = a1_*rho_*k_/max(a1_*omega_, b1_*F23()*sqrt(S2)); mut_.correctBoundaryConditions(); // Re-calculate thermal diffusivity diff --git a/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.H b/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.H index 6cc40f6b683e81e8955566e3803b5462597c9298..222242a74cb2fe26e851b53688f564a80c9affe5 100644 --- a/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.H +++ b/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.H @@ -38,7 +38,7 @@ Description Nov. 2001 \endverbatim - with the addition of the F3 term for rough walls from + with the addition of the optional F3 term for rough walls from \verbatim Hellsten, A. "Some Improvements in Menter’s k-omega-SST turbulence model" @@ -80,6 +80,7 @@ Description a1 0.31; b1 1.0; c1 10.0; + F3 no; } \endverbatim @@ -138,6 +139,8 @@ protected: dimensionedScalar b1_; dimensionedScalar c1_; + Switch F3_; + //- Wall distance // Note: different to wall distance in parent RASModel @@ -156,6 +159,7 @@ protected: tmp<volScalarField> F1(const volScalarField& CDkOmega) const; tmp<volScalarField> F2() const; tmp<volScalarField> F3() const; + tmp<volScalarField> F23() const; tmp<volScalarField> blend ( diff --git a/src/turbulenceModels/incompressible/RAS/kOmegaSST/kOmegaSST.C b/src/turbulenceModels/incompressible/RAS/kOmegaSST/kOmegaSST.C index cff0d48146077e69f3a7af4a18a303b671d33588..290a4c73a862c3b4cdb22f6ae76928472b444c53 100644 --- a/src/turbulenceModels/incompressible/RAS/kOmegaSST/kOmegaSST.C +++ b/src/turbulenceModels/incompressible/RAS/kOmegaSST/kOmegaSST.C @@ -98,6 +98,19 @@ tmp<volScalarField> kOmegaSST::F3() const } +tmp<volScalarField> kOmegaSST::F23() const +{ + tmp<volScalarField> f23(F2()); + + if (F3_) + { + f23() *= F3(); + } + + return f23; +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // kOmegaSST::kOmegaSST @@ -219,6 +232,15 @@ kOmegaSST::kOmegaSST 10.0 ) ), + F3_ + ( + Switch::lookupOrAddToDict + ( + "F3", + coeffDict_, + false + ) + ), y_(mesh_), @@ -268,7 +290,7 @@ kOmegaSST::kOmegaSST / max ( a1_*omega_, - b1_*F2()*F3()*sqrt(2.0)*mag(symm(fvc::grad(U_))) + b1_*F23()*sqrt(2.0)*mag(symm(fvc::grad(U_))) ) ); nut_.correctBoundaryConditions(); @@ -362,6 +384,7 @@ bool kOmegaSST::read() a1_.readIfPresent(coeffDict()); b1_.readIfPresent(coeffDict()); c1_.readIfPresent(coeffDict()); + F3_.readIfPresent("F3", coeffDict()); return true; } @@ -439,7 +462,7 @@ void kOmegaSST::correct() // Re-calculate viscosity - nut_ = a1_*k_/max(a1_*omega_, b1_*F2()*F3()*sqrt(S2)); + nut_ = a1_*k_/max(a1_*omega_, b1_*F23()*sqrt(S2)); nut_.correctBoundaryConditions(); } diff --git a/src/turbulenceModels/incompressible/RAS/kOmegaSST/kOmegaSST.H b/src/turbulenceModels/incompressible/RAS/kOmegaSST/kOmegaSST.H index 88feb733bc1814cd9479cbb38c2fe9dd7432495d..07d16d273210a68d882aeaaebef81b02d5b6683c 100644 --- a/src/turbulenceModels/incompressible/RAS/kOmegaSST/kOmegaSST.H +++ b/src/turbulenceModels/incompressible/RAS/kOmegaSST/kOmegaSST.H @@ -36,7 +36,7 @@ Description Nov. 2001. \endverbatim - with the addition of the F3 term for rough walls from + with the addition of the optional F3 term for rough walls from \verbatim Hellsten, A. "Some Improvements in Menter’s k-omega-SST turbulence model" @@ -77,6 +77,7 @@ Description a1 0.31; b1 1.0; c1 10.0; + F3 no; } \endverbatim @@ -132,6 +133,9 @@ protected: dimensionedScalar b1_; dimensionedScalar c1_; + Switch F3_; + + //- Wall distance field // Note: different to wall distance in parent RASModel wallDist y_; @@ -148,6 +152,7 @@ protected: tmp<volScalarField> F1(const volScalarField& CDkOmega) const; tmp<volScalarField> F2() const; tmp<volScalarField> F3() const; + tmp<volScalarField> F23() const; tmp<volScalarField> blend (