in mass formulation p_neumann has to be rho_fr times larger than before rho_fr = 1000 1e-5 * rho_fr = 1e-2