porosityModel, sampledIsoSurfaceSurface: Removed unused multi cellZone support

This reduces the code complexity and allows for a future plan to make zones
dynamic rather than a static fixed indexed list.

src/OpenFOAM/meshes/polyMesh/polyMeshIO.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -306,128 +306,69 @@ Foam::polyMesh::readUpdateState Foam::polyMesh::readUpdate()
         geometricD_ = Zero;
         solutionD_ = Zero;
 
-        // Zones
+        // pointZones
-        meshPointZones newPointZones
+        {
-        (
+            meshPointZones newPointZones
-            IOobject
             (
-                "pointZones",
+                IOobject
-                facesInst,
+                (
-                meshSubDir,
+                    "pointZones",
-                *this,
+                    facesInst,
-                IOobject::READ_IF_PRESENT,
+                    meshSubDir,
-                IOobject::NO_WRITE,
+                    *this,
-                false
+                    IOobject::READ_IF_PRESENT,
-            ),
+                    IOobject::NO_WRITE,
-            *this
+                    false
-        );
+                ),
-
+                *this
-        label oldSize = pointZones_.size();
-
-        if (newPointZones.size() <= pointZones_.size())
-        {
-            pointZones_.setSize(newPointZones.size());
-        }
-
-        // Reset existing ones
-        forAll(pointZones_, czI)
-        {
-            pointZones_[czI] = newPointZones[czI];
-        }
-
-        // Extend with extra ones
-        pointZones_.setSize(newPointZones.size());
-
-        for (label czI = oldSize; czI < newPointZones.size(); czI++)
-        {
-            pointZones_.set(czI, newPointZones[czI].clone(pointZones_));
-        }
-
-        pointZones_.instance() = facesInst;
-
-
-        meshFaceZones newFaceZones
-        (
-            IOobject
-            (
-                "faceZones",
-                facesInst,
-                meshSubDir,
-                *this,
-                IOobject::READ_IF_PRESENT,
-                IOobject::NO_WRITE,
-                false
-            ),
-            *this
-        );
-
-        oldSize = faceZones_.size();
-
-        if (newFaceZones.size() <= faceZones_.size())
-        {
-            faceZones_.setSize(newFaceZones.size());
-        }
-
-        // Reset existing ones
-        forAll(faceZones_, fzI)
-        {
-            faceZones_[fzI].resetAddressing
-            (
-                newFaceZones[fzI],
-                newFaceZones[fzI].flipMap()
             );
+
+            pointZones_.swap(newPointZones);
+            pointZones_.instance() = facesInst;
         }
 
-        // Extend with extra ones
+        // faceZones
-        faceZones_.setSize(newFaceZones.size());
-
-        for (label fzI = oldSize; fzI < newFaceZones.size(); fzI++)
         {
-            faceZones_.set(fzI, newFaceZones[fzI].clone(faceZones_));
+            meshFaceZones newFaceZones
-        }
-
-        faceZones_.instance() = facesInst;
-
-
-        meshCellZones newCellZones
-        (
-            IOobject
             (
-                "cellZones",
+                IOobject
-                facesInst,
+                (
-                meshSubDir,
+                    "faceZones",
-                *this,
+                    facesInst,
-                IOobject::READ_IF_PRESENT,
+                    meshSubDir,
-                IOobject::NO_WRITE,
+                    *this,
-                false
+                    IOobject::READ_IF_PRESENT,
-            ),
+                    IOobject::NO_WRITE,
-            *this
+                    false
-        );
+                ),
+                *this
+            );
 
-        oldSize = cellZones_.size();
+            faceZones_.swap(newFaceZones);
-
+            faceZones_.instance() = facesInst;
-        if (newCellZones.size() <= cellZones_.size())
-        {
-            cellZones_.setSize(newCellZones.size());
         }
 
-        // Reset existing ones
+        // cellZones
-        forAll(cellZones_, czI)
         {
-            cellZones_[czI] = newCellZones[czI];
+            meshCellZones newCellZones
+            (
+                IOobject
+                (
+                    "cellZones",
+                    facesInst,
+                    meshSubDir,
+                    *this,
+                    IOobject::READ_IF_PRESENT,
+                    IOobject::NO_WRITE,
+                    false
+                ),
+                *this
+            );
+
+            cellZones_.swap(newCellZones);
+            cellZones_.instance() = facesInst;
         }
 
-        // Extend with extra ones
-        cellZones_.setSize(newCellZones.size());
-
-        for (label czI = oldSize; czI < newCellZones.size(); czI++)
-        {
-            cellZones_.set(czI, newCellZones[czI].clone(cellZones_));
-        }
-
-        cellZones_.instance() = facesInst;
-
-
         // Re-read tet base points
         tetBasePtIsPtr_ = readTetBasePtIs();
 

src/OpenFOAM/meshes/polyMesh/zones/MeshZones/MeshZones.C

@@ -270,35 +270,24 @@ Foam::wordList Foam::MeshZones<ZoneType, MeshType>::names() const
 
 
 template<class ZoneType, class MeshType>
-Foam::labelList Foam::MeshZones<ZoneType, MeshType>::findIndices
+bool Foam::MeshZones<ZoneType, MeshType>::found
 (
-    const wordRe& key
+    const word& zoneName
 ) const
 {
-    labelList indices;
+    if (zoneName != word::null)
-
-    if (!key.empty())
     {
-        if (key.isPattern())
+        forAll(*this, i)
         {
-            indices = findStrings(key, this->names());
+            if (zoneName == operator[](i).name())
-        }
-        else
-        {
-            indices.setSize(this->size());
-            label nFound = 0;
-            forAll(*this, i)
             {
-                if (key == operator[](i).name())
+                return true;
-                {
-                    indices[nFound++] = i;
-                }
             }
-            indices.setSize(nFound);
         }
     }
 
-    return indices;
+    // Not found
+    return false;
 }
 
 
@@ -337,6 +326,39 @@ Foam::label Foam::MeshZones<ZoneType, MeshType>::findIndex
 }
 
 
