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
|