diff --git a/src/dynamicMesh/fvMeshDistribute/fvMeshDistribute.C b/src/dynamicMesh/fvMeshDistribute/fvMeshDistribute.C
index b28fa8a3a6a62c40f14e257165f354a7d90951ad..26d5582ac77c2a8c83940d4b1a495a80b82723d5 100644
--- a/src/dynamicMesh/fvMeshDistribute/fvMeshDistribute.C
+++ b/src/dynamicMesh/fvMeshDistribute/fvMeshDistribute.C
@@ -639,9 +639,14 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::repatch
     // shared points (see mergeSharedPoints below). So temporarily points
     // and edges do not match!
 
+    // TBD: temporarily unset mesh moving to avoid problems in meshflux
+    //      mapping. To be investigated.
+    const bool oldMoving = mesh_.moving(false);
     autoPtr<mapPolyMesh> mapPtr = meshMod.changeMesh(mesh_, false, true);
+    mesh_.moving(oldMoving);
     mapPolyMesh& map = *mapPtr;
 
+
     // Update fields. No inflation, parallel sync.
     mesh_.updateMesh(map);
 
@@ -766,7 +771,12 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::mergeSharedPoints
     fvMeshAdder::mergePoints(mesh_, pointToMaster, meshMod);
 
     // Change the mesh (no inflation). Note: parallel comms allowed.
+
+    // TBD: temporarily unset mesh moving to avoid problems in meshflux
+    //      mapping. To be investigated.
+    const bool oldMoving = mesh_.moving(false);
     autoPtr<mapPolyMesh> mapPtr = meshMod.changeMesh(mesh_, false, true);
+    mesh_.moving(oldMoving);
     mapPolyMesh& map = *mapPtr;
 
     // Update fields. No inflation, parallel sync.
@@ -1402,7 +1412,12 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::doRemoveCells
     saveInternalFields(tFlds);
 
     // Change the mesh. No inflation. Note: no parallel comms allowed.
+
+    // TBD: temporarily unset mesh moving to avoid problems in meshflux
+    //      mapping. To be investigated.
+    const bool oldMoving = mesh_.moving(false);
     autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh_, false, false);
+    mesh_.moving(oldMoving);
 
     // Update fields
     mesh_.updateMesh(map());
@@ -2825,6 +2840,10 @@ Foam::autoPtr<Foam::mapDistributePolyMesh> Foam::fvMeshDistribute::distribute
     const label nOldInternalFaces = mesh_.nInternalFaces();
     const labelList oldFaceOwner(mesh_.faceOwner());
 
+    // TBD: temporarily unset mesh moving to avoid problems in meshflux
+    //      mapping. To be investigated.
+    const bool oldMoving = mesh_.moving(false);
+
     fvMeshAdder::add
     (
         Pstream::myProcNo(),    // index of mesh to modify (== mesh_)
@@ -2842,6 +2861,9 @@ Foam::autoPtr<Foam::mapDistributePolyMesh> Foam::fvMeshDistribute::distribute
         constructPointMap
     );
 