+template<class ZoneType, class MeshType>
+Foam::labelList Foam::MeshZones<ZoneType, MeshType>::findIndices
+(
+    const wordRe& key
+) const
+{
+    labelList indices;
+
+    if (!key.empty())
+    {
+        if (key.isPattern())
+        {
+            indices = findStrings(key, this->names());
+        }
+        else
+        {
+            indices.setSize(this->size());
+            label nFound = 0;
+            forAll(*this, i)
+            {
+                if (key == operator[](i).name())
+                {
+                    indices[nFound++] = i;
+                }
+            }
+            indices.setSize(nFound);
+        }
+    }
+
+    return indices;
+}
+
+
 template<class ZoneType, class MeshType>
 Foam::label Foam::MeshZones<ZoneType, MeshType>::findZoneID
 (
@@ -384,6 +406,33 @@ Foam::PackedBoolList Foam::MeshZones<ZoneType, MeshType>::findMatching
 }
 
 
+template<class ZoneType, class MeshType>
+void Foam::MeshZones<ZoneType, MeshType>::append
+(
+    const word& zoneName,
+    const labelList& cells
+) const
+{
+    MeshZones<ZoneType, MeshType>& zones =
+        const_cast<MeshZones<ZoneType, MeshType>&>(*this);
+
+    const label zoneID = zones.size();
+    zones.setSize(zoneID + 1);
+
+    zones.set
+    (
+        zoneID,
+        new ZoneType
+        (
+            zoneName,
+            cells,
+            zoneID,
+            *this
+        )
+    );
+}
+
+
 template<class ZoneType, class MeshType>
 void Foam::MeshZones<ZoneType, MeshType>::clearAddressing()
 {

src/OpenFOAM/meshes/polyMesh/zones/MeshZones/MeshZones.H

@@ -143,18 +143,31 @@ public:
         //- Return a list of zone names
         wordList names() const;
 
-        //- Return zone indices for all matches
+        //- Return true if the given zoneName is present
-        labelList findIndices(const wordRe&) const;
+        bool found(const word& zoneName) const;
 
         //- Return zone index for the first match, return -1 if not found
         label findIndex(const wordRe&) const;
 
+        //- Return zone indices for all matches
+        labelList findIndices(const wordRe&) const;
+
         //- Find zone index given a name
         label findZoneID(const word& zoneName) const;
 
         //- Mark cells that match the zone specification
         PackedBoolList findMatching(const wordRe&) const;
 
+        using PtrList<ZoneType>::append;
+
+        //- Append a zone
+        //  Temporary function pending the rewrite of zones
+        void append
+        (
+            const word& zoneName,
+            const labelList& cells
+        ) const;
+
         //- Clear addressing
         void clearAddressing();
 

src/OpenFOAM/meshes/polyMesh/zones/faceZone/faceZone.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -518,6 +518,13 @@ bool Foam::faceZone::checkParallelSync(const bool report) const
 }
 
 
+void Foam::faceZone::swap(faceZone& fz)
+{
+    zone::swap(fz);
+    flipMap_.swap(fz.flipMap_);
+}
+
+
 void Foam::faceZone::movePoints(const pointField& p)
 {
     if (patchPtr_)
@@ -526,6 +533,7 @@ void Foam::faceZone::movePoints(const pointField& p)
     }
 }
 
+
 void Foam::faceZone::write(Ostream& os) const
 {
     os  << nl << name()

src/OpenFOAM/meshes/polyMesh/zones/faceZone/faceZone.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -288,6 +288,9 @@ public:
         //  true if in error.
         virtual bool checkParallelSync(const bool report = false) const;
 
+        //- Swap two faceZones
+        virtual void swap(faceZone&);
+
         //- Correct patch after moving points
         virtual void movePoints(const pointField&);
 

src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.C

@@ -47,20 +47,14 @@ void kEpsilonLopesdaCosta<BasicMomentumTransportModel>::setPorosityCoefficient
 {
     if (pm.dict().found(C.name()))
     {
-        const labelList& cellZoneIDs = pm.cellZoneIDs();
-
         const scalar Cpm = pm.dict().lookup<scalar>(C.name());
 
-        forAll(cellZoneIDs, zonei)
+        const labelList& cells = this->mesh_.cellZones()[pm.zoneName()];
-        {
-            const labelList& cells =
-                this->mesh_.cellZones()[cellZoneIDs[zonei]];
 
-            forAll(cells, i)
+        forAll(cells, i)
-            {
+        {
-                const label celli = cells[i];
+            const label celli = cells[i];
-                C[celli] = Cpm;
+            C[celli] = Cpm;
-            }
         }
     }
 }
@@ -75,21 +69,16 @@ void kEpsilonLopesdaCosta<BasicMomentumTransportModel>::setCdAv
 {
     if (pm.dict().found(C.name()))
     {
-        const labelList& cellZoneIDs = pm.cellZoneIDs();
         const scalarField& Av = pm.Av();
 
         const scalar Cpm = pm.dict().lookup<scalar>(C.name());
 
-        forAll(cellZoneIDs, zonei)
+        const labelList& cells = this->mesh_.cellZones()[pm.zoneName()];
-        {
-            const labelList& cells =
-                this->mesh_.cellZones()[cellZoneIDs[zonei]];
 
-            forAll(cells, i)
+        forAll(cells, i)
-            {
+        {
-                const label celli = cells[i];
+            const label celli = cells[i];
-                C[celli] = Cpm*Av[celli];
+            C[celli] = Cpm*Av[celli];
-            }
         }
     }
 }

src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C

@@ -258,33 +258,15 @@ Foam::porosityModels::powerLawLopesdaCostaZone::powerLawLopesdaCostaZone
     }
 
     // Create the normalised height field
-    scalarField zNorm(zBottom/(zBottom + zTop));
+    const scalarField zNorm(zBottom/(zBottom + zTop));
 
     // Create the porosity surface area per unit volume zone field
     Av_ = AvFunc->value(zNorm);
 
     // Create the porous region cellZone and add to the mesh cellZones
-
+    if (!mesh.cellZones().found(zoneName_))
-    meshCellZones& cellZones = const_cast<meshCellZones&>(mesh.cellZones());
-
-    label zoneID = cellZones.findZoneID(zoneName_);
-
-    if (zoneID == -1)
     {
-        zoneID = cellZones.size();
+        mesh.cellZones().append(zoneName_, porousCells);
-        cellZones.setSize(zoneID + 1);
-
-        cellZones.set
-        (
-            zoneID,
-            new cellZone
-            (
-                zoneName_,
-                porousCells,
-                zoneID,
-                cellZones
-            )
-        );
     }
     else
     {

src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCostaTemplates.C

@@ -38,18 +38,16 @@ void Foam::porosityModels::powerLawLopesdaCosta::apply
 {
     const scalar C1m1b2 = (C1_ - 1.0)/2.0;
 
-    forAll(cellZoneIDs_, zonei)
+    const labelList& cells =
+        mesh_.cellZones()[powerLawLopesdaCostaZone::zoneName_];
+
+    forAll(cells, i)
     {
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zonei]];
+        const label celli = cells[i];
 
-        forAll(cells, i)
+        Udiag[celli] +=
-        {
+            V[celli]*rho[celli]
-            const label celli = cells[i];
+           *Cd_*Av_[i]*pow(magSqr(U[celli]), C1m1b2);
-
-            Udiag[celli] +=
-                V[celli]*rho[celli]
-               *Cd_*Av_[i]*pow(magSqr(U[celli]), C1m1b2);
-        }
     }
 }
 
@@ -64,22 +62,20 @@ void Foam::porosityModels::powerLawLopesdaCosta::apply
 {
     const scalar C1m1b2 = (C1_ - 1.0)/2.0;
 
-    forAll(cellZoneIDs_, zonei)
+    const labelList& cells =
+        mesh_.cellZones()[powerLawLopesdaCostaZone::zoneName_];
+
+    forAll(cells, i)
     {
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zonei]];
+        const label celli = cells[i];
 
-        forAll(cells, i)
+        AU[celli] =
-        {
+            AU[celli]
-            const label celli = cells[i];
+          + I
-
+           *(
-            AU[celli] =
+               0.5*rho[celli]*Cd_*Av_[i]
-                AU[celli]
+              *pow(magSqr(U[celli]), C1m1b2)
-              + I
+            );
-               *(
-                   0.5*rho[celli]*Cd_*Av_[i]
-                  *pow(magSqr(U[celli]), C1m1b2)
-                );
-        }
     }
 }
 

