From 6a41ec1da4b46eecb6d217a25cc3216d3a76cfee Mon Sep 17 00:00:00 2001 From: eginer Date: Mon, 21 Oct 2019 19:19:26 +0200 Subject: [PATCH] casscf really good --- src/casscf/NEED | 2 +- src/casscf/casscf.irp.f | 85 +--------------- src/casscf/change_bitmasks.irp.f | 14 --- src/casscf/cipsi_routines.irp.f | 75 --------------- src/casscf/cisd_routine.irp.f | 85 ---------------- src/casscf/cisdtq_routine.irp.f | 47 --------- src/casscf/densities.irp.f | 2 +- src/casscf/gradient.irp.f | 2 - src/casscf/h_apply.irp.f | 18 ---- src/casscf/hessian.irp.f | 71 ++++---------- src/generators_fluid/NEED | 1 - src/generators_fluid/README.rst | 0 src/generators_fluid/extract_cas.irp.f | 23 ----- src/generators_fluid/generators.irp.f | 101 -------------------- src/generators_fluid/generators_cas.irp.f | 69 ------------- src/generators_fluid/generators_hf.irp.f | 51 ---------- src/generators_fluid/generators_hf_sd.irp.f | 80 ---------------- 17 files changed, 26 insertions(+), 700 deletions(-) delete mode 100644 src/casscf/change_bitmasks.irp.f delete mode 100644 src/casscf/cipsi_routines.irp.f delete mode 100644 src/generators_fluid/NEED delete mode 100644 src/generators_fluid/README.rst delete mode 100644 src/generators_fluid/extract_cas.irp.f delete mode 100644 src/generators_fluid/generators.irp.f delete mode 100644 src/generators_fluid/generators_cas.irp.f delete mode 100644 src/generators_fluid/generators_hf.irp.f delete mode 100644 src/generators_fluid/generators_hf_sd.irp.f diff --git a/src/casscf/NEED b/src/casscf/NEED index b992ff71..d9da718e 100644 --- a/src/casscf/NEED +++ b/src/casscf/NEED @@ -1,4 +1,4 @@ cipsi selectors_full -generators_fluid +generators_cas two_body_rdm diff --git a/src/casscf/casscf.irp.f b/src/casscf/casscf.irp.f index 8b5a365f..8fe77fcc 100644 --- a/src/casscf/casscf.irp.f +++ b/src/casscf/casscf.irp.f @@ -4,86 +4,9 @@ program casscf ! TODO : Put the documentation of the program here END_DOC no_vvvv_integrals = .True. - SOFT_TOUCH no_vvvv_integrals - threshold_davidson = 1.d-7 - 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 + pt2_max = 0.02 + SOFT_TOUCH no_vvvv_integrals pt2_max + call run end subroutine run @@ -107,7 +30,7 @@ subroutine run call write_double(6,energy_improvement, 'Predicted energy improvement') 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 call save_mos diff --git a/src/casscf/change_bitmasks.irp.f b/src/casscf/change_bitmasks.irp.f deleted file mode 100644 index cad6ec38..00000000 --- a/src/casscf/change_bitmasks.irp.f +++ /dev/null @@ -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 - diff --git a/src/casscf/cipsi_routines.irp.f b/src/casscf/cipsi_routines.irp.f deleted file mode 100644 index 272a7116..00000000 --- a/src/casscf/cipsi_routines.irp.f +++ /dev/null @@ -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 diff --git a/src/casscf/cisd_routine.irp.f b/src/casscf/cisd_routine.irp.f index a4cbfcfb..e69de29b 100644 --- a/src/casscf/cisd_routine.irp.f +++ b/src/casscf/cisd_routine.irp.f @@ -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 diff --git a/src/casscf/cisdtq_routine.irp.f b/src/casscf/cisdtq_routine.irp.f index 0479d462..e69de29b 100644 --- a/src/casscf/cisdtq_routine.irp.f +++ b/src/casscf/cisdtq_routine.irp.f @@ -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 - diff --git a/src/casscf/densities.irp.f b/src/casscf/densities.irp.f index 292067b4..3d1ff0f9 100644 --- a/src/casscf/densities.irp.f +++ b/src/casscf/densities.irp.f @@ -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) do t = 1, n_act_orb 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) = act_two_rdm_spin_trace_mo(t,v,u,x) enddo enddo enddo diff --git a/src/casscf/gradient.irp.f b/src/casscf/gradient.irp.f index f00bc7c8..6bf8b93b 100644 --- a/src/casscf/gradient.irp.f +++ b/src/casscf/gradient.irp.f @@ -171,11 +171,9 @@ BEGIN_PROVIDER [real*8, gradvec2, (nMonoEx)] norm_grad+=gradvec2(indx)*gradvec2(indx) end do norm_grad=sqrt(norm_grad) -! if (bavard) then write(6,*) write(6,*) ' Norm of the orbital gradient (via D, P and integrals): ', norm_grad write(6,*) -! endif END_PROVIDER diff --git a/src/casscf/h_apply.irp.f b/src/casscf/h_apply.irp.f index 6fcb2900..e69de29b 100644 --- a/src/casscf/h_apply.irp.f +++ b/src/casscf/h_apply.irp.f @@ -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 - diff --git a/src/casscf/hessian.irp.f b/src/casscf/hessian.irp.f index 06aed6ef..52be1b76 100644 --- a/src/casscf/hessian.irp.f +++ b/src/casscf/hessian.irp.f @@ -536,9 +536,6 @@ real*8 function hessmat_taub(t,a,u,b) integer :: v3,x3 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) aa=list_virt(a) if (t == u) then @@ -548,87 +545,59 @@ real*8 function hessmat_taub(t,a,u,b) t2=0.D0 t3=0.D0 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 vv=list_act(v) v3=v+n_core_inact_orb do x=1,n_act_orb xx=list_act(x) x3=x+n_core_inact_orb - t2+=(P0tuvx_no(t,x,v,t)+P0tuvx_no(t,x,t,v))* & - bielec_pxxq_no(aa,x3,v3,aa) - end do - end do - do y=1,n_act_orb - do x=1,n_act_orb - xx=list_act(x) - do v=1,n_act_orb - t3-=P0tuvx_no(t,v,x,y)*bielecCI_no(t,v,y,xx) + t2+=2.D0*(P0tuvx_no(t,t,v,x)*bielec_pqxx_no(aa,aa,v3,x3) & + +(P0tuvx_no(t,x,v,t)+P0tuvx_no(t,x,t,v))* & + bielec_pxxq_no(aa,x3,v3,aa)) + do y=1,n_act_orb + t3-=2.D0*P0tuvx_no(t,v,x,y)*bielecCI_no(t,v,y,xx) end do end do end do - term=t1+2.d0*(t2+t3) + term=t1+t2+t3 else bb=list_virt(b) ! ta/tb b/=a - term=0.5d0*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 + term=occnum(tt)*Fipq(aa,bb) do v=1,n_act_orb vv=list_act(v) v3=v+n_core_inact_orb do x=1,n_act_orb xx=list_act(x) x3=x+n_core_inact_orb - term= term + (P0tuvx_no(t,x,v,t)+P0tuvx_no(t,x,t,v)) & - *bielec_pxxq_no(aa,x3,v3,bb) + term+=2.D0*(P0tuvx_no(t,t,v,x)*bielec_pqxx_no(aa,bb,v3,x3) & + +(P0tuvx_no(t,x,v,t)+P0tuvx_no(t,x,t,v)) & + *bielec_pxxq_no(aa,x3,v3,bb)) end do end do - term += term end if else ! ta/ub t/=u uu=list_act(u) 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 - do x=1,n_act_orb - do v=1,n_act_orb - term += P0tuvx_no_t(u,v,x)*bielec_pqxx_no_2(v,x) - term += bielec_pxxq_no_2(x,v) * (P0tuvx_no_t(x,v,u)+P0tuvx_no_t(x,u,v)) + do v=1,n_act_orb + vv=list_act(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 - term = 6.d0*term if (a.eq.b) then term-=0.5D0*(occnum(tt)*Fipq(uu,tt)+occnum(uu)*Fipq(tt,uu)) do v=1,n_act_orb do y=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) end do end do diff --git a/src/generators_fluid/NEED b/src/generators_fluid/NEED deleted file mode 100644 index d3d4d2c7..00000000 --- a/src/generators_fluid/NEED +++ /dev/null @@ -1 +0,0 @@ -determinants diff --git a/src/generators_fluid/README.rst b/src/generators_fluid/README.rst deleted file mode 100644 index e69de29b..00000000 diff --git a/src/generators_fluid/extract_cas.irp.f b/src/generators_fluid/extract_cas.irp.f deleted file mode 100644 index 9cdaf27f..00000000 --- a/src/generators_fluid/extract_cas.irp.f +++ /dev/null @@ -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 diff --git a/src/generators_fluid/generators.irp.f b/src/generators_fluid/generators.irp.f deleted file mode 100644 index 153ab605..00000000 --- a/src/generators_fluid/generators.irp.f +++ /dev/null @@ -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 - diff --git a/src/generators_fluid/generators_cas.irp.f b/src/generators_fluid/generators_cas.irp.f deleted file mode 100644 index b6d83e0a..00000000 --- a/src/generators_fluid/generators_cas.irp.f +++ /dev/null @@ -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 - diff --git a/src/generators_fluid/generators_hf.irp.f b/src/generators_fluid/generators_hf.irp.f deleted file mode 100644 index d4d2e728..00000000 --- a/src/generators_fluid/generators_hf.irp.f +++ /dev/null @@ -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 diff --git a/src/generators_fluid/generators_hf_sd.irp.f b/src/generators_fluid/generators_hf_sd.irp.f deleted file mode 100644 index 9c13a5a0..00000000 --- a/src/generators_fluid/generators_hf_sd.irp.f +++ /dev/null @@ -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 -