From 4bfe42b305e25029a7a17b32edcfa241c24f2962 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Joachim=20Sch=C3=B6berl?= <joachim.schoeberl@tuwien.ac.at>
Date: Wed, 7 Aug 2019 23:32:14 +0200
Subject: [PATCH] more PointIndex::BASE 0/1 fixes

---
 libsrc/meshing/improve2.cpp | 13 ++++++----
 libsrc/meshing/improve3.cpp |  4 +--
 libsrc/meshing/meshtype.cpp | 20 +++++++-------
 libsrc/occ/occgenmesh.cpp   | 52 ++++++++++++++++++++++---------------
 4 files changed, 50 insertions(+), 39 deletions(-)

diff --git a/libsrc/meshing/improve2.cpp b/libsrc/meshing/improve2.cpp
index 582839e3..0dd3c311 100644
--- a/libsrc/meshing/improve2.cpp
+++ b/libsrc/meshing/improve2.cpp
@@ -562,10 +562,13 @@ namespace netgen
 	    PointIndex pi1 = elem[j];
 	    PointIndex pi2 = elem[(j+1) % 3];
 
+            /*
 	    if (pi1 < PointIndex::BASE || 
 		pi2 < PointIndex::BASE)
 	      continue;
-
+            */
+            if (!pi1.IsValid() || !pi2.IsValid())
+	      continue;              
 	    /*
 	      INDEX_2 i2(pi1, pi2);
 	      i2.Sort();
@@ -597,8 +600,8 @@ namespace netgen
 
 	    double loch = mesh.GetH (mesh[pi1]);
 
-	    INDEX_2 si2 (pi1, pi2);
-	    si2.Sort();
+	    // INDEX_2 si2 (pi1, pi2);
+	    // si2.Sort();
 
 	    /*	  
 	      if (edgetested.Used (si2))
@@ -737,8 +740,8 @@ namespace netgen
 	      
 		Element2d *el1p(NULL);
 		int l = 0;
-		while(mesh[elementsonnode[pi1][l]].IsDeleted() && l<elementsonnode.EntrySize(pi1)) l++;
-		if(l<elementsonnode.EntrySize(pi1))
+		while(mesh[elementsonnode[pi1][l]].IsDeleted() && l<elementsonnode[pi1].Size()) l++;
+		if(l<elementsonnode[pi1].Size())
 		  el1p = &mesh[elementsonnode[pi1][l]];
 		else
 		  cerr << "OOPS!" << endl;
diff --git a/libsrc/meshing/improve3.cpp b/libsrc/meshing/improve3.cpp
index 8d6482d0..5774caf9 100644
--- a/libsrc/meshing/improve3.cpp
+++ b/libsrc/meshing/improve3.cpp
@@ -748,7 +748,7 @@ void MeshOptimize3d :: SwapImprove (Mesh & mesh, OPTIMIZEGOAL goal,
 		  el31[3] = pi4;
 		}
 	      
-	      pi5 = 0;
+	      pi5.Invalidate();
 	      for (int k = 0; k < 3; k++)   // JS, 201212
 		{
 		  const Element & elemk = mesh[hasbothpoints[k]];
@@ -2315,7 +2315,7 @@ void MeshOptimize3d :: SwapImprove2 (Mesh & mesh, OPTIMIZEGOAL goal)
 {
   static Timer t("MeshOptimize3d::SwapImprove2"); RegionTimer reg(t);
   
-  PointIndex pi1(0), pi2(0), pi3(0), pi4(0), pi5(0);
+  PointIndex pi1, pi2, pi3, pi4, pi5;
   Element el21(TET), el22(TET), el31(TET), el32(TET), el33(TET);
 
   int cnt = 0;
diff --git a/libsrc/meshing/meshtype.cpp b/libsrc/meshing/meshtype.cpp
index 3d818647..95127c7c 100644
--- a/libsrc/meshing/meshtype.cpp
+++ b/libsrc/meshing/meshtype.cpp
@@ -185,7 +185,7 @@ namespace netgen
   { 
     for (int i = 0; i < ELEMENT2D_MAXPOINTS; i++)
       {
-	pnum[i] = 0;
+        pnum[i].Invalidate();
 	geominfo[i].trignum = 0;
       }
     np = anp;
@@ -210,7 +210,7 @@ namespace netgen
   { 
     for (int i = 0; i < ELEMENT2D_MAXPOINTS; i++)
       {
-	pnum[i] = 0;
+        pnum[i].Invalidate();
 	geominfo[i].trignum = 0;
       }
 
@@ -235,9 +235,9 @@ namespace netgen
     pnum[2] = pi3;
     np = 3;
     typ = TRIG;
-    pnum[3] = 0;
-    pnum[4] = 0;
-    pnum[5] = 0;
+    
+    for (int i = 0; i < ELEMENT2D_MAXPOINTS; i++)
+        pnum[i].Invalidate();
   
     for (int i = 0; i < ELEMENT2D_MAXPOINTS; i++)
       geominfo[i].trignum = 0;
@@ -1000,9 +1000,8 @@ namespace netgen
   Element :: Element (int anp)
   {
     np = anp;
-    int i;
-    for (i = 0; i < ELEMENT_MAXPOINTS; i++)
-      pnum[i] = 0;
+    for (int i = 0; i < ELEMENT_MAXPOINTS; i++)
+        pnum[i].Invalidate();
     index = 0;
     flags.marked = 1;
     flags.badel = 0;
@@ -1051,9 +1050,8 @@ namespace netgen
   {
     SetType (type);
 
-    int i;
-    for (i = 0; i < ELEMENT_MAXPOINTS; i++)
-      pnum[i] = 0;
+    for (int i = 0; i < ELEMENT_MAXPOINTS; i++)
+        pnum[i].Invalidate();
     index = 0;
     flags.marked = 1;
     flags.badel = 0;
diff --git a/libsrc/occ/occgenmesh.cpp b/libsrc/occ/occgenmesh.cpp
index d80f35e0..1f330201 100644
--- a/libsrc/occ/occgenmesh.cpp
+++ b/libsrc/occ/occgenmesh.cpp
@@ -331,7 +331,8 @@ namespace netgen
 
         bool exists = 0;
         if (merge_solids)
-          for (PointIndex pi = 1; pi <= mesh.GetNP(); pi++)
+          // for (PointIndex pi = 1; pi <= mesh.GetNP(); pi++)
+          for (PointIndex pi : Range(mesh.Points()))
             if ( Dist2 (mesh[pi], Point<3>(mp)) < eps*eps)
               {
                 exists = 1;
@@ -346,7 +347,9 @@ namespace netgen
     (*testout) << "different vertices = " << mesh.GetNP() << endl;
 
 
-    int first_ep = mesh.GetNP()+1;
+    // int first_ep = mesh.GetNP()+1;
+    PointIndex first_ep = mesh.Points().End();
+    auto vertexrange = mesh.Points().Range();
 
     NgArray<int> face2solid[2];
     for (int i = 0; i<2; i++)
@@ -437,7 +440,6 @@ namespace netgen
               {
                 curr++;
                 (*testout) << "edge nr " << curr << endl;
-
                 multithread.percent = 100 * curr / double (total);
                 if (multithread.terminate) return;
 
@@ -473,54 +475,62 @@ namespace netgen
                 NgArray <double> params;
 
                 DivideEdge (edge, mp, params, mesh, mparam);
- 
-                NgArray <int> pnums;
+
+                NgArray <PointIndex> pnums;
                 pnums.SetSize (mp.Size()+2);
 
                 if (!merge_solids)
                   {
-                    pnums[0] = geom.vmap.FindIndex (TopExp::FirstVertex (edge));
-                    pnums[pnums.Size()-1] = geom.vmap.FindIndex (TopExp::LastVertex (edge));
+                    pnums[0] = geom.vmap.FindIndex (TopExp::FirstVertex (edge)) + PointIndex::BASE-1;
+                    pnums[pnums.Size()-1] = geom.vmap.FindIndex (TopExp::LastVertex (edge)) + PointIndex::BASE-1;
                   }
                 else
                   {
                     Point<3> fp = occ2ng (BRep_Tool::Pnt (TopExp::FirstVertex (edge)));
                     Point<3> lp = occ2ng (BRep_Tool::Pnt (TopExp::LastVertex (edge)));
 
-                    pnums[0] = -1;
-                    pnums.Last() = -1;
-                    for (PointIndex pi = 1; pi < first_ep; pi++)
+                    pnums[0] = PointIndex::INVALID;
+                    pnums.Last() = PointIndex::INVALID;
+                    // for (PointIndex pi = 1; pi < first_ep; pi++)
+                    for (PointIndex pi : vertexrange)
                       {
                         if (Dist2 (mesh[pi], fp) < eps*eps) pnums[0] = pi;
                         if (Dist2 (mesh[pi], lp) < eps*eps) pnums.Last() = pi;
                       }
                   }
 
-
                 for (int i = 1; i <= mp.Size(); i++)
                   {
                     bool exists = 0;
                     tsearch.Start();
-                    int j;
+                    PointIndex j;
+                    /*
                     for (j = first_ep; j <= mesh.GetNP(); j++)
                       if ((mesh.Point(j)-Point<3>(mp[i-1])).Length() < eps)
                         {
                           exists = 1;
                           break;
                         }
+                    */
+                    for (j = first_ep; j < mesh.Points().End(); j++)
+                      if ((mesh.Point(j)-Point<3>(mp[i-1])).Length() < eps)
+                        {
+                          exists = 1;
+                          break;
+                        }
+                      
                     tsearch.Stop();
                     
                     if (exists)
                       pnums[i] = j;
                     else
                       {
-                        mesh.AddPoint (mp[i-1]);
+                        pnums[i] = mesh.AddPoint (mp[i-1]);
                         (*testout) << "add meshpoint " << mp[i-1] << endl;
-                        pnums[i] = mesh.GetNP();
+                        // pnums[i] = mesh.GetNP();
                       }
                   }
                 (*testout) << "NP = " << mesh.GetNP() << endl;
