mirror of
https://github.com/QuantumPackage/qp2.git
synced 2024-12-23 12:55:37 +01:00
fixed casscf
This commit is contained in:
eginer
parent
e3cbde9ac8
commit
9c1932eb04
@ -1,4 +1,4 @@
|
|||||||
cipsi
|
cipsi
|
||||||
selectors_full
|
selectors_full
|
||||||
generators_fluid
|
generators_cas
|
||||||
two_body_rdm
|
two_body_rdm
|
||||||
|
|||||||
@ -4,85 +4,8 @@ program casscf
|
|||||||
! TODO : Put the documentation of the program here
|
! TODO : Put the documentation of the program here
|
||||||
END_DOC
|
END_DOC
|
||||||
no_vvvv_integrals = .True.
|
no_vvvv_integrals = .True.
|
||||||
SOFT_TOUCH no_vvvv_integrals
|
pt2_max = 0.02
|
||||||
threshold_davidson = 1.d-7
|
SOFT_TOUCH no_vvvv_integrals pt2_max
|
||||||
touch threshold_davidson
|
|
||||||
if(cisd_guess)then
|
|
||||||
logical :: converged
|
|
||||||
integer :: iteration
|
|
||||||
double precision :: energy
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'USING A CISD WAVE FUNCTION AS GUESS FOR THE MCSCF WF'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
converged = .False.
|
|
||||||
iteration = 0
|
|
||||||
generators_type = "HF"
|
|
||||||
touch generators_type
|
|
||||||
read_wf = .False.
|
|
||||||
touch read_wf
|
|
||||||
logical :: do_cisdtq
|
|
||||||
do_cisdtq = .True.
|
|
||||||
double precision :: thr
|
|
||||||
thr = 5.d-3
|
|
||||||
do while (.not.converged)
|
|
||||||
call cisd_scf_iteration(converged,iteration,energy,thr)
|
|
||||||
if(HF_index.ne.1.and.iteration.gt.0)then
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'The HF determinant is not the dominant determinant in the CISD WF ...'
|
|
||||||
print*,'Therefore we skip the CISD WF ..'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
do_cisdtq = .False.
|
|
||||||
exit
|
|
||||||
endif
|
|
||||||
if(iteration.gt.15.and..not.converged)then
|
|
||||||
print*,'It seems that the orbital optimization for the CISD WAVE FUNCTION CANNOT CONVERGE ...'
|
|
||||||
print*,'Passing to CISDTQ WAVE FUNCTION'
|
|
||||||
exit
|
|
||||||
endif
|
|
||||||
enddo
|
|
||||||
if(do_cisdtq)then
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'SWITCHING WITH A CISDTQ WAVE FUNCTION AS GUESS FOR THE MCSCF WF'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
converged = .False.
|
|
||||||
iteration = 0
|
|
||||||
read_wf = .False.
|
|
||||||
touch read_wf
|
|
||||||
pt2_max = 0.01d0
|
|
||||||
touch pt2_max
|
|
||||||
energy = 0.d0
|
|
||||||
do while (.not.converged)
|
|
||||||
call cisdtq_scf_iteration(converged,iteration,energy,thr)
|
|
||||||
if(HF_index.ne.1.and.iteration.gt.0)then
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'The HF determinant is not the dominant determinant in the CISDTQ WF ...'
|
|
||||||
print*,'Therefore we skip the CISDTQ WF ..'
|
|
||||||
print*,'*******************************'
|
|
||||||
print*,'*******************************'
|
|
||||||
exit
|
|
||||||
endif
|
|
||||||
if(iteration.gt.15.and..not.converged)then
|
|
||||||
print*,'It seems that the orbital optimization for the CISDTQ WAVE FUNCTION CANNOT CONVERGE ...'
|
|
||||||
print*,'Passing to CISDTQ WAVE FUNCTION'
|
|
||||||
exit
|
|
||||||
endif
|
|
||||||
enddo
|
|
||||||
endif
|
|
||||||
endif
|
|
||||||
read_wf = .False.
|
|
||||||
touch read_wf
|
|
||||||
pt2_max = 0.0d0
|
|
||||||
touch pt2_max
|
|
||||||
! call run_cipsi_scf
|
|
||||||
call run
|
call run
|
||||||
end
|
end
|
||||||
|
|
||||||
@ -107,7 +30,7 @@ subroutine run
|
|||||||
call write_double(6,energy_improvement, 'Predicted energy improvement')
|
call write_double(6,energy_improvement, 'Predicted energy improvement')
|
||||||
|
|
||||||
converged = dabs(energy_improvement) < thresh_scf
|
converged = dabs(energy_improvement) < thresh_scf
|
||||||
! pt2_max = dabs(energy_improvement / pt2_relative_error)
|
pt2_max = dabs(energy_improvement / pt2_relative_error)
|
||||||
|
|
||||||
mo_coef = NewOrbs
|
mo_coef = NewOrbs
|
||||||
call save_mos
|
call save_mos
|
||||||
|
|||||||
@ -1,14 +0,0 @@
|
|||||||
subroutine only_act_bitmask
|
|
||||||
implicit none
|
|
||||||
integer :: i,j,k
|
|
||||||
do k = 1, N_generators_bitmask
|
|
||||||
do j = 1, 6
|
|
||||||
do i = 1, N_int
|
|
||||||
generators_bitmask(i,1,j,k) = act_bitmask(i,1)
|
|
||||||
generators_bitmask(i,2,j,k) = act_bitmask(i,2)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
touch generators_bitmask
|
|
||||||
end
|
|
||||||
|
|
||||||
@ -1,75 +0,0 @@
|
|||||||
subroutine run_cipsi_scf
|
|
||||||
implicit none
|
|
||||||
double precision :: energy_old, energy, extrap,extrap_old,pt2_max_begin
|
|
||||||
logical :: converged
|
|
||||||
integer :: iteration
|
|
||||||
print*,'*********************************'
|
|
||||||
print*,'*********************************'
|
|
||||||
print*,' DOING THE CIPSI-SCF '
|
|
||||||
print*,'*********************************'
|
|
||||||
converged = .False.
|
|
||||||
pt2_max_begin = pt2_max
|
|
||||||
energy = 0.d0
|
|
||||||
extrap = 0.d0
|
|
||||||
mo_label = "MCSCF"
|
|
||||||
iteration = 1
|
|
||||||
threshold_davidson = 1.d-09
|
|
||||||
touch threshold_davidson
|
|
||||||
do while (.not.converged)
|
|
||||||
print*,''
|
|
||||||
call write_int(6,iteration,'CI STEP OF THE ITERATION = ')
|
|
||||||
call write_double(6,pt2_max,'PT2 MAX = ')
|
|
||||||
!call cisd_guess_wf
|
|
||||||
generators_type = "CAS"
|
|
||||||
touch generators_type
|
|
||||||
call run_stochastic_cipsi
|
|
||||||
call change_orb_cipsi(converged,iteration,energy)
|
|
||||||
if(iteration.gt.n_it_scf_max.and..not.converged)then
|
|
||||||
print*,'It seems that the orbital optimization for the CISDTQ WAVE FUNCTION CANNOT CONVERGE ...'
|
|
||||||
print*,'The required delta E was :',thresh_scf
|
|
||||||
print*,'The obtained delta E was :',extrap - extrap_old
|
|
||||||
print*,'After ',iteration,'iterations ...'
|
|
||||||
print*,'Getting out of the SCF loop ...'
|
|
||||||
exit
|
|
||||||
endif
|
|
||||||
iteration += 1
|
|
||||||
enddo
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
||||||
subroutine change_orb_cipsi(converged,iteration,energy)
|
|
||||||
implicit none
|
|
||||||
double precision :: energy_old, extrap,extrap_old,pt2_max_begin
|
|
||||||
double precision, intent(inout):: energy
|
|
||||||
logical, intent(out) :: converged
|
|
||||||
integer, intent(in) :: iteration
|
|
||||||
extrap_old = energy
|
|
||||||
energy = eone+etwo+ecore
|
|
||||||
extrap = extrapolated_energy(2,1)
|
|
||||||
|
|
||||||
call write_time(6)
|
|
||||||
call write_int(6,iteration,'CAS-SCF iteration')
|
|
||||||
call write_double(6,energy,'CAS-SCF variational energy')
|
|
||||||
call write_double(6,extrap,'CAS-SCF extrapolated energy')
|
|
||||||
call write_double(6,extrap - extrap_old,'Change in extrapolated energy')
|
|
||||||
energy = extrap
|
|
||||||
call write_double(6,energy_improvement, 'Predicted energy improvement')
|
|
||||||
|
|
||||||
converged = dabs(extrap - extrap_old) < thresh_scf
|
|
||||||
pt2_max = dabs(extrap - extrap_old) * 10.d0
|
|
||||||
pt2_max = min(pt2_max,1.d-2)
|
|
||||||
pt2_max = max(pt2_max,1.d-10)
|
|
||||||
if(N_det.gt.10**6)then
|
|
||||||
pt2_max = max(pt2_max,1.d-2)
|
|
||||||
endif
|
|
||||||
|
|
||||||
mo_coef = NewOrbs
|
|
||||||
call save_mos
|
|
||||||
call map_deinit(mo_integrals_map)
|
|
||||||
N_det = N_det/2
|
|
||||||
psi_det = psi_det_sorted
|
|
||||||
psi_coef = psi_coef_sorted
|
|
||||||
read_wf = .True.
|
|
||||||
FREE mo_integrals_map mo_two_e_integrals_in_map
|
|
||||||
SOFT_TOUCH mo_coef N_det pt2_max psi_det psi_coef
|
|
||||||
end
|
|
||||||
@ -1,85 +0,0 @@
|
|||||||
subroutine cisd_scf_iteration(converged,iteration,energy,thr)
|
|
||||||
implicit none
|
|
||||||
double precision, intent(in) :: thr
|
|
||||||
logical, intent(out) :: converged
|
|
||||||
integer, intent(inout) :: iteration
|
|
||||||
double precision, intent(out) :: energy
|
|
||||||
converged = .False.
|
|
||||||
call only_act_bitmask
|
|
||||||
N_det = N_det_generators
|
|
||||||
psi_coef = psi_coef_generators
|
|
||||||
psi_det = psi_det_generators
|
|
||||||
touch N_det psi_coef psi_det
|
|
||||||
call run_cisd
|
|
||||||
call change_orb_cisd(converged,iteration,energy,thr)
|
|
||||||
end
|
|
||||||
|
|
||||||
|
|
||||||
subroutine cisd_guess_wf
|
|
||||||
implicit none
|
|
||||||
call only_act_bitmask
|
|
||||||
N_det = N_det_generators
|
|
||||||
psi_coef = psi_coef_generators
|
|
||||||
psi_det = psi_det_generators
|
|
||||||
touch N_det psi_coef psi_det
|
|
||||||
generators_type = "HF"
|
|
||||||
touch generators_type
|
|
||||||
call run_cisd
|
|
||||||
touch N_det psi_coef psi_det psi_coef_sorted psi_det_sorted psi_det_sorted_order psi_average_norm_contrib_sorted
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
subroutine change_orb_cisd(converged,iteration,energy,thr)
|
|
||||||
implicit none
|
|
||||||
double precision, intent(in) :: thr
|
|
||||||
logical, intent(inout) :: converged
|
|
||||||
integer, intent(inout) :: iteration
|
|
||||||
double precision, intent(inout) :: energy
|
|
||||||
double precision :: energy_old
|
|
||||||
energy_old = energy
|
|
||||||
|
|
||||||
energy = eone+etwo+ecore
|
|
||||||
|
|
||||||
call write_time(6)
|
|
||||||
call write_int(6,iteration,'CISD-SCF iteration')
|
|
||||||
call write_double(6,energy,'CISD-SCF energy')
|
|
||||||
call write_double(6,energy_improvement, 'Predicted energy improvement')
|
|
||||||
converged = dabs(energy_improvement) < thr
|
|
||||||
|
|
||||||
mo_coef = NewOrbs
|
|
||||||
call save_mos
|
|
||||||
call map_deinit(mo_integrals_map)
|
|
||||||
FREE mo_integrals_map mo_two_e_integrals_in_map
|
|
||||||
iteration += 1
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
||||||
subroutine run_cisd
|
|
||||||
implicit none
|
|
||||||
integer :: i
|
|
||||||
|
|
||||||
if(pseudo_sym)then
|
|
||||||
call H_apply_cisd_sym
|
|
||||||
else
|
|
||||||
call H_apply_cisd
|
|
||||||
endif
|
|
||||||
print *, 'N_det = ', N_det
|
|
||||||
print*,'******************************'
|
|
||||||
print *, 'Energies of the states:'
|
|
||||||
do i = 1,N_states
|
|
||||||
print *, i, CI_energy(i)
|
|
||||||
enddo
|
|
||||||
if (N_states > 1) then
|
|
||||||
print*,'******************************'
|
|
||||||
print*,'Excitation energies '
|
|
||||||
do i = 2, N_states
|
|
||||||
print*, i ,CI_energy(i) - CI_energy(1)
|
|
||||||
enddo
|
|
||||||
endif
|
|
||||||
psi_coef = ci_eigenvectors
|
|
||||||
SOFT_TOUCH psi_coef
|
|
||||||
call save_wavefunction
|
|
||||||
|
|
||||||
end
|
|
||||||
@ -1,47 +0,0 @@
|
|||||||
subroutine cisdtq_scf_iteration(converged,iteration,energy,thr)
|
|
||||||
implicit none
|
|
||||||
double precision, intent(in) :: thr
|
|
||||||
logical, intent(out) :: converged
|
|
||||||
integer, intent(inout) :: iteration
|
|
||||||
double precision, intent(inout) :: energy
|
|
||||||
converged = .False.
|
|
||||||
call only_act_bitmask
|
|
||||||
generators_type = "HF_SD"
|
|
||||||
threshold_generators = 0.99d0
|
|
||||||
touch threshold_generators
|
|
||||||
touch generators_type
|
|
||||||
selection_factor = 5
|
|
||||||
touch selection_factor
|
|
||||||
call run_stochastic_cipsi
|
|
||||||
call change_orb_cisdtq(converged,iteration,energy,thr)
|
|
||||||
end
|
|
||||||
|
|
||||||
subroutine change_orb_cisdtq(converged,iteration,energy,thr)
|
|
||||||
implicit none
|
|
||||||
double precision, intent(in) :: thr
|
|
||||||
logical, intent(inout) :: converged
|
|
||||||
integer, intent(inout) :: iteration
|
|
||||||
double precision, intent(inout) :: energy
|
|
||||||
double precision :: extrap,extrap_old,pt2_max_begin
|
|
||||||
extrap_old = energy
|
|
||||||
extrap = extrapolated_energy(2,1)
|
|
||||||
energy = extrap
|
|
||||||
|
|
||||||
call write_time(6)
|
|
||||||
call write_int(6,iteration,'CISDTQ-SCF iteration')
|
|
||||||
call write_double(6,energy,'CISDTQ-SCF variational energy')
|
|
||||||
call write_double(6,extrap,'CISDTQ-SCF extrapolated energy')
|
|
||||||
call write_double(6,extrap - extrap_old,'Change in extrapolated energy')
|
|
||||||
|
|
||||||
converged = dabs(extrap - extrap_old) < thr
|
|
||||||
pt2_max = dabs(extrap - extrap_old) * 10.d0
|
|
||||||
pt2_max = max(pt2_max,1.d-10)
|
|
||||||
|
|
||||||
mo_coef = NewOrbs
|
|
||||||
call save_mos
|
|
||||||
call map_deinit(mo_integrals_map)
|
|
||||||
FREE mo_integrals_map mo_two_e_integrals_in_map
|
|
||||||
iteration += 1
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
||||||
@ -56,8 +56,8 @@ BEGIN_PROVIDER [real*8, P0tuvx, (n_act_orb,n_act_orb,n_act_orb,n_act_orb) ]
|
|||||||
uu = list_act(u)
|
uu = list_act(u)
|
||||||
do t = 1, n_act_orb
|
do t = 1, n_act_orb
|
||||||
tt = list_act(t)
|
tt = list_act(t)
|
||||||
! P0tuvx(t,u,v,x) = state_av_act_two_rdm_openmp_spin_trace_mo(t,v,u,x)
|
|
||||||
P0tuvx(t,u,v,x) = state_av_act_two_rdm_spin_trace_mo(t,v,u,x)
|
P0tuvx(t,u,v,x) = state_av_act_two_rdm_spin_trace_mo(t,v,u,x)
|
||||||
|
! P0tuvx(t,u,v,x) = act_two_rdm_spin_trace_mo(t,v,u,x)
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
|||||||
@ -171,11 +171,9 @@ BEGIN_PROVIDER [real*8, gradvec2, (nMonoEx)]
|
|||||||
norm_grad+=gradvec2(indx)*gradvec2(indx)
|
norm_grad+=gradvec2(indx)*gradvec2(indx)
|
||||||
end do
|
end do
|
||||||
norm_grad=sqrt(norm_grad)
|
norm_grad=sqrt(norm_grad)
|
||||||
! if (bavard) then
|
|
||||||
write(6,*)
|
write(6,*)
|
||||||
write(6,*) ' Norm of the orbital gradient (via D, P and integrals): ', norm_grad
|
write(6,*) ' Norm of the orbital gradient (via D, P and integrals): ', norm_grad
|
||||||
write(6,*)
|
write(6,*)
|
||||||
! endif
|
|
||||||
|
|
||||||
END_PROVIDER
|
END_PROVIDER
|
||||||
|
|
||||||
|
|||||||
@ -1,18 +0,0 @@
|
|||||||
! Generates subroutine H_apply_cisd
|
|
||||||
! ----------------------------------
|
|
||||||
|
|
||||||
BEGIN_SHELL [ /usr/bin/env python2 ]
|
|
||||||
from generate_h_apply import H_apply
|
|
||||||
H = H_apply("cisd",do_double_exc=True)
|
|
||||||
print H
|
|
||||||
|
|
||||||
from generate_h_apply import H_apply
|
|
||||||
H = H_apply("cisdtq",do_double_exc=True)
|
|
||||||
H.set_selection_pt2("epstein_nesbet_2x2")
|
|
||||||
print H
|
|
||||||
|
|
||||||
H = H_apply("cisd_sym",do_double_exc=True)
|
|
||||||
H.filter_only_connected_to_hf()
|
|
||||||
print H
|
|
||||||
END_SHELL
|
|
||||||
|
|
||||||
@ -536,9 +536,6 @@ real*8 function hessmat_taub(t,a,u,b)
|
|||||||
integer :: v3,x3
|
integer :: v3,x3
|
||||||
real*8 :: term,t1,t2,t3
|
real*8 :: term,t1,t2,t3
|
||||||
|
|
||||||
double precision,allocatable :: P0tuvx_no_t(:,:,:)
|
|
||||||
double precision :: bielec_pqxx_no_2(n_act_orb,n_act_orb)
|
|
||||||
double precision :: bielec_pxxq_no_2(n_act_orb,n_act_orb)
|
|
||||||
tt=list_act(t)
|
tt=list_act(t)
|
||||||
aa=list_virt(a)
|
aa=list_virt(a)
|
||||||
if (t == u) then
|
if (t == u) then
|
||||||
@ -548,87 +545,59 @@ real*8 function hessmat_taub(t,a,u,b)
|
|||||||
t2=0.D0
|
t2=0.D0
|
||||||
t3=0.D0
|
t3=0.D0
|
||||||
t1-=occnum(tt)*Fipq(tt,tt)
|
t1-=occnum(tt)*Fipq(tt,tt)
|
||||||
do x=1,n_act_orb
|
|
||||||
xx=list_act(x)
|
|
||||||
x3=x+n_core_inact_orb
|
|
||||||
do v=1,n_act_orb
|
|
||||||
vv=list_act(v)
|
|
||||||
v3=v+n_core_inact_orb
|
|
||||||
t2+=P0tuvx_no(t,t,v,x)*bielec_pqxx_no(aa,aa,v3,x3)
|
|
||||||
end do
|
|
||||||
end do
|
|
||||||
do v=1,n_act_orb
|
do v=1,n_act_orb
|
||||||
vv=list_act(v)
|
vv=list_act(v)
|
||||||
v3=v+n_core_inact_orb
|
v3=v+n_core_inact_orb
|
||||||
do x=1,n_act_orb
|
do x=1,n_act_orb
|
||||||
xx=list_act(x)
|
xx=list_act(x)
|
||||||
x3=x+n_core_inact_orb
|
x3=x+n_core_inact_orb
|
||||||
t2+=(P0tuvx_no(t,x,v,t)+P0tuvx_no(t,x,t,v))* &
|
t2+=2.D0*(P0tuvx_no(t,t,v,x)*bielec_pqxx_no(aa,aa,v3,x3) &
|
||||||
bielec_pxxq_no(aa,x3,v3,aa)
|
+(P0tuvx_no(t,x,v,t)+P0tuvx_no(t,x,t,v))* &
|
||||||
end do
|
bielec_pxxq_no(aa,x3,v3,aa))
|
||||||
end do
|
|
||||||
do y=1,n_act_orb
|
do y=1,n_act_orb
|
||||||
do x=1,n_act_orb
|
t3-=2.D0*P0tuvx_no(t,v,x,y)*bielecCI_no(t,v,y,xx)
|
||||||
xx=list_act(x)
|
|
||||||
do v=1,n_act_orb
|
|
||||||
t3-=P0tuvx_no(t,v,x,y)*bielecCI_no(t,v,y,xx)
|
|
||||||
end do
|
end do
|
||||||
end do
|
end do
|
||||||
end do
|
end do
|
||||||
term=t1+2.d0*(t2+t3)
|
term=t1+t2+t3
|
||||||
else
|
else
|
||||||
bb=list_virt(b)
|
bb=list_virt(b)
|
||||||
! ta/tb b/=a
|
! ta/tb b/=a
|
||||||
term=0.5d0*occnum(tt)*Fipq(aa,bb)
|
term=occnum(tt)*Fipq(aa,bb)
|
||||||
do x=1,n_act_orb
|
|
||||||
xx=list_act(x)
|
|
||||||
x3=x+n_core_inact_orb
|
|
||||||
do v=1,n_act_orb
|
|
||||||
vv=list_act(v)
|
|
||||||
v3=v+n_core_inact_orb
|
|
||||||
term = term + P0tuvx_no(t,t,v,x)*bielec_pqxx_no(aa,bb,v3,x3)
|
|
||||||
end do
|
|
||||||
end do
|
|
||||||
do v=1,n_act_orb
|
do v=1,n_act_orb
|
||||||
vv=list_act(v)
|
vv=list_act(v)
|
||||||
v3=v+n_core_inact_orb
|
v3=v+n_core_inact_orb
|
||||||
do x=1,n_act_orb
|
do x=1,n_act_orb
|
||||||
xx=list_act(x)
|
xx=list_act(x)
|
||||||
x3=x+n_core_inact_orb
|
x3=x+n_core_inact_orb
|
||||||
term= term + (P0tuvx_no(t,x,v,t)+P0tuvx_no(t,x,t,v)) &
|
term+=2.D0*(P0tuvx_no(t,t,v,x)*bielec_pqxx_no(aa,bb,v3,x3) &
|
||||||
*bielec_pxxq_no(aa,x3,v3,bb)
|
+(P0tuvx_no(t,x,v,t)+P0tuvx_no(t,x,t,v)) &
|
||||||
|
*bielec_pxxq_no(aa,x3,v3,bb))
|
||||||
end do
|
end do
|
||||||
end do
|
end do
|
||||||
term += term
|
|
||||||
end if
|
end if
|
||||||
else
|
else
|
||||||
! ta/ub t/=u
|
! ta/ub t/=u
|
||||||
uu=list_act(u)
|
uu=list_act(u)
|
||||||
bb=list_virt(b)
|
bb=list_virt(b)
|
||||||
allocate(P0tuvx_no_t(n_act_orb,n_act_orb,n_act_orb))
|
|
||||||
P0tuvx_no_t(:,:,:) = P0tuvx_no(t,:,:,:)
|
|
||||||
do x=1,n_act_orb
|
|
||||||
x3=x+n_core_inact_orb
|
|
||||||
do v=1,n_act_orb
|
|
||||||
v3=v+n_core_inact_orb
|
|
||||||
bielec_pqxx_no_2(v,x) = bielec_pqxx_no(aa,bb,v3,x3)
|
|
||||||
bielec_pxxq_no_2(v,x) = bielec_pxxq_no(aa,v3,x3,bb)
|
|
||||||
end do
|
|
||||||
end do
|
|
||||||
term=0.D0
|
term=0.D0
|
||||||
do x=1,n_act_orb
|
|
||||||
do v=1,n_act_orb
|
do v=1,n_act_orb
|
||||||
term += P0tuvx_no_t(u,v,x)*bielec_pqxx_no_2(v,x)
|
vv=list_act(v)
|
||||||
term += bielec_pxxq_no_2(x,v) * (P0tuvx_no_t(x,v,u)+P0tuvx_no_t(x,u,v))
|
v3=v+n_core_inact_orb
|
||||||
|
do x=1,n_act_orb
|
||||||
|
xx=list_act(x)
|
||||||
|
x3=x+n_core_inact_orb
|
||||||
|
term+=2.D0*(P0tuvx_no(t,u,v,x)*bielec_pqxx_no(aa,bb,v3,x3) &
|
||||||
|
+(P0tuvx_no(t,x,v,u)+P0tuvx_no(t,x,u,v)) &
|
||||||
|
*bielec_pxxq_no(aa,x3,v3,bb))
|
||||||
end do
|
end do
|
||||||
end do
|
end do
|
||||||
term = 6.d0*term
|
|
||||||
if (a.eq.b) then
|
if (a.eq.b) then
|
||||||
term-=0.5D0*(occnum(tt)*Fipq(uu,tt)+occnum(uu)*Fipq(tt,uu))
|
term-=0.5D0*(occnum(tt)*Fipq(uu,tt)+occnum(uu)*Fipq(tt,uu))
|
||||||
do v=1,n_act_orb
|
do v=1,n_act_orb
|
||||||
do y=1,n_act_orb
|
do y=1,n_act_orb
|
||||||
do x=1,n_act_orb
|
do x=1,n_act_orb
|
||||||
term-=P0tuvx_no_t(v,x,y)*bielecCI_no(x,y,v,uu)
|
term-=P0tuvx_no(t,v,x,y)*bielecCI_no(x,y,v,uu)
|
||||||
term-=P0tuvx_no(u,v,x,y)*bielecCI_no(x,y,v,tt)
|
term-=P0tuvx_no(u,v,x,y)*bielecCI_no(x,y,v,tt)
|
||||||
end do
|
end do
|
||||||
end do
|
end do
|
||||||
|
|||||||
@ -1 +0,0 @@
|
|||||||
determinants
|
|
||||||
@ -1,23 +0,0 @@
|
|||||||
subroutine extract_cas
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! Replaces the total wave function by the normalized projection on the CAS.
|
|
||||||
END_DOC
|
|
||||||
|
|
||||||
integer :: i,j,k
|
|
||||||
do k=1,N_states
|
|
||||||
do j=1,N_det_generators
|
|
||||||
psi_coef(j,k) = psi_coef_generators(j,k)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
|
|
||||||
do j=1,N_det_generators
|
|
||||||
do k=1,N_int
|
|
||||||
psi_det(k,1,j) = psi_det_generators(k,1,j)
|
|
||||||
psi_det(k,2,j) = psi_det_generators(k,2,j)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
N_det = N_det_generators
|
|
||||||
|
|
||||||
SOFT_TOUCH N_det psi_det psi_coef
|
|
||||||
end
|
|
||||||
@ -1,101 +0,0 @@
|
|||||||
use bitmasks
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ character*(32), generators_type]
|
|
||||||
implicit none
|
|
||||||
generators_type = trim("CAS")
|
|
||||||
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer, N_det_generators ]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! Number of generator detetrminants
|
|
||||||
END_DOC
|
|
||||||
if(generators_type == "CAS")then
|
|
||||||
N_det_generators = N_det_generators_CAS
|
|
||||||
else if (generators_type == "HF")then
|
|
||||||
N_det_generators = N_det_generators_HF
|
|
||||||
else if (generators_type == "HF_SD")then
|
|
||||||
N_det_generators = N_det_generators_HF_SD
|
|
||||||
endif
|
|
||||||
N_det_generators = max(N_det_generators,1)
|
|
||||||
call write_int(6,N_det_generators,'Number of generators')
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators, (N_int,2,psi_det_size) ]
|
|
||||||
&BEGIN_PROVIDER [ double precision, psi_coef_generators, (psi_det_size,N_states) ]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! For Single reference wave functions, the generator is the
|
|
||||||
! Hartree-Fock determinant
|
|
||||||
END_DOC
|
|
||||||
|
|
||||||
if(generators_type == "CAS")then
|
|
||||||
psi_det_generators(1:N_int,1:2,1:N_det_generators_CAS) = psi_det_generators_CAS(1:N_int,1:2,1:N_det_generators_CAS)
|
|
||||||
psi_coef_generators(1:N_det_generators_CAS,1:N_states) = psi_coef_generators_CAS(1:N_det_generators_CAS,1:N_states)
|
|
||||||
else if (generators_type == "HF")then
|
|
||||||
psi_det_generators(1:N_int,1:2,1:N_det_generators_HF) = psi_det_generators_HF(1:N_int,1:2,1:N_det_generators_HF)
|
|
||||||
psi_coef_generators(1:N_det_generators_HF,1:N_states) = psi_coef_generators_HF(1:N_det_generators_HF,1:N_states)
|
|
||||||
else if (generators_type == "HF_SD")then
|
|
||||||
psi_det_generators(1:N_int,1:2,1:N_det_generators_HF_SD) = psi_det_generators_HF_SD(1:N_int,1:2,1:N_det_generators_HF_SD)
|
|
||||||
psi_coef_generators(1:N_det_generators_HF_SD,1:N_states) = psi_coef_generators_HF_SD(1:N_det_generators_HF_SD,1:N_states)
|
|
||||||
endif
|
|
||||||
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer(bit_kind), psi_det_sorted_gen, (N_int,2,psi_det_size) ]
|
|
||||||
&BEGIN_PROVIDER [ double precision, psi_coef_sorted_gen, (psi_det_size,N_states) ]
|
|
||||||
&BEGIN_PROVIDER [ integer, psi_det_sorted_gen_order, (psi_det_size) ]
|
|
||||||
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! For Single reference wave functions, the generator is the
|
|
||||||
! Hartree-Fock determinant
|
|
||||||
END_DOC
|
|
||||||
if(generators_type == "CAS")then
|
|
||||||
psi_det_sorted_gen = psi_det_sorted_gen_CAS
|
|
||||||
psi_coef_sorted_gen = psi_coef_sorted_gen_CAS
|
|
||||||
psi_det_sorted_gen_order = psi_det_sorted_gen_CAS_order
|
|
||||||
else if(generators_type == "HF")then
|
|
||||||
psi_det_sorted_gen = 0_bit_kind
|
|
||||||
psi_coef_sorted_gen = 0.d0
|
|
||||||
psi_det_sorted_gen_order = 0
|
|
||||||
else if(generators_type == "HF_SD")then
|
|
||||||
psi_det_sorted_gen = psi_det_sorted_gen_HF_SD
|
|
||||||
psi_coef_sorted_gen = psi_coef_sorted_gen_HF_SD
|
|
||||||
psi_det_sorted_gen_order = psi_det_sorted_gen_HF_SD_order
|
|
||||||
endif
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [integer, degree_max_generators]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! Max degree of excitation (respect to HF) of the generators
|
|
||||||
END_DOC
|
|
||||||
integer :: i,degree
|
|
||||||
degree_max_generators = 0
|
|
||||||
do i = 1, N_det_generators
|
|
||||||
call get_excitation_degree(HF_bitmask,psi_det_generators(1,1,i),degree,N_int)
|
|
||||||
if(degree .gt. degree_max_generators)then
|
|
||||||
degree_max_generators = degree
|
|
||||||
endif
|
|
||||||
enddo
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer, size_select_max]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! Size of the select_max array
|
|
||||||
END_DOC
|
|
||||||
size_select_max = 10000
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ double precision, select_max, (size_select_max) ]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! Memo to skip useless selectors
|
|
||||||
END_DOC
|
|
||||||
select_max = huge(1.d0)
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
@ -1,69 +0,0 @@
|
|||||||
use bitmasks
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer, N_det_generators_CAS ]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! Number of generator detetrminants
|
|
||||||
END_DOC
|
|
||||||
integer :: i,k,l
|
|
||||||
logical :: good
|
|
||||||
integer, external :: number_of_holes,number_of_particles
|
|
||||||
call write_time(6)
|
|
||||||
N_det_generators_CAS = 0
|
|
||||||
do i=1,N_det
|
|
||||||
good = ( number_of_holes(psi_det_sorted(1,1,i)) ==0).and.(number_of_particles(psi_det_sorted(1,1,i))==0 )
|
|
||||||
if (good) then
|
|
||||||
N_det_generators_CAS += 1
|
|
||||||
endif
|
|
||||||
enddo
|
|
||||||
N_det_generators_CAS = max(N_det_generators_CAS,1)
|
|
||||||
call write_int(6,N_det_generators_CAS,'Number of generators_CAS')
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators_CAS, (N_int,2,psi_det_size) ]
|
|
||||||
&BEGIN_PROVIDER [ double precision, psi_coef_generators_CAS, (psi_det_size,N_states) ]
|
|
||||||
&BEGIN_PROVIDER [ integer(bit_kind), psi_det_sorted_gen_CAS, (N_int,2,psi_det_size) ]
|
|
||||||
&BEGIN_PROVIDER [ double precision, psi_coef_sorted_gen_CAS, (psi_det_size,N_states) ]
|
|
||||||
&BEGIN_PROVIDER [ integer, psi_det_sorted_gen_CAS_order, (psi_det_size) ]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! For Single reference wave functions, the gen_CASerator is the
|
|
||||||
! Hartree-Fock determinant
|
|
||||||
END_DOC
|
|
||||||
integer :: i, k, l, m
|
|
||||||
logical :: good
|
|
||||||
integer, external :: number_of_holes,number_of_particles
|
|
||||||
integer, allocatable :: nongen_CAS(:)
|
|
||||||
integer :: inongen_CAS
|
|
||||||
|
|
||||||
allocate(nongen_CAS(N_det))
|
|
||||||
|
|
||||||
inongen_CAS = 0
|
|
||||||
m=0
|
|
||||||
do i=1,N_det
|
|
||||||
good = ( number_of_holes(psi_det_sorted(1,1,i)) ==0).and.(number_of_particles(psi_det_sorted(1,1,i))==0 )
|
|
||||||
if (good) then
|
|
||||||
m = m+1
|
|
||||||
psi_det_sorted_gen_CAS_order(i) = m
|
|
||||||
do k=1,N_int
|
|
||||||
psi_det_generators_CAS(k,1,m) = psi_det_sorted(k,1,i)
|
|
||||||
psi_det_generators_CAS(k,2,m) = psi_det_sorted(k,2,i)
|
|
||||||
enddo
|
|
||||||
psi_coef_generators_CAS(m,:) = psi_coef_sorted(i,:)
|
|
||||||
else
|
|
||||||
inongen_CAS += 1
|
|
||||||
nongen_CAS(inongen_CAS) = i
|
|
||||||
endif
|
|
||||||
enddo
|
|
||||||
ASSERT (m == N_det_generators_CAS)
|
|
||||||
|
|
||||||
psi_det_sorted_gen_CAS(:,:,:N_det_generators_CAS) = psi_det_generators_CAS(:,:,:N_det_generators_CAS)
|
|
||||||
psi_coef_sorted_gen_CAS(:N_det_generators_CAS, :) = psi_coef_generators_CAS(:N_det_generators_CAS, :)
|
|
||||||
do i=1,inongen_CAS
|
|
||||||
psi_det_sorted_gen_CAS_order(nongen_CAS(i)) = N_det_generators_CAS+i
|
|
||||||
psi_det_sorted_gen_CAS(:,:,N_det_generators_CAS+i) = psi_det_sorted(:,:,nongen_CAS(i))
|
|
||||||
psi_coef_sorted_gen_CAS(N_det_generators_CAS+i, :) = psi_coef_sorted(nongen_CAS(i),:)
|
|
||||||
end do
|
|
||||||
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
@ -1,51 +0,0 @@
|
|||||||
|
|
||||||
use bitmasks
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer, N_det_generators_HF ]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! For Single reference wave functions, the number of generators is 1 : the
|
|
||||||
! Hartree-Fock determinant
|
|
||||||
END_DOC
|
|
||||||
N_det_generators_HF = 1
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators_HF, (N_int,2,psi_det_size) ]
|
|
||||||
&BEGIN_PROVIDER [ double precision, psi_coef_generators_HF, (psi_det_size,N_states) ]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! For Single reference wave functions, the generator is the
|
|
||||||
! Hartree-Fock determinant
|
|
||||||
END_DOC
|
|
||||||
psi_det_generators_HF = 0_bit_kind
|
|
||||||
integer :: i,j
|
|
||||||
integer :: degree
|
|
||||||
|
|
||||||
do i=1,N_int
|
|
||||||
psi_det_generators_HF(i,1,1) = HF_bitmask(i,1)
|
|
||||||
psi_det_generators_HF(i,2,1) = HF_bitmask(i,2)
|
|
||||||
enddo
|
|
||||||
|
|
||||||
do j=1,N_det
|
|
||||||
call get_excitation_degree(HF_bitmask,psi_det(1,1,j),degree,N_int)
|
|
||||||
if (degree == 0) then
|
|
||||||
exit
|
|
||||||
endif
|
|
||||||
end do
|
|
||||||
|
|
||||||
psi_det_generators_HF(:,:,1) = psi_det(:,:,j)
|
|
||||||
psi_coef_generators_HF(1,:) = 1.d0
|
|
||||||
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer , HF_index ]
|
|
||||||
implicit none
|
|
||||||
integer :: j,degree
|
|
||||||
do j=1,N_det
|
|
||||||
call get_excitation_degree(HF_bitmask,psi_det_sorted(1,1,j),degree,N_int)
|
|
||||||
if (degree == 0) then
|
|
||||||
HF_index = j
|
|
||||||
exit
|
|
||||||
endif
|
|
||||||
end do
|
|
||||||
END_PROVIDER
|
|
||||||
@ -1,80 +0,0 @@
|
|||||||
|
|
||||||
use bitmasks
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer, N_det_generators_HF_SD ]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! For Single reference wave functions, the number of generators is 1 : the
|
|
||||||
! Hartree-Fock determinant
|
|
||||||
END_DOC
|
|
||||||
N_det_generators_HF_SD = 0
|
|
||||||
integer :: i,degree
|
|
||||||
double precision :: thr
|
|
||||||
double precision :: accu
|
|
||||||
accu = 0.d0
|
|
||||||
thr = threshold_generators
|
|
||||||
do i = 1, N_det
|
|
||||||
call get_excitation_degree(HF_bitmask,psi_det_sorted(1,1,i),degree,N_int)
|
|
||||||
if(degree.le.2.and. accu .le. thr )then
|
|
||||||
accu += psi_coef_sorted(i,1)**2
|
|
||||||
N_det_generators_HF_SD += 1
|
|
||||||
endif
|
|
||||||
enddo
|
|
||||||
!print*,''
|
|
||||||
!print*,'N_det_generators_HF_SD = ',N_det_generators_HF_SD
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators_HF_SD, (N_int,2,psi_det_size) ]
|
|
||||||
&BEGIN_PROVIDER [ double precision, psi_coef_generators_HF_SD, (psi_det_size,N_states) ]
|
|
||||||
&BEGIN_PROVIDER [ integer(bit_kind), psi_det_sorted_gen_HF_SD, (N_int,2,psi_det_size) ]
|
|
||||||
&BEGIN_PROVIDER [ double precision, psi_coef_sorted_gen_HF_SD, (psi_det_size,N_states) ]
|
|
||||||
&BEGIN_PROVIDER [ integer, psi_det_sorted_gen_HF_SD_order, (psi_det_size) ]
|
|
||||||
implicit none
|
|
||||||
BEGIN_DOC
|
|
||||||
! For Single reference wave functions, the generator is the
|
|
||||||
! Hartree-Fock determinant
|
|
||||||
END_DOC
|
|
||||||
psi_det_generators_HF_SD = 0_bit_kind
|
|
||||||
integer :: i,j,k
|
|
||||||
integer :: degree
|
|
||||||
double precision :: thr
|
|
||||||
double precision :: accu
|
|
||||||
integer, allocatable :: nongen(:)
|
|
||||||
integer :: inongen
|
|
||||||
|
|
||||||
allocate(nongen(N_det))
|
|
||||||
|
|
||||||
thr = threshold_generators
|
|
||||||
|
|
||||||
accu = 0.d0
|
|
||||||
k = 0
|
|
||||||
inongen = 0
|
|
||||||
do j=1,N_det
|
|
||||||
call get_excitation_degree(HF_bitmask,psi_det_sorted(1,1,j),degree,N_int)
|
|
||||||
if(degree.le.2.and. accu.le.thr )then
|
|
||||||
accu += psi_coef_sorted(j,1)**2
|
|
||||||
k += 1
|
|
||||||
psi_det_sorted_gen_HF_SD_order(j) = k
|
|
||||||
do i = 1, N_int
|
|
||||||
psi_det_generators_HF_SD(i,1,k) = psi_det_sorted(i,1,j)
|
|
||||||
psi_det_generators_HF_SD(i,2,k) = psi_det_sorted(i,2,j)
|
|
||||||
enddo
|
|
||||||
do i = 1, N_states
|
|
||||||
psi_coef_generators_HF_SD(k,i) = psi_coef_sorted(j,i)
|
|
||||||
enddo
|
|
||||||
else
|
|
||||||
inongen += 1
|
|
||||||
nongen(inongen) = j
|
|
||||||
endif
|
|
||||||
end do
|
|
||||||
|
|
||||||
psi_det_sorted_gen_HF_SD(:,:,:N_det_generators_HF_SD) = psi_det_generators_HF_SD(:,:,:N_det_generators_HF_SD)
|
|
||||||
psi_coef_sorted_gen_HF_SD(:N_det_generators_HF_SD, :) = psi_coef_generators_HF_SD(:N_det_generators_HF_SD, :)
|
|
||||||
do i=1,inongen
|
|
||||||
psi_det_sorted_gen_HF_SD_order(nongen(i)) = N_det_generators_HF_SD+i
|
|
||||||
psi_det_sorted_gen_HF_SD(:,:,N_det_generators_HF_SD+i) = psi_det_sorted(:,:,nongen(i))
|
|
||||||
psi_coef_sorted_gen_HF_SD(N_det_generators_HF_SD+i, :) = psi_coef_sorted(nongen(i),:)
|
|
||||||
end do
|
|
||||||
|
|
||||||
END_PROVIDER
|
|
||||||
|
|
||||||
Loading…
Reference in New Issue
Block a user