Ended magnetic field line simulation

charge/charge.c

@@ -9,20 +9,94 @@
  * 
  */
 #include "../draw/draw.h"
-#include "../vector/vector.h"
+#include "../vec2/vec2.h"
 #include "../utilities/utils.h"
 #include "charge.h"
 
 bool compute_e(charge_t c, vec2 p, double eps, vec2 *e)
 {
+    vec2 qP = vec2_sub(p, c.pos);
+    double norm = vec2_norm(qP);
+    double bigE = K * c.q / pow(norm, 2);    // E = k * (Q / r^2)
+    *e = vec2_mul(vec2_div(qP, norm), bigE); // e = E * qP / ||qP||
+    return (norm > eps);
 }
 
 bool compute_total_normalized_e(charge_t *charges, int num_charges, vec2 p, double eps, vec2 *e)
 {
+    vec2 sum;
+    vec2 tmp_e;
+
+    for (int i = 0; i < num_charges; i++)
+    {
+        if (!compute_e(charges[i], p, eps, &tmp_e))
+            return false;
+        sum = vec2_add(sum, tmp_e);
+    }
+
+    *e = sum;
+    return true;
+}
+
+bool compute_p_next(charge_t *charges, int num_charges, double eps, double delta, vec2 *p, bool opposed)
+{
+    vec2 v;
+    bool can_draw = compute_total_normalized_e(charges, num_charges, *p, eps, &v);
+
+    v = vec2_div(vec2_mul(v, delta), vec2_norm(v));
+
+    if (opposed)
+        *p = vec2_sub(*p, v); // P + deltax * E / ||E||
+    else
+        *p = vec2_add(*p, v); // P - deltax * E / ||E||
+
+    return can_draw;
+}
+
+bool pos_contrain_in_universe(vec2 pos, double x0, double x1, double y0, double y1)
+{
+    return (pos.x <= x1 && pos.x >= x0 && pos.y <= y1 && pos.y >= y0);
+}
+
+bool pos_not_too_close(vec2 p, charge_t *c, int num_charges)
+{
+    for (int i = 0; i < num_charges; i++)
+    {
+        double norm = vec2_norm(vec2_sub(p, c[i].pos));
+        if (norm < EPS)
+            return false;
+    }
+    return true;
+}
+
+bool field_line_segment(struct gfx_context_t *ctxt, charge_t *charges, int num_charges, double dx, double x0, double x1, double y0, double y1, vec2 old, vec2 p, bool side)
+{
+    while (pos_contrain_in_universe(p, x0, x1, y0, y1) && pos_not_too_close(p, charges, num_charges))
+    {
+        coordinates_t old_coord = position_to_coordinates(ctxt->width, ctxt->height, x0, x1, y0, y1, old);
+        coordinates_t pi_coord = position_to_coordinates(ctxt->width, ctxt->height, x0, x1, y0, y1, p);
+
+        gfx_draw_line(ctxt, pi_coord, old_coord, 0xFFFFFF);
+
+        old = p; // Save current point
+        compute_p_next(charges, num_charges, EPS, dx, &p, side);
+    }
+
+    return true;
 }
 
-bool draw_field_line(struct gfx_context_t *ctxt, charge_t *charges, int num_charges, double dx, vec2 pos0, double x0, double x1, double y0, double y1)
+void draw_field_line(struct gfx_context_t *ctxt, charge_t *charges, int num_charges, double dx, vec2 pos0, double x0, double x1, double y0, double y1)
 {
+    vec2 p = pos0;
+    vec2 old = p;
+
+    bool ended = false;
+
+    for (int side = 0; side <= 1; side++)
+    {
+        old = pos0; // Reset to start from the same point for both sides
+        field_line_segment(ctxt, charges, num_charges, dx, x0, x1, y0, y1, old, p, (bool)side);
+    }
 }
 
 void draw_charges(struct gfx_context_t *ctxt, charge_t *charges, int num_charges, double x0, double x1, double y0, double y1)
@@ -34,8 +108,7 @@ void draw_charges(struct gfx_context_t *ctxt, charge_t *charges, int num_charges
     coordinates_t ch0;
     coordinates_t ch1;
 
-    int radius = 25;
-    int chargePadding = 10;
+    int chargePadding = CHARGE_RADIUS / 2;
 
     // Draw charges
     for (int i = 0; i < num_charges; i++)
@@ -53,7 +126,7 @@ void draw_charges(struct gfx_context_t *ctxt, charge_t *charges, int num_charges
             gfx_draw_line(ctxt, coordinates_create(chCoord.row - chargePadding, chCoord.column), coordinates_create(chCoord.row + chargePadding, chCoord.column), COLOR_RED);
         }
 
-        gfx_draw_circle(ctxt, chCoord, radius, charges[i].q < 0 ? COLOR_BLUE : COLOR_RED);
+        gfx_draw_circle(ctxt, chCoord, CHARGE_RADIUS, charges[i].q < 0 ? COLOR_BLUE : COLOR_RED);
     }
 }
 
@@ -67,6 +140,6 @@ void generate_points(vec2 points[], int nb_points)
         x = rand_one();
         y = rand_one();
 
-        points[i] = init_vector(x, y);
+        points[i] = vec2_create(x, y);
     }
 }
\ No newline at end of file

charge/charge.h

@@ -11,10 +11,15 @@
 #ifndef _CHARGE_H_
 #define _CHARGE_H_
 
