From 2ecf720f45933ac56a3520b238c14fb351ee74cf Mon Sep 17 00:00:00 2001 From: Cory Bloor Date: Wed, 26 Feb 2020 23:13:26 -0700 Subject: [PATCH] Extract acceleration function --- main.c | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/main.c b/main.c index 739e0a2..9605d90 100644 --- a/main.c +++ b/main.c @@ -701,6 +701,11 @@ void fmadd(double a[Y][X], double b, double c[Y][X]) { } } +// Acceleration due to pressure: -grad(p) / density +float accel(float delta_p) { + return -invf(k_density * k_side_length) * delta_p; +} + void project(float dt, float u[Y][X], float v[Y][X], float uout[Y][X], float vout[Y][X]) { // Using -div(u) = laplacian(p) * dt / density. // Solve Ap = b, with A = laplacian * dx², @@ -779,7 +784,7 @@ void project(float dt, float u[Y][X], float v[Y][X], float uout[Y][X], float vou if (u_property(g_solid, v2i(x,y))) { uout[y][x] = 0.f; } else if (u_property(g_fluid, v2i(x,y))) { - uout[y][x] = u[y][x] - invf(k_density*k_side_length) * dt * (p[y][x+1] - p[y][x]); + uout[y][x] = u[y][x] + accel(p[y][x+1] - p[y][x]) * dt; } else { // air uout[y][x] = 0.f; } @@ -792,7 +797,7 @@ void project(float dt, float u[Y][X], float v[Y][X], float uout[Y][X], float vou if (v_property(g_solid, v2i(x,y))) { vout[y][x] = 0.f; } else if (v_property(g_fluid, v2i(x,y))) { - vout[y][x] = v[y][x] - invf(k_density*k_side_length) * dt * (p[y+1][x] - p[y][x]); + vout[y][x] = v[y][x] + accel(p[y+1][x] - p[y][x]) * dt; } else { // air vout[y][x] = 0.f; }