diff --git a/apps/gpu/DrivenCavityMultiGPU/DrivenCavityMultiGPU.cpp b/apps/gpu/DrivenCavityMultiGPU/DrivenCavityMultiGPU.cpp
index a353f97a990f9456b1a43b44ffa6d2f72e93821f..f2e75eed402dd387be607cc2203c7d29587b22de 100755
--- a/apps/gpu/DrivenCavityMultiGPU/DrivenCavityMultiGPU.cpp
+++ b/apps/gpu/DrivenCavityMultiGPU/DrivenCavityMultiGPU.cpp
@@ -135,14 +135,14 @@ void runVirtualFluids(const vf::basics::ConfigurationFile &config)
     const real zSplit = 0.0;
 
     if (numberOfProcesses == 2) {
-        gridBuilderFacade->addDomainSplit(zSplit, MultipleGridBuilderFacade::CoordDirection::z);
+        gridBuilderFacade->addDomainSplit(zSplit, Axis::z);
     } else if (numberOfProcesses == 4) {
-        gridBuilderFacade->addDomainSplit(xSplit, MultipleGridBuilderFacade::CoordDirection::x);
-        gridBuilderFacade->addDomainSplit(zSplit, MultipleGridBuilderFacade::CoordDirection::z);
+        gridBuilderFacade->addDomainSplit(xSplit, Axis::x);
+        gridBuilderFacade->addDomainSplit(zSplit, Axis::z);
     } else if (numberOfProcesses == 8) {
-        gridBuilderFacade->addDomainSplit(xSplit, MultipleGridBuilderFacade::CoordDirection::x);
-        gridBuilderFacade->addDomainSplit(ySplit, MultipleGridBuilderFacade::CoordDirection::y);
-        gridBuilderFacade->addDomainSplit(zSplit, MultipleGridBuilderFacade::CoordDirection::z);
+        gridBuilderFacade->addDomainSplit(xSplit, Axis::x);
+        gridBuilderFacade->addDomainSplit(ySplit, Axis::y);
+        gridBuilderFacade->addDomainSplit(zSplit, Axis::z);
     }
 
     // create grids
diff --git a/apps/gpu/SphereMultiGPU/SphereMultiGPU.cpp b/apps/gpu/SphereMultiGPU/SphereMultiGPU.cpp
index d6b6e0f68ab79ae30e326f88bac687302e78cec7..2f83f616ceac0c76577fafc2917e0c2167a078ea 100755
--- a/apps/gpu/SphereMultiGPU/SphereMultiGPU.cpp
+++ b/apps/gpu/SphereMultiGPU/SphereMultiGPU.cpp
@@ -136,14 +136,14 @@ void runVirtualFluids(const vf::basics::ConfigurationFile& config)
     const real zSplit = 0.0;
 
     if (numberOfProcesses == 2) {
-        gridBuilderFacade->addDomainSplit(zSplit, MultipleGridBuilderFacade::CoordDirection::z);
+        gridBuilderFacade->addDomainSplit(zSplit, Axis::z);
     } else if (numberOfProcesses == 4) {
-        gridBuilderFacade->addDomainSplit(xSplit, MultipleGridBuilderFacade::CoordDirection::y);
-        gridBuilderFacade->addDomainSplit(zSplit, MultipleGridBuilderFacade::CoordDirection::z);
+        gridBuilderFacade->addDomainSplit(xSplit, Axis::y);
+        gridBuilderFacade->addDomainSplit(zSplit, Axis::z);
     } else if (numberOfProcesses == 8) {
-        gridBuilderFacade->addDomainSplit(xSplit, MultipleGridBuilderFacade::CoordDirection::x);
-        gridBuilderFacade->addDomainSplit(ySplit, MultipleGridBuilderFacade::CoordDirection::y);
-        gridBuilderFacade->addDomainSplit(zSplit, MultipleGridBuilderFacade::CoordDirection::z);
+        gridBuilderFacade->addDomainSplit(xSplit, Axis::x);
+        gridBuilderFacade->addDomainSplit(ySplit, Axis::y);
+        gridBuilderFacade->addDomainSplit(zSplit, Axis::z);
     }
 
     // create grids
diff --git a/src/basics/geometry3d/Axis.h b/src/basics/geometry3d/Axis.h
index 4bd924db197f35c73a45420031d490f0be957020..9462b827cd14b99d71265989ef87a86ca45f3b2e 100644
--- a/src/basics/geometry3d/Axis.h
+++ b/src/basics/geometry3d/Axis.h
@@ -11,13 +11,17 @@ enum Axis {
     z = 2,
 };
 
+namespace axis
+{
+
 const std::map<Axis, std::array<double, 3>> unitVectors{ { x, { 1, 0, 0 } },
                                                          { y, { 0, 1, 0 } },
                                                          { z, { 0, 0, 1 } } };
 
-namespace axis
-{
 std::string to_string(Axis axis);
-}
+
+const std::array<Axis, 3> allAxes = { Axis::x, Axis::y, Axis::z };
+
+} // namespace axis
 
 #endif
diff --git a/src/gpu/GridGenerator/geometries/Cylinder/Cylinder.cpp b/src/gpu/GridGenerator/geometries/Cylinder/Cylinder.cpp
index a82e8e01f71f315d9eba3d69f268bfe268ebce43..0d4eed024acf5a095fc1701b51eb7bac64bbf89c 100644
--- a/src/gpu/GridGenerator/geometries/Cylinder/Cylinder.cpp
+++ b/src/gpu/GridGenerator/geometries/Cylinder/Cylinder.cpp
@@ -1,6 +1,8 @@
 #include "Cylinder.h"
 #include <numeric>
 
