mirror of
https://gitlab.com/scemama/qp_plugins_scemama.git
synced 2024-06-02 03:15:25 +02:00
200 lines
5.0 KiB
Fortran
200 lines
5.0 KiB
Fortran
|
|
program invFSVD_det
|
|
|
|
BEGIN_DOC
|
|
! build CI matrix from SVD vectors
|
|
END_DOC
|
|
|
|
implicit none
|
|
|
|
read_wf = .True.
|
|
TOUCH read_wf
|
|
|
|
call run()
|
|
|
|
end
|
|
|
|
|
|
subroutine run
|
|
|
|
implicit none
|
|
|
|
integer :: mm, nn, i_state, low_rank
|
|
integer :: i, j, k, l
|
|
double precision :: tmp1, eps_det
|
|
double precision, allocatable :: U_SVD(:,:), D_SVD(:), V_SVD(:,:), A_SVD(:,:), US(:,:)
|
|
double precision, allocatable :: newpsi_rows(:), newpsi_columns(:), newpsi_values(:)
|
|
double precision :: t0, t1, t2
|
|
|
|
call CPU_TIME(t0)
|
|
|
|
i_state = 1
|
|
mm = n_det_alpha_unique
|
|
nn = n_det_beta_unique
|
|
print *, ' matrix dimensions:', mm,'x',nn
|
|
print *, ' N det:', N_det
|
|
|
|
call ezfio_get_spindeterminants_n_svd_coefs(low_rank)
|
|
print *, ' SVD rank:', low_rank
|
|
|
|
allocate( U_SVD(mm,low_rank), D_SVD(low_rank), V_SVD(nn,low_rank) )
|
|
call ezfio_get_spindeterminants_psi_svd_alpha(U_SVD)
|
|
call ezfio_get_spindeterminants_psi_svd_beta (V_SVD)
|
|
call ezfio_get_spindeterminants_psi_svd_coefs(D_SVD)
|
|
print *, ' read EZFIO SVD vectors OK'
|
|
|
|
! US = U x S :
|
|
allocate( US(mm,low_rank) )
|
|
US(:,:) = 0.d0
|
|
do i = 1, mm
|
|
do l = 1, low_rank
|
|
US(i,l) = U_SVD(i,l) * D_SVD(l)
|
|
enddo
|
|
enddo
|
|
print *, ' U x D: ok'
|
|
deallocate( U_SVD , D_SVD )
|
|
|
|
! A_TSVD = US x V.T
|
|
allocate( A_SVD(mm,nn) )
|
|
A_SVD(:,:) = 0.d0
|
|
call dgemm( 'N', 'T', mm, nn, low_rank, 1.d0 &
|
|
, US , size(US ,1) &
|
|
, V_SVD , size(V_SVD,1) &
|
|
, 0.d0, A_SVD, size(A_SVD,1) )
|
|
print *, ' U x D x Vt: ok'
|
|
deallocatE(US)
|
|
|
|
call ezfio_set_spindeterminants_n_states(N_states)
|
|
call ezfio_set_spindeterminants_n_det_alpha(n_det_alpha_unique)
|
|
call ezfio_set_spindeterminants_n_det_beta(n_det_beta_unique)
|
|
|
|
print *, 'EZFIO set n_det_alpha and n_det_beta'
|
|
|
|
! ---------------
|
|
eps_det = 1d-13
|
|
print *, ' det coeff thresh for consideration =', eps_det
|
|
! ---------------
|
|
|
|
|
|
k = 0
|
|
do j = 1, nn
|
|
do i = 1, mm
|
|
tmp1 = A_SVD(i,j)
|
|
if( dabs(tmp1) .lt. (eps_det) ) cycle
|
|
k = k + 1
|
|
enddo
|
|
enddo
|
|
print *, ' non zero elements = ', k
|
|
|
|
|
|
! call generate_all_alpha_beta_det_products
|
|
! call save_wavefunction
|
|
|
|
call all_ab_det(eps_det, A_SVD)
|
|
call update_SVDWF(eps_det, k, A_SVD)
|
|
|
|
deallocate( A_SVD )
|
|
|
|
|
|
call CPU_TIME(t2)
|
|
print *, ''
|
|
print *, ' end after'
|
|
print *, (t2-t0)/60.d0, 'min'
|
|
|
|
end
|
|
|
|
|
|
|
|
subroutine update_SVDWF(eps_det, n_ab, A_SVD)
|
|
|
|
implicit none
|
|
|
|
integer , intent(in) :: n_ab
|
|
double precision, intent(in) :: eps_det
|
|
double precision, intent(in) :: A_SVD(n_det_alpha_unique,n_det_beta_unique)
|
|
|
|
integer :: i, j, l
|
|
integer(bit_kind), allocatable :: tmp_det(:,:,:)
|
|
double precision, allocatable :: tmp_coef(:)
|
|
|
|
|
|
PROVIDE H_apply_buffer_allocated
|
|
|
|
allocate( tmp_coef(n_ab) )
|
|
allocate (tmp_det(N_int,2,n_ab))
|
|
|
|
l = 0
|
|
do j = 1, N_det_beta_unique
|
|
do i = 1, N_det_alpha_unique
|
|
|
|
if( (dabs(A_SVD(i,j)) .ge. eps_det) .and. (l.le.n_ab) ) then
|
|
l = l + 1
|
|
tmp_det(1:N_int,1,l) = psi_det_alpha_unique(1:N_int,i)
|
|
tmp_det(1:N_int,2,l) = psi_det_beta_unique (1:N_int,j)
|
|
tmp_coef(l) = A_SVD(i,j)
|
|
endif
|
|
|
|
enddo
|
|
enddo
|
|
print *, ' l_max =', l
|
|
print *, ' n_ab =', n_ab
|
|
|
|
!call save_wavefunction_general(ndet,nstates ,psidet ,dim_psicoef,psicoef )
|
|
call save_wavefunction_general(l ,n_states,tmp_det,l ,tmp_coef)
|
|
|
|
deallocate(tmp_coef)
|
|
deallocate(tmp_det)
|
|
|
|
end subroutine update_SVDWF
|
|
|
|
|
|
|
|
|
|
|
|
subroutine all_ab_det(eps_det, A_SVD)
|
|
|
|
implicit none
|
|
|
|
double precision, intent(in) :: eps_det
|
|
double precision, intent(in) :: A_SVD(n_det_alpha_unique,n_det_beta_unique)
|
|
|
|
integer :: i, j, k, l
|
|
integer :: iproc
|
|
integer(bit_kind), allocatable :: tmp_det(:,:,:)
|
|
|
|
integer, external :: get_index_in_psi_det_sorted_bit
|
|
logical, external :: is_in_wavefunction
|
|
|
|
PROVIDE H_apply_buffer_allocated
|
|
|
|
!$OMP PARALLEL DEFAULT(NONE) SHARED(psi_coef_sorted_bit, N_det_beta_unique, &
|
|
!$OMP N_det_alpha_unique, N_int, psi_det_alpha_unique, psi_det_beta_unique, &
|
|
!$OMP eps_det, A_SVD ) &
|
|
!$OMP PRIVATE(i,j,k,l,tmp_det,iproc)
|
|
!$ iproc = omp_get_thread_num()
|
|
allocate (tmp_det(N_int,2,N_det_alpha_unique))
|
|
!$OMP DO SCHEDULE(static,8)
|
|
do j=1,N_det_beta_unique
|
|
l = 1
|
|
do i=1,N_det_alpha_unique
|
|
do k=1,N_int
|
|
tmp_det(k,1,l) = psi_det_alpha_unique(k,i)
|
|
tmp_det(k,2,l) = psi_det_beta_unique (k,j)
|
|
enddo
|
|
if( (.not.is_in_wavefunction(tmp_det(1,1,l),N_int)) .and. (dabs(A_SVD(i,j)) .ge. eps_det) ) then
|
|
l = l + 1
|
|
endif
|
|
enddo
|
|
call fill_H_apply_buffer_no_selection(l-1, tmp_det, N_int, iproc)
|
|
enddo
|
|
!$OMP END DO
|
|
deallocate(tmp_det)
|
|
!$OMP END PARALLEL
|
|
|
|
call copy_H_apply_buffer_to_wf
|
|
SOFT_TOUCH psi_det psi_coef N_det N_det_beta_unique N_det_alpha_unique psi_det_alpha_unique psi_det_beta_unique
|
|
|
|
print *, ' new n_det =', N_det
|
|
|
|
end subroutine all_ab_det()
|