The development of a time domain boundary element method for axisymmetric quasi-static poroelasticity is discussed. This new formulation, for the complete Biot consolidation theory, has the distinct advantage of being written exclusively in terms of boundary variables. Thus, no volume discretization is required, and the approach is ideally suited for geotechnical problems involving media of infinite extent. In the presentation, the required axisymmetric integral equations and kernel functions are first developed from the corresponding three-dimensional theory. In particular, emphasis is placed on the analytical and numerical treatment of the kernels. This is followed by an overview of the numerical implementation, and a demonstration of its merits via the consideration of several examples.