src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2012-2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -54,8 +54,6 @@ Foam::porosityModels::DarcyForchheimer::DarcyForchheimer
     porosityModel(name, modelType, mesh, dict, cellZoneName),
     dXYZ_("d", dimless/sqr(dimLength), coeffs_),
     fXYZ_("f", dimless/dimLength, coeffs_),
-    D_(cellZoneIDs_.size()),
-    F_(cellZoneIDs_.size()),
     rhoName_(coeffs_.lookupOrDefault<word>("rho", "rho")),
     muName_(coeffs_.lookupOrDefault<word>("mu", "mu")),
     nuName_(coeffs_.lookupOrDefault<word>("nu", "nu"))
@@ -79,58 +77,52 @@ void Foam::porosityModels::DarcyForchheimer::calcTransformModelData()
 {
     if (coordSys_.R().uniform())
     {
-        forAll(cellZoneIDs_, zoneI)
+        D_.setSize(1);
-        {
+        F_.setSize(1);
-            D_[zoneI].setSize(1);
-            F_[zoneI].setSize(1);
 
-            D_[zoneI][0] = Zero;
+        D_[0] = Zero;
-            D_[zoneI][0].xx() = dXYZ_.value().x();
+        D_[0].xx() = dXYZ_.value().x();
-            D_[zoneI][0].yy() = dXYZ_.value().y();
+        D_[0].yy() = dXYZ_.value().y();
-            D_[zoneI][0].zz() = dXYZ_.value().z();
+        D_[0].zz() = dXYZ_.value().z();
 
-            D_[zoneI][0] = coordSys_.R().transform(Zero, D_[zoneI][0]);
+        D_[0] = coordSys_.R().transform(Zero, D_[0]);
 
-            // leading 0.5 is from 1/2*rho
+        // leading 0.5 is from 1/2*rho
-            F_[zoneI][0] = Zero;
+        F_[0] = Zero;
-            F_[zoneI][0].xx() = 0.5*fXYZ_.value().x();
+        F_[0].xx() = 0.5*fXYZ_.value().x();
-            F_[zoneI][0].yy() = 0.5*fXYZ_.value().y();
+        F_[0].yy() = 0.5*fXYZ_.value().y();
-            F_[zoneI][0].zz() = 0.5*fXYZ_.value().z();
+        F_[0].zz() = 0.5*fXYZ_.value().z();
 
-            F_[zoneI][0] = coordSys_.R().transform(Zero, F_[zoneI][0]);
+        F_[0] = coordSys_.R().transform(Zero, F_[0]);
-        }
     }
     else
     {
-        forAll(cellZoneIDs_, zoneI)
+        const labelList& cells = mesh_.cellZones()[zoneName_];
+
+        D_.setSize(cells.size());
+        F_.setSize(cells.size());
+
+        forAll(cells, i)
         {
-            const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+            D_[i] = Zero;
+            D_[i].xx() = dXYZ_.value().x();
+            D_[i].yy() = dXYZ_.value().y();
+            D_[i].zz() = dXYZ_.value().z();
 
-            D_[zoneI].setSize(cells.size());
+            // leading 0.5 is from 1/2*rho
-            F_[zoneI].setSize(cells.size());
+            F_[i] = Zero;
-
+            F_[i].xx() = 0.5*fXYZ_.value().x();
-            forAll(cells, i)
+            F_[i].yy() = 0.5*fXYZ_.value().y();
-            {
+            F_[i].zz() = 0.5*fXYZ_.value().z();
-                D_[zoneI][i] = Zero;
-                D_[zoneI][i].xx() = dXYZ_.value().x();
-                D_[zoneI][i].yy() = dXYZ_.value().y();
-                D_[zoneI][i].zz() = dXYZ_.value().z();
-
-                // leading 0.5 is from 1/2*rho
-                F_[zoneI][i] = Zero;
-                F_[zoneI][i].xx() = 0.5*fXYZ_.value().x();
-                F_[zoneI][i].yy() = 0.5*fXYZ_.value().y();
-                F_[zoneI][i].zz() = 0.5*fXYZ_.value().z();
-            }
-
-            const coordinateRotation& R = coordSys_.R
-            (
-                UIndirectList<vector>(mesh_.C(), cells)()
-            );
-
-            D_[zoneI] = R.transform(D_[zoneI]);
-            F_[zoneI] = R.transform(F_[zoneI]);
         }
+
+        const coordinateRotation& R = coordSys_.R
+        (
+            UIndirectList<vector>(mesh_.C(), cells)()
+        );
+
+        D_ = R.transform(D_);
+        F_ = R.transform(F_);
     }
 
     if (debug && (mesh_.time().writeTime() || mesh_.time().timeIndex() == 0))
@@ -162,8 +154,8 @@ void Foam::porosityModels::DarcyForchheimer::calcTransformModelData()
             dimensionedTensor(fXYZ_.dimensions(), Zero)
         );
 
