Merge branch '7-update-of-omm_callback-in-libomm' into 'master'

Resolve "Update of omm_callback in libOMM"

Closes #7

See merge request ElectronicStructureLibrary/omm-bundle!15

MatrixSwitch/examples/example_pdcsr_pddbc.F90

@@ -127,7 +127,7 @@ program example_pdcsr_pddbc
   if (mpi_rank==0) then
      do ii=1,SIZE(row_blk_sizes)
         do jj=1,SIZE(col_blk_sizes)
-           call m_get_element(C,ii,jj,myblock,found_block)
+           call m_get_element(C,ii,jj,myblock,found=found_block)
            if (found_block) then
               print *, "block found at (",ii,",",jj,") on rank 0. Values= ", myblock
               do ll=1, block_size

MatrixSwitch/src/MatrixSwitch_ops.F90

@@ -33,7 +33,7 @@ module MatrixSwitch_ops
 #ifdef HAVE_MPI
   character(1), save :: ms_lap_order !< Ordering of the BLACS process grid.
 
-  logical, save :: ms_mpi_use2D !< Is 2D MPI grid used by default?
+  logical, save :: ms_mpi_use2D !< Is 2D process grid used by default?
   integer, save :: ms_mpi_comm !< MPI communicator.
   integer, save :: ms_mpi_size !< Number of MPI processes.
   integer, save :: ms_mpi_rank !< MPI process rank.
@@ -41,14 +41,14 @@ module MatrixSwitch_ops
   integer, save :: ms_lap_npcol !< Number of columns in the BLACS process grid.
   integer, save :: ms_lap_bs_def !< Default block size.
   integer, save :: ms_lap_bs_num !< Number of block size exceptions.
-  integer, save :: ms_lap_icontxt !< BLACS context handle used by MatrixSwitch for 2D case.
+  integer, save :: ms_lap_icontxt !< BLACS context handle used by MatrixSwitch.
   integer, save :: ms_lap_icontxt_1D !< BLACS context handle used by MatrixSwitch for 1D case.
   integer, allocatable, save :: ms_lap_bs_list(:,:) !< Block size exception list.
 
 #ifdef HAVE_DBCSR
   logical, save :: ms_dbcsr_init = .false. !< check if dbcsr is set.
-  integer, save :: ms_dbcsr_group = mpi_comm_null !< Group communication for 2D case
-  integer, save :: ms_dbcsr_brd_group = mpi_comm_null !< Group communication for 1D case
+  integer, save :: ms_dbcsr_group = mpi_comm_null !< Group communication.
+  integer, save :: ms_dbcsr_group_1D = mpi_comm_null !< Group communication for 1D case.
 #endif
 
 #endif
@@ -75,14 +75,14 @@ module MatrixSwitch_ops
      logical :: iaux3_is_allocated=.false. !< Is iaux3 directly allocated or it is a pointer?
      logical :: iaux4_is_allocated=.false. !< Is iaux4 directly allocated or it is a pointer?
      logical :: dval_is_allocated=.false. !< Is dval directly allocated or it is a pointer?
-     logical :: csr_dval_is_allocated=.false. !< Is dval directly allocated or it is a pointer?
+     logical :: csr_dval_is_allocated=.false. !< Is csr_dval directly allocated or it is a pointer?
      logical :: zval_is_allocated=.false. !< Is zval directly allocated or it is a pointer?
-     logical :: use2D !< Is the 2D MPI grid used?
+     logical :: use2D !< Is the 2D process grid used?
 
      integer :: dim1 !< Row dimension size of the matrix.
      integer :: dim2 !< Column dimension size of the matrix.
-     integer :: csr_nrows=0 !< Local number of csr rows.
-     integer :: csr_nze=0 !< Local number of csr nonempty elements.
+     integer :: csr_nrows=0 !< Local number of rows of the csr matrix.
+     integer :: csr_nze=0 !< Local number of nonempty elements of the csr matrix.
      integer :: blk_size1=0 !< Block size for rows
      integer :: blk_size2=0 !< Block size for columns
      integer, pointer :: iaux1(:) => null() !< Auxiliary information for certain storage formats.

