Fix AXK implementation and extend it to SOSEX (cp2k#3468)

hfp · Jun 12, 2024 · 578b4e2 · 578b4e2
1 parent 014179d
commit 578b4e2
Show file tree

Hide file tree

Showing 47 changed files with 2,230 additions and 1,445 deletions.
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
@@ -710,7 +710,7 @@ list(
   restraint.F
   ri_environment_methods.F
   rmsd.F
-  rpa_axk.F
+  rpa_exchange.F
   rpa_communication.F
   rpa_grad.F
   rpa_gw.F

diff --git a/src/aobasis/basis_set_types.F b/src/aobasis/basis_set_types.F
@@ -45,7 +45,8 @@ MODULE basis_set_types
                                               nsoset
    USE orbital_symbols,                 ONLY: cgf_symbol,&
                                               sgf_symbol
-   USE orbital_transformation_matrices, ONLY: orbtramat
+   USE orbital_transformation_matrices, ONLY: init_spherical_harmonics,&
+                                              orbtramat
    USE sto_ng,                          ONLY: get_sto_ng
    USE string_utilities,                ONLY: integer_to_string,&
                                               remove_word,&
@@ -822,6 +823,15 @@ SUBROUTINE init_cphi_and_sphi(gto_basis_set)
 
       gto_basis_set%sphi = 0.0_dp
       IF (n .GT. 0) THEN
+         lmax = -1
+         ! Ensure proper setup of orbtramat
+         DO iset = 1, gto_basis_set%nset
+            DO ishell = 1, gto_basis_set%nshell(iset)
+               lmax = MAX(lmax, gto_basis_set%l(ishell, iset))
+            END DO
+         END DO
+         CALL init_spherical_harmonics(lmax, -1)
+
          DO iset = 1, gto_basis_set%nset
             DO ishell = 1, gto_basis_set%nshell(iset)
                l = gto_basis_set%l(ishell, iset)

diff --git a/src/bse_util.F b/src/bse_util.F
@@ -20,8 +20,6 @@ MODULE bse_util
                                               cp_fm_struct_type
    USE cp_fm_types,                     ONLY: cp_fm_create,&
                                               cp_fm_get_info,&
-                                              cp_fm_indxg2l,&
-                                              cp_fm_indxg2p,&
                                               cp_fm_release,&
                                               cp_fm_set_all,&
                                               cp_fm_to_fm_submat_general,&
@@ -160,9 +158,9 @@ SUBROUTINE fm_general_add_bse(fm_out, fm_in, beta, nrow_secidx_in, ncol_secidx_i
       CHARACTER(LEN=*), PARAMETER :: routineN = 'fm_general_add_bse'
 
       INTEGER :: col_idx_loc, dummy, handle, handle2, i_entry_rec, idx_col_out, idx_row_out, ii, &
-         iproc, jj, ncol_block_in, ncol_block_out, ncol_local_in, ncol_local_out, npcol, nprocs, &
-         nprow, nrow_block_in, nrow_block_out, nrow_local_in, nrow_local_out, proc_send, &
-         row_idx_loc, send_pcol, send_prow
+         iproc, jj, ncol_block_in, ncol_block_out, ncol_local_in, ncol_local_out, nprocs, &
+         nrow_block_in, nrow_block_out, nrow_local_in, nrow_local_out, proc_send, row_idx_loc, &
+         send_pcol, send_prow
       INTEGER, ALLOCATABLE, DIMENSION(:)                 :: entry_counter, num_entries_rec, &
                                                             num_entries_send
       INTEGER, DIMENSION(4)                              :: indices_in
@@ -187,9 +185,6 @@ SUBROUTINE fm_general_add_bse(fm_out, fm_in, beta, nrow_secidx_in, ncol_secidx_i
                           nrow_block=nrow_block_out, &
                           ncol_block=ncol_block_out)
 
-      nprow = fm_out%matrix_struct%context%num_pe(1)
-      npcol = fm_out%matrix_struct%context%num_pe(2)
-
       ALLOCATE (num_entries_rec(0:para_env_out%num_pe - 1))
       ALLOCATE (num_entries_send(0:para_env_out%num_pe - 1))
 
@@ -237,11 +232,8 @@ SUBROUTINE fm_general_add_bse(fm_out, fm_in, beta, nrow_secidx_in, ncol_secidx_i
             idx_row_out = indices_in(reordering(2)) + (indices_in(reordering(1)) - 1)*nrow_secidx_out
             idx_col_out = indices_in(reordering(4)) + (indices_in(reordering(3)) - 1)*ncol_secidx_out
 
-            send_prow = cp_fm_indxg2p(idx_row_out, nrow_block_out, dummy, &
-                                      fm_out%matrix_struct%first_p_pos(1), nprow)
-
-            send_pcol = cp_fm_indxg2p(idx_col_out, ncol_block_out, dummy, &
-                                      fm_out%matrix_struct%first_p_pos(2), npcol)
+            send_prow = fm_out%matrix_struct%g2p_row(idx_row_out)
+            send_pcol = fm_out%matrix_struct%g2p_col(idx_col_out)
 
             proc_send = fm_out%matrix_struct%context%blacs2mpi(send_prow, send_pcol)
 
@@ -320,11 +312,8 @@ SUBROUTINE fm_general_add_bse(fm_out, fm_in, beta, nrow_secidx_in, ncol_secidx_i
             idx_row_out = indices_in(reordering(2)) + (indices_in(reordering(1)) - 1)*nrow_secidx_out
             idx_col_out = indices_in(reordering(4)) + (indices_in(reordering(3)) - 1)*ncol_secidx_out
 
-            send_prow = cp_fm_indxg2p(idx_row_out, nrow_block_out, dummy, &
-                                      fm_out%matrix_struct%first_p_pos(1), nprow)
-
-            send_pcol = cp_fm_indxg2p(idx_col_out, ncol_block_out, dummy, &
-                                      fm_out%matrix_struct%first_p_pos(2), npcol)
+            send_prow = fm_out%matrix_struct%g2p_row(idx_row_out)
+            send_pcol = fm_out%matrix_struct%g2p_col(idx_col_out)
 
             proc_send = fm_out%matrix_struct%context%blacs2mpi(send_prow, send_pcol)
             entry_counter(proc_send) = entry_counter(proc_send) + 1
@@ -362,15 +351,12 @@ SUBROUTINE fm_general_add_bse(fm_out, fm_in, beta, nrow_secidx_in, ncol_secidx_i
       nprocs = para_env_out%num_pe
 
 !$OMP PARALLEL DO DEFAULT(NONE) &
-!$OMP SHARED(fm_out, nprocs, nrow_block_out, ncol_block_out, &
-!$OMP          num_entries_rec, buffer_rec, beta, dummy, nprow, npcol) &
+!$OMP SHARED(fm_out, nprocs, num_entries_rec, buffer_rec, beta) &
 !$OMP PRIVATE(iproc, i_entry_rec, ii, jj)
       DO iproc = 0, nprocs - 1
          DO i_entry_rec = 1, num_entries_rec(iproc)
-            ii = cp_fm_indxg2l(buffer_rec(iproc)%indx(i_entry_rec, 1), nrow_block_out, &
-                               dummy, dummy, nprow)
-            jj = cp_fm_indxg2l(buffer_rec(iproc)%indx(i_entry_rec, 2), ncol_block_out, &
-                               dummy, dummy, npcol)
+            ii = fm_out%matrix_struct%g2l_row(buffer_rec(iproc)%indx(i_entry_rec, 1))
+            jj = fm_out%matrix_struct%g2l_col(buffer_rec(iproc)%indx(i_entry_rec, 2))
 
             fm_out%local_data(ii, jj) = fm_out%local_data(ii, jj) + beta*buffer_rec(iproc)%msg(i_entry_rec)
          END DO
@@ -430,8 +416,8 @@ SUBROUTINE truncate_fm(fm_out, fm_in, ncol_in, &
 
       INTEGER :: col_idx_loc, dummy, handle, handle2, i_entry_rec, idx_col_first, idx_col_in, &
          idx_col_out, idx_col_sec, idx_row_in, ii, iproc, jj, ncol_block_in, ncol_block_out, &
-         ncol_local_in, ncol_local_out, npcol, nprocs, nprow, nrow_block_in, nrow_block_out, &
-         nrow_local_in, nrow_local_out, proc_send, row_idx_loc, send_pcol, send_prow
+         ncol_local_in, ncol_local_out, nprocs, nrow_block_in, nrow_block_out, nrow_local_in, &
+         nrow_local_out, proc_send, row_idx_loc, send_pcol, send_prow
       INTEGER, ALLOCATABLE, DIMENSION(:)                 :: entry_counter, num_entries_rec, &
                                                             num_entries_send
       INTEGER, DIMENSION(:), POINTER                     :: col_indices_in, col_indices_out, &
@@ -465,9 +451,6 @@ SUBROUTINE truncate_fm(fm_out, fm_in, ncol_in, &
                           nrow_block=nrow_block_out, &
                           ncol_block=ncol_block_out)
 
-      nprow = fm_out%matrix_struct%context%num_pe(1)
-      npcol = fm_out%matrix_struct%context%num_pe(2)
-
       ALLOCATE (num_entries_rec(0:para_env_out%num_pe - 1))
       ALLOCATE (num_entries_send(0:para_env_out%num_pe - 1))
 
@@ -529,11 +512,8 @@ SUBROUTINE truncate_fm(fm_out, fm_in, ncol_in, &
          DO row_idx_loc = 1, nrow_local_in
             idx_row_in = row_indices_in(row_idx_loc)
 
-            send_prow = cp_fm_indxg2p(idx_row_in, nrow_block_out, dummy, &
-                                      fm_out%matrix_struct%first_p_pos(1), nprow)
-
-            send_pcol = cp_fm_indxg2p(idx_col_out, ncol_block_out, dummy, &
-                                      fm_out%matrix_struct%first_p_pos(2), npcol)
+            send_prow = fm_out%matrix_struct%g2p_row(idx_row_in)
+            send_pcol = fm_out%matrix_struct%g2p_col(idx_col_out)
 
             proc_send = fm_out%matrix_struct%context%blacs2mpi(send_prow, send_pcol)
 
@@ -628,11 +608,9 @@ SUBROUTINE truncate_fm(fm_out, fm_in, ncol_in, &
          DO row_idx_loc = 1, nrow_local_in
             idx_row_in = row_indices_in(row_idx_loc)
 
-            send_prow = cp_fm_indxg2p(idx_row_in, nrow_block_out, dummy, &
-                                      fm_out%matrix_struct%first_p_pos(1), nprow)
+            send_prow = fm_out%matrix_struct%g2p_row(idx_row_in)
 
-            send_pcol = cp_fm_indxg2p(idx_col_out, ncol_block_out, dummy, &
-                                      fm_out%matrix_struct%first_p_pos(2), npcol)
+            send_pcol = fm_out%matrix_struct%g2p_col(idx_col_out)
 
             proc_send = fm_out%matrix_struct%context%blacs2mpi(send_prow, send_pcol)
             entry_counter(proc_send) = entry_counter(proc_send) + 1
@@ -671,15 +649,12 @@ SUBROUTINE truncate_fm(fm_out, fm_in, ncol_in, &
       nprocs = para_env_out%num_pe
 
 !$OMP PARALLEL DO DEFAULT(NONE) &
-!$OMP SHARED(fm_out, nprocs, nrow_block_out, ncol_block_out, &
-!$OMP          num_entries_rec, buffer_rec, nprow, npcol, dummy) &
+!$OMP SHARED(fm_out, nprocs, num_entries_rec, buffer_rec) &
 !$OMP PRIVATE(iproc, i_entry_rec, ii, jj)
       DO iproc = 0, nprocs - 1
          DO i_entry_rec = 1, num_entries_rec(iproc)
-            ii = cp_fm_indxg2l(buffer_rec(iproc)%indx(i_entry_rec, 1), nrow_block_out, &
-                               dummy, dummy, nprow)
-            jj = cp_fm_indxg2l(buffer_rec(iproc)%indx(i_entry_rec, 2), ncol_block_out, &
-                               dummy, dummy, npcol)
+            ii = fm_out%matrix_struct%g2l_row(buffer_rec(iproc)%indx(i_entry_rec, 1))
+            jj = fm_out%matrix_struct%g2l_col(buffer_rec(iproc)%indx(i_entry_rec, 2))
 
             fm_out%local_data(ii, jj) = fm_out%local_data(ii, jj) + buffer_rec(iproc)%msg(i_entry_rec)
          END DO

diff --git a/src/common/bibliography.F b/src/common/bibliography.F
@@ -91,7 +91,7 @@ MODULE bibliography
                             Rycroft2009, Thomas2015, Brehm2018, Brehm2020, Shigeta2001, Heinecke2016, &
                             Brehm2021, Bussy2021a, Bussy2021b, Ditler2021, Ditler2022, Mattiat2019, &
                             Mattiat2022, Belleflamme2023, Knizia2013, Musaelian2023, Eriksen2020, &
-                            Bussy2023, Wang2018, Zeng2023, Graml2024, Solca2024
+                            Bussy2023, Wang2018, Zeng2023, Graml2024, Solca2024, Freeman1977, Gruneis2009
 
 CONTAINS
 
@@ -5010,6 +5010,46 @@ SUBROUTINE add_all_references()
                          "ER"), &
                          DOI="10.1007/978-3-031-61763-8_13")
 
+      CALL add_reference(key=Freeman1977, ISI_record=s2a( &
+                         "TY JOUR", &
+                         "PB American Physical Society", &
+                         "ID 10.1103/PhysRevB.15.5512", &
+                         "DO 10.1103/PhysRevB.15.5512", &
+                         "TI Coupled-cluster expansion applied to the electron gas: Inclusion of ring and exchange effects", &
+                         "PY 1977", &
+                         "UR https://link.aps.org/doi/10.1103/PhysRevB.15.5512", &
+                         "JF Physical Review B", &
+                         "JA Phys. Rev. B", &
+                         "J1 PRB", &
+                         "VL 15", &
+                         "IS 12", &
+                         "SP 5512", &
+                         "EP 5521", &
+                         "AU Freeman, David L.", &
+                         "ER"), &
+                         DOI="10.1103/PhysRevB.15.5512")
+
+      CALL add_reference(key=Gruneis2009, ISI_record=s2a( &
+                         "TY JOUR", &
+                         "AU Grüneis, Andreas", &
+                         "   Marsman, Martijn", &
+                         "   Harl, Judith", &
+                         "   Schimka, Laurids", &
+                         "   Kresse, Georg", &
+                         "TI Making the random phase approximation to electronic correlation accurate", &
+                         "PY 2009", &
+                         "Y1 2009/10/21", &
+                         "JO The Journal of Chemical Physics", &
+                         "JA J. Chem. Phys.", &
+                         "VL 131", &
+                         "IS 15", &
+                         "SP 154115", &
+                         "SN 0021-9606", &
+                         "Y2 5/1/2024", &
+                         "UR https://doi.org/10.1063/1.3250347", &
+                         "ER"), &
+                         DOI="10.1063/1.3250347")
+
    END SUBROUTINE add_all_references
 
 END MODULE bibliography
diff --git a/src/fm/cp_cfm_basic_linalg.F b/src/fm/cp_cfm_basic_linalg.F
@@ -84,9 +84,8 @@ SUBROUTINE cp_cfm_det(matrix_a, det_a)
       COMPLEX(KIND=dp), DIMENSION(:), POINTER  :: diag
 
 #if defined(__SCALAPACK)
-      INTEGER, EXTERNAL                        :: indxl2g
-      INTEGER                                  :: myprow, nprow, npcol, nrow_local, nrow_block, irow_local, &
-                                                  mypcol, ncol_local, ncol_block, icol_local, j
+      INTEGER                                  :: myprow, nprow, npcol, nrow_local, irow_local, &
+                                                  mypcol, ncol_local, icol_local, j
       INTEGER, DIMENSION(9)                    :: desca
 #endif
 
@@ -111,21 +110,19 @@ SUBROUTINE cp_cfm_det(matrix_a, det_a)
       nprow = matrix_lu%matrix_struct%context%num_pe(1)
       npcol = matrix_lu%matrix_struct%context%num_pe(2)
       nrow_local = matrix_lu%matrix_struct%nrow_locals(myprow)
-      nrow_block = matrix_lu%matrix_struct%nrow_block
-      ncol_block = matrix_lu%matrix_struct%ncol_block
       ncol_local = matrix_lu%matrix_struct%ncol_locals(mypcol)
 
       DO irow_local = 1, nrow_local
-         i = indxl2g(irow_local, nrow_block, myprow, matrix_lu%matrix_struct%first_p_pos(1), nprow)
+         i = matrix_lu%matrix_struct%row_indices(irow_local)
          DO icol_local = 1, ncol_local
-            j = indxl2g(icol_local, ncol_block, mypcol, matrix_lu%matrix_struct%first_p_pos(2), npcol)
+            j = matrix_lu%matrix_struct%col_indices(icol_local)
             IF (i == j) diag(i) = matrix_lu%local_data(irow_local, icol_local)
          END DO
       END DO
       CALL matrix_lu%matrix_struct%para_env%sum(diag)
       determinant = PRODUCT(diag)
       DO irow_local = 1, nrow_local
-         i = indxl2g(irow_local, nrow_block, myprow, matrix_lu%matrix_struct%first_p_pos(1), nprow)
+         i = matrix_lu%matrix_struct%row_indices(irow_local)
          IF (ipivot(irow_local) /= i) P = P + 1
       END DO
       CALL matrix_lu%matrix_struct%para_env%sum(P)

diff --git a/src/fm/cp_fm_basic_linalg.F b/src/fm/cp_fm_basic_linalg.F
@@ -105,7 +105,6 @@ SUBROUTINE cp_fm_det(matrix_a, det_a)
       REAL(KIND=dp), DIMENSION(:), POINTER     :: diag
 
 #if defined(__SCALAPACK)
-      INTEGER, EXTERNAL                        :: indxl2g
       INTEGER                                  :: myprow, nprow, npcol, nrow_local, nrow_block, irow_local
       INTEGER, DIMENSION(9)                    :: desca
 #endif
@@ -134,7 +133,7 @@ SUBROUTINE cp_fm_det(matrix_a, det_a)
       nrow_local = matrix_lu%matrix_struct%nrow_locals(myprow)
       nrow_block = matrix_lu%matrix_struct%nrow_block
       DO irow_local = 1, nrow_local
-         i = indxl2g(irow_local, nrow_block, myprow, matrix_lu%matrix_struct%first_p_pos(1), nprow)
+         i = matrix_lu%matrix_struct%row_indices(irow_local)
          IF (ipivot(irow_local) /= i) P = P + 1
       END DO
       CALL matrix_lu%matrix_struct%para_env%sum(P)
@@ -1305,7 +1304,6 @@ SUBROUTINE cp_fm_upper_to_full(matrix, work)
                                                   ncol_block, ncol_local, &
                                                   nrow_block, nrow_local
       INTEGER, DIMENSION(9)                    :: desca, descc
-      INTEGER, EXTERNAL                        :: indxl2g
       REAL(KIND=dp), DIMENSION(:, :), POINTER   :: c
       REAL(KIND=sp), DIMENSION(:, :), POINTER   :: c_sp
 #endif
@@ -1342,11 +1340,9 @@ SUBROUTINE cp_fm_upper_to_full(matrix, work)
       descc(:) = matrix%matrix_struct%descriptor(:)
 
       DO icol_local = 1, ncol_local
-         icol_global = indxl2g(icol_local, ncol_block, mypcol, &
-                               matrix%matrix_struct%first_p_pos(2), npcol)
+         icol_global = matrix%matrix_struct%col_indices(icol_local)
          DO irow_local = 1, nrow_local
-            irow_global = indxl2g(irow_local, nrow_block, myprow, &
-                                  matrix%matrix_struct%first_p_pos(1), nprow)
+            irow_global = matrix%matrix_struct%row_indices(irow_local)
             IF (irow_global > icol_global) THEN
                IF (matrix%use_sp) THEN
                   c_sp(irow_local, icol_local) = 0.0_sp