diff --git a/src/turbulenceModels/incompressible/RAS/Make/files b/src/turbulenceModels/incompressible/RAS/Make/files
index 51a3bbfa1..9eba12e34 100644
--- a/src/turbulenceModels/incompressible/RAS/Make/files
+++ b/src/turbulenceModels/incompressible/RAS/Make/files
@@ -18,6 +18,9 @@ NonlinearKEShih/NonlinearKEShih.C
LienLeschzinerLowRe/LienLeschzinerLowRe.C
LamBremhorstKE/LamBremhorstKE.C
+coupledKEpsilon/coupledKEpsilon.C
+coupledKOmegaSST/coupledKOmegaSST.C
+
/* Wall functions */
wallFunctions = derivedFvPatchFields/wallFunctions
diff --git a/src/turbulenceModels/incompressible/RAS/coupledKEpsilon/coupledKEpsilon.C b/src/turbulenceModels/incompressible/RAS/coupledKEpsilon/coupledKEpsilon.C
new file mode 100644
index 000000000..5940f6d05
--- /dev/null
+++ b/src/turbulenceModels/incompressible/RAS/coupledKEpsilon/coupledKEpsilon.C
@@ -0,0 +1,333 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | foam-extend: Open Source CFD
+ \\ / O peration |
+ \\ / A nd | For copyright notice see file Copyright
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+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 "coupledKEpsilon.H"
+#include "fvBlockMatrix.H"
+#include "addToRunTimeSelectionTable.H"
+
+#include "backwardsCompatibilityWallFunctions.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace incompressible
+{
+namespace RASModels
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(coupledKEpsilon, 0);
+addToRunTimeSelectionTable(RASModel, coupledKEpsilon, dictionary);
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+coupledKEpsilon::coupledKEpsilon
+(
+ const volVectorField& U,
+ const surfaceScalarField& phi,
+ transportModel& lamTransportModel
+)
+:
+ RASModel(typeName, U, phi, lamTransportModel),
+
+ Cmu_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "Cmu",
+ coeffDict_,
+ 0.09
+ )
+ ),
+ C1_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "C1",
+ coeffDict_,
+ 1.44
+ )
+ ),
+ C2_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "C2",
+ coeffDict_,
+ 1.92
+ )
+ ),
+ sigmaEps_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "sigmaEps",
+ coeffDict_,
+ 1.3
+ )
+ ),
+
+ k_
+ (
+ IOobject
+ (
+ "k",
+ runTime_.timeName(),
+ U_.db(),
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE
+ ),
+ autoCreateK("k", mesh_)
+ ),
+ epsilon_
+ (
+ IOobject
+ (
+ "epsilon",
+ runTime_.timeName(),
+ U_.db(),
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE
+ ),
+ autoCreateEpsilon("epsilon", mesh_)
+ ),
+ nut_
+ (
+ IOobject
+ (
+ "nut",
+ runTime_.timeName(),
+ U_.db(),
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE
+ ),
+ autoCreateNut("nut", mesh_)
+ ),
+ kEpsilon_
+ (
+ IOobject
+ (
+ "kEpsilon",
+ runTime_.timeName(),
+ U_.db(),
+ IOobject::NO_READ,
+ IOobject::NO_WRITE
+ ),
+ mesh_,
+ dimensionedVector2("zero", dimless, vector2::zero)
+ )
+{
+ nut_ = Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
+ nut_ = min(nut_, nuRatio()*nu());
+ nut_.correctBoundaryConditions();
+
+ printCoeffs();
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+tmp coupledKEpsilon::R() const
+{
+ return tmp
+ (
+ 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 coupledKEpsilon::devReff() const
+{
+ return tmp
+ (
+ new volSymmTensorField
+ (
+ IOobject
+ (
+ "devRhoReff",
+ runTime_.timeName(),
+ U_.db(),
+ IOobject::NO_READ,
+ IOobject::NO_WRITE
+ ),
+ -nuEff()*dev(twoSymm(fvc::grad(U_)))
+ )
+ );
+}
+
+
+tmp coupledKEpsilon::divDevReff(volVectorField& U) const
+{
+ return
+ (
+ - fvm::laplacian(nuEff(), U)
+ - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
+ );
+}
+
+
+bool coupledKEpsilon::read()
+{
+ if (RASModel::read())
+ {
+ Cmu_.readIfPresent(coeffDict());
+ C1_.readIfPresent(coeffDict());
+ C2_.readIfPresent(coeffDict());
+ sigmaEps_.readIfPresent(coeffDict());
+
+ return true;
+ }
+ else
+ {
+ return false;
+ }
+}
+
+
+void coupledKEpsilon::correct()
+{
+ // Bound in case of topological change
+ // HJ, 22/Aug/2007
+ if (mesh_.changing())
+ {
+ bound(k_, k0_);
+ bound(epsilon_, epsilon0_);
+ }
+
+ RASModel::correct();
+
+ if (!turbulence_)
+ {
+ return;
+ }
+
+ // Make coupled matrix
+ fvBlockMatrix keEqn(kEpsilon_);
+
+ volScalarField G("RASModel::G", nut_*2*magSqr(symm(fvc::grad(U_))));
+
+ // Update epsilon and G at the wall
+ epsilon_.boundaryField().updateCoeffs();
+
+ // Dissipation equation
+ {
+ fvScalarMatrix epsEqn
+ (
+ fvm::ddt(epsilon_)
+ + fvm::div(phi_, epsilon_)
+ + fvm::SuSp(-fvc::div(phi_), epsilon_)
+ - fvm::laplacian(DepsilonEff(), epsilon_)
+ + fvm::Sp(C2_*epsilon_/k_, epsilon_)
+ ==
+ C1_*G*epsilon_/k_
+ );
+
+ epsEqn.relax();
+ epsEqn.completeAssembly();
+
+ keEqn.insertEquation(1, epsEqn);
+
+ // Add coupling term:
+ // G_epsilon = C1*Cmu*(symm(grad(U))) k
+ // but with wall function corrections: must be calculated from G
+ // HJ, 27/Apr/2015
+ volScalarField coupling
+ (
+ "coupling",
+ -C1_*G*epsilon_/sqr(k_)
+ );
+ scalarField& couplingIn = coupling.internalField();
+
+ // Eliminate coupling in wall function cells
+ labelList eliminated = epsEqn.eliminatedEqns().toc();
+
+ forAll (eliminated, cellI)
+ {
+ couplingIn[eliminated[cellI]] *= 0;
+ }
+
+ // Insert coupling
+ keEqn.insertEquationCoupling(1, 0, coupling);
+ }
+
+ // Turbulent kinetic energy equation
+ {
+ fvScalarMatrix kEqn
+ (
+ fvm::ddt(k_)
+ + fvm::div(phi_, k_)
+ + fvm::SuSp(-fvc::div(phi_), k_)
+ - fvm::laplacian(DkEff(), k_)
+ + fvm::SuSp((epsilon_ - G)/k_, k_)
+ );
+
+ kEqn.relax();
+
+ keEqn.insertEquation(0, kEqn);
+ }
+
+ // Update source coupling: coupling terms eliminated from source
+ keEqn.updateSourceCoupling();
+
+ keEqn.solve();
+
+ // Retrieve solution
+ keEqn.retrieveSolution(0, k_.internalField());
+ keEqn.retrieveSolution(1, epsilon_.internalField());
+
+ k_.correctBoundaryConditions();
+ epsilon_.correctBoundaryConditions();
+
+ bound(epsilon_, epsilon0_);
+ bound(k_, k0_);
+
+ // Re-calculate viscosity
+ nut_ = Cmu_*sqr(k_)/epsilon_;
+ nut_ = min(nut_, nuRatio()*nu());
+ nut_.correctBoundaryConditions();
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace RASModels
+} // End namespace incompressible
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/turbulenceModels/incompressible/RAS/coupledKEpsilon/coupledKEpsilon.H b/src/turbulenceModels/incompressible/RAS/coupledKEpsilon/coupledKEpsilon.H
new file mode 100644
index 000000000..42444564f
--- /dev/null
+++ b/src/turbulenceModels/incompressible/RAS/coupledKEpsilon/coupledKEpsilon.H
@@ -0,0 +1,178 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | foam-extend: Open Source CFD
+ \\ / O peration |
+ \\ / A nd | For copyright notice see file Copyright
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+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::coupledKEpsilon
+
+Description
+ Standard k-epsilon turbulence model for incompressible flows using
+ a block-coupled solution.
+
+ The default model coefficients correspond to the following:
+ @verbatim
+ coupledKEpsilonCoeffs
+ {
+ Cmu 0.09;
+ C1 1.44;
+ C2 1.92;
+ sigmaEps 1.3;
+ }
+ @endverbatim
+
+Author
+ Hrvoje Jasak, Wikki Ltd. All rights reserved.
+
+SourceFiles
+ coupledKEpsilon.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef coupledKEpsilon_H
+#define coupledKEpsilon_H
+
+#include "RASModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace incompressible
+{
+namespace RASModels
+{
+
+/*---------------------------------------------------------------------------*\
+ Class coupledKEpsilon Declaration
+\*---------------------------------------------------------------------------*/
+
+class coupledKEpsilon
+:
+ public RASModel
+{
+ // Private data
+
+ // Model coefficients
+
+ dimensionedScalar Cmu_;
+ dimensionedScalar C1_;
+ dimensionedScalar C2_;
+ dimensionedScalar sigmaEps_;
+
+
+ // Fields
+
+ volScalarField k_;
+ volScalarField epsilon_;
+ volScalarField nut_;
+
+ // Block vector field for k-epsilon
+ volVector2Field kEpsilon_;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("coupledKEpsilon");
+
+ // Constructors
+
+ //- Construct from components
+ coupledKEpsilon
+ (
+ const volVectorField& U,
+ const surfaceScalarField& phi,
+ transportModel& transport
+ );
+
+
+ //- Destructor
+ virtual ~coupledKEpsilon()
+ {}
+
+
+ // Member Functions
+
+ //- Return the turbulence viscosity
+ virtual tmp nut() const
+ {
+ return nut_;
+ }
+
+ //- Return the effective diffusivity for k
+ tmp DkEff() const
+ {
+ return tmp
+ (
+ new volScalarField("DkEff", nut_ + nu())
+ );
+ }
+
+ //- Return the effective diffusivity for epsilon
+ tmp DepsilonEff() const
+ {
+ return tmp
+ (
+ new volScalarField("DepsilonEff", nut_/sigmaEps_ + nu())
+ );
+ }
+
+ //- Return the turbulence kinetic energy
+ virtual tmp k() const
+ {
+ return k_;
+ }
+
+ //- Return the turbulence kinetic energy dissipation rate
+ virtual tmp epsilon() const
+ {
+ return epsilon_;
+ }
+
+ //- Return the Reynolds stress tensor
+ virtual tmp R() const;
+
+ //- Return the effective stress tensor including the laminar stress
+ virtual tmp devReff() const;
+
+ //- Return the source term for the momentum equation
+ virtual tmp 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
+} // End namespace incompressible
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/turbulenceModels/incompressible/RAS/coupledKOmegaSST/coupledKOmegaSST.C b/src/turbulenceModels/incompressible/RAS/coupledKOmegaSST/coupledKOmegaSST.C
new file mode 100644
index 000000000..4b86a39fc
--- /dev/null
+++ b/src/turbulenceModels/incompressible/RAS/coupledKOmegaSST/coupledKOmegaSST.C
@@ -0,0 +1,469 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | foam-extend: Open Source CFD
+ \\ / O peration |
+ \\ / A nd | For copyright notice see file Copyright
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+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 "coupledKOmegaSST.H"
+#include "fvBlockMatrix.H"
+#include "addToRunTimeSelectionTable.H"
+
+#include "backwardsCompatibilityWallFunctions.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace incompressible
+{
+namespace RASModels
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(coupledKOmegaSST, 0);
+addToRunTimeSelectionTable(RASModel, coupledKOmegaSST, dictionary);
+
+// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
+
+tmp coupledKOmegaSST::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 coupledKOmegaSST::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));
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+coupledKOmegaSST::coupledKOmegaSST
+(
+ const volVectorField& U,
+ const surfaceScalarField& phi,
+ transportModel& lamTransportModel
+)
+:
+ RASModel(typeName, U, phi, lamTransportModel),
+
+ alphaK1_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "alphaK1",
+ coeffDict_,
+ 0.85034
+ )
+ ),
+ alphaK2_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "alphaK2",
+ coeffDict_,
+ 1.0
+ )
+ ),
+ alphaOmega1_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "alphaOmega1",
+ coeffDict_,
+ 0.5
+ )
+ ),
+ alphaOmega2_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "alphaOmega2",
+ coeffDict_,
+ 0.85616
+ )
+ ),
+ gamma1_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "gamma1",
+ coeffDict_,
+ 0.5532
+ )
+ ),
+ gamma2_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "gamma2",
+ coeffDict_,
+ 0.4403
+ )
+ ),
+ beta1_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "beta1",
+ coeffDict_,
+ 0.075
+ )
+ ),
+ beta2_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "beta2",
+ coeffDict_,
+ 0.0828
+ )
+ ),
+ betaStar_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "betaStar",
+ coeffDict_,
+ 0.09
+ )
+ ),
+ a1_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "a1",
+ coeffDict_,
+ 0.31
+ )
+ ),
+ c1_
+ (
+ dimensioned::lookupOrAddToDict
+ (
+ "c1",
+ coeffDict_,
+ 10.0
+ )
+ ),
+
+ 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(omega_, omega0_);
+
+ nut_ = a1_*k_/max(a1_*omega_, F2()*sqrt(2.0)*mag(symm(fvc::grad(U_))));
+ nut_.correctBoundaryConditions();
+
+ printCoeffs();
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+tmp coupledKOmegaSST::R() const
+{
+ return tmp
+ (
+ 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 coupledKOmegaSST::devReff() const
+{
+ return tmp
+ (
+ new volSymmTensorField
+ (
+ IOobject
+ (
+ "devRhoReff",
+ runTime_.timeName(),
+ U_.db(),
+ IOobject::NO_READ,
+ IOobject::NO_WRITE
+ ),
+ -nuEff()*dev(twoSymm(fvc::grad(U_)))
+ )
+ );
+}
+
+
+tmp coupledKOmegaSST::divDevReff(volVectorField& U) const
+{
+ return
+ (
+ - fvm::laplacian(nuEff(), U)
+ - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
+ );
+}
+
+
+bool coupledKOmegaSST::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());
+ c1_.readIfPresent(coeffDict());
+
+ return true;
+ }
+ else
+ {
+ return false;
+ }
+}
+
+
+void coupledKOmegaSST::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();
+ }
+
+ volScalarField S2 = magSqr(symm(fvc::grad(U_)));
+ volScalarField G("RASModel::G", nut_*2*S2);
+
+ // Make coupled matrix
+ fvBlockMatrix kOmegaEqn(kOmega_);
+
+ // 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);
+
+ Info<< "Coupled k-omega" << endl;
+
+ // Turbulent frequency equation
+ {
+ fvScalarMatrix omegaEqn
+ (
+ fvm::ddt(omega_)
+ + fvm::div(phi_, omega_)
+ + fvm::SuSp(-fvc::div(phi_), omega_)
+ - fvm::laplacian(DomegaEff(F1), omega_)
+ + fvm::SuSp
+ (
+ beta(F1)*omega_
+ + (F1 - scalar(1))*CDkOmega/omega_,
+ omega_
+ )
+ ==
+ gamma(F1)*2*S2
+ );
+
+ omegaEqn.relax();
+ omegaEqn.completeAssembly();
+
+ kOmegaEqn.insertEquation(1, omegaEqn);
+
+ // Add coupling term
+ volScalarField coupling
+ (
+ "coupling",
+ -gamma(F1)*2*S2/k_
+ );
+ scalarField& couplingIn = coupling.internalField();
+
+ // Eliminate coupling in wall function cells
+ labelList eliminated = omegaEqn.eliminatedEqns().toc();
+
+ forAll (eliminated, cellI)
+ {
+ couplingIn[eliminated[cellI]] *= 0;
+ }
+
+ // Insert coupling
+ kOmegaEqn.insertEquationCoupling(1, 0, coupling);
+ }
+
+ // Turbulent kinetic energy equation
+ {
+ fvScalarMatrix kEqn
+ (
+ fvm::ddt(k_)
+ + fvm::div(phi_, k_)
+ + fvm::SuSp(-fvc::div(phi_), k_)
+ - fvm::laplacian(DkEff(F1), k_)
+ + fvm::SuSp
+ (
+ betaStar_*omega_
+ - min(G/k_, c1_*betaStar_*omega_),
+ k_
+ )
+ );
+
+ kEqn.relax();
+
+ kOmegaEqn.insertEquation(0, kEqn);
+ }
+
+ // Update source coupling: coupling terms eliminated from source
+ kOmegaEqn.updateSourceCoupling();
+
+ kOmegaEqn.solve();
+
+ // Retrieve solution
+ kOmegaEqn.retrieveSolution(0, k_.internalField());
+ kOmegaEqn.retrieveSolution(1, omega_.internalField());
+
+ k_.correctBoundaryConditions();
+ omega_.correctBoundaryConditions();
+
+ bound(k_, k0_);
+ bound(omega_, omega0_);
+
+ // Re-calculate viscosity
+ nut_ = a1_*k_/max(a1_*omega_, F2()*sqrt(2*S2));
+ nut_ = min(nut_, nuRatio()*nu());
+ nut_.correctBoundaryConditions();
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace RASModels
+} // End namespace incompressible
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/turbulenceModels/incompressible/RAS/coupledKOmegaSST/coupledKOmegaSST.H b/src/turbulenceModels/incompressible/RAS/coupledKOmegaSST/coupledKOmegaSST.H
new file mode 100644
index 000000000..dd5b0d575
--- /dev/null
+++ b/src/turbulenceModels/incompressible/RAS/coupledKOmegaSST/coupledKOmegaSST.H
@@ -0,0 +1,295 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | foam-extend: Open Source CFD
+ \\ / O peration |
+ \\ / A nd | For copyright notice see file Copyright
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+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::coupledKOmegaSST
+
+Description
+ Implementation of the k-omega-SST turbulence model for incompressible
+ flows using a block-coupled solution.
+
+ 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
+
+ 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
+ coupledKOmegaSSTCoeffs
+ {
+ alphaK1 0.85034;
+ alphaK2 1.0;
+ alphaOmega1 0.5;
+ alphaOmega2 0.85616;
+ beta1 0.075;
+ beta2 0.0828;
+ betaStar 0.09;
+ gamma1 0.5532;
+ gamma2 0.4403;
+ a1 0.31;
+ c1 10.0;
+ }
+ @endverbatim
+
+Author
+ Hrvoje Jasak, Wikki Ltd. All rights reserved.
+
+SourceFiles
+ coupledKOmegaSST.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef coupledKOmegaSST_H
+#define coupledKOmegaSST_H
+
+#include "RASModel.H"
+#include "wallDist.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace incompressible
+{
+namespace RASModels
+{
+
+/*---------------------------------------------------------------------------*\
+ Class kOmega Declaration
+\*---------------------------------------------------------------------------*/
+
+class coupledKOmegaSST
+:
+ 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 c1_;
+
+
+ //- 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 F1(const volScalarField& CDkOmega) const;
+ tmp F2() const;
+
+ tmp blend
+ (
+ const volScalarField& F1,
+ const dimensionedScalar& psi1,
+ const dimensionedScalar& psi2
+ ) const
+ {
+ return F1*(psi1 - psi2) + psi2;
+ }
+
+ tmp alphaK
+ (
+ const volScalarField& F1
+ ) const
+ {
+ return blend(F1, alphaK1_, alphaK2_);
+ }
+
+ tmp alphaOmega
+ (
+ const volScalarField& F1
+ ) const
+ {
+ return blend(F1, alphaOmega1_, alphaOmega2_);
+ }
+
+ tmp beta
+ (
+ const volScalarField& F1
+ ) const
+ {
+ return blend(F1, beta1_, beta2_);
+ }
+
+ tmp gamma
+ (
+ const volScalarField& F1
+ ) const
+ {
+ return blend(F1, gamma1_, gamma2_);
+ }
+
+
+public:
+
+ //- Runtime type information
+ TypeName("coupledKOmegaSST");
+
+
+ // Constructors
+
+ //- Construct from components
+ coupledKOmegaSST
+ (
+ const volVectorField& U,
+ const surfaceScalarField& phi,
+ transportModel& transport
+ );
+
+
+ //- Destructor
+
+ virtual ~coupledKOmegaSST()
+ {}
+
+
+ // Member Functions
+
+ //- Return the turbulence viscosity
+ virtual tmp nut() const
+ {
+ return nut_;
+ }
+
+ //- Return the effective diffusivity for k
+ tmp DkEff(const volScalarField& F1) const
+ {
+ return tmp
+ (
+ new volScalarField("DkEff", alphaK(F1)*nut_ + nu())
+ );
+ }
+
+ //- Return the effective diffusivity for omega
+ tmp DomegaEff(const volScalarField& F1) const
+ {
+ return tmp
+ (
+ new volScalarField("DomegaEff", alphaOmega(F1)*nut_ + nu())
+ );
+ }
+
+ //- Return the turbulence kinetic energy
+ virtual tmp k() const
+ {
+ return k_;
+ }
+
+ //- Return the turbulence specific dissipation rate
+ virtual tmp omega() const
+ {
+ return omega_;
+ }
+
+ //- Return the turbulence kinetic energy dissipation rate
+ virtual tmp epsilon() const
+ {
+ return tmp
+ (
+ new volScalarField
+ (
+ IOobject
+ (
+ "epsilon",
+ mesh_.time().timeName(),
+ mesh_
+ ),
+ betaStar_*k_*omega_,
+ omega_.boundaryField().types()
+ )
+ );
+ }
+
+ //- Return the Reynolds stress tensor
+ virtual tmp R() const;
+
+ //- Return the effective stress tensor including the laminar stress
+ virtual tmp devReff() const;
+
+ //- Return the source term for the momentum equation
+ virtual tmp 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
+
+// ************************************************************************* //