exaFOAM/foam-extend4.1-coherent-io
Archived
1
1
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity

New version of block-coupled k-omega SST model. Robert Keser

Browse source
This commit is contained in:
Hrvoje Jasak 2016-03-16 17:54:03 +00:00
parent 2d5052d714
commit e2f568cd84
2 changed files with 826 additions and 0 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

552
src/turbulenceModels/incompressible/RAS/blockKOmegaSST/blockKOmegaSST.C Normal file
View file

@ -0,0 +1,552 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "blockKOmegaSST.H"
#include "fvBlockMatrix.H"
#include "addToRunTimeSelectionTable.H"
#include "backwardsCompatibilityWallFunctions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace incompressible
{
namespace RASModels
{
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
defineTypeNameAndDebug(blockKOmegaSST, 0);
addToRunTimeSelectionTable(RASModel, blockKOmegaSST, dictionary);
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
tmp<volScalarField> blockKOmegaSST::F1(const volScalarField& CDkOmega) const
{
volScalarField CDkOmegaPlus = max
(
CDkOmega,
dimensionedScalar("1.0e-10", dimless/sqr(dimTime), 1.0e-10)
);
volScalarField arg1 = min
(
min
(
max
(
(scalar(1)/betaStar_)*sqrt(k_)/(omega_*y_),
scalar(500)*nu()/(sqr(y_)*omega_)
),
(4*alphaOmega2_)*k_/(CDkOmegaPlus*sqr(y_))
),
scalar(10)
);
return tanh(pow4(arg1));
}
tmp<volScalarField> blockKOmegaSST::F2() const
{
volScalarField arg2 = min
(
max
(
(scalar(2)/betaStar_)*sqrt(k_)/(omega_*y_),
scalar(500)*nu()/(sqr(y_)*omega_)
),
scalar(100)
);
return tanh(sqr(arg2));
}
tmp<volScalarField> blockKOmegaSST::F3() const
{
tmp<volScalarField> arg3 = min
(
150*nu()/(omega_*sqr(y_)),
scalar(10)
);
return 1 - tanh(pow4(arg3));
}
tmp<volScalarField> blockKOmegaSST::F23() const
{
tmp<volScalarField> f23(F2());
if (F3_)
{
f23() *= F3();
}
return f23;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
blockKOmegaSST::blockKOmegaSST
(
const volVectorField& U,
const surfaceScalarField& phi,
transportModel& lamTransportModel
)
:
RASModel(typeName, U, phi, lamTransportModel),
alphaK1_
(
dimensionedScalar::lookupOrAddToDict
(
"alphaK1",
coeffDict_,
0.85034
)
),
alphaK2_
(
dimensionedScalar::lookupOrAddToDict
(
"alphaK2",
coeffDict_,
1.0
)
),
alphaOmega1_
(
dimensionedScalar::lookupOrAddToDict
(
"alphaOmega1",
coeffDict_,
0.5
)
),
alphaOmega2_
(
dimensionedScalar::lookupOrAddToDict
(
"alphaOmega2",
coeffDict_,
0.856
)
),
gamma1_
(
dimensionedScalar::lookupOrAddToDict
(
"gamma1",
coeffDict_,
5.0/9.0
)
),
gamma2_
(
dimensionedScalar::lookupOrAddToDict
(
"gamma2",
coeffDict_,
0.44
)
),
beta1_
(
dimensionedScalar::lookupOrAddToDict
(
"beta1",
coeffDict_,
0.075
)
),
beta2_
(
dimensionedScalar::lookupOrAddToDict
(
"beta2",
coeffDict_,
0.0828
)
),
betaStar_
(
dimensionedScalar::lookupOrAddToDict
(
"betaStar",
coeffDict_,
0.09
)
),
a1_
(
dimensionedScalar::lookupOrAddToDict
(
"a1",
coeffDict_,
0.31
)
),
b1_
(
dimensionedScalar::lookupOrAddToDict
(
"b1",
coeffDict_,
1.0
)
),
c1_
(
dimensionedScalar::lookupOrAddToDict
(
"c1",
coeffDict_,
10.0
)
),
F3_
(
Switch::lookupOrAddToDict
(
"F3",
coeffDict_,
false
)
),
y_(mesh_),
k_
(
IOobject
(
"k",
runTime_.timeName(),
U_.db(),
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
autoCreateK("k", mesh_, U_.db())
),
omega_
(
IOobject
(
"omega",
runTime_.timeName(),
U_.db(),
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
autoCreateOmega("omega", mesh_, U_.db())
),
nut_
(
IOobject
(
"nut",
runTime_.timeName(),
U_.db(),
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
autoCreateNut("nut", mesh_, U_.db())
),
kOmega_
(
IOobject
(
"kOmega",
runTime_.timeName(),
U_.db(),
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh_,
dimensionedVector2("zero", dimless, vector2::zero)
)
{
bound(k_, k0_);
bound(omega_, omega0_);
nut_ =
(
a1_*k_/
max
(
a1_*omega_,
b1_*F23()*sqrt(2.0)*mag(symm(fvc::grad(U_)))
)
);
nut_.correctBoundaryConditions();
printCoeffs();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
tmp<volSymmTensorField> blockKOmegaSST::R() const
{
return tmp<volSymmTensorField>
(
new volSymmTensorField
(
IOobject
(
"R",
runTime_.timeName(),
U_.db(),
IOobject::NO_READ,
IOobject::NO_WRITE
),
((2.0/3.0)*I)*k_ - nut_*twoSymm(fvc::grad(U_)),
k_.boundaryField().types()
)
);
}
tmp<volSymmTensorField> blockKOmegaSST::devReff() const
{
return tmp<volSymmTensorField>
(
new volSymmTensorField
(
IOobject
(
"devRhoReff",
runTime_.timeName(),
U_.db(),
IOobject::NO_READ,
IOobject::NO_WRITE
),
-nuEff()*dev(twoSymm(fvc::grad(U_)))
)
);
}
tmp<fvVectorMatrix> blockKOmegaSST::divDevReff(volVectorField& U) const
{
return
(
- fvm::laplacian(nuEff(), U)
- fvc::div(nuEff()*dev(fvc::grad(U)().T()))
);
}
bool blockKOmegaSST::read()
{
if (RASModel::read())
{
alphaK1_.readIfPresent(coeffDict());
alphaK2_.readIfPresent(coeffDict());
alphaOmega1_.readIfPresent(coeffDict());
alphaOmega2_.readIfPresent(coeffDict());
gamma1_.readIfPresent(coeffDict());
gamma2_.readIfPresent(coeffDict());
beta1_.readIfPresent(coeffDict());
beta2_.readIfPresent(coeffDict());
betaStar_.readIfPresent(coeffDict());
a1_.readIfPresent(coeffDict());
b1_.readIfPresent(coeffDict());
c1_.readIfPresent(coeffDict());
F3_.readIfPresent("F3", coeffDict());
return true;
}
else
{
return false;
}
}
void blockKOmegaSST::correct()
{
// Bound in case of topological change
// HJ, 22/Aug/2007
if (mesh_.changing())
{
bound(k_, k0_);
bound(omega_, omega0_);
}
RASModel::correct();
if (!turbulence_)
{
return;
}
if (mesh_.changing())
{
y_.correct();
}
const volScalarField S2(2*magSqr(symm(fvc::grad(U_))));
volScalarField G("RASModel::G", nut_*S2);
// Make coupled matrix
fvBlockMatrix<vector2> kOmegaEqn(kOmega_);
// Update omega and G at the wall
omega_.boundaryField().updateCoeffs();
const volScalarField CDkOmega
(
(2*alphaOmega2_)*(fvc::grad(k_) & fvc::grad(omega_))/omega_
);
const volScalarField F1(this->F1(CDkOmega));
// Record coupling equations when omega is set. Same equations
// will be eliminated from k coupling as well due to wall functions
// HJ, 13/Nov/2015
labelList eliminated;
// Turbulent frequency equation
{
fvScalarMatrix omegaEqn
(
fvm::ddt(omega_)
+ fvm::div(phi_, omega_)
+ fvm::SuSp(-fvc::div(phi_), omega_)
- fvm::laplacian(DomegaEff(F1), omega_)
==
gamma(F1)
*min
(
S2,
(c1_/a1_)*betaStar_*omega_*max(a1_*omega_, b1_*F23()*sqrt(S2))
)
- fvm::Sp(2*beta(F1)*omega_, omega_)
+ beta(F1)*sqr(omega_)
- fvm::SuSp
(
(F1 - scalar(1))*CDkOmega/omega_,
omega_
)
);
omegaEqn.relax();
omegaEqn.completeAssembly();
kOmegaEqn.insertEquation(1, omegaEqn);
// Add coupling term
volScalarField couplingOmega
(
"couplingOmega",
-pos((scalar(1) - F1)*CDkOmega)*(scalar(1) - F1)*CDkOmega/k_
);
scalarField& couplingOmegaIn = couplingOmega.internalField();
// Collect cell indices where wall functions are active from the
// omega wall functions
eliminated = omegaEqn.eliminatedEqns().toc();
// Eliminate coupling in wall function cells for omega
forAll (eliminated, cellI)
{
couplingOmegaIn[eliminated[cellI]] *= 0;
}
// Insert coupling
kOmegaEqn.insertEquationCoupling(1, 0, couplingOmega);
}
// Turbulent kinetic energy equation
{
fvScalarMatrix kEqn
(
fvm::ddt(k_)
+ fvm::div(phi_, k_)
+ fvm::SuSp(-fvc::div(phi_), k_)
- fvm::laplacian(DkEff(F1), k_)
==
min(G, c1_*betaStar_*k_*omega_)
- fvm::Sp(betaStar_*omega_, k_)
);
kEqn.relax();
kEqn.completeAssembly();
kOmegaEqn.insertEquation(0, kEqn);
// Add coupling term
volScalarField couplingK
(
"couplingK",
-pos(G - c1_*betaStar_*k_*omega_)*c1_*betaStar_*k_
);
scalarField& couplingKIn = couplingK.internalField();
// Eliminate coupling in wall function cells for k
forAll (eliminated, cellI)
{
couplingKIn[eliminated[cellI]] *= 0;
}
// Insert coupling
kOmegaEqn.insertEquationCoupling(0, 1, couplingK);
}
// Update source coupling: coupling terms eliminated from source
kOmegaEqn.updateSourceCoupling();
kOmegaEqn.solve();
// Retrieve solution
kOmegaEqn.retrieveSolution(0, k_.internalField());
kOmegaEqn.retrieveSolution(1, omega_.internalField());
bound(k_, k0_);
bound(omega_, omega0_);
k_.correctBoundaryConditions();
omega_.correctBoundaryConditions();
// Re-calculate viscosity
// Fixed sqrt(2) error. HJ, 10/Jun/2015
nut_ = a1_*k_/max(a1_*omega_, b1_*F23()*sqrt(S2));
nut_ = min(nut_, nuRatio()*nu());
nut_.correctBoundaryConditions();
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace RASModels
} // End namespace incompressible
} // End namespace Foam
// ************************************************************************* //

274
src/turbulenceModels/incompressible/RAS/blockKOmegaSST/blockKOmegaSST.H Normal file
View file

@ -0,0 +1,274 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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 <http://www.gnu.org/licenses/>.
Class
Foam::incompressible::RASModels::blockKOmegaSST
Description
Implementation of the k-omega-SST turbulence model for incompressible
flows using a block-coupled solution.
The default model coefficients correspond to the following:
@verbatim
blockKOmegaSSTCoeffs
{
alphaK1 0.85;
alphaK2 1.0;
alphaOmega1 0.5;
alphaOmega2 0.856;
beta1 0.075;
beta2 0.0828;
betaStar 0.09;
gamma1 5/9;
gamma2 0.44;
a1 0.31;
b1 1.0;
c1 10.0;
F3 no;
}
@endverbatim
SourceFiles
blockKOmegaSST.C
\*---------------------------------------------------------------------------*/
#ifndef blockKOmegaSST_H
#define blockKOmegaSST_H
#include "RASModel.H"
#include "wallDist.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace incompressible
{
namespace RASModels
{
/*---------------------------------------------------------------------------*\
Class kOmega Declaration
\*---------------------------------------------------------------------------*/
class blockKOmegaSST
:
public RASModel
{
// Private data
// Model coefficients
dimensionedScalar alphaK1_;
dimensionedScalar alphaK2_;
dimensionedScalar alphaOmega1_;
dimensionedScalar alphaOmega2_;
dimensionedScalar gamma1_;
dimensionedScalar gamma2_;
dimensionedScalar beta1_;
dimensionedScalar beta2_;
dimensionedScalar betaStar_;
dimensionedScalar a1_;
dimensionedScalar b1_;
dimensionedScalar c1_;
Switch F3_;
//- Wall distance field
// Note: different to wall distance in parent RASModel
wallDist y_;
// Fields
volScalarField k_;
volScalarField omega_;
volScalarField nut_;
// Block vector field for k-omega
volVector2Field kOmega_;
// 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:
//- Runtime type information
TypeName("blockKOmegaSST");
// Constructors
//- Construct from components
blockKOmegaSST
(
const volVectorField& U,
const surfaceScalarField& phi,
transportModel& transport
);
//- Destructor
virtual ~blockKOmegaSST()
{}
// Member Functions
//- Return the turbulence viscosity
virtual tmp<volScalarField> nut() const
{
return nut_;
}
//- Return the effective diffusivity for k
tmp<volScalarField> DkEff(const volScalarField& F1) const
{
return tmp<volScalarField>
(
new volScalarField("DkEff", alphaK(F1)*nut_ + nu())
);
}
//- Return the effective diffusivity for omega
tmp<volScalarField> DomegaEff(const volScalarField& F1) const
{
return tmp<volScalarField>
(
new volScalarField("DomegaEff", alphaOmega(F1)*nut_ + nu())
);
}
//- Return the turbulence kinetic energy
virtual tmp<volScalarField> k() const
{
return k_;
}
//- Return the turbulence specific dissipation rate
virtual tmp<volScalarField> omega() const
{
return omega_;
}
//- Return the turbulence kinetic energy dissipation rate
virtual tmp<volScalarField> epsilon() const
{
return tmp<volScalarField>
(
new volScalarField
(
IOobject
(
"epsilon",
mesh_.time().timeName(),
mesh_
),
betaStar_*k_*omega_,
omega_.boundaryField().types()
)
);
}
//- Return the Reynolds stress tensor
virtual tmp<volSymmTensorField> R() const;
//- Return the effective stress tensor including the laminar stress
virtual tmp<volSymmTensorField> devReff() const;
//- Return the source term for the momentum equation
virtual tmp<fvVectorMatrix> divDevReff(volVectorField& U) const;
//- Solve the turbulence equations and correct the turbulence viscosity
virtual void correct();
//- Read RASProperties dictionary
virtual bool read();
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace RASModels
} // namespace incompressible
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //