Navigation Menu

Skip to content

Commit

Permalink
Completed boundaryField() -> boundaryFieldRef()
Browse files Browse the repository at this point in the history
Resolves bug-report http://www.openfoam.org/mantisbt/view.php?id=1938

Because C++ does not support overloading based on the return-type there
is a problem defining both const and non-const member functions which
are resolved based on the const-ness of the object for which they are
called rather than the intent of the programmer declared via the
const-ness of the returned type.  The issue for the "boundaryField()"
member function is that the non-const version increments the
event-counter and checks the state of the stored old-time fields in case
the returned value is altered whereas the const version has no
side-effects and simply returns the reference.  If the the non-const
function is called within the patch-loop the event-counter may overflow.
To resolve this it in necessary to avoid calling the non-const form of
"boundaryField()" if the results is not altered and cache the reference
outside the patch-loop when mutation of the patch fields is needed.

The most straight forward way of resolving this problem is to name the
const and non-const forms of the member functions differently e.g. the
non-const form could be named:

    mutableBoundaryField()
    mutBoundaryField()
    nonConstBoundaryField()
    boundaryFieldRef()

Given that in C++ a reference is non-const unless specified as const:
"T&" vs "const T&" the logical convention would be

    boundaryFieldRef()
    boundaryFieldConstRef()

and given that the const form which is more commonly used is it could
simply be named "boundaryField()" then the logical convention is

    GeometricBoundaryField& boundaryFieldRef();

    inline const GeometricBoundaryField& boundaryField() const;

This is also consistent with the new "tmp" class for which non-const
access to the stored object is obtained using the ".ref()" member function.

This new convention for non-const access to the components of
GeometricField will be applied to "dimensionedInternalField()" and "internalField()" in the
future, i.e. "dimensionedInternalFieldRef()" and "internalFieldRef()".
  • Loading branch information
Henry Weller committed Apr 25, 2016
1 parent a8bf4be commit a4e2afa
Show file tree
Hide file tree
Showing 34 changed files with 148 additions and 93 deletions.
Expand Up @@ -130,7 +130,7 @@ void PDRkEpsilon::correct()
tgradU.clear();

// Update espsilon and G at the wall
epsilon_.boundaryField().updateCoeffs();
epsilon_.boundaryFieldRef().updateCoeffs();

// Add the blockage generation term so that it is included consistently
// in both the k and epsilon equations
Expand Down Expand Up @@ -163,7 +163,7 @@ void PDRkEpsilon::correct()

epsEqn.ref().relax();

epsEqn.ref().boundaryManipulate(epsilon_.boundaryField());
epsEqn.ref().boundaryManipulate(epsilon_.boundaryFieldRef());

solve(epsEqn);
bound(epsilon_, epsilonMin_);
Expand Down
Expand Up @@ -119,9 +119,11 @@ Foam::tmp<Foam::volScalarField> Foam::XiEqModels::SCOPEXiEq::XiEq() const
}
}

volScalarField::GeometricBoundaryField& xieqBf = xieq.boundaryFieldRef();

forAll(xieq.boundaryField(), patchi)
{
scalarField& xieqp = xieq.boundaryField()[patchi];
scalarField& xieqp = xieqBf[patchi];
const scalarField& Kp = K.boundaryField()[patchi];
const scalarField& Map = Ma.boundaryField()[patchi];
const scalarField& upBySup = upBySu.boundaryField()[patchi];
Expand Down
Expand Up @@ -266,9 +266,11 @@ Foam::tmp<Foam::volScalarField> Foam::laminarFlameSpeedModels::SCOPE::Su0pTphi
Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi);
}

forAll(Su0.boundaryField(), patchi)
volScalarField::GeometricBoundaryField& Su0Bf = Su0.boundaryFieldRef();

forAll(Su0Bf, patchi)
{
scalarField& Su0p = Su0.boundaryField()[patchi];
scalarField& Su0p = Su0Bf[patchi];
const scalarField& pp = p.boundaryField()[patchi];
const scalarField& Tup = Tu.boundaryField()[patchi];

Expand Down Expand Up @@ -313,9 +315,11 @@ Foam::tmp<Foam::volScalarField> Foam::laminarFlameSpeedModels::SCOPE::Su0pTphi
Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi[celli]);
}

forAll(Su0.boundaryField(), patchi)
volScalarField::GeometricBoundaryField& Su0Bf = Su0.boundaryFieldRef();

forAll(Su0Bf, patchi)
{
scalarField& Su0p = Su0.boundaryField()[patchi];
scalarField& Su0p = Su0Bf[patchi];
const scalarField& pp = p.boundaryField()[patchi];
const scalarField& Tup = Tu.boundaryField()[patchi];
const scalarField& phip = phi.boundaryField()[patchi];
Expand Down Expand Up @@ -365,9 +369,11 @@ Foam::tmp<Foam::volScalarField> Foam::laminarFlameSpeedModels::SCOPE::Ma
ma[celli] = Ma(phi[celli]);
}

forAll(ma.boundaryField(), patchi)
volScalarField::GeometricBoundaryField& maBf = ma.boundaryFieldRef();

forAll(maBf, patchi)
{
scalarField& map = ma.boundaryField()[patchi];
scalarField& map = maBf[patchi];
const scalarField& phip = phi.boundaryField()[patchi];

forAll(map, facei)
Expand Down
Expand Up @@ -189,7 +189,7 @@ int main(int argc, char *argv[])
rhoU.dimensionedInternalField()
/rho.dimensionedInternalField();
U.correctBoundaryConditions();
rhoU.boundaryField() == rho.boundaryField()*U.boundaryField();
rhoU.boundaryFieldRef() == rho.boundaryField()*U.boundaryField();

if (!inviscid)
{
Expand Down Expand Up @@ -223,7 +223,7 @@ int main(int argc, char *argv[])
e = rhoE/rho - 0.5*magSqr(U);
e.correctBoundaryConditions();
thermo.correct();
rhoE.boundaryField() ==
rhoE.boundaryFieldRef() ==
rho.boundaryField()*
(
e.boundaryField() + 0.5*magSqr(U.boundaryField())
Expand All @@ -244,7 +244,7 @@ int main(int argc, char *argv[])
rho.dimensionedInternalField()
/psi.dimensionedInternalField();
p.correctBoundaryConditions();
rho.boundaryField() == psi.boundaryField()*p.boundaryField();
rho.boundaryFieldRef() == psi.boundaryField()*p.boundaryField();

turbulence->correct();

Expand Down
Expand Up @@ -182,7 +182,7 @@ int main(int argc, char *argv[])
rhoU.dimensionedInternalField()
/rho.dimensionedInternalField();
U.correctBoundaryConditions();
rhoU.boundaryField() == rho.boundaryField()*U.boundaryField();
rhoU.boundaryFieldRef() == rho.boundaryField()*U.boundaryField();

if (!inviscid)
{
Expand Down Expand Up @@ -216,7 +216,7 @@ int main(int argc, char *argv[])
e = rhoE/rho - 0.5*magSqr(U);
e.correctBoundaryConditions();
thermo.correct();
rhoE.boundaryField() ==
rhoE.boundaryFieldRef() ==
rho.boundaryField()*
(
e.boundaryField() + 0.5*magSqr(U.boundaryField())
Expand All @@ -237,7 +237,7 @@ int main(int argc, char *argv[])
rho.dimensionedInternalField()
/psi.dimensionedInternalField();
p.correctBoundaryConditions();
rho.boundaryField() == psi.boundaryField()*p.boundaryField();
rho.boundaryFieldRef() == psi.boundaryField()*p.boundaryField();

turbulence->correct();

Expand Down
Expand Up @@ -923,7 +923,7 @@ Foam::tmp<Foam::volScalarField> Foam::multiphaseMixtureThermo::K
{
tmp<surfaceVectorField> tnHatfv = nHatfv(alpha1, alpha2);

correctContactAngle(alpha1, alpha2, tnHatfv.ref().boundaryField());
correctContactAngle(alpha1, alpha2, tnHatfv.ref().boundaryFieldRef());

// Simple expression for curvature
return -fvc::div(tnHatfv & mesh_.Sf());
Expand Down
5 changes: 4 additions & 1 deletion applications/solvers/multiphase/interFoam/alphaEqn.H
Expand Up @@ -54,11 +54,14 @@
phic += (mixture.cAlpha()*icAlpha)*fvc::interpolate(mag(U));
}

surfaceScalarField::GeometricBoundaryField& phicBf =
phic.boundaryFieldRef();

// Do not compress interface at non-coupled boundary faces
// (inlets, outlets etc.)
forAll(phic.boundaryField(), patchi)
{
fvsPatchScalarField& phicp = phic.boundaryField()[patchi];
fvsPatchScalarField& phicp = phicBf[patchi];

if (!phicp.coupled())
{
Expand Down
Expand Up @@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
Expand Down Expand Up @@ -134,7 +134,7 @@ void Foam::threePhaseInterfaceProperties::calculateK()
// Face unit interface normal
surfaceVectorField nHatfv(gradAlphaf/(mag(gradAlphaf) + deltaN_));

correctContactAngle(nHatfv.boundaryField());
correctContactAngle(nHatfv.boundaryFieldRef());

// Face unit interface normal flux
nHatf_ = nHatfv & Sf;
Expand Down
Expand Up @@ -124,11 +124,13 @@ void Foam::multiphaseSystem::solveAlphas()
);
}

surfaceScalarField::GeometricBoundaryField& alphaPhiCorrBf =
alphaPhiCorr.boundaryFieldRef();

// Ensure that the flux at inflow BCs is preserved
forAll(alphaPhiCorr.boundaryField(), patchi)
forAll(alphaPhiCorrBf, patchi)
{
fvsPatchScalarField& alphaPhiCorrp =
alphaPhiCorr.boundaryField()[patchi];
fvsPatchScalarField& alphaPhiCorrp = alphaPhiCorrBf[patchi];

if (!alphaPhiCorrp.coupled())
{
Expand Down Expand Up @@ -372,7 +374,7 @@ Foam::tmp<Foam::volScalarField> Foam::multiphaseSystem::K
{
tmp<surfaceVectorField> tnHatfv = nHatfv(phase1, phase2);

correctContactAngle(phase1, phase2, tnHatfv.ref().boundaryField());
correctContactAngle(phase1, phase2, tnHatfv.ref().boundaryFieldRef());

// Simple expression for curvature
return -fvc::div(tnHatfv & mesh_.Sf());
Expand Down Expand Up @@ -666,6 +668,9 @@ Foam::tmp<Foam::volVectorField> Foam::multiphaseSystem::Svm
}
}

volVectorField::GeometricBoundaryField& SvmBf =
tSvm.ref().boundaryFieldRef();

// Remove virtual mass at fixed-flux boundaries
forAll(phase.phi().boundaryField(), patchi)
{
Expand All @@ -677,7 +682,7 @@ Foam::tmp<Foam::volVectorField> Foam::multiphaseSystem::Svm
)
)
{
tSvm.ref().boundaryField()[patchi] = Zero;
SvmBf[patchi] = Zero;
}
}

Expand Down Expand Up @@ -713,6 +718,8 @@ Foam::multiphaseSystem::dragCoeffs() const
)
).ptr();

volScalarField::GeometricBoundaryField& Kbf = Kptr->boundaryFieldRef();

// Remove drag at fixed-flux boundaries
forAll(dm.phase1().phi().boundaryField(), patchi)
{
Expand All @@ -724,7 +731,7 @@ Foam::multiphaseSystem::dragCoeffs() const
)
)
{
Kptr->boundaryField()[patchi] = 0.0;
Kbf[patchi] = 0.0;
}
}

