From 3bca66892c26687b68c37187c0ca5be878f5145c Mon Sep 17 00:00:00 2001
From: Thomas Witkowski <thomas.witkowski@gmx.de>
Date: Thu, 27 May 2010 12:57:16 +0000
Subject: [PATCH] Removed allowFirstRef from ProblemVec.

---
 AMDiS/src/BasisFunction.h |   3 +-
 AMDiS/src/ElInfo.cc       |   4 +-
 AMDiS/src/ElInfo2d.cc     | 130 +++++++++++++++++++++++++++++++-------
 AMDiS/src/ElInfo2d.h      |  18 ++++--
 AMDiS/src/ElInfo3d.cc     |   4 +-
 AMDiS/src/ProblemVec.cc   |  29 ++-------
 AMDiS/src/ProblemVec.h    |  10 ---
 7 files changed, 130 insertions(+), 68 deletions(-)

diff --git a/AMDiS/src/BasisFunction.h b/AMDiS/src/BasisFunction.h
index 6203bb6c..d8038de6 100644
--- a/AMDiS/src/BasisFunction.h
+++ b/AMDiS/src/BasisFunction.h
@@ -197,8 +197,7 @@ namespace AMDiS {
     /// WorldVector<double> valued interpol function.
     virtual const WorldVector<double>* interpol(const ElInfo *el_info, int no, 
 						const int *b_no,
-						AbstractFunction<WorldVector<double>,
-						             WorldVector<double> > *f, 
+						AbstractFunction<WorldVector<double>, WorldVector<double> > *f, 
 						WorldVector<double> *vec) = 0;
 
     /// Returns the i-th local basis function
diff --git a/AMDiS/src/ElInfo.cc b/AMDiS/src/ElInfo.cc
index 748cb1e6..567f7f74 100644
--- a/AMDiS/src/ElInfo.cc
+++ b/AMDiS/src/ElInfo.cc
@@ -14,9 +14,9 @@
 
 namespace AMDiS {
 
-  std::vector<std::map<unsigned long, mtl::dense2D<double> > > ElInfo::subElemMatrices(4);
+  std::vector<std::map<unsigned long, mtl::dense2D<double> > > ElInfo::subElemMatrices(5);
 
-  std::vector<std::map<unsigned long, mtl::dense2D<double> > > ElInfo::subElemGradMatrices(4);
+  std::vector<std::map<unsigned long, mtl::dense2D<double> > > ElInfo::subElemGradMatrices(5);
 
   ElInfo::ElInfo(Mesh *aMesh) 
     : mesh(aMesh),
diff --git a/AMDiS/src/ElInfo2d.cc b/AMDiS/src/ElInfo2d.cc
index 1b89ab5a..60c9f204 100644
--- a/AMDiS/src/ElInfo2d.cc
+++ b/AMDiS/src/ElInfo2d.cc
@@ -14,12 +14,6 @@
 
 namespace AMDiS {
 
-  double ElInfo2d::mat_d1_val[3][3] = {{1.0, 0.0, 0.0}, 
-				       {0.0, 1.0, 0.0}, 
-				       {0.0, 0.0, 1.0}};
-  mtl::dense2D<double> ElInfo2d::mat_d1(mat_d1_val);
-
-  
   double ElInfo2d::mat_d1_left_val[3][3] = {{0.0, 1.0, 0.5}, 
   					    {0.0, 0.0, 0.5},
   					    {1.0, 0.0, 0.0}}; 
@@ -33,15 +27,6 @@ namespace AMDiS {
 
 
 
-
-  double ElInfo2d::mat_d2_val[6][6] = {{1.0, 0.0, 0.0, 0.0, 0.0, 0.0}, 
-				       {0.0, 1.0, 0.0, 0.0, 0.0, 0.0}, 
-				       {0.0, 0.0, 1.0, 0.0, 0.0, 0.0},
-				       {0.0, 0.0, 0.0, 1.0, 0.0, 0.0},
-				       {0.0, 0.0, 0.0, 0.0, 1.0, 0.0},
-				       {0.0, 0.0, 0.0, 0.0, 0.0, 1.0}};
-  mtl::dense2D<double> ElInfo2d::mat_d2(mat_d2_val);
-
   double ElInfo2d::mat_d2_left_val[6][6] = {{0.0, 1.0, 0.0, 0.375, -0.125, 0.0}, 
 					    {0.0, 0.0, 0.0, -0.125, -0.125, 0.0}, 
 					    {1.0, 0.0, 0.0, 0.0, 0.0, 0.0},
@@ -59,6 +44,68 @@ namespace AMDiS {
   mtl::dense2D<double> ElInfo2d::mat_d2_right(mat_d2_right_val);
 
 
+
+  double ElInfo2d::mat_d3_left_val[10][10] = {{0.0, 1.0, 0.0, 0.3125, 0.0, 0.0, 0.0625, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0625, 0.0, 0.0, 0.0625, 0.0, 0.0, 0.0625},
+					      {1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.5, 1.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.375},
+					      {0.0, 0.0, 0.5625, 0.9375, 1.0, 0.0, 0.0, 0.0, 0.0, 0.1875},
+					      {0.0, 0.0, 0.5625, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.5, 0.0, 0.0, 0.75}};
+  mtl::dense2D<double> ElInfo2d::mat_d3_left(mat_d3_left_val);
+
+  double ElInfo2d::mat_d3_right_val[10][10] = {{0.0, 0.0, 0.0, 0.0625, 0.0, 0.0, 0.0625, 0.0, 0.0, 0.0625},
+					       {1.0, 0.0, 0.0, 0.0625, 0.0, 0.0, 0.3125, 0.0, 0.0, 0.0},
+					       {0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.375},
+					       {0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.5625, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.5625, 0.0, 0.0, 1.0, 0.9375, 0.0, 0.0, 0.1875},
+					       {0.0, 0.0, 0.0, 0.5, 1.0, 0.0, 0.0, 0.0, 0.0, 0.75}};
+  mtl::dense2D<double> ElInfo2d::mat_d3_right(mat_d3_right_val);
+
+
+
+  double ElInfo2d::mat_d4_left_val[15][15] = {{0.0, 1.0, 0.0, 0.273437, 0.0, 0.0, 0.023437, 0.0, 0.0, 0.0, 0.0, 0.0, 0.023437, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.023437, 0.023437, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.1875, 0.0, 0.0, 0.0, 0.125, 0.0625, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.5, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.375, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.1875, 0.0, 0.0, 1.0, 0.0, 0.3125, 0.0},
+					      {0.0, 0.0, 0.0, 1.093750, 1.0, 0.468750, 0.0, 0.0, 0.031250, 0.0, 0.0, 0.0, 0.0, 0.156250, 0.0},
+					      {0.0, 0.0, 1.0, 0.0, 0.0, 0.703125, 0.140625, 0.0, 0.015625, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.218750, 0.0, 0.0, 0.0, 0.0, 0.031250, 0.0, 0.0, 0.0, 0.093750, 0.156250, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.5625, 0.0, 0.0, 0.0, 0.0, 0.0, 0.375, 0.9375, 1.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.5625, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					      {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.75, 0.0, 0.0, 0.0, 0.75, 0.0, 0.0}};  
+  mtl::dense2D<double> ElInfo2d::mat_d4_left(mat_d4_left_val);
+
+  double ElInfo2d::mat_d4_right_val[15][15] = {{0.0, 0.0, 0.0, 0.0, 0.0, 0.023437, 0.023437, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {1.0, 0.0, 0.0, 0.0, 0.0, 0.023437, 0.0, 0.0, 0.273437, 0.0, 0.0, 0.0, 0.0, 0.023437, 0.0},
+					       {0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.1875, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.3125, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.375, 0.0},
+					       {0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.1875, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0625, 0.125, 0.0},
+					       {0.0, 0.0, 0.0, 0.031250, 0.0, 0.0, 0.0, 0.0, 0.218750, 0.0, 0.0, 0.0, 0.156250, 0.093750, 0.0},
+					       {0.0, 0.0, 1.0, 0.015625, 0.0, 0.140625, 0.703125, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.031250, 0.0, 0.0, 0.468750, 1.0, 1.093750, 0.0, 0.0, 0.0, 0.156250, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.0, 0.0, 0.5625, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
+					       {0.0, 0.0, 0.0, 0.0, 0.0, 0.5625, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.9375, 0.375, 1.0},
+					       {0.0, 0.0, 0.0, 0.75, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.75, 0.0}};  
+  mtl::dense2D<double> ElInfo2d::mat_d4_right(mat_d4_right_val);
+
+
   ElInfo2d::ElInfo2d(Mesh *aMesh) 
     : ElInfo(aMesh) 
   {
@@ -723,9 +770,9 @@ namespace AMDiS {
       switch (degree) {
       case 1:
 	{
-	  dense2D<double> mat(mat_d1);
-	  dense2D<double> tmpMat(num_rows(mat), num_rows(mat));
-	  
+	  dense2D<double> mat(3, 3), tmpMat(3, 3);
+	  mat = 1;
+
 	  for (int i = 0; i < refinementPathLength; i++) {
 	    if (refinementPath & (1 << i)) {
 	      tmpMat = mat_d1_right * mat;
@@ -741,9 +788,9 @@ namespace AMDiS {
 	break;
       case 2:
 	{
-	  dense2D<double> mat(mat_d2);
-	  dense2D<double> tmpMat(num_rows(mat), num_rows(mat));
-	  
+	  dense2D<double> mat(6, 6), tmpMat(6, 6);
+	  mat = 1;
+
 	  for (int i = 0; i < refinementPathLength; i++) {
 	    if (refinementPath & (1 << i)) {
 	      tmpMat = mat_d2_right * mat;
@@ -757,6 +804,43 @@ namespace AMDiS {
 	  subElemMatrices[2][refinementPath] = mat;  
 	}
 	break;
+      case 3:
+	{
+	  dense2D<double> mat(10, 10), tmpMat(10, 10);
+	  mat = 1;
+
+	  for (int i = 0; i < refinementPathLength; i++) {
+	    if (refinementPath & (1 << i)) {
+	      tmpMat = mat_d3_right * mat;
+	      mat = tmpMat;
+	    } else  {
+	      tmpMat = mat_d3_left * mat;
+	      mat = tmpMat;
+	    }
+	  }
+
+	  subElemMatrices[3][refinementPath] = mat;  
+	}
+	break;
+      case 4:
+	{
+	  dense2D<double> mat(15, 15), tmpMat(15, 15);
+	  mat = 1;
+
+	  for (int i = 0; i < refinementPathLength; i++) {
+	    if (refinementPath & (1 << i)) {
+	      tmpMat = mat_d4_right * mat;
+	      mat = tmpMat;
+	    } else  {
+	      tmpMat = mat_d4_left * mat;
+	      mat = tmpMat;
+	    }
+	  }
+
+	  subElemMatrices[4][refinementPath] = mat;  
+	}
+	break;
+
       default:
 	ERROR_EXIT("Not supported for basis function degree: %d\n", degree);
       }
@@ -775,8 +859,8 @@ namespace AMDiS {
     using namespace mtl;
 
     if (subElemGradMatrices[degree].count(refinementPath) == 0) {
-      dense2D<double> mat(mat_d1);
-      dense2D<double> tmpMat(num_rows(mat), num_rows(mat));
+      dense2D<double> mat(3, 3), tmpMat(3, 3);
+      mat = 1;
 
       double test_left[3][3] = {{0.0, 0.0, 0.5},
 				{-0.5, -0.5, 0.0},
diff --git a/AMDiS/src/ElInfo2d.h b/AMDiS/src/ElInfo2d.h
index 54f50651..09cb0f88 100644
--- a/AMDiS/src/ElInfo2d.h
+++ b/AMDiS/src/ElInfo2d.h
@@ -66,23 +66,29 @@ namespace AMDiS {
     /// Temp vectors for function \ref calcGrdLambda.
     WorldVector<double> *e1, *e2, *normal;
 
-    static double mat_d1_val[3][3];
-    static mtl::dense2D<double> mat_d1;
-
     static double mat_d1_left_val[3][3];
     static mtl::dense2D<double> mat_d1_left;
 
     static double mat_d1_right_val[3][3];
     static mtl::dense2D<double> mat_d1_right;
 
-    static double mat_d2_val[6][6];
-    static mtl::dense2D<double> mat_d2;
-
     static double mat_d2_left_val[6][6];
     static mtl::dense2D<double> mat_d2_left;
 
     static double mat_d2_right_val[6][6];
     static mtl::dense2D<double> mat_d2_right;
+
+    static double mat_d3_left_val[10][10];
+    static mtl::dense2D<double> mat_d3_left;
+
+    static double mat_d3_right_val[10][10];
+    static mtl::dense2D<double> mat_d3_right;
+
+    static double mat_d4_left_val[15][15];
+    static mtl::dense2D<double> mat_d4_left;
+
+    static double mat_d4_right_val[15][15];
+    static mtl::dense2D<double> mat_d4_right;
   };
 
 }
diff --git a/AMDiS/src/ElInfo3d.cc b/AMDiS/src/ElInfo3d.cc
index ed0c5865..bdb97fd8 100644
--- a/AMDiS/src/ElInfo3d.cc
+++ b/AMDiS/src/ElInfo3d.cc
@@ -759,8 +759,8 @@ namespace AMDiS {
 	    if (refinementPath & (1 << i)) {
 	      if ((level + i) % 3 == 0)
 		tmpMat = mat_d4_l0_right * mat;
-// 	      else
-// 		tmpMat = mat_d4_l12_right * mat;
+ 	      else
+ 		tmpMat = mat_d4_l12_right * mat;
 
 	      mat = tmpMat;
 	    } else  {
diff --git a/AMDiS/src/ProblemVec.cc b/AMDiS/src/ProblemVec.cc
index f2c28870..dd8f3598 100644
--- a/AMDiS/src/ProblemVec.cc
+++ b/AMDiS/src/ProblemVec.cc
@@ -560,10 +560,6 @@ namespace AMDiS {
   {
     FUNCNAME("ProblemVec::markElements()");
 
-    // to enforce albert-like behavior: refinement even if space tolerance
-    // here is reached already because of time adaption
-    allowFirstRefinement();
-
     TEST_EXIT_DBG(static_cast<unsigned int>(nComponents) == marker.size())
       ("Wrong number of markers!\n");
 
@@ -611,18 +607,11 @@ namespace AMDiS {
   {
     FUNCNAME("ProblemVec::oneIteration()");
 
-    if (allowFirstRef) {
-      for (int i = 0; i < nComponents; i++)
-	adaptInfo->allowRefinement(true, i);
-      
-      allowFirstRef = false;
-    } else {
-      for (int i = 0; i < nComponents; i++)
- 	if (adaptInfo->spaceToleranceReached(i))
- 	  adaptInfo->allowRefinement(false, i);
- 	else
- 	  adaptInfo->allowRefinement(true, i);	
-    }
+    for (int i = 0; i < nComponents; i++)
+      if (adaptInfo->spaceToleranceReached(i))
+	adaptInfo->allowRefinement(false, i);
+      else
+	adaptInfo->allowRefinement(true, i);	
 
     return StandardProblemIteration::oneIteration(adaptInfo, toDo);
   }
@@ -1532,17 +1521,13 @@ namespace AMDiS {
     }
 
     ElementFileWriter::writeFile(vec, this->getFESpace(comp), name);
-//     ElementFileWriter fw(name, this->getFESpace(comp), vec);
-//     fw.writeFiles(adaptInfo, true);    
   }
 
 
   void ProblemVec::serialize(std::ostream &out) 
   {
     FUNCNAME("ProblemVec::serialize()");
-
-    SerUtil::serialize(out, allowFirstRef);
-    
+   
     for (int i = 0; i < static_cast<int>(meshes.size()); i++)
       meshes[i]->serialize(out);
 
@@ -1554,8 +1539,6 @@ namespace AMDiS {
   {
     FUNCNAME("ProblemVec::deserialize()");
 
-    SerUtil::deserialize(in, allowFirstRef);
-
     for (int i = 0; i < static_cast<int>(meshes.size()); i++)
       meshes[i]->deserialize(in);
 
diff --git a/AMDiS/src/ProblemVec.h b/AMDiS/src/ProblemVec.h
index 49968702..9daf183b 100644
--- a/AMDiS/src/ProblemVec.h
+++ b/AMDiS/src/ProblemVec.h
@@ -69,7 +69,6 @@ namespace AMDiS {
 	systemMatrix(NULL),
 	useGetBound(true),
 	info(10),
-	allowFirstRef(false),
 	computeExactError(false),
 	boundaryConditionSet(false),
 	writeAsmInfo(false)
@@ -259,12 +258,6 @@ namespace AMDiS {
     /// Adds periodic boundary conditions.
     virtual void addPeriodicBC(BoundaryType type, int row, int col);
 
-    /// Implementation of ProblemStatBase::allowFirstRefinement().
-    inline void allowFirstRefinement() 
-    {
-      allowFirstRef = true;
-    }
-
     /** \brief
      * This function assembles a DOFMatrix and a DOFVector for the case,
      * the meshes from row and col FE-space are equal.
@@ -594,9 +587,6 @@ namespace AMDiS {
     /// Info level.
     int info;
 
-    /// Allows one refinement although the adapt tolerance is reached already.
-    bool allowFirstRef;
-
     /** \brief
      * This vectors stores pointers to functions defining the exact solution of
      * the problem. This may be used to compute the real error of the computed
-- 
GitLab