Refactoring of using periodic DOF mapping in parallel computations.

AMDiS/src/AMDiS.h

@@ -81,7 +81,6 @@
 #include "Operator.h"
 #include "OperatorTerm.h"
 #include "Parametric.h"
-#include "PeriodicMap.h"
 #include "PeriodicBC.h"
 #include "ProblemStat.h"
 #include "ProblemInstat.h"

AMDiS/src/Mesh.cc

@@ -34,7 +34,6 @@
 #include "CoarseningManager.h"
 #include "DOFIterator.h"
 #include "VertexVector.h"
-#include "PeriodicMap.h"
 #include "Projection.h"
 #include "ElInfoStack.h"
 #include "Serializer.h"

AMDiS/src/PeriodicMap.h

-// ============================================================================
-// ==                                                                        ==
-// == AMDiS - Adaptive multidimensional simulations                          ==
-// ==                                                                        ==
-// ==  http://www.amdis-fem.org                                              ==
-// ==                                                                        ==
-// ============================================================================
-//
-// Software License for AMDiS
-//
-// Copyright (c) 2010 Dresden University of Technology 
-// All rights reserved.
-// Authors: Simon Vey, Thomas Witkowski et al.
-//
-// This file is part of AMDiS
-//
-// See also license.opensource.txt in the distribution.
-
-
-
-/** \file PeriodicMap.h */
-
-#ifndef AMDIS_PERIODICMAP_H
-#define AMDIS_PERIODICMAP_H
-
-#include "Global.h"
-
-#include <map>
-
-namespace AMDiS {
-
-  class PeriodicMap
-  {
-  public:
-    void setEntry(DegreeOfFreedom key, DegreeOfFreedom entry) {
-      // no trivial entries!
-      if (key == entry) 
-	return;
-
-      // if a key equal to entry exists ...
-      if (getEntry(entry) >= 0) {
-	if (getEntry(entry) == key) {
-	  return;
-	}
-	// ... let it be equal entries
-	setEntry(key, getEntry(entry));
-	return;
-      }
-
-      // replace entries equal to key
-      std::map<DegreeOfFreedom, DegreeOfFreedom>::iterator it;
-      for (it = periodicMap.begin(); it != periodicMap.end(); ++it) {
-	if (it->second == key) {
-	  it->second = entry;
-	}
-      }
-    
-      // if key exists already
-      if (getEntry(key) >= 0) {
-	// insert new entry with old entry as key
-	setEntry(getEntry(key), entry);
-      }
-
-      // set entry
-      periodicMap[key] = entry;
-    };
-
-    DegreeOfFreedom getEntry(DegreeOfFreedom key) {
-      std::map<DegreeOfFreedom, DegreeOfFreedom>::iterator it;
-      it = periodicMap.find(key);
-      if (it != periodicMap.end()) {
-	return it->second;
-      } else {
-	return -1;
-      }
-    };
-
-  
-
-  protected:
-    std::map<DegreeOfFreedom, DegreeOfFreedom> periodicMap;
-  };
-
-}
-
-#endif

AMDiS/src/io/MacroReader.cc

@@ -10,6 +10,10 @@
 // See also license.opensource.txt in the distribution.
 
 
+#include <string.h>
+#include <map>
+#include <iostream>
+#include <fstream>
 #include "MacroReader.h"
 #include "MacroWriter.h"
 #include "MacroElement.h"
@@ -17,17 +21,12 @@
 #include "Boundary.h"
 #include "FiniteElemSpace.h"
 #include "Mesh.h"
-#include <string.h>
 #include "FixVec.h"
-#include "PeriodicMap.h"
 #include "ElInfo.h"
 #include "Initfile.h"
 #include "DOFIterator.h"
 #include "LeafData.h"
 #include "VertexVector.h"
