diff --git a/src/main/c/src/fluid/queue/boundsolver.c b/src/main/c/src/fluid/queue/boundsolver.c index c19ce574..47202b95 100644 --- a/src/main/c/src/fluid/queue/boundsolver.c +++ b/src/main/c/src/fluid/queue/boundsolver.c @@ -10,19 +10,22 @@ /** * Define as 1 to source values from surrounding chunks */ -#define USE_BOUNDS 0 +#define USE_BOUNDS 1 -static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ +static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal, float * adjacencyMatrix){ int i, j; int neighborIndex; + float adjacentVal; //x+ face neighborIndex = CK(2,1,1); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / ((DIM-2)*(DIM-2)); for(i = 1; i < DIM-1; i++){ for(j = 1; j < DIM-1; j++){ - arrays[CENTER_LOC][IX(DIM-1, i, j)] = arrays[neighborIndex][IX( 1, i, j)]; + arrays[CENTER_LOC][IX(DIM-1, i, j)] = arrays[neighborIndex][IX( 1, i, j)] + adjacentVal; } } } else { @@ -36,9 +39,11 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x- face neighborIndex = CK(0,1,1); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / ((DIM-2)*(DIM-2)); for(i = 1; i < DIM-1; i++){ for(j = 1; j < DIM-1; j++){ - arrays[CENTER_LOC][IX(0, i, j)] = arrays[neighborIndex][IX(DIM-2, i, j)]; + arrays[CENTER_LOC][IX(0, i, j)] = arrays[neighborIndex][IX(DIM-2, i, j)] + adjacentVal; } } } else { @@ -52,9 +57,11 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //y+ face neighborIndex = CK(1,2,1); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / ((DIM-2)*(DIM-2)); for(i = 1; i < DIM-1; i++){ for(j = 1; j < DIM-1; j++){ - arrays[CENTER_LOC][IX( i, DIM-1, j)] = arrays[neighborIndex][IX( i, 1, j)]; + arrays[CENTER_LOC][IX( i, DIM-1, j)] = arrays[neighborIndex][IX( i, 1, j)] + adjacentVal; } } } else { @@ -68,9 +75,11 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //y- face neighborIndex = CK(1,0,1); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / ((DIM-2)*(DIM-2)); for(i = 1; i < DIM-1; i++){ for(j = 1; j < DIM-1; j++){ - arrays[CENTER_LOC][IX(i, 0, j)] = arrays[neighborIndex][IX( i, DIM-2, j)]; + arrays[CENTER_LOC][IX(i, 0, j)] = arrays[neighborIndex][IX( i, DIM-2, j)] + adjacentVal; } } } else { @@ -84,9 +93,11 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //z+ face neighborIndex = CK(1,1,2); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / ((DIM-2)*(DIM-2)); for(i = 1; i < DIM-1; i++){ for(j = 1; j < DIM-1; j++){ - arrays[CENTER_LOC][IX( i, j, DIM-1)] = arrays[neighborIndex][IX( i, j, 1)]; + arrays[CENTER_LOC][IX( i, j, DIM-1)] = arrays[neighborIndex][IX( i, j, 1)] + adjacentVal; } } } else { @@ -100,9 +111,11 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //z- face neighborIndex = CK(1,1,0); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / ((DIM-2)*(DIM-2)); for(i = 1; i < DIM-1; i++){ for(j = 1; j < DIM-1; j++){ - arrays[CENTER_LOC][IX(i, j, 0)] = arrays[neighborIndex][IX( i, j, DIM-2)]; + arrays[CENTER_LOC][IX(i, j, 0)] = arrays[neighborIndex][IX( i, j, DIM-2)] + adjacentVal; } } } else { @@ -121,8 +134,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x+ y+ edge neighborIndex = CK(2,2,1); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX(DIM-1, DIM-1, i)] = arrays[neighborIndex][IX( 1, 1, i)]; + arrays[CENTER_LOC][IX(DIM-1, DIM-1, i)] = arrays[neighborIndex][IX( 1, 1, i)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -133,8 +148,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x+ y- edge neighborIndex = CK(2,0,1); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX(DIM-1, 0, i)] = arrays[neighborIndex][IX( 1, DIM-2, i)]; + arrays[CENTER_LOC][IX(DIM-1, 0, i)] = arrays[neighborIndex][IX( 1, DIM-2, i)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -145,8 +162,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x- y+ edge neighborIndex = CK(0,2,1); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX( 0, DIM-1, i)] = arrays[neighborIndex][IX( DIM-2, 1, i)]; + arrays[CENTER_LOC][IX( 0, DIM-1, i)] = arrays[neighborIndex][IX( DIM-2, 1, i)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -157,8 +176,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x- y- edge neighborIndex = CK(0,0,1); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX( 0, 0, i)] = arrays[neighborIndex][IX( DIM-2, DIM-2, i)]; + arrays[CENTER_LOC][IX( 0, 0, i)] = arrays[neighborIndex][IX( DIM-2, DIM-2, i)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -177,8 +198,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x+ z+ edge neighborIndex = CK(2,1,2); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX(DIM-1, i, DIM-1)] = arrays[neighborIndex][IX( 1, i, 1)]; + arrays[CENTER_LOC][IX(DIM-1, i, DIM-1)] = arrays[neighborIndex][IX( 1, i, 1)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -189,8 +212,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x+ z- edge neighborIndex = CK(2,1,0); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX(DIM-1, i, 0)] = arrays[neighborIndex][IX( 1, i, DIM-2)]; + arrays[CENTER_LOC][IX(DIM-1, i, 0)] = arrays[neighborIndex][IX( 1, i, DIM-2)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -201,8 +226,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x- z+ edge neighborIndex = CK(0,1,2); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX( 0, i, DIM-1)] = arrays[neighborIndex][IX( DIM-2, i, 1)]; + arrays[CENTER_LOC][IX( 0, i, DIM-1)] = arrays[neighborIndex][IX( DIM-2, i, 1)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -213,8 +240,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x- z- edge neighborIndex = CK(0,1,0); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX( 0, i, 0)] = arrays[neighborIndex][IX( DIM-2, i, DIM-2)]; + arrays[CENTER_LOC][IX( 0, i, 0)] = arrays[neighborIndex][IX( DIM-2, i, DIM-2)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -232,8 +261,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //y+ z+ edge neighborIndex = CK(1,2,2); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX( i, DIM-1, DIM-1)] = arrays[neighborIndex][IX( i, 1, 1)]; + arrays[CENTER_LOC][IX( i, DIM-1, DIM-1)] = arrays[neighborIndex][IX( i, 1, 1)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -244,8 +275,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //y+ z- edge neighborIndex = CK(1,2,0); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX( i,DIM-1, 0)] = arrays[neighborIndex][IX( i, 1, DIM-2)]; + arrays[CENTER_LOC][IX( i,DIM-1, 0)] = arrays[neighborIndex][IX( i, 1, DIM-2)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -256,8 +289,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //y- z+ edge neighborIndex = CK(1,0,2); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX( i, 0, DIM-1)] = arrays[neighborIndex][IX( i, DIM-2, 1)]; + arrays[CENTER_LOC][IX( i, 0, DIM-1)] = arrays[neighborIndex][IX( i, DIM-2, 1)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -268,8 +303,10 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //y- z- edge neighborIndex = CK(1,0,0); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + adjacentVal = adjacentVal / (DIM-2); for(i = 1; i < DIM-1; i++){ - arrays[CENTER_LOC][IX( i, 0, 0)] = arrays[neighborIndex][IX( i, DIM-2, DIM-2)]; + arrays[CENTER_LOC][IX( i, 0, 0)] = arrays[neighborIndex][IX( i, DIM-2, DIM-2)] + adjacentVal; } } else { for(i = 1; i < DIM; i++){ @@ -285,7 +322,8 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x+ y+ z+ corner neighborIndex = CK(2,2,2); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ - arrays[CENTER_LOC][IX(DIM-1, DIM-1, DIM-1)] = arrays[neighborIndex][IX( 1, 1, 1)]; + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + arrays[CENTER_LOC][IX(DIM-1, DIM-1, DIM-1)] = arrays[neighborIndex][IX( 1, 1, 1)] + adjacentVal; } else { arrays[CENTER_LOC][IX(DIM-1, DIM-1, DIM-1)] = fillVal; } @@ -293,7 +331,8 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x+ y+ z- corner neighborIndex = CK(2,2,0); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ - arrays[CENTER_LOC][IX(DIM-1, DIM-1, 0)] = arrays[neighborIndex][IX( 1, 1, DIM-2)]; + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + arrays[CENTER_LOC][IX(DIM-1, DIM-1, 0)] = arrays[neighborIndex][IX( 1, 1, DIM-2)] + adjacentVal; } else { arrays[CENTER_LOC][IX(DIM-1, DIM-1, 0)] = fillVal; } @@ -303,7 +342,8 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x+ y- z+ corner neighborIndex = CK(2,0,2); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ - arrays[CENTER_LOC][IX(DIM-1, 0, DIM-1)] = arrays[neighborIndex][IX( 1, DIM-2, 1)]; + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + arrays[CENTER_LOC][IX(DIM-1, 0, DIM-1)] = arrays[neighborIndex][IX( 1, DIM-2, 1)] + adjacentVal; } else { arrays[CENTER_LOC][IX(DIM-1, 0, DIM-1)] = fillVal; } @@ -311,7 +351,8 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x+ y- z- corner neighborIndex = CK(2,0,0); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ - arrays[CENTER_LOC][IX(DIM-1, 0, 0)] = arrays[neighborIndex][IX( 1, DIM-2, DIM-2)]; + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + arrays[CENTER_LOC][IX(DIM-1, 0, 0)] = arrays[neighborIndex][IX( 1, DIM-2, DIM-2)] + adjacentVal; } else { arrays[CENTER_LOC][IX(DIM-1, 0, 0)] = fillVal; } @@ -321,7 +362,8 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x- y+ z+ corner neighborIndex = CK(0,2,2); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ - arrays[CENTER_LOC][IX(0, DIM-1, DIM-1)] = arrays[neighborIndex][IX( DIM-2, 1, 1)]; + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + arrays[CENTER_LOC][IX(0, DIM-1, DIM-1)] = arrays[neighborIndex][IX( DIM-2, 1, 1)] + adjacentVal; } else { arrays[CENTER_LOC][IX(0, DIM-1, DIM-1)] = fillVal; } @@ -329,7 +371,8 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x- y+ z- corner neighborIndex = CK(0,2,0); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ - arrays[CENTER_LOC][IX(0, DIM-1, 0)] = arrays[neighborIndex][IX( DIM-2, 1, DIM-2)]; + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + arrays[CENTER_LOC][IX(0, DIM-1, 0)] = arrays[neighborIndex][IX( DIM-2, 1, DIM-2)] + adjacentVal; } else { arrays[CENTER_LOC][IX(0, DIM-1, 0)] = fillVal; } @@ -339,7 +382,8 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x- y- z+ corner neighborIndex = CK(0,0,2); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ - arrays[CENTER_LOC][IX(0, 0, DIM-1)] = arrays[neighborIndex][IX( DIM-2, DIM-2, 1)]; + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + arrays[CENTER_LOC][IX(0, 0, DIM-1)] = arrays[neighborIndex][IX( DIM-2, DIM-2, 1)] + adjacentVal; } else { arrays[CENTER_LOC][IX(0, 0, DIM-1)] = fillVal; } @@ -347,7 +391,8 @@ static inline void fluid_solve_bounds_checker(float ** arrays, float fillVal){ //x- y- z- corner neighborIndex = CK(0,0,0); if(USE_BOUNDS && arrays[neighborIndex] != NULL){ - arrays[CENTER_LOC][IX(0, 0, 0)] = arrays[neighborIndex][IX( DIM-2, DIM-2, DIM-2)]; + adjacentVal = adjacencyMatrix == NULL ? 0 : adjacencyMatrix[neighborIndex]; + arrays[CENTER_LOC][IX(0, 0, 0)] = arrays[neighborIndex][IX( DIM-2, DIM-2, DIM-2)] + adjacentVal; } else { arrays[CENTER_LOC][IX(0, 0, 0)] = fillVal; } @@ -366,17 +411,17 @@ LIBRARY_API void fluid_solve_bounds(int numReadIn, Chunk ** chunkViewC, Environm int neighborIndex; for(chunkIndex = 0; chunkIndex < numReadIn; chunkIndex++){ Chunk * current = chunkViewC[chunkIndex]; - fluid_solve_bounds_checker(current->d,0); - fluid_solve_bounds_checker(current->d0,0); - fluid_solve_bounds_checker(current->u,0); - fluid_solve_bounds_checker(current->v,0); - fluid_solve_bounds_checker(current->w,0); - fluid_solve_bounds_checker(current->u0,0); - fluid_solve_bounds_checker(current->v0,0); - fluid_solve_bounds_checker(current->w0,0); - fluid_solve_bounds_checker(current->bounds,BOUND_MAX_VALUE); - fluid_solve_bounds_checker(current->divergenceCache,0); - fluid_solve_bounds_checker(current->pressureCache,0); + fluid_solve_bounds_checker(current->d,0,current->incomingDensity); + fluid_solve_bounds_checker(current->d0,0,NULL); + fluid_solve_bounds_checker(current->u,0,NULL); + fluid_solve_bounds_checker(current->v,0,NULL); + fluid_solve_bounds_checker(current->w,0,NULL); + fluid_solve_bounds_checker(current->u0,0,NULL); + fluid_solve_bounds_checker(current->v0,0,NULL); + fluid_solve_bounds_checker(current->w0,0,NULL); + fluid_solve_bounds_checker(current->bounds,BOUND_MAX_VALUE,NULL); + fluid_solve_bounds_checker(current->divergenceCache,0,NULL); + fluid_solve_bounds_checker(current->pressureCache,0,current->incomingPressure); } } diff --git a/src/main/c/src/fluid/sim/pressurecell/pressure.c b/src/main/c/src/fluid/sim/pressurecell/pressure.c index 7beeb70a..82712f4b 100644 --- a/src/main/c/src/fluid/sim/pressurecell/pressure.c +++ b/src/main/c/src/fluid/sim/pressurecell/pressure.c @@ -7,7 +7,7 @@ /** * Divergence level we warn at */ -#define DIVERGENCE_WARNING_VALUE 100.0f +#define DIVERGENCE_WARNING_VALUE 300.0f /** * Approximates the pressure for this chunk