From c03132b4f94c69949e2f4e809008102c3a179b42 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Thu, 24 Apr 2014 12:31:15 +0200 Subject: [PATCH] Minor changes in Hartree-Fock --- src/BiInts/ao_bi_integrals.irp.f | 1 + src/Hartree_Fock/Fock_matrix_mo.irp.f | 3 +- src/Hartree_Fock/HF_density_matrix_ao.irp.f | 87 +++++++++++++++++++++ src/Hartree_Fock/hartree_fock.ezfio_config | 2 + src/Hartree_Fock/options.irp.f | 39 ++++++++- src/Makefile.common | 2 +- src/MonoInts/ao_mono_ints.irp.f | 14 ++++ src/MonoInts/mo_mono_ints.irp.f | 4 +- 8 files changed, 146 insertions(+), 6 deletions(-) create mode 100644 src/Hartree_Fock/HF_density_matrix_ao.irp.f diff --git a/src/BiInts/ao_bi_integrals.irp.f b/src/BiInts/ao_bi_integrals.irp.f index 904122d6..fbd6fb98 100644 --- a/src/BiInts/ao_bi_integrals.irp.f +++ b/src/BiInts/ao_bi_integrals.irp.f @@ -194,6 +194,7 @@ BEGIN_PROVIDER [ logical, ao_bielec_integrals_in_map ] use map_module BEGIN_DOC ! Map of Atomic integrals + ! i(r1) j(r2) 1/r12 k(r1) l(r2) END_DOC integer :: i,j,k,l diff --git a/src/Hartree_Fock/Fock_matrix_mo.irp.f b/src/Hartree_Fock/Fock_matrix_mo.irp.f index 38b19362..d9e05510 100644 --- a/src/Hartree_Fock/Fock_matrix_mo.irp.f +++ b/src/Hartree_Fock/Fock_matrix_mo.irp.f @@ -1,4 +1,4 @@ - BEGIN_PROVIDER [ double precision, Fock_matrix_mo, (mo_tot_num,mo_tot_num) ] + BEGIN_PROVIDER [ double precision, Fock_matrix_mo, (mo_tot_num_align,mo_tot_num) ] &BEGIN_PROVIDER [ double precision, Fock_matrix_diag_mo, (mo_tot_num)] implicit none BEGIN_DOC @@ -17,7 +17,6 @@ ! END_DOC integer :: i,j,n - double precision :: get_mo_bielec_integral if (elec_alpha_num == elec_beta_num) then Fock_matrix_mo = Fock_matrix_alpha_mo else diff --git a/src/Hartree_Fock/HF_density_matrix_ao.irp.f b/src/Hartree_Fock/HF_density_matrix_ao.irp.f new file mode 100644 index 00000000..b0d58344 --- /dev/null +++ b/src/Hartree_Fock/HF_density_matrix_ao.irp.f @@ -0,0 +1,87 @@ + BEGIN_PROVIDER [ double precision, HF_density_matrix_ao_alpha, (ao_num_align,ao_num) ] +&BEGIN_PROVIDER [ double precision, HF_density_matrix_ao_beta, (ao_num_align,ao_num) ] + implicit none + BEGIN_DOC + ! Alpha and Beta density matrix in the AO basis + END_DOC + integer :: i,j,k,l1,l2 + integer, allocatable :: mo_occ(:,:) + + allocate ( mo_occ(elec_alpha_num,2) ) + call bitstring_to_list( HF_bitmask(1,1), mo_occ(1,1), j, N_int) + ASSERT ( j==elec_alpha_num ) + + call bitstring_to_list( HF_bitmask(1,2), mo_occ(1,2), j, N_int) + ASSERT ( j==elec_beta_num ) + + do j=1,ao_num + !DIR$ VECTOR ALIGNED + do i=1,ao_num_align + HF_density_matrix_ao_alpha(i,j) = 0.d0 + HF_density_matrix_ao_beta (i,j) = 0.d0 + enddo + do k=1,elec_beta_num + l1 = mo_occ(k,1) + l2 = mo_occ(k,2) + !DIR$ VECTOR ALIGNED + do i=1,ao_num + HF_density_matrix_ao_alpha(i,j) = HF_density_matrix_ao_alpha(i,j) +& + mo_coef(i,l1) * mo_coef(j,l1) + HF_density_matrix_ao_beta (i,j) = HF_density_matrix_ao_beta (i,j) +& + mo_coef(i,l2) * mo_coef(j,l2) + enddo + enddo + do k=elec_beta_num+1,elec_alpha_num + l1 = mo_occ(k,1) + !DIR$ VECTOR ALIGNED + do i=1,ao_num + HF_density_matrix_ao_alpha(i,j) = HF_density_matrix_ao_alpha(i,j) +& + mo_coef(i,l1) * mo_coef(j,l1) + enddo + enddo + enddo + deallocate(mo_occ) +END_PROVIDER + +BEGIN_PROVIDER [ double precision, HF_density_matrix_ao, (ao_num_align,ao_num) ] + implicit none + BEGIN_DOC + ! Density matrix in the AO basis + END_DOC + integer :: i,j,k,l1,l2 + integer, allocatable :: mo_occ(:,:) + + allocate ( mo_occ(elec_alpha_num,2) ) + call bitstring_to_list( HF_bitmask(1,1), mo_occ(1,1), j, N_int) + ASSERT ( j==elec_alpha_num ) + + call bitstring_to_list( HF_bitmask(1,2), mo_occ(1,2), j, N_int) + ASSERT ( j==elec_beta_num ) + + do j=1,ao_num + !DIR$ VECTOR ALIGNED + do i=1,ao_num_align + HF_density_matrix_ao(i,j) = 0.d0 + enddo + do k=1,elec_beta_num + l1 = mo_occ(k,1) + l2 = mo_occ(k,2) + !DIR$ VECTOR ALIGNED + do i=1,ao_num + HF_density_matrix_ao(i,j) = HF_density_matrix_ao(i,j) + & + mo_coef(i,l1) * mo_coef(j,l1) + & + mo_coef(i,l2) * mo_coef(j,l2) + enddo + enddo + do k=elec_beta_num+1,elec_alpha_num + l1 = mo_occ(k,1) + !DIR$ VECTOR ALIGNED + do i=1,ao_num + HF_density_matrix_ao(i,j) = HF_density_matrix_ao(i,j) + & + mo_coef(i,l1) * mo_coef(j,l1) + enddo + enddo + enddo + deallocate(mo_occ) +END_PROVIDER + diff --git a/src/Hartree_Fock/hartree_fock.ezfio_config b/src/Hartree_Fock/hartree_fock.ezfio_config index bbb383f9..134486f4 100644 --- a/src/Hartree_Fock/hartree_fock.ezfio_config +++ b/src/Hartree_Fock/hartree_fock.ezfio_config @@ -1,4 +1,6 @@ hartree_fock thresh_scf double precision n_it_scf_max integer + direct logical + diis logical diff --git a/src/Hartree_Fock/options.irp.f b/src/Hartree_Fock/options.irp.f index b3770418..303ce89d 100644 --- a/src/Hartree_Fock/options.irp.f +++ b/src/Hartree_Fock/options.irp.f @@ -16,7 +16,7 @@ BEGIN_PROVIDER [ double precision,thresh_SCF ] END_PROVIDER -BEGIN_PROVIDER [ integer ,n_it_scf_max] +BEGIN_PROVIDER [ integer, n_it_scf_max] implicit none BEGIN_DOC ! Maximum number of SCF iterations @@ -34,3 +34,40 @@ BEGIN_PROVIDER [ integer ,n_it_scf_max] END_PROVIDER + +BEGIN_PROVIDER [ logical, do_direct_SCF ] + implicit none + BEGIN_DOC +! If True, compute integrals on the fly + END_DOC + + logical :: has + PROVIDE ezfio_filename + call ezfio_has_Hartree_Fock_direct(has) + if (has) then + call ezfio_get_Hartree_Fock_direct(do_direct_SCF) + else + do_direct_SCF = .False. + call ezfio_set_Hartree_Fock_direct(do_direct_SCF) + endif + +END_PROVIDER + +BEGIN_PROVIDER [ logical, do_DIIS ] + implicit none + BEGIN_DOC +! If True, compute integrals on the fly + END_DOC + + logical :: has + PROVIDE ezfio_filename + call ezfio_has_Hartree_Fock_DIIS(has) + if (has) then + call ezfio_get_Hartree_Fock_DIIS(do_DIIS) + else + do_DIIS = .False. + call ezfio_set_Hartree_Fock_DIIS(do_DIIS) + endif + +END_PROVIDER + diff --git a/src/Makefile.common b/src/Makefile.common index 4c2b2b2b..2dad85e2 100644 --- a/src/Makefile.common +++ b/src/Makefile.common @@ -113,7 +113,7 @@ clean_links: endif LIB+=$(EZFIO) $(MKL) -IRPF90+=$(patsubst %, -I %, $(INCLUDE_DIRS)) +IRPF90+=$(patsubst %, -I %, $(INCLUDE_DIRS)) $(IRPF90_FLAGS) irpf90.make: $(filter-out IRPF90_temp/%, $(wildcard */*.irp.f)) $(wildcard *.irp.f) $(wildcard *.inc.f) Makefile $(EZFIO) NEEDED_MODULES $(IRPF90) diff --git a/src/MonoInts/ao_mono_ints.irp.f b/src/MonoInts/ao_mono_ints.irp.f index f77924d3..9a135dda 100644 --- a/src/MonoInts/ao_mono_ints.irp.f +++ b/src/MonoInts/ao_mono_ints.irp.f @@ -119,3 +119,17 @@ BEGIN_PROVIDER [ double precision, ao_overlap_abs,(ao_num_align,ao_num) ] !$OMP END PARALLEL DO END_PROVIDER +BEGIN_PROVIDER [ double precision, ao_mono_elec_integral,(ao_num_align,ao_num)] + implicit none + integer :: i,j,n,l + BEGIN_DOC + ! array of the mono electronic hamiltonian on the AOs basis + ! : sum of the kinetic and nuclear electronic potential + END_DOC + do j = 1, ao_num + do i = 1, ao_num + ao_mono_elec_integral(i,j) = ao_nucl_elec_integral(i,j) + ao_kinetic_integral(i,j) + enddo + enddo +END_PROVIDER + diff --git a/src/MonoInts/mo_mono_ints.irp.f b/src/MonoInts/mo_mono_ints.irp.f index 956c6454..e466c9cc 100644 --- a/src/MonoInts/mo_mono_ints.irp.f +++ b/src/MonoInts/mo_mono_ints.irp.f @@ -39,8 +39,8 @@ BEGIN_PROVIDER [ double precision, mo_mono_elec_integral,(mo_tot_num_align,mo_to ! array of the mono electronic hamiltonian on the MOs basis ! : sum of the kinetic and nuclear electronic potential END_DOC - do i = 1, mo_tot_num - do j = 1, mo_tot_num + do j = 1, mo_tot_num + do i = 1, mo_tot_num mo_mono_elec_integral(i,j) = mo_nucl_elec_integral(i,j) + mo_kinetic_integral(i,j) enddo enddo