Line data    Source code

       1             : !--------------------------------------------------------------------------------------------------!
       2             : !   CP2K: A general program to perform molecular dynamics simulations                              !
       3             : !   Copyright 2000-2024 CP2K developers group <https://cp2k.org>                                   !
       4             : !                                                                                                  !
       5             : !   SPDX-License-Identifier: GPL-2.0-or-later                                                      !
       6             : !--------------------------------------------------------------------------------------------------!
       7             : 
       8             : ! **************************************************************************************************
       9             : !> \brief   Interface to (sca)lapack for the Cholesky based procedures
      10             : !> \author  VW
      11             : !> \date    2009-11-09
      12             : !> \version 0.8
      13             : !>
      14             : !> <b>Modification history:</b>
      15             : !> - Created 2009-11-09
      16             : ! **************************************************************************************************
      17             : MODULE cp_dbcsr_diag
      18             : 
      19             :    USE cp_blacs_env,                    ONLY: cp_blacs_env_type
      20             :    USE cp_cfm_diag,                     ONLY: cp_cfm_heevd
      21             :    USE cp_cfm_types,                    ONLY: cp_cfm_create,&
      22             :                                               cp_cfm_release,&
      23             :                                               cp_cfm_type
      24             :    USE cp_dbcsr_operations,             ONLY: copy_cfm_to_dbcsr,&
      25             :                                               copy_dbcsr_to_cfm,&
      26             :                                               copy_dbcsr_to_fm,&
      27             :                                               copy_fm_to_dbcsr
      28             :    USE cp_fm_diag,                      ONLY: choose_eigv_solver,&
      29             :                                               cp_fm_power,&
      30             :                                               cp_fm_syevx
      31             :    USE cp_fm_struct,                    ONLY: cp_fm_struct_create,&
      32             :                                               cp_fm_struct_release,&
      33             :                                               cp_fm_struct_type
      34             :    USE cp_fm_types,                     ONLY: cp_fm_create,&
      35             :                                               cp_fm_release,&
      36             :                                               cp_fm_type
      37             :    USE dbcsr_api,                       ONLY: dbcsr_get_info,&
      38             :                                               dbcsr_type
      39             :    USE kinds,                           ONLY: dp
      40             :    USE message_passing,                 ONLY: mp_para_env_type
      41             : #include "base/base_uses.f90"
      42             : 
      43             :    IMPLICIT NONE
      44             : 
      45             :    PRIVATE
      46             : 
      47             :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_dbcsr_diag'
      48             : 
      49             :    ! Public subroutines
      50             : 
      51             :    PUBLIC :: cp_dbcsr_syevd, &
      52             :              cp_dbcsr_syevx, &
      53             :              cp_dbcsr_heevd, &
      54             :              cp_dbcsr_power
      55             : 
      56             : CONTAINS
      57             : 
      58             : ! **************************************************************************************************
      59             : !> \brief ...
      60             : !> \param matrix ...
      61             : !> \param eigenvectors ...
      62             : !> \param eigenvalues ...
      63             : !> \param para_env ...
      64             : !> \param blacs_env ...
      65             : ! **************************************************************************************************
      66       45199 :    SUBROUTINE cp_dbcsr_syevd(matrix, eigenvectors, eigenvalues, para_env, blacs_env)
      67             : 
      68             :       ! Computes all eigenvalues and vectors of a real symmetric matrix
      69             :       ! should be quite a bit faster than syevx for that case
      70             :       ! especially in parallel with thightly clustered evals
      71             :       ! needs more workspace in the worst case, but much better distributed
      72             : 
      73             :       TYPE(dbcsr_type)                                   :: matrix, eigenvectors
      74             :       REAL(KIND=dp), DIMENSION(:), INTENT(OUT)           :: eigenvalues
      75             :       TYPE(mp_para_env_type), POINTER                    :: para_env
      76             :       TYPE(cp_blacs_env_type), POINTER                   :: blacs_env
      77             : 
      78             :       CHARACTER(len=*), PARAMETER                        :: routineN = 'cp_dbcsr_syevd'
      79             : 
      80             :       INTEGER                                            :: handle, nfullrows_total
      81             :       TYPE(cp_fm_struct_type), POINTER                   :: fm_struct
      82             :       TYPE(cp_fm_type)                                   :: fm_eigenvectors, fm_matrix
      83             : 
      84       45199 :       CALL timeset(routineN, handle)
      85             : 
      86       45199 :       NULLIFY (fm_struct)
      87       45199 :       CALL dbcsr_get_info(matrix, nfullrows_total=nfullrows_total)
      88             : 
      89             :       CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
      90       45199 :                                ncol_global=nfullrows_total, para_env=para_env)
      91       45199 :       CALL cp_fm_create(fm_matrix, fm_struct, name="fm_matrix")
      92       45199 :       CALL cp_fm_create(fm_eigenvectors, fm_struct, name="fm_eigenvectors")
      93       45199 :       CALL cp_fm_struct_release(fm_struct)
      94             : 
      95       45199 :       CALL copy_dbcsr_to_fm(matrix, fm_matrix)
      96             : 
      97       45199 :       CALL choose_eigv_solver(fm_matrix, fm_eigenvectors, eigenvalues)
      98             : 
      99       45199 :       CALL copy_fm_to_dbcsr(fm_eigenvectors, eigenvectors)
     100             : 
     101       45199 :       CALL cp_fm_release(fm_matrix)
     102       45199 :       CALL cp_fm_release(fm_eigenvectors)
     103             : 
     104       45199 :       CALL timestop(handle)
     105             : 
     106       45199 :    END SUBROUTINE cp_dbcsr_syevd
     107             : 
     108             : ! **************************************************************************************************
     109             : !> \brief   compute eigenvalues and optionally eigenvectors of a real symmetric matrix using scalapack.
     110             : !>          If eigenvectors are required this routine will replicate a full matrix on each CPU...
     111             : !>          if more than a handful of vectors are needed, use cp_dbcsr_syevd instead
     112             : !> \param matrix ...
     113             : !> \param eigenvectors ...
     114             : !> \param eigenvalues ...
     115             : !> \param neig ...
     116             : !> \param work_syevx ...
     117             : !> \param para_env ...
     118             : !> \param blacs_env ...
     119             : !> \par     matrix is supposed to be in upper triangular form, and overwritten by this routine
     120             : !>          neig   is the number of vectors needed (default all)
     121             : !>          work_syevx evec calculation only, is the fraction of the working buffer allowed (1.0 use full buffer)
     122             : !>                     reducing this saves time, but might cause the routine to fail
     123             : ! **************************************************************************************************
     124           0 :    SUBROUTINE cp_dbcsr_syevx(matrix, eigenvectors, eigenvalues, neig, work_syevx, &
     125             :                              para_env, blacs_env)
     126             : 
     127             :       ! Diagonalise the symmetric n by n matrix using the LAPACK library.
     128             : 
     129             :       TYPE(dbcsr_type), POINTER                          :: matrix
     130             :       TYPE(dbcsr_type), OPTIONAL, POINTER                :: eigenvectors
     131             :       REAL(KIND=dp), DIMENSION(:), INTENT(OUT)           :: eigenvalues
     132             :       INTEGER, INTENT(IN), OPTIONAL                      :: neig
     133             :       REAL(KIND=dp), INTENT(IN), OPTIONAL                :: work_syevx
     134             :       TYPE(mp_para_env_type), POINTER                    :: para_env
     135             :       TYPE(cp_blacs_env_type), POINTER                   :: blacs_env
     136             : 
     137             :       CHARACTER(LEN=*), PARAMETER                        :: routineN = 'cp_dbcsr_syevx'
     138             : 
     139             :       INTEGER                                            :: handle, n, neig_local
     140             :       TYPE(cp_fm_struct_type), POINTER                   :: fm_struct
     141             :       TYPE(cp_fm_type)                                   :: fm_eigenvectors, fm_matrix
     142             : 
     143           0 :       CALL timeset(routineN, handle)
     144             : 
     145             :       ! by default all
     146           0 :       CALL dbcsr_get_info(matrix, nfullrows_total=n)
     147           0 :       neig_local = n
     148           0 :       IF (PRESENT(neig)) neig_local = neig
     149           0 :       IF (neig_local == 0) RETURN
     150             : 
     151           0 :       NULLIFY (fm_struct)
     152             :       CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=n, &
     153           0 :                                ncol_global=n, para_env=para_env)
     154           0 :       CALL cp_fm_create(fm_matrix, fm_struct, name="fm_matrix")
     155             : 
     156           0 :       CALL copy_dbcsr_to_fm(matrix, fm_matrix)
     157             : 
     158           0 :       IF (PRESENT(eigenvectors)) THEN
     159           0 :          CALL cp_fm_create(fm_eigenvectors, fm_struct, name="fm_eigenvectors")
     160           0 :          CALL cp_fm_syevx(fm_matrix, fm_eigenvectors, eigenvalues, neig, work_syevx)
     161           0 :          CALL copy_fm_to_dbcsr(fm_eigenvectors, eigenvectors)
     162           0 :          CALL cp_fm_release(fm_eigenvectors)
     163             :       ELSE
     164           0 :          CALL cp_fm_syevx(fm_matrix, eigenvalues=eigenvalues, neig=neig, work_syevx=work_syevx)
     165             :       END IF
     166             : 
     167           0 :       CALL cp_fm_struct_release(fm_struct)
     168           0 :       CALL cp_fm_release(fm_matrix)
     169             : 
     170           0 :       CALL timestop(handle)
     171             : 
     172           0 :    END SUBROUTINE cp_dbcsr_syevx
     173             : 
     174             : ! **************************************************************************************************
     175             : !> \brief ...
     176             : !> \param matrix ...
     177             : !> \param eigenvectors ...
     178             : !> \param eigenvalues ...
     179             : !> \param para_env ...
     180             : !> \param blacs_env ...
     181             : ! **************************************************************************************************
     182        2570 :    SUBROUTINE cp_dbcsr_heevd(matrix, eigenvectors, eigenvalues, para_env, blacs_env)
     183             : 
     184             :       TYPE(dbcsr_type)                                   :: matrix
     185             :       TYPE(dbcsr_type), OPTIONAL, POINTER                :: eigenvectors
     186             :       REAL(KIND=dp), DIMENSION(:), INTENT(OUT)           :: eigenvalues
     187             :       TYPE(mp_para_env_type), POINTER                    :: para_env
     188             :       TYPE(cp_blacs_env_type), POINTER                   :: blacs_env
     189             : 
     190             :       CHARACTER(len=*), PARAMETER                        :: routineN = 'cp_dbcsr_heevd'
     191             : 
     192             :       INTEGER                                            :: handle, nfullrows_total
     193             :       TYPE(cp_cfm_type)                                  :: fm_eigenvectors, fm_matrix
     194             :       TYPE(cp_fm_struct_type), POINTER                   :: fm_struct
     195             : 
     196        2570 :       CALL timeset(routineN, handle)
     197             : 
     198        2570 :       NULLIFY (fm_struct)
     199        2570 :       CALL dbcsr_get_info(matrix, nfullrows_total=nfullrows_total)
     200             : 
     201             :       CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
     202        2570 :                                ncol_global=nfullrows_total, para_env=para_env)
     203        2570 :       CALL cp_cfm_create(fm_matrix, fm_struct, name="fm_matrix")
     204        2570 :       CALL cp_cfm_create(fm_eigenvectors, fm_struct, name="fm_eigenvectors")
     205        2570 :       CALL cp_fm_struct_release(fm_struct)
     206             : 
     207        2570 :       CALL copy_dbcsr_to_cfm(matrix, fm_matrix)
     208             : 
     209        2570 :       CALL cp_cfm_heevd(fm_matrix, fm_eigenvectors, eigenvalues)
     210             : 
     211        2570 :       CALL copy_cfm_to_dbcsr(fm_eigenvectors, eigenvectors)
     212             : 
     213        2570 :       CALL cp_cfm_release(fm_matrix)
     214        2570 :       CALL cp_cfm_release(fm_eigenvectors)
     215             : 
     216        2570 :       CALL timestop(handle)
     217             : 
     218        2570 :    END SUBROUTINE cp_dbcsr_heevd
     219             : 
     220             : ! **************************************************************************************************
     221             : !> \brief ...
     222             : !> \param matrix ...
     223             : !> \param exponent ...
     224             : !> \param threshold ...
     225             : !> \param n_dependent ...
     226             : !> \param para_env ...
     227             : !> \param blacs_env ...
     228             : !> \param verbose ...
     229             : !> \param eigenvectors ...
     230             : !> \param eigenvalues ...
     231             : ! **************************************************************************************************
     232         164 :    SUBROUTINE cp_dbcsr_power(matrix, exponent, threshold, n_dependent, para_env, blacs_env, verbose, eigenvectors, eigenvalues)
     233             :       TYPE(dbcsr_type), INTENT(INOUT)                    :: matrix
     234             :       REAL(dp), INTENT(IN)                               :: exponent, threshold
     235             :       INTEGER, INTENT(OUT)                               :: n_dependent
     236             :       TYPE(mp_para_env_type), POINTER                    :: para_env
     237             :       TYPE(cp_blacs_env_type), POINTER                   :: blacs_env
     238             :       LOGICAL, INTENT(IN), OPTIONAL                      :: verbose
     239             :       TYPE(dbcsr_type), INTENT(INOUT), OPTIONAL          :: eigenvectors
     240             :       REAL(KIND=dp), DIMENSION(2), INTENT(OUT), OPTIONAL :: eigenvalues
     241             : 
     242             :       CHARACTER(len=*), PARAMETER                        :: routineN = 'cp_dbcsr_power'
     243             : 
     244             :       INTEGER                                            :: handle, nfullrows_total
     245             :       REAL(KIND=dp), DIMENSION(2)                        :: eigenvalues_prv
     246             :       TYPE(cp_fm_struct_type), POINTER                   :: fm_struct
     247             :       TYPE(cp_fm_type)                                   :: fm_eigenvectors, fm_matrix
     248             : 
     249          82 :       CALL timeset(routineN, handle)
     250             : 
     251          82 :       NULLIFY (fm_struct)
     252          82 :       CALL dbcsr_get_info(matrix, nfullrows_total=nfullrows_total)
     253             : 
     254             :       CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
     255          82 :                                ncol_global=nfullrows_total, para_env=para_env)
     256          82 :       CALL cp_fm_create(fm_matrix, fm_struct, name="fm_matrix")
     257          82 :       CALL cp_fm_create(fm_eigenvectors, fm_struct, name="fm_eigenvectors")
     258          82 :       CALL cp_fm_struct_release(fm_struct)
     259             : 
     260          82 :       CALL copy_dbcsr_to_fm(matrix, fm_matrix)
     261             : 
     262          82 :       CALL cp_fm_power(fm_matrix, fm_eigenvectors, exponent, threshold, n_dependent, verbose, eigenvalues_prv)
     263             : 
     264          82 :       CALL copy_fm_to_dbcsr(fm_matrix, matrix)
     265          82 :       CALL cp_fm_release(fm_matrix)
     266             : 
     267          82 :       IF (PRESENT(eigenvalues)) eigenvalues(:) = eigenvalues_prv
     268          82 :       IF (PRESENT(eigenvectors)) CALL copy_fm_to_dbcsr(fm_eigenvectors, eigenvectors)
     269             : 
     270          82 :       CALL cp_fm_release(fm_eigenvectors)
     271             : 
     272          82 :       CALL timestop(handle)
     273             : 
     274          82 :    END SUBROUTINE
     275             : 
     276             : END MODULE cp_dbcsr_diag