10
0
mirror of https://github.com/LCPQ/quantum_package synced 2024-12-23 04:43:50 +01:00

Corrected Hartree-Fock MOs when not converged

This commit is contained in:
Anthony Scemama 2015-11-13 10:51:03 +01:00
parent 27eb074c26
commit e033596481
4 changed files with 32 additions and 41 deletions

View File

@ -48,10 +48,7 @@ subroutine run
E0 = HF_energy E0 = HF_energy
thresh_SCF = 1.d-10
call damping_SCF
mo_label = "Canonical" mo_label = "Canonical"
TOUCH mo_label mo_coef call damping_SCF
call save_mos
end end

View File

@ -86,7 +86,7 @@ subroutine damping_SCF
if ((E_half > E).and.(E_new < E)) then if ((E_half > E).and.(E_new < E)) then
lambda = 1.d0 lambda = 1.d0
exit exit
else if ((E_half > E).and.(lambda > 5.d-2)) then else if ((E_half > E).and.(lambda > 5.d-4)) then
lambda = 0.5d0 * lambda lambda = 0.5d0 * lambda
E_new = E_half E_new = E_half
else else

View File

@ -102,9 +102,11 @@ subroutine tamiser(key, idx, no, n, Nint, N_key)
k = no k = no
j = 2*k j = 2*k
do while(j <= n) do while(j <= n)
if(j < n .and. det_inf(key(:,:,j), key(:,:,j+1), Nint)) then if(j < n) then
if (det_inf(key(:,:,j), key(:,:,j+1), Nint)) then
j = j+1 j = j+1
endif endif
endif
if(det_inf(key(:,:,k), key(:,:,j), Nint)) then if(det_inf(key(:,:,k), key(:,:,j), Nint)) then
tmp(:,:) = key(:,:,k) tmp(:,:) = key(:,:,k)
key(:,:,k) = key(:,:,j) key(:,:,k) = key(:,:,j)
@ -113,7 +115,7 @@ subroutine tamiser(key, idx, no, n, Nint, N_key)
idx(k) = idx(j) idx(k) = idx(j)
idx(j) = tmpidx idx(j) = tmpidx
k = j k = j
j = 2*k j = k+k
else else
return return
endif endif

View File

@ -9,13 +9,19 @@ program first_guess
double precision :: E double precision :: E
integer, allocatable :: kept(:) integer, allocatable :: kept(:)
integer :: nelec_kept(2) integer :: nelec_kept(2)
character :: occ_char, keep_char
PROVIDE H_apply_buffer_allocated PROVIDE H_apply_buffer_allocated psi_det
allocate (orb_energy(mo_tot_num), kept(0:mo_tot_num)) allocate (orb_energy(mo_tot_num), kept(0:mo_tot_num))
nelec_kept(1:2) = 0 nelec_kept(1:2) = 0
kept(0) = 0 kept(0) = 0
print *, 'Orbital energies'
print *, '================'
print *, ''
do i=1,mo_tot_num do i=1,mo_tot_num
keep_char = ' '
occ_char = '-'
orb_energy(i) = mo_mono_elec_integral(i,i) orb_energy(i) = mo_mono_elec_integral(i,i)
do j=1,elec_beta_num do j=1,elec_beta_num
if (i==j) cycle if (i==j) cycle
@ -24,8 +30,9 @@ program first_guess
do j=1,elec_alpha_num do j=1,elec_alpha_num
orb_energy(i) += mo_bielec_integral_jj(i,j) orb_energy(i) += mo_bielec_integral_jj(i,j)
enddo enddo
if ( (orb_energy(i) > -1.d0).and.(orb_energy(i) < .5d0) ) then if ( (orb_energy(i) > -.5d0).and.(orb_energy(i) < .1d0) ) then
kept(0) += 1 kept(0) += 1
keep_char = 'X'
kept( kept(0) ) = i kept( kept(0) ) = i
if (i <= elec_beta_num) then if (i <= elec_beta_num) then
nelec_kept(2) += 1 nelec_kept(2) += 1
@ -34,8 +41,17 @@ program first_guess
nelec_kept(1) += 1 nelec_kept(1) += 1
endif endif
endif endif
if (i <= elec_alpha_num) then
if (i <= elec_beta_num) then
occ_char = '#'
else
occ_char = '+'
endif
endif
print '(I4, 3X, A, 3X, F10.6, 3X, A)', i, occ_char, orb_energy(i), keep_char
enddo enddo
integer, allocatable :: list (:,:) integer, allocatable :: list (:,:)
integer(bit_kind), allocatable :: string(:,:) integer(bit_kind), allocatable :: string(:,:)
allocate ( list(N_int*bit_kind_size,2), string(N_int,2) ) allocate ( list(N_int*bit_kind_size,2), string(N_int,2) )
@ -50,34 +66,10 @@ program first_guess
N_det_beta_unique = 1 N_det_beta_unique = 1
integer :: i1,i2,i3,i4,i5,i6,i7,i8,i9 integer :: i1,i2,i3,i4,i5,i6,i7,i8,i9
! select case (nelec_kept(1))
!
! case(0)
! continue
!
! case(1)
! do i1=kept(0),1,-1
! list(elec_alpha_num-0,1) = kept(i1)
! call list_to_bitstring( string(1,1), list(1,1), elec_alpha_num, N_int)
! N_det_alpha_unique += 1
! psi_det_alpha_unique(:,N_det_alpha_unique) = string(:,1)
! enddo
!
! case(2)
! do i1=kept(0),1,-1
! list(elec_alpha_num-0,1) = kept(i1)
! call list_to_bitstring( string(1,1), list(1,1), elec_alpha_num, N_int)
! N_det_alpha_unique += 1
! psi_det_alpha_unique(:,N_det_alpha_unique) = string(:,1)
! do i2=i1-1,1,-1
! list(elec_alpha_num-1,1) = kept(i2)
! call list_to_bitstring( string(1,1), list(1,1), elec_alpha_num, N_int)
! N_det_alpha_unique += 1
! psi_det_alpha_unique(:,N_det_alpha_unique) = string(:,1)
! enddo
! enddo
psi_det_size = kept(0)**(nelec_kept(1)+nelec_kept(2)) psi_det_size = kept(0)**(nelec_kept(1)+nelec_kept(2))
print *, kept(0), nelec_kept(:)
call write_int(6,psi_det_size,'psi_det_size')
TOUCH psi_det_size TOUCH psi_det_size
BEGIN_SHELL [ /usr/bin/python ] BEGIN_SHELL [ /usr/bin/python ]