anchor.diff:  
   !--- Excitation generation probabilities ---
   
  -    pure function calc_pgen_ueg_no_renorm(sys, max_na) result(pgen)
  +    pure function calc_pgen_ueg_no_renorm(sys, max_na, ij_spin) result(pgen)
   
           ! NOTE: this assumes that the corresponding excitation (i,j)->(a,b) is
  -        ! actually allowed.  Cases where this is not true must be explcitly
  +        ! actually allowed.  Cases where this is not true must be explicitly
           ! handled.
   
           ! In:
  
positive.code: 
              ! 3. Return generation probability and conecting matrix element.
  
              pgen = calc_pgen_ueg_no_renorm(sys, max_na)
  
              call find_excitation_permutation2(cdet%f, connection)
              hmatel = slater_condon2_ueg_excit(sys, connection%from_orb(1), connection%from_orb(2), &
  
negative.code:             ! Ok, b is now completely defined.
              kb = ij_k - basis_fns(a)%l
  
              ! Implementation detail: b must be down spin if possible as a is up spin if possible.
              select case(ij_spin)
              case(2)
                  b = ueg_basis_index(kb, 1)
              case(0)
                  b = ueg_basis_index(kb, -1)
              case(-2)
                  b = ueg_basis_index(kb, -1)
              end select
  
              ! Excitation is forbidden if b is already occupied!
  