MatrixSwitch/src/MatrixSwitch_wrapper.F90

@@ -939,7 +939,7 @@ contains
   !============================================================================!
   !> @brief Wrapper to get matrix element (real version).
   !============================================================================!
-  subroutine m_dget_element(C,i,j,alpha,found,label)
+  subroutine m_dget_element(C,i,j,alpha,label,found)
     implicit none
 
     !**** INPUT ***********************************!
@@ -958,14 +958,14 @@ contains
 
     !**********************************************!
 
-    call m_get_element_orig(ms_matrices(ms_lookup(C)),i,j,alpha,found,label)
+    call m_get_element_orig(ms_matrices(ms_lookup(C)),i,j,alpha,label,found)
 
   end subroutine m_dget_element
 
   !============================================================================!
   !> @brief Wrapper to get matrix element (complex version).
   !============================================================================!
-  subroutine m_zget_element(C,i,j,alpha,found,label)
+  subroutine m_zget_element(C,i,j,alpha,label,found)
     implicit none
 
     !**** INPUT ***********************************!
@@ -984,14 +984,14 @@ contains
 
     !**********************************************!
 
-    call m_get_element_orig(ms_matrices(ms_lookup(C)),i,j,alpha,found,label)
+    call m_get_element_orig(ms_matrices(ms_lookup(C)),i,j,alpha,label,found)
 
   end subroutine m_zget_element
 
   !============================================================================!
   !> @brief Wrapper to get block (real version).
   !============================================================================!
-  subroutine m_dget_block(C,i,j,alpha,found,label)
+  subroutine m_dget_block(C,i,j,alpha,label,found)
     implicit none
 
     !**** INPUT ***********************************!
@@ -1010,7 +1010,7 @@ contains
 
     !**********************************************!
 
-    call m_get_element_orig(ms_matrices(ms_lookup(C)),i,j,alpha,found,label)
+    call m_get_element_orig(ms_matrices(ms_lookup(C)),i,j,alpha,label,found)
 
   end subroutine m_dget_block
 

libOMM/src/Makefile.manual

@@ -24,7 +24,7 @@ clean :
 	rm -f *.mod
 
 .F90.o : 
-	$(FORTRAN) $(OPTS) $(FPPFLAGS) -c $< -I$(MSWINCLUDEPATH) -I$(PSPINCLUDEPATH) -o $@
+	$(FORTRAN) $(OPTS) $(FPPFLAGS) -c $< -I$(MSWINCLUDEPATH) -I$(PSPINCLUDEPATH) $(DBCSRINC) -o $@
 
 .f90.o : 
-	$(FORTRAN) $(OPTS) -c $< -I$(MSWINCLUDEPATH) -I$(PSPINCLUDEPATH) -o $@
+	$(FORTRAN) $(OPTS) -c $< -I$(MSWINCLUDEPATH) -I$(PSPINCLUDEPATH) $(DBCSRINC) -o $@

libOMM/src/omm_ops.F90

@@ -77,6 +77,7 @@ subroutine calc_HW_callback(H,C,HW,HC,m_operation)
 
   interface
     subroutine H(C,HC)
+      ! The Hamiltonian should be already shifted by eta (if != 0) 
       use MatrixSwitch
       implicit none
       type(matrix), intent(in) :: C ! WF coeffs. matrix
@@ -162,10 +163,10 @@ subroutine calc_coeff(HW,SW,HWd,SWd,HWdd,SWdd,SWdH,coeff,m_operation)
 
   type(matrix), intent(inout) :: HW ! hamiltonian matrix in WF basis
   type(matrix), intent(inout) :: SW ! overlap matrix in WF basis
-  type(matrix), intent(inout) :: HWd ! g^T*h*c
-  type(matrix), intent(inout) :: SWd ! g^T*s*c
-  type(matrix), intent(inout) :: HWdd ! g^T*h*g
-  type(matrix), intent(inout) :: SWdd ! g^T*h*g
+  type(matrix), intent(inout) :: HWd ! g^dag*h*c
+  type(matrix), intent(inout) :: SWd ! g^dag*s*c
+  type(matrix), intent(inout) :: HWdd ! g^dag*h*g
+  type(matrix), intent(inout) :: SWdd ! g^dag*h*g
   type(matrix), intent(inout) :: SWdH ! work matrix
 
   !**** LOCAL ***********************************!
