diff --git a/applications/solvers/multiphase/reactingMultiphaseEulerFoam/pUf/pEqn.H b/applications/solvers/multiphase/reactingMultiphaseEulerFoam/pUf/pEqn.H
index 08d97031b044fd0379c63b333fde4638e4cdc742..b66aabffce6e986fb1db584af0f7772c02eea69b 100644
--- a/applications/solvers/multiphase/reactingMultiphaseEulerFoam/pUf/pEqn.H
+++ b/applications/solvers/multiphase/reactingMultiphaseEulerFoam/pUf/pEqn.H
@@ -272,10 +272,7 @@ while (pimple.correct())
                     ).ptr()
                 );
 
-                deleteDemandDrivenData
-                (
-                    pEqnComps[phasei].faceFluxCorrectionPtr()
-                );
+                pEqnComps[phasei].faceFluxCorrectionPtr(nullptr);
 
                 pEqnComps[phasei].relax();
             }
diff --git a/src/finiteArea/faMatrices/faMatrix/faMatrix.C b/src/finiteArea/faMatrices/faMatrix/faMatrix.C
index 9db099316f26eca8238d8143cccc9e8caec46eb1..ac90e7101280522ca012ac3d104faa5784bc8f17 100644
--- a/src/finiteArea/faMatrices/faMatrix/faMatrix.C
+++ b/src/finiteArea/faMatrices/faMatrix/faMatrix.C
@@ -190,8 +190,7 @@ Foam::faMatrix<Type>::faMatrix
     dimensions_(ds),
     source_(psi.size(), Zero),
     internalCoeffs_(psi.mesh().boundary().size()),
-    boundaryCoeffs_(psi.mesh().boundary().size()),
-    faceFluxCorrectionPtr_(nullptr)
+    boundaryCoeffs_(psi.mesh().boundary().size())
 {
     DebugInFunction
         << "constructing faMatrix<Type> for field " << psi_.name()
@@ -231,19 +230,17 @@ Foam::faMatrix<Type>::faMatrix(const faMatrix<Type>& fam)
     dimensions_(fam.dimensions_),
     source_(fam.source_),
     internalCoeffs_(fam.internalCoeffs_),
-    boundaryCoeffs_(fam.boundaryCoeffs_),
-    faceFluxCorrectionPtr_(nullptr)
+    boundaryCoeffs_(fam.boundaryCoeffs_)
 {
     DebugInFunction
         << "Copying faMatrix<Type> for field " << psi_.name() << endl;
 
     if (fam.faceFluxCorrectionPtr_)
     {
-        faceFluxCorrectionPtr_ =
-            new GeometricField<Type, faePatchField, edgeMesh>
-            (
-                *(fam.faceFluxCorrectionPtr_)
-            );
+        faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
+        (
+            *(fam.faceFluxCorrectionPtr_)
+        );
     }
 }
 
@@ -256,8 +253,7 @@ Foam::faMatrix<Type>::faMatrix(const tmp<faMatrix<Type>>& tmat)
     dimensions_(tmat().dimensions_),
     source_(tmat.constCast().source_, tmat.movable()),
     internalCoeffs_(tmat.constCast().internalCoeffs_, tmat.movable()),
-    boundaryCoeffs_(tmat.constCast().boundaryCoeffs_, tmat.movable()),
-    faceFluxCorrectionPtr_(nullptr)
+    boundaryCoeffs_(tmat.constCast().boundaryCoeffs_, tmat.movable())
 {
     DebugInFunction
         << "Copy/Move faMatrix<Type> for field " << psi_.name() << endl;
@@ -266,16 +262,15 @@ Foam::faMatrix<Type>::faMatrix(const tmp<faMatrix<Type>>& tmat)
     {
         if (tmat.movable())
         {
-            faceFluxCorrectionPtr_ = tmat().faceFluxCorrectionPtr_;
-            tmat().faceFluxCorrectionPtr_ = nullptr;
+            faceFluxCorrectionPtr_ =
+                std::move(tmat.constCast().faceFluxCorrectionPtr_);
         }
-        else
+        else if (tmat().faceFluxCorrectionPtr_)
         {
-            faceFluxCorrectionPtr_ =
-                new GeometricField<Type, faePatchField, edgeMesh>
-                (
-                    *(tmat().faceFluxCorrectionPtr_)
-                );
+            faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
+            (
+                *(tmat().faceFluxCorrectionPtr_)
+            );
         }
     }
 
@@ -290,8 +285,6 @@ Foam::faMatrix<Type>::~faMatrix()
 {
     DebugInFunction
         << "Destroying faMatrix<Type> for field " << psi_.name() << endl;
-
-    deleteDemandDrivenData(faceFluxCorrectionPtr_);
 }
 
 
@@ -788,9 +781,10 @@ void Foam::faMatrix<Type>::operator=(const faMatrix<Type>& famv)
     }
     else if (famv.faceFluxCorrectionPtr_)
     {
-        faceFluxCorrectionPtr_ =
-            new GeometricField<Type, faePatchField, edgeMesh>
-            (*famv.faceFluxCorrectionPtr_);
+        faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
+        (
+            *famv.faceFluxCorrectionPtr_
+        );
     }
 }
 
