diff --git a/src/main/c/includes/fluid/queue/chunk.h b/src/main/c/includes/fluid/queue/chunk.h
index 987e8a39..fead4f85 100644
--- a/src/main/c/includes/fluid/queue/chunk.h
+++ b/src/main/c/includes/fluid/queue/chunk.h
@@ -84,6 +84,11 @@ typedef struct {
      */
     float * wTempCache;
 
+    /**
+     * Temp cache for storing pressure during current iteration
+     */
+    float * pressureTempCache;
+
     /**
      * The bitmask which tracks valid neighbors
      */
diff --git a/src/main/c/includes/fluid/sim/pressurecell/pressure.h b/src/main/c/includes/fluid/sim/pressurecell/pressure.h
index 8d0b6dce..cd594afb 100644
--- a/src/main/c/includes/fluid/sim/pressurecell/pressure.h
+++ b/src/main/c/includes/fluid/sim/pressurecell/pressure.h
@@ -12,5 +12,10 @@
 */
 LIBRARY_API void pressurecell_approximate_pressure(Environment * environment, Chunk * chunk);
 
+/**
+ * Approximates the divergence for this chunk
+*/
+LIBRARY_API void pressurecell_approximate_divergence(Environment * environment, Chunk * chunk);
+
 
 #endif
\ No newline at end of file
diff --git a/src/main/c/includes/fluid/sim/pressurecell/pressurecell.h b/src/main/c/includes/fluid/sim/pressurecell/pressurecell.h
index 2c787cb8..c3ed98f2 100644
--- a/src/main/c/includes/fluid/sim/pressurecell/pressurecell.h
+++ b/src/main/c/includes/fluid/sim/pressurecell/pressurecell.h
@@ -21,7 +21,7 @@
  * @param environment The environment data
  * @param timestep The timestep to simulate by
  */
-LIBRARY_API void fluid_grid2_simulate(int numChunks, Chunk ** passedInChunks, Environment * environment, float timestep);
+LIBRARY_API void fluid_pressurecell_simulate(int numChunks, Chunk ** passedInChunks, Environment * environment, float timestep);
 
 
 
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 9dce7894..39426586 100644
--- a/src/main/c/includes/fluid/sim/pressurecell/solver_consts.h
+++ b/src/main/c/includes/fluid/sim/pressurecell/solver_consts.h
@@ -22,5 +22,10 @@
  */
 #define FLUID_PRESSURECELL_VISCOSITY_CONSTANT 0.0001f
 
+/**
+ * Amount that density contributes to the pressure
+ */
+#define FLUID_PRESSURECELL_DENSITY_CONST 0.001f
+
 
 #endif
\ No newline at end of file
diff --git a/src/main/c/src/fluid/env/environment.c b/src/main/c/src/fluid/env/environment.c
index 83f7cde1..8d09b422 100644
--- a/src/main/c/src/fluid/env/environment.c
+++ b/src/main/c/src/fluid/env/environment.c
@@ -17,6 +17,7 @@ LIBRARY_API Environment * fluid_environment_create(){
     rVal->queue.cellularQueue = NULL;
     rVal->queue.gridQueue = NULL;
     rVal->queue.grid2Queue = NULL;
+    rVal->consts.dt = 0.02f;
 
     //allocate arrays
     fluid_environment_allocate_arrays(rVal);
diff --git a/src/main/c/src/fluid/queue/chunk.c b/src/main/c/src/fluid/queue/chunk.c
index e96f5be2..bfe6d07d 100644
--- a/src/main/c/src/fluid/queue/chunk.c
+++ b/src/main/c/src/fluid/queue/chunk.c
@@ -12,5 +12,6 @@ LIBRARY_API Chunk * chunk_create(){
     rVal->uTempCache = (float *)calloc(1,DIM*DIM*DIM*sizeof(float));
     rVal->vTempCache = (float *)calloc(1,DIM*DIM*DIM*sizeof(float));
     rVal->wTempCache = (float *)calloc(1,DIM*DIM*DIM*sizeof(float));
+    rVal->pressureTempCache = (float *)calloc(1,DIM*DIM*DIM*sizeof(float));
     return rVal;
 }
diff --git a/src/main/c/src/fluid/sim/pressurecell/density.c b/src/main/c/src/fluid/sim/pressurecell/density.c
index 50409c33..b3f80219 100644
--- a/src/main/c/src/fluid/sim/pressurecell/density.c
+++ b/src/main/c/src/fluid/sim/pressurecell/density.c
@@ -46,5 +46,71 @@ LIBRARY_API void pressurecell_diffuse_density(Environment * environment, Chunk *
  * Advects the density of this chunk
 */
 LIBRARY_API void pressurecell_advect_density(Environment * environment, Chunk * chunk){
+    int x, y, z;
+    float * densityArr = chunk->d[CENTER_LOC];
+    float * densityTemp = chunk->dTempCache;
+    float * u = chunk->u[CENTER_LOC];
+    float * v = chunk->v[CENTER_LOC];
+    float * w = chunk->w[CENTER_LOC];
+    int x0, x1, y0, y1, z0, z1;
+    float xp, yp, zp;
+    float s0, s1, t0, t1, u0, u1;
+    float interpolated;
+    for(y = 1; y < DIM-1; y++){
+        //TODO: eventually skip y levels if there is no density to advect
+        for(z = 1; z < DIM-1; z++){
+            for(x = 1; x < DIM-1; x++){
+                //calculate the real (float) position we are at
+                xp = x - u[IX(x,y,z)] * environment->consts.dt;
+                yp = y - v[IX(x,y,z)] * environment->consts.dt;
+                zp = z - w[IX(x,y,z)] * environment->consts.dt;
 
+                //clamp to border
+                x0 = xp;
+                y0 = yp;
+                z0 = zp;
+                //find far border
+                x1 = x0 + 1;
+                y1 = y0 + 1;
+                z1 = z0 + 1;
+
+                //calculate the percentage to sample from each side of border
+                s1 = xp-x0;
+                s0 = 1-s1;
+
+                t1 = yp-y0;
+                t0 = 1-t1;
+
+                u1 = zp-z0;
+                u0 = 1-u1;
+
+
+                if(
+                    x0 < 0 || y0 < 0 || z0 < 0 ||
+                    x0 > DIM-1 || y0 > DIM-1 || z0 > DIM-1 ||
+                    x1 < 0 || y1 < 0 || z1 < 0 ||
+                    x1 > DIM-1 || y1 > DIM-1 || z1 > DIM-1
+                ){
+                    printf("advect dens: %d %d %d   %d %d %d ---  %f %f %f\n", x0, y0, z0, x1, y1, z1, x, y, z);
+                    fflush(stdout);
+                }
+
+                interpolated = 
+                s0*(
+                    t0*u0*densityTemp[IX(x0,y0,z0)]+
+                    t1*u0*densityTemp[IX(x0,y1,z0)]+
+                    t0*u1*densityTemp[IX(x0,y0,z1)]+
+                    t1*u1*densityTemp[IX(x0,y1,z1)]
+                )+
+                s1*(
+                    t0*u0*densityTemp[IX(x1,y0,z0)]+
+                    t1*u0*densityTemp[IX(x1,y1,z0)]+
+                    t0*u1*densityTemp[IX(x1,y0,z1)]+
+                    t1*u1*densityTemp[IX(x1,y1,z1)]
+                );
+
+                densityArr[IX(x,y,z)] = interpolated;
+            }
+        }
+    }
 }
\ No newline at end of file
diff --git a/src/main/c/src/fluid/sim/pressurecell/pressure.c b/src/main/c/src/fluid/sim/pressurecell/pressure.c
index 34645150..4b0a2ade 100644
--- a/src/main/c/src/fluid/sim/pressurecell/pressure.c
+++ b/src/main/c/src/fluid/sim/pressurecell/pressure.c
@@ -1,6 +1,7 @@
 
 
 #include "fluid/sim/pressurecell/pressure.h"
+#include "fluid/sim/pressurecell/solver_consts.h"
 
 
 
@@ -8,6 +9,61 @@
  * Approximates the pressure for this chunk
 */
 LIBRARY_API void pressurecell_approximate_pressure(Environment * environment, Chunk * chunk){
-
+    int x, y, z;
+    //values stored across frames
+    float * presureCache = chunk->pressureCache[CENTER_LOC];
+    float * divArr = chunk->divergenceCache[CENTER_LOC];
+    //temporary caches
+    float * pressureTemp = chunk->pressureTempCache;
+    float gridSpacing = 1;
+    float du, dv, dw;
+    float newPressure;
+    for(z = 1; z < DIM-1; z++){
+        for(y = 1; y < DIM-1; y++){
+            for(x = 1; x < DIM-1; x++){
+                //compute the new pressure value
+                newPressure = 
+                    (divArr[IX(x+1,y,z)] - divArr[IX(x-1,y,z)]) / (2.0f * gridSpacing) +
+                    (divArr[IX(x+1,y,z)] - divArr[IX(x-1,y,z)]) / (2.0f * gridSpacing) +
+                    (divArr[IX(x+1,y,z)] - divArr[IX(x-1,y,z)]) / (2.0f * gridSpacing)
+                ;
+                newPressure = newPressure / (gridSpacing * 2);
+                pressureTemp[IX(x,y,z)] = newPressure;
+            }
+        }
+    }
+}
+
+/**
+ * Approximates the divergence for this chunk
+*/
+LIBRARY_API void pressurecell_approximate_divergence(Environment * environment, Chunk * chunk){
+    int x, y, z;
+    //values stored across frames
+    float * uArr = chunk->u[CENTER_LOC];
+    float * vArr = chunk->v[CENTER_LOC];
+    float * wArr = chunk->w[CENTER_LOC];
+    float * divArr = chunk->divergenceCache[CENTER_LOC];
+    float * presureCache = chunk->pressureCache[CENTER_LOC];
+    //temporary caches
+    float * pressureTemp = chunk->pressureTempCache;
+    float gridSpacing = 1;
+    float du, dv, dw;
+    float newDivergence;
+    for(z = 1; z < DIM-1; z++){
+        for(y = 1; y < DIM-1; y++){
+            for(x = 1; x < DIM-1; x++){
+                //compute divergence
+                du = (uArr[IX(x+1,y,z)] - uArr[IX(x-1,y,z)]) / (2.0f * gridSpacing);
+                dv = (vArr[IX(x+1,y,z)] - vArr[IX(x-1,y,z)]) / (2.0f * gridSpacing);
+                dw = (wArr[IX(x+1,y,z)] - wArr[IX(x-1,y,z)]) / (2.0f * gridSpacing);
+                newDivergence = du+dv+dw;
+                divArr[IX(x,y,z)] = newDivergence;
+                
+                //store pressure value from this frame
+                presureCache[IX(x,y,z)] = pressureTemp[IX(x,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 7b57a877..705346c3 100644
--- a/src/main/c/src/fluid/sim/pressurecell/pressurecell.c
+++ b/src/main/c/src/fluid/sim/pressurecell/pressurecell.c
@@ -9,12 +9,10 @@
 #include "fluid/env/utilities.h"
 #include "fluid/queue/chunkmask.h"
 #include "fluid/queue/chunk.h"
-#include "fluid/sim/grid2/flux.h"
-#include "fluid/sim/grid2/grid2.h"
-#include "fluid/sim/grid2/solver_consts.h"
-#include "fluid/sim/grid2/velocity.h"
-#include "fluid/sim/grid2/density.h"
-#include "fluid/sim/grid2/utilities.h"
+#include "fluid/sim/pressurecell/density.h"
+#include "fluid/sim/pressurecell/pressure.h"
+#include "fluid/sim/pressurecell/solver_consts.h"
+#include "fluid/sim/pressurecell/velocity.h"
 
 /**
  * Used for storing timings
@@ -33,11 +31,65 @@ LIBRARY_API void fluid_pressurecell_simulate(
     //This is the section of non-parallel code
     //
 
-    //update ODE solver data
-    environment->state.grid2.diffuseData.dt = timestep;
+    //
+    //This is the section of parallel code
+    //
+
+    //
+    // Velocity phase
+    //
+
     for(int i = 0; i < numChunks; i++){
         Chunk * currentChunk = chunks[i];
-        fluid_grid2_update_ghost_flux(environment,currentChunk);
+        pressurecell_add_velocity(environment,currentChunk);
+    }
+
+    for(int i = 0; i < numChunks; i++){
+        Chunk * currentChunk = chunks[i];
+        pressurecell_diffuse_velocity(environment,currentChunk);
+    }
+
+    for(int i = 0; i < numChunks; i++){
+        Chunk * currentChunk = chunks[i];
+        pressurecell_advect_velocity(environment,currentChunk);
+    }
+
+    for(int i = 0; i < numChunks; i++){
+        Chunk * currentChunk = chunks[i];
+        pressurecell_approximate_pressure(environment,currentChunk);
+    }
+
+    for(int i = 0; i < numChunks; i++){
+        Chunk * currentChunk = chunks[i];
+        pressurecell_interpolate_velocity(environment,currentChunk);
+    }
+
+    for(int i = 0; i < numChunks; i++){
+        Chunk * currentChunk = chunks[i];
+        pressurecell_approximate_divergence(environment,currentChunk);
+    }
+
+
+
+
+    //
+    // Density phase
+    //
+
+
+    for(int i = 0; i < numChunks; i++){
+        Chunk * currentChunk = chunks[i];
+        pressurecell_add_density(environment,currentChunk);
+    }
+
+    for(int i = 0; i < numChunks; i++){
+        Chunk * currentChunk = chunks[i];
+        pressurecell_diffuse_density(environment,currentChunk);
+    }
+
+    for(int i = 0; i < numChunks; i++){
+        Chunk * currentChunk = chunks[i];
+        pressurecell_advect_density(environment,currentChunk);
     }
 }
 
diff --git a/src/main/c/src/fluid/sim/pressurecell/velocity.c b/src/main/c/src/fluid/sim/pressurecell/velocity.c
index 29fce49b..689f216a 100644
--- a/src/main/c/src/fluid/sim/pressurecell/velocity.c
+++ b/src/main/c/src/fluid/sim/pressurecell/velocity.c
@@ -80,12 +80,143 @@ LIBRARY_API void pressurecell_diffuse_velocity(Environment * environment, Chunk
  * Advects the velocity of this chunk
 */
 LIBRARY_API void pressurecell_advect_velocity(Environment * environment, Chunk * chunk){
+    int x, y, z;
+    float * uArr = chunk->u0[CENTER_LOC];
+    float * uTemp = chunk->uTempCache;
+    float * vArr = chunk->v0[CENTER_LOC];
+    float * vTemp = chunk->vTempCache;
+    float * wArr = chunk->w0[CENTER_LOC];
+    float * wTemp = chunk->wTempCache;
+    int x0, x1, y0, y1, z0, z1;
+    float xp, yp, zp;
+    float s0, s1, t0, t1, u0, u1;
+    float interpolated;
+    for(y = 1; y < DIM-1; y++){
+        //TODO: eventually skip y levels if there is no density to advect
+        for(z = 1; z < DIM-1; z++){
+            for(x = 1; x < DIM-1; x++){
+                //calculate the real (float) position we are at
+                xp = x - uTemp[IX(x,y,z)] * environment->consts.dt;
+                yp = y - vTemp[IX(x,y,z)] * environment->consts.dt;
+                zp = z - wTemp[IX(x,y,z)] * environment->consts.dt;
 
+                //clamp to border
+                x0 = xp;
+                y0 = yp;
+                z0 = zp;
+                //find far border
+                x1 = x0 + 1;
+                y1 = y0 + 1;
+                z1 = z0 + 1;
+
+                //calculate the percentage to sample from each side of border
+                s1 = xp-x0;
+                s0 = 1-s1;
+
+                t1 = yp-y0;
+                t0 = 1-t1;
+
+                u1 = zp-z0;
+                u0 = 1-u1;
+
+
+                if(
+                    x0 < 0 || y0 < 0 || z0 < 0 ||
+                    x0 > DIM-1 || y0 > DIM-1 || z0 > DIM-1 ||
+                    x1 < 0 || y1 < 0 || z1 < 0 ||
+                    x1 > DIM-1 || y1 > DIM-1 || z1 > DIM-1
+                ){
+                    printf("advect dens: %d %d %d   %d %d %d ---  %f %f %f\n", x0, y0, z0, x1, y1, z1, x, y, z);
+                    fflush(stdout);
+                }
+
+                interpolated = 
+                s0*(
+                    t0*u0*uTemp[IX(x0,y0,z0)]+
+                    t1*u0*uTemp[IX(x0,y1,z0)]+
+                    t0*u1*uTemp[IX(x0,y0,z1)]+
+                    t1*u1*uTemp[IX(x0,y1,z1)]
+                )+
+                s1*(
+                    t0*u0*uTemp[IX(x1,y0,z0)]+
+                    t1*u0*uTemp[IX(x1,y1,z0)]+
+                    t0*u1*uTemp[IX(x1,y0,z1)]+
+                    t1*u1*uTemp[IX(x1,y1,z1)]
+                );
+                uArr[IX(x,y,z)] = interpolated;
+
+                interpolated = 
+                s0*(
+                    t0*u0*vTemp[IX(x0,y0,z0)]+
+                    t1*u0*vTemp[IX(x0,y1,z0)]+
+                    t0*u1*vTemp[IX(x0,y0,z1)]+
+                    t1*u1*vTemp[IX(x0,y1,z1)]
+                )+
+                s1*(
+                    t0*u0*vTemp[IX(x1,y0,z0)]+
+                    t1*u0*vTemp[IX(x1,y1,z0)]+
+                    t0*u1*vTemp[IX(x1,y0,z1)]+
+                    t1*u1*vTemp[IX(x1,y1,z1)]
+                );
+                vArr[IX(x,y,z)] = interpolated;
+
+                interpolated = 
+                s0*(
+                    t0*u0*wTemp[IX(x0,y0,z0)]+
+                    t1*u0*wTemp[IX(x0,y1,z0)]+
+                    t0*u1*wTemp[IX(x0,y0,z1)]+
+                    t1*u1*wTemp[IX(x0,y1,z1)]
+                )+
+                s1*(
+                    t0*u0*wTemp[IX(x1,y0,z0)]+
+                    t1*u0*wTemp[IX(x1,y1,z0)]+
+                    t0*u1*wTemp[IX(x1,y0,z1)]+
+                    t1*u1*wTemp[IX(x1,y1,z1)]
+                );
+                wArr[IX(x,y,z)] = interpolated;
+            }
+        }
+    }
 }
 
 /**
  * Interpolates between the advected velocity and the previous frame's velocity by the pressure divergence amount
 */
 LIBRARY_API void pressurecell_interpolate_velocity(Environment * environment, Chunk * chunk){
+    int x, y, z;
+    float * presureCache = chunk->pressureCache[CENTER_LOC];
+    float * uArr = chunk->u[CENTER_LOC];
+    float * vArr = chunk->v[CENTER_LOC];
+    float * wArr = chunk->w[CENTER_LOC];
+    float * uTemp = chunk->uTempCache;
+    float * vTemp = chunk->vTempCache;
+    float * wTemp = chunk->wTempCache;
+    //temporary caches
+    float * pressureTemp = chunk->pressureTempCache;
+    float gridSpacing = 1;
+    float du, dv, dw;
+    float pressureDivergence;
+    for(z = 1; z < DIM-1; z++){
+        for(y = 1; y < DIM-1; y++){
+            for(x = 1; x < DIM-1; x++){
+                //compute the new pressure value
+                pressureDivergence = 
+                    (presureCache[IX(x+1,y,z)] - presureCache[IX(x-1,y,z)]) / (2.0f * gridSpacing) +
+                    (presureCache[IX(x+1,y,z)] - presureCache[IX(x-1,y,z)]) / (2.0f * gridSpacing) +
+                    (presureCache[IX(x+1,y,z)] - presureCache[IX(x-1,y,z)]) / (2.0f * gridSpacing)
+                ;
+                pressureDivergence = pressureDivergence / (gridSpacing * 2);
+                
+                //project the pressure gradient onto the velocity field
+                uTemp[IX(x,y,z)] = uTemp[IX(x,y,z)] - pressureDivergence;
+                vTemp[IX(x,y,z)] = vTemp[IX(x,y,z)] - pressureDivergence;
+                wTemp[IX(x,y,z)] = wTemp[IX(x,y,z)] - pressureDivergence;
 
+                //map the new velocity field onto the old one
+                uArr[IX(x,y,z)] = uTemp[IX(x,y,z)];
+                vArr[IX(x,y,z)] = vTemp[IX(x,y,z)];
+                wArr[IX(x,y,z)] = wTemp[IX(x,y,z)];
+            }
+        }
+    }
 }
\ No newline at end of file
diff --git a/src/test/c/fluid/sim/pressurecell/advection_tests.c b/src/test/c/fluid/sim/pressurecell/advection_tests.c
new file mode 100644
index 00000000..34cc1e42
--- /dev/null
+++ b/src/test/c/fluid/sim/pressurecell/advection_tests.c
@@ -0,0 +1,119 @@
+#include <math.h>
+
+#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/pressurecell/density.h"
+#include "fluid/sim/pressurecell/velocity.h"
+#include "fluid/sim/pressurecell/solver_consts.h"
+#include "math/ode/multigrid.h"
+#include "../../../util/chunk_test_utils.h"
+#include "../../../util/test.h"
+
+/**
+ * Error margin for tests
+ */
+#define FLUID_PRESSURE_CELL_ERROR_MARGIN 0.00001f
+
+/**
+ * Number of chunks
+ */
+#define CHUNK_DIM 4
+
+/**
+ * Testing advecting values
+ */
+int fluid_sim_pressurecell_advection_test1(){
+    printf("fluid_sim_pressurecell_advection_test1\n");
+    int rVal = 0;
+    Environment * env = fluid_environment_create();
+    Chunk ** queue = NULL;
+    queue = createChunkGrid(env,CHUNK_DIM,CHUNK_DIM,CHUNK_DIM);
+    int chunkCount = arrlen(queue);
+
+
+
+    //setup chunk values
+    Chunk * currentChunk = queue[0];
+    currentChunk->dTempCache[IX(4,4,4)] = MAX_FLUID_VALUE;
+    currentChunk->u[CENTER_LOC][IX(4,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+    currentChunk->u[CENTER_LOC][IX(5,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+
+    //actually simulate
+    pressurecell_advect_density(env,currentChunk);
+
+    //test the result
+    float expected, actual;
+
+    //
+    // cell that originall had values
+    //
+    expected = MAX_FLUID_VALUE - FLUID_PRESSURECELL_MAX_VELOCITY * env->consts.dt * MAX_FLUID_VALUE;
+    actual = currentChunk->d[CENTER_LOC][IX(4,4,4)];
+    if(fabs(expected - actual) > FLUID_PRESSURE_CELL_ERROR_MARGIN){
+        rVal += assertEqualsFloat(expected,actual,"Failed to advect density correctly (4,4,4)!  expected: %f    actual: %f  \n");
+    }
+
+
+    //
+    // neighbors
+    //
+    expected = FLUID_PRESSURECELL_MAX_VELOCITY * env->consts.dt * MAX_FLUID_VALUE;
+    actual = currentChunk->d[CENTER_LOC][IX(5,4,4)];
+    if(fabs(expected - actual) > FLUID_PRESSURE_CELL_ERROR_MARGIN){
+        rVal += assertEqualsFloat(expected,actual,"Failed to advect density correctly (5,4,4)!  expected: %f    actual: %f  \n");
+    }
+
+    return rVal;
+}
+
+/**
+ * Testing advecting values
+ */
+int fluid_sim_pressurecell_advection_test2(){
+    printf("fluid_sim_pressurecell_advection_test2\n");
+    int rVal = 0;
+    Environment * env = fluid_environment_create();
+    Chunk ** queue = NULL;
+    queue = createChunkGrid(env,CHUNK_DIM,CHUNK_DIM,CHUNK_DIM);
+    int chunkCount = arrlen(queue);
+
+
+
+    //setup chunk values
+    Chunk * currentChunk = queue[0];
+    currentChunk->uTempCache[IX(4,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+    currentChunk->u[CENTER_LOC][IX(4,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+
+    //actually simulate
+    pressurecell_advect_velocity(env,currentChunk);
+
+    //test the result
+    float expected, actual;
+    //
+    // cell that originall had values
+    //
+    expected = FLUID_PRESSURECELL_MAX_VELOCITY - FLUID_PRESSURECELL_MAX_VELOCITY * env->consts.dt;
+    actual = currentChunk->u0[CENTER_LOC][IX(4,4,4)];
+    if(fabs(expected - actual) > FLUID_PRESSURE_CELL_ERROR_MARGIN){
+        rVal += assertEqualsFloat(expected,actual,"Failed to diffuse velocity correctly (4,4,4)!  expected: %f    actual: %f  \n");
+    }
+
+    return rVal;
+}
+
+/**
+ * Testing advecting values
+ */
+int fluid_sim_pressurecell_advection_tests(int argc, char **argv){
+    int rVal = 0;
+
+    rVal += fluid_sim_pressurecell_advection_test1();
+    rVal += fluid_sim_pressurecell_advection_test2();
+
+    return rVal;
+}
\ No newline at end of file
diff --git a/src/test/c/fluid/sim/pressurecell/pressure_tests.c b/src/test/c/fluid/sim/pressurecell/pressure_tests.c
new file mode 100644
index 00000000..b2cee159
--- /dev/null
+++ b/src/test/c/fluid/sim/pressurecell/pressure_tests.c
@@ -0,0 +1,75 @@
+#include <math.h>
+
+#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/pressurecell/pressure.h"
+#include "fluid/sim/pressurecell/solver_consts.h"
+#include "math/ode/multigrid.h"
+#include "../../../util/chunk_test_utils.h"
+#include "../../../util/test.h"
+
+/**
+ * Error margin for tests
+ */
+#define FLUID_PRESSURE_CELL_ERROR_MARGIN 0.00001f
+
+/**
+ * Number of chunks
+ */
+#define CHUNK_DIM 4
+
+/**
+ * Testing pressure values
+ */
+int fluid_sim_pressurecell_pressure_test1(){
+    printf("fluid_sim_pressurecell_pressure_test1\n");
+    int rVal = 0;
+    Environment * env = fluid_environment_create();
+    Chunk ** queue = NULL;
+    queue = createChunkGrid(env,CHUNK_DIM,CHUNK_DIM,CHUNK_DIM);
+    int chunkCount = arrlen(queue);
+
+
+
+    //setup chunk values
+    float deltaDensity = 0.01f;
+    Chunk * currentChunk = queue[0];
+    currentChunk->d[CENTER_LOC][IX(4,4,4)] = MAX_FLUID_VALUE - deltaDensity;
+    currentChunk->d0[CENTER_LOC][IX(4,4,4)] = deltaDensity;
+    currentChunk->u0[CENTER_LOC][IX(4,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+    currentChunk->v0[CENTER_LOC][IX(4,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+    currentChunk->w0[CENTER_LOC][IX(4,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+
+    //actually simulate
+    pressurecell_approximate_pressure(env,currentChunk);
+
+    //test the result
+    float expected, actual;
+
+    //
+    // cell that originall had values
+    //
+    // expected = MAX_FLUID_VALUE - FLUID_PRESSURECELL_MAX_VELOCITY * env->consts.dt * MAX_FLUID_VALUE;
+    // actual = currentChunk->pressureCache[CENTER_LOC][IX(4,4,4)];
+    // if(fabs(expected - actual) > FLUID_PRESSURE_CELL_ERROR_MARGIN){
+    //     rVal += assertEqualsFloat(expected,actual,"Failed to advect density correctly (4,4,4)!  expected: %f    actual: %f  \n");
+    // }
+
+    return rVal;
+}
+
+/**
+ * Testing pressure values
+ */
+int fluid_sim_pressurecell_pressure_tests(int argc, char **argv){
+    int rVal = 0;
+
+    rVal += fluid_sim_pressurecell_pressure_test1();
+
+    return rVal;
+}
\ No newline at end of file
diff --git a/src/test/c/fluid/sim/pressurecell/sim_e2e_tests.c b/src/test/c/fluid/sim/pressurecell/sim_e2e_tests.c
new file mode 100644
index 00000000..a21398e2
--- /dev/null
+++ b/src/test/c/fluid/sim/pressurecell/sim_e2e_tests.c
@@ -0,0 +1,60 @@
+#include <math.h>
+
+#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/pressurecell/pressure.h"
+#include "fluid/sim/pressurecell/pressurecell.h"
+#include "fluid/sim/pressurecell/solver_consts.h"
+#include "math/ode/multigrid.h"
+#include "../../../util/chunk_test_utils.h"
+#include "../../../util/test.h"
+
+
+/**
+ * Number of chunks
+ */
+#define CHUNK_DIM 4
+
+/**
+ * Testing pressure values
+ */
+int fluid_sim_pressurecell_sim_e2e_test1(){
+    printf("fluid_sim_pressurecell_sim_e2e_test1\n");
+    int rVal = 0;
+    Environment * env = fluid_environment_create();
+    Chunk ** queue = NULL;
+    queue = createChunkGrid(env,CHUNK_DIM,CHUNK_DIM,CHUNK_DIM);
+    int chunkCount = arrlen(queue);
+
+
+
+    //setup chunk values
+    float deltaDensity = 0.01f;
+    Chunk * currentChunk = queue[0];
+    currentChunk->d[CENTER_LOC][IX(4,4,4)] = MAX_FLUID_VALUE - deltaDensity;
+    currentChunk->d0[CENTER_LOC][IX(4,4,4)] = deltaDensity;
+    currentChunk->u0[CENTER_LOC][IX(4,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+    currentChunk->v0[CENTER_LOC][IX(4,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+    currentChunk->w0[CENTER_LOC][IX(4,4,4)] = FLUID_PRESSURECELL_MAX_VELOCITY;
+
+    //actually simulate
+    fluid_pressurecell_simulate(chunkCount,queue,env,env->consts.dt);
+
+    return rVal;
+}
+
+/**
+ * Testing pressure values
+ */
+int fluid_sim_pressurecell_sim_e2e_tests(int argc, char **argv){
+    int rVal = 0;
+
+    rVal += fluid_sim_pressurecell_sim_e2e_test1();
+
+    return rVal;
+}
\ No newline at end of file