PTEROSOR/QuantumPackage
10
0
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-11-03 12:43:48 +01:00
Code Releases Activity

Merge branch 'dev-stable' of github.com:QuantumPackage/qp2 into dev-stable

Browse Source
This commit is contained in:
Anthony Scemama 2023-10-13 15:10:55 +02:00
commit 4257185f12
4 changed files with 176 additions and 9 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
src/fci/fci.irp.f
View File

@ -42,7 +42,7 @@ program fci
write(json_unit,json_array_open_fmt) 'fci'
double precision, allocatable :: Ev(:),PT2(:)
allocate(Ev(N_states), PT2(N_state))
allocate(Ev(N_states), PT2(N_states))
if (do_pt2) then
call run_stochastic_cipsi(Ev,PT2)
else

82
src/non_hermit_dav/biorthog.irp.f
View File

@ -270,7 +270,7 @@ subroutine non_hrmt_bieig(n, A, thr_d, thr_nd, leigvec, reigvec, n_real_eigv, ei
integer, intent(out) :: n_real_eigv
double precision, intent(out) :: reigvec(n,n), leigvec(n,n), eigval(n)
integer :: i, j
integer :: i, j,k
integer :: n_good
double precision :: thr, thr_cut, thr_diag, thr_norm
double precision :: accu_d, accu_nd
@ -278,6 +278,8 @@ subroutine non_hrmt_bieig(n, A, thr_d, thr_nd, leigvec, reigvec, n_real_eigv, ei
integer, allocatable :: list_good(:), iorder(:)
double precision, allocatable :: WR(:), WI(:), VL(:,:), VR(:,:)
double precision, allocatable :: S(:,:)
double precision, allocatable :: phi_1_tilde(:),phi_2_tilde(:),chi_1_tilde(:),chi_2_tilde(:)
allocate(phi_1_tilde(n),phi_2_tilde(n),chi_1_tilde(n),chi_2_tilde(n))
! -------------------------------------------------------------------------------------
@ -301,11 +303,78 @@ subroutine non_hrmt_bieig(n, A, thr_d, thr_nd, leigvec, reigvec, n_real_eigv, ei
call lapack_diag_non_sym(n, A, WR, WI, VL, VR)
!call lapack_diag_non_sym_new(n, A, WR, WI, VL, VR)
!print *, ' '
!print *, ' eigenvalues'
!do i = 1, n
! write(*, '(1000(F16.10,X))') WR(i), WI(i)
!enddo
print *, ' '
print *, ' eigenvalues'
i = 1
do while(i .le. n)
write(*, '(I3,X,1000(F16.10,X))')i, WR(i), WI(i)
if(.false.)then
if(WI(i).ne.0.d0)then
print*,'*****************'
print*,'WARNING ! IMAGINARY EIGENVALUES !!!'
write(*, '(1000(F16.10,X))') WR(i), WI(i+1)
! phi = VR(:,i), psi = VR(:,i+1), |Phi_i> = phi + j psi , |Phi_i+1> = phi - j psi
! chi = VL(:,i), xhi = VL(:,i+1), |Chi_i> = chi + j xhi , |Chi_i+1> = chi - j xhi
!
accu_chi_phi = 0.d0
accu_xhi_psi = 0.d0
accu_chi_psi = 0.d0
accu_xhi_phi = 0.d0
double precision :: accu_chi_phi, accu_xhi_psi, accu_chi_psi, accu_xhi_phi
double precision :: mat_ovlp(2,2),eigval_tmp(2),eigvec(2,2),mat_ovlp_orig(2,2)
do j = 1, n
accu_chi_phi += VL(j,i) * VR(j,i)
accu_xhi_psi += VL(j,i+1) * VR(j,i+1)
accu_chi_psi += VL(j,i) * VR(j,i+1)
accu_xhi_phi += VL(j,i+1) * VR(j,i)
enddo
mat_ovlp_orig(1,1) = accu_chi_phi
mat_ovlp_orig(2,1) = accu_xhi_phi
mat_ovlp_orig(1,2) = accu_chi_psi
mat_ovlp_orig(2,2) = accu_xhi_psi
print*,'old overlap matrix '
write(*,'(100(F16.10,X))')mat_ovlp_orig(1:2,1)
write(*,'(100(F16.10,X))')mat_ovlp_orig(1:2,2)
mat_ovlp(1,1) = accu_xhi_phi
mat_ovlp(2,1) = accu_chi_phi
mat_ovlp(1,2) = accu_xhi_psi
mat_ovlp(2,2) = accu_chi_psi
!print*,'accu_chi_phi = ',accu_chi_phi
!print*,'accu_xhi_psi = ',accu_xhi_psi
!print*,'accu_chi_psi = ',accu_chi_psi
!print*,'accu_xhi_phi = ',accu_xhi_phi
print*,'new overlap matrix '
write(*,'(100(F16.10,X))')mat_ovlp(1:2,1)
write(*,'(100(F16.10,X))')mat_ovlp(1:2,2)
call lapack_diag(eigval_tmp,eigvec,mat_ovlp,2,2)
print*,'eigval_tmp(1) = ',eigval_tmp(1)
print*,'eigvec(1) = ',eigvec(1:2,1)
print*,'eigval_tmp(2) = ',eigval_tmp(2)
print*,'eigvec(2) = ',eigvec(1:2,2)
print*,'*****************'
phi_1_tilde = 0.d0
phi_2_tilde = 0.d0
chi_1_tilde = 0.d0
chi_2_tilde = 0.d0
do j = 1, n
phi_1_tilde(j) += VR(j,i) * eigvec(1,1) + VR(j,i+1) * eigvec(2,1)
phi_2_tilde(j) += VR(j,i) * eigvec(1,2) + VR(j,i+1) * eigvec(2,2)
chi_1_tilde(j) += VL(j,i+1) * eigvec(1,1) + VL(j,i) * eigvec(2,1)
chi_2_tilde(j) += VL(j,i+1) * eigvec(1,2) + VL(j,i) * eigvec(2,2)
enddo
VR(1:n,i) = phi_1_tilde(1:n)
VR(1:n,i+1) = phi_2_tilde(1:n)
! Vl(1:n,i) = -chi_1_tilde(1:n)
! Vl(1:n,i+1) = chi_2_tilde(1:n)
i+=1
endif
endif
i+=1
enddo
!print *, ' right eigenvect bef'
!do i = 1, n
! write(*, '(1000(F16.10,X))') VR(:,i)
@ -429,6 +498,7 @@ subroutine non_hrmt_bieig(n, A, thr_d, thr_nd, leigvec, reigvec, n_real_eigv, ei
! call impose_orthog_degen_eigvec(n, eigval, reigvec)
! call impose_orthog_degen_eigvec(n, eigval, leigvec)
call reorder_degen_eigvec(n, eigval, leigvec, reigvec)
call impose_biorthog_degen_eigvec(n, eigval, leigvec, reigvec)

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

@ -1944,6 +1944,96 @@ subroutine check_orthog(n, m, V, accu_d, accu_nd, S)
end subroutine check_orthog
! ---
subroutine reorder_degen_eigvec(n, e0, L0, R0)
implicit none
integer, intent(in) :: n
double precision, intent(in) :: e0(n)
double precision, intent(inout) :: L0(n,n), R0(n,n)
logical :: complex_root
integer :: i, j, k, m
double precision :: ei, ej, de, de_thr
double precision :: accu_d, accu_nd
integer, allocatable :: deg_num(:)
double precision, allocatable :: L(:,:), R(:,:), S(:,:), S_inv_half(:,:)
! ---
allocate( deg_num(n) )
do i = 1, n
deg_num(i) = 1
enddo
de_thr = thr_degen_tc
do i = 1, n-1
ei = e0(i)
! already considered in degen vectors
if(deg_num(i).eq.0) cycle
do j = i+1, n
ej = e0(j)
de = dabs(ei - ej)
if(de .lt. de_thr) then
deg_num(i) = deg_num(i) + 1
deg_num(j) = 0
endif
enddo
enddo
do i = 1, n
if(deg_num(i) .gt. 1) then
print *, ' degen on', i, deg_num(i), e0(i)
endif
enddo
! ---
do i = 1, n
m = deg_num(i)
if(m .gt. 1) then
allocate(L(n,m))
allocate(R(n,m),S(m,m))
do j = 1, m
L(1:n,j) = L0(1:n,i+j-1)
R(1:n,j) = R0(1:n,i+j-1)
enddo
call dgemm( 'T', 'N', m, m, n, 1.d0 &
, L, size(L, 1), R, size(R, 1) &
, 0.d0, S, size(S, 1) )
print*,'Overlap matrix '
accu_nd = 0.D0
do j = 1, m
write(*,'(100(F16.10,X))')S(1:m,j)
do k = 1, m
if(j==k)cycle
accu_nd += dabs(S(j,k))
enddo
enddo
print*,'accu_nd = ',accu_nd
! if(accu_nd .gt.1.d-10)then
! stop
! endif
do j = 1, m
L0(1:n,i+j-1) = L(1:n,j)
R0(1:n,i+j-1) = R(1:n,j)
enddo
deallocate(L, R,S)
endif
enddo
end subroutine reorder_degen_eigvec
subroutine impose_biorthog_degen_eigvec(n, e0, L0, R0)
@ -2030,7 +2120,7 @@ subroutine impose_biorthog_degen_eigvec(n, e0, L0, R0)
!call bi_ortho_s_inv_half(m, L, R, S_inv_half)
!deallocate(S, S_inv_half)
! call impose_biorthog_svd(n, m, L, R)
call impose_biorthog_svd(n, m, L, R)
! call impose_biorthog_inverse(n, m, L, R)
!call impose_biorthog_qr(n, m, thr_d, thr_nd, L, R)
@ -2046,7 +2136,7 @@ subroutine impose_biorthog_degen_eigvec(n, e0, L0, R0)
endif
enddo
call impose_biorthog_inverse(n, n, L0, R0)
! call impose_biorthog_inverse(n, n, L0, R0)
end subroutine impose_biorthog_degen_eigvec

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

@ -38,6 +38,13 @@
thr_d = 1.d-6
thr_nd = 1.d-6
thr_deg = 1.d-3
do i = 1, mo_num
do j = 1, mo_num
if(dabs(dm_tmp(j,i)).lt.thr_d)then
dm_tmp(j,i) = 0.d0
endif
enddo
enddo
! if(n_core_orb.ne.0)then
! call diag_mat_per_fock_degen_core( fock_diag, dm_tmp, list_core, n_core_orb, mo_num, thr_d, thr_nd, thr_deg &
! , natorb_tc_leigvec_mo, natorb_tc_reigvec_mo, natorb_tc_eigval)