+using namespace axis;
+
 Cylinder::Cylinder(double centerX, double centerY, double centerZ, double radius, double height, Axis rotationalAxis)
     : center({ centerX, centerY, centerZ }), radius(radius), height(height), rotationalAxis(rotationalAxis)
 {
diff --git a/src/gpu/GridGenerator/grid/MultipleGridBuilderFacade.cpp b/src/gpu/GridGenerator/grid/MultipleGridBuilderFacade.cpp
index 3e18525b6fe7ca8bc4ffaf42f20e48a22a6ea5a3..fdb54761ce3f40f11b3ee9f3d6033aec409c07f5 100644
--- a/src/gpu/GridGenerator/grid/MultipleGridBuilderFacade.cpp
+++ b/src/gpu/GridGenerator/grid/MultipleGridBuilderFacade.cpp
@@ -1,4 +1,5 @@
 #include <algorithm>
+#include <climits>
 #include <iostream>
 #include <iterator>
 #include <stdexcept>
@@ -6,10 +7,12 @@
 #include <utility>
 
 #include "MultipleGridBuilderFacade.h"
-#include "grid/GridBuilder/MultipleGridBuilder.h"
-#include "grid/GridDimensions.h"
 #include "geometries/BoundingBox/BoundingBox.h"
 #include "geometries/Object.h"
+#include "grid/GridBuilder/MultipleGridBuilder.h"
+#include "grid/GridDimensions.h"
+
+using namespace CommunicationDirections;
 
 MultipleGridBuilderFacade::MultipleGridBuilderFacade(SPtr<MultipleGridBuilder> gridBuilder,
                                                      SPtr<GridDimensions> gridDimensions,
@@ -32,12 +35,13 @@ void MultipleGridBuilderFacade::createGrids(uint generatePart)
     createGridsHasBeenCalled = true;
 
     this->calculateNumberOfSubdomains();
-    this->numberOfGridsTotal = this->numberGridsX * this->numberGridsY * this->numberGridsZ;
+    this->numberOfSubdomainsTotal =
+        this->numberOfSubdomains[Axis::x] * this->numberOfSubdomains[Axis::y] * this->numberOfSubdomains[Axis::z];
 
-    if (numberOfGridsTotal > 1 && !this->overlapOfSubdomains)
+    if (numberOfSubdomainsTotal > 1 && !this->overlapOfSubdomains)
         throw std::runtime_error("OverlapOfSubdomains in MultipleGridBuilderFacade is NaN.");
 
-    if (generatePart >= numberOfGridsTotal)
+    if (generatePart >= numberOfSubdomainsTotal)
         throw std::runtime_error("Invalid id for subdomain: It is greater or equal to numberOfSubdomains");
 
     this->sortSplitLocations();
@@ -61,9 +65,9 @@ void MultipleGridBuilderFacade::createGrids(uint generatePart)
 
 void MultipleGridBuilderFacade::calculateNumberOfSubdomains()
 {
-    this->numberGridsX = (uint)(this->xSplits.size() + 1);
-    this->numberGridsY = (uint)(this->ySplits.size() + 1);
-    this->numberGridsZ = (uint)(this->zSplits.size() + 1);
+    this->numberOfSubdomains[Axis::x] = (uint)(this->xSplits.size() + 1);
+    this->numberOfSubdomains[Axis::y] = (uint)(this->ySplits.size() + 1);
+    this->numberOfSubdomains[Axis::z] = (uint)(this->zSplits.size() + 1);
 }
 
 void MultipleGridBuilderFacade::sortSplitLocations()
@@ -76,12 +80,12 @@ void MultipleGridBuilderFacade::sortSplitLocations()
 
 void MultipleGridBuilderFacade::calculatedIndexOfPart(uint generatePart)
 {
-    this->xIndex = this->getX3D(generatePart);
-    this->yIndex = this->getY3D(generatePart);
-    this->zIndex = this->getZ3D(generatePart);
+    this->index.at(Axis::x) = this->getX3D(generatePart);
+    this->index.at(Axis::y) = this->getY3D(generatePart);
+    this->index.at(Axis::z) = this->getZ3D(generatePart);
 }
 
-void MultipleGridBuilderFacade::checkSplitLocations(const std::vector<real> &splits, real lowerBound, real upperBound) const
+void MultipleGridBuilderFacade::checkSplitLocations(const std::vector<real>& splits, real lowerBound, real upperBound) const
 {
     if (splits.empty()) return;
 
@@ -124,12 +128,12 @@ void MultipleGridBuilderFacade::configureSubDomainGrids()
     std::vector<real> zValues = { this->gridDimensionsDomain->minZ, this->gridDimensionsDomain->maxZ };
     zValues.insert(std::prev(zValues.end()), this->zSplits.begin(), this->zSplits.end());
 
-    real xMinCoarseGrid = xValues[xIndex];
-    real yMinCoarseGrid = yValues[yIndex];
-    real zMinCoarseGrid = zValues[zIndex];
-    real xMaxCoarseGrid = xValues[xIndex + 1];
-    real yMaxCoarseGrid = yValues[yIndex + 1];
-    real zMaxCoarseGrid = zValues[zIndex + 1];
+    real xMinCoarseGrid = xValues[index.at(Axis::x)];
+    real yMinCoarseGrid = yValues[index.at(Axis::y)];
+    real zMinCoarseGrid = zValues[index.at(Axis::z)];
+    real xMaxCoarseGrid = xValues[index.at(Axis::x) + 1];
+    real yMaxCoarseGrid = yValues[index.at(Axis::y) + 1];
+    real zMaxCoarseGrid = zValues[index.at(Axis::z) + 1];
 
     // add overlap
     xMinCoarseGrid -= (hasNeighbors[CommunicationDirections::MX]) ? overlapOfSubdomains.value() : 0;
@@ -148,17 +152,17 @@ void MultipleGridBuilderFacade::configureSubDomainGrids()
 
     // set subdomain boxes
     // subdomain boxes are only needed on multiple gpus
-    if ((numberGridsX * numberGridsY * numberGridsZ) > 1) {
-        gridBuilder->setSubDomainBox(std::make_shared<BoundingBox>(xValues[xIndex], xValues[xIndex + 1],
-                                                                   yValues[yIndex], yValues[yIndex + 1],
-                                                                   zValues[zIndex], zValues[zIndex + 1]));
+    if ((numberOfSubdomains[Axis::x] * numberOfSubdomains[Axis::y] * numberOfSubdomains[Axis::z]) > 1) {
+        gridBuilder->setSubDomainBox(std::make_shared<BoundingBox>(xValues[index.at(Axis::x)], xValues[index.at(Axis::x) + 1],
+                                                                   yValues[index.at(Axis::y)], yValues[index.at(Axis::y) + 1],
+                                                                   zValues[index.at(Axis::z)], zValues[index.at(Axis::z) + 1]));
     }
 }
 
 void MultipleGridBuilderFacade::setUpCommunicationNeighbors()
 {
     // Communication is only needed on multiple gpus
-    if (numberOfGridsTotal == 1) return;
+    if (numberOfSubdomainsTotal == 1) return;
 
     if (hasNeighbors.empty())
         throw std::runtime_error("checkForNeighbors() has to be called befor calling setUpCommunicationNeighbors()");
@@ -169,22 +173,28 @@ void MultipleGridBuilderFacade::setUpCommunicationNeighbors()
 
             switch (direction) {
                 case CommunicationDirections::MX:
-                    gridBuilder->setCommunicationProcess(direction, getIndex1D(xIndex - 1, yIndex, zIndex));
+                    gridBuilder->setCommunicationProcess(
+                        direction, getIndex1D(index.at(Axis::x) - 1, index.at(Axis::y), index.at(Axis::z)));
                     break;
                 case CommunicationDirections::MY:
-                    gridBuilder->setCommunicationProcess(direction, getIndex1D(xIndex, yIndex - 1, zIndex));
+                    gridBuilder->setCommunicationProcess(
+                        direction, getIndex1D(index.at(Axis::x), index.at(Axis::y) - 1, index.at(Axis::z)));
                     break;
                 case CommunicationDirections::MZ:
-                    gridBuilder->setCommunicationProcess(direction, getIndex1D(xIndex, yIndex, zIndex - 1));
+                    gridBuilder->setCommunicationProcess(
+                        direction, getIndex1D(index.at(Axis::x), index.at(Axis::y), index.at(Axis::z) - 1));
                     break;
                 case CommunicationDirections::PX:
-                    gridBuilder->setCommunicationProcess(direction, getIndex1D(xIndex + 1, yIndex, zIndex));
+                    gridBuilder->setCommunicationProcess(
+                        direction, getIndex1D(index.at(Axis::x) + 1, index.at(Axis::y), index.at(Axis::z)));
                     break;
                 case CommunicationDirections::PY:
-                    gridBuilder->setCommunicationProcess(direction, getIndex1D(xIndex, yIndex + 1, zIndex));
+                    gridBuilder->setCommunicationProcess(
+                        direction, getIndex1D(index.at(Axis::x), index.at(Axis::y) + 1, index.at(Axis::z)));
                     break;
                 case CommunicationDirections::PZ:
-                    gridBuilder->setCommunicationProcess(direction, getIndex1D(xIndex, yIndex, zIndex + 1));
+                    gridBuilder->setCommunicationProcess(
+                        direction, getIndex1D(index.at(Axis::x), index.at(Axis::y), index.at(Axis::z) + 1));
                     break;
             }
         }
@@ -193,24 +203,24 @@ void MultipleGridBuilderFacade::setUpCommunicationNeighbors()
 
 void MultipleGridBuilderFacade::checkForNeighbors()
 {
-    hasNeighbors[CommunicationDirections::MX] = (xIndex > 0);
-    hasNeighbors[CommunicationDirections::MY] = (yIndex > 0);
-    hasNeighbors[CommunicationDirections::MZ] = (zIndex > 0);
-    hasNeighbors[CommunicationDirections::PX] = (xIndex < numberGridsX - 1);
-    hasNeighbors[CommunicationDirections::PY] = (yIndex < numberGridsY - 1);
-    hasNeighbors[CommunicationDirections::PZ] = (zIndex < numberGridsZ - 1);
+    hasNeighbors[CommunicationDirections::MX] = (index.at(Axis::x) > 0);
+    hasNeighbors[CommunicationDirections::MY] = (index.at(Axis::y) > 0);
+    hasNeighbors[CommunicationDirections::MZ] = (index.at(Axis::z) > 0);
+    hasNeighbors[CommunicationDirections::PX] = (index.at(Axis::x) < numberOfSubdomains[Axis::x] - 1);
+    hasNeighbors[CommunicationDirections::PY] = (index.at(Axis::y) < numberOfSubdomains[Axis::y] - 1);
+    hasNeighbors[CommunicationDirections::PZ] = (index.at(Axis::z) < numberOfSubdomains[Axis::z] - 1);
 }
 
 void MultipleGridBuilderFacade::addFineGridsToGridBuilder()
 {
-    for (auto const &grid : fineGrids) {
+    for (auto const& grid : fineGrids) {
         gridBuilder->addGrid(grid.first, grid.second);
     }
 }
 
 void MultipleGridBuilderFacade::addGeometriesToGridBuilder()
 {
-    for (auto const &geometry : geometries) {
+    for (auto const& geometry : geometries) {
         gridBuilder->addGeometry(geometry);
     }
 }
@@ -220,7 +230,7 @@ void MultipleGridBuilderFacade::setOverlapOfSubdomains(real overlap)
     this->overlapOfSubdomains = overlap;
 }
 
-void MultipleGridBuilderFacade::addDomainSplit(real coordinate, MultipleGridBuilderFacade::CoordDirection direction)
+void MultipleGridBuilderFacade::addDomainSplit(real coordinate, Axis direction)
 {
     if (this->createGridsHasBeenCalled)
         throw std::runtime_error("MultipleGridBuilderFacade::addSplit() should be called before createGrids().");
@@ -262,17 +272,17 @@ void MultipleGridBuilderFacade::setNumberOfLayersForRefinement(uint numberOfLaye
 
 uint MultipleGridBuilderFacade::getX3D(uint index1D) const
 {
-    return index1D % numberGridsX;
+    return index1D % numberOfSubdomains[Axis::x];
 }
 
 uint MultipleGridBuilderFacade::getY3D(uint index1D) const
 {
-    return (index1D / numberGridsX) % numberGridsY;
+    return (index1D / numberOfSubdomains[Axis::x]) % numberOfSubdomains[Axis::y];
 }
 
 uint MultipleGridBuilderFacade::getZ3D(uint index1D) const
 {
-    return index1D / (numberGridsX * numberGridsY);
+    return index1D / (numberOfSubdomains[Axis::x] * numberOfSubdomains[Axis::y]);
 }
 
 std::array<uint, 3> MultipleGridBuilderFacade::convertToIndices3D(uint index1D) const
@@ -285,16 +295,21 @@ std::array<uint, 3> MultipleGridBuilderFacade::convertToIndices3D(uint index1D)
 
 uint MultipleGridBuilderFacade::getIndex1D(uint xIndex, uint yIndex, uint zIndex) const
 {
-    return xIndex + yIndex * numberGridsX + zIndex * numberGridsX * numberGridsY;
+    return xIndex + yIndex * numberOfSubdomains[Axis::x] + zIndex * numberOfSubdomains[Axis::x] * numberOfSubdomains[Axis::y];
+}
+
+uint MultipleGridBuilderFacade::getIndex1D(const std::array<uint, 3>& index3D) const
+{
+    return getIndex1D(index3D[Axis::x], index3D[Axis::y], index3D[Axis::z]);
 }
 
-void MultipleGridBuilderFacade::setSlipBoundaryCondition(SideType sideType, real normalX, real normalY, real normalZ)
+void MultipleGridBuilderFacade::setSlipBoundaryCondition(SideType sideType, real normalX, real normalY, real normalZ) const
 {
     setBoundaryCondition(sideType, [&]() { gridBuilder->setSlipBoundaryCondition(sideType, normalX, normalY, normalZ); });
 }
 
 void MultipleGridBuilderFacade::setStressBoundaryCondition(SideType sideType, real normalX, real normalY, real normalZ,
-                                                           uint samplingOffset, real z0, real dx)
+                                                           uint samplingOffset, real z0, real dx) const
 {
     setBoundaryCondition(sideType, [&]() {
         gridBuilder->setStressBoundaryCondition(sideType, normalX, normalY, normalZ, samplingOffset, z0, dx);
@@ -303,7 +318,8 @@ void MultipleGridBuilderFacade::setStressBoundaryCondition(SideType sideType, re
 
 void MultipleGridBuilderFacade::setPrecursorBoundaryCondition(SideType sideType, SPtr<FileCollection> fileCollection,
                                                               int timeStepsBetweenReads, real velocityX, real velocityY,
-                                                              real velocityZ, std::vector<uint> fileLevelToGridLevelMap)
+                                                              real velocityZ,
+                                                              std::vector<uint> fileLevelToGridLevelMap) const
 {
     setBoundaryCondition(sideType, [&]() {
         gridBuilder->setPrecursorBoundaryCondition(sideType, fileCollection, timeStepsBetweenReads, velocityX, velocityY,
@@ -311,76 +327,83 @@ void MultipleGridBuilderFacade::setPrecursorBoundaryCondition(SideType sideType,
     });
 }
 
-void MultipleGridBuilderFacade::setVelocityBoundaryCondition(SideType sideType, real vx, real vy, real vz)
+void MultipleGridBuilderFacade::setVelocityBoundaryCondition(SideType sideType, real vx, real vy, real vz) const
 {
     setBoundaryCondition(sideType, [&]() { gridBuilder->setVelocityBoundaryCondition(sideType, vx, vy, vz); });
 }
 
-void MultipleGridBuilderFacade::setPressureBoundaryCondition(SideType sideType, real rho)
+void MultipleGridBuilderFacade::setPressureBoundaryCondition(SideType sideType, real rho) const
 {
     setBoundaryCondition(sideType, [&]() { gridBuilder->setPressureBoundaryCondition(sideType, rho); });
 }
 
-void MultipleGridBuilderFacade::setNoSlipBoundaryCondition(SideType sideType)
+void MultipleGridBuilderFacade::setNoSlipBoundaryCondition(SideType sideType) const
 {
     setBoundaryCondition(sideType, [&]() { gridBuilder->setNoSlipBoundaryCondition(sideType); });
 }
 
-void MultipleGridBuilderFacade::setPeriodicBoundaryCondition(bool periodic_X, bool periodic_Y, bool periodic_Z)
+bool MultipleGridBuilderFacade::isFinalSubdomainInDirection(CommunicationDirection direction) const
 {
-    bool localPeriodicityX = false;
-    bool localPeriodicityY = false;
-    bool localPeriodicityZ = false;
+    return !hasNeighbors.at(direction);
+}
 
-    if (periodic_X) {
-        if (numberGridsX == 1) {
-            localPeriodicityX = true;
-        }
-        if (numberGridsX > 1 && !hasNeighbors[CommunicationDirections::MX]) {
-            // set last grid in x-direction as communication neighbor
-            gridBuilder->findCommunicationIndices(CommunicationDirections::MX);
-            gridBuilder->setCommunicationProcess(CommunicationDirections::MX, getIndex1D(numberGridsX - 1, yIndex, zIndex));
-        } else if (numberGridsX > 1 && !hasNeighbors[CommunicationDirections::PX]) {
-            // set first grid in x-direction as communication neighbor
-            gridBuilder->findCommunicationIndices(CommunicationDirections::PX);
-            gridBuilder->setCommunicationProcess(CommunicationDirections::PX, getIndex1D(0, yIndex, zIndex));
-        }
+uint MultipleGridBuilderFacade::getIndexOfFinalSubdomainInDirection(CommunicationDirection direction) const
+{
+    std::array<uint, 3> resultIndex3D = index;
+    const Axis axis = communicationDirectionToAxes.at(direction);
+
+    if (isNegative(direction)) {
+        resultIndex3D[axis] = 0; // first subdomain index in direction
+        return getIndex1D(resultIndex3D);
+    }
+    if (isPositive(direction)) {
+        resultIndex3D[axis] = numberOfSubdomains[axis] - 1; // last subdomain index in direction
+        return getIndex1D(resultIndex3D);
     }
+    return UINT_MAX;
+}
 
-    if (periodic_Y) {
-        if (numberGridsY == 1) {
-            localPeriodicityY = true;
-        }
-        if (numberGridsY > 1 && !hasNeighbors[CommunicationDirections::MY]) {
-            // set last grid in x-direction as communication neighbor
-            gridBuilder->findCommunicationIndices(CommunicationDirections::MY);
-            gridBuilder->setCommunicationProcess(CommunicationDirections::MY, getIndex1D(xIndex, numberGridsY - 1, zIndex));
-        } else if (numberGridsY > 1 && !hasNeighbors[CommunicationDirections::PY]) {
-            // set first grid in x-direction as communication neighbor
-            gridBuilder->findCommunicationIndices(CommunicationDirections::PY);
-            gridBuilder->setCommunicationProcess(CommunicationDirections::PY, getIndex1D(xIndex, 0, zIndex));
+void MultipleGridBuilderFacade::setPeriodicBoundaryCondition(bool periodic_X, bool periodic_Y, bool periodic_Z) const
+{
+    setPeriodicBoundaryCondition({ periodic_X, periodic_Y, periodic_Z });
+}
+
+void MultipleGridBuilderFacade::setPeriodicBoundaryCondition(const std::array<bool, 3>& periodicity) const
+{
+
+    std::array<bool, 3> localPeriodicity = { false, false, false };
+
+    for (const auto coordDirection : axis::allAxes) {
+        if (!periodicity[coordDirection]) {
+            continue;
         }
-    }
 
-    if (periodic_Z) {
-        if (numberGridsZ == 1) {
-            localPeriodicityZ = true;
+        // only one grid in direction --> set local periodicity
+        if (numberOfSubdomains[coordDirection] == 1) {
+            localPeriodicity[coordDirection] = true;
+            continue;
         }
-        if (numberGridsZ > 1 && !hasNeighbors[CommunicationDirections::MZ]) {
-            // set last grid in x-direction as communication neighbor
-            gridBuilder->findCommunicationIndices(CommunicationDirections::MZ);
-            gridBuilder->setCommunicationProcess(CommunicationDirections::MZ, getIndex1D(xIndex, yIndex, numberGridsZ - 1));
-        } else if (numberGridsZ > 1 && !hasNeighbors[CommunicationDirections::PZ]) {
-            // set first grid in x-direction as communication neighbor
-            gridBuilder->findCommunicationIndices(CommunicationDirections::PZ);
-            gridBuilder->setCommunicationProcess(CommunicationDirections::PZ, getIndex1D(xIndex, yIndex, 0));
+
+        // non-local periodicity --> set communication neighbors
+        const CommunicationDirection negativeDirection = getNegativeDirectionAlongAxis(coordDirection);
+        const CommunicationDirection positiveDirection = getPositiveDirectionAlongAxis(coordDirection);
+
+        if (isFinalSubdomainInDirection(negativeDirection)) {
+            // set final grid in positive direction as communication neighbor
+            gridBuilder->findCommunicationIndices(negativeDirection);
+            gridBuilder->setCommunicationProcess(negativeDirection, getIndexOfFinalSubdomainInDirection(positiveDirection));
+        } else if (isFinalSubdomainInDirection(positiveDirection)) {
+            // set final grid in negative direction as communication neighbor
+            gridBuilder->findCommunicationIndices(positiveDirection);
+            gridBuilder->setCommunicationProcess(positiveDirection, getIndexOfFinalSubdomainInDirection(negativeDirection));
         }
     }
 
-    gridBuilder->setPeriodicBoundaryCondition(localPeriodicityX, localPeriodicityY, localPeriodicityZ);
+    gridBuilder->setPeriodicBoundaryCondition(localPeriodicity[Axis::x], localPeriodicity[Axis::y],
+                                              localPeriodicity[Axis::z]);
 }
 
-SPtr<MultipleGridBuilder> MultipleGridBuilderFacade::getGridBuilder()
+SPtr<MultipleGridBuilder> MultipleGridBuilderFacade::getGridBuilder() const
 {
     return gridBuilder;
 }
\ No newline at end of file
diff --git a/src/gpu/GridGenerator/grid/MultipleGridBuilderFacade.h b/src/gpu/GridGenerator/grid/MultipleGridBuilderFacade.h
index ae139d32c79143cf2a8b5e3970729853e6db3e0c..a36d8e10a0e15b7126b3ed3e87f831cdb0fbafb9 100644
--- a/src/gpu/GridGenerator/grid/MultipleGridBuilderFacade.h
+++ b/src/gpu/GridGenerator/grid/MultipleGridBuilderFacade.h
@@ -45,6 +45,7 @@
 
 #include <basics/DataTypes.h>
 #include <basics/constants/NumericConstants.h>
+#include <basics/geometry3d/Axis.h>
 
 #include "grid/BoundaryConditions/Side.h"
 #include "utilities/communication.h"
@@ -85,8 +86,6 @@ using namespace vf::basics::constant;
 class MultipleGridBuilderFacade
 {
 public:
-    enum CoordDirection { x, y, z };
-
     MultipleGridBuilderFacade(SPtr<MultipleGridBuilder> gridBuilder, SPtr<GridDimensions> gridDimensions,
                               std::optional<real> overlapOfSubdomains = std::nullopt);
 
@@ -94,7 +93,7 @@ public:
 
     //! \brief split the domain in the passed direction at the passed coordinate
     //! \details multiple splits can be added to a domain
-    void addDomainSplit(real coordinate, MultipleGridBuilderFacade::CoordDirection direction);
+    void addDomainSplit(real coordinate, Axis direction);
 
     //! \brief set the overlap of the subdomains
     void setOverlapOfSubdomains(real overlap);
@@ -116,18 +115,19 @@ public:
     void createGrids(uint generatePart);
 
     // Boundary conditions, call after createGrids()
-    void setSlipBoundaryCondition(SideType sideType, real normalX, real normalY, real normalZ);
+    void setSlipBoundaryCondition(SideType sideType, real normalX, real normalY, real normalZ) const;
     void setStressBoundaryCondition(SideType sideType, real normalX, real normalY, real normalZ, uint samplingOffset,
-                                    real z0, real dx);
-    void setVelocityBoundaryCondition(SideType sideType, real vx, real vy, real vz);
-    void setPressureBoundaryCondition(SideType sideType, real rho);
-    void setNoSlipBoundaryCondition(SideType sideType);
-    void setPeriodicBoundaryCondition(bool periodic_X, bool periodic_Y, bool periodic_Z);
+                                    real z0, real dx) const;
+    void setVelocityBoundaryCondition(SideType sideType, real vx, real vy, real vz) const;
+    void setPressureBoundaryCondition(SideType sideType, real rho) const;
+    void setNoSlipBoundaryCondition(SideType sideType) const;
+    void setPeriodicBoundaryCondition(bool periodic_X, bool periodic_Y, bool periodic_Z) const;
+    void setPeriodicBoundaryCondition(const std::array<bool, 3>& periodicity) const;
     void setPrecursorBoundaryCondition(SideType sideType, SPtr<FileCollection> fileCollection, int timeStepsBetweenReads,
                                        real velocityX = c0o1, real velocityY = c0o1, real velocityZ = c0o1,
-                                       std::vector<uint> fileLevelToGridLevelMap = {});
+                                       std::vector<uint> fileLevelToGridLevelMap = {}) const;
 
-    SPtr<MultipleGridBuilder> getGridBuilder();
+    SPtr<MultipleGridBuilder> getGridBuilder() const;
 
 protected:
     // index calculations
@@ -136,6 +136,7 @@ protected:
     uint getY3D(uint index1D) const;
     uint getZ3D(uint index1D) const;
     uint getIndex1D(uint xIndex, uint yIndex, uint zIndex) const;
+    uint getIndex1D(const std::array<uint, 3>& index3D) const;
 
 private:
     //! \brief calculate the number of subdomains in all coordinate directions
@@ -150,24 +151,30 @@ private:
     //! \brief for each direction, calculate if the current subdomain has a neighbor in this direction
     void checkForNeighbors();
 
-    //! \brief set up coarse grids and subdomain boxes for all grids
+    //! \brief set up coarse grids and subdomain boxes for all subdomains
     void configureSubDomainGrids();
 
     //! \brief set up the communication to neighboring subdomains
     void setUpCommunicationNeighbors();
 
     //! \brief check if all locations for domain splits are inside the grid bounds and there are no duplicates.
-    void checkSplitLocations(const std::vector<real> &splits, real lowerBound, real upperBound) const;
+    void checkSplitLocations(const std::vector<real>& splits, real lowerBound, real upperBound) const;
 
-    //! \brief add fine grids to the gridBuilder
+    //! \brief add fine grids to the grid builder
     void addFineGridsToGridBuilder();
 
     //! \brief add geometries to the grid builder
     void addGeometriesToGridBuilder();
 
+    //! \brief check whether a subdomain is the last one in a direction
+    bool isFinalSubdomainInDirection(CommunicationDirections::CommunicationDirection direction) const;
+
+    //! \brief get 1D index of the final subdomain in a direction, in the other directions it has the same position as the current subdomain
+    uint getIndexOfFinalSubdomainInDirection(CommunicationDirections::CommunicationDirection direction) const;
+
     //! \brief call the grid builder's setter for a boundary condition
     template <typename function>
-    void setBoundaryCondition(SideType sideType, function boundaryConditionFunction)
+    void setBoundaryCondition(SideType sideType, function boundaryConditionFunction) const
     {
         if (!createGridsHasBeenCalled) {
             throw std::runtime_error(
@@ -175,16 +182,14 @@ private:
         }
 
         if (sideType == SideType::GEOMETRY ||
-            !hasNeighbors[static_cast<CommunicationDirections::CommunicationDirection>(static_cast<int>(sideType))]) {
+            !hasNeighbors.at(static_cast<CommunicationDirections::CommunicationDirection>(static_cast<int>(sideType)))) {
             boundaryConditionFunction();
         }
     }
 
 protected:
-    //! \brief total number of subdomains (computed)
-    uint numberGridsX;
-    uint numberGridsY;
-    uint numberGridsZ;
+    //! \brief total number of subdomains in each coordinate direction (computed)
+    std::array<uint, 3> numberOfSubdomains;
 
 private:
     const SPtr<MultipleGridBuilder> gridBuilder;
@@ -192,7 +197,7 @@ private:
     // basic grid dimension (set in constructor)
     const SPtr<GridDimensions> gridDimensionsDomain;
 
-    uint numberOfGridsTotal;
+    uint numberOfSubdomainsTotal;
 
     //! \brief coordinates, signifying where the domain is split into subdomains (must be set in a setter)
     std::vector<real> xSplits;
@@ -203,9 +208,7 @@ private:
     std::optional<real> overlapOfSubdomains = std::nullopt;
 
     //! \brief index of the current subdomains in relation to all subdomains (computed)
-    uint xIndex;
-    uint yIndex;
-    uint zIndex;
+    std::array<uint, 3> index;
 
     //! \brief hasNeighbors, indicates if the current subdomains has a neighbor in a specific direction (computed)
     //! \details use the enum CommunciationDirection to access the data
diff --git a/src/gpu/GridGenerator/grid/MultipleGridBuilderFacadeTest.cpp b/src/gpu/GridGenerator/grid/MultipleGridBuilderFacadeTest.cpp
index a6160eebd335d83943108875fc31991f9a69e1f4..20f058f97389761ac41b70e8884a834db657129f 100644
--- a/src/gpu/GridGenerator/grid/MultipleGridBuilderFacadeTest.cpp
+++ b/src/gpu/GridGenerator/grid/MultipleGridBuilderFacadeTest.cpp
@@ -26,27 +26,29 @@ class MultipleGridBuilderFacadeTest : public MultipleGridBuilderFacade
 
     void setNumberOfGrids(uint x, uint y, uint z)
     {
-        this->numberGridsX = x;
-        this->numberGridsY = y;
-        this->numberGridsZ = z;
+        numberOfSubdomains = { x, y, z };
     };
 
-    uint getX3D(const uint index1D)
+    uint getX3D(uint index1D)
     {
         return MultipleGridBuilderFacade::getX3D(index1D);
     }
-    uint getY3D(const uint index1D)
+    uint getY3D(uint index1D)
     {
         return MultipleGridBuilderFacade::getY3D(index1D);
     }
-    uint getZ3D(const uint index1D)
+    uint getZ3D(uint index1D)
     {
         return MultipleGridBuilderFacade::getZ3D(index1D);
     }
-    uint getIndex1D(const uint xIndex, const uint yIndex, const uint zIndex)
+    uint getIndex1D(uint xIndex, uint yIndex, uint zIndex)
     {
         return MultipleGridBuilderFacade::getIndex1D(xIndex, yIndex, zIndex);
     }
+    uint getIndex1D(const std::array<uint, 3>& index3D)
+    {
+        return MultipleGridBuilderFacade::getIndex1D(index3D);
+    }
 };
 
 class MockMultipleGridBuilder : public MultipleGridBuilder
@@ -138,16 +140,20 @@ TEST(MultipleGridBuilderFacadeTest, transformComponentsTo1DCoordinate)
     sut.setNumberOfGrids(2, 3, 4);
 
     EXPECT_THAT(sut.getIndex1D(0, 0, 0), testing::Eq(0));
+    EXPECT_THAT(sut.getIndex1D({0, 0, 0}), testing::Eq(0));
     EXPECT_THAT(sut.getIndex1D(1, 2, 1), testing::Eq(11));
+    EXPECT_THAT(sut.getIndex1D({1, 2, 1}), testing::Eq(11));
     EXPECT_THAT(sut.getIndex1D(0, 0, 2), testing::Eq(12));
+    EXPECT_THAT(sut.getIndex1D({0, 0, 2}), testing::Eq(12));
     EXPECT_THAT(sut.getIndex1D(1, 2, 3), testing::Eq(23));
+    EXPECT_THAT(sut.getIndex1D({1, 2, 3}), testing::Eq(23));
 }
 
 TEST(MultipleGridBuilderFacadeTest, noOverlapOnMultiGpu_Throws)
 {
     auto gridBuilder = std::make_shared<MultipleGridBuilder>();
     auto sut = MultipleGridBuilderFacadeTest(gridBuilder, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0);
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
+    sut.addDomainSplit(1.0, Axis::x);
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
 }
 
@@ -206,12 +212,12 @@ protected:
     void createNewSut()
     {
         sut = std::make_unique<MultipleGridBuilderFacadeTest>(mockGridBuilder, std::make_unique<GridDimensions>(0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1), 0.1);
-        sut->addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
-        sut->addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::y);
-        sut->addDomainSplit(2.0, MultipleGridBuilderFacade::CoordDirection::y);
-        sut->addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::z);
-        sut->addDomainSplit(2.0, MultipleGridBuilderFacade::CoordDirection::z);
-        sut->addDomainSplit(3.0, MultipleGridBuilderFacade::CoordDirection::z);
+        sut->addDomainSplit(1.0, Axis::x);
+        sut->addDomainSplit(1.0, Axis::y);
+        sut->addDomainSplit(2.0, Axis::y);
+        sut->addDomainSplit(1.0, Axis::z);
+        sut->addDomainSplit(2.0, Axis::z);
+        sut->addDomainSplit(3.0, Axis::z);
     }
 };
 
@@ -241,7 +247,7 @@ TEST_F(MultipleGridBuilderFacadeTest_24subdomains, createGridsMultiGpu)
 TEST(MultipleGridBuilderFacadeTest, xSplitToLarge)
 {
     auto sut = MultipleGridBuilderFacadeTest(nullptr, 0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1, 0.1);
-    sut.addDomainSplit(10.0, MultipleGridBuilderFacade::CoordDirection::x); // xSplit > maxX
+    sut.addDomainSplit(10.0, Axis::x); // xSplit > maxX
 
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
     EXPECT_THROW(sut.createGrids(1), std::runtime_error);
@@ -250,7 +256,7 @@ TEST(MultipleGridBuilderFacadeTest, xSplitToLarge)
 TEST(MultipleGridBuilderFacadeTest, xSplitToSmall)
 {
     auto sut = MultipleGridBuilderFacadeTest(nullptr, 0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1, 0.1);
-    sut.addDomainSplit(-1.0, MultipleGridBuilderFacade::CoordDirection::x); // xSplit < minX
+    sut.addDomainSplit(-1.0, Axis::x); // xSplit < minX
 
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
     EXPECT_THROW(sut.createGrids(1), std::runtime_error);
@@ -259,8 +265,8 @@ TEST(MultipleGridBuilderFacadeTest, xSplitToSmall)
 TEST(MultipleGridBuilderFacadeTest, ySplitToLarge)
 {
     auto sut = MultipleGridBuilderFacadeTest(nullptr, 0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1, 0.1);
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::y);  // valid ySplit
-    sut.addDomainSplit(10.0, MultipleGridBuilderFacade::CoordDirection::y); // ySplit > maxY
+    sut.addDomainSplit(1.0, Axis::y);  // valid ySplit
+    sut.addDomainSplit(10.0, Axis::y); // ySplit > maxY
 
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
     EXPECT_THROW(sut.createGrids(1), std::runtime_error);
@@ -269,8 +275,8 @@ TEST(MultipleGridBuilderFacadeTest, ySplitToLarge)
 TEST(MultipleGridBuilderFacadeTest, ySplitToSmall)
 {
     auto sut = MultipleGridBuilderFacadeTest(nullptr, 0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1, 0.1);
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::y);
-    sut.addDomainSplit(-1.0, MultipleGridBuilderFacade::CoordDirection::y); // ySplit < minY
+    sut.addDomainSplit(1.0, Axis::y);
+    sut.addDomainSplit(-1.0, Axis::y); // ySplit < minY
 
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
     EXPECT_THROW(sut.createGrids(1), std::runtime_error);
@@ -279,9 +285,9 @@ TEST(MultipleGridBuilderFacadeTest, ySplitToSmall)
 TEST(MultipleGridBuilderFacadeTest, zSplitToLarge)
 {
     auto sut = MultipleGridBuilderFacadeTest(nullptr, 0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1, 0.1);
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::z);
-    sut.addDomainSplit(10.0, MultipleGridBuilderFacade::CoordDirection::z); // zSplit > maxZ
-    sut.addDomainSplit(2.0, MultipleGridBuilderFacade::CoordDirection::z);
+    sut.addDomainSplit(1.0, Axis::z);
+    sut.addDomainSplit(10.0, Axis::z); // zSplit > maxZ
+    sut.addDomainSplit(2.0, Axis::z);
 
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
     EXPECT_THROW(sut.createGrids(1), std::runtime_error);
@@ -290,9 +296,9 @@ TEST(MultipleGridBuilderFacadeTest, zSplitToLarge)
 TEST(MultipleGridBuilderFacadeTest, zSplitToSmall)
 {
     auto sut = MultipleGridBuilderFacadeTest(nullptr, 0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1, 0.1);
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::z);
-    sut.addDomainSplit(-1.0, MultipleGridBuilderFacade::CoordDirection::z); // zSplit < minZ
-    sut.addDomainSplit(2.0, MultipleGridBuilderFacade::CoordDirection::z);
+    sut.addDomainSplit(1.0, Axis::z);
+    sut.addDomainSplit(-1.0, Axis::z); // zSplit < minZ
+    sut.addDomainSplit(2.0, Axis::z);
 
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
     EXPECT_THROW(sut.createGrids(1), std::runtime_error);
@@ -301,9 +307,9 @@ TEST(MultipleGridBuilderFacadeTest, zSplitToSmall)
 TEST(MultipleGridBuilderFacadeTest, sameSplitTwiceY)
 {
     auto sut = MultipleGridBuilderFacadeTest(nullptr, 0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1, 0.1);
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::y);
-    sut.addDomainSplit(2.0, MultipleGridBuilderFacade::CoordDirection::y);
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::y);
+    sut.addDomainSplit(1.0, Axis::y);
+    sut.addDomainSplit(2.0, Axis::y);
+    sut.addDomainSplit(1.0, Axis::y);
 
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
     EXPECT_THROW(sut.createGrids(1), std::runtime_error);
@@ -312,8 +318,8 @@ TEST(MultipleGridBuilderFacadeTest, sameSplitTwiceY)
 TEST(MultipleGridBuilderFacadeTest, sameSplitTwiceZ)
 {
     auto sut = MultipleGridBuilderFacadeTest(nullptr, 0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1);
-    sut.addDomainSplit(0.9, MultipleGridBuilderFacade::CoordDirection::z);
-    sut.addDomainSplit(0.9, MultipleGridBuilderFacade::CoordDirection::z);
+    sut.addDomainSplit(0.9, Axis::z);
+    sut.addDomainSplit(0.9, Axis::z);
 
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
     EXPECT_THROW(sut.createGrids(1), std::runtime_error);
@@ -322,8 +328,8 @@ TEST(MultipleGridBuilderFacadeTest, sameSplitTwiceZ)
 TEST(MultipleGridBuilderFacadeTest, sameSplitTwiceX)
 {
     auto sut = MultipleGridBuilderFacadeTest(nullptr, 0.0, 2.0, 0.0, 3.0, 0.0, 4.0, 0.1, 0.1);
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
+    sut.addDomainSplit(1.0, Axis::x);
+    sut.addDomainSplit(1.0, Axis::x);
 
     EXPECT_THROW(sut.createGrids(0), std::runtime_error);
     EXPECT_THROW(sut.createGrids(1), std::runtime_error);
@@ -628,6 +634,7 @@ TEST_F(MultipleGridBuilderFacadeTest_24subdomains, periodicAllDirectionsMultiGPU
     bool periodic_Y = true;
     bool periodic_Z = true;
 
+    // no local periodicity, periodicity is realized by setting up inter-gpu communication
     EXPECT_CALL(*mockGridBuilder, setPeriodicBoundaryCondition(periodic_X, periodic_Y, periodic_Z)).Times(0);
     EXPECT_CALL(*mockGridBuilder, setPeriodicBoundaryCondition(false, false, false)).Times(24);
     for (int i = 0; i < 24; i++) {
@@ -744,7 +751,7 @@ TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, periodicAllDirectionsSing
 
 TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, periodicXY2GPUs)
 {
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
+    sut.addDomainSplit(1.0, Axis::x);
     sut.setOverlapOfSubdomains(0.1);
 
     EXPECT_CALL(*mockGridBuilder, findCommunicationIndices(CommunicationDirections::PX, testing::_));
@@ -758,7 +765,7 @@ TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, periodicXY2GPUs)
 
 TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, periodicYZ2GPUs)
 {
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
+    sut.addDomainSplit(1.0, Axis::x);
     sut.setOverlapOfSubdomains(0.1);
 
     EXPECT_CALL(*mockGridBuilder, findCommunicationIndices(CommunicationDirections::PX, testing::_));
@@ -809,7 +816,7 @@ TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, addFineGrid_createGridCal
 
 TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, addFineGrid_createGridCallsFineGridsInProcess0)
 {
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
+    sut.addDomainSplit(1.0, Axis::x);
     sut.setOverlapOfSubdomains(0.1);
     std::shared_ptr<Object> fineGrid = std::make_shared<Sphere>(0.0, 0.0, 0.0, 10.0);
     sut.addFineGrid(fineGrid, 1);
@@ -821,7 +828,7 @@ TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, addFineGrid_createGridCal
 
 TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, addFineGrid_createGridCallsFineGridsInProcess1)
 {
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
+    sut.addDomainSplit(1.0, Axis::x);
     sut.setOverlapOfSubdomains(0.1);
     std::shared_ptr<Object> fineGrid = std::make_shared<Sphere>(0.0, 0.0, 0.0, 10.0);
     sut.addFineGrid(fineGrid, 1);
@@ -849,7 +856,7 @@ TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, noFineGrid_createGrid_doe
 
 TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, addGeometry_createGridCallsAddGeometryFunctionOfGridBuilderProcess0)
 {
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
+    sut.addDomainSplit(1.0, Axis::x);
     sut.setOverlapOfSubdomains(0.1);
     std::shared_ptr<Object> geometry = std::make_shared<Sphere>(0.0, 0.0, 0.0, 10.0);
     sut.addGeometry(geometry);
@@ -860,7 +867,7 @@ TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, addGeometry_createGridCal
 
 TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, addGeometry_createGridCallsAddGeometryFunctionOfGridBuilderProcess1)
 {
-    sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x);
+    sut.addDomainSplit(1.0, Axis::x);
     sut.setOverlapOfSubdomains(0.1);
     std::shared_ptr<Object> geometry = std::make_shared<Sphere>(0.0, 0.0, 0.0, 10.0);
     sut.addGeometry(geometry);
@@ -893,7 +900,7 @@ TEST_F(MultipleGridBuilderFacadeTest_CreateMockAndSut, addSplit_calledAfterCreat
 {
     sut.createGrids(0);
     std::shared_ptr<Object> fineGrid = std::make_shared<Cuboid>(-0.25, -0.25, -0.25, 0.25, 0.25, 0.25);
-    EXPECT_THROW(sut.addDomainSplit(1.0, MultipleGridBuilderFacade::CoordDirection::x), std::runtime_error);
+    EXPECT_THROW(sut.addDomainSplit(1.0, Axis::x), std::runtime_error);
 }
 
 TEST(MultipleGridBuilderFacadeTest, createGrids_createGridsMoreThanOnceForSamePart_throws)
diff --git a/src/gpu/GridGenerator/utilities/communciationTest.cpp b/src/gpu/GridGenerator/utilities/communciationTest.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..eaecde97fab81c9ead136444179ec8ec7eb2d6ab
--- /dev/null
+++ b/src/gpu/GridGenerator/utilities/communciationTest.cpp
@@ -0,0 +1,39 @@
+#include <gmock/gmock.h>
+
+#include "communication.h"
+
+using namespace CommunicationDirections;
+
+TEST(communicationTest, isNegative)
+{
+    EXPECT_TRUE(isNegative(CommunicationDirection::MX));
+    EXPECT_TRUE(isNegative(CommunicationDirection::MY));
+    EXPECT_TRUE(isNegative(CommunicationDirection::MZ));
+    EXPECT_FALSE(isNegative(CommunicationDirection::PX));
+    EXPECT_FALSE(isNegative(CommunicationDirection::PY));
+    EXPECT_FALSE(isNegative(CommunicationDirection::PZ));
+}
+
+TEST(communicationTest, isPositive)
+{
+    EXPECT_TRUE(isPositive(CommunicationDirection::PX));
+    EXPECT_TRUE(isPositive(CommunicationDirection::PY));
+    EXPECT_TRUE(isPositive(CommunicationDirection::PZ));
+    EXPECT_FALSE(isPositive(CommunicationDirection::MX));
+    EXPECT_FALSE(isPositive(CommunicationDirection::MY));
+    EXPECT_FALSE(isPositive(CommunicationDirection::MZ));
+}
+
+TEST(communicationTest, getNegativeDirectionAlongAxis)
+{
+    EXPECT_THAT(getNegativeDirectionAlongAxis(Axis::x), CommunicationDirection::MX);
+    EXPECT_THAT(getNegativeDirectionAlongAxis(Axis::y), CommunicationDirection::MY);
+    EXPECT_THAT(getNegativeDirectionAlongAxis(Axis::z), CommunicationDirection::MZ);
+}
+
+TEST(communicationTest, getPositiveDirectionAlongAxis)
+{
+    EXPECT_THAT(getPositiveDirectionAlongAxis(Axis::x), CommunicationDirection::PX);
+    EXPECT_THAT(getPositiveDirectionAlongAxis(Axis::y), CommunicationDirection::PY);
+    EXPECT_THAT(getPositiveDirectionAlongAxis(Axis::z), CommunicationDirection::PZ);
+}
\ No newline at end of file
diff --git a/src/gpu/GridGenerator/utilities/communication.cpp b/src/gpu/GridGenerator/utilities/communication.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..cc9a52be979576fcde136e3a39aeea64e89067f0
--- /dev/null
+++ b/src/gpu/GridGenerator/utilities/communication.cpp
@@ -0,0 +1,49 @@
+#include "communication.h"
+
+using namespace CommunicationDirections;
+
+bool CommunicationDirections::isNegative(CommunicationDirection direction)
+{
+    return direction == CommunicationDirection::MX || direction == CommunicationDirection::MY ||
+           direction == CommunicationDirection::MZ;
+}
+
+bool CommunicationDirections::isPositive(CommunicationDirection direction)
+{
+    return direction == CommunicationDirection::PX || direction == CommunicationDirection::PY ||
+           direction == CommunicationDirection::PZ;
+}
+
+CommunicationDirection CommunicationDirections::getNegativeDirectionAlongAxis(Axis axis)
+{
+    switch (axis) {
+        case Axis::x:
+            return MX;
+            break;
+        case Axis::y:
+            return MY;
+            break;
+        case Axis::z:
+            return MZ;
+            break;
+        default:
+            throw std::runtime_error("Unknown coordinate direction" + axis::to_string(axis));
+    }
+}
+
+CommunicationDirection CommunicationDirections::getPositiveDirectionAlongAxis(Axis axis)
+{
+    switch (axis) {
+        case Axis::x:
+            return PX;
+            break;
+        case Axis::y:
+            return PY;
+            break;
+        case Axis::z:
+            return PZ;
+            break;
+        default:
+            throw std::runtime_error("Unknown coordinate direction" + axis::to_string(axis));
+    }
+}
\ No newline at end of file
diff --git a/src/gpu/GridGenerator/utilities/communication.h b/src/gpu/GridGenerator/utilities/communication.h
index f8c89f09dc3d80791cfeee0420964bfb1a463a0d..fd1b1feefc9ffc8e16f7e4a8f1ba94bd15d3c890 100644
--- a/src/gpu/GridGenerator/utilities/communication.h
+++ b/src/gpu/GridGenerator/utilities/communication.h
@@ -26,25 +26,37 @@
 //  You should have received a copy of the GNU General Public License along
 //  with VirtualFluids (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>.
 //
-//! \file Communication.h
-//! \ingroup utilities
 //! \author Soeren Peters, Stephan Lenz
 //=======================================================================================
 #ifndef Communication_H
 #define Communication_H
 
+#include <basics/geometry3d/Axis.h>
+
 #include "grid/BoundaryConditions/Side.h"
 
 // has to have the same order as SideType in Side.h
-namespace CommunicationDirections {
-    enum CommunicationDirection{
-        MX = static_cast<int>(SideType::MX),
-        PX = static_cast<int>(SideType::PX),
-        MY = static_cast<int>(SideType::MY),
-        PY = static_cast<int>(SideType::PY),
-        MZ = static_cast<int>(SideType::MZ),
-        PZ = static_cast<int>(SideType::PZ)
-    };
-}
+namespace CommunicationDirections
+{
+enum CommunicationDirection {
+    MX = static_cast<int>(SideType::MX),
+    PX = static_cast<int>(SideType::PX),
+    MY = static_cast<int>(SideType::MY),
+    PY = static_cast<int>(SideType::PY),
+    MZ = static_cast<int>(SideType::MZ),
+    PZ = static_cast<int>(SideType::PZ)
+};
+
+bool isNegative(CommunicationDirection direction);
+bool isPositive(CommunicationDirection direction);
+
+const std::map<CommunicationDirection, Axis> communicationDirectionToAxes { { MX, Axis::x }, { PX, Axis::x },
+                                                                            { MY, Axis::y }, { PY, Axis::y },
+                                                                            { MZ, Axis::z }, { PZ, Axis::z } };
+
+CommunicationDirection getNegativeDirectionAlongAxis(Axis axis);
+CommunicationDirection getPositiveDirectionAlongAxis(Axis axis);
+
+} // namespace CommunicationDirections
 
 #endif // Communication_H
\ No newline at end of file