From 2d759b106a0470fd9d6efc08f0320467352ffce9 Mon Sep 17 00:00:00 2001
From: Jan Eitzinger <jan@moebiusband.org>
Date: Sun, 11 Feb 2024 21:47:53 +0100
Subject: [PATCH] Debug 2D mpi version. Does not yet work v2 upwards.

---
 BasicSolver/2D-mpi/Makefile       |   2 +
 BasicSolver/2D-mpi/config.mk      |   5 +-
 BasicSolver/2D-mpi/dcavity.par    |  10 +-
 BasicSolver/2D-mpi/src/comm-v1.c  |  80 +--------
 BasicSolver/2D-mpi/src/comm-v2.c  | 257 +++++++++++-----------------
 BasicSolver/2D-mpi/src/comm-v3.c  | 274 +++++++++++-------------------
 BasicSolver/2D-mpi/src/comm.c     | 123 ++++++++++++++
 BasicSolver/2D-mpi/src/comm.h     |  10 +-
 BasicSolver/2D-mpi/src/main.c     |  42 ++++-
 BasicSolver/2D-mpi/src/progress.c |  41 ++---
 BasicSolver/2D-mpi/src/solver.c   |  18 +-
 BasicSolver/2D-mpi/velocity.png   | Bin 9306 -> 0 bytes
 BasicSolver/2D-seq/Makefile       |   2 +
 BasicSolver/2D-seq/dcavity.par    |   2 +-
 BasicSolver/3D-mpi/src/progress.c |   1 -
 15 files changed, 406 insertions(+), 461 deletions(-)
 create mode 100644 BasicSolver/2D-mpi/src/comm.c
 delete mode 100644 BasicSolver/2D-mpi/velocity.png

diff --git a/BasicSolver/2D-mpi/Makefile b/BasicSolver/2D-mpi/Makefile
index 12ce373..3d972dc 100644
--- a/BasicSolver/2D-mpi/Makefile
+++ b/BasicSolver/2D-mpi/Makefile
@@ -48,6 +48,8 @@ clean:
 distclean: clean
 	$(info ===>  DIST CLEAN)
 	@rm -f $(TARGET)
+	@rm -f *.dat
+	@rm -f *.png
 
 info:
 	$(info $(CFLAGS))
diff --git a/BasicSolver/2D-mpi/config.mk b/BasicSolver/2D-mpi/config.mk
index 01c6298..f9ec1d1 100644
--- a/BasicSolver/2D-mpi/config.mk
+++ b/BasicSolver/2D-mpi/config.mk
@@ -2,11 +2,12 @@
 TAG ?= CLANG
 ENABLE_MPI ?= true
 ENABLE_OPENMP ?= false
-COMM_TYPE ?= v3
+COMM_TYPE ?= v2
 
 #Feature options
 OPTIONS +=  -DARRAY_ALIGNMENT=64
-#OPTIONS +=  -DVERBOSE
+OPTIONS +=  -DVERBOSE
+# OPTIONS +=  -DTEST
 #OPTIONS +=  -DVERBOSE_AFFINITY
 #OPTIONS +=  -DVERBOSE_DATASIZE
 #OPTIONS +=  -DVERBOSE_TIMER
diff --git a/BasicSolver/2D-mpi/dcavity.par b/BasicSolver/2D-mpi/dcavity.par
index b4013d6..2c81044 100644
--- a/BasicSolver/2D-mpi/dcavity.par
+++ b/BasicSolver/2D-mpi/dcavity.par
@@ -15,7 +15,7 @@ bcRight    1			#
 gx    0.0			# Body forces (e.g. gravity)
 gy    0.0			#
 
-re    10.0		    # Reynolds number
+re    100.0		    # Reynolds number
 
 u_init    0.0		# initial value for velocity in x-direction
 v_init    0.0		# initial value for velocity in y-direction
@@ -26,13 +26,13 @@ p_init    0.0		# initial value for pressure
 
 xlength    1.0		# domain size in x-direction
 ylength    1.0		# domain size in y-direction
-imax       100		# number of interior cells in x-direction
-jmax       100		# number of interior cells in y-direction
+imax       80		# number of interior cells in x-direction
+jmax       80		# number of interior cells in y-direction
 
 # Time Data:
 # ---------
 
-te      5.0		# final time
+te      10.0		# final time
 dt     0.02	    # time stepsize
 tau     0.5		# safety factor for time stepsize control (<0 constant delt)
 
@@ -41,6 +41,6 @@ tau     0.5		# safety factor for time stepsize control (<0 constant delt)
 
 itermax  1000		# maximal number of pressure iteration in one time step
 eps      0.001		# stopping tolerance for pressure iteration
-omg      1.7		# relaxation parameter for SOR iteration
+omg      1.9		# relaxation parameter for SOR iteration
 gamma    0.9		# upwind differencing factor gamma
 #===============================================================================
diff --git a/BasicSolver/2D-mpi/src/comm-v1.c b/BasicSolver/2D-mpi/src/comm-v1.c
index 8dfd8c0..8f0ea3e 100644
--- a/BasicSolver/2D-mpi/src/comm-v1.c
+++ b/BasicSolver/2D-mpi/src/comm-v1.c
@@ -1,21 +1,15 @@
 /*
- * Copyright (C)  NHR@FAU, University Erlangen-Nuremberg.
+ * Copyright (C) NHR@FAU, University Erlangen-Nuremberg.
  * All rights reserved. This file is part of nusif-solver.
  * Use of this source code is governed by a MIT style
  * license that can be found in the LICENSE file.
  */
-#include <stdio.h>
 #include <stdlib.h>
 
 #include "comm.h"
 
-#if defined(_MPI)
+#ifdef _MPI
 // subroutines local to this module