+    mesh_.moving(oldMoving);
+
+
     if (debug)
     {
         Pout<< nl << "ADDED REMOTE MESHES:" << endl;
diff --git a/src/finiteVolume/fvMesh/fvGeometryScheme/fvGeometryScheme/fvGeometryScheme.C b/src/finiteVolume/fvMesh/fvGeometryScheme/fvGeometryScheme/fvGeometryScheme.C
index c703ad8aa20a947d48c038ab764e1eb7e2c9935e..88686584c0addce5a37494f6e303fab22639ca28 100644
--- a/src/finiteVolume/fvMesh/fvGeometryScheme/fvGeometryScheme/fvGeometryScheme.C
+++ b/src/finiteVolume/fvMesh/fvGeometryScheme/fvGeometryScheme/fvGeometryScheme.C
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020-2021 OpenCFD Ltd.
+    Copyright (C) 2020-2022 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -63,29 +63,33 @@ bool Foam::fvGeometryScheme::setMeshPhi() const
     const faceList& faces = mesh_.faces();
     const scalar rdt = 1.0/mesh_.time().deltaTValue();
 
-    auto& meshPhi = const_cast<fvMesh&>(mesh_).setPhi();
-    auto& meshPhii = meshPhi.primitiveFieldRef();
-    forAll(meshPhii, facei)
+    auto tmeshPhi(const_cast<fvMesh&>(mesh_).setPhi());
+    if (tmeshPhi)
     {
-        const face& f = faces[facei];
-        meshPhii[facei] = f.sweptVol(oldPoints, currPoints)*rdt;
-    }
-
-    auto& meshPhiBf = meshPhi.boundaryFieldRef();
-    for (auto& meshPhip : meshPhiBf)
-    {
-        if (!meshPhip.size())
+        auto& meshPhi = tmeshPhi.ref();
+        auto& meshPhii = meshPhi.primitiveFieldRef();
+        forAll(meshPhii, facei)
         {
-            // Empty patches
-            continue;
+            const face& f = faces[facei];
+            meshPhii[facei] = f.sweptVol(oldPoints, currPoints)*rdt;
         }
 
-        const auto& pp = meshPhip.patch().patch();
-
-        forAll(pp, facei)
+        auto& meshPhiBf = meshPhi.boundaryFieldRef();
+        for (auto& meshPhip : meshPhiBf)
         {
-            const face& f = pp[facei];
-            meshPhip[facei] = f.sweptVol(oldPoints, currPoints)*rdt;
+            if (!meshPhip.size())
+            {
+                // Empty patches
+                continue;
+            }
+
+            const auto& pp = meshPhip.patch().patch();
+
+            forAll(pp, facei)
+            {
+                const face& f = pp[facei];
+                meshPhip[facei] = f.sweptVol(oldPoints, currPoints)*rdt;
+            }
         }
     }
 
diff --git a/src/finiteVolume/fvMesh/fvGeometryScheme/solidBody/solidBodyFvGeometryScheme.C b/src/finiteVolume/fvMesh/fvGeometryScheme/solidBody/solidBodyFvGeometryScheme.C
index 7c48fb9fb6ffa694f41198e4bb1dc39d2abdbe24..74a8e66c63e9c34767b22413f5c0a379875f9cf0 100644
--- a/src/finiteVolume/fvMesh/fvGeometryScheme/solidBody/solidBodyFvGeometryScheme.C
+++ b/src/finiteVolume/fvMesh/fvGeometryScheme/solidBody/solidBodyFvGeometryScheme.C
@@ -220,40 +220,45 @@ void Foam::solidBodyFvGeometryScheme::movePoints()
 
         const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
         const faceList& faces = mesh_.faces();
-        const scalar rdt = 1.0/mesh_.time().deltaTValue();
 
-        // Set the mesh flux
-        auto& meshPhi = const_cast<fvMesh&>(mesh_).setPhi();
-        auto& meshPhii = meshPhi.primitiveFieldRef();
-        auto& meshPhiBf = meshPhi.boundaryFieldRef();
+        auto tmeshPhi(const_cast<fvMesh&>(mesh_).setPhi());
+        if (tmeshPhi)
+        {
+            const scalar rdt = 1.0/mesh_.time().deltaTValue();
 
-        //meshPhi == dimensionedScalar(dimVolume/dimTime, Zero);
-        meshPhii = Zero;
-        meshPhiBf == Zero;
+            // Set the mesh flux
+            auto& meshPhi = tmeshPhi.ref();
+            auto& meshPhii = meshPhi.primitiveFieldRef();
+            auto& meshPhiBf = meshPhi.boundaryFieldRef();
 
-        forAll(changedFaceIDs_, i)
-        {
-            const face& f = faces[changedFaceIDs_[i]];
+            //meshPhi == dimensionedScalar(dimVolume/dimTime, Zero);
+            meshPhii = Zero;
+            meshPhiBf == Zero;
 
-            if (changedPatchIDs_[i] == -1)
-            {
-                const label facei = changedFaceIDs_[i];
-                meshPhii[facei] = f.sweptVol(oldPoints, currPoints)*rdt;
-            }
-            else
+            forAll(changedFaceIDs_, i)
             {
-                const label patchi = changedPatchIDs_[i];
-                const polyPatch& pp = pbm[patchi];
+                const face& f = faces[changedFaceIDs_[i]];
 
-                if (isA<emptyPolyPatch>(pp))
+                if (changedPatchIDs_[i] == -1)
                 {
-                    continue;
+                    const label facei = changedFaceIDs_[i];
+                    meshPhii[facei] = f.sweptVol(oldPoints, currPoints)*rdt;
                 }
+                else
+                {
+                    const label patchi = changedPatchIDs_[i];
+                    const polyPatch& pp = pbm[patchi];
+
+                    if (isA<emptyPolyPatch>(pp))
+                    {
+                        continue;
+                    }
 
-                const label patchFacei = changedFaceIDs_[i] - pp.start();
+                    const label patchFacei = changedFaceIDs_[i] - pp.start();
 
-                meshPhiBf[patchi][patchFacei] =
-                    f.sweptVol(oldPoints, currPoints)*rdt;
+                    meshPhiBf[patchi][patchFacei] =
+                        f.sweptVol(oldPoints, currPoints)*rdt;
+                }
             }
         }
 
diff --git a/src/finiteVolume/fvMesh/fvMesh.H b/src/finiteVolume/fvMesh/fvMesh.H
index defd0f228efdde0ec78d56d480a8ace002de8466..58a4d3d8a8f93af276dfba0980163bd4639fe32e 100644
--- a/src/finiteVolume/fvMesh/fvMesh.H
+++ b/src/finiteVolume/fvMesh/fvMesh.H
@@ -482,8 +482,8 @@ public:
             //- these fvPatches.
             void removeFvBoundary();
 
-            //- Return cell face motion fluxes
-            surfaceScalarField& setPhi();
+            //- Return cell face motion fluxes (or null)
+            refPtr<surfaceScalarField> setPhi();
 
             //- Return old-time cell volumes
             DimensionedField<scalar, volMesh>& setV0();
diff --git a/src/finiteVolume/fvMesh/fvMeshGeometry.C b/src/finiteVolume/fvMesh/fvMeshGeometry.C
index bfe6d43193025c0ec8f7d7f7981db1fd6e99159a..853ce036e391ab351b19eb57dcbe196f35e66236 100644
--- a/src/finiteVolume/fvMesh/fvMeshGeometry.C
+++ b/src/finiteVolume/fvMesh/fvMeshGeometry.C
@@ -415,16 +415,19 @@ const Foam::surfaceScalarField& Foam::fvMesh::phi() const
 }
 
 
-Foam::surfaceScalarField& Foam::fvMesh::setPhi()
+Foam::refPtr<Foam::surfaceScalarField> Foam::fvMesh::setPhi()
 {
     if (!phiPtr_)
     {
-        FatalErrorInFunction
-            << "mesh flux field does not exist, is the mesh actually moving?"
-            << abort(FatalError);
+        return nullptr;
+    }
+    else
+    {
+        // Return non-const reference
+        refPtr<surfaceScalarField> p;
+        p.ref(*phiPtr_);
+        return p;
     }
-
-    return *phiPtr_;
 }
 
 
diff --git a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C
index 966cdb26a3fd2a78fbbd4bd449caa65552cb82f6..b573d3a2b2667a6610cfef6cf526f959ee3d3d2c 100644
--- a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C
+++ b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicACMI/cyclicACMIFvPatch.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2013-2016 OpenFOAM Foundation
-    Copyright (C) 2019 OpenCFD Ltd.
+    Copyright (C) 2019,2022 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -345,7 +345,7 @@ void Foam::cyclicACMIFvPatch::movePoints()
         const labelListList& newTgtAddr = AMI().tgtAddress();
 
         const fvMesh& mesh = boundaryMesh().mesh();
-        surfaceScalarField& meshPhi = const_cast<fvMesh&>(mesh).setPhi();
+        surfaceScalarField& meshPhi = const_cast<fvMesh&>(mesh).setPhi().ref();
         surfaceScalarField::Boundary& meshPhiBf = meshPhi.boundaryFieldRef();
 
         // Note: phip and phiNonOverlap will be different sizes if new faces
diff --git a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C
index 86d45f1637ca2efbc0508f412776d42a6f71ddf7..5b74ec455da57b200b3f3a620d6f7360d947b933 100644
--- a/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C
+++ b/src/finiteVolume/fvMesh/fvPatches/constraint/cyclicAMI/cyclicAMIFvPatch.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
-    Copyright (C) 2019-2020 OpenCFD Ltd.
+    Copyright (C) 2019-2020,2022 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -267,7 +267,7 @@ void Foam::cyclicAMIFvPatch::movePoints()
     // src = src + mapped(tgt) and tgt = tgt + mapped(src)?
 
     const fvMesh& mesh = boundaryMesh().mesh();
-    surfaceScalarField& meshPhi = const_cast<fvMesh&>(mesh).setPhi();
+    surfaceScalarField& meshPhi = const_cast<fvMesh&>(mesh).setPhi().ref();
     surfaceScalarField::Boundary& meshPhiBf = meshPhi.boundaryFieldRef();
 
     if (cyclicAMIPolyPatch_.owner())