This work presents a new multilayer nonhydrostatic formulation for surface water waves. The new governing equations define velocities and pressure at an arbitrary location of a vertical layer and only contain spatial derivatives of maximum second order. Stoke-type Fourier and shoaling analyses are carried out to scrutinize the mathematical properties of the new formulation, subsequently optimizing the representative interface and the location to define variables in each layer to improve model accuracy. Following the analysis, the one-layer model exhibits accurate linear and nonlinear characteristics up to kd = I, demonstrating similar solution accuracy to the existing second-order Boussinesq-type models. The two-layer model with optimized coefficients can maintain its linear and nonlinear accuracy up to kd = 4I, which boasts of better solution accuracy a larger application range than most existing fourth-order Boussinesq model and two-layer Boussinesq models. The three-layer model presents accurate linear and nonlinear characteristics up to kd = 10Ï, effectively removing any shallow water limitation. The current multilayer nonhydrostatic water wave model does not predefine the vertical flow structures, and more accurate vertical velocity distributions can be obtained by considering the velocity profiles in coefficient optimization