The Heisenberg spin Hamiltonian, Hex, with spin-exchange coupling constants (Ji,j) between n effective spins (ESs) generated in a polynuclear transition-metal complex was derived from a multiple grand canonical (MGC) spin density functional theory combining a set of 2n-1 mutually independent ES arrangement (ESA) states defined by the principal GC internal energy functional (UUHFD) by using Dirac's spin-dependent electron-spin permutation operator. The variation principle for minimizing the grand potential in each ESA state yielded independently an orthonormal set of self-consistent natural LCAO-MO's. In all the ESA states, three sets of the expected values of UUHFD, Hex, and S² tot (Stot is the total spin operator) were combined to formulate n ES densities, and each set of the exchange-correlation density overlap integral ratios {SiES, jES} between iES and jES, and thereby the mean isotropic spin-exchange coupling constants {Ji,j}. The derived formulas were benchmark-tested and demonstrated the best quantitative agreement with 11 experimental results. [ABSTRACT FROM AUTHOR]