1
0
mirror of https://gitlab.com/scemama/qp_plugins_scemama.git synced 2024-12-23 04:43:38 +01:00
qp_plugins_scemama/devel/svdwf/TSVD_err.irp.f

117 lines
2.7 KiB
FortranFixed
Raw Normal View History

2021-11-02 16:18:07 +01:00
program TSVD_err
BEGIN_DOC
! study precision variation with truncated SVD
END_DOC
implicit none
read_wf = .True.
TOUCH read_wf
call run()
end
subroutine run
implicit none
integer :: mm, nn, i_state, low_rank, lrank_min, lrank_max
integer :: i, j, k, l
double precision :: norm_psi, err_verif, err_tmp
double precision, allocatable :: U_FSVD(:,:), D_FSVD(:), Vt_FSVD(:,:), A_FSVD(:,:)
double precision, allocatable :: US(:,:), A_TSVD(:,:)
double precision :: tmp_acc
double precision, allocatable :: D2_FSVD(:)
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
allocate( A_FSVD(mm,nn), U_FSVD(mm,mm), D_FSVD(min(mm,nn)), Vt_FSVD(nn,nn) )
norm_psi = 0.d0
A_FSVD(:,:) = 0.d0
do k = 1, N_det
i = psi_bilinear_matrix_rows(k)
j = psi_bilinear_matrix_columns(k)
A_FSVD(i,j) = psi_bilinear_matrix_values(k,i_state)
norm_psi += psi_bilinear_matrix_values(k,i_state) * psi_bilinear_matrix_values(k,i_state)
enddo
call svd_s( A_FSVD, size(A_FSVD,1), U_FSVD, size(U_FSVD,1), D_FSVD, Vt_FSVD, size(Vt_FSVD,1), mm, nn)
print *, ' --- Full SVD: ok after --- '
call CPU_TIME(t1)
print *, (t1-t0)/60.d0, 'min'
! ! US = U x SV :
! allocate( US(mm,min(mm,nn)) )
! US(:,:) = 0.d0
! do i = 1, mm
! do l = 1, min(mm,nn)
! US(i,l) = U_FSVD(i,l) * D_FSVD(l)
! enddo
! enddo
!
! allocate( A_TSVD(mm,nn) )
! A_TSVD(:,:) = 0.d0
! lrank_min = 1
! lrank_max = min(mm,nn)
! do low_rank = lrank_min, lrank_max, 1
!
! ! A_TSVD = US x Vt_FSVD
! call dgemm( 'N', 'N', mm, nn, low_rank, 1.d0 &
! , US , size(US ,1) &
! , Vt_FSVD, size(Vt_FSVD,1) &
! , 0.d0, A_TSVD, size(A_TSVD,1) )
! ! |A_FSVD-A_TSVD|
! err_verif = 0.d0
! do j = 1, nn
! do i = 1, mm
! err_tmp = A_FSVD(i,j) - A_TSVD(i,j)
! err_verif = err_verif + err_tmp * err_tmp
! enddo
! enddo
! print*, low_rank, 100.d0 * dsqrt(err_verif/norm_psi)
! enddo
! deallocate( A_TSVD , US )
deallocate( A_FSVD, U_FSVD, Vt_FSVD )
allocate( D2_FSVD(min(mm,nn)) )
do i = 1, min(mm,nn)
tmp_acc = 0.d0
do j = 1, i
tmp_acc = tmp_acc + D_FSVD(j) * D_FSVD(j)
enddo
D2_FSVD(i) = tmp_acc
enddo
open(unit=11, file='d_svd.txt',action='write')
do i = 1, min(mm,nn)
write(11, *) D_FSVD(i), D2_FSVD(i), 1.d0-D2_FSVD(i)
enddo
close(11)
deallocate( D2_FSVD )
deallocate( D_FSVD )
call CPU_TIME(t2)
print *, ''
print *, ' end after'
print *, (t2-t0)/60.d0, 'min'
end