From 3f9271f2bed34597b8ee39550a3485f96eacb576 Mon Sep 17 00:00:00 2001 From: Vuko Vukcevic Date: Mon, 29 May 2017 13:58:48 +0200 Subject: [PATCH] Additional wall functions: omega/nutMEWTWallFunctions Modified enhanced wall treatment by Sutalo, wall functions are sensitive to convection and pressure gradient effects. Author: Filip Sutalo, Merge: Vuko Vukcevic --- .../incompressible/RAS/Make/files | 2 + .../nutCWTWallFunctionFvPatchScalarField.H | 2 +- .../nutMEWTWallFunctionFvPatchScalarField.C | 384 ++++++++++++++++++ .../nutMEWTWallFunctionFvPatchScalarField.H | 221 ++++++++++ .../omegaCWTWallFunctionFvPatchScalarField.C | 6 +- .../omegaCWTWallFunctionFvPatchScalarField.H | 2 +- .../omegaMEWTWallFunctionFvPatchScalarField.C | 345 ++++++++++++++++ .../omegaMEWTWallFunctionFvPatchScalarField.H | 205 ++++++++++ 8 files changed, 1163 insertions(+), 4 deletions(-) create mode 100644 src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutMEWTWallFunction/nutMEWTWallFunctionFvPatchScalarField.C create mode 100644 src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutMEWTWallFunction/nutMEWTWallFunctionFvPatchScalarField.H create mode 100644 src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaMEWTWallFunction/omegaMEWTWallFunctionFvPatchScalarField.C create mode 100644 src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaMEWTWallFunction/omegaMEWTWallFunctionFvPatchScalarField.H diff --git a/src/turbulenceModels/incompressible/RAS/Make/files b/src/turbulenceModels/incompressible/RAS/Make/files index b0f1a0944..5d6c61b2c 100644 --- a/src/turbulenceModels/incompressible/RAS/Make/files +++ b/src/turbulenceModels/incompressible/RAS/Make/files @@ -33,6 +33,7 @@ $(nutWallFunctions)/nutSpalartAllmarasStandardWallFunction/nutSpalartAllmarasSta $(nutWallFunctions)/nutSpalartAllmarasStandardRoughWallFunction/nutSpalartAllmarasStandardRoughWallFunctionFvPatchScalarField.C $(nutWallFunctions)/nutLowReWallFunction/nutLowReWallFunctionFvPatchScalarField.C $(nutWallFunctions)/nutCWTWallFunction/nutCWTWallFunctionFvPatchScalarField.C +$(nutWallFunctions)/nutMEWTWallFunction/nutMEWTWallFunctionFvPatchScalarField.C epsilonWallFunctions = $(wallFunctions)/epsilonWallFunctions $(epsilonWallFunctions)/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C @@ -40,6 +41,7 @@ $(epsilonWallFunctions)/epsilonWallFunction/epsilonWallFunctionFvPatchScalarFiel omegaWallFunctions = $(wallFunctions)/omegaWallFunctions $(omegaWallFunctions)/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C $(omegaWallFunctions)/omegaCWTWallFunction/omegaCWTWallFunctionFvPatchScalarField.C +$(omegaWallFunctions)/omegaMEWTWallFunction/omegaMEWTWallFunctionFvPatchScalarField.C kqRWallFunctions = $(wallFunctions)/kqRWallFunctions $(kqRWallFunctions)/kqRWallFunction/kqRWallFunctionFvPatchFields.C diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutCWTWallFunction/nutCWTWallFunctionFvPatchScalarField.H b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutCWTWallFunction/nutCWTWallFunctionFvPatchScalarField.H index 14880b673..a4eb35ccc 100644 --- a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutCWTWallFunction/nutCWTWallFunctionFvPatchScalarField.H +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutCWTWallFunction/nutCWTWallFunctionFvPatchScalarField.H @@ -118,7 +118,7 @@ protected: public: //- Runtime type information - TypeName("nutCWTWallTreatment"); + TypeName("nutCWTWallFunction"); // Constructors diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutMEWTWallFunction/nutMEWTWallFunctionFvPatchScalarField.C b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutMEWTWallFunction/nutMEWTWallFunctionFvPatchScalarField.C new file mode 100644 index 000000000..5fe8aa937 --- /dev/null +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutMEWTWallFunction/nutMEWTWallFunctionFvPatchScalarField.C @@ -0,0 +1,384 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 4.0 + \\ / 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 "nutMEWTWallFunctionFvPatchScalarField.H" +#include "RASModel.H" +#include "fvPatchFieldMapper.H" +#include "volFields.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +const Foam::debug::tolerancesSwitch +Foam::incompressible::RASModels:: +nutMEWTWallFunctionFvPatchScalarField::dimlessAFactorTol_ +( + "dimlessAFactorMEWTTolerance", + 1e-6 +); + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace incompressible +{ +namespace RASModels +{ + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +void nutMEWTWallFunctionFvPatchScalarField::checkType() +{ + if (!patch().isWall()) + { + FatalErrorIn("nutMEWTWallFunctionFvPatchScalarField::checkType()") + << "Invalid wall function specification" << nl + << " Patch type for patch " << patch().name() + << " must be wall" << nl + << " Current patch type is " << patch().type() << nl << endl + << abort(FatalError); + } +} + + +scalar nutMEWTWallFunctionFvPatchScalarField::calcYPlusLam +( + const scalar kappa, + const scalar E +) const +{ + scalar ypl = 11.0; + + for (int i = 0; i < 10; i++) + { + ypl = log(E*ypl)/kappa; + } + + return ypl; +} + + +tmp nutMEWTWallFunctionFvPatchScalarField::calcNut() const +{ + const label patchI = patch().index(); + + const RASModel& rasModel = db().lookupObject("RASProperties"); + const scalarField& y = rasModel.y()[patchI]; + const tmp tk = rasModel.k(); + const volScalarField& k = tk(); + const scalarField& nuw = rasModel.nu().boundaryField()[patchI]; + + const scalar Cmu25 = pow(Cmu_, 0.25); + + // Get normals + const vectorField n = patch().nf(); + + // Patch velocity field at this wall + const fvPatchVectorField& Uw = + lookupPatchField(UName_); + + const scalarField magGradUw = mag(Uw.snGrad()); + const vectorField UwIn = Uw.patchInternalField(); + + // Patch internal velocity field tangential to the wall + const vectorField UwInTang = UwIn - (UwIn & n)*n; + const scalarField magUwInTang = mag(UwInTang); + + // Calculate tangential direction for patch cells + const vectorField tDir = UwInTang/magUwInTang; + + // Wall-velocity vector field tangential to the wall + const vectorField UwTang = Uw - (Uw & n)*n; + const scalarField magUwTang = mag(UwTang); + + + // Pressure terms + const volScalarField& p = + this->dimensionedInternalField().mesh().lookupObject + < + volScalarField + >(pName_); + + // Pressure gradient + const volVectorField gradp = fvc::grad(p); + + // Pressure gradient in wall adjacent cell + const vectorField gradPIn = + gradp.boundaryField()[this->patch().index()].patchInternalField(); + + // Pressure gradient projected on the wall parallel velocity + const scalarField gradpTang= gradPIn & tDir; + + + // Convective terms + const volVectorField& U = + this->dimensionedInternalField().mesh().lookupObject + < + volVectorField + >(UName_); + + const surfaceScalarField& phi = + this->dimensionedInternalField().mesh().lookupObject + < + surfaceScalarField + >("phi"); + + const volVectorField convection = fvc::div(phi, U); + + const vectorField convectionIn = + convection.boundaryField()[this->patch().index()].patchInternalField(); + + // Convection term projected on the wall parallel velocity + const scalarField convectionTang = convectionIn & tDir; + + // Needed to calculate yPlus + const scalarField& eddyVis = + lookupPatchField(nutName_); + + tmp tnutw(new scalarField(patch().size(), SMALL)); + scalarField& nutw = tnutw(); + + // Get face cells + const unallocLabelList& fc = patch().faceCells(); + + forAll(nutw, faceI) + { + const label faceCellI = fc[faceI]; + const scalar uStar = Cmu25*sqrt(k[faceCellI]); + + // Calculate yPlus + const scalar yPlus = + sqrt + ( + (eddyVis[faceI]+nuw[faceI])*magGradUw[faceI] + )* + y[faceI]/ + (nuw[faceI] + SMALL); + + // Relative tangential velocity + const scalar magUrel = magUwInTang[faceI] - magUwTang[faceI]; + + // Dimless A factor + scalar A = nuw[faceI]* + (gradpTang[faceI] + convectionTang[faceI])/(pow(uStar, 3) + SMALL); + + // Numerical stabilisation of the A factor + if (A < SMALL) + { + A += dimlessAFactorTol_(); + } + + // Helper variables + const scalar S1 = sqrt(max(SMALL, 1.0 + A*yPlus)); + const scalar p1 = sqrt(max(SMALL, 1.0 +6.0*A)); + const scalar uPlusT = + (1.0/kappa_)*log(6.0*E_) + -(1.0/kappa_)*(2.0*p1 + log(mag(p1 - 1.0)) + log(p1 + 1.0)); + + // Dimless velocity in log layer + const scalar uLogPlus = + (1.0/kappa_)*(2.0*S1 + log(mag(S1 - 1.0)) + log(S1 + 1.0)) + uPlusT; + + // Friction velocity in viscous sublayer + const scalar uTauVis = sqrt(nuw[faceI]*magUrel/(y[faceI] + SMALL)); + + // Friction velocity in log layer + const scalar uTauLog = magUrel/(uLogPlus + SMALL); + + // Kader blending for friction velocity + const scalar gamma = -0.01*pow(yPlus, 4)/(1.0 + 5.0*yPlus); + const scalar uTau = uTauVis*exp(gamma) + uTauLog*exp(1.0/gamma); + + // Need to limit nutw for stability reasons since at some point, yPlus + // is 0 and Kader blending is not defined + nutw[faceI] = + max + ( + SMALL, + sqr(uTau)/(magGradUw[faceI] + SMALL) - nuw[faceI] + ); + } + + return tnutw; +} + + +void nutMEWTWallFunctionFvPatchScalarField::writeLocalEntries(Ostream& os) const +{ + writeEntryIfDifferent(os, "U", "U", UName_); + writeEntryIfDifferent(os, "p", "p", pName_); + writeEntryIfDifferent(os, "nut", "nut", nutName_); + os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl; + os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl; + os.writeKeyword("E") << E_ << token::END_STATEMENT << nl; +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +nutMEWTWallFunctionFvPatchScalarField::nutMEWTWallFunctionFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF +) +: + fixedValueFvPatchScalarField(p, iF), + UName_("U"), + pName_("p"), + nutName_("nut"), + Cmu_(0.09), + kappa_(0.41), + E_(9.8), + yPlusLam_(calcYPlusLam(kappa_, E_)) +{ + checkType(); +} + + +nutMEWTWallFunctionFvPatchScalarField::nutMEWTWallFunctionFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF, + const dictionary& dict +) +: + fixedValueFvPatchScalarField(p, iF, dict), + UName_(dict.lookupOrDefault("U", "U")), + pName_(dict.lookupOrDefault("p", "p")), + nutName_(dict.lookupOrDefault("nut", "nut")), + Cmu_(dict.lookupOrDefault("Cmu", 0.09)), + kappa_(dict.lookupOrDefault("kappa", 0.41)), + E_(dict.lookupOrDefault("E", 9.8)), + yPlusLam_(calcYPlusLam(kappa_, E_)) +{ + checkType(); +} + + +nutMEWTWallFunctionFvPatchScalarField::nutMEWTWallFunctionFvPatchScalarField +( + const nutMEWTWallFunctionFvPatchScalarField& ptf, + const fvPatch& p, + const DimensionedField& iF, + const fvPatchFieldMapper& mapper +) +: + fixedValueFvPatchScalarField(ptf, p, iF, mapper), + UName_(ptf.UName_), + pName_(ptf.pName_), + nutName_(ptf.nutName_), + Cmu_(ptf.Cmu_), + kappa_(ptf.kappa_), + E_(ptf.E_), + yPlusLam_(ptf.yPlusLam_) +{ + checkType(); +} + + +nutMEWTWallFunctionFvPatchScalarField::nutMEWTWallFunctionFvPatchScalarField +( + const nutMEWTWallFunctionFvPatchScalarField& wfpsf +) +: + fixedValueFvPatchScalarField(wfpsf), + UName_(wfpsf.UName_), + pName_(wfpsf.pName_), + nutName_(wfpsf.nutName_), + Cmu_(wfpsf.Cmu_), + kappa_(wfpsf.kappa_), + E_(wfpsf.E_), + yPlusLam_(wfpsf.yPlusLam_) +{ + checkType(); +} + + +nutMEWTWallFunctionFvPatchScalarField::nutMEWTWallFunctionFvPatchScalarField +( + const nutMEWTWallFunctionFvPatchScalarField& wfpsf, + const DimensionedField& iF +) +: + fixedValueFvPatchScalarField(wfpsf, iF), + UName_(wfpsf.UName_), + pName_(wfpsf.pName_), + nutName_(wfpsf.nutName_), + Cmu_(wfpsf.Cmu_), + kappa_(wfpsf.kappa_), + E_(wfpsf.E_), + yPlusLam_(wfpsf.yPlusLam_) +{ + checkType(); +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void nutMEWTWallFunctionFvPatchScalarField::updateCoeffs() +{ + operator==(calcNut()); + + fixedValueFvPatchScalarField::updateCoeffs(); +} + + +tmp nutMEWTWallFunctionFvPatchScalarField::yPlus() const +{ + const label patchI = patch().index(); + + const RASModel& rasModel = db().lookupObject("RASProperties"); + const scalarField& y = rasModel.y()[patchI]; + + const tmp tk = rasModel.k(); + const volScalarField& k = tk(); + const scalarField kwc = k.boundaryField()[patchI].patchInternalField(); + const scalarField& nuw = rasModel.nu().boundaryField()[patchI]; + + return pow(Cmu_, 0.25)*y*sqrt(kwc)/nuw; +} + + +void nutMEWTWallFunctionFvPatchScalarField::write(Ostream& os) const +{ + fvPatchField::write(os); + writeLocalEntries(os); + writeEntry("value", os); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +makePatchTypeField(fvPatchScalarField, nutMEWTWallFunctionFvPatchScalarField); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace incompressible +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutMEWTWallFunction/nutMEWTWallFunctionFvPatchScalarField.H b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutMEWTWallFunction/nutMEWTWallFunctionFvPatchScalarField.H new file mode 100644 index 000000000..12bb3a15c --- /dev/null +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutMEWTWallFunction/nutMEWTWallFunctionFvPatchScalarField.H @@ -0,0 +1,221 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 3.2 + \\ / 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::nutMEWTWallFunctionFvPatchScalarField + +Description + Improved boundary condition for turbulent (kinematic) viscosity, taking into + account convection and pressure gradient terms. + + Reference (bibtex entry): + + @article{sutaloThesis2017, + author = {\v{S}utalo, F.} + title = {{Evaluation of Variants of Enhanced Wall Treatment Wall + Function in Turbulent Flow Simulations}}, + school = {Faculty of Mechanical Engineering and Naval Architecture, + University of Zagreb}, + year = {2017}, + } + +SourceFiles + nutMEWTWallFunctionFvPatchScalarField.C + +Author + Filip Sutalo, FMENA Zagreb. All rights reserved. + +\*---------------------------------------------------------------------------*/ + +#ifndef nutMEWTWallFunctionFvPatchScalarField_H +#define nutMEWTWallFunctionFvPatchScalarField_H + +#include "fixedValueFvPatchFields.H" +#include "tolerancesSwitch.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace incompressible +{ +namespace RASModels +{ + +/*---------------------------------------------------------------------------*\ + Class nutMEWTWallFunctionFvPatchScalarField Declaration +\*---------------------------------------------------------------------------*/ + +class nutMEWTWallFunctionFvPatchScalarField +: + public fixedValueFvPatchScalarField +{ + + // Private static data member + + //- Tolerance for the A dimensionless parameter + static const debug::tolerancesSwitch dimlessAFactorTol_; + + +protected: + + // Protected data + + //- Name of velocity field + word UName_; + + //- Name of pressure field + word pName_; + + //- Name of eddy viscosity field + word nutName_; + + //- Name of omega field + word omegaName_; + + //- Cmu coefficient + scalar Cmu_; + + //- Von Karman constant + scalar kappa_; + + //- E coefficient + scalar E_; + + //- Y+ at the edge of the laminar sublayer + scalar yPlusLam_; + + + // Protected member functions + + //- Check the type of the patch + virtual void checkType(); + + //- Calculate the Y+ at the edge of the laminar sublayer + virtual scalar calcYPlusLam(const scalar kappa, const scalar E) const; + + //- Calculate the turbulence viscosity + virtual tmp calcNut() const; + + //- Write local wall function variables + virtual void writeLocalEntries(Ostream&) const; + + +public: + + //- Runtime type information + TypeName("nutMEWTWallFunction"); + + + // Constructors + + //- Construct from patch and internal field + nutMEWTWallFunctionFvPatchScalarField + ( + const fvPatch&, + const DimensionedField& + ); + + //- Construct from patch, internal field and dictionary + nutMEWTWallFunctionFvPatchScalarField + ( + const fvPatch&, + const DimensionedField&, + const dictionary& + ); + + //- Construct by mapping given + // nutMEWTWallFunctionFvPatchScalarField + // onto a new patch + nutMEWTWallFunctionFvPatchScalarField + ( + const nutMEWTWallFunctionFvPatchScalarField&, + const fvPatch&, + const DimensionedField&, + const fvPatchFieldMapper& + ); + + //- Construct as copy + nutMEWTWallFunctionFvPatchScalarField + ( + const nutMEWTWallFunctionFvPatchScalarField& + ); + + //- Construct and return a clone + virtual tmp clone() const + { + return tmp + ( + new nutMEWTWallFunctionFvPatchScalarField(*this) + ); + } + + //- Construct as copy setting internal field reference + nutMEWTWallFunctionFvPatchScalarField + ( + const nutMEWTWallFunctionFvPatchScalarField&, + const DimensionedField& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp clone + ( + const DimensionedField& iF + ) const + { + return tmp + ( + new nutMEWTWallFunctionFvPatchScalarField(*this, iF) + ); + } + + + // Member functions + + // Evaluation functions + + //- Calculate and return the yPlus at the boundary + virtual tmp yPlus() const; + + //- Update the coefficients associated with the patch field + virtual void updateCoeffs(); + + + // I-O + + //- Write + virtual void write(Ostream&) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace incompressible +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaCWTWallFunction/omegaCWTWallFunctionFvPatchScalarField.C b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaCWTWallFunction/omegaCWTWallFunctionFvPatchScalarField.C index 9cacb7eb0..d48f8a34d 100644 --- a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaCWTWallFunction/omegaCWTWallFunctionFvPatchScalarField.C +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaCWTWallFunction/omegaCWTWallFunctionFvPatchScalarField.C @@ -131,7 +131,7 @@ omegaCWTWallFunctionFvPatchScalarField::omegaCWTWallFunctionFvPatchScalarField ) : fixedInternalValueFvPatchField(owfpsf), - pName_(owfpsf.pName_), + pName_(owfpsf.pName_), UName_(owfpsf.UName_), kName_(owfpsf.kName_), GName_(owfpsf.GName_), @@ -232,6 +232,7 @@ void omegaCWTWallFunctionFvPatchScalarField::updateCoeffs() // Magnitude of the surface normal gradient velocity const scalarField magGradUw = mag(Uw.snGrad()); + // Pressure effects const volScalarField& p = this->dimensionedInternalField().mesh().lookupObject @@ -244,9 +245,10 @@ void omegaCWTWallFunctionFvPatchScalarField::updateCoeffs() const vectorField gradPIn = gradp.boundaryField()[this->patch().index()].patchInternalField(); - // Pressure gradient projected on the wall parallel velocity + // Pressure gradient projected on the wall parallel velocity direction const scalarField gradpTang= gradPIn & tDir; + // Convective terms const volVectorField& U = this->dimensionedInternalField().mesh().lookupObject diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaCWTWallFunction/omegaCWTWallFunctionFvPatchScalarField.H b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaCWTWallFunction/omegaCWTWallFunctionFvPatchScalarField.H index 47bca6c6c..85b823b7d 100644 --- a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaCWTWallFunction/omegaCWTWallFunctionFvPatchScalarField.H +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaCWTWallFunction/omegaCWTWallFunctionFvPatchScalarField.H @@ -113,7 +113,7 @@ class omegaCWTWallFunctionFvPatchScalarField public: //- Runtime type information - TypeName("omegaCWTWallTreatment"); + TypeName("omegaCWTWallFunction"); // Constructors diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaMEWTWallFunction/omegaMEWTWallFunctionFvPatchScalarField.C b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaMEWTWallFunction/omegaMEWTWallFunctionFvPatchScalarField.C new file mode 100644 index 000000000..932a5a079 --- /dev/null +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaMEWTWallFunction/omegaMEWTWallFunctionFvPatchScalarField.C @@ -0,0 +1,345 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 3.2 + \\ / 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 "omegaMEWTWallFunctionFvPatchScalarField.H" +#include "RASModel.H" +#include "fvPatchFieldMapper.H" +#include "volFields.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace incompressible +{ +namespace RASModels +{ + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void omegaMEWTWallFunctionFvPatchScalarField::checkType() +{ + if (!patch().isWall()) + { + FatalErrorIn("omegaMEWTWallFunctionFvPatchScalarField::checkType()") + << "Invalid wall function specification" << nl + << " Patch type for patch " << patch().name() + << " must be wall" << nl + << " Current patch type is " << patch().type() << nl << endl + << abort(FatalError); + } +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +omegaMEWTWallFunctionFvPatchScalarField::omegaMEWTWallFunctionFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF +) +: + fixedInternalValueFvPatchField(p, iF), + pName_("p"), + UName_("U"), + kName_("k"), + GName_("RASModel::G"), + nuName_("nu"), + nutName_("nut"), + Cmu_(0.09), + kappa_(0.41), + E_(9.8), + beta1_(0.075) +{ + checkType(); +} + + +omegaMEWTWallFunctionFvPatchScalarField::omegaMEWTWallFunctionFvPatchScalarField +( + const omegaMEWTWallFunctionFvPatchScalarField& ptf, + const fvPatch& p, + const DimensionedField& iF, + const fvPatchFieldMapper& mapper +) +: + fixedInternalValueFvPatchField(ptf, p, iF, mapper), + pName_(ptf.pName_), + UName_(ptf.UName_), + kName_(ptf.kName_), + GName_(ptf.GName_), + nuName_(ptf.nuName_), + nutName_(ptf.nutName_), + Cmu_(ptf.Cmu_), + kappa_(ptf.kappa_), + E_(ptf.E_), + beta1_(ptf.beta1_) +{ + checkType(); +} + + +omegaMEWTWallFunctionFvPatchScalarField::omegaMEWTWallFunctionFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF, + const dictionary& dict +) +: + fixedInternalValueFvPatchField(p, iF, dict), + pName_(dict.lookupOrDefault("p", "p")), + UName_(dict.lookupOrDefault("U", "U")), + kName_(dict.lookupOrDefault("k", "k")), + GName_(dict.lookupOrDefault("G", "RASModel::G")), + nuName_(dict.lookupOrDefault("nu", "nu")), + nutName_(dict.lookupOrDefault("nut", "nut")), + Cmu_(dict.lookupOrDefault("Cmu", 0.09)), + kappa_(dict.lookupOrDefault("kappa", 0.41)), + E_(dict.lookupOrDefault("E", 9.8)), + beta1_(dict.lookupOrDefault("beta1", 0.075)) +{ + checkType(); +} + + +omegaMEWTWallFunctionFvPatchScalarField::omegaMEWTWallFunctionFvPatchScalarField +( + const omegaMEWTWallFunctionFvPatchScalarField& owfpsf +) +: + fixedInternalValueFvPatchField(owfpsf), + pName_(owfpsf.pName_), + UName_(owfpsf.UName_), + kName_(owfpsf.kName_), + GName_(owfpsf.GName_), + nuName_(owfpsf.nuName_), + nutName_(owfpsf.nutName_), + Cmu_(owfpsf.Cmu_), + kappa_(owfpsf.kappa_), + E_(owfpsf.E_), + beta1_(owfpsf.beta1_) +{ + checkType(); +} + + +omegaMEWTWallFunctionFvPatchScalarField::omegaMEWTWallFunctionFvPatchScalarField +( + const omegaMEWTWallFunctionFvPatchScalarField& owfpsf, + const DimensionedField& iF +) +: + fixedInternalValueFvPatchField(owfpsf, iF), + pName_(owfpsf.pName_), + UName_(owfpsf.UName_), + kName_(owfpsf.kName_), + GName_(owfpsf.GName_), + nuName_(owfpsf.nuName_), + nutName_(owfpsf.nutName_), + Cmu_(owfpsf.Cmu_), + kappa_(owfpsf.kappa_), + E_(owfpsf.E_), + beta1_(owfpsf.beta1_) +{ + checkType(); +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void omegaMEWTWallFunctionFvPatchScalarField::updateCoeffs() +{ + if (updated()) + { + return; + } + + // If G field is not present, execute zero gradient evaluation + // HJ, 20/Mar/2011 + if (!db().foundObject(GName_)) + { + InfoIn("void omegaMEWTWallFunctionFvPatchScalarField::updateCoeffs()") + << "Cannot access " << GName_ << " field for patch " + << patch().name() << ". Evaluating as zeroGradient" + << endl; + + fvPatchScalarField::updateCoeffs(); + zeroGradientFvPatchScalarField::evaluate(); + + return; + } + + const RASModel& rasModel = db().lookupObject("RASProperties"); + const scalarField& y = rasModel.y()[patch().index()]; + + volScalarField& G = const_cast + (db().lookupObject(GName_)); + + // Note: omega is now a refValue and set in fixedInternalValueFvPatchField + // HJ, 3/Aug/2011 + scalarField& omega = refValue(); + + const scalarField& nuw = + lookupPatchField(nuName_); + + const scalarField& nutw = + lookupPatchField(nutName_); + + // Velocity field at the patch + const fvPatchVectorField& Uw = + lookupPatchField(UName_); + + // Patch normals + const vectorField n = patch().nf(); + + // Velocity patch internal field + const vectorField UwIn = Uw.patchInternalField(); + + // Velocity vector tangential to the wall + const vectorField UwInTang = UwIn - (UwIn & n)*n; + const scalarField magUwInTang = mag(UwInTang); + + // Calculate tangential direction for patch cells + const vectorField tDir = UwInTang/magUwInTang; + + // Magnitude of the surface normal gradient velocity + const scalarField magGradUw = mag(Uw.snGrad()); + + + // Pressure effects + const volScalarField& p = + this->dimensionedInternalField().mesh().lookupObject + (pName_); + + // Pressure gradient + const volVectorField gradp = fvc::grad(p); + + // Pressure gradient in wall adjacent cell + const vectorField gradPIn = + gradp.boundaryField()[this->patch().index()].patchInternalField(); + + // Pressure gradient projected on the wall parallel velocity direction + const scalarField gradpTang= gradPIn & tDir; + + + // Convective terms + const volVectorField& U = + this->dimensionedInternalField().mesh().lookupObject + < + volVectorField + >(UName_); + + const surfaceScalarField& phi = + this->dimensionedInternalField().mesh().lookupObject + < + surfaceScalarField + >("phi"); + + const volVectorField convection = fvc::div(phi, U); + + const vectorField convectionIn = + convection.boundaryField()[this->patch().index()].patchInternalField(); + + // Convection term projected on the wall parallel velocity + const scalarField convectionTang = convectionIn & tDir; + + // Get face cells + const unallocLabelList& fc = patch().faceCells(); + + // Set omega and G + forAll(nutw, faceI) + { + const label faceCellI = fc[faceI]; + + // Tangential stress at the wall + const scalar tauw = (nuw[faceI] + nutw[faceI])*magGradUw[faceI]; + + const scalar yPlus = sqrt(tauw)*y[faceI]/nuw[faceI]; + + // Velocity gradient for viscous sublayer + const scalar dudyVis= magGradUw[faceI]; + + // Velocity gradient for log layer + const scalar dudyLog = + sqrt + ( + max + ( + SMALL, + (gradpTang[faceI] + convectionTang[faceI])*y[faceI] +tauw + ) + )/(kappa_*y[faceI]); + + // Kader blending for velocity gradient + const scalar gamma = -0.01*pow(yPlus, 4)/(1.0 + 5.0*yPlus); + const scalar dudy = dudyVis*exp(gamma) + dudyLog*exp(1.0/gamma); + + const scalar omegaVis = 6.0*nuw[faceI]/(beta1_*sqr(y[faceI])); + const scalar omegaLog = dudyLog/(sqrt(Cmu_)); + + // Menter blend for omega + omega[faceI] = sqrt(sqr(omegaVis) + sqr(omegaLog)); + + // Production term + G[faceCellI]= tauw*dudy; + } + + // TODO: perform averaging for cells sharing more than one boundary face + + fixedInternalValueFvPatchField::updateCoeffs(); +} + +void omegaMEWTWallFunctionFvPatchScalarField::write(Ostream& os) const +{ + fixedInternalValueFvPatchField::write(os); + writeEntryIfDifferent(os, "U", "U", UName_); + writeEntryIfDifferent(os, "p", "p", pName_); + writeEntryIfDifferent(os, "k", "k", kName_); + writeEntryIfDifferent(os, "G", "RASModel::G", GName_); + writeEntryIfDifferent(os, "nu", "nu", nuName_); + writeEntryIfDifferent(os, "nut", "nut", nutName_); + os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl; + os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl; + os.writeKeyword("E") << E_ << token::END_STATEMENT << nl; + os.writeKeyword("beta1") << beta1_ << token::END_STATEMENT << nl; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +makePatchTypeField +( + fvPatchScalarField, + omegaMEWTWallFunctionFvPatchScalarField +); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace incompressible +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaMEWTWallFunction/omegaMEWTWallFunctionFvPatchScalarField.H b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaMEWTWallFunction/omegaMEWTWallFunctionFvPatchScalarField.H new file mode 100644 index 000000000..edee0bd10 --- /dev/null +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaMEWTWallFunction/omegaMEWTWallFunctionFvPatchScalarField.H @@ -0,0 +1,205 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 3.2 + \\ / 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::omegaMEWTWallFunctionFvPatchScalarField + +Description + Improved boundary condition for specific dissipation, taking into account + convection and pressure gradient terms. + + Reference (bibtex entry): + + @article{sutaloThesis2017, + author = {\v{S}utalo, F.} + title = {{Evaluation of Variants of Enhanced Wall Treatment Wall + Function in Turbulent Flow Simulations}}, + school = {Faculty of Mechanical Engineering and Naval Architecture, + University of Zagreb}, + year = {2017}, + } + +SourceFiles + omegaMEWTWallFunctionFvPatchScalarField.C + +Author + Filip Sutalo, FMENA Zagreb. All rights reserved. + +\*---------------------------------------------------------------------------*/ + +#ifndef omegaMEWTWallFunctionFvPatchScalarField_H +#define omegaMEWTWallFunctionFvPatchScalarField_H + +#include "fixedInternalValueFvPatchFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace incompressible +{ +namespace RASModels +{ + +/*---------------------------------------------------------------------------*\ + Class omegaMEWTWallFunctionFvPatchScalarField Declaration +\*---------------------------------------------------------------------------*/ + +class omegaMEWTWallFunctionFvPatchScalarField +: + public fixedInternalValueFvPatchScalarField +{ + // Private data + + //- Name of pressure field + word pName_; + + //- Name of velocity field + word UName_; + + //- Name of turbulence kinetic energy field + word kName_; + + //- Name of turbulence generation field + word GName_; + + //- Name of laminar viscosity field + word nuName_; + + //- Name of turbulent viscosity field + word nutName_; + + //- Cmu coefficient + scalar Cmu_; + + //- Von Karman constant + scalar kappa_; + + //- E coefficient + scalar E_; + + //- beta1 coefficient + scalar beta1_; + + + // Private member functions + + //- Check the type of the patch + void checkType(); + + +public: + + //- Runtime type information + TypeName("omegaMEWTWallFunction"); + + + // Constructors + + //- Construct from patch and internal field + omegaMEWTWallFunctionFvPatchScalarField + ( + const fvPatch&, + const DimensionedField& + ); + + //- Construct from patch, internal field and dictionary + omegaMEWTWallFunctionFvPatchScalarField + ( + const fvPatch&, + const DimensionedField&, + const dictionary& + ); + + //- Construct by mapping given + // omegaMEWTWallFunctionFvPatchScalarField + // onto a new patch + omegaMEWTWallFunctionFvPatchScalarField + ( + const omegaMEWTWallFunctionFvPatchScalarField&, + const fvPatch&, + const DimensionedField&, + const fvPatchFieldMapper& + ); + + //- Construct as copy + omegaMEWTWallFunctionFvPatchScalarField + ( + const omegaMEWTWallFunctionFvPatchScalarField& + ); + + //- Construct and return a clone + virtual tmp clone() const + { + return tmp + ( + new omegaMEWTWallFunctionFvPatchScalarField(*this) + ); + } + + //- Construct as copy setting internal field reference + omegaMEWTWallFunctionFvPatchScalarField + ( + const omegaMEWTWallFunctionFvPatchScalarField&, + const DimensionedField& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp clone + ( + const DimensionedField& iF + ) const + { + return tmp + ( + new omegaMEWTWallFunctionFvPatchScalarField(*this, iF) + ); + } + + + // Member functions + + // Evaluation functions + + //- Update the coefficients associated with the patch field + virtual void updateCoeffs(); + + + // I-O + + //- Write + void write(Ostream&) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace incompressible +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* //