Corrected OpenMP bugs

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():

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
+
+

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
+ 

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

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 <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
+* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
+* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
+* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
+* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
+* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
+* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
+

src/MOGuess/localize.irp.f

@@ -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
-