ex5: finish pthreads

exercise5/handout_pthreads/wave_2d_pthread.c

@@ -10,7 +10,7 @@
 // TASK: T1a
 // Include the pthreads library
 // BEGIN: T1a
-;
+#include <pthread.h>
 // END: T1a
 
 // Option to change numerical precision
@@ -20,7 +20,8 @@ typedef double real_t;
 // TASK: T1b
 // Pthread management
 // BEGIN: T1b
-int_t n_threads = 1;
+int_t n_threads = 4;
+pthread_barrier_t barrier;
 // END: T1b
 
 // Performance measurement
@@ -86,7 +87,8 @@ void time_step(int_t thread_id) {
   // BEGIN: T3
   for (int_t i = 0; i < N; i += 1)
     for (int_t j = 0; j < N; j++)
-      U_nxt(i, j) = -U_prv(i, j) + 2.0 * U(i, j) + (dt * dt * c * c) / (h * h) * (U(i - 1, j) + U(i + 1, j) + U(i, j - 1) + U(i, j + 1) - 4.0 * U(i, j));
+      if (j % n_threads == thread_id)
+        U_nxt(i, j) = -U_prv(i, j) + 2.0 * U(i, j) + (dt * dt * c * c) / (h * h) * (U(i - 1, j) + U(i + 1, j) + U(i, j - 1) + U(i, j + 1) - 4.0 * U(i, j));
   // END: T3
 }
 
@@ -95,10 +97,14 @@ void time_step(int_t thread_id) {
 void boundary_condition(int_t thread_id) {
   // BEGIN: T4
   for (int_t i = 0; i < N; i += 1) {
+    if (i % n_threads != thread_id)
+      continue;
     U(i, -1) = U(i, 1);
     U(i, N) = U(i, N - 2);
   }
   for (int_t j = 0; j < N; j += 1) {
+    if (j % n_threads != thread_id)
+      continue;
     U(-1, j) = U(1, j);
     U(N, j) = U(N - 2, j);
   }
@@ -121,16 +127,21 @@ void *simulate(void *id) {
   // BEGIN: T5
   // Go through each time step
   for (int_t iteration = 0; iteration <= max_iteration; iteration++) {
-    if ((iteration % snapshot_freq) == 0) {
-      domain_save(iteration / snapshot_freq);
-    }
+    int_t tid = *((int_t *)id);
+    if (tid == 0)
+      if ((iteration % snapshot_freq) == 0)
+        domain_save(iteration / snapshot_freq);
 
     // Derive step t+1 from steps t and t-1
-    boundary_condition(0);
-    time_step(0);
+    boundary_condition(tid);
+    pthread_barrier_wait(&barrier);
+    time_step(tid);
+    pthread_barrier_wait(&barrier);
 
     // Rotate the time step buffers
-    move_buffer_window();
+    if (tid == 0)
+      move_buffer_window();
+    pthread_barrier_wait(&barrier);
   }
   // END: T5
 }
@@ -151,7 +162,7 @@ int main(int argc, char **argv) {
   // TASK: T1c
   // Initialise pthreads
   // BEGIN: T1c
-  ;
+  pthread_barrier_init(&barrier, NULL, n_threads);
   // END: T1b
 
   // Set up the initial state of the domain
@@ -163,7 +174,16 @@ int main(int argc, char **argv) {
   // TASK: T2
   // Run the integration loop
   // BEGIN: T2
-  simulate(NULL);
+  pthread_t threads[n_threads];
+  int_t tids[n_threads];
+
+  for (int i = 0; i < n_threads; i++) {
+    tids[i] = i;
+    pthread_create(&threads[tids[i]], NULL, &simulate, &tids[i]);
+  }
+
+  for (int i = 0; i < n_threads; i++)
+    pthread_join(threads[tids[i]], NULL);
   // END: T2
 
   // Report how long we spent in the integration stage
@@ -178,7 +198,7 @@ int main(int argc, char **argv) {
   // TASK: T1d
   // Finalise pthreads
   // BEGIN: T1d
-  ;
+  pthread_barrier_destroy(&barrier);
   // END: T1d
 
   exit(EXIT_SUCCESS);