(git:a145afa)
Loading...
Searching...
No Matches
hfx_ace_methods.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 Adaptively Compressed Exchange (ACE) operator for HFX.
10!> Reference: Lin, J. Chem. Theory Comput. 2016, 12, 5, 2242-2249
11!>
12!> Algorithm (per spin):
13!>
14!> BUILD (first call, or every rebuild_frequency steps):
15!> 1. Full HFX: ks_matrix = K_HFX + H_core, energy%ex = E_x(full)
16!> 2. K_AO = ks_matrix - H_core (nao x nao, negative semidefinite)
17!> 3. C_occ = first nocc columns of mo_coeff, redistributed to a layout
18!> compatible with K_AO so that PDGEMM works correctly
19!> 4. xi = K_AO * C_occ (nao x nocc)
20!> 5. M = C_occ^T * xi (nocc x nocc, negative definite)
21!> 6. -M = U^T U via Cholesky (U upper triangular, stored in M_fm)
22!> 7. W = xi * U^{-1} (nao x nocc, the ACE projector)
23!> 8. Apply (see below) to update ks_matrix and energy%ex
24!>
25!> APPLY (all other steps):
26!> ks_matrix = H_core - W * W^T
27!> E_x = -0.5 * Tr[W^T * P * W]
28!>
29!> Diagnostics (controlled by DBG_STALE / DBG_EXACT_EX module flags):
30!>
31!> DIAG A projector staleness (cheap, runs every APPLY step)
32!> Computes ||W^T C_occ^current||_F / ||W^T C_occ^BUILD||_F.
33!> Ratio = 1 → W still accurate. Ratio -> 0 → W is stale.
34!>
35!> DIAG B exact vs ACE exchange energy (expensive: one full HFX per APPLY)
36!> Calls full HFX with just_energy=.TRUE. to get E_x^exact[P^k] and
37!> compares to E_x^ACE[P^k]. Growing |delta| confirms stale W.
38!> ACE ks_matrix and energy%ex are restored after the diagnostic.
39!>
40!> \author Ritama Kar
41! **************************************************************************************************
42
43MODULE hfx_ace_methods
44
45 USE admm_types, ONLY: admm_type,&
46 get_admm_env
47 USE bibliography, ONLY: lin2016ace,&
48 cite_reference
49 USE cp_blacs_env, ONLY: cp_blacs_env_type
50 USE cp_control_types, ONLY: dft_control_type
51 USE cp_dbcsr_api, ONLY: dbcsr_add,&
52 dbcsr_copy,&
53 dbcsr_create,&
54 dbcsr_p_type,&
55 dbcsr_release,&
56 dbcsr_set,&
57 dbcsr_type
58 USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
59 cp_dbcsr_plus_fm_fm_t
60 USE cp_fm_basic_linalg, ONLY: cp_fm_scale,&
61 cp_fm_trace,&
62 cp_fm_triangular_multiply
63 USE cp_fm_cholesky, ONLY: cp_fm_cholesky_decompose
64 USE cp_fm_struct, ONLY: cp_fm_struct_create,&
65 cp_fm_struct_release,&
66 cp_fm_struct_type
67 USE cp_fm_types, ONLY: cp_fm_create,&
68 cp_fm_get_info,&
69 cp_fm_release,&
70 cp_fm_to_fm,&
71 cp_fm_type
72 USE cp_log_handling, ONLY: cp_get_default_logger,&
73 cp_logger_get_default_io_unit,&
74 cp_logger_type
75 USE hfx_admm_utils, ONLY: hfx_ks_matrix
76 USE hfx_types, ONLY: hfx_type
77 USE input_section_types, ONLY: section_vals_type
78 USE kinds, ONLY: dp
79 USE message_passing, ONLY: mp_para_env_type
80 USE parallel_gemm_api, ONLY: parallel_gemm
81 USE pw_types, ONLY: pw_r3d_rs_type
82 USE qs_energy_types, ONLY: qs_energy_type
83 USE qs_environment_types, ONLY: get_qs_env,&
84 qs_environment_type
85 USE qs_mo_types, ONLY: get_mo_set,&
86 mo_set_type
87 USE qs_rho_types, ONLY: qs_rho_get,&
88 qs_rho_type
89 USE scf_control_types, ONLY: scf_control_type
90#include "./base/base_uses.f90"
91
92 IMPLICIT NONE
93 PRIVATE
94
95 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'hfx_ace_methods'
96
97 ! -----------------------------------------------------------------------
98 ! Module-level state: persists across SCF steps within one run.
99 !
100 ! ace_W(1,ispin) ACE projector W, shape nao x nocc
101 ! ace_is_built .FALSE. until at least one successful BUILD
102 ! ace_step_counter counts calls since last BUILD
103 ! ace_W_ref_norm ||W^T C_occ^BUILD||_F stored at BUILD for DIAG A
104 ! -----------------------------------------------------------------------
105 TYPE(cp_fm_type), ALLOCATABLE, SAVE :: ace_W(:, :)
106 LOGICAL, SAVE :: ace_is_built = .false.
107 INTEGER, SAVE :: ace_step_counter = 0
108 REAL(dp), SAVE :: ace_W_ref_norm = 0.0_dp
109 INTEGER, SAVE :: ace_geo_step = 0 ! NEW
110
111 ! -----------------------------------------------------------------------
112 ! Debug / diagnostic flags — set all .FALSE. for production.
113 !
114 ! DBG_ROUTING BUILD/APPLY/DEFER decisions, counters
115 ! DBG_BUILD norms during BUILD
116 ! DBG_ENERGY E_x(ACE) at every step; DIAG C on BUILD steps
117 ! DBG_STALE DIAG A: staleness ratio at every APPLY step (cheap)
118 ! DBG_EXACT_EX DIAG B: full HFX energy at every APPLY step (expensive)
119 ! -----------------------------------------------------------------------
120 LOGICAL, PARAMETER, PRIVATE :: DBG_ROUTING = .false.
121 LOGICAL, PARAMETER, PRIVATE :: DBG_BUILD = .false.
122 LOGICAL, PARAMETER, PRIVATE :: DBG_ENERGY = .false.
123 LOGICAL, PARAMETER, PRIVATE :: DBG_STALE = .false.
124 LOGICAL, PARAMETER, PRIVATE :: DBG_EXACT_EX = .false.
125
126 LOGICAL, SAVE :: ace_dynamic_mode = .false.
127 ! Set to .TRUE. by hfx_ace_set_dynamic_mode before geo_opt/MD starts.
128 ! Stays .FALSE. for ENERGY/ENERGY_FORCE single-point runs.
129
130 PUBLIC :: hfx_ace_ks_matrix, hfx_ace_release, hfx_ace_set_dynamic_mode
131
132CONTAINS
133
134! **************************************************************************************************
135!> \brief Main ACE entry point, replacing hfx_ks_matrix in qs_ks_methods.
136!> \param qs_env ...
137!> \param ks_matrix ...
138!> \param rho ...
139!> \param energy ...
140!> \param calculate_forces ...
141!> \param just_energy ...
142!> \param v_rspace_new ...
143!> \param v_tau_rspace ...
144!> \param ace_rebuild_frequency ...
145!> \param ext_xc_section ...
146! **************************************************************************************************
147 SUBROUTINE hfx_ace_ks_matrix(qs_env, ks_matrix, rho, energy, &
148 calculate_forces, just_energy, &
149 v_rspace_new, v_tau_rspace, &
150 ace_rebuild_frequency, ext_xc_section)
151
152 TYPE(qs_environment_type), POINTER :: qs_env
153 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: ks_matrix
154 TYPE(qs_rho_type), POINTER :: rho
155 TYPE(qs_energy_type), POINTER :: energy
156 LOGICAL, INTENT(IN) :: calculate_forces, just_energy
157 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: v_rspace_new, v_tau_rspace
158 INTEGER, INTENT(IN) :: ace_rebuild_frequency
159 TYPE(section_vals_type), OPTIONAL, POINTER :: ext_xc_section
160
161 CHARACTER(LEN=*), PARAMETER :: routinen = 'hfx_ace_ks_matrix'
162
163 INTEGER :: handle, iw, n_rep_hf, nspins, &
164 rebuild_freq
165 LOGICAL :: ace_built_now, rebuild_ace
166 REAL(dp) :: ex_ace
167 TYPE(cp_logger_type), POINTER :: logger
168 TYPE(dft_control_type), POINTER :: dft_control
169 TYPE(hfx_type), DIMENSION(:, :), POINTER :: x_data
170 TYPE(scf_control_type), POINTER :: scf_control
171
172 CALL timeset(routinen, handle)
173
174 CALL cite_reference(lin2016ace)
175 NULLIFY (logger, dft_control, x_data, scf_control)
176
177 logger => cp_get_default_logger()
178 iw = cp_logger_get_default_io_unit(logger)
179
180 CALL get_qs_env(qs_env, x_data=x_data, dft_control=dft_control)
181 n_rep_hf = SIZE(x_data, 1)
182 nspins = dft_control%nspins
183
184 IF (n_rep_hf /= 1) cpabort("ACE: only one &HF section is supported.")
185 IF (dft_control%nimages /= 1) THEN
186 cpabort("ACE: k-points / multiple images are not implemented.")
187 END IF
188
189 ! ACE requires explicit MO coefficients (C_occ) which are only available
190 ! with diagonalization-based SCF. OT never constructs mo_coeff during
191 ! the SCF, so the ACE projector build loop would silently get garbage.
192 CALL get_qs_env(qs_env, scf_control=scf_control)
193 IF (scf_control%use_ot) THEN
194 cpabort("ACE: OT doesn't work, use diagonalization-based SCF.")
195 END IF
196
197 rebuild_freq = max(1, ace_rebuild_frequency)
198
199 ! ------------------------------------------------------------------
200 ! Bypass A: energy-only call
201 ! ------------------------------------------------------------------
202 IF (just_energy) THEN
203 IF (dbg_routing .AND. iw > 0) THEN
204 WRITE (iw, '(T2,A)') 'ACE | just_energy=T: full HFX (no matrix update)'
205 END IF
206 CALL hfx_call(qs_env, ks_matrix, rho, energy, &
207 calculate_forces, just_energy, &
208 v_rspace_new, v_tau_rspace, ext_xc_section)
209 CALL timestop(handle)
210 RETURN
211 END IF
212
213 ! ------------------------------------------------------------------
214 ! Bypass B: ionic forces requested
215 ! ------------------------------------------------------------------
216 IF (calculate_forces) THEN
217 IF (dbg_routing .AND. iw > 0) THEN
218 WRITE (iw, '(T2,A)') 'ACE | calculate_forces=T: full HFX for exact forces'
219 END IF
220 CALL hfx_call(qs_env, ks_matrix, rho, energy, &
221 calculate_forces, just_energy, &
222 v_rspace_new, v_tau_rspace, ext_xc_section)
223 ace_is_built = .false.
224 ace_step_counter = 0
225 ace_geo_step = ace_geo_step + 1 ! NEW: step 0 done, ACE active from now
226 CALL timestop(handle)
227 RETURN
228 END IF
229
230 ! ------------------------------------------------------------------
231 ! Bypass C: first geometry step.
232 !
233 ! The ATOMIC initial guess gives C_occ far from self-consistency,
234 ! which would produce an inaccurate W. Running full HFX for the
235 ! entire first geometry step ensures that wavefunction extrapolation
236 ! delivers a near-converged C_occ to geometry step 1, making the
237 ! ACE BUILD there accurate from the first application.
238 ! ------------------------------------------------------------------
239 IF (ace_geo_step == 0 .AND. ace_dynamic_mode) THEN
240 IF (iw > 0) WRITE (iw, '(T2,A)') &
241 'ACE | geo_step=0 (MD/GEO_OPT): full HFX for reference wavefunction'
242 CALL hfx_call(qs_env, ks_matrix, rho, energy, &
243 .false., just_energy, &
244 v_rspace_new, v_tau_rspace, ext_xc_section)
245 CALL timestop(handle)
246 RETURN
247 END IF
248
249 ! ------------------------------------------------------------------
250 ! Rebuild decision
251 ! ------------------------------------------------------------------
252 rebuild_ace = (.NOT. ace_is_built) .OR. &
253 (mod(ace_step_counter, rebuild_freq) == 0)
254
255 IF (dbg_routing .AND. iw > 0) THEN
256 WRITE (iw, '(/,T2,A)') repeat('-', 56)
257 WRITE (iw, '(T2,A)') 'ACE | hfx_ace_ks_matrix'
258 WRITE (iw, '(T4,A,L1)') 'ace_is_built = ', ace_is_built
259 WRITE (iw, '(T4,A,L1)') 'rebuild_ace = ', rebuild_ace
260 WRITE (iw, '(T4,A,I6)') 'step_counter = ', ace_step_counter
261 WRITE (iw, '(T4,A,I6)') 'rebuild_freq = ', rebuild_freq
262 WRITE (iw, '(T4,A,I4)') 'nspins = ', nspins
263 WRITE (iw, '(T4,A)') merge('-> BUILD', '-> APPLY', rebuild_ace)
264 WRITE (iw, '(T2,A)') repeat('-', 56)
265 END IF
266
267 IF (rebuild_ace) THEN
268
269 ace_built_now = .false.
270 CALL hfx_ace_build_projector(qs_env, ks_matrix, rho, energy, &
271 just_energy, &
272 v_rspace_new, v_tau_rspace, &
273 nspins, iw, ace_built_now, &
274 ext_xc_section)
275 IF (ace_built_now) THEN
276 ace_is_built = .true.
277 ace_step_counter = 1
278 IF (dbg_routing .AND. iw > 0) THEN
279 WRITE (iw, '(T4,A)') 'ACE | W built. Projector live from next step.'
280 END IF
281 ELSE
282 ace_is_built = .false.
283 ace_step_counter = 0
284 IF (dbg_routing .AND. iw > 0) THEN
285 WRITE (iw, '(T4,A)') 'ACE | Build deferred (C_occ=0). Full HFX in ks_matrix.'
286 END IF
287 END IF
288
289 ELSE
290
291 CALL hfx_ace_apply_projector(qs_env, ks_matrix, rho, energy, nspins, iw)
292 ace_step_counter = ace_step_counter + 1
293
294 ! ----------------------------------------------------------------
295 ! DIAGNOSTIC B: compare E_x^ACE[P^k] with E_x^exact[P^k].
296 !
297 ! Calls full HFX (just_energy=.TRUE.) to get the exact exchange
298 ! energy at the current ACE-converging density P^k. Compares
299 ! with E_x^ACE[P^k] already stored in energy%ex.
300 !
301 ! Growing |delta| over the SCF confirms the root cause: W was
302 ! built from C_occ^(step 1), which is far from the converged
303 ! C_occ, so K_ACE = -WW^T no longer represents K_x accurately.
304 !
305 ! After the comparison, hfx_ace_apply_projector is called a
306 ! second time to restore ks_matrix and energy%ex to ACE values
307 ! so the SCF continues correctly. Cost: +2 full HFX per step.
308 ! ----------------------------------------------------------------
309 IF (dbg_exact_ex) THEN
310 ex_ace = energy%ex
311
312 CALL hfx_call(qs_env, ks_matrix, rho, energy, &
313 .false., .true., &
314 v_rspace_new, v_tau_rspace, ext_xc_section)
315 IF (iw > 0) THEN
316 WRITE (iw, '(/,T2,A)') repeat('-', 56)
317 WRITE (iw, '(T2,A,I6)') 'ACE DIAG B | ace_step_counter = ', ace_step_counter
318 WRITE (iw, '(T4,A,F20.10)') 'E_x(exact, P^k) = ', energy%ex
319 WRITE (iw, '(T4,A,F20.10)') 'E_x(ACE, P^k) = ', ex_ace
320 WRITE (iw, '(T4,A,ES12.4)') '|delta| = ', abs(ex_ace - energy%ex)
321 WRITE (iw, '(T4,A)') &
322 '|delta|->0 on BUILD step; growth confirms stale projector'
323 WRITE (iw, '(T2,A)') repeat('-', 56)
324 END IF
325
326 ! Restore ACE ks_matrix and energy%ex
327 CALL hfx_ace_apply_projector(qs_env, ks_matrix, rho, energy, nspins, iw)
328 END IF
329
330 END IF
331
332 IF (dbg_routing .AND. iw > 0) THEN
333 WRITE (iw, '(T4,A,F20.10)') 'energy%ex on exit = ', energy%ex
334 WRITE (iw, '(T4,A,I6)') 'step_counter now = ', ace_step_counter
335 END IF
336
337 CALL timestop(handle)
338
339 END SUBROUTINE hfx_ace_ks_matrix
340
341! **************************************************************************************************
342!> \brief Build the ACE projector W.
343!> \param qs_env ...
344!> \param ks_matrix ...
345!> \param rho ...
346!> \param energy ...
347!> \param just_energy ...
348!> \param v_rspace_new ...
349!> \param v_tau_rspace ...
350!> \param nspins ...
351!> \param iw ...
352!> \param build_succeeded ...
353!> \param ext_xc_section ...
354! **************************************************************************************************
355 SUBROUTINE hfx_ace_build_projector(qs_env, ks_matrix, rho, energy, &
356 just_energy, &
357 v_rspace_new, v_tau_rspace, &
358 nspins, iw, build_succeeded, &
359 ext_xc_section)
360
361 TYPE(qs_environment_type), POINTER :: qs_env
362 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: ks_matrix
363 TYPE(qs_rho_type), POINTER :: rho
364 TYPE(qs_energy_type), POINTER :: energy
365 LOGICAL, INTENT(IN) :: just_energy
366 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: v_rspace_new, v_tau_rspace
367 INTEGER, INTENT(IN) :: nspins, iw
368 LOGICAL, INTENT(OUT) :: build_succeeded
369 TYPE(section_vals_type), OPTIONAL, POINTER :: ext_xc_section
370
371 CHARACTER(LEN=*), PARAMETER :: routinen = 'hfx_ace_build_projector'
372
373 INTEGER :: handle, info_chol, ispin, nao, nmo, nocc
374 LOGICAL :: do_admm
375 REAL(dp) :: ehfx_full, frob
376 REAL(dp), DIMENSION(:), POINTER :: occ_nums
377 TYPE(admm_type), POINTER :: admm_env
378 TYPE(cp_blacs_env_type), POINTER :: blacs_env
379 TYPE(cp_fm_struct_type), POINTER :: fmstruct
380 TYPE(cp_fm_type) :: a_ref_fm, c_occ_fm, k_fm, m_fm, xi_fm
381 TYPE(cp_fm_type), POINTER :: mo_coeff
382 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: ks_aux_fit
383 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_h
384 TYPE(dbcsr_type) :: k_ao_dbcsr
385 TYPE(dft_control_type), POINTER :: dft_control
386 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos, mos_for_ace
387 TYPE(mp_para_env_type), POINTER :: para_env
388
389! A_ref_fm: temporary nocc x nocc scratch for DIAG A reference norm
390
391 CALL timeset(routinen, handle)
392 NULLIFY (blacs_env, para_env, mos, mo_coeff, matrix_h, occ_nums, fmstruct)
393 NULLIFY (dft_control, admm_env, mos_for_ace, ks_aux_fit)
394
395 build_succeeded = .false.
396
397 CALL get_qs_env(qs_env, blacs_env=blacs_env, para_env=para_env, &
398 mos=mos, matrix_h_kp=matrix_h, &
399 dft_control=dft_control)
400
401 do_admm = dft_control%do_admm
402
403 IF (do_admm) THEN
404 CALL get_qs_env(qs_env, admm_env=admm_env)
405 CALL get_admm_env(admm_env, &
406 matrix_ks_aux_fit=ks_aux_fit, &
407 mos_aux_fit=mos_for_ace)
408 ELSE
409 mos_for_ace => mos
410 END IF
411
412 ! Step 1: full HFX
413 CALL hfx_call(qs_env, ks_matrix, rho, energy, &
414 .false., just_energy, &
415 v_rspace_new, v_tau_rspace, ext_xc_section)
416 ehfx_full = energy%ex
417
418 IF (dbg_build .AND. iw > 0) THEN
419 WRITE (iw, '(/,T2,A,F20.10)') 'ACE BUILD | E_x(full HFX) = ', ehfx_full
420 END IF
421
422 ! Allocate / resize module storage
423 IF (ALLOCATED(ace_w)) THEN
424 IF (SIZE(ace_w, 2) /= nspins) CALL hfx_ace_release()
425 END IF
426 IF (.NOT. ALLOCATED(ace_w)) ALLOCATE (ace_w(1, nspins))
427
428 ! Reset reference norm; accumulated per spin in the loop below
429 ace_w_ref_norm = 0.0_dp
430
431 ! ----------------------------------------------------------------
432 ! Per-spin build loop
433 ! ----------------------------------------------------------------
434 DO ispin = 1, nspins
435
436 IF (mos_for_ace(ispin)%use_mo_coeff_b) THEN
437 CALL copy_dbcsr_to_fm(mos_for_ace(ispin)%mo_coeff_b, &
438 mos_for_ace(ispin)%mo_coeff)
439 END IF
440
441 CALL get_mo_set(mos_for_ace(ispin), mo_coeff=mo_coeff, &
442 nao=nao, nmo=nmo, homo=nocc, &
443 occupation_numbers=occ_nums)
444
445 IF (nocc <= 0) cpabort("ACE: homo <= 0.")
446 IF (nocc > nmo) cpabort("ACE: homo > nmo.")
447 cpassert(ASSOCIATED(mo_coeff))
448
449 CALL cp_fm_trace(mo_coeff, mo_coeff, frob)
450
451 IF (dbg_build .AND. iw > 0) THEN
452 WRITE (iw, '(/,T2,A,I4)') 'ACE BUILD | ispin = ', ispin
453 WRITE (iw, '(T4,A,I8)') 'nao = ', nao
454 WRITE (iw, '(T4,A,I8)') 'nmo = ', nmo
455 WRITE (iw, '(T4,A,I8)') 'nocc = ', nocc
456 WRITE (iw, '(T4,A,L1)') 'use_mo_coeff_b = ', mos_for_ace(ispin)%use_mo_coeff_b
457 WRITE (iw, '(T4,A,ES12.4)') '||mo_coeff||_F = ', sqrt(max(frob, 0.0_dp))
458 END IF
459
460 IF (frob < 1.0e-20_dp) THEN
461 IF (dbg_build .AND. iw > 0) THEN
462 WRITE (iw, '(T4,A)') 'mo_coeff=0: build deferred to next step.'
463 END IF
464 CALL timestop(handle)
465 RETURN
466 END IF
467
468 ! Step 2: K_AO
469 IF (do_admm) THEN
470 CALL dbcsr_create(k_ao_dbcsr, template=ks_aux_fit(ispin)%matrix, &
471 name="K_ACE_aux")
472 CALL dbcsr_copy(k_ao_dbcsr, ks_aux_fit(ispin)%matrix)
473 ELSE
474 CALL dbcsr_create(k_ao_dbcsr, template=ks_matrix(ispin, 1)%matrix, &
475 name="K_AO")
476 CALL dbcsr_copy(k_ao_dbcsr, ks_matrix(ispin, 1)%matrix)
477 CALL dbcsr_add(k_ao_dbcsr, matrix_h(1, 1)%matrix, 1.0_dp, -1.0_dp)
478 END IF
479
480 NULLIFY (fmstruct)
481 CALL cp_fm_struct_create(fmstruct, context=blacs_env, para_env=para_env, &
482 nrow_global=nao, ncol_global=nao)
483 CALL cp_fm_create(k_fm, fmstruct, name="K_dense")
484 CALL cp_fm_struct_release(fmstruct)
485 CALL copy_dbcsr_to_fm(k_ao_dbcsr, k_fm)
486 CALL dbcsr_release(k_ao_dbcsr)
487
488 IF (dbg_build .AND. iw > 0) THEN
489 CALL cp_fm_trace(k_fm, k_fm, frob)
490 WRITE (iw, '(T4,A,ES12.4)') '||K_AO||_F = ', sqrt(max(frob, 0.0_dp))
491 END IF
492
493 ! Step 3: C_occ
494 NULLIFY (fmstruct)
495 CALL cp_fm_struct_create(fmstruct, context=blacs_env, para_env=para_env, &
496 nrow_global=nao, ncol_global=nocc)
497 CALL cp_fm_create(c_occ_fm, fmstruct, name="C_occ")
498 CALL cp_fm_create(xi_fm, fmstruct, name="xi")
499 CALL cp_fm_struct_release(fmstruct)
500
501 CALL cp_fm_to_fm(mo_coeff, c_occ_fm)
502
503 ! Step 4: xi = K_AO * C_occ
504 CALL parallel_gemm('N', 'N', nao, nocc, nao, &
505 1.0_dp, k_fm, c_occ_fm, 0.0_dp, xi_fm)
506 CALL cp_fm_release(k_fm)
507
508 IF (dbg_build .AND. iw > 0) THEN
509 CALL cp_fm_trace(xi_fm, xi_fm, frob)
510 WRITE (iw, '(T4,A,ES12.4)') '||xi||_F = ', sqrt(max(frob, 0.0_dp))
511 END IF
512
513 ! Step 5: M = C_occ^T * xi
514 NULLIFY (fmstruct)
515 CALL cp_fm_struct_create(fmstruct, context=blacs_env, para_env=para_env, &
516 nrow_global=nocc, ncol_global=nocc)
517 CALL cp_fm_create(m_fm, fmstruct, name="M")
518 CALL cp_fm_struct_release(fmstruct)
519
520 CALL parallel_gemm('T', 'N', nocc, nocc, nao, &
521 1.0_dp, c_occ_fm, xi_fm, 0.0_dp, m_fm)
522 CALL cp_fm_release(c_occ_fm)
523
524 IF (dbg_build .AND. iw > 0) THEN
525 CALL cp_fm_trace(m_fm, m_fm, frob)
526 WRITE (iw, '(T4,A,ES12.4)') '||M||_F = ', sqrt(max(frob, 0.0_dp))
527 END IF
528
529 ! Step 6: Cholesky of -M = U^T U
530 CALL cp_fm_scale(-1.0_dp, m_fm)
531 CALL cp_fm_cholesky_decompose(m_fm, n=nocc, info_out=info_chol)
532
533 IF (info_chol /= 0) THEN
534 IF (iw > 0) THEN
535 WRITE (iw, '(T4,A,I6)') 'ACE | Cholesky failed, info = ', info_chol
536 WRITE (iw, '(T4,A,F20.10)') 'ACE | E_x(full) = ', ehfx_full
537 WRITE (iw, '(T4,A,I8,A,I8)') 'ACE | nao=', nao, ' nocc=', nocc
538 END IF
539 cpabort("ACE: Cholesky of -M failed (not positive definite).")
540 END IF
541
542 IF (dbg_build .AND. iw > 0) THEN
543 WRITE (iw, '(T4,A)') 'Cholesky OK (info=0).'
544 END IF
545
546 ! Step 7: W = xi * U^{-1}
547 IF (ASSOCIATED(ace_w(1, ispin)%matrix_struct)) THEN
548 CALL cp_fm_release(ace_w(1, ispin))
549 END IF
550
551 CALL cp_fm_create(ace_w(1, ispin), xi_fm%matrix_struct, name="W_ACE")
552 CALL cp_fm_to_fm(xi_fm, ace_w(1, ispin))
553
554 CALL cp_fm_triangular_multiply(m_fm, ace_w(1, ispin), &
555 side='R', uplo_tr='U', &
556 transpose_tr=.false., &
557 invert_tr=.true., &
558 n_rows=nao, n_cols=nocc, &
559 alpha=1.0_dp)
560
561 CALL cp_fm_release(xi_fm)
562 CALL cp_fm_release(m_fm)
563
564 IF (dbg_build .AND. iw > 0) THEN
565 CALL cp_fm_trace(ace_w(1, ispin), ace_w(1, ispin), frob)
566 WRITE (iw, '(T4,A,I4,A,2I8,A,ES12.4)') &
567 'W spin=', ispin, ' shape=', nao, nocc, &
568 ' ||W||_F=', sqrt(max(frob, 0.0_dp))
569 END IF
570
571 ! ----------------------------------------------------------------
572 ! DIAG A reference: compute ||W^T C_occ^BUILD||_F for this spin.
573 !
574 ! C_occ_fm was released after step 5, but mo_coeff is still valid
575 ! (it is a pointer into mos_for_ace, not allocated here).
576 ! We re-create C_occ_fm from mo_coeff.
577 !
578 ! Theory: W^T C_occ^BUILD = U^{-T}(-U^T U) = -U → norm = ||U||_F
579 ! (computed directly rather than storing U).
580 ! ----------------------------------------------------------------
581 IF (dbg_stale) THEN
582 NULLIFY (fmstruct)
583 CALL cp_fm_struct_create(fmstruct, context=blacs_env, &
584 para_env=para_env, &
585 nrow_global=nao, ncol_global=nocc)
586 CALL cp_fm_create(c_occ_fm, fmstruct, name="C_occ_ref_diag")
587 CALL cp_fm_struct_release(fmstruct)
588 CALL cp_fm_to_fm(mo_coeff, c_occ_fm)
589
590 NULLIFY (fmstruct)
591 CALL cp_fm_struct_create(fmstruct, context=blacs_env, &
592 para_env=para_env, &
593 nrow_global=nocc, ncol_global=nocc)
594 CALL cp_fm_create(a_ref_fm, fmstruct, name="WtC_ref")
595 CALL cp_fm_struct_release(fmstruct)
596
597 CALL parallel_gemm('T', 'N', nocc, nocc, nao, &
598 1.0_dp, ace_w(1, ispin), c_occ_fm, 0.0_dp, a_ref_fm)
599 CALL cp_fm_trace(a_ref_fm, a_ref_fm, frob)
600 ace_w_ref_norm = ace_w_ref_norm + sqrt(max(frob, 0.0_dp))
601
602 CALL cp_fm_release(c_occ_fm)
603 CALL cp_fm_release(a_ref_fm)
604 END IF
605
606 END DO ! ispin
607
608 IF (dbg_stale .AND. iw > 0) THEN
609 WRITE (iw, '(/,T2,A)') repeat('-', 56)
610 WRITE (iw, '(T2,A)') 'ACE DIAG A | Reference norm stored at BUILD'
611 WRITE (iw, '(T4,A,ES12.4)') &
612 '||W^T C_occ^BUILD||_F (sum over spins) = ', ace_w_ref_norm
613 WRITE (iw, '(T4,A)') &
614 'Staleness ratio = 1.0 at BUILD step; decreasing means W is becoming stale'
615 WRITE (iw, '(T2,A)') repeat('-', 56)
616 END IF
617
618 build_succeeded = .true.
619
620 ! Step 8: apply immediately (DIAG C ratio printed via ehfx_full_ref)
621 CALL hfx_ace_apply_projector(qs_env, ks_matrix, rho, energy, &
622 nspins, iw, ehfx_full_ref=ehfx_full)
623
624 CALL timestop(handle)
625
626 END SUBROUTINE hfx_ace_build_projector
627
628! **************************************************************************************************
629!> \brief Apply the stored ACE projector.
630!> \param qs_env ...
631!> \param ks_matrix ...
632!> \param rho ...
633!> \param energy ...
634!> \param nspins ...
635!> \param iw ...
636!> \param ehfx_full_ref ...
637! **************************************************************************************************
638 SUBROUTINE hfx_ace_apply_projector(qs_env, ks_matrix, rho, energy, &
639 nspins, iw, ehfx_full_ref)
640
641 TYPE(qs_environment_type), POINTER :: qs_env
642 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: ks_matrix
643 TYPE(qs_rho_type), POINTER :: rho
644 TYPE(qs_energy_type), POINTER :: energy
645 INTEGER, INTENT(IN) :: nspins, iw
646 REAL(dp), INTENT(IN), OPTIONAL :: ehfx_full_ref
647
648 CHARACTER(LEN=*), PARAMETER :: routinen = 'hfx_ace_apply_projector'
649
650 INTEGER :: handle, ispin, nao, nao_d, nmo_d, nocc, &
651 nocc_d
652 LOGICAL :: do_admm
653 REAL(dp) :: ehfx_ace, frob_a, stale_norm, trace_val
654 TYPE(admm_type), POINTER :: admm_env
655 TYPE(cp_blacs_env_type), POINTER :: blacs_env
656 TYPE(cp_fm_struct_type), POINTER :: fmstruct, fmstruct_diag
657 TYPE(cp_fm_type) :: a_diag_fm, c_occ_diag, p_fm, pw_fm
658 TYPE(cp_fm_type), POINTER :: mo_coeff_diag
659 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: ks_aux_fit, ks_aux_fit_hfx
660 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_h, rho_ao
661 TYPE(dft_control_type), POINTER :: dft_control
662 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos_aux_diag, mos_diag
663 TYPE(mp_para_env_type), POINTER :: para_env
664 TYPE(qs_rho_type), POINTER :: rho_aux_fit
665
666! ------------------------------------------------------------------
667! DIAG A local variables
668! stale_norm ||W^T C_occ^current||_F summed over spins
669! frob_A scratch for cp_fm_trace
670! C_occ_diag current C_occ redistributed to W layout
671! A_diag_fm nocc x nocc overlap W^T C_occ^current
672! mos_diag primary mos (non-ADMM path)
673! mos_aux_diag auxiliary mos (ADMM path)
674! mo_coeff_diag pointer to the relevant mo_coeff
675! nao_d, nmo_d, nocc_d dimensions from get_mo_set
676! ------------------------------------------------------------------
677
678 CALL timeset(routinen, handle)
679 NULLIFY (blacs_env, para_env, matrix_h, rho_ao, fmstruct)
680 NULLIFY (dft_control, admm_env, ks_aux_fit, ks_aux_fit_hfx, rho_aux_fit)
681 NULLIFY (mos_diag, mos_aux_diag, mo_coeff_diag, fmstruct_diag)
682 cpassert(ALLOCATED(ace_w))
683
684 CALL get_qs_env(qs_env, blacs_env=blacs_env, para_env=para_env, &
685 matrix_h_kp=matrix_h, dft_control=dft_control)
686 do_admm = dft_control%do_admm
687
688 IF (do_admm) THEN
689 CALL get_qs_env(qs_env, admm_env=admm_env)
690 CALL get_admm_env(admm_env, &
691 matrix_ks_aux_fit=ks_aux_fit, &
692 matrix_ks_aux_fit_hfx=ks_aux_fit_hfx, &
693 rho_aux_fit=rho_aux_fit)
694 CALL qs_rho_get(rho_aux_fit, rho_ao_kp=rho_ao)
695 ELSE
696 CALL qs_rho_get(rho, rho_ao_kp=rho_ao)
697 END IF
698
699 DO ispin = 1, nspins
700 CALL cp_fm_get_info(ace_w(1, ispin), nrow_global=nao, ncol_global=nocc)
701
702 IF (do_admm) THEN
703 CALL dbcsr_set(ks_matrix(ispin, 1)%matrix, 0.0_dp)
704 CALL dbcsr_add(ks_matrix(ispin, 1)%matrix, &
705 matrix_h(1, 1)%matrix, 1.0_dp, 1.0_dp)
706 CALL dbcsr_set(ks_aux_fit(ispin)%matrix, 0.0_dp)
707 CALL cp_dbcsr_plus_fm_fm_t( &
708 sparse_matrix=ks_aux_fit(ispin)%matrix, &
709 matrix_v=ace_w(1, ispin), &
710 ncol=nocc, &
711 alpha=-1.0_dp, &
712 keep_sparsity=.true.)
713 CALL dbcsr_add(ks_aux_fit_hfx(ispin)%matrix, &
714 ks_aux_fit(ispin)%matrix, 0.0_dp, 1.0_dp)
715 ELSE
716 CALL dbcsr_set(ks_matrix(ispin, 1)%matrix, 0.0_dp)
717 CALL cp_dbcsr_plus_fm_fm_t( &
718 sparse_matrix=ks_matrix(ispin, 1)%matrix, &
719 matrix_v=ace_w(1, ispin), &
720 ncol=nocc, &
721 alpha=-1.0_dp, &
722 keep_sparsity=.true.)
723 CALL dbcsr_add(ks_matrix(ispin, 1)%matrix, &
724 matrix_h(1, 1)%matrix, 1.0_dp, 1.0_dp)
725 END IF
726 END DO
727
728 ! Exchange energy: E_x = -0.5 * sum_spin Tr[ W^T * P * W ]
729 ehfx_ace = 0.0_dp
730 DO ispin = 1, nspins
731 CALL cp_fm_get_info(ace_w(1, ispin), nrow_global=nao, ncol_global=nocc)
732
733 NULLIFY (fmstruct)
734 CALL cp_fm_struct_create(fmstruct, context=blacs_env, para_env=para_env, &
735 nrow_global=nao, ncol_global=nao)
736 CALL cp_fm_create(p_fm, fmstruct, name="P_dense")
737 CALL cp_fm_struct_release(fmstruct)
738 CALL copy_dbcsr_to_fm(rho_ao(ispin, 1)%matrix, p_fm)
739
740 CALL cp_fm_create(pw_fm, ace_w(1, ispin)%matrix_struct, name="PW")
741 CALL parallel_gemm('N', 'N', nao, nocc, nao, &
742 1.0_dp, p_fm, ace_w(1, ispin), 0.0_dp, pw_fm)
743 CALL cp_fm_trace(ace_w(1, ispin), pw_fm, trace_val)
744 ehfx_ace = ehfx_ace - 0.5_dp*trace_val
745
746 CALL cp_fm_release(p_fm)
747 CALL cp_fm_release(pw_fm)
748
749 IF (dbg_energy .AND. iw > 0) THEN
750 WRITE (iw, '(T4,A,I4,A,F20.10)') &
751 'ispin=', ispin, ' E_x(ACE) += ', -0.5_dp*trace_val
752 END IF
753 END DO
754
755 energy%ex = ehfx_ace
756
757 ! DIAG C: BUILD-step consistency check (printed when ehfx_full_ref present)
758 IF (dbg_energy .AND. iw > 0) THEN
759 WRITE (iw, '(T2,A,F20.10)') 'ACE | E_x(ACE) = ', ehfx_ace
760 IF (PRESENT(ehfx_full_ref)) THEN
761 WRITE (iw, '(T2,A,F20.10)') 'ACE | E_x(full) = ', ehfx_full_ref
762 WRITE (iw, '(T2,A,ES12.4)') 'ACE | |delta| = ', abs(ehfx_ace - ehfx_full_ref)
763 WRITE (iw, '(T2,A)') '(|delta| should be ~0 on BUILD steps; small is good)'
764 END IF
765 END IF
766
767 ! ----------------------------------------------------------------
768 ! DIAG A: projector staleness check.
769 !
770 ! Computes ||W^T C_occ^current||_F (summed over spins) and divides
771 ! by ace_W_ref_norm = ||W^T C_occ^BUILD||_F stored at BUILD time.
772 !
773 ! staleness_ratio:
774 ! 1.0 → C_occ hasn't changed since BUILD; projector is fresh
775 ! < 1 → C_occ has rotated; how much depends on the SCF dynamics
776 ! → 0 → C_occ is orthogonal to the BUILD-time span; W is useless
777 !
778 ! For non-ADMM: C_occ comes from primary mos (nao_orb x nocc).
779 ! For ADMM: C_occ comes from mos_aux_fit (nao_aux x nocc_aux),
780 ! consistent with ace_W dimensions.
781 ! ----------------------------------------------------------------
782 IF (dbg_stale .AND. ace_w_ref_norm > 0.0_dp) THEN
783 stale_norm = 0.0_dp
784 CALL get_qs_env(qs_env, mos=mos_diag)
785
786 DO ispin = 1, nspins
787 CALL cp_fm_get_info(ace_w(1, ispin), nrow_global=nao_d, ncol_global=nocc_d)
788
789 IF (do_admm) THEN
790 CALL get_admm_env(admm_env, mos_aux_fit=mos_aux_diag)
791 IF (mos_aux_diag(ispin)%use_mo_coeff_b) THEN
792 CALL copy_dbcsr_to_fm(mos_aux_diag(ispin)%mo_coeff_b, &
793 mos_aux_diag(ispin)%mo_coeff)
794 END IF
795 CALL get_mo_set(mos_aux_diag(ispin), mo_coeff=mo_coeff_diag, &
796 nao=nao_d, nmo=nmo_d, homo=nocc_d)
797 ELSE
798 IF (mos_diag(ispin)%use_mo_coeff_b) THEN
799 CALL copy_dbcsr_to_fm(mos_diag(ispin)%mo_coeff_b, &
800 mos_diag(ispin)%mo_coeff)
801 END IF
802 CALL get_mo_set(mos_diag(ispin), mo_coeff=mo_coeff_diag, &
803 nao=nao_d, nmo=nmo_d, homo=nocc_d)
804 END IF
805
806 NULLIFY (fmstruct_diag)
807 CALL cp_fm_struct_create(fmstruct_diag, context=blacs_env, &
808 para_env=para_env, &
809 nrow_global=nao_d, ncol_global=nocc_d)
810 CALL cp_fm_create(c_occ_diag, fmstruct_diag, name="C_stale")
811 CALL cp_fm_struct_release(fmstruct_diag)
812 CALL cp_fm_to_fm(mo_coeff_diag, c_occ_diag)
813
814 NULLIFY (fmstruct_diag)
815 CALL cp_fm_struct_create(fmstruct_diag, context=blacs_env, &
816 para_env=para_env, &
817 nrow_global=nocc_d, ncol_global=nocc_d)
818 CALL cp_fm_create(a_diag_fm, fmstruct_diag, name="WtC_stale")
819 CALL cp_fm_struct_release(fmstruct_diag)
820
821 CALL parallel_gemm('T', 'N', nocc_d, nocc_d, nao_d, &
822 1.0_dp, ace_w(1, ispin), c_occ_diag, 0.0_dp, a_diag_fm)
823 CALL cp_fm_trace(a_diag_fm, a_diag_fm, frob_a)
824 stale_norm = stale_norm + sqrt(max(frob_a, 0.0_dp))
825
826 CALL cp_fm_release(c_occ_diag)
827 CALL cp_fm_release(a_diag_fm)
828 END DO
829 IF (iw > 0) THEN
830 WRITE (iw, '(/,T2,A)') repeat('-', 56)
831 WRITE (iw, '(T2,A,I6)') 'ACE DIAG A | ace_step_counter = ', ace_step_counter
832 WRITE (iw, '(T4,A,ES12.4)') '||W^T C_occ^current||_F = ', stale_norm
833 WRITE (iw, '(T4,A,ES12.4)') '||W^T C_occ^BUILD||_F (ref) = ', ace_w_ref_norm
834 WRITE (iw, '(T4,A,F10.6)') 'staleness ratio (1=fresh, 0=stale) = ', &
835 stale_norm/max(ace_w_ref_norm, 1.0e-30_dp)
836 WRITE (iw, '(T2,A)') repeat('-', 56)
837 END IF
838 END IF
839
840 CALL timestop(handle)
841
842 END SUBROUTINE hfx_ace_apply_projector
843
844! **************************************************************************************************
845!> \brief Release all ACE storage and reset state flags.
846!> \param iw_opt ...
847! **************************************************************************************************
848 SUBROUTINE hfx_ace_release(iw_opt)
849
850 INTEGER, INTENT(IN), OPTIONAL :: iw_opt
851
852 INTEGER :: i, iw, j
853
854 iw = -1
855 IF (PRESENT(iw_opt)) iw = iw_opt
856
857 IF (ALLOCATED(ace_w)) THEN
858 DO j = 1, SIZE(ace_w, 2)
859 DO i = 1, SIZE(ace_w, 1)
860 IF (ASSOCIATED(ace_w(i, j)%matrix_struct)) CALL cp_fm_release(ace_w(i, j))
861 END DO
862 END DO
863 DEALLOCATE (ace_w)
864 END IF
865
866 ace_is_built = .false.
867 ace_step_counter = 0
868 ace_w_ref_norm = 0.0_dp
869 ace_geo_step = 0
870 ace_dynamic_mode = .false. ! ADDED: reset dynamic mode on release, so it must be explicitly re-enabled for GEO_OPT/MD runs
871
872 IF (iw > 0) WRITE (iw, '(T2,A)') 'ACE | storage released, counters reset'
873
874 END SUBROUTINE hfx_ace_release
875
876 ! **************************************************************************************************
877 !> \brief Mark this run as dynamic (GEO_OPT/MD) so Bypass C fires for geo step 0.
878 !> Call this once from the geo_opt or MD driver before the first SCF.
879 !> \param is_dynamic .TRUE. for GEO_OPT/MD, .FALSE. to reset.
880 ! **************************************************************************************************
881! **************************************************************************************************
882!> \brief ...
883!> \param is_dynamic ...
884! **************************************************************************************************
885 SUBROUTINE hfx_ace_set_dynamic_mode(is_dynamic)
886 LOGICAL, INTENT(IN) :: is_dynamic
887
888 ace_dynamic_mode = is_dynamic
889 END SUBROUTINE hfx_ace_set_dynamic_mode
890
891 ! **************************************************************************************************
892 !> \brief Private helper: call hfx_ks_matrix with or without ext_xc_section.
893 !> \param qs_env ...
894 !> \param ks_matrix ...
895 !> \param rho ...
896 !> \param energy ...
897 !> \param calculate_forces ...
898 !> \param just_energy ...
899 !> \param v_rspace_new ...
900 !> \param v_tau_rspace ...
901 !> \param ext_xc_section ...
902 ! **************************************************************************************************
903! **************************************************************************************************
904!> \brief ...
905!> \param qs_env ...
906!> \param ks_matrix ...
907!> \param rho ...
908!> \param energy ...
909!> \param calculate_forces ...
910!> \param just_energy ...
911!> \param v_rspace_new ...
912!> \param v_tau_rspace ...
913!> \param ext_xc_section ...
914! **************************************************************************************************
915 SUBROUTINE hfx_call(qs_env, ks_matrix, rho, energy, &
916 calculate_forces, just_energy, &
917 v_rspace_new, v_tau_rspace, ext_xc_section)
918
919 TYPE(qs_environment_type), POINTER :: qs_env
920 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: ks_matrix
921 TYPE(qs_rho_type), POINTER :: rho
922 TYPE(qs_energy_type), POINTER :: energy
923 LOGICAL, INTENT(IN) :: calculate_forces, just_energy
924 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: v_rspace_new, v_tau_rspace
925 TYPE(section_vals_type), OPTIONAL, POINTER :: ext_xc_section
926
927 IF (PRESENT(ext_xc_section)) THEN
928 CALL hfx_ks_matrix(qs_env, ks_matrix, rho, energy, &
929 calculate_forces, just_energy, &
930 v_rspace_new, v_tau_rspace, &
931 ext_xc_section=ext_xc_section)
932 ELSE
933 CALL hfx_ks_matrix(qs_env, ks_matrix, rho, energy, &
934 calculate_forces, just_energy, &
935 v_rspace_new, v_tau_rspace)
936 END IF
937
938 END SUBROUTINE hfx_call
939
940END MODULE hfx_ace_methods
cp_dbcsr_api::dbcsr_create
Definition cp_dbcsr_api.F:194
cp_fm_basic_linalg::cp_fm_trace
Definition cp_fm_basic_linalg.F:78
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:88
cp_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:83
parallel_gemm_api::parallel_gemm
Definition parallel_gemm_api.F:39
admm_types
Types and set/get functions for auxiliary density matrix methods.
Definition admm_types.F:15
admm_types::get_admm_env
subroutine, public get_admm_env(admm_env, mo_derivs_aux_fit, mos_aux_fit, sab_aux_fit, sab_aux_fit_asymm, sab_aux_fit_vs_orb, matrix_s_aux_fit, matrix_s_aux_fit_kp, matrix_s_aux_fit_vs_orb, matrix_s_aux_fit_vs_orb_kp, task_list_aux_fit, matrix_ks_aux_fit, matrix_ks_aux_fit_kp, matrix_ks_aux_fit_im, matrix_ks_aux_fit_dft, matrix_ks_aux_fit_hfx, matrix_ks_aux_fit_dft_kp, matrix_ks_aux_fit_hfx_kp, rho_aux_fit, rho_aux_fit_buffer, admm_dm)
Get routine for the ADMM env.
Definition admm_types.F:599
bibliography
collects all references to literature in CP2K as new algorithms / method are included from literature...
Definition bibliography.F:21
bibliography::lin2016ace
integer, save, public lin2016ace
Definition bibliography.F:36
cp_blacs_env
methods related to the blacs parallel environment
Definition cp_blacs_env.F:15
cp_control_types
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
Definition cp_control_types.F:12
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_release
subroutine, public dbcsr_release(matrix)
...
Definition cp_dbcsr_api.F:1142
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_operations
DBCSR operations in CP2K.
Definition cp_dbcsr_operations.F:18
cp_dbcsr_operations::copy_dbcsr_to_fm
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
Definition cp_dbcsr_operations.F:212
cp_dbcsr_operations::cp_dbcsr_plus_fm_fm_t
subroutine, public cp_dbcsr_plus_fm_fm_t(sparse_matrix, matrix_v, matrix_g, ncol, alpha, keep_sparsity, symmetry_mode)
performs the multiplication sparse_matrix+dense_mat*dens_mat^T if matrix_g is not explicitly given,...
Definition cp_dbcsr_operations.F:693
cp_fm_basic_linalg
Basic linear algebra operations for full matrices.
Definition cp_fm_basic_linalg.F:14
cp_fm_basic_linalg::cp_fm_scale
subroutine, public cp_fm_scale(alpha, matrix_a)
scales a matrix matrix_a = alpha * matrix_b
Definition cp_fm_basic_linalg.F:1502
cp_fm_basic_linalg::cp_fm_triangular_multiply
subroutine, public cp_fm_triangular_multiply(triangular_matrix, matrix_b, side, transpose_tr, invert_tr, uplo_tr, unit_diag_tr, n_rows, n_cols, alpha)
multiplies in place by a triangular matrix: matrix_b = alpha op(triangular_matrix) matrix_b or (if si...
Definition cp_fm_basic_linalg.F:1413
cp_fm_cholesky
various cholesky decomposition related routines
Definition cp_fm_cholesky.F:14
cp_fm_cholesky::cp_fm_cholesky_decompose
subroutine, public cp_fm_cholesky_decompose(matrix, n, info_out)
used to replace a symmetric positive def. matrix M with its cholesky decomposition U: M = U^T * U,...
Definition cp_fm_cholesky.F:61
cp_fm_struct
represent the structure of a full matrix
Definition cp_fm_struct.F:14
cp_fm_struct::cp_fm_struct_create
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
Definition cp_fm_struct.F:132
cp_fm_struct::cp_fm_struct_release
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition cp_fm_struct.F:355
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_get_info
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
Definition cp_fm_types.F:1038
cp_fm_types::cp_fm_create
subroutine, public cp_fm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
creates a new full matrix with the given structure
Definition cp_fm_types.F:165
cp_log_handling
various routines to log and control the output. The idea is that decisions about where to log should ...
Definition cp_log_handling.F:41
cp_log_handling::cp_logger_get_default_io_unit
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
Definition cp_log_handling.F:523
cp_log_handling::cp_get_default_logger
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
Definition cp_log_handling.F:234
hfx_ace_methods
Adaptively Compressed Exchange (ACE) operator for HFX. Reference: Lin, J. Chem. Theory Comput....
Definition hfx_ace_methods.F:43
hfx_ace_methods::hfx_ace_set_dynamic_mode
subroutine, public hfx_ace_set_dynamic_mode(is_dynamic)
Mark this run as dynamic (GEO_OPT/MD) so Bypass C fires for geo step 0. Call this once from the geo_o...
Definition hfx_ace_methods.F:886
hfx_ace_methods::hfx_ace_release
subroutine, public hfx_ace_release(iw_opt)
Release all ACE storage and reset state flags.
Definition hfx_ace_methods.F:849
hfx_ace_methods::hfx_ace_ks_matrix
subroutine, public hfx_ace_ks_matrix(qs_env, ks_matrix, rho, energy, calculate_forces, just_energy, v_rspace_new, v_tau_rspace, ace_rebuild_frequency, ext_xc_section)
Main ACE entry point, replacing hfx_ks_matrix in qs_ks_methods.
Definition hfx_ace_methods.F:151
hfx_admm_utils
Utilities for hfx and admm methods.
Definition hfx_admm_utils.F:17
hfx_admm_utils::hfx_ks_matrix
subroutine, public hfx_ks_matrix(qs_env, matrix_ks, rho, energy, calculate_forces, just_energy, v_rspace_new, v_tau_rspace, ext_xc_section)
Add the hfx contributions to the Hamiltonian.
Definition hfx_admm_utils.F:1146
hfx_types
Types and set/get functions for HFX.
Definition hfx_types.F:16
input_section_types
objects that represent the structure of input sections and the data contained in an input section
Definition input_section_types.F:15
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:23
parallel_gemm_api
basic linear algebra operations for full matrixes
Definition parallel_gemm_api.F:14
pw_types
Definition pw_types.F:24
qs_energy_types
Definition qs_energy_types.F:14
qs_environment_types
Definition qs_environment_types.F:14
qs_environment_types::get_qs_env
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, mimic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_vhxc, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, xcint_weights, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
Definition qs_environment_types.F:534
qs_mo_types
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
qs_mo_types::get_mo_set
subroutine, public get_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, mo_coeff, mo_coeff_b, uniform_occupation, kts, mu, flexible_electron_count)
Get the components of a MO set data structure.
Definition qs_mo_types.F:401
qs_rho_types
superstucture that hold various representations of the density and keeps track of which ones are vali...
Definition qs_rho_types.F:18
qs_rho_types::qs_rho_get
subroutine, public qs_rho_get(rho_struct, rho_ao, rho_ao_im, rho_ao_kp, rho_ao_im_kp, rho_r, drho_r, rho_g, drho_g, tau_r, tau_g, rho_r_valid, drho_r_valid, rho_g_valid, drho_g_valid, tau_r_valid, tau_g_valid, tot_rho_r, tot_rho_g, rho_r_sccs, soft_valid, complex_rho_ao)
returns info about the density described by this object. If some representation is not available an e...
Definition qs_rho_types.F:229
scf_control_types
parameters that control an scf iteration
Definition scf_control_types.F:17
admm_types::admm_type
stores some data used in wavefunction fitting
Definition admm_types.F:120
cp_blacs_env::cp_blacs_env_type
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
Definition cp_blacs_env.F:53
cp_control_types::dft_control_type
Definition cp_control_types.F:757
cp_dbcsr_api::dbcsr_p_type
Definition cp_dbcsr_api.F:172
cp_dbcsr_api::dbcsr_type
Definition cp_dbcsr_api.F:176
cp_fm_struct::cp_fm_struct_type
keeps the information about the structure of a full matrix
Definition cp_fm_struct.F:82
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:114
cp_log_handling::cp_logger_type
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Definition cp_log_handling.F:140
hfx_types::hfx_type
stores some data used in construction of Kohn-Sham matrix
Definition hfx_types.F:514
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:127
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:746
pw_types::pw_r3d_rs_type
Definition pw_types.F:62
qs_energy_types::qs_energy_type
Definition qs_energy_types.F:25
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:222
qs_mo_types::mo_set_type
Definition qs_mo_types.F:47
qs_rho_types::qs_rho_type
keeps the density in various representations, keeping track of which ones are valid.
Definition qs_rho_types.F:62
scf_control_types::scf_control_type
Definition scf_control_types.F:126