From 537a10945ad71dad90450ad7a75ffd2df3f5d71f Mon Sep 17 00:00:00 2001
From: Joachim Schoeberl <joachim.schoeberl@tuwien.ac.at>
Date: Thu, 8 Jan 2015 21:18:22 +0000
Subject: [PATCH] boundary layer

---
 libsrc/meshing/boundarylayer.cpp | 23 ++++++++++++++++++++++-
 libsrc/meshing/boundarylayer.hpp |  3 ++-
 2 files changed, 24 insertions(+), 2 deletions(-)

diff --git a/libsrc/meshing/boundarylayer.cpp b/libsrc/meshing/boundarylayer.cpp
index 48dbae00..1ca3589e 100644
--- a/libsrc/meshing/boundarylayer.cpp
+++ b/libsrc/meshing/boundarylayer.cpp
@@ -195,6 +195,7 @@ namespace netgen
          // surface elements being added during the process
          int np = mesh.GetNP();
          int nse = mesh.GetNSE();
+         int ne = mesh.GetNE();
 
          // Safety measure to ensure no issues with mesh 
          // consistency
@@ -509,7 +510,7 @@ namespace netgen
                  Element el(types[classify]);
                  for (int i = 0; i < 6; i++)
                    el[i] = nums[vertices[classify][i]];
-                 el.SetIndex(blp.matnr);
+                 el.SetIndex(blp.new_matnr);
                  // cout << "el = " << el << endl;
                  if (classify != 0)
                    mesh.AddVolumeElement(el);
@@ -537,6 +538,26 @@ namespace netgen
                   }
               }
            }
+         for (int i = 1; i <= ne; i++)
+           {
+             Element & el = mesh.VolumeElement(i);
+             if(el.GetIndex() != blp.bulk_matnr)
+              {
+                for (int j = 1; j <= el.GetNP(); j++)
+                  {
+                    // Check (Doublecheck) if the corresponding point has a 
+                    // copy available for remapping
+                    if (mapto.Get(el.PNum(j)))
+                      {
+                        // Map the surface elements to the new points
+                        el.PNum(j) = mapto.Get(el.PNum(j));
+                      }
+                  }
+              }
+           }
+
+
+
          
          // Lock all the prism points so that the rest of the mesh can be 
          // optimised without invalidating the entire mesh
diff --git a/libsrc/meshing/boundarylayer.hpp b/libsrc/meshing/boundarylayer.hpp
index 73f3de5c..cbad16af 100644
--- a/libsrc/meshing/boundarylayer.hpp
+++ b/libsrc/meshing/boundarylayer.hpp
@@ -14,7 +14,8 @@ public:
   // parameters by Philippose ..
   Array<int> surfid;
   int prismlayers = 1;
-  int matnr = 1;
+  int bulk_matnr = 1;
+  int new_matnr = 1;
   double hfirst = 0.01;
   double growthfactor = 1;
   bool optimize = true;