+/**
+ * @file nbody.c
+ * @purpose Implementations for N-Body simulator
+ * @author Sam Moore
+ * @date August 2012
+ */
+
+#include <assert.h> // To make sure I don't do anything stupid.
+#include <math.h> //For maths and stuff
+#include <stdio.h> //For writing data to files
+#include "nbody.h" //Declarations of functions and structures in this file
+
+// --- Body functions --- //
+
+/**
+ * @function Body_Init
+ * @purpose Initialise a Body object with mass, position and velocity
+ * @param body - Target Body
+ * @param m - The mass of the Body
+ * @param x[...] - The position vector of the body
+ * @param v[...] - The velocity vector of the body
+ */
+void Body_Init(Body * body, float m, float x[DIMENSIONS], float v[DIMENSIONS])
+{
+ assert(body != NULL);
+ body->m = m;
+ for (unsigned i = 0; i < DIMENSIONS; ++i)
+ {
+ body->x[i] = x[i];
+ body->v[i] = v[i];
+ body->F[i] = 0.00;
+ }
+ body->exists = true;
+}
+
+/**
+ * @function Body_Destroy
+ * @purpose Destroy a Body
+ * @param body - Target Body
+ */
+void Body_Destroy(Body * body)
+{
+ assert(body != NULL);
+ body->exists = false;
+}
+
+/**
+ * @function Body_Step
+ * @purpose Perform a single step, updating target Body
+ * @param body - Target Body
+ * @param system - System that the Body is part of (used for universal constants, etc)
+ */
+void Body_Step(Body * body, System * system)
+{
+ assert(body != NULL && system != NULL && body->exists);
+
+ //fprintf(stderr, "DEBUG - Body %p\t%d\n", (void*)(body), body->exists);
+ for (unsigned i = 0; i < DIMENSIONS; ++i)
+ {
+ body->x[i] += body->v[i] * system->dt; //Update the position vector
+ body->v[i] += (body->F[i] / body->m) * system->dt; //Update the velocity vector
+ body->F[i] = 0.00; //Reset the net force vector
+ }
+
+ //Do collision checking
+
+}
+
+/**
+ * @function Body_Force
+ * @purpose Update the net force on Body A due to Body B
+ * @param A - Body to update; subject of force from B
+ * @param B - Body causing force on A; NOT CHANGED BY THIS FUNCTION
+ * @param system - System which Body A is a part of (used for universal constants, etc)
+ */
+void Body_Force(Body * A, Body * B, System * system)
+{
+
+ assert(A != NULL && B != NULL && system != NULL && A->exists && B->exists);
+ if (A == B)
+ return;
+
+ float r2 = 0.00;
+ for (unsigned i = 0; i < DIMENSIONS; ++i)
+ r2 += powf((B->x[i] - A->x[i]), 2);
+
+ float k = system->G * (A->m * B->m) / powf(r2, 1.50);
+ for (unsigned i = 0; i < DIMENSIONS; ++i)
+ A->F[i] += k * (B->x[i] - A->x[i]);
+
+}
+
+/**
+ * @function Body_Radius
+ * @purpose Calculate the radius of a body
+ * @param body - Target Body
+ * @returns The radius of body
+ */
+float Body_Radius(Body * body)
+{
+ assert(body != NULL && body->exists);
+ float radius = 5.0;
+ return radius;
+}
+
+/**
+ * @function Body_WriteData
+ * @purpose Write data about the Body to a file
+ * @param body - The target Body
+ * @param system - System to which the Body belongs
+ * @param file - File to write data to
+ */
+void Body_WriteData(Body * body, System * system, FILE * file)
+{
+ assert(body != NULL && system != NULL && file != NULL && body->exists);
+ fprintf(file,"%d\t%p\t%f", system->step, (void*)(body), body->m);
+ for (unsigned i = 0; i < DIMENSIONS; ++i)
+ fprintf(file,"\t%f", body->x[i]);
+ for (unsigned i = 0; i < DIMENSIONS; ++i)
+ fprintf(file,"\t%f", body->v[i]);
+ for (unsigned i = 0; i < DIMENSIONS; ++i)
+ fprintf(file,"\t%f", body->F[i]);
+
+ fprintf(file,"\n");
+
+}
+
+// --- System functions --- //
+
+/**
+ * @function System_Init
+ * @purpose Initialise a System object
+ * @param system - Target System
+ * @param x - The origin of the System (currently used for offset of graphics)
+ * @param dt - Time step
+ * @param G - Universal Gravitational Constant
+ * @param nMax - Number of Body objects this System will contain (DEFAULT = 10)
+ */
+void System_Init(System * system, float x[DIMENSIONS], float dt, float G, unsigned nMax)
+{
+ assert(system != NULL);
+ for (unsigned i = 0; i < DIMENSIONS; ++i)
+ system->x[i] = x[i];
+
+
+ system->dt = dt;
+ system->G = G;
+ system->nMax = nMax;
+ system->n = 0;
+
+ system->bodies = (Body*)(calloc(sizeof(Body), nMax)); //Allocate space for nMax bodies
+ for (unsigned a = 0; a < nMax; ++a)
+ (system->bodies+a)->exists = false;
+}
+
+/**
+ * @function System_Destroy
+ * @purpose Destroy all Bodies in a system
+ * @param system - Target system
+ */
+
+void System_Destroy(System * system)
+{
+ assert(system != NULL);
+ if (system->bodies != NULL)
+ free(system->bodies);
+ system->nMax = 0;
+}
+
+/**
+ * @function System_Step
+ * @purpose Compute single step for all Bodies in System
+ * @param system - Target system
+ */
+void System_Step(System * system)
+{
+ assert(system != NULL);
+
+ // Update forces
+ for (unsigned a = 0; a < system->nMax; ++a)
+ {
+ if (system->bodies[a].exists)
+ {
+ for (unsigned b = 0; b < system->nMax; ++b)
+ {
+ //fprintf(stderr, "DEBUG1 - %p\t%p\t%d\t%d\n", (void*)(system->bodies+a), (void*)(system->bodies+b), system->bodies[a].exists, system->bodies[b].exists);
+ if (system->bodies[b].exists)
+ {
+ //fprintf(stderr, "DEBUG2 - %p\t%p\t%d\t%d\n", (void*)(system->bodies+a), (void*)(system->bodies+b), system->bodies[a].exists, system->bodies[b].exists);
+ Body_Force(&(system->bodies[a]), &(system->bodies[b]), system); //Update force on bodies[a] due to bodies[b]
+ //fprintf(stderr, "DEBUG3 - %p\t%p\t%d\t%d\n", (void*)(system->bodies+a), (void*)(system->bodies+b), system->bodies[a].exists, system->bodies[b].exists);
+ }
+ }
+ }
+ }
+
+ // Perform steps
+ for (unsigned a = 0; a < system->nMax; ++a)
+ {
+ if ((system->bodies+a)->exists)
+ Body_Step(system->bodies+a, system);
+ }
+}
+
+/**
+ * @function System_AddBody
+ * @purpose Add a Body to a System. Resizes the System's bodies array if necessary.
+ * @param system - Target System
+ * @param m - Mass of the new Body
+ * @param x[...] - position vector for the new Body
+ * @param v[...] - velocity vector for the new Body
+ * @returns Pointer to new Body
+ */
+
+Body * System_AddBody(System * system, float m, float x[DIMENSIONS], float v[DIMENSIONS])
+{
+ assert(system != NULL);
+
+ unsigned a = 0;
+ for (a = 0; a < system->nMax; ++a)
+ {
+ if (!(system->bodies+a)->exists)
+ break;
+ }
+
+ if (a >= system->nMax)
+ {
+ system->nMax = (system->nMax+1) * 2;
+ system->bodies = (Body*)(realloc(system->bodies, sizeof(Body) * system->nMax));
+ for (unsigned b = a; b < system->nMax; ++b)
+ (system->bodies+b)->exists = false;
+ }
+
+ Body_Init(system->bodies+a, m, x, v);
+ system->n += 1;
+
+ return (system->bodies+a);
+
+}
+
+/**
+ * @function System_WriteData
+ * @purpose Write all information about System to a file
+ * @param system - Target System
+ * @param file - Target FILE (must be opened for writing)
+ */
+void System_WriteData(System * system, FILE * file)
+{
+ assert(system != NULL && file != NULL);
+ for (unsigned a = 0; a < system->nMax; ++a)
+ {
+ if ((system->bodies+a)->exists)
+ Body_WriteData(system->bodies+a, system, file);
+ }
+ fprintf(file, "\n");
+
+}
+
+// --- Graphical functions --- //
+
+#ifdef _GRAPHICS_H
+
+/**
+ * @function Body_Draw
+ * @purpose Draw an individual Body
+ * @param body - Target Body
+ * @param system - System to which the Body belongs
+ */
+
+void Body_Draw(Body * body, System * system)
+{
+ assert(body != NULL && system != NULL && body->exists);
+ int pixel_pos[DIMENSIONS];
+ for (unsigned i = 0; i < DIMENSIONS; ++i)
+ pixel_pos[i] = (int)(system->x[i] - body->x[i]);
+
+ float rgb[3] = {0.0f, 0.0f, 0.0f};
+ if (body->m >= 100.0)
+ {
+ rgb[0] = 0.0f; rgb[1] = 0.0f; rgb[2] = 1.0f;
+ }
+ if (body->m >= 10.0)
+ {
+ rgb[0] = 0.5f; rgb[1] = 0.5f; rgb[2] = 1.0f;
+ }
+ else if (body->m >= 1.0)
+ {
+ rgb[0] = 0.8f; rgb[1] = 0.2f; rgb[2] = 0.1f;
+ }
+ else if (body->m >= 0.1)
+ {
+ rgb[0] = 0.0f; rgb[1] = 0.5f; rgb[2] = 0.0f;
+ }
+ else if (body->m >= 0.01)
+ {
+ rgb[0] = 0.5f; rgb[1] = 0.5f; rgb[2] = 0.5f;
+ }
+
+
+ Graphics_Circle(pixel_pos,Body_Radius(body), rgb[0], rgb[1], rgb[2]);
+
+ //fprintf(stdout, "DEBUG - Body_Draw - Called\n");
+}
+
+/**
+ * @function System_Draw
+ * @purpose Draw a System of Bodies
+ * @param system - Target System to draw
+ */
+void System_Draw(System * system)
+{
+ assert(system != NULL);
+ for (unsigned a = 0; a < system->nMax; ++a)
+ {
+ if ((system->bodies+a)->exists)
+ Body_Draw(system->bodies+a, system);
+ }
+}
+
+/**
+ * @function Process_Events
+ * @purpose Handle any SDL events recieved.
+ */
+void Process_Events()
+{
+ SDL_Event event;
+ while (SDL_PollEvent(&event))
+ {
+ switch(event.type)
+ {
+ case SDL_QUIT:
+ exit(EXIT_SUCCESS);
+ break;
+ }
+ }
+ //SDL_Delay(1);
+}
+
+#endif //_GRAPHICS_H
+
+//EOF
+
+