-        UIndirectList<tensor>(Dout, mesh_.cellZones()[cellZoneIDs_[0]]) = D_[0];
+        UIndirectList<tensor>(Dout, mesh_.cellZones()[zoneName_]) = D_;
-        UIndirectList<tensor>(Fout, mesh_.cellZones()[cellZoneIDs_[0]]) = F_[0];
+        UIndirectList<tensor>(Fout, mesh_.cellZones()[zoneName_]) = F_;
 
         Dout.write();
         Fout.write();

src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H

@@ -79,10 +79,10 @@ class DarcyForchheimer
         dimensionedVector fXYZ_;
 
         //- Darcy coefficient - converted from dXYZ [1/m^2]
-        List<tensorField> D_;
+        tensorField D_;
 
         //- Forchheimer coefficient - converted from fXYZ [1/m]
-        List<tensorField> F_;
+        tensorField F_;
 
         //- Name of density field
         word rhoName_;

src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2012-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -36,25 +36,19 @@ void Foam::porosityModels::DarcyForchheimer::apply
     const vectorField& U
 ) const
 {
-    forAll(cellZoneIDs_, zoneI)
+    const labelList& cells = mesh_.cellZones()[zoneName_];
+
+    forAll(cells, i)
     {
-        const tensorField& dZones = D_[zoneI];
+        const label celli = cells[i];
-        const tensorField& fZones = F_[zoneI];
+        const label j = this->fieldIndex(i);
+        const tensor Cd =
+            mu[celli]*D_[j] + (rho[celli]*mag(U[celli]))*F_[j];
 
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+        const scalar isoCd = tr(Cd);
 
-        forAll(cells, i)
+        Udiag[celli] += V[celli]*isoCd;
-        {
+        Usource[celli] -= V[celli]*((Cd - I*isoCd) & U[celli]);
-            const label celli = cells[i];
-            const label j = this->fieldIndex(i);
-            const tensor Cd =
-                mu[celli]*dZones[j] + (rho[celli]*mag(U[celli]))*fZones[j];
-
-            const scalar isoCd = tr(Cd);
-
-            Udiag[celli] += V[celli]*isoCd;
-            Usource[celli] -= V[celli]*((Cd - I*isoCd) & U[celli]);
-        }
     }
 }
 
@@ -68,22 +62,16 @@ void Foam::porosityModels::DarcyForchheimer::apply
     const vectorField& U
 ) const
 {
-    forAll(cellZoneIDs_, zoneI)
+    const labelList& cells = mesh_.cellZones()[zoneName_];
+
+    forAll(cells, i)
     {
-        const tensorField& dZones = D_[zoneI];
+        const label celli = cells[i];
-        const tensorField& fZones = F_[zoneI];
+        const label j = this->fieldIndex(i);
+        const tensor D = D_[j];
+        const tensor F = F_[j];
 
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+        AU[celli] += mu[celli]*D + (rho[celli]*mag(U[celli]))*F;
-
-        forAll(cells, i)
-        {
-            const label celli = cells[i];
-            const label j = this->fieldIndex(i);
-            const tensor D = dZones[j];
-            const tensor F = fZones[j];
-
-            AU[celli] += mu[celli]*D + (rho[celli]*mag(U[celli]))*F;
-        }
     }
 }
 

src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2012-2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -50,23 +50,17 @@ void Foam::porosityModels::fixedCoeff::apply
     const scalar rho
 ) const
 {
-    forAll(cellZoneIDs_, zoneI)
+    const labelList& cells = mesh_.cellZones()[zoneName_];
+
+    forAll(cells, i)
     {
-        const tensorField& alphaZones = alpha_[zoneI];
+        const label celli = cells[i];
-        const tensorField& betaZones = beta_[zoneI];
+        const label j = fieldIndex(i);
+        const tensor Cd = rho*(alpha_[j] + beta_[j]*mag(U[celli]));
+        const scalar isoCd = tr(Cd);
 
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+        Udiag[celli] += V[celli]*isoCd;
-
+        Usource[celli] -= V[celli]*((Cd - I*isoCd) & U[celli]);
-        forAll(cells, i)
-        {
-            const label celli = cells[i];
-            const label j = fieldIndex(i);
-            const tensor Cd = rho*(alphaZones[j] + betaZones[j]*mag(U[celli]));
-            const scalar isoCd = tr(Cd);
-
-            Udiag[celli] += V[celli]*isoCd;
-            Usource[celli] -= V[celli]*((Cd - I*isoCd) & U[celli]);
-        }
     }
 }
 
