AbstractThis paper is concerned with the numerical approximation of initial‐boundary‐value problems of a three‐parameter family of Bona–Smith systems, derived as a model for the propagation of surface waves under a physical Boussinesq regime. The work proposed here is focused on the corresponding problem with Dirichlet boundary conditions and its approximation in space with spectral methods based on Jacobi polynomials, which are defined from the orthogonality with respect to some weighted inner product. Well‐posedness of the problem on the corresponding weighted Sobolev spaces is first analyzed and existence and uniqueness of solution, locally in time, are proved. Then, the spectral Galerkin semidiscrete scheme and some detailed comments on its implementation are introduced. The existence of numerical solution and error estimates on those weighted Sobolev spaces are established. Finally, the choice of the time integrator to complete the full discretization takes care of different stability issues that may be relevant when approximating the semidiscrete system. Some numerical experiments illustrate the results.