In this work, the two-center Dirac equation is solvednumerically using an extension of an adapted B-spline basis set methodpreviously implemented in relativistic atomic calculations (Fischer, C. F.;Zatsarinny, O. Comput. Phys. Commun. 2009, 180, 879). The robustnessof the chosen numerical method, which avoids the appearance ofspurious states common in other approaches, allows us to investigatemolecular photoionization within a relativistic framework by simplyadapting those methods already available in the nonrelativistic case(Brosolo, M.; Decleva, P. Chem. Phys. 1992, 159, 185; Brosolo, M.;Decleva, P.; Lisini, A. Mol. Opt. Phys. 1992, 25, 3345). First, lightdiatomic molecules (i.e., H2+ and HeH2+) are investigated with thepurpose of testing the validity and efficiency of the method. Then, a seriesof one-electron molecular hydrides (i.e., HF9+, HCl17+ and HI53+) isexplored by computing the total photoionization cross sections, asymmetry β-parameters and partial phase shifts. The presentmethodology can be easily extended to treat N-electron molecules following previous approaches in nonrelativistic calculations(Plesiat, E.; Decleva, P.; Martin, F. Phys. Chem. Chem. Phys. 2012, 14, 10853). The inclusion of a second photon can be alsoaccomplished just like in atomic investigations aiming at reproducing pump−probe experiments capable to extract thephotoionization time-delays (Vinbladh, J.; Dahlstrom, J. M.; Lindroth, E. Phys. Rev A 2019, 100, 043424; Vinblach, J.; Dahlstrom, J.M.; Lindroth, E. Atoms 2022, 10, 80).