From 6caf6e9c3716d61c84aa21cc04f119c7e85e0a9e Mon Sep 17 00:00:00 2001
From: Thomas Witkowski <thomas.witkowski@gmx.de>
Date: Tue, 12 May 2009 07:56:55 +0000
Subject: [PATCH] Simplification of assembling on different meshes.

---
 AMDiS/src/AMDiS_fwd.h             |   5 ++
 AMDiS/src/Assembler.cc            |  24 +++---
 AMDiS/src/FirstOrderAssembler.cc  | 103 -------------------------
 AMDiS/src/FirstOrderAssembler.h   |  59 +--------------
 AMDiS/src/SecondOrderAssembler.cc | 120 ++---------------------------
 AMDiS/src/SecondOrderAssembler.h  |  27 -------
 AMDiS/src/SubAssembler.h          |   6 --
 AMDiS/src/ZeroOrderAssembler.cc   | 122 ------------------------------
 AMDiS/src/ZeroOrderAssembler.h    |  27 +------
 9 files changed, 26 insertions(+), 467 deletions(-)

diff --git a/AMDiS/src/AMDiS_fwd.h b/AMDiS/src/AMDiS_fwd.h
index 70b3400e..2e1440a0 100644
--- a/AMDiS/src/AMDiS_fwd.h
+++ b/AMDiS/src/AMDiS_fwd.h
@@ -69,6 +69,11 @@ namespace AMDiS {
   class Projection;
   class PreconditionerScal;
   class Quadrature; 
+  class Q00PsiPhi;
+  class Q0Psi;
+  class Q10PsiPhi;
+  class Q01PsiPhi;
+  class Q1Psi;
   class RCNeighbourList;
   class RefinementManager;
   class RobinBC;
diff --git a/AMDiS/src/Assembler.cc b/AMDiS/src/Assembler.cc
index f8b164f9..50a4f276 100644
--- a/AMDiS/src/Assembler.cc
+++ b/AMDiS/src/Assembler.cc
@@ -109,19 +109,23 @@ namespace AMDiS {
   
     ElementMatrix& mat = rememberElMat ? elementMatrix : userMat;
 
+    TEST_EXIT(zeroOrderAssembler)("Not yet implemented for not zero order assembler!\n");
+
     if (secondOrderAssembler)
-      secondOrderAssembler->calculateElementMatrix(rowElInfo, colElInfo, 
-						   smallElInfo, largeElInfo, mat);   
+      secondOrderAssembler->calculateElementMatrix(smallElInfo, mat);
     if (firstOrderAssemblerGrdPsi)
-      firstOrderAssemblerGrdPsi->calculateElementMatrix(rowElInfo, colElInfo, 
-							smallElInfo, largeElInfo, mat);    
+      firstOrderAssemblerGrdPsi->calculateElementMatrix(smallElInfo, mat);
     if (firstOrderAssemblerGrdPhi)
-      firstOrderAssemblerGrdPhi->calculateElementMatrix(rowElInfo, colElInfo, 
-							smallElInfo, largeElInfo, mat);    
-    if (zeroOrderAssembler) 
-      zeroOrderAssembler->calculateElementMatrix(rowElInfo, colElInfo, 
-						 smallElInfo, largeElInfo, mat);
-  
+      firstOrderAssemblerGrdPhi->calculateElementMatrix(smallElInfo, mat);
+    if (zeroOrderAssembler)
+      zeroOrderAssembler->calculateElementMatrix(smallElInfo, mat);
+
+    const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
+    ElementMatrix tmpMat(nRow, nRow);
+    tmpMat = m * mat;
+    mat = tmpMat;
+
+
     if (rememberElMat)
       userMat += factor * elementMatrix;
   }
diff --git a/AMDiS/src/FirstOrderAssembler.cc b/AMDiS/src/FirstOrderAssembler.cc
index 03bb3c2c..7df4bcb2 100644
--- a/AMDiS/src/FirstOrderAssembler.cc
+++ b/AMDiS/src/FirstOrderAssembler.cc
@@ -147,59 +147,6 @@ namespace AMDiS {
     }
   }
 
