diff --git a/src/main/c/includes/fluid/sim/pressurecell/solver_consts.h b/src/main/c/includes/fluid/sim/pressurecell/solver_consts.h index 467880b7..86de6de6 100644 --- a/src/main/c/includes/fluid/sim/pressurecell/solver_consts.h +++ b/src/main/c/includes/fluid/sim/pressurecell/solver_consts.h @@ -104,7 +104,7 @@ /** * Enables clamping small density values to 0 */ -#define FLUID_PRESSURECELL_ENABLE_CLAMP_MIN_DENSITY 0 +#define FLUID_PRESSURECELL_ENABLE_CLAMP_MIN_DENSITY 1 /** * Enables renormalizing the velocity field to a max value of 1 during projection diff --git a/src/main/c/src/fluid/sim/pressurecell/pressure.c b/src/main/c/src/fluid/sim/pressurecell/pressure.c index 78cb31e0..59d11a35 100644 --- a/src/main/c/src/fluid/sim/pressurecell/pressure.c +++ b/src/main/c/src/fluid/sim/pressurecell/pressure.c @@ -4,7 +4,10 @@ #include "fluid/sim/pressurecell/solver_consts.h" #include "math/ode/multigrid_parallel.h" - +/** + * Divergence level we warn at + */ +#define DIVERGENCE_WARNING_VALUE 100.0f /** * Approximates the pressure for this chunk @@ -85,7 +88,7 @@ LIBRARY_API void pressurecell_approximate_pressure(Environment * environment, Ch phi0[IX(x,y,z)] = divCache[IX(x,y,z)]; pressureTemp[IX(x,y,z)] = pressureCache[IX(x,y,z)] * FLUID_PRESSURECELL_PRESSURE_BACKDOWN_FACTOR; // pressureTemp[IX(x,y,z)] = 0; - if(divCache[IX(x,y,z)] > 3){ + if(divCache[IX(x,y,z)] > DIVERGENCE_WARNING_VALUE){ printf("invalid divergence!\n"); printf("%f \n", divCache[IX(x,y,z)]); printf("\n"); @@ -294,7 +297,7 @@ LIBRARY_API void pressurecell_approximate_divergence(Environment * environment, // } newDivergence = (du+dv+dw) * (-0.5f * FLUID_PRESSURECELL_SPACING); // divArr[IX(x,y,z)] = divArr[IX(x,y,z)] + newDivergence - FLUID_PRESSURECELL_RESIDUAL_MULTIPLIER * divArr[IX(x,y,z)] + outflowDiv; - if(newDivergence > FLUID_PRESSURECELL_MAX_DIVERGENCE || newDivergence < -FLUID_PRESSURECELL_MAX_DIVERGENCE){ + if(newDivergence > DIVERGENCE_WARNING_VALUE || newDivergence < -DIVERGENCE_WARNING_VALUE){ printf("Invalid divergence! \n"); printf("%f \n",newDivergence); printf("%f %f \n", uArr[IX(x+1,y,z)], uArr[IX(x-1,y,z)]); diff --git a/src/main/c/src/fluid/sim/pressurecell/pressurecell.c b/src/main/c/src/fluid/sim/pressurecell/pressurecell.c index 3378910b..450c518c 100644 --- a/src/main/c/src/fluid/sim/pressurecell/pressurecell.c +++ b/src/main/c/src/fluid/sim/pressurecell/pressurecell.c @@ -143,6 +143,7 @@ LIBRARY_API void fluid_pressurecell_simulate( Chunk * currentChunk = chunks[i]; pressurecell_advect_density(environment,currentChunk); if(FLUID_PRESSURECELL_ENABLE_RECAPTURE){ + //d->dTemp->d fluid_pressurecell_recapture_density(environment,currentChunk); } } @@ -153,6 +154,7 @@ LIBRARY_API void fluid_pressurecell_simulate( // for(int i = 0; i < numChunks; i++){ Chunk * currentChunk = chunks[i]; + //d->d fluid_pressurecell_normalize_chunk(environment,currentChunk); } @@ -171,6 +173,7 @@ LIBRARY_API void fluid_pressurecell_simulate( // for(int i = 0; i < numChunks; i++){ Chunk * currentChunk = chunks[i]; + //uTemp->u pressurecell_copy_for_next_frame(environment,currentChunk); fluid_pressurecell_clearArr(currentChunk->d0[CENTER_LOC]); fluid_pressurecell_clearArr(currentChunk->u0[CENTER_LOC]); diff --git a/src/main/c/src/fluid/sim/pressurecell/velocity.c b/src/main/c/src/fluid/sim/pressurecell/velocity.c index 594fc707..a1c71394 100644 --- a/src/main/c/src/fluid/sim/pressurecell/velocity.c +++ b/src/main/c/src/fluid/sim/pressurecell/velocity.c @@ -6,6 +6,11 @@ #include "math/ode/multigrid_parallel.h" +/** + * Velocity magnitude at which we warn + */ +#define WARNING_VELOCITY_MAGNITUDE 1000 + /** * Adds velocity from the delta buffer to this chunk */ @@ -45,7 +50,7 @@ LIBRARY_API void pressurecell_add_velocity(Environment * environment, Chunk * ch uTemp[IX(x,y,z)] = uArr[IX(x,y,z)] + uSourceArr[IX(x,y,z)]; vTemp[IX(x,y,z)] = vArr[IX(x,y,z)] + vSourceArr[IX(x,y,z)]; wTemp[IX(x,y,z)] = wArr[IX(x,y,z)] + wSourceArr[IX(x,y,z)]; - if(vTemp[IX(x,y,z)] > 5){ + if(vTemp[IX(x,y,z)] > WARNING_VELOCITY_MAGNITUDE){ printf("Invalid add velocity!\n"); printf("%f %f %f \n", vTemp[IX(x,y,z)], vArr[IX(x,y,z)], vSourceArr[IX(x,y,z)] ); printf("\n"); @@ -179,7 +184,7 @@ LIBRARY_API void pressurecell_diffuse_velocity(Environment * environment, Chunk ) ) * a ; - if(vArr[IX(x,y,z)] > 5){ + if(vArr[IX(x,y,z)] > WARNING_VELOCITY_MAGNITUDE){ printf("Invalid diffuse!\n"); printf("%f \n", vArr[IX(x,y,z)]); printf("%f\n", vTemp[IX(x,y,z)]); @@ -242,6 +247,7 @@ LIBRARY_API void pressurecell_advect_velocity(Environment * environment, Chunk * float xp, yp, zp; float s0, s1, t0, t1, u0, u1; float interpolatedU, interpolatedV, interpolatedW; + float vecU, vecV, vecW; float magnitude, maxMagnitude; float interpConst = environment->consts.dt / (FLUID_PRESSURECELL_SPACING * FLUID_PRESSURECELL_SPACING); for(y = 1; y < DIM-1; y++){ @@ -258,10 +264,31 @@ LIBRARY_API void pressurecell_advect_velocity(Environment * environment, Chunk * for(y = 1; y < DIM-1; y++){ for(z = 1; z < DIM-1; z++){ for(x = 1; x < DIM-1; x++){ + //figure how far we're advecting + vecU = uArr[IX(x,y,z)] * interpConst; + vecV = vArr[IX(x,y,z)] * interpConst; + vecW = wArr[IX(x,y,z)] * interpConst; + if(vecU > 0.999f){ + vecU = 0.999f; + } else if(vecU < -0.999f){ + vecU = -0.999f; + } + if(vecV > 0.999f){ + vecV = 0.999f; + } else if(vecV < -0.999f){ + vecV = -0.999f; + } + if(vecW > 0.999f){ + vecW = 0.999f; + } else if(vecW < -0.999f){ + vecW = -0.999f; + } + + //calculate the real (float) position we are at - xp = x - uArr[IX(x,y,z)] * interpConst; - yp = y - vArr[IX(x,y,z)] * interpConst; - zp = z - wArr[IX(x,y,z)] * interpConst; + xp = x - vecU; + yp = y - vecV; + zp = z - vecW; //clamp to border x0 = xp; @@ -380,7 +407,7 @@ LIBRARY_API void pressurecell_advect_velocity(Environment * environment, Chunk * // interpolatedW = interpolatedW / magnitude; // } - if(magnitude > 10){ + if(magnitude > WARNING_VELOCITY_MAGNITUDE){ printf("advect invalid set: %f %f %f \n", uArr[IX(x,y,z)], vArr[IX(x,y,z)], wArr[IX(x,y,z)]); int a = 1; int b = 0; @@ -498,7 +525,7 @@ LIBRARY_API double pressurecell_project_velocity(Environment * environment, Chun maxMagnitude = magnitude; } - if(magnitude != magnitude || magnitude > 10){ + if(magnitude != magnitude || magnitude > WARNING_VELOCITY_MAGNITUDE){ printf("invalid magnitude! %f\n", magnitude); printf("%f %f %f\n", pressureDifferenceX, pressureDifferenceY, pressureDifferenceZ); printf("%f %f %f\n", uArr[IX(x,y,z)], vArr[IX(x,y,z)], wArr[IX(x,y,z)]);