/**
* @file nbody.c
- * @author Sam Moore (20503628) 2012
- * @purpose N-Body simulator - Definition of simulation functions; single threaded version
+ * @purpose Implentation of multi-threaded N-Body simulator using pthreads
+ * @author Sam Moore (20503628) - 2012
*/
-#include <stdlib.h>
-#include <stdio.h>
-#include <math.h>
-#include <time.h>
-#include <string.h>
-#include <stdbool.h>
+#include "nbody.h" // Declarations
+#include "../single-thread/nbody.c" // Include all functions from the single threaded version
-#include "nbody.h"
+#include "graphics.h" // For declaration of Graphics_Run only
+
+// --- Variable declarations --- //
+
+pthread_t compute_thread; // The thread responsible for computations; it spawns worker threads
+
+pthread_t * worker_thread = NULL; //Array of worker threads responsible for Force and Position updates
+System * sub_system = NULL; //Array of Systems used to divide up the main "universe" System for worker threads
+System * nested_sub_system = NULL; //Array of Systems for division of "universe" for the nested worker threads
+
+pthread_mutex_t mutex_runstate; // Mutex around the runstate
+
+
+
+pthread_attr_t attr; //thread attribute for the workers.
+
+Barrier force_barrier; // I laughed at this variable name. A bit sad really.
+Barrier position_barrier;
+
+
+Barrier graphics_barrier;
-/**
- * Prints the body on screen
- */
-void Body_Print(Body * a)
-{
- printf("Body %p M=%f X=%f Y=%f Z=%f Fx=%f Fy=%f Fz=%f Vx=%f Vy=%f Vz=%f\n",
- (void*)a, a->mass, a->x[0], a->x[1], a->x[2], a->F[0], a->F[1], a->F[2], a->v[0], a->v[1], a->v[2]);
-}
/**
- * Computing forces
+ * @function Compute_Thread
+ * @purpose Thread - Continuously computes steps for a system of bodies. Seperate from graphics functions.
+ * Spawns worker threads to divide up computation.
+ * @param arg - Can be cast to the System for which computations are performed (ie: &universe)
*/
-void Body_Force(Body * a, System * s)
+void * Compute_Thread(void * arg)
{
- double distance;
- double con;
- double gd;
- for (unsigned i = 0; i < DIMENSIONS; ++i)
- a->F[i] = 0;
+ System * s = (System*)(arg); //cast argument to a System*
+
+
+ // If no number of threads provided, use the default value, unless someone changed that to a stupid value
+ if (options.num_threads <= 0)
+ options.num_threads = (DEFAULT_WORKING_THREADS > 1) ? DEFAULT_WORKING_THREADS : 1;
+
+ if (options.nested_threads <= 0)
+ options.nested_threads = 1;
+
+ // Do a sanity check; there is no point spawning more threads than bodies.
+ if (options.num_threads > s->N)
+ {
+ fprintf(stderr,
+ "Warning: Using %u threads instead of %u specified, because there are only %u bodies to simulate!\n",
+ s->N, options.num_threads, s->N);
+ options.num_threads = s->N;
+ }
+
+ pthread_attr_init(&attr);
+ pthread_attr_setdetachstate(&attr, PTHREAD_CREATE_DETACHED); //Needs to be detached, so that memory can be reused.
+
+ if (options.num_threads > 1) // Allocate worker threads and sub systems, as long as there would be more than 1
+ {
+ worker_thread = Allocate_Threads(options.num_threads);
+ sub_system = Split_System(&universe, options.num_threads);
+
+ if (options.nested_threads > 1)
+ nested_sub_system = Split_System(&universe, options.nested_threads);
+
+ #ifdef PERSISTENT_THREADS
+ for (unsigned i = 0; i < options.num_threads; ++i)
+ {
+ if (pthread_create(worker_thread+i, & attr, Worker_Thread, (void*)(sub_system+i)) != 0)
+ {
+ perror("In compute thread, couldn't create worker thread");
+ QuitProgram(true);
+ pthread_exit(NULL);
+ }
+ }
+ #endif //PERSISTENT_THREADS
+
+
+ }
+ #ifdef PERSISTENT_THREADS
+ else
+ {
+ while (!ExitCondition())
+ {
+ if (options.verbosity != 0 && universe.steps % options.verbosity == 1)
+ DisplayStatistics();
+
+ // Just do everything in this thread
+ System_Forces(s, s);
+ System_Positions(s);
+ }
+ QuitProgram(false);
+ pthread_exit(NULL);
+ }
- for (unsigned index = 0; index < s->N; ++index)
+ #else
+ // The main computation loop
+ while (!ExitCondition())
{
- Body * b = s->body+index;
- if (b == a)
+ if (options.verbosity != 0 && universe.steps % options.verbosity == 1)
+ DisplayStatistics();
+
+ if (options.num_threads <= 1)
+ {
+ // Just do everything in this thread
+ System_Forces(s, s);
+ System_Positions(s);
continue;
+ }
+
- distance = 0.0;
- for (unsigned i = 0; i < DIMENSIONS; ++i)
- distance += square(b->x[i] - a->x[i]);
- distance = sqrt(distance);
- con = G * a->mass * b->mass / square(distance);
- gd = con / distance;
- for (unsigned i = 0; i < DIMENSIONS; ++i)
- a->F[i] += gd * (b->x[i] - a->x[i]);
+
+ //Compute forces
+ for (unsigned i = 0; i < options.num_threads; ++i)
+ {
+ if (pthread_create(worker_thread+i, &attr, Force_Thread, (void*)(sub_system+i)) != 0)
+ {
+ perror("In compute thread, couldn't create worker thread (force)");
+ QuitProgram(true);
+ pthread_exit(NULL);
+ }
+
+ }
+
+
+ Barrier_Wait(&force_barrier);
+
+ //All the forces are now computed
+
+ if (options.draw_graphics && options.pedantic_graphics)
+ {
+ Barrier_Wait(&graphics_barrier);
+ }
+
+ //Compute positions
+ for (unsigned i = 0; i < options.num_threads; ++i)
+ {
+ if (pthread_create(worker_thread+i, &attr, Position_Thread, (void*)(sub_system+i)) != 0)
+ {
+ perror("In compute thread, couldn't create worker thread (position)");
+ QuitProgram(true);
+ pthread_exit(NULL);
+ }
+ }
+
+ //Wait for positions to be computed
+ Barrier_Wait(&position_barrier);
+
+
+ //Update number of steps computed
+ universe.steps += 1;
+
+
+
+
+
+
}
+ QuitProgram(false);
+ #endif //PERSISTENT_THREADS
+ return NULL;
}
/**
- * Compute velocities
+ * @function BeforeDraw
+ * @purpose Called in graphics thread before the draw loop
+ * When --pedantic-graphics enabled, will wait for position computations to finish before drawing
+ * Otherwise does nothing
*/
-void Body_Velocity(Body * a)
+void BeforeDraw()
{
- for (unsigned i = 0; i < DIMENSIONS; ++i)
- a->v[i] += a->F[i] / a->mass * DELTA_T;
+ if (options.verbosity != 0 && universe.steps % options.verbosity == 0)
+ DisplayStatistics();
+ //printf("BEFORE DRAW\n");
+ if (!options.pedantic_graphics)
+ return;
+
+ Barrier_Wait(&position_barrier);
+
+
+ Barrier_Enter(&graphics_barrier);
}
/**
- * Compute positions
+ * @function AfterDraw
+ * @purpose Called in graphics thread after the draw loop
+ * When --pedantic-graphics is supplied, will signal computation thread that drawing is finished
+ * So that positions can be safely altered
+ * Otherwise does nothing
*/
-void Body_Position(Body * a)
+void AfterDraw()
{
- for (unsigned i = 0; i < DIMENSIONS; ++i)
- a->x[i] += a->v[i] * DELTA_T;
+ universe.steps += 1;
+ if (!options.pedantic_graphics)
+ return;
+ Barrier_Leave(&graphics_barrier);
+
}
+#ifdef PERSISTENT_THREADS
+
/**
- * Main compute function
+ * @function Worker_Thread
+ * @purpose Thread - Calculate stuff
*/
-void System_Compute(System * s)
+void * Worker_Thread(void * arg)
{
-// clock_t start, finish;
+ System * s = (System*)(arg);
-// start = clock();
- for (unsigned i = 0; i < s->N; ++i)
+ pthread_t * nested_workers = NULL;
+ System_ForcePair * system_pairs = NULL;
+ System * nested_position = NULL;
+
+ Barrier nested_barrier;
+ Barrier_Init(&nested_barrier);
+
+ printf("options.nested_threads == %d\n", (int)(options.nested_threads));
+
+ if (options.nested_threads != 1)
{
- Body_Force(s->body+i, s);
- Body_Velocity(s->body+i);
- Body_Position(s->body+i);
+
+ system_pairs = (System_ForcePair*)(calloc(options.nested_threads, sizeof(System_ForcePair)));
+ if (system_pairs == NULL) // Handle tedious error cases
+ {
+ perror("Couldn't allocate array of system pairs");
+ QuitProgram(true);
+ pthread_exit(NULL);
+ }
+ nested_workers = Allocate_Threads(options.nested_threads);
+ nested_position =
+
+ for (unsigned i = 0; i < options.nested_threads; ++i)
+ {
+ system_pairs[i].A = s;
+ system_pairs[i].B = nested_sub_system+i;
+ }
}
+ while (!ExitCondition())
+ {
+
+ if (options.nested_threads == 1)
+ {
+ Barrier_Enter(&force_barrier);
+ System_Forces(s, &universe);
+ Barrier_Leave(&force_barrier);
+ }
+ else
+ {
+ for (unsigned i = 0; i < options.nested_threads; ++i)
+ {
+ if (pthread_create(nested_workers+i, &attr, Force_Thread, (void*)(system_pairs+i)) != 0)
+ {
+ perror("In worker thread, couldn't create nested worker thread (force)");
+ QuitProgram(true);
+ free(nested_workers);
+ pthread_exit(NULL);
+ }
+ }
+ }
+ //printf("Computed forces for %p\n", arg);
+ Barrier_Wait(&force_barrier);
+ //printf("Computed ALL forces\n");
+
+ if (options.draw_graphics && options.pedantic_graphics)
+ Barrier_Wait(&graphics_barrier);
+
+ Barrier_Enter(&position_barrier);
+
+ System_Positions(s);
+
+ Barrier_Leave(&position_barrier);
+ //printf("Computed positions for %p\n", arg);
+ Barrier_Wait(&position_barrier);
+ //printf("Computed ALL positions\n");
+ }
+ printf("Worker thread exits\n");
+ QuitProgram(false);
+ pthread_exit(NULL);
}
+#endif //PERSISTENT_THREADS
+/**
+ * @function Force_Thread
+ * @purpose Thread - Calculates the forces on objects in a System
+ * @param s - Can be cast to the System*
+ */
+void * Force_Thread(void * s)
+{
+ System_ForcePair * pair = (System_ForcePair*)s;
+ Barrier_Enter(&force_barrier);
+
+ System_Forces(pair->A, pair->B); //Simple wrapper
+
+ Barrier_Leave(&force_barrier);
+ return NULL;
+}
+
+
+/**
+ * @function Position_Thread
+ * @purpose Thread - Calculates the positions of objects in a System
+ * @param s - Can be cast to the System*
+ */
+void * Position_Thread(void * s)
+{
+ Barrier_Enter(&position_barrier);
+ System_Positions((System*)s); // Simple wrapper
+ Barrier_Leave(&position_barrier);
+ return NULL;
+}
-/*
- * This function reads an input file. You can change it if you choose a
- * different file format
+/**
+ * @function QuitProgram
+ * @purpose This function can either be called by the main thread in order to signal other threads
+ * that it wants to exit. The main thread then calls pthread_join before exiting.
+ * It can also be called by a child thread to request the main thread to exit.
+ * It is only used this way if there is an unrecovarable error (ie: Can't allocate memory in a child thread)
*/
-#define LINE_SIZE BUFSIZ
-void System_Init(System * s, char *fileName)
+void QuitProgram(bool error)
{
- char line[LINE_SIZE];
- char * token;
- FILE * file;
+ if (runstate == QUIT || runstate == QUIT_ERROR)
+ return; //Don't do anything if already quitting
+ pthread_mutex_lock(&mutex_runstate); // aquire mutex
+ if (error) // set the runstate
+ runstate = QUIT_ERROR;
+ else
+ runstate = QUIT;
+ pthread_mutex_unlock(&mutex_runstate); //release mutex
+}
- file = fopen(fileName, "rt");
- s->N = atoi(fgets(line, LINE_SIZE, file));
- s->body = (Body*) calloc((size_t)s->N, sizeof(Body));
+/**
+ * @function Thread_Cleanup
+ * @purpose Will be called in the main thread when exit() is called
+ * Automatically tells all other threads to quit (if they haven't already been told)
+ * Then waits for them to finish.
+ * Also frees memory associated with the worker threads.
+ */
+void Thread_Cleanup(void)
+{
+ if (runstate == RUN) // If this is true, as far as child threads are concerned, the simulation is still running
+ QuitProgram(false); // So call QuitProgram which will set runstate, and cause child threads to exit
+ pthread_join(compute_thread, NULL);
+ free(worker_thread);
+ free(sub_system);
+}
- //printf("----------------------Initial field-------------------------------\n");
- for (unsigned i = 0; i < s->N; ++i)
+/**
+ * @function Simulation_Run
+ * @purpose Initialise and start the simulation. Will be called in the main thread.
+ * Replaces the single-threaded macro that does nothing, and sets up the compute thread
+ * @param argc - Number of arguments - Passed to Graphics_Run if needed
+ * @param argv - Argument strings - Passed to Graphics_Run if needed
+ */
+void Simulation_Run(int argc, char ** argv)
+{
+ atexit(Thread_Cleanup);
+
+ Barrier_Init(&force_barrier);
+ Barrier_Init(&position_barrier);
+ Barrier_Init(&graphics_barrier);
+
+
+ if (options.draw_graphics)
{
- if (fgets(line, LINE_SIZE, file) != NULL)
+ // The graphics are enabled, so create a thread to do computations
+ // Graphics are done in the main loop
+ //printf("Graphics are enabled\n");
+ #ifdef PERSISTENT_THREADS
+ Compute_Thread((void*)(&universe));
+ #else
+ if (pthread_create(&compute_thread, NULL, Compute_Thread, (void*)&universe) != 0)
{
- Body * a = s->body+i;
- token = strtok(line, ",; ");
- a->mass = atof(token);
-
- for (unsigned j = 0; j < DIMENSIONS; ++j)
- {
- token = strtok(NULL, ",; ");
- a->x[j] = atof(token);
- }
- for (unsigned j = 0; j < DIMENSIONS; ++j)
- {
- token = strtok(NULL, ",; ");
- a->v[j] = atof(token);
- }
- //Body_Print(a);
+ perror("Error creating compute thread");
+ exit(EXIT_FAILURE);
}
+ #endif //PERSISTENT_THREADS
+ //printf("Run compute thread\n");
+ Graphics_Run(argc, argv);
}
-
- //printf("--------------------------------------------------------------\n");
+ else
- fclose(file);
+ Compute_Thread((void*)(&universe)); // Graphics are disabled, so do computations in the main thread
}
-/**
- * Cleans up the universe by freeing all memory
- */
-void Universe_Cleanup()
+
+
+void Barrier_Init(Barrier * b)
{
- free(universe.body);
- pthread_exit(NULL);
+ b->threads_busy = 0;
+}
+
+void Barrier_Enter(Barrier * b)
+{
+ pthread_mutex_lock(&(b->mutex));
+ b->threads_busy += 1;
+ pthread_mutex_unlock(&(b->mutex));
+}
+
+void Barrier_Leave(Barrier * b)
+{
+ pthread_mutex_lock(&(b->mutex));
+ b->threads_busy -= 1;
+ if (b->threads_busy == 0)
+ {
+ pthread_cond_signal(&(b->threads_done_cv));
+ }
+ pthread_mutex_unlock(&(b->mutex));
+}
+
+void Barrier_Wait(Barrier * b)
+{
+ pthread_mutex_lock(&(b->mutex));
+ while (b->threads_busy > 0)
+ pthread_cond_wait(&(b->threads_done_cv), &(b->mutex));
+ pthread_mutex_unlock(&(b->mutex));
}
/**
- * Thread - Continuously computes steps for a system of bodies
+ * @function Split_System
+ * @purpose Helper to divide one system into an array of systems
+ * Each sub system will have N = (s->N / n) bodies in it
+ * @param s - The original system (typically &universe)
+ * @param n - The number of sub systems in the array
+ *
+ * WARNING: It is the caller's responsibility to free() the returned array
*/
-void * Compute_Thread(void * s)
+System * Split_System(System * s, unsigned n)
{
- while (true)
+ System * result = (System*)(calloc(n, sizeof(System)));
+ if (result == NULL)
{
- if (terminate) pthread_exit(NULL);
- System_Compute((System*)s);
+ perror("Couldn't create array of sub systems");
+ QuitProgram(true);
+ pthread_exit(NULL);
}
-}
+ unsigned n_per_system = (s->N) / n;
+ unsigned remainder = (s->N) % n;
+ for (unsigned i = 0; i < n; ++i)
+ {
+ result[i].N = n_per_system;
+ if (i == n-1)
+ result[i].N += remainder;
+ result[i].body = (s->body) + (n_per_system * i);
+ result[i].steps = 0;
+ }
+ return result;
+}
+/**
+ * @function Allocate_Threads
+ * @purpose Helper function to allocate an array of pthread_t objects
+ * Handles all the pointless, er, "important" error checking that should be done
+ * @param n - Number of threads in the array
+ *
+ * WARNING: Remember to free() the array!!!
+ */
+pthread_t * Allocate_Threads(unsigned n)
+{
+ pthread_t * result = (pthread_t*)(calloc(n, sizeof(pthread_t)));
+ if (result == NULL)
+ {
+ perror("Unable to allocate memory for threads");
+ QuitProgram(true);
+ pthread_exit(NULL);
+ }
+ return result;
+}
git.ucc.asn.au / matches / honours.git / blobdiff
UCC git Repository :: git.ucc.asn.au