     lbm_vars *h_vars = (lbm_vars*)malloc(sizeof(lbm_vars));
     initObstacles(h_vars);
     initVelocity(h_vars);
-    initRho(h_vars);
     
     // Initialization of the populations at equilibrium with the given velocity.
     for (int y = 0; y < NY; y++) {
         for (int x = 0; x < NX; x++) {
-            h_equilibrium(h_vars->fin[x][y], h_vars->rho[x][y], h_vars->vel[x][y]);
+            h_equilibrium(h_vars->fin[x][y], 1.0, h_vars->vel[x][y]);
         }
     }
     

                s2 += d_vars->fin[0][y][d_consts.col[1][i]];
                s3 += d_vars->fin[0][y][d_consts.col[2][i]];
            }
            d_vars->rho[0][y] = 1./(1 - u0) * (s2 + 2*s3);
        }

        // Compute equilibrium
        double feq[9];
        d_equilibrium(feq, d_vars->rho[x][y], u0, u1);

        if (x == 0) {
            for (size_t i = 0, f = d_consts.col[0][i]; i < 3; f = d_consts.col[0][++i]) {

    }
}

static void initCol(size_t* col, ssize_t v0)
{
    for (int f = 0, i = 0; f < 9 && i < 3; f++) {
