course/semester2/pprog/assignment1/openmp/nbody.c~

   1 /**
   2  * @file nbody.c
   3  * @purpose OpenMP version of N-Body simulator, implementation
   4  * @author Sam Moore (20503628) - 2012
   5  */
   6
   7 #include "nbody.h"
   8 #include "../single-thread/nbody.c" // Include original code
   9 #include "graphics.h" // For Graphics_Run function only
  10
  11
  12 #ifdef OVER_ENGINEERED
  13 System * sub_system = NULL;
  14 #endif //OVER_ENGINEERED
  15
  16 /**
  17  * @function Simulation_Run
  18  * @purpose Start the simulation
  19  * @param argc, argv - Passed to Graphics_Run if needed, otherwise unused
  20  */
  21 void Simulation_Run(int argc, char ** argv)
  22 {
  23
  24         if (options.num_threads == 0)
  25                 options.num_threads = omp_get_max_threads();
  26
  27         if (options.draw_graphics)
  28         {
  29                 Graphics_Run(argc, argv);
  30                 if (options.pedantic_graphics == false)
  31                 {
  32                         options.num_threads += 1;
  33                         omp_set_nested(1);
  34                         if (omp_get_nested() == 0)
  35                         {
  36                                 fprintf(stderr, "Couldn't enable nested parallelism\n.");
  37                                 exit(EXIT_FAILURE);
  38                         }
  39
  40
  41                         #pragma omp parallel sections num_threads(2)
  42                         {
  43                                 #pragma omp section
  44                                 {
  45                                         Compute();
  46                                 }
  47                                 #pragma omp section
  48                                 {
  49                                         glutMainLoop();
  50                                 }
  51                         }
  52                         return;
  53                 }
  54         }
  55         #ifdef OVER_ENGINEERED
  56         Compute();
  57         #else
  58         Compute2();
  59         #endif //OVER_ENGINEERED
  60
  61
  62
  63 }
  64
  65
  66 void Compute(void)
  67 {
  68         omp_set_num_threads(options.num_threads);
  69         if (omp_get_max_threads() != options.num_threads)
  70         {
  71                 fprintf(stderr, "Couldn't use %d threads!\n", options.num_threads);
  72                 exit(EXIT_FAILURE);
  73         }
  74
  75         #ifdef OVER_ENGINEERED
  76         sub_system = Split_System(&universe, options.num_threads); // I could use omp for loops instead, but I like this way
  77         #endif //OVER_ENGINEERED
  78
  79         #pragma omp parallel
  80         {
  81                 #ifdef OVER_ENGINEERED
  82                 int index = omp_get_thread_num();
  83                 #endif //OVER_ENGINEERED
  84                 while (!ExitCondition())
  85                 {
  86
  87                         #ifdef OVER_ENGINEERED
  88                                 System_Forces(sub_system+index, &universe);
  89                                 #pragma omp barrier // Look! An explicit barrier!
  90                                 System_Positions(sub_system+index);
  91                                 #pragma omp barrier
  92                         #else
  93                                 #pragma omp for
  94                                 for (unsigned a = 0; a < universe.N; ++a)
  95                                 {
  96                                         for (unsigned i = 0; i < DIMENSIONS; ++i)
  97                                                 (universe.body+a)->F[i] = 0;
  98                                         //printf("Thread %d computes force for body %d\n", omp_get_thread_num(), a);
  99                                         omp_set_num_threads(options.nested_threads);
 100                                         #pragma omp parallel for
 101                                         for (unsigned b = 0; b < universe.N; ++b)
 102                                         {
 103                                                 if (b != a)
 104                                                         Body_Force(universe.body+a, universe.body+b);
 105                                         }
 106                                 }
 107
 108                                 #pragma omp for
 109                                 for (unsigned a = 0; a < universe.N; ++a)
 110                                 {
 111                                         //printf("Thread %d computes position for body %d\n", omp_get_thread_num(), a);
 112                                         Body_Velocity(universe.body+a);
 113                                         Body_Position(universe.body+a);
 114                                 }
 115
 116                         #endif //OVER_ENGINEERED
 117
 118                         #pragma omp master
 119                         {
 120                                 // The if statement has a cost, but if can be removed easily in an optimised version.
 121                                 // Only the main thread is allowed to mess with glut, because glut is evil (not thread safe).
 122                                 // Although this function may be nested, that will only happen if the if would be false.
 123                                 // ("if the if would be false"; I am awesome at explaining code)
 124                                 if (options.draw_graphics && options.pedantic_graphics)
 125                                 {
 126                                         glutMainLoopEvent();
 127                                         Graphics_Display();
 128                                 }
 129                         }
 130                         #pragma omp single
 131                         {
 132                                 // Only one thread should execute this, and we don't care which one.
 133                                 // So much easier than the pthread version!
 134                                 universe.steps += 1;
 135                                 if (options.verbosity != 0 && universe.steps % options.verbosity == 0)
 136                                         DisplayStatistics();
 137                         }
 138
 139                 }
 140         }
 141 }
 142