PTEROSOR/QuantumPackage
10
0
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-11-14 18:13:51 +01:00
Code Releases Activity

adapted the tc_natorb with block diagonalization

Browse Source
This commit is contained in:
eginer 2022-10-29 20:07:33 +02:00
parent c7dd091b63
commit 1e2f8ab44f
5 changed files with 33 additions and 15 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

19
src/non_hermit_dav/lapack_diag_non_hermit.irp.f
View File

@ -930,11 +930,15 @@ subroutine check_EIGVEC(n, m, A, eigval, leigvec, reigvec, thr_diag, thr_norm, s
tmp_abs = tmp_abs + tmp
V_nrm = V_nrm + U_nrm
print *, j, tmp, U_nrm
write(*,'(I4,X,(100(F25.16,X)))')j,eigval(j), tmp, U_nrm
enddo
tmp_rel = tmp_abs / tmp_nrm
if(tmp_abs.lt.10.d-10)then
tmp_rel = thr_diag/10.d0
else
tmp_rel = tmp_abs / tmp_nrm
endif
tmp_dif = dabs(V_nrm - dble(m))
if( stop_ifnot .and. ((tmp_rel .gt. thr_diag) .or. (tmp_dif .gt. thr_norm)) ) then
@ -969,13 +973,20 @@ subroutine check_EIGVEC(n, m, A, eigval, leigvec, reigvec, thr_diag, thr_norm, s
tmp_abs = tmp_abs + tmp
V_nrm = V_nrm + U_nrm
print *, j, tmp, U_nrm
write(*,'(I4,X,(100(F25.16,X)))')j,eigval(j), tmp, U_nrm
enddo
tmp_rel = tmp_abs / tmp_nrm
if(tmp_abs.lt.10.d-10)then
tmp_rel = thr_diag/10.d0
else
tmp_rel = tmp_abs / tmp_nrm
endif
if( stop_ifnot .and. ((tmp_rel .gt. thr_diag) .or. (tmp_dif .gt. thr_norm)) ) then
print *, ' error in left-eigenvectors'
print *, ' err tol = ',thr_diag, thr_norm
print *, '(tmp_rel .gt. thr_diag) = ',(tmp_rel .gt. thr_diag)
print *, '(tmp_dif .gt. thr_norm) = ',(tmp_dif .gt. thr_norm)
print *, ' err estim = ', tmp_abs, tmp_rel
print *, ' CR norm = ', V_nrm
stop

16
src/tc_bi_ortho/tc_natorb.irp.f
View File

@ -8,17 +8,23 @@
! natorb_tc_leigvec_mo : LEFT eigenvectors of the ground state transition matrix (equivalent of natural orbitals)
! natorb_tc_eigval : eigenvalues of the ground state transition matrix (equivalent of the occupation numbers). WARNINING :: can be negative !!
END_DOC
double precision, allocatable :: dm_tmp(:,:)
double precision, allocatable :: dm_tmp(:,:),fock_diag(:)
double precision :: thr_deg
integer :: i,j,k,n_real
allocate( dm_tmp(mo_num,mo_num))
allocate( dm_tmp(mo_num,mo_num),fock_diag(mo_num))
dm_tmp(:,:) = -tc_transition_matrix(:,:,1,1)
print*,'dm_tmp'
do i = 1, mo_num
fock_diag(i) = fock_matrix_tc_mo_tot(i,i)
write(*,'(100(F16.10,X))')-dm_tmp(:,i)
enddo
call non_hrmt_bieig( mo_num, dm_tmp&
, natorb_tc_leigvec_mo, natorb_tc_reigvec_mo&
, n_real, natorb_tc_eigval )
thr_deg = thr_degen_tc
call diag_mat_per_fock_degen(fock_diag,dm_tmp,mo_num,thr_deg,&
natorb_tc_leigvec_mo,natorb_tc_reigvec_mo,&
natorb_tc_eigval)
! call non_hrmt_bieig( mo_num, dm_tmp&
! , natorb_tc_leigvec_mo, natorb_tc_reigvec_mo&
! , n_real, natorb_tc_eigval )
double precision :: accu
accu = 0.d0
do i = 1, n_real

2
src/tc_scf/minimize_tc_angles.irp.f
View File

@ -6,7 +6,7 @@ program print_angles
! my_n_pt_r_grid = 10 ! small grid for quick debug
! my_n_pt_a_grid = 14 ! small grid for quick debug
touch my_n_pt_r_grid my_n_pt_a_grid
call sort_by_tc_fock
! call sort_by_tc_fock
call minimize_tc_orb_angles
end

1
src/tools/NEED
View File

@ -2,3 +2,4 @@ fci
mo_two_e_erf_ints
aux_quantities
hartree_fock
dft_utils_in_r

10
src/utils/block_diag_degen.irp.f
View File

@ -29,10 +29,10 @@ subroutine diag_mat_per_fock_degen(fock_diag,mat_ref,n,thr_deg,leigvec,reigvec,e
reigvec_unsrtd = 0.d0
eigval_unsrtd = 0.d0
allocate(list_degen(mo_num,0:mo_num))
allocate(list_degen(n,0:n))
! obtain degeneracies
call give_degen_full_list(fock_diag,mo_num,thr_deg,list_degen,n_degen_list)
call give_degen_full_list(fock_diag,n,thr_deg,list_degen,n_degen_list)
allocate(iorder(n_degen_list),list_degen_sorted(n_degen_list))
do i =1, n_degen_list
n_degen = list_degen(i,0)
@ -102,9 +102,9 @@ subroutine diag_mat_per_fock_degen(fock_diag,mat_ref,n,thr_deg,leigvec,reigvec,e
deallocate(mat_tmp,list_same_degen)
deallocate(eigval_tmp,leigvec_tmp,reigvec_tmp)
enddo
if(icount.ne.n)then
print*,'pb !! (icount.ne.n)'
print*,'icount,n',icount,n
if(icount_eigval.ne.n)then
print*,'pb !! (icount_eigval.ne.n)'
print*,'icount_eigval,n',icount_eigval,n
stop
endif