diff --git a/src/finiteVolume/fvMesh/zoneDistribute/zoneDistribute.C b/src/finiteVolume/fvMesh/zoneDistribute/zoneDistribute.C
index 0aa301237d..b52943570c 100644
--- a/src/finiteVolume/fvMesh/zoneDistribute/zoneDistribute.C
+++ b/src/finiteVolume/fvMesh/zoneDistribute/zoneDistribute.C
@@ -27,6 +27,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "zoneDistribute.H"
+#include "processorPolyPatch.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -44,10 +45,48 @@ Foam::zoneDistribute::zoneDistribute(const fvMesh& mesh)
     stencil_(zoneCPCStencil::New(mesh)),
     globalNumbering_(stencil_.globalNumbering()),
     send_(UPstream::nProcs()),
-    pBufs_(UPstream::commsTypes::nonBlocking)
+    pBufs_(UPstream::commsTypes::nonBlocking),
+    cyclicBoundaryCells_(mesh.nCells(), false)
 {
     // Don't clear storage on persistent buffer
     pBufs_.allowClearRecv(false);
+
+    // Loop over boundary patches and store cells with a face on a cyclic patch
+    bool hasCyclicPatches = false;
+    forAll(mesh.boundaryMesh(), patchi)
+    {
+        const cyclicPolyPatch* cpp =
+            isA<cyclicPolyPatch>(mesh.boundaryMesh()[patchi]);
+
+        if (cpp)
+        {
+            cyclicBoundaryCells_.set(cpp->faceCells());
+            hasCyclicPatches = true;
+        }
+    }
+
+    // Populate cyclicCentres_
+    if(hasCyclicPatches)
+    {
+        // Make a boolList from the bitSet
+        boolList isCyclicCell(mesh.nCells(), false);
+
+        forAll(cyclicBoundaryCells_, celli)
+        {
+            if (cyclicBoundaryCells_.test(celli))
+            {
+                isCyclicCell[celli] = true;
+            }
+        }
+
+        // Use getFields to get map of cell centres across processor boundaries
+        setUpCommforZone(isCyclicCell, true);
+
+        cyclicCentres_.reset
+        (
+            new Map<vectorField>(getFields(isCyclicCell, mesh_.C()))
+        );
+    }
 }
 
 
@@ -141,5 +180,226 @@ void Foam::zoneDistribute::setUpCommforZone
     }
 }
 
