LCOV - code coverage report
Current view: top level - src - spectral_localizer_flatten.F (source / functions) Coverage Total Hit
Test: CP2K Regtests (git:9ccafda) Lines: 95.7 % 92 88
Test Date: 2026-10-15 00:55:44 Functions: 100.0 % 2 2

            Line data    Source code
       1              : !--------------------------------------------------------------------------------------------------!
       2              : !   CP2K: A general program to perform molecular dynamics simulations                              !
       3              : !   Copyright 2000-2026 CP2K developers group <https://cp2k.org>                                   !
       4              : !                                                                                                  !
       5              : !   SPDX-License-Identifier: GPL-2.0-or-later                                                      !
       6              : !--------------------------------------------------------------------------------------------------!
       7              : 
       8              : ! **************************************************************************************************
       9              : !> \brief Metric-covariant spectral flattening with existing distributed matrix iterations.
      10              : !> \author Thomas D. Kuehne, tkuehne@cp2k.org
      11              : ! **************************************************************************************************
      12              : MODULE spectral_localizer_flatten
      13              :    USE cp_dbcsr_api,                    ONLY: dbcsr_add,&
      14              :                                               dbcsr_create,&
      15              :                                               dbcsr_get_info,&
      16              :                                               dbcsr_multiply,&
      17              :                                               dbcsr_release,&
      18              :                                               dbcsr_type
      19              :    USE cp_dbcsr_contrib,                ONLY: dbcsr_add_on_diag,&
      20              :                                               dbcsr_frobenius_norm
      21              :    USE ieee_arithmetic,                 ONLY: ieee_is_finite
      22              :    USE iterate_matrix,                  ONLY: invert_Hotelling,&
      23              :                                               matrix_sign_Newton_Schulz
      24              :    USE kinds,                           ONLY: dp
      25              :    USE message_passing,                 ONLY: mp_para_env_type
      26              :    USE spectral_localizer,              ONLY: localizer_matrix_evaluate
      27              :    USE spectral_localizer_dbcsr,        ONLY: localizer_dbcsr_check,&
      28              :                                               localizer_dbcsr_dense,&
      29              :                                               localizer_dbcsr_from_real,&
      30              :                                               localizer_dbcsr_gap_pencil,&
      31              :                                               localizer_dbcsr_realify,&
      32              :                                               localizer_dbcsr_release,&
      33              :                                               localizer_matrix_type
      34              : 
      35              :    IMPLICIT NONE
      36              :    PRIVATE
      37              :    PUBLIC :: localizer_flatten_dbcsr, localizer_flatten_checked
      38              : 
      39              : CONTAINS
      40              : 
      41              : ! **************************************************************************************************
      42              : !> \brief Establish a positive AO metric and electronic gap before iterating the matrix sign.
      43              : !> \param h Covariant physical Hamiltonian, unchanged
      44              : !> \param s Covariant physical metric, unchanged
      45              : !> \param energy Reference energy
      46              : !> \param scale Positive flattened-band energy
      47              : !> \param tolerance Matrix-sign residual and Hermiticity tolerance
      48              : !> \param metric_tolerance Minimum overlap eigenvalue ratio
      49              : !> \param gap_tolerance Absolute minimum original electronic gap
      50              : !> \param para_env Communicator
      51              : !> \param flattened Newly initialized shifted Hamiltonian only on success
      52              : !> \param gap Original electronic gap lower and upper bounds, not the flattened localizer gap
      53              : !> \param metric_ratio Overlap eigenvalue ratio
      54              : !> \param residuals Four matrix-sign residuals
      55              : !> \param status 0 success, 1 unresolved input gap, -1 input, -2 metric, -4 solver,
      56              : !>        -6 matrix-sign iteration or structure failure
      57              : ! **************************************************************************************************
      58           38 :    SUBROUTINE localizer_flatten_checked(h, s, energy, scale, tolerance, metric_tolerance, gap_tolerance, &
      59              :                                         para_env, flattened, gap, metric_ratio, residuals, status)
      60              :       TYPE(localizer_matrix_type), INTENT(INOUT)         :: h, s
      61              :       REAL(dp), INTENT(IN)                               :: energy, scale, tolerance, &
      62              :                                                             metric_tolerance, gap_tolerance
      63              :       TYPE(mp_para_env_type), POINTER                    :: para_env
      64              :       TYPE(localizer_matrix_type), INTENT(OUT)           :: flattened
      65              :       REAL(dp), INTENT(OUT)                              :: gap(2), metric_ratio, residuals(4)
      66              :       INTEGER, INTENT(OUT)                               :: status
      67              : 
      68           38 :       COMPLEX(dp), ALLOCATABLE                           :: dense(:, :), dense_metric(:, :)
      69              :       INTEGER                                            :: invariant, n, ns
      70              :       LOGICAL                                            :: valid
      71              :       TYPE(localizer_matrix_type)                        :: metric, pencil
      72              : 
      73           38 :       status = -1
      74           38 :       gap = 0.0_dp
      75           38 :       metric_ratio = 0.0_dp
      76          190 :       residuals = HUGE(1.0_dp)
      77          228 :       IF (.NOT. ALL(ieee_is_finite([energy, scale, tolerance, metric_tolerance, gap_tolerance]))) RETURN
      78           38 :       IF (scale <= 0.0_dp .OR. tolerance <= 0.0_dp .OR. tolerance >= 1.0_dp) RETURN
      79           38 :       IF (metric_tolerance <= 0.0_dp .OR. metric_tolerance >= 1.0_dp .OR. gap_tolerance < 0.0_dp) RETURN
      80           38 :       CALL dbcsr_get_info(h%real_part, nfullrows_total=n)
      81           38 :       CALL dbcsr_get_info(s%real_part, nfullrows_total=ns)
      82           38 :       IF (n < 1 .OR. n /= ns) RETURN
      83           38 :       CALL localizer_dbcsr_check(h, para_env, tolerance, valid)
      84           38 :       IF (.NOT. valid) RETURN
      85           38 :       CALL localizer_dbcsr_check(s, para_env, tolerance, valid)
      86           38 :       IF (.NOT. valid) RETURN
      87           38 :       CALL localizer_dbcsr_gap_pencil(h, s, energy, pencil, metric)
      88           38 :       CALL localizer_dbcsr_dense(pencil, para_env, dense)
      89           38 :       CALL localizer_dbcsr_dense(metric, para_env, dense_metric)
      90           38 :       IF (para_env%is_source()) THEN
      91              :          CALL localizer_matrix_evaluate(dense, dense_metric(:n, :n), .FALSE., tolerance, metric_tolerance, &
      92           19 :                                         gap_tolerance, invariant, gap(1), metric_ratio, status)
      93           19 :          IF (status > 1) status = -4
      94           19 :          gap(2) = gap(1)
      95              :       END IF
      96           38 :       CALL localizer_dbcsr_release(pencil)
      97           38 :       CALL localizer_dbcsr_release(metric)
      98           38 :       CALL para_env%bcast(status, para_env%source)
      99           38 :       CALL para_env%bcast(gap, para_env%source)
     100           38 :       CALL para_env%bcast(metric_ratio, para_env%source)
     101           38 :       IF (status /= 0) RETURN
     102           30 :       CALL para_env%bcast(invariant, para_env%source)
     103           30 :       IF (invariant /= 0) THEN
     104            0 :          status = -4
     105            0 :          RETURN
     106              :       END IF
     107           30 :       CALL localizer_flatten_dbcsr(h, s, energy, scale, tolerance, para_env, flattened, residuals, status)
     108           30 :       IF (status /= 0) status = -6
     109           76 :    END SUBROUTINE localizer_flatten_checked
     110              : 
     111              : ! **************************************************************************************************
     112              : !> \brief Form scale*S*sign(S^-1*(H-energy*S)) without dense AO matrices or S^-1/2.
     113              : !> \param h Covariant Hermitian Hamiltonian; unchanged
     114              : !> \param s Positive definite AO metric; must be validated by the caller
     115              : !> \param energy Reference energy, separated from the generalized spectrum at input accuracy
     116              : !> \param scale Positive energy scale of the flattened bands
     117              : !> \param tolerance Residual tolerance, not a certificate of the original Hamiltonian gap
     118              : !> \param para_env DBCSR communicator
     119              : !> \param flattened Newly initialized covariant operator only on success
     120              : !> \param residuals Inverse, sign-square, complex-structure and Hermiticity residuals
     121              : !> \param status 0 success, -1 invalid input, 1 unconverged sign/inverse, 2 broken structure
     122              : !> \note Caller must establish a positive metric and a resolved electronic gap independently.
     123              : !>       Matrix-sign convergence alone cannot certify the occupied subspace near a gap closing.
     124              : ! **************************************************************************************************
     125          334 :    SUBROUTINE localizer_flatten_dbcsr(h, s, energy, scale, tolerance, para_env, flattened, residuals, status)
     126              :       TYPE(localizer_matrix_type), INTENT(INOUT)         :: h, s
     127              :       REAL(dp), INTENT(IN)                               :: energy, scale, tolerance
     128              :       TYPE(mp_para_env_type), POINTER                    :: para_env
     129              :       TYPE(localizer_matrix_type), INTENT(OUT)           :: flattened
     130              :       REAL(dp), INTENT(OUT)                              :: residuals(4)
     131              :       INTEGER, INTENT(OUT)                               :: status
     132              : 
     133              :       INTEGER                                            :: n, ns
     134              :       LOGICAL                                            :: valid
     135              :       REAL(dp)                                           :: filter, norm
     136              :       TYPE(dbcsr_type)                                   :: a, b, inverse, sign_matrix, trial, work
     137              : 
     138           52 :       status = -1
     139          260 :       residuals = HUGE(1.0_dp)
     140          202 :       IF (.NOT. ALL(ieee_is_finite([energy, scale, tolerance]))) RETURN
     141           50 :       IF (scale <= 0.0_dp .OR. tolerance <= 0.0_dp .OR. tolerance >= 1.0_dp) RETURN
     142           48 :       CALL dbcsr_get_info(h%real_part, nfullrows_total=n)
     143           48 :       CALL dbcsr_get_info(s%real_part, nfullrows_total=ns)
     144           48 :       IF (n < 1 .OR. n /= ns) RETURN
     145           48 :       CALL localizer_dbcsr_check(h, para_env, tolerance, valid)
     146           48 :       IF (.NOT. valid) RETURN
     147           46 :       CALL localizer_dbcsr_check(s, para_env, tolerance, valid)
     148           46 :       IF (.NOT. valid) RETURN
     149           46 :       CALL localizer_dbcsr_realify(h, a)
     150           46 :       CALL localizer_dbcsr_realify(s, b)
     151           46 :       CALL dbcsr_add(a, b, 1.0_dp, -energy)
     152           46 :       CALL dbcsr_create(inverse, template=b)
     153           46 :       CALL dbcsr_create(sign_matrix, template=b)
     154           46 :       CALL dbcsr_create(trial, template=b)
     155           46 :       CALL dbcsr_create(work, template=b)
     156           46 :       filter = MIN(1.e-16_dp, 0.01_dp*tolerance**2)
     157          126 :       DO
     158           46 :          norm = dbcsr_frobenius_norm(a)
     159           46 :          status = 1
     160           46 :          IF (.NOT. ieee_is_finite(norm) .OR. norm == 0.0_dp) EXIT
     161              :          CALL invert_Hotelling(inverse, b, 0.1_dp*tolerance, filter_eps=filter, &
     162           44 :                                accelerator_order=0, silent=.TRUE.)
     163           44 :          CALL dbcsr_multiply('N', 'N', 1.0_dp, b, inverse, 0.0_dp, work)
     164           44 :          CALL dbcsr_add_on_diag(work, -1.0_dp)
     165           44 :          residuals(1) = dbcsr_frobenius_norm(work)
     166           44 :          IF (.NOT. ieee_is_finite(residuals(1)) .OR. residuals(1) > tolerance) EXIT
     167           44 :          CALL dbcsr_multiply('N', 'N', 1.0_dp, inverse, a, 0.0_dp, trial)
     168           44 :          CALL matrix_sign_Newton_Schulz(sign_matrix, trial, filter, sign_order=2, iounit=-1)
     169           44 :          CALL dbcsr_multiply('N', 'N', 1.0_dp, sign_matrix, sign_matrix, 0.0_dp, work)
     170           44 :          CALL dbcsr_add_on_diag(work, -1.0_dp)
     171           44 :          residuals(2) = dbcsr_frobenius_norm(work)
     172           44 :          IF (.NOT. ieee_is_finite(residuals(2)) .OR. residuals(2) > tolerance) EXIT
     173           42 :          CALL dbcsr_multiply('N', 'N', scale, b, sign_matrix, 0.0_dp, trial)
     174           42 :          CALL localizer_dbcsr_from_real(trial, h, flattened)
     175           42 :          CALL dbcsr_release(work)
     176           42 :          CALL localizer_dbcsr_realify(flattened, work)
     177           42 :          norm = MAX(dbcsr_frobenius_norm(trial), TINY(1.0_dp))
     178           42 :          CALL dbcsr_add(work, trial, 1.0_dp, -1.0_dp)
     179           42 :          residuals(3) = dbcsr_frobenius_norm(work)/norm
     180           42 :          CALL localizer_dbcsr_check(flattened, para_env, tolerance, valid, residuals(4))
     181           42 :          status = 2
     182          420 :          IF (.NOT. valid .OR. .NOT. ALL(ieee_is_finite(residuals)) .OR. MAXVAL(residuals) > tolerance) THEN
     183            0 :             CALL localizer_dbcsr_release(flattened)
     184            0 :             EXIT
     185              :          END IF
     186           42 :          status = 0
     187           46 :          EXIT
     188              :       END DO
     189           46 :       CALL dbcsr_release(a)
     190           46 :       CALL dbcsr_release(b)
     191           46 :       CALL dbcsr_release(inverse)
     192           46 :       CALL dbcsr_release(sign_matrix)
     193           46 :       CALL dbcsr_release(trial)
     194           46 :       CALL dbcsr_release(work)
     195           52 :    END SUBROUTINE localizer_flatten_dbcsr
     196              : 
     197              : END MODULE spectral_localizer_flatten
        

Generated by: LCOV version 2.0-1