The charge renormalization procedure for the calculation of the correlation energy of atoms utilizing the analytically known large-D limit solutions for the exact and Hartree-Fock equations is extended to diatomic molecules. This procedure is based on the variation of the nuclear charge, Z, and internuclear distance, R, of the Hartree-Fock Hamiltonian such that the Hartree-Fock energy will be significantly closer to the exact energy. We calculate to first order in delta Z the leading contribution to the correlation energy by changing the nuclear charge to some renormalized nuclear charge, Z(i)(R)-->Z(i)+delta Z(i). To first order in delta Z, this leads to an approximate expression, E(corr)(Z(a),Z(b),R)=(partial derivative E(HF)/partial derivative Z(a))delta Z(a) + (partial derivative E(HF)/partial derivative Z(b))delta Z(b), for the correlation energy based on the charge renormalization parameter delta Z, which is fixed systematically from the large-D limit. The theory is applied to diatomic molecules. Near the equilibrium, we are predicting the correlation energy typically with 80% or greater accuracy in a completely self-consistent and systematic way with no additional cost to the Hartree-Fock calculation. An improved approach to estimating the correlation energy for all R is outlined.