+Foam::List<Foam::label> Foam::zoneDistribute::getCyclicPatches
+(
+    const label celli,
+    const label globalIdx
+) const
+{
+    // Initialise cyclic patch label list
+    List<label> patches(0);
+
+    // If celli is not on a cyclic boundary, return the empty list
+    if (!cyclicBoundaryCells_.test(celli))
+    {
+        return patches;
+    }
+
+    const polyBoundaryMesh& bMesh = mesh_.boundaryMesh();
+
+    // So celli belongs to at least one cyclic patch.
+    // Let us figure out which.
+    if (globalNumbering_.isLocal(globalIdx)) // celli and globalIdx on same proc
+    {
+        // Making list of cyclic patches to which celli belongs
+        List<label> celliCyclicPatches;
+        forAll(bMesh, patchi)
+        {
+            if (isA<cyclicPolyPatch>(bMesh[patchi]))
+            {
+                // Note: Probably not efficient due to use of found(celli) but
+                // typically only used for very cells (interface cells and their
+                // point neighbours on cyclic boundaries).
+                if (bMesh[patchi].faceCells().found(celli))
+                {
+                    celliCyclicPatches.append(patchi);
+                }
+            }
+        }
+
+        // Get all local point neighbor cells of celli, i.e. all point
+        // neighbours that are not on the other side of a cyclic patch.
+        List<label> localPointNeiCells(0);
+        const labelList& cellPoints = mesh_.cellPoints()[celli];
+
+        for (const label cellPoint : cellPoints)
+        {
+            const labelList& pointKCells = mesh_.pointCells()[cellPoint];
+
+            for (const label pointKCell : pointKCells)
+            {
+                if (!localPointNeiCells.found(pointKCell))
+                {
+                    localPointNeiCells.append(pointKCell);
+                }
+            }
+        }
+
+        // Since globalIdx is a global cell index obtained from the point
+        // neighbour list, stencil[celli], all cells in this that are not in
+        // localPointNeiCells must be cyclic neighbour cells.
+        const label localIdx = globalNumbering_.toLocal(globalIdx);
+        if (!localPointNeiCells.found(localIdx))
+        {
+            for (const label patchi : celliCyclicPatches)
+            {
+                // find the corresponding cyclic neighbor patch ID
+                const cyclicPolyPatch& cpp =
+                    static_cast<const cyclicPolyPatch&>(bMesh[patchi]);
+
+                const label neiPatch = cpp.neighbPatchID();
+
+                // Check if the cell globalIdx is on neiPatch.
+                // If it is, append neiPatch to list of
+                if (bMesh[neiPatch].faceCells().found(localIdx))
+                {
+                    patches.append(neiPatch);
+                    // Here it may be possible to append patchi and do:
+                    //
+                    //    cpp.transformPosition()
+                    //
+                    // instead of
+                    //
+                    //    cpp.neighbPatch().transformPosition() in getPosition()
+                }
+            }
+        }
+    }
+    else // celli and globalIdx on differet processors
+    {
+        // Note: The following is needed if a celli is located at the interface
+        // (plicRDF), on a cyclic patch and a processor patch. In this case
+        // globalIdx may be on a different processor, but requires the
+        // transformation from a cyclic patch on the processor of celli.
+
+        const List<label>& faces = mesh_.cells()[celli];
+
+        // Loop over all faces of celli and find cyclic patches
+        for (const label facei : faces)
+        {
+            if (mesh_.isInternalFace(facei)) continue;
+
+            const label patchi = bMesh.whichPatch(facei);
+            const cyclicPolyPatch* cpp = isA<cyclicPolyPatch>(bMesh[patchi]);
+
+            if (cpp)
+            {
+                // Get the neighbor cell across the cyclic face
+                const label cycNeiPatch = cpp->neighbPatchID();
+                const label cycNeiCell =
+                    bMesh[cycNeiPatch].faceCells()[facei - cpp->start()];
+                const List<label>& cycNeiCellFaces = mesh_.cells()[cycNeiCell];
+
+                // Loop over all the faces of the neighbor cell
+                for (const label cycNeiCellFace : cycNeiCellFaces)
+                {
+                    if (mesh_.isInternalFace(cycNeiCellFace))
+                    {
+                        continue;
+                    }
+
+                    // Check if the neighbor cell has processor patches
+                    const label neiPatch = bMesh.whichPatch(cycNeiCellFace);
+                    const processorPolyPatch* ppp =
+                        isA<processorPolyPatch>(bMesh[neiPatch]);
+
+                    if (ppp)
+                    {
+                        // Avoid duplicate entries
+                        if (patches.found(cycNeiPatch))
+                        {
+                            continue;
+                        }
+
+                        // Since we can not access any information on globalIdx
+                        // we use the cell centre map from the stencil to
+                        // identify the cell.
+                        const label localFaceID = cycNeiCellFace - ppp->start();
+                        const vector neiCentre =
+                            ppp->neighbFaceCellCentres()[localFaceID];
+
+                        forAll(cyclicCentres_()[celli], k)
+                        {
+                            if
+                            (
+                                (
+                                    mag(cyclicCentres_()[celli][k] - neiCentre)
+                                        < 100*SMALL
+                                )
+                                && (stencil_[celli][k] == globalIdx)
+                            )
+                            {
+                                patches.append(cycNeiPatch);
+                                // Here an alternative might be to append patchi
+                                // and do:
+                                //     cpp.transformPosition()
+                                // instead of
+                                //     cpp.neighbPatch().transformPosition()
+                                // in getValue()
+                            }
+                        }
+                    }
+                }
+            }
+        }
+    }
+
+    return patches;
+}
+
+
+Foam::vector Foam::zoneDistribute::getPosition
+(
+    const VolumeField<vector>& positions,
+    const Map<vector>& valuesFromOtherProc,
+    const label gblIdx,
+    const List<label> cyclicPatchID
+) const
+{
+    // Position vector, possibly from other processor, to be returned
+    vector position(getValue(positions, valuesFromOtherProc, gblIdx));
+
+    // Dealing with position transformation across cyclic patches.
+    // If no transformation is required (most cases), cyclicPatchID is empty
+    forAll(cyclicPatchID, i)
+    {
+        const label patchi = cyclicPatchID[i];
+
+        const cyclicPolyPatch& cpp =
+            static_cast<const cyclicPolyPatch&>
+                (
+                    positions.mesh().boundaryMesh()[patchi]
+                );
+
+        if (cpp.transform() != coupledPolyPatch::transformType::ROTATIONAL)
+        {
+            cpp.neighbPatch().transformPosition(position, 0);
+        }
+        else if (globalNumbering_.isLocal(gblIdx))
+        {
+            const label localIdx = globalNumbering_.toLocal(gblIdx);
+
+            for (const label facei : mesh_.cells()[localIdx])
+            {
+                if (mesh_.boundaryMesh().whichPatch(facei) == cyclicPatchID[i])
+                {
+                    cpp.neighbPatch().transformPosition(position, facei);
+                    continue;
+                }
+            }
+        }
+        else
+        {
+            FatalErrorInFunction
+                << "Rotational cyclic patches are not supported in parallel.\n"
+                << "Try to decompose the domain so that the rotational cyclic"
+                << "patch is not split in between processors."
+                << exit(FatalError);
+        }
+    }
+
+    return position;
+}
+
 
 // ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/zoneDistribute/zoneDistribute.H b/src/finiteVolume/fvMesh/zoneDistribute/zoneDistribute.H
