Added first implementation of flux computation on cells for 2D and 3D.

Rigorous testing still required. There is a tricky corner case in 3D where only 2 facets are involved in the flux, need to think about how to deal with this.