OMP version

src/impl/omp.c

@@ -0,0 +1,47 @@
+/*
+ * OpenMP version.
+ */
+
+#include <stdio.h>
+#include <math.h>
+#include <omp.h>
+#include "../config.h"
+#include "../utils.h"
+
+double *compute_jacobi(int n, double init_value, double threshold, borders b, int *iterations) {
+  double *x;
+  double max_diff, new_x;
+  int i, j;
+  int nb = n + 2; // n plus the border
+
+  /* Initialize boundary regions */
+  x = create_sa_matrix(n + 2, n + 2);
+  for (i = 1; i <= n; i++) {
+    x[IDX(nb, 0, i)] = b.north;
+    x[IDX(nb, n + 1, i)] = b.south;
+    x[IDX(nb, i, 0)] = b.west;
+    x[IDX(nb, i, n + 1)] = b.east;
+  }
+  /* Initialize the rest of the matrix */
+  for (i = 1; i <= n; i++) {
+    for (j = 1; j <= n; j++) {
+      x[IDX(nb, i, j)] = init_value;
+    }
+  }
+  /* Iterative refinement of x until values converge */
+  omp_set_num_threads(4);
+  *iterations = 0;
+  do {
+    max_diff = 0;
+    #pragma omp parallel for schedule(static, 20) reduction (max:max_diff) private(new_x, j) firstprivate(n, nb) shared(x)
+    for (i = 1; i <= n; i++) {
+      for (j = 1; j <= n; j++) {
+        new_x = 0.25 * (x[IDX(nb, i - 1, j)] + x[IDX(nb, i, j + 1)] + x[IDX(nb, i + 1, j)] + x[IDX(nb, i, j - 1)]);
+        max_diff = (double) fmax(max_diff, fabs(new_x - x[IDX(nb, i, j)]));
+        x[IDX(nb, i, j)] = new_x;
+      }
+    }
+    (*iterations)++;
+  } while (max_diff > threshold);
+  return x;
+}

src/main.c

@@ -35,7 +35,9 @@ int main(int argc, char* argv[]) {
 
   printf("clock time: %fs\n", (double)(endtime - starttime) / CLOCKS_PER_SEC);
   printf("Iterations: %d\n", iterations);
-  print_sa_matrix(x, n + 2, n + 2);
+  if (n < 10) {
+    print_sa_matrix(x, n + 2, n + 2);
+  }
   destroy_sa_matrix(x);
 
   return 0;

src/main_mpi.c

@@ -61,7 +61,9 @@ int main(int argc, char* argv[]) {
     endwtime = MPI_Wtime();
     printf("Wall clock time: %fs\n", endwtime - startwtime);
     printf("Iterations: %d\n", iterations);
-    /* print_sa_matrix(x, n + 2, n + 2); */
+    if (n < 10) {
+      print_sa_matrix(x, n + 2, n + 2);
+    }
   }
 
   destroy_sa_matrix(x);