From 3b3f7a7de9209a22868058f6f11168aee7be4235 Mon Sep 17 00:00:00 2001
From: Anthony Scemama <scemama@irsamc.ups-tlse.fr>
Date: Sun, 5 Dec 2021 22:45:21 +0100
Subject: [PATCH] Added truncate_wf

---
 RELEASE_NOTES.org           |  1 +
 src/tools/truncate_wf.irp.f | 98 +++++++++++++++++++++++++++++++++++++
 2 files changed, 99 insertions(+)
 create mode 100644 src/tools/truncate_wf.irp.f

diff --git a/RELEASE_NOTES.org b/RELEASE_NOTES.org
index 98830f3f..586bf431 100644
--- a/RELEASE_NOTES.org
+++ b/RELEASE_NOTES.org
@@ -58,6 +58,7 @@
       symmetry in matrices
     - qp_export_as_tgz exports also plugin codes
     - Added a basis module containing basis set information
+    - Added qp_run truncate_wf
 
 *** Code
 
diff --git a/src/tools/truncate_wf.irp.f b/src/tools/truncate_wf.irp.f
new file mode 100644
index 00000000..6c66c8ec
--- /dev/null
+++ b/src/tools/truncate_wf.irp.f
@@ -0,0 +1,98 @@
+program truncate_wf
+ implicit none
+ BEGIN_DOC
+! Truncate the wave function
+ END_DOC
+ read_wf = .True.
+ if (s2_eig) then
+  call routine_s2
+ else
+  call routine
+ endif
+end
+
+subroutine routine
+ implicit none
+ integer :: ndet_max
+ print*, 'Max number of determinants ?'
+ read(5,*) ndet_max
+ integer(bit_kind), allocatable :: psi_det_tmp(:,:,:)
+ double precision, allocatable :: psi_coef_tmp(:,:)
+ allocate(psi_det_tmp(N_int,2,ndet_max),psi_coef_tmp(ndet_max, N_states))
+ 
+ integer :: i,j
+ double precision :: accu(N_states)
+ accu = 0.d0
+ do i = 1, ndet_max
+  do j = 1, N_int
+   psi_det_tmp(j,1,i) = psi_det_sorted(j,1,i)
+   psi_det_tmp(j,2,i) = psi_det_sorted(j,2,i)
+  enddo
+  do j = 1, N_states
+   psi_coef_tmp(i,j) = psi_coef_sorted(i,j)
+   accu(j) += psi_coef_tmp(i,j) **2
+  enddo
+ enddo
+ do j = 1, N_states
+  accu(j) = 1.d0/dsqrt(accu(j))
+ enddo
+ do j = 1, N_states
+  do i = 1, ndet_max
+   psi_coef_tmp(i,j) = psi_coef_tmp(i,j) * accu(j)
+  enddo
+ enddo
+
+ call save_wavefunction_general(ndet_max,N_states,psi_det_tmp,size(psi_coef_tmp,1),psi_coef_tmp)
+ 
+end
+
+subroutine routine_s2
+ implicit none
+ integer :: ndet_max
+ double precision :: wmin
+ integer(bit_kind), allocatable :: psi_det_tmp(:,:,:)
+ double precision, allocatable :: psi_coef_tmp(:,:)
+ integer :: i,j,k
+ double precision :: accu(N_states)
+
+ print *,  'Weights of the CFG'
+ do i=1,N_det
+   print *,  i, real(weight_configuration(det_to_configuration(i),:)), real(sum(weight_configuration(det_to_configuration(i),:)))
+ enddo
+ print*, 'Min weight of the configuration?'
+ read(5,*) wmin
+
+ ndet_max = 0
+ do i=1,N_det
+   if (maxval(weight_configuration( det_to_configuration(i),:)) < wmin) cycle
+   ndet_max = ndet_max+1
+ enddo
+
+ allocate(psi_det_tmp(N_int,2,ndet_max),psi_coef_tmp(ndet_max, N_states))
+ 
+ accu = 0.d0
+ k=0
+ do i = 1, N_det
+  if (maxval(weight_configuration( det_to_configuration(i),:)) < wmin) cycle
+  k = k+1
+  do j = 1, N_int
+   psi_det_tmp(j,1,k) = psi_det(j,1,i)
+   psi_det_tmp(j,2,k) = psi_det(j,2,i)
+  enddo
+  do j = 1, N_states
+   psi_coef_tmp(k,j) = psi_coef(i,j)
+   accu(j) += psi_coef_tmp(k,j)**2
+  enddo
+ enddo
+ do j = 1, N_states
+  accu(j) = 1.d0/dsqrt(accu(j))
+ enddo
+ do j = 1, N_states
+  do i = 1, ndet_max
+   psi_coef_tmp(i,j) = psi_coef_tmp(i,j) * accu(j)
+  enddo
+ enddo
+
+ call save_wavefunction_general(ndet_max,N_states,psi_det_tmp,size(psi_coef_tmp,1),psi_coef_tmp)
+ 
+end