We present a new implementation for computing spin–orbit couplings (SOCs) within a time-dependent density-functional theory (TD-DFT) framework in the standard spin-conserving formulation as well in the spin–flip variant (SF-TD-DFT). This approach employs the Breit–Pauli Hamiltonian and Wigner–Eckart’s theorem applied to the reduced one-particle transition density matrices, together with the spin–orbit mean-field treatment of the two-electron contributions. We use a state-interaction procedure and compute the SOC matrix elements using zero-order non-relativistic states. Benchmark calculations using several closed-shell organic molecules, diradicals, and a single-molecule magnet illustrate the efficiency of the SOC protocol. The results for organic molecules (described by standard TD-DFT) show that SOCs are insensitive to the choice of the functional or basis sets, as long as the states of the same characters are compared. In contrast, the SF-TD-DFT results for small diradicals (CH 2 , [Formula: see text], SiH 2 , and [Formula: see text]) show strong functional dependence. The spin-reversal energy barrier in a Fe(III) single-molecule magnet computed using non-collinear SF-TD-DFT (PBE0, ωPBEh/cc-pVDZ) agrees well with the experimental estimate.