diff --git a/experiment/newton/Makefile b/experiment/newton/Makefile index d43118c11726a0697cbb3bc4a8c9a6e399baed30..adec1dbfd2678a7fce0a8715aee9aab3d2503979 100644 --- a/experiment/newton/Makefile +++ b/experiment/newton/Makefile @@ -1,6 +1,7 @@ C_GCC=gcc C_CLANG=clang -CFLAGS=-Wall -Wextra -std=gnu11 -g -fsanitize=leak -fsanitize=undefined -fsanitize=address +# CFLAGS=-Wall -Wextra -std=gnu11 -g -fsanitize=leak -fsanitize=undefined -fsanitize=address +CFLAGS=-Wall -Wextra -std=gnu11 -O3 LIBS=-lSDL2 -lm diff --git a/experiment/newton/newton.c b/experiment/newton/newton.c index 777b61f4f72c0a9df5686fde1c3f55cc22a69ac6..b766ffe59ff6417154c4a878da06bd32e8d11a07 100644 --- a/experiment/newton/newton.c +++ b/experiment/newton/newton.c @@ -1,95 +1,95 @@ +#include "gfx.h" +#include <math.h> +#include <stdbool.h> #include <stdio.h> #include <stdlib.h> -#include <stdbool.h> -#include <math.h> #include <time.h> -#include "gfx.h" #define GRAV 9.8 #define MASS 0.1 /// @param context graphical context to use. -static void render(struct gfx_context_t *context, double pos) -{ - int x = (int)(pos * 1000); - int y = context->height / 2; - uint32_t color = COLOR_WHITE; // 8-bit per color channel - for (int i = 0; i < 5; ++i) - { - for (int j = 0; j < 5; ++j) - { - gfx_putpixel(context, x + i, y + j, color); - } +static void render(struct gfx_context_t *context, double pos) { + int x = (int)(pos * 1000); + int y = context->height / 2; + uint32_t color = COLOR_WHITE; // 8-bit per color channel + for (int i = 0; i < 5; ++i) { + for (int j = 0; j < 5; ++j) { + gfx_putpixel(context, x + i, y + j, color); } - gfx_present(context); + } + gfx_present(context); } -double verlet(double p, double p_1, double a, double dt) -{ - double pos = 2.0 * p - p_1 + dt * dt * a; - return pos; +double verlet(double p, double p_1, double a, double dt) { + double pos = 2.0 * p - p_1 + dt * dt * a; + return pos; } -int main(int argc, char **argv) -{ - double force = 0.0; - if (argc == 2) - { - force = atof(argv[1]); - } - else - { - printf("Error: missing force parameter.\n\ - Normal use is: ./newton <force>\n\ +int main(int argc, char **argv) { + double force = 0.0; + bool debug = false; + if (argc == 2) { + force = atof(argv[1]); + } else if (argc == 3) { + force = atof(argv[1]); + debug = (0 == strcmp("-d", argv[2])); + } else { + printf("Error: missing force parameter.\n\ + Normal use is: ./newton <force> -d\n\ + where -d is for debug mode and is optional\n\ Example: ./newton 1.0\n"); - return EXIT_FAILURE; - } - int width = 1900, height = 50; - struct gfx_context_t *ctxt = gfx_create("Newton's Law", width, height); - if (!ctxt) - { - fprintf(stderr, "Graphics initialization failed!\n"); - return EXIT_FAILURE; - } + return EXIT_FAILURE; + } + int width = 1900, height = 50; + struct gfx_context_t *ctxt = gfx_create("Newton's Law", width, height); + if (!ctxt) { + fprintf(stderr, "Graphics initialization failed!\n"); + return EXIT_FAILURE; + } - double p = 0.0; - double p_1 = 0.0; - double p_2 = 0.0; - double dt = 0.02; + double p = 0.0; + double p_1 = 0.0; + double p_2 = 0.0; + double dt = 0.02; - struct timespec time = {0, dt * 1000000000}; + struct timespec time = {0, dt * 1e9}; - int num = sqrt(width / 1000.0 / 0.5 / (force / MASS)) / dt; - printf("%d\n", num); + int num = sqrt(width / 1000.0 / 0.5 / (force / MASS)) / dt; + printf("%d\n", num); - gfx_clear(ctxt, COLOR_BLACK); + gfx_clear(ctxt, COLOR_BLACK); - double avg_acc = 0.0; - for (int i = 0; i < num; ++i) - { - double vel = (p - p_1) / dt; - double vel_1 = (p_1 - p_2) / dt; - double acc = (vel - vel_1) / dt; - if (i > 0) - { - printf("t = %f, pos = %f, vel = %f, acc = %f\n", i * dt, p, vel, acc); - avg_acc += acc; - } - render(ctxt, p); - nanosleep(&time, NULL); - double p_tmp = p; - double noise = ((rand() % 1000) - 500) / 100000000.0; - p = verlet(p, p_1, force / MASS, dt) + noise; // add some noise of the order of - p_2 = p_1; - p_1 = p_tmp; + double avg_acc = 0.0; + for (int i = 0; i < num; ++i) { + double vel = (p - p_1) / dt; + double vel_1 = (p_1 - p_2) / dt; + double acc = (vel - vel_1) / dt; + if (i > 0) { + if (debug) { + printf("t = %f, pos = %f, vel = %f, acc = %f\n", i * dt, p, vel, acc); + } else { + printf("t = %f, pos = %f\n", i * dt, p); + } + avg_acc += acc; } + render(ctxt, p); + nanosleep(&time, NULL); + double p_tmp = p; + double noise = ((rand() % 1000) - 500) / 100000000.0; + p = verlet(p, p_1, force / MASS, dt) + + noise; // add some noise of the order of + p_2 = p_1; + p_1 = p_tmp; + } + if (debug) { printf("Simulation ended: avg acceleration = %f\n", avg_acc / (num - 1)); + } - while (gfx_keypressed() != SDLK_ESCAPE) - { - } + while (gfx_keypressed() != SDLK_ESCAPE) { + } - gfx_destroy(ctxt); + gfx_destroy(ctxt); - return EXIT_SUCCESS; + return EXIT_SUCCESS; } \ No newline at end of file