@@ -78,23 +72,16 @@ void Foam::porosityModels::fixedCoeff::apply
     const scalar rho
 ) const
 {
+    const labelList& cells = mesh_.cellZones()[zoneName_];
 
-    forAll(cellZoneIDs_, zoneI)
+    forAll(cells, i)
     {
-        const tensorField& alphaZones = alpha_[zoneI];
+        const label celli = cells[i];
-        const tensorField& betaZones = beta_[zoneI];
+        const label j = fieldIndex(i);
+        const tensor alpha = alpha_[j];
+        const tensor beta = beta_[j];
 
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+        AU[celli] += rho*(alpha + beta*mag(U[celli]));
-
-        forAll(cells, i)
-        {
-            const label celli = cells[i];
-            const label j = fieldIndex(i);
-            const tensor alpha = alphaZones[j];
-            const tensor beta = betaZones[j];
-
-            AU[celli] += rho*(alpha + beta*mag(U[celli]));
-        }
     }
 }
 
@@ -112,9 +99,7 @@ Foam::porosityModels::fixedCoeff::fixedCoeff
 :
     porosityModel(name, modelType, mesh, dict, cellZoneName),
     alphaXYZ_("alpha", dimless/dimTime, coeffs_),
-    betaXYZ_("beta", dimless/dimLength, coeffs_),
+    betaXYZ_("beta", dimless/dimLength, coeffs_)
-    alpha_(cellZoneIDs_.size()),
-    beta_(cellZoneIDs_.size())
 {
     adjustNegativeResistance(alphaXYZ_);
     adjustNegativeResistance(betaXYZ_);
@@ -135,54 +120,48 @@ void Foam::porosityModels::fixedCoeff::calcTransformModelData()
 {
     if (coordSys_.R().uniform())
     {
-        forAll(cellZoneIDs_, zoneI)
+        alpha_.setSize(1);
-        {
+        beta_.setSize(1);
-            alpha_[zoneI].setSize(1);
-            beta_[zoneI].setSize(1);
 
-            alpha_[zoneI][0] = Zero;
+        alpha_[0] = Zero;
-            alpha_[zoneI][0].xx() = alphaXYZ_.value().x();
+        alpha_[0].xx() = alphaXYZ_.value().x();
-            alpha_[zoneI][0].yy() = alphaXYZ_.value().y();
+        alpha_[0].yy() = alphaXYZ_.value().y();
-            alpha_[zoneI][0].zz() = alphaXYZ_.value().z();
+        alpha_[0].zz() = alphaXYZ_.value().z();
-            alpha_[zoneI][0] = coordSys_.R().transform(Zero, alpha_[zoneI][0]);
+        alpha_[0] = coordSys_.R().transform(Zero, alpha_[0]);
 
-            beta_[zoneI][0] = Zero;
+        beta_[0] = Zero;
-            beta_[zoneI][0].xx() = betaXYZ_.value().x();
+        beta_[0].xx() = betaXYZ_.value().x();
-            beta_[zoneI][0].yy() = betaXYZ_.value().y();
+        beta_[0].yy() = betaXYZ_.value().y();
-            beta_[zoneI][0].zz() = betaXYZ_.value().z();
+        beta_[0].zz() = betaXYZ_.value().z();
-            beta_[zoneI][0] = coordSys_.R().transform(Zero, beta_[zoneI][0]);
+        beta_[0] = coordSys_.R().transform(Zero, beta_[0]);
-        }
     }
     else
     {
-        forAll(cellZoneIDs_, zoneI)
+        const labelList& cells = mesh_.cellZones()[zoneName_];
+
+        alpha_.setSize(cells.size());
+        beta_.setSize(cells.size());
+
+        forAll(cells, i)
         {
-            const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+            alpha_[i] = Zero;
+            alpha_[i].xx() = alphaXYZ_.value().x();
+            alpha_[i].yy() = alphaXYZ_.value().y();
+            alpha_[i].zz() = alphaXYZ_.value().z();
 
-            alpha_[zoneI].setSize(cells.size());
+            beta_[i] = Zero;
-            beta_[zoneI].setSize(cells.size());
+            beta_[i].xx() = betaXYZ_.value().x();
-
+            beta_[i].yy() = betaXYZ_.value().y();
-            forAll(cells, i)
+            beta_[i].zz() = betaXYZ_.value().z();
-            {
-                alpha_[zoneI][i] = Zero;
-                alpha_[zoneI][i].xx() = alphaXYZ_.value().x();
-                alpha_[zoneI][i].yy() = alphaXYZ_.value().y();
-                alpha_[zoneI][i].zz() = alphaXYZ_.value().z();
-
-                beta_[zoneI][i] = Zero;
-                beta_[zoneI][i].xx() = betaXYZ_.value().x();
-                beta_[zoneI][i].yy() = betaXYZ_.value().y();
-                beta_[zoneI][i].zz() = betaXYZ_.value().z();
-            }
-
-            const coordinateRotation& R = coordSys_.R
-            (
-                UIndirectList<vector>(mesh_.C(), cells)()
-            );
-
-            alpha_[zoneI] = R.transform(alpha_[zoneI]);
-            beta_[zoneI] = R.transform(beta_[zoneI]);
         }
+
+        const coordinateRotation& R = coordSys_.R
+        (
+            UIndirectList<vector>(mesh_.C(), cells)()
+        );
+
+        alpha_ = R.transform(alpha_);
+        beta_ = R.transform(beta_);
     }
 }
 

src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2012-2021 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -69,10 +69,10 @@ class fixedCoeff
         dimensionedVector betaXYZ_;
 
         //- Model alpha coefficient - converted from alphaXYZ [1/s]
-        List<tensorField> alpha_;
+        tensorField alpha_;
 
         //- Model beta coefficient - converted from betaXYZ [1/m]
-        List<tensorField> beta_;
+        tensorField beta_;
 
 
     // Private Member Functions

src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2012-2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -100,23 +100,17 @@ Foam::porosityModel::porosityModel
     mesh_(mesh),
     dict_(dict),
     coeffs_(dict.optionalSubDict(modelType + "Coeffs")),
-    zoneName_(cellZoneName),
+    zoneName_
-    cellZoneIDs_(),
+    (
+        cellZoneName != word::null
+      ? cellZoneName
+      : dict_.lookup<word>("cellZone")
+    ),
     coordSys_(coordinateSystem::New(mesh, coeffs_))
 {
-    if (zoneName_ == word::null)
-    {
-        dict_.lookup("cellZone") >> zoneName_;
-    }
-
-    cellZoneIDs_ = mesh_.cellZones().findIndices(zoneName_);
-
     Info<< "    creating porous zone: " << zoneName_ << endl;
 
-    bool foundZone = !cellZoneIDs_.empty();
+    if (mesh_.cellZones().findIndex(zoneName_) == -1)
-    reduce(foundZone, orOp<bool>());
-
-    if (!foundZone && Pstream::master())
     {
         FatalErrorInFunction
             << "cannot find porous cellZone " << zoneName_
@@ -155,11 +149,7 @@ Foam::tmp<Foam::vectorField> Foam::porosityModel::porosityModel::force
     const_cast<porosityModel&>(*this).transformModelData();
 
     tmp<vectorField> tforce(new vectorField(U.size(), Zero));
-
+    this->calcForce(U, rho, mu, tforce.ref());
-    if (!cellZoneIDs_.empty())
-    {
-        this->calcForce(U, rho, mu, tforce.ref());
-    }
 
     return tforce;
 }
@@ -167,11 +157,6 @@ Foam::tmp<Foam::vectorField> Foam::porosityModel::porosityModel::force
 
 void Foam::porosityModel::addResistance(fvVectorMatrix& UEqn)
 {
-    if (cellZoneIDs_.empty())
-    {
-        return;
-    }
-
     transformModelData();
     this->correct(UEqn);
 }
@@ -184,19 +169,14 @@ void Foam::porosityModel::addResistance
     bool correctAUprocBC
 )
 {
-    if (cellZoneIDs_.empty())
-    {
-        return;
-    }
-
     transformModelData();
     this->correct(UEqn, AU);
 
     if (correctAUprocBC)
     {
-        // Correct the boundary conditions of the tensorial diagonal to ensure
+        // Correct the boundary conditions of the tensorial diagonal to
-        // processor boundaries are correctly handled when AU^-1 is interpolated
+        // ensure processor boundaries are correctly handled when AU^-1 is
-        // for the pressure equation.
+        // interpolated for the pressure equation.
         AU.correctBoundaryConditions();
     }
 }
@@ -213,7 +193,6 @@ bool Foam::porosityModel::read(const dictionary& dict)
     coeffs_ = dict.optionalSubDict(type() + "Coeffs");
 
     dict.lookup("cellZone") >> zoneName_;
-    cellZoneIDs_ = mesh_.cellZones().findIndices(zoneName_);
 
     return true;
 }

