A new three-dimensional, non-hydrostatic free surface flow model is presented. For simulating water wave motions over uneven bottoms, the model employs an explicit project method on a Cartesian the staggered gird system to solve the complete three-dimensional Navier–Stokes equations. A bi-conjugated gradient method with a pre-conditioning procedure is used to solve the resulting matrix system. The model is capable of resolving non-hydrostatic pressure by incorporating the integral method of the top-layer pressure treatment, and predicting wave propagation and interaction over irregular bottom by including a partial bottom-cell treatment. Four examples of surface wave propagation are used to demonstrate the capability of the model. Using a small of vertical layers (e.g. 2–3 layers), it is shown that the model could effectively and accurately resolve wave shoaling, non-linearity, dispersion, fission, refraction, and diffraction phenomena.