index dce29d9f86..575fcbd51e 100644
--- a/src/finiteVolume/fvMesh/zoneDistribute/zoneDistribute.H
+++ b/src/finiteVolume/fvMesh/zoneDistribute/zoneDistribute.H
@@ -115,6 +115,15 @@ class zoneDistribute
         //- Persistent set of exchange buffers
         PstreamBuffers pBufs_;
 
+        //- Cell labels for all cells with a face on a cyclic boundary
+        bitSet cyclicBoundaryCells_;
+
+        //- Holds for each cell on cyclic patch the centres of the cells in its
+        //  point neighbour stencil.
+        //  Used in getCyclicPatches to identify neigbhour patch ID of point
+        //  neighbours on other side of a processorPolyPatch.
+        autoPtr<Map<vectorField>> cyclicCentres_;
+
 
     // Private Member Functions
 
@@ -175,6 +184,19 @@ public:
             return globalNumbering_;
         }
 
+        //- Finds and returns list of all cyclic patch labels to which celli's
+        //  point neighbour cell, globalIdx, belongs. celli and globalIdx touch
+        //  in at least one point on these patces. globalIdx typically belongs
+        //  to stencil_[celli]. The returned label list is used to transform
+        //  positions across cyclic boundaries e.g. to be able to calculate
+        //  distances between cell centres and interface centres in plicRDF
+        //  across such boundaries.
+        List<label> getCyclicPatches
+        (
+            const label celli,
+            const label globalIdx
+        ) const;
+
         //- Gives patchNumber and patchFaceNumber for a given
         //- Geometric volume field
         template<typename Type>
@@ -185,6 +207,14 @@ public:
             const label gblIdx
         ) const;
 
+        vector getPosition
+        (
+            const VolumeField<vector>& positions,
+            const Map<vector>& valuesFromOtherProc,
+            const label gblIdx,
+            const List<label> cyclicPatchID = List<label>()
+        ) const;
+
         //- Returns stencil and provides a Map with globalNumbering
         template<typename Type>
         Map<Field<Type>> getFields
@@ -193,6 +223,17 @@ public:
             const VolumeField<Type>& phi
         );
 
+        //- Returns stencil and provides a Map with globalNumbering
+        //  Note: Not used currently (v2412) but needed for future surface
+        //  tension modelling.
+        template<typename Type>
+        Map<Field<Type>> getPositionFields
+        (
+            const boolList& zone,
+            const VolumeField<Type>& phi,
+            const bool& checkTransformation = false
+        );
+
         //- Returns stencil and provides a Map with globalNumbering
         template<typename Type>
         Map<Type> getDatafromOtherProc
