A simplified scheme for relativistic density functional computations in the zeroth-order regular approximation (ZORA) to the Dirac equation is presented. The potential function in the kinetic energy operator is approximated by the potential generated from the superposition of the charge of the constituting atoms. The transition state method and ESA are adopted in bonding energy calculations to eliminate the gauge dependence error of the ZORA method. A two-step procedure is adopted to solve the ZORA equation: the scalar relativistic ZORA equation is first solved self-consistently, then the spin-orbit interaction is incorporated into the computation. The spin-orbit coupling matrix is made real by adopting the symmetry functions of irreducible representations of relevant double groups as basis sets and properly choosing their phase to avoid the complex arithmetic. The calculated results for several molecules show that this simplified scheme can be satisfactorily used for theoretical studies of compounds containing fairly heavy elements.