From 9b969cd4bf0503e6355cbbeb7a6705b56705f7c9 Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@esi-group.com>
Date: Mon, 14 May 2018 17:01:54 +0100
Subject: [PATCH] ENH: overset movement for patches (issue #2)

BUG: general patch movement was inconsistent
---
 src/catalyst/volMesh/foamVtkFvMeshAdaptor.H   |   6 +
 .../foamVtkFvMeshAdaptorFieldTemplates.C      |  79 ++++-----
 .../volMesh/foamVtkFvMeshAdaptorGeom.C        | 166 +++++++++++++-----
 3 files changed, 167 insertions(+), 84 deletions(-)

diff --git a/src/catalyst/volMesh/foamVtkFvMeshAdaptor.H b/src/catalyst/volMesh/foamVtkFvMeshAdaptor.H
index 025f14a..9ebf492 100644
--- a/src/catalyst/volMesh/foamVtkFvMeshAdaptor.H
+++ b/src/catalyst/volMesh/foamVtkFvMeshAdaptor.H
@@ -184,6 +184,12 @@ private:
         //  There will be several for groups, but only one for regular patches.
         void convertGeometryPatches();
 
+        //- Ghost cells - internal mesh
+        void applyGhostingInternal(const labelUList& types);
+
+        //- Ghost cells - patches
+        void applyGhostingPatches(const labelUList& types);
+
         //- Ghost cells to affect the visibility of the geometry.
         //  Currently only used for overset.
         void applyGhosting();
diff --git a/src/catalyst/volMesh/foamVtkFvMeshAdaptorFieldTemplates.C b/src/catalyst/volMesh/foamVtkFvMeshAdaptorFieldTemplates.C
index 4914144..57717e8 100644
--- a/src/catalyst/volMesh/foamVtkFvMeshAdaptorFieldTemplates.C
+++ b/src/catalyst/volMesh/foamVtkFvMeshAdaptorFieldTemplates.C
@@ -83,18 +83,13 @@ void Foam::vtk::fvMeshAdaptor::convertVolField
         if (!iter.found() || !iter.object().dataset)
         {
             // Should not happen, but for safety require a vtk geometry
+            Pout<<"Cache miss for VTP patch " << longName << endl;
             continue;
         }
 
         foamVtpData& vtpData = iter.object();
         auto dataset = vtpData.dataset;
 
-        const labelList& patchIds = vtpData.additionalIds();
-        if (patchIds.empty())
-        {
-            continue;
-        }
-
         // This is slightly roundabout, but we deal with groups and with
         // single patches.
         // For groups (spanning several patches) it is fairly messy to
@@ -111,53 +106,48 @@ void Foam::vtk::fvMeshAdaptor::convertVolField
         );
 
         vtkSmartPointer<vtkFloatArray> pdata;
-        const bool allowPdata = (interpField && patchIds.size() == 1);
 
