Relativistic spin-orbit density functional theory (DFT) methods have been implemented in the molecular Gaussian DFT and pseudopotential planewave DFT modules of the NWChem electronic-structure program. The Gaussian basis set implementation is based upon the zeroth-order regular approximation (ZORA) while the planewave implementation uses spin-orbit pseudopotentials that are directly generated from the atomic Dirac-Kohn-Sham wave functions or atomic ZORA-Kohn-Sham wave functions. Compared to solving the full Dirac equation these methods are computationally efficient but robust enough for a realistic description of relativistic effects such as spin-orbit splitting, molecular orbital hybridization, and core effects. Both methods have been applied to a variety of small molecules, including I2, IF, HI, Br2, Bi2, AuH, and Au2, using various exchange-correlation functionals. Our results are in good agreement with experiment and previously reported calculations.