In this paper we will calculate the effect of spin-orbit coupling on properties of closed shell molecules, using the zero-order regular approximation to the Dirac equation. Results are obtained using density functionals including density gradient corrections. Close agreement with experiment is obtained for the calculated molecular properties of a number of heavy element diatomic molecules.