diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kOmegaSSTSato/kOmegaSSTSato.C b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kOmegaSSTSato/kOmegaSSTSato.C new file mode 100644 index 0000000000000000000000000000000000000000..3ab51aeb2607d1f4d7200773f7300719ea0a7e37 --- /dev/null +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kOmegaSSTSato/kOmegaSSTSato.C @@ -0,0 +1,171 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2014 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 "kOmegaSSTSato.H" +#include "addToRunTimeSelectionTable.H" +#include "twoPhaseSystem.H" +#include "dragModel.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace RASModels +{ + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template<class BasicTurbulenceModel> +kOmegaSSTSato<BasicTurbulenceModel>::kOmegaSSTSato +( + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName, + const word& type +) +: + kOmegaSST<BasicTurbulenceModel> + ( + alpha, + rho, + U, + alphaRhoPhi, + phi, + transport, + propertiesName, + type + ), + + gasTurbulencePtr_(NULL), + + Cmub_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "Cmub", + this->coeffDict_, + 0.6 + ) + ) +{ + if (type == typeName) + { + correctNut(); + this->printCoeffs(type); + } +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template<class BasicTurbulenceModel> +bool kOmegaSSTSato<BasicTurbulenceModel>::read() +{ + if (kOmegaSST<BasicTurbulenceModel>::read()) + { + Cmub_.readIfPresent(this->coeffDict()); + + return true; + } + else + { + return false; + } +} + +template<class BasicTurbulenceModel> +const PhaseCompressibleTurbulenceModel +< + typename BasicTurbulenceModel::transportModel +>& +kOmegaSSTSato<BasicTurbulenceModel>::gasTurbulence() const +{ + if (!gasTurbulencePtr_) + { + const volVectorField& U = this->U_; + + const transportModel& liquid = this->transport(); + const twoPhaseSystem& fluid = liquid.fluid(); + const transportModel& gas = fluid.otherPhase(liquid); + + gasTurbulencePtr_ = + &U.db() + .lookupObject<PhaseCompressibleTurbulenceModel<transportModel> > + ( + IOobject::groupName + ( + turbulenceModel::propertiesName, + gas.name() + ) + ); + } + + return *gasTurbulencePtr_; +} + + +template<class BasicTurbulenceModel> +void kOmegaSSTSato<BasicTurbulenceModel>::correctNut() +{ + const PhaseCompressibleTurbulenceModel<transportModel>& gasTurbulence = + this->gasTurbulence(); + + volScalarField yPlus + ( + pow(this->betaStar_, 0.25)*this->y_*sqrt(this->k_)/this->nu() + ); + + this->nut_ = + this->a1_*this->k_ + /max + ( + this->a1_*this->omega_, + this->F23()*sqrt(2.0)*mag(symm(fvc::grad(this->U_))) + ) + + sqr(1 - exp(-yPlus/16.0)) + *Cmub_*gasTurbulence.transport().d()*gasTurbulence.alpha() + *(mag(this->U_ - gasTurbulence.U())); + + this->nut_.correctBoundaryConditions(); +} + + +template<class BasicTurbulenceModel> +void kOmegaSSTSato<BasicTurbulenceModel>::correct() +{ + kOmegaSST<BasicTurbulenceModel>::correct(); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace Foam + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kOmegaSSTSato/kOmegaSSTSato.H b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kOmegaSSTSato/kOmegaSSTSato.H new file mode 100644 index 0000000000000000000000000000000000000000..f8118d284a02ce7b9211411d65c0fa9d41ed80a5 --- /dev/null +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kOmegaSSTSato/kOmegaSSTSato.H @@ -0,0 +1,217 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2014 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::RASModels::kOmegaSSTSato + +Group + grpRASTurbulence + +Description + Implementation of the k-omega-SST turbulence model for dispersed bubbly + flows with Sato (1981) bubble induced turbulent viscosity model. + + Bubble induced turbulent viscosity model described in: + \verbatim + Sato, Y., Sadatomi, M. + "Momentum and heat transfer in two-phase bubble flow - I, Theory" + International Journal of Multiphase FLow 7, pp. 167-177, 1981. + \endverbatim + + Turbulence model described in: + \verbatim + Menter, F., Esch, T. + "Elements of Industrial Heat Transfer Prediction" + 16th Brazilian Congress of Mechanical Engineering (COBEM), + Nov. 2001 + \endverbatim + + 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" + 29th AIAA Fluid Dynamics Conference, + AIAA-98-2554, + June 1998. + \endverbatim + + Note that this implementation is written in terms of alpha diffusion + coefficients rather than the more traditional sigma (alpha = 1/sigma) so + that the blending can be applied to all coefficuients in a consistent + manner. The paper suggests that sigma is blended but this would not be + consistent with the blending of the k-epsilon and k-omega models. + + Also note that the error in the last term of equation (2) relating to + sigma has been corrected. + + Wall-functions are applied in this implementation by using equations (14) + to specify the near-wall omega as appropriate. + + The blending functions (15) and (16) are not currently used because of the + uncertainty in their origin, range of applicability and that is y+ becomes + sufficiently small blending u_tau in this manner clearly becomes nonsense. + + The default model coefficients correspond to the following: + \verbatim + kOmegaSSTCoeffs + { + alphaK1 0.85034; + alphaK2 1.0; + alphaOmega1 0.5; + alphaOmega2 0.85616; + Prt 1.0; // only for compressible + beta1 0.075; + beta2 0.0828; + betaStar 0.09; + gamma1 0.5532; + gamma2 0.4403; + a1 0.31; + b1 1.0; + c1 10.0; + F3 no; + Cmub 0.6; + } + \endverbatim + +SourceFiles + kOmegaSST.C + +SourceFiles + kOmegaSSTSato.C + +\*---------------------------------------------------------------------------*/ + +#ifndef kOmegaSSTSato_H +#define kOmegaSSTSato_H + +#include "kOmegaSST.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace RASModels +{ + +/*---------------------------------------------------------------------------*\ + Class kOmegaSSTSato Declaration +\*---------------------------------------------------------------------------*/ + +template<class BasicTurbulenceModel> +class kOmegaSSTSato +: + public kOmegaSST<BasicTurbulenceModel> +{ + // Private data + + mutable const PhaseCompressibleTurbulenceModel + < + typename BasicTurbulenceModel::transportModel + > *gasTurbulencePtr_; + + + // Private Member Functions + + //- Return the turbulence model for the gas phase + const PhaseCompressibleTurbulenceModel + < + typename BasicTurbulenceModel::transportModel + >& + gasTurbulence() const; + + + // Disallow default bitwise copy construct and assignment + kOmegaSSTSato(const kOmegaSSTSato&); + kOmegaSSTSato& operator=(const kOmegaSSTSato&); + + +protected: + + // Protected data + + // Model coefficients + dimensionedScalar Cmub_; + + + // Protected Member Functions + + virtual void correctNut(); + +public: + + typedef typename BasicTurbulenceModel::alphaField alphaField; + typedef typename BasicTurbulenceModel::rhoField rhoField; + typedef typename BasicTurbulenceModel::transportModel transportModel; + + + //- Runtime type information + TypeName("kOmegaSSTSato"); + + + // Constructors + + //- Construct from components + kOmegaSSTSato + ( + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName = turbulenceModel::propertiesName, + const word& type = typeName + ); + + + //- Destructor + virtual ~kOmegaSSTSato() + {} + + + // Member Functions + + //- Read model coefficients if they have changed + virtual bool read(); + + //- Solve the turbulence equations and correct the turbulence viscosity + virtual void correct(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "kOmegaSSTSato.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/TurbulenceModels/compressible/compressibleTurbulenceModels.C b/src/TurbulenceModels/compressible/compressibleTurbulenceModels.C index 81ed6aeafccb134e510b34d9802f6ba0219f43db..054dee8884593e0e0bc0b35b6918a235e443662a 100644 --- a/src/TurbulenceModels/compressible/compressibleTurbulenceModels.C +++ b/src/TurbulenceModels/compressible/compressibleTurbulenceModels.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2013-2014 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -55,6 +55,9 @@ makeRASModel(kEpsilon); #include "buoyantKEpsilon.H" makeRASModel(buoyantKEpsilon); +#include "kOmegaSST.H" +makeRASModel(kOmegaSST); + #include "Smagorinsky.H" makeLESModel(Smagorinsky); diff --git a/src/TurbulenceModels/incompressible/incompressibleTurbulenceModels.C b/src/TurbulenceModels/incompressible/incompressibleTurbulenceModels.C index 50a813841e360f81304595ad6ad3f0e4cc6ce0ba..361c6f4b75b173d48b07e8d7409ed529ca8c054d 100644 --- a/src/TurbulenceModels/incompressible/incompressibleTurbulenceModels.C +++ b/src/TurbulenceModels/incompressible/incompressibleTurbulenceModels.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2013-2014 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -52,6 +52,9 @@ makeBaseTurbulenceModel #include "kEpsilon.H" makeRASModel(kEpsilon); +#include "kOmegaSST.H" +makeRASModel(kOmegaSST); + #include "Smagorinsky.H" makeLESModel(Smagorinsky); diff --git a/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSST/kOmegaSST.C b/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSST/kOmegaSST.C new file mode 100644 index 0000000000000000000000000000000000000000..e952c1f48d7d04e5987027bfb0fc829611331f67 --- /dev/null +++ b/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSST/kOmegaSST.C @@ -0,0 +1,464 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2014 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 "kOmegaSST.H" +#include "bound.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace RASModels +{ + +// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * // + +template<class BasicTurbulenceModel> +tmp<volScalarField> kOmegaSST<BasicTurbulenceModel>::kOmegaSST::F1 +( + const volScalarField& CDkOmega +) const +{ + tmp<volScalarField> CDkOmegaPlus = max + ( + CDkOmega, + dimensionedScalar("1.0e-10", dimless/sqr(dimTime), 1.0e-10) + ); + + tmp<volScalarField> arg1 = min + ( + min + ( + max + ( + (scalar(1)/betaStar_)*sqrt(k_)/(omega_*y_), + scalar(500)*(this->mu()/this->rho_)/(sqr(y_)*omega_) + ), + (4*alphaOmega2_)*k_/(CDkOmegaPlus*sqr(y_)) + ), + scalar(10) + ); + + return tanh(pow4(arg1)); +} + +template<class BasicTurbulenceModel> +tmp<volScalarField> kOmegaSST<BasicTurbulenceModel>::kOmegaSST::F2() const +{ + tmp<volScalarField> arg2 = min + ( + max + ( + (scalar(2)/betaStar_)*sqrt(k_)/(omega_*y_), + scalar(500)*(this->mu()/this->rho_)/(sqr(y_)*omega_) + ), + scalar(100) + ); + + return tanh(sqr(arg2)); +} + +template<class BasicTurbulenceModel> +tmp<volScalarField> kOmegaSST<BasicTurbulenceModel>::kOmegaSST::F3() const +{ + tmp<volScalarField> arg3 = min + ( + 150*(this->mu()/this->rho_)/(omega_*sqr(y_)), + scalar(10) + ); + + return 1 - tanh(pow4(arg3)); +} + +template<class BasicTurbulenceModel> +tmp<volScalarField> kOmegaSST<BasicTurbulenceModel>::kOmegaSST::F23() const +{ + tmp<volScalarField> f23(F2()); + + if (F3_) + { + f23() *= F3(); + } + + return f23; +} + + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template<class BasicTurbulenceModel> +kOmegaSST<BasicTurbulenceModel>::kOmegaSST +( + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName, + const word& type +) +: + eddyViscosity<RASModel<BasicTurbulenceModel> > + ( + type, + alpha, + rho, + U, + alphaRhoPhi, + phi, + transport, + propertiesName + ), + + alphaK1_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "alphaK1", + this->coeffDict_, + 0.85034 + ) + ), + alphaK2_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "alphaK2", + this->coeffDict_, + 1.0 + ) + ), + alphaOmega1_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "alphaOmega1", + this->coeffDict_, + 0.5 + ) + ), + alphaOmega2_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "alphaOmega2", + this->coeffDict_, + 0.85616 + ) + ), + Prt_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "Prt", + this->coeffDict_, + 1.0 + ) + ), + gamma1_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "gamma1", + this->coeffDict_, + 0.5532 + ) + ), + gamma2_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "gamma2", + this->coeffDict_, + 0.4403 + ) + ), + beta1_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "beta1", + this->coeffDict_, + 0.075 + ) + ), + beta2_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "beta2", + this->coeffDict_, + 0.0828 + ) + ), + betaStar_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "betaStar", + this->coeffDict_, + 0.09 + ) + ), + a1_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "a1", + this->coeffDict_, + 0.31 + ) + ), + b1_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "b1", + this->coeffDict_, + 1.0 + ) + ), + c1_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "c1", + this->coeffDict_, + 10.0 + ) + ), + F3_ + ( + Switch::lookupOrAddToDict + ( + "F3", + this->coeffDict_, + false + ) + ), + + y_(this->mesh_), + + k_ + ( + IOobject + ( + IOobject::groupName("k", U.group()), + this->runTime_.timeName(), + this->mesh_, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + this->mesh_ + ), + omega_ + ( + IOobject + ( + IOobject::groupName("omega", U.group()), + this->runTime_.timeName(), + this->mesh_, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + this->mesh_ + ) +{ + + bound(k_, this->kMin_); + bound(omega_, this->omegaMin_); + + this->nut_ = + ( + a1_*k_ + /max + ( + a1_*omega_, + b1_*F23()*sqrt(2.0)*mag(symm(fvc::grad(this->U_))) + ) + ); + this->nut_.correctBoundaryConditions(); +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template<class BasicTurbulenceModel> +bool kOmegaSST<BasicTurbulenceModel>::read() +{ + if (eddyViscosity<RASModel<BasicTurbulenceModel> >::read()) + { + alphaK1_.readIfPresent(this->coeffDict()); + alphaK2_.readIfPresent(this->coeffDict()); + alphaOmega1_.readIfPresent(this->coeffDict()); + alphaOmega2_.readIfPresent(this->coeffDict()); + beta1_.readIfPresent(this->coeffDict()); + beta2_.readIfPresent(this->coeffDict()); + betaStar_.readIfPresent(this->coeffDict()); + gamma1_.readIfPresent(this->coeffDict()); + gamma2_.readIfPresent(this->coeffDict()); + a1_.readIfPresent(this->coeffDict()); + b1_.readIfPresent(this->coeffDict()); + c1_.readIfPresent(this->coeffDict()); + F3_.readIfPresent("F3", this->coeffDict()); + + return true; + } + else + { + return false; + } +} + + +template<class BasicTurbulenceModel> +void kOmegaSST<BasicTurbulenceModel>::correctNut() +{ + this->nut_ = + a1_*k_/max(a1_*omega_, F23()*sqrt(2.0)*mag(symm(fvc::grad(this->U_)))); + this->nut_.correctBoundaryConditions(); +} + + +template<class BasicTurbulenceModel> +tmp<fvScalarMatrix> kOmegaSST<BasicTurbulenceModel>::kSource() const +{ + return tmp<fvScalarMatrix> + ( + new fvScalarMatrix + ( + k_, + dimVolume*this->rho_.dimensions()*k_.dimensions()/dimTime + ) + ); +} + + +template<class BasicTurbulenceModel> +tmp<fvScalarMatrix> kOmegaSST<BasicTurbulenceModel>::omegaSource() const +{ + return tmp<fvScalarMatrix> + ( + new fvScalarMatrix + ( + omega_, + dimVolume*this->rho_.dimensions()*omega_.dimensions()/dimTime + ) + ); +} + + +template<class BasicTurbulenceModel> +void kOmegaSST<BasicTurbulenceModel>::correct() +{ + if (!this->turbulence_) + { + return; + } + + // Local references + const alphaField& alpha = this->alpha_; + const rhoField& rho = this->rho_; + const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_; + const volVectorField& U = this->U_; + volScalarField& nut = this->nut_; + + eddyViscosity<RASModel<BasicTurbulenceModel> >::correct(); + + volScalarField divU(fvc::div(fvc::absolute(this->phi(), U))); + + tmp<volTensorField> tgradU = fvc::grad(U); + volScalarField S2(2*magSqr(symm(tgradU()))); + volScalarField GbyMu((tgradU() && dev(twoSymm(tgradU())))); + volScalarField G(this->GName(), nut*GbyMu); + tgradU.clear(); + + // Update omega and G at the wall + omega_.boundaryField().updateCoeffs(); + + volScalarField CDkOmega + ( + (2*alphaOmega2_)*(fvc::grad(k_) & fvc::grad(omega_))/omega_ + ); + + volScalarField F1(this->F1(CDkOmega)); + volScalarField rhoGammaF1(rho*gamma(F1)); + + // Turbulent frequency equation + tmp<fvScalarMatrix> omegaEqn + ( + fvm::ddt(alpha, rho, omega_) + + fvm::div(alphaRhoPhi, omega_) + - fvm::Sp(fvc::ddt(alpha, rho) + fvc::div(alphaRhoPhi), omega_) + - fvm::laplacian(alpha*rho*DomegaEff(F1), omega_) + == + alpha*rhoGammaF1*GbyMu + - fvm::SuSp((2.0/3.0)*alpha*rhoGammaF1*divU, omega_) + - fvm::Sp(alpha*rho*beta(F1)*omega_, omega_) + - fvm::SuSp + ( + alpha*rho*(F1 - scalar(1))*CDkOmega/omega_, + omega_ + ) + + omegaSource() + ); + + omegaEqn().relax(); + + omegaEqn().boundaryManipulate(omega_.boundaryField()); + + solve(omegaEqn); + bound(omega_, this->omegaMin_); + + // Turbulent kinetic energy equation + tmp<fvScalarMatrix> kEqn + ( + fvm::ddt(alpha, rho, k_) + + fvm::div(alphaRhoPhi, k_) + - fvm::Sp(fvc::ddt(alpha, rho) + fvc::div(alphaRhoPhi), k_) + - fvm::laplacian(alpha*rho*DkEff(F1), k_) + == + min(alpha*rho*G, (c1_*betaStar_)*alpha*rho*k_*omega_) + - fvm::SuSp((2.0/3.0)*alpha*rho*divU, k_) + - fvm::Sp(alpha*rho*betaStar_*omega_, k_) + + kSource() + ); + + kEqn().relax(); + solve(kEqn); + bound(k_, this->kMin_); + + correctNut(); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSST/kOmegaSST.H b/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSST/kOmegaSST.H new file mode 100644 index 0000000000000000000000000000000000000000..85dee85ffe0f64f23760502ad1b1e58bf978a6fe --- /dev/null +++ b/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSST/kOmegaSST.H @@ -0,0 +1,314 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2014 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::RASModels::kOmegaSST + +Group + grpRASTurbulence + +Description + Implementation of the k-omega-SST turbulence model for compressible flows. + + Turbulence model described in: + \verbatim + Menter, F., Esch, T. + "Elements of Industrial Heat Transfer Prediction" + 16th Brazilian Congress of Mechanical Engineering (COBEM), + Nov. 2001 + \endverbatim + + 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" + 29th AIAA Fluid Dynamics Conference, + AIAA-98-2554, + June 1998. + \endverbatim + + Note that this implementation is written in terms of alpha diffusion + coefficients rather than the more traditional sigma (alpha = 1/sigma) so + that the blending can be applied to all coefficuients in a consistent + manner. The paper suggests that sigma is blended but this would not be + consistent with the blending of the k-epsilon and k-omega models. + + Also note that the error in the last term of equation (2) relating to + sigma has been corrected. + + Wall-functions are applied in this implementation by using equations (14) + to specify the near-wall omega as appropriate. + + The blending functions (15) and (16) are not currently used because of the + uncertainty in their origin, range of applicability and that if y+ becomes + sufficiently small blending u_tau in this manner clearly becomes nonsense. + + The default model coefficients correspond to the following: + \verbatim + kOmegaSSTCoeffs + { + alphaK1 0.85034; + alphaK2 1.0; + alphaOmega1 0.5; + alphaOmega2 0.85616; + Prt 1.0; // only for compressible + beta1 0.075; + beta2 0.0828; + betaStar 0.09; + gamma1 0.5532; + gamma2 0.4403; + a1 0.31; + b1 1.0; + c1 10.0; + F3 no; + } + \endverbatim + +SourceFiles + kOmegaSST.C + +\*---------------------------------------------------------------------------*/ + +#ifndef kOmegaSST_H +#define kOmegaSST_H + +#include "RASModel.H" +#include "eddyViscosity.H" +#include "wallDist.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace RASModels +{ + +/*---------------------------------------------------------------------------*\ + Class kOmegaSST Declaration +\*---------------------------------------------------------------------------*/ + +template<class BasicTurbulenceModel> +class kOmegaSST +: + public eddyViscosity<RASModel<BasicTurbulenceModel> > +{ + // Private Member Functions + + // Disallow default bitwise copy construct and assignment + kOmegaSST(const kOmegaSST&); + kOmegaSST& operator=(const kOmegaSST&); + + +protected: + + // Protected data + + // Model coefficients + + dimensionedScalar alphaK1_; + dimensionedScalar alphaK2_; + + dimensionedScalar alphaOmega1_; + dimensionedScalar alphaOmega2_; + + dimensionedScalar Prt_; + + dimensionedScalar gamma1_; + dimensionedScalar gamma2_; + + dimensionedScalar beta1_; + dimensionedScalar beta2_; + + dimensionedScalar betaStar_; + + dimensionedScalar a1_; + dimensionedScalar b1_; + dimensionedScalar c1_; + + Switch F3_; + + + //- Wall distance + // Note: different to wall distance in parent RASModel + wallDist y_; + + // Fields + + volScalarField k_; + volScalarField omega_; + + + // Protected Member Functions + + virtual void correctNut(); + virtual tmp<fvScalarMatrix> kSource() const; + virtual tmp<fvScalarMatrix> omegaSource() const; + + + // Private Member Functions + + tmp<volScalarField> F1(const volScalarField& CDkOmega) const; + tmp<volScalarField> F2() const; + tmp<volScalarField> F3() const; + tmp<volScalarField> F23() const; + + tmp<volScalarField> blend + ( + const volScalarField& F1, + const dimensionedScalar& psi1, + const dimensionedScalar& psi2 + ) const + { + return F1*(psi1 - psi2) + psi2; + } + + tmp<volScalarField> alphaK(const volScalarField& F1) const + { + return blend(F1, alphaK1_, alphaK2_); + } + + tmp<volScalarField> alphaOmega(const volScalarField& F1) const + { + return blend(F1, alphaOmega1_, alphaOmega2_); + } + + tmp<volScalarField> beta(const volScalarField& F1) const + { + return blend(F1, beta1_, beta2_); + } + + tmp<volScalarField> gamma(const volScalarField& F1) const + { + return blend(F1, gamma1_, gamma2_); + } + + +public: + + typedef typename BasicTurbulenceModel::alphaField alphaField; + typedef typename BasicTurbulenceModel::rhoField rhoField; + typedef typename BasicTurbulenceModel::transportModel transportModel; + + + //- Runtime type information + TypeName("kOmegaSST"); + + + // Constructors + + //- Construct from components + kOmegaSST + ( + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName = turbulenceModel::propertiesName, + const word& type = typeName + ); + + + //- Destructor + virtual ~kOmegaSST() + {} + + + // Member Functions + + //- Re-read model coefficients if they have changed + virtual bool read(); + + //- Return the effective diffusivity for k + tmp<volScalarField> DkEff(const volScalarField& F1) const + { + return tmp<volScalarField> + ( + new volScalarField("DkEff", alphaK(F1)*this->nut_ + this->nu()) + ); + } + + //- Return the effective diffusivity for omega + tmp<volScalarField> DomegaEff(const volScalarField& F1) const + { + return tmp<volScalarField> + ( + new volScalarField + ( + "DomegaEff", + alphaOmega(F1)*this->nut_ + this->nu() + ) + ); + } + + //- Return the turbulence kinetic energy + virtual tmp<volScalarField> k() const + { + return k_; + } + + //- Return the turbulence kinetic energy dissipation rate + virtual tmp<volScalarField> epsilon() const + { + return tmp<volScalarField> + ( + new volScalarField + ( + IOobject + ( + "epsilon", + this->mesh_.time().timeName(), + this->mesh_ + ), + betaStar_*k_*omega_, + omega_.boundaryField().types() + ) + ); + } + + //- Return the turbulence kinetic energy dissipation rate + virtual tmp<volScalarField> omega() const + { + return omega_; + } + + //- Solve the turbulence equations and correct the turbulence viscosity + virtual void correct(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // +#ifdef NoRepository +# include "kOmegaSST.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // +#endif + +// ************************************************************************* // diff --git a/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.C b/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.C index e1e69bfe5eb277f5ed731a8ca7b2dac833f849ea..999676d4f118d1b6fb8755951d591c8aa9a9d497 100644 --- a/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.C +++ b/src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.C @@ -410,7 +410,7 @@ void kOmegaSST::correct() // Re-calculate viscosity mut_ = a1_*rho_*k_ - /max(a1_*omega_, F2()*sqrt(2.0)*mag(symm(fvc::grad(U_)))); + /max(a1_*omega_, F23()*sqrt(2.0)*mag(symm(fvc::grad(U_)))); mut_.correctBoundaryConditions(); // Re-calculate thermal diffusivity