-#include <map>
-#include <iostream>
-#include <fstream>
 
 namespace AMDiS {
 
@@ -64,7 +63,7 @@ namespace AMDiS {
       int mode = -1; // 0: drop dofs, 1: associate dofs
       int result;
       BoundaryType boundaryType;
-      PeriodicMap periodicMap;
+      MacroReader::PeriodicMap periodicMap;
 
       fscanf(file, "%*s %d", &n);
       fscanf(file, "%*s %*s %*s %*s %*s %*s %*s %*s %*s %*s %*s");

AMDiS/src/io/MacroReader.h

@@ -41,6 +41,53 @@ namespace AMDiS {
    */
   class MacroReader
   {
+    // Is used to read periodic macros
+    class PeriodicMap
+    {
+    public:
+      void setEntry(DegreeOfFreedom key, DegreeOfFreedom entry) 
+      {
+	// no trivial entries!
+	if (key == entry) 
+	  return;
+	
+	// if a key equal to entry exists ...
+	if (getEntry(entry) >= 0) {
+	  if (getEntry(entry) == key)
+	    return;
+	  
+	  // ... let it be equal entries
+	  setEntry(key, getEntry(entry));
+	  return;
+	}
+	
+	// replace entries equal to key
+	for (map<DegreeOfFreedom, DegreeOfFreedom>::iterator it = 
+	       periodicMap.begin(); it != periodicMap.end(); ++it)
+	  if (it->second == key)
+	    it->second = entry;
+	
+	// if key exists already, insert new entry with old entry as key
+	if (getEntry(key) >= 0)
+	  setEntry(getEntry(key), entry);
+	
+	// set entry
+	periodicMap[key] = entry;
+      }
+            
+      DegreeOfFreedom getEntry(DegreeOfFreedom key) 
+      {
+	map<DegreeOfFreedom, DegreeOfFreedom>::iterator it = periodicMap.find(key);
+	if (it != periodicMap.end())
+	  return it->second;
+	return -1;
+      }
+    
+    protected:
+      map<DegreeOfFreedom, DegreeOfFreedom> periodicMap;
+    };
+
+
   public:
     /// Creates a Mesh by reading the macro file with the given filename.
     static MacroInfo* readMacro(std::string filename, 

AMDiS/src/parallel/InteriorBoundary.h

@@ -95,6 +95,8 @@ namespace AMDiS {
     ExcludeList excludedSubstructures;
   };
 
+
+
   /** \brief 
    * Defines one atomic part of the boundary, i.e., two boundary objects where the
    * boundary goes through.
@@ -118,6 +120,8 @@ namespace AMDiS {
     BoundaryType type;
   };
 
+
+
   /** \brief
    * Defines the interior boundary, i.e. a bound within the domain. It is used for
    * the classical domain decomposition parallelization.

AMDiS/src/parallel/MeshDistributor.cc

@@ -2014,8 +2014,7 @@ namespace AMDiS {
 
     // Clear all periodic DOF mappings calculated before. We do it from scratch.
     periodicDofs.clear();
-    periodicDofMap.clear();
-    periodicDofAssociations.clear();
+    periodicMap.clear();
 
     for (unsigned int i = 0; i < feSpaces.size(); i++)
       createPeriodicMap(feSpaces[i]);
@@ -2064,10 +2063,8 @@ namespace AMDiS {
 	    DegreeOfFreedom globalDof1 = 
 	      dofFeData[feSpace].mapLocalGlobalDofs[*(dofs1[j])];
 
-	    if (periodicDofAssociations[feSpace][globalDof0].count(type) == 0) {
-	      periodicDofMap[feSpace][type][globalDof0] = globalDof1;
-	      periodicDofAssociations[feSpace][globalDof0].insert(type);
-	    }
+	    if (!periodicMap.isPeriodicOnBound(feSpace, type, globalDof0))
+	      periodicMap.add(feSpace, type, globalDof0, globalDof1);
 	  }
 	}
 
@@ -2121,10 +2118,8 @@ namespace AMDiS {
 
 	// Check if this global DOF with the corresponding boundary type was
 	// not added before by another periodic boundary from other rank.
-	if (periodicDofAssociations[feSpace][globalDofIndex].count(type) == 0) {
-	  periodicDofMap[feSpace][type][globalDofIndex] = mapGlobalDofIndex;
-	  periodicDofAssociations[feSpace][globalDofIndex].insert(type);
-	}
+	if (!periodicMap.isPeriodicOnBound(feSpace, type, globalDofIndex))
+	  periodicMap.add(feSpace, type, globalDofIndex, mapGlobalDofIndex);
       }
     }
 
@@ -2151,17 +2146,15 @@ namespace AMDiS {
 	  DegreeOfFreedom globalDof = 
 	    dofFeData[feSpace].mapLocalGlobalDofs[*dofs[i]];
 
-	  TEST_EXIT_DBG(periodicDofAssociations[feSpace].count(globalDof))
-	    ("Should hot happen!\n");      
+	  std::set<BoundaryType>& assoc = 
+	    periodicMap.getAssociations(feSpace, globalDof);
+	  TEST_EXIT_DBG(assoc.size() > 0)("Should not happen!\n");
 
-	  for (std::set<BoundaryType>::iterator perAscIt = periodicDofAssociations[feSpace][globalDof].begin();
-	       perAscIt != periodicDofAssociations[feSpace][globalDof].end(); ++perAscIt)
-	    if (*perAscIt >= -3) {
-	      TEST_EXIT_DBG(periodicDofMap[feSpace][*perAscIt].count(globalDof) == 1)
-		("Should not happen!\n");
+	  for (std::set<BoundaryType>::iterator perAscIt =  assoc.begin();
+	       perAscIt != assoc.end(); ++perAscIt)
+	    if (*perAscIt >= -3)
 	      perObjMap[*perAscIt][globalDof] = 
-		periodicDofMap[feSpace][*perAscIt][globalDof];
-	    }
+		periodicMap.map(feSpace, *perAscIt, globalDof);
 	}
       }
 
@@ -2172,19 +2165,12 @@ namespace AMDiS {
     stdMpi2.startCommunication();
 
     for (map<int, PeriodicDofMap>::iterator it = stdMpi2.getRecvData().begin();
-	 it != stdMpi2.getRecvData().end(); ++it) {
+	 it != stdMpi2.getRecvData().end(); ++it)
       for (PeriodicDofMap::iterator perIt = it->second.begin();
-	   perIt != it->second.end(); ++perIt) {
+	   perIt != it->second.end(); ++perIt)
 	for (DofMapping::iterator dofIt = perIt->second.begin();
-	     dofIt != perIt->second.end(); ++dofIt) {
-	  TEST_EXIT_DBG(periodicDofMap[feSpace][perIt->first].count(dofIt->second) == 0 ||
-			periodicDofMap[feSpace][perIt->first][dofIt->second] == dofIt->first)
-	    ("Should not happen!\n");
-
-	  periodicDofMap[feSpace][perIt->first][dofIt->second] = dofIt->first;
-	}
-      }
-    }
+	     dofIt != perIt->second.end(); ++dofIt)
+	  periodicMap.add(feSpace, perIt->first, dofIt->second, dofIt->first);
   }
 
 
@@ -2275,11 +2261,7 @@ namespace AMDiS {
       SerUtil::serialize(out, dofFeData[feSpaces[i]].mapLocalDofIndex);
     }
 
-    for (unsigned int i = 0; i < nFeSpace; i++)
-      serialize(out, periodicDofMap[feSpaces[i]]);
-
-    for (unsigned int i = 0; i < nFeSpace; i++)
-      serialize(out, periodicDofAssociations[feSpaces[i]]);
+    periodicMap.serialize(out, feSpaces);
 
     SerUtil::serialize(out, macroElementNeighbours);
 
@@ -2352,12 +2334,7 @@ namespace AMDiS {
       SerUtil::deserialize(in, dofFeData[feSpaces[i]].mapLocalDofIndex);
     }
 
-
-    for (unsigned int i = 0; i < nFeSpace; i++)
-      deserialize(in, periodicDofMap[feSpaces[i]]);
-
-    for (unsigned int i = 0; i < nFeSpace; i++)
-      deserialize(in, periodicDofAssociations[feSpaces[i]]);
+    periodicMap.deserialize(in, feSpaces);
 
     SerUtil::deserialize(in, macroElementNeighbours);
 
@@ -2379,73 +2356,4 @@ namespace AMDiS {
     deserialized = true;
   }
 
-
-  void MeshDistributor::serialize(ostream &out, PeriodicDofMap &data)
-  {
-    int mapSize = data.size();
-    SerUtil::serialize(out, mapSize);
-
-    for (PeriodicDofMap::iterator it = data.begin(); it != data.end(); ++it) {
-      int type = it->first;
-      DofMapping dofMap = it->second;
-
-      SerUtil::serialize(out, type);
-      SerUtil::serialize(out, dofMap);
-    }
-  }
-
-
-  void MeshDistributor::serialize(ostream &out, map<int, std::set<int> >& data)
-  {
-    int mapSize = data.size();
-    SerUtil::serialize(out, mapSize);
-
-    for (map<int, std::set<int> >::iterator it = data.begin(); it != data.end(); ++it) {
-      int dof = it->first;
-      std::set<int> typeSet = it->second;
-
-      SerUtil::serialize(out, dof);
-      SerUtil::serialize(out, typeSet);
-    }
-  }
-
-
-  void MeshDistributor::deserialize(istream &in, PeriodicDofMap &data)
-  {
-    data.clear();
-
-    int mapSize = 0;
-    SerUtil::deserialize(in, mapSize);
-
-    for (int i = 0; i < mapSize; i++) {
-      int type;
-      DofMapping dofMap;
-
-      SerUtil::deserialize(in, type);
-      SerUtil::deserialize(in, dofMap);
-
-      data[type] = dofMap;
-    }
-  }
-
-
-  void MeshDistributor::deserialize(istream &in, 
-				    map<int, std::set<int> >& data)
-  {
-    data.clear();
-
-    int mapSize = 0;
-    SerUtil::deserialize(in, mapSize);
-
-    for (int i = 0; i < mapSize; i++) {
-      int dof;
-      std::set<int> typeSet;
-
-      SerUtil::deserialize(in, dof);
-      SerUtil::deserialize(in, typeSet);
-
-      data[dof] = typeSet;
-    }
-  }
-
 }

AMDiS/src/parallel/MeshDistributor.h

@@ -30,6 +30,7 @@
 #include "parallel/ParallelTypes.h"
 #include "parallel/MeshPartitioner.h"
 #include "parallel/InteriorBoundary.h"
+#include "parallel/PeriodicMap.h"
 #include "parallel/StdMpi.h"
 #include "AMDiS_fwd.h"
 #include "Global.h"
@@ -277,58 +278,10 @@ namespace AMDiS {
       return dofFeData[feSpace].mapLocalDofIndex[dof];
     }
 
-    /// Returns the periodic mapping for all boundary DOFs in rank.
-    inline PeriodicDofMap& getPeriodicMapping(const FiniteElemSpace *feSpace)
+    /// Returns the periodic mapping handler, \ref periodicMap.
+    inline PeriodicMap& getPeriodicMap()
     {
-      return periodicDofMap[feSpace];
-    }
-
-    /// Returns for a global DOF index of a given FE space its periodic mapping
-    /// for a given  boundary type.
-    inline int getPeriodicMapping(const FiniteElemSpace *feSpace,
-				  BoundaryType type,
-				  int globalDofIndex) 
-    {
-      FUNCNAME("MeshDistributor::getPeriodicMapping()");
-
-      TEST_EXIT_DBG(periodicDofMap.count(feSpace))("Should not happen!\n");
-      TEST_EXIT_DBG(periodicDofMap[feSpace][type].count(globalDofIndex) == 1)
-	("There is no periodic association for global DOF %d for boundary type %d!\n",
-	 globalDofIndex, type);
-
-      return periodicDofMap[feSpace][type][globalDofIndex];
-    }
-
-    /// For a given global DOF index, this function returns the set of periodic
-    /// associations, i.e., the boundary types the DOF is associated to, for 
-    /// this DOF.
-    inline std::set<BoundaryType>& getPerDofAssociations(const FiniteElemSpace* feSpace,
-							 int globalDofIndex)
-    {
-      FUNCNAME("MeshDistributor::getPerDofAssociations()");
-
-      TEST_EXIT_DBG(periodicDofAssociations.count(feSpace))
-	("Should not happen!\n");
-      TEST_EXIT_DBG(periodicDofAssociations[feSpace].count(globalDofIndex)) 
- 	("Should not happen!\n"); 
-
-      return periodicDofAssociations[feSpace][globalDofIndex];
-    }
-
-    /// Returns true, if the DOF (global index) is a periodic DOF.
-    inline bool isPeriodicDof(const FiniteElemSpace *feSpace, int globalDofIndex)
-    {
-      return (periodicDofAssociations[feSpace].count(globalDofIndex) > 0 &&
-	       periodicDofAssociations[feSpace][globalDofIndex].size() > 0);
-    }
-
-    /// Returns true, if the DOF (global index) of a given FE space is a 
-    /// periodic DOF for the given boundary type.
-    inline bool isPeriodicDof(const FiniteElemSpace *feSpace,
-			      BoundaryType type,
-			      int globalDofIndex)
-    {
-      return (periodicDofMap[feSpace][type].count(globalDofIndex) > 0);
+      return periodicMap;
     }
 
     DofComm& getSendDofs()
@@ -556,16 +509,6 @@ namespace AMDiS {
 		     map<int, map<const FiniteElemSpace*, DofContainer> > &data,
 		     map<const FiniteElemSpace*, map<int, const DegreeOfFreedom*> > &dofMap);
 
-    /// Writes a periodic dof mapping to an output stream.
-    void serialize(ostream &out, PeriodicDofMap &data);
-
-    void serialize(ostream &out, map<int, std::set<int> >& data);
-
-    /// Reads a periodic dof mapping from an input stream.
-    void deserialize(istream &in, PeriodicDofMap &data);
-
-    void deserialize(istream &in, map<int, std::set<int> >& data);
-
     /// Writes a mapping from dof pointers to some values to an output stream.
     template<typename T>
     void serialize(ostream &out, map<const DegreeOfFreedom*, T> &data)
@@ -710,24 +653,8 @@ namespace AMDiS {
      */
     DofComm periodicDofs;
 
-    /** \brief
-     * If periodic boundaries are used, this map stores, for each periodic 
-     * boundary type, for all DOFs in rank's partition (that are on periodic 
-     * boundaries), the corresponding mapped periodic DOFs. The mapping is 
-     * defined by using global DOF indices.
-     */
-    PeriodicDofMapFeSpace periodicDofMap;
-    
-    /** \brief
-     * If periodic boundaries are used, this map stores to each periodic DOF in 
-     * rank's partition the set of periodic boundaries the DOF is associated to.
-     * In 2D, most DOFs are only on one periodic boundary. Only, e.g., in a box 
-     * with all boundaries being periodic, the four corners are associated by 
-     * two different boundaries.
-     */
-    map<const FiniteElemSpace*, map<DegreeOfFreedom, std::set<BoundaryType> > > periodicDofAssociations;
+    PeriodicMap periodicMap;
 
-    
     /// This set of values must be interchanged between ranks when the mesh is 
     /// repartitioned.
     vector<DOFVector<double>*> interchangeVectors;

AMDiS/src/parallel/ParallelDebug.cc

@@ -157,9 +157,11 @@ namespace AMDiS {
 
     // === 1. check: All periodic DOFs should have at least a correct number ===
     // === of periodic associations.                                         ===
-   
-    for (map<int, std::set<BoundaryType> >::iterator it = pdb.periodicDofAssociations[feSpace].begin();
-	 it != pdb.periodicDofAssociations[feSpace].end(); ++it) {
+
+    PeriodicMap &perMap = pdb.getPeriodicMap();
+    for (map<int, std::set<BoundaryType> >::iterator it = 
+	   perMap.periodicDofAssociations[feSpace].begin();
+	 it != perMap.periodicDofAssociations[feSpace].end(); ++it) {
       WorldVector<double> c;
       pdb.mesh->getDofIndexCoords(it->first, pdb.feSpaces[0], c);
       int nAssoc = it->second.size();
@@ -178,7 +180,7 @@ namespace AMDiS {
       for (int i = 1; i < pdb.mpiSize; i++)
 	stdMpi.recv(i);
     } else {
-      stdMpi.send(0, pdb.periodicDofMap[feSpace]);
+      stdMpi.send(0, perMap.periodicDofMap[feSpace]);
     }
 
     stdMpi.startCommunication();
@@ -190,7 +192,7 @@ namespace AMDiS {
     if (pdb.mpiRank == 0) {
       // Stores to each rank the periodic DOF mappings of this rank.
       map<int, PeriodicDofMap> rankToMaps;
-      PeriodicDofMap dofMap = pdb.periodicDofMap[feSpace];
+      PeriodicDofMap dofMap = perMap.periodicDofMap[feSpace];
       rankToMaps[0] = dofMap;
 
       for (int i = 1; i < pdb.mpiSize; i++) {
@@ -678,6 +680,8 @@ namespace AMDiS {
     ERROR_EXIT("Function must be rewritten!\n");
 
 #if 0
+    PeriodicMap &perMap = pdb.getPeriodicMap();
+
     const FiniteElemSpace* feSpace = pdb.feSpaces[0];
 
     typedef map<DegreeOfFreedom, DegreeOfFreedom> DofMapping;
@@ -686,8 +690,8 @@ namespace AMDiS {
     if (rank == -1 || pdb.mpiRank == rank) {
       cout << "====== DOF MAP PERIODIC ====== " << endl;
 
-      for (PeriodicDofMap::iterator it = pdb.periodicDofMap.begin();
-	   it != pdb.periodicDofMap.end(); ++it) {
+      for (PeriodicDofMap::iterator it = perMap.periodicDofMap.begin();
+	   it != perMap.periodicDofMap.end(); ++it) {
 	cout << "DOF MAP " << it->first << ": ";
 	for (std::set<DegreeOfFreedom>::iterator dofit = it->second.begin();
 	     dofit != it->second.end(); ++dofit)

AMDiS/src/parallel/ParallelTypes.h

@@ -62,12 +62,6 @@ namespace AMDiS {
   
   typedef map<const DegreeOfFreedom*, DegreeOfFreedom> DofIndexMap;
   
-  /// Mapps a boundar type, i.e., a boundary identifier index, to a periodic 
-  /// DOF mapping.
-  typedef map<BoundaryType, DofMapping> PeriodicDofMap;
-
-  typedef map<const FiniteElemSpace*, PeriodicDofMap> PeriodicDofMapFeSpace;
-  
   typedef vector<MeshStructure> MeshCodeVec;
 }
 

AMDiS/src/parallel/PeriodicMap.cc

+//
+// Software License for AMDiS
+//
+// Copyright (c) 2010 Dresden University of Technology 
+// All rights reserved.
+// Authors: Simon Vey, Thomas Witkowski et al.
+//
+// This file is part of AMDiS
+//
+// See also license.opensource.txt in the distribution.
+
+
+#include <fstream>
+#include "parallel/PeriodicMap.h"
+
+namespace AMDiS {
+
+  void PeriodicMap::serialize(ostream &out, 
+			      vector<const FiniteElemSpace*> feSpaces)
+  {
+    FUNCNAME("PeriodicMap::serialize()");
+    
+    int nFeSpace = static_cast<int>(feSpaces.size());
+    
+    for (int i = 0; i < nFeSpace; i++)
+      serialize(out, periodicDofMap[feSpaces[i]]);
+    
+    for (int i = 0; i < nFeSpace; i++)
+      serialize(out, periodicDofAssociations[feSpaces[i]]);
+  }
+
+
+  void PeriodicMap::deserialize(istream &in, 
+				vector<const FiniteElemSpace*> feSpaces)
+  {      
+    FUNCNAME("PeriodicMap::deserialize()");
+    
+    int nFeSpace = static_cast<int>(feSpaces.size());
+    
+    for (int i = 0; i < nFeSpace; i++)
+      deserialize(in, periodicDofMap[feSpaces[i]]);
+    
+    for (int i = 0; i < nFeSpace; i++)
+      deserialize(in, periodicDofAssociations[feSpaces[i]]);
+  }
+  
+  
+  void PeriodicMap::serialize(ostream &out, PeriodicDofMap &data)
+  {
+    int mapSize = data.size();
+    SerUtil::serialize(out, mapSize);
+    
+    for (PeriodicDofMap::iterator it = data.begin(); it != data.end(); ++it) {
+      int type = it->first;
+      DofMapping dofMap = it->second;
+      
+      SerUtil::serialize(out, type);
+      SerUtil::serialize(out, dofMap);
+    }
+  }
+  
+
+  void PeriodicMap::serialize(ostream &out, std::map<int, std::set<int> >& data)
+  {