-  void Stand10::calculateElementMatrix(const ElInfo *rowElInfo, 
-				       const ElInfo *colElInfo,
-				       const ElInfo *smallElInfo,
-				       const ElInfo *largeElInfo,
-				       ElementMatrix& mat)
-  {
-    FUNCNAME("Stand10::calculateElementMatrix()");
-
-    TEST_EXIT((nRow <= 3) && (nCol <= 3))("not yet!\n");
-
-    DimVec<double> grdPsi(dim, DEFAULT_VALUE, 0.0);
-
-    const BasisFunction *psi = owner->getRowFESpace()->getBasisFcts();
-    const BasisFunction *phi = owner->getColFESpace()->getBasisFcts();
-    const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
-
-    int nPoints = quadrature->getNumPoints();
-    int myRank = omp_get_thread_num();
-    VectorOfFixVecs<DimVec<double> > &Lb = tmpLb[myRank];
-    Lb.resize(nPoints);
-
-    for (int iq = 0; iq < nPoints; iq++)
-      Lb[iq].set(0.0);
-
-    for (int i = 0; i < static_cast<int>(terms[myRank].size()); i++) {
-      (static_cast<FirstOrderTerm*>((terms[myRank][i])))->
-	getLb(smallElInfo, nPoints, Lb);
-    }
-
-    for (int iq = 0; iq < nPoints; iq++)
-      Lb[iq] *= smallElInfo->getDet();
-
-    DimVec<double> val1(dim, DEFAULT_VALUE, 0.0);
-
-    for (int i = 0; i < nRow; i++) {
-      for (int j = 0; j < nCol; j++) {	
-	for (int iq = 0; iq < nPoints; iq++) {
-	  
-	  val1.set(0.0);
-
-	  for (int k = 0; k < nCol; k++) {
-	    (*(psi->getGrdPhi(k)))(quadrature->getLambda(iq), grdPsi);
-	    grdPsi *= m[i][k];
-	    val1 += grdPsi;
-	  }
-
-	  mat[i][j] += quadrature->getWeight(iq) * (Lb[iq] * val1) *
-	    (*(phi->getPhi(j)))(quadrature->getLambda(iq));
-	}
-      }
-    }
-  }
-
   void Stand10::calculateElementVector(const ElInfo *elInfo, ElementVector& vec)
   {
     DimVec<double> grdPsi(dim, DEFAULT_VALUE, 0.0);
@@ -365,7 +312,6 @@ namespace AMDiS {
     }
   }
 
-
   Stand01::Stand01(Operator *op, Assembler *assembler, Quadrature *quad) 
     : FirstOrderAssembler(op, assembler, quad, false, GRD_PHI)
   {}
@@ -404,55 +350,6 @@ namespace AMDiS {
     } 
   }
 
-  void Stand01::calculateElementMatrix(const ElInfo *rowElInfo, 
-				       const ElInfo *colElInfo,
-				       const ElInfo *smallElInfo,
-				       const ElInfo *largeElInfo,
-				       ElementMatrix& mat)
-  {
-    FUNCNAME("Stand01::calculateElementMatrix()");
-
-    TEST_EXIT((nRow <= 3) && (nCol <= 3))("not yet!\n");
-
-    DimVec<double> grdPhi(dim, DEFAULT_VALUE, 0.0);
-
-    const BasisFunction *psi = owner->getRowFESpace()->getBasisFcts();
-    const BasisFunction *phi = owner->getColFESpace()->getBasisFcts();
-    const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
-
-    int nPoints = quadrature->getNumPoints();
-    int myRank = omp_get_thread_num();
-    VectorOfFixVecs<DimVec<double> > &Lb = tmpLb[myRank];
-    Lb.resize(nPoints);
-
-    for (int iq = 0; iq < nPoints; iq++)
-      Lb[iq].set(0.0);
-
-    for (int i = 0; i < static_cast<int>(terms[myRank].size()); i++)
-      (static_cast<FirstOrderTerm*>((terms[myRank][i])))->
-	getLb(smallElInfo, nPoints, Lb);
-
-    for (int iq = 0; iq < nPoints; iq++)
-      Lb[iq] *= smallElInfo->getDet();
-
-    for (int i = 0; i < nRow; i++) {
-      for (int j = 0; j < nCol; j++) {	
-	for (int iq = 0; iq < nPoints; iq++) {
-
-	  double val0 = 0.0;
-
-	  for (int k = 0; k < nCol; k++)
-	    val0 = m[i][k] * (*(psi->getPhi(k)))(quadrature->getLambda(iq));
-
-	  (*(phi->getGrdPhi(j)))(quadrature->getLambda(iq), grdPhi);
-
-	  mat[i][j] += quadrature->getWeight(iq) * ((Lb[iq] * val0) * grdPhi);
-
-	}
-      }
-    }
-  }
-
   Quad01::Quad01(Operator *op, Assembler *assembler, Quadrature *quad) 
     : FirstOrderAssembler(op, assembler, quad, true, GRD_PHI)
   {}
