(git:08825b9)
Loading...
Searching...
No Matches
spectral_localizer_flatten.F
Go to the documentation of this file.
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! **************************************************************************************************
13 USE cp_dbcsr_api, ONLY: dbcsr_add,&
21 USE ieee_arithmetic, ONLY: ieee_is_finite
24 USE kinds, ONLY: dp
34
35 IMPLICIT NONE
36 PRIVATE
38
39CONTAINS
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 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 COMPLEX(dp), ALLOCATABLE :: dense(:, :), dense_metric(:, :)
69 INTEGER :: invariant, n, ns
70 LOGICAL :: valid
71 TYPE(localizer_matrix_type) :: metric, pencil
72
73 status = -1
74 gap = 0.0_dp
75 metric_ratio = 0.0_dp
76 residuals = huge(1.0_dp)
77 IF (.NOT. all(ieee_is_finite([energy, scale, tolerance, metric_tolerance, gap_tolerance]))) RETURN
78 IF (scale <= 0.0_dp .OR. tolerance <= 0.0_dp .OR. tolerance >= 1.0_dp) RETURN
79 IF (metric_tolerance <= 0.0_dp .OR. metric_tolerance >= 1.0_dp .OR. gap_tolerance < 0.0_dp) RETURN
80 CALL dbcsr_get_info(h%real_part, nfullrows_total=n)
81 CALL dbcsr_get_info(s%real_part, nfullrows_total=ns)
82 IF (n < 1 .OR. n /= ns) RETURN
83 CALL localizer_dbcsr_check(h, para_env, tolerance, valid)
84 IF (.NOT. valid) RETURN
85 CALL localizer_dbcsr_check(s, para_env, tolerance, valid)
86 IF (.NOT. valid) RETURN
87 CALL localizer_dbcsr_gap_pencil(h, s, energy, pencil, metric)
88 CALL localizer_dbcsr_dense(pencil, para_env, dense)
89 CALL localizer_dbcsr_dense(metric, para_env, dense_metric)
90 IF (para_env%is_source()) THEN
91 CALL localizer_matrix_evaluate(dense, dense_metric(:n, :n), .false., tolerance, metric_tolerance, &
92 gap_tolerance, invariant, gap(1), metric_ratio, status)
93 IF (status > 1) status = -4
94 gap(2) = gap(1)
95 END IF
96 CALL localizer_dbcsr_release(pencil)
97 CALL localizer_dbcsr_release(metric)
98 CALL para_env%bcast(status, para_env%source)
99 CALL para_env%bcast(gap, para_env%source)
100 CALL para_env%bcast(metric_ratio, para_env%source)
101 IF (status /= 0) RETURN
102 CALL para_env%bcast(invariant, para_env%source)
103 IF (invariant /= 0) THEN
104 status = -4
105 RETURN
106 END IF
107 CALL localizer_flatten_dbcsr(h, s, energy, scale, tolerance, para_env, flattened, residuals, status)
108 IF (status /= 0) status = -6
109 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 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 status = -1
139 residuals = huge(1.0_dp)
140 IF (.NOT. all(ieee_is_finite([energy, scale, tolerance]))) RETURN
141 IF (scale <= 0.0_dp .OR. tolerance <= 0.0_dp .OR. tolerance >= 1.0_dp) RETURN
142 CALL dbcsr_get_info(h%real_part, nfullrows_total=n)
143 CALL dbcsr_get_info(s%real_part, nfullrows_total=ns)
144 IF (n < 1 .OR. n /= ns) RETURN
145 CALL localizer_dbcsr_check(h, para_env, tolerance, valid)
146 IF (.NOT. valid) RETURN
147 CALL localizer_dbcsr_check(s, para_env, tolerance, valid)
148 IF (.NOT. valid) RETURN
149 CALL localizer_dbcsr_realify(h, a)
150 CALL localizer_dbcsr_realify(s, b)
151 CALL dbcsr_add(a, b, 1.0_dp, -energy)
152 CALL dbcsr_create(inverse, template=b)
153 CALL dbcsr_create(sign_matrix, template=b)
154 CALL dbcsr_create(trial, template=b)
155 CALL dbcsr_create(work, template=b)
156 filter = min(1.e-16_dp, 0.01_dp*tolerance**2)
157 DO
158 norm = dbcsr_frobenius_norm(a)
159 status = 1
160 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 accelerator_order=0, silent=.true.)
163 CALL dbcsr_multiply('N', 'N', 1.0_dp, b, inverse, 0.0_dp, work)
164 CALL dbcsr_add_on_diag(work, -1.0_dp)
165 residuals(1) = dbcsr_frobenius_norm(work)
166 IF (.NOT. ieee_is_finite(residuals(1)) .OR. residuals(1) > tolerance) EXIT
167 CALL dbcsr_multiply('N', 'N', 1.0_dp, inverse, a, 0.0_dp, trial)
168 CALL matrix_sign_newton_schulz(sign_matrix, trial, filter, sign_order=2, iounit=-1)
169 CALL dbcsr_multiply('N', 'N', 1.0_dp, sign_matrix, sign_matrix, 0.0_dp, work)
170 CALL dbcsr_add_on_diag(work, -1.0_dp)
171 residuals(2) = dbcsr_frobenius_norm(work)
172 IF (.NOT. ieee_is_finite(residuals(2)) .OR. residuals(2) > tolerance) EXIT
173 CALL dbcsr_multiply('N', 'N', scale, b, sign_matrix, 0.0_dp, trial)
174 CALL localizer_dbcsr_from_real(trial, h, flattened)
175 CALL dbcsr_release(work)
176 CALL localizer_dbcsr_realify(flattened, work)
177 norm = max(dbcsr_frobenius_norm(trial), tiny(1.0_dp))
178 CALL dbcsr_add(work, trial, 1.0_dp, -1.0_dp)
179 residuals(3) = dbcsr_frobenius_norm(work)/norm
180 CALL localizer_dbcsr_check(flattened, para_env, tolerance, valid, residuals(4))
181 status = 2
182 IF (.NOT. valid .OR. .NOT. all(ieee_is_finite(residuals)) .OR. maxval(residuals) > tolerance) THEN
183 CALL localizer_dbcsr_release(flattened)
184 EXIT
185 END IF
186 status = 0
187 EXIT
188 END DO
189 CALL dbcsr_release(a)
190 CALL dbcsr_release(b)
191 CALL dbcsr_release(inverse)
192 CALL dbcsr_release(sign_matrix)
193 CALL dbcsr_release(trial)
194 CALL dbcsr_release(work)
195 END SUBROUTINE localizer_flatten_dbcsr
196
subroutine, public dbcsr_multiply(transa, transb, alpha, matrix_a, matrix_b, beta, matrix_c, first_row, last_row, first_column, last_column, first_k, last_k, retain_sparsity, filter_eps, flop)
...
subroutine, public dbcsr_get_info(matrix, nblkrows_total, nblkcols_total, nfullrows_total, nfullcols_total, nblkrows_local, nblkcols_local, nfullrows_local, nfullcols_local, my_prow, my_pcol, local_rows, local_cols, proc_row_dist, proc_col_dist, row_blk_size, col_blk_size, row_blk_offset, col_blk_offset, distribution, name, matrix_type, group)
...
subroutine, public dbcsr_release(matrix)
...
subroutine, public dbcsr_add(matrix_a, matrix_b, alpha_scalar, beta_scalar)
...
subroutine, public dbcsr_add_on_diag(matrix, alpha)
Adds the given scalar to the diagonal of the matrix. Reserves any missing diagonal blocks.
real(dp) function, public dbcsr_frobenius_norm(matrix)
Compute the frobenius norm of a dbcsr matrix.
Routines useful for iterative matrix calculations.
subroutine, public matrix_sign_newton_schulz(matrix_sign, matrix, threshold, sign_order, iounit)
compute the sign a matrix using Newton-Schulz iterations
subroutine, public invert_hotelling(matrix_inverse, matrix, threshold, use_inv_as_guess, norm_convergence, filter_eps, accelerator_order, max_iter_lanczos, eps_lanczos, silent)
invert a symmetric positive definite matrix by Hotelling's method explicit symmetrization makes this ...
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Interface to the message passing library MPI.
Distributed block-sparse finite localizers, retaining the covariant AO metric.
subroutine, public localizer_dbcsr_from_real(realified, template, matrix)
Extract the two left quadrants of R(A); does not average away structure errors.
subroutine, public localizer_dbcsr_dense(matrix, para_env, dense)
Collect the size-limited dense oracle; sparse solvers must not call this routine.
subroutine, public localizer_dbcsr_gap_pencil(h, s, energy, pencil, metric)
Pair H-E*S with its negative to reuse the metric-aware localizer gap solvers.
subroutine, public localizer_dbcsr_release(matrix)
Release both real components of a complex sparse operator.
subroutine, public localizer_dbcsr_check(matrix, para_env, tolerance, valid, relative_error)
Check both components before a sparse solver discards the upper triangle.
subroutine, public localizer_dbcsr_realify(matrix, realified)
Distributed real matrix R(A)=[Re(A),-Im(A);Im(A),Re(A)].
Metric-covariant spectral flattening with existing distributed matrix iterations.
subroutine, public localizer_flatten_checked(h, s, energy, scale, tolerance, metric_tolerance, gap_tolerance, para_env, flattened, gap, metric_ratio, residuals, status)
Establish a positive AO metric and electronic gap before iterating the matrix sign.
subroutine, public localizer_flatten_dbcsr(h, s, energy, scale, tolerance, para_env, flattened, residuals, status)
Form scale*S*sign(S^-1*(H-energy*S)) without dense AO matrices or S^-1/2.
Dense reference spectral localizers in a nonorthogonal finite basis.
subroutine, public localizer_matrix_evaluate(localizer, s, aii, matrix_tol, metric_tol, gap_tol, invariant, gap, metric_ratio, status)
Dense reference evaluation of an assembled localizer with metric diag(S,S).
stores all the informations relevant to an mpi environment