@@ -835,8 +829,7 @@ void Foam::faMatrix<Type>::operator+=(const faMatrix<Type>& famv)
     }
     else if (famv.faceFluxCorrectionPtr_)
     {
-        faceFluxCorrectionPtr_ = new
-        GeometricField<Type, faePatchField, edgeMesh>
+        faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
         (
             *famv.faceFluxCorrectionPtr_
         );
@@ -869,9 +862,10 @@ void Foam::faMatrix<Type>::operator-=(const faMatrix<Type>& famv)
     }
     else if (famv.faceFluxCorrectionPtr_)
     {
-        faceFluxCorrectionPtr_ =
-            new GeometricField<Type, faePatchField, edgeMesh>
-            (-*famv.faceFluxCorrectionPtr_);
+        faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
+        (
+            -*famv.faceFluxCorrectionPtr_
+        );
     }
 }
 
diff --git a/src/finiteArea/faMatrices/faMatrix/faMatrix.H b/src/finiteArea/faMatrices/faMatrix/faMatrix.H
index dd56ad68dbf4bf4085530827a00b50e3ac018040..0b4681537c04c95e2a455c36f334c360bb16a2bd 100644
--- a/src/finiteArea/faMatrices/faMatrix/faMatrix.H
+++ b/src/finiteArea/faMatrices/faMatrix/faMatrix.H
@@ -151,7 +151,7 @@ private:
         FieldField<Field, Type> boundaryCoeffs_;
 
         //- Face flux field for non-orthogonal correction
-        mutable faceFluxFieldType* faceFluxCorrectionPtr_;
+        mutable std::unique_ptr<faceFluxFieldType> faceFluxCorrectionPtr_;
 
 
 protected:
@@ -353,8 +353,7 @@ public:
             }
 
             //- Declare return type of the faceFluxCorrectionPtr() function
-            typedef GeometricField<Type, faePatchField, edgeMesh>
-                *faceFluxFieldPtrType;
+            typedef std::unique_ptr<faceFluxFieldType> faceFluxFieldPtrType;
 
             //- Return pointer to face-flux non-orthogonal correction field
             faceFluxFieldPtrType& faceFluxCorrectionPtr()
@@ -362,6 +361,12 @@ public:
                 return faceFluxCorrectionPtr_;
             }
 
