Abstract | The present work describes a procedure for the numerical evaluation of the classical integral-transform solution of the transient elastodynamic point-load (axisymmetric) Lamb’s problem. This solution involves integrals of rapidly oscillatory functions over semi-infinite intervals and inversion of one-sided (time) Laplace transforms. These features introduce difficulties for a numerical treatment and constitute a challenging problem in trying to obtain results for quantities (e.g. displacements) in the interior of the half-space. To deal with the oscillatory integrands, which in addition may take very large values (pseudo-pole behavior) at certain points, we follow the concept of Longman’s method but using as accelerator in the summation procedure a modified Epsilon algorithm instead of the standard Euler’s transformation. Also, an adaptive procedure using the Gauss 32-point rule is introduced to integrate in the vicinity of the pseudo-pole. The numerical Laplace-transform inversion is based on the robust Fourier-series technique of Dubner/Abate-Crump-Durbin. Extensive results are given for sub-surface displacements, whereas the limit-case results for the surface displacements compare very favorably with previous exact results.