I had a nice day but am still hungry....

AMDiS/AMDISConfig.cmake.in

@@ -60,6 +60,7 @@ endif(Boost_FOUND)
 
 set(AMDIS_NEED_ZOLTAN @ENABLE_ZOLTAN@)
 set(AMDIS_HAS_PARALLEL_DOMAIN @ENABLE_PARALLEL_DOMAIN@)
+set(AMDIS_OPENMP @ENABLE_OPENMP@)
 set(AMDIS_NEED_UMFPACK @ENABLE_UMFPACK@)
 set(AMDIS_NEED_PNG @ENABLE_PNG@)
 set(AMDIS_NEED_BDDCML @ENABLE_BDDCML@)
@@ -72,6 +73,20 @@ set(AMDIS_VERSION @CurrentRevision@)
 set(AMDIS_MAJOR_VERSION @AMDIS_MAJOR@)
 set(AMDIS_MINOR_VERSION @AMDIS_MINOR@)
 
+if(CMAKE_BUILD_TYPE STREQUAL "")
+	set(CMAKE_BUILD_TYPE "Release")
+endif()
+
+if(AMDIS_OPENMP)
+  if(CMAKE_CXX_COMPILER MATCHES ".*icpc")
+        list(APPEND AMDIS_COMPILEFLAGS "-openmp")
+        list(APPEND AMDIS_LIBRARIES "-openmp")
+  else()
+    list(APPEND AMDIS_COMPILEFLAGS "-fopenmp")
+        list(APPEND AMDIS_LIBRARIES "-fopenmp")
+  endif()
+endif(AMDIS_OPENMP)
+
 if(AMDIS_NEED_UMFPACK)
 	set(AMDIS_UMFPACK_PATH @UMFPACK_PATH@)
 	list(APPEND AMDIS_INCLUDE_DIRS ${AMDIS_UMFPACK_PATH})

AMDiS/src/Debug.cc

@@ -66,13 +66,14 @@ namespace AMDiS {
 #endif
 
     
-    void writeDofIndexMesh(const FiniteElemSpace *feSpace)
+    void writeDofIndexMesh(const FiniteElemSpace *feSpace,
+			   string filename)
     {
       DOFVector<double> tmp(feSpace, "tmp");
       DOFIterator<double> it(&tmp, USED_DOFS);
       for (it.reset(); !it.end(); ++it)
 	*it = it.getDOFIndex();
-      VtkWriter::writeFile(tmp, "dofindex.vtu");
+      VtkWriter::writeFile(tmp, filename);
     }
 
 

AMDiS/src/Debug.h

@@ -68,9 +68,11 @@ namespace AMDiS {
      * Create a vtu file with name 'dofindex.vtu'. All nodes in the mesh are colored
      * by the global DOF index.
      *
-     * \param[in]  feSpace  The FE space to be used.
+     * \param[in]  feSpace   The FE space to be used.
+     * \param[in]  filename  Name of the VTU file
      */
-    void writeDofIndexMesh(const FiniteElemSpace *feSpace);
+    void writeDofIndexMesh(const FiniteElemSpace *feSpace, 
+			   string filename = "dofindex.vtu");
 
     void colorEdgeInMesh(const FiniteElemSpace *feSpace,
 			 Element *el, 

AMDiS/src/FirstOrderTerm.cc

@@ -638,17 +638,23 @@ namespace AMDiS {
     }
   }
 
-  void Vec3FctAtQP_FOT::eval( int nPoints,
-			      const mtl::dense_vector<double>& uhAtQP,
-			      const mtl::dense_vector<WorldVector<double> >& grdUhAtQP,
-			      const mtl::dense_vector<WorldMatrix<double> >& D2UhAtQP,
-			      mtl::dense_vector<double>& result,
-			      double fac) 
+  void Vec3FctAtQP_FOT::eval(int nPoints,
+			     const mtl::dense_vector<double>& uhAtQP,
+			     const mtl::dense_vector<WorldVector<double> >& grdUhAtQP,
+			     const mtl::dense_vector<WorldMatrix<double> >& D2UhAtQP,
+			     mtl::dense_vector<double>& result,
+			     double fac) 
   {
-    if (num_rows(grdUhAtQP) > 0)
+    if (num_rows(grdUhAtQP) == 0)
+      return;
+
+    if (bOne > -1) {
+      ERROR_EXIT("Not yet implemented!\n");
+    } else {
       for (int iq = 0; iq < nPoints; iq++)
 	result[iq] += fac * (*f)(vec1AtQPs[iq], vec2AtQPs[iq] ,vec3AtQPs[iq]) * 
 	  ((*b) * grdUhAtQP[iq]);
+    }
   }
 
 

AMDiS/src/ProblemStat.cc

@@ -864,6 +864,12 @@ namespace AMDiS {
 	if (assembleMatrix && matrix->getBoundaryManager())
 	  matrix->getBoundaryManager()->initMatrix(matrix);
 
+	// This OpenMP barrier is required for the case of periodc boundary 
+	// conditions. In this case, we may have a data race between 
+	// fillBoundaryCondition and exitMatrix, where both make use of the
+	// vertex vector associated.
+#pragma omp barrier
+
 	// The simplest case: either the right hand side has no operaters, no aux
 	// fe spaces, or all aux fe spaces are equal to the row and col fe space.
 	assembleOnOneMesh(componentSpaces[rowComponent],
@@ -872,12 +878,14 @@ namespace AMDiS {
 			  ((rowComponent == colComponent) && asmVector) ? 
 			  rhs->getDOFVector(rowComponent) : 
 			  NULL);
+
+#pragma omp barrier
       
 	assembledMatrix[rowComponent][colComponent] = true;
 
 	if (assembleMatrix)
 	  matrix->finishInsertion();
-
+	
  	if (assembleMatrix && matrix->getBoundaryManager())
  	  matrix->getBoundaryManager()->exitMatrix(matrix);	
 

AMDiS/src/SecondOrderTerm.h

@@ -66,17 +66,16 @@ namespace AMDiS {
 	      bool symm,
 	      double factor) const;
 
-    /** \brief
-     * Evaluation of \f$ \Lambda \cdot A \cdot \Lambda^t\f$ for \f$ A \f$
-     * the matrix having a ONE in the position \f$ (K,L) \f$
-     * and ZEROS in all other positions.
-     */
+    /// Evaluation of \f$ \Lambda \cdot A \cdot \Lambda^t\f$ for \f$ A \f$
+    /// the matrix having a ONE in the position \f$ (K,L) \f$
+    /// and ZEROS in all other positions.
     static void lalt_kl(const DimVec<WorldVector<double> >& Lambda,
 			int k, int l,
 			mtl::dense2D<double>& LALt,
 			double factor);
 
-    /// Evaluation of \f$ \Lambda \cdot A \cdot \Lambda^t\f$ for A equal to the identity.
+    /// Evaluation of \f$ \Lambda \cdot A \cdot \Lambda^t\f$ for A equal to 
+    /// the identity.
     inline void l1lt(const DimVec<WorldVector<double> >& Lambda, 
 		     mtl::dense2D<double>& LALt,
 		     double factor) const