Merge branch 'develop' of develop.openfoam.com:Development/OpenFOAM-plus into develop

This commit is contained in:
mattijs
2016-11-02 14:48:25 +00:00

View File

@ -57,57 +57,37 @@ static const scalar Cmu(0.09);
static const scalar kappa(0.41); static const scalar kappa(0.41);
Foam::tmp<Foam::volVectorField> createSimplifiedU(const volVectorField& U) template<class Type>
{ void correctProcessorPatches
tmp<volVectorField> tU (
( GeometricField<Type, fvPatchField, volMesh>& vf
new volVectorField )
(
IOobject
(
"Udash",
U.mesh().time().timeName(),
U.mesh(),
IOobject::NO_READ
),
U.mesh(),
dimensionedVector("0", dimVelocity, vector::zero),
zeroGradientFvPatchField<vector>::typeName
)
);
// Assign the internal value to the original field
tU.ref() = U;
return tU;
}
void correctProcessorPatches(volScalarField& vf)
{ {
if (!Pstream::parRun()) if (!Pstream::parRun())
{ {
return; return;
} }
typedef GeometricField<Type, fvPatchField, volMesh> volFieldType;
// Not possible to use correctBoundaryConditions on fields as they may // Not possible to use correctBoundaryConditions on fields as they may
// use local info as opposed to the constraint values employed here, // use local info as opposed to the constraint values employed here,
// but still need to update processor patches // but still need to update processor patches
volScalarField::Boundary& bf = vf.boundaryFieldRef(); typename volFieldType::Boundary& bf = vf.boundaryFieldRef();
forAll(bf, patchI) forAll(bf, patchi)
{ {
if (isA<processorFvPatchField<scalar>>(bf[patchI])) if (isA<processorFvPatchField<Type>>(bf[patchi]))
{ {
bf[patchI].initEvaluate(); bf[patchi].initEvaluate();
} }
} }
forAll(bf, patchI) forAll(bf, patchi)
{ {
if (isA<processorFvPatchField<scalar>>(bf[patchI])) if (isA<processorFvPatchField<Type>>(bf[patchi]))
{ {
bf[patchI].evaluate(); bf[patchi].evaluate();
} }
} }
} }
@ -142,7 +122,7 @@ void blendField
// - operation may use inconsistent fields wrt these local // - operation may use inconsistent fields wrt these local
// manipulations // manipulations
//fld.correctBoundaryConditions(); //fld.correctBoundaryConditions();
correctProcessorPatches(fld); correctProcessorPatches<scalar>(fld);
Info<< "Writing " << fieldName << nl << endl; Info<< "Writing " << fieldName << nl << endl;
fld.write(); fld.write();
@ -181,7 +161,7 @@ void calcOmegaField
// - operation may use inconsistent fields wrt these local // - operation may use inconsistent fields wrt these local
// manipulations // manipulations
// omega.correctBoundaryConditions(); // omega.correctBoundaryConditions();
correctProcessorPatches(omega); correctProcessorPatches<scalar>(omega);
Info<< "Writing omega\n" << endl; Info<< "Writing omega\n" << endl;
omega.write(); omega.write();
@ -215,7 +195,7 @@ void setField
// - operation may use inconsistent fields wrt these local // - operation may use inconsistent fields wrt these local
// manipulations // manipulations
// fld.correctBoundaryConditions(); // fld.correctBoundaryConditions();
correctProcessorPatches(fld); correctProcessorPatches<scalar>(fld);
Info<< "Writing " << fieldName << nl << endl; Info<< "Writing " << fieldName << nl << endl;
fld.write(); fld.write();
@ -267,7 +247,8 @@ tmp<volScalarField> calcNut
// case for the Templated Turbulence models. The call to correct // case for the Templated Turbulence models. The call to correct
// below will evolve the turbulence model equations and update nut, // below will evolve the turbulence model equations and update nut,
// whereas only nut update is required. Need to revisit. // whereas only nut update is required. Need to revisit.
turbulence->correct(); // turbulence->correct();
turbulence->correctEnergyTransport();
return tmp<volScalarField>(new volScalarField(turbulence->nut())); return tmp<volScalarField>(new volScalarField(turbulence->nut()));
} }
@ -290,8 +271,8 @@ tmp<volScalarField> calcNut
// Note: in previous versions of the code, nut was initialised on // Note: in previous versions of the code, nut was initialised on
// construction of the turbulence model. This is no longer the // construction of the turbulence model. This is no longer the
// case for the Templated Turbulence models. The call to correct // case for the Templated Turbulence models. The call to correct
// below will evolve the turbulence model equations and update nut, // below will evolve the turbulence model equations and update
// whereas only nut update is required. Need to revisit. // nut, whereas only nut update is required. Need to revisit.
turbulence->correct(); turbulence->correct();
return tmp<volScalarField>(new volScalarField(turbulence->nut())); return tmp<volScalarField>(new volScalarField(turbulence->nut()));
@ -347,38 +328,32 @@ int main(int argc, char *argv[])
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Create a copy of the U field where BCs are simplified to zeroGradient
// to enable boundary condition update without requiring other fields,
// e.g. phi. We can then call correctBoundaryConditions on this field to
// enable appropriate behaviour for processor patches.
volVectorField Udash(createSimplifiedU(U));
// Modify velocity by applying a 1/7th power law boundary-layer // Modify velocity by applying a 1/7th power law boundary-layer
// u/U0 = (y/ybl)^(1/7) // u/U0 = (y/ybl)^(1/7)
// assumes U0 is the same as the current cell velocity // assumes U0 is the same as the current cell velocity
Info<< "Setting boundary layer velocity" << nl << endl; Info<< "Setting boundary layer velocity" << nl << endl;
scalar yblv = ybl.value(); scalar yblv = ybl.value();
forAll(Udash, celli) forAll(U, celli)
{ {
if (y[celli] <= yblv) if (y[celli] <= yblv)
{ {
mask[celli] = 1; mask[celli] = 1;
Udash[celli] *= ::pow(y[celli]/yblv, (1.0/7.0)); U[celli] *= ::pow(y[celli]/yblv, (1.0/7.0));
} }
} }
mask.correctBoundaryConditions(); mask.correctBoundaryConditions();
Udash.correctBoundaryConditions(); correctProcessorPatches<vector>(U);
// Retrieve nut from turbulence model // Retrieve nut from turbulence model
volScalarField nut(calcNut(mesh, Udash)); volScalarField nut(calcNut(mesh, U));
// Blend nut using boundary layer profile // Blend nut using boundary layer profile
volScalarField S("S", mag(dev(symm(fvc::grad(Udash))))); volScalarField S("S", mag(dev(symm(fvc::grad(U)))));
nut = (1 - mask)*nut + mask*sqr(kappa*min(y, ybl))*::sqrt(2)*S; nut = (1 - mask)*nut + mask*sqr(kappa*min(y, ybl))*::sqrt(2)*S;
// Do not correct BC - wall functions will 'undo' manipulation above // Do not correct BC - wall functions will 'undo' manipulation above
// by using nut from turbulence model // by using nut from turbulence model
correctProcessorPatches(nut); correctProcessorPatches<scalar>(nut);
nut.write(); nut.write();
// Boundary layer turbulence kinetic energy // Boundary layer turbulence kinetic energy
@ -395,10 +370,8 @@ int main(int argc, char *argv[])
calcOmegaField(mesh, mask, kBL, epsilonBL); calcOmegaField(mesh, mask, kBL, epsilonBL);
setField(mesh, "nuTilda", nut); setField(mesh, "nuTilda", nut);
// Write the updated U field
// Copy internal field Udash into U before writing
Info<< "Writing U\n" << endl; Info<< "Writing U\n" << endl;
U = Udash;
U.write(); U.write();
Info<< nl << "ExecutionTime = " << runTime.elapsedCpuTime() << " s" Info<< nl << "ExecutionTime = " << runTime.elapsedCpuTime() << " s"