fvMesh: Generalised the handing of old-time fields during mesh motion and topology change

Topology change occurs before the time-increment and hence the oldest time
field (old-time in the case of 1st order time schemes, old-old-time in the case
of 2nd-order time schemes) is not actually needed as it is replaced by the
current time-field after time-increment so there is no purpose to mapping this
field.  However, it is necessary to keep track of the existence of the
oldest-time field to ensure the correct number of old-time fields are cached for
the time-scheme.  This development allows fvMesh to delete the redundant
oldest-time fields in such a manner that GeometricField can reinstate them
correctly after time-increment which is more efficient and more reliable than
attempting to map them and done previously.

Additionally fvMesh movement, which occurs after time-increment, now ensure all
old-time fields are up-to-date before NCC stitcher mapping so that both fields
and their old-time values are mapped consistently.  This removes the need for
old-time field caching calls in MapGeometricFields, fvMeshAdder and
fvMeshStitcher, thus simplifying the code and improving maintainability.

applications/test/nullObject/Test-nullObject.C

@@ -40,10 +40,34 @@ int main()
     }
 
 
+    // Test const pointer and const reference to the nullObject
+
+    const SimpleClass* constPtrToNull(NullObjectConstPtr<SimpleClass>());
+    const SimpleClass& constRefToNull(NullObjectRef<SimpleClass>());
+
+    if (isNull(constPtrToNull))
+    {
+        Info<< "Pass: constPtrToNull is null" << endl;
+    }
+    else
+    {
+        Info<< "FAIL: constPtrToNull is not null" << endl;
+    }
+
+    if (isNull(constRefToNull))
+    {
+        Info<< "Pass: constRefToNull is null" << endl;
+    }
+    else
+    {
+        Info<< "FAIL: constRefToNull is not null" << endl;
+    }
+
+
     // Test pointer and reference to the nullObject
 
-    const SimpleClass* ptrToNull(NullObjectPtr<SimpleClass>());
+    SimpleClass* ptrToNull(NullObjectPtr<SimpleClass>());
-    const SimpleClass& refToNull(*ptrToNull);
+    SimpleClass& refToNull(*ptrToNull);
 
     if (isNull(ptrToNull))
     {

src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.C

@@ -1068,6 +1068,8 @@ template<class Type, template<class> class PatchField, class GeoMesh>
 void Foam::GeometricField<Type, PatchField, GeoMesh>::storeOldTime() const
 {
     if (field0Ptr_)
+    {
+        if (notNull(field0Ptr_))
         {
             field0Ptr_->storeOldTime();
 
@@ -1086,6 +1088,12 @@ void Foam::GeometricField<Type, PatchField, GeoMesh>::storeOldTime() const
                 field0Ptr_->writeOpt() = this->writeOpt();
             }
         }
+        else
+        {
+            // Reinstate old-time field
+            oldTime();
+        }
+    }
 }
 
 