+            //- Set pointer to face-flux non-orthogonal correction field
+            void faceFluxCorrectionPtr(faceFluxFieldType* flux)
+            {
+                faceFluxCorrectionPtr_.reset(flux);
+            }
+
             //- True if face-flux non-orthogonal correction field exists
             bool hasFaceFluxCorrection() const noexcept
             {
diff --git a/src/finiteArea/finiteArea/laplacianSchemes/gaussFaLaplacianScheme/gaussFaLaplacianScheme.C b/src/finiteArea/finiteArea/laplacianSchemes/gaussFaLaplacianScheme/gaussFaLaplacianScheme.C
index 46e7d4bb4f3502bb177523dc985d30905cc2b5d1..b5236b2369e9787d3ef4f1d3bcc829751914c503 100644
--- a/src/finiteArea/finiteArea/laplacianSchemes/gaussFaLaplacianScheme/gaussFaLaplacianScheme.C
+++ b/src/finiteArea/finiteArea/laplacianSchemes/gaussFaLaplacianScheme/gaussFaLaplacianScheme.C
@@ -82,8 +82,10 @@ gaussLaplacianScheme<Type>::famLaplacian
     {
         if (this->mesh().fluxRequired(vf.name()))
         {
-            fam.faceFluxCorrectionPtr() = new
-            GeometricField<Type, faePatchField, edgeMesh>
+            fam.faceFluxCorrectionPtr() = std::make_unique
+            <
+                GeometricField<Type, faePatchField, edgeMesh>
+            >
             (
                 gammaMagSf*this->tlnGradScheme_().correction(vf)
             );
diff --git a/src/finiteVolume/finiteVolume/laplacianSchemes/gaussLaplacianScheme/gaussLaplacianScheme.C b/src/finiteVolume/finiteVolume/laplacianSchemes/gaussLaplacianScheme/gaussLaplacianScheme.C
index 20f20e80e1b4df8b11b79b96804819de65f402b7..27f420f72eeba9c31b0b1c91e72ca8afb0c16ca0 100644
--- a/src/finiteVolume/finiteVolume/laplacianSchemes/gaussLaplacianScheme/gaussLaplacianScheme.C
+++ b/src/finiteVolume/finiteVolume/laplacianSchemes/gaussLaplacianScheme/gaussLaplacianScheme.C
@@ -194,7 +194,7 @@ gaussLaplacianScheme<Type, GType>::fvmLaplacian
 
     if (mesh.fluxRequired(vf.name()))
     {
-        fvm.faceFluxCorrectionPtr() = tfaceFluxCorrection.ptr();
+        fvm.faceFluxCorrectionPtr(tfaceFluxCorrection.ptr());
     }
 
     return tfvm;
diff --git a/src/finiteVolume/finiteVolume/laplacianSchemes/gaussLaplacianScheme/gaussLaplacianSchemes.C b/src/finiteVolume/finiteVolume/laplacianSchemes/gaussLaplacianScheme/gaussLaplacianSchemes.C
index 1201bad1fe2f66ae24587653e2da1b8d49d94a6b..7d19f5b2656c2fb8c48e39e67aa6cb57466516f6 100644
--- a/src/finiteVolume/finiteVolume/laplacianSchemes/gaussLaplacianScheme/gaussLaplacianSchemes.C
+++ b/src/finiteVolume/finiteVolume/laplacianSchemes/gaussLaplacianScheme/gaussLaplacianSchemes.C
@@ -61,8 +61,10 @@ Foam::fv::gaussLaplacianScheme<Foam::Type, Foam::scalar>::fvmLaplacian         \
     {                                                                          \
         if (mesh.fluxRequired(vf.name()))                                      \
         {                                                                      \
-            fvm.faceFluxCorrectionPtr() = new                                  \
-            GeometricField<Type, fvsPatchField, surfaceMesh>                   \
+            fvm.faceFluxCorrectionPtr() = std::make_unique                     \
+            <                                                                  \
+                GeometricField<Type, fvsPatchField, surfaceMesh>               \
+            >                                                                  \
             (                                                                  \
                 gammaMagSf*this->tsnGradScheme_().correction(vf)               \
             );                                                                 \
diff --git a/src/finiteVolume/finiteVolume/laplacianSchemes/relaxedNonOrthoGaussLaplacianScheme/relaxedNonOrthoGaussLaplacianScheme.C b/src/finiteVolume/finiteVolume/laplacianSchemes/relaxedNonOrthoGaussLaplacianScheme/relaxedNonOrthoGaussLaplacianScheme.C
index ef0efdb299fe2026b77b40934784c78cc711ef12..44e554ad400a705d89ec7e3409296bbc805de612 100644
--- a/src/finiteVolume/finiteVolume/laplacianSchemes/relaxedNonOrthoGaussLaplacianScheme/relaxedNonOrthoGaussLaplacianScheme.C
+++ b/src/finiteVolume/finiteVolume/laplacianSchemes/relaxedNonOrthoGaussLaplacianScheme/relaxedNonOrthoGaussLaplacianScheme.C
@@ -222,7 +222,7 @@ relaxedNonOrthoGaussLaplacianScheme<Type, GType>::fvmLaplacian
 
     if (mesh.fluxRequired(vf.name()))
     {
-        fvm.faceFluxCorrectionPtr() = trelaxedCorrection.ptr();
+        fvm.faceFluxCorrectionPtr(trelaxedCorrection.ptr());
     }
 
     return tfvm;
diff --git a/src/finiteVolume/finiteVolume/laplacianSchemes/relaxedNonOrthoGaussLaplacianScheme/relaxedNonOrthoGaussLaplacianSchemes.C b/src/finiteVolume/finiteVolume/laplacianSchemes/relaxedNonOrthoGaussLaplacianScheme/relaxedNonOrthoGaussLaplacianSchemes.C
index ce32fd20fbbb7f9b4ed3f2494365640e9b1db887..662e909dd5095c1b86164aa9966a11b07e7b8c31 100644
--- a/src/finiteVolume/finiteVolume/laplacianSchemes/relaxedNonOrthoGaussLaplacianScheme/relaxedNonOrthoGaussLaplacianSchemes.C
+++ b/src/finiteVolume/finiteVolume/laplacianSchemes/relaxedNonOrthoGaussLaplacianScheme/relaxedNonOrthoGaussLaplacianSchemes.C
@@ -98,7 +98,7 @@ fvmLaplacian                                                                   \
                                                                                \
         if (mesh.fluxRequired(vf.name()))                                      \
         {                                                                      \
-            fvm.faceFluxCorrectionPtr() = trelaxedCorrection.ptr();            \
+            fvm.faceFluxCorrectionPtr(trelaxedCorrection.ptr());               \
         }                                                                      \
     }                                                                          \
                                                                                \
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
index bf2ab42635a6bee80179197d9f820be4fd526238..254436784f3e1c4e492f68a249b335a9b81dc0f1 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
@@ -369,8 +369,7 @@ Foam::fvMatrix<Type>::fvMatrix
     dimensions_(ds),
     source_(psi.size(), Zero),
     internalCoeffs_(psi.mesh().boundary().size()),
-    boundaryCoeffs_(psi.mesh().boundary().size()),
-    faceFluxCorrectionPtr_(nullptr)
+    boundaryCoeffs_(psi.mesh().boundary().size())
 {
     DebugInFunction
         << "Constructing fvMatrix<Type> for field " << psi_.name() << endl;
@@ -410,19 +409,17 @@ Foam::fvMatrix<Type>::fvMatrix(const fvMatrix<Type>& fvm)
     dimensions_(fvm.dimensions_),
     source_(fvm.source_),
     internalCoeffs_(fvm.internalCoeffs_),
-    boundaryCoeffs_(fvm.boundaryCoeffs_),
-    faceFluxCorrectionPtr_(nullptr)
+    boundaryCoeffs_(fvm.boundaryCoeffs_)
 {
     DebugInFunction
         << "Copying fvMatrix<Type> for field " << psi_.name() << endl;
 
     if (fvm.faceFluxCorrectionPtr_)
     {
-        faceFluxCorrectionPtr_ =
-            new GeometricField<Type, fvsPatchField, surfaceMesh>
-            (
-                *(fvm.faceFluxCorrectionPtr_)
-            );
+        faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
+        (
+            *(fvm.faceFluxCorrectionPtr_)
+        );
     }
 }
 
@@ -438,8 +435,7 @@ Foam::fvMatrix<Type>::fvMatrix(const tmp<fvMatrix<Type>>& tmat)
     dimensions_(tmat().dimensions_),
     source_(tmat.constCast().source_, tmat.movable()),
     internalCoeffs_(tmat.constCast().internalCoeffs_, tmat.movable()),
-    boundaryCoeffs_(tmat.constCast().boundaryCoeffs_, tmat.movable()),
-    faceFluxCorrectionPtr_(nullptr)
+    boundaryCoeffs_(tmat.constCast().boundaryCoeffs_, tmat.movable())
 {
     DebugInFunction
         << "Copy/move fvMatrix<Type> for field " << psi_.name() << endl;
@@ -448,16 +444,15 @@ Foam::fvMatrix<Type>::fvMatrix(const tmp<fvMatrix<Type>>& tmat)
     {
         if (tmat.movable())
         {
-            faceFluxCorrectionPtr_ = tmat().faceFluxCorrectionPtr_;
-            tmat().faceFluxCorrectionPtr_ = nullptr;
+            faceFluxCorrectionPtr_ =
+                std::move(tmat.constCast().faceFluxCorrectionPtr_);
         }
-        else
+        else if (tmat().faceFluxCorrectionPtr_)
         {
-            faceFluxCorrectionPtr_ =
-                new GeometricField<Type, fvsPatchField, surfaceMesh>
-                (
-                    *(tmat().faceFluxCorrectionPtr_)
-                );
+            faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
+            (
+                *(tmat().faceFluxCorrectionPtr_)
+            );
         }
     }
 
@@ -473,7 +468,6 @@ Foam::fvMatrix<Type>::~fvMatrix()
     DebugInFunction
         << "Destroying fvMatrix<Type> for field " << psi_.name() << endl;
 
-    deleteDemandDrivenData(faceFluxCorrectionPtr_);
     subMatrices_.clear();
 }
 
@@ -1577,9 +1571,10 @@ void Foam::fvMatrix<Type>::operator=(const fvMatrix<Type>& fvmv)
     }
     else if (fvmv.faceFluxCorrectionPtr_)
     {
-        faceFluxCorrectionPtr_ =
-            new GeometricField<Type, fvsPatchField, surfaceMesh>
-        (*fvmv.faceFluxCorrectionPtr_);
+        faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
+        (
+            *fvmv.faceFluxCorrectionPtr_
+        );
     }
 
     useImplicit_ = fvmv.useImplicit_;
@@ -1631,8 +1626,7 @@ void Foam::fvMatrix<Type>::operator+=(const fvMatrix<Type>& fvmv)
     }
     else if (fvmv.faceFluxCorrectionPtr_)
     {
-        faceFluxCorrectionPtr_ = new
-        GeometricField<Type, fvsPatchField, surfaceMesh>
+        faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
         (
             *fvmv.faceFluxCorrectionPtr_
         );
@@ -1669,9 +1663,10 @@ void Foam::fvMatrix<Type>::operator-=(const fvMatrix<Type>& fvmv)
     }
     else if (fvmv.faceFluxCorrectionPtr_)
     {
-        faceFluxCorrectionPtr_ =
-            new GeometricField<Type, fvsPatchField, surfaceMesh>
-        (-*fvmv.faceFluxCorrectionPtr_);
+        faceFluxCorrectionPtr_ = std::make_unique<faceFluxFieldType>
+        (
+            -*fvmv.faceFluxCorrectionPtr_
+        );
     }
 }
 
@@ -2002,7 +1997,7 @@ Foam::tmp<Foam::fvMatrix<Type>> Foam::correction
 
     // Delete the faceFluxCorrection from the correction matrix
     // as it does not have a clear meaning or purpose
-    deleteDemandDrivenData(tAcorr.ref().faceFluxCorrectionPtr());
+    tAcorr.ref().faceFluxCorrectionPtr(nullptr);
 
     return tAcorr;
 }
@@ -2018,7 +2013,7 @@ Foam::tmp<Foam::fvMatrix<Type>> Foam::correction
 
     // Delete the faceFluxCorrection from the correction matrix
     // as it does not have a clear meaning or purpose
-    deleteDemandDrivenData(tAcorr.ref().faceFluxCorrectionPtr());
+    tAcorr.ref().faceFluxCorrectionPtr(nullptr);
 
     return tAcorr;
 }
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H
index c00236f7c9f8a5da6539e4fefb304872df54a5f9..78f18d22c3e97b8dc76dc7efc3928b50320a7f70 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H
@@ -172,7 +172,7 @@ private:
         FieldField<Field, Type> boundaryCoeffs_;
 
         //- Face flux field for non-orthogonal correction
-        mutable faceFluxFieldType* faceFluxCorrectionPtr_;
+        mutable std::unique_ptr<faceFluxFieldType> faceFluxCorrectionPtr_;
 
 
 protected:
@@ -494,8 +494,7 @@ public:
             }
 
             //- Declare return type of the faceFluxCorrectionPtr() function
-            typedef GeometricField<Type, fvsPatchField, surfaceMesh>
-                *faceFluxFieldPtrType;
+            typedef std::unique_ptr<faceFluxFieldType> faceFluxFieldPtrType;
 
             //- Return pointer to face-flux non-orthogonal correction field
             faceFluxFieldPtrType& faceFluxCorrectionPtr()
@@ -503,6 +502,12 @@ public:
                 return faceFluxCorrectionPtr_;
             }
 
+            //- Set pointer to face-flux non-orthogonal correction field
+            void faceFluxCorrectionPtr(faceFluxFieldType* flux)
+            {
+                faceFluxCorrectionPtr_.reset(flux);
+            }
+
             //- True if face-flux non-orthogonal correction field exists
             bool hasFaceFluxCorrection() const noexcept
             {