anchor.diff:    logical :: ok
     integer :: s1, s2, p1, p2, ib, istate
     integer(bit_kind) :: mask(N_int, 2), det(N_int, 2)
  -  double precision :: e_pert, delta_E, val, Hii, max_e_pert
  +  double precision :: e_pert, delta_E, val, Hii, max_e_pert, tmp
     double precision, external :: diag_H_mat_elem_fock
     
     
  
positive.code:     max_e_pert = 0d0
      
      do istate=1,N_states
        val = vect(istate, p1)
        delta_E = E0(istate) - Hii
        if (delta_E < 0.d0) then
          e_pert = 0.5d0 * (-dsqrt(delta_E * delta_E + 4.d0 * val * val) - delta_E)
        else
          e_pert = 0.5d0 * ( dsqrt(delta_E * delta_E + 4.d0 * val * val) - delta_E)
        endif
        pt2(istate) += e_pert
        if(dabs(e_pert) > dabs(max_e_pert)) max_e_pert = e_pert
      end do
  
negative.code:       
        do istate=1,N_states
          delta_E = E0(istate) - Hii
          val = mat(istate, p1, p2) + mat(istate, p1, p2) 
          tmp = dsqrt(delta_E * delta_E + val * val)
          if (delta_E < 0.d0) then
            tmp = -tmp
          endif
          e_pert = 0.5d0 * ( tmp - delta_E)
          pt2(istate) = pt2(istate) + e_pert
          max_e_pert = min(e_pert,max_e_pert)
        end do
        
        if(dabs(max_e_pert) > buf%mini) then
  