@@ -215,7 +216,7 @@ end subroutine calc_coeff
 
 !================================================!
 ! matrix factorization                           !
-! C = U^T*U                                      !
+! C = U^dag*U                                      !
 !================================================!
 subroutine m_factorize(C,label)
   implicit none

libOMM/src/omm_params.F90

@@ -10,15 +10,15 @@ public
 
 !**** PARAMS ************************************!
 
-integer, parameter :: dp=selected_real_kind(15,300)
+integer, parameter :: dp = selected_real_kind(15,300)
 
-real(dp), parameter :: Pi=3.141592653589793238462643383279502884197_dp
+real(dp), parameter :: Pi = 3.141592653589793238462643383279502884197_dp
 
-complex(dp), parameter :: cmplx_1=(1.0_dp,0.0_dp)
-complex(dp), parameter :: cmplx_i=(0.0_dp,1.0_dp)
-complex(dp), parameter :: cmplx_0=(0.0_dp,0.0_dp)
+complex(dp), parameter :: cmplx_1 = (1.0_dp, 0.0_dp)
+complex(dp), parameter :: cmplx_i = (0.0_dp, 1.0_dp)
+complex(dp), parameter :: cmplx_0 = (0.0_dp, 0.0_dp)
 
-logical, save :: ms_scalapack_running=.false.
+logical, save :: ms_scalapack_running = .false.
 
 integer, save :: log_unit
 

libOMM/src/omm_rand.F90

@@ -28,11 +28,11 @@ contains
         character(10) :: system_time
         real(dp) :: rtime
 #ifdef NORAND
-        seed=123456
+        seed = 123456
 #else
         call date_and_time(time=system_time)
         read (system_time,*) rtime
-        seed = int(rtime*100.0_dp)
+        seed = int(rtime * 100.0_dp)
 #endif
 #ifdef HAVE_MPI
         seed = seed * (ms_mpi_rank + 1)
@@ -55,8 +55,8 @@ contains
         integer, parameter :: c = 12345
         integer(i64), parameter :: m = 2_i64**31
 
-        x = int(mod(a*x+c,m))
-        r = real(x,dp)/m
+        x = int(mod(a * x + c, m))
+        r = real(x, dp) / m
 
     end subroutine omm_bsd_lcg
 

libOMM/src/omm_wrapper.F90

@@ -2,8 +2,8 @@
 #include "config.h"
 #endif
 
-subroutine omm_wrapper(m,n,n_occ,H,S,new_S,e_min,D_min,calc_ED,eta,C_min,init_C,T,scale_T,flavour,np,ip,cg_tol,&
-                       long_out,dealloc,m_storage,m_operation)
+subroutine omm_wrapper(m,n,n_occ,H,S,new_S,e_min,D_min,calc_ED,eta,C_min,init_C,T,scale_T,&
+                       flavour,np,ip,cg_tol,long_out,dealloc,m_storage,m_operation)
 
   use omm_params, only : dp
   use MatrixSwitch_wrapper
