VolumeFromBoundary surviving for edges

Author: JSchoeberl

comp/meshaccess.cpp

@@ -2435,30 +2435,35 @@ namespace ngcomp
           vol_cf->Evaluate(ir, values);
           return;
         }
-      const MeshAccess & ma = *static_cast<const MeshAccess*> (trafo.GetMesh());
-      int facet = ma.GetElFacets(ei)[0];
-      ArrayMem<int,2> elnums;
-      ma.GetFacetElements (facet, elnums);
 
-      for (auto elnr : elnums)
+      const MeshAccess & ma = *static_cast<const MeshAccess*> (trafo.GetMesh());
+      if (ei.IsBoundary())
         {
-          ElementId vei(VOL, elnr);
-          int locfacnr = ma.GetElFacets(vei).Pos(facet);
-
-          ElementTransformation & vol_trafo = ma.GetTrafo (vei, lh);
-          if (!vol_cf->DefinedOn(vol_trafo)) continue;
-
-          Facet2ElementTrafo f2el(vol_trafo.GetElementType(), ma.GetElVertices(vei));
-          Array<int> surfvnums { ma.GetElVertices(ei) };
-          Facet2SurfaceElementTrafo f2sel(trafo.GetElementType(), surfvnums);
-          auto & ir_ref = f2sel.Inverse(ir.IR(), lh);
-          auto & ir_vol = f2el(locfacnr, ir_ref, lh);
-          auto & mir_vol = vol_trafo(ir_vol, lh);
-          mir_vol.ComputeNormalsAndMeasure (vol_trafo.GetElementType(), locfacnr);
-          vol_cf -> Evaluate (mir_vol, values);
-          return;
+          int facet = ma.GetElFacets(ei)[0];
+          ArrayMem<int,2> elnums;
+          ma.GetFacetElements (facet, elnums);
+          
+          for (auto elnr : elnums)
+            {
+              ElementId vei(VOL, elnr);
+              int locfacnr = ma.GetElFacets(vei).Pos(facet);
+              
+              ElementTransformation & vol_trafo = ma.GetTrafo (vei, lh);
+              if (!vol_cf->DefinedOn(vol_trafo)) continue;
+              
+              Facet2ElementTrafo f2el(vol_trafo.GetElementType(), ma.GetElVertices(vei));
+              Array<int> surfvnums { ma.GetElVertices(ei) };
+              Facet2SurfaceElementTrafo f2sel(trafo.GetElementType(), surfvnums);
+              auto & ir_ref = f2sel.Inverse(ir.IR(), lh);
+              auto & ir_vol = f2el(locfacnr, ir_ref, lh);
+              auto & mir_vol = vol_trafo(ir_vol, lh);
+              mir_vol.ComputeNormalsAndMeasure (vol_trafo.GetElementType(), locfacnr);
+              vol_cf -> Evaluate (mir_vol, values);
+              return;
+            }
         }
-
+      
+      values.AddSize(Dimension(), ir.Size()) = T(0.0);
     }
 
     template <typename MIR, typename T, ORDERING ORD>