src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2012-2021 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -81,11 +81,8 @@ protected:
         //- Model coefficients dictionary
         dictionary coeffs_;
 
-        //- Name(s) of cell-zone
+        //- Name of cellZone
-        wordRe zoneName_;
+        word zoneName_;
-
-        //- Cell zone IDs
-        labelList cellZoneIDs_;
 
         //- Local co-ordinate system
         coordinateSystem coordSys_;
@@ -205,8 +202,8 @@ public:
         //- Return const access to the porosity model name
         inline const word& name() const;
 
-        //- Return const access to the cell zone IDs
+        //- Return const access to the cell zone name
-        inline const labelList& cellZoneIDs() const;
+        inline const word& zoneName() const;
 
         //- Return dictionary used for model construction
         const dictionary& dict() const;

src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2012-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,9 +35,9 @@ inline const Foam::dictionary& Foam::porosityModel::dict() const
 }
 
 
-inline const Foam::labelList& Foam::porosityModel::cellZoneIDs() const
+inline const Foam::word& Foam::porosityModel::zoneName() const
 {
-    return cellZoneIDs_;
+    return zoneName_;
 }
 
 

src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLawTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2012-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -34,20 +34,17 @@ void Foam::porosityModels::powerLaw::apply
     const vectorField& U
 ) const
 {
+    const labelList& cells = mesh_.cellZones()[zoneName_];
+
     const scalar C0 = C0_;
     const scalar C1m1b2 = (C1_ - 1.0)/2.0;
 
-    forAll(cellZoneIDs_, zoneI)
+    forAll(cells, i)
     {
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+        const label celli = cells[i];
 
-        forAll(cells, i)
+        Udiag[celli] +=
-        {
+            V[celli]*rho[celli]*C0*pow(magSqr(U[celli]), C1m1b2);
-            const label celli = cells[i];
-
-            Udiag[celli] +=
-                V[celli]*rho[celli]*C0*pow(magSqr(U[celli]), C1m1b2);
-        }
     }
 }
 
@@ -60,20 +57,17 @@ void Foam::porosityModels::powerLaw::apply
     const vectorField& U
 ) const
 {
+    const labelList& cells = mesh_.cellZones()[zoneName_];
+
     const scalar C0 = C0_;
     const scalar C1m1b2 = (C1_ - 1.0)/2.0;
 
-    forAll(cellZoneIDs_, zoneI)
+    forAll(cells, i)
     {
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+        const label celli = cells[i];
 
-        forAll(cells, i)
+        AU[celli] =
-        {
+            AU[celli] + I*(rho[celli]*C0*pow(magSqr(U[celli]), C1m1b2));
-            const label celli = cells[i];
-
-            AU[celli] =
-                AU[celli] + I*(rho[celli]*C0*pow(magSqr(U[celli]), C1m1b2));
-        }
     }
 }
 

src/finiteVolume/cfdTools/general/porosityModel/solidification/solidificationTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2017-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2017-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -38,21 +38,18 @@ void Foam::porosityModels::solidification::apply
     const volVectorField& U
 ) const
 {
+    const labelList& cells = mesh_.cellZones()[zoneName_];
+
     const volScalarField& T = mesh_.lookupObject<volScalarField>
     (
         IOobject::groupName(TName_, U.group())
     );
 
-    forAll(cellZoneIDs_, zoneI)
+    forAll(cells, i)
     {
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+        const label celli = cells[i];
-
+        Udiag[celli] +=
-        forAll(cells, i)
+            V[celli]*alpha[celli]*rho[celli]*D_->value(T[celli]);
-        {
-            const label celli = cells[i];
-            Udiag[celli] +=
-                V[celli]*alpha[celli]*rho[celli]*D_->value(T[celli]);
-        }
     }
 }
 