+#define K 8.988e9
+#define E 1.602e-19
+#define CHARGE_RADIUS 10
+#define EPS 0.025
+
 #include <stdlib.h>
 #include <stdbool.h>
 #include "../utilities/utils.h"
-#include "../vector/vector.h"
+#include "../vec2/vec2.h"
 #include "../gfx/gfx.h"
 
 /**
@@ -42,23 +47,52 @@ bool compute_total_normalized_e(charge_t *charges, int num_charges, vec2 p,
                                 double eps, vec2 *e);
 
 /**
- * @brief Compute and then draw all the points belonging to a field line,
- *        starting from pos0
+ * @brief Compute the next particle
  * 
- * @param ctxt Context to draw on
  * @param charges All the charges
  * @param num_charges Number of charges in action
- * @param dx Space between computed points on the field line
- * @param pos0 First post to start the computing of the field line
- * @param x0 Minal x of the universe
- * @param x1 Maximal x of the universe
- * @param y0 Minal y of the universe
- * @param y1 Maximal y of the universe
- * @return bool false if pos0 not a valid pos (e.g. too close to a charge)
+ * @param eps Minimal gap between field and charge
+ * @param delta Distance between two points in the field line
+ * @param p Current point
+ * @param opposed Is the next point towards charge
+ * @return true 
+ * @return false 
+ */
+bool compute_p_next(charge_t *charges, int num_charges, double eps, double delta, vec2 *p, bool opposed);
+
+/**
+ * @brief Draw a segment of a field line between two given points
+ * 
+ * @param ctxt Context to draw on
+ * @param charges All the charges
+ * @param num_charges Number of charges
+ * @param dx Distance between next point and current point
+ * @param x0 Minimal X
+ * @param x1 Maximal X
+ * @param y0 Minimal Y
+ * @param y1 Maximal Y
+ * @param old Old pos
+ * @param p Next pos
+ * @param side Side of the line
+ * @return true 
+ * @return false 
+ */
+bool field_line_segment(struct gfx_context_t *ctxt, charge_t *charges, int num_charges, double dx, double x0, double x1, double y0, double y1, vec2 old, vec2 p, bool side);
+
+/**
+ * @brief Draw a field line
+ * 
+ * @param ctxt Context to draw on
+ * @param charges All the charges in action
+ * @param num_charges Number of charges
+ * @param dx Distance between  current and next point to draw line
+ * @param pos0 Line starting point (will draw on both side of the point)
+ * @param x0 Minimal X
+ * @param x1 Maximal X
+ * @param y0 Minimal Y
+ * @param y1 Maximal Y
  */
-bool draw_field_line(struct gfx_context_t *ctxt, charge_t *charges,
-                     int num_charges, double dx, vec2 pos0, double x0,
-                     double x1, double y0, double y1);
+void draw_field_line(struct gfx_context_t *ctxt, charge_t *charges, int num_charges, double dx, vec2 pos0, double x0, double x1, double y0, double y1);
 
 /**
  * @brief Draw all charges

main.c

@@ -21,8 +21,9 @@
 #define x1 1 // Maximal x of the universe
 #define y0 0 // Minimal y of the universe
 #define y1 1 // Maximal y of the universe
-#define NB_CHARGES 2
-#define NB_POINTS 100
+#define NB_CHARGES 3
+#define NB_POINTS 200
+#define DELTA (1 / sqrt(pow(WINDOW_WIDTH, 2) + pow(WINDOW_HEIGHT, 2)))
 
 int main(void)
 {
@@ -36,17 +37,21 @@ int main(void)
         return EXIT_FAILURE;
     }
 
+    // Generate all the charges
     charge_t charges[NB_CHARGES];
-    charges[0] = charge_create(10, init_vector(0.25, 0.5));
-    charges[1] = charge_create(-10, init_vector(0.75, 0.5));
+    for (int i = 0; i < NB_CHARGES; i++)
+    {
+        charges[i] = charge_create(rand_one() / E * (rand_one() > 0.5 ? 1 : -1), vec2_create(rand_one(), rand_one()));
+    }
     draw_charges(ctxt, charges, NB_CHARGES, x0, x1, y0, y1);
 
+    // Generate all the points
     vec2 points[NB_POINTS];
     generate_points(points, NB_POINTS);
     for (int i = 0; i < NB_POINTS; i++)
     {
         coordinates_t c = position_to_coordinates(WINDOW_WIDTH, WINDOW_HEIGHT, x0, x1, y0, y1, points[i]);
-        gfx_putpixel(ctxt, c.column, c.row, COLOR_WHITE);
+        draw_field_line(ctxt, charges, NB_CHARGES, DELTA, points[i], x0, x1, y0, y1);
     }
 
     // GFX Draw loop

utilities/utils.c

 #include <math.h>
 #include <stdlib.h>
-#include "../vector/vector.h"
+#include "../vec2/vec2.h"
 #include "utils.h"
 
 // Transform a position in the univers [x0,y0]x[x1,y1] to a screen position

utilities/utils.h

 #ifndef _UTILS_H_
 #define _UTILS_H_
 
-#define K 8.988e9;
-#define E 1.602e-19;
-
 #include <stdint.h>
-#include "../vector/vector.h"
+#include "../vec2/vec2.h"
 #include "../draw/draw.h"
 
 typedef struct _charge_t