A continuous and reliable velocity field is required for the simulation of mass and heat transport in an aquifer system. In this paper, a cubic-spline technique is proposed to calculate nodal Darcian velocities after the flow field is obtained. With this technique, a continuous velocity field is produced, and mass conservation is satisfied to an acceptable degree. Its requirements for computer storage are very small, and the calculation speed is very fast; it can be implemented in parallel. A synthetic diffusion problem is used to examine the accuracy of this technique. The results agree very well with the analytical solution. Application to another problem shows that this technique can be applied to any heterogeneous and anisotropic aquifer. The results of the two sample problems are also compared with the Galerkin finite element method <Yeh, 1981> and the dual mesh method <Batu, 1984>. It is observed that its agreement with the Galerkin method is much better than the dual mesh method.