-static int sizeOfRank(int rank, int size, int N)
-{
-    return N / size + ((N % size > rank) ? 1 : 0);
-}
-
 static int sum(int* sizes, int position)
 {
     int sum = 0;
@@ -49,20 +43,9 @@ static void gatherArray(
 #endif // defined _MPI
 
 // exported subroutines
-void commReduction(double* v, int op)
-{
-#if defined(_MPI)
-    if (op == MAX) {
-        MPI_Allreduce(MPI_IN_PLACE, v, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD);
-    } else if (op == SUM) {
-        MPI_Allreduce(MPI_IN_PLACE, v, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
-    }
-#endif
-}
-
 int commIsBoundary(Comm* c, int direction)
 {
-#if defined(_MPI)
+#ifdef _MPI
     switch (direction) {
     case LEFT:
         return 1;
@@ -84,7 +67,7 @@ int commIsBoundary(Comm* c, int direction)
 
 void commExchange(Comm* c, double* grid)
 {
-#if defined(_MPI)
+#ifdef _MPI
     MPI_Request requests[4] = { MPI_REQUEST_NULL,
         MPI_REQUEST_NULL,
         MPI_REQUEST_NULL,
@@ -116,7 +99,7 @@ void commExchange(Comm* c, double* grid)
 
 void commShift(Comm* c, double* f, double* g)
 {
-#if defined(_MPI)
+#ifdef _MPI
     MPI_Request requests[2] = { MPI_REQUEST_NULL, MPI_REQUEST_NULL };
 
     /* shift G */
@@ -143,7 +126,6 @@ void commShift(Comm* c, double* f, double* g)
 #endif
 }
 
-// TODO Merge with seq
 void commCollectResult(Comm* c,
     double* ug,
     double* vg,
@@ -154,6 +136,7 @@ void commCollectResult(Comm* c,
     int jmax,
     int imax)
 {
+#ifdef _MPI
     int *rcvCounts, *displs;
     int cnt = c->jmaxLocal * (imax + 2);
 
@@ -183,49 +166,12 @@ void commCollectResult(Comm* c,
     gatherArray(c, cnt, rcvCounts, displs, p, pg);
     gatherArray(c, cnt, rcvCounts, displs, u, ug);
     gatherArray(c, cnt, rcvCounts, displs, v, vg);
-}
-
-void commPrintConfig(Comm* c)
-{
-#if defined(_MPI)
-    fflush(stdout);
-    MPI_Barrier(MPI_COMM_WORLD);
-    if (commIsMaster(c)) {
-        printf("Communication setup:\n");
-    }
-
-    for (int i = 0; i < c->size; i++) {
-        if (i == c->rank) {
-            printf("\tRank %d of %d\n", c->rank, c->size);
-            printf("\tNeighbours (bottom, top, left, right): %d %d, %d, %d\n",
-                c->neighbours[BOTTOM],
-                c->neighbours[TOP],
-                c->neighbours[LEFT],
-                c->neighbours[RIGHT]);
-            printf("\tCoordinates (j,i) %d %d\n", c->coords[JDIM], c->coords[IDIM]);
-            printf("\tLocal domain size (j,i) %dx%d\n", c->jmaxLocal, c->imaxLocal);
-            fflush(stdout);
-        }
-    }
-    MPI_Barrier(MPI_COMM_WORLD);
-#endif
-}
-
-void commInit(Comm* c, int argc, char** argv)
-{
-#if defined(_MPI)
-    MPI_Init(&argc, &argv);
-    MPI_Comm_rank(MPI_COMM_WORLD, &(c->rank));
-    MPI_Comm_size(MPI_COMM_WORLD, &(c->size));
-#else
-    c->rank = 0;
-    c->size = 1;
 #endif
 }
 
 void commPartition(Comm* c, int jmax, int imax)
 {
-#if defined(_MPI)
+#ifdef _MPI
     c->imaxLocal = imax;
     c->jmaxLocal = sizeOfRank(c->rank, c->size, jmax);
 #else
@@ -233,15 +179,3 @@ void commPartition(Comm* c, int jmax, int imax)
     c->jmaxLocal = jmax;
 #endif
 }
-
-void commFinalize(Comm* c)
-{
-#if defined(_MPI)
-    for (int i = 0; i < NDIRS; i++) {
-        MPI_Type_free(&c->sbufferTypes[i]);
-        MPI_Type_free(&c->rbufferTypes[i]);
-    }
-
-    MPI_Finalize();
-#endif
-}
diff --git a/BasicSolver/2D-mpi/src/comm-v2.c b/BasicSolver/2D-mpi/src/comm-v2.c
index 42abb66..79b6687 100644
--- a/BasicSolver/2D-mpi/src/comm-v2.c
+++ b/BasicSolver/2D-mpi/src/comm-v2.c
@@ -1,42 +1,28 @@
 /*
- * Copyright (C)  NHR@FAU, University Erlangen-Nuremberg.
+ * Copyright (C) NHR@FAU, University Erlangen-Nuremberg.
  * All rights reserved. This file is part of nusif-solver.
  * Use of this source code is governed by a MIT style
  * license that can be found in the LICENSE file.
  */
-#include <stddef.h>
 #include <stdio.h>
 #include <stdlib.h>
 
 #include "comm.h"
 
-#if defined(_MPI)
+#ifdef _MPI
 // subroutines local to this module
-static int sizeOfRank(int rank, int size, int N)
+static int sum(int* sizes, int position)
 {
-    return N / size + ((N % size > rank) ? 1 : 0);
+    int sum = 0;
+
+    for (int i = 0; i < position; i++) {
+        sum += sizes[i];
+    }
+
+    return sum;
 }
 
-static void setupCommunication(Comm* c, int direction, int layer)
-{
-    MPI_Datatype type;
-    size_t dblsize = sizeof(double);
-    int imaxLocal  = c->imaxLocal;
-    int jmaxLocal  = c->jmaxLocal;
-    int sizes[NDIMS];
-    int subSizes[NDIMS];
-    int starts[NDIMS];
-    int offset = 0;
-}
-
-static void assembleResult(Comm* c,
-    double* src,
-    double* dst,
-    int imaxLocal[],
-    int jmaxLocal[],
-    int offset[],
-    int jmax,
-    int imax)
+static void assembleResult(Comm* c, double* src, double* dst, int jmax, int imax)
 {
     MPI_Request* requests;
     int numRequests = 1;
@@ -49,11 +35,27 @@ static void assembleResult(Comm* c,
 
     requests = (MPI_Request*)malloc(numRequests * sizeof(MPI_Request));
 
-    /* all ranks send their bulk array */
+    /* all ranks send their bulk array, but including external boundary layer */
     MPI_Datatype bulkType;
-    int oldSizes[NDIMS] = { c->jmaxLocal + 2, c->imaxLocal + 2 };
-    int newSizes[NDIMS] = { c->jmaxLocal, c->imaxLocal };
+    int oldSizes[NDIMS] = { c->imaxLocal + 2, c->jmaxLocal + 2 };
+    int newSizes[NDIMS] = { c->imaxLocal, c->jmaxLocal };
     int starts[NDIMS]   = { 1, 1 };
+
+    if (commIsBoundary(c, LEFT)) {
+        newSizes[IDIM] += 1;
+        starts[IDIM] = 0;
+    }
+    if (commIsBoundary(c, RIGHT)) {
+        newSizes[IDIM] += 1;
+    }
+    if (commIsBoundary(c, BOTTOM)) {
+        newSizes[JDIM] += 1;
+        starts[JDIM] = 0;
+    }
+    if (commIsBoundary(c, TOP)) {
+        newSizes[JDIM] += 1;
+    }
+
     MPI_Type_create_subarray(NDIMS,
         oldSizes,
         newSizes,
@@ -62,16 +64,23 @@ static void assembleResult(Comm* c,
         MPI_DOUBLE,
         &bulkType);
     MPI_Type_commit(&bulkType);
-
     MPI_Isend(src, 1, bulkType, 0, 0, c->comm, &requests[0]);
 
+    int newSizesI[c->size];
+    int newSizesJ[c->size];
+    MPI_Gather(&newSizes[IDIM], 1, MPI_INT, newSizesI, 1, MPI_INT, 0, MPI_COMM_WORLD);
+    MPI_Gather(&newSizes[JDIM], 1, MPI_INT, newSizesJ, 1, MPI_INT, 0, MPI_COMM_WORLD);
+
     /* rank 0 assembles the subdomains */
     if (c->rank == 0) {
         for (int i = 0; i < c->size; i++) {
             MPI_Datatype domainType;
-            int oldSizes[NDIMS] = { jmax, imax };
-            int newSizes[NDIMS] = { jmaxLocal[i], imaxLocal[i] };
-            int starts[NDIMS]   = { offset[i * NDIMS + JDIM], offset[i * NDIMS + IDIM] };
+            int oldSizes[NDIMS] = { imax + 2, jmax + 2 };
+            int newSizes[NDIMS] = { newSizesI[i], newSizesJ[i] };
+            int coords[NDIMS];
+            MPI_Cart_coords(c->comm, i, NDIMS, coords);
+            int starts[NDIMS] = { sum(newSizesI, coords[IDIM]),
+                sum(newSizesJ, coords[JDIM]) };
             MPI_Type_create_subarray(NDIMS,
                 oldSizes,
                 newSizes,
@@ -87,34 +96,12 @@ static void assembleResult(Comm* c,
 
     MPI_Waitall(numRequests, requests, MPI_STATUSES_IGNORE);
 }
-
-static int sum(int* sizes, int position)
-{
-    int sum = 0;
-
-    for (int i = 0; i < position; i++) {
-        sum += sizes[i];
-    }
-
-    return sum;
-}
 #endif // defined _MPI
 
 // exported subroutines
-void commReduction(double* v, int op)
-{
-#if defined(_MPI)
-    if (op == MAX) {
-        MPI_Allreduce(MPI_IN_PLACE, v, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD);
-    } else if (op == SUM) {
-        MPI_Allreduce(MPI_IN_PLACE, v, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
-    }
-#endif
-}
-
 int commIsBoundary(Comm* c, int direction)
 {
-#if defined(_MPI)
+#ifdef _MPI
     switch (direction) {
     case LEFT:
         return c->coords[IDIM] == 0;
@@ -136,63 +123,51 @@ int commIsBoundary(Comm* c, int direction)
 
 void commExchange(Comm* c, double* grid)
 {
-#if defined(_MPI)
-    double* buf[8];
+#ifdef _MPI
+    double* sbuf[NDIRS];
+    double* rbuf[NDIRS];
     MPI_Request requests[8];
     for (int i = 0; i < 8; i++)
         requests[i] = MPI_REQUEST_NULL;
 
-    buf[0] = grid + 1;                                           // recv bottom
-    buf[1] = grid + (c->imaxLocal + 2) + 1;                      // send bottom
-    buf[2] = grid + (c->jmaxLocal + 1) * (c->imaxLocal + 2) + 1; // recv top
-    buf[3] = grid + (c->jmaxLocal) * (c->imaxLocal + 2) + 1;     // send top
-    buf[4] = grid + (c->imaxLocal + 2);                          // recv left
-    buf[5] = grid + (c->imaxLocal + 2) + 1;                      // send left
-    buf[6] = grid + (c->imaxLocal + 2) + (c->imaxLocal + 1);     // recv right
-    buf[7] = grid + (c->imaxLocal + 2) + (c->imaxLocal);         // send right
+    rbuf[LEFT]   = grid + (c->imaxLocal + 2);
+    rbuf[RIGHT]  = grid + (c->imaxLocal + 2) + (c->imaxLocal + 1);
+    rbuf[BOTTOM] = grid + 1;
+    rbuf[TOP]    = grid + (c->jmaxLocal + 1) * (c->imaxLocal + 2) + 1;
+    sbuf[LEFT]   = grid + (c->imaxLocal + 2) + 1;
+    sbuf[RIGHT]  = grid + (c->imaxLocal + 2) + (c->imaxLocal);
+    sbuf[BOTTOM] = grid + (c->imaxLocal + 2) + 1;
+    sbuf[TOP]    = grid + (c->jmaxLocal) * (c->imaxLocal + 2) + 1;
+    // buf[0] = grid + (c->imaxLocal + 2);                          // recv left
+    // buf[1] = grid + (c->imaxLocal + 2) + 1;                      // send left
+    // buf[2] = grid + (c->imaxLocal + 2) + (c->imaxLocal + 1);     // recv right
+    // buf[3] = grid + (c->imaxLocal + 2) + (c->imaxLocal);         // send right
+    // buf[4] = grid + 1;                                           // recv bottom
+    // buf[5] = grid + (c->imaxLocal + 2) + 1;                      // send bottom
+    // buf[6] = grid + (c->jmaxLocal + 1) * (c->imaxLocal + 2) + 1; // recv top
+    // buf[7] = grid + (c->jmaxLocal) * (c->imaxLocal + 2) + 1;     // send top
 
-    for (int i = 2; i < 4; i++) {
+    for (int i = 0; i < NDIRS; i++) {
         int tag = 0;
         if (c->neighbours[i] != MPI_PROC_NULL) {
+            // printf("DEBUG: Rank %d - SendRecv with %d\n", c->rank, c->neighbours[i]);
             tag = c->neighbours[i];
         }
-        /* exchange ghost cells with bottom/top neighbor */
-        MPI_Irecv(buf[i * 2],
+        MPI_Irecv(rbuf[i],
             1,
-            c->rbufferTypes[0],
+            c->bufferTypes[i],
             c->neighbours[i],
             tag,
             c->comm,
             &requests[i * 2]);
-        MPI_Isend(buf[(i * 2) + 1],
+        MPI_Isend(sbuf[i],
             1,
-            c->rbufferTypes[0],
+            c->bufferTypes[i],
             c->neighbours[i],
             c->rank,
             c->comm,
             &requests[i * 2 + 1]);
     }
-    for (int i = 0; i < 2; i++) {
-        int tag = 0;
-        if (c->neighbours[i] != MPI_PROC_NULL) {
-            tag = c->neighbours[i];
-        }
-        /* exchange ghost cells with left/right neighbor */
-        MPI_Irecv(buf[i * 2 + 4],
-            1,
-            c->sbufferTypes[0],
-            c->neighbours[i],
-            tag,
-            c->comm,
-            &requests[i * 2 + 4]);
-        MPI_Isend(buf[i * 2 + 5],
-            1,
-            c->sbufferTypes[0],
-            c->neighbours[i],
-            c->rank,
-            c->comm,
-            &requests[(i * 2) + 5]);
-    }
 
     MPI_Waitall(8, requests, MPI_STATUSES_IGNORE);
 #endif
@@ -200,7 +175,7 @@ void commExchange(Comm* c, double* grid)
 
 void commShift(Comm* c, double* f, double* g)
 {
-#if defined(_MPI)
+#ifdef _MPI
     MPI_Request requests[4] = { MPI_REQUEST_NULL,
         MPI_REQUEST_NULL,
         MPI_REQUEST_NULL,
@@ -211,7 +186,7 @@ void commShift(Comm* c, double* f, double* g)
     /* receive ghost cells from bottom neighbor */
     MPI_Irecv(buf,
         1,
-        c->rbufferTypes[0],
+        c->bufferTypes[BOTTOM],
         c->neighbours[BOTTOM],
         0,
         c->comm,
@@ -219,16 +194,28 @@ void commShift(Comm* c, double* f, double* g)
 
     buf = g + (c->jmaxLocal) * (c->imaxLocal + 2) + 1;
     /* send ghost cells to top neighbor */
-    MPI_Isend(buf, 1, c->rbufferTypes[0], c->neighbours[TOP], 0, c->comm, &requests[1]);
+    MPI_Isend(buf, 1, c->bufferTypes[TOP], c->neighbours[TOP], 0, c->comm, &requests[1]);
 
     /* shift F */
     buf = f + (c->imaxLocal + 2);
     /* receive ghost cells from left neighbor */
-    MPI_Irecv(buf, 1, c->sbufferTypes[0], c->neighbours[LEFT], 1, c->comm, &requests[2]);
+    MPI_Irecv(buf,
+        1,
+        c->bufferTypes[LEFT],
+        c->neighbours[LEFT],
+        1,
+        c->comm,
+        &requests[2]);
 
-    buf = f + (c->imaxLocal + 2) + (c->imaxLocal);
+    buf = f + (c->imaxLocal + 2) + (c->imaxLocal + 1);
     /* send ghost cells to right neighbor */
-    MPI_Isend(buf, 1, c->sbufferTypes[0], c->neighbours[RIGHT], 1, c->comm, &requests[3]);
+    MPI_Isend(buf,
+        1,
+        c->bufferTypes[RIGHT],
+        c->neighbours[RIGHT],
+        1,
+        c->comm,
+        &requests[3]);
 
     MPI_Waitall(4, requests, MPI_STATUSES_IGNORE);
 #endif
@@ -244,6 +231,8 @@ void commCollectResult(Comm* c,
     int jmax,
     int imax)
 {
+#ifdef _MPI
+
     int offset[c->size * NDIMS];
     int imaxLocal[c->size];
     int jmaxLocal[c->size];
@@ -270,56 +259,19 @@ void commCollectResult(Comm* c,
     }
 
     /* collect P */
-    assembleResult(c, p, pg, imaxLocal, jmaxLocal, offset, jmax, imax);
+    assembleResult(c, p, pg, jmax, imax);
 
     /* collect U */
-    assembleResult(c, u, ug, imaxLocal, jmaxLocal, offset, jmax, imax);
+    assembleResult(c, u, ug, jmax, imax);
 
     /* collect V */
-    assembleResult(c, v, vg, imaxLocal, jmaxLocal, offset, jmax, imax);
-}
-
-void commPrintConfig(Comm* c)
-{
-#if defined(_MPI)
-    fflush(stdout);
-    MPI_Barrier(MPI_COMM_WORLD);
-    if (commIsMaster(c)) {
-        printf("Communication setup:\n");
-    }
-
-    for (int i = 0; i < c->size; i++) {
-        if (i == c->rank) {
-            printf("\tRank %d of %d\n", c->rank, c->size);
-            printf("\tNeighbours (bottom, top, left, right): %d %d, %d, %d\n",
-                c->neighbours[BOTTOM],
-                c->neighbours[TOP],
-                c->neighbours[LEFT],
-                c->neighbours[RIGHT]);
-            printf("\tCoordinates (j,i) %d %d\n", c->coords[JDIM], c->coords[IDIM]);
-            printf("\tLocal domain size (j,i) %dx%d\n", c->jmaxLocal, c->imaxLocal);
-            fflush(stdout);
-        }
-    }
-    MPI_Barrier(MPI_COMM_WORLD);
-#endif
-}
-
-void commInit(Comm* c, int argc, char** argv)
-{
-#if defined(_MPI)
-    MPI_Init(&argc, &argv);
-    MPI_Comm_rank(MPI_COMM_WORLD, &(c->rank));
-    MPI_Comm_size(MPI_COMM_WORLD, &(c->size));
-#else
-    c->rank = 0;
-    c->size = 1;
+    assembleResult(c, v, vg, jmax, imax);
 #endif
 }
 
 void commPartition(Comm* c, int jmax, int imax)
 {
-#if defined(_MPI)
+#ifdef _MPI
     int dims[NDIMS]    = { 0, 0 };
     int periods[NDIMS] = { 0, 0 };
     MPI_Dims_create(c->size, NDIMS, dims);
@@ -331,25 +283,20 @@ void commPartition(Comm* c, int jmax, int imax)
     c->imaxLocal = sizeOfRank(c->rank, dims[IDIM], imax);
     c->jmaxLocal = sizeOfRank(c->rank, dims[JDIM], jmax);
 
-    MPI_Type_contiguous(c->imaxLocal, MPI_DOUBLE, &c->rbufferTypes[0]);
-    MPI_Type_commit(&c->rbufferTypes[0]);
+    MPI_Datatype jBufferType;
+    MPI_Type_contiguous(c->imaxLocal, MPI_DOUBLE, &jBufferType);
+    MPI_Type_commit(&jBufferType);
 
-    MPI_Type_vector(c->jmaxLocal, 1, c->imaxLocal + 2, MPI_DOUBLE, &c->sbufferTypes[0]);
-    MPI_Type_commit(&c->sbufferTypes[0]);
+    MPI_Datatype iBufferType;
+    MPI_Type_vector(c->jmaxLocal, 1, c->imaxLocal + 2, MPI_DOUBLE, &iBufferType);
+    MPI_Type_commit(&iBufferType);
+
+    c->bufferTypes[LEFT]   = iBufferType;
+    c->bufferTypes[RIGHT]  = iBufferType;
+    c->bufferTypes[BOTTOM] = jBufferType;
+    c->bufferTypes[TOP]    = jBufferType;
 #else
     c->imaxLocal = imax;
     c->jmaxLocal = jmax;
 #endif
 }
-
-void commFinalize(Comm* c)
-{
-#if defined(_MPI)
-    for (int i = 0; i < NDIRS; i++) {
-        MPI_Type_free(&c->sbufferTypes[i]);
-        MPI_Type_free(&c->rbufferTypes[i]);
-    }
-
-    MPI_Finalize();
-#endif
-}
diff --git a/BasicSolver/2D-mpi/src/comm-v3.c b/BasicSolver/2D-mpi/src/comm-v3.c
index b39cc98..e2e9d8c 100644
--- a/BasicSolver/2D-mpi/src/comm-v3.c
+++ b/BasicSolver/2D-mpi/src/comm-v3.c
@@ -1,5 +1,5 @@
 /*
- * Copyright (C)  NHR@FAU, University Erlangen-Nuremberg.
+ * Copyright (C) NHR@FAU, University Erlangen-Nuremberg.
  * All rights reserved. This file is part of nusif-solver.
  * Use of this source code is governed by a MIT style
  * license that can be found in the LICENSE file.
@@ -10,82 +10,20 @@
 
 #include "comm.h"
 
-#if defined(_MPI)
+#ifdef _MPI
 // subroutines local to this module
-static int sizeOfRank(int rank, int size, int N)
+static int sum(int* sizes, int position)
 {
-    return N / size + ((N % size > rank) ? 1 : 0);
+    int sum = 0;
+
+    for (int i = 0; i < position; i++) {
+        sum += sizes[i];
+    }
+
+    return sum;
 }
 
-static void setupCommunication(Comm* c, int direction, int layer)
-{
-    MPI_Datatype type;
-    size_t dblsize = sizeof(double);
-    int imaxLocal  = c->imaxLocal;
-    int jmaxLocal  = c->jmaxLocal;
-    int sizes[NDIMS];
-    int subSizes[NDIMS];
-    int starts[NDIMS];
-    int offset = 0;
-
-    sizes[IDIM] = imaxLocal + 2;
-    sizes[JDIM] = jmaxLocal + 2;
-
-    if (layer == HALO) {
-        offset = 1;
-    }
-
-    switch (direction) {
-    case LEFT:
-        subSizes[IDIM] = 1;
-        subSizes[JDIM] = jmaxLocal;
-        starts[IDIM]   = 1 - offset;
-        starts[JDIM]   = 1;
-        break;
-    case RIGHT:
-        subSizes[IDIM] = 1;
-        subSizes[JDIM] = jmaxLocal;
-        starts[IDIM]   = imaxLocal + offset;
-        starts[JDIM]   = 1;
-        break;
-    case BOTTOM:
-        subSizes[IDIM] = imaxLocal;
-        subSizes[JDIM] = 1;
-        starts[IDIM]   = 1;
-        starts[JDIM]   = 1 - offset;
-        break;
-    case TOP:
-        subSizes[IDIM] = imaxLocal;
-        subSizes[JDIM] = 1;
-        starts[IDIM]   = 1;
-        starts[JDIM]   = jmaxLocal + offset;
-        break;
-    }
-
-    MPI_Type_create_subarray(NDIMS,
-        sizes,
-        subSizes,
-        starts,
-        MPI_ORDER_C,
-        MPI_DOUBLE,
-        &type);
-    MPI_Type_commit(&type);
-
-    if (layer == HALO) {
-        c->rbufferTypes[direction] = type;
-    } else if (layer == BULK) {
-        c->sbufferTypes[direction] = type;
-    }
-}
-
-static void assembleResult(Comm* c,
-    double* src,
-    double* dst,
-    int imaxLocal[],
-    int jmaxLocal[],
-    int offset[],
-    int jmax,
-    int imax)
+static void assembleResult(Comm* c, double* src, double* dst, int jmax, int imax)
 {
     MPI_Request* requests;
     int numRequests = 1;
@@ -98,11 +36,27 @@ static void assembleResult(Comm* c,
 
     requests = (MPI_Request*)malloc(numRequests * sizeof(MPI_Request));
 
-    /* all ranks send their bulk array */
+    /* all ranks send their bulk array, but including external boundary layer */
     MPI_Datatype bulkType;
     int oldSizes[NDIMS] = { c->jmaxLocal + 2, c->imaxLocal + 2 };
     int newSizes[NDIMS] = { c->jmaxLocal, c->imaxLocal };
     int starts[NDIMS]   = { 1, 1 };
+
+    if (commIsBoundary(c, LEFT)) {
+        newSizes[IDIM] += 1;
+        starts[IDIM] = 0;
+    }
+    if (commIsBoundary(c, RIGHT)) {
+        newSizes[IDIM] += 1;
+    }
+    if (commIsBoundary(c, BOTTOM)) {
+        newSizes[JDIM] += 1;
+        starts[JDIM] = 0;
+    }
+    if (commIsBoundary(c, TOP)) {
+        newSizes[JDIM] += 1;
+    }
+
     MPI_Type_create_subarray(NDIMS,
         oldSizes,
         newSizes,
@@ -111,16 +65,23 @@ static void assembleResult(Comm* c,
         MPI_DOUBLE,
         &bulkType);
     MPI_Type_commit(&bulkType);
-
     MPI_Isend(src, 1, bulkType, 0, 0, c->comm, &requests[0]);
 
+    int newSizesI[c->size];
+    int newSizesJ[c->size];
+    MPI_Gather(&newSizes[IDIM], 1, MPI_INT, newSizesI, 1, MPI_INT, 0, MPI_COMM_WORLD);
+    MPI_Gather(&newSizes[JDIM], 1, MPI_INT, newSizesJ, 1, MPI_INT, 0, MPI_COMM_WORLD);
+
     /* rank 0 assembles the subdomains */
     if (c->rank == 0) {
         for (int i = 0; i < c->size; i++) {
             MPI_Datatype domainType;
-            int oldSizes[NDIMS] = { jmax, imax };
-            int newSizes[NDIMS] = { jmaxLocal[i], imaxLocal[i] };
-            int starts[NDIMS]   = { offset[i * NDIMS + JDIM], offset[i * NDIMS + IDIM] };
+            int oldSizes[NDIMS] = { jmax + 2, imax + 2 };
+            int newSizes[NDIMS] = { newSizesJ[i], newSizesI[i] };
+            int coords[NDIMS];
+            MPI_Cart_coords(c->comm, i, NDIMS, coords);
+            int starts[NDIMS] = { sum(newSizesJ, coords[JDIM]),
+                sum(newSizesI, coords[IDIM]) };
             MPI_Type_create_subarray(NDIMS,
                 oldSizes,
                 newSizes,
@@ -136,34 +97,12 @@ static void assembleResult(Comm* c,
 
     MPI_Waitall(numRequests, requests, MPI_STATUSES_IGNORE);
 }
-
-static int sum(int* sizes, int position)
-{
-    int sum = 0;
-
-    for (int i = 0; i < position; i++) {
-        sum += sizes[i];
-    }
-
-    return sum;
-}
 #endif // defined _MPI
 
 // exported subroutines
-void commReduction(double* v, int op)
-{
-#if defined(_MPI)
-    if (op == MAX) {
-        MPI_Allreduce(MPI_IN_PLACE, v, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD);
-    } else if (op == SUM) {
-        MPI_Allreduce(MPI_IN_PLACE, v, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
-    }
-#endif
-}
-
 int commIsBoundary(Comm* c, int direction)
 {
-#if defined(_MPI)
+#ifdef _MPI
     switch (direction) {
     case LEFT:
         return c->coords[IDIM] == 0;
@@ -185,25 +124,24 @@ int commIsBoundary(Comm* c, int direction)
 
 void commExchange(Comm* c, double* grid)
 {
-#if defined(_MPI)
-    int counts[NDIRS]      = { 1, 1, 1, 1 };
-    MPI_Aint displs[NDIRS] = { 0, 0, 0, 0 };
+#ifdef _MPI
+    int counts[NDIRS] = { 1, 1, 1, 1 };
 
     MPI_Neighbor_alltoallw(grid,
         counts,
-        displs,
-        c->sbufferTypes,
+        c->sdispls,
+        c->bufferTypes,
         grid,
         counts,
-        displs,
-        c->rbufferTypes,
+        c->rdispls,
+        c->bufferTypes,
         c->comm);
 #endif
 }
 
 void commShift(Comm* c, double* f, double* g)
 {
-#if defined(_MPI)
+#ifdef _MPI
     MPI_Request requests[4] = { MPI_REQUEST_NULL,
         MPI_REQUEST_NULL,
         MPI_REQUEST_NULL,
@@ -211,25 +149,35 @@ void commShift(Comm* c, double* f, double* g)
 
     /* shift G */
     /* receive ghost cells from bottom neighbor */
-    MPI_Irecv(g,
+    double* buf = g + 1;
+    MPI_Irecv(buf,
         1,
-        c->rbufferTypes[BOTTOM],
+        c->bufferTypes[BOTTOM],
         c->neighbours[BOTTOM],
         0,
         c->comm,
         &requests[0]);
 
     /* send ghost cells to top neighbor */
-    MPI_Isend(g, 1, c->sbufferTypes[TOP], c->neighbours[TOP], 0, c->comm, &requests[1]);
+    buf = g + (c->jmaxLocal) * (c->imaxLocal + 2) + 1;
+    MPI_Isend(buf, 1, c->bufferTypes[TOP], c->neighbours[TOP], 0, c->comm, &requests[1]);
 
     /* shift F */
     /* receive ghost cells from left neighbor */
-    MPI_Irecv(f, 1, c->rbufferTypes[LEFT], c->neighbours[LEFT], 1, c->comm, &requests[2]);
+    buf = f + (c->imaxLocal + 2);
+    MPI_Irecv(buf,
+        1,
+        c->bufferTypes[LEFT],
+        c->neighbours[LEFT],
+        1,
+        c->comm,
+        &requests[2]);
 
     /* send ghost cells to right neighbor */
-    MPI_Isend(f,
+    buf = f + (c->imaxLocal + 2) + (c->imaxLocal);
+    MPI_Isend(buf,
         1,
-        c->sbufferTypes[RIGHT],
+        c->bufferTypes[RIGHT],
         c->neighbours[RIGHT],
         1,
         c->comm,
@@ -239,7 +187,6 @@ void commShift(Comm* c, double* f, double* g)
 #endif
 }
 
-// TODO Merge with seq
 void commCollectResult(Comm* c,
     double* ug,
     double* vg,
@@ -250,6 +197,8 @@ void commCollectResult(Comm* c,
     int jmax,
     int imax)
 {
+#ifdef _MPI
+
     int offset[c->size * NDIMS];
     int imaxLocal[c->size];
     int jmaxLocal[c->size];
@@ -276,56 +225,19 @@ void commCollectResult(Comm* c,
     }
 
     /* collect P */
-    assembleResult(c, p, pg, imaxLocal, jmaxLocal, offset, jmax, imax);
+    assembleResult(c, p, pg, jmax, imax);
 
     /* collect U */
-    assembleResult(c, u, ug, imaxLocal, jmaxLocal, offset, jmax, imax);
+    assembleResult(c, u, ug, jmax, imax);
 
     /* collect V */
-    assembleResult(c, v, vg, imaxLocal, jmaxLocal, offset, jmax, imax);
-}
-
-void commPrintConfig(Comm* c)
-{
-#if defined(_MPI)
-    fflush(stdout);
-    MPI_Barrier(MPI_COMM_WORLD);
-    if (commIsMaster(c)) {
-        printf("Communication setup:\n");
-    }
-
-    for (int i = 0; i < c->size; i++) {
-        if (i == c->rank) {
-            printf("\tRank %d of %d\n", c->rank, c->size);
-            printf("\tNeighbours (bottom, top, left, right): %d %d, %d, %d\n",
-                c->neighbours[BOTTOM],
-                c->neighbours[TOP],
-                c->neighbours[LEFT],
-                c->neighbours[RIGHT]);
-            printf("\tCoordinates (j,i) %d %d\n", c->coords[JDIM], c->coords[IDIM]);
-            printf("\tLocal domain size (j,i) %dx%d\n", c->jmaxLocal, c->imaxLocal);
-            fflush(stdout);
-        }
-    }
-    MPI_Barrier(MPI_COMM_WORLD);
-#endif
-}
-
-void commInit(Comm* c, int argc, char** argv)
-{
-#if defined(_MPI)
-    MPI_Init(&argc, &argv);
-    MPI_Comm_rank(MPI_COMM_WORLD, &(c->rank));
-    MPI_Comm_size(MPI_COMM_WORLD, &(c->size));
-#else
-    c->rank = 0;
-    c->size = 1;
+    assembleResult(c, v, vg, jmax, imax);
 #endif
 }
 
 void commPartition(Comm* c, int jmax, int imax)
 {
-#if defined(_MPI)
+#ifdef _MPI
     int dims[NDIMS]    = { 0, 0 };
     int periods[NDIMS] = { 0, 0 };
     MPI_Dims_create(c->size, NDIMS, dims);
@@ -337,29 +249,35 @@ void commPartition(Comm* c, int jmax, int imax)
     c->imaxLocal = sizeOfRank(c->rank, dims[IDIM], imax);
     c->jmaxLocal = sizeOfRank(c->rank, dims[JDIM], jmax);
 
-    // setup buffer types for communication
-    setupCommunication(c, LEFT, BULK);
-    setupCommunication(c, LEFT, HALO);
-    setupCommunication(c, RIGHT, BULK);
-    setupCommunication(c, RIGHT, HALO);
-    setupCommunication(c, BOTTOM, BULK);
-    setupCommunication(c, BOTTOM, HALO);
-    setupCommunication(c, TOP, BULK);
-    setupCommunication(c, TOP, HALO);
+    MPI_Datatype jBufferType;
+    MPI_Type_contiguous(c->imaxLocal, MPI_DOUBLE, &jBufferType);
+    MPI_Type_commit(&jBufferType);
+
+    MPI_Datatype iBufferType;
+    MPI_Type_vector(c->jmaxLocal, 1, c->imaxLocal + 2, MPI_DOUBLE, &iBufferType);
+    MPI_Type_commit(&iBufferType);
+
+    // in the order of the dimensions i->0, j->1
+    // first negative direction, then positive direction
+    size_t dblsize         = sizeof(double);
+    int imaxLocal          = c->imaxLocal;
+    int jmaxLocal          = c->jmaxLocal;
+    c->bufferTypes[LEFT]   = iBufferType;
+    c->bufferTypes[RIGHT]  = iBufferType;
+    c->bufferTypes[BOTTOM] = jBufferType;
+    c->bufferTypes[TOP]    = jBufferType;
+
+    c->sdispls[LEFT]   = ((imaxLocal + 2) + 1) * dblsize;             // send left
+    c->sdispls[RIGHT]  = ((imaxLocal + 2) + imaxLocal) * dblsize;     // send right
+    c->sdispls[BOTTOM] = ((imaxLocal + 2) + 1) * dblsize;             // send bottom
+    c->sdispls[TOP]    = (jmaxLocal * (imaxLocal + 2) + 1) * dblsize; // send top
+
+    c->rdispls[LEFT]   = (imaxLocal + 2) * dblsize;                         // recv left
+    c->rdispls[RIGHT]  = ((imaxLocal + 2) + (imaxLocal + 1)) * dblsize;     // recv right
+    c->rdispls[BOTTOM] = 1 * dblsize;                                       // recv bottom
+    c->rdispls[TOP]    = ((jmaxLocal + 1) * (imaxLocal + 2) + 1) * dblsize; // recv top
 #else
     c->imaxLocal = imax;
     c->jmaxLocal = jmax;
 #endif
 }
-
-void commFinalize(Comm* c)
-{
-#if defined(_MPI)
-    for (int i = 0; i < NDIRS; i++) {
-        MPI_Type_free(&c->sbufferTypes[i]);
-        MPI_Type_free(&c->rbufferTypes[i]);
-    }
-
-    MPI_Finalize();
-#endif
-}
diff --git a/BasicSolver/2D-mpi/src/comm.c b/BasicSolver/2D-mpi/src/comm.c
new file mode 100644
index 0000000..5568aa4
--- /dev/null
+++ b/BasicSolver/2D-mpi/src/comm.c
@@ -0,0 +1,123 @@
+/*
+ * Copyright (C)  NHR@FAU, University Erlangen-Nuremberg.
+ * All rights reserved. This file is part of nusif-solver.
+ * Use of this source code is governed by a MIT style
+ * license that can be found in the LICENSE file.
+ */
+#include <stdio.h>
+#include <stdlib.h>
+
+#include "comm.h"
+
+// subroutines local to this module
+int sizeOfRank(int rank, int size, int N)
+{
+    return N / size + ((N % size > rank) ? 1 : 0);
+}
+
+void commReduction(double* v, int op)
+{
+#ifdef _MPI
+    if (op == MAX) {
+        MPI_Allreduce(MPI_IN_PLACE, v, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD);
+    } else if (op == SUM) {
+        MPI_Allreduce(MPI_IN_PLACE, v, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
+    }
+#endif
+}
+
+void commPrintConfig(Comm* c)
+{
+#ifdef _MPI
+    fflush(stdout);
+    MPI_Barrier(MPI_COMM_WORLD);
+    if (commIsMaster(c)) {
+        printf("Communication setup:\n");
+    }
+
+    for (int i = 0; i < c->size; i++) {
+        if (i == c->rank) {
+            printf("\tRank %d of %d\n", c->rank, c->size);
+            printf("\tNeighbours (bottom, top, left, right): %d %d, %d, %d\n",
+                c->neighbours[BOTTOM],
+                c->neighbours[TOP],
+                c->neighbours[LEFT],
+                c->neighbours[RIGHT]);
+            printf("\tCoordinates (i,j) %d %d\n", c->coords[IDIM], c->coords[JDIM]);
+            printf("\tLocal domain size (i,j) %dx%d\n", c->imaxLocal, c->jmaxLocal);
+            fflush(stdout);
+        }
+    }
+    MPI_Barrier(MPI_COMM_WORLD);
+#endif
+}
+
+void commInit(Comm* c, int argc, char** argv)
+{
+#ifdef _MPI
+    MPI_Init(&argc, &argv);
+    MPI_Comm_rank(MPI_COMM_WORLD, &(c->rank));
+    MPI_Comm_size(MPI_COMM_WORLD, &(c->size));
+#else
+    c->rank = 0;
+    c->size = 1;
+#endif
+}
+
+void commTestInit(Comm* c, double* p, double* f, double* g)
+{
+    int imax = c->imaxLocal;
+    int jmax = c->jmaxLocal;
+    int rank = c->rank;
+
+    for (int j = 0; j < jmax + 2; j++) {
+        for (int i = 0; i < imax + 2; i++) {
+            p[j * (imax + 2) + i] = rank;
+            f[j * (imax + 2) + i] = rank;
+            g[j * (imax + 2) + i] = rank;
+        }
+    }
+}
+
+static void testWriteFile(char* filename, double* grid, int imax, int jmax)
+{
+    FILE* fp = fopen(filename, "w");
+
+    if (fp == NULL) {
+        printf("Error!\n");
+        exit(EXIT_FAILURE);
+    }
+
+    for (int j = 0; j < jmax + 2; j++) {
+        for (int i = 0; i < imax + 2; i++) {
+            fprintf(fp, "%f ", grid[j * (imax + 2) + i]);
+        }
+        fprintf(fp, "\n");
+    }
+
+    fclose(fp);
+}
+
+void commTestWrite(Comm* c, double* p, double* f, double* g)
+{
+    int imax = c->imaxLocal;
+    int jmax = c->jmaxLocal;
+    int rank = c->rank;
+
+    char filename[30];
+    snprintf(filename, 30, "ptest-%d.dat", rank);
+    testWriteFile(filename, p, imax, jmax);
+
+    snprintf(filename, 30, "ftest-%d.dat", rank);
+    testWriteFile(filename, f, imax, jmax);
+
+    snprintf(filename, 30, "gtest-%d.dat", rank);
+    testWriteFile(filename, g, imax, jmax);
+}
+
+void commFinalize(Comm* c)
+{
+#ifdef _MPI
+    MPI_Finalize();
+#endif
+}
diff --git a/BasicSolver/2D-mpi/src/comm.h b/BasicSolver/2D-mpi/src/comm.h
index 5c69857..53864c8 100644
--- a/BasicSolver/2D-mpi/src/comm.h
+++ b/BasicSolver/2D-mpi/src/comm.h
@@ -11,7 +11,7 @@
 #endif
 
 enum direction { LEFT = 0, RIGHT, BOTTOM, TOP, NDIRS };
-enum dimension { JDIM = 0, IDIM, NDIMS };
+enum dimension { IDIM = 0, JDIM, NDIMS };
 enum layer { HALO = 0, BULK };
 enum op { MAX = 0, SUM };
 
@@ -20,15 +20,19 @@ typedef struct {
     int size;
 #if defined(_MPI)
     MPI_Comm comm;
-    MPI_Datatype sbufferTypes[NDIRS];
-    MPI_Datatype rbufferTypes[NDIRS];
+    MPI_Datatype bufferTypes[NDIRS];
+    MPI_Aint sdispls[NDIRS];
+    MPI_Aint rdispls[NDIRS];
 #endif
     int neighbours[NDIRS];
     int coords[NDIMS], dims[NDIMS];
     int imaxLocal, jmaxLocal;
 } Comm;
 
+extern int sizeOfRank(int rank, int size, int N);
 extern void commInit(Comm* c, int argc, char** argv);
+extern void commTestInit(Comm* c, double* p, double* f, double* g);
+extern void commTestWrite(Comm* c, double* p, double* f, double* g);
 extern void commFinalize(Comm* c);
 extern void commPartition(Comm* c, int jmax, int imax);
 extern void commPrintConfig(Comm*);
diff --git a/BasicSolver/2D-mpi/src/main.c b/BasicSolver/2D-mpi/src/main.c
index 69a9fce..6e2ffda 100644
--- a/BasicSolver/2D-mpi/src/main.c
+++ b/BasicSolver/2D-mpi/src/main.c
@@ -1,5 +1,5 @@
 /*
- * Copyright (C)  NHR@FAU, University Erlangen-Nuremberg.
+ * Copyright (C) NHR@FAU, University Erlangen-Nuremberg.
  * All rights reserved.
  * Use of this source code is governed by a MIT-style
  * license that can be found in the LICENSE file.
@@ -15,6 +15,26 @@
 #include "solver.h"
 #include "timing.h"
 
+static void writeResults(Solver* s)
+{
+#ifdef _MPI
+    size_t bytesize = (s->imax + 2) * (s->jmax + 2) * sizeof(double);
+
+    double* ug = allocate(64, bytesize);
+    double* vg = allocate(64, bytesize);
+    double* pg = allocate(64, bytesize);
+
+    commCollectResult(&s->comm, ug, vg, pg, s->u, s->v, s->p, s->jmax, s->imax);
+    writeResult(s, ug, vg, pg);
+
+    free(ug);
+    free(vg);
+    free(pg);
+#else
+    writeResult(s, s->u, s->v, s->p);
+#endif
+}
+
 int main(int argc, char** argv)
 {
     int rank;
@@ -35,7 +55,18 @@ int main(int argc, char** argv)
     if (commIsMaster(&s.comm)) {
         printParameter(&p);
     }
+
     initSolver(&s, &p);
+#ifdef TEST
+    commPrintConfig(&s.comm);
+    commTestInit(&s.comm, s.p, s.f, s.g);
+    commExchange(&s.comm, s.p);
+    commShift(&s.comm, s.f, s.g);
+    commTestWrite(&s.comm, s.p, s.f, s.g);
+    writeResults(&s);
+    commFinalize(&s.comm);
+    exit(EXIT_SUCCESS);
+#endif
 #ifndef VERBOSE
     initProgress(s.te);
 #endif
@@ -70,15 +101,8 @@ int main(int argc, char** argv)
     if (commIsMaster(&s.comm)) {
         printf("Solution took %.2fs\n", timeStop - timeStart);
     }
-    size_t bytesize = s.imax * s.jmax * sizeof(double);
-
-    double* ug = allocate(64, bytesize);
-    double* vg = allocate(64, bytesize);
-    double* pg = allocate(64, bytesize);
-
-    commCollectResult(&s.comm, ug, vg, pg, s.u, s.v, s.p, s.jmax, s.imax);
-    writeResult(&s, ug, vg, pg);
 
+    writeResults(&s);
     commFinalize(&s.comm);
     return EXIT_SUCCESS;
 }
diff --git a/BasicSolver/2D-mpi/src/progress.c b/BasicSolver/2D-mpi/src/progress.c
index 517c2f8..be98fce 100644
--- a/BasicSolver/2D-mpi/src/progress.c
+++ b/BasicSolver/2D-mpi/src/progress.c
@@ -7,54 +7,45 @@
 #include <math.h>
 #include <mpi.h>
 #include <stdio.h>
-#include <stdlib.h>
 #include <string.h>
 
 #include "progress.h"
 
 static double _end;
 static int _current;
-static int _rank = -1;
 
 void initProgress(double end)
 {
-    MPI_Comm_rank(MPI_COMM_WORLD, &_rank);
     _end     = end;
     _current = 0;
 
-    if (_rank == 0) {
-        printf("[          ]");
-        fflush(stdout);
-    }
+    printf("[          ]");
+    fflush(stdout);
 }
 
 void printProgress(double current)
 {
-    if (_rank == 0) {
-        int new = (int)rint((current / _end) * 10.0);
+    int new = (int)rint((current / _end) * 10.0);
 
-        if (new > _current) {
-            char progress[11];
-            _current    = new;
-            progress[0] = 0;
+    if (new > _current) {
+        char progress[11];
+        _current    = new;
+        progress[0] = 0;
 
-            for (int i = 0; i < 10; i++) {
-                if (i < _current) {
-                    sprintf(progress + strlen(progress), "#");
-                } else {
-                    sprintf(progress + strlen(progress), " ");
-                }
+        for (int i = 0; i < 10; i++) {
+            if (i < _current) {
+                sprintf(progress + strlen(progress), "#");
+            } else {
+                sprintf(progress + strlen(progress), " ");
             }
-            printf("\r[%s]", progress);
         }
-        fflush(stdout);
+        printf("\r[%s]", progress);
     }
+    fflush(stdout);
 }
 
 void stopProgress()
 {
-    if (_rank == 0) {
-        printf("\n");
-        fflush(stdout);
-    }
+    printf("\n");
+    fflush(stdout);
 }
diff --git a/BasicSolver/2D-mpi/src/solver.c b/BasicSolver/2D-mpi/src/solver.c
index 96cb64d..83d4cd4 100644
--- a/BasicSolver/2D-mpi/src/solver.c
+++ b/BasicSolver/2D-mpi/src/solver.c
@@ -511,9 +511,9 @@ void writeResult(Solver* s, double* u, double* v, double* p)
         exit(EXIT_FAILURE);
     }
 
-    for (int j = 1; j < jmax; j++) {
+    for (int j = 1; j <= jmax; j++) {
         y = (double)(j - 0.5) * dy;
-        for (int i = 1; i < imax; i++) {
+        for (int i = 1; i <= imax; i++) {
             x = (double)(i - 0.5) * dx;
             fprintf(fp, "%.2f %.2f %f\n", x, y, p[j * (imax) + i]);
         }
@@ -529,14 +529,14 @@ void writeResult(Solver* s, double* u, double* v, double* p)
         exit(EXIT_FAILURE);
     }
 
-    for (int j = 1; j < jmax; j++) {
+    for (int j = 1; j <= jmax; j++) {
         y = dy * (j - 0.5);
-        for (int i = 1; i < imax; i++) {
-            x            = dx * (i - 0.5);
-            double vel_u = (u[j * (imax) + i] + u[j * (imax) + (i - 1)]) / 2.0;
-            double vel_v = (v[j * (imax) + i] + v[(j - 1) * (imax) + i]) / 2.0;
-            double len   = sqrt((vel_u * vel_u) + (vel_v * vel_v));
-            fprintf(fp, "%.2f %.2f %f %f %f\n", x, y, vel_u, vel_v, len);
+        for (int i = 1; i <= imax; i++) {
+            x           = dx * (i - 0.5);
+            double velU = (u[j * (imax + 2) + i] + u[j * (imax + 2) + (i - 1)]) / 2.0;
+            double velV = (v[j * (imax + 2) + i] + v[(j - 1) * (imax + 2) + i]) / 2.0;
+            double len  = sqrt((velU * velU) + (velV * velV));
+            fprintf(fp, "%.2f %.2f %f %f %f\n", x, y, velU, velV, len);
         }
     }
 
diff --git a/BasicSolver/2D-mpi/velocity.png b/BasicSolver/2D-mpi/velocity.png
deleted file mode 100644
index a8c9f05c0921deb610bbfa718d22e88f394d794f..0000000000000000000000000000000000000000
GIT binary patch
literal 0
HcmV?d00001

literal 9306
zcma)C3tW<C`+lrU%n~Y_X`-^WwRO@m)6~StTI(~jt<3V6tz}wys5~GDd=6!srT*7S
zO$}S6Tc-A{6ipPvJd~mpDN98oCuAz1Xn+du|3%BKzCMH>@}8dKbzk>&-_P^<<L+JC
z7cDSf003CD<Et;e0l)wOz+8Zs3*RyLy7qzQV)xFy+aL&n@8sp>X%002P@ozB0kOEo
z4Ky}_&qD=JjWwXS)leYYgOn5i+$aE|KxF92wbiTPs&3C7!u8C}&F$>$JUu;+963TJ
zle4n2Zr!?7Q&ZF2+}zjK2X`UexG`OvUi7C>>MeZ+L1KO)RCc#SEQWXh?hPD<f5{)F
ze7HWvXX-k|%q$Wf3IH5h%K*>}5dtlbAd_<tdAF28fwqBh%E|r<@lxTLD#o$$WTUaM
zt~`o0bPb3T#lHY9KaZcj+s*r!4aEddOiU<l@Rve4SVK7oQ8eRsgYw|-K`1X0s;PnC
z`9L)s#w95v%`4;Fh3;}{ZW1B4JjlAbx*DFH30QyzCM$vEr{L3dVB<#MvKf521$b=*
z2i(9R960tRIK2ZT>;%cXK<d{ZV-KM31Nr+wv9DfffF3hQuQpik{<nIq-|4l5=yinZ
z^&Hdt>%^QF5pxEk=15P@k)59NDrSx%Ztetr?j&KZDq*ggG!LT8gZ@1a`d(k1q_0ZW
zpZH<^_>c4DDF$O_3`SB7#A%3^=?LM^hR@F#^3N~$^McXi45Q|(h4mK~)?_cL$Tj}s
zlJRvq^2%jY_OFZ2U0M9o)g_7fOYsFJk-seq`Tc`{>!y2eEdN)L+2_UP4katr-LzQw
zAGFymOXO`U{ZeaC_90-b1mz!rici3uRe<>^s9FuG*8tW!P_rJ~wFT@Apw<rDvj=sX
zLA@hrZ~~1k;Jzz(fCWvTgJvAyYy&OZLF*3ia3^@Q3q0NpxZi*$dqJBIXg>gW{@~9*
z@bn<)_!e{?0bNG{{}|{F2i=jN=M;Dr0|a>RoCtc!;IBl`mjr}A08t7MrGXbegZ}g2
zWfmC70fTv9=n4=Q0P%GoDF)I%z;GEDxdUWXVDv5+tpj5Zz^fJ@e+(4uK-mG5Jz%^S
zOuPV-gFrP5RAWFr4l@fZmDzT|Ha3j>>cA-gjPGgw^v<s_UIBnL*zv{Zd*gD3_(_Gw
zRhGQrDj|J;irMmZ`jLZge<HmXzWvGh(>|(2nJbZXh7>8gct&Eu(^Dd`YGN=!%%o>+
zy(Q-^g*$!wV=+h@db}N2^GKN)u51zK3TlotSQrXl)QUX0pHcxaJ!*L}z8^<y56?+L
zQPOUi39qF|w`2!HAHffEqG^$qF<$P78Ux>GgeWS0i92w8d>vPACMxJe8H0OAf4|#6
zDvzc~`<`XXGd3r*SKHIV`gU{n;a6*pMIG)6C>~#qy2S5V=`UPP7ws;xxc+pt7yQWL
z`gv3tiVc@5eSj;H?KN(DaM{-|*a&mu<Zt^5U7k2i7ixQsvGD5ocW1QKtIL*h4`tDo
zy@lSEVYvlxDdV=VVU+7t@0Li-`SJ79_x^=gvu{i7W#27(Z+TyW+k`8Wq0~Kr`;=Q0
z@fbVGa;_@-T0md>%Wj84*I~6CB6wTlLhJe%CSh2G3a5FgeYTroh}9jdx}na>i5pEl
zWsIKGTE{BZVORf+J56nPq&sTwn?b`J@{*6K%QOvEB8>&prE?15(oe<6eggSX|IhqB
z<oXJ@td-*Gzc4otS$g;W^pxAEM~m8glBMd2h?8>05Qh=G!UkL4R+6}kecn)0WT)4}
z5OKT!XZwn4yK|yIbOzny1jQCK#TIzP78sbzYlnSlLKjwc?Xvoz$arM0fG!jj=5|B3
zWZz$Hfa)eEQKJ6u2kU#GS0;tZzFcN*cxJ0a%}-#;u-*P6@plJG5)-PgnzJPdk`$lC
z?4Q!j+<q2K?5AOo{h6(4TP;_VTE)m*(ByyD4AWD0bXINM6o1EbGt1>u{2jJhLE211
zSgN<jiXE$yBAvd7Wbv8IN*vyyA$F4$wp-|d>7tiNWA0-LWtB$Vq-BF12pTTBgmw1M
z=Y*b#Z;U$K736rs&%S~YfAHDY4Z<ABGqmIqdPOwp-teEf8O|h$d%3f}YADN=JK*h=
z+8U+Yi3lo|6Vs%l8yF6fzxw;dNK|(5dZ)wTkeE}+TJQ=b4>QEBz)j)%#g?rt_{N7a
zG&7VC_C$uR<Mx^MIq}nT*AvdR4<tP34PZxE7)zvf7>fqJdWat|Hf)6H;=4TFB1vW}
zA-(juyNix*X<3yZX<JonJ!Rary&GaCYRN76BGjcn;$H9uR_00E<Cv57ds=-QH-|YN
zm&?>LOin?q!P7vRl<$m>UdrZ2ig>xK!Y);iK4ENLn`Bv+ifK5-vZnCNw5=C=EZrub
zZ)Ml=8##Uvfs*FjlJO*zWl0n>BsL*U+Jx!q8hAGGfDd`!TxWBDS-X{uT)c5^<H<%|
zKR!F61Q(~W6etPdjGmP)t~JV6>3=CkCYw)rNRp+*y92VFNa$Q4DwW|i!dxI2>ity|
z8?9P$wv#}0@ClVBEpRPOmXhU5zDbA8o;jq{Cp3hdb>&p@`5|Y;bh`Z=ay6sWY?tAR
zme$)(Ev%$2_NPc47E3Ew)q7%?2fK}F&Y2CSM~?^kDnDRzlndyj_7@Yzmu#L|wnQuF
zgY@-5Uc>4GS1M)R!V>7}Xngx`sK5nTdQ|-HjY?YPkfb^_va8A^__Rx%q>m8aEBlVy
zpR}leH^*{~{A8&^c5ifQ9YI-f)Cc2Ibxv+WkDE9_5Xn4Sk%X3#jN|zt8Oe#-V7l4Z
z6q`}>HG{W!f#6!Y8KA2^`!!m26B;?+?v7aFegup{kG`AWrY%vWGfAi^VnoJoOYmv%
z@jlbpEg4%mWeeboQ?RC~2&gwdIXoY@?_aucF2LM?y^w-xr#6E;dn$A*BRvs5bweD3
z8@8s}W9$HmeG%>z5$xHx@O8~1%%oKn)#Ri1HS*7zM;8gml8>)izqgU8?kVjkH6ze{
z!IHN}^3O)5yT6?*Wxbj*<69&7#~j>_As`PMjQ-4Y{=1RsIlw^B2~$4#zkryo`OusV
zJ@NSJbL~K+=?cV)H$beA#aJYV%3|(D#zb8kYXy1#hYEapT<#?ck6R<SSy}aXJOn)a
z3L!Bs;E_$-DM#JocsCO!wqxdf)&2pS+|%@}$>Hu)*}#fg+cOeYNZg&l&PxWAw3Sqb
zabIVUzikAqX^HE-g}@gfK+D3**9y3s0-4I?cHpQrH8vng5=s~w`H@iaJ9u`P$T^{u
z2O&z>u}j$zaVEg^7gYC{xNoVvoAiqrAb#zKfC6j=Y#L)6tmqb}#w`W5u3Lv%B3ey-
ze4&NF+)!|yv&p#kG^x0G19hTbv63fjT2n-x56)6=e)v9aXb4uOonWwWBDAw~#yDjD
ztj6)t(KxKZH^zCF>hBP2n|6Z1ETvD2&HrNQ-=213Rv^4|WGVTvb|7XvtyjCJ!S0HE
zsUr}<vjSl?`j)#&Q-ni~3%SA(c0;ya&8oGQY&_vbvzp)FG_YA(*)hd`n5^-s1~sv_
z;}8v@q>e?h&W4V$i|t&;0On^@u;GHS;7BHtoZ87qD}8`6?rd@}b(#<OrAbnW{G?2f
zKdcbN5GVGpStwwlMJYy9ly&`Mj0;v+PN*a<PaYnx_{_*T=O{V){ze16JXRepF1~qK
z<zVR5N<ZLZDB?#RaL3I1;Xn1LYL|m|KJ+m#Z-3Rb4^(&Xx-PK%<&`GVrB-6$<wxg+
zI+mP@&MQSDDSyC;R#jm2<l(GyZAlNCR)KHUq1o9PUOE2JAv-$x#L4_?lErhVzD`uy
z{OamFl4F-mS8k<UY>t9|JU?prs#bITuS(+(J*po{zeUKo!`lfh0S>QL!BhPJDZs&j
zWQHwdhIY10D>|&p&e|CdYZ_}7x0zu}(X4Fg(ba9f)XtU}UFK<LJggt&A{~J!oE->V
z-KI}F5HkVFX6!A-zZ+yV0qT&BK$yPKKxsY6*|>KTE5J4<z&}(4tN2G$pUs^rE)y-}
z871?F&bMRdQoeMLYn?16xg$`XaNH=1fZ(OheRwe<BH+^#2YY%}S|mqx_^ai}jdLkB
zOWB0{U=k!`qSt#UtnGkeBZm2HtSs`y-*JKwBdNwIj=vKrd=L%{G}BPPa54f`o`@2@
z>-u=xU8WJPH48z<?IdaI*dOXdrsSEeWp_6x!{#PdW!m#<Cm>p2m>aiZ54&J<apm;-
z;5cX{FLip<`{6|IpBc!EwqrZBOB&2{JM9LZ&Rt+hql>jlnkEDJT}K0ZY0p6Z_SRpu
zlN1he4ceLhHxTW!0&zu0Agube0x@j|4M`)llN1heNbOAj8;F~;0%7om>C^6_*fufb
zc!o>F3wk-=*-}^snds!uo&+C<>>+hrz4G~#M#mx2<Yh6nUHE+=vYVB0eh9roYtd|u
zt(Dq>ksC6$+23+Ed<&)|p{k2m9_pg5KocY{k4Y>sZr88j9CtNfW7pU;ZE7}+D48q=
zxBXnG1q~KHYEL%54GnkDh@l3Rtte(CZE^AY%MqL&wrBL&@Ko#dC-2b+6xx>#4(QWq
zDC_Nh2r}6QsW_aK-LQ=ghnvR77=KmMDseR2`$>2}xVb!4;&{Ju17E31f}Qi*VN|hd
z`wjlb?@@SGiFT-F6uwHkxnYIx;cDm1jKYhx^B#5{iH`XIIYc{WW)$A}ZF9e83Bsz<
z5r`AB1JR))5ErxqF{AJf+LJrjdB}To1Y*S-0gzr&*|G@794pJ~aqt=}3Duoh5KNK_
z;g!IHmGF%pQGqLvv5+t}!U*HG3BxQkgHH`HZnxp&BRSCWS8E^4l)4A@Vvc}^?udwh
zWGTFI&G?x3gBiGuV*5E@i18QeQyi=g%?InK()L_bcuAmZl*E1+aQUUAVx1{_07bfT
z4oC%44^z2y*~krZDIc7`OpaB#eDQf8&hd<wS0$fF{G*u|b(WT?=VpXqZoXu>!L=$*
zx}2>#U4HPQO^FP^!D*#t6XE?U#^_ni6*P<J&+zu5c5}hJEzvIeGg?RX*3n#R-k9s%
zs6L~0tXV{V#)1202Vxe{p8?{rb|7YWoAGuo{_Y$!yp7Tkh?LnK_#M$-_o|i{eIM5X
zo?Ti-jqEnAy-k)sYH?aLDTa4uOV<0g3V<SYABr7ykseb!`RY+jLVol{%U@lB0zT0s
z{m7Z5zpNNlog>`}C+UwqU<W;H%0)VyZNM^mJGU=&2h~}9)50DBacmd%-tNC$Vt|>4
zC`&|LvNTD6JlR}X$iBU>xFgCAA_@efn|>i+i1?b*XfZYNMC)<=XN3vyYPczy<)MdF
z3c@YTw>~ZYbfb9joO#EL1g(;5efKx6f;Wa8{zF7Wz}|DRnqF=UOZsso%F<jhk~l7i
zV?;8t)4VQZReDaQ)^l-pGHm)$?v&+g5Q5IM{mw~|SoC8%4thGjho@^`@+SAuz#ph1
z{HJE+A+|tA{d_tr59w_>>Swc#JPetghw(ZB!JU<d*is#V_)I$xGnu59c2Ut}l2>&E
zV%3{WQc#3<UdUAs4G9LHJ0Wf7QZ~ECoyZXu@NHSuTL;76sH)HATB71qS9PE@6+EI!
zA8fLOx>8F_5z}FS8rU9`9_6I2AUMC!aKD22dJm~2$LbCD*L#$!&A_H5>@{+C#R#Rp
z@b%791^(#h;dc`KJ&1@)Yfqg-vm;_GfNLfyumZ2VcV9_*K2-_c79taEkh{~PCdOgj
zN3ve?-lq6<2;g<UUeLzH5vrLjbKxL%7ParayNYlw*`Pg_)Hug+9Xavs?3}oxBPYJr
zkrQd!m%tiM{GcNzXtO65`*j2Yr(GpxRu@*YEP>%=)?FQeKx)qvH9*w74Fu__CbcKy
zpRNvZ+`EPrR`y`e7R|BrCzP00(R;O}ergY>?z^;!Ym@C3*V|SjOIuzSNk+ABZRK50
z+NL*zQjDzr`(2qn#)1}{S(*AW)I5zG-tfLM*Vj?8c^MN{Vq(d$-1>fnNIdfFzkB<b
z8Wc&|J^C`6`rc>O)Se%4rM`Wa+h4y0{9vJA6O(u<vJJzWc=##U;Co2=$SRxzUmZ3B
zFG4jZBZgBaE}NvU^vA$M?sZ(8Wcv+0H*9ZzHm&knn%<lDnq9taOLs&dZzTpv;-Y_3
z#aN`Wcu@y@Xcv(*^L$RCffN0?Q1Y0xf6}{<S1|q)y?Qt*+||3F<Qel}Xkcq()NlPt
zJA0d$yek8(Pe|rGxkCb3;lGZQEj^;oxQG;WOI`*vbAsEY1LzXt(x>7y2%CwG9*gMh
zz8q7_W*`~l3}KWK6-W<T=Jm4L^@$*8xzH;>o?SukyLmui@O1z>3_B<+;5-hhb*A~&
zdbZtKhVBc#@sh^rvN*Dc)jcc@IU1;dew4G0i|8L$W;)X@=GQ@C{G71<6_2iTvid#9
z7caRgGr|JM+*YjhN#C*-9_`{TsKKw6Cq23d9%6_$mNeQDd;wC<lVfp4J&*9|xS@Ld
zupcdZPZ5U7qjlU6FqU_zJWBl(`GWD?qJ~5H-$Cak@EPlyPJ7~FuM4a;qGK!>Q}>P7
z;)v3wpWphBB=F-BXcrm1%*NorE{ia`!3>D^SN0yXIsxxLp1l(uLp|SdEwN`n-7{JW
z&+$xltyL3ARct&qsdyMxU3)vxk3%~p_c2-Dry1MN-KvG_RMwNZICOZoDXE1!c6k>A
z*+i4J@FKd`CwRn1{}$C{!IG5@ie7eygbn*8XXA#D<VU8^5ye<WgeT)yuVdx#C=17z
zxb8&6dSIC`^7aIWwSMeE9H)RB5%?X#6#9g_I%0<eeW>4=>OgaHwdsqL;JwM*G=_m}
zIO|FZ`o)707_SE0m7tTN-IaCp1$d~bY_w>i@ABq??O_>60pF{Q-P2m(`b0UT`hy{f
zh4NcdEc!T{?;9H!qNMCo6P}c)<qRXJN0QqYVQ-bXs)=}i=<)>;VxWIdP+3w0AEpRi
zso&n@nmUG{Ul%bPnv$h0@N~#$D-r@9^nYR0v$hLzE}dL8P^TG62;YXjzcyo3cz;!>
z5GM4f=6?Zpr(oBQk2=Cn?gkan6nk?e_;voG&ekU*v4{xcv6CUk1OM98Dpu$Es`FpA
z{ggUlyKk&Az6X885gKS;P<_>h=8H3kFfQeD8bT;*)yK02Y<d>B$H~T=$=8J+r=6bo
z-Mm31kc>ZX)vIuM!koleIlvU@)gN;pb47Kmc*%`t`(%O}W7}`!Mp(Y){T(Zp=sWn#
zRIL0&zfJF0hU6$r`P`R^m3J;cQTuuH>%Ct4{l@Ll4m967mUWO-*a_^Rq0_$wFK2i_
pEJ@qyaQhj5edD)NN6Hb%s&{|M_jPS|=4)>3*tYA7V)uiI{|_I`UIG9B

diff --git a/BasicSolver/2D-seq/Makefile b/BasicSolver/2D-seq/Makefile
index 84c5eff..aab8407 100644
--- a/BasicSolver/2D-seq/Makefile
+++ b/BasicSolver/2D-seq/Makefile
@@ -47,6 +47,8 @@ clean:
 distclean: clean
 	$(info ===>  DIST CLEAN)
 	@rm -f $(TARGET)
+	@rm -f *.dat
+	@rm -f *.png
 
 info:
 	$(info $(CFLAGS))
diff --git a/BasicSolver/2D-seq/dcavity.par b/BasicSolver/2D-seq/dcavity.par
index 4241393..3ba6a1c 100644
--- a/BasicSolver/2D-seq/dcavity.par
+++ b/BasicSolver/2D-seq/dcavity.par
@@ -15,7 +15,7 @@ bcRight    1			#
 gx    0.0			# Body forces (e.g. gravity)
 gy    0.0			#
 
-re    10.0		    # Reynolds number
+re    100.0		    # Reynolds number
 
 u_init    0.0		# initial value for velocity in x-direction
 v_init    0.0		# initial value for velocity in y-direction
diff --git a/BasicSolver/3D-mpi/src/progress.c b/BasicSolver/3D-mpi/src/progress.c
index d5393c4..f03dc02 100644
--- a/BasicSolver/3D-mpi/src/progress.c
+++ b/BasicSolver/3D-mpi/src/progress.c
@@ -6,7 +6,6 @@
  */
 #include <math.h>
 #include <stdio.h>
-#include <stdlib.h>
 #include <string.h>
 
 #include "progress.h"