From a4e2afa4b39e2f33f7b66051844cccf0296be520 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Mon, 25 Apr 2016 16:16:05 +0100 Subject: [PATCH] Completed boundaryField() -> boundaryFieldRef() 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()". --- .../turbulence/PDRkEpsilon/PDRkEpsilon.C | 4 +-- .../XiModels/XiEqModels/SCOPEXiEq/SCOPEXiEq.C | 4 ++- .../SCOPE/SCOPELaminarFlameSpeed.C | 18 ++++++++---- .../rhoCentralDyMFoam/rhoCentralDyMFoam.C | 6 ++-- .../rhoCentralFoam/rhoCentralFoam.C | 6 ++-- .../multiphaseMixtureThermo.C | 2 +- .../solvers/multiphase/interFoam/alphaEqn.H | 5 +++- .../threePhaseInterfaceProperties.C | 4 +-- .../multiphaseSystem/multiphaseSystem.C | 19 +++++++++---- .../multiphase/multiphaseEulerFoam/pEqn.H | 2 +- .../multiphaseMixture/multiphaseMixture.C | 2 +- .../saturationModels/polynomial/polynomial.C | 4 ++- .../wallLubricationModel.C | 4 ++- .../BlendedInterfacialModel.C | 5 +++- .../ThermalPhaseChangePhaseSystem.C | 4 +-- .../multiphaseSystem/multiphaseSystem.C | 8 ++++-- .../reactingMultiphaseEulerFoam/pU/pEqn.H | 7 +++-- .../reactingTwoPhaseEulerFoam/pU/pEqn.H | 28 ++++++++++++------- .../reactingTwoPhaseEulerFoam/pUf/pEqn.H | 2 +- .../JohnsonJacksonSchaefferFrictionalStress.C | 4 ++- .../Schaeffer/SchaefferFrictionalStress.C | 4 ++- .../kineticTheoryModel/kineticTheoryModel.C | 2 +- .../phasePressureModel/phasePressureModel.C | 4 +-- .../twoPhaseSystem/twoPhaseSystem.C | 8 ++++-- .../multiphase/twoPhaseEulerFoam/pU/pEqn.H | 28 ++++++++++++------- .../multiphase/twoPhaseEulerFoam/pUf/pEqn.H | 2 +- .../kineticTheoryModel/kineticTheoryModel.C | 2 +- .../phasePressureModel/phasePressureModel.C | 4 +-- .../BlendedInterfacialModel.C | 5 +++- .../twoPhaseSystem/twoPhaseSystem.C | 10 ++++--- .../GeometricField/GeometricField.H | 26 ++++++++--------- .../RAS/mixtureKEpsilon/mixtureKEpsilon.C | 2 +- .../constrainPressure/constrainPressure.C | 2 +- .../fvMatrices/solvers/MULES/MULESTemplates.C | 4 +-- 34 files changed, 148 insertions(+), 93 deletions(-) diff --git a/applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.C b/applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.C index afe97f7f7..4a2a35635 100644 --- a/applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.C +++ b/applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.C @@ -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 @@ -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_); diff --git a/applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/SCOPEXiEq/SCOPEXiEq.C b/applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/SCOPEXiEq/SCOPEXiEq.C index 0208818ee..6ca24555b 100644 --- a/applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/SCOPEXiEq/SCOPEXiEq.C +++ b/applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/SCOPEXiEq/SCOPEXiEq.C @@ -119,9 +119,11 @@ Foam::tmp 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]; diff --git a/applications/solvers/combustion/PDRFoam/laminarFlameSpeed/SCOPE/SCOPELaminarFlameSpeed.C b/applications/solvers/combustion/PDRFoam/laminarFlameSpeed/SCOPE/SCOPELaminarFlameSpeed.C index fd07b0b3c..17f92a85f 100644 --- a/applications/solvers/combustion/PDRFoam/laminarFlameSpeed/SCOPE/SCOPELaminarFlameSpeed.C +++ b/applications/solvers/combustion/PDRFoam/laminarFlameSpeed/SCOPE/SCOPELaminarFlameSpeed.C @@ -266,9 +266,11 @@ Foam::tmp 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]; @@ -313,9 +315,11 @@ Foam::tmp 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]; @@ -365,9 +369,11 @@ Foam::tmp 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) diff --git a/applications/solvers/compressible/rhoCentralFoam/rhoCentralDyMFoam/rhoCentralDyMFoam.C b/applications/solvers/compressible/rhoCentralFoam/rhoCentralDyMFoam/rhoCentralDyMFoam.C index 3b43b5dce..fb0fbc707 100644 --- a/applications/solvers/compressible/rhoCentralFoam/rhoCentralDyMFoam/rhoCentralDyMFoam.C +++ b/applications/solvers/compressible/rhoCentralFoam/rhoCentralDyMFoam/rhoCentralDyMFoam.C @@ -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) { @@ -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()) @@ -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(); diff --git a/applications/solvers/compressible/rhoCentralFoam/rhoCentralFoam.C b/applications/solvers/compressible/rhoCentralFoam/rhoCentralFoam.C index 3754632d3..083a1b1b5 100644 --- a/applications/solvers/compressible/rhoCentralFoam/rhoCentralFoam.C +++ b/applications/solvers/compressible/rhoCentralFoam/rhoCentralFoam.C @@ -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) { @@ -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()) @@ -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(); diff --git a/applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/multiphaseMixtureThermo.C b/applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/multiphaseMixtureThermo.C index 14bbbe870..db4b5972e 100644 --- a/applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/multiphaseMixtureThermo.C +++ b/applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/multiphaseMixtureThermo.C @@ -923,7 +923,7 @@ Foam::tmp Foam::multiphaseMixtureThermo::K { tmp 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()); diff --git a/applications/solvers/multiphase/interFoam/alphaEqn.H b/applications/solvers/multiphase/interFoam/alphaEqn.H index 9d6dd2dc1..f75f4d0f5 100644 --- a/applications/solvers/multiphase/interFoam/alphaEqn.H +++ b/applications/solvers/multiphase/interFoam/alphaEqn.H @@ -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()) { diff --git a/applications/solvers/multiphase/interFoam/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.C b/applications/solvers/multiphase/interFoam/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.C index 8343ceb79..d5ee2f25f 100644 --- a/applications/solvers/multiphase/interFoam/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.C +++ b/applications/solvers/multiphase/interFoam/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.C @@ -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 @@ -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; diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C b/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C index 3f514904e..62e2a97a3 100644 --- a/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C +++ b/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C @@ -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()) { @@ -372,7 +374,7 @@ Foam::tmp Foam::multiphaseSystem::K { tmp 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()); @@ -666,6 +668,9 @@ Foam::tmp Foam::multiphaseSystem::Svm } } + volVectorField::GeometricBoundaryField& SvmBf = + tSvm.ref().boundaryFieldRef(); + // Remove virtual mass at fixed-flux boundaries forAll(phase.phi().boundaryField(), patchi) { @@ -677,7 +682,7 @@ Foam::tmp Foam::multiphaseSystem::Svm ) ) { - tSvm.ref().boundaryField()[patchi] = Zero; + SvmBf[patchi] = Zero; } } @@ -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) { @@ -724,7 +731,7 @@ Foam::multiphaseSystem::dragCoeffs() const ) ) { - Kptr->boundaryField()[patchi] = 0.0; + Kbf[patchi] = 0.0; } } diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H b/applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H index bbcc36f4e..f65fdf8de 100644 --- a/applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H +++ b/applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H @@ -203,7 +203,7 @@ setSnGrad ( - p_rgh.boundaryField(), + p_rgh.boundaryFieldRef(), ( phiHbyA.boundaryField() - MRF.relative(phib) )/(mesh.magSf().boundaryField()*rAUf.boundaryField()) diff --git a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C index ea08ba298..f1fed0103 100644 --- a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C +++ b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C @@ -523,7 +523,7 @@ Foam::tmp Foam::multiphaseMixture::K { tmp 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()); diff --git a/applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/polynomial/polynomial.C b/applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/polynomial/polynomial.C index a2fbcd0bb..886da7c20 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/polynomial/polynomial.C +++ b/applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/polynomial/polynomial.C @@ -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) diff --git a/applications/solvers/multiphase/reactingEulerFoam/interfacialModels/wallLubricationModels/wallLubricationModel/wallLubricationModel.C b/applications/solvers/multiphase/reactingEulerFoam/interfacialModels/wallLubricationModels/wallLubricationModel/wallLubricationModel.C index a2a95cce0..784f323e0 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/interfacialModels/wallLubricationModels/wallLubricationModel/wallLubricationModel.C +++ b/applications/solvers/multiphase/reactingEulerFoam/interfacialModels/wallLubricationModels/wallLubricationModel/wallLubricationModel.C @@ -50,11 +50,13 @@ Foam::tmp Foam::wallLubricationModel::zeroGradWalls volVectorField& Fi = tFi.ref(); const fvPatchList& patches = Fi.mesh().boundary(); + volVectorField::GeometricBoundaryField& FiBf = Fi.boundaryFieldRef(); + forAll(patches, patchi) { if (isA(patches[patchi])) { - fvPatchVectorField& Fiw = Fi.boundaryField()[patchi]; + fvPatchVectorField& Fiw = FiBf[patchi]; Fiw = Fiw.patchInternalField(); } } diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.C index f216859b0..a19b995bc 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.C @@ -36,6 +36,9 @@ void Foam::BlendedInterfacialModel::correctFixedFluxBCs GeometricField& field ) const { + typename GeometricField::GeometricBoundaryField& fieldBf = + field.boundaryFieldRef(); + forAll(phase1_.phi()().boundaryField(), patchi) { if @@ -46,7 +49,7 @@ void Foam::BlendedInterfacialModel::correctFixedFluxBCs ) ) { - field.boundaryField()[patchi] = Zero; + fieldBf[patchi] = Zero; } } } diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C index ff1e0b099..8ba5a4a50 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C @@ -353,9 +353,9 @@ void Foam::ThermalPhaseChangePhaseSystem::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 diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C index b81e174db..6504f943a 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C @@ -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()) { @@ -477,7 +479,7 @@ Foam::tmp Foam::multiphaseSystem::K { tmp 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()); diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H index 05a049b4c..e078414a2 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H @@ -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 @@ -191,7 +194,7 @@ while (pimple.correct()) || isA(mesh.boundary()[patchi]) ) { - phiCorrCoeff.boundaryField()[patchi] = 0; + phiCorrCoeffBf[patchi] = 0; } } @@ -281,7 +284,7 @@ while (pimple.correct()) setSnGrad ( - p_rgh.boundaryField(), + p_rgh.boundaryFieldRef(), ( phiHbyA.boundaryField() - phib )/(mesh.magSf().boundaryField()*rAUf.boundaryField()) diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/pEqn.H index ba312ab24..008134c21 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/pEqn.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/pEqn.H @@ -149,17 +149,25 @@ while (pimple.correct()) surfaceScalarField phiCorrCoeff1(pos(alphaf1Bar - 0.99)); surfaceScalarField phiCorrCoeff2(pos(0.01 - alphaf1Bar)); - forAll(mesh.boundary(), patchi) { - // Set ddtPhiCorr to 0 on non-coupled boundaries - if - ( - !mesh.boundary()[patchi].coupled() - || isA(mesh.boundary()[patchi]) - ) + surfaceScalarField::GeometricBoundaryField& phiCorrCoeff1Bf = + phiCorrCoeff1.boundaryFieldRef(); + + surfaceScalarField::GeometricBoundaryField& phiCorrCoeff2Bf = + phiCorrCoeff2.boundaryFieldRef(); + + forAll(mesh.boundary(), patchi) { - phiCorrCoeff1.boundaryField()[patchi] = 0; - phiCorrCoeff2.boundaryField()[patchi] = 0; + // Set ddtPhiCorr to 0 on non-coupled boundaries + if + ( + !mesh.boundary()[patchi].coupled() + || isA(mesh.boundary()[patchi]) + ) + { + phiCorrCoeff1Bf[patchi] = 0; + phiCorrCoeff2Bf[patchi] = 0; + } } } @@ -215,7 +223,7 @@ while (pimple.correct()) // Update the fixedFluxPressure BCs to ensure flux consistency setSnGrad ( - p_rgh.boundaryField(), + p_rgh.boundaryFieldRef(), ( phiHbyA.boundaryField() - ( diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/pEqn.H index 444f2de81..01c120d48 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/pEqn.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/pEqn.H @@ -201,7 +201,7 @@ while (pimple.correct()) // Update the fixedFluxPressure BCs to ensure flux consistency setSnGrad ( - p_rgh.boundaryField(), + p_rgh.boundaryFieldRef(), ( phiHbyA.boundaryField() - ( diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/frictionalStressModel/JohnsonJacksonSchaeffer/JohnsonJacksonSchaefferFrictionalStress.C b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/frictionalStressModel/JohnsonJacksonSchaeffer/JohnsonJacksonSchaefferFrictionalStress.C index 91b49b30e..bb9683c39 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/frictionalStressModel/JohnsonJacksonSchaeffer/JohnsonJacksonSchaefferFrictionalStress.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/frictionalStressModel/JohnsonJacksonSchaeffer/JohnsonJacksonSchaefferFrictionalStress.C @@ -163,11 +163,13 @@ JohnsonJacksonSchaeffer::nu const fvPatchList& patches = phase.mesh().boundary(); const volVectorField& U = phase.U(); + volScalarField::GeometricBoundaryField& nufBf = nuf.boundaryFieldRef(); + forAll(patches, patchi) { if (!patches[patchi].coupled()) { - nuf.boundaryField()[patchi] = + nufBf[patchi] = ( pf.boundaryField()[patchi]*sin(phi_.value()) /( diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.C b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.C index 8c072802d..0853ec845 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.C @@ -152,11 +152,13 @@ Foam::kineticTheoryModels::frictionalStressModels::Schaeffer::nu const fvPatchList& patches = phase.mesh().boundary(); const volVectorField& U = phase.U(); + volScalarField::GeometricBoundaryField& nufBf = nuf.boundaryFieldRef(); + forAll(patches, patchi) { if (!patches[patchi].coupled()) { - nuf.boundaryField()[patchi] = + nufBf[patchi] = ( pf.boundaryField()[patchi]*sin(phi_.value()) /( diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C index ff3be20ab..d0a5b023d 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C @@ -295,7 +295,7 @@ Foam::RASModels::kineticTheoryModel::pPrime() const ); volScalarField::GeometricBoundaryField& bpPrime = - tpPrime.ref().boundaryField(); + tpPrime.ref().boundaryFieldRef(); forAll(bpPrime, patchi) { diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/phasePressureModel/phasePressureModel.C b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/phasePressureModel/phasePressureModel.C index 9d3b53c78..0edda3819 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/phasePressureModel/phasePressureModel.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/phasePressureModel/phasePressureModel.C @@ -165,7 +165,7 @@ Foam::RASModels::phasePressureModel::pPrime() const ); volScalarField::GeometricBoundaryField& bpPrime = - tpPrime.ref().boundaryField(); + tpPrime.ref().boundaryFieldRef(); forAll(bpPrime, patchi) { @@ -193,7 +193,7 @@ Foam::RASModels::phasePressureModel::pPrimef() const ); surfaceScalarField::GeometricBoundaryField& bpPrime = - tpPrime.ref().boundaryField(); + tpPrime.ref().boundaryFieldRef(); forAll(bpPrime, patchi) { diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C index 597d15e7f..fc86ede74 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C @@ -316,11 +316,13 @@ void Foam::twoPhaseSystem::solve() ) ); + surfaceScalarField::GeometricBoundaryField& alphaPhic1Bf = + alphaPhic1.boundaryFieldRef(); + // Ensure that the flux at inflow BCs is preserved - forAll(alphaPhic1.boundaryField(), patchi) + forAll(alphaPhic1Bf, patchi) { - fvsPatchScalarField& alphaPhic1p = - alphaPhic1.boundaryField()[patchi]; + fvsPatchScalarField& alphaPhic1p = alphaPhic1Bf[patchi]; if (!alphaPhic1p.coupled()) { diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/pU/pEqn.H b/applications/solvers/multiphase/twoPhaseEulerFoam/pU/pEqn.H index b41887774..eb3eeb562 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/pU/pEqn.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/pU/pEqn.H @@ -146,17 +146,25 @@ while (pimple.correct()) surfaceScalarField phiCorrCoeff1(pos(alphaf1Bar - 0.99)); surfaceScalarField phiCorrCoeff2(pos(0.01 - alphaf1Bar)); - forAll(mesh.boundary(), patchi) { - // Set ddtPhiCorr to 0 on non-coupled boundaries - if - ( - !mesh.boundary()[patchi].coupled() - || isA(mesh.boundary()[patchi]) - ) + surfaceScalarField::GeometricBoundaryField& phiCorrCoeff1Bf = + phiCorrCoeff1.boundaryFieldRef(); + + surfaceScalarField::GeometricBoundaryField& phiCorrCoeff2Bf = + phiCorrCoeff2.boundaryFieldRef(); + + forAll(mesh.boundary(), patchi) { - phiCorrCoeff1.boundaryField()[patchi] = 0; - phiCorrCoeff2.boundaryField()[patchi] = 0; + // Set ddtPhiCorr to 0 on non-coupled boundaries + if + ( + !mesh.boundary()[patchi].coupled() + || isA(mesh.boundary()[patchi]) + ) + { + phiCorrCoeff1Bf[patchi] = 0; + phiCorrCoeff2Bf[patchi] = 0; + } } } @@ -212,7 +220,7 @@ while (pimple.correct()) // Update the fixedFluxPressure BCs to ensure flux consistency setSnGrad ( - p_rgh.boundaryField(), + p_rgh.boundaryFieldRef(), ( phiHbyA.boundaryField() - ( diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/pUf/pEqn.H b/applications/solvers/multiphase/twoPhaseEulerFoam/pUf/pEqn.H index 6855b6fb9..bb78b2799 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/pUf/pEqn.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/pUf/pEqn.H @@ -200,7 +200,7 @@ while (pimple.correct()) // Update the fixedFluxPressure BCs to ensure flux consistency setSnGrad ( - p_rgh.boundaryField(), + p_rgh.boundaryFieldRef(), ( phiHbyA.boundaryField() - ( diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C index 6a5a474b3..8d915bc39 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C @@ -279,7 +279,7 @@ Foam::RASModels::kineticTheoryModel::pPrime() const ); volScalarField::GeometricBoundaryField& bpPrime = - tpPrime.ref().boundaryField(); + tpPrime.ref().boundaryFieldRef(); forAll(bpPrime, patchi) { diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/phasePressureModel/phasePressureModel.C b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/phasePressureModel/phasePressureModel.C index 0c5906d90..3245dec4d 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/phasePressureModel/phasePressureModel.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseCompressibleTurbulenceModels/phasePressureModel/phasePressureModel.C @@ -171,7 +171,7 @@ Foam::RASModels::phasePressureModel::pPrime() const ); volScalarField::GeometricBoundaryField& bpPrime = - tpPrime.ref().boundaryField(); + tpPrime.ref().boundaryFieldRef(); forAll(bpPrime, patchi) { @@ -199,7 +199,7 @@ Foam::RASModels::phasePressureModel::pPrimef() const ); surfaceScalarField::GeometricBoundaryField& bpPrime = - tpPrime.ref().boundaryField(); + tpPrime.ref().boundaryFieldRef(); forAll(bpPrime, patchi) { diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/BlendedInterfacialModel/BlendedInterfacialModel.C b/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/BlendedInterfacialModel/BlendedInterfacialModel.C index ecf674e6e..a3e187e7e 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/BlendedInterfacialModel/BlendedInterfacialModel.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/BlendedInterfacialModel/BlendedInterfacialModel.C @@ -36,6 +36,9 @@ void Foam::BlendedInterfacialModel::correctFixedFluxBCs GeometricField& field ) const { + typename GeometricField::GeometricBoundaryField& fieldBf = + field.boundaryFieldRef(); + forAll(pair_.phase1().phi().boundaryField(), patchi) { if @@ -46,7 +49,7 @@ void Foam::BlendedInterfacialModel::correctFixedFluxBCs ) ) { - field.boundaryField()[patchi] = Zero; + fieldBf[patchi] = Zero; } } } diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C b/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C index 76af4a310..a8e810f73 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2013-2015 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2013-2016 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -446,11 +446,13 @@ void Foam::twoPhaseSystem::solve() ) ); + surfaceScalarField::GeometricBoundaryField& alphaPhic1Bf = + alphaPhic1.boundaryFieldRef(); + // Ensure that the flux at inflow BCs is preserved - forAll(alphaPhic1.boundaryField(), patchi) + forAll(alphaPhic1Bf, patchi) { - fvsPatchScalarField& alphaPhic1p = - alphaPhic1.boundaryField()[patchi]; + fvsPatchScalarField& alphaPhic1p = alphaPhic1Bf[patchi]; if (!alphaPhic1p.coupled()) { diff --git a/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.H b/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.H index 2f1b71f7f..1d2c3e825 100644 --- a/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.H +++ b/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.H @@ -429,30 +429,28 @@ public: // Member Functions - //- Return dimensioned internal field + //- Return a reference to the dimensioned internal field + // Note: this increments the event counter and checks the + // old-time fields; avoid in loops. DimensionedInternalField& dimensionedInternalField(); - //- Return dimensioned internal field + //- Return a const-reference to the dimensioned internal field inline const DimensionedInternalField& dimensionedInternalField() const; - //- Return internal field + //- Return a reference to the internal field + // Note: this increments the event counter and checks the + // old-time fields; avoid in loops. InternalField& internalField(); - //- Return internal field + //- Return a const-reference to the internal field inline const InternalField& internalField() const; - //- Return reference to GeometricBoundaryField + //- Return a reference to the boundary field + // Note: this increments the event counter and checks the + // old-time fields; avoid in loops. GeometricBoundaryField& boundaryFieldRef(); - //- Return reference to GeometricBoundaryField - #ifndef BOUNDARY_FIELD_REF - GeometricBoundaryField& boundaryField() - { - return boundaryFieldRef(); - } - #endif - - //- Return reference to GeometricBoundaryField for const field + //- Return const-reference to the boundary field inline const GeometricBoundaryField& boundaryField() const; //- Return the time index of the field diff --git a/src/TurbulenceModels/phaseCompressible/RAS/mixtureKEpsilon/mixtureKEpsilon.C b/src/TurbulenceModels/phaseCompressible/RAS/mixtureKEpsilon/mixtureKEpsilon.C index 32f448224..01b55287f 100644 --- a/src/TurbulenceModels/phaseCompressible/RAS/mixtureKEpsilon/mixtureKEpsilon.C +++ b/src/TurbulenceModels/phaseCompressible/RAS/mixtureKEpsilon/mixtureKEpsilon.C @@ -196,7 +196,7 @@ void mixtureKEpsilon::correctInletOutlet const volScalarField& refVsf ) const { - volScalarField::GeometricBoundaryField& bf = vsf.boundaryField(); + volScalarField::GeometricBoundaryField& bf = vsf.boundaryFieldRef(); const volScalarField::GeometricBoundaryField& refBf = refVsf.boundaryField(); diff --git a/src/finiteVolume/cfdTools/general/constrainPressure/constrainPressure.C b/src/finiteVolume/cfdTools/general/constrainPressure/constrainPressure.C index 760a33546..e9f1b40ec 100644 --- a/src/finiteVolume/cfdTools/general/constrainPressure/constrainPressure.C +++ b/src/finiteVolume/cfdTools/general/constrainPressure/constrainPressure.C @@ -44,7 +44,7 @@ void Foam::constrainPressure { const fvMesh& mesh = p.mesh(); - volScalarField::GeometricBoundaryField& pBf = p.boundaryField(); + volScalarField::GeometricBoundaryField& pBf = p.boundaryFieldRef(); const volVectorField::GeometricBoundaryField& UBf = U.boundaryField(); const surfaceScalarField::GeometricBoundaryField& phiHbyABf = diff --git a/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C b/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C index 6673e5d98..b18852a1a 100644 --- a/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C +++ b/src/finiteVolume/fvMatrices/solvers/MULES/MULESTemplates.C @@ -609,7 +609,7 @@ void Foam::MULES::limitSum(SurfaceScalarFieldList& phiPsiCorrs) limitSum(phiPsiCorrsInternal); } - surfaceScalarField::GeometricBoundaryField& bfld = + const surfaceScalarField::GeometricBoundaryField& bfld = phiPsiCorrs[0].boundaryField(); forAll(bfld, patchi) @@ -622,7 +622,7 @@ void Foam::MULES::limitSum(SurfaceScalarFieldList& phiPsiCorrs) phiPsiCorrsPatch.set ( phasei, - &phiPsiCorrs[phasei].boundaryField()[patchi] + &phiPsiCorrs[phasei].boundaryFieldRef()[patchi] ); }