47 Cmu_(dimensioned<scalar>::lookupOrAddToDict(
51 C1_(dimensioned<scalar>::lookupOrAddToDict(
55 C2_(dimensioned<scalar>::lookupOrAddToDict(
59 C3_(dimensioned<scalar>::lookupOrAddToDict(
63 sigmak_(dimensioned<scalar>::lookupOrAddToDict(
67 sigmaEps_(dimensioned<scalar>::lookupOrAddToDict(
72 epsilon_(const_cast<volScalarField&>(
73 mesh_.thisDb().lookupObject<volScalarField>(
"epsilon"))),
77 mesh.time().timeName(),
82 dimensionedScalar(
"epsilonRes", dimensionSet(0, 0, 0, 0, 0, 0, 0), 0.0),
83 zeroGradientFvPatchField<scalar>::typeName),
84 k_(const_cast<volScalarField&>(
85 mesh_.thisDb().lookupObject<volScalarField>(
"k"))),
89 mesh.time().timeName(),
94 dimensionedScalar(
"kRes", dimensionSet(0, 0, 0, 0, 0, 0, 0), 0.0),
95 zeroGradientFvPatchField<scalar>::typeName),
100 mesh.time().timeName(),
102 IOobject::READ_IF_PRESENT,
103 IOobject::AUTO_WRITE),
105 dimensionedScalar(
"betaFIK", dimensionSet(0, 0, 0, 0, 0, 0, 0), 1.0),
110 mesh.time().timeName(),
112 IOobject::READ_IF_PRESENT,
113 IOobject::AUTO_WRITE),
115 dimensionedScalar(
"betaFIEpsilon", dimensionSet(0, 0, 0, 0, 0, 0, 0), 1.0),
122 epsilonRes_.dimensions().reset(dimensionSet(0, 2, -4, 0, 0, 0, 0));
123 kRes_.dimensions().reset(dimensionSet(0, 2, -3, 0, 0, 0, 0));
128 epsilonRes_.dimensions().reset(dimensionSet(1, -1, -4, 0, 0, 0, 0));
129 kRes_.dimensions().reset(dimensionSet(1, -1, -3, 0, 0, 0, 0));
133 label nWallFaces = 0;
136 if (
epsilon_.boundaryField()[patchI].type() ==
"epsilonWallFunction"
137 and
epsilon_.boundaryField()[patchI].size() > 0)
151 tmp<volTensorField> tgradU = fvc::grad(
U_);
152 GPtr_.reset(
new volScalarField::Internal(
154 nut_.v() * (dev(twoSymm(tgradU().v())) && tgradU().v())));
162 return tmp<fvScalarMatrix>(
171 return tmp<fvScalarMatrix>(
180 return tmp<volScalarField>(
188 return tmp<volScalarField>(
222 word stateName = modelStates[idxI];
223 if (stateName ==
"nut")
225 modelStates[idxI] =
"epsilon";
226 modelStates.append(
"k");
241 nut_.correctBoundaryConditions();
258 k_.correctBoundaryConditions();
283 epsilon_.correctBoundaryConditions();
297 if (
epsilon_.boundaryField()[patchI].type() ==
"epsilonWallFunction"
298 and
epsilon_.boundaryField()[patchI].size() > 0)
300 const UList<label>& faceCells =
mesh_.boundaryMesh()[patchI].faceCells();
322 if (
epsilon_.boundaryField()[patchI].type() ==
"epsilonWallFunction"
323 &&
epsilon_.boundaryField()[patchI].size() > 0)
325 const UList<label>& faceCells =
mesh_.boundaryMesh()[patchI].faceCells();
388 forAll(stateCon[idxI], idxJ)
390 word conStateName = stateCon[idxI][idxJ];
391 if (conStateName ==
"nut")
393 stateCon[idxI][idxJ] =
"epsilon";
394 stateCon[idxI].append(
"k");
440 if (
mesh_.thisDb().foundObject<volScalarField>(
"p"))
444 else if (
mesh_.thisDb().foundObject<volScalarField>(
"p_rgh"))
451 "Neither p nor p_rgh was found in mesh.thisDb()!"
452 "addModelResidualCon failed to setup turbulence residuals!")
462 {
"U",
"epsilon",
"k",
"phi"},
463 {
"U",
"epsilon",
"k"},
464 {
"U",
"epsilon",
"k"}
469 {
"U",
"epsilon",
"k",
"phi"},
470 {
"U",
"epsilon",
"k"},
471 {
"U",
"epsilon",
"k"}
479 {
"U",
"T", pName,
"epsilon",
"k",
"phi"},
480 {
"U",
"T", pName,
"epsilon",
"k"},
481 {
"U",
"T", pName,
"epsilon",
"k"}
486 {
"U",
"T", pName,
"epsilon",
"k",
"phi"},
487 {
"U",
"T", pName,
"epsilon",
"k"},
488 {
"U",
"T", pName,
"epsilon",
"k"}
506 dictionary dummyOptions;
507 dummyOptions.set(
"printToScreen", printToScreen);
536 word divKScheme =
"div(phi,k)";
537 word divEpsilonScheme =
"div(phi,epsilon)";
541 label printToScreen = options.lookupOrDefault(
"printToScreen", 0);
545 isPC = options.getLabel(
"isPC");
549 divKScheme =
"div(pc)";
550 divEpsilonScheme =
"div(pc)";
554 volScalarField
rho = this->
rho();
559 volScalarField::Internal divU(
560 fvc::div(fvc::absolute(
phi_ / fvc::interpolate(
rho),
U_))().v());
562 tmp<volTensorField> tgradU = fvc::grad(
U_);
563 volScalarField::Internal& G =
const_cast<volScalarField::Internal&
>(
GPtr_());
564 G =
nut_.v() * (dev(twoSymm(tgradU().v())) && tgradU().v());
570 epsilon_.boundaryFieldRef().updateCoeffs();
579 tmp<fvScalarMatrix> epsEqn(
588 epsEqn.ref().relax();
590 epsEqn.ref().boundaryManipulate(
epsilon_.boundaryFieldRef());
596 SolverPerformance<scalar> solverEpsilon =
solve(epsEqn);
613 tmp<fvScalarMatrix> kEqn(
618 - fvm::SuSp((2.0 / 3.0) *
phase_() *
rho() * divU,
k_)
629 SolverPerformance<scalar> solverK =
solve(kEqn);
659 if (varName !=
"k" && varName !=
"epsilon")
662 "varName not valid. It has to be k or epsilon")
666 volScalarField
rho = this->
rho();
671 volScalarField::Internal divU(
672 fvc::div(fvc::absolute(
phi_ / fvc::interpolate(
rho),
U_))().v());
674 tmp<volTensorField> tgradU = fvc::grad(
U_);
675 volScalarField::Internal& G =
const_cast<volScalarField::Internal&
>(
GPtr_());
676 G =
nut_.v() * (dev(twoSymm(tgradU().v())) && tgradU().v());
682 if (varName ==
"epsilon")
685 fvScalarMatrix epsEqn(
700 upper = epsEqn.upper();
701 lower = epsEqn.lower();
703 else if (varName ==
"k")
710 - fvm::SuSp((2.0 / 3.0) *
phase_() *
rho() * divU,
k_)
717 upper = kEqn.upper();
718 lower = kEqn.lower();
728 tmp<volTensorField> tgradU = fvc::grad(
U_);
729 volScalarField::Internal& G =
const_cast<volScalarField::Internal&
>(
GPtr_());
730 G =
nut_.v() * (dev(twoSymm(tgradU().v())) && tgradU().v());
732 volScalarField
rho = this->
rho();
734 volScalarField::Internal P =
phase_() *
rho() * G;
739 PoD[cellI] = P[cellI] / (
D[cellI] + P[cellI] + 1e-16);
750 tmp<volTensorField> tgradU = fvc::grad(
U_);
751 volScalarField::Internal& G =
const_cast<volScalarField::Internal&
>(
GPtr_());
752 G =
nut_.v() * (dev(twoSymm(tgradU().v())) && tgradU().v());
754 volScalarField
rho = this->
rho();
756 volScalarField::Internal P =
phase_() *
rho() * G;
761 CoP[cellI] = C[cellI] / (P[cellI] + C[cellI] + 1e-16);