Introducing trexion convention

bin/qp_convert_output_to_ezfio

@@ -44,6 +44,7 @@ def write_ezfio(res, filename):
 
     res.clean_uncontractions()
     ezfio.set_file(filename)
+    ezfio.set_ezfio_files_ezfio_convention(20250211)
 
     #  _
     # |_ |  _   _ _|_ ._ _  ._   _
@@ -172,6 +173,7 @@ def write_ezfio(res, filename):
     # ~#~#~#~#~ #
 
     ezfio.set_ao_basis_ao_coef(coef)
+    ezfio.set_basis_ao_normalized(True)
     ezfio.set_ao_basis_ao_expo(expo)
     ezfio.set_ao_basis_ao_basis("Read by resultsFile")
 

external/ezfio

@@ -1 +1 @@
-Subproject commit dba01c4fe0ff7b84c5ecfb1c7c77ec68781311b3
+Subproject commit d02132ea79217c16fd24242e8f8b8a6c3ff68091

ocaml/qp_create_ezfio.ml

@@ -265,6 +265,7 @@ let run ?o b au c d m p cart xyz_file =
   let write_file () =
       (* Create EZFIO *)
       Ezfio.set_file ezfio_file;
+      Ezfio.set_files_ezfio_ezfio_convention 20250211;
 
       (* Write Pseudo *)
       let pseudo =

scripts/qp_import_trexio.py

@@ -84,6 +84,7 @@ def write_ezfio(trexio_filename, filename):
 
     ezfio.set_file(filename)
     ezfio.set_trexio_trexio_file(trexio_filename)
+    ezfio.set_ezfio_files_ezfio_convention(20250211)
 
     print("Nuclei\t\t...\t", end=' ')
 
@@ -315,8 +316,8 @@ def write_ezfio(trexio_filename, filename):
                 power_x.append(x)
                 power_y.append(y)
                 power_z.append(z)
-                coefficient.append(coef[i])
+                coefficient.append(list(coef[i]))
-                exponent.append(expo[i])
+                exponent.append(list(expo[i]))
                 num_prim.append(num_prim0[i])
 
         assert (len(coefficient) == ao_num)
@@ -326,15 +327,15 @@ def write_ezfio(trexio_filename, filename):
 
         prim_num_max = max( [ len(x) for x in coefficient ] )
 
-        ao_normalization = trexio.read_ao_normalization(trexio_file_cart)
-        for i, coef in enumerate(coefficient):
-           for j in range(len(coef)):
-              coef[j] *= ao_normalization[i]
-
         for i in range(ao_num):
             coefficient[i] += [0. for j in range(len(coefficient[i]), prim_num_max)]
             exponent   [i] += [0. for j in range(len(exponent[i]), prim_num_max)]
 
+        ao_normalization = trexio.read_ao_normalization(trexio_file_cart)
+        for i in range(ao_num):
+           for j in range(prim_num_max):
+               coefficient[i][j] *= ao_normalization[i]
+
         coefficient = reduce(lambda x, y: x + y, coefficient, [])
         exponent    = reduce(lambda x, y: x + y, exponent   , [])
 
@@ -345,6 +346,7 @@ def write_ezfio(trexio_filename, filename):
                 coef.append(coefficient[j])
                 expo.append(exponent[j])
 
+
         ezfio.set_ao_basis_ao_coef(coef)
         ezfio.set_ao_basis_ao_expo(expo)
 

src/ao_basis/aos.irp.f

@@ -53,42 +53,84 @@ END_PROVIDER
   C_A(3) = 0.d0
   ao_coef_normalized = 0.d0
 
-  do i=1,ao_num
+  if (primitives_normalized) then
 
-    powA(1) = ao_power(i,1)
+    if (ezfio_convention >= 20250211) then
-    powA(2) = ao_power(i,2)
+      ! Same primitive normalization factors for all AOs of the same shell, or read from trexio file
-    powA(3) = ao_power(i,3)
 
-    ! Normalization of the primitives
+      do i=1,ao_num
-    if (primitives_normalized) then
+        k=1
-      do j=1,ao_prim_num(i)
+        do while (k<=prim_num .and. shell_index(k) /= ao_shell(i))
-        call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j), &
+          k = k+1
-           powA,powA,overlap_x,overlap_y,overlap_z,norm,nz)
+        end do
-        ao_coef_normalized(i,j) = ao_coef(i,j)/dsqrt(norm)
+        do j=1,ao_prim_num(i)
+          ao_coef_normalized(i,j) = ao_coef(i,j)*prim_normalization_factor(k+j-1)
+        enddo
       enddo
