For each monoanionic phosphorylated amino acid, we performed two TorsionDrive calculations that scanned the backbone (φ, ψ) and side chain (χ 1 , χ 2 ) dihedral angles respectively, each on a 2-D grid with 15 • resolution. We observed that some of these minimizations led to proton transfer from a backbone N-H group to the phosphate, and rectified this by constraining the two N-H bond lengths to 1.01 Å. These minimizations produce two 24×24 grids with an associated geometry at each point. Because this procedure involves several times as many energy minimizations as grid points (i.e. thousands per dipeptide), they are carried out using a relatively inexpensive level of theory, B97-D3/6-31G*. The lowest-energy structures on each grid point are then re-minimized with their constraints at the RI-MP2/augcc-pVDZ level of theory, 57 followed by a RI-MP2/aug-ccPVTZ gradient calculation and a RI-MP2/aug-cc-pV(T,Q)Z single-point energy calculation to estimate the MP2/CBS energy using Helgaker's two-point extrapolation. Because the TorsionDrive procedure was highly costly, we took a less expensive approach for the dianionic residues by removing the hydrogen on the phosphate group from each structure in the monoanionic grids and running the constrained minimizations independently; single-point energy and gradient calculations were then performed following the same procedure as the monoanionic case.