Improved ARH file reader and writer.

AMDiS/src/Initfile.cc

@@ -87,14 +87,10 @@ namespace AMDiS {
         std::string paramName = variableReplacement(InitfileInternal::trim(parser.name));
         std::string paramValue = variableReplacement(InitfileInternal::trim(parser.value));
         operator[](paramName) = paramValue;
-#ifndef DEBUG
         int info = 0;
         get("parameter information", info, 0);
         if (info >= 2)
-#endif
-        {
-          std::cout << "read [" << paramName << " => " << paramValue << "]" << std::endl;
-        }
+          std::cout << "read [" << paramName << " => " << paramValue << "]\n";
       } else if (pos0 != std::string::npos &&
 		 sw[pos0] == '#'  &&
 		 static_cast<size_t>(sw.find("#include")) == pos0) {

AMDiS/src/io/ArhReader.cc

@@ -27,7 +27,9 @@ namespace AMDiS {
   void ArhReader::read(string filename, Mesh *mesh,
 		       DOFVector<double>* vec0, 
 		       DOFVector<double>* vec1,
-		       DOFVector<double>* vec2)
+		       DOFVector<double>* vec2,
+		       bool writeParallel,
+		       int nProcs)
   {
     vector<DOFVector<double>*> vecs(0);
     if (vec0 != NULL)
@@ -37,15 +39,50 @@ namespace AMDiS {
     if (vec2 != NULL)
       vecs.push_back(vec2);
 
-    ArhReader::read(filename, mesh, vecs);
+    ArhReader::read(filename, mesh, vecs, writeParallel, nProcs);
   }
 
 
   void ArhReader::read(string filename, Mesh *mesh, 
-		       vector<DOFVector<double>*> vecs)
+		       vector<DOFVector<double>*> vecs,
+		       bool writeParallel,
+		       int nProcs)
   {
     FUNCNAME("ArhReader::read()");
 
+#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
+    if (writeParallel) {
+      using boost::lexical_cast;
+      int sPos = filename.find(".arh");
+      TEST_EXIT(sPos >= 0)("Failed to find file postfix!\n");
+      string name = filename.substr(0, sPos);      
+
+      if (nProcs == -1) {     
+	string procFilename = name + "-p" + lexical_cast<string>(MPI::COMM_WORLD.Get_rank()) + "-.arh";
+	readFile(procFilename, mesh, vecs);
+      } else {
+	for (int i = 0; i < nProcs; i++) {
+	  string procFilename = name + "-p" + lexical_cast<string>(i) + "-.arh";
+	  readFile(procFilename, mesh, vecs);
+	}
+      }
+    } else {
+      readFile(filename, mesh, vecs);
+    }
+#else
+    readFile(filename, mesh, vecs);
+#endif
+
+    MSG("ARH file read from: %s\n", filename.c_str());
+  }
+
+  
+  void ArhReader::readFile(string filename, 
+			   Mesh *mesh,
+			   vector<DOFVector<double>*> vecs)
+  {
+    FUNCNAME("ArhReader::readFile()");
+
     // === Get set of all macro elements in mesh. ===
     std::set<int> macroInMesh;
     for (std::deque<MacroElement*>::iterator it = mesh->getMacroElements().begin();
@@ -65,6 +102,7 @@ namespace AMDiS {
       ERROR_EXIT("Should not happen!\n");
     }
 
+
     ifstream file;
     file.open(filename.c_str(), ios::in | ios::binary);
 
@@ -112,7 +150,7 @@ namespace AMDiS {
 
     file.close();
 
-    MSG("ARH file read from: %s\n", filename.c_str());
+    delete refManager;
   }
 
 

AMDiS/src/io/ArhReader.h

@@ -27,20 +27,31 @@
 
 namespace AMDiS {
 
+  using namespace std;
+
   class ArhReader
   {
   public:
-    static void read(std::string filename, Mesh *mesh,
+    static void read(string filename, Mesh *mesh,
 		     DOFVector<double>* vec0 = NULL, 
 		     DOFVector<double>* vec1 = NULL,
-		     DOFVector<double>* vec2 = NULL);
+		     DOFVector<double>* vec2 = NULL,		     
+		     bool writeParallel = true,
+		     int nProcs = -1);
 
-    static void read(std::string filename, Mesh *mesh, 
-		     std::vector<DOFVector<double>*> vecs);
+    static void read(string filename, Mesh *mesh, 
+		     vector<DOFVector<double>*> vecs,
+		     bool writeParallel = true,
+		     int nProcs = -1);
 
   private:
+
+    static void readFile(string filename, 
+			 Mesh *mesh,
+			 vector<DOFVector<double>*> vecs);
+
     static void setDofValues(int macroElIndex, Mesh *mesh, 
-			     std::vector<double>& values, DOFVector<double>* vec);
+			     vector<double>& values, DOFVector<double>* vec);
   };
 }
 

AMDiS/src/io/ArhWriter.cc

@@ -42,10 +42,22 @@ namespace AMDiS {
 
 
   void ArhWriter::write(string filename, Mesh *mesh, 
-			vector<DOFVector<double>*> vecs)
+			vector<DOFVector<double>*> vecs,
+			bool writeParallel)
   {
     FUNCNAME("ArhWriter::write()");
 
+#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
+    if (writeParallel) {
+      using boost::lexical_cast;
+      
+      int sPos = filename.find(".arh");
+      TEST_EXIT(sPos >= 0)("Failed to find file postfix!\n");
+      string name = filename.substr(0, sPos);      
+      filename = name + "-p" + lexical_cast<string>(MPI::COMM_WORLD.Get_rank()) + "-.arh";
+    }
+#endif
+
     ofstream file;
     file.open(filename.c_str(), ios::out | ios::binary | ios::trunc);
 

AMDiS/src/io/ArhWriter.h

@@ -39,7 +39,8 @@ namespace AMDiS {
 		      DOFVector<double>* vec2 = NULL);
 
     static void write(std::string filename, Mesh *mesh, 
-		      std::vector<DOFVector<double>*> vecs);
+		      std::vector<DOFVector<double>*> vecs,
+		      bool writeParallel = true);
 
   protected:
     static void writeMacroElement(std::ofstream &file, 

AMDiS/src/io/VtkWriter.hh

@@ -143,7 +143,7 @@ namespace AMDiS {
     DOFVector<list<WorldVector<double> > >::Iterator coordIt(dofCoords, USED_DOFS);
 
     file << scientific;
-    file.precision(5);
+    file.precision(15);
          
     // Write the values for all vertex DOFs.
     for (intPointIt.reset(), valueIt.reset(), coordIt.reset();

AMDiS/src/parallel/MeshDistributor.cc

@@ -857,6 +857,14 @@ namespace AMDiS {
 
     // === Print imbalance factor. ===
 
+    printImbalanceFactor();
+  }
+
+  
+  void MeshDistributor::printImbalanceFactor()
+  {
+    FUNCNAME("MeshDistributor::printImbalanceFactor()");
+
     vector<int> nDofsInRank(mpiSize);
     int nDofs = mesh->getDofAdmin(0).getUsedDofs();
     mpiComm.Gather(&nDofs, 1, MPI_INT, &(nDofsInRank[0]), 1, MPI_INT, 0);
@@ -1328,6 +1336,26 @@ namespace AMDiS {
     check3dValidMesh();
 
     MSG("Mesh repartitioning needed %.5f seconds\n", MPI::Wtime() - timePoint);
+
+    
+    // === Print DOF information to screen. ===
+
+    nDofs = mesh->getDofAdmin(0).getUsedDofs();
+    mpiComm.Gather(&nDofs, 1, MPI_INT, &(nDofsInRank[0]), 1, MPI_INT, 0);
+
+    if (mpiRank == 0) {
+      int nOverallDofs = 0;
+      int minDofs = numeric_limits<int>::max();
+      int maxDofs = numeric_limits<int>::min();
+      for (int i = 0; i < mpiSize; i++) {
+	nOverallDofs += nDofsInRank[i];
+	minDofs = std::min(minDofs, nDofsInRank[i]);
+	maxDofs = std::max(maxDofs, nDofsInRank[i]);
+      }      
+     
+      MSG("Overall DOFs: %d    Min DOFs: %d    Max DOFs: %d\n", 
+	  nOverallDofs, minDofs, maxDofs);
+    }
   }
 
 

AMDiS/src/parallel/MeshDistributor.h

@@ -92,6 +92,16 @@ namespace AMDiS {
      */
     void checkMeshChange(bool tryRepartition = true);
 
+    /** \brief
+     * Checks if is required to repartition the mesh. If this is the case, a new
+     * partition will be created and the mesh will be redistributed between the
+     * ranks.
+     */
+    void repartitionMesh();
+    
+    /// Calculates the imbalancing factor and prints it to screen.
+    void printImbalanceFactor();
+
     /** \brief
      * Test, if the mesh consists of macro elements only. The mesh partitioning 
      * of the parallelization works for macro meshes only and would fail, if the 
@@ -380,13 +390,6 @@ namespace AMDiS {
      */
     bool checkAndAdaptBoundary(RankToBoundMap &allBound);
   
-    /** \brief
-     * Checks if is required to repartition the mesh. If this is the case, a new
-     * partition will be created and the mesh will be redistributed between the
-     * ranks.
-     */
-    void repartitionMesh();
-
     /// Sets \ref isRankDof to all matrices and rhs vectors in a given 
     /// stationary problem.
     void setRankDofs(ProblemStatSeq *probStat);

AMDiS/src/parallel/PetscSolverGlobalMatrix.cc

@@ -118,7 +118,8 @@ namespace AMDiS {
   {
     FUNCNAME("PetscSolverGlobalMatrix::fillPetscRhs()");
 
-    TEST_EXIT_DBG(vec)("NO DOF vector defined!\n");
+    TEST_EXIT_DBG(vec)("No DOF vector defined!\n");
+    TEST_EXIT_DBG(meshDistributor)("No mesh distributor defined!\n");
 
     int nComponents = vec->getSize();
     int nRankRows = meshDistributor->getNumberRankDofs() * nComponents;

AMDiS/src/time/RosenbrockAdaptInstationary.cc

@@ -56,7 +56,8 @@ namespace AMDiS {
       tauGamma(1.0),
       minusTauGamma(-1.0),
       invTauGamma(1.0),
-      minusInvTauGamma(-1.0)
+      minusInvTauGamma(-1.0),
+      dbgTimestepStudy(false)
   {
     FUNCNAME("RosenbrockAdaptInstationary::RosenbrockAdaptInstationary()");
     initConstructor(&problemStat);