-        label coffset = 0; // The write offset into cell-data
-        for (label patchId : patchIds)
-        {
-            const fvPatchField<Type>& ptf = fld.boundaryField()[patchId];
+        const fvPatchField<Type>& ptf = fld.boundaryField()[patchId];
 
-            if
+        if
+        (
+            isType<emptyFvPatchField<Type>>(ptf)
+         ||
             (
-                isType<emptyFvPatchField<Type>>(ptf)
-             ||
-                (
-                    extrapPatch
-                 && !polyPatch::constraintType(patches[patchId].type())
-                )
+                extrapPatch
+             && !polyPatch::constraintType(patches[patchId].type())
             )
-            {
-                fvPatch p(ptf.patch().patch(), mesh.boundary());
+        )
+        {
+            fvPatch p(ptf.patch().patch(), mesh.boundary());
 
-                tmp<Field<Type>> tpptf
-                (
-                    fvPatchField<Type>(p, fld).patchInternalField()
-                );
+            tmp<Field<Type>> tpptf
+            (
+                fvPatchField<Type>(p, fld).patchInternalField()
+            );
 
-                coffset += transcribeFloatData(cdata, tpptf(), coffset);
+            transcribeFloatData(cdata, tpptf());
 
-                if (allowPdata && patchId < patchInterpList.size())
-                {
-                    pdata = vtk::Tools::convertFieldToVTK
-                    (
-                        fld.name(),
-                        patchInterpList[patchId].faceToPointInterpolate(tpptf)()
-                    );
-                }
+            if (interpField && patchId < patchInterpList.size())
+            {
+                pdata = vtk::Tools::convertFieldToVTK
+                (
+                    fld.name(),
+                    patchInterpList[patchId].faceToPointInterpolate(tpptf)()
+                );
             }
-            else
+        }
+        else
+        {
+            transcribeFloatData(cdata, ptf);
+
+            if (interpField && patchId < patchInterpList.size())
             {
-                coffset += transcribeFloatData(cdata, ptf, coffset);
-
-                if (allowPdata && patchId < patchInterpList.size())
-                {
-                    pdata = vtk::Tools::convertFieldToVTK
-                    (
-                        fld.name(),
-                        patchInterpList[patchId].faceToPointInterpolate(ptf)()
-                    );
-                }
+                pdata = vtk::Tools::convertFieldToVTK
+                (
+                    fld.name(),
+                    patchInterpList[patchId].faceToPointInterpolate(ptf)()
+                );
             }
         }
 
@@ -214,6 +204,7 @@ void Foam::vtk::fvMeshAdaptor::convertVolFieldInternal
     if (!iter.found() || !iter.object().dataset)
     {
         // Should not happen, but for safety require a vtk geometry
+        Pout<<"Cache miss for VTU " << longName << endl;
         return;
     }
     foamVtuData& vtuData = iter.object();
diff --git a/src/catalyst/volMesh/foamVtkFvMeshAdaptorGeom.C b/src/catalyst/volMesh/foamVtkFvMeshAdaptorGeom.C
index e74d6e9..7867c12 100644
--- a/src/catalyst/volMesh/foamVtkFvMeshAdaptorGeom.C
+++ b/src/catalyst/volMesh/foamVtkFvMeshAdaptorGeom.C
@@ -56,7 +56,7 @@ void Foam::vtk::fvMeshAdaptor::convertGeometryInternal()
             {
                 Info<< "reuse " << longName << nl;
             }
-            vtuData.reuse(); // reuse
+            vtuData.reuse();  // No movement - simply reuse
             return;
         }
         else if (meshState_ == polyMesh::POINTS_MOVED)
@@ -87,6 +87,8 @@ void Foam::vtk::fvMeshAdaptor::convertGeometryInternal()
 
 void Foam::vtk::fvMeshAdaptor::convertGeometryPatches()
 {
+    // PATCHES
+
     const polyBoundaryMesh& patches = mesh_.boundaryMesh();
     const label npatches = this->nPatches();
 
@@ -102,28 +104,22 @@ void Foam::vtk::fvMeshAdaptor::convertGeometryPatches()
         {
             if (meshState_ == polyMesh::UNCHANGED)
             {
-                // Without movement is easy.
                 if (debug)
                 {
                     Info<< "reuse " << longName << nl;
                 }
-                vtpData.reuse();
+                vtpData.reuse();  // No movement - simply reuse
                 continue;
             }
             else if (meshState_ == polyMesh::POINTS_MOVED)
             {
-                // Point movement on single patch is OK
-
-                const labelList& patchIds = vtpData.additionalIds();
-                if (patchIds.size() == 1)
+                // Point movement on single patch
+                if (debug)
                 {
-                    vtkgeom = vtpData.getCopy();
-                    vtkgeom->SetPoints
-                    (
-                        vtk::Tools::Patch::points(patches[patchIds[0]])
-                    );
-                    continue;
+                    Info<< "move points " << longName << nl;
                 }
+                vtkgeom = vtpData.getCopy();
+                vtkgeom->SetPoints(vtk::Tools::Patch::points(pp));
             }
         }
 
@@ -135,13 +131,13 @@ void Foam::vtk::fvMeshAdaptor::convertGeometryPatches()
                 << longName << endl;
         }
 
-        // Store good patch id as additionalIds
-        vtpData.additionalIds() = {patchId};
+        // Unused information
+        vtpData.additionalIds().clear();
 
         // This is somewhat inconsistent, since we currently only have
         // normal (non-grouped) patches but this may change in the future.
 
-        vtkgeom = vtk::Tools::Patch::mesh(patches[patchId]);
+        vtkgeom = vtk::Tools::Patch::mesh(pp);
 
         if (vtkgeom)
         {
@@ -156,22 +152,11 @@ void Foam::vtk::fvMeshAdaptor::convertGeometryPatches()
 }
 
 
-// These need to be rebuild here, since the component mesh may just have point
-// motion without topology changes.
-void Foam::vtk::fvMeshAdaptor::applyGhosting()
+void Foam::vtk::fvMeshAdaptor::applyGhostingInternal(const labelUList& types)
 {
-    if (!usingVolume())
-    {
-        return;
-    }
-
-    const auto* stencilPtr =
-        mesh_.lookupObjectPtr<cellCellStencilObject>
-        (
-            cellCellStencilObject::typeName
-        );
+    // MESH = "internal"
 
-    if (!stencilPtr)
+    if (types.empty() || !usingVolume())
     {
         return;
     }
@@ -182,12 +167,12 @@ void Foam::vtk::fvMeshAdaptor::applyGhosting()
     if (!iter.found() || !iter.object().dataset)
     {
         // Should not happen, but for safety require a vtk geometry
+        Pout<<"Cache miss for VTU " << longName << endl;
         return;
     }
     foamVtuData& vtuData = iter.object();
     auto dataset = vtuData.dataset;
 
-    const auto& stencil = *stencilPtr;
     const labelUList& cellMap = vtuData.cellMap();
 
     auto vtkgcell = dataset->GetCellGhostArray();
@@ -202,38 +187,139 @@ void Foam::vtk::fvMeshAdaptor::applyGhosting()
     //     vtkgpoint = dataset->AllocatePointGhostArray();
     // }
 
-    UList<uint8_t> gcell =
-        vtk::Tools::asUList(vtkgcell, cellMap.size());
+    UList<uint8_t> gcell = vtk::Tools::asUList(vtkgcell, cellMap.size());
 
     vtkgcell->FillValue(0); // Initialize to zero
 
-    const labelUList& types = stencil.cellTypes();
-
     forAll(cellMap, i)
     {
         const label cellType = types[cellMap[i]];
 
-        if (cellType == cellCellStencil::INTERPOLATED)
+        if (cellType == cellCellStencil::HOLE)
         {
+            // Need duplicate (not just HIDDENCELL) for it to be properly
+            // recognized
             gcell[i] |=
             (
                 vtkDataSetAttributes::DUPLICATECELL
+              | vtkDataSetAttributes::HIDDENCELL
             );
         }
-        else if (cellType == cellCellStencil::HOLE)
+        #if 0
+        // No special treatment for INTERPOLATED.
+        // This causes duplicate/overlapping values, but avoids holes
+        // in the results
+        else if (cellType == cellCellStencil::INTERPOLATED)
         {
-            // Need duplicate (not just HIDDENCELL) for it to be properly
-            // recognized
             gcell[i] |=
             (
                 vtkDataSetAttributes::DUPLICATECELL
-              | vtkDataSetAttributes::HIDDENCELL
             );
         }
+        #endif
     }
 
     dataset->GetCellData()->AddArray(vtkgcell);
 }
 
 
+void Foam::vtk::fvMeshAdaptor::applyGhostingPatches(const labelUList& types)
+{
+    // PATCHES
+
+    if (types.empty())
+    {
+        return;
+    }
+
+    const polyBoundaryMesh& patches = mesh_.boundaryMesh();
+    const label npatches = this->nPatches();
+
+    for (label patchId=0; patchId < npatches; ++patchId)
+    {
+        const polyPatch& pp = patches[patchId];
+        const word& longName = pp.name();
+
+        auto iter = cachedVtp_.find(longName);
+        if (!iter.found() || !iter.object().dataset)
+        {
+            // Should not happen, but for safety require a vtk geometry
+            Pout<<"Cache miss for VTP patch " << longName << endl;
+            continue;
+        }
+
+        foamVtpData& vtpData = iter.object();
+        auto dataset = vtpData.dataset;
+
+        auto vtkgcell = dataset->GetCellGhostArray();
+        if (!vtkgcell)
+        {
+            vtkgcell = dataset->AllocateCellGhostArray();
+        }
+
+
+        // Determine face ghosting based on interior cells
+        const labelUList& bCells = pp.faceCells();
+
+        const label len = bCells.size();
+
+        UList<uint8_t> gcell = vtk::Tools::asUList(vtkgcell, len);
+
+        vtkgcell->FillValue(0); // Initialize to zero
+
+        for (label i=0; i < len; ++i)
+        {
+            const label celli = bCells[i];
+
+            const label cellType = types[celli];
+
+            if (cellType == cellCellStencil::HOLE)
+            {
+                // Need duplicate (not just HIDDENCELL) for it to be properly
+                // recognized
+                gcell[i] |=
+                (
+                    vtkDataSetAttributes::DUPLICATECELL
+                  | vtkDataSetAttributes::HIDDENCELL
+                );
+            }
+            #if 0
+            // No special treatment for INTERPOLATED.
+            // This causes duplicate/overlapping values, but avoids holes
+            // in the results
+            else if (cellType == cellCellStencil::INTERPOLATED)
+            {
+                gcell[i] |=
+                (
+                    vtkDataSetAttributes::DUPLICATECELL
+                );
+            }
+            #endif
+        }
+
+        dataset->GetCellData()->AddArray(vtkgcell);
+    }
+}
+
+
+// These need to be rebuild here, since the component mesh may just have point
+// motion without topology changes.
+void Foam::vtk::fvMeshAdaptor::applyGhosting()
+{
+    const auto* stencilPtr =
+        mesh_.lookupObjectPtr<cellCellStencilObject>
+        (
+            cellCellStencilObject::typeName
+        );
+
+    if (stencilPtr)
+    {
+        const labelUList& types = stencilPtr->cellTypes();
+
+        applyGhostingInternal(types);
+        applyGhostingPatches(types);
+    }
+}
+
+
 // ************************************************************************* //
-- 
GitLab