diff --git a/sim/CMakeLists.txt b/sim/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..2ebd04ca3698fb274934fa51755544e8f6dbe9b2
--- /dev/null
+++ b/sim/CMakeLists.txt
@@ -0,0 +1,10 @@
+cmake_minimum_required(VERSION 3.10)
+
+project(main LANGUAGES C)
+
+
+link_libraries( -lGL -lglfw -lm -ggdb -lgraphene-1.0 -g -lpthread -O2)
+
+add_executable(main ${PROJECT_SOURCE_DIR}/src/main.c ${PROJECT_SOURCE_DIR}/src/glad.c)
+include_directories( PRIVATE include)
+include_directories(PRIVATE /usr/lib64/graphene-1.0/include) #beware: not system agnostic
diff --git a/sim/include/importData.h b/sim/include/importData.h
new file mode 100644
index 0000000000000000000000000000000000000000..b9d0806be8a6647aa36fb404fa933f20e6f292a3
--- /dev/null
+++ b/sim/include/importData.h
@@ -0,0 +1,85 @@
+#ifndef SIMDATA
+#define SIMDATA
+typedef struct {
+ float x;
+ float y;
+ float z;
+} Vector;
+
+typedef struct {
+ long noOfBodies;
+ long noOfRecords;
+ float *allocator;
+ float *traceAllocator;
+ char *filename;
+ long size;
+ bool run;
+ int speed;
+ bool trace;
+ unsigned int *indices;
+ unsigned int *count;
+ unsigned int noOfTracepoints;
+} SimData;
+
+void importData(SimData *sim) {
+
+ char *charptr1, *charptr2;
+ FILE *output = fopen(sim -> filename, "r");
+
+ //get the file size to make sure we have enough space to read it!
+ struct stat filestats;
+ stat(sim->filename, &filestats);
+ long size = filestats.st_size; //since we read line by line this is large enough!
+ //buffer is freed at the end of this scope
+ char* buffer = (char*)malloc(size);
+
+ //random line ignore
+ getline(&buffer, &size, output);
+ printf("Vertex File opened: %s\n", buffer);
+
+ //this line has information about the number of Record Positions
+ getline(&buffer, &size, output);
+ charptr1 = &buffer[0];
+ sim->noOfRecords = strtol(charptr1, &charptr2, 0);
+ printf("number of Recorded Points: :%ld\n", sim -> noOfRecords);
+
+ //this line tells us how many bodies are part of the system
+ getline(&buffer, &size, output);
+ charptr1 = &buffer[0];
+ sim -> noOfBodies = strtol(charptr1, &charptr2, 0);
+ printf("number of Bodies:::%ld\n", sim -> noOfBodies);
+
+ //now we can use the noOfBodies and the NoOfRecords to allocate enough memory to
+ //hold all the positions in memory
+ //pointer is owned by caller and will be freed when the program terminates
+ sim -> size = sim -> noOfRecords * sim -> noOfBodies*3;
+ sim -> allocator = malloc(size*sizeof(float));
+ if (!sim->allocator) {
+ exit(1);
+ }
+ sim -> traceAllocator = malloc(size*sizeof(float));
+ if (!sim->traceAllocator) {
+ exit(1);
+ }
+ charptr1 = &buffer[0];
+
+ for (long record = 0; record < sim->noOfRecords; record++) {
+ long offset = record * 3*sim -> noOfBodies;
+ long offset_trace = sim->noOfRecords;
+
+ getline(&buffer, &size, output);
+ charptr1 = &buffer[0];
+ for (long body = 0; body <sim -> noOfBodies*3; body+=3) {
+ sim->allocator[offset + body] = strtof(charptr1, &charptr2);
+ sim->allocator[offset + body+1] = strtof(charptr2, &charptr2);
+ sim->allocator[offset + body+2] = strtof(charptr2, &charptr2);
+ sim->traceAllocator[body*offset_trace +record*3] = sim->allocator[offset + body];
+ sim->traceAllocator[body*offset_trace +1+record*3] = sim->allocator[offset + body+1];
+ sim->traceAllocator[body*offset_trace +2+record*3] = sim->allocator[offset + body+2];
+ charptr1 = charptr2;
+ }
+
+ }
+ free(buffer);
+}
+#endif
diff --git a/sim/include/main.h b/sim/include/main.h
new file mode 100644
index 0000000000000000000000000000000000000000..3e47dd471c22664ba076fe158fe6fc7efab692e2
--- /dev/null
+++ b/sim/include/main.h
@@ -0,0 +1,52 @@
+#ifndef MAIN
+#define MAIN
+
+#include "importData.h"
+#include "Camera.h"
+
+enum textureFormat {png, jpg};
+struct texture{
+ char* name;
+ enum textureFormat format;
+ bool flip;
+};
+
+void processInput(
+ GLFWwindow *window,
+ int key,
+ int scancode,
+ int action,
+ int mods
+);
+
+void updateCamera();
+void mouse_callback(GLFWwindow* window, double xpos, double ypos);
+void framebuffer_size_callback(GLFWwindow* window, int width, int height);
+void compileAndAttachShaders( unsigned int Program,char vertexSource[], char fragmentSource[]);
+
+void loadTexture(struct texture Tex) {
+
+ int width, height, nrChannels;
+ stbi_set_flip_vertically_on_load(Tex.flip);
+
+ unsigned char* data = stbi_load(Tex.name, &width, &height, &nrChannels, 0);
+ if (data) {
+
+ switch (Tex.format) {
+ case jpg:
+ glTexImage2D(GL_TEXTURE_2D, 0, GL_RGB, width, height, 0, GL_RGB, GL_UNSIGNED_BYTE, data);
+ break;
+ case png:
+ glTexImage2D(GL_TEXTURE_2D, 0, GL_RGB, width, height, 0, GL_RGBA, GL_UNSIGNED_BYTE, data);
+ break;
+ }
+
+ glGenerateMipmap(GL_TEXTURE_2D);
+ } else {
+ printf("Failed to load texture");
+ exit(EXIT_FAILURE);
+ }
+ stbi_image_free(data);
+
+}
+#endif
diff --git a/sim/include/openGLProgram.h b/sim/include/openGLProgram.h
new file mode 100644
index 0000000000000000000000000000000000000000..58ec78e0c730a9a96013899fbfe1540b446414eb
--- /dev/null
+++ b/sim/include/openGLProgram.h
@@ -0,0 +1,86 @@
+#include <glad/glad.h>
+
+#ifndef openGLProgram
+#define openGLProgram
+void compileAndAttachShaders( unsigned int Program, char vertexSource[], char fragmentSource[]) {
+
+ unsigned int vertexShader, fragmentShader;
+ char vertexbuffer[2024], fragmentbuffer[2024], infoLog[512];
+ int success;
+ long charsRead;
+
+ FILE* SourceFile = fopen(vertexSource, "r");
+ if (SourceFile) {
+ charsRead = fread(vertexbuffer,sizeof(char), sizeof(vertexbuffer)/sizeof(char), SourceFile);
+ vertexbuffer[charsRead] = '\0';
+
+ const char* vertexptr = &vertexbuffer[0];
+
+ fclose(SourceFile);
+ if ( charsRead > 0&& charsRead < sizeof(vertexbuffer)/sizeof(char)) {
+ vertexShader = glCreateShader(GL_VERTEX_SHADER);
+ glShaderSource(vertexShader, 1, &vertexptr, NULL);
+ } else {
+ printf("Failed to open %s as vertex Source File\n", vertexSource);
+ exit(EXIT_FAILURE);
+ }
+ } else {
+ printf("Failed to open file: %s\n", vertexSource);
+ exit(EXIT_FAILURE);
+ }
+ glCompileShader(vertexShader);
+
+ glGetShaderiv(vertexShader, GL_COMPILE_STATUS, &success);
+ if (!success) {
+ glGetShaderInfoLog(vertexShader, 512, NULL, infoLog);
+ printf("ERROR::SHADER::VERTEX::COMPILATION_FAILED\n");
+ printf("%s\n", infoLog);
+ printf("\n\n%s", vertexbuffer);
+ return exit(EXIT_FAILURE);
+ }
+
+
+ SourceFile = fopen(fragmentSource, "r");
+ if (SourceFile) {
+ charsRead = fread(fragmentbuffer,sizeof(char), sizeof(fragmentbuffer)/sizeof(char), SourceFile);
+ fragmentbuffer[charsRead] = '\0';
+ const char* fragmentptr = &fragmentbuffer[0];
+ fclose(SourceFile);
+ if (charsRead > 0 && charsRead < sizeof(fragmentbuffer)/sizeof(char)) {
+ fragmentShader = glCreateShader(GL_FRAGMENT_SHADER);
+ glShaderSource(fragmentShader, 1, &fragmentptr, NULL);
+ } else {
+ printf("Failed to open %s as vertex Source File\n", fragmentSource);
+ exit(EXIT_FAILURE);
+ }
+ } else {
+ printf("Failed to open file: %s\n", fragmentSource);
+ exit(EXIT_FAILURE);
+ }
+ glCompileShader(fragmentShader);
+
+ glGetShaderiv(fragmentShader, GL_COMPILE_STATUS, &success);
+ if (!success) {
+ glGetShaderInfoLog(fragmentShader, 512, NULL, infoLog);
+ printf("ERROR::SHADER::FRAGMENT::COMPILATION_FAILED\n");
+ printf("%s\n", infoLog);
+ printf("\n\n%s", fragmentbuffer);
+ return exit(EXIT_FAILURE);
+ }
+
+ glAttachShader(Program, vertexShader);
+ glAttachShader(Program, fragmentShader);
+ glLinkProgram(Program);
+
+ glGetProgramiv(Program, GL_LINK_STATUS, &success);
+ if (!success) {
+ glGetProgramInfoLog(Program, 512, NULL, infoLog);
+ printf("ERROR::GL::Program::LINKING_FAILES\n");
+ printf("%s\n", infoLog);
+ return exit(EXIT_FAILURE);
+ }
+ glDeleteShader(vertexShader);
+ glDeleteShader(fragmentShader);
+}
+
+#endif
diff --git a/sim/shaders/fragmentshader.fs b/sim/shaders/fragmentshader.fs
new file mode 100644
index 0000000000000000000000000000000000000000..614339e1838644dca8bd63f50e757f60b0c9bca3
--- /dev/null
+++ b/sim/shaders/fragmentshader.fs
@@ -0,0 +1,39 @@
+#version 450 core
+
+
+out vec4 FragColor;
+uniform bool Trace;
+
+void point(in vec2 gl_PointCoord, out vec4 FragColor) {
+ vec2 center = gl_PointCoord-vec2(0.5,0.5);
+ float radius = dot(center, center);
+ if (radius > 0.25)
+ {
+ gl_FragDepth = 0.999;
+ discard;
+ } if ( radius <0.13)
+ {
+ FragColor = vec4(1.0,1.0,1.0,1.0);
+ } else
+ {
+ float brightness = -4 * (radius-0.13) +0.8;
+ FragColor = vec4(brightness, brightness, brightness, 0.5);
+ }
+
+}
+void trace(in vec2 gl_PointCoord, out vec4 FragColor) {
+
+
+ FragColor = vec4(0.84,0.36,0.055,1.0);
+
+}
+
+void main()
+{
+ gl_FragDepth = gl_FragCoord.z;
+ if (Trace) {
+ trace(gl_PointCoord, FragColor);
+ } else {
+ point(gl_PointCoord, FragColor);
+ }
+}
diff --git a/sim/shaders/gui.fs b/sim/shaders/gui.fs
new file mode 100644
index 0000000000000000000000000000000000000000..cfda982046d4bbffe944ae28648e85b43c6edd19
--- /dev/null
+++ b/sim/shaders/gui.fs
@@ -0,0 +1,11 @@
+#version 330 core
+
+out vec4 FragColor;
+in vec4 vertexColor;
+in float depth;
+
+
+void main()
+{
+ FragColor = vertexColor;
+}
diff --git a/sim/shaders/gui.vs b/sim/shaders/gui.vs
new file mode 100644
index 0000000000000000000000000000000000000000..48517ee31605b1548a260a278aca0d2001b453de
--- /dev/null
+++ b/sim/shaders/gui.vs
@@ -0,0 +1,14 @@
+#version 330 core
+layout (location = 0) in vec3 aPos;
+out vec4 vertexColor;
+
+uniform mat4 view;
+
+void main()
+{
+ gl_Position = vec4(aPos, 1.0);
+ gl_PointSize = 30;
+ vertexColor = vec4(0.45,0.85,0.88,1.0);
+
+
+}
diff --git a/sim/shaders/vertexshader.vs b/sim/shaders/vertexshader.vs
new file mode 100644
index 0000000000000000000000000000000000000000..798d93c270ffd4b1ff9c9924fdabba88c9774602
--- /dev/null
+++ b/sim/shaders/vertexshader.vs
@@ -0,0 +1,22 @@
+#version 330 core
+layout (location = 1) in vec3 aPos;
+//layout (location = 16) in float
+//out vec4 vertexColor;
+
+
+uniform mat4 view;
+uniform mat4 proj;
+uniform float worldScale;
+uniform float pointScale;
+void main()
+{
+ mat4 scaling = mat4(worldScale);
+ scaling[3][3] = 1.0f;
+ vec4 World = proj*inverse(view)*scaling*vec4(aPos, 1.0);
+ gl_PointSize = 1/(World.z)*pointScale;
+
+ gl_Position = World;
+ //vertexColor = vec4(1,1,1,0.5);
+
+
+}
diff --git a/sim/src/main.c b/sim/src/main.c
new file mode 100644
index 0000000000000000000000000000000000000000..5e956dca06f384c09dea45108fae9fbb5e739fe0
--- /dev/null
+++ b/sim/src/main.c
@@ -0,0 +1,440 @@
+//open gl related libraries
+#include <glad/glad.h>
+#include <GLFW/glfw3.h>
+
+//libc imports
+#include <stdio.h>
+#include <stdbool.h>
+#include <math.h>
+#include <sys/stat.h>
+#include <threads.h>
+
+//image library for importing pictures
+#define STB_IMAGE_IMPLEMENTATION
+#include "stb_image.h"
+
+//linear algebra library
+#include <graphene-1.0/graphene.h>
+
+//project specific files
+#include "main.h"
+#include "importData.h"
+#include "openGLProgram.h"
+
+
+float pause_indicator[] = {
+ //first triangle
+
+ 0.90f, 0.98f, -1.0, // top right
+ 0.94f, 0.8f, -1.0, // bottom right
+ 0.90f, 0.8f, -1.0,
+
+ 0.90f, 0.98f, -1.0, // top right
+ 0.94f, 0.98f, -1.0, // bottom right
+ 0.94f, 0.8f, -1.0,
+
+ 0.95f, 0.98f, -1.0, // top right
+ 0.99f, 0.8f, -1.0, // bottom right
+ 0.95f, 0.8f, -1.0,
+
+ 0.95f, 0.98f, -1.0, // top right
+ 0.99f, 0.98f, -1.0, // bottom right
+ 0.99f, 0.8f, -1.0,
+
+
+};
+
+float play_indicator[] = {
+ 0.9f, 1.0f, -1.0, // top right
+ 1.0f, 0.9f, -1.0, // bottom right
+ 0.9f, 0.8f, -1.0, // bottom left
+};
+
+unsigned int indices[] = {
+ 0, 1, 3, //first triangle
+ 1, 2, 3 //second triangle
+};
+
+float deltaTime = 0.0f; // Time between current frame and last frame
+float lastFrame = 0.0f;
+float lastX = 0, lastY = 0;
+
+graphene_vec3_t direction;
+Camera cam;
+SimData sim;
+GLFWwindow *window;
+
+
+void GLAPIENTRY
+MessageCallback( GLenum source,
+ GLenum type,
+ GLuint id,
+ GLenum severity,
+ GLsizei length,
+ const GLchar* message,
+ const void* userParam ) {
+ fprintf( stderr, "GL CALLBACK: %s type = 0x%x, severity = 0x%x, message = %s\n",
+ ( type == GL_DEBUG_TYPE_ERROR ? "** GL ERROR **" : "" ),
+ type, severity, message );
+}
+
+typedef struct {
+ float distance;
+ unsigned int index;
+} TransparentPoint;
+
+int main(int argc, char *argv[]) {
+
+ glfwInit();
+
+ char *filename = argv[1];
+
+ sim.filename = filename;
+ sim.allocator = NULL;
+ sim.run = false;
+ sim.speed = 1;
+ sim.noOfTracepoints = 1000;
+ sim.trace = false;
+ importData(&sim);
+ TransparentPoint *root = malloc(sizeof(TransparentPoint) * (sim.noOfBodies));
+ if (!root) {
+ exit(1);
+ }
+
+
+
+ //init glfw and configure the enum that is the state of glfw
+ glfwWindowHint(GLFW_CONTEXT_VERSION_MAJOR, 4);
+ glfwWindowHint(GLFW_CONTEXT_VERSION_MINOR, 5);
+ glfwWindowHint(GLFW_OPENGL_PROFILE, GLFW_OPENGL_CORE_PROFILE);
+ glfwWindowHint(GLFW_SAMPLES,1);
+ glfwWindowHint(GLFW_CLIENT_API, GLFW_OPENGL_API);
+ glfwWindowHint(GLFW_CONTEXT_CREATION_API, GLFW_NATIVE_CONTEXT_API);
+
+ glfwWindowHint(GLFW_DEPTH_BITS, 24);
+ window = glfwCreateWindow(1000,1000,"Visualize", NULL, NULL);
+ if (window == NULL) {
+ printf("Failed to create GLFW window");
+ glfwTerminate();
+ return -1;
+ }
+ glfwSetKeyCallback(window, processInput);
+ glfwSetInputMode(window, GLFW_CURSOR, GLFW_CURSOR_DISABLED);
+
+
+ glfwMakeContextCurrent(window);
+ if (!gladLoadGLLoader((GLADloadproc)glfwGetProcAddress)) {
+ printf("failed to initialize GLAD");
+ return -1;
+ }
+
+ glViewport(0, 0, 1000, 1000);
+ glfwSetFramebufferSizeCallback(window, framebuffer_size_callback);
+ glEnable ( GL_DEBUG_OUTPUT );
+ glDebugMessageCallback( MessageCallback, 0 );
+
+ //to change things related to GL_POINT Primitive
+ glEnable(GL_PROGRAM_POINT_SIZE);
+
+
+ glEnable(GL_LINE_SMOOTH);
+ glLineWidth(3.0);
+ glEnable(GL_BLEND);
+ glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
+ unsigned int shaderProgram;
+ shaderProgram = glCreateProgram();
+
+ //Prepare the shaders
+ compileAndAttachShaders( shaderProgram, "./shaders/vertexshader.vs", "./shaders/fragmentshader.fs");
+ glUseProgram(shaderProgram);
+ glEnable(GL_DEPTH_TEST);
+ glDepthFunc(GL_LESS);
+ unsigned int VBO, EBO, VAO;
+ glGenVertexArrays(1,&VAO);
+ glGenBuffers(1,&EBO);
+ glGenBuffers(1, &VBO);
+
+ //after binding the VAO can be used
+ glBindVertexArray(VAO);
+
+ glBindBuffer(GL_ARRAY_BUFFER, VBO);
+ glBufferData(GL_ARRAY_BUFFER, sim.size*sizeof(float), sim.allocator, GL_DYNAMIC_DRAW);
+
+ glVertexAttribPointer(1,3,GL_FLOAT, GL_FALSE, 0, (void*)0);
+ glEnableVertexAttribArray(1);
+ unsigned int VAO_trace, VBO_trace;
+ glGenVertexArrays(1,&VAO_trace);
+ glGenBuffers(1,&VBO_trace);
+
+ //after binding the VAO can be used
+ glBindVertexArray(VAO_trace);
+
+ glBindBuffer(GL_ARRAY_BUFFER, VBO_trace);
+ int offset = 1;
+ glBufferData(GL_ARRAY_BUFFER, sim.size*sizeof(float), sim.traceAllocator, GL_DYNAMIC_DRAW);
+ glVertexAttribPointer(1,3,GL_FLOAT, GL_FALSE,offset *3* sizeof(float), (void*)0);
+ glEnableVertexAttribArray(1);
+
+
+ float arr[16];
+ graphene_matrix_t models;
+ graphene_vec3_t vector;
+ graphene_matrix_t projection_mat;
+ graphene_matrix_t* projection = graphene_matrix_init_perspective(&projection_mat, 45,1,0.1,10000.0);
+ graphene_matrix_to_float(projection, arr);
+
+ glUniformMatrix4fv(glGetUniformLocation( shaderProgram, "proj"), 1, GL_FALSE, arr);
+
+ int t = 0;
+ graphene_point3d_t point;
+
+ graphene_vec3_t* eye =graphene_vec3_init(&cam.position, 0.0, 0.0, 3);
+
+ CameraInit(&cam, shaderProgram);
+ glfwSetCursorPos(window, 0, 0);
+ glUniform1f(glGetUniformLocation(cam.program, "pointScale"), cam.pointScale);
+
+ unsigned int gui;
+ gui = glCreateProgram();
+
+ compileAndAttachShaders( gui, "./shaders/gui.vs", "./shaders/gui.fs");
+ glUseProgram(gui);
+ unsigned int VBO_play, VAO_play, VBO_pause, VAO_pause;
+ glGenVertexArrays(1,&VAO_play);
+ glGenVertexArrays(1,&VAO_pause);
+ glGenBuffers(1,&VBO_play);
+ glGenBuffers(1,&VBO_pause);
+
+ glBindVertexArray(VAO_play);
+ glBindBuffer(GL_ARRAY_BUFFER,VBO_play);
+ glBufferData(GL_ARRAY_BUFFER, 9*sizeof(float), play_indicator, GL_STATIC_DRAW);
+ glVertexAttribPointer(0,3,GL_FLOAT, GL_FALSE, 0, (void*)0);
+ glEnableVertexAttribArray(0);
+ glBindBuffer(GL_ARRAY_BUFFER, 0);
+
+ glBindVertexArray(VAO_pause);
+ glBindBuffer(GL_ARRAY_BUFFER,VBO_pause);
+ glBufferData(GL_ARRAY_BUFFER,36*sizeof(float), pause_indicator, GL_STATIC_DRAW);
+ glVertexAttribPointer(0,3,GL_FLOAT, GL_FALSE, 0, (void*)0);
+ glEnableVertexAttribArray(0);
+ glBindBuffer(GL_ARRAY_BUFFER, 0);
+
+ graphene_vec4_t zvec;
+ graphene_vec4_t vertex;
+ int trace = sim.noOfTracepoints;
+ while (!glfwWindowShouldClose(window) && t < sim.noOfRecords) {
+ //empty screen
+ glDepthMask(GL_TRUE);
+ glClear(GL_COLOR_BUFFER_BIT);
+ glClear(GL_DEPTH_BUFFER_BIT);
+ glClearColor(0.0f,0.0f,0.0f,0.0f);
+
+ //render bodies
+ glUseProgram(shaderProgram);
+ glfwPollEvents();
+ updateCamera();
+ graphene_matrix_get_row(&cam.view, 2, &zvec);
+
+ for (int k = 0; k < sim.noOfBodies; k++) {
+ graphene_vec4_t *tmp = graphene_vec4_init(
+ &vertex,
+ sim.traceAllocator[t+k*sim.noOfRecords*3],
+ sim.traceAllocator[t+k*sim.noOfRecords*3+1],
+ sim.traceAllocator[t+k*sim.noOfRecords*3+2],
+ 1.0f);
+
+ root[k].distance = graphene_vec4_dot(&vertex, &zvec);
+ root[k].index = k;
+ }
+ glBindVertexArray(VAO_trace);
+ glUniform1f(glGetUniformLocation(cam.program, "Trace"), 0);
+ for (int body = 0; body < sim.noOfBodies; body++) {
+ glDrawArrays(GL_POINTS, t+body*sim.noOfRecords,1);
+ }
+ int trace = sim.noOfTracepoints;
+ if (sim.trace == true) {
+ glBindVertexArray(VAO_trace);
+ glUniform1f(glGetUniformLocation(cam.program, "Trace"), 1);
+
+ for (int j = 0; j <sim.noOfBodies; j++) {
+ glDrawArrays(
+ GL_LINE_STRIP,
+ sim.noOfRecords*j/1+ fmax(0, t-trace),
+ fmin(t/1, trace)
+ );
+ }
+ }
+
+ //render gui
+ glUseProgram(gui);
+ glBindVertexArray(VAO_play);
+ if (sim.run) {
+ t+= sim.speed;
+ glBindVertexArray(VAO_play);
+
+ } else {
+ glBindVertexArray(VAO_pause);
+ sim.run = 0;
+ }
+ if (t == 0 && sim.speed < 0) {
+ sim.run = 0;
+ sim.speed = 1;
+ }
+ glDrawArrays(GL_TRIANGLES, 0,3+9*(1-sim.run));
+
+ //swap
+ glfwSwapBuffers(window);
+ }
+ glfwTerminate();
+ free(root);
+ return 0;
+}
+
+void processInput(
+ GLFWwindow *window,
+ int key,
+ int scancode,
+ int action,
+ int mods
+ ) {
+
+ if (key ==GLFW_KEY_SPACE && action == GLFW_PRESS) {
+ sim.run = -1*(sim.run -1);
+ }
+
+ else if (key == GLFW_KEY_COMMA && action == GLFW_PRESS) {
+ sim.speed -= 1;
+ }
+ else if (key == GLFW_KEY_PERIOD && action == GLFW_PRESS) {
+ sim.speed += 1;
+ }
+ else if (key == GLFW_KEY_C && action == GLFW_PRESS) {
+ sim.trace = (sim.trace-1)*-1;
+ }
+ else if (key == GLFW_KEY_M && action == GLFW_PRESS) {
+ sim.noOfTracepoints =(int)sim.noOfTracepoints* 2;
+ }
+ else if (key == GLFW_KEY_N && action == GLFW_PRESS) {
+ sim.noOfTracepoints = (int)sim.noOfTracepoints/2;
+ }
+
+
+}
+
+void updateCamera () {
+ //acquire the lock to write to the fields of the Camera
+ float arr[16];
+ double xpos, ypos;
+ graphene_vec3_t tmp;
+
+ float currentFrame = glfwGetTime();
+ deltaTime = currentFrame - lastFrame;
+ printf("Frametime : %lf\n", deltaTime);
+ lastFrame = currentFrame;
+
+ float speed = cam.speed*2;
+ glfwGetCursorPos(window, &xpos, &ypos);
+ glfwSetCursorPos(window, 0,0 );
+
+ cam.phi += cam.sensitivity *(float)(-xpos);
+ cam.theta += cam.sensitivity * (float)(ypos );
+ graphene_vec3_t *left = graphene_vec3_init(
+ &cam.left,
+ cos(cam.phi),
+ 0,
+ -sin(cam.phi)
+ );
+
+ graphene_vec3_cross(&cam.left, &cam.direction, &cam.up);
+ graphene_vec3_negate(&cam.up, &cam.up);
+
+
+
+ graphene_vec3_t *direction = graphene_vec3_init(
+ &cam.direction,
+ sin(cam.phi) * cos(cam.theta),
+ sin(cam.theta),
+ cos(cam.phi) * cos(cam.theta)
+
+ );
+
+ if (glfwGetKey(window, GLFW_KEY_ESCAPE) == GLFW_PRESS)
+ glfwSetWindowShouldClose(window, true);
+
+ if (glfwGetKey(window, GLFW_KEY_W) == GLFW_PRESS) {
+ graphene_vec3_scale( &cam.direction,(speed *1.2)* deltaTime, &tmp);
+ graphene_vec3_add(&cam.position, &tmp, &cam.position);
+
+
+ }
+ if (glfwGetKey(window, GLFW_KEY_S) == GLFW_PRESS) {
+ graphene_vec3_scale( &cam.direction, (speed *1.2)* deltaTime, &tmp);
+ graphene_vec3_subtract(&cam.position, &tmp, &cam.position);
+
+ }
+ if (glfwGetKey(window, GLFW_KEY_A) == GLFW_PRESS) {
+ graphene_vec3_scale( &cam.left, speed * deltaTime, &tmp);
+ graphene_vec3_add(&cam.position, &tmp, &cam.position);
+ }
+ if (glfwGetKey(window, GLFW_KEY_D) == GLFW_PRESS) {
+ graphene_vec3_scale(&cam.left, speed * deltaTime, &tmp);
+ graphene_vec3_subtract(&cam.position, &tmp, &cam.position);
+ }
+ if (glfwGetKey(window, GLFW_KEY_LEFT_CONTROL) == GLFW_PRESS) {
+ graphene_vec3_scale( &cam.up, speed * deltaTime, &tmp);
+ graphene_vec3_add(&cam.position, &tmp, &cam.position);
+ }
+ if (glfwGetKey(window, GLFW_KEY_LEFT_SHIFT) == GLFW_PRESS) {
+ graphene_vec3_scale( &cam.up, speed * deltaTime, &tmp);
+ graphene_vec3_subtract(&cam.position, &tmp, &cam.position);
+ }
+ if (glfwGetKey(window, GLFW_KEY_P) == GLFW_PRESS) {
+ cam.pointScale *= 1.05;
+ glUniform1f(glGetUniformLocation(cam.program, "pointScale"), cam.pointScale);
+ }
+ if (glfwGetKey(window, GLFW_KEY_O) == GLFW_PRESS) {
+ cam.pointScale /= 1.05;
+ glUniform1f(glGetUniformLocation(cam.program, "pointScale"), cam.pointScale);
+ }
+ if (glfwGetKey(window, GLFW_KEY_RIGHT) == GLFW_PRESS) {
+ cam.speed *=1.1f;
+ }
+ if (glfwGetKey(window, GLFW_KEY_LEFT) == GLFW_PRESS) {
+ cam.speed *=0.9f;
+ }
+ if (glfwGetKey(window, GLFW_KEY_L) == GLFW_PRESS) {
+ cam.worldScale *= 1.1f;
+ }
+ if (glfwGetKey(window, GLFW_KEY_K) == GLFW_PRESS) {
+ cam.worldScale /= 1.1f;
+ }
+ graphene_matrix_t worldScale;
+ graphene_matrix_t *worldScale_pointer = graphene_matrix_init_scale(
+ &worldScale,
+ cam.worldScale,
+ cam.worldScale,
+ cam.worldScale
+ );
+ graphene_matrix_to_float(&worldScale, arr);
+ glUniform1f(glGetUniformLocation(cam.program, "worldScale"), cam.worldScale);
+
+
+
+ graphene_vec3_subtract(&cam.position, &cam.direction, &tmp);
+
+ graphene_matrix_t* view = graphene_matrix_init_look_at(
+ &cam.view,
+ &cam.position,
+ &tmp,
+ &cam.up
+ );
+
+
+ graphene_matrix_to_float(&cam.view, arr);
+ glUniformMatrix4fv(glGetUniformLocation( cam.program, "view"), 1, GL_FALSE, arr);
+}
+
+void framebuffer_size_callback(GLFWwindow* window, int width, int height) {
+
+ glViewport(0,0,1800,1800);
+}
diff --git a/src/.gitkeep b/src/.gitkeep
deleted file mode 100644
index e69de29bb2d1d6434b8b29ae775ad8c2e48c5391..0000000000000000000000000000000000000000
diff --git a/src/jupiter.c b/src/jupiter.c
deleted file mode 100644
index 7cbe00aa5c19e9012139a9f4c1806224eac3d08e..0000000000000000000000000000000000000000
--- a/src/jupiter.c
+++ /dev/null
@@ -1,166 +0,0 @@
-#include <stdio.h>
-#include <stdlib.h>
-#include <math.h>
-
-#define pi 3.141592653589793
-#define solar_mass 1
-#define G (4 * pi * pi)
-#define days_per_year 365.24
-#define JD 5.2
-#define VJ0 2.75535903023
-
-
-typedef struct {
- double x[3];
- double v[3];
- double a[3];
- double mass;
-} Planet;
-
-//Unit scaling: Mass in solar mass, distance in AU, velocity in AU/year, and constant G in AU^3/solar mass*year^2
-
-Planet Planets[] = {
- {
- {0,0,0},
- {0,0,0},
- {0,0,0},
- solar_mass
- },
- {
- {JD, 0, 0},
- {0, VJ0, 0},
- {0,0,0},
- 9.54791938424326609e-04 * solar_mass
- },
- {
- {-JD, 0, 0},
- {0, -VJ0, 0},
- {0, 0, 0},
- 2.85885980666130812e-04 * solar_mass
- }
-};
-static inline double distance( const double x1[3], const double x2[3]) {
- const double x = x1[0] - x2[0];
- const double y = x1[1] - x2[1];
- const double z = x1[2] - x2[2];
-
- return sqrt(x*x + y*y + z*z + 0.001);
-
-}
-
-void update(int steps){
- long run = 0;
- const int noOfBodies = sizeof(Planets)/sizeof(Planet);
- const double dt = 0.01;
-
- while (run < steps) {
- for (int l = 0; l < noOfBodies; l++) {
- Planets[l].v[0] += 0.5 * Planets[l].a[0]* dt;
- Planets[l].v[1] += 0.5 * Planets[l].a[1]* dt;
- Planets[l].v[2] += 0.5 * Planets[l].a[2]* dt;
- }
- for (int k = 0; k < noOfBodies; k++) {
- Planets[k].x[0] += dt * Planets[k].v[0];
- Planets[k].x[1] += dt * Planets[k].v[1];
- Planets[k].x[2] += dt * Planets[k].v[2];
- Planets[k].a[0] = 0;
- Planets[k].a[1] = 0;
- Planets[k].a[2] = 0;
- }
-
- for (int i = 0; i < noOfBodies; i++) {
- for (int j = i+1; j < noOfBodies; j++) {
- double dist = distance(Planets[i].x, Planets[j].x);
- double a = G/pow(dist, 3);
-
- Planets[i].a[0] += a * Planets[j].mass * ( Planets[j].x[0] - Planets[i].x[0] );
- Planets[i].a[1] += a * Planets[j].mass * ( Planets[j].x[1] - Planets[i].x[1] );
- Planets[i].a[2] += a * Planets[j].mass * ( Planets[j].x[2] - Planets[i].x[2] );
- Planets[j].a[0] += a * Planets[i].mass * ( Planets[i].x[0] - Planets[j].x[0] );
- Planets[j].a[1] += a * Planets[i].mass * ( Planets[i].x[1] - Planets[j].x[1] );
- Planets[j].a[2] += a * Planets[i].mass * ( Planets[i].x[2] - Planets[j].x[2] );
- }
- }
-
- for (int m = 0; m < noOfBodies; m++) {
- Planets[m].v[0] += 0.5 * Planets[m].a[0]* dt;
- Planets[m].v[1] += 0.5 * Planets[m].a[1]* dt;
- Planets[m].v[2] += 0.5 * Planets[m].a[2]* dt;
- }
-
- run++;}
-}
-
-int main() {
- FILE *starting = fopen("starting.txt", "w+");
- if (starting == NULL) {
- printf("HELP\n COULDNT OPEN FILE\nSIGSEVV");
- } else {
- printf("Output File opened\n");
- }
- const int noOfBodies = sizeof(Planets)/sizeof(Planet);
- double s =0;
- double ds=0.001;
-
- for(s=1.5; s<=3; s+=ds){
-
- Planets[2].x[0]=-s*JD;
- Planets[2].x[1]=0;
- Planets[2].x[2]=0;
- Planets[2].v[0]=0;
- Planets[2].v[1]=-sqrt(1/s)*VJ0;
- Planets[2].v[2]=0;
- for(int i=0; i<3;i++){
- for(int k=0; k<3; k++){
- Planets[k].a[i]=0;
- }
- }
- for(int i=0; i<3; i++){
- Planets[0].x[i]=0;
- Planets[0].v[i]=0;
- }
- Planets[1].x[0]=JD;
- Planets[1].x[1]=0;
- Planets[1].x[2]=0;
- Planets[1].v[0]=0;
- Planets[1].v[1]=VJ0;
- Planets[1].v[2]=0;
- update(10000);
- double r10= sqrt(Planets[2].x[1]*Planets[2].x[1]+Planets[2].x[0]*Planets[2].x[0]);
- update(10000);
- double r1= sqrt(Planets[2].x[1]*Planets[2].x[1]+Planets[2].x[0]*Planets[2].x[0]);
-
- Planets[2].x[0]=(-s*JD)-0.01;
- Planets[2].x[1]=0;
- Planets[2].x[2]=0;
- Planets[2].v[0]=0;
- Planets[2].v[1]=-sqrt(G/((s*JD)-0.01));
- Planets[2].v[2]=0;
- for(int i=0; i<3;i++){
- for(int k=0; k<3; k++){
- Planets[k].a[i]=0;
- }
- }
- for(int i=0; i<3; i++){
- Planets[0].x[i]=0;
- Planets[0].v[i]=0;
- }
- Planets[1].x[0]=JD;
- Planets[1].x[1]=0;
- Planets[1].x[2]=0;
- Planets[1].v[0]=0;
- Planets[1].v[1]=VJ0;
- Planets[1].v[2]=0;
- update(10000);
- double r20= sqrt(Planets[2].x[1]*Planets[2].x[1]+Planets[2].x[0]*Planets[2].x[0]);
- update(10000);
- double r2=sqrt(Planets[2].x[1]*Planets[2].x[1]+Planets[2].x[0]*Planets[2].x[0]);
- double d0=fabs(r10-r20);
- double d=fabs(r2-r1);
- double l= log(d/d0);
- fprintf(starting, "%lf %lf\n", r10, l);
-
- }
- fclose(starting);
-}
-
diff --git a/src/nbody.c b/src/nbody.c
deleted file mode 100644
index 60de895d661aba1d518ab77d04bd0c3e0835dc25..0000000000000000000000000000000000000000
--- a/src/nbody.c
+++ /dev/null
@@ -1,186 +0,0 @@
-#include <stdio.h>
-#include <stdlib.h>
-#include <math.h>
-
-#define pi 3.141592653589793
-#define solar_mass 1
-#define G (4 * pi * pi)
-#define days_per_year 365.24
-
-//compile options gcc -O2 main.c -Wall -lm
-
-typedef struct {
- double x[3];
- double v[3];
- double a[3];
- double mass;
-} Planet;
-
-//Unit scaling: Mass in solar mass, distance in AU, velocity in AU/year, and constant G in AU^3/solar mass*year^2
-
-Planet Planets[] = {
- {
- {0,0,0},
- {0,0,0},
- {0,0,0},
- solar_mass
- },
- {
- {4.84143144246472090e+00,
- -1.16032004402742839e+00,
- -1.03622044471123109e-01},
- {1.66007664274403694e-03 * days_per_year,
- 7.69901118419740425e-03 * days_per_year,
- -6.90460016972063023e-05 * days_per_year},
- {0,0,0},
- 9.54791938424326609e-04 * solar_mass
- },
- {
- {8.34336671824457987e+00,
- 4.12479856412430479e+00,
- -4.03523417114321381e-01},
- {-2.76742510726862411e-03 * days_per_year,
- 4.99852801234917238e-03 * days_per_year,
- 2.30417297573763929e-05 * days_per_year},
- 2.85885980666130812e-04 * solar_mass
- },
- {
- {1.28943695621391310e+01,
- -1.51111514016986312e+01,
- -2.23307578892655734e-01},
- {2.96460137564761618e-03 * days_per_year,
- 2.37847173959480950e-03 * days_per_year,
- -2.96589568540237556e-05 * days_per_year},
- 4.36624404335156298e-05 * solar_mass
- },
- {
- {1.53796971148509165e+01,
- -2.59193146099879641e+01,
- 1.79258772950371181e-01},
- {2.68067772490389322e-03 * days_per_year,
- 1.62824170038242295e-03 * days_per_year,
- -9.51592254519715870e-05 * days_per_year},
- 5.15138902046611451e-05 * solar_mass
- }
-};
-
-//TODO add softening factor
-static inline double distance( const double x1[3], const double x2[3]) {
- const double x = x1[0] - x2[0];
- const double y = x1[1] - x2[1];
- const double z = x1[2] - x2[2];
-
- return sqrt(x*x + y*y + z*z + 0.001);
-
-}
-double energy(int n){
- double e;
- int i=0, j=0, k=0;
- e = 0.0;
-
- for(i=0; i<n; i++){
- for(k=0; k<3; k++){
- e += 0.5 * Planets[i].mass * Planets[i].v[k] * Planets[i].v[k];
- }
- for(j=i+1; j < n; j++) {
- double dist = distance(Planets[i].x, Planets[j].x);
- e -= (G * Planets[i].mass * Planets[j].mass) / dist;
- }
- }
- return e;
-}
-
-static inline double radius(const double r[3]){
- const double x = r[0];
- const double y = r[1];
- const double z = r[2];
- return sqrt(x*x+y*y+z*z);
-}
-
-static inline double velocity(const double v[3]){
- const double vx = v[0];
- const double vy = v[1];
- const double vz = v[2];
- return sqrt(vx*vx+vy*vy+vz*vz);
-}
-
-double an_mom(int n){
- double L;
- int i=0;
- L=0.0;
-
- for(i=0; i<n; i++){
- double v = velocity(Planets[i].v);
- double r = radius(Planets[i].x);
- L+= Planets[i].mass * v * r;
- }
- return L;
-}
-
-
-int main(int argc, char * argv[]) {
-
- const long steps = atoi(argv[1]);
- FILE *output = fopen("output.xyz", "w+");
- if (output == NULL) {
- printf("HELP\n COULDNT OPEN FILE\nSIGSEVV"); //TODO stop here
- } else {
- printf("Output File opened\n");
- }
- //this could also be done via the preprocessor
- const int noOfBodies = sizeof(Planets)/sizeof(Planet);
- printf("size: %d\n", noOfBodies);
- const double dt = 0.01;
-
- printf("energy before: %.9f\n", energy(noOfBodies));
- printf("total angular momentum before: %.9f\n", an_mom(noOfBodies));
-
- long run = 0;
-
- while (run < steps) {
- if (run %(long int) 10== 0){
- fprintf(output, "%d\n%d.000\n", noOfBodies, run+1);
- }
- for (int l = 0; l < noOfBodies; l++) {
- Planets[l].v[0] += 0.5 * Planets[l].a[0]* dt;
- Planets[l].v[1] += 0.5 * Planets[l].a[1]* dt;
- Planets[l].v[2] += 0.5 * Planets[l].a[2]* dt;
- if (run %(long int) 10== 0) {
- fprintf(output, "%d %lf %lf %lf\n", l, Planets[l].x[0], Planets[l].x[1], Planets[l].x[2]);
- }
- }
- for (int k = 0; k < noOfBodies; k++) {
- Planets[k].x[0] += dt * Planets[k].v[0];
- Planets[k].x[1] += dt * Planets[k].v[1];
- Planets[k].x[2] += dt * Planets[k].v[2];
- Planets[k].a[0] = 0;
- Planets[k].a[1] = 0;
- Planets[k].a[2] = 0;
- }
-
- for (int i = 0; i < noOfBodies; i++) {
- for (int j = i+1; j < noOfBodies; j++) {
- double dist = distance(Planets[i].x, Planets[j].x);
- double a = G/pow(dist, 3);
-
- Planets[i].a[0] += a * Planets[j].mass * ( Planets[j].x[0] - Planets[i].x[0] );
- Planets[i].a[1] += a * Planets[j].mass * ( Planets[j].x[1] - Planets[i].x[1] );
- Planets[i].a[2] += a * Planets[j].mass * ( Planets[j].x[2] - Planets[i].x[2] );
- Planets[j].a[0] += a * Planets[i].mass * ( Planets[i].x[0] - Planets[j].x[0] );
- Planets[j].a[1] += a * Planets[i].mass * ( Planets[i].x[1] - Planets[j].x[1] );
- Planets[j].a[2] += a * Planets[i].mass * ( Planets[i].x[2] - Planets[j].x[2] );
- }
- }
-
- for (int m = 0; m < noOfBodies; m++) {
- Planets[m].v[0] += 0.5 * Planets[m].a[0]* dt;
- Planets[m].v[1] += 0.5 * Planets[m].a[1]* dt;
- Planets[m].v[2] += 0.5 * Planets[m].a[2]* dt;
- }
- run++;
- }
- printf("energy after: %.9f\n", energy(noOfBodies));
- printf("total angular momentum after: %.9f\n", an_mom(noOfBodies));
- fclose(output);
-}
-