parent faces

libsrc/include/nginterface_v2.hpp

@@ -345,6 +345,8 @@ namespace netgen
     void Curve (int order);
     int GetCurveOrder ();
 
+    void EnableTable (string name, bool set);
+
     void Refine (NG_REFINEMENT_TYPE reftype,
                  void (*taskmanager)(function<void(int,int)>) = &DummyTaskManager2,
                  void (*tracer)(string, bool) = &DummyTracer2);

libsrc/interface/nginterface_v2.cpp

@@ -1137,6 +1137,13 @@ namespace netgen
     return mesh->GetCurvedElements().GetOrder();
   }
 
+  void Ngx_Mesh :: EnableTable (string name, bool set)
+  {
+    mesh->GetTopology().EnableTable (name, set);
+  }
+
+
+  
   template <>
   DLL_HEADER void Ngx_Mesh :: SetRefinementFlag<2> (size_t elnr, bool flag)
   {

libsrc/meshing/python_mesh.cpp

@@ -1129,15 +1129,20 @@ project_boundaries : Optional[str] = None
 
     .def ("EnableTable", [] (Mesh & self, string name, bool set)
           {
+            const_cast<MeshTopology&>(self.GetTopology()).EnableTable(name, set);
+            /*
             if (name == "edges")
               const_cast<MeshTopology&>(self.GetTopology()).SetBuildEdges(set);
             else if (name == "faces")
               const_cast<MeshTopology&>(self.GetTopology()).SetBuildFaces(set);
             else if (name == "parentedges")
-              const_cast<MeshTopology&>(self.GetTopology()).SetBuildHierarchy(set);
+              const_cast<MeshTopology&>(self.GetTopology()).SetBuildParentEdges(set);
+            else if (name == "parentfaces")
+              const_cast<MeshTopology&>(self.GetTopology()).SetBuildParentFaces(set);
 	    else
 	      throw Exception ("noting known about table "+name +"\n"
-			       "knwon are 'edges', 'faces', 'parentedges'");
+			       "knwon are 'edges', 'faces', 'parentedges', 'parentfaces'");
+            */
           },
           py::arg("name"), py::arg("set")=true)
     

libsrc/meshing/topology.cpp

@@ -55,6 +55,24 @@ namespace netgen
   bool MeshTopology :: NeedsUpdate() const
   { return (timestamp <= mesh->GetTimeStamp()); }
 
+  
+  void MeshTopology :: EnableTable (string name, bool set)
+  {
+    if (name == "edges")
+      SetBuildEdges(set);
+    else if (name == "faces")
+      SetBuildFaces(set);
+    else if (name == "parentedges")
+      SetBuildParentEdges(set);
+    else if (name == "parentfaces")
+      SetBuildParentFaces(set);
+    else
+      throw Exception ("noting known about table "+name +"\n"
+                       "knwon are 'edges', 'faces', 'parentedges', 'parentfaces'");
+  }
+
+
+  
   template <typename FUNC>
   void LoopOverEdges (const Mesh & mesh, MeshTopology & top, PointIndex v,
                       FUNC func)
@@ -679,7 +697,7 @@ namespace netgen
            } );
 
 
-	if (build_hierarchy)
+	if (build_parent_edges)
 	  {
 	    static Timer t("build_hierarchy"); RegionTimer reg(t);
 	    cnt = 0;
@@ -1438,9 +1456,11 @@ namespace netgen
         // NgProfiler::StopTimer (timer2c);
 
 
-        if (build_hierarchy)
+        if (build_parent_faces)
           {
             // tets only
+            cout << "build face hierarchy:" << endl;
+            cout << "f2v = " << face2vert << endl;
             
             ngcore::ClosedHashTable<INT<3>, int> v2f(nv);
             for (auto i : Range(face2vert))
@@ -1459,11 +1479,16 @@ namespace netgen
             for (auto i : Range(nfa))
               {
                 INT<3,PointIndex> f3(face2vert[i][0], face2vert[i][1], face2vert[i][2]);
-                PointIndex vmax = Max(f3);
-                if (vmax >= mesh->mlbetweennodes.Size()+PointIndex::BASE)
-                  continue;
 
-                auto parents = mesh->mlbetweennodes[vmax];
+                // find a vertex, such that one of its parent is a trig vertex
+
+                for (int k = 0; k < 3; k++)
+                  {
+                    PointIndex v3 = f3[k];
+                    
+                    if (v3 >= mesh->mlbetweennodes.Size()+PointIndex::BASE)
+                      continue;
+                    auto parents = mesh->mlbetweennodes[v3];
 
                     // is face part of one parent face (boundary-face) ?
                     for (int j = 0; j < 2; j++)
@@ -1474,23 +1499,29 @@ namespace netgen
                             PointIndex v1 = parents[1-j];
                             
                             // the third one, on the tip
-                        PointIndex v2 = f3[0]+f3[1]+f3[2] - v0 - v1;
+                            PointIndex v2 = f3[0]+f3[1]+f3[2] - v0 - v3;
-                        
-                        int classnr = 0;
-                        if (v0 > v1) { Swap (v0, v1); classnr += 1; }
-                        if (v1 > v2) { Swap (v1, v2); classnr += 2; }
-                        if (v0 > v1) { Swap (v0, v1); classnr += 4; }
 
                             INT<3> parentverts(v0, v1, v2);
+                            parentverts.Sort();
+
+                            
+                            int classnr = 0;
+                            if (v2 > v3) { Swap (v2, v3); classnr += 1; }                            
+                            if (v0 > v1) { Swap (v0, v1); classnr += 2; }
+                            if (v1 > v2) { Swap (v1, v2); classnr += 4; }
+                            if (v0 > v1) { Swap (v0, v1); classnr += 8; }
+
+                            if (v2f.Used(parentverts))
+                              {
                                 int pafacenr = v2f[parentverts];
                                 cout << "parent-face = " << pafacenr << endl;
                                 parent_faces[i] = { classnr, { pafacenr, -1, -1, -1 } };
                               }
                           }
                       }
-                   
                   }
-        
+              }
+          }
       }
     
 

libsrc/meshing/topology.hpp

@@ -45,7 +45,8 @@ class MeshTopology
   const Mesh * mesh;
   bool buildedges;
   bool buildfaces;
-  bool build_hierarchy = false; // may be changed to default = false
+  bool build_parent_edges = false; // may be changed to default = false
+  bool build_parent_faces = false; // may be changed to default = false
 
   NgArray<INDEX_2> edge2vert;
   NgArray<INDEX_4> face2vert;
@@ -80,14 +81,14 @@ public:
   { buildedges = be; }
   void SetBuildFaces (bool bf)
   { buildfaces = bf; }
-  void SetBuildHierarchy (bool bh) { build_hierarchy = bh; }
+  void SetBuildParentEdges (bool bh) { build_parent_edges = bh; }
+  void SetBuildParentFaces (bool bh) { build_parent_faces = bh; }
 
-  bool HasEdges () const
+  void EnableTable (string name, bool set);
-  { return buildedges; }
+
-  bool HasFaces () const
+  bool HasEdges () const  { return buildedges; }
-  { return buildfaces; }
+  bool HasFaces () const  { return buildfaces; }
-  bool HasParentEdges () const
+  bool HasParentEdges () const { return build_parent_edges; }
-  { return build_hierarchy; }
 
   void Update(NgTaskManager tm = &DummyTaskManager, NgTracer tracer = &DummyTracer);
   bool NeedsUpdate() const;