Bug fixed for post-processing meshes with topological changes

applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3FoamFaceField.H

@@ -53,6 +53,7 @@ void Foam::vtkPV3Foam::convertFaceField
 {
     const label nComp = pTraits<Type>::nComponents;
     const label nInternalFaces = mesh.nInternalFaces();
+    const label nFaces = mesh.nFaces();
     const labelList& faceOwner = mesh.faceOwner();
     const labelList& faceNeigh = mesh.faceNeighbour();
 
@@ -77,22 +78,31 @@ void Foam::vtkPV3Foam::convertFaceField
     // for boundary faces: owner
     forAll(faceLabels, faceI)
     {
+        // Bug fix: filter inactive faces.  HJ, 21/Mar/2011
         const label faceNo = faceLabels[faceI];
+
         if (faceNo < nInternalFaces)
         {
             Type t = 0.5*(tf[faceOwner[faceNo]] + tf[faceNeigh[faceNo]]);
 
-            for (direction d=0; d<nComp; d++)
+            for (direction d = 0; d < nComp; d++)
+            {
+                vec[d] = component(t, d);
+            }
+        }
+        else if (faceNo < nFaces)
+        {
+            const Type& t = tf[faceOwner[faceNo]];
+            for (direction d = 0; d < nComp; d++)
             {
                 vec[d] = component(t, d);
             }
         }
         else
         {
-            const Type& t = tf[faceOwner[faceNo]];
+            for (direction d = 0; d < nComp; d++)
-            for (direction d=0; d<nComp; d++)
             {
-                vec[d] = component(t, d);
+                vec[d] = 0;
             }
         }
         vtkOpenFOAMTupleRemap<Type>(vec);
@@ -124,6 +134,7 @@ void Foam::vtkPV3Foam::convertFaceField
 {
     const label nComp = pTraits<Type>::nComponents;
     const label nInternalFaces = mesh.nInternalFaces();
+    const label nFaces = mesh.nFaces();
     const labelList& faceOwner = mesh.faceOwner();
     const labelList& faceNeigh = mesh.faceNeighbour();
 
@@ -151,21 +162,29 @@ void Foam::vtkPV3Foam::convertFaceField
     {
         const label faceNo = iter.key();
 
+        // Bug fix: filter inactive faces.  HJ, 21/Mar/2011
         if (faceNo < nInternalFaces)
         {
             Type t = 0.5*(tf[faceOwner[faceNo]] + tf[faceNeigh[faceNo]]);
 
-            for (direction d=0; d<nComp; d++)
+            for (direction d = 0; d < nComp; d++)
+            {
+                vec[d] = component(t, d);
+            }
+        }
+        else if (faceNo < nFaces)
+        {
+            const Type& t = tf[faceOwner[faceNo]];
+            for (direction d = 0; d < nComp; d++)
             {
                 vec[d] = component(t, d);
             }
         }
         else
         {
-            const Type& t = tf[faceOwner[faceNo]];
+            for (direction d = 0; d < nComp; d++)
-            for (direction d=0; d<nComp; d++)
             {
-                vec[d] = component(t, d);
+                vec[d] = 0;
             }
         }
         vtkOpenFOAMTupleRemap<Type>(vec);

applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3FoamMeshZone.C

@@ -54,12 +54,22 @@ vtkPolyData* Foam::vtkPV3Foam::faceZoneVTKMesh
 
     // Construct primitivePatch of faces in faceZone
 
-    const faceList& meshFaces = mesh.faces();
+    const faceList& meshFaces = mesh.allFaces();
     faceList patchFaces(faceLabels.size());
+    label npf = 0;
+
+    // Filter faces that are not in live mesh
+    // Bug fix.  HJ, 21/Mar/2011
     forAll(faceLabels, faceI)
     {
-        patchFaces[faceI] = meshFaces[faceLabels[faceI]];
+        if (faceLabels[faceI] < mesh.nFaces())
+        {
+            patchFaces[npf] = meshFaces[faceLabels[faceI]];
+            npf++;
         }
+    }
+    patchFaces.setSize(npf);
+
     primitiveFacePatch p(patchFaces, mesh.points());
 
 
@@ -69,7 +79,7 @@ vtkPolyData* Foam::vtkPV3Foam::faceZoneVTKMesh
     const pointField& points = p.localPoints();
 
     vtkPoints* vtkpoints = vtkPoints::New();
-    vtkpoints->Allocate( points.size() );
+    vtkpoints->Allocate(points.size());
     forAll(points, i)
     {
         vtkPV3FoamInsertNextPoint(vtkpoints, points[i]);
@@ -124,14 +134,28 @@ vtkPolyData* Foam::vtkPV3Foam::pointZoneVTKMesh
         printMemory();
     }
 
-    const pointField& meshPoints = mesh.points();
+    const pointField& meshPoints = mesh.allPoints();
+
+    // Filter point labels to include only live points
+    labelList pl(pointLabels.size());
+    label npl = 0;
+
+    forAll (pointLabels, pointI)
+    {
+        if (pointLabels[pointI] < mesh.nPoints())
+        {
+            pl[npl] = pointLabels[pointI];
+            npl++;
+        }
+    }
+    pl.setSize(npl);
 
     vtkPoints *vtkpoints = vtkPoints::New();
-    vtkpoints->Allocate( pointLabels.size() );
+    vtkpoints->Allocate( pl.size());
 
     forAll(pointLabels, pointI)
     {
-        vtkPV3FoamInsertNextPoint(vtkpoints, meshPoints[pointLabels[pointI]]);
+        vtkPV3FoamInsertNextPoint(vtkpoints, meshPoints[pl[pointI]]);
     }
 
     vtkmesh->SetPoints(vtkpoints);