-
                 //(*testout) << pnums[pnums.Size()-1] << endl;
 
                 for (int i = 1; i <= mp.Size()+1; i++)
@@ -623,7 +633,7 @@ namespace netgen
 
     double starttime = GetTime();
 
-    NgArray<int> glob2loc(noldp);
+    NgArray<int, PointIndex::BASE> glob2loc(noldp);
 
     //int projecttype = PARAMETERSPACE;
 
@@ -696,12 +706,12 @@ namespace netgen
                   {
                     for (int j = 1; j <= 2; j++)
                       {
-                        int pi = (j == 1) ? seg[0] : seg[1];
-                        if (!glob2loc.Get(pi))
+                        PointIndex pi = (j == 1) ? seg[0] : seg[1];
+                        if (glob2loc[pi] == 0)
                           {
                             meshing.AddPoint (mesh.Point(pi), pi);
                             cntp++;
-                            glob2loc.Elem(pi) = cntp;
+                            glob2loc[pi] = cntp;
                           }
                       }
                   }
@@ -719,7 +729,7 @@ namespace netgen
                     gi1.u = seg.epgeominfo[1].u;
                     gi1.v = seg.epgeominfo[1].v;
 
-                    meshing.AddBoundaryElement (glob2loc.Get(seg[0]), glob2loc.Get(seg[1]), gi0, gi1);
+                    meshing.AddBoundaryElement (glob2loc[seg[0]], glob2loc[seg[1]], gi0, gi1);
                     //(*testout) << gi0.u << " " << gi0.v << endl;
                     //(*testout) << gi1.u << " " << gi1.v << endl;
                   }
@@ -770,7 +780,7 @@ namespace netgen
 
                         if (locpnum[j] == 0)
                           {
-                            int pi = (j == 0) ? seg[0] : seg[1];
+                            PointIndex pi = (j == 0) ? seg[0] : seg[1];
                             meshing.AddPoint (mesh.Point(pi), pi);
 
                             gis.SetSize (gis.Size()+1);