(git:5e7fe52)
Loading...
Searching...
No Matches
qs_scf_subspace.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 History management and matrix combination for ADIIS.
10! **************************************************************************************************
11MODULE qs_scf_subspace
12
13 USE cp_dbcsr_api, ONLY: dbcsr_add,&
14 dbcsr_copy,&
15 dbcsr_create,&
16 dbcsr_p_type,&
17 dbcsr_set
18 USE cp_dbcsr_contrib, ONLY: dbcsr_dot
19 USE ieee_arithmetic, ONLY: ieee_is_finite
20 USE kinds, ONLY: dp
21 USE qs_scf_subspace_math, ONLY: qs_scf_subspace_build_adiis_model,&
22 qs_scf_subspace_fifo_slot,&
23 simplex_qp_nonfinite_input,&
24 simplex_qp_pairwise_stationary,&
25 simplex_qp_success,&
26 simplex_quadratic_minimize
27 USE qs_scf_subspace_types, ONLY: qs_scf_subspace_buffer_clear,&
28 qs_scf_subspace_buffer_type
29#include "./base/base_uses.f90"
30
31 IMPLICIT NONE
32
33 PRIVATE
34
35 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_scf_subspace'
36
37 INTEGER, PARAMETER, PUBLIC :: scf_subspace_invalid_history = 100
38
39 PUBLIC :: qs_scf_subspace_build, qs_scf_subspace_push, qs_scf_subspace_restart
40
41CONTAINS
42
43! **************************************************************************************************
44!> \brief Restart the ADIIS history from the current strictly paired P, F[P] state.
45!> \param buffer Persistent SCF subspace history.
46!> \param fock Current raw F[P], indexed by spin and real-space cell.
47!> \param density Current input P that generated fock, with matching indices.
48!> \param state_energy ...
49!> \param restarted Whether the current state was stored successfully.
50! **************************************************************************************************
51 SUBROUTINE qs_scf_subspace_restart(buffer, fock, density, state_energy, restarted)
52
53 TYPE(qs_scf_subspace_buffer_type), INTENT(INOUT) :: buffer
54 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: fock, density
55 REAL(kind=dp), INTENT(IN) :: state_energy
56 LOGICAL, INTENT(OUT) :: restarted
57
58 restarted = .false.
59 CALL qs_scf_subspace_buffer_clear(buffer)
60 CALL qs_scf_subspace_push(buffer, fock, density, state_energy, restarted)
61 IF (.NOT. restarted) RETURN
62
63 buffer%coefficients = 0.0_dp
64 buffer%coefficients(1) = 1.0_dp
65 buffer%last_status = simplex_qp_success
66 buffer%last_restart = .true.
67 buffer%use_combined_fock = .false.
68 buffer%last_objective = 0.0_dp
69 buffer%last_old_fock_weight = 0.0_dp
70
71 END SUBROUTINE qs_scf_subspace_restart
72
73! **************************************************************************************************
74!> \brief Append one accepted, evaluated P,F[P] state to the ADIIS history.
75!>
76!> The SCF driver owns state acceptance. Trial endpoints and arbitrary initial guesses
77!> must not call this routine.
78!> \param buffer Persistent SCF subspace history.
79!> \param fock Current raw F[P], indexed by spin and real-space cell.
80!> \param density Current input P that generated fock, with matching indices.
81!> \param state_energy ...
82!> \param pushed Whether the state was stored successfully.
83! **************************************************************************************************
84 SUBROUTINE qs_scf_subspace_push(buffer, fock, density, state_energy, pushed)
85
86 TYPE(qs_scf_subspace_buffer_type), INTENT(INOUT) :: buffer
87 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: fock, density
88 REAL(kind=dp), INTENT(IN) :: state_energy
89 LOGICAL, INTENT(OUT) :: pushed
90
91 LOGICAL :: storage_valid
92
93 pushed = .false.
94 buffer%use_combined_fock = .false.
95 buffer%last_restart = .false.
96 buffer%last_objective = 0.0_dp
97 buffer%last_old_fock_weight = 0.0_dp
98 IF (ALLOCATED(buffer%coefficients)) buffer%coefficients = 0.0_dp
99
100 IF (buffer%nbuffer < 1) THEN
101 buffer%last_status = scf_subspace_invalid_history
102 RETURN
103 END IF
104 IF (buffer%ncall < 0 .OR. buffer%nstored < 0 .OR. buffer%nstored > buffer%nbuffer .OR. &
105 buffer%nstored /= min(buffer%ncall, buffer%nbuffer)) THEN
106 buffer%last_status = scf_subspace_invalid_history
107 RETURN
108 END IF
109
110 CALL ensure_matrix_storage(buffer, fock, density, storage_valid)
111 IF (.NOT. storage_valid) THEN
112 buffer%last_status = scf_subspace_invalid_history
113 RETURN
114 END IF
115
116 CALL push_paired_state(buffer, fock, density, state_energy, pushed)
117 IF (.NOT. pushed) THEN
118 buffer%last_status = simplex_qp_nonfinite_input
119 RETURN
120 END IF
121 buffer%last_status = simplex_qp_success
122
123 END SUBROUTINE qs_scf_subspace_push
124
125! **************************************************************************************************
126!> \brief Solve the ADIIS model from previously accepted history and form an effective KS matrix.
127!> \param buffer Persistent SCF subspace history.
128! **************************************************************************************************
129 SUBROUTINE qs_scf_subspace_build(buffer)
130
131 TYPE(qs_scf_subspace_buffer_type), INTENT(INOUT) :: buffer
132
133 INTEGER :: i, j, m, physical
134 INTEGER, ALLOCATABLE, DIMENSION(:) :: slots
135 LOGICAL :: model_valid, solver_usable
136 REAL(kind=dp) :: coefficient_sum, model_delta, &
137 model_tolerance
138 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: difference, linear
139 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: chronological_pf, hessian
140
141 buffer%use_combined_fock = .false.
142 buffer%last_restart = .false.
143 buffer%last_objective = 0.0_dp
144 buffer%last_old_fock_weight = 0.0_dp
145 IF (ALLOCATED(buffer%coefficients)) buffer%coefficients = 0.0_dp
146
147 IF (buffer%nbuffer < 1 .OR. buffer%nstored < 1) THEN
148 buffer%last_status = scf_subspace_invalid_history
149 RETURN
150 END IF
151 IF (buffer%ncall < 0 .OR. buffer%nstored < 0 .OR. buffer%nstored > buffer%nbuffer .OR. &
152 buffer%nstored /= min(buffer%ncall, buffer%nbuffer)) THEN
153 buffer%last_status = scf_subspace_invalid_history
154 RETURN
155 END IF
156
157 IF (.NOT. ASSOCIATED(buffer%fock) .OR. .NOT. ASSOCIATED(buffer%density) .OR. &
158 .NOT. ASSOCIATED(buffer%combined_fock)) THEN
159 buffer%last_status = scf_subspace_invalid_history
160 RETURN
161 END IF
162 IF (.NOT. ALLOCATED(buffer%generation) .OR. .NOT. ALLOCATED(buffer%pf_metric) .OR. &
163 .NOT. ALLOCATED(buffer%coefficients) .OR. .NOT. ALLOCATED(buffer%state_energy)) THEN
164 buffer%last_status = scf_subspace_invalid_history
165 RETURN
166 END IF
167 IF (SIZE(buffer%generation) /= buffer%nbuffer .OR. &
168 SIZE(buffer%pf_metric, 1) /= buffer%nbuffer .OR. SIZE(buffer%pf_metric, 2) /= buffer%nbuffer .OR. &
169 SIZE(buffer%coefficients) /= buffer%nbuffer .OR. SIZE(buffer%state_energy) /= buffer%nbuffer .OR. &
170 count(buffer%generation > 0) /= buffer%nstored .OR. any(buffer%generation < 0)) THEN
171 buffer%last_status = scf_subspace_invalid_history
172 RETURN
173 END IF
174
175 m = buffer%nstored
176 ALLOCATE (slots(m), chronological_pf(m, m), hessian(m, m), linear(m))
177 CALL chronological_slots(buffer, slots)
178
179 DO i = 1, m
180 DO j = 1, m
181 chronological_pf(i, j) = buffer%pf_metric(slots(i), slots(j))
182 END DO
183 END DO
184
185 CALL qs_scf_subspace_build_adiis_model(chronological_pf, m, hessian, linear, model_valid)
186
187 buffer%coefficients = 0.0_dp
188 IF (.NOT. model_valid) THEN
189 CALL restart_newest_after_failure(buffer, scf_subspace_invalid_history)
190 RETURN
191 END IF
192
193 CALL simplex_quadratic_minimize(hessian, linear, buffer%coefficients(1:m), &
194 buffer%last_objective, buffer%last_status, preferred_index=m)
195
196 solver_usable = buffer%last_status == simplex_qp_success .OR. &
197 buffer%last_status == simplex_qp_pairwise_stationary
198 IF (.NOT. solver_usable) THEN
199 CALL restart_newest_after_failure(buffer, buffer%last_status)
200 RETURN
201 END IF
202 IF (.NOT. ieee_is_finite(buffer%last_objective) .OR. &
203 .NOT. all(ieee_is_finite(buffer%coefficients(1:m)))) THEN
204 CALL restart_newest_after_failure(buffer, simplex_qp_nonfinite_input)
205 RETURN
206 END IF
207
208 coefficient_sum = sum(buffer%coefficients(1:m))
209 IF (minval(buffer%coefficients(1:m)) < -100.0_dp*epsilon(1.0_dp) .OR. &
210 abs(coefficient_sum - 1.0_dp) > 1000.0_dp*epsilon(1.0_dp)) THEN
211 CALL restart_newest_after_failure(buffer, scf_subspace_invalid_history)
212 RETURN
213 END IF
214
215 ALLOCATE (difference(m))
216 difference(:) = buffer%coefficients(1:m)
217 difference(m) = difference(m) - 1.0_dp
218 model_delta = 0.5_dp*dot_product(difference, matmul(hessian, difference)) + &
219 dot_product(difference, hessian(:, m) + linear)
220 model_tolerance = 1000.0_dp*epsilon(1.0_dp)*real(m, kind=dp)* &
221 max(1.0_dp, maxval(abs(hessian)), maxval(abs(linear)))
222 IF (.NOT. ieee_is_finite(model_delta) .OR. model_delta > model_tolerance) THEN
223 IF (ieee_is_finite(model_delta)) THEN
224 CALL restart_newest_after_failure(buffer, scf_subspace_invalid_history)
225 ELSE
226 CALL restart_newest_after_failure(buffer, simplex_qp_nonfinite_input)
227 END IF
228 RETURN
229 END IF
230
231 buffer%last_old_fock_weight = max(0.0_dp, 1.0_dp - buffer%coefficients(m))
232 IF (buffer%last_old_fock_weight <= 100.0_dp*epsilon(1.0_dp)) RETURN
233
234 CALL zero_fock_combination(buffer)
235 DO j = 1, m
236 physical = slots(j)
237 IF (buffer%coefficients(j) <= 0.0_dp) cycle
238 CALL add_fock_to_combination(buffer, physical, buffer%coefficients(j))
239 END DO
240 buffer%use_combined_fock = .true.
241
242 END SUBROUTINE qs_scf_subspace_build
243
244! **************************************************************************************************
245!> \brief Retain only the newest accepted state after an unusable ADIIS model.
246!> \param buffer Persistent SCF subspace history.
247!> \param failure_status Status code that caused the restart.
248! **************************************************************************************************
249 SUBROUTINE restart_newest_after_failure(buffer, failure_status)
250
251 TYPE(qs_scf_subspace_buffer_type), INTENT(INOUT) :: buffer
252 INTEGER, INTENT(IN) :: failure_status
253
254 INTEGER :: icell, ispin, newest
255 REAL(kind=dp) :: contribution, newest_energy, newest_pf
256
257 IF (buffer%nstored < 1) THEN
258 buffer%last_status = failure_status
259 RETURN
260 END IF
261 IF (count(buffer%generation > 0) /= buffer%nstored) THEN
262 buffer%last_status = failure_status
263 RETURN
264 END IF
265 newest = maxloc(buffer%generation, dim=1)
266 newest_energy = buffer%state_energy(newest)
267 newest_pf = 0.0_dp
268 DO icell = 1, SIZE(buffer%fock, 3)
269 DO ispin = 1, SIZE(buffer%fock, 2)
270 CALL dbcsr_dot(buffer%density(newest, ispin, icell)%matrix, &
271 buffer%fock(newest, ispin, icell)%matrix, contribution)
272 newest_pf = newest_pf + contribution
273 END DO
274 END DO
275 IF (.NOT. ieee_is_finite(newest_energy) .OR. .NOT. ieee_is_finite(newest_pf)) THEN
276 buffer%last_status = failure_status
277 RETURN
278 END IF
279
280 buffer%ncall = 1
281 buffer%nstored = 1
282 buffer%generation = 0
283 buffer%generation(newest) = 1
284 buffer%state_energy = huge(1.0_dp)
285 buffer%state_energy(newest) = newest_energy
286 buffer%pf_metric = 0.0_dp
287 buffer%pf_metric(newest, newest) = newest_pf
288 buffer%coefficients = 0.0_dp
289 buffer%coefficients(1) = 1.0_dp
290 buffer%last_status = failure_status
291 buffer%last_restart = .true.
292 buffer%use_combined_fock = .false.
293 buffer%diis_state_valid = .false.
294 buffer%last_objective = 0.0_dp
295 buffer%last_old_fock_weight = 0.0_dp
296 buffer%diis_weight = 0.0_dp
297
298 END SUBROUTINE restart_newest_after_failure
299
300! **************************************************************************************************
301!> \brief Allocate matrix and scalar storage lazily from current matrix templates.
302!> \param buffer Persistent SCF subspace history.
303!> \param fock Current Fock-matrix templates.
304!> \param density Current density-matrix templates.
305!> \param valid Whether dimensions and existing allocation match.
306! **************************************************************************************************
307 SUBROUTINE ensure_matrix_storage(buffer, fock, density, valid)
308
309 TYPE(qs_scf_subspace_buffer_type), INTENT(INOUT) :: buffer
310 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: fock, density
311 LOGICAL, INTENT(OUT) :: valid
312
313 INTEGER :: icell, islot, ispin, ncell, nspin
314
315 valid = ASSOCIATED(fock) .AND. ASSOCIATED(density)
316 IF (.NOT. valid) RETURN
317 nspin = SIZE(fock, 1)
318 ncell = SIZE(fock, 2)
319 valid = nspin > 0 .AND. ncell > 0 .AND. &
320 SIZE(density, 1) == nspin .AND. SIZE(density, 2) == ncell
321 IF (.NOT. valid) RETURN
322 DO icell = 1, ncell
323 DO ispin = 1, nspin
324 valid = ASSOCIATED(fock(ispin, icell)%matrix) .AND. &
325 ASSOCIATED(density(ispin, icell)%matrix)
326 IF (.NOT. valid) RETURN
327 END DO
328 END DO
329
330 IF (ASSOCIATED(buffer%fock)) THEN
331 valid = SIZE(buffer%fock, 1) == buffer%nbuffer .AND. &
332 SIZE(buffer%fock, 2) == nspin .AND. SIZE(buffer%fock, 3) == ncell
333 IF (.NOT. valid) RETURN
334 valid = ASSOCIATED(buffer%density)
335 IF (.NOT. valid) RETURN
336 valid = SIZE(buffer%density, 1) == buffer%nbuffer .AND. &
337 SIZE(buffer%density, 2) == nspin .AND. SIZE(buffer%density, 3) == ncell
338 IF (.NOT. valid) RETURN
339 valid = ASSOCIATED(buffer%combined_fock)
340 IF (.NOT. valid) RETURN
341 valid = SIZE(buffer%combined_fock, 1) == nspin .AND. &
342 SIZE(buffer%combined_fock, 2) == ncell
343 IF (.NOT. valid) RETURN
344 valid = ALLOCATED(buffer%generation) .AND. ALLOCATED(buffer%pf_metric) .AND. &
345 ALLOCATED(buffer%coefficients) .AND. ALLOCATED(buffer%state_energy)
346 IF (.NOT. valid) RETURN
347 valid = SIZE(buffer%generation) == buffer%nbuffer .AND. &
348 SIZE(buffer%pf_metric, 1) == buffer%nbuffer .AND. &
349 SIZE(buffer%pf_metric, 2) == buffer%nbuffer .AND. &
350 SIZE(buffer%coefficients) == buffer%nbuffer .AND. &
351 SIZE(buffer%state_energy) == buffer%nbuffer
352 IF (.NOT. valid) RETURN
353 DO icell = 1, ncell
354 DO ispin = 1, nspin
355 valid = ASSOCIATED(buffer%combined_fock(ispin, icell)%matrix)
356 IF (.NOT. valid) RETURN
357 DO islot = 1, buffer%nbuffer
358 valid = ASSOCIATED(buffer%fock(islot, ispin, icell)%matrix) .AND. &
359 ASSOCIATED(buffer%density(islot, ispin, icell)%matrix)
360 IF (.NOT. valid) RETURN
361 END DO
362 END DO
363 END DO
364 RETURN
365 END IF
366
367 valid = .NOT. ASSOCIATED(buffer%density) .AND. &
368 .NOT. ASSOCIATED(buffer%combined_fock) .AND. .NOT. ALLOCATED(buffer%generation) .AND. &
369 .NOT. ALLOCATED(buffer%pf_metric) .AND. .NOT. ALLOCATED(buffer%coefficients) .AND. &
370 .NOT. ALLOCATED(buffer%state_energy)
371 IF (.NOT. valid) RETURN
372
373 ALLOCATE (buffer%fock(buffer%nbuffer, nspin, ncell))
374 ALLOCATE (buffer%density(buffer%nbuffer, nspin, ncell))
375 ALLOCATE (buffer%combined_fock(nspin, ncell))
376 ALLOCATE (buffer%generation(buffer%nbuffer), buffer%state_energy(buffer%nbuffer))
377 ALLOCATE (buffer%pf_metric(buffer%nbuffer, buffer%nbuffer), buffer%coefficients(buffer%nbuffer))
378 buffer%generation = 0
379 buffer%state_energy = huge(1.0_dp)
380 buffer%pf_metric = 0.0_dp
381 buffer%coefficients = 0.0_dp
382
383 DO icell = 1, ncell
384 DO ispin = 1, nspin
385 DO islot = 1, buffer%nbuffer
386 ALLOCATE (buffer%fock(islot, ispin, icell)%matrix)
387 CALL dbcsr_create(buffer%fock(islot, ispin, icell)%matrix, &
388 template=fock(ispin, icell)%matrix)
389 ALLOCATE (buffer%density(islot, ispin, icell)%matrix)
390 CALL dbcsr_create(buffer%density(islot, ispin, icell)%matrix, &
391 template=density(ispin, icell)%matrix)
392 END DO
393 ALLOCATE (buffer%combined_fock(ispin, icell)%matrix)
394 CALL dbcsr_create(buffer%combined_fock(ispin, icell)%matrix, &
395 template=fock(ispin, icell)%matrix)
396 END DO
397 END DO
398
399 END SUBROUTINE ensure_matrix_storage
400
401! **************************************************************************************************
402!> \brief Atomically append one strictly paired P, F[P] state.
403!> \param buffer Persistent SCF subspace history.
404!> \param fock Current raw F[P].
405!> \param density Current P.
406!> \param state_energy ...
407!> \param pushed Whether all scalar checks passed and the state was committed.
408! **************************************************************************************************
409 SUBROUTINE push_paired_state(buffer, fock, density, state_energy, pushed)
410
411 TYPE(qs_scf_subspace_buffer_type), INTENT(INOUT) :: buffer
412 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: fock, density
413 REAL(kind=dp), INTENT(IN) :: state_energy
414 LOGICAL, INTENT(OUT) :: pushed
415
416 INTEGER :: icell, islot, ispin, target_slot
417 LOGICAL :: retained
418 REAL(kind=dp) :: current_current
419 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: current_old, old_current
420
421 pushed = .false.
422 IF (.NOT. ieee_is_finite(state_energy)) RETURN
423
424 ! Keep the baseline deterministic and auditable: fill empty slots in physical
425 ! order, then evict the oldest accepted state. Energy remains metadata.
426 target_slot = qs_scf_subspace_fifo_slot(buffer%generation)
427 IF (target_slot <= 0) RETURN
428
429 ALLOCATE (current_old(buffer%nbuffer), old_current(buffer%nbuffer))
430 current_old = 0.0_dp
431 old_current = 0.0_dp
432
433 CALL paired_matrix_dot(density, fock, current_current)
434 IF (.NOT. ieee_is_finite(current_current)) RETURN
435
436 DO islot = 1, buffer%nbuffer
437 retained = buffer%generation(islot) > 0 .AND. islot /= target_slot
438 IF (.NOT. retained) cycle
439 CALL history_current_dot(buffer, islot, fock, density, old_current(islot), current_old(islot))
440 IF (.NOT. ieee_is_finite(old_current(islot)) .OR. &
441 .NOT. ieee_is_finite(current_old(islot))) RETURN
442 END DO
443
444 DO icell = 1, SIZE(fock, 2)
445 DO ispin = 1, SIZE(fock, 1)
446 CALL dbcsr_copy(buffer%fock(target_slot, ispin, icell)%matrix, fock(ispin, icell)%matrix)
447 CALL dbcsr_copy(buffer%density(target_slot, ispin, icell)%matrix, density(ispin, icell)%matrix)
448 END DO
449 END DO
450
451 buffer%pf_metric(target_slot, :) = 0.0_dp
452 buffer%pf_metric(:, target_slot) = 0.0_dp
453 buffer%pf_metric(target_slot, target_slot) = current_current
454 DO islot = 1, buffer%nbuffer
455 retained = buffer%generation(islot) > 0 .AND. islot /= target_slot
456 IF (.NOT. retained) cycle
457 buffer%pf_metric(islot, target_slot) = old_current(islot)
458 buffer%pf_metric(target_slot, islot) = current_old(islot)
459 END DO
460
461 buffer%ncall = buffer%ncall + 1
462 buffer%generation(target_slot) = buffer%ncall
463 buffer%state_energy(target_slot) = state_energy
464 buffer%nstored = min(buffer%nstored + 1, buffer%nbuffer)
465 pushed = .true.
466
467 END SUBROUTINE push_paired_state
468
469! **************************************************************************************************
470!> \brief Compute Tr(P F), summed across every spin and real-space cell.
471!> \param density Density matrices.
472!> \param fock Fock matrices.
473!> \param value Summed trace contraction.
474! **************************************************************************************************
475 SUBROUTINE paired_matrix_dot(density, fock, value)
476
477 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: density, fock
478 REAL(kind=dp), INTENT(OUT) :: value
479
480 INTEGER :: icell, ispin
481 REAL(kind=dp) :: contribution
482
483 value = 0.0_dp
484 DO icell = 1, SIZE(fock, 2)
485 DO ispin = 1, SIZE(fock, 1)
486 CALL dbcsr_dot(density(ispin, icell)%matrix, fock(ispin, icell)%matrix, contribution)
487 value = value + contribution
488 END DO
489 END DO
490
491 END SUBROUTINE paired_matrix_dot
492
493! **************************************************************************************************
494!> \brief Compute both cross contractions between one old and the current state.
495!> \param buffer Persistent SCF subspace history.
496!> \param islot Physical old-history slot.
497!> \param current_fock Current F[P].
498!> \param current_density Current P.
499!> \param old_current Tr(P_old F_current).
500!> \param current_old Tr(P_current F_old).
501! **************************************************************************************************
502 SUBROUTINE history_current_dot(buffer, islot, current_fock, current_density, old_current, current_old)
503
504 TYPE(qs_scf_subspace_buffer_type), INTENT(IN) :: buffer
505 INTEGER, INTENT(IN) :: islot
506 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: current_fock, current_density
507 REAL(kind=dp), INTENT(OUT) :: old_current, current_old
508
509 INTEGER :: icell, ispin
510 REAL(kind=dp) :: contribution
511
512 old_current = 0.0_dp
513 current_old = 0.0_dp
514 DO icell = 1, SIZE(current_fock, 2)
515 DO ispin = 1, SIZE(current_fock, 1)
516 CALL dbcsr_dot(buffer%density(islot, ispin, icell)%matrix, &
517 current_fock(ispin, icell)%matrix, contribution)
518 old_current = old_current + contribution
519 CALL dbcsr_dot(current_density(ispin, icell)%matrix, &
520 buffer%fock(islot, ispin, icell)%matrix, contribution)
521 current_old = current_old + contribution
522 END DO
523 END DO
524
525 END SUBROUTINE history_current_dot
526
527! **************************************************************************************************
528!> \brief Return physical slots ordered from oldest to newest.
529!> \param buffer Persistent SCF subspace history.
530!> \param slots Chronologically ordered physical slot indices.
531! **************************************************************************************************
532 PURE SUBROUTINE chronological_slots(buffer, slots)
533
534 TYPE(qs_scf_subspace_buffer_type), INTENT(IN) :: buffer
535 INTEGER, DIMENSION(:), INTENT(OUT) :: slots
536
537 INTEGER :: i, islot, j, tmp
538
539 j = 0
540 DO islot = 1, buffer%nbuffer
541 IF (buffer%generation(islot) <= 0) cycle
542 j = j + 1
543 IF (j <= SIZE(slots)) slots(j) = islot
544 END DO
545
546 ! The ADIIS model uses the most recent state as its expansion point, so
547 ! return physical slots ordered by insertion generation with the newest
548 ! state last. The history size is small, making insertion sort sufficient.
549 DO i = 2, SIZE(slots)
550 tmp = slots(i)
551 j = i - 1
552 DO WHILE (j >= 1)
553 IF (buffer%generation(slots(j)) <= buffer%generation(tmp)) EXIT
554 slots(j + 1) = slots(j)
555 j = j - 1
556 END DO
557 slots(j + 1) = tmp
558 END DO
559
560 END SUBROUTINE chronological_slots
561
562! **************************************************************************************************
563!> \brief Zero the effective KS-matrix combination.
564!> \param buffer Persistent SCF subspace history and combination target.
565! **************************************************************************************************
566 SUBROUTINE zero_fock_combination(buffer)
567
568 TYPE(qs_scf_subspace_buffer_type), INTENT(INOUT) :: buffer
569
570 INTEGER :: icell, ispin
571
572 DO icell = 1, SIZE(buffer%combined_fock, 2)
573 DO ispin = 1, SIZE(buffer%combined_fock, 1)
574 CALL dbcsr_set(buffer%combined_fock(ispin, icell)%matrix, 0.0_dp)
575 END DO
576 END DO
577
578 END SUBROUTINE zero_fock_combination
579
580! **************************************************************************************************
581!> \brief Add one physical history slot to the effective KS-matrix combination.
582!> \param buffer Persistent SCF subspace history and combination target.
583!> \param physical Physical history slot.
584!> \param coefficient Weight of this history matrix.
585! **************************************************************************************************
586 SUBROUTINE add_fock_to_combination(buffer, physical, coefficient)
587
588 TYPE(qs_scf_subspace_buffer_type), INTENT(INOUT) :: buffer
589 INTEGER, INTENT(IN) :: physical
590 REAL(kind=dp), INTENT(IN) :: coefficient
591
592 INTEGER :: icell, ispin
593
594 DO icell = 1, SIZE(buffer%combined_fock, 2)
595 DO ispin = 1, SIZE(buffer%combined_fock, 1)
596 CALL dbcsr_add(buffer%combined_fock(ispin, icell)%matrix, &
597 buffer%fock(physical, ispin, icell)%matrix, &
598 alpha_scalar=1.0_dp, beta_scalar=coefficient)
599 END DO
600 END DO
601
602 END SUBROUTINE add_fock_to_combination
603
604END MODULE qs_scf_subspace
cp_dbcsr_api::dbcsr_create
Definition cp_dbcsr_api.F:194
cp_dbcsr_api
Definition cp_dbcsr_api.F:8
cp_dbcsr_api::dbcsr_copy
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
Definition cp_dbcsr_api.F:371
cp_dbcsr_api::dbcsr_set
subroutine, public dbcsr_set(matrix, alpha)
...
Definition cp_dbcsr_api.F:1204
cp_dbcsr_api::dbcsr_add
subroutine, public dbcsr_add(matrix_a, matrix_b, alpha_scalar, beta_scalar)
...
Definition cp_dbcsr_api.F:254
cp_dbcsr_contrib
Definition cp_dbcsr_contrib.F:8
cp_dbcsr_contrib::dbcsr_dot
subroutine, public dbcsr_dot(matrix_a, matrix_b, trace)
Computes the dot product of two matrices, also known as the trace of their matrix product.
Definition cp_dbcsr_contrib.F:367
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
qs_scf_subspace_math
Small, matrix-free mathematical kernels used by ADIIS.
Definition qs_scf_subspace_math.F:18
qs_scf_subspace_math::simplex_qp_success
integer, parameter, public simplex_qp_success
Definition qs_scf_subspace_math.F:27
qs_scf_subspace_math::qs_scf_subspace_build_adiis_model
subroutine, public qs_scf_subspace_build_adiis_model(pf_metric, newest, hessian, linear, valid)
Construct the canonical ADIIS model relative to the newest history entry.
Definition qs_scf_subspace_math.F:158
qs_scf_subspace_math::simplex_qp_nonfinite_input
integer, parameter, public simplex_qp_nonfinite_input
Definition qs_scf_subspace_math.F:29
qs_scf_subspace_math::simplex_quadratic_minimize
subroutine, public simplex_quadratic_minimize(hessian, linear, coeff, objective, status, preferred_index)
Minimize a quadratic model over the probability simplex.
Definition qs_scf_subspace_math.F:61
qs_scf_subspace_math::qs_scf_subspace_fifo_slot
pure integer function, public qs_scf_subspace_fifo_slot(generation)
Select the next physical slot for a strict FIFO history.
Definition qs_scf_subspace_math.F:207
qs_scf_subspace_math::simplex_qp_pairwise_stationary
integer, parameter, public simplex_qp_pairwise_stationary
Definition qs_scf_subspace_math.F:31
qs_scf_subspace_types
Data types for the ADIIS SCF subspace accelerator.
Definition qs_scf_subspace_types.F:11
qs_scf_subspace_types::qs_scf_subspace_buffer_clear
pure subroutine, public qs_scf_subspace_buffer_clear(buffer)
Clear logical history while retaining allocated matrix storage.
Definition qs_scf_subspace_types.F:87
qs_scf_subspace
History management and matrix combination for ADIIS.
Definition qs_scf_subspace.F:11
qs_scf_subspace::scf_subspace_invalid_history
integer, parameter, public scf_subspace_invalid_history
Definition qs_scf_subspace.F:37
qs_scf_subspace::qs_scf_subspace_build
subroutine, public qs_scf_subspace_build(buffer)
Solve the ADIIS model from previously accepted history and form an effective KS matrix.
Definition qs_scf_subspace.F:130
qs_scf_subspace::qs_scf_subspace_restart
subroutine, public qs_scf_subspace_restart(buffer, fock, density, state_energy, restarted)
Restart the ADIIS history from the current strictly paired P, F[P] state.
Definition qs_scf_subspace.F:52
qs_scf_subspace::qs_scf_subspace_push
subroutine, public qs_scf_subspace_push(buffer, fock, density, state_energy, pushed)
Append one accepted, evaluated P,F[P] state to the ADIIS history.
Definition qs_scf_subspace.F:85
cp_dbcsr_api::dbcsr_p_type
Definition cp_dbcsr_api.F:172
qs_scf_subspace_types::qs_scf_subspace_buffer_type
History buffer holding strictly paired P and F[P] SCF states.
Definition qs_scf_subspace_types.F:36