diff --git a/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/foamToTetDualMesh.C b/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/foamToTetDualMesh.C index 72b8877c7a..6dcc567e19 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/foamToTetDualMesh.C +++ b/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/foamToTetDualMesh.C @@ -106,11 +106,10 @@ void ReadAndMapFields } else if (index < 0) { - label facei = -index-1; - label bFacei = facei - mesh.nInternalFaces(); - if (bFacei >= 0) + const label facei = -index-1; + const label patchi = mesh.boundaryMesh().patchID(facei); + if (patchi >= 0) { - label patchi = mesh.boundaryMesh().patchID()[bFacei]; label localFacei = mesh.boundaryMesh()[patchi].whichFace ( facei diff --git a/applications/utilities/preProcessing/setFields/setFields.C b/applications/utilities/preProcessing/setFields/setFields.C index a5c984cfad..2865c0e1df 100644 --- a/applications/utilities/preProcessing/setFields/setFields.C +++ b/applications/utilities/preProcessing/setFields/setFields.C @@ -391,7 +391,7 @@ bool setFaceFieldType } else { - label patchi = mesh.boundaryMesh().patchID()[bFacei]; + const label patchi = mesh.boundaryMesh().patchID()[bFacei]; allBoundaryValues[bFacei] = fieldValue; ++nChanged[patchi]; diff --git a/src/OpenFOAM/meshes/polyMesh/polyBoundaryMesh/polyBoundaryMesh.C b/src/OpenFOAM/meshes/polyMesh/polyBoundaryMesh/polyBoundaryMesh.C index 2f2ce8a568..4dca521071 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyBoundaryMesh/polyBoundaryMesh.C +++ b/src/OpenFOAM/meshes/polyMesh/polyBoundaryMesh/polyBoundaryMesh.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2018-2022 OpenCFD Ltd. + Copyright (C) 2018-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -450,6 +450,31 @@ const Foam::labelList& Foam::polyBoundaryMesh::patchID() const } +Foam::label Foam::polyBoundaryMesh::patchID(const label meshFacei) const +{ + const label bndFacei = (meshFacei - mesh_.nInternalFaces()); + + return + ( + (bndFacei >= 0 && bndFacei < mesh_.nBoundaryFaces()) + ? this->patchID()[bndFacei] + : -1 + ); +} + + +Foam::labelList +Foam::polyBoundaryMesh::patchID(const labelUList& meshFaceIndices) const +{ + labelList output(meshFaceIndices.size()); + forAll(meshFaceIndices, i) + { + output[i] = patchID(meshFaceIndices[i]); + } + return output; +} + + const Foam::HashTable& Foam::polyBoundaryMesh::groupPatchIDs() const { @@ -789,22 +814,22 @@ Foam::label Foam::polyBoundaryMesh::findPatchID Foam::labelPair -Foam::polyBoundaryMesh::whichPatchFace(const label faceIndex) const +Foam::polyBoundaryMesh::whichPatchFace(const label meshFacei) const { - if (faceIndex < mesh().nInternalFaces()) + if (meshFacei < mesh().nInternalFaces()) { // Internal face: return (-1, meshFace) - return labelPair(-1, faceIndex); + return labelPair(-1, meshFacei); } - else if (faceIndex >= mesh().nFaces()) + else if (meshFacei >= mesh().nFaces()) { // Bounds error: abort FatalErrorInFunction - << "Face " << faceIndex + << "Face " << meshFacei << " out of bounds. Number of geometric faces " << mesh().nFaces() << abort(FatalError); - return labelPair(-1, faceIndex); + return labelPair(-1, meshFacei); } @@ -812,44 +837,47 @@ Foam::polyBoundaryMesh::whichPatchFace(const label faceIndex) const // Find out which patch index and local patch face the specified // mesh face belongs to by comparing label with patch start labels. + + // TBD: use patchIDPtr_ if it exists? + const polyPatchList& patches = *this; const label patchi = findLower ( patches, - faceIndex, + meshFacei, 0, // Must include the start in the comparison [](const polyPatch& p, label val) { return (p.start() <= val); } ); - if (patchi < 0 || !patches[patchi].range().found(faceIndex)) + if (patchi < 0 || !patches[patchi].range().contains(meshFacei)) { // If not in any of above, it is trouble! FatalErrorInFunction - << "Face " << faceIndex << " not found in any of the patches " + << "Face " << meshFacei << " not found in any of the patches " << flatOutput(names()) << nl << "The patches appear to be inconsistent with the mesh :" << " internalFaces:" << mesh().nInternalFaces() << " total number of faces:" << mesh().nFaces() << abort(FatalError); - return labelPair(-1, faceIndex); + return labelPair(-1, meshFacei); } // (patch, local face index) - return labelPair(patchi, faceIndex - patches[patchi].start()); + return labelPair(patchi, meshFacei - patches[patchi].start()); } Foam::labelPairList -Foam::polyBoundaryMesh::whichPatchFace(const labelUList& faceIndices) const +Foam::polyBoundaryMesh::whichPatchFace(const labelUList& meshFaceIndices) const { - labelPairList output(faceIndices.size()); - forAll(faceIndices, i) + labelPairList output(meshFaceIndices.size()); + forAll(meshFaceIndices, i) { - output[i] = whichPatchFace(faceIndices[i]); + output[i] = whichPatchFace(meshFaceIndices[i]); } return output; } diff --git a/src/OpenFOAM/meshes/polyMesh/polyBoundaryMesh/polyBoundaryMesh.H b/src/OpenFOAM/meshes/polyMesh/polyBoundaryMesh/polyBoundaryMesh.H index 8787136c47..fa2b321a8c 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyBoundaryMesh/polyBoundaryMesh.H +++ b/src/OpenFOAM/meshes/polyMesh/polyBoundaryMesh/polyBoundaryMesh.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2018-2022 OpenCFD Ltd. + Copyright (C) 2018-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -247,20 +247,31 @@ public: //- Lookup mesh face index and return (patchi, patchFacei) tuple //- or (-1, meshFacei) for internal faces - labelPair whichPatchFace(const label faceIndex) const; + labelPair whichPatchFace(const label meshFacei) const; //- Lookup mesh face indices and return (patchi, patchFacei) tuples - labelPairList whichPatchFace(const labelUList& faceIndices) const; + labelPairList whichPatchFace(const labelUList& meshFaceIndices) const; - //- Return patch index for a given mesh face index - label whichPatch(const label faceIndex) const + //- Return patch index for a given mesh face index. + //- Uses binary search. + label whichPatch(const label meshFacei) const { - return whichPatchFace(faceIndex).first(); + return whichPatchFace(meshFacei).first(); } //- Per boundary face label the patch index const labelList& patchID() const; + //- Return patch index for a given mesh face index. + //- Uses direct lookup into patchID() list. + //- Returns -1 for internal or out-of-range faces + label patchID(const label meshFacei) const; + + //- Lookup mesh face indices and return patch indices. + //- Uses direct lookup into patchID() list. + //- Returns values of -1 for internal or out-of-range faces + labelList patchID(const labelUList& meshFaceIndices) const; + //- The patch indices per patch group const HashTable& groupPatchIDs() const; diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshTetDecomposition/polyMeshTetDecomposition.C b/src/OpenFOAM/meshes/polyMesh/polyMeshTetDecomposition/polyMeshTetDecomposition.C index 717fde2a0e..36a98b990b 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyMeshTetDecomposition/polyMeshTetDecomposition.C +++ b/src/OpenFOAM/meshes/polyMesh/polyMeshTetDecomposition/polyMeshTetDecomposition.C @@ -268,14 +268,14 @@ Foam::labelList Foam::polyMeshTetDecomposition::findFaceBasePts fI++, bFI++ ) { - label patchi = mesh.boundaryMesh().patchID()[bFI]; + const label patchi = mesh.boundaryMesh().patchID()[bFI]; if (patches[patchi].coupled()) { - const coupledPolyPatch& pp = + const coupledPolyPatch& cpp = refCast(patches[patchi]); - if (pp.owner()) + if (cpp.owner()) { boundaryFaceTetBasePtIs[bFI] = findSharedBasePoint ( @@ -340,11 +340,11 @@ Foam::labelList Foam::polyMeshTetDecomposition::findFaceBasePts continue; } - label patchi = mesh.boundaryMesh().patchID()[bFI]; + const label patchi = mesh.boundaryMesh().patchID()[bFI]; if (patches[patchi].coupled()) { - const coupledPolyPatch& pp = + const coupledPolyPatch& cpp = refCast(patches[patchi]); // Calculated base points on coupled faces are those of @@ -369,7 +369,7 @@ Foam::labelList Foam::polyMeshTetDecomposition::findFaceBasePts // +: 6 - 1 = 5 // #: 6 - 2 = 4 - if (!pp.owner()) + if (!cpp.owner()) { bFTetBasePtI = mesh.faces()[fI].size() - bFTetBasePtI; } @@ -476,7 +476,7 @@ bool Foam::polyMeshTetDecomposition::checkFaceTets } else { - label patchi = patches.patchID()[facei - mesh.nInternalFaces()]; + const label patchi = patches.patchID(facei); if (patches[patchi].coupled()) { diff --git a/src/finiteArea/faMesh/faMeshTopology.C b/src/finiteArea/faMesh/faMeshTopology.C index f91f3ec23f..daa5ac9131 100644 --- a/src/finiteArea/faMesh/faMeshTopology.C +++ b/src/finiteArea/faMesh/faMeshTopology.C @@ -138,10 +138,7 @@ Foam::faMesh::getBoundaryEdgeConnections() const const label patchFacei = edgeFaces[0]; const label meshFacei = faceLabels_[patchFacei]; - const label bndFacei = (meshFacei - mesh().nInternalFaces()); - - /// const label patchId = pbm.whichPatch(meshFacei); - const label patchId = pbm.patchID()[bndFacei]; + const label patchId = pbm.patchID(meshFacei); // Primary bookkeeping { diff --git a/src/finiteVolume/fields/fvPatchFields/derived/pressurePIDControlInletVelocity/pressurePIDControlInletVelocityFvPatchVectorFieldTemplates.C b/src/finiteVolume/fields/fvPatchFields/derived/pressurePIDControlInletVelocity/pressurePIDControlInletVelocityFvPatchVectorFieldTemplates.C index 7c30a117c6..b594bd28b5 100644 --- a/src/finiteVolume/fields/fvPatchFields/derived/pressurePIDControlInletVelocity/pressurePIDControlInletVelocityFvPatchVectorFieldTemplates.C +++ b/src/finiteVolume/fields/fvPatchFields/derived/pressurePIDControlInletVelocity/pressurePIDControlInletVelocityFvPatchVectorFieldTemplates.C @@ -39,7 +39,8 @@ void Foam::pressurePIDControlInletVelocityFvPatchVectorField::faceZoneAverage Type& average ) const { - const fvMesh& mesh(patch().boundaryMesh().mesh()); + const fvMesh& mesh = patch().boundaryMesh().mesh(); + const auto& pbm = mesh.boundaryMesh(); bitSet isMasterFace(syncTools::getInternalOrMasterFaces(mesh)); @@ -48,33 +49,30 @@ void Foam::pressurePIDControlInletVelocityFvPatchVectorField::faceZoneAverage area = 0; average = Type(Zero); - forAll(zone, faceI) + for (const label meshFacei : zone) { - const label f(zone[faceI]); - - if (mesh.isInternalFace(f)) + if (mesh.isInternalFace(meshFacei)) { - const scalar da(mesh.magSf()[f]); + const scalar da = mesh.magSf()[meshFacei]; area += da; - average += da*field[f]; + average += da*field[meshFacei]; } - else if (isMasterFace[f]) + else if (isMasterFace[meshFacei]) { - const label bf(f-mesh.nInternalFaces()); - const label patchID = mesh.boundaryMesh().patchID()[bf]; - const label lf(mesh.boundaryMesh()[patchID].whichFace(f)); - const scalar da(mesh.magSf().boundaryField()[patchID][lf]); + const label patchi = pbm.patchID(meshFacei); + const label patchFacei = pbm[patchi].whichFace(meshFacei); + const scalar da = mesh.magSf().boundaryField()[patchi][patchFacei]; area += da; - average += da*field.boundaryField()[patchID][lf]; + average += da*field.boundaryField()[patchi][patchFacei]; } } reduce(area, sumOp()); reduce(average, sumOp()); - average /= area; + average /= (area + VSMALL); } diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.C b/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.C index f5dd0c9519..354c7c2706 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.C +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.C @@ -50,13 +50,13 @@ Type Foam::interpolationCellPatchConstrained::interpolate const label facei ) const { - if (facei >= 0 && facei >= this->psi_.mesh().nInternalFaces()) - { - // Use boundary value - const polyBoundaryMesh& pbm = this->psi_.mesh().boundaryMesh(); - label patchi = pbm.patchID()[facei-this->psi_.mesh().nInternalFaces()]; - label patchFacei = pbm[patchi].whichFace(facei); + const auto& pbm = this->psi_.mesh().boundaryMesh(); + const label patchi = pbm.patchID(facei); + if (patchi >= 0) + { + // Boundary value + const label patchFacei = pbm[patchi].whichFace(facei); return this->psi_.boundaryField()[patchi][patchFacei]; } else diff --git a/src/fvOptions/sources/derived/directionalPressureGradientExplicitSource/directionalPressureGradientExplicitSource.C b/src/fvOptions/sources/derived/directionalPressureGradientExplicitSource/directionalPressureGradientExplicitSource.C index df48d608ef..d436078a18 100644 --- a/src/fvOptions/sources/derived/directionalPressureGradientExplicitSource/directionalPressureGradientExplicitSource.C +++ b/src/fvOptions/sources/derived/directionalPressureGradientExplicitSource/directionalPressureGradientExplicitSource.C @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2015-2021 OpenCFD Ltd. + Copyright (C) 2015-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -412,7 +412,7 @@ void Foam::fv::directionalPressureGradientExplicitSource::correct } else if (fZone.flipMap()[i]) { - label patchI = pbm.patchID()[faceI-mesh_.nInternalFaces()]; + const label patchI = pbm.patchID(faceI); label localFaceI = pbm[patchI].whichFace(faceI); scalar w = mesh_.magSf().boundaryField()[patchI][localFaceI]; diff --git a/src/lagrangian/basic/InteractionLists/InteractionLists.C b/src/lagrangian/basic/InteractionLists/InteractionLists.C index bd3ccbbd56..28609e9891 100644 --- a/src/lagrangian/basic/InteractionLists/InteractionLists.C +++ b/src/lagrangian/basic/InteractionLists/InteractionLists.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2019-2022 OpenCFD Ltd. + Copyright (C) 2019-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -522,10 +522,7 @@ void Foam::InteractionLists::buildInteractionLists() const face& f = mesh_.faces()[wallFaceIndex]; - label patchi = mesh_.boundaryMesh().patchID() - [ - wallFaceIndex - mesh_.nInternalFaces() - ]; + const label patchi = mesh_.boundaryMesh().patchID(wallFaceIndex); referredWallFaces_[rWFI] = referredWallFace ( @@ -945,14 +942,10 @@ void Foam::InteractionLists::prepareWallDataToRefer() globalTransforms.transformIndex(wfiat) ); - label patchi = mesh_.boundaryMesh().patchID() - [ - wallFaceIndex - mesh_.nInternalFaces() - ]; + const label patchi = mesh_.boundaryMesh().patchID(wallFaceIndex); - label patchFacei = - wallFaceIndex - - mesh_.boundaryMesh()[patchi].start(); + const label patchFacei = + mesh_.boundaryMesh()[patchi].whichFace(wallFaceIndex); // Need to transform velocity when tensor transforms are // supported diff --git a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/FacePostProcessing/FacePostProcessing.C b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/FacePostProcessing/FacePostProcessing.C index e358155e50..1240702ee8 100644 --- a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/FacePostProcessing/FacePostProcessing.C +++ b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/FacePostProcessing/FacePostProcessing.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2016-2022 OpenCFD Ltd. + Copyright (C) 2016-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -275,10 +275,12 @@ Foam::FacePostProcessing::FacePostProcessing const faceZoneMesh& fzm = owner.mesh().faceZones(); const surfaceScalarField& magSf = owner.mesh().magSf(); const polyBoundaryMesh& pbm = owner.mesh().boundaryMesh(); + forAll(faceZoneNames, i) { const word& zoneName = faceZoneNames[i]; label zoneI = fzm.findZoneID(zoneName); + if (zoneI != -1) { zoneIDs.append(zoneI); @@ -291,17 +293,15 @@ Foam::FacePostProcessing::FacePostProcessing Info<< " " << zoneName << " faces: " << nFaces << nl; scalar totArea = 0.0; - forAll(fz, j) + for (const label facei : fz) { - label facei = fz[j]; if (facei < owner.mesh().nInternalFaces()) { - totArea += magSf[fz[j]]; + totArea += magSf[facei]; } else { - label bFacei = facei - owner.mesh().nInternalFaces(); - label patchi = pbm.patchID()[bFacei]; + const label patchi = pbm.patchID(facei); const polyPatch& pp = pbm[patchi]; if diff --git a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/RemoveParcels/RemoveParcels.C b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/RemoveParcels/RemoveParcels.C index dd92c540d5..3dad8bb247 100644 --- a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/RemoveParcels/RemoveParcels.C +++ b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/RemoveParcels/RemoveParcels.C @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2020 OpenCFD Ltd. + Copyright (C) 2020-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -195,7 +195,7 @@ Foam::RemoveParcels::RemoveParcels Info<< " " << zoneName << " faces: " << nFaces << nl; scalar totArea = 0.0; - for (label facei : fz) + for (const label facei : fz) { if (facei < owner.mesh().nInternalFaces()) { @@ -203,8 +203,7 @@ Foam::RemoveParcels::RemoveParcels } else { - label bFacei = facei - owner.mesh().nInternalFaces(); - label patchi = pbm.patchID()[bFacei]; + const label patchi = pbm.patchID(facei); const polyPatch& pp = pbm[patchi]; if diff --git a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyLayerDriver.C b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyLayerDriver.C index f2a038ba65..fb84d4622b 100644 --- a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyLayerDriver.C +++ b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyLayerDriver.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2015 OpenFOAM Foundation - Copyright (C) 2015-2022 OpenCFD Ltd. + Copyright (C) 2015-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -843,7 +843,7 @@ void Foam::snappyLayerDriver::handleWarpedFaces { const face& f = pp[i]; label faceI = pp.addressing()[i]; - label patchI = patches.patchID()[faceI-mesh.nInternalFaces()]; + label patchI = patches.patchID(faceI); // It is hard to calculate some length scale if not in relative // mode so disable this check. diff --git a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.C b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.C index d0c67cef8e..c56eb54fad 100644 --- a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.C +++ b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2015 OpenFOAM Foundation - Copyright (C) 2015-2022 OpenCFD Ltd. + Copyright (C) 2015-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -1029,23 +1029,14 @@ Foam::label Foam::snappyRefineDriver::danglingCellRefine forAll(cells, celli) { - const cell& cFaces = cells[celli]; - label nIntFaces = 0; - forAll(cFaces, i) + for (const label meshFacei : cells[celli]) { - label bFacei = cFaces[i]-mesh.nInternalFaces(); - if (bFacei < 0) + const label patchi = pbm.patchID(meshFacei); + + if (patchi < 0 || pbm[patchi].coupled()) { - nIntFaces++; - } - else - { - label patchi = pbm.patchID()[bFacei]; - if (pbm[patchi].coupled()) - { - nIntFaces++; - } + ++nIntFaces; } } diff --git a/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.C b/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.C index ea5106b4cf..b56b6700ae 100644 --- a/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.C +++ b/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.C @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2017-2021 OpenCFD Ltd. + Copyright (C) 2017-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -211,7 +211,7 @@ Foam::label Foam::voxelMeshSearch::searchProcPatch const pointField& cellCentres = mesh_.cellCentres(); const polyBoundaryMesh& bMeshes = mesh_.boundaryMesh(); - label patchi = bMeshes.patchID()[faceID-mesh_.nInternalFaces()]; + const label patchi = bMeshes.patchID(faceID); const polyPatch& bMeshPatch = bMeshes[patchi]; if (!isA(bMeshPatch)) diff --git a/src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvection.C b/src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvection.C index 27b94e3db3..d9b3f75e15 100644 --- a/src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvection.C +++ b/src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvection.C @@ -368,19 +368,18 @@ void Foam::isoAdvection::timeIntegratedFlux() const surfaceScalarField::Boundary& phib = phi_.boundaryField(); const surfaceScalarField::Boundary& magSfb = mesh_.magSf().boundaryField(); surfaceScalarField::Boundary& dVfb = dVf_.boundaryFieldRef(); - const label nInternalFaces = mesh_.nInternalFaces(); // Loop through boundary surface faces forAll(bsFaces_, i) { // Get boundary face index (in the global list) const label facei = bsFaces_[i]; - const label patchi = boundaryMesh.patchID()[facei - nInternalFaces]; - const label start = boundaryMesh[patchi].start(); + const label patchi = boundaryMesh.patchID(facei); if (!phib[patchi].empty()) { - const label patchFacei = facei - start; + const label patchFacei = boundaryMesh[patchi].whichFace(facei); + const scalar phiP = phib[patchi][patchFacei]; if (phiP >= 0) @@ -601,7 +600,7 @@ void Foam::isoAdvection::checkIfOnProcPatch(const label facei) if (!mesh_.isInternalFace(facei)) { const polyBoundaryMesh& pbm = mesh_.boundaryMesh(); - const label patchi = pbm.patchID()[facei - mesh_.nInternalFaces()]; + const label patchi = pbm.patchID(facei); if (isA(pbm[patchi]) && !pbm[patchi].empty()) { diff --git a/src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvectionTemplates.C b/src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvectionTemplates.C index cf2d4c0d4a..e51b4c9246 100644 --- a/src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvectionTemplates.C +++ b/src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvectionTemplates.C @@ -50,7 +50,7 @@ Type Foam::isoAdvection::faceValue const polyBoundaryMesh& pbm = mesh_.boundaryMesh(); // Boundary face. Find out which face of which patch - const label patchi = pbm.patchID()[facei - mesh_.nInternalFaces()]; + const label patchi = pbm.patchID(facei); if (patchi < 0 || patchi >= pbm.size()) { @@ -89,7 +89,7 @@ void Foam::isoAdvection::setFaceValue const polyBoundaryMesh& pbm = mesh_.boundaryMesh(); // Boundary face. Find out which face of which patch - const label patchi = pbm.patchID()[facei - mesh_.nInternalFaces()]; + const label patchi = pbm.patchID(facei); if (patchi < 0 || patchi >= pbm.size()) { @@ -106,7 +106,6 @@ void Foam::isoAdvection::setFaceValue } const label patchFacei = pp.whichFace(facei); - f.boundaryFieldRef()[patchi][patchFacei] = value; } }