+
     else
+        ! GAMESS convention for primitive factors
+
+      do i=1,ao_num
+        powA(1) = ao_power(i,1)
+        powA(2) = ao_power(i,2)
+        powA(3) = ao_power(i,3)
+
+        do j=1,ao_prim_num(i)
+          call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j), &
+             powA,powA,overlap_x,overlap_y,overlap_z,norm,nz)
+          ao_coef_normalized(i,j) = ao_coef(i,j)/dsqrt(norm)
+        enddo
+      enddo
+
+    endif
+
+  else
+
+    do i=1,ao_num
       do j=1,ao_prim_num(i)
         ao_coef_normalized(i,j) = ao_coef(i,j)
       enddo
-    endif
-
-    ! Normalization of the contracted basis functions
-    norm = 0.d0
-    do j=1,ao_prim_num(i)
-      do k=1,ao_prim_num(i)
-        call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,k),powA,powA,overlap_x,overlap_y,overlap_z,c,nz)
-        norm = norm+c*ao_coef_normalized(i,j)*ao_coef_normalized(i,k)
-      enddo
     enddo
-    ao_coef_normalization_factor(i) = 1.d0/dsqrt(norm)
 
-    if (ao_normalized) then
+  endif
-      do j=1,ao_prim_num(i)
+
-        ao_coef_normalized(i,j) = ao_coef_normalized(i,j) * ao_coef_normalization_factor(i)
+  double precision, allocatable :: self_overlap(:)
+  allocate(self_overlap(ao_num))
+
+  do i=1,ao_num
+    powA(1) = ao_power(i,1)
+    powA(2) = ao_power(i,2)
+    powA(3) = ao_power(i,3)
+    self_overlap(i) = 0.d0
+    do j=1,ao_prim_num(i)
+      do k=1,j-1
+        call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,k),powA,powA,overlap_x,overlap_y,overlap_z,c,nz)
+        self_overlap(i) = self_overlap(i) + 2.d0*c*ao_coef_normalized(i,j)*ao_coef_normalized(i,k)
       enddo
-    else
+      call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j),powA,powA,overlap_x,overlap_y,overlap_z,c,nz)
+      self_overlap(i) = self_overlap(i) +c*ao_coef_normalized(i,j)*ao_coef_normalized(i,j)
+    enddo
+  enddo
+
+  if (ao_normalized) then
+
+    do i=1,ao_num
+      ao_coef_normalization_factor(i) = 1.d0/dsqrt(self_overlap(i))
+    enddo
+
+  else
+
+    do i=1,ao_num
       ao_coef_normalization_factor(i) = 1.d0
-    endif
+    enddo
+
+  endif
+
+  do i=1,ao_num
+    do j=1,ao_prim_num(i)
+      ao_coef_normalized(i,j) = ao_coef_normalized(i,j) * ao_coef_normalization_factor(i)
+    enddo
   enddo
 
 END_PROVIDER

src/basis/EZFIO.cfg

@@ -84,4 +84,3 @@ type: logical
 doc: If true, normalize the basis functions
 interface: ezfio, provider, ocaml
 default: false
-

src/basis/basis.irp.f

