Introduced stack for subElemMatrices.

AMDiS/src/ElInfo.cc

@@ -14,6 +14,8 @@
 
 namespace AMDiS {
 
+  std::map<unsigned long, mtl::dense2D<double> > ElInfo::subElemMatrices;
+
   ElInfo::ElInfo(Mesh *aMesh) 
     : mesh(aMesh),
       element(NULL),

AMDiS/src/ElInfo.h

@@ -544,6 +544,8 @@ namespace AMDiS {
     mtl::dense2D<double> subElemCoordsMat_so;
 
   public:
+    static std::map<unsigned long, mtl::dense2D<double> > subElemMatrices;
+
     /** \brief 
      * child_vertex[el_type][child][i] = father's local vertex index of new 
      * vertex i. 4 stands for the newly generated vertex .

AMDiS/src/ElInfo2d.cc

@@ -697,18 +697,25 @@ namespace AMDiS {
     switch (degree) {
     case 1:
       {
-      mat = mat_d1;
-      dense2D<double> tmpMat(num_rows(mat), num_rows(mat));
-
-      for (int i = 0; i < refinementPathLength; i++) {
-	if (refinementPath & (1 << i)) {
-	  tmpMat = mat_d1_right * mat;
-	  mat = tmpMat;
-	} else  {
-	  tmpMat = mat_d1_left * mat;
-	  mat = tmpMat;
+	if (subElemMatrices.count(refinementPath) > 0) {
+	  mat = subElemMatrices[refinementPath];
+	  return;
 	}
-      }
+
+	mat = mat_d1;
+	dense2D<double> tmpMat(num_rows(mat), num_rows(mat));
+	
+	for (int i = 0; i < refinementPathLength; i++) {
+	  if (refinementPath & (1 << i)) {
+	    tmpMat = mat_d1_right * mat;
+	    mat = tmpMat;
+	  } else  {
+	    tmpMat = mat_d1_left * mat;
+	    mat = tmpMat;
+	  }
+	}
+
+	subElemMatrices[refinementPath] = mat;
       }
       break;
     default:

AMDiS/src/ProblemVec.cc

@@ -627,6 +627,8 @@ namespace AMDiS {
 
     //    printOpenmpTraverseInfo(this, true);
 
+    //    std::cout << "ElInfo = " << ElInfo::subElemMatrices.size() << std::endl;
+
     clock_t first = clock();
 #ifdef _OPENMP
     double wtime = omp_get_wtime();