diff --git a/src/density.irp.f b/src/density.irp.f
index 9361f57..45ed3e0 100644
--- a/src/density.irp.f
+++ b/src/density.irp.f
@@ -94,6 +94,7 @@ BEGIN_PROVIDER [ real, density_beta_value_p ]
  integer :: exc(4)
  PROVIDE det
  PROVIDE elec_beta_num
+
  do k=1,det_num
   do l=1,det_num
    exc(1) = abs(det_exc(k,l,1))
diff --git a/src/det.irp.f b/src/det.irp.f
index 7248cc8..07a4b52 100644
--- a/src/det.irp.f
+++ b/src/det.irp.f
@@ -33,7 +33,7 @@ BEGIN_PROVIDER  [ integer, det, (elec_alpha_num-mo_closed_num,det_num,2) ]
 END_PROVIDER
 
 
-BEGIN_PROVIDER [ integer, det_exc, (det_num, det_num, 2) ]
+BEGIN_PROVIDER [ integer*1, det_exc, (det_num, det_num, 2) ]
  implicit none
  BEGIN_DOC  
 ! Degree of excitation between two determinants. Indices are alpha, beta
@@ -141,7 +141,7 @@ subroutine get_single_excitation(k,l,m,n,p)
  do i=1,elec_num_2(p)-mo_closed_num
   found = .False.
   do j=1,elec_num_2(p)-mo_closed_num
-   if (det(i,k,p) == det(i,l,p)) then
+   if (det(j,k,p) == det(i,l,p)) then
     found = .True.
     exit
    endif
@@ -188,7 +188,7 @@ subroutine get_double_excitation(k,l,m,n,r,s,p)
  do i=1,elec_num_2(p)-mo_closed_num
   found = .False.
   do j=1,elec_num_2(p)-mo_closed_num
-   if (det(i,k,p) == det(i,l,p)) then
+   if (det(j,k,p) == det(i,l,p)) then
     found = .True.
     exit
    endif
diff --git a/src/eplf_function.irp.f b/src/eplf_function.irp.f
index 61b37b8..905f19e 100644
--- a/src/eplf_function.irp.f
+++ b/src/eplf_function.irp.f
@@ -38,7 +38,7 @@ BEGIN_PROVIDER [ double precision, mo_eplf_integral_matrix, (mo_num,mo_num) ]
  integer :: i, j, k, l
  double precision :: t
  PROVIDE ao_eplf_integral_matrix
- PROVIDE mo_coef
+ PROVIDE mo_coef mo_coef_transp
  do i=1,mo_num
   do j=i,mo_num
     mo_eplf_integral_matrix(j,i) = 0.d0
diff --git a/src/mo_point.irp.f b/src/mo_point.irp.f
index b36c1c0..cc80ce8 100644
--- a/src/mo_point.irp.f
+++ b/src/mo_point.irp.f
@@ -16,6 +16,23 @@ BEGIN_PROVIDER [ real, mo_value_p, (mo_num) ]
 
 END_PROVIDER
 
+BEGIN_PROVIDER [ real, mo_value_prod_p, (mo_num,mo_num) ]
+ implicit none
+
+ BEGIN_DOC
+! Products of 2 molecular orbitals
+ END_DOC
+
+ integer :: i, j, k
+
+ do j=1,mo_num
+   do k=1,mo_num
+    mo_value_prod_p(k,j) = mo_value_p(k)*mo_value_p(j)
+   enddo
+ enddo
+
+END_PROVIDER
+
 BEGIN_PROVIDER [ real, mo_grad_p, (mo_num,3) ]
  implicit none