diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutURoughWallFunction/nutURoughWallFunctionFvPatchScalarField.C b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutURoughWallFunction/nutURoughWallFunctionFvPatchScalarField.C new file mode 100644 index 000000000..15af4c0b5 --- /dev/null +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutURoughWallFunction/nutURoughWallFunctionFvPatchScalarField.C @@ -0,0 +1,308 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 4.1 + \\ / A nd | Web: http://www.foam-extend.org + \\/ M anipulation | For copyright notice see file Copyright +------------------------------------------------------------------------------- +License + This file is part of foam-extend. + + foam-extend 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. + + foam-extend 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 foam-extend. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "nutURoughWallFunctionFvPatchScalarField.H" +#include "RASModel.H" +#include "fvPatchFieldMapper.H" +#include "volFields.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace incompressible +{ +namespace RASModels +{ + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +tmp nutURoughWallFunctionFvPatchScalarField::calcNut() const +{ + const label patchI = patch().index(); + + const turbulenceModel& turbModel = + db().lookupObject("turbulenceModel"); + + const scalarField& y = turbModel.y()[patchI]; + const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchI]; + const scalarField& nuw = turbModel.nu().boundaryField()[patchI]; + + // The flow velocity at the adjacent cell centre + const scalarField magUp(mag(Uw.patchInternalField() - Uw)); + + tmp tyPlus = calcYPlus(magUp); + const scalarField& yPlus = tyPlus(); + + tmp tnutw(new scalarField(patch().size(), 0.0)); + scalarField& nutw = tnutw(); + + forAll(yPlus, faceI) + { + if (yPlus[faceI] > yPlusLam_) + { + const scalar Re = magUp[faceI]*y[faceI]/nuw[faceI] + ROOTVSMALL; + nutw[faceI] = nuw[faceI]*(sqr(yPlus[faceI])/Re - 1); + } + } + + return tnutw; +} + + +tmp nutURoughWallFunctionFvPatchScalarField::calcYPlus +( + const scalarField& magUp +) const +{ + const label patchI = patch().index(); + + const turbulenceModel& turbModel = + db().lookupObject("turbulenceModel"); + + const scalarField& y = turbModel.y()[patchI]; + const scalarField& nuw = turbModel.nu().boundaryField()[patchI]; + + tmp tyPlus(new scalarField(patch().size(), 0.0)); + scalarField& yPlus = tyPlus(); + + if (roughnessHeight_ > 0.0) + { + // Rough Walls + const scalar c_1 = 1/(90 - 2.25) + roughnessConstant_; + static const scalar c_2 = 2.25/(90 - 2.25); + static const scalar c_3 = 2.0*atan(1.0)/log(90/2.25); + static const scalar c_4 = c_3*log(2.25); + + //if (KsPlusBasedOnYPlus_) + { + // If KsPlus is based on YPlus the extra term added to the law + // of the wall will depend on yPlus + forAll(yPlus, faceI) + { + const scalar magUpara = magUp[faceI]; + const scalar Re = magUpara*y[faceI]/nuw[faceI]; + const scalar kappaRe = kappa_*Re; + + scalar yp = yPlusLam_; + const scalar ryPlusLam = 1.0/yp; + + label iter = 0; + scalar yPlusLast = 0.0; + scalar dKsPlusdYPlus = roughnessHeight_/y[faceI]; + + // Additional tuning parameter - nominally = 1 + dKsPlusdYPlus *= roughnessFactor_; + + do + { + yPlusLast = yp; + + // The non-dimensional roughness height + scalar KsPlus = yp*dKsPlusdYPlus; + + // The extra term in the law-of-the-wall + scalar G = 0.0; + + scalar yPlusGPrime = 0.0; + + if (KsPlus >= 90) + { + const scalar t_1 = 1 + roughnessConstant_*KsPlus; + G = log(t_1); + yPlusGPrime = roughnessConstant_*KsPlus/t_1; + } + else if (KsPlus > 2.25) + { + const scalar t_1 = c_1*KsPlus - c_2; + const scalar t_2 = c_3*log(KsPlus) - c_4; + const scalar sint_2 = sin(t_2); + const scalar logt_1 = log(t_1); + G = logt_1*sint_2; + yPlusGPrime = + (c_1*sint_2*KsPlus/t_1) + (c_3*logt_1*cos(t_2)); + } + + scalar denom = 1.0 + log(E_*yp) - G - yPlusGPrime; + if (mag(denom) > VSMALL) + { + yp = (kappaRe + yp*(1 - yPlusGPrime))/denom; + } + } while + ( + mag(ryPlusLam*(yp - yPlusLast)) > 0.0001 + && ++iter < 10 + && yp > VSMALL + ); + + yPlus[faceI] = max(0.0, yp); + } + } + } + else + { + // Smooth Walls + forAll(yPlus, faceI) + { + const scalar magUpara = magUp[faceI]; + const scalar Re = magUpara*y[faceI]/nuw[faceI]; + const scalar kappaRe = kappa_*Re; + + scalar yp = yPlusLam_; + const scalar ryPlusLam = 1.0/yp; + + label iter = 0; + scalar yPlusLast = 0.0; + + do + { + yPlusLast = yp; + yp = (kappaRe + yp)/(1.0 + log(E_*yp)); + + } while (mag(ryPlusLam*(yp - yPlusLast)) > 0.0001 && ++iter < 10); + + yPlus[faceI] = max(0.0, yp); + } + } + + return tyPlus; +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +nutURoughWallFunctionFvPatchScalarField::nutURoughWallFunctionFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF +) +: + nutWallFunctionFvPatchScalarField(p, iF), + roughnessHeight_(0.0), + roughnessConstant_(0.0), + roughnessFactor_(0.0) +{} + + +nutURoughWallFunctionFvPatchScalarField::nutURoughWallFunctionFvPatchScalarField +( + const nutURoughWallFunctionFvPatchScalarField& ptf, + const fvPatch& p, + const DimensionedField& iF, + const fvPatchFieldMapper& mapper +) +: + nutWallFunctionFvPatchScalarField(ptf, p, iF, mapper), + roughnessHeight_(ptf.roughnessHeight_), + roughnessConstant_(ptf.roughnessConstant_), + roughnessFactor_(ptf.roughnessFactor_) +{} + + +nutURoughWallFunctionFvPatchScalarField::nutURoughWallFunctionFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF, + const dictionary& dict +) +: + nutWallFunctionFvPatchScalarField(p, iF, dict), + roughnessHeight_(readScalar(dict.lookup("roughnessHeight"))), + roughnessConstant_(readScalar(dict.lookup("roughnessConstant"))), + roughnessFactor_(readScalar(dict.lookup("roughnessFactor"))) +{} + + +nutURoughWallFunctionFvPatchScalarField::nutURoughWallFunctionFvPatchScalarField +( + const nutURoughWallFunctionFvPatchScalarField& rwfpsf +) +: + nutWallFunctionFvPatchScalarField(rwfpsf), + roughnessHeight_(rwfpsf.roughnessHeight_), + roughnessConstant_(rwfpsf.roughnessConstant_), + roughnessFactor_(rwfpsf.roughnessFactor_) +{} + + +nutURoughWallFunctionFvPatchScalarField::nutURoughWallFunctionFvPatchScalarField +( + const nutURoughWallFunctionFvPatchScalarField& rwfpsf, + const DimensionedField& iF +) +: + nutWallFunctionFvPatchScalarField(rwfpsf, iF), + roughnessHeight_(rwfpsf.roughnessHeight_), + roughnessConstant_(rwfpsf.roughnessConstant_), + roughnessFactor_(rwfpsf.roughnessFactor_) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +tmp nutURoughWallFunctionFvPatchScalarField::yPlus() const +{ + const label patchI = patch().index(); + + const turbulenceModel& turbModel = + db().lookupObject("turbulenceModel"); + + const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchI]; + tmp magUp = mag(Uw.patchInternalField() - Uw); + + return calcYPlus(magUp()); +} + + +void nutURoughWallFunctionFvPatchScalarField::write(Ostream& os) const +{ + fvPatchField::write(os); + writeLocalEntries(os); + os.writeKeyword("roughnessHeight") + << roughnessHeight_ << token::END_STATEMENT << nl; + os.writeKeyword("roughnessConstant") + << roughnessConstant_ << token::END_STATEMENT << nl; + os.writeKeyword("roughnessFactor") + << roughnessFactor_ << token::END_STATEMENT << nl; + writeEntry("value", os); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +makePatchTypeField +( + fvPatchScalarField, + nutURoughWallFunctionFvPatchScalarField +); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace incompressible +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutURoughWallFunction/nutURoughWallFunctionFvPatchScalarField.H b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutURoughWallFunction/nutURoughWallFunctionFvPatchScalarField.H new file mode 100644 index 000000000..c308a2d6b --- /dev/null +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutURoughWallFunction/nutURoughWallFunctionFvPatchScalarField.H @@ -0,0 +1,240 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 4.1 + \\ / A nd | Web: http://www.foam-extend.org + \\/ M anipulation | For copyright notice see file Copyright +------------------------------------------------------------------------------- +License + This file is part of foam-extend. + + foam-extend 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. + + foam-extend 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 foam-extend. If not, see . + +Class + Foam::incompressible::RASModels::nutURoughWallFunctionFvPatchScalarField + +Description + This boundary condition provides a turbulent kinematic viscosity condition + when using wall functions for rough walls, based on velocity. + +Usage + \table + Property | Description | Required | Default value + roughnessHeight | roughness height | yes | + roughnessConstant | roughness constanr | yes | + roughnessFactor | scaling factor | yes | + \endtable + + Example of the boundary condition specification: + \verbatim + + { + type nutURoughWallFunction; + roughnessHeight 1e-5; + roughnessConstant 0.5; + roughnessFactor 1; + } + \endverbatim + +See also + Foam::nutWallFunctionFvPatchScalarField + +SourceFiles + nutURoughWallFunctionFvPatchScalarField.C + +\*---------------------------------------------------------------------------*/ + +#ifndef nutURoughWallFunctionFvPatchScalarField_H +#define nutURoughWallFunctionFvPatchScalarField_H + +#include "nutWallFunctionFvPatchScalarField.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace incompressible +{ +namespace RASModels +{ + +/*---------------------------------------------------------------------------*\ + Class nutURoughWallFunctionFvPatchScalarField Declaration +\*---------------------------------------------------------------------------*/ + +class nutURoughWallFunctionFvPatchScalarField +: + public nutWallFunctionFvPatchScalarField +{ + // Private data + + // Roughness model parameters + + //- Height + scalar roughnessHeight_; + + //- Constant + scalar roughnessConstant_; + + //- Scale factor + scalar roughnessFactor_; + + + // Protected Member Functions + + //- Calculate yPLus + virtual tmp calcYPlus(const scalarField& magUp) const; + + //- Calculate the turbulence viscosity + virtual tmp calcNut() const; + + +public: + + //- Runtime type information + TypeName("nutURoughWallFunction"); + + + // Constructors + + //- Construct from patch and internal field + nutURoughWallFunctionFvPatchScalarField + ( + const fvPatch&, + const DimensionedField& + ); + + //- Construct from patch, internal field and dictionary + nutURoughWallFunctionFvPatchScalarField + ( + const fvPatch&, + const DimensionedField&, + const dictionary& + ); + + //- Construct by mapping given + // nutURoughWallFunctionFvPatchScalarField + // onto a new patch + nutURoughWallFunctionFvPatchScalarField + ( + const nutURoughWallFunctionFvPatchScalarField&, + const fvPatch&, + const DimensionedField&, + const fvPatchFieldMapper& + ); + + //- Construct as copy + nutURoughWallFunctionFvPatchScalarField + ( + const nutURoughWallFunctionFvPatchScalarField& + ); + + //- Construct and return a clone + virtual tmp clone() const + { + return tmp + ( + new nutURoughWallFunctionFvPatchScalarField(*this) + ); + } + + //- Construct as copy setting internal field reference + nutURoughWallFunctionFvPatchScalarField + ( + const nutURoughWallFunctionFvPatchScalarField&, + const DimensionedField& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp clone + ( + const DimensionedField& iF + ) const + { + return tmp + ( + new nutURoughWallFunctionFvPatchScalarField(*this, iF) + ); + } + + + // Member functions + + // Access + + //- Return the roughness height + scalar roughnessHeight() const + { + return roughnessHeight_; + } + + //- Return reference to the roughness height to allow adjustment + scalar& roughnessHeight() + { + return roughnessHeight_; + } + + + //- Return the roughness constant scale + scalar roughnessConstant() const + { + return roughnessConstant_; + } + + //- Return reference to the roughness constant to allow adjustment + scalar& roughnessConstant() + { + return roughnessConstant_; + } + + //- Return the roughness scale factor + scalar roughnessFactor() const + { + return roughnessFactor_; + } + + //- Return reference to the roughness scale factor to allow + // adjustment + scalar& roughnessFactor() + { + return roughnessFactor_; + } + + + // I-O + + // Evaluation functions + + //- Calculate and return the yPlus at the boundary + virtual tmp yPlus() const; + + + // I-O + + //- Write + virtual void write(Ostream& os) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace incompressible +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* //