Equations are proposed for computing from ab initio wave functions quantities like [S-A.S-B], which appear in the Heisenberg model Hamiltonian of magnetism. These equations are based on projection operators derived from Lowdin orthogonalization. They result in local spin operators S-A which obey the definition of angular momentum operators and commute with each other. These equations are evaluated for several typical closed and open shell molecules. For closed shells in the single Slater determinant approximation, [S-A.S-B] is -3/8 of the bond-order and [S(A)2] is +3/8 of the total number of bonds to center A. For open shells there are additional contributions from the unpaired electrons. In favorable cases, these additional terms have the value assumed as the whole answer in the usual applications of the Heisenberg Hamiltonian.