course/semester2/pprog/assignment1/single-thread/nbody.c

   1 /**
   2  * @file nbody.c
   3  * @author Sam Moore (20503628) 2012
   4  * @purpose N-Body simulator - Definition of simulation functions; single threaded version
   5  */
   6
   7 #include <stdlib.h>
   8 #include <stdio.h>
   9 #include <math.h>
  10 #include <time.h>
  11 #include <string.h>
  12
  13 #include "nbody.h"
  14
  15 /**
  16  * Prints the body on screen
  17  */
  18 void Body_Print(Body * a)
  19 {
  20         printf("Body %p M=%f X=%f Y=%f Z=%f Fx=%f Fy=%f Fz=%f Vx=%f Vy=%f Vz=%f\n",
  21            (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]);
  22 }
  23
  24 /**
  25  * Computing forces
  26  */
  27 void Body_Force(Body * a, System * s)
  28 {
  29         double distance;
  30         double con;
  31         double gd;
  32
  33         for (unsigned i = 0; i < DIMENSIONS; ++i)
  34                 a->F[i] = 0;
  35
  36         for (unsigned index = 0; index < s->N; ++index)
  37         {
  38                 Body * b = s->body+index;
  39                 if (b == a)
  40                         continue;
  41
  42                 distance = 0.0;
  43                 for (unsigned i = 0; i < DIMENSIONS; ++i)
  44                         distance += square(b->x[i] - a->x[i]);
  45                 distance = sqrt(distance);
  46                 con = G * a->mass * b->mass / square(distance);
  47                 gd = con / distance;
  48                 for (unsigned i = 0; i < DIMENSIONS; ++i)
  49                         a->F[i] += gd * (b->x[i] - a->x[i]);
  50         }
  51 }
  52
  53 /**
  54  * Compute velocities
  55  */
  56 void Body_Velocity(Body * a)
  57 {
  58         for (unsigned i = 0; i < DIMENSIONS; ++i)
  59                 a->v[i] += a->F[i] / a->mass * DELTA_T;
  60 }
  61
  62 /**
  63  * Compute positions
  64  */
  65 void Body_Position(Body * a)
  66 {
  67         for (unsigned i = 0; i < DIMENSIONS; ++i)
  68                 a->x[i] += a->v[i] * DELTA_T;
  69 }
  70
  71 /**
  72  * Main compute function
  73  */
  74 void System_Compute(System * s)
  75 {
  76 //      clock_t start, finish;
  77
  78 //      start = clock();
  79
  80         for (unsigned i = 0; i < s->N; ++i)
  81         {
  82                 Body_Force(s->body+i, s);
  83                 Body_Velocity(s->body+i);
  84                 Body_Position(s->body+i);
  85         }
  86
  87 }
  88
  89
  90
  91
  92 /*
  93  * This function reads an input file. You can change it if you choose a
  94  * different file format
  95  */
  96 #define LINE_SIZE BUFSIZ
  97 void System_Init(System * s, char *fileName)
  98 {
  99         char line[LINE_SIZE];
 100         char * token;
 101         FILE * file;
 102
 103         file = fopen(fileName, "rt");
 104         s->N = atoi(fgets(line, LINE_SIZE, file));
 105         s->body = (Body*) calloc((size_t)s->N, sizeof(Body));
 106
 107         //printf("----------------------Initial field-------------------------------\n");
 108
 109         for (unsigned i = 0; i < s->N; ++i)
 110         {
 111                 if (fgets(line, LINE_SIZE, file) != NULL)
 112                 {
 113                         Body * a = s->body+i;
 114                         token = strtok(line, ",; ");
 115                         a->mass = atof(token);
 116
 117                         for (unsigned j = 0; j < DIMENSIONS; ++j)
 118                         {
 119                                 token = strtok(NULL, ",; ");
 120                                 a->x[j] = atof(token);
 121                         }
 122                         for (unsigned j = 0; j < DIMENSIONS; ++j)
 123                         {
 124                                 token = strtok(NULL, ",; ");
 125                                 a->v[j] = atof(token);
 126                         }
 127                         //Body_Print(a);
 128                 }
 129         }
 130
 131         //printf("--------------------------------------------------------------\n");
 132
 133         fclose(file);
 134 }
 135
 136 /**
 137  * Cleans up the universe by freeing all memory
 138  */
 139 void Universe_Cleanup()
 140 {
 141         free(universe.body);
 142 }
 143