diff --git a/AMDiS/src/FirstOrderAssembler.h b/AMDiS/src/FirstOrderAssembler.h
index 4b1cac19..5f32e6c1 100644
--- a/AMDiS/src/FirstOrderAssembler.h
+++ b/AMDiS/src/FirstOrderAssembler.h
@@ -1,14 +1,11 @@
 #ifndef AMDIS_FIRST_ORDER_ASSEMBLER_H
 #define AMDIS_FIRST_ORDER_ASSEMBLER_H
 
+#include "AMDiS_fwd.h"
 #include "SubAssembler.h"
 
 namespace AMDiS {
 
-  class Q10PsiPhi;
-  class Q01PsiPhi;
-  class Q1Psi;
-
   /**
    * \ingroup Assembler
    * 
@@ -79,13 +76,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    ///
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat);
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo *, ElementVector&);
   };
@@ -106,13 +96,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    /// 
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat);
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo*, ElementVector&) {
       ERROR_EXIT("should not be called\n");
@@ -135,16 +118,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    ///
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat) 
-    {
-      ERROR_EXIT("CEM not yet\n");
-    }
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo *, ElementVector&);
   };
@@ -164,16 +137,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    ///
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat) 
-    {
-      ERROR_EXIT("CEM not yet\n");
-    }
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo*, ElementVector&) {
       ERROR_EXIT("should not be called\n");
@@ -196,16 +159,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    ///
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat) 
-    {
-      ERROR_EXIT("CEM not yet\n");
-    }
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo*, ElementVector&);
 
@@ -235,16 +188,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    ///
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat)
-    {
-      ERROR_EXIT("CEM not yet\n");
-    }
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo*, ElementVector&) {
       ERROR_EXIT("should not be called\n");
diff --git a/AMDiS/src/SecondOrderAssembler.cc b/AMDiS/src/SecondOrderAssembler.cc
index eea1e692..603b9f19 100644
--- a/AMDiS/src/SecondOrderAssembler.cc
+++ b/AMDiS/src/SecondOrderAssembler.cc
@@ -70,11 +70,11 @@ namespace AMDiS {
     if (!optimized) {
       newAssembler = NEW Stand2(op, assembler, quad);
     } else {
-	if (pwConst) {
-	    newAssembler = NEW Pre2(op, assembler, quad);
-	} else {
-	    newAssembler = NEW Quad2(op, assembler, quad);
-	}
+      if (pwConst) {
+	newAssembler = NEW Pre2(op, assembler, quad);
+      } else {
+	newAssembler = NEW Quad2(op, assembler, quad);
+      }
     }
 
     subAssemblers->push_back(newAssembler);
@@ -321,114 +321,4 @@ namespace AMDiS {
     delete [] LALt;
   }
 
-  void Stand2::calculateElementMatrix(const ElInfo *rowElInfo, 
-				      const ElInfo *colElInfo,
-				      const ElInfo *smallElInfo,
-				      const ElInfo *largeElInfo,
-				      ElementMatrix& mat)
-  {
-    FUNCNAME("Stand2::calculateElementMatrix()");
-
-    TEST_EXIT((nRow <= 3) && (nCol <= 3))("not yet!\n");
-
-    bool p = false;
-    /*
-    if (smallElInfo->getLevel() == largeElInfo->getLevel() + 1) {
-      p = true;
-      std::cout << "--------------------------------------" << std::endl;
-      std::cout << "i-th child = " << smallElInfo->getIChild() << std::endl;
-      std::cout << "Mesh row = " << owner->getRowFESpace()->getMesh() << std::endl;
-      std::cout << "Mesh col = " << owner->getColFESpace()->getMesh() << std::endl;
-      std::cout << "Mesh main = " << rowElInfo->getMesh() << std::endl;
-      std::cout << "Mesh aux = " << colElInfo->getMesh() << std::endl;
-    }
-    */
-
-    DimVec<double> grdPsi(dim, NO_INIT);
-    VectorOfFixVecs<DimVec<double> > grdPhi(dim, nCol, NO_INIT);
-
-    const BasisFunction *psi = owner->getRowFESpace()->getBasisFcts();
-    const BasisFunction *phi = owner->getColFESpace()->getBasisFcts();
-    const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
-
-    int nPoints = quadrature->getNumPoints();
-
-    DimMat<double> **LALt = NEW DimMat<double>*[nPoints];
-    for (int iq = 0; iq < nPoints; iq++) {
-      LALt[iq] = NEW DimMat<double>(dim, NO_INIT);
-      LALt[iq]->set(0.0);
-    }
-
-    int myRank = omp_get_thread_num();
-
-    if (p) {
-      std::cout << "Terms = " << terms[myRank].size() << std::endl;
-    }
-
-    for (int i = 0; i < static_cast<int>(terms[myRank].size()); i++) {
-      (static_cast<SecondOrderTerm*>(terms[myRank][i]))
-	->getLALt(smallElInfo, nPoints, LALt);
-    }
-
-    if (p) {
-      for (int i = 0; i < nPoints; i++) {
-	std::cout << "LALt for iq = " << i << std::endl;
-	LALt[i]->print();
-      }
-    }
-
-    DimVec<double> val1(dim, DEFAULT_VALUE, 0.0);
-
-    if (p)
-      std::cout << "Multiply LALt with det = " << smallElInfo->getDet() << std::endl;
-
-    for (int iq = 0; iq < nPoints; iq++) {
-      (*LALt[iq]) *= smallElInfo->getDet();
-    }
-     
-    /*
-    for (int i = 0; i < nRow; i++) {
-      for (int j = 0; j < nCol; j++) {	
-	for (int iq = 0; iq < nPoints; iq++) {
-
-	  val1.set(0.0);
-
-	  for (int k = 0; k < nCol; k++) {
-	    (*(phi->getGrdPhi(k)))(quadrature->getLambda(iq), grdPsi);
-	    grdPsi *= m[j][k];
-	    val1 += grdPsi;
-	  }
-
-	  (*(psi->getGrdPhi(i)))(quadrature->getLambda(iq), grdPsi);
-
-	  mat[i][j] += quadrature->getWeight(iq) *
-	      (grdPsi * ((*LALt[iq]) * val1));
-	}
-      }
-      }*/
-
-      for (int iq = 0; iq < nPoints; iq++) {
-	for (int i = 0; i < nCol; i++) {
-	  (*(phi->getGrdPhi(i)))(quadrature->getLambda(iq), grdPhi[i]);
-	  if (p) {
-	    std::cout << "grdPhi[" << i << "]" << std::endl;
-	    grdPhi[i].print();
-	  }
-	}
-
-	for (int i = 0; i < nRow; i++) {
-	  (*(psi->getGrdPhi(i)))(quadrature->getLambda(iq), grdPsi);
-	  for (int j = 0; j < nCol; j++) {
-	    mat[i][j] += quadrature->getWeight(iq) * 
-	      (grdPsi * ((*LALt[iq]) * grdPhi[j]));
-	  }
-	}
-      }
-
-    for (int iq = 0; iq < nPoints; iq++) 
-      delete LALt[iq];
-
-    delete [] LALt;
-  }
-
 }