@@ -17,23 +17,29 @@ subroutine omm_wrapper(m,n,n_occ,H,S,new_S,e_min,D_min,calc_ED,eta,C_min,init_C,
   !**** INPUT ***********************************!
 
   character(5), intent(in) :: m_storage ! label identifying the MatrixSwitch storage format
-  character(3), intent(in), optional :: m_operation ! label identifying the MatrixSwitch implementation of the operations to use
+  character(3), intent(in), optional :: m_operation ! label identifying the MatrixSwitch 
+                                                    ! implementation of the operations to use
   character(*), intent(in) :: H ! Hamiltonian matrix
-  character(*), intent(in) :: S ! overlap matrix (or its Cholesky-factorized upper triangular matrix)
-  character(*), intent(in) :: D_min ! density (or energy-weighted density) matrix
-  character(*), intent(in) :: C_min ! WF coeffs.
+  character(*), intent(in) :: S ! overlap matrix (or its Cholesky-factorized upper 
+                                ! triangular matrix)
+  character(*), intent(in) :: D_min ! density (or energy-density) matrix
+  character(*), intent(in) :: C_min ! Hermitian conjugate of the WF coeffs., C^dag
   character(*), intent(in) :: T ! kinetic energy matrix
 
 #ifdef CBIND
   logical(c_bool), intent(in) :: new_S ! is the S matrix new for this value of ip?
-  logical(c_bool), intent(in) :: calc_ED ! calculate the energy-weighted density matrix from the existing WF coeffs.?
-  logical(c_bool), intent(in) :: init_C ! have the WF coeffs. been initialized or modified externally to libOMM?
+  logical(c_bool), intent(in) :: calc_ED ! calculate the energy-density matrix from 
+                                         ! the existing WF coeffs.?
+  logical(c_bool), intent(in) :: init_C ! have the WF coeffs. been initialized or 
+                                        ! modified externally to libOMM?
   logical(c_bool), intent(in) :: long_out ! print detailed output?
   logical(c_bool), intent(in) :: dealloc ! deallocate all internal matrices?
 #else
   logical, intent(in) :: new_S ! is the S matrix new for this value of ip?
-  logical, intent(in) :: calc_ED ! calculate the energy-weighted density matrix from the existing WF coeffs.?
-  logical, intent(in) :: init_C ! have the WF coeffs. been initialized or modified externally to libOMM?
+  logical, intent(in) :: calc_ED ! calculate the energy-density matrix 
+                                 ! from the existing WF coeffs.?
+  logical, intent(in) :: init_C ! have the WF coeffs. been initialized or 
+                                ! modified externally to libOMM?
   logical, intent(in) :: long_out ! print detailed output?
   logical, intent(in) :: dealloc ! deallocate all internal matrices?
 #endif
@@ -50,7 +56,8 @@ subroutine omm_wrapper(m,n,n_occ,H,S,new_S,e_min,D_min,calc_ED,eta,C_min,init_C,
   integer, intent(in) :: n_occ ! number of occupied states
 
   real(dp), intent(in) :: eta ! eigenspectrum shift parameter
-  real(dp), intent(in) :: cg_tol ! convergence tolerance of CG minimization (if negative, default of 1.0d-9 is used)
+  real(dp), intent(in) :: cg_tol ! convergence tolerance of CG minimization 
+                                 ! (if negative, default of 1.0d-9 is used)
   real(dp), intent(in) :: scale_T ! kinetic energy scale for the preconditioning
 
   !**** OUTPUT **********************************!
@@ -67,11 +74,11 @@ subroutine omm_wrapper(m,n,n_occ,H,S,new_S,e_min,D_min,calc_ED,eta,C_min,init_C,
 
   !**********************************************!
 
-  new_S_conv=new_S
-  calc_ED_conv=calc_ED
-  init_C_conv=init_C
-  long_out_conv=long_out
-  dealloc_conv=dealloc
+  new_S_conv = new_S
+  calc_ED_conv = calc_ED
+  init_C_conv = init_C
+  long_out_conv = long_out
+  dealloc_conv = dealloc
 
   call omm(m, &
            n, &

tomato/src/Makefile.manual

@@ -24,7 +24,7 @@ clean :
 	rm -f *.mod
 
 .F90.o : 
-	$(FORTRAN) $(OPTS) $(FPPFLAGS) -c $< -I$(MSWINCLUDEPATH) -I$(PSPINCLUDEPATH) -o $@
+	$(FORTRAN) $(OPTS) $(FPPFLAGS) -c $< -I$(MSWINCLUDEPATH) -I$(PSPINCLUDEPATH) $(DBCSRINC) -o $@
 
 .f90.o : 
-	$(FORTRAN) $(OPTS) -c $< -I$(MSWINCLUDEPATH) -I$(PSPINCLUDEPATH) -o $@
+	$(FORTRAN) $(OPTS) -c $< -I$(MSWINCLUDEPATH) -I$(PSPINCLUDEPATH) $(DBCSRINC) -o $@