This paper presents analytical formulations for systematically deriving the solutions of Biot’s poroelasticity in saturated multilayered media of either full-space or halfspace extents. The number of the saturated multilayer media is either $n+2$ for full-space extent or $n+1$ for halfspace extent, where $n$ is a positive or zero integer. The applied loadings include the internal forces and liquid source for full-space and both internal and external loadings for halfspace region with eight cases of four boundary conditions. The mathematical tools for the formulations are classical and include the two-dimensional Fourier transform, the Hankel transform, Laplace transform as well as linear algebra. The solutions are expressed in matrix forms and each matrix is explicitly expressed with clear physical meaning and well-defined elements. The matrix solutions in the Fourier and Laplace transform domains are axially symmetric about the vertical axis. The internal and boundary conditions can be four-dimensional and the matrix solutions in the physical domain are also four-dimensional.
