diff --git a/src/ao_one_e_ints/pseudopot.f90 b/src/ao_one_e_ints/pseudopot.f90
index aaa456ba..6c8e5c83 100644
--- a/src/ao_one_e_ints/pseudopot.f90
+++ b/src/ao_one_e_ints/pseudopot.f90
@@ -1862,6 +1862,7 @@ double precision function int_prod_bessel(l,gam,n,m,a,b,arg)
       qk = dble(q)
       two_qkmp1 = 2.d0*(qk+mk)+1.d0
       do k=0,q-1
+        ! possible FPE here. To be checked
         s_q_k = two_qkmp1*qk*inverses(k)*s_q_k
 !        if (s_q_k < 1.d-32) then
 !          s_q_k = 0.d0
diff --git a/src/cipsi/selection.irp.f b/src/cipsi/selection.irp.f
index 071d9294..c325020f 100644
--- a/src/cipsi/selection.irp.f
+++ b/src/cipsi/selection.irp.f
@@ -780,24 +780,23 @@ subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_d
         endif
       enddo
 
-      ! Gram-Schmidt using input overlap matrix
-      do istate=1,N_states
-        do jstate=1,istate-1
-          if ( (pt2_overlap(jstate,istate) == 0.d0).or.(pt2_overlap(jstate,jstate) == 0.d0) ) cycle
-          coef(istate) = coef(istate) - pt2_overlap(jstate,istate)/pt2_overlap(jstate,jstate) * coef(jstate)
-        enddo
-      enddo
+!      ! Gram-Schmidt using input overlap matrix
+!      do istate=1,N_states
+!        do jstate=1,istate-1
+!          if ( (pt2_overlap(jstate,istate) == 0.d0).or.(pt2_overlap(jstate,jstate) == 0.d0) ) cycle
+!          coef(istate) = coef(istate) - pt2_overlap(jstate,istate)/pt2_overlap(jstate,jstate) * coef(jstate)
+!        enddo
+!      enddo
 
       do istate=1, N_states
-        do jstate=1,N_states
-          pt2_data % overlap(jstate,istate) += coef(jstate) * coef(istate)
-        enddo
-      enddo
 
-      do istate=1,N_states
         alpha_h_psi = mat(istate, p1, p2)
         e_pert = coef(istate) * alpha_h_psi
 
+        do jstate=1,N_states
+          pt2_data % overlap(jstate,istate) += coef(jstate) * alpha_h_psi
+        enddo
+
         pt2_data % variance(istate) += alpha_h_psi * alpha_h_psi
         pt2_data % pt2(istate)      += e_pert