diff --git a/src/finiteVolume/fvMesh/zoneDistribute/zoneDistributeI.H b/src/finiteVolume/fvMesh/zoneDistribute/zoneDistributeI.H
index 9be368d10f..e94beb1719 100644
--- a/src/finiteVolume/fvMesh/zoneDistribute/zoneDistributeI.H
+++ b/src/finiteVolume/fvMesh/zoneDistribute/zoneDistributeI.H
@@ -128,7 +128,59 @@ Foam::Map<Foam::Field<Type>> Foam::zoneDistribute::getFields
 
             for (const label gblIdx : stencil_[celli])
             {
-                tmpField.append(getValue(phi,neiValues,gblIdx));
+                tmpField.append(getValue(phi, neiValues, gblIdx));
+            }
+
+            stencilWithValues.emplace(celli, tmpField);
+        }
+    }
+
+    return stencilWithValues;
+}
+
+template<typename Type>
+Foam::Map<Foam::Field<Type>> Foam::zoneDistribute::getPositionFields
+(
+    const boolList& zone,
+    const VolumeField<Type>& phi,
+    const bool& checkTransformation
+)
+{
+    if (zone.size() != phi.size())
+    {
+        FatalErrorInFunction
+            << "size of zone: " << zone.size()
+            << "size of phi:" <<  phi.size()
+            << "do not match. Did the mesh change?"
+            << exit(FatalError);
+    }
+
+
+    // Get values from other proc
+    Map<Type> neiValues = getDatafromOtherProc(zone, phi);
+
+    Map<Field<Type>> stencilWithValues;
+
+    DynamicField<Type> tmpField(128);
+
+    forAll(zone, celli)
+    {
+        if (zone[celli])
+        {
+            tmpField.clear();
+
+            for (const label gblIdx : stencil_[celli])
+            {
+                List<label> cyclicPatches(0);
+                if(checkTransformation)
+                {
+                    cyclicPatches = getCyclicPatches(celli, gblIdx);
+                }
+
+                tmpField.append
+                (
+                    getPosition(phi, neiValues, gblIdx, cyclicPatches)
+                );
             }
 
             stencilWithValues.emplace(celli, tmpField);
diff --git a/src/transportModels/geometricVoF/reconstructedDistanceFunction/reconstructedDistanceFunction.C b/src/transportModels/geometricVoF/reconstructedDistanceFunction/reconstructedDistanceFunction.C
index da0c0a0b47..5b0cd84e63 100644
--- a/src/transportModels/geometricVoF/reconstructedDistanceFunction/reconstructedDistanceFunction.C
+++ b/src/transportModels/geometricVoF/reconstructedDistanceFunction/reconstructedDistanceFunction.C
@@ -282,9 +282,18 @@ const Foam::volScalarField&  Foam::reconstructedDistanceFunction::constructRDF
                     if (mag(n) != 0)
                     {
                         n /= mag(n);
-                        vector c = distribute.getValue(centre,mapCentres,gblIdx);
-                        vector distanceToIntSeg = (c - p);
-                        scalar distToSurf = distanceToIntSeg & (n);
+                        vector c
+                        (
+                            distribute.getPosition
+                            (
+                                centre,
+                                mapCentres,
+                                gblIdx,
+                                distribute.getCyclicPatches(celli, gblIdx)
+                            )
+                        );
+                        vector distanceToIntSeg(c - p);
+                        scalar distToSurf = distanceToIntSeg & n;
                         scalar weight = 0;
 
                         if (mag(distanceToIntSeg) != 0)
@@ -336,10 +345,18 @@ const Foam::volScalarField&  Foam::reconstructedDistanceFunction::constructRDF
                         if (mag(n) != 0)
                         {
                             n /= mag(n);
-                            vector c =
-                                distribute.getValue(centre, mapCentres, gblIdx);
-                            vector distanceToIntSeg = (c - p);
-                            scalar distToSurf = distanceToIntSeg & (n);
+                            vector c
+                            (
+                                distribute.getPosition
+                                (
+                                    centre,
+                                    mapCentres,
+                                    gblIdx,
+                                    distribute.getCyclicPatches(pCellI, gblIdx)
+                                )
+                            );
+                            vector distanceToIntSeg(c - p);
+                            scalar distToSurf = distanceToIntSeg & n;
                             scalar weight = 0;
 
                             if (mag(distanceToIntSeg) != 0)
diff --git a/src/transportModels/geometricVoF/reconstructionSchemes/plicSchemes/plicRDF/plicRDF.C b/src/transportModels/geometricVoF/reconstructionSchemes/plicSchemes/plicRDF/plicRDF.C
index 62db3416f8..8d591f2371 100644
--- a/src/transportModels/geometricVoF/reconstructionSchemes/plicSchemes/plicRDF/plicRDF.C
+++ b/src/transportModels/geometricVoF/reconstructionSchemes/plicSchemes/plicRDF/plicRDF.C
@@ -96,8 +96,16 @@ void Foam::reconstruction::plicRDF::interpolateNormal()
             if (mag(n) != 0)
             {
                 n /= mag(n);
-                vector centre =
-                    exchangeFields.getValue(centre_, mapCentre, gblIdx);
+                vector centre
+                (
+                    exchangeFields.getPosition
+                    (
+                        centre_,
+                        mapCentre,
+                        gblIdx,
+                        exchangeFields.getCyclicPatches(celli, gblIdx)
+                    )
+                );
                 vector distanceToIntSeg = (tensor::I- n*n) & (p - centre);
                 estimatedNormal += n /max(mag(distanceToIntSeg), SMALL);
                 weight += 1/max(mag(distanceToIntSeg), SMALL);
@@ -146,7 +154,11 @@ void Foam::reconstruction::plicRDF::interpolateNormal()
                 const label gblIdx = stencil[celli][i];
                 cellCentre.append
                 (
-                    exchangeFields.getValue(mesh_.C(), mapCC, gblIdx)
+                    exchangeFields.getPosition
+                    (
+                        mesh_.C(), mapCC, gblIdx,
+                        exchangeFields.getCyclicPatches(celli, gblIdx)
+                    )
                 );
                 alphaValues.append
                 (
@@ -193,7 +205,11 @@ void Foam::reconstruction::plicRDF::gradSurf(const volScalarField& phi)
         {
             cellCentre.append
             (
-                exchangeFields.getValue(mesh_.C(), mapCC, gblIdx)
+                exchangeFields.getPosition
+                (
+                    mesh_.C(), mapCC, gblIdx,
+                    exchangeFields.getCyclicPatches(celli, gblIdx)
+                )
             );
             phiValues.append
             (
diff --git a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/Allclean b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/Allclean
index fb1f384730..c13797ec72 100755
--- a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/Allclean
+++ b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/Allclean
@@ -5,4 +5,6 @@ cd "${0%/*}" || exit                                # Run from this directory
 
 cleanCase0
 
+rm -f Phase1
+
 #------------------------------------------------------------------------------
diff --git a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/Allrun b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/Allrun
index 321d8db8b1..8070263534 100755
--- a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/Allrun
+++ b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/Allrun
@@ -5,10 +5,20 @@ cd "${0%/*}" || exit                                # Run from this directory
 
 restore0Dir
 
+touch case.foam
+
 runApplication blockMesh
 
 runApplication setAlphaField
 
-runApplication $(getApplication)
+runApplication decomposePar
+
+runParallel $(getApplication)
+
+grep "Phase-1" log.$(getApplication) | cut -d' ' -f5 > Phase1
+
+echo "Volume conservation:"
+echo "Initial Phase-1: " $(head -1 Phase1)
+echo "Final Phase-1  : " $(tail -1 Phase1)
 
 #------------------------------------------------------------------------------
diff --git a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/controlDict b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/controlDict
index bfd00f8afb..beccd14af8 100644
--- a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/controlDict
+++ b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/controlDict
@@ -35,7 +35,7 @@ purgeWrite      0;
 
 writeFormat     ascii;
 
-writePrecision  6;
+writePrecision  12;
 
 writeCompression off;
 
diff --git a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/decomposeParDict b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/decomposeParDict
index 6c6f4b95f8..0caa1f4b9e 100644
--- a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/decomposeParDict
+++ b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/decomposeParDict
@@ -24,5 +24,13 @@ coeffs
     n   ( 2 1 2 );
 }
 
+constraints
+{
+    preservePatches
+    {
+        type            preservePatches;
+        patches         (left right top bottom);
+    }
+}
 
 // ************************************************************************* //
diff --git a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/fvSolution b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/fvSolution
index f242d11157..9a18a72b82 100644
--- a/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/fvSolution
+++ b/tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/fvSolution
@@ -22,10 +22,10 @@ solvers
     {
         isoFaceTol      1e-6;
         surfCellTol     1e-6;
-        snapTol         1e-10;
+        snapTol         0;
         nAlphaBounds    3;
-        clip            true;
-        reconstructionScheme isoAlpha; // isoAlpha
+        clip            false;
+        reconstructionScheme plicRDF;
 
         nAlphaSubCycles 1;
         cAlpha          1; // Read but not used by interIsoFoam