diff --git a/scripts/update_README.py b/scripts/update_README.py index 00e038b7..fc769db7 100755 --- a/scripts/update_README.py +++ b/scripts/update_README.py @@ -95,18 +95,14 @@ def update_needed(data): return data -import subprocess +import os def git_add(): """Executes: git add README.rst if git is present on the machine.""" command = "git add "+README - - try: - subprocess.call(command.split()) - except OSError: - pass + os.system(command+" &> /dev/null") def main(): diff --git a/src/Hartree_Fock/Fock_matrix_mo.irp.f b/src/Hartree_Fock/Fock_matrix_mo.irp.f index 9b5c193d..04d2f1cc 100644 --- a/src/Hartree_Fock/Fock_matrix_mo.irp.f +++ b/src/Hartree_Fock/Fock_matrix_mo.irp.f @@ -94,7 +94,7 @@ END_PROVIDER double precision :: ao_bielec_integral !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: ao_ints_idx, ao_ints_val if (do_direct_integrals) then - PROVIDE all_utils + PROVIDE all_utils ao_overlap_abs ao_integrals_threshold !$OMP PARALLEL DEFAULT(NONE) & !$OMP PRIVATE(i,j,l,k1,k,integral) & !$OMP SHARED(ao_num,Fock_matrix_alpha_ao,ao_mono_elec_integral,& @@ -107,6 +107,7 @@ END_PROVIDER Fock_matrix_beta_ao (i,j) = ao_mono_elec_integral(i,j) do l=1,ao_num do k=1,ao_num + PROVIDE HF_density_matrix_ao_alpha HF_density_matrix_ao_beta if ((abs(HF_density_matrix_ao_alpha(k,l)) > 1.d-9).or. & (abs(HF_density_matrix_ao_beta (k,l)) > 1.d-9)) then integral = (HF_density_matrix_ao_alpha(k,l)+HF_density_matrix_ao_beta (k,l)) * ao_bielec_integral(k,l,i,j) @@ -207,3 +208,12 @@ BEGIN_PROVIDER [ double precision, Fock_matrix_beta_mo, (mo_tot_num_align,mo_tot deallocate(T) END_PROVIDER +BEGIN_PROVIDER [ double precision, HF_energy ] + implicit none + BEGIN_DOC + ! Hartree-Fock energy + END_DOC + HF_energy = nuclear_repulsion + ref_bitmask_energy +END_PROVIDER + + diff --git a/src/Hartree_Fock/diagonalize_fock.irp.f b/src/Hartree_Fock/diagonalize_fock.irp.f index 00f3275b..c34a1415 100644 --- a/src/Hartree_Fock/diagonalize_fock.irp.f +++ b/src/Hartree_Fock/diagonalize_fock.irp.f @@ -1,20 +1,20 @@ -subroutine diagonalize_fock() - implicit none - - double precision, allocatable :: mo_coef_new(:,:), R(:,:),eigvalues(:) - !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: mo_coef_new, R - - allocate(R(size(Fock_matrix_mo,1),mo_tot_num)) - allocate(mo_coef_new(ao_num_align,mo_tot_num),eigvalues(mo_tot_num)) - mo_coef_new = mo_coef - - call lapack_diag(eigvalues,R,Fock_matrix_mo,size(Fock_matrix_mo,1),mo_tot_num) - - call dgemm('N','N',ao_num,mo_tot_num,mo_tot_num,1.d0,mo_coef_new,size(mo_coef_new,1),R,size(R,1),0.d0,mo_coef,size(mo_coef,1)) - deallocate(mo_coef_new,R,eigvalues) - - mo_label = "Canonical" - SOFT_TOUCH mo_coef mo_label - call clear_mo_map -end - + BEGIN_PROVIDER [ double precision, diagonal_Fock_matrix_mo, (mo_tot_num) ] +&BEGIN_PROVIDER [ double precision, eigenvectors_Fock_matrix_mo, (ao_num_align,mo_tot_num) ] + implicit none + BEGIN_DOC + ! Diagonal Fock matrix in the MO basis + END_DOC + + double precision, allocatable :: R(:,:) + !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: R + + allocate(R(mo_tot_num_align,mo_tot_num)) + + call lapack_diag(diagonal_Fock_matrix_mo,R,Fock_matrix_mo,size(Fock_matrix_mo,1),& + mo_tot_num) + + call dgemm('N','N',ao_num,mo_tot_num,mo_tot_num,1.d0,mo_coef,size(mo_coef,1),& + R,size(R,1),0.d0,eigenvectors_Fock_matrix_mo,size(eigenvectors_Fock_matrix_mo,1)) + deallocate(R) +END_PROVIDER + diff --git a/src/Hartree_Fock/mo_SCF_iterations.irp.f b/src/Hartree_Fock/mo_SCF_iterations.irp.f index c8f3d113..61969f0b 100644 --- a/src/Hartree_Fock/mo_SCF_iterations.irp.f +++ b/src/Hartree_Fock/mo_SCF_iterations.irp.f @@ -5,16 +5,20 @@ program scf_iteration double precision :: E0 integer :: i_it - E0 = ref_bitmask_energy + nuclear_repulsion + E0 = HF_energy i_it = 0 n_it_scf_max = 100 SCF_energy_before = huge(1.d0) SCF_energy_after = E0 print *, E0 - do while (dabs(SCF_energy_before - SCF_energy_after) > thresh_SCF) + mo_label = "Canonical" + thresh_SCF = 1.d-10 + do while (i_it < 10 .or. dabs(SCF_energy_before - SCF_energy_after) > thresh_SCF) SCF_energy_before = SCF_energy_after - call diagonalize_fock() - SCF_energy_after = ref_bitmask_energy + nuclear_repulsion + mo_coef = eigenvectors_Fock_matrix_mo + TOUCH mo_coef mo_label + call clear_mo_map + SCF_energy_after = HF_energy print*,SCF_energy_after i_it +=1 if(i_it > n_it_scf_max)exit @@ -28,6 +32,6 @@ program scf_iteration endif mo_label = "Canonical" TOUCH mo_label mo_coef -! call save_mos + call save_mos end diff --git a/src/MOGuess/README.rst b/src/MOGuess/README.rst index 14ade90f..7a81c93d 100644 --- a/src/MOGuess/README.rst +++ b/src/MOGuess/README.rst @@ -2,3 +2,17 @@ MOGuess Module ============== +Needed Modules +============== + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +* `AOs `_ +* `Ezfio_files `_ +* `MOs `_ +* `Nuclei `_ +* `Output `_ +* `Utils `_ +* `MonoInts `_ + diff --git a/src/MOGuess/localize.irp.f b/src/MOGuess/localize.irp.f deleted file mode 100644 index 53790e02..00000000 --- a/src/MOGuess/localize.irp.f +++ /dev/null @@ -1,114 +0,0 @@ -!TODO Ecrire un cholesky avec bitmask - - -subroutine localize_mos(mask, nint) - implicit none - use bitmasks - integer, intent(in) :: nint - integer(bit_kind), intent(in) :: mask(nint) - integer :: i,j,k,l - double precision, allocatable :: DM(:,:) - double precision, allocatable :: mo_coef_new(:,:), R(:,:) - integer :: n - !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: DM, mo_coef_new, R - integer :: rank - integer, parameter :: n_core = 2 - - allocate(R(mo_tot_num,mo_tot_num)) - allocate(DM(ao_num_align,ao_num)) - allocate(mo_coef_new(ao_num_align,mo_tot_num)) - n = ao_num - mo_coef_new = mo_coef - -BEGIN_TEMPLATE - DM = 0.d0 - if ($START < $END) then - do k=$START, $END - do j=1,n - !DEC$ VECTOR ALIGNED - do i=1,n - DM(i,j) += mo_coef_new(i,k)*mo_coef_new(j,k) - enddo - enddo - enddo - call cholesky_mo(n,$END-$START+1,DM,mo_coef_new(1,$START),size(mo_coef_new,1),-1.d0,rank) - endif -SUBST [ START, END ] - 1 ; n_core ;; -END_TEMPLATE - - deallocate(DM) - call find_rotation(mo_coef,ao_num_align,mo_coef_new,ao_num,R,mo_tot_num) - mo_coef = mo_coef_new - deallocate(mo_coef_new) - - double precision,allocatable :: mo_energy_new(:) - integer, allocatable :: iorder(:) - allocate(mo_energy_new(mo_tot_num),iorder(mo_tot_num)) - - do i=1,mo_tot_num - iorder(i) = i - mo_energy_new(i) = 0.d0 - do k=1,mo_tot_num - mo_energy_new(i) += R(k,i)*R(k,i)*mo_energy(k) - enddo - enddo - mo_energy = mo_energy_new - call dsort(mo_energy(1),iorder(1),n_core) - allocate (mo_coef_new(ao_num_align,mo_tot_num)) - mo_coef_new = mo_coef - do j=1,mo_tot_num - do i=1,ao_num - mo_coef(i,j) = mo_coef_new(i,iorder(j)) - enddo - enddo - deallocate (mo_coef_new,R) - deallocate(mo_energy_new,iorder) - mo_label = 'localized' - - SOFT_TOUCH mo_coef mo_energy mo_label -end - - - - - -subroutine cholesky_mo(n,m,P,C,LDC,tol_in,rank) - implicit none - integer, intent(in) :: n,m, LDC - double precision, intent(in) :: P(LDC,n) - double precision, intent(out) :: C(LDC,m) - double precision, intent(in) :: tol_in - integer, intent(out) :: rank - - integer :: info - integer :: i,k - integer :: ipiv(n) - double precision:: tol - double precision, allocatable :: W(:,:), work(:) - !DEC$ ATTRIBUTES ALIGN: 32 :: W - !DEC$ ATTRIBUTES ALIGN: 32 :: work - !DEC$ ATTRIBUTES ALIGN: 32 :: ipiv - - allocate(W(LDC,n),work(2*n)) - tol=tol_in - - info = 0 - do i=1,n - do k=1,i - W(i,k) = P(i,k) - enddo - do k=i+1,n - W(i,k) = 0. - enddo - enddo - call DPSTRF('L', n, W, LDC, ipiv, rank, tol, work, info ) - do i=1,n - do k=1,min(m,rank) - C(ipiv(i),k) = W(i,k) - enddo - enddo - - deallocate(W,work) -end -