From ccde47d2c07e5b5c9a6b5fb38fc38e407d1938d7 Mon Sep 17 00:00:00 2001
From: Christopher Lackner <christopher.lackner@gmx.at>
Date: Sun, 16 Oct 2016 09:45:16 +0200
Subject: [PATCH] added multielementtransformation

---
 libsrc/interface/nginterface_v2.cpp | 63 +++++++++++++++++++++++++++++
 libsrc/meshing/curvedelems.cpp      |  1 +
 2 files changed, 64 insertions(+)

diff --git a/libsrc/interface/nginterface_v2.cpp b/libsrc/interface/nginterface_v2.cpp
index d3e8fc83..d7548127 100644
--- a/libsrc/interface/nginterface_v2.cpp
+++ b/libsrc/interface/nginterface_v2.cpp
@@ -535,6 +535,17 @@ namespace netgen
     if (dxdxi) dxdxi[0] = dx(0);
   }
 
+  template <> DLL_HEADER void Ngx_Mesh ::
+  ElementTransformation<0,2> (int elnr,
+			      const double *xi,
+			      double * x,
+			      double * dxdxi) const
+  {
+    PointIndex pnum = mesh->pointelements[elnr].pnum;
+    if (x)
+      for (int i = 0; i< 2; i++) x[i] = (*mesh)[pnum](i);
+  }
+
 
   template <> DLL_HEADER void Ngx_Mesh :: 
   ElementTransformation<0,1> (int elnr, 
@@ -581,6 +592,15 @@ namespace netgen
     mesh->GetCurvedElements().CalcMultiPointSurfaceTransformation<3> (elnr, npts, xi, sxi, x, sx, dxdxi, sdxdxi);
   }
 
+  template <> DLL_HEADER void Ngx_Mesh ::
+  MultiElementTransformation<1,3> (int elnr, int npts,
+				   const double * xi, size_t sxi,
+                                   double * x, size_t sx,
+                                   double * dxdxi, size_t sdxdxi) const
+  {
+    mesh->GetCurvedElements().CalcMultiPointSegmentTransformation<3> (elnr, npts, xi, sxi, x, sx, dxdxi, sdxdxi);
+  }
+
   template <> DLL_HEADER void Ngx_Mesh :: 
   MultiElementTransformation<1,2> (int elnr, int npts,
                                    const double * xi, size_t sxi,
@@ -600,6 +620,16 @@ namespace netgen
       ElementTransformation<1,1> (elnr, xi + i*sxi, x+i*sx, dxdxi+i*sdxdxi);
   }
 
+  template <> DLL_HEADER void Ngx_Mesh ::
+  MultiElementTransformation<0,2> (int elnr, int npts,
+                                   const double * xi, size_t sxi,
+                                   double * x, size_t sx,
+                                   double * dxdxi, size_t sdxdxi) const
+  {
+    for (int i = 0; i < npts; i++)
+      ElementTransformation<0,2> (elnr, xi + i*sxi, x+i*sx, dxdxi+i*sdxdxi);
+  }
+
 
   template <> DLL_HEADER void Ngx_Mesh :: 
   MultiElementTransformation<0,1> (int elnr, int npts,
@@ -698,6 +728,15 @@ namespace netgen
     */
   }
 
+  template<> DLL_HEADER void Ngx_Mesh ::
+  MultiElementTransformation<0,2> (int elnr, int npts,
+				   const __m256d *xi, size_t sxi,
+				   __m256d * x, size_t sx,
+				   __m256d * dxdxi, size_t sdxdxi) const
+  {
+    cout << "MultiElementtransformation<0,2> simd not implemented" << endl;
+  }
+  
   template<> DLL_HEADER void Ngx_Mesh :: 
   MultiElementTransformation<0,1> (int elnr, int npts,
                                    const __m256d * xi, size_t sxi,
@@ -707,6 +746,30 @@ namespace netgen
     cout << "multi-eltrafo simd called, 0,1,simd" << endl;
   }
 
+  template<> DLL_HEADER void Ngx_Mesh :: 
+  MultiElementTransformation<1,3> (int elnr, int npts,
+                                   const __m256d * xi, size_t sxi,
+                                   __m256d * x, size_t sx,
+                                   __m256d * dxdxi, size_t sdxdxi) const
+  {
+    double hxi[4][1];
+    double hx[4][3];
+    double hdxdxi[4][3];
+    for (int j = 0; j<4;j++)
+      hxi[j][0] = ((double*)&(xi[0]))[j];
+    MultiElementTransformation<1,3> (elnr, 4, &hxi[0][0], 1, &hx[0][0], 2, &hdxdxi[0][0],4);
+    for(int j=0; j<4; j++)
+      for(int k=0; k<3; k++)
+	((double*)&(x[k]))[j] = hx[j][k];
+    for(int j=0; j< 4; j++)
+      for (int k = 0; k<3; k++)
+	((double*) & (dxdxi[k]))[j] = hdxdxi[j][k];
+
+    xi += sxi;
+    x += sx;
+    dxdxi += sdxdxi;
+  }
+  
   template<> DLL_HEADER void Ngx_Mesh :: 
   MultiElementTransformation<1,2> (int elnr, int npts,
                                    const __m256d * xi, size_t sxi,
diff --git a/libsrc/meshing/curvedelems.cpp b/libsrc/meshing/curvedelems.cpp
index 6f05ed86..14b0d4df 100644
--- a/libsrc/meshing/curvedelems.cpp
+++ b/libsrc/meshing/curvedelems.cpp
@@ -3686,6 +3686,7 @@ namespace netgen
       }
   }
 
+
   template void CurvedElements :: 
   CalcMultiPointSegmentTransformation<2> (SegmentIndex elnr, int npts,
 					  const double * xi, size_t sxi,