Support for Apple M1

basiclinalg/ngblas.cpp

@@ -3792,7 +3792,7 @@ namespace ngbla
     return ssum;
   }
 
-#elif defined __SSE__
+#elif defined NETGEN_ARCH_AMD64
 
   double MatKernelMaskedScalAB (size_t n,
 				double * pa, size_t da,
@@ -3853,21 +3853,20 @@ namespace ngbla
 				double * pb, size_t db,
 				const BitArray & ba)
   {
-    double sum = 0;
     double vhsum[8] = { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 };
     int i(0);
-    for ( ; i+8 < fa.Size(); i += 8)
+    for ( ; i+8 <= n; i += 8)
       {
 	for (int j = 0; j < 8; j++)
 	  {
-	    double hprod = fa(i+j)*fb(i+j);
+	    double hprod = pa[i+j]*pb[i+j];
 	    if (ba.Test(i+j))
 	      vhsum[j] += hprod;
 	  }
       }
-    for ( ; i < fa.Size(); i++)
+    for ( ; i < n; i++)
       if (ba.Test(i))
-	sum += fa(i)*fb(i);
+	vhsum[0] += pa[i]*pb[i];
     for (int j = 1; j < 8; j++)
       vhsum[0] += vhsum[j];
     return vhsum[0];

basiclinalg/ngblas.hpp

@@ -3,7 +3,7 @@
 
 // optimized matrix kernels
 
-#ifdef __clang__
+#if defined(__clang__ ) && defined(NETGEN_ARCH_AMD64)
 #define REGCALL __regcall
 #else
 #define REGCALL

comp/meshaccess.cpp

@@ -2497,13 +2497,13 @@ namespace ngcomp
   atomic<size_t> ProgressOutput :: cnt;
   thread_local size_t ProgressOutput :: thd_cnt = 0;
   // thread_local double ProgressOutput :: thd_prev_time = WallTime();
-  thread_local size_t ProgressOutput :: thd_prev_time = __rdtsc();
-  size_t tsc_wait = 0.05*2.7e9; // rough 
+  thread_local size_t ProgressOutput :: thd_prev_time = GetTimeCounter();
+  size_t tsc_wait = 0.05*(1.0/seconds_per_tick);
   void ProgressOutput :: Update ()
   {
     thd_cnt++;
     // double time = WallTime();
-    size_t time = __rdtsc();
+    size_t time = GetTimeCounter();
     // if (time > thd_prev_time+0.05)
     if (time > thd_prev_time+tsc_wait)
       {

fem/intrule.cpp

@@ -3277,8 +3277,14 @@ namespace ngfem
     nip = ir.GetNIP();
     this->size = (ir.Size()+SIMD<IntegrationPoint>::Size()-1) / SIMD<IntegrationPoint>::Size();
 
+#ifdef NETGEN_ARCH_AMD64
     this -> mem_to_delete = (SIMD<IntegrationPoint>*)
       _mm_malloc(this->size*sizeof(SIMD<IntegrationPoint>), SIMD<double>::Size()*sizeof(double));
+#else // NETGEN_ARCH_AMD64
+    this -> mem_to_delete = (SIMD<IntegrationPoint>*)
+      malloc(this->size*sizeof(SIMD<IntegrationPoint>));
+#endif // NETGEN_ARCH_AMD64
+
     this->data = this->mem_to_delete;
 
     dimension = ir.Dim();

fem/intrule.hpp

@@ -2283,7 +2283,13 @@ namespace ngfem
     SIMD_IntegrationRule (int nip, LocalHeap & lh);
     NGS_DLL_HEADER ~SIMD_IntegrationRule ()
     {
-      if (mem_to_delete) _mm_free(mem_to_delete);
+#ifdef NETGEN_ARCH_AMD64
+      if (mem_to_delete)
+          _mm_free(mem_to_delete);
+#else // NETGEN_ARCH_AMD64
+      if (mem_to_delete)
+          free(mem_to_delete);
+#endif // NETGEN_ARCH_AMD64
       mem_to_delete = nullptr;
     }
 

include/ngs_stdcpp_include.hpp

@@ -137,7 +137,6 @@
 #endif
 
 
-#include <immintrin.h>
 
 
 #ifndef __assume

linalg/sparsematrix_impl.hpp

@@ -21,6 +21,7 @@ namespace ngla
   void SparseMatrixTM<TM> ::
   PrefetchRow (int rownr) const
   {
+#ifdef NETGEN_ARCH_AMD64
 #ifdef __GNUC__
     size_t fi = firsti[rownr], fin = firsti[rownr+1];
     // int * pi = &colnr[fi], * pin = &colnr[fin];
@@ -38,6 +39,7 @@ namespace ngla
         vi += 64/sizeof(double);
       }
 #endif
+#endif // NETGEN_ARCH_AMD64
     ;
   }
 

ngstd/autodiffdiff.hpp

@@ -456,12 +456,12 @@ inline AutoDiffDiff<D, SCAL> sqrt (const AutoDiffDiff<D, SCAL> & x)
   AutoDiffDiff<D, SCAL> res;
   res.Value() = sqrt(x.Value());
   for (int j = 0; j < D; j++)
-    res.DValue(j) = IfZero(x.DValue(j),0.,0.5 / res.Value() * x.DValue(j));
+    res.DValue(j) = IfZero(x.DValue(j),SCAL{0.},0.5 / res.Value() * x.DValue(j));
 
 
   for (int i = 0; i < D; i++)
     for (int j = 0; j < D; j++)
-      res.DDValue(i,j) = IfZero(x.DDValue(i,j)+x.DValue(i) * x.DValue(j),0.,0.5/res.Value() * x.DDValue(i,j) - 0.25 / (x.Value()*res.Value()) * x.DValue(i) * x.DValue(j));
+      res.DDValue(i,j) = IfZero(x.DDValue(i,j)+x.DValue(i) * x.DValue(j),SCAL{0.},0.5/res.Value() * x.DDValue(i,j) - 0.25 / (x.Value()*res.Value()) * x.DValue(i) * x.DValue(j));
 
   return res;
 }