  The paper is concerned with the mathematical theory and numerical
approximation of systems of partial differential equations (pde) of hyperbolic,
pseudo-parabolic type. Some mathematical properties of the
initial-boundary-value problem (ibvp) with Dirichlet boundary conditions are
first studied. They include the weak formulation, well-posedness and existence
of traveling wave solutions connecting two states, when the equations are
considered as a variant of a conservation law. Then, the numerical
approximation consists of a spectral approximation in space based on Legendre
polynomials along with a temporal discretization with strong stability
preserving (SSP) property. The convergence of the semidiscrete approximation is
proved under suitable regularity conditions on the data. The choice of the
temporal discretization is justified in order to guarantee the stability of the
full discretization when dealing with nonsmooth initial conditions. A
computational study explores the performance of the fully discrete scheme with
regular and nonregular data.
