The studies of the dynamic behaviors of nonlinear models arising in ocean engineering play a significant role in our daily activities. In this study, we investigate the coupled Boussinesq equation which arises in the shallow water waves for two-layered fluid flow. The modified exp $$(-\varphi (\zeta ))$$ -expansion function method is utilized in reaching the solutions to this equation such as the topological kink-type soliton and singular soliton solutions. The interesting 2D and 3D graphics of the obtained analytical solutions in this study are presented. Via one of the reported analytical solutions, the finite forward difference method is used in obtaining the approximate numerical and exact solutions to this equation. The Fourier–Von Neumann analysis is used in checking the stability of the used numerical method with the studied model. The $$L_{2}$$ and $$L_{\infty }$$ error norms are computed. We finally present a comprehensive conclusion to this study.