As an alternative to domain discretization methods, the boundary element method (BEM) provides a powerful tool for the calculation of dynamic structural response in frequency and time domain. Field equations of motion and boundary conditions are cast into boundary integral equations (BIE), which are discretized only on the boundary. Fundamental solutions are used as weighting functions in the BIE which fulfil the Sommerfeld radiation condition, i.e., the energy radiation into a surrounding medium is modelled correctly. Therefore, infinite and semi-infinite domains can be effectively treated by the method. The soil represents such a semi-infinite domain in soil-structure-interaction problems. The response to vibratory loads superimposed to static pre-loads can often be calculated by linear viscoelastic constitutive equations. Conventional viscoelastic constitutive equations can be generalized by taking fractional order time derivatives into account. In the present paper two time domain BEM approaches including generalized viscoelastic behaviour are compared with the Laplace domain BEM approach and subsequent numerical inverse transformation. One of the presented time domain approaches uses an analytical integration of the elastodynamic BIE in a time step. Viscoelastic constitutive properties are introduced after Laplace transformation by means of an elastic–viscoelastic correspondence principle. The transient response is obtained by inverse transformation in each time step. The other time domain approach is based on the so-called `convolution quadrature method'. In this formulation, the convolution integral in the BIE is numerically approximated by a quadrature formula whose weights are determined by the same Laplace transformed fundamental solutions used in the first method and a linear multistep method. A numerical study of wave propagation problems in 3-d viscoelastic continuum is performed for comparing the three BEM formulations.