A three-dimensional numerical model based on the full Navier-Stokes equations (NSE) inσ-coordinate is developed in this study. The σ-coordinate transformation is first introduced to map the irregular physical domain with the wavy free surface and uneven bottom to the regular computational domain with the shape of a rectangular prism. Using the chain rule of partial differentiation, a new set of governing equations is derived in the σ-coordinate from the original NSE defined in the Cartesian coordinate. The operator splitting method (Li and Yu, Int. J. Num. Meth. Fluids 1996; 23: 485-501), which splits the solution procedure into the advection, diffusion, and propagation steps, is used to solve the modified NSE. The model is first tested for mass and energy conservation as well as mesh convergence by using an example of water sloshing in a confined tank. Excellent agreements between numerical results and analytical solutions are obtained. The model is then used to simulate two-and three-dimensional solitary waves propagating in constant depth. Very good agreements between numerical results and analytical solutions are obtained for both free surface displacements and velocities. Finally, a more realistic case of periodic wave train passing through a submerged breakwater is simulated. Comparisons between numerical results and experimental data are promising. The model is proven to be an accurate tool for consequent studies of wave-structure interaction. Copyright