A three-dimensional (3D) gradient-augmented level set (GALS) two-phase flow model with a pretreated reinitialization procedure is developed to simulate violent sloshing in a cuboid tank. Based on a two-dimensional (2D) GALS method, 3D Hermite, and 3D Lagrange polynomial schemes are derived to interpolate the level set function and the velocity field at arbitrary positions over a cell, respectively. A reinitialization procedure is performed on a 3D narrow band to treat the strongly distorted interface and improve computational efficiency. In addition, an identification-correction technique is proposed and incorporated into the reinitialization procedure to treat the tiny droplet which can distort the free surface shape, even lead to computation failure. To validate the accuracy of the present GALS method and the effectiveness of the proposed identification-correction technique, a 3D velocity advection case is first simulated. The present method is validated to have better mass conservation property than the classical level set and original GALS methods. Also, distorted and thin interfaces are well captured on all grid resolutions by the present GALS method. Then, sloshing under coupled surge and sway excitation, sloshing under rotational excitation are simulated. Good agreements are obtained when the present wave and pressure results are compared with the experimental and numerical results. In addition, the highly nonlinear free surface is observed, and the relationship between the excitation frequency and the impulsive pressure is investigated.