@@ -66,21 +63,18 @@ void Foam::porosityModels::solidification::apply
     const volVectorField& U
 ) const
 {
+    const labelList& cells = mesh_.cellZones()[zoneName_];
+
     const volScalarField& T = mesh_.lookupObject<volScalarField>
     (
         IOobject::groupName(TName_, U.group())
     );
 
-    forAll(cellZoneIDs_, zoneI)
+    forAll(cells, i)
     {
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+        const label celli = cells[i];
-
+        AU[celli] +=
-        forAll(cells, i)
+            tensor::I*alpha[celli]*rho[celli]*D_->value(T[celli]);
-        {
-            const label celli = cells[i];
-            AU[celli] +=
-                tensor::I*alpha[celli]*rho[celli]*D_->value(T[celli]);
-        }
     }
 }
 

src/functionObjects/forces/forces/forces.C

@@ -991,21 +991,14 @@ void Foam::functionObjects::forces::calcForcesMoment()
 
             const vectorField fPTot(pm.force(U, rho, mu));
 
-            const labelList& cellZoneIDs = pm.cellZoneIDs();
+            const cellZone& cZone = mesh_.cellZones()[pm.zoneName()];
+            const vectorField d(mesh_.C(), cZone);
+            const vectorField fP(fPTot, cZone);
+            const vectorField Md(d - coordSys_.origin());
 
-            forAll(cellZoneIDs, i)
+            const vectorField fDummy(Md.size(), Zero);
-            {
-                const label zoneI = cellZoneIDs[i];
-                const cellZone& cZone = mesh_.cellZones()[zoneI];
 
-                const vectorField d(mesh_.C(), cZone);
+            applyBins(Md, fDummy, fDummy, fP, d);
-                const vectorField fP(fPTot, cZone);
-                const vectorField Md(d - coordSys_.origin());
-
-                const vectorField fDummy(Md.size(), Zero);
-
-                applyBins(Md, fDummy, fDummy, fP, d);
-            }
         }
     }
 

src/fvModels/derived/porosityForce/porosityForce.C

@@ -24,9 +24,8 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "porosityForce.H"
-#include "fvMesh.H"
-#include "fvMatrices.H"
 #include "porosityModel.H"
+#include "fvMatrices.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -51,6 +50,20 @@ namespace fv
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+void Foam::fv::porosityForce::reset()
+{
+    porosityPtr_.reset
+    (
+        porosityModel::New
+        (
+            name(),
+            mesh(),
+            coeffs()
+        ).ptr()
+    );
+}
+
+
 void Foam::fv::porosityForce::readCoeffs()
 {
     if (coeffs().found("UNames"))
@@ -62,20 +75,10 @@ void Foam::fv::porosityForce::readCoeffs()
         UNames_ = wordList(1, coeffs().lookupOrDefault<word>("U", "U"));
     }
 
-    porosityPtr_.reset
+    reset();
-    (
-        porosityModel::New
-        (
-            name(),
-            mesh(),
-            coeffs(),
-            set_.cellSetName()
-        ).ptr()
-    );
 }
 
 
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::fv::porosityForce::porosityForce
@@ -87,7 +90,6 @@ Foam::fv::porosityForce::porosityForce
 )
 :
     fvModel(name, modelType, mesh, dict),
-    set_(mesh, coeffs()),
     UNames_(),
     porosityPtr_(nullptr)
 {
@@ -144,20 +146,20 @@ void Foam::fv::porosityForce::addSup
 
 bool Foam::fv::porosityForce::movePoints()
 {
-    set_.movePoints();
+    // Currently there is no mechanism to update the porous media orientation
     return true;
 }
 
 
 void Foam::fv::porosityForce::topoChange(const polyTopoChangeMap& map)
 {
-    set_.topoChange(map);
+    reset();
 }
 
 
 void Foam::fv::porosityForce::mapMesh(const polyMeshMap& map)
 {
-    set_.mapMesh(map);
+    reset();
 }
 
 
@@ -166,7 +168,7 @@ void Foam::fv::porosityForce::distribute
     const polyDistributionMap& map
 )
 {
-    set_.distribute(map);
+    reset();
 }
 
 