Expand Down
2 changes: 1 addition & 1 deletion applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H
Expand Up @@ -203,7 +203,7 @@

setSnGrad<fixedFluxPressureFvPatchScalarField>
(
p_rgh.boundaryField(),
p_rgh.boundaryFieldRef(),
(
phiHbyA.boundaryField() - MRF.relative(phib)
)/(mesh.magSf().boundaryField()*rAUf.boundaryField())
Expand Down
Expand Up @@ -523,7 +523,7 @@ Foam::tmp<Foam::volScalarField> Foam::multiphaseMixture::K
{
tmp<surfaceVectorField> tnHatfv = nHatfv(alpha1, alpha2);

correctContactAngle(alpha1, alpha2, tnHatfv.ref().boundaryField());
correctContactAngle(alpha1, alpha2, tnHatfv.ref().boundaryFieldRef());

// Simple expression for curvature
return -fvc::div(tnHatfv & mesh_.Sf());
Expand Down
Expand Up @@ -118,9 +118,11 @@ Foam::saturationModels::polynomial::Tsat
Tsat[celli] = C_.value(p[celli]);
}

volScalarField::GeometricBoundaryField& TsatBf = Tsat.boundaryFieldRef();

forAll(Tsat.boundaryField(), patchi)
{
scalarField& Tsatp = Tsat.boundaryField()[patchi];
scalarField& Tsatp = TsatBf[patchi];
const scalarField& pp = p.boundaryField()[patchi];

forAll(Tsatp, facei)
Expand Down
Expand Up @@ -50,11 +50,13 @@ Foam::tmp<Foam::volVectorField> Foam::wallLubricationModel::zeroGradWalls
volVectorField& Fi = tFi.ref();
const fvPatchList& patches = Fi.mesh().boundary();

volVectorField::GeometricBoundaryField& FiBf = Fi.boundaryFieldRef();

forAll(patches, patchi)
{
if (isA<wallFvPatch>(patches[patchi]))
{
fvPatchVectorField& Fiw = Fi.boundaryField()[patchi];
fvPatchVectorField& Fiw = FiBf[patchi];
Fiw = Fiw.patchInternalField();
}
}
Expand Down
Expand Up @@ -36,6 +36,9 @@ void Foam::BlendedInterfacialModel<ModelType>::correctFixedFluxBCs
GeometricField& field
) const
{
typename GeometricField::GeometricBoundaryField& fieldBf =
field.boundaryFieldRef();

forAll(phase1_.phi()().boundaryField(), patchi)
{
if
Expand All @@ -46,7 +49,7 @@ void Foam::BlendedInterfacialModel<ModelType>::correctFixedFluxBCs
)
)
{
field.boundaryField()[patchi] = Zero;
fieldBf[patchi] = Zero;
}
}
}
Expand Down
Expand Up @@ -353,9 +353,9 @@ void Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::correctThermo()