@@ -1093,9 +1101,16 @@ template<class Type, template<class> class PatchField, class GeoMesh>
 Foam::label Foam::GeometricField<Type, PatchField, GeoMesh>::nOldTimes() const
 {
     if (field0Ptr_)
+    {
+        if (isNull(field0Ptr_))
+        {
+            return 1;
+        }
+        else
         {
             return field0Ptr_->nOldTimes() + 1;
         }
+    }
     else
     {
         return 0;
@@ -1107,8 +1122,12 @@ template<class Type, template<class> class PatchField, class GeoMesh>
 const Foam::GeometricField<Type, PatchField, GeoMesh>&
 Foam::GeometricField<Type, PatchField, GeoMesh>::oldTime() const
 {
-    if (!field0Ptr_)
+    if (!field0Ptr_ || isNull(field0Ptr_))
     {
+        // Set field0Ptr_ null to ensure the old-time field constructor
+        // does not construct the old-old-time field
+        field0Ptr_ = nullptr;
+
         field0Ptr_ = new GeometricField<Type, PatchField, GeoMesh>
         (
             IOobject
@@ -1218,11 +1237,46 @@ Foam::GeometricField<Type, PatchField, GeoMesh>::prevIter() const
 template<class Type, template<class> class PatchField, class GeoMesh>
 void Foam::GeometricField<Type, PatchField, GeoMesh>::clearOldTimes()
 {
+    if (field0Ptr_ && notNull(field0Ptr_))
+    {
         deleteDemandDrivenData(field0Ptr_);
+    }
+
     deleteDemandDrivenData(fieldPrevIterPtr_);
 }
 
 
+template<class Type, template<class> class PatchField, class GeoMesh>
+void Foam::GeometricField<Type, PatchField, GeoMesh>::nullOldestTime()
+{
+    // Check that the field is not an old-time field
+    if (!isOldTime())
+    {
+        // Search for the oldest old-time field and set to nullObject
+        nullOldestTimeFound();
+    }
+}
+
+
+template<class Type, template<class> class PatchField, class GeoMesh>
+void Foam::GeometricField<Type, PatchField, GeoMesh>::nullOldestTimeFound()
+{
+    if (field0Ptr_ && notNull(field0Ptr_))
+    {
+        if (field0Ptr_->field0Ptr_)
+        {
+            field0Ptr_->nullOldestTimeFound();
+        }
+        else
+        {
+            deleteDemandDrivenData(field0Ptr_);
+            field0Ptr_ =
+                NullObjectPtr<GeometricField<Type, PatchField, GeoMesh>>();
+        }
+    }
+}
+
+
 template<class Type, template<class> class PatchField, class GeoMesh>
 void Foam::GeometricField<Type, PatchField, GeoMesh>::
 correctBoundaryConditions()

src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.H

@@ -130,6 +130,12 @@ private:
         //- Read the field - create the field dictionary on-the-fly
         void readFields();
 
+        //- Store the old-time field
+        void storeOldTime() const;
+
+        //- Find oldest time and set to nullObject
+        void nullOldestTimeFound();
+
 
 public:
 
@@ -407,9 +413,6 @@ public:
         //- Store the old-time fields
         void storeOldTimes() const;
 
-        //- Store the old-time field
-        void storeOldTime() const;
-
         //- Return the number of old time fields stored
         label nOldTimes() const;
 
@@ -439,6 +442,14 @@ public:
         //- Delete old time and previous iteration fields
         void clearOldTimes();
 
+        //- Set oldest time field pointer to nullObjectPtr
+        //  to remove the oldest time field while maintaining a tag
+        //  so that it is reinstated on the next call to oldTime().
+        //
+        //  This is used to temporarily remove the oldest-time field
+        //  prior to mesh topology change.
+        void nullOldestTime();
+
         //- Correct boundary field
         void correctBoundaryConditions();
 

src/OpenFOAM/fields/GeometricFields/GeometricField/MapGeometricFields.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-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -78,31 +78,6 @@ void MapGeometricFields
             lookupClass<GeometricField<Type, PatchField, GeoMesh>>()
     );
 
-    // It is necessary to enforce that all old-time fields are stored
-    // before the mapping is performed.  Otherwise, if the
-    // old-time-level field is mapped before the field itself, sizes
-    // will not match.
-
-    for
-    (
-        typename HashTable<const GeometricField<Type, PatchField, GeoMesh>*>::
-            iterator fieldIter = fields.begin();
-        fieldIter != fields.end();
-        ++fieldIter
-    )
-    {
-        GeometricField<Type, PatchField, GeoMesh>& field =
-            const_cast<GeometricField<Type, PatchField, GeoMesh>&>
-            (*fieldIter());
-
-        // Note: check can be removed once pointFields are actually stored on
-        //      the pointMesh instead of now on the polyMesh!
-        if (&field.mesh() == &mapper.mesh())
-        {
-            field.storeOldTimes();
-        }
-    }
-
     for
     (
         typename HashTable<const GeometricField<Type, PatchField, GeoMesh>*>::

src/OpenFOAM/primitives/nullObject/nullObject.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2014-2019 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2014-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -68,13 +68,17 @@ public:
 extern const NullObject* nullObjectPtr;
 
 
-//- Return reference to the nullObject of type T
+//- Return const reference to the nullObject of type T
 template<class T>
 inline const T& NullObjectRef();
 
+//- Return const pointer to the nullObject of type T
+template<class T>
+inline const T* NullObjectConstPtr();
+
 //- Return pointer to the nullObject of type T
 template<class T>
-inline const T* NullObjectPtr();
+inline T* NullObjectPtr();
 
 //- Return reference to the nullObject of type T
 template<class T>
@@ -99,6 +103,15 @@ template<class T>
 inline bool notNull(const T* t);
 
 
+//- Return true if t is a pointer to the nullObject of type T
+template<class T>
+inline bool isNull(T* t);
+
+//- Return true if t is not a pointer to the nullObject of type T
+template<class T>
+inline bool notNull(T* t);
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam

src/OpenFOAM/primitives/nullObject/nullObjectI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2014-2019 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2014-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -36,11 +36,17 @@ inline T&& Foam::NullObjectMove()
 }
 
 template<class T>
-inline const T* Foam::NullObjectPtr()
+inline const T* Foam::NullObjectConstPtr()
 {
     return reinterpret_cast<const T*>(nullObjectPtr);
 }
 
+template<class T>
+inline T* Foam::NullObjectPtr()
+{
+    return const_cast<T*>(reinterpret_cast<const T*>(nullObjectPtr));
+}
+
 
 template<class T>
 inline bool Foam::isNull(const T& t)
@@ -68,4 +74,17 @@ inline bool Foam::notNull(const T* t)
 }
 
 
+template<class T>
+inline bool Foam::isNull(T* t)
+{
+    return t == NullObjectPtr<T>();
+}
+
+template<class T>
+inline bool Foam::notNull(T* t)
+{
+    return t != NullObjectPtr<T>();
+}
+
+
 // ************************************************************************* //

src/dynamicMesh/fvMeshAdder/fvMeshAdderTemplates.C

@@ -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-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -270,30 +270,6 @@ void Foam::fvMeshAdder::MapVolFields
         ()
     );
 
-    // It is necessary to enforce that all old-time fields are stored
-    // before the mapping is performed.  Otherwise, if the
-    // old-time-level field is mapped before the field itself, sizes
-    // will not match.
-
-    for
-    (
-        typename HashTable<const GeometricField<Type, fvPatchField, volMesh>*>::
-            iterator fieldIter = fields.begin();
-        fieldIter != fields.end();
-        ++fieldIter
-    )
-    {
-        if (debug)
-        {
-            Pout<< "MapVolFields : Storing old time for " << fieldIter()->name()
-                << endl;
-        }
-
-        const_cast<GeometricField<Type, fvPatchField, volMesh>*>(fieldIter())
-            ->storeOldTimes();
-    }
-
-
     for
     (
         typename HashTable<const GeometricField<Type, fvPatchField, volMesh>*>::
@@ -589,29 +565,6 @@ void Foam::fvMeshAdder::MapSurfaceFields
         meshToAdd.objectRegistry::lookupClass<fldType>()
     );
 
-    // It is necessary to enforce that all old-time fields are stored
-    // before the mapping is performed.  Otherwise, if the
-    // old-time-level field is mapped before the field itself, sizes
-    // will not match.
-
-    for
-    (
-        typename HashTable<const fldType*>::
-            iterator fieldIter = fields.begin();
-        fieldIter != fields.end();
-        ++fieldIter
-    )
-    {
-        if (debug)
-        {
-            Pout<< "MapSurfaceFields : Storing old time for "
-                << fieldIter()->name() << endl;
-        }
-
-        const_cast<fldType*>(fieldIter())->storeOldTimes();
-    }
-
-
     for
     (
         typename HashTable<const fldType*>::
@@ -908,29 +861,6 @@ void Foam::fvMeshAdder::MapPointFields
     HashTable<const fldType*> fields(mesh.thisDb().lookupClass<fldType>());
     HashTable<const fldType*> fieldsToAdd(meshToAdd.lookupClass<fldType>());
 
-    // It is necessary to enforce that all old-time fields are stored
-    // before the mapping is performed.  Otherwise, if the
-    // old-time-level field is mapped before the field itself, sizes
-    // will not match.
-
-    for
-    (
-        typename HashTable<const fldType*>::
-            iterator fieldIter = fields.begin();
-        fieldIter != fields.end();
-        ++fieldIter
-    )
-    {
-        if (debug)
-        {
-            Pout<< "MapPointFields : Storing old time for "
-                << fieldIter()->name() << endl;
-        }
-
-        const_cast<fldType*>(fieldIter())->storeOldTimes();
-    }
-
-
     for
     (
         typename HashTable<const fldType*>::

src/finiteVolume/fvMesh/fvMesh.C

@@ -254,6 +254,22 @@ void Foam::fvMesh::storeOldVol(const scalarField& V)
 }
 
 
+void Foam::fvMesh::storeOldTimeFields()
+{
+    storeOldTimeFields<PointField>();
+    storeOldTimeFields<VolField>();
+    storeOldTimeFields<SurfaceField>();
+}
+
+
+void Foam::fvMesh::nullOldestTimeFields()
+{
+    nullOldestTimeFields<PointField>();
+    nullOldestTimeFields<VolField>();
+    nullOldestTimeFields<SurfaceField>();
+}
+
+
 void Foam::fvMesh::clearOut()
 {
     clearGeom();
@@ -584,6 +600,16 @@ bool Foam::fvMesh::update()
 {
     if (!conformal()) stitcher_->disconnect(true, true);
 
+    if
+    (
+        stitcher_->stitches()
+     || topoChanger_->dynamic()
+     || distributor_->dynamic()
+    )
+    {
+        nullOldestTimeFields();
+    }
+
     const bool hasV00 = V00Ptr_;
     deleteDemandDrivenData(V00Ptr_);
 
@@ -625,6 +651,15 @@ bool Foam::fvMesh::move()
 {
     if (!conformal()) stitcher_->disconnect(true, true);
 
+    if (curTimeIndex_ < time().timeIndex() && stitcher_->stitches())
+    {
+        // Store all old-time fields. If we don't do this then we risk
+        // triggering a store in the middle of mapping and potentially
+        // overwriting a mapped old-time field with a not-yet-mapped new-time
+        // field.
+        storeOldTimeFields();
+    }
+
     // Do not set moving false
     // Once the mesh starts moving it is considered to be moving
     // for the rest of the run

src/finiteVolume/fvMesh/fvMesh.H

@@ -231,6 +231,36 @@ private:
             //- Preserve old volume(s)
             void storeOldVol(const scalarField&);
 
+            //- Return the list of current fields of type GeoField<Type>
+            //  i.e. fields that are not old-time or old-old-time etc.
+            template<class Type, template<class> class GeoField>
+            UPtrList<GeoField<Type>> curFields();
+
+            //- Store old times fields for all fields of a given type and class
+            template<class Type, template<class> class GeoField>
+            void storeOldTimeFields();
+
+            //- Store old times fields for all fields of a given class
+            template<template<class> class GeoField>
+            void storeOldTimeFields();
+
+            //- Store old times fields for all fields
+            //  Used prior to field mapping following mesh-motion
+            //  and NCC stitching
+            void storeOldTimeFields();
+
+            //- Null oldest time field for all fields of a given type and class
+            template<class Type, template<class> class GeoField>
+            void nullOldestTimeFields();
+
+            //- Null oldest time field for all fields of a given class
+            template<template<class> class GeoField>
+            void nullOldestTimeFields();
+
+            //- Null oldest time field for all fields
+            //  Used prior to field mapping following mesh topology change
+            void nullOldestTimeFields();
+
 
        // Make topological data
 

src/finiteVolume/fvMesh/fvMeshStitchers/fvMeshStitcher/fvMeshStitcher.C

@@ -1369,12 +1369,6 @@ bool Foam::fvMeshStitcher::connect
     // Make the mesh non-conformal
     mesh_.unconform(polyFacesBf, Sf, Cf);
 
-    // Store all old-time fields. If we don't do this then we risk triggering a
-    // store in the middle of mapping and potentially overwriting a mapped
-    // old-time field with a not-yet-mapped new-time field.
-    storeOldTimeFields<VolField>();
-    storeOldTimeFields<SurfaceField>();
-
     // Resize all the affected patch fields
     resizePatchFields<VolField>();
     resizePatchFields<SurfaceField>();

src/finiteVolume/fvMesh/fvMeshStitchers/fvMeshStitcher/fvMeshStitcher.H

@@ -174,14 +174,6 @@ private:
 
         // Field Mapping
 
-            //- Store old times for all fields of a given type and class
-            template<class Type, template<class> class GeoField>
-            void storeOldTimeFields();
-
-            //- Store old times for all fields of a given class
-            template<template<class> class GeoField>
-            void storeOldTimeFields();
-
             //- Resize the patch fields of a given type and class to match the
             //  mesh
             template<class Type, template<class> class GeoField>

src/finiteVolume/fvMesh/fvMeshStitchers/fvMeshStitcher/fvMeshStitcherTemplates.C

@@ -40,27 +40,6 @@ Description
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-template<class Type, template<class> class GeoField>
-void Foam::fvMeshStitcher::storeOldTimeFields()
-{
-    HashTable<GeoField<Type>*> fields(mesh_.lookupClass<GeoField<Type>>());
-    forAllIter(typename HashTable<GeoField<Type>*>, fields, iter)
-    {
-        iter()->storeOldTimes();
-    }
-}
-
-
-template<template<class> class GeoField>
-void Foam::fvMeshStitcher::storeOldTimeFields()
-{
-    #define StoreOldTimeFields(Type, nullArg) \
-        storeOldTimeFields<Type, GeoField>();
-    FOR_ALL_FIELD_TYPES(StoreOldTimeFields);
-    #undef StoreOldTimeFields
-}
-
-
 template<class Type, template<class> class GeoField>
 void Foam::fvMeshStitcher::resizePatchFields()
 {

src/finiteVolume/fvMesh/fvMeshTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2015-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -42,4 +42,68 @@ typename Foam::pTraits<Type>::labelType Foam::fvMesh::validComponents() const
 }
 
 
+template<class Type, template<class> class GeoField>
+Foam::UPtrList<GeoField<Type>> Foam::fvMesh::curFields()
+{
+    HashTable<GeoField<Type>*> fields(lookupClass<GeoField<Type>>());
+    UPtrList<GeoField<Type>> curFields(fields.size());
+
+    label i = 0;
+    forAllIter(typename HashTable<GeoField<Type>*>, fields, iter)
+    {
+        if (!iter()->isOldTime())
+        {
+            curFields.set(i++, iter());
+        }
+    }
+    curFields.setSize(i);
+
+    return curFields;
+}
+
+
+template<class Type, template<class> class GeoField>
+void Foam::fvMesh::storeOldTimeFields()
+{
+    UPtrList<GeoField<Type>> curFields(this->curFields<Type, GeoField>());
+
+    forAll(curFields, i)
+    {
+        curFields[i].storeOldTimes();
+    }
+}
+
+
+template<template<class> class GeoField>
+void Foam::fvMesh::storeOldTimeFields()
+{
+    #define StoreOldTimeFields(Type, nullArg) \
+        storeOldTimeFields<Type, GeoField>();
+    FOR_ALL_FIELD_TYPES(StoreOldTimeFields);
+    #undef StoreOldTimeFields
+}
+
+
+template<class Type, template<class> class GeoField>
+void Foam::fvMesh::nullOldestTimeFields()
+{
+    UPtrList<GeoField<Type>> curFields(this->curFields<Type, GeoField>());
+
+    forAll(curFields, i)
+    {
+        curFields[i].nullOldestTime();
+    }
+}
+
+
+template<template<class> class GeoField>
+void Foam::fvMesh::nullOldestTimeFields()
+{
+    #define nullOldestTimeFields(Type, nullArg) \
+        nullOldestTimeFields<Type, GeoField>();
+    FOR_ALL_FIELD_TYPES(nullOldestTimeFields);
+    #undef nullOldestTimeFields
+}
+
+
 // ************************************************************************* //