@@ -174,7 +176,6 @@ bool Foam::fv::porosityForce::read(const dictionary& dict)
 {
     if (fvModel::read(dict))
     {
-        set_.read(coeffs());
         readCoeffs();
         return true;
     }

src/fvModels/derived/porosityForce/porosityForce.H

@@ -66,7 +66,6 @@ SourceFiles
 #define porosityForce_H
 
 #include "fvModel.H"
-#include "fvCellSet.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -78,7 +77,6 @@ class porosityModel;
 namespace fv
 {
 
-
 /*---------------------------------------------------------------------------*\
                    Class porosityForce Declaration
 \*---------------------------------------------------------------------------*/
@@ -89,9 +87,6 @@ class porosityForce
 {
     // Private data
 
-        //- The set of cells the fvConstraint applies to
-        fvCellSet set_;
-
         //- Names of the velocity fields
         wordList UNames_;
 
@@ -101,6 +96,10 @@ class porosityForce
 
     // Private Member Functions
 
+        //- Reset the porosity model
+        //  during construction, re-reading and mesh-change
+        void reset();
+
         //- Non-virtual read
         void readCoeffs();
 

src/fvModels/interRegion/interRegionPorosityForce/interRegionPorosityForce.C

@@ -88,17 +88,8 @@ Foam::fv::interRegionPorosityForce::interRegionPorosityForce
 
     const word zoneName(name + ":porous");
 
-    const meshCellZones& cellZones = mesh.cellZones();
+    if (!mesh.cellZones().found(zoneName))
-    label zoneID = cellZones.findZoneID(zoneName);
-
-    if (zoneID == -1)
     {
-        meshCellZones& cz = const_cast<meshCellZones&>(cellZones);
-
-        zoneID = cz.size();
-
-        cz.setSize(zoneID + 1);
-
         // Scan the porous region filter for all cells containing porosity
         labelList porousCells(mesh.nCells());
 
@@ -112,19 +103,7 @@ Foam::fv::interRegionPorosityForce::interRegionPorosityForce
         }
         porousCells.setSize(i);
 
-        cz.set
+        mesh.cellZones().append(zoneName, porousCells);
-        (
-            zoneID,
-            new cellZone
-            (
-                zoneName,
-                porousCells,
-                zoneID,
-                cellZones
-            )
-        );
-
-        cz.clearAddressing();
     }
     else
     {

src/meshTools/cutPoly/cutPolyIsoSurface.C

@@ -44,7 +44,7 @@ Foam::cutPolyIsoSurface::cutPolyIsoSurface
     const polyMesh& mesh,
     const scalarField& pAlphas,
     const scalar isoAlpha,
-    const labelList& zoneIDs
+    const word& cellZoneName
 )
 :
     points_(),
@@ -89,14 +89,13 @@ Foam::cutPolyIsoSurface::cutPolyIsoSurface
 
         nCutCells += !cellCuts[celli].empty();
     };
-    if (!isNull<labelList>(zoneIDs))
+
+    if (cellZoneName != word::null)
     {
-        forAll(zoneIDs, i)
+        const labelList& zoneCells = mesh.cellZones()[cellZoneName];
+        forAll(zoneCells, zoneCelli)
         {
-            forAll(mesh.cellZones()[zoneIDs[i]], zoneCelli)
+            cutCell(zoneCells[zoneCelli]);
-            {
-                cutCell(mesh.cellZones()[zoneIDs[i]][zoneCelli]);
-            }
         }
     }
     else

src/meshTools/cutPoly/cutPolyIsoSurface.H

@@ -81,7 +81,7 @@ public:
             const polyMesh& mesh,
             const scalarField& pAlphas,
             const scalar isoAlpha,
-            const labelList& zoneIDs = NullObjectRef<labelList>()
+            const word& cellZoneName = word::null
         );
 
         //- Construct by combining a list of iso-surfaces

src/meshTools/sets/faceSources/zoneToFace/zoneToFace.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -25,7 +25,7 @@ Class
     Foam::zoneToFace
 
 Description
-    A topoSetSource to select faces based on faceZone.
+    A topoSetSource to select faces from a faceZone.
 
 SourceFiles
     zoneToFace.C
@@ -51,11 +51,9 @@ class zoneToFace
 :
     public topoSetSource
 {
-
     // Private Data
 
-
+        //- Name/regular expression of faceZone
-        //- Name/regular expression of cellZone
         wordRe zoneName_;
 
 
@@ -102,7 +100,6 @@ public:
             const topoSetSource::setAction action,
             topoSet&
         ) const;
-
 };
 
 

src/sampling/sampledSurface/sampledCutPlane/sampledCutPlane.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -66,7 +66,7 @@ Foam::sampledSurfaces::cutPlane::calcIsoSurf() const
     // Construct an iso-surface at the given distance
     return autoPtr<cutPolyIsoSurface>
     (
-        new cutPolyIsoSurface(mesh(), pointDistance, 0, zoneIDs())
+        new cutPolyIsoSurface(mesh(), pointDistance, 0, zoneName())
     );
 }
 

src/sampling/sampledSurface/sampledDistanceSurface/sampledDistanceSurface.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -90,7 +90,7 @@ Foam::sampledSurfaces::distanceSurface::calcIsoSurf() const
     // Construct an iso-surface at the given distance
     return autoPtr<cutPolyIsoSurface>
     (
-        new cutPolyIsoSurface(mesh(), pointDistance, distance_, zoneIDs())
+        new cutPolyIsoSurface(mesh(), pointDistance, distance_, zoneName())
     );
 }
 

src/sampling/sampledSurface/sampledIsoSurface/sampledIsoSurface.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -60,7 +60,7 @@ Foam::sampledSurfaces::isoSurface::calcIsoSurf() const
         isoSurfs.set
         (
             i,
-            new cutPolyIsoSurface(mesh(), pField, isoValues_[i], zoneIDs())
+            new cutPolyIsoSurface(mesh(), pField, isoValues_[i], zoneName())
         );
     }
 

src/sampling/sampledSurface/sampledIsoSurfaceSurface/sampledIsoSurfaceSurface.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2022-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -45,12 +45,11 @@ Foam::sampledSurfaces::sampledIsoSurfaceSurface::sampledIsoSurfaceSurface
 )
 :
     sampledSurface(name, mesh, dict),
-    zoneName_(dict.lookupOrDefault("zone", wordRe::null)),
+    zoneName_(dict.lookupOrDefault("zone", word::null)),
-    zoneIDs_(mesh.cellZones().findIndices(zoneName_)),
     isoSurfPtr_(nullptr),
     isoSurfTimeIndex_(-1)
 {
-    if (zoneName_ != wordRe::null && zoneIDs_.empty())
+    if (zoneName_ != word::null && !mesh.cellZones().found(zoneName_))
     {
         WarningInFunction
             << "Cell zone " << zoneName_

src/sampling/sampledSurface/sampledIsoSurfaceSurface/sampledIsoSurfaceSurface.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2022-2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -57,11 +57,8 @@ private:
 
     // Private Data
 
-        //- If restricted to zones, name of this zone or a regular expression
+        //- Optional name of the cellZone
-        const wordRe zoneName_;
+        const word zoneName_;
-
-        //- If restricted to zones, the indices of the zones
-        const labelList zoneIDs_;
 
         //- Constructed iso surface
         mutable autoPtr<cutPolyIsoSurface> isoSurfPtr_;
@@ -91,10 +88,10 @@ protected:
 
     // Protected Member Functions
 
-        //- Access the zone indices
+        //- Access the zoneName
-        inline const labelList& zoneIDs() const
+        inline const word& zoneName() const
         {
-            return zoneIDs_.empty() ? NullObjectRef<labelList>() : zoneIDs_;
+            return zoneName_;
         }
 
         //- Access the time index