@@ -1,72 +1,3 @@
-BEGIN_PROVIDER [ double precision, shell_normalization_factor , (shell_num) ]
-  implicit none
-  BEGIN_DOC
-  ! Number of primitives per |AO|
-  END_DOC
-
-  logical                        :: has
-  PROVIDE ezfio_filename
-  if (.not.ao_normalized) then
-    shell_normalization_factor = 1.d0
-    return
-  endif
-
-  if (mpi_master) then
-    if (size(shell_normalization_factor) == 0) return
-
-    call ezfio_has_basis_shell_normalization_factor(has)
-    if (has) then
-      write(6,'(A)') '.. >>>>> [ IO READ: shell_normalization_factor ] <<<<< ..'
-      call ezfio_get_basis_shell_normalization_factor(shell_normalization_factor)
-    else
-
-      double precision               :: norm,overlap_x,overlap_y,overlap_z,C_A(3), c
-      integer                        :: l, powA(3), nz
-      integer                        :: i,j,k
-      nz=100
-      C_A(1) = 0.d0
-      C_A(2) = 0.d0
-      C_A(3) = 0.d0
-
-      do i=1,shell_num
-
-        powA(1) = shell_ang_mom(i)
-        powA(2) = 0
-        powA(3) = 0
-
-        norm = 0.d0
-        do k=1, prim_num
-          if (shell_index(k) /= i) cycle
-          do j=1, prim_num
-            if (shell_index(j) /= i) cycle
-            call overlap_gaussian_xyz(C_A,C_A,prim_expo(j),prim_expo(k), &
-              powA,powA,overlap_x,overlap_y,overlap_z,c,nz)
-            norm = norm+c*prim_coef(j)*prim_coef(k) * prim_normalization_factor(j) * prim_normalization_factor(k)
-          enddo
-        enddo
-        shell_normalization_factor(i) = 1.d0/dsqrt(norm)
-      enddo
-
-    endif
-  endif
-  IRP_IF MPI_DEBUG
-  print *,  irp_here, mpi_rank
-  call MPI_BARRIER(MPI_COMM_WORLD, ierr)
-  IRP_ENDIF
-  IRP_IF MPI
-  include 'mpif.h'
-  integer                        :: ierr
-  call MPI_BCAST( shell_normalization_factor, (shell_num), MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr)
-  if (ierr /= MPI_SUCCESS) then
-    stop 'Unable to read shell_normalization_factor with MPI'
-  endif
-  IRP_ENDIF
-
-  call write_time(6)
-
-END_PROVIDER
-
-
 BEGIN_PROVIDER [ double precision, prim_normalization_factor , (prim_num) ]
   implicit none
   BEGIN_DOC

src/trexio/export_trexio_routines.irp.f

@@ -15,6 +15,8 @@ subroutine export_trexio(update,full_path)
 
   integer, external :: getunitandopen
 
+  integer :: i,j,l
+
   if (full_path) then
     fp = trexio_filename
     call system('realpath '//trim(fp)//' > '//trim(fp)//'.tmp')
@@ -271,7 +273,7 @@ subroutine export_trexio(update,full_path)
      call trexio_assert(rc, TREXIO_SUCCESS)
 
      allocate(factor(shell_num))
-     factor(1:shell_num) = shell_normalization_factor(1:shell_num)
+     factor(1:shell_num) = 1.d0
      rc = trexio_write_basis_shell_factor(f(1), factor)
      call trexio_assert(rc, TREXIO_SUCCESS)
 
@@ -312,22 +314,27 @@ subroutine export_trexio(update,full_path)
     rc = trexio_write_ao_shell(f(1), ao_shell)
     call trexio_assert(rc, TREXIO_SUCCESS)
 
-    integer :: i, pow0(3), powA(3), j, l, nz
+    if (ezfio_convention >= 20250211) then
-    double precision :: normA, norm0, C_A(3), overlap_x, overlap_z, overlap_y, c
+      rc = trexio_write_ao_normalization(f(1), ao_coef_normalization_factor)
-    nz=100
+      print *, ao_coef_normalization_factor(:)
+    else
+      integer :: pow0(3), powA(3), nz
+      double precision :: normA, norm0, C_A(3), overlap_x, overlap_z, overlap_y, c
+      nz=100
 
-    C_A(1) = 0.d0
+      C_A(1) = 0.d0
-    C_A(2) = 0.d0
+      C_A(2) = 0.d0
-    C_A(3) = 0.d0
+      C_A(3) = 0.d0
 
-    allocate(factor(ao_num))
+      allocate(factor(ao_num))
-    do i=1,ao_num
+      do i=1,ao_num
-      l = ao_first_of_shell(ao_shell(i))
+        l = ao_first_of_shell(ao_shell(i))
-      factor(i) = (ao_coef_normalized(i,1)+tiny(1.d0))/(ao_coef_normalized(l,1)+tiny(1.d0))
+        factor(i) = (ao_coef_normalized(i,1)+tiny(1.d0))/(ao_coef_normalized(l,1)+tiny(1.d0))
-    enddo
+      enddo
-    rc = trexio_write_ao_normalization(f(1), factor)
+      rc = trexio_write_ao_normalization(f(1), factor)
+      deallocate(factor)
+    endif
     call trexio_assert(rc, TREXIO_SUCCESS)
-    deallocate(factor)
 
   endif