Complete benchmark rovibrational energy linelists calculated for the primordial polar molecules of the universe, namely HD + , HD, and the HeH + isotopologues, with accuracy up to 10 −2 cm −1 for low-lying states, are presented. To allow for these calculations to be performed, new high-accuracy potential energy curves, which include the diagonal Born-Oppenheimer adiabatic corrections and the leading relativistic corrections, are determined. Also, a new approach for calculating non-adiabatic corrections involving an effective vibrational nuclear mass obtained based on the atoms-in-molecules theory is employed. The vibrational and rotational masses are taken as being different and dependent on the nuclear distance. Accurate dipole moment curves are calculated and used to generate lists of Einstein A-coefficients. The energy linelists and the sets of Einstein A-coefficients for HD are upgrades of previous calculations including quasibound states, while for HD + and HeH + and its isotopologues the present results represent significant improvement over the previous calculations. The results obtained here suggest that, with the inclusion of the non-adiabatic corrections, the accuracy limit at least for low-lying states might have been reached. Thus, further progress should involve accounting for even smaller effects such as the quantumelectrodynamics corrections. The present results represent the state-of-the-art of theoretical spectroscopy of the primordial polar molecules.