diff --git a/AMDiS/src/SecondOrderAssembler.h b/AMDiS/src/SecondOrderAssembler.h
index fed24e08..15808461 100644
--- a/AMDiS/src/SecondOrderAssembler.h
+++ b/AMDiS/src/SecondOrderAssembler.h
@@ -62,13 +62,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    ///
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat);
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo *, ElementVector&) {
       ERROR_EXIT("should not be called\n");
@@ -92,16 +85,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    ///
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat) 
-    {
-      ERROR_EXIT("CEM not yet\n");
-    }
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo *, ElementVector&) {
       ERROR_EXIT("should not be called\n");
@@ -133,16 +116,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    ///
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat) 
-    {
-      ERROR_EXIT("Not yet implemented!\n");
-    }
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo *, ElementVector&) {
       ERROR_EXIT("should not be called\n");
diff --git a/AMDiS/src/SubAssembler.h b/AMDiS/src/SubAssembler.h
index 4eaf222d..1cadb91b 100644
--- a/AMDiS/src/SubAssembler.h
+++ b/AMDiS/src/SubAssembler.h
@@ -49,12 +49,6 @@ namespace AMDiS {
     virtual void calculateElementMatrix(const ElInfo *elInfo,
 					ElementMatrix& mat) = 0;
 
-    virtual void calculateElementMatrix(const ElInfo *rowElInfo,
-					const ElInfo *colElInfo,
-					const ElInfo *smallElInfo,
-					const ElInfo *largeElInfo,
-					ElementMatrix& mat) = 0;
-
     /** \brief
      * Calculates the element vector for elInfo and adds it to vec. Memory
      * for vec must be provided by the caller.
diff --git a/AMDiS/src/ZeroOrderAssembler.cc b/AMDiS/src/ZeroOrderAssembler.cc
index 8d040afe..b93e7b73 100644
--- a/AMDiS/src/ZeroOrderAssembler.cc
+++ b/AMDiS/src/ZeroOrderAssembler.cc
@@ -137,72 +137,6 @@ namespace AMDiS {
     }
   }
 
-  void StandardZOA::calculateElementMatrix(const ElInfo *rowElInfo,
-					   const ElInfo *colElInfo,
-					   const ElInfo *smallElInfo,
-					   const ElInfo *largeElInfo,
-					   ElementMatrix& mat) 
-  {
-    FUNCNAME("StandardZOA::calculateElementMatrix()");
-
-    TEST_EXIT((nRow <= 3) && (nCol <= 3))("not yet!\n");
-
-    calculateElementMatrix(smallElInfo, mat);
-
-    const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
-
-    ElementMatrix tmpMat(nRow, nRow);
-
-    for (int i = 0; i < nRow; i++) {
-      for (int j = 0; j < nRow; j++) {
-	tmpMat[i][j] = 0.0;
-
-	for (int k = 0; k < nRow; k++) {
-	  tmpMat[i][j] += m[i][j] * mat[j][i];
-	}
-      }
-    }
-
-    for (int i = 0; i < nRow; i++) {
-      for (int j = 0; j < nRow; j++) {
-	mat[i][j] = tmpMat[i][j];
-      }
-    }
-
-    /*
-    const BasisFunction *psi = owner->getRowFESpace()->getBasisFcts();
-    const BasisFunction *phi = owner->getColFESpace()->getBasisFcts();
-    const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
-
-    int myRank = omp_get_thread_num();
-    int nPoints = quadrature->getNumPoints();
-    std::vector<double> c(nPoints, 0.0);
-
-    std::vector<OperatorTerm*>::iterator termIt;
-    for (termIt = terms[myRank].begin(); termIt != terms[myRank].end(); ++termIt)
-      (static_cast<ZeroOrderTerm*>((*termIt)))->getC(smallElInfo, nPoints, c);
-
-    for (int i = 0; i < nRow; i++) {
-      for (int j = 0; j < nCol; j++) {
-
-	for (int iq = 0; iq < nPoints; iq++) {
-	  double val0 = c[iq] * smallElInfo->getDet() * quadrature->getWeight(iq) * 
-	    (*(phi->getPhi(i)))(quadrature->getLambda(iq));
-
-	  double val1 = 0.0;
-	  for (int k = 0; k < nCol; k++)
-	    val1 += m[j][k] * (*(psi->getPhi(k)))(quadrature->getLambda(iq));
-
-	  val0 *= val1;
-
-	  mat[i][j] += val0;
-	}
-      }
-    }
-    */
-
-  }
-
   void StandardZOA::calculateElementVector(const ElInfo *elInfo, ElementVector& vec)
   {
     int nPoints = quadrature->getNumPoints();
@@ -285,62 +219,6 @@ namespace AMDiS {
     }
   }
 
-  void FastQuadZOA::calculateElementMatrix(const ElInfo *rowElInfo,
-					   const ElInfo *colElInfo,
-					   const ElInfo *smallElInfo,
-					   const ElInfo *largeElInfo,
-					   ElementMatrix& mat) 
-  {
-    FUNCNAME("FastQuadZOA::calculateElementMatrix()");
-
-    ERROR_EXIT("CHECK FIRST, IF IT WILL REALLY WORK!\n");
-
-    int nPoints = quadrature->getNumPoints();
-    int myRank = omp_get_thread_num();
-    const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
-
-    if (firstCall) {
-#ifdef _OPENMP
-#pragma omp critical
-#endif 
-      {
-	const BasisFunction *basFcts = owner->getRowFESpace()->getBasisFcts();
-	psiFast = updateFastQuadrature(psiFast, basFcts, INIT_PHI);
-	basFcts = owner->getColFESpace()->getBasisFcts();
-	phiFast = updateFastQuadrature(phiFast, basFcts, INIT_PHI);
-	firstCall = false;
-      }
-    }
-
-    std::vector<double> c(nPoints, 0.0);
-    std::vector<OperatorTerm*>::iterator termIt;
-    for (termIt = terms[myRank].begin(); termIt != terms[myRank].end(); ++termIt) {
-      (static_cast<ZeroOrderTerm*>((*termIt)))->getC(rowElInfo, nPoints, c);
-    }
-
-    for (int iq = 0; iq < nPoints; iq++) {
-      c[iq] *= smallElInfo->getDet();
-      
-      const double *psi = psiFast->getPhi(iq);
-      const double *phi = phiFast->getPhi(iq);
-
-      for (int i = 0; i < nCol; i++) { 
-	for (int j = 0; j < nRow; j++) { 
-	  double val = quadrature->getWeight(iq) * c[iq] * psi[i];
-
-	  double tmpval = 0.0;
-	  for (int k = 0; k < nCol; k++) {
-	    tmpval += m[j][k] * phi[k];
-	  }
-	  val *= tmpval;
-
-	  mat[j][i] += val;
-	}
-      }
-
-    }    
-  }
-
   void FastQuadZOA::calculateElementVector(const ElInfo *elInfo, ElementVector& vec)
   {
     int myRank = omp_get_thread_num();
diff --git a/AMDiS/src/ZeroOrderAssembler.h b/AMDiS/src/ZeroOrderAssembler.h
index 99bf4578..0b3e63fa 100644
--- a/AMDiS/src/ZeroOrderAssembler.h
+++ b/AMDiS/src/ZeroOrderAssembler.h
@@ -2,13 +2,11 @@
 #define AMDIS_ZERO_ORDER_ASSEMBLER_H
 
 #include <vector>
+#include "AMDiS_fwd.h"
 #include "SubAssembler.h"
 
 namespace AMDiS {
 
-  class Q00PsiPhi;
-  class Q0Psi;
-
   /**
    * \ingroup Assembler
    * 
@@ -65,12 +63,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat);
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo *elInfo, ElementVector& vec);
   };
@@ -91,12 +83,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat);
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo *elInfo, ElementVector& vec);
   };
@@ -117,17 +103,6 @@ namespace AMDiS {
     /// Implements SubAssembler::calculateElementMatrix().
     void calculateElementMatrix(const ElInfo *elInfo, ElementMatrix& mat);
 
-    void calculateElementMatrix(const ElInfo *rowElInfo,
-				const ElInfo *colElInfo,
-				const ElInfo *smallElInfo,
-				const ElInfo *largeElInfo,
-				ElementMatrix& mat) 
-    {
-      FUNCNAME("PrecalcZOA::calculateElementMatrix()");
-      
-      ERROR_EXIT("CEM not yet\n");
-    }
-
     /// Implements SubAssembler::calculateElementVector().
     void calculateElementVector(const ElInfo *elInfo, ElementVector& vec);
 
-- 
GitLab