Variational Dirac-Coulomb explicitly correlated computations for atoms and molecules

The Dirac-Coulomb equation with positive-energy projection is solved using explicitly correlated Gaussian functions. The algorithm and computational procedure aims for a parts-per-billion convergence of the energy to provide a starting point for further comparison and further developments in relation with high-resolution atomic and molecular spectroscopy. Besides a detailed discussion of the implementation of the fundamental spinor structure, permutation and point-group symmetries, various options for the positive-energy projection procedure are presented. The no-pair Dirac-Coulomb energy converged to a parts-per-billion precision is compared with perturbative results for atomic and molecular systems with small nuclear charge numbers. The subsequent paper [Paper II: D. Ferenc, P. Jeszenszki, and E. M��tyus (2022)] describes the implementation of the Breit interaction in this framework.