We have developed a higher order perfectly matched layer (PML) formulation to improve the absorption performance for finite-difference time-domain seismic modeling. First, we outlined a new unsplit “correction” approach, which allowed for traditional, first-order PMLs to be added directly to existing codes in a straightforward manner. Then, using this framework, we constructed a PML formulation that can be used to construct higher order PMLs of arbitrary order. The greater number of degrees of freedom associated with the higher order PML allow for enhanced flexibility of the PML stretching functions, thus potentially facilitating enhanced absorption performance. We found that the new approach can offer increased elastodynamic absorption, particularly for evanescent waves. We also discovered that the extra degrees of freedom associated with the higher order PML required careful optimization if enhanced absorption was to be achieved. Furthermore, these extra degrees of freedom increased the computational requirements in comparison with first-order schemes. We reached our formulations using one compact equation thus increasing the ease of implementation. Additionally, the formulations are based on a recursive integration approach that reduce PML memory requirements, and do not require special consideration for corner regions. We tested the new formulations to determine their ability to absorb body waves and surface waves. We also tested standard staggered grid stencils and rotated staggered grid stencils.