// Limit the H[12] boundary field to avoid /0
const scalar HLimit = 1e-4;
H1.boundaryField() =
H1.boundaryFieldRef() =
max(H1.boundaryField(), phase1.boundaryField()*HLimit);
H2.boundaryField() =
H2.boundaryFieldRef() =
max(H2.boundaryField(), phase2.boundaryField()*HLimit);

volScalarField mDotL
Expand Down
Expand Up @@ -134,11 +134,13 @@ void Foam::multiphaseSystem::solveAlphas()
);
}

surfaceScalarField::GeometricBoundaryField& alphaPhiCorrBf =
alphaPhiCorr.boundaryFieldRef();

// Ensure that the flux at inflow BCs is preserved
forAll(alphaPhiCorr.boundaryField(), patchi)
{
fvsPatchScalarField& alphaPhiCorrp =
alphaPhiCorr.boundaryField()[patchi];
fvsPatchScalarField& alphaPhiCorrp = alphaPhiCorrBf[patchi];

if (!alphaPhiCorrp.coupled())
{
Expand Down Expand Up @@ -477,7 +479,7 @@ Foam::tmp<Foam::volScalarField> Foam::multiphaseSystem::K
{
tmp<surfaceVectorField> tnHatfv = nHatfv(phase1, phase2);

correctContactAngle(phase1, phase2, tnHatfv.ref().boundaryField());
correctContactAngle(phase1, phase2, tnHatfv.ref().boundaryFieldRef());

// Simple expression for curvature
return -fvc::div(tnHatfv & mesh_.Sf());
Expand Down
Expand Up @@ -182,6 +182,9 @@ while (pimple.correct())
);
surfaceScalarField phiCorrCoeff(pos(alphafBar - 0.99));

surfaceScalarField::GeometricBoundaryField& phiCorrCoeffBf =
phiCorrCoeff.boundaryFieldRef();

forAll(mesh.boundary(), patchi)
{
// Set ddtPhiCorr to 0 on non-coupled boundaries
Expand All @@ -191,7 +194,7 @@ while (pimple.correct())
|| isA<cyclicAMIFvPatch>(mesh.boundary()[patchi])
)
{
phiCorrCoeff.boundaryField()[patchi] = 0;
phiCorrCoeffBf[patchi] = 0;
}
}

Expand Down Expand Up @@ -281,7 +284,7 @@ while (pimple.correct())

setSnGrad<fixedFluxPressureFvPatchScalarField>
(
p_rgh.boundaryField(),
p_rgh.boundaryFieldRef(),
(
phiHbyA.boundaryField() - phib
)/(mesh.magSf().boundaryField()*rAUf.boundaryField())
Expand Down

0 comments on commit a4e2afa

Please sign in to comment.