(git:15c1bfc)
Loading...
Searching...
No Matches
qs_core_hamiltonian.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Calculation of the core Hamiltonian integral matrix <a|H|b> over
10!> Cartesian Gaussian-type functions.
11!>
12!> <a|H|b> = <a|T|b> + <a|V|b>
13!>
14!> Kinetic energy:
15!>
16!> <a|T|b> = <a|-nabla**2/2|b>
17!> \_______________/
18!> |
19!> kinetic
20!>
21!> Nuclear potential energy:
22!>
23!> a) Allelectron calculation:
24!>
25!> erfc(r)
26!> <a|V|b> = -Z*<a|---------|b>
27!> r
28!>
29!> 1 - erf(r)
30!> = -Z*<a|------------|b>
31!> r
32!>
33!> 1 erf(r)
34!> = -Z*(<a|---|b> - <a|--------|b>)
35!> r r
36!>
37!> 1
38!> = -Z*(<a|---|b> - N*<ab||c>)
39!> r
40!>
41!> -Z
42!> = <a|---|b> + Z*N*<ab||c>
43!> r
44!> \_______/ \_____/
45!> | |
46!> nuclear coulomb
47!>
48!> b) Pseudopotential calculation (Goedecker, Teter and Hutter; GTH):
49!>
50!> <a|V|b> = <a|(V(local) + V(non-local))|b>
51!>
52!> = <a|(V(local)|b> + <a|V(non-local))|b>
53!>
54!> <a|V(local)|b> = <a|-Z(eff)*erf(SQRT(2)*alpha*r)/r +
55!> (C1 + C2*(alpha*r)**2 + C3*(alpha*r)**4 +
56!> C4*(alpha*r)**6)*exp(-(alpha*r)**2/2))|b>
57!>
58!> <a|V(non-local)|b> = <a|p(l,i)>*h(i,j)*<p(l,j)|b>
59!> \par Literature
60!> S. Goedecker, M. Teter and J. Hutter, Phys. Rev. B 54, 1703 (1996)
61!> C. Hartwigsen, S. Goedecker and J. Hutter, Phys. Rev. B 58, 3641 (1998)
62!> M. Krack and M. Parrinello, Phys. Chem. Chem. Phys. 2, 2105 (2000)
63!> S. Obara and A. Saika, J. Chem. Phys. 84, 3963 (1986)
64!> \par History
65!> - Joost VandeVondele (April 2003) : added LSD forces
66!> - Non-redundant calculation of the non-local part of the GTH PP
67!> (22.05.2003,MK)
68!> - New parallelization scheme (27.06.2003,MK)
69!> - OpenMP version (07.12.2003,JGH)
70!> - Binary search loop for VPPNL operators (09.01.2004,JGH,MK)
71!> - Refactoring of pseudopotential and nuclear attraction integrals (25.02.2009,JGH)
72!> - General refactoring (01.10.2010,JGH)
73!> - Refactoring related to the new kinetic energy and overlap routines (07.2014,JGH)
74!> - k-point functionality (07.2015,JGH)
75!> \author Matthias Krack (14.09.2000,21.03.02)
76! **************************************************************************************************
81 USE cp_dbcsr_api, ONLY: dbcsr_add,&
86 dbcsr_set,&
88 dbcsr_type_antisymmetric
96 USE cp_output_handling, ONLY: cp_p_file,&
105 USE kinds, ONLY: default_string_length,&
106 dp
110 USE qs_condnum, ONLY: overlap_condnum
117 USE qs_kind_types, ONLY: qs_kind_type
118 USE qs_ks_types, ONLY: get_ks_env,&
123 USE qs_oce_types, ONLY: allocate_oce_set,&
127 USE qs_rho_types, ONLY: qs_rho_get,&
129 USE virial_types, ONLY: virial_type
130#include "./base/base_uses.f90"
131
132 IMPLICIT NONE
133
134 PRIVATE
135
136 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_core_hamiltonian'
137
140
141CONTAINS
142
143! **************************************************************************************************
144!> \brief Cosntruction of the QS Core Hamiltonian Matrix
145!> \param qs_env ...
146!> \param calculate_forces ...
147!> \author Creation (11.03.2002,MK)
148!> Non-redundant calculation of the non-local part of the GTH PP (22.05.2003,MK)
149!> New parallelization scheme (27.06.2003,MK)
150! **************************************************************************************************
151 SUBROUTINE build_core_hamiltonian_matrix(qs_env, calculate_forces)
152
153 TYPE(qs_environment_type), POINTER :: qs_env
154 LOGICAL, INTENT(IN) :: calculate_forces
155
156 CHARACTER(LEN=*), PARAMETER :: routinen = 'build_core_hamiltonian_matrix'
157
158 INTEGER :: handle, img, iw, nder, nders, nimages, &
159 nkind
160 LOGICAL :: h_is_complex, norml1, norml2, ofdft, &
161 use_arnoldi, use_virial
162 REAL(kind=dp) :: eps_filter, eps_fit
163 REAL(kind=dp), DIMENSION(2) :: condnum
164 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
165 TYPE(cp_blacs_env_type), POINTER :: blacs_env
166 TYPE(cp_logger_type), POINTER :: logger
167 TYPE(dbcsr_distribution_type), POINTER :: dbcsr_dist
168 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_h, matrix_p, matrix_s, matrix_t, &
169 matrix_w
170 TYPE(dft_control_type), POINTER :: dft_control
171 TYPE(kg_environment_type), POINTER :: kg_env
172 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
173 POINTER :: sab_orb, sap_oce
174 TYPE(oce_matrix_type), POINTER :: oce
175 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
176 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
177 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
178 TYPE(qs_ks_env_type), POINTER :: ks_env
179 TYPE(qs_rho_type), POINTER :: rho
180 TYPE(virial_type), POINTER :: virial
181
182 IF (calculate_forces) THEN
183 CALL timeset(routinen//"_forces", handle)
184 ELSE
185 CALL timeset(routinen, handle)
186 END IF
187
188 NULLIFY (logger)
189 logger => cp_get_default_logger()
190
191 NULLIFY (dft_control)
192 CALL get_qs_env(qs_env=qs_env, dft_control=dft_control)
193
194 ! is this a orbital-free method calculation
195 ofdft = dft_control%qs_control%ofgpw
196
197 nimages = dft_control%nimages
198 IF (ofdft) THEN
199 cpassert(nimages == 1)
200 END IF
201
202 nders = 0
203 IF (calculate_forces) THEN
204 nder = 1
205 ELSE
206 IF (cp_print_key_should_output(logger%iter_info, qs_env%input, &
207 "DFT%PRINT%AO_MATRICES/DERIVATIVES") /= 0) THEN
208 nder = 1
209 ELSE
210 nder = 0
211 END IF
212 END IF
213
214 IF ((cp_print_key_should_output(logger%iter_info, qs_env%input, &
215 "DFT%PRINT%AO_MATRICES/OVERLAP") /= 0 .AND. &
216 btest(cp_print_key_should_output(logger%iter_info, qs_env%input, &
217 "DFT%PRINT%AO_MATRICES/DERIVATIVES"), cp_p_file))) THEN
218 nders = 1
219 END IF
220
221 ! the delta pulse in the periodic case needs the momentum operator,
222 ! which is equivalent to the derivative of the overlap matrix
223 IF (ASSOCIATED(dft_control%rtp_control)) THEN
224 IF (dft_control%rtp_control%apply_delta_pulse .AND. &
225 dft_control%rtp_control%periodic) THEN
226 nders = 1
227 END IF
228 END IF
229
230 IF (dft_control%tddfpt2_control%enabled) THEN
231 nders = 1
232 IF (dft_control%do_admm) THEN
233 IF (dft_control%admm_control%purification_method /= do_admm_purify_none) &
234 CALL cp_abort(__location__, &
235 "Only purification method NONE is possible with TDDFT at the moment")
236 END IF
237 END IF
238
239 ! filter for new matrices
240 eps_filter = dft_control%qs_control%eps_filter_matrix
241 !
242 NULLIFY (ks_env)
243 CALL get_qs_env(qs_env=qs_env, ks_env=ks_env)
244 NULLIFY (matrix_s, matrix_t)
245 CALL get_qs_env(qs_env=qs_env, kinetic_kp=matrix_t, matrix_s_kp=matrix_s)
246 NULLIFY (sab_orb)
247 CALL get_qs_env(qs_env=qs_env, sab_orb=sab_orb)
248 NULLIFY (rho, force, matrix_p, matrix_w)
249 IF (calculate_forces) THEN
250 CALL get_qs_env(qs_env=qs_env, force=force, matrix_w_kp=matrix_w)
251 CALL get_qs_env(qs_env=qs_env, rho=rho)
252 CALL qs_rho_get(rho, rho_ao_kp=matrix_p)
253 ! *** If LSD, then combine alpha density and beta density to
254 ! *** total density: alpha <- alpha + beta and
255 ! *** spin density: beta <- alpha - beta
256 ! (since all things can be computed based on the sum of these matrices anyway)
257 ! (matrix_p is restored at the end of the run, matrix_w is left in its modified state
258 ! (as it should not be needed afterwards)
259 IF (SIZE(matrix_p, 1) == 2) THEN
260 DO img = 1, nimages
261 CALL dbcsr_add(matrix_p(1, img)%matrix, matrix_p(2, img)%matrix, &
262 alpha_scalar=1.0_dp, beta_scalar=1.0_dp)
263 CALL dbcsr_add(matrix_p(2, img)%matrix, matrix_p(1, img)%matrix, &
264 alpha_scalar=-2.0_dp, beta_scalar=1.0_dp)
265 CALL dbcsr_add(matrix_w(1, img)%matrix, matrix_w(2, img)%matrix, &
266 alpha_scalar=1.0_dp, beta_scalar=1.0_dp)
267 END DO
268 END IF
269 ELSE
270 NULLIFY (matrix_p, matrix_w)
271 END IF
272
273 ! S matrix
274 CALL build_overlap_matrix(ks_env, nderivative=nders, matrixkp_s=matrix_s, &
275 matrix_name="OVERLAP MATRIX", &
276 basis_type_a="ORB", &
277 basis_type_b="ORB", &
278 sab_nl=sab_orb, calculate_forces=calculate_forces, &
279 matrixkp_p=matrix_w)
280
281 IF (calculate_forces) THEN
282 ! *** If LSD, then recover alpha density and beta density ***
283 ! *** from the total density (1) and the spin density (2) ***
284 ! *** The W matrix is neglected, since it will be destroyed ***
285 ! *** in the calling force routine after leaving this routine ***
286 IF (SIZE(matrix_p, 1) == 2) THEN
287 DO img = 1, nimages
288 CALL dbcsr_add(matrix_p(1, img)%matrix, matrix_p(2, img)%matrix, &
289 alpha_scalar=0.5_dp, beta_scalar=0.5_dp)
290 CALL dbcsr_add(matrix_p(2, img)%matrix, matrix_p(1, img)%matrix, &
291 alpha_scalar=-1.0_dp, beta_scalar=1.0_dp)
292 END DO
293 END IF
294 END IF
295
296 ! T matrix
297 CALL kinetic_energy_matrix(qs_env, matrixkp_t=matrix_t, &
298 matrix_p=matrix_p, &
299 matrix_name="KINETIC ENERGY MATRIX", &
300 basis_type="ORB", &
301 sab_orb=sab_orb, &
302 calculate_forces=calculate_forces, &
303 eps_filter=eps_filter)
304
305 ! (Re-)allocate H matrix based on overlap matrix
306 CALL get_ks_env(ks_env, complex_ks=h_is_complex)
307 CALL qs_matrix_h_allocate(qs_env, matrix_s(1, 1)%matrix, is_complex=h_is_complex)
308
309 NULLIFY (matrix_h)
310 CALL get_qs_env(qs_env, matrix_h_kp=matrix_h)
311
312 IF (.NOT. ofdft) THEN
313 DO img = 1, nimages
314 CALL dbcsr_copy(matrix_h(1, img)%matrix, matrix_t(1, img)%matrix, &
315 keep_sparsity=.true., name="CORE HAMILTONIAN MATRIX")
316 END DO
317 END IF
318
319 NULLIFY (qs_kind_set, atomic_kind_set, particle_set)
320 CALL get_qs_env(qs_env=qs_env, qs_kind_set=qs_kind_set, atomic_kind_set=atomic_kind_set, &
321 particle_set=particle_set)
322
323 ! *** core and pseudopotentials
324 CALL core_matrices(qs_env, matrix_h, matrix_p, calculate_forces, nder)
325
326 ! *** CNEO nuclear V_core
327 CALL cneo_core_matrices(qs_env, calculate_forces, nder)
328
329 ! *** GAPW one-center-expansion (oce) matrices
330 NULLIFY (sap_oce)
331 CALL get_qs_env(qs_env=qs_env, sap_oce=sap_oce)
332 NULLIFY (oce)
333 IF (dft_control%qs_control%gapw .OR. dft_control%qs_control%gapw_xc) THEN
334 CALL get_qs_env(qs_env=qs_env, oce=oce)
335 CALL create_oce_set(oce)
336 nkind = SIZE(atomic_kind_set)
337 CALL allocate_oce_set(oce, nkind)
338 eps_fit = dft_control%qs_control%gapw_control%eps_fit
339 IF (ASSOCIATED(sap_oce)) &
340 CALL build_oce_matrices(oce%intac, calculate_forces, nder, qs_kind_set, particle_set, &
341 sap_oce, eps_fit)
342 END IF
343
344 ! *** KG atomic potentials for nonadditive kinetic energy
345 IF (dft_control%qs_control%do_kg) THEN
346 IF (qs_env%kg_env%tnadd_method == kg_tnadd_atomic) THEN
347 CALL get_qs_env(qs_env=qs_env, kg_env=kg_env, virial=virial, dbcsr_dist=dbcsr_dist)
348 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
349 CALL build_tnadd_mat(kg_env, matrix_p, force, virial, calculate_forces, use_virial, &
350 qs_kind_set, atomic_kind_set, particle_set, sab_orb, dbcsr_dist)
351 END IF
352 END IF
353
354 ! *** Put the core Hamiltonian matrix in the QS environment ***
355 CALL set_qs_env(qs_env, oce=oce)
356 CALL set_ks_env(ks_env, matrix_s_kp=matrix_s, kinetic_kp=matrix_t, matrix_h_kp=matrix_h)
357
358 ! *** Print matrices if requested
359 CALL dump_info_core_hamiltonian(qs_env, calculate_forces)
360
361 ! *** Overlap condition number
362 IF (.NOT. calculate_forces) THEN
363 IF (cp_print_key_should_output(logger%iter_info, qs_env%input, &
364 "DFT%PRINT%OVERLAP_CONDITION") .NE. 0) THEN
365 iw = cp_print_key_unit_nr(logger, qs_env%input, "DFT%PRINT%OVERLAP_CONDITION", &
366 extension=".Log")
367 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%OVERLAP_CONDITION%1-NORM", l_val=norml1)
368 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%OVERLAP_CONDITION%DIAGONALIZATION", l_val=norml2)
369 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%OVERLAP_CONDITION%ARNOLDI", l_val=use_arnoldi)
370 CALL get_qs_env(qs_env=qs_env, blacs_env=blacs_env)
371 CALL overlap_condnum(matrix_s, condnum, iw, norml1, norml2, use_arnoldi, blacs_env)
372 END IF
373 END IF
374
375 CALL timestop(handle)
376
377 END SUBROUTINE build_core_hamiltonian_matrix
378
379! **************************************************************************************************
380!> \brief Possibly prints matrices after the construction of the Core
381!> Hamiltonian Matrix
382!> \param qs_env ...
383!> \param calculate_forces ...
384! **************************************************************************************************
385 SUBROUTINE dump_info_core_hamiltonian(qs_env, calculate_forces)
386 TYPE(qs_environment_type), POINTER :: qs_env
387 LOGICAL, INTENT(IN) :: calculate_forces
388
389 CHARACTER(LEN=*), PARAMETER :: routinen = 'dump_info_core_hamiltonian'
390
391 INTEGER :: after, handle, i, ic, iw, output_unit
392 LOGICAL :: omit_headers
393 TYPE(cp_logger_type), POINTER :: logger
394 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_v
395 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrixkp_h, matrixkp_s, matrixkp_t
396 TYPE(mp_para_env_type), POINTER :: para_env
397
398 CALL timeset(routinen, handle)
399
400 NULLIFY (logger, matrix_v, para_env)
401 logger => cp_get_default_logger()
402 CALL get_qs_env(qs_env, para_env=para_env)
403
404 ! Print the distribution of the overlap matrix blocks
405 ! this duplicates causes duplicate printing at the force calc
406 IF (.NOT. calculate_forces) THEN
407 IF (btest(cp_print_key_should_output(logger%iter_info, &
408 qs_env%input, "PRINT%DISTRIBUTION"), cp_p_file)) THEN
409 output_unit = cp_print_key_unit_nr(logger, qs_env%input, "PRINT%DISTRIBUTION", &
410 extension=".distribution")
411 CALL get_qs_env(qs_env, matrix_s_kp=matrixkp_s)
412 CALL cp_dbcsr_write_matrix_dist(matrixkp_s(1, 1)%matrix, output_unit, para_env)
413 CALL cp_print_key_finished_output(output_unit, logger, qs_env%input, "PRINT%DISTRIBUTION")
414 END IF
415 END IF
416
417 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%AO_MATRICES%OMIT_HEADERS", l_val=omit_headers)
418 ! Print the overlap integral matrix, if requested
419 IF (btest(cp_print_key_should_output(logger%iter_info, &
420 qs_env%input, "DFT%PRINT%AO_MATRICES/OVERLAP"), cp_p_file)) THEN
421 iw = cp_print_key_unit_nr(logger, qs_env%input, "DFT%PRINT%AO_MATRICES/OVERLAP", &
422 extension=".Log")
423 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%AO_MATRICES%NDIGITS", i_val=after)
424 after = min(max(after, 1), 16)
425 CALL get_qs_env(qs_env, matrix_s_kp=matrixkp_s)
426 IF (ASSOCIATED(matrixkp_s)) THEN
427 DO ic = 1, SIZE(matrixkp_s, 2)
428 CALL cp_dbcsr_write_sparse_matrix(matrixkp_s(1, ic)%matrix, 4, after, qs_env, para_env, &
429 output_unit=iw, omit_headers=omit_headers)
430 END DO
431 IF (btest(cp_print_key_should_output(logger%iter_info, qs_env%input, &
432 "DFT%PRINT%AO_MATRICES/DERIVATIVES"), cp_p_file)) THEN
433 DO ic = 1, SIZE(matrixkp_s, 2)
434 DO i = 2, SIZE(matrixkp_s, 1)
435 CALL cp_dbcsr_write_sparse_matrix(matrixkp_s(i, ic)%matrix, 4, after, qs_env, para_env, &
436 output_unit=iw, omit_headers=omit_headers)
437 END DO
438 END DO
439 END IF
440 END IF
441 CALL cp_print_key_finished_output(iw, logger, qs_env%input, &
442 "DFT%PRINT%AO_MATRICES/OVERLAP")
443 END IF
444
445 ! Print the kinetic energy integral matrix, if requested
446 IF (btest(cp_print_key_should_output(logger%iter_info, &
447 qs_env%input, "DFT%PRINT%AO_MATRICES/KINETIC_ENERGY"), cp_p_file)) THEN
448 iw = cp_print_key_unit_nr(logger, qs_env%input, "DFT%PRINT%AO_MATRICES/KINETIC_ENERGY", &
449 extension=".Log")
450 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%AO_MATRICES%NDIGITS", i_val=after)
451 after = min(max(after, 1), 16)
452 CALL get_qs_env(qs_env, kinetic_kp=matrixkp_t)
453 IF (ASSOCIATED(matrixkp_t)) THEN
454 DO ic = 1, SIZE(matrixkp_t, 2)
455 CALL cp_dbcsr_write_sparse_matrix(matrixkp_t(1, ic)%matrix, 4, after, qs_env, para_env, &
456 output_unit=iw, omit_headers=omit_headers)
457 END DO
458 END IF
459 CALL cp_print_key_finished_output(iw, logger, qs_env%input, &
460 "DFT%PRINT%AO_MATRICES/KINETIC_ENERGY")
461 END IF
462
463 ! Print the potential energy matrix, if requested
464 IF (btest(cp_print_key_should_output(logger%iter_info, &
465 qs_env%input, "DFT%PRINT%AO_MATRICES/POTENTIAL_ENERGY"), cp_p_file)) THEN
466 iw = cp_print_key_unit_nr(logger, qs_env%input, "DFT%PRINT%AO_MATRICES/POTENTIAL_ENERGY", &
467 extension=".Log")
468 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%AO_MATRICES%NDIGITS", i_val=after)
469 after = min(max(after, 1), 16)
470 CALL get_qs_env(qs_env, matrix_h_kp=matrixkp_h, kinetic_kp=matrixkp_t)
471 IF (ASSOCIATED(matrixkp_h)) THEN
472 IF (SIZE(matrixkp_h, 2) == 1) THEN
473 CALL dbcsr_allocate_matrix_set(matrix_v, 1)
474 ALLOCATE (matrix_v(1)%matrix)
475 CALL dbcsr_copy(matrix_v(1)%matrix, matrixkp_h(1, 1)%matrix, name="POTENTIAL ENERGY MATRIX")
476 CALL dbcsr_add(matrix_v(1)%matrix, matrixkp_t(1, 1)%matrix, &
477 alpha_scalar=1.0_dp, beta_scalar=-1.0_dp)
478 CALL cp_dbcsr_write_sparse_matrix(matrix_v(1)%matrix, 4, after, qs_env, &
479 para_env, output_unit=iw, omit_headers=omit_headers)
480 CALL dbcsr_deallocate_matrix_set(matrix_v)
481 ELSE
482 cpwarn("Printing of potential energy matrix not implemented for k-points")
483 END IF
484 END IF
485 CALL cp_print_key_finished_output(iw, logger, qs_env%input, &
486 "DFT%PRINT%AO_MATRICES/POTENTIAL_ENERGY")
487 END IF
488
489 ! Print the core Hamiltonian matrix, if requested
490 IF (btest(cp_print_key_should_output(logger%iter_info, &
491 qs_env%input, "DFT%PRINT%AO_MATRICES/CORE_HAMILTONIAN"), cp_p_file)) THEN
492 iw = cp_print_key_unit_nr(logger, qs_env%input, "DFT%PRINT%AO_MATRICES/CORE_HAMILTONIAN", &
493 extension=".Log")
494 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%AO_MATRICES%NDIGITS", i_val=after)
495 after = min(max(after, 1), 16)
496 CALL get_qs_env(qs_env, matrix_h_kp=matrixkp_h)
497 IF (ASSOCIATED(matrixkp_h)) THEN
498 DO ic = 1, SIZE(matrixkp_h, 2)
499 CALL cp_dbcsr_write_sparse_matrix(matrixkp_h(1, ic)%matrix, 4, after, qs_env, para_env, &
500 output_unit=iw, omit_headers=omit_headers)
501 END DO
502 END IF
503 CALL cp_print_key_finished_output(iw, logger, qs_env%input, &
504 "DFT%PRINT%AO_MATRICES/CORE_HAMILTONIAN")
505 END IF
506
507 CALL timestop(handle)
508
509 END SUBROUTINE dump_info_core_hamiltonian
510
511! **************************************************************************************************
512!> \brief (Re-)allocate matrix_h based on the template (typically the overlap matrix)
513!> \param qs_env ...
514!> \param template ...
515!> \param is_complex ...
516! **************************************************************************************************
517 SUBROUTINE qs_matrix_h_allocate(qs_env, template, is_complex)
518 TYPE(qs_environment_type) :: qs_env
519 TYPE(dbcsr_type), INTENT(in) :: template
520 LOGICAL, INTENT(in) :: is_complex
521
522 CHARACTER(LEN=default_string_length) :: headline
523 INTEGER :: img, nimages
524 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_h, matrix_h_im
525 TYPE(dft_control_type), POINTER :: dft_control
526 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
527 POINTER :: sab_orb
528 TYPE(qs_ks_env_type), POINTER :: ks_env
529
530 NULLIFY (matrix_h, matrix_h_im, sab_orb, dft_control, ks_env)
531 CALL get_qs_env(qs_env=qs_env, &
532 matrix_h_kp=matrix_h, &
533 matrix_h_im_kp=matrix_h_im, &
534 sab_orb=sab_orb, &
535 dft_control=dft_control, &
536 ks_env=ks_env)
537
538 nimages = dft_control%nimages
539 CALL dbcsr_allocate_matrix_set(matrix_h, 1, nimages)
540 headline = "CORE HAMILTONIAN MATRIX"
541 DO img = 1, nimages
542 ALLOCATE (matrix_h(1, img)%matrix)
543 CALL dbcsr_create(matrix_h(1, img)%matrix, name=trim(headline), template=template)
544 CALL cp_dbcsr_alloc_block_from_nbl(matrix_h(1, img)%matrix, sab_orb)
545 CALL dbcsr_set(matrix_h(1, img)%matrix, 0.0_dp)
546 END DO
547 CALL set_ks_env(ks_env, matrix_h_kp=matrix_h)
548
549 IF (is_complex) THEN
550 headline = "IMAGINARY PART OF CORE HAMILTONIAN MATRIX"
551 CALL dbcsr_allocate_matrix_set(matrix_h_im, 1, nimages)
552 DO img = 1, nimages
553 ALLOCATE (matrix_h_im(1, img)%matrix)
554 CALL dbcsr_create(matrix_h_im(1, img)%matrix, name=trim(headline), template=template, &
555 matrix_type=dbcsr_type_antisymmetric)
556 CALL cp_dbcsr_alloc_block_from_nbl(matrix_h_im(1, img)%matrix, sab_orb)
557 CALL dbcsr_set(matrix_h_im(1, img)%matrix, 0.0_dp)
558 END DO
559 CALL set_ks_env(ks_env, matrix_h_im_kp=matrix_h_im)
560 END IF
561
562 END SUBROUTINE qs_matrix_h_allocate
563
564! **************************************************************************************************
565!> \brief (Re-)allocates matrix_h_im from matrix_h
566!> \param qs_env ...
567! **************************************************************************************************
569 TYPE(qs_environment_type) :: qs_env
570
571 CHARACTER(LEN=default_string_length) :: headline
572 INTEGER :: image, nimages
573 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_h, matrix_h_im
574 TYPE(dbcsr_type), POINTER :: template
575 TYPE(dft_control_type), POINTER :: dft_control
576 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
577 POINTER :: sab_orb
578 TYPE(qs_ks_env_type), POINTER :: ks_env
579
580 NULLIFY (matrix_h_im, matrix_h, dft_control, template, sab_orb, ks_env)
581
582 CALL get_qs_env(qs_env, &
583 matrix_h_im_kp=matrix_h_im, &
584 matrix_h_kp=matrix_h, &
585 dft_control=dft_control, &
586 sab_orb=sab_orb, &
587 ks_env=ks_env)
588
589 nimages = dft_control%nimages
590
591 cpassert(nimages .EQ. SIZE(matrix_h, 2))
592
593 CALL dbcsr_allocate_matrix_set(matrix_h_im, 1, nimages)
594
595 DO image = 1, nimages
596 headline = "IMAGINARY CORE HAMILTONIAN MATRIX"
597 ALLOCATE (matrix_h_im(1, image)%matrix)
598 template => matrix_h(1, image)%matrix ! base on real part, but anti-symmetric
599 CALL dbcsr_create(matrix=matrix_h_im(1, image)%matrix, template=template, &
600 name=trim(headline), matrix_type=dbcsr_type_antisymmetric)
601 CALL cp_dbcsr_alloc_block_from_nbl(matrix_h_im(1, image)%matrix, sab_orb)
602 CALL dbcsr_set(matrix_h_im(1, image)%matrix, 0.0_dp)
603 END DO
604 CALL set_ks_env(ks_env, matrix_h_im_kp=matrix_h_im)
605
607
608END MODULE qs_core_hamiltonian
Define the atomic kind types and their sub types.
methods related to the blacs parallel environment
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
subroutine, public dbcsr_set(matrix, alpha)
...
subroutine, public dbcsr_add(matrix_a, matrix_b, alpha_scalar, beta_scalar)
...
DBCSR operations in CP2K.
DBCSR output in CP2K.
subroutine, public cp_dbcsr_write_matrix_dist(matrix, output_unit, para_env)
Print the distribution of a sparse matrix.
subroutine, public cp_dbcsr_write_sparse_matrix(sparse_matrix, before, after, qs_env, para_env, first_row, last_row, first_col, last_col, scale, output_unit, omit_headers)
...
various routines to log and control the output. The idea is that decisions about where to log should ...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
integer, parameter, public cp_p_file
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public do_admm_purify_none
integer, parameter, public kg_tnadd_atomic
objects that represent the structure of input sections and the data contained in an input section
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Types needed for a Kim-Gordon-like partitioning into molecular subunits.
Calculation of the local potential contribution of the nonadditive kinetic energy <a|V(local)|b> = <a...
subroutine, public build_tnadd_mat(kg_env, matrix_p, force, virial, calculate_forces, use_virial, qs_kind_set, atomic_kind_set, particle_set, sab_orb, dbcsr_dist)
...
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public default_string_length
Definition kinds.F:57
Interface to the message passing library MPI.
Define the data structure for the particle information.
A collection of functions used by CNEO-DFT (see J. Chem. Theory Comput. 2025, 21, 16,...
subroutine, public cneo_core_matrices(qs_env, calculate_forces, nder)
...
Calculation of overlap matrix condition numbers.
Definition qs_condnum.F:13
subroutine, public overlap_condnum(matrixkp_s, condnum, iunit, norml1, norml2, use_arnoldi, blacs_env)
Calculation of the overlap matrix Condition Number.
Definition qs_condnum.F:66
Calculation of the core Hamiltonian integral matrix <a|H|b> over Cartesian Gaussian-type functions.
subroutine, public build_core_hamiltonian_matrix(qs_env, calculate_forces)
Cosntruction of the QS Core Hamiltonian Matrix.
subroutine, public dump_info_core_hamiltonian(qs_env, calculate_forces)
Possibly prints matrices after the construction of the Core Hamiltonian Matrix.
subroutine, public qs_matrix_h_allocate_imag_from_real(qs_env)
(Re-)allocates matrix_h_im from matrix_h
Calculation of the core Hamiltonian integral matrix <a|H|b> over Cartesian Gaussian-type functions.
subroutine, public kinetic_energy_matrix(qs_env, matrixkp_t, matrix_t, matrix_p, matrix_name, calculate_forces, nderivative, sab_orb, eps_filter, basis_type, debug_forces, debug_stress)
Calculate kinetic energy matrix and possible relativistic correction.
subroutine, public core_matrices(qs_env, matrix_h, matrix_p, calculate_forces, nder, ec_env, dcdr_env, ec_env_matrices, basis_type, debug_forces, debug_stress, atcore)
...
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, 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, 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, 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.
subroutine, public set_qs_env(qs_env, super_cell, mos, qmmm, qmmm_periodic, ewald_env, ewald_pw, mpools, rho_external, external_vxc, mask, scf_control, rel_control, qs_charges, ks_env, ks_qmmm_env, wf_history, scf_env, active_space, input, oce, rho_atom_set, rho0_atom_set, rho0_mpole, run_rtp, rtp, rhoz_set, rhoz_tot, ecoul_1c, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, efield, rhoz_cneo_set, linres_control, xas_env, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, ls_scf_env, do_transport, transport_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, harris_env, gcp_env, mp2_env, bs_env, kg_env, force, kpoints, wanniercentres, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Set the QUICKSTEP environment.
Define the quickstep kind type and their sub types.
subroutine, public set_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, vppl, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, kpoints, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, task_list, task_list_soft, subsys, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env)
...
subroutine, public get_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, rho, rho_xc, vppl, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, task_list, task_list_soft, kpoints, do_kpoints, atomic_kind_set, qs_kind_set, cell, cell_ref, use_ref_cell, particle_set, energy, force, local_particles, local_molecules, molecule_kind_set, molecule_set, subsys, cp_subsys, virial, results, atprop, nkind, natom, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env, nelectron_total, nelectron_spin)
...
Define the neighbor list data types and the corresponding functionality.
Routines for the construction of the coefficients for the expansion of the atomic densities rho1_hard...
subroutine, public build_oce_matrices(intac, calculate_forces, nder, qs_kind_set, particle_set, sap_oce, eps_fit)
Set up the sparse matrix for the coefficients of one center expansions This routine uses the same log...
subroutine, public allocate_oce_set(oce_set, nkind)
Allocate and initialize the matrix set of oce coefficients.
subroutine, public create_oce_set(oce_set)
...
Calculation of overlap matrix, its derivatives and forces.
Definition qs_overlap.F:19
subroutine, public build_overlap_matrix(ks_env, matrix_s, matrixkp_s, matrix_name, nderivative, basis_type_a, basis_type_b, sab_nl, calculate_forces, matrix_p, matrixkp_p)
Calculation of the overlap matrix over Cartesian Gaussian functions.
Definition qs_overlap.F:120
superstucture that hold various representations of the density and keeps track of which ones are vali...
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...
Provides all information about an atomic kind.
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Contains all the info needed for KG runs...
stores all the informations relevant to an mpi environment
Provides all information about a quickstep kind.
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
keeps the density in various representations, keeping track of which ones are valid.