mirror of
https://github.com/LCPQ/quantum_package
synced 2024-12-23 04:43:50 +01:00
Working on mrpt
This commit is contained in:
parent
7769ea536c
commit
a592143744
@ -94,7 +94,6 @@ subroutine davidson_diag_hjj_mrcc(dets_in,u_in,H_jj,energies,dim_in,sze,N_st,N_s
|
|||||||
double precision, allocatable :: overlap(:,:)
|
double precision, allocatable :: overlap(:,:)
|
||||||
double precision :: u_dot_v, u_dot_u
|
double precision :: u_dot_v, u_dot_u
|
||||||
|
|
||||||
integer, allocatable :: kl_pairs(:,:)
|
|
||||||
integer :: k_pairs, kl
|
integer :: k_pairs, kl
|
||||||
|
|
||||||
integer :: iter2
|
integer :: iter2
|
||||||
@ -144,7 +143,6 @@ subroutine davidson_diag_hjj_mrcc(dets_in,u_in,H_jj,energies,dim_in,sze,N_st,N_s
|
|||||||
sze_8 = align_double(sze)
|
sze_8 = align_double(sze)
|
||||||
|
|
||||||
allocate( &
|
allocate( &
|
||||||
kl_pairs(2,N_st_diag*(N_st_diag+1)/2), &
|
|
||||||
W(sze_8,N_st_diag,davidson_sze_max), &
|
W(sze_8,N_st_diag,davidson_sze_max), &
|
||||||
U(sze_8,N_st_diag,davidson_sze_max), &
|
U(sze_8,N_st_diag,davidson_sze_max), &
|
||||||
R(sze_8,N_st_diag), &
|
R(sze_8,N_st_diag), &
|
||||||
@ -360,7 +358,6 @@ subroutine davidson_diag_hjj_mrcc(dets_in,u_in,H_jj,energies,dim_in,sze,N_st,N_s
|
|||||||
call write_time(iunit)
|
call write_time(iunit)
|
||||||
|
|
||||||
deallocate ( &
|
deallocate ( &
|
||||||
kl_pairs, &
|
|
||||||
W, residual_norm, &
|
W, residual_norm, &
|
||||||
U, overlap, &
|
U, overlap, &
|
||||||
R, c, &
|
R, c, &
|
||||||
@ -649,7 +646,6 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
double precision, allocatable :: overlap(:,:)
|
double precision, allocatable :: overlap(:,:)
|
||||||
double precision :: u_dot_v, u_dot_u
|
double precision :: u_dot_v, u_dot_u
|
||||||
|
|
||||||
integer, allocatable :: kl_pairs(:,:)
|
|
||||||
integer :: k_pairs, kl
|
integer :: k_pairs, kl
|
||||||
|
|
||||||
integer :: iter2
|
integer :: iter2
|
||||||
@ -661,7 +657,7 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
character*(16384) :: write_buffer
|
character*(16384) :: write_buffer
|
||||||
double precision :: to_print(3,N_st)
|
double precision :: to_print(3,N_st)
|
||||||
double precision :: cpu, wall
|
double precision :: cpu, wall
|
||||||
integer :: shift, shift2
|
integer :: shift, shift2, itermax
|
||||||
include 'constants.include.F'
|
include 'constants.include.F'
|
||||||
|
|
||||||
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: U, W, R, S, y, h, lambda
|
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: U, W, R, S, y, h, lambda
|
||||||
@ -711,21 +707,28 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
delta = 0.d0
|
delta = 0.d0
|
||||||
endif
|
endif
|
||||||
|
|
||||||
|
itermax = min(davidson_sze_max, sze/N_st_diag)
|
||||||
allocate( &
|
allocate( &
|
||||||
kl_pairs(2,N_st_diag*(N_st_diag+1)/2), &
|
W(sze_8,N_st_diag*itermax), &
|
||||||
W(sze_8,N_st_diag*davidson_sze_max), &
|
U(sze_8,N_st_diag*itermax), &
|
||||||
U(sze_8,N_st_diag*davidson_sze_max), &
|
S(sze_8,N_st_diag*itermax), &
|
||||||
R(sze_8,N_st_diag), &
|
h(N_st_diag*itermax,N_st_diag*itermax), &
|
||||||
S(sze_8,N_st_diag*davidson_sze_max), &
|
y(N_st_diag*itermax,N_st_diag*itermax), &
|
||||||
h(N_st_diag*davidson_sze_max,N_st_diag*davidson_sze_max), &
|
s_(N_st_diag*itermax,N_st_diag*itermax), &
|
||||||
y(N_st_diag*davidson_sze_max,N_st_diag*davidson_sze_max), &
|
s_tmp(N_st_diag*itermax,N_st_diag*itermax), &
|
||||||
s_(N_st_diag*davidson_sze_max,N_st_diag*davidson_sze_max), &
|
|
||||||
s_tmp(N_st_diag*davidson_sze_max,N_st_diag*davidson_sze_max), &
|
|
||||||
residual_norm(N_st_diag), &
|
residual_norm(N_st_diag), &
|
||||||
overlap(N_st_diag,N_st_diag), &
|
c(N_st_diag*itermax), &
|
||||||
c(N_st_diag*davidson_sze_max), &
|
s2(N_st_diag*itermax), &
|
||||||
s2(N_st_diag*davidson_sze_max), &
|
lambda(N_st_diag*itermax))
|
||||||
lambda(N_st_diag*davidson_sze_max))
|
|
||||||
|
h = 0.d0
|
||||||
|
s_ = 0.d0
|
||||||
|
s_tmp = 0.d0
|
||||||
|
U = 0.d0
|
||||||
|
W = 0.d0
|
||||||
|
S = 0.d0
|
||||||
|
y = 0.d0
|
||||||
|
|
||||||
|
|
||||||
ASSERT (N_st > 0)
|
ASSERT (N_st > 0)
|
||||||
ASSERT (N_st_diag >= N_st)
|
ASSERT (N_st_diag >= N_st)
|
||||||
@ -738,25 +741,25 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
|
|
||||||
converged = .False.
|
converged = .False.
|
||||||
|
|
||||||
do k=1,N_st
|
|
||||||
call normalize(u_in(1,k),sze)
|
|
||||||
enddo
|
|
||||||
|
|
||||||
do k=N_st+1,N_st_diag
|
|
||||||
do i=1,sze
|
|
||||||
double precision :: r1, r2
|
double precision :: r1, r2
|
||||||
|
do k=N_st+1,N_st_diag-2,2
|
||||||
|
do i=1,sze
|
||||||
call random_number(r1)
|
call random_number(r1)
|
||||||
call random_number(r2)
|
call random_number(r2)
|
||||||
u_in(i,k) = dsqrt(-2.d0*dlog(r1))*dcos(dtwo_pi*r2)
|
r1 = dsqrt(-2.d0*dlog(r1))
|
||||||
|
r2 = dtwo_pi*r2
|
||||||
|
u_in(i,k) = r1*dcos(r2)
|
||||||
|
u_in(i,k+1) = r1*dsin(r2)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
do k=N_st_diag-1,N_st_diag
|
||||||
|
do i=1,sze
|
||||||
|
call random_number(r1)
|
||||||
|
call random_number(r2)
|
||||||
|
r1 = dsqrt(-2.d0*dlog(r1))
|
||||||
|
r2 = dtwo_pi*r2
|
||||||
|
u_in(i,k) = r1*dcos(r2)
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
! Gram-Schmidt
|
|
||||||
! ------------
|
|
||||||
call dgemv('T',sze,k-1,1.d0,u_in,size(u_in,1), &
|
|
||||||
u_in(1,k),1,0.d0,c,1)
|
|
||||||
call dgemv('N',sze,k-1,-1.d0,u_in,size(u_in,1), &
|
|
||||||
c,1,1.d0,u_in(1,k),1)
|
|
||||||
call normalize(u_in(1,k),sze)
|
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
|
|
||||||
@ -773,11 +776,11 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
shift = N_st_diag*(iter-1)
|
shift = N_st_diag*(iter-1)
|
||||||
shift2 = N_st_diag*iter
|
shift2 = N_st_diag*iter
|
||||||
|
|
||||||
|
call ortho_qr(U,size(U,1),sze,shift2)
|
||||||
|
|
||||||
! Compute |W_k> = \sum_i |i><i|H|u_k>
|
! Compute |W_k> = \sum_i |i><i|H|u_k>
|
||||||
! -----------------------------------------
|
! -----------------------------------------
|
||||||
|
|
||||||
|
|
||||||
call H_S2_u_0_mrcc_nstates(W(1,shift+1),S(1,shift+1),U(1,shift+1),H_jj,S2_jj,sze,dets_in,Nint,&
|
call H_S2_u_0_mrcc_nstates(W(1,shift+1),S(1,shift+1),U(1,shift+1),H_jj,S2_jj,sze,dets_in,Nint,&
|
||||||
istate,N_st_diag,sze_8)
|
istate,N_st_diag,sze_8)
|
||||||
|
|
||||||
@ -786,19 +789,6 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
! -------------------------------------------
|
! -------------------------------------------
|
||||||
|
|
||||||
|
|
||||||
! do l=1,N_st_diag
|
|
||||||
! do k=1,N_st_diag
|
|
||||||
! do iter2=1,iter-1
|
|
||||||
! h(k,iter2,l,iter) = u_dot_v(U(1,k,iter2),W(1,l,iter),sze)
|
|
||||||
! h(k,iter,l,iter2) = h(k,iter2,l,iter)
|
|
||||||
! enddo
|
|
||||||
! enddo
|
|
||||||
! do k=1,l
|
|
||||||
! h(k,iter,l,iter) = u_dot_v(U(1,k,iter),W(1,l,iter),sze)
|
|
||||||
! h(l,iter,k,iter) = h(k,iter,l,iter)
|
|
||||||
! enddo
|
|
||||||
! enddo
|
|
||||||
|
|
||||||
call dgemm('T','N', shift2, N_st_diag, sze, &
|
call dgemm('T','N', shift2, N_st_diag, sze, &
|
||||||
1.d0, U, size(U,1), W(1,shift+1), size(W,1), &
|
1.d0, U, size(U,1), W(1,shift+1), size(W,1), &
|
||||||
0.d0, h(1,shift+1), size(h,1))
|
0.d0, h(1,shift+1), size(h,1))
|
||||||
@ -829,7 +819,7 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
if (s2_eig) then
|
if (s2_eig) then
|
||||||
logical :: state_ok(N_st_diag*davidson_sze_max)
|
logical :: state_ok(N_st_diag*davidson_sze_max)
|
||||||
do k=1,shift2
|
do k=1,shift2
|
||||||
state_ok(k) = (dabs(s2(k)-expected_s2) < 0.3d0)
|
state_ok(k) = (dabs(s2(k)-expected_s2) < 0.6d0)
|
||||||
enddo
|
enddo
|
||||||
do k=1,shift2
|
do k=1,shift2
|
||||||
if (.not. state_ok(k)) then
|
if (.not. state_ok(k)) then
|
||||||
@ -851,22 +841,6 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
! Express eigenvectors of h in the determinant basis
|
! Express eigenvectors of h in the determinant basis
|
||||||
! --------------------------------------------------
|
! --------------------------------------------------
|
||||||
|
|
||||||
! do k=1,N_st_diag
|
|
||||||
! do i=1,sze
|
|
||||||
! U(i,shift2+k) = 0.d0
|
|
||||||
! W(i,shift2+k) = 0.d0
|
|
||||||
! S(i,shift2+k) = 0.d0
|
|
||||||
! enddo
|
|
||||||
! do l=1,N_st_diag*iter
|
|
||||||
! do i=1,sze
|
|
||||||
! U(i,shift2+k) = U(i,shift2+k) + U(i,l)*y(l,k)
|
|
||||||
! W(i,shift2+k) = W(i,shift2+k) + W(i,l)*y(l,k)
|
|
||||||
! S(i,shift2+k) = S(i,shift2+k) + S(i,l)*y(l,k)
|
|
||||||
! enddo
|
|
||||||
! enddo
|
|
||||||
! enddo
|
|
||||||
!
|
|
||||||
!
|
|
||||||
call dgemm('N','N', sze, N_st_diag, shift2, &
|
call dgemm('N','N', sze, N_st_diag, shift2, &
|
||||||
1.d0, U, size(U,1), y, size(y,1), 0.d0, U(1,shift2+1), size(U,1))
|
1.d0, U, size(U,1), y, size(y,1), 0.d0, U(1,shift2+1), size(U,1))
|
||||||
call dgemm('N','N', sze, N_st_diag, shift2, &
|
call dgemm('N','N', sze, N_st_diag, shift2, &
|
||||||
@ -877,80 +851,36 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
! Compute residual vector
|
! Compute residual vector
|
||||||
! -----------------------
|
! -----------------------
|
||||||
|
|
||||||
! do k=1,N_st_diag
|
|
||||||
! print *, s2(k)
|
|
||||||
! s2(k) = u_dot_v(U(1,shift2+k), S(1,shift2+k), sze) + S_z2_Sz
|
|
||||||
! print *, s2(k)
|
|
||||||
! print *, ''
|
|
||||||
! pause
|
|
||||||
! enddo
|
|
||||||
do k=1,N_st_diag
|
do k=1,N_st_diag
|
||||||
do i=1,sze
|
do i=1,sze
|
||||||
R(i,k) = (lambda(k) * U(i,shift2+k) - W(i,shift2+k) ) &
|
U(i,shift2+k) = (lambda(k) * U(i,shift2+k) - W(i,shift2+k) ) &
|
||||||
* (1.d0 + s2(k) * U(i,shift2+k) - S(i,shift2+k) - S_z2_Sz)
|
* (1.d0 + s2(k) * U(i,shift2+k) - S(i,shift2+k) - S_z2_Sz &
|
||||||
|
)/max(H_jj(i) - lambda (k),1.d-2)
|
||||||
enddo
|
enddo
|
||||||
if (k <= N_st) then
|
if (k <= N_st) then
|
||||||
residual_norm(k) = u_dot_u(R(1,k),sze)
|
residual_norm(k) = u_dot_u(U(1,shift2+k),sze)
|
||||||
to_print(1,k) = lambda(k) + nuclear_repulsion
|
to_print(1,k) = lambda(k) + nuclear_repulsion
|
||||||
to_print(2,k) = s2(k)
|
to_print(2,k) = s2(k)
|
||||||
to_print(3,k) = residual_norm(k)
|
to_print(3,k) = residual_norm(k)
|
||||||
if (residual_norm(k) > 1.e9) then
|
|
||||||
stop 'Davidson failed'
|
|
||||||
endif
|
|
||||||
endif
|
endif
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
write(iunit,'(X,I3,X,100(X,F16.10,X,F11.6,X,E11.3))') iter, to_print(:,1:N_st)
|
write(iunit,'(X,I3,X,100(X,F16.10,X,F11.6,X,E11.3))') iter, to_print(1:3,1:N_st)
|
||||||
call davidson_converged(lambda,residual_norm,wall,iter,cpu,N_st,converged)
|
call davidson_converged(lambda,residual_norm,wall,iter,cpu,N_st,converged)
|
||||||
|
do k=1,N_st
|
||||||
|
if (residual_norm(k) > 1.e8) then
|
||||||
|
print *, ''
|
||||||
|
stop 'Davidson failed'
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
if (converged) then
|
if (converged) then
|
||||||
exit
|
exit
|
||||||
endif
|
endif
|
||||||
|
|
||||||
! Davidson step
|
|
||||||
! -------------
|
|
||||||
|
|
||||||
do k=1,N_st_diag
|
|
||||||
do i=1,sze
|
|
||||||
U(i,shift2+k) = - R(i,k)/max(H_jj(i) - lambda(k),1.d-2)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
|
|
||||||
! Gram-Schmidt
|
|
||||||
! ------------
|
|
||||||
|
|
||||||
do k=1,N_st_diag
|
|
||||||
|
|
||||||
! do l=1,N_st_diag*iter
|
|
||||||
! c(1) = u_dot_v(U(1,shift2+k),U(1,l),sze)
|
|
||||||
! do i=1,sze
|
|
||||||
! U(i,k,iter+1) = U(i,shift2+k) - c(1) * U(i,l)
|
|
||||||
! enddo
|
|
||||||
! enddo
|
|
||||||
!
|
|
||||||
call dgemv('T',sze,N_st_diag*iter,1.d0,U,size(U,1), &
|
|
||||||
U(1,shift2+k),1,0.d0,c,1)
|
|
||||||
call dgemv('N',sze,N_st_diag*iter,-1.d0,U,size(U,1), &
|
|
||||||
c,1,1.d0,U(1,shift2+k),1)
|
|
||||||
!
|
|
||||||
! do l=1,k-1
|
|
||||||
! c(1) = u_dot_v(U(1,shift2+k),U(1,shift2+l),sze)
|
|
||||||
! do i=1,sze
|
|
||||||
! U(i,k,iter+1) = U(i,shift2+k) - c(1) * U(i,shift2+l)
|
|
||||||
! enddo
|
|
||||||
! enddo
|
|
||||||
!
|
|
||||||
call dgemv('T',sze,k-1,1.d0,U(1,shift2+1),size(U,1), &
|
|
||||||
U(1,shift2+k),1,0.d0,c,1)
|
|
||||||
call dgemv('N',sze,k-1,-1.d0,U(1,shift2+1),size(U,1), &
|
|
||||||
c,1,1.d0,U(1,shift2+k),1)
|
|
||||||
|
|
||||||
call normalize( U(1,shift2+k), sze )
|
|
||||||
enddo
|
|
||||||
|
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
if (.not.converged) then
|
if (.not.converged) then
|
||||||
iter = davidson_sze_max-1
|
iter = itermax-1
|
||||||
endif
|
endif
|
||||||
|
|
||||||
! Re-contract to u_in
|
! Re-contract to u_in
|
||||||
@ -960,15 +890,6 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
energies(k) = lambda(k)
|
energies(k) = lambda(k)
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
! do k=1,N_st_diag
|
|
||||||
! do i=1,sze
|
|
||||||
! do l=1,iter*N_st_diag
|
|
||||||
! u_in(i,k) += U(i,l)*y(l,k)
|
|
||||||
! enddo
|
|
||||||
! enddo
|
|
||||||
! enddo
|
|
||||||
! enddo
|
|
||||||
|
|
||||||
call dgemm('N','N', sze, N_st_diag, N_st_diag*iter, 1.d0, &
|
call dgemm('N','N', sze, N_st_diag, N_st_diag*iter, 1.d0, &
|
||||||
U, size(U,1), y, size(y,1), 0.d0, u_in, size(u_in,1))
|
U, size(U,1), y, size(y,1), 0.d0, u_in, size(u_in,1))
|
||||||
|
|
||||||
@ -983,7 +904,6 @@ subroutine davidson_diag_hjj_sjj_mrcc(dets_in,u_in,H_jj,S2_jj,energies,dim_in,sz
|
|||||||
call write_time(iunit)
|
call write_time(iunit)
|
||||||
|
|
||||||
deallocate ( &
|
deallocate ( &
|
||||||
kl_pairs, &
|
|
||||||
W, residual_norm, &
|
W, residual_norm, &
|
||||||
U, overlap, &
|
U, overlap, &
|
||||||
R, c, S, &
|
R, c, S, &
|
||||||
|
|||||||
35
plugins/MRPT_Utils/give_2h2p.irp.f
Normal file
35
plugins/MRPT_Utils/give_2h2p.irp.f
Normal file
@ -0,0 +1,35 @@
|
|||||||
|
subroutine give_2h2p(contrib_2h2p)
|
||||||
|
implicit none
|
||||||
|
double precision, intent(out) :: contrib_2h2p(N_states)
|
||||||
|
integer :: i,j,k,l,m
|
||||||
|
integer :: iorb,jorb,korb,lorb
|
||||||
|
|
||||||
|
double precision :: get_mo_bielec_integral
|
||||||
|
double precision :: direct_int,exchange_int
|
||||||
|
double precision :: numerator,denominator(N_states)
|
||||||
|
|
||||||
|
contrib_2h2p = 0.d0
|
||||||
|
do i = 1, n_inact_orb
|
||||||
|
iorb = list_inact(i)
|
||||||
|
do j = 1, n_inact_orb
|
||||||
|
jorb = list_inact(j)
|
||||||
|
do k = 1, n_virt_orb
|
||||||
|
korb = list_virt(k)
|
||||||
|
do l = 1, n_virt_orb
|
||||||
|
lorb = list_virt(l)
|
||||||
|
direct_int = get_mo_bielec_integral(iorb,jorb,korb,lorb ,mo_integrals_map)
|
||||||
|
exchange_int = get_mo_bielec_integral(iorb,jorb,lorb,korb ,mo_integrals_map)
|
||||||
|
numerator = 3.d0 * direct_int*direct_int + exchange_int*exchange_int -2.d0 * exchange_int * direct_int
|
||||||
|
do m = 1, N_states
|
||||||
|
denominator(m) = fock_core_inactive_total_spin_trace(iorb,m) + fock_core_inactive_total_spin_trace(jorb,m) &
|
||||||
|
-fock_virt_total_spin_trace(korb,m) - fock_virt_total_spin_trace(lorb,m)
|
||||||
|
contrib_2h2p(m) += numerator / denominator(m)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
contrib_2h2p = contrib_2h2p*0.5d0
|
||||||
|
|
||||||
|
end
|
||||||
|
|
||||||
@ -267,18 +267,18 @@ END_PROVIDER
|
|||||||
allocate (eigenvectors(size(H_matrix_all_dets,1),N_det))
|
allocate (eigenvectors(size(H_matrix_all_dets,1),N_det))
|
||||||
allocate (eigenvalues(N_det))
|
allocate (eigenvalues(N_det))
|
||||||
call lapack_diag(eigenvalues,eigenvectors, &
|
call lapack_diag(eigenvalues,eigenvectors, &
|
||||||
Hmatrix_dressed_pt2_new_symmetrized(1,1,1),N_det,N_det)
|
H_matrix_all_dets,size(H_matrix_all_dets,1),N_det)
|
||||||
CI_electronic_dressed_pt2_new_energy(:) = 0.d0
|
CI_electronic_energy(:) = 0.d0
|
||||||
if (s2_eig) then
|
if (s2_eig) then
|
||||||
i_state = 0
|
i_state = 0
|
||||||
allocate (s2_eigvalues(N_det))
|
allocate (s2_eigvalues(N_det))
|
||||||
allocate(index_good_state_array(N_det),good_state_array(N_det))
|
allocate(index_good_state_array(N_det),good_state_array(N_det))
|
||||||
good_state_array = .False.
|
good_state_array = .False.
|
||||||
|
call u_0_S2_u_0(s2_eigvalues,eigenvectors,N_det,psi_det,N_int,&
|
||||||
|
N_det,size(eigenvectors,1))
|
||||||
do j=1,N_det
|
do j=1,N_det
|
||||||
call get_s2_u0(psi_det,eigenvectors(1,j),N_det,size(eigenvectors,1),s2)
|
|
||||||
s2_eigvalues(j) = s2
|
|
||||||
! Select at least n_states states with S^2 values closed to "expected_s2"
|
! Select at least n_states states with S^2 values closed to "expected_s2"
|
||||||
if(dabs(s2-expected_s2).le.0.3d0)then
|
if(dabs(s2_eigvalues(j)-expected_s2).le.0.5d0)then
|
||||||
i_state +=1
|
i_state +=1
|
||||||
index_good_state_array(i_state) = j
|
index_good_state_array(i_state) = j
|
||||||
good_state_array(j) = .True.
|
good_state_array(j) = .True.
|
||||||
@ -291,10 +291,10 @@ END_PROVIDER
|
|||||||
! Fill the first "i_state" states that have a correct S^2 value
|
! Fill the first "i_state" states that have a correct S^2 value
|
||||||
do j = 1, i_state
|
do j = 1, i_state
|
||||||
do i=1,N_det
|
do i=1,N_det
|
||||||
CI_dressed_pt2_new_eigenvectors(i,j) = eigenvectors(i,index_good_state_array(j))
|
CI_eigenvectors(i,j) = eigenvectors(i,index_good_state_array(j))
|
||||||
enddo
|
enddo
|
||||||
CI_electronic_dressed_pt2_new_energy(j) = eigenvalues(index_good_state_array(j))
|
CI_electronic_energy(j) = eigenvalues(index_good_state_array(j))
|
||||||
CI_dressed_pt2_new_eigenvectors_s2(j) = s2_eigvalues(index_good_state_array(j))
|
CI_eigenvectors_s2(j) = s2_eigvalues(index_good_state_array(j))
|
||||||
enddo
|
enddo
|
||||||
i_other_state = 0
|
i_other_state = 0
|
||||||
do j = 1, N_det
|
do j = 1, N_det
|
||||||
@ -303,16 +303,13 @@ END_PROVIDER
|
|||||||
if(i_state+i_other_state.gt.n_states_diag)then
|
if(i_state+i_other_state.gt.n_states_diag)then
|
||||||
exit
|
exit
|
||||||
endif
|
endif
|
||||||
call get_s2_u0(psi_det,eigenvectors(1,j),N_det,size(eigenvectors,1),s2)
|
|
||||||
do i=1,N_det
|
do i=1,N_det
|
||||||
CI_dressed_pt2_new_eigenvectors(i,i_state+i_other_state) = eigenvectors(i,j)
|
CI_eigenvectors(i,i_state+i_other_state) = eigenvectors(i,j)
|
||||||
enddo
|
enddo
|
||||||
CI_electronic_dressed_pt2_new_energy(i_state+i_other_state) = eigenvalues(j)
|
CI_electronic_energy(i_state+i_other_state) = eigenvalues(j)
|
||||||
CI_dressed_pt2_new_eigenvectors_s2(i_state+i_other_state) = s2
|
CI_eigenvectors_s2(i_state+i_other_state) = s2_eigvalues(i_state+i_other_state)
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
deallocate(index_good_state_array,good_state_array)
|
|
||||||
|
|
||||||
else
|
else
|
||||||
print*,''
|
print*,''
|
||||||
print*,'!!!!!!!! WARNING !!!!!!!!!'
|
print*,'!!!!!!!! WARNING !!!!!!!!!'
|
||||||
@ -320,27 +317,28 @@ END_PROVIDER
|
|||||||
print*,' and the ',N_states_diag,'states requested'
|
print*,' and the ',N_states_diag,'states requested'
|
||||||
print*,' We did not find any state with S^2 values close to ',expected_s2
|
print*,' We did not find any state with S^2 values close to ',expected_s2
|
||||||
print*,' We will then set the first N_states eigenvectors of the H matrix'
|
print*,' We will then set the first N_states eigenvectors of the H matrix'
|
||||||
print*,' as the CI_dressed_pt2_new_eigenvectors'
|
print*,' as the CI_eigenvectors'
|
||||||
print*,' You should consider more states and maybe ask for diagonalize_s2 to be .True. or just enlarge the CI space'
|
print*,' You should consider more states and maybe ask for s2_eig to be .True. or just enlarge the CI space'
|
||||||
print*,''
|
print*,''
|
||||||
do j=1,min(N_states_diag,N_det)
|
do j=1,min(N_states_diag,N_det)
|
||||||
do i=1,N_det
|
do i=1,N_det
|
||||||
CI_dressed_pt2_new_eigenvectors(i,j) = eigenvectors(i,j)
|
CI_eigenvectors(i,j) = eigenvectors(i,j)
|
||||||
enddo
|
enddo
|
||||||
CI_electronic_dressed_pt2_new_energy(j) = eigenvalues(j)
|
CI_electronic_energy(j) = eigenvalues(j)
|
||||||
CI_dressed_pt2_new_eigenvectors_s2(j) = s2_eigvalues(j)
|
CI_eigenvectors_s2(j) = s2_eigvalues(j)
|
||||||
enddo
|
enddo
|
||||||
endif
|
endif
|
||||||
|
deallocate(index_good_state_array,good_state_array)
|
||||||
deallocate(s2_eigvalues)
|
deallocate(s2_eigvalues)
|
||||||
else
|
else
|
||||||
|
call u_0_S2_u_0(CI_eigenvectors_s2,eigenvectors,N_det,psi_det,N_int,&
|
||||||
|
min(N_det,N_states_diag),size(eigenvectors,1))
|
||||||
! Select the "N_states_diag" states of lowest energy
|
! Select the "N_states_diag" states of lowest energy
|
||||||
do j=1,min(N_det,N_states_diag)
|
do j=1,min(N_det,N_states_diag)
|
||||||
call get_s2_u0(psi_det,eigenvectors(1,j),N_det,N_det,s2)
|
|
||||||
do i=1,N_det
|
do i=1,N_det
|
||||||
CI_dressed_pt2_new_eigenvectors(i,j) = eigenvectors(i,j)
|
CI_eigenvectors(i,j) = eigenvectors(i,j)
|
||||||
enddo
|
enddo
|
||||||
CI_electronic_dressed_pt2_new_energy(j) = eigenvalues(j)
|
CI_electronic_energy(j) = eigenvalues(j)
|
||||||
CI_dressed_pt2_new_eigenvectors_s2(j) = s2
|
|
||||||
enddo
|
enddo
|
||||||
endif
|
endif
|
||||||
deallocate(eigenvectors,eigenvalues)
|
deallocate(eigenvectors,eigenvalues)
|
||||||
|
|||||||
Loading…
Reference in New Issue
Block a user