diff --git a/AMDiS/src/BoundaryObject.h b/AMDiS/src/BoundaryObject.h
index 11d7ceda201b618e9dc5907d9c7e91f9d934e764..37f55b542cb4790484f52e773151cf3a40a07d21 100644
--- a/AMDiS/src/BoundaryObject.h
+++ b/AMDiS/src/BoundaryObject.h
@@ -104,7 +104,8 @@ namespace AMDiS {
    */
   struct AtomicBoundary {
     AtomicBoundary()
-      : type(INTERIOR)
+    : type(INTERIOR),
+      maxLevel(0)
     {}
 
     bool operator==(const AtomicBoundary& other) const;
@@ -119,6 +120,8 @@ namespace AMDiS {
     /// boundaries. Till now it is used only for periodic boundaries, which are
     /// also handles as interior boundaries.
     BoundaryType type;
+
+    int maxLevel;
   };
 
 
diff --git a/AMDiS/src/parallel/ElementObjectDatabase.cc b/AMDiS/src/parallel/ElementObjectDatabase.cc
index f169e8fe7a7ad7f95c434af5375c1968ae47fe00..1554531a4afb81af6a2efcf71aad6c200b245bee 100644
--- a/AMDiS/src/parallel/ElementObjectDatabase.cc
+++ b/AMDiS/src/parallel/ElementObjectDatabase.cc
@@ -388,6 +388,10 @@ namespace AMDiS {
     vertexInRank.clear();
     for (map<DegreeOfFreedom, vector<ElementObjectData> >::iterator it = vertexElements.begin();
 	 it != vertexElements.end(); ++it) {
+
+      vector<std::set<int> > ranksInLevel;
+      ranksInLevel.resize(nLevel);
+
       for (vector<ElementObjectData>::iterator it2 = it->second.begin(); 
 	   it2 != it->second.end(); ++it2) {
 	int elementInRank = macroElementRankMap[it2->elIndex];
@@ -400,7 +404,18 @@ namespace AMDiS {
 	  vertexOwner[it->first][level] = 
 	    std::max(vertexOwner[it->first][level], levelId);
 	}
+
+	ranksInLevel[0].insert(elementInRank);
+	for (int level = 1; level < nLevel; level++)
+	  ranksInLevel[level].insert(levelData.getLevelId(level, elementInRank));
       }
+
+      int maxLevel = 0;
+      for (int level = 1; level < nLevel; level++)
+	if (ranksInLevel[level].size() > 1)
+	  maxLevel = level;
+
+      vertexMaxLevel[it->first] = maxLevel;
     }
     
 
@@ -408,6 +423,10 @@ namespace AMDiS {
     edgeInRank.clear();
     for (map<DofEdge, vector<ElementObjectData> >::iterator it = edgeElements.begin();
 	 it != edgeElements.end(); ++it) {
+
+      vector<std::set<int> > ranksInLevel;
+      ranksInLevel.resize(nLevel);
+
       for (vector<ElementObjectData>::iterator it2 = it->second.begin(); 
 	   it2 != it->second.end(); ++it2) {
 	int elementInRank = macroElementRankMap[it2->elIndex];
@@ -420,7 +439,18 @@ namespace AMDiS {
 	  edgeOwner[it->first][level] = 
 	    std::max(edgeOwner[it->first][level], levelId);
 	}
+
+	ranksInLevel[0].insert(elementInRank);
+	for (int level = 1; level < nLevel; level++)
+	  ranksInLevel[level].insert(levelData.getLevelId(level, elementInRank));
       }
+
+      int maxLevel = 0;
+      for (int level = 1; level < nLevel; level++)
+	if (ranksInLevel[level].size() > 1)
+	  maxLevel = level;
+
+      edgeMaxLevel[it->first] = maxLevel;
     }
     
 
@@ -428,6 +458,10 @@ namespace AMDiS {
     faceInRank.clear();
     for (map<DofFace, vector<ElementObjectData> >::iterator it = faceElements.begin();
 	 it != faceElements.end(); ++it) {
+
+      vector<std::set<int> > ranksInLevel;
+      ranksInLevel.resize(nLevel);
+
       for (vector<ElementObjectData>::iterator it2 = it->second.begin(); 
 	   it2 != it->second.end(); ++it2) {
 	int elementInRank = macroElementRankMap[it2->elIndex];
@@ -440,7 +474,18 @@ namespace AMDiS {
 	  faceOwner[it->first][level] = 
 	    std::max(faceOwner[it->first][level], levelId);
 	}
+
+	ranksInLevel[0].insert(elementInRank);
+	for (int level = 1; level < nLevel; level++)
+	  ranksInLevel[level].insert(levelData.getLevelId(level, elementInRank));
       }
+
+      int maxLevel = 0;
+      for (int level = 1; level < nLevel; level++)
+	if (ranksInLevel[level].size() > 1)
+	  maxLevel = level;
+
+      faceMaxLevel[it->first] = maxLevel;
     }
   }
 
diff --git a/AMDiS/src/parallel/ElementObjectDatabase.h b/AMDiS/src/parallel/ElementObjectDatabase.h
index a7ecb5989c2a9c6d7f39b2955cde10dcd92dace7..96164cb1ec37146b824891a327f34e2ebf88f3b2 100644
--- a/AMDiS/src/parallel/ElementObjectDatabase.h
+++ b/AMDiS/src/parallel/ElementObjectDatabase.h
@@ -279,6 +279,26 @@ namespace AMDiS {
       }
     }
 
+    /// Returns the rank owner of the current iterator position.
+    int getIterateMaxLevel()
+    {
+      switch (iterGeoPos) {
+      case VERTEX:
+	return vertexMaxLevel[vertexIter->first];
+	break;
+      case EDGE:
+	return edgeMaxLevel[edgeIter->first];
+	break;
+      case FACE:
+	return faceMaxLevel[faceIter->first];
+	break;
+      default:
+	ERROR_EXIT("Should not happen!\n");
+
+	// Will never be reached, just to avoid compiler warnings.
+	return -1;
+      }
+    }
 
     /// Returns the rank owner of a vertex DOF.
     int getOwner(DegreeOfFreedom vertex, int level)
@@ -554,6 +574,13 @@ namespace AMDiS {
     /// Defines for all faces the rank that ownes this face.
     map<DofFace, LevelRank> faceOwner;
 
+    
+    map<DegreeOfFreedom, int> vertexMaxLevel;
+
+    map<DofEdge, int> edgeMaxLevel;
+
+    map<DofFace, int> faceMaxLevel;
+
 
     /// Defines to each vertex DOF a map that maps to each rank number the element
     /// objects that have this vertex DOF in common.
diff --git a/AMDiS/src/parallel/InteriorBoundary.cc b/AMDiS/src/parallel/InteriorBoundary.cc
index 62b13aa6faf0ab12a77a0b5d27d7a214cbfda4fa..3dc9816788d3b969658d7d52a011f705c1f93061 100644
--- a/AMDiS/src/parallel/InteriorBoundary.cc
+++ b/AMDiS/src/parallel/InteriorBoundary.cc
@@ -18,11 +18,11 @@
 
 namespace AMDiS {
 
-  AtomicBoundary& InteriorBoundary::getNewAtomic(int level, int rank)
+  AtomicBoundary& InteriorBoundary::getNewAtomic(int rank)
   {
-    int size = boundary[level][rank].size();
-    boundary[level][rank].resize(size + 1);
-    return boundary[level][rank][size];
+    int size = boundary[rank].size();
+    boundary[rank].resize(size + 1);
+    return boundary[rank][size];
   }
 
 
@@ -34,20 +34,20 @@ namespace AMDiS {
       return false;
 
     for (unsigned int level = 0; level < boundary.size(); level++) {
-      for (RankToBoundMap::const_iterator it = boundary[level].begin(); 
-	   it != boundary[level].end(); ++it) {
-	if (other2.boundary[level].count(it->first) == 0)
+      for (RankToBoundMap::const_iterator it = boundary.begin(); 
+	   it != boundary.end(); ++it) {
+	if (other2.boundary.count(it->first) == 0)
 	  return false;     
 	
-	if (other2.boundary[level][it->first].size() != it->second.size())
+	if (other2.boundary[it->first].size() != it->second.size())
 	  return false;
 	
 	for (unsigned int i = 0; i < it->second.size(); i++) {
 	  std::vector<AtomicBoundary>::iterator bIt = 
-	    find(other2.boundary[level][it->first].begin(),
-		 other2.boundary[level][it->first].end(), it->second[i]);
+	    find(other2.boundary[it->first].begin(),
+		 other2.boundary[it->first].end(), it->second[i]);
 	  
-	  if (bIt == other2.boundary[level][it->first].end())
+	  if (bIt == other2.boundary[it->first].end())
 	    return false;	
 	}      
       }
@@ -57,14 +57,11 @@ namespace AMDiS {
   }
 
 
-  void InteriorBoundary::reset(int level)
+  void InteriorBoundary::clear()
   { 
-    FUNCNAME("InteriorBoundary::reset()");
-
-    nLevel = level;
+    FUNCNAME("InteriorBoundary::clear()");
 
     boundary.clear();
-    boundary.resize(nLevel);
   }
 
 
diff --git a/AMDiS/src/parallel/InteriorBoundary.h b/AMDiS/src/parallel/InteriorBoundary.h
index bdc02b773648b55be2e199ff5ac2d5075512752a..ab830e41ed80d4ac2e0272ca3f4c8bc9d7e7c6ac 100644
--- a/AMDiS/src/parallel/InteriorBoundary.h
+++ b/AMDiS/src/parallel/InteriorBoundary.h
@@ -45,49 +45,36 @@ namespace AMDiS {
     /// Iterator for the interior boundary object.
     class iterator {      
     public:
-      iterator(InteriorBoundary &b)
+      iterator(InteriorBoundary &b, int traverseLevel = 0)
 	: bound(b),
-	  level(0)
+	  level(traverseLevel)
       {
 	reset();
       }
 
-      iterator(InteriorBoundary &b, int level)
-	: bound(b),
-	  level(level)
-      {
-	TEST_EXIT_DBG(level < bound.boundary.size())
-	  ("Should not happen!\n");
-
-	reset();
-      }
-
       /// Set the iterator to the first position.
       void reset()
       {
-	mapIt = bound.boundary[level].begin();
+	mapIt = bound.boundary.begin();
 	nextNonempty();
-
-	if (mapIt != bound.boundary[level].end())
-	  vecIt = mapIt->second.begin();	
       }
 
       /// Test if iterator is at the final position.
       bool end() const
       {
-	return (mapIt == bound.boundary[level].end());
+	return (mapIt == bound.boundary.end());
       }
 
       /// Move iterator to the next position.
       void operator++()
       {
-	++vecIt;
+	do {
+	  ++vecIt;
+	} while (vecIt->maxLevel < level && vecIt != mapIt->second.end());
+
 	if (vecIt == mapIt->second.end()) {
 	  ++mapIt;
 	  nextNonempty();
-
-	  if (mapIt != bound.boundary[level].end())
-	    vecIt = mapIt->second.begin();	  
 	}	
       }
 
@@ -105,9 +92,6 @@ namespace AMDiS {
       {
 	++mapIt;
 	nextNonempty();
-
-	if (mapIt != bound.boundary[level].end())
-	  vecIt = mapIt->second.begin();	
       }
 
       inline int getRank() 
@@ -119,14 +103,33 @@ namespace AMDiS {
 
       inline void nextNonempty()
       {
-	if (mapIt == bound.boundary[level].end())
-	  return;
-
-	while (mapIt->second.size() == 0) {
+	do {
+	  // Return, we are at the end.
+	  if (mapIt == bound.boundary.end())
+	    return;
+	  
+	  // Search for the next non empty boundary map.
+	  while (mapIt->second.size() == 0) {
+	    ++mapIt;
+	    if (mapIt == bound.boundary.end())
+	      return;	  
+	  }
+	  
+	  vecIt = mapIt->second.begin();
+	  
+	  // Search for the next atomic boundary on the mesh level
+	  while (vecIt->maxLevel < level && vecIt != mapIt->second.end())
+	    ++vecIt;
+	  
+	  // If vector iterator is not at the end, we have found one and
+	  // can return.
+	  if (vecIt != mapIt->second.end())
+	    return;
+	  
+	  // In this case, no boundary on the given level is found, continue
+	  // with next rank.
 	  ++mapIt;
-	  if (mapIt == bound.boundary[level].end())
-	    return;	  
-	}
+	} while (true);
       }
 
     protected:
@@ -140,15 +143,9 @@ namespace AMDiS {
     };
 
   public:
-    InteriorBoundary(int l = 1) 
-      : nLevel(l)
-    {
-      boundary.resize(nLevel);
-    }
+    void clear();
 
-    void reset(int nLevel);
-
-    AtomicBoundary& getNewAtomic(int level, int rank);
+    AtomicBoundary& getNewAtomic(int rank);
 
     /// Writes this object to a file.
     void serialize(ostream &out);
@@ -166,10 +163,7 @@ namespace AMDiS {
     void deserializeExcludeList(istream &in, ExcludeList &list);    
 
   public:
-    vector<RankToBoundMap> boundary;    
-
-  protected:
-    int nLevel;
+    RankToBoundMap boundary;    
   };
 }
 
diff --git a/AMDiS/src/parallel/MeshDistributor.cc b/AMDiS/src/parallel/MeshDistributor.cc
index ee86021131c18282f24e499f0bc0f8425b3eae86..ea6040ea03267eed17dbbbe324883bc9ce4a936f 100644
--- a/AMDiS/src/parallel/MeshDistributor.cc
+++ b/AMDiS/src/parallel/MeshDistributor.cc
@@ -1576,19 +1576,9 @@ namespace AMDiS {
 
     // === Clear all relevant data structures. ===
 
-    int nLevel = levelData.getLevelNumber();
-    rankIntBoundary.reset(nLevel);
-    otherIntBoundary.reset(nLevel);
-    periodicBoundary.reset(nLevel);
-
-    for (int level = 0; level < nLevel; level++)
-      createBoundaryData(level);  
-  }
-
-
-  void MeshDistributor::createBoundaryData(int level)
-  {
-    FUNCNAME("MeshDistributor::createBoundaryData()");
+    rankIntBoundary.clear();
+    otherIntBoundary.clear();
+    periodicBoundary.clear();
 
     // === Create interior boundary data structure. ===
 
@@ -1600,10 +1590,12 @@ namespace AMDiS {
 	if (!(objData.count(mpiRank) && objData.size() > 1))
 	  continue;
 
-	int owner = elObjDb.getIterateOwner(level);
+	int owner = elObjDb.getIterateOwner(0);
 	ElementObjectData& rankBoundEl = objData[mpiRank];
 	
-	AtomicBoundary bound;	    	    
+	AtomicBoundary bound;
+	bound.maxLevel = elObjDb.getIterateMaxLevel();
+
 	bound.rankObj.el = macroElIndexMap[rankBoundEl.elIndex];
 	bound.rankObj.elIndex = rankBoundEl.elIndex;
 	bound.rankObj.elType = macroElIndexTypeMap[rankBoundEl.elIndex];
@@ -1634,7 +1626,7 @@ namespace AMDiS {
 	    
 	    bound.type = INTERIOR;
 	    
-	    AtomicBoundary& b = rankIntBoundary.getNewAtomic(level, it2->first);
+	    AtomicBoundary& b = rankIntBoundary.getNewAtomic(it2->first);
 	    b = bound;
 	    if (geoIndex == EDGE)
 	      b.neighObj.reverseMode = 
@@ -1658,7 +1650,7 @@ namespace AMDiS {
 	  
 	  bound.type = INTERIOR;
 	  
-	  AtomicBoundary& b = otherIntBoundary.getNewAtomic(level, owner);
+	  AtomicBoundary& b = otherIntBoundary.getNewAtomic(owner);
 	  b = bound;	    
 	  if (geoIndex == EDGE)
 	    b.rankObj.reverseMode =
@@ -1707,11 +1699,11 @@ namespace AMDiS {
 	      bound.rankObj.ithObj == 1 ||
 	      bound.rankObj.elIndex == 78 &&
 	      bound.rankObj.ithObj == 2) {
-	    AtomicBoundary& b = periodicBoundary.getNewAtomic(level, otherElementRank);
+	    AtomicBoundary& b = periodicBoundary.getNewAtomic(otherElementRank);
 	    b = bound;
 	  }
 	} else {
-	  AtomicBoundary& b = periodicBoundary.getNewAtomic(level, otherElementRank);
+	  AtomicBoundary& b = periodicBoundary.getNewAtomic(otherElementRank);
 	  b = bound;	    
 	}
       }
@@ -1746,7 +1738,7 @@ namespace AMDiS {
 	
 	bound.type = it->second;
 	
-	AtomicBoundary& b = periodicBoundary.getNewAtomic(level, otherElementRank);
+	AtomicBoundary& b = periodicBoundary.getNewAtomic(otherElementRank);
 	b = bound;
      
 	if (mpiRank > otherElementRank)
@@ -1792,7 +1784,7 @@ namespace AMDiS {
 	
 	bound.type = it->second;
 	
-	AtomicBoundary& b = periodicBoundary.getNewAtomic(level, otherElementRank);
+	AtomicBoundary& b = periodicBoundary.getNewAtomic(otherElementRank);
 	b = bound;
      
 	if (mpiRank > otherElementRank)
@@ -1811,8 +1803,8 @@ namespace AMDiS {
     // === share the bounday.                                                 ===
 
     StdMpi<vector<AtomicBoundary> > stdMpi(mpiComm);
-    stdMpi.send(rankIntBoundary.boundary[level]);
-    stdMpi.recv(otherIntBoundary.boundary[level]);
+    stdMpi.send(rankIntBoundary.boundary);
+    stdMpi.recv(otherIntBoundary.boundary);
     stdMpi.startCommunication();
 
 
@@ -1821,8 +1813,8 @@ namespace AMDiS {
     // === the same order. If not, the atomic boundaries are swaped to the    ===
     // === correct order.                                                     ===
 
-    for (RankToBoundMap::iterator rankIt = otherIntBoundary.boundary[level].begin();
-	 rankIt != otherIntBoundary.boundary[level].end(); ++rankIt) {
+    for (RankToBoundMap::iterator rankIt = otherIntBoundary.boundary.begin();
+	 rankIt != otherIntBoundary.boundary.end(); ++rankIt) {
 
       // === We have received from rank "rankIt->first" the ordered list of   ===
       // === element indices. Now, we have to sort the corresponding list in  ===
@@ -1856,29 +1848,28 @@ namespace AMDiS {
 
     // === Do the same for the periodic boundaries. ===
 
-    if (periodicBoundary.boundary[level].size() > 0) {
+    if (periodicBoundary.boundary.size() > 0) {
       stdMpi.clear();
 
-      InteriorBoundary sendBounds(levelData.getLevelNumber());
-      InteriorBoundary recvBounds(levelData.getLevelNumber());
-      for (RankToBoundMap::iterator rankIt = periodicBoundary.boundary[level].begin();
-	   rankIt != periodicBoundary.boundary[level].end(); ++rankIt) {
+      InteriorBoundary sendBounds, recvBounds;
+      for (RankToBoundMap::iterator rankIt = periodicBoundary.boundary.begin();
+	   rankIt != periodicBoundary.boundary.end(); ++rankIt) {
 
 	if (rankIt->first == mpiRank)
 	  continue;
 
 	if (rankIt->first < mpiRank)
-	  sendBounds.boundary[level][rankIt->first] = rankIt->second;
+	  sendBounds.boundary[rankIt->first] = rankIt->second;
 	else
-	  recvBounds.boundary[level][rankIt->first] = rankIt->second;	
+	  recvBounds.boundary[rankIt->first] = rankIt->second;	
       }
 
-      stdMpi.send(sendBounds.boundary[level]);
-      stdMpi.recv(recvBounds.boundary[level]);
+      stdMpi.send(sendBounds.boundary);
+      stdMpi.recv(recvBounds.boundary);
       stdMpi.startCommunication();
 
-      for (RankToBoundMap::iterator rankIt = periodicBoundary.boundary[level].begin();
-	   rankIt != periodicBoundary.boundary[level].end(); ++rankIt) {
+      for (RankToBoundMap::iterator rankIt = periodicBoundary.boundary.begin();
+	   rankIt != periodicBoundary.boundary.end(); ++rankIt) {
 
  	if (rankIt->first <= mpiRank)
  	  continue;
@@ -1889,12 +1880,12 @@ namespace AMDiS {
 	  BoundaryObject &recvNeighObj = 
 	    stdMpi.getRecvData()[rankIt->first][j].neighObj;
 
-	  if (periodicBoundary.boundary[level][rankIt->first][j].neighObj != recvRankObj ||
-	      periodicBoundary.boundary[level][rankIt->first][j].rankObj != recvNeighObj) {
+	  if (periodicBoundary.boundary[rankIt->first][j].neighObj != recvRankObj ||
+	      periodicBoundary.boundary[rankIt->first][j].rankObj != recvNeighObj) {
 	    unsigned int k = j + 1;	    
 	    for (; k < rankIt->second.size(); k++)
-	      if (periodicBoundary.boundary[level][rankIt->first][k].neighObj == recvRankObj &&
-		  periodicBoundary.boundary[level][rankIt->first][k].rankObj == recvNeighObj)
+	      if (periodicBoundary.boundary[rankIt->first][k].neighObj == recvRankObj &&
+		  periodicBoundary.boundary[rankIt->first][k].rankObj == recvNeighObj)
 		break;
 	    
 	    // The element must always be found, because the list is just in 
@@ -1908,7 +1899,7 @@ namespace AMDiS {
 	  } 
 	}
       }     
-    } // periodicBoundary.boundary[level].size() > 0
+    } // periodicBoundary.boundary.size() > 0
   }
 
 
@@ -2112,8 +2103,6 @@ namespace AMDiS {
   {
     FUNCNAME("MeshDistributor::createPeriodicMap()");
 
-    TEST_EXIT(levelData.getLevelNumber() == 1)("Not yet implemented for multilevel!\n");
-
     // Clear all periodic DOF mappings calculated before. We do it from scratch.
     periodicDofs.init(levelData.getLevelNumber());
     periodicMap.clear();
@@ -2122,7 +2111,7 @@ namespace AMDiS {
     // periodicMap must be still cleared before: if we do repartitioning and
     // there were periodic boundaries in subdomain before and after repartitioning
     // there are no more periodic boundaries.
-    if (periodicBoundary.boundary[0].size() == 0)
+    if (periodicBoundary.boundary.size() == 0)
       return;
 
     for (unsigned int i = 0; i < feSpaces.size(); i++)
@@ -2134,8 +2123,6 @@ namespace AMDiS {
   {
     FUNCNAME("MeshDistributor::createPeriodicMap()");
 
-    TEST_EXIT(levelData.getLevelNumber() == 1)("No yet implemented for multilevel stuff!\n");
-
     StdMpi<vector<int> > stdMpi(mpiComm, false);
 
     // === Each rank traverse its periodic boundaries and sends the DOF      ===
@@ -2143,8 +2130,8 @@ namespace AMDiS {
 
     map<int, vector<int> > rankToDofType;
 
-    for (RankToBoundMap::iterator it = periodicBoundary.boundary[0].begin();
-	 it != periodicBoundary.boundary[0].end(); ++it) {
+    for (RankToBoundMap::iterator it = periodicBoundary.boundary.begin();
+	 it != periodicBoundary.boundary.end(); ++it) {
 
       if (it->first == mpiRank) {
 	// Here we have a periodic boundary within rank's subdomain. So we can
@@ -2212,8 +2199,8 @@ namespace AMDiS {
     // === DOFs from the other ranks.                                        ===
 
 
-    for (RankToBoundMap::iterator it = periodicBoundary.boundary[0].begin();
-	 it != periodicBoundary.boundary[0].end(); ++it) {
+    for (RankToBoundMap::iterator it = periodicBoundary.boundary.begin();
+	 it != periodicBoundary.boundary.end(); ++it) {
       DofContainer& dofs = periodicDofs.getDofContainer(it->first, feSpace);
       vector<int>& types = rankToDofType[it->first];
 
@@ -2236,8 +2223,8 @@ namespace AMDiS {
     StdMpi<PeriodicDofMap> stdMpi2(mpiComm);
 
 
-    for (RankToBoundMap::iterator it = periodicBoundary.boundary[0].begin();
-	 it != periodicBoundary.boundary[0].end(); ++it) {
+    for (RankToBoundMap::iterator it = periodicBoundary.boundary.begin();
+	 it != periodicBoundary.boundary.end(); ++it) {
       if (it->first == mpiRank)
 	continue;
 
diff --git a/AMDiS/src/parallel/MeshDistributor.h b/AMDiS/src/parallel/MeshDistributor.h
index 0f584268b8f3a2bcde9776bd7f4d4ab4ad25d140..7aeb3037deeabb38bf54f5063e4a8ac7ebb4b369 100644
--- a/AMDiS/src/parallel/MeshDistributor.h
+++ b/AMDiS/src/parallel/MeshDistributor.h
@@ -297,8 +297,6 @@ namespace AMDiS {
 
     void createBoundaryData();
 
-    void createBoundaryData(int level);
-
     void createBoundaryDofs();
 
     void createBoundaryDofs(const FiniteElemSpace *feSpace, int level);
diff --git a/AMDiS/src/parallel/ParallelDebug.cc b/AMDiS/src/parallel/ParallelDebug.cc
index fbcac563997bff65bdae71cf7a25ca07fcdd2bfb..2cef40f53eef19c89f3006feffba1059f23ca31a 100644
--- a/AMDiS/src/parallel/ParallelDebug.cc
+++ b/AMDiS/src/parallel/ParallelDebug.cc
@@ -39,8 +39,8 @@ namespace AMDiS {
 
     // === Send rank's boundary information. ===
 
-    for (RankToBoundMap::iterator rankIt = pdb.rankIntBoundary.boundary[0].begin();
-	 rankIt != pdb.rankIntBoundary.boundary[0].end(); ++rankIt) {
+    for (RankToBoundMap::iterator rankIt = pdb.rankIntBoundary.boundary.begin();
+	 rankIt != pdb.rankIntBoundary.boundary.end(); ++rankIt) {
 
       int nSendInt = rankIt->second.size();
       int* buffer = new int[nSendInt];
@@ -56,8 +56,8 @@ namespace AMDiS {
 
     // === Receive information from other ranks about the interior boundaries. ====
 
-    for (RankToBoundMap::iterator rankIt = pdb.otherIntBoundary.boundary[0].begin();
-	 rankIt != pdb.otherIntBoundary.boundary[0].end(); ++rankIt) {
+    for (RankToBoundMap::iterator rankIt = pdb.otherIntBoundary.boundary.begin();
+	 rankIt != pdb.otherIntBoundary.boundary.end(); ++rankIt) {
       int nRecvInt = rankIt->second.size();
       int *buffer = new int[nRecvInt];
       recvBuffers.push_back(buffer);
@@ -69,8 +69,8 @@ namespace AMDiS {
 
     // === To the last, do the same of periodic boundaries. ===
 
-    for (RankToBoundMap::iterator rankIt = pdb.periodicBoundary.boundary[0].begin();       
-	 rankIt != pdb.periodicBoundary.boundary[0].end(); ++rankIt) {
+    for (RankToBoundMap::iterator rankIt = pdb.periodicBoundary.boundary.begin();       
+	 rankIt != pdb.periodicBoundary.boundary.end(); ++rankIt) {
       if (rankIt->first == pdb.mpiRank)
 	continue;
 
@@ -103,8 +103,8 @@ namespace AMDiS {
     // === and after this the periodic ones.                                   ===
 
     int bufCounter = 0;
-    for (RankToBoundMap::iterator rankIt = pdb.otherIntBoundary.boundary[0].begin();
-	 rankIt != pdb.otherIntBoundary.boundary[0].end(); ++rankIt) {
+    for (RankToBoundMap::iterator rankIt = pdb.otherIntBoundary.boundary.begin();
+	 rankIt != pdb.otherIntBoundary.boundary.end(); ++rankIt) {
 
       TEST_EXIT(rankIt->second.size() == 
 		pdb.otherIntBoundary.boundary[rankIt->first].size())
@@ -112,7 +112,7 @@ namespace AMDiS {
 
       for (unsigned int i = 0; i < rankIt->second.size(); i++) {
 	int elIndex1 = recvBuffers[bufCounter][i];
-	int elIndex2 = pdb.otherIntBoundary.boundary[0][rankIt->first][i].neighObj.elIndex;
+	int elIndex2 = pdb.otherIntBoundary.boundary[rankIt->first][i].neighObj.elIndex;
 
 	TEST_EXIT(elIndex1 == elIndex2)("Wrong element index at interior boundary!\n");
       }
@@ -121,18 +121,18 @@ namespace AMDiS {
     }
 
 
-    for (RankToBoundMap::iterator rankIt = pdb.periodicBoundary.boundary[0].begin();
-	 rankIt != pdb.periodicBoundary.boundary[0].end(); ++rankIt) {
+    for (RankToBoundMap::iterator rankIt = pdb.periodicBoundary.boundary.begin();
+	 rankIt != pdb.periodicBoundary.boundary.end(); ++rankIt) {
       if (rankIt->first == pdb.mpiRank)
 	continue;
 
       for (unsigned int i = 0; i < rankIt->second.size(); i++) {
 	int elIndex1 = recvBuffers[bufCounter][i];
-	int elIndex2 = pdb.periodicBoundary.boundary[0][rankIt->first][i].neighObj.elIndex;
+	int elIndex2 = pdb.periodicBoundary.boundary[rankIt->first][i].neighObj.elIndex;
 
 	TEST_EXIT(elIndex1 == elIndex2)
 	  ("Wrong element index at periodic boundary el %d with rank %d: %d %d\n", 
-	   pdb.periodicBoundary.boundary[0][rankIt->first][i].rankObj.elIndex,
+	   pdb.periodicBoundary.boundary[rankIt->first][i].rankObj.elIndex,
 	   rankIt->first, elIndex1, elIndex2);
       }
 
@@ -256,8 +256,8 @@ namespace AMDiS {
     RankToCoords sendCoords;
     map<int, vector<BoundaryType> > rankToDofType;
 
-    for (InteriorBoundary::RankToBoundMap::iterator it = pdb.periodicBoundary.boundary[0].begin();
-	 it != pdb.periodicBoundary.boundary[0].end(); ++it) {
+    for (InteriorBoundary::RankToBoundMap::iterator it = pdb.periodicBoundary.boundary.begin();
+	 it != pdb.periodicBoundary.boundary.end(); ++it) {
       if (it->first == pdb.mpiRank)
 	continue;
 
diff --git a/AMDiS/src/parallel/ParallelDofMapping.cc b/AMDiS/src/parallel/ParallelDofMapping.cc
index d981fee3782399be4c94a02d898b9e9282e7c97e..6690c6044e8d5fc972dc9502a6700722929cbbf9 100644
--- a/AMDiS/src/parallel/ParallelDofMapping.cc
+++ b/AMDiS/src/parallel/ParallelDofMapping.cc
@@ -52,7 +52,7 @@ namespace AMDiS {
       
       nOverallDofs[i] = 0;
       rStartDofs[i] = 0;
-      mpi::getDofNumbering(levelData->getMpiComm(i), 
+      mpi::getDofNumbering(levelData->getMpiComm(0), 
 			   nRankDofs[i], rStartDofs[i], nOverallDofs[i]);
       
       // === If required, compute also the global indices. ===
@@ -86,7 +86,6 @@ namespace AMDiS {
     // === Send all global indices of DOFs that are owned by the rank to all ===
     // === other ranks that also include this DOF.                           ===
 
-    //    StdMpi<vector<int> > stdMpi(levelData->getMpiComm(level));
     StdMpi<vector<int> > stdMpi(levelData->getMpiComm(0));
 
     for (DofComm::Iterator it(*sendDofs, level, feSpace); 
@@ -350,7 +349,6 @@ namespace AMDiS {
       // === Communicate the matrix indices for all DOFs that are on some ===
       // === interior boundaries.                                         ===
 
-      //      StdMpi<vector<DegreeOfFreedom> > stdMpi(levelData->getMpiComm(level));
       StdMpi<vector<DegreeOfFreedom> > stdMpi(levelData->getMpiComm(0));
       for (DofComm::Iterator it(*sendDofs, level, feSpaces[i]); 
 	   !it.end(); it.nextRank()) {
diff --git a/AMDiS/src/parallel/PetscSolverFeti.cc b/AMDiS/src/parallel/PetscSolverFeti.cc
index 9705970ce56d52a085a33033d59f0e00fb331b82..aad8af990c9976db2200de464b5ce23e656a9bcc 100644
--- a/AMDiS/src/parallel/PetscSolverFeti.cc
+++ b/AMDiS/src/parallel/PetscSolverFeti.cc
@@ -308,15 +308,6 @@ namespace AMDiS {
     if (fetiPreconditioner != FETI_NONE)
       interiorDofMap.update();
 
-    for (unsigned int i = 0; i < meshDistributor->getFeSpaces().size(); i++) {
-      const FiniteElemSpace *feSpace = meshDistributor->getFeSpace(i);
-
-      MSG("TEST DUALS %d %d %d\n",
-	  dualDofMap[feSpace].nLocalDofs[meshLevel], 
-	  dualDofMap[feSpace].nRankDofs[meshLevel], 
-	  dualDofMap[feSpace].nOverallDofs[meshLevel]);
-    }	
-
     for (unsigned int i = 0; i < meshDistributor->getFeSpaces().size(); i++) {
       const FiniteElemSpace *feSpace = meshDistributor->getFeSpace(i);
       createLagrange(feSpace);
@@ -361,9 +352,10 @@ namespace AMDiS {
 
     /// Set of DOF indices that are considered to be primal variables.
     
-    DofIndexSet primals;
     DofContainerSet& vertices = 
       meshDistributor->getBoundaryDofInfo(feSpace, meshLevel).geoDofs[VERTEX];
+
+    DofIndexSet primals;
     for (DofContainerSet::iterator it = vertices.begin(); 
 	 it != vertices.end(); ++it)
       primals.insert(**it);
diff --git a/AMDiS/src/parallel/PetscSolverFeti.h b/AMDiS/src/parallel/PetscSolverFeti.h
index e1dfc73121f82296fdee7560fe25846006ea2d02..a33ac6f9a5f65052569ff3f742e14203f812db04 100644
--- a/AMDiS/src/parallel/PetscSolverFeti.h
+++ b/AMDiS/src/parallel/PetscSolverFeti.h
@@ -80,22 +80,22 @@ namespace AMDiS {
 
     int getNumberOfPrimals()
     {
-      return primalDofMap.getOverallDofs(0);
+      return primalDofMap.getOverallDofs(meshLevel);
     }
 
     int getNumberOfRankPrimals()
     {
-      return primalDofMap.getRankDofs(0);
+      return primalDofMap.getRankDofs(meshLevel);
     }
 
     int getNumberOfDuals()
     {
-      return dualDofMap.getOverallDofs(0);
+      return dualDofMap.getOverallDofs(meshLevel);
     }
 
     int getNumberOfRankDuals()
     {
-      return dualDofMap.getRankDofs(0);
+      return dualDofMap.getRankDofs(meshLevel);
     }
 
   protected:
diff --git a/test/mpi/src/test0001.cc b/test/mpi/src/test0001.cc
index ad599746d862be0bc943d87b41d2d785efbb3117..b833be448cd8f84984077944e87a158827771c2a 100644
--- a/test/mpi/src/test0001.cc
+++ b/test/mpi/src/test0001.cc
@@ -45,58 +45,64 @@ BOOST_AUTO_TEST_CASE(amdis_mpi_simple_partitioning)
     neighTest.insert(1); neighTest.insert(4); neighTest.insert(5);
     break;
   case 1:
-    neighTest.insert(0); neighTest.insert(2); neighTest.insert(5); neighTest.insert(6);
+    neighTest.insert(0); neighTest.insert(2); neighTest.insert(4); neighTest.insert(5); neighTest.insert(6);
     break;
   case 2:
-    neighTest.insert(1); neighTest.insert(3); neighTest.insert(6); neighTest.insert(7);
+    neighTest.insert(1); neighTest.insert(3); neighTest.insert(5); neighTest.insert(6); neighTest.insert(7);
     break;
   case 3:
-    neighTest.insert(2); neighTest.insert(7);
+    neighTest.insert(2); neighTest.insert(6); neighTest.insert(7);
     break;
   case 4:
-    neighTest.insert(0); neighTest.insert(5); neighTest.insert(8); neighTest.insert(9);
+    neighTest.insert(0); neighTest.insert(1); neighTest.insert(5); neighTest.insert(8); neighTest.insert(9);
     break;
   case 5:
-    neighTest.insert(0); neighTest.insert(1); neighTest.insert(4); neighTest.insert(6);
-    neighTest.insert(9); neighTest.insert(10);
+    neighTest.insert(0); neighTest.insert(1); neighTest.insert(2); 
+    neighTest.insert(4); neighTest.insert(6);
+    neighTest.insert(8); neighTest.insert(9); neighTest.insert(10);
     break;
   case 6:
-    neighTest.insert(1); neighTest.insert(2); neighTest.insert(5); neighTest.insert(7);
-    neighTest.insert(10); neighTest.insert(11);
+    neighTest.insert(1); neighTest.insert(2); neighTest.insert(3);
+    neighTest.insert(5); neighTest.insert(7);
+    neighTest.insert(9); neighTest.insert(10); neighTest.insert(11);
     break;
   case 7:
-    neighTest.insert(2); neighTest.insert(3); neighTest.insert(6); neighTest.insert(11);
+    neighTest.insert(2); neighTest.insert(3); neighTest.insert(6); neighTest.insert(10); neighTest.insert(11);
     break;
   case 8:
-    neighTest.insert(4); neighTest.insert(9); neighTest.insert(12); neighTest.insert(13);
+    neighTest.insert(4); neighTest.insert(5); neighTest.insert(9); neighTest.insert(12); neighTest.insert(13);
     break;
   case 9:
-    neighTest.insert(4); neighTest.insert(5); neighTest.insert(8); neighTest.insert(10);
-    neighTest.insert(13); neighTest.insert(14);
+    neighTest.insert(4); neighTest.insert(5); neighTest.insert(6); 
+    neighTest.insert(8); neighTest.insert(10);
+    neighTest.insert(12); neighTest.insert(13); neighTest.insert(14);
     break;
   case 10:
-    neighTest.insert(5); neighTest.insert(6); neighTest.insert(9); neighTest.insert(11);
-    neighTest.insert(14); neighTest.insert(15);
+    neighTest.insert(5); neighTest.insert(6); neighTest.insert(7);
+    neighTest.insert(9); neighTest.insert(11);
+    neighTest.insert(13); neighTest.insert(14); neighTest.insert(15);
     break;
   case 11:
-    neighTest.insert(6); neighTest.insert(7); neighTest.insert(10); neighTest.insert(15);
+    neighTest.insert(6); neighTest.insert(7); neighTest.insert(10); neighTest.insert(14); neighTest.insert(15);
     break;
   case 12:
-    neighTest.insert(8); neighTest.insert(13);
+    neighTest.insert(8); neighTest.insert(9); neighTest.insert(13);
     break;
   case 13:
-    neighTest.insert(8); neighTest.insert(9); neighTest.insert(12); neighTest.insert(14);
+    neighTest.insert(8); neighTest.insert(9); neighTest.insert(10); neighTest.insert(12); neighTest.insert(14);
     break;
   case 14:
-    neighTest.insert(9); neighTest.insert(10); neighTest.insert(13); neighTest.insert(15);
+    neighTest.insert(9); neighTest.insert(10); neighTest.insert(11); neighTest.insert(13); neighTest.insert(15);
     break;
   case 15:
     neighTest.insert(10); neighTest.insert(11); neighTest.insert(14);
     break;
   } 
 
-  BOOST_REQUIRE(levelNeighbours == neighTest);
-
+  if (levelNeighbours != neighTest) {
+    cout << "RANK " << MPI::COMM_WORLD.Get_rank() << ": " << levelNeighbours.size() << " " << neighTest.size() << "\n";
+    BOOST_ERROR("Neighbours are not set correctly on level 0!\n");
+  }
 
 
   levelRanks = md.getLevelRanks(1);
@@ -116,7 +122,9 @@ BOOST_AUTO_TEST_CASE(amdis_mpi_simple_partitioning)
     break;
   }
 
-  BOOST_REQUIRE(levelRanks == rankTest);
+  if (levelRanks != rankTest) {
+    BOOST_ERROR("Wrong ranks in level 1 subdomain!\n");
+  }
 
   levelNeighbours = md.getLevelNeighbours(1);
   neighTest.clear();
@@ -126,41 +134,43 @@ BOOST_AUTO_TEST_CASE(amdis_mpi_simple_partitioning)
     neighTest.insert(2); neighTest.insert(6);
     break;
   case 2:
-    neighTest.insert(1); 
+    neighTest.insert(1);  neighTest.insert(5); 
     break;
   case 4:
     neighTest.insert(8); neighTest.insert(9);
     break;
   case 5:
-    neighTest.insert(6); neighTest.insert(9); neighTest.insert(10);
+    neighTest.insert(2); neighTest.insert(6); neighTest.insert(8); neighTest.insert(9); neighTest.insert(10);
     break;
   case 6:
-    neighTest.insert(1); neighTest.insert(5); neighTest.insert(10); neighTest.insert(11);
+    neighTest.insert(1); neighTest.insert(5); neighTest.insert(9); neighTest.insert(10); neighTest.insert(11);
     break;
   case 7:
-    neighTest.insert(11);
+    neighTest.insert(10); neighTest.insert(11);
     break;
   case 8:
-    neighTest.insert(4);
+    neighTest.insert(4); neighTest.insert(5);
     break;
   case 9:
-    neighTest.insert(4); neighTest.insert(5); neighTest.insert(10); neighTest.insert(14);
+    neighTest.insert(4); neighTest.insert(5); neighTest.insert(6); neighTest.insert(10); neighTest.insert(14);
     break;
   case 10:
-    neighTest.insert(5); neighTest.insert(6); neighTest.insert(9); 
+    neighTest.insert(5); neighTest.insert(6); neighTest.insert(7); neighTest.insert(9); neighTest.insert(13);
     break;
   case 11:
     neighTest.insert(6); neighTest.insert(7);
     break;
   case 13:
-    neighTest.insert(14);
+    neighTest.insert(10); neighTest.insert(14);
     break;
   case 14:
     neighTest.insert(9); neighTest.insert(13);
     break;
   }
 
-  BOOST_REQUIRE(levelNeighbours == neighTest);
+  if (levelNeighbours != neighTest) {
+    BOOST_ERROR("Neighbours are not set correctly on level 1!\n");
+  }
 }
 
 
diff --git a/test/mpi/src/test0003.cc b/test/mpi/src/test0003.cc
index 64845688f969945cef220b1d7ccf32bbaaa60bf8..2706d109d873d61d01da6f27ba1c66e111d9210e 100644
--- a/test/mpi/src/test0003.cc
+++ b/test/mpi/src/test0003.cc
@@ -1,12 +1,14 @@
 #define BOOST_TEST_DYN_LINK
-#define BOOST_TEST_MODULE 0002
+#define BOOST_TEST_MODULE 0003
 #define BOOST_TEST_NO_MAIN
 
 #include <boost/test/unit_test.hpp>
+#include <boost/lexical_cast.hpp>
 #include <AMDiS.h>
 
 using namespace AMDiS;
 using namespace std;
+using boost::lexical_cast;
 
 /*
   - 
@@ -34,8 +36,13 @@ BOOST_AUTO_TEST_CASE(amdis_mpi_feti_multilevel)
   feti.initialize(feSpaces);
   feti.createFetiData();
 
-  BOOST_REQUIRE(feti.getNumberOfPrimals() == 9);
-  BOOST_REQUIRE(feti.getNumberOfDuals() == 16);
+  if (feti.getNumberOfPrimals() != 9) {
+    BOOST_ERROR("Number of primals is " + lexical_cast<string>(feti.getNumberOfPrimals()) + " and not 9!\n");
+  }
+
+  if (feti.getNumberOfDuals() != 16) {
+    BOOST_ERROR("Number of duals is " + lexical_cast<string>(feti.getNumberOfDuals()) + " and not 16!\n");
+  }
 }