From 506cced6958c10eb35a5b22ee63f9ade3462bd25 Mon Sep 17 00:00:00 2001 From: austin Date: Mon, 9 Dec 2024 19:38:10 -0500 Subject: [PATCH] density advection tests --- .vscode/settings.json | 3 +- src/main/c/includes/fluid/sim/grid2/density.h | 2 +- src/main/c/src/fluid/sim/grid2/density.c | 48 ++- src/main/c/src/fluid/sim/grid2/fluidsim.c | 53 +--- .../fluid/sim/grid2/density_advection_tests.c | 285 ++++++++++++++++++ 5 files changed, 345 insertions(+), 46 deletions(-) create mode 100644 src/test/c/fluid/sim/grid2/density_advection_tests.c diff --git a/.vscode/settings.json b/.vscode/settings.json index 32fb9d47..fdd3769d 100644 --- a/.vscode/settings.json +++ b/.vscode/settings.json @@ -41,6 +41,7 @@ "mathutils.h": "c", "velocity.h": "c", "density.h": "c", - "grid2.h": "c" + "grid2.h": "c", + "math.h": "c" } } \ No newline at end of file diff --git a/src/main/c/includes/fluid/sim/grid2/density.h b/src/main/c/includes/fluid/sim/grid2/density.h index 52b6400a..a5d1bd37 100644 --- a/src/main/c/includes/fluid/sim/grid2/density.h +++ b/src/main/c/includes/fluid/sim/grid2/density.h @@ -45,7 +45,7 @@ LIBRARY_API void fluid_grid2_solveDiffuseDensity( /** * Advects the density based on the vectors */ -void fluid_grid2_advectDensity(float ** d, float ** d0, float ** ur, float ** vr, float ** wr, float dt); +LIBRARY_API void fluid_grid2_advectDensity(float ** d, float ** d0, float ** ur, float ** vr, float ** wr, float dt); diff --git a/src/main/c/src/fluid/sim/grid2/density.c b/src/main/c/src/fluid/sim/grid2/density.c index 531f3e32..a8d3c24a 100644 --- a/src/main/c/src/fluid/sim/grid2/density.c +++ b/src/main/c/src/fluid/sim/grid2/density.c @@ -85,7 +85,7 @@ LIBRARY_API void fluid_grid2_solveDiffuseDensity( /** * Advects the density based on the vectors */ -void fluid_grid2_advectDensity(float ** d, float ** d0, float ** ur, float ** vr, float ** wr, float dt){ +LIBRARY_API void fluid_grid2_advectDensity(float ** d, float ** d0, float ** ur, float ** vr, float ** wr, float dt){ int i, j, k, i0, j0, k0, i1, j1, k1; int m,n,o; float x, y, z, s0, t0, s1, t1, u1, u0, dtx,dty,dtz; @@ -214,6 +214,7 @@ void fluid_grid2_advectDensity(float ** d, float ** d0, float ** ur, float ** vr if(k1 >= DIM){ k1 = DIM - 1; } + // if(k1 < 0){ // k1 = 0; // } @@ -230,6 +231,51 @@ void fluid_grid2_advectDensity(float ** d, float ** d0, float ** ur, float ** vr t0*u1*center_d0[IX(i1,j0,k1)]+ t1*u1*center_d0[IX(i1,j1,k1)] ); + + // if(i == 3 && j == 2 && k == 2){ + // printf("density at <%d,%d,%d> \n",i,j,k); + // printf("sample point precise: %.2f %.2f %.2f\n",x,y,z); + // printf("sample box range: <%d,%d,%d> -> <%d,%d,%d> \n",i0,j0,k0,i1,j1,k1); + // printf("sample values: %.2f %.2f %.2f %.2f \n", + // center_d0[IX(i0,j0,k0)], + // center_d0[IX(i0,j1,k0)], + // center_d0[IX(i0,j0,k1)], + // center_d0[IX(i0,j1,k1)] + // ); + // printf("sample values: %.2f %.2f %.2f %.2f \n", + // center_d0[IX(i1,j0,k0)], + // center_d0[IX(i1,j1,k0)], + // center_d0[IX(i1,j0,k1)], + // center_d0[IX(i1,j1,k1)] + // ); + // //print ints + // // printf("i0: %d\n",i0); + // // printf("j0: %d\n",j0); + // // printf("k0: %d\n",k0); + // // printf("i1: %d\n",i1); + // // printf("j1: %d\n",j1); + // // printf("k1: %d\n",k1); + // // printf("m: %d\n",m); + // // printf("n: %d\n",n); + // // printf("o: %d\n",o); + + // //print floats + // // printf("x: %f\n",x); + // // printf("y: %f\n",y); + // // printf("z: %f\n",z); + + // // printf("t0: %f\n",s0); + // // printf("s0: %f\n",t0); + // // printf("t1: %f\n",s1); + // // printf("s1: %f\n",t1); + // // printf("u0: %f\n",u0); + // // printf("u1: %f\n",u1); + + // // printf("dtx: %f\n",dtx); + // // printf("dty: %f\n",dty); + // // printf("dtz: %f\n",dtz); + // printf("\n"); + // } } } } diff --git a/src/main/c/src/fluid/sim/grid2/fluidsim.c b/src/main/c/src/fluid/sim/grid2/fluidsim.c index b42474b2..22ddab4d 100644 --- a/src/main/c/src/fluid/sim/grid2/fluidsim.c +++ b/src/main/c/src/fluid/sim/grid2/fluidsim.c @@ -74,22 +74,9 @@ void fluid_grid2_simulate( //swap all vector fields //swap vector fields - float * tmpArr; - for(int j = 0; j < 27; j++){ - tmpArr = currentChunk->u[j]; - currentChunk->u[j] = currentChunk->u0[j]; - currentChunk->u0[j] = tmpArr; - } - for(int j = 0; j < 27; j++){ - tmpArr = currentChunk->v[j]; - currentChunk->v[j] = currentChunk->v0[j]; - currentChunk->v0[j] = tmpArr; - } - for(int j = 0; j < 27; j++){ - tmpArr = currentChunk->w[j]; - currentChunk->w[j] = currentChunk->w0[j]; - currentChunk->w0[j] = tmpArr; - } + fluid_grid2_flip_arrays(currentChunk->u,currentChunk->u0); + fluid_grid2_flip_arrays(currentChunk->v,currentChunk->v0); + fluid_grid2_flip_arrays(currentChunk->w,currentChunk->w0); //solve vector diffusion for(int l = 0; l < FLUID_GRID2_VECTOR_DIFFUSE_TIMES; l++){ //solve vector diffusion @@ -128,22 +115,10 @@ void fluid_grid2_simulate( //swap all vector fields //swap vector fields + fluid_grid2_flip_arrays(currentChunk->u,currentChunk->u0); + fluid_grid2_flip_arrays(currentChunk->v,currentChunk->v0); + fluid_grid2_flip_arrays(currentChunk->w,currentChunk->w0); - for(int j = 0; j < 27; j++){ - tmpArr = currentChunk->u[j]; - currentChunk->u[j] = currentChunk->u0[j]; - currentChunk->u0[j] = tmpArr; - } - for(int j = 0; j < 27; j++){ - tmpArr = currentChunk->v[j]; - currentChunk->v[j] = currentChunk->v0[j]; - currentChunk->v0[j] = tmpArr; - } - for(int j = 0; j < 27; j++){ - tmpArr = currentChunk->w[j]; - currentChunk->w[j] = currentChunk->w0[j]; - currentChunk->w0[j] = tmpArr; - } //advect vectors across boundaries //advect fluid_grid2_advectVectors(currentChunk->u,currentChunk->v,currentChunk->w,currentChunk->u0,currentChunk->v0,currentChunk->w0,timestep); @@ -195,15 +170,10 @@ void fluid_grid2_simulate( for(int i = 0; i < numChunks; i++){ Chunk * currentChunk = chunks[i]; fluid_grid2_addDensity(environment,currentChunk->d,currentChunk->d0,timestep); + //swap all density arrays //swap vector fields - - float * tmpArr; - for(int j = 0; j < 27; j++){ - tmpArr = currentChunk->d[j]; - currentChunk->d[j] = currentChunk->d0[j]; - currentChunk->d0[j] = tmpArr; - } + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); //diffuse density for(int l = 0; l < FLUID_GRID2_LINEARSOLVERTIMES; l++){ fluid_grid2_solveDiffuseDensity(currentChunk->d,currentChunk->d0,FLUID_GRID2_DIFFUSION_CONSTANT,timestep); @@ -211,11 +181,8 @@ void fluid_grid2_simulate( } //swap all density arrays //swap vector fields - for(int j = 0; j < 27; j++){ - tmpArr = currentChunk->d[j]; - currentChunk->d[j] = currentChunk->d0[j]; - currentChunk->d0[j] = tmpArr; - } + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + //advect density fluid_grid2_advectDensity(currentChunk->d,currentChunk->d0,currentChunk->u,currentChunk->v,currentChunk->w,timestep); } diff --git a/src/test/c/fluid/sim/grid2/density_advection_tests.c b/src/test/c/fluid/sim/grid2/density_advection_tests.c new file mode 100644 index 00000000..475a0aa1 --- /dev/null +++ b/src/test/c/fluid/sim/grid2/density_advection_tests.c @@ -0,0 +1,285 @@ +#include + +#include "stb/stb_ds.h" + +#include "fluid/queue/boundsolver.h" +#include "fluid/queue/chunkmask.h" +#include "fluid/queue/chunk.h" +#include "fluid/env/environment.h" +#include "fluid/env/utilities.h" +#include "fluid/sim/grid2/density.h" +#include "fluid/sim/grid2/solver_consts.h" +#include "fluid/sim/grid2/utilities.h" +#include "../../../util/chunk_test_utils.h" +#include "../../../util/test.h" + + +/** + * Amount to place in density advection tests + */ +#define FLUID_GRID2_DENSITY_ADVECTION_TESTS_PLACEMENT_VAL 1.0f + +/** + * Center of the advection cell + */ +#define FLUID_GRID2_DENSITY_ADVECTION_CELL_CENTER 24 + +/** + * Creates a convection cell for testing advection + */ +void fluid_sim_grid2_density_advection_setup_convection_cell(Chunk ** queue){ + int chunkCount = arrlen(queue); + int realX, realY, realZ; + int worldX, worldY, worldZ; + for(int chunkIndex = 0; chunkIndex < chunkCount; chunkIndex++){ + Chunk * chunk = queue[chunkIndex]; + worldX = chunk->x; + worldY = chunk->y; + worldZ = chunk->z; + for(int x = 0; x < DIM; x++){ + for(int y = 0; y < DIM; y++){ + for(int z = 0; z < DIM; z++){ + double angle = ((x - FLUID_GRID2_DENSITY_ADVECTION_CELL_CENTER),(y - FLUID_GRID2_DENSITY_ADVECTION_CELL_CENTER)); + if(x == 0){ + realX = DIM-2 + (CHUNK_SPACING * (worldX - 1)); + + } else if(x == DIM-1){ + realX = 1 + (CHUNK_SPACING * (worldX + 1)); + + } else { + realX = x + (CHUNK_SPACING * worldX); + } + if(y == 0){ + realY = DIM-2 + (CHUNK_SPACING * (worldY - 1)); + + } else if(y == DIM-1){ + realY = 1 + (CHUNK_SPACING * (worldY + 1)); + + } else { + realY = y + (CHUNK_SPACING * worldY); + } + if(z == 0){ + realZ = DIM-2 + (CHUNK_SPACING * (worldZ - 1)); + + } else if(z == DIM-1){ + realZ = 1 + (CHUNK_SPACING * (worldZ + 1)); + + } else { + realZ = z + (CHUNK_SPACING * worldZ); + } + chunk->u[CENTER_LOC][IX(x,y,z)] = (float)sin(angle); + chunk->v[CENTER_LOC][IX(x,y,z)] = (float)cos(angle); + chunk->w[CENTER_LOC][IX(x,y,z)] = 0; + } + } + } + } +} + +/** + * Testing density diffusion + */ +int fluid_sim_grid2_density_advection_test1(){ + printf("fluid_sim_grid2_density_advection_test1\n"); + int rVal = 0; + Environment * env = fluid_environment_create(); + Chunk ** queue = NULL; + queue = createChunkGrid(env,3,3,3); + + + + //setup chunk values + Chunk * currentChunk = queue[0]; + currentChunk->d[CENTER_LOC][IX(2,2,2)] = MAX_FLUID_VALUE; + float beforeSum = chunk_queue_sum_density(queue); + chunk_fill_real(queue[0]->u[CENTER_LOC],FLUID_GRID2_DENSITY_ADVECTION_TESTS_PLACEMENT_VAL); + + //actually advect + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + fluid_grid2_advectDensity(currentChunk->d,currentChunk->d0,currentChunk->u,currentChunk->v,currentChunk->w,FLUID_GRID2_SIM_STEP); + + //sum the result + float afterSum = chunk_queue_sum_density(queue); + + if(fabs(beforeSum - afterSum) > FLUID_GRID2_REALLY_SMALL_VALUE){ + rVal += assertEqualsFloat(beforeSum,afterSum,"Density advection step changed density sum! %f %f \n"); + } + + return rVal; +} + +/** + * Testing density diffusion + */ +int fluid_sim_grid2_density_advection_test2(){ + printf("fluid_sim_grid2_density_advection_test2\n"); + int rVal = 0; + Environment * env = fluid_environment_create(); + Chunk ** queue = NULL; + queue = createChunkGrid(env,3,3,3); + + + + //setup chunk values + Chunk * currentChunk = queue[0]; + currentChunk->d[CENTER_LOC][IX(2,2,2)] = MAX_FLUID_VALUE; + float beforeSum = chunk_queue_sum_density(queue); + fluid_sim_grid2_density_advection_setup_convection_cell(queue); + + //actually simulate + int frameCount = 50; + for(int frame = 0; frame < frameCount; frame++){ + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + fluid_grid2_advectDensity(currentChunk->d,currentChunk->d0,currentChunk->u,currentChunk->v,currentChunk->w,FLUID_GRID2_SIM_STEP); + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + } + + //test the result + float afterSum = chunk_queue_sum_density(queue); + if(fabs(beforeSum - afterSum) > FLUID_GRID2_REALLY_SMALL_VALUE){ + rVal += assertEqualsFloat(beforeSum,afterSum,"Density advection step changed density sum! %f %f \n"); + } + + return rVal; +} + +/** + * Testing density diffusion + */ +int fluid_sim_grid2_density_advection_test3(){ + printf("fluid_sim_grid2_density_advection_test3\n"); + int rVal = 0; + Environment * env = fluid_environment_create(); + Chunk ** queue = NULL; + queue = createChunkGrid(env,3,3,3); + + + + + //setup chunk values + Chunk * currentChunk = queue[0]; + currentChunk->d[CENTER_LOC][IX(2,2,2)] = MAX_FLUID_VALUE; + float beforeSum = chunk_queue_sum_density(queue); + fluid_sim_grid2_density_advection_setup_convection_cell(queue); + + //actually simulate + int frameCount = 400; + for(int frame = 0; frame < frameCount; frame++){ + int chunkCount = arrlen(queue); + for(int chunkIndex = 0; chunkIndex < 1; chunkIndex++){ + currentChunk = queue[chunkIndex]; + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + fluid_grid2_advectDensity(currentChunk->d,currentChunk->d0,currentChunk->u,currentChunk->v,currentChunk->w,FLUID_GRID2_SIM_STEP); + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + } + } + + //test the result + float afterSum = chunk_queue_sum_density(queue); + if(fabs(beforeSum - afterSum) > FLUID_GRID2_REALLY_SMALL_VALUE){ + rVal += assertEqualsFloat(beforeSum,afterSum,"Density advection step changed density sum! %f %f \n"); + } + + return rVal; +} + +/** + * Testing density diffusion + */ +int fluid_sim_grid2_density_advection_test4(){ + printf("fluid_sim_grid2_density_advection_test4\n"); + int rVal = 0; + Environment * env = fluid_environment_create(); + Chunk ** queue = NULL; + queue = createChunkGrid(env,3,3,3); + + + + + //setup chunk values + Chunk * currentChunk = queue[0]; + currentChunk->d[CENTER_LOC][IX(2,2,2)] = MAX_FLUID_VALUE; + currentChunk->d[CENTER_LOC][IX(2,4,2)] = MAX_FLUID_VALUE; + currentChunk->d[CENTER_LOC][IX(2,2,7)] = MAX_FLUID_VALUE; + currentChunk->d[CENTER_LOC][IX(12,2,2)] = MAX_FLUID_VALUE; + float beforeSum = chunk_queue_sum_density(queue); + fluid_sim_grid2_density_advection_setup_convection_cell(queue); + + //actually simulate + int frameCount = 400; + for(int frame = 0; frame < frameCount; frame++){ + int chunkCount = arrlen(queue); + for(int chunkIndex = 0; chunkIndex < 1; chunkIndex++){ + currentChunk = queue[chunkIndex]; + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + fluid_grid2_advectDensity(currentChunk->d,currentChunk->d0,currentChunk->u,currentChunk->v,currentChunk->w,FLUID_GRID2_SIM_STEP); + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + } + } + + //test the result + float afterSum = chunk_queue_sum_density(queue); + if(fabs(beforeSum - afterSum) > FLUID_GRID2_REALLY_SMALL_VALUE){ + rVal += assertEqualsFloat(beforeSum,afterSum,"Density advection step changed density sum! %f %f \n"); + } + + return rVal; +} + +/** + * Testing density diffusion + */ +int fluid_sim_grid2_density_advection_test5(){ + printf("fluid_sim_grid2_density_advection_test4\n"); + int rVal = 0; + Environment * env = fluid_environment_create(); + Chunk ** queue = NULL; + queue = createChunkGrid(env,3,3,3); + + + + + //setup chunk values + Chunk * currentChunk = queue[3 * 3 + 3]; + chunk_fill_real(currentChunk->d[CENTER_LOC],MAX_FLUID_VALUE); + float beforeSum = chunk_queue_sum_density(queue); + fluid_sim_grid2_density_advection_setup_convection_cell(queue); + + //actually simulate + int frameCount = 400; + for(int frame = 0; frame < frameCount; frame++){ + int chunkCount = arrlen(queue); + for(int chunkIndex = 0; chunkIndex < 1; chunkIndex++){ + currentChunk = queue[chunkIndex]; + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + fluid_grid2_advectDensity(currentChunk->d,currentChunk->d0,currentChunk->u,currentChunk->v,currentChunk->w,FLUID_GRID2_SIM_STEP); + fluid_grid2_flip_arrays(currentChunk->d,currentChunk->d0); + } + } + + //test the result + float afterSum = chunk_queue_sum_density(queue); + if(fabs(beforeSum - afterSum) > FLUID_GRID2_REALLY_SMALL_VALUE){ + rVal += assertEqualsFloat(beforeSum,afterSum,"Density advection step changed density sum! %f %f \n"); + } + + return rVal; +} + + + + +/** + * Testing density diffusion + */ +int fluid_sim_grid2_density_advection_tests(int argc, char **argv){ + int rVal = 0; + + rVal += fluid_sim_grid2_density_advection_test1(); + rVal += fluid_sim_grid2_density_advection_test2(); + rVal += fluid_sim_grid2_density_advection_test3(); + rVal += fluid_sim_grid2_density_advection_test4(); + rVal += fluid_sim_grid2_density_advection_test5(); + + return rVal; +} \ No newline at end of file