(git:5c1df3d)
Loading...
Searching...
No Matches
rt_propagation_output.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 Routine for the real time propagation output.
10!> \author Florian Schiffmann (02.09)
11! **************************************************************************************************
12
13MODULE rt_propagation_output
14 USE atomic_kind_types, ONLY: atomic_kind_type
15 USE cell_types, ONLY: cell_type
16 USE cp_control_types, ONLY: dft_control_type,&
17 rtp_control_type
18 USE cp_dbcsr_api, ONLY: &
19 dbcsr_add, dbcsr_binary_write, dbcsr_convert_offsets_to_sizes, dbcsr_copy, dbcsr_create, &
20 dbcsr_deallocate_matrix, dbcsr_desymmetrize, dbcsr_distribution_type, dbcsr_filter, &
21 dbcsr_get_info, dbcsr_get_occupation, dbcsr_init_p, dbcsr_iterator_blocks_left, &
22 dbcsr_iterator_next_block, dbcsr_iterator_start, dbcsr_iterator_stop, dbcsr_iterator_type, &
23 dbcsr_multiply, dbcsr_p_type, dbcsr_scale, dbcsr_set, dbcsr_type, &
24 dbcsr_type_antisymmetric, dbcsr_type_no_symmetry, dbcsr_type_symmetric
25 USE cp_dbcsr_contrib, ONLY: dbcsr_checksum,&
26 dbcsr_trace
27 USE cp_dbcsr_cp2k_link, ONLY: cp_dbcsr_alloc_block_from_nbl
28 USE cp_dbcsr_operations, ONLY: cp_dbcsr_sm_fm_multiply,&
29 dbcsr_allocate_matrix_set,&
30 dbcsr_deallocate_matrix_set
31 USE cp_fm_basic_linalg, ONLY: cp_fm_scale_and_add
32 USE cp_fm_struct, ONLY: cp_fm_struct_create,&
33 cp_fm_struct_double,&
34 cp_fm_struct_release,&
35 cp_fm_struct_type
36 USE cp_fm_types, ONLY: cp_fm_create,&
37 cp_fm_get_info,&
38 cp_fm_release,&
39 cp_fm_set_all,&
40 cp_fm_type
41 USE cp_log_handling, ONLY: cp_get_default_logger,&
42 cp_logger_get_default_io_unit,&
43 cp_logger_get_default_unit_nr,&
44 cp_logger_type,&
45 cp_to_string
46 USE cp_output_handling, ONLY: cp_iter_string,&
47 cp_p_file,&
48 cp_print_key_finished_output,&
49 cp_print_key_should_output,&
50 cp_print_key_unit_nr,&
51 cp_printkey_is_on
52 USE cp_realspace_grid_cube, ONLY: cp_pw_to_cube
53 USE efield_utils, ONLY: make_field
54 USE greenx_interface, ONLY: greenx_refine_ft
55 USE input_constants, ONLY: ehrenfest,&
56 real_time_propagation
57 USE input_section_types, ONLY: section_get_ivals,&
58 section_vals_get_subs_vals,&
59 section_vals_type
60 USE kahan_sum, ONLY: accurate_sum
61 USE kinds, ONLY: default_path_length,&
62 dp
63 USE machine, ONLY: m_flush
64 USE mathconstants, ONLY: twopi
65 USE message_passing, ONLY: mp_comm_type
66 USE parallel_gemm_api, ONLY: parallel_gemm
67 USE particle_list_types, ONLY: particle_list_type
68 USE particle_methods, ONLY: get_particle_set
69 USE particle_types, ONLY: particle_type
70 USE physcon, ONLY: evolt,&
71 femtoseconds
72 USE pw_env_types, ONLY: pw_env_get,&
73 pw_env_type
74 USE pw_methods, ONLY: pw_zero
75 USE pw_pool_types, ONLY: pw_pool_type
76 USE pw_types, ONLY: pw_c1d_gs_type,&
77 pw_r3d_rs_type
78 USE qs_energy_types, ONLY: qs_energy_type
79 USE qs_environment_types, ONLY: get_qs_env,&
80 qs_environment_type
81 USE qs_kind_types, ONLY: get_qs_kind_set,&
82 qs_kind_type
83 USE qs_linres_current, ONLY: calculate_jrho_resp
84 USE qs_linres_types, ONLY: current_env_type
85 USE qs_mo_io, ONLY: write_rt_mos_to_restart
86 USE qs_moments, ONLY: build_local_moment_matrix
87 USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type
88 USE qs_operators_ao, ONLY: build_lin_mom_matrix
89 USE qs_rho_types, ONLY: qs_rho_get,&
90 qs_rho_type
91 USE qs_scf_post_gpw, ONLY: qs_scf_post_moments,&
92 write_mo_dependent_results,&
93 write_mo_free_results
94 USE qs_scf_post_tb, ONLY: scf_post_calculation_tb
95 USE qs_scf_types, ONLY: qs_scf_env_type
96 USE qs_subsys_types, ONLY: qs_subsys_get,&
97 qs_subsys_type
98 USE rt_projection_mo_utils, ONLY: compute_and_write_proj_mo
99 USE rt_propagation_ft, ONLY: multi_fft
100 USE rt_propagation_types, ONLY: get_rtp,&
101 rt_prop_type
102 USE rt_propagation_utils, ONLY: write_rtp_mo_cubes,&
103 write_rtp_mos_to_output_unit
104 USE rt_propagation_velocity_gauge, ONLY: velocity_gauge_nl_commutator
105#include "../base/base_uses.f90"
106
107 IMPLICIT NONE
108
109 PRIVATE
110
111 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'rt_propagation_output'
112
113 PUBLIC :: rt_prop_output, &
114 rt_convergence, &
115 rt_convergence_density, &
116 report_density_occupation, &
117 print_moments, &
118 calc_local_moment, &
119 print_ft, &
120 print_rt_file
121
122 INTEGER, PARAMETER, PUBLIC :: rt_file_comp_both = 0, &
123 rt_file_comp_real = 1, &
124 rt_file_comp_imag = 2
125
126CONTAINS
127
128! **************************************************************************************************
129!> \brief ...
130!> \param qs_env ...
131!> \param run_type ...
132!> \param delta_iter ...
133!> \param used_time ...
134! **************************************************************************************************
135 SUBROUTINE rt_prop_output(qs_env, run_type, delta_iter, used_time)
136 TYPE(qs_environment_type), POINTER :: qs_env
137 INTEGER, INTENT(in) :: run_type
138 REAL(dp), INTENT(in), OPTIONAL :: delta_iter, used_time
139
140 INTEGER :: i, n_electrons, n_proj, natom, nspin, &
141 output_unit, spin, unit_nr
142 INTEGER, ALLOCATABLE, DIMENSION(:) :: first_sgf, last_sgf
143 INTEGER, DIMENSION(2) :: nelectron_spin
144 INTEGER, DIMENSION(:), POINTER :: row_blk_sizes
145 LOGICAL :: new_file
146 REAL(dp) :: orthonormality, strace, tot_rho_r, trace
147 REAL(kind=dp), DIMENSION(3) :: field, reference_point, vec_pot
148 REAL(kind=dp), DIMENSION(:), POINTER :: qs_tot_rho_r
149 REAL(kind=dp), DIMENSION(:, :), POINTER :: j_int
150 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
151 TYPE(cell_type), POINTER :: cell
152 TYPE(cp_fm_type), DIMENSION(:), POINTER :: mos_new
153 TYPE(cp_logger_type), POINTER :: logger
154 TYPE(dbcsr_distribution_type), POINTER :: dbcsr_dist
155 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_cim, matrix_cre, matrix_s, p_im, &
156 p_re, p_xyz, rho_new
157 TYPE(dbcsr_type), POINTER :: tmp_ao
158 TYPE(dft_control_type), POINTER :: dft_control
159 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
160 POINTER :: sab_all, sab_orb, sap_ppnl
161 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
162 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
163 TYPE(qs_rho_type), POINTER :: rho
164 TYPE(rt_prop_type), POINTER :: rtp
165 TYPE(section_vals_type), POINTER :: dft_section, input, rtp_section
166
167 NULLIFY (logger, dft_control)
168
169 logger => cp_get_default_logger()
170 CALL get_qs_env(qs_env, &
171 rtp=rtp, &
172 matrix_s=matrix_s, &
173 input=input, &
174 rho=rho, &
175 particle_set=particle_set, &
176 atomic_kind_set=atomic_kind_set, &
177 qs_kind_set=qs_kind_set, &
178 dft_control=dft_control, sab_all=sab_all, sab_orb=sab_orb, &
179 sap_ppnl=sap_ppnl, cell=cell, &
180 dbcsr_dist=dbcsr_dist, nelectron_spin=nelectron_spin)
181
182 rtp_section => section_vals_get_subs_vals(input, "DFT%REAL_TIME_PROPAGATION")
183
184 CALL get_qs_kind_set(qs_kind_set, nelectron=n_electrons)
185 n_electrons = n_electrons - dft_control%charge
186
187 CALL qs_rho_get(rho_struct=rho, tot_rho_r=qs_tot_rho_r)
188
189 tot_rho_r = accurate_sum(qs_tot_rho_r)
190
191 output_unit = cp_print_key_unit_nr(logger, rtp_section, "PRINT%PROGRAM_RUN_INFO", &
192 extension=".scfLog")
193
194 IF (output_unit > 0) THEN
195 WRITE (output_unit, fmt="(/,(T3,A,T40,I5))") &
196 "Information at iteration step:", rtp%iter
197 WRITE (unit=output_unit, fmt="((T3,A,T41,2F20.10))") &
198 "Total electronic density (r-space): ", &
199 tot_rho_r, &
200 tot_rho_r + &
201 REAL(n_electrons, dp)
202 WRITE (unit=output_unit, fmt="((T3,A,T59,F22.14))") &
203 "Total energy:", rtp%energy_new
204 IF (run_type == ehrenfest) THEN
205 WRITE (unit=output_unit, fmt="((T3,A,T61,F20.14))") &
206 "Energy difference to previous iteration step:", rtp%energy_new - rtp%energy_old
207 END IF
208 IF (run_type == real_time_propagation) THEN
209 WRITE (unit=output_unit, fmt="((T3,A,T61,F20.14))") &
210 "Energy difference to initial state:", rtp%energy_new - rtp%energy_old
211 END IF
212 IF (PRESENT(delta_iter)) THEN
213 WRITE (unit=output_unit, fmt="((T3,A,T61,E20.6))") &
214 "Convergence:", delta_iter
215 END IF
216 IF (rtp%converged) THEN
217 IF (run_type == real_time_propagation) THEN
218 WRITE (unit=output_unit, fmt="((T3,A,T61,F12.2))") &
219 "Time needed for propagation:", used_time
220 END IF
221 WRITE (unit=output_unit, fmt="(/,(T3,A,3X,F16.14))") &
222 "CONVERGENCE REACHED", rtp%energy_new - rtp%energy_old
223 END IF
224 END IF
225
226 IF (rtp%converged) THEN
227 IF (.NOT. rtp%linear_scaling) THEN
228 CALL get_rtp(rtp=rtp, mos_new=mos_new)
229 CALL rt_calculate_orthonormality(orthonormality, &
230 mos_new, matrix_s(1)%matrix)
231 IF (output_unit > 0) THEN
232 WRITE (output_unit, fmt="(/,(T3,A,T60,F20.10))") &
233 "Max deviation from orthonormalization:", orthonormality
234 END IF
235 END IF
236 END IF
237
238 IF (output_unit > 0) THEN
239 CALL m_flush(output_unit)
240 END IF
241 CALL cp_print_key_finished_output(output_unit, logger, rtp_section, &
242 "PRINT%PROGRAM_RUN_INFO")
243
244 IF (rtp%converged) THEN
245 dft_section => section_vals_get_subs_vals(input, "DFT")
246 IF (btest(cp_print_key_should_output(logger%iter_info, &
247 dft_section, "REAL_TIME_PROPAGATION%PRINT%FIELD"), cp_p_file)) THEN
248 CALL print_field_applied(qs_env, dft_section)
249 END IF
250 CALL make_moment(qs_env)
251 IF (btest(cp_print_key_should_output(logger%iter_info, &
252 dft_section, "REAL_TIME_PROPAGATION%PRINT%E_CONSTITUENTS"), cp_p_file)) THEN
253 CALL print_rtp_energy_components(qs_env, dft_section)
254 END IF
255 IF (.NOT. dft_control%qs_control%dftb) THEN
256 CALL write_available_results(qs_env=qs_env, rtp=rtp)
257 END IF
258
259 IF (rtp%linear_scaling) THEN
260 CALL get_rtp(rtp=rtp, rho_new=rho_new)
261
262 ! Probably have to rebuild the moment matrix, since atoms can also move, in principle
263 IF (dft_control%rtp_control%save_local_moments) THEN
264 ! Save the field value
265 IF (dft_control%apply_efield_field) THEN
266 CALL make_field(dft_control, field, qs_env%sim_step, qs_env%sim_time)
267 rtp%fields(:, rtp%istep + rtp%i_start + 1) = cmplx(field(:), 0.0, kind=dp)
268 END IF
269 IF (.NOT. dft_control%rtp_control%fixed_ions) THEN
270 CALL build_local_moment_matrix(qs_env, rtp%local_moments, 1, reference_point)
271 END IF
272 ! TODO : Is symmetric rho possible?
273 ! Spin + complex parts
274 ! Extensions setup
275 CALL calc_local_moment(rtp%local_moments, rho_new, &
276 rtp%local_moments_work, rtp%moments(:, :, rtp%istep + rtp%i_start + 1))
277 ! Time 1 is zero (start) time
278 rtp%times(rtp%istep + rtp%i_start + 1) = qs_env%sim_time
279 output_unit = cp_logger_get_default_io_unit(logger)
280 CALL print_moments(section_vals_get_subs_vals(rtp_section, "PRINT%MOMENTS"), output_unit, &
281 rtp%moments(:, :, rtp%istep + rtp%i_start + 1), qs_env%sim_time, rtp%track_imag_density)
282 END IF
283
284 IF (btest(cp_print_key_should_output(logger%iter_info, &
285 dft_section, "REAL_TIME_PROPAGATION%PRINT%RESTART"), cp_p_file)) THEN
286 CALL write_rt_p_to_restart(rho_new, .false.)
287 END IF
288 IF (btest(cp_print_key_should_output(logger%iter_info, &
289 dft_section, "REAL_TIME_PROPAGATION%PRINT%RESTART_HISTORY"), cp_p_file)) THEN
290 CALL write_rt_p_to_restart(rho_new, .true.)
291 END IF
292 IF (.NOT. dft_control%qs_control%dftb) THEN
293 !Not sure if these things could also work with dftb or not
294 IF (btest(cp_print_key_should_output(logger%iter_info, &
295 dft_section, "REAL_TIME_PROPAGATION%PRINT%CURRENT"), cp_p_file)) THEN
296 DO spin = 1, SIZE(rho_new)/2
297 CALL rt_current(qs_env, rho_new(2*spin)%matrix, dft_section, spin, SIZE(rho_new)/2)
298 END DO
299 END IF
300 END IF
301 ELSE
302 CALL get_rtp(rtp=rtp, mos_new=mos_new)
303 IF (.NOT. dft_control%qs_control%dftb .AND. .NOT. dft_control%qs_control%xtb) THEN
304 IF (rtp%track_imag_density) THEN
305 NULLIFY (p_im, p_xyz)
306 CALL dbcsr_allocate_matrix_set(p_xyz, 3)
307
308! Linear momentum operator
309! prepare for allocation
310 natom = SIZE(particle_set, 1)
311 ALLOCATE (first_sgf(natom))
312 ALLOCATE (last_sgf(natom))
313 CALL get_particle_set(particle_set, qs_kind_set, &
314 first_sgf=first_sgf, &
315 last_sgf=last_sgf)
316 ALLOCATE (row_blk_sizes(natom))
317 CALL dbcsr_convert_offsets_to_sizes(first_sgf, row_blk_sizes, last_sgf)
318 DEALLOCATE (first_sgf)
319 DEALLOCATE (last_sgf)
320
321 ALLOCATE (p_xyz(1)%matrix)
322 CALL dbcsr_create(matrix=p_xyz(1)%matrix, &
323 name="p_xyz", &
324 dist=dbcsr_dist, matrix_type=dbcsr_type_antisymmetric, &
325 row_blk_size=row_blk_sizes, col_blk_size=row_blk_sizes, &
326 mutable_work=.true.)
327 CALL cp_dbcsr_alloc_block_from_nbl(p_xyz(1)%matrix, sab_orb)
328 CALL dbcsr_set(p_xyz(1)%matrix, 0.0_dp)
329 DO i = 2, 3
330 ALLOCATE (p_xyz(i)%matrix)
331 CALL dbcsr_copy(p_xyz(i)%matrix, p_xyz(1)%matrix, "p_xyz-"//trim(adjustl(cp_to_string(i))))
332 CALL dbcsr_set(p_xyz(i)%matrix, 0.0_dp)
333 END DO
334 CALL build_lin_mom_matrix(qs_env, p_xyz)
335 DEALLOCATE (row_blk_sizes)
336
337 nspin = SIZE(mos_new)/2
338 CALL qs_rho_get(rho, rho_ao_im=p_im)
339 ALLOCATE (j_int(nspin, 3))
340 j_int = 0.0_dp
341
342 NULLIFY (tmp_ao)
343 CALL dbcsr_init_p(tmp_ao)
344 CALL dbcsr_create(tmp_ao, template=matrix_s(1)%matrix, matrix_type=dbcsr_type_no_symmetry, name="tmp")
345 CALL cp_dbcsr_alloc_block_from_nbl(tmp_ao, sab_all)
346 CALL dbcsr_set(tmp_ao, 0.0_dp)
347
348 DO i = 1, 3
349 strace = 0.0_dp
350 DO spin = 1, nspin
351 CALL dbcsr_set(tmp_ao, 0.0_dp)
352 CALL dbcsr_multiply("T", "N", 1.0_dp, p_im(spin)%matrix, p_xyz(i)%matrix, &
353 0.0_dp, tmp_ao)
354 CALL dbcsr_trace(tmp_ao, trace)
355 strace = strace + trace
356! dft_control%rtp_control%vec_pot(1)
357 j_int(spin, i) = -trace + dft_control%rtp_control%vec_pot(i)*nelectron_spin(spin)
358!! j_int(spin, i) = strace
359 END DO
360 END DO
361
362 ! Nonlocal pseudopotential correction to the integrated current:
363 ! the velocity operator contains i*[Vnl, r], which with C_d = [Vnl~, r_d]
364 ! contributes j_nl(d) = -Tr[P_re^T C_im(d)] - Tr[P_im^T C_re(d)].
365 ! Absent in all-electron calculations (no nonlocal projectors).
366 IF (ASSOCIATED(sap_ppnl)) THEN
367 ! the commutator is taken with the gauge-transformed Vnl only if
368 ! the propagation Hamiltonian uses it
369 vec_pot = 0.0_dp
370 IF (dft_control%rtp_control%velocity_gauge .AND. &
371 dft_control%rtp_control%nl_gauge_transform) THEN
372 vec_pot = dft_control%rtp_control%vec_pot
373 END IF
374
375 NULLIFY (matrix_cre, matrix_cim)
376 CALL dbcsr_allocate_matrix_set(matrix_cre, 3)
377 CALL dbcsr_allocate_matrix_set(matrix_cim, 3)
378 DO i = 1, 3
379 CALL dbcsr_init_p(matrix_cre(i)%matrix)
380 CALL dbcsr_create(matrix_cre(i)%matrix, template=matrix_s(1)%matrix, &
381 matrix_type=dbcsr_type_antisymmetric, &
382 name="nl commutator real part")
383 CALL cp_dbcsr_alloc_block_from_nbl(matrix_cre(i)%matrix, sab_orb)
384 CALL dbcsr_set(matrix_cre(i)%matrix, 0.0_dp)
385 CALL dbcsr_init_p(matrix_cim(i)%matrix)
386 CALL dbcsr_create(matrix_cim(i)%matrix, template=matrix_s(1)%matrix, &
387 matrix_type=dbcsr_type_symmetric, &
388 name="nl commutator imaginary part")
389 CALL cp_dbcsr_alloc_block_from_nbl(matrix_cim(i)%matrix, sab_orb)
390 CALL dbcsr_set(matrix_cim(i)%matrix, 0.0_dp)
391 END DO
392
393 CALL velocity_gauge_nl_commutator(qs_env, matrix_cre, matrix_cim, vec_pot)
394
395 CALL qs_rho_get(rho, rho_ao=p_re)
396 DO i = 1, 3
397 DO spin = 1, nspin
398 CALL dbcsr_set(tmp_ao, 0.0_dp)
399 CALL dbcsr_multiply("T", "N", 1.0_dp, p_re(spin)%matrix, matrix_cim(i)%matrix, &
400 0.0_dp, tmp_ao)
401 CALL dbcsr_trace(tmp_ao, trace)
402 j_int(spin, i) = j_int(spin, i) - trace
403 CALL dbcsr_set(tmp_ao, 0.0_dp)
404 CALL dbcsr_multiply("T", "N", 1.0_dp, p_im(spin)%matrix, matrix_cre(i)%matrix, &
405 0.0_dp, tmp_ao)
406 CALL dbcsr_trace(tmp_ao, trace)
407 j_int(spin, i) = j_int(spin, i) - trace
408 END DO
409 END DO
410
411 CALL dbcsr_deallocate_matrix_set(matrix_cre)
412 CALL dbcsr_deallocate_matrix_set(matrix_cim)
413 END IF
414
415 IF (btest(cp_print_key_should_output(logger%iter_info, &
416 dft_section, "REAL_TIME_PROPAGATION%PRINT%CURRENT_INT"), cp_p_file)) THEN
417
418 output_unit = cp_logger_get_default_io_unit(logger)
419 unit_nr = cp_print_key_unit_nr(logger, dft_section, &
420 "REAL_TIME_PROPAGATION%PRINT%CURRENT_INT", extension=".dat", is_new_file=new_file)
421
422 IF (output_unit > 0) THEN
423 WRITE (unit=output_unit, fmt="(T2,A,E24.16)") &
424 "RTP_CURRENT| CheckSum j_int=", sqrt(sum(j_int**2))
425 IF (new_file) THEN
426 IF (nspin == 2) THEN
427 WRITE (unit=unit_nr, fmt='("#",5X,A,4X,A,2X,A,2(10X,A),4X,A,2(10X,A))') &
428 "Step Nr.", "Time[fs]", "ALPHA jint[X]", "jint[Y]", "jint[Z]", &
429 "BETA jint[X]", "jint[Y]", "jint[Z]"
430 ELSE
431 WRITE (unit=unit_nr, fmt='("#",5X,A,4X,A,8X,A,2(10X,A))') "Step Nr.", "Time[fs]", &
432 "jint[X]", "jint[Y]", "jint[Z]"
433 END IF
434 END IF
435
436 IF (nspin == 2) THEN
437 WRITE (unit=unit_nr, fmt="(I10,F16.6,6(F16.8,1X))") qs_env%sim_step, qs_env%sim_time*femtoseconds, &
438 j_int(1, 1:3), j_int(2, 1:3)
439 ELSE
440 WRITE (unit=unit_nr, fmt="(I10,F16.6,3(F16.8,1X))") qs_env%sim_step, qs_env%sim_time*femtoseconds, &
441 j_int(1, 1:3)
442 END IF
443 END IF
444 CALL cp_print_key_finished_output(unit_nr, logger, dft_section, &
445 "REAL_TIME_PROPAGATION%PRINT%CURRENT_INT")
446 END IF
447 DEALLOCATE (j_int)
448
449 IF (btest(cp_print_key_should_output(logger%iter_info, &
450 dft_section, "REAL_TIME_PROPAGATION%PRINT%CURRENT"), cp_p_file)) THEN
451 DO spin = 1, nspin
452 CALL rt_current(qs_env, p_im(spin)%matrix, dft_section, spin, nspin)
453 END DO
454 END IF
455 CALL dbcsr_deallocate_matrix(tmp_ao)
456 CALL dbcsr_deallocate_matrix_set(p_xyz)
457 END IF
458
459! projection of molecular orbitals
460 IF (dft_control%rtp_control%is_proj_mo) THEN
461 DO n_proj = 1, SIZE(dft_control%rtp_control%proj_mo_list)
462 CALL compute_and_write_proj_mo(qs_env, mos_new, &
463 dft_control%rtp_control%proj_mo_list(n_proj)%proj_mo, n_proj)
464 END DO
465 END IF
466 END IF
467 CALL write_rt_mos_to_restart(qs_env%mos, mos_new, particle_set, &
468 dft_section, qs_kind_set)
469 END IF
470 END IF
471
472 rtp%energy_old = rtp%energy_new
473
474 IF (.NOT. rtp%converged .AND. rtp%iter >= dft_control%rtp_control%max_iter) THEN
475 CALL cp_abort(__location__, "EMD did not converge, either increase MAX_ITER "// &
476 "or use a smaller TIMESTEP")
477 END IF
478
479 END SUBROUTINE rt_prop_output
480
481! **************************************************************************************************
482!> \brief computes the effective orthonormality of a set of mos given an s-matrix
483!> orthonormality is the max deviation from unity of the C^T S C
484!> \param orthonormality ...
485!> \param mos_new ...
486!> \param matrix_s ...
487!> \author Florian Schiffmann (02.09)
488! **************************************************************************************************
489 SUBROUTINE rt_calculate_orthonormality(orthonormality, mos_new, matrix_s)
490 REAL(kind=dp), INTENT(out) :: orthonormality
491 TYPE(cp_fm_type), DIMENSION(:), POINTER :: mos_new
492 TYPE(dbcsr_type), OPTIONAL, POINTER :: matrix_s
493
494 CHARACTER(len=*), PARAMETER :: routinen = 'rt_calculate_orthonormality'
495
496 INTEGER :: handle, i, im, ispin, j, k, n, &
497 ncol_local, nrow_local, nspin, re
498 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
499 REAL(kind=dp) :: alpha, max_alpha, max_beta
500 TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
501 TYPE(cp_fm_type) :: overlap_re, svec_im, svec_re
502
503 NULLIFY (tmp_fm_struct)
504
505 CALL timeset(routinen, handle)
506
507 nspin = SIZE(mos_new)/2
508 max_alpha = 0.0_dp
509 max_beta = 0.0_dp
510 DO ispin = 1, nspin
511 re = ispin*2 - 1
512 im = ispin*2
513 ! get S*C
514 CALL cp_fm_create(svec_re, mos_new(im)%matrix_struct)
515 CALL cp_fm_create(svec_im, mos_new(im)%matrix_struct)
516 CALL cp_fm_get_info(mos_new(im), &
517 nrow_global=n, ncol_global=k)
518 CALL cp_dbcsr_sm_fm_multiply(matrix_s, mos_new(re), &
519 svec_re, k)
520 CALL cp_dbcsr_sm_fm_multiply(matrix_s, mos_new(im), &
521 svec_im, k)
522
523 ! get C^T (S*C)
524 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=k, ncol_global=k, &
525 para_env=mos_new(re)%matrix_struct%para_env, &
526 context=mos_new(re)%matrix_struct%context)
527 CALL cp_fm_create(overlap_re, tmp_fm_struct)
528
529 CALL cp_fm_struct_release(tmp_fm_struct)
530
531 CALL parallel_gemm('T', 'N', k, k, n, 1.0_dp, mos_new(re), &
532 svec_re, 0.0_dp, overlap_re)
533 CALL parallel_gemm('T', 'N', k, k, n, 1.0_dp, mos_new(im), &
534 svec_im, 1.0_dp, overlap_re)
535
536 CALL cp_fm_release(svec_re)
537 CALL cp_fm_release(svec_im)
538
539 CALL cp_fm_get_info(overlap_re, nrow_local=nrow_local, ncol_local=ncol_local, &
540 row_indices=row_indices, col_indices=col_indices)
541 DO i = 1, nrow_local
542 DO j = 1, ncol_local
543 alpha = overlap_re%local_data(i, j)
544 IF (row_indices(i) == col_indices(j)) alpha = alpha - 1.0_dp
545 max_alpha = max(max_alpha, abs(alpha))
546 END DO
547 END DO
548 CALL cp_fm_release(overlap_re)
549 END DO
550 CALL mos_new(1)%matrix_struct%para_env%max(max_alpha)
551 CALL mos_new(1)%matrix_struct%para_env%max(max_beta)
552 orthonormality = max_alpha
553
554 CALL timestop(handle)
555
556 END SUBROUTINE rt_calculate_orthonormality
557
558! **************************************************************************************************
559!> \brief computes the convergence criterion for RTP and EMD
560!> \param rtp ...
561!> \param matrix_s Overlap matrix without the derivatives
562!> \param delta_mos ...
563!> \param delta_eps ...
564!> \author Florian Schiffmann (02.09)
565! **************************************************************************************************
566
567 SUBROUTINE rt_convergence(rtp, matrix_s, delta_mos, delta_eps)
568 TYPE(rt_prop_type), POINTER :: rtp
569 TYPE(dbcsr_type), POINTER :: matrix_s
570 TYPE(cp_fm_type), DIMENSION(:), INTENT(IN) :: delta_mos
571 REAL(dp), INTENT(out) :: delta_eps
572
573 CHARACTER(len=*), PARAMETER :: routinen = 'rt_convergence'
574 REAL(kind=dp), PARAMETER :: one = 1.0_dp, zero = 0.0_dp
575
576 INTEGER :: handle, i, icol, im, ispin, j, lcol, &
577 lrow, nao, newdim, nmo, nspin, re
578 LOGICAL :: double_col, double_row
579 REAL(kind=dp) :: alpha, max_alpha
580 TYPE(cp_fm_struct_type), POINTER :: newstruct, newstruct1, tmp_fm_struct
581 TYPE(cp_fm_type) :: work, work1, work2
582 TYPE(cp_fm_type), DIMENSION(:), POINTER :: mos_new
583
584 NULLIFY (tmp_fm_struct)
585
586 CALL timeset(routinen, handle)
587
588 CALL get_rtp(rtp=rtp, mos_new=mos_new)
589
590 nspin = SIZE(delta_mos)/2
591 max_alpha = 0.0_dp
592
593 DO i = 1, SIZE(mos_new)
594 CALL cp_fm_scale_and_add(-one, delta_mos(i), one, mos_new(i))
595 END DO
596
597 DO ispin = 1, nspin
598 re = ispin*2 - 1
599 im = ispin*2
600
601 double_col = .true.
602 double_row = .false.
603 CALL cp_fm_struct_double(newstruct, &
604 delta_mos(re)%matrix_struct, &
605 delta_mos(re)%matrix_struct%context, &
606 double_col, &
607 double_row)
608
609 CALL cp_fm_create(work, matrix_struct=newstruct)
610 CALL cp_fm_create(work1, matrix_struct=newstruct)
611
612 CALL cp_fm_get_info(delta_mos(re), ncol_local=lcol, ncol_global=nmo, &
613 nrow_global=nao)
614 CALL cp_fm_get_info(work, ncol_global=newdim)
615
616 CALL cp_fm_set_all(work, zero, zero)
617
618 DO icol = 1, lcol
619 work%local_data(:, icol) = delta_mos(re)%local_data(:, icol)
620 work%local_data(:, icol + lcol) = delta_mos(im)%local_data(:, icol)
621 END DO
622
623 CALL cp_dbcsr_sm_fm_multiply(matrix_s, work, work1, ncol=newdim)
624
625 CALL cp_fm_release(work)
626
627 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=nmo, ncol_global=nmo, &
628 para_env=delta_mos(re)%matrix_struct%para_env, &
629 context=delta_mos(re)%matrix_struct%context)
630 CALL cp_fm_struct_double(newstruct1, &
631 tmp_fm_struct, &
632 delta_mos(re)%matrix_struct%context, &
633 double_col, &
634 double_row)
635
636 CALL cp_fm_create(work, matrix_struct=newstruct1)
637 CALL cp_fm_create(work2, matrix_struct=newstruct1)
638
639 CALL parallel_gemm("T", "N", nmo, newdim, nao, one, delta_mos(re), &
640 work1, zero, work)
641
642 CALL parallel_gemm("T", "N", nmo, newdim, nao, one, delta_mos(im), &
643 work1, zero, work2)
644
645 CALL cp_fm_get_info(work, nrow_local=lrow)
646 DO i = 1, lrow
647 DO j = 1, lcol
648 alpha = sqrt((work%local_data(i, j) + work2%local_data(i, j + lcol))**2 + &
649 (work%local_data(i, j + lcol) - work2%local_data(i, j))**2)
650 max_alpha = max(max_alpha, abs(alpha))
651 END DO
652 END DO
653
654 CALL cp_fm_release(work)
655 CALL cp_fm_release(work1)
656 CALL cp_fm_release(work2)
657 CALL cp_fm_struct_release(tmp_fm_struct)
658 CALL cp_fm_struct_release(newstruct)
659 CALL cp_fm_struct_release(newstruct1)
660
661 END DO
662
663 CALL delta_mos(1)%matrix_struct%para_env%max(max_alpha)
664 delta_eps = sqrt(max_alpha)
665
666 CALL timestop(handle)
667
668 END SUBROUTINE rt_convergence
669
670! **************************************************************************************************
671!> \brief computes the convergence criterion for RTP and EMD based on the density matrix
672!> \param rtp ...
673!> \param delta_P ...
674!> \param delta_eps ...
675!> \author Samuel Andermatt (02.14)
676! **************************************************************************************************
677
678 SUBROUTINE rt_convergence_density(rtp, delta_P, delta_eps)
679
680 TYPE(rt_prop_type), POINTER :: rtp
681 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: delta_p
682 REAL(dp), INTENT(out) :: delta_eps
683
684 CHARACTER(len=*), PARAMETER :: routinen = 'rt_convergence_density'
685 REAL(kind=dp), PARAMETER :: one = 1.0_dp, zero = 0.0_dp
686
687 INTEGER :: col_atom, handle, i, ispin, row_atom
688 REAL(dp) :: alpha, max_alpha
689 REAL(dp), DIMENSION(:, :), POINTER :: block_values
690 TYPE(dbcsr_iterator_type) :: iter
691 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_new
692 TYPE(dbcsr_type), POINTER :: tmp
693 TYPE(mp_comm_type) :: group
694
695 CALL timeset(routinen, handle)
696
697 CALL get_rtp(rtp=rtp, rho_new=rho_new)
698
699 DO i = 1, SIZE(rho_new)
700 CALL dbcsr_add(delta_p(i)%matrix, rho_new(i)%matrix, one, -one)
701 END DO
702 !get the maximum value of delta_P
703 DO i = 1, SIZE(delta_p)
704 !square all entries of both matrices
705 CALL dbcsr_iterator_start(iter, delta_p(i)%matrix)
706 DO WHILE (dbcsr_iterator_blocks_left(iter))
707 CALL dbcsr_iterator_next_block(iter, row_atom, col_atom, block_values)
708 block_values = block_values*block_values
709 END DO
710 CALL dbcsr_iterator_stop(iter)
711 END DO
712 NULLIFY (tmp)
713 ALLOCATE (tmp)
714 CALL dbcsr_create(tmp, template=delta_p(1)%matrix, matrix_type="N")
715 DO ispin = 1, SIZE(delta_p)/2
716 CALL dbcsr_desymmetrize(delta_p(2*ispin - 1)%matrix, tmp)
717 CALL dbcsr_add(delta_p(2*ispin)%matrix, tmp, one, one)
718 END DO
719 !the absolute values are now in the even entries of delta_P
720 max_alpha = zero
721 DO ispin = 1, SIZE(delta_p)/2
722 CALL dbcsr_iterator_start(iter, delta_p(2*ispin)%matrix)
723 DO WHILE (dbcsr_iterator_blocks_left(iter))
724 CALL dbcsr_iterator_next_block(iter, row_atom, col_atom, block_values)
725 alpha = maxval(block_values)
726 IF (alpha > max_alpha) max_alpha = alpha
727 END DO
728 CALL dbcsr_iterator_stop(iter)
729 END DO
730 CALL dbcsr_get_info(delta_p(1)%matrix, group=group)
731 CALL group%max(max_alpha)
732 delta_eps = sqrt(max_alpha)
733 CALL dbcsr_deallocate_matrix(tmp)
734 CALL timestop(handle)
735
736 END SUBROUTINE rt_convergence_density
737
738! **************************************************************************************************
739!> \brief interface to qs_moments. Does only work for nonperiodic dipole
740!> \param qs_env ...
741!> \author Florian Schiffmann (02.09)
742! **************************************************************************************************
743
744 SUBROUTINE make_moment(qs_env)
745
746 TYPE(qs_environment_type), POINTER :: qs_env
747
748 CHARACTER(len=*), PARAMETER :: routinen = 'make_moment'
749
750 INTEGER :: handle, output_unit
751 TYPE(cp_logger_type), POINTER :: logger
752 TYPE(dft_control_type), POINTER :: dft_control
753
754 CALL timeset(routinen, handle)
755
756 NULLIFY (dft_control)
757
758 logger => cp_get_default_logger()
759 output_unit = cp_logger_get_default_io_unit(logger)
760 CALL get_qs_env(qs_env, dft_control=dft_control)
761 IF (dft_control%qs_control%dftb) THEN
762 CALL scf_post_calculation_tb(qs_env, "DFTB", .false.)
763 ELSE IF (dft_control%qs_control%xtb) THEN
764 CALL scf_post_calculation_tb(qs_env, "xTB", .false.)
765 ELSE
766 CALL qs_scf_post_moments(qs_env%input, logger, qs_env, output_unit)
767 END IF
768 CALL timestop(handle)
769
770 END SUBROUTINE make_moment
771
772! **************************************************************************************************
773!> \brief Reports the sparsity pattern of the complex density matrix
774!> \param filter_eps ...
775!> \param rho ...
776!> \author Samuel Andermatt (09.14)
777! **************************************************************************************************
778
779 SUBROUTINE report_density_occupation(filter_eps, rho)
780
781 REAL(kind=dp) :: filter_eps
782 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho
783
784 CHARACTER(len=*), PARAMETER :: routinen = 'report_density_occupation'
785
786 INTEGER :: handle, i, im, ispin, re, unit_nr
787 REAL(kind=dp) :: eps, occ
788 TYPE(cp_logger_type), POINTER :: logger
789 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: tmp
790
791 CALL timeset(routinen, handle)
792
793 logger => cp_get_default_logger()
794 unit_nr = cp_logger_get_default_io_unit(logger)
795 NULLIFY (tmp)
796 CALL dbcsr_allocate_matrix_set(tmp, SIZE(rho))
797 DO i = 1, SIZE(rho)
798 CALL dbcsr_init_p(tmp(i)%matrix)
799 CALL dbcsr_create(tmp(i)%matrix, template=rho(i)%matrix)
800 CALL dbcsr_copy(tmp(i)%matrix, rho(i)%matrix)
801 END DO
802 DO ispin = 1, SIZE(rho)/2
803 re = 2*ispin - 1
804 im = 2*ispin
805 eps = max(filter_eps, 1.0e-11_dp)
806 DO WHILE (eps < 1.1_dp)
807 CALL dbcsr_filter(tmp(re)%matrix, eps)
808 occ = dbcsr_get_occupation(tmp(re)%matrix)
809 IF (unit_nr > 0) WRITE (unit_nr, fmt="((T3,A,I1,A,F15.12,A,T61,F20.10))") "Occupation of rho spin ", &
810 ispin, " eps ", eps, " real: ", occ
811 eps = eps*10
812 END DO
813 eps = max(filter_eps, 1.0e-11_dp)
814 DO WHILE (eps < 1.1_dp)
815 CALL dbcsr_filter(tmp(im)%matrix, eps)
816 occ = dbcsr_get_occupation(tmp(im)%matrix)
817 IF (unit_nr > 0) WRITE (unit_nr, fmt="((T3,A,I1,A,F15.12,A,T61,F20.10))") "Occupation of rho spin ", &
818 ispin, " eps ", eps, " imag: ", occ
819 eps = eps*10.0_dp
820 END DO
821 END DO
822 CALL dbcsr_deallocate_matrix_set(tmp)
823 CALL timestop(handle)
824
825 END SUBROUTINE report_density_occupation
826
827! **************************************************************************************************
828!> \brief Writes the density matrix and the atomic positions to a restart file
829!> \param rho_new ...
830!> \param history ...
831!> \author Samuel Andermatt (09.14)
832! **************************************************************************************************
833
834 SUBROUTINE write_rt_p_to_restart(rho_new, history)
835
836 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_new
837 LOGICAL :: history
838
839 CHARACTER(LEN=*), PARAMETER :: routinen = 'write_rt_p_to_restart'
840
841 CHARACTER(LEN=default_path_length) :: file_name, project_name
842 INTEGER :: handle, im, ispin, re, unit_nr
843 REAL(kind=dp) :: cs_pos
844 TYPE(cp_logger_type), POINTER :: logger
845
846 CALL timeset(routinen, handle)
847 logger => cp_get_default_logger()
848 IF (logger%para_env%is_source()) THEN
849 unit_nr = cp_logger_get_default_unit_nr(logger, local=.true.)
850 ELSE
851 unit_nr = -1
852 END IF
853
854 project_name = logger%iter_info%project_name
855 DO ispin = 1, SIZE(rho_new)/2
856 re = 2*ispin - 1
857 im = 2*ispin
858 IF (history) THEN
859 WRITE (file_name, '(A,I0,A)') &
860 trim(project_name)//"_LS_DM_SPIN_RE", ispin, "_"//trim(cp_iter_string(logger%iter_info))//"_RESTART.dm"
861 ELSE
862 WRITE (file_name, '(A,I0,A)') trim(project_name)//"_LS_DM_SPIN_RE", ispin, "_RESTART.dm"
863 END IF
864 cs_pos = dbcsr_checksum(rho_new(re)%matrix, pos=.true.)
865 IF (unit_nr > 0) THEN
866 WRITE (unit_nr, '(T2,A,E20.8)') "Writing restart DM "//trim(file_name)//" with checksum: ", cs_pos
867 END IF
868 CALL dbcsr_binary_write(rho_new(re)%matrix, file_name)
869 IF (history) THEN
870 WRITE (file_name, '(A,I0,A)') &
871 trim(project_name)//"_LS_DM_SPIN_IM", ispin, "_"//trim(cp_iter_string(logger%iter_info))//"_RESTART.dm"
872 ELSE
873 WRITE (file_name, '(A,I0,A)') trim(project_name)//"_LS_DM_SPIN_IM", ispin, "_RESTART.dm"
874 END IF
875 cs_pos = dbcsr_checksum(rho_new(im)%matrix, pos=.true.)
876 IF (unit_nr > 0) THEN
877 WRITE (unit_nr, '(T2,A,E20.8)') "Writing restart DM "//trim(file_name)//" with checksum: ", cs_pos
878 END IF
879 CALL dbcsr_binary_write(rho_new(im)%matrix, file_name)
880 END DO
881
882 CALL timestop(handle)
883
884 END SUBROUTINE write_rt_p_to_restart
885
886! **************************************************************************************************
887!> \brief Collocation of the current and printing of it in a cube file
888!> \param qs_env ...
889!> \param P_im ...
890!> \param dft_section ...
891!> \param spin ...
892!> \param nspin ...
893!> \author Samuel Andermatt (06.15)
894! **************************************************************************************************
895 SUBROUTINE rt_current(qs_env, P_im, dft_section, spin, nspin)
896 TYPE(qs_environment_type), POINTER :: qs_env
897 TYPE(dbcsr_type), POINTER :: p_im
898 TYPE(section_vals_type), POINTER :: dft_section
899 INTEGER :: spin, nspin
900
901 CHARACTER(len=*), PARAMETER :: routinen = 'rt_current'
902
903 CHARACTER(len=1) :: char_spin
904 CHARACTER(len=14) :: ext
905 CHARACTER(len=2) :: sdir
906 INTEGER :: dir, handle, print_unit
907 INTEGER, DIMENSION(:), POINTER :: stride
908 LOGICAL :: mpi_io
909 TYPE(cp_logger_type), POINTER :: logger
910 TYPE(current_env_type) :: current_env
911 TYPE(dbcsr_type), POINTER :: tmp, zero
912 TYPE(particle_list_type), POINTER :: particles
913 TYPE(pw_c1d_gs_type) :: gs
914 TYPE(pw_env_type), POINTER :: pw_env
915 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
916 TYPE(pw_r3d_rs_type) :: rs
917 TYPE(qs_subsys_type), POINTER :: subsys
918
919 CALL timeset(routinen, handle)
920
921 logger => cp_get_default_logger()
922 CALL get_qs_env(qs_env=qs_env, subsys=subsys, pw_env=pw_env)
923 CALL qs_subsys_get(subsys, particles=particles)
924 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool)
925
926 NULLIFY (zero, tmp)
927 ALLOCATE (zero, tmp)
928 CALL dbcsr_create(zero, template=p_im)
929 CALL dbcsr_copy(zero, p_im)
930 CALL dbcsr_set(zero, 0.0_dp)
931 CALL dbcsr_create(tmp, template=p_im)
932 CALL dbcsr_copy(tmp, p_im)
933 IF (nspin == 1) THEN
934 CALL dbcsr_scale(tmp, 0.5_dp)
935 END IF
936 current_env%gauge = -1
937 current_env%gauge_init = .false.
938 CALL auxbas_pw_pool%create_pw(rs)
939 CALL auxbas_pw_pool%create_pw(gs)
940
941 NULLIFY (stride)
942 ALLOCATE (stride(3))
943
944 DO dir = 1, 3
945
946 CALL pw_zero(rs)
947 CALL pw_zero(gs)
948
949 CALL calculate_jrho_resp(zero, tmp, zero, zero, dir, dir, rs, gs, qs_env, current_env, retain_rsgrid=.true.)
950
951 stride = section_get_ivals(dft_section, "REAL_TIME_PROPAGATION%PRINT%CURRENT%STRIDE")
952
953 IF (dir == 1) THEN
954 sdir = "-x"
955 ELSE IF (dir == 2) THEN
956 sdir = "-y"
957 ELSE
958 sdir = "-z"
959 END IF
960 WRITE (char_spin, "(I1)") spin
961
962 ext = "-SPIN-"//char_spin//sdir//".cube"
963 mpi_io = .true.
964 print_unit = cp_print_key_unit_nr(logger, dft_section, "REAL_TIME_PROPAGATION%PRINT%CURRENT", &
965 extension=ext, file_status="REPLACE", file_action="WRITE", &
966 log_filename=.false., mpi_io=mpi_io)
967
968 CALL cp_pw_to_cube(rs, print_unit, "EMD current", particles=particles, stride=stride, &
969 mpi_io=mpi_io)
970
971 CALL cp_print_key_finished_output(print_unit, logger, dft_section, "REAL_TIME_PROPAGATION%PRINT%CURRENT", &
972 mpi_io=mpi_io)
973
974 END DO
975
976 CALL auxbas_pw_pool%give_back_pw(rs)
977 CALL auxbas_pw_pool%give_back_pw(gs)
978
979 CALL dbcsr_deallocate_matrix(zero)
980 CALL dbcsr_deallocate_matrix(tmp)
981
982 DEALLOCATE (stride)
983
984 CALL timestop(handle)
985
986 END SUBROUTINE rt_current
987
988! **************************************************************************************************
989!> \brief Interface routine to trigger writing of results available from normal
990!> SCF. Can write MO-dependent and MO free results (needed for call from
991!> the linear scaling code)
992!> Update: trigger also some of prints for time-dependent runs
993!> \param qs_env ...
994!> \param rtp ...
995!> \par History
996!> 2022-11 Update [Guillaume Le Breton]
997! **************************************************************************************************
998 SUBROUTINE write_available_results(qs_env, rtp)
999 TYPE(qs_environment_type), POINTER :: qs_env
1000 TYPE(rt_prop_type), POINTER :: rtp
1001
1002 CHARACTER(len=*), PARAMETER :: routinen = 'write_available_results'
1003
1004 INTEGER :: handle
1005 TYPE(qs_scf_env_type), POINTER :: scf_env
1006
1007 CALL timeset(routinen, handle)
1008
1009 CALL get_qs_env(qs_env, scf_env=scf_env)
1010 IF (rtp%linear_scaling) THEN
1011 CALL write_mo_free_results(qs_env)
1012 ELSE
1013 CALL write_mo_free_results(qs_env)
1014 CALL write_mo_dependent_results(qs_env, scf_env)
1015 ! Time-dependent MO print
1016 CALL write_rtp_mos_to_output_unit(qs_env, rtp)
1017 CALL write_rtp_mo_cubes(qs_env, rtp)
1018 END IF
1019
1020 CALL timestop(handle)
1021
1022 END SUBROUTINE write_available_results
1023
1024! **************************************************************************************************
1025!> \brief Print the field applied to the system. Either the electric
1026!> field or the vector potential depending on the gauge used
1027!> \param qs_env ...
1028!> \param dft_section ...
1029!> \par History
1030!> 2023-01 Created [Guillaume Le Breton]
1031! **************************************************************************************************
1032 SUBROUTINE print_field_applied(qs_env, dft_section)
1033 TYPE(qs_environment_type), POINTER :: qs_env
1034 TYPE(section_vals_type), POINTER :: dft_section
1035
1036 CHARACTER(LEN=3), DIMENSION(3) :: rlab
1037 CHARACTER(LEN=default_path_length) :: filename
1038 INTEGER :: i, i_step, output_unit, unit_nr
1039 LOGICAL :: new_file
1040 REAL(kind=dp) :: field(3)
1041 TYPE(cp_logger_type), POINTER :: logger
1042 TYPE(dft_control_type), POINTER :: dft_control
1043 TYPE(rt_prop_type), POINTER :: rtp
1044
1045 NULLIFY (dft_control)
1046
1047 logger => cp_get_default_logger()
1048 output_unit = cp_logger_get_default_io_unit(logger)
1049
1050 CALL get_qs_env(qs_env, dft_control=dft_control, rtp=rtp)
1051
1052 i_step = rtp%istep
1053
1054 unit_nr = cp_print_key_unit_nr(logger, dft_section, &
1055 "REAL_TIME_PROPAGATION%PRINT%FIELD", extension=".dat", is_new_file=new_file)
1056
1057 IF (output_unit > 0) THEN
1058 rlab = [CHARACTER(LEN=3) :: "X", "Y", "Z"]
1059 IF (unit_nr /= output_unit) THEN
1060 INQUIRE (unit=unit_nr, name=filename)
1061 WRITE (unit=output_unit, fmt="(/,T2,A,2(/,T3,A),/)") &
1062 "FIELD", "The field applied is written to the file:", &
1063 trim(filename)
1064 ELSE
1065 WRITE (unit=output_unit, fmt="(/,T2,A)") "FIELD APPLIED [a.u.]"
1066 WRITE (unit=output_unit, fmt="(T5,3(A,A,E16.8,1X))") &
1067 (trim(rlab(i)), "=", dft_control%rtp_control%field(i), i=1, 3)
1068 END IF
1069
1070 IF (new_file) THEN
1071 IF (dft_control%apply_efield_field) THEN
1072 WRITE (unit=unit_nr, fmt='("#",5X,A,8X,A,3(6X,A))') "Step Nr.", "Time[fs]", " Field X", " Field Y", " Field Z"
1073 ELSE IF (dft_control%apply_vector_potential) THEN
1074 WRITE (unit=unit_nr, fmt='("#",5X,A,8X,A,6(6X,A))') "Step Nr.", "Time[fs]", " Field X", " Field Y", " Field Z", &
1075 " Vec. Pot. X", " Vec. Pot. Y", " Vec. Pot. Z"
1076 END IF
1077 END IF
1078
1079 field = 0.0_dp
1080 IF (dft_control%apply_efield_field) THEN
1081 CALL make_field(dft_control, field, qs_env%sim_step, qs_env%sim_time)
1082 WRITE (unit=unit_nr, fmt="(I10,F16.6,3(F16.8,1X))") qs_env%sim_step, qs_env%sim_time*femtoseconds, &
1083 field(1), field(2), field(3)
1084! DO i=1,3
1085! IF (ABS(field(i))< 10E-10) field(i) = 0.0_dp
1086! END IF
1087 ELSE IF (dft_control%apply_vector_potential) THEN
1088 WRITE (unit=unit_nr, fmt="(I10,F16.6,6(F16.8,1X))") qs_env%sim_step, qs_env%sim_time*femtoseconds, &
1089 dft_control%rtp_control%field(1), dft_control%rtp_control%field(2), dft_control%rtp_control%field(3), &
1090 dft_control%rtp_control%vec_pot(1), dft_control%rtp_control%vec_pot(2), dft_control%rtp_control%vec_pot(3)
1091 END IF
1092
1093 END IF
1094
1095 CALL cp_print_key_finished_output(unit_nr, logger, dft_section, &
1096 "REAL_TIME_PROPAGATION%PRINT%FIELD")
1097
1098 END SUBROUTINE print_field_applied
1099
1100! **************************************************************************************************
1101!> \brief Print the components of the total energy used in an RTP calculation
1102!> \param qs_env ...
1103!> \param dft_section ...
1104!> \par History
1105!> 2024-02 Created [ANB]
1106! **************************************************************************************************
1107 SUBROUTINE print_rtp_energy_components(qs_env, dft_section)
1108 TYPE(qs_environment_type), POINTER :: qs_env
1109 TYPE(section_vals_type), POINTER :: dft_section
1110
1111 CHARACTER(LEN=default_path_length) :: filename
1112 INTEGER :: i_step, output_unit, unit_nr
1113 LOGICAL :: new_file
1114 TYPE(cp_logger_type), POINTER :: logger
1115 TYPE(dft_control_type), POINTER :: dft_control
1116 TYPE(qs_energy_type), POINTER :: energy
1117 TYPE(rt_prop_type), POINTER :: rtp
1118
1119 NULLIFY (dft_control, energy, rtp)
1120
1121 logger => cp_get_default_logger()
1122 output_unit = cp_logger_get_default_io_unit(logger)
1123
1124 CALL get_qs_env(qs_env, dft_control=dft_control, rtp=rtp, energy=energy)
1125 i_step = rtp%istep
1126
1127 unit_nr = cp_print_key_unit_nr(logger, dft_section, &
1128 "REAL_TIME_PROPAGATION%PRINT%E_CONSTITUENTS", extension=".ener", &
1129 file_action="WRITE", is_new_file=new_file)
1130
1131 IF (output_unit > 0) THEN
1132 IF (unit_nr /= output_unit) THEN
1133 INQUIRE (unit=unit_nr, name=filename)
1134 WRITE (unit=output_unit, fmt="(/,T2,A,2(/,T3,A),/)") &
1135 "ENERGY_CONSTITUENTS", "Total Energy constituents written to file:", &
1136 trim(filename)
1137 ELSE
1138 WRITE (unit=output_unit, fmt="(/,T2,A)") "ENERGY_CONSTITUENTS"
1139 END IF
1140
1141 IF (new_file) THEN
1142 ! NOTE that these are not all terms contributing to the total energy for RTP, only a selection of those
1143 ! most significant / impactful. Therefore the printed components likely will not add up to the total energy.
1144 WRITE (unit=unit_nr, fmt='("#",5X,A,8X,A,10(6X,A))') "Step Nr.", "Time[fs]", &
1145 "Total ener.[a.u.]", "core[a.u.] ", " overlap [a.u.]", "hartree[a.u.]", " exc. [a.u.] ", &
1146 " hartree 1c[a.u.]", "exc 1c[a.u.] ", "exc admm[a.u.]", "exc 1c admm[a.u.]", "efield LG"
1147
1148 END IF
1149 WRITE (unit=unit_nr, fmt="(I10,F20.6,10(F20.9))") &
1150 qs_env%sim_step, qs_env%sim_time*femtoseconds, &
1151 energy%total, energy%core, energy%core_overlap, energy%hartree, energy%exc, &
1152 energy%hartree_1c, energy%exc1, energy%exc_aux_fit, energy%exc1_aux_fit, energy%efield_core
1153
1154 END IF
1155
1156 CALL cp_print_key_finished_output(unit_nr, logger, dft_section, &
1157 "REAL_TIME_PROPAGATION%PRINT%E_CONSTITUENTS")
1158
1159 END SUBROUTINE print_rtp_energy_components
1160
1161! **************************************************************************************************
1162!> \brief Print the dipole moments into a file
1163!> \param moments_section Section of the input defining the file/stream to print the moments to
1164!> \param info_unit Unit where standard output from the program is written - for add. identifiers
1165!> \param moments Actual moment values (for specific time step)
1166!> \param time Current simulation time
1167!> \param imag_opt Whether to calculate the imaginary part
1168!> \param append_opt ...
1169!> \par History
1170!> 10.2025 Created [Marek]
1171! **************************************************************************************************
1172 SUBROUTINE print_moments(moments_section, info_unit, moments, time, imag_opt, append_opt)
1173 TYPE(section_vals_type), POINTER :: moments_section
1174 INTEGER :: info_unit
1175 COMPLEX(kind=dp), DIMENSION(:, :) :: moments
1176 REAL(kind=dp), OPTIONAL :: time
1177 LOGICAL, OPTIONAL :: imag_opt, append_opt
1178
1179 CHARACTER(len=14), DIMENSION(4) :: file_extensions
1180 CHARACTER(len=21) :: prefix
1181 COMPLEX(kind=dp), DIMENSION(3, 1) :: moment_t
1182 INTEGER :: i, ndir, nspin, print_unit
1183 LOGICAL :: append, imaginary
1184 TYPE(cp_logger_type), POINTER :: logger
1185
1186! Index 1 : spin, Index 2 : direction
1187
1188 nspin = SIZE(moments, 1)
1189 ndir = SIZE(moments, 2)
1190
1191 IF (nspin < 1) cpabort("Zero spin index size in print moments!")
1192 IF (ndir < 1) cpabort("Zero direction index size in print moments!")
1193
1194 imaginary = .true.
1195 IF (PRESENT(imag_opt)) imaginary = imag_opt
1196
1197 append = .true.
1198 IF (PRESENT(append_opt)) append = append_opt
1199
1200 ! Get the program run info unit and target unit
1201 ! If these are the same (most likely the case of __STD_OUT__), add
1202 ! extra identifier to the printed output
1203 file_extensions(1) = "_SPIN_A_RE.dat"
1204 file_extensions(2) = "_SPIN_A_IM.dat"
1205 file_extensions(3) = "_SPIN_B_RE.dat"
1206 file_extensions(4) = "_SPIN_B_IM.dat"
1207 logger => cp_get_default_logger()
1208 DO i = 1, nspin
1209 moment_t(:, 1) = moments(i, :)
1210 ! Real part
1211 print_unit = cp_print_key_unit_nr(logger, moments_section, extension=file_extensions(2*i - 1))
1212 IF (print_unit == info_unit) THEN
1213 ! print with prefix
1214 prefix = " MOMENTS_TRACE_RE|"
1215 IF (append) THEN
1216 ! Print without headers
1217 CALL print_rt_file(print_unit, xvals=[time], yvals=moment_t, &
1218 prefix=prefix, prefix_format="(A18)", &
1219 xscale_opt=femtoseconds, comp_opt=rt_file_comp_real)
1220 ELSE
1221 ! Print with headers
1222 CALL print_rt_file(print_unit, &
1223 headers=["# Time [fs]", " re(mom_t) x [at.u.]", &
1224 " re(mom_t) y [at.u.]", " re(mom_t) z [at.u.]"], &
1225 xvals=[time], yvals=moment_t, &
1226 prefix=prefix, prefix_format="(A18)", &
1227 xscale_opt=femtoseconds, comp_opt=rt_file_comp_real)
1228 END IF
1229 ELSE
1230 ! Print without prefix
1231 IF (append) THEN
1232 ! Print without headers
1233 CALL print_rt_file(print_unit, xvals=[time], yvals=moment_t, &
1234 xscale_opt=femtoseconds, comp_opt=rt_file_comp_real)
1235 ELSE
1236 ! Print with headers
1237 CALL print_rt_file(print_unit, &
1238 headers=["# Time [fs]", " re(mom_t) x [at.u.]", &
1239 " re(mom_t) y [at.u.]", " re(mom_t) z [at.u.]"], &
1240 xvals=[time], yvals=moment_t, &
1241 xscale_opt=femtoseconds, comp_opt=rt_file_comp_real)
1242 END IF
1243 END IF
1244 CALL cp_print_key_finished_output(print_unit, logger, moments_section)
1245 ! Same for imaginary part
1246 IF (imaginary) THEN
1247 print_unit = cp_print_key_unit_nr(logger, moments_section, extension=file_extensions(2*i))
1248 IF (print_unit == info_unit) THEN
1249 ! print with prefix
1250 prefix = " MOMENTS_TRACE_IM|"
1251 IF (append) THEN
1252 ! Print without headers
1253 CALL print_rt_file(print_unit, xvals=[time], yvals=moment_t, &
1254 prefix=prefix, prefix_format="(A18)", &
1255 xscale_opt=femtoseconds, comp_opt=rt_file_comp_imag)
1256 ELSE
1257 ! Print with headers
1258 CALL print_rt_file(print_unit, &
1259 headers=["# Time [fs]", " im(mom_t) x [at.u.]", &
1260 " im(mom_t) y [at.u.]", " im(mom_t) z [at.u.]"], &
1261 xvals=[time], yvals=moment_t, &
1262 prefix=prefix, prefix_format="(A18)", &
1263 xscale_opt=femtoseconds, comp_opt=rt_file_comp_imag)
1264 END IF
1265 ELSE
1266 ! Print without prefix
1267 IF (append) THEN
1268 ! Print without headers
1269 CALL print_rt_file(print_unit, xvals=[time], yvals=moment_t, &
1270 xscale_opt=femtoseconds, comp_opt=rt_file_comp_imag)
1271 ELSE
1272 ! Print with headers
1273 CALL print_rt_file(print_unit, &
1274 headers=["# Time [fs]", " im(mom_t) x [at.u.]", &
1275 " im(mom_t) y [at.u.]", " im(mom_t) z [at.u.]"], &
1276 xvals=[time], yvals=moment_t, &
1277 xscale_opt=femtoseconds, comp_opt=rt_file_comp_imag)
1278 END IF
1279 END IF
1280 CALL cp_print_key_finished_output(print_unit, logger, moments_section)
1281 END IF
1282 END DO
1283
1284 END SUBROUTINE print_moments
1285
1286! **************************************************************************************************
1287!> \brief Calculate the values of real/imaginary parts of moments in all directions
1288!> \param moment_matrices Local matrix representations of dipole (position) operator
1289!> \param density_matrices Density matrices (spin and real+complex parts)
1290!> \param work Extra dbcsr matrix for work
1291!> \param moment Resulting moments (spin and direction)
1292!> \param imag_opt Whether to calculate the imaginary part of the moment
1293!> \par History
1294!> 10.2025 Created [Marek]
1295! **************************************************************************************************
1296 SUBROUTINE calc_local_moment(moment_matrices, density_matrices, work, moment, imag_opt)
1297 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: moment_matrices, density_matrices
1298 TYPE(dbcsr_type) :: work
1299 COMPLEX(kind=dp), DIMENSION(:, :) :: moment
1300 LOGICAL, OPTIONAL :: imag_opt
1301
1302 INTEGER :: i, k, nspin
1303 LOGICAL :: imag
1304 REAL(kind=dp) :: real_moment
1305
1306 imag = .false.
1307 IF (PRESENT(imag_opt)) imag = imag_opt
1308 nspin = SIZE(density_matrices)/2
1309
1310 DO i = 1, nspin
1311 DO k = 1, 3
1312 CALL dbcsr_multiply("N", "N", -1.0_dp, &
1313 density_matrices(2*i - 1)%matrix, moment_matrices(k)%matrix, &
1314 0.0_dp, work)
1315 CALL dbcsr_trace(work, real_moment)
1316 moment(i, k) = cmplx(real_moment, 0.0, kind=dp)
1317 IF (imag) THEN
1318 CALL dbcsr_multiply("N", "N", -1.0_dp, &
1319 density_matrices(2*i)%matrix, moment_matrices(k)%matrix, &
1320 0.0_dp, work)
1321 CALL dbcsr_trace(work, real_moment)
1322 moment(i, k) = moment(i, k) + cmplx(0.0, real_moment, kind=dp)
1323 END IF
1324 END DO
1325 END DO
1326
1327 END SUBROUTINE calc_local_moment
1328
1329! **************************************************************************************************
1330!> \brief Calculate and print the Fourier transforms + polarizabilites from moment trace
1331!> \param rtp_section The RTP input section (needed to access PRINT configurations)
1332!> \param moments Moment trace
1333!> \param times Corresponding times
1334!> \param fields Corresponding fields
1335!> \param rtc rt_control_type that includes metadata
1336!> \param info_opt ...
1337!> \param cell If present, used to change the delta peak representation to be in units of reciprocal lattice
1338!> \par History
1339!> 10.2025 Created [Marek]
1340! **************************************************************************************************
1341 SUBROUTINE print_ft(rtp_section, moments, times, fields, rtc, info_opt, cell)
1342 TYPE(section_vals_type), POINTER :: rtp_section
1343 COMPLEX(kind=dp), DIMENSION(:, :, :), POINTER :: moments
1344 REAL(kind=dp), DIMENSION(:), POINTER :: times
1345 COMPLEX(kind=dp), DIMENSION(:, :), POINTER :: fields
1346 TYPE(rtp_control_type), POINTER :: rtc
1347 INTEGER, OPTIONAL :: info_opt
1348 TYPE(cell_type), OPTIONAL, POINTER :: cell
1349
1350 CHARACTER(len=*), PARAMETER :: routinen = 'print_ft'
1351
1352 CHARACTER(len=11), DIMENSION(2) :: file_extensions
1353 CHARACTER(len=20), ALLOCATABLE, DIMENSION(:) :: headers
1354 CHARACTER(len=21) :: prefix
1355 CHARACTER(len=5) :: prefix_format
1356 COMPLEX(kind=dp), ALLOCATABLE, DIMENSION(:) :: omegas_complex, omegas_pade
1357 COMPLEX(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: field_results, field_results_pade, &
1358 pol_results, pol_results_pade, pol_results_pade_spin_total, pol_results_spin_total, &
1359 results, results_pade, results_pade_spin_total, results_spin_total, value_series
1360 INTEGER :: ft_unit, handle, i, idx_omega_zero, &
1361 info_unit, k, k_static, n, n_elems, &
1362 n_pade, nspin
1363 LOGICAL :: do_moments_ft, do_polarizability
1364 REAL(kind=dp) :: damping, t0
1365 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: omegas, omegas_pade_real
1366 REAL(kind=dp), DIMENSION(3) :: delta_vec
1367 TYPE(cp_logger_type), POINTER :: logger
1368 TYPE(section_vals_type), POINTER :: moment_ft_section, pol_section
1369
1370 CALL timeset(routinen, handle)
1371! For results, using spin * direction for first index, e.g. for nspin = 2
1372! results(1,:) = (spin=1 and direction=1,:),
1373! results(5,:) = (spin=2 and direction=2,:)
1374
1375 logger => cp_get_default_logger()
1376
1377 moment_ft_section => section_vals_get_subs_vals(rtp_section, "PRINT%MOMENTS_FT")
1378 pol_section => section_vals_get_subs_vals(rtp_section, "PRINT%POLARIZABILITY")
1379
1380 nspin = SIZE(moments, 1)
1381 n = SIZE(times)
1382 n_elems = SIZE(rtc%print_pol_elements, 1)
1383
1384 info_unit = -1
1385 IF (PRESENT(info_opt)) info_unit = info_opt
1386
1387 ! NOTE : Allows for at most 2 spin species
1388 file_extensions(1) = "_SPIN_A.dat"
1389 file_extensions(2) = "_SPIN_B.dat"
1390
1391 ! Determine whether MOMENTS_FT and/or polarizability needs to be calculated
1392 do_moments_ft = cp_printkey_is_on(logger%iter_info, moment_ft_section)
1393 do_polarizability = cp_printkey_is_on(logger%iter_info, pol_section)
1394 do_polarizability = do_polarizability .AND. (n_elems > 0)
1395
1396 damping = rtc%ft_damping
1397 t0 = rtc%ft_t0
1398
1399 ! Determine field ft if polarizability required
1400 IF (do_polarizability) THEN
1401 ALLOCATE (field_results(3, n))
1402 IF (rtc%apply_delta_pulse) THEN
1403 ! Constant real FT
1404 IF (PRESENT(cell)) THEN
1405 delta_vec(:) = (real(rtc%delta_pulse_direction(1), kind=dp)*cell%h_inv(1, :) + &
1406 REAL(rtc%delta_pulse_direction(2), kind=dp)*cell%h_inv(2, :) + &
1407 REAL(rtc%delta_pulse_direction(3), kind=dp)*cell%h_inv(3, :)) &
1408 *twopi*rtc%delta_pulse_scale
1409 ELSE
1410 delta_vec(:) = real(rtc%delta_pulse_direction(:), kind=dp)*rtc%delta_pulse_scale
1411 END IF
1412 DO k = 1, 3
1413 field_results(k, :) = cmplx(delta_vec(k), 0.0, kind=dp)
1414 END DO
1415 ELSE
1416 ! Do explicit FT
1417 CALL multi_fft(times, fields, field_results, &
1418 damping_opt=damping, t0_opt=t0, subtract_initial_opt=.true.)
1419 END IF
1420 END IF
1421
1422 IF (do_moments_ft .OR. do_polarizability) THEN
1423 ! We need to transform at least the moments
1424 ! NOTE : Might be able to save some memory by only doing FT of actually
1425 ! required moments, but for now, doing FT of all moment directions
1426 ALLOCATE (results(3*nspin, n))
1427 ALLOCATE (omegas(n))
1428 ALLOCATE (value_series(3*nspin, n))
1429 DO i = 1, nspin
1430 DO k = 1, 3
1431 value_series(3*(i - 1) + k, :) = moments(i, k, :)
1432 END DO
1433 END DO
1434 ! TODO : Choose whether the initial subtraction is applied in &FT section?
1435 CALL multi_fft(times, value_series, results, omegas, &
1436 damping_opt=damping, t0_opt=t0, subtract_initial_opt=.true.)
1437 DEALLOCATE (value_series)
1438 DO i = 1, nspin
1439 ! Output to FT file, if needed
1440 ft_unit = cp_print_key_unit_nr(logger, moment_ft_section, extension=file_extensions(i), &
1441 file_form="FORMATTED", file_position="REWIND")
1442 ! Print header
1443 IF (ft_unit > 0) THEN
1444 ALLOCATE (headers(7))
1445 headers(2) = " x,real [at.u.]"
1446 headers(3) = " x,imag [at.u.]"
1447 headers(4) = " y,real [at.u.]"
1448 headers(5) = " y,imag [at.u.]"
1449 headers(6) = " z,real [at.u.]"
1450 headers(7) = " z,imag [at.u.]"
1451 IF (info_unit == ft_unit) THEN
1452 headers(1) = "# Energy [eV]"
1453 prefix = " MOMENTS_FT|"
1454 prefix_format = "(A12)"
1455 CALL print_rt_file(ft_unit, headers, omegas, results(3*(i - 1) + 1:3*(i - 1) + 3, :), &
1456 prefix, prefix_format, evolt)
1457 ELSE
1458 headers(1) = "# omega [at.u.]"
1459 CALL print_rt_file(ft_unit, headers, omegas, results(3*(i - 1) + 1:3*(i - 1) + 3, :))
1460 END IF
1461 DEALLOCATE (headers)
1462 END IF
1463 CALL cp_print_key_finished_output(ft_unit, logger, moment_ft_section)
1464 END DO
1465 ! Spin-summed total moments FT (open shell only; inert for nspin=1)
1466 IF (nspin > 1) THEN
1467 ALLOCATE (results_spin_total(3, n))
1468 results_spin_total(:, :) = (0.0_dp, 0.0_dp)
1469 DO i = 1, nspin
1470 DO k = 1, 3
1471 results_spin_total(k, :) = results_spin_total(k, :) + results(3*(i - 1) + k, :)
1472 END DO
1473 END DO
1474 ft_unit = cp_print_key_unit_nr(logger, moment_ft_section, extension="_SPIN_TOTAL.dat", &
1475 file_form="FORMATTED", file_position="REWIND")
1476 IF (ft_unit > 0) THEN
1477 ALLOCATE (headers(7))
1478 headers(2) = " x,real [at.u.]"
1479 headers(3) = " x,imag [at.u.]"
1480 headers(4) = " y,real [at.u.]"
1481 headers(5) = " y,imag [at.u.]"
1482 headers(6) = " z,real [at.u.]"
1483 headers(7) = " z,imag [at.u.]"
1484 IF (info_unit == ft_unit) THEN
1485 headers(1) = "# Energy [eV]"
1486 prefix = " MOMENTS_FT|"
1487 prefix_format = "(A12)"
1488 CALL print_rt_file(ft_unit, headers, omegas, results_spin_total, &
1489 prefix, prefix_format, evolt)
1490 ELSE
1491 headers(1) = "# omega [at.u.]"
1492 CALL print_rt_file(ft_unit, headers, omegas, results_spin_total)
1493 END IF
1494 DEALLOCATE (headers)
1495 END IF
1496 CALL cp_print_key_finished_output(ft_unit, logger, moment_ft_section)
1497 DEALLOCATE (results_spin_total)
1498 END IF
1499 END IF
1500
1501 IF (rtc%pade_requested .AND. (do_moments_ft .OR. do_polarizability)) THEN
1502 ALLOCATE (omegas_complex(SIZE(omegas)))
1503 omegas_complex(:) = cmplx(omegas(:), 0.0, kind=dp)
1504 n_pade = int((rtc%pade_e_max - rtc%pade_e_min)/rtc%pade_e_step)
1505 ALLOCATE (omegas_pade(n_pade))
1506 ALLOCATE (omegas_pade_real(n_pade))
1507 ! Construct omegas_pade and omegas_complex
1508 DO i = 1, n_pade
1509 omegas_pade_real(i) = (i - 1)*rtc%pade_e_step + rtc%pade_e_min
1510 omegas_pade(i) = cmplx(omegas_pade_real(i), 0.0, kind=dp)
1511 END DO
1512 ALLOCATE (results_pade(nspin*3, n_pade), source=cmplx(0.0, 0.0, kind=dp))
1513 DO i = 1, nspin
1514 DO k = 1, 3
1515 CALL greenx_refine_ft(rtc%pade_fit_e_min, rtc%pade_fit_e_max, omegas_complex, results(3*(i - 1) + k, :), &
1516 omegas_pade, results_pade(3*(i - 1) + k, :))
1517 END DO
1518 ! Print to a file
1519 ft_unit = cp_print_key_unit_nr(logger, moment_ft_section, extension="_PADE"//file_extensions(i), &
1520 file_form="FORMATTED", file_position="REWIND")
1521 IF (ft_unit > 0) THEN
1522 ALLOCATE (headers(7))
1523 headers(2) = " x,real,pade [at.u.]"
1524 headers(3) = " x,imag,pade [at.u.]"
1525 headers(4) = " y,real,pade [at.u.]"
1526 headers(5) = " y,imag,pade [at.u.]"
1527 headers(6) = " z,real,pade [at.u.]"
1528 headers(7) = " z,imag,pade [at.u.]"
1529 IF (info_unit == ft_unit) THEN
1530 headers(1) = "# Energy [eV]"
1531 prefix = " MOMENTS_FT_PADE|"
1532 prefix_format = "(A17)"
1533 CALL print_rt_file(ft_unit, headers, omegas_pade_real, results_pade(3*(i - 1) + 1:3*(i - 1) + 3, :), &
1534 prefix, prefix_format, evolt)
1535 ELSE
1536 headers(1) = "# omega [at.u.]"
1537 CALL print_rt_file(ft_unit, headers, omegas_pade_real, results_pade(3*(i - 1) + 1:3*(i - 1) + 3, :))
1538 END IF
1539 DEALLOCATE (headers)
1540 END IF
1541 END DO
1542 ! Spin-summed total moments-FT Padé (open shell only; inert for nspin=1)
1543 IF (nspin > 1) THEN
1544 ALLOCATE (results_pade_spin_total(3, n_pade))
1545 results_pade_spin_total(:, :) = (0.0_dp, 0.0_dp)
1546 DO i = 1, nspin
1547 DO k = 1, 3
1548 results_pade_spin_total(k, :) = results_pade_spin_total(k, :) + results_pade(3*(i - 1) + k, :)
1549 END DO
1550 END DO
1551 ft_unit = cp_print_key_unit_nr(logger, moment_ft_section, extension="_PADE_SPIN_TOTAL.dat", &
1552 file_form="FORMATTED", file_position="REWIND")
1553 IF (ft_unit > 0) THEN
1554 ALLOCATE (headers(7))
1555 headers(2) = " x,real,pade [at.u.]"
1556 headers(3) = " x,imag,pade [at.u.]"
1557 headers(4) = " y,real,pade [at.u.]"
1558 headers(5) = " y,imag,pade [at.u.]"
1559 headers(6) = " z,real,pade [at.u.]"
1560 headers(7) = " z,imag,pade [at.u.]"
1561 IF (info_unit == ft_unit) THEN
1562 headers(1) = "# Energy [eV]"
1563 prefix = " MOMENTS_FT_PADE|"
1564 prefix_format = "(A17)"
1565 CALL print_rt_file(ft_unit, headers, omegas_pade_real, results_pade_spin_total, &
1566 prefix, prefix_format, evolt)
1567 ELSE
1568 headers(1) = "# omega [at.u.]"
1569 CALL print_rt_file(ft_unit, headers, omegas_pade_real, results_pade_spin_total)
1570 END IF
1571 DEALLOCATE (headers)
1572 END IF
1573 DEALLOCATE (results_pade_spin_total)
1574 END IF
1575 END IF
1576
1577 IF (do_polarizability) THEN
1578 ! get the polarizability elements, as required
1579 ALLOCATE (pol_results(n_elems, n))
1580 DO i = 1, nspin
1581 DO k = 1, n_elems
1582 ! NOTE - field is regularized to small value
1583 pol_results(k, :) = results(3*(i - 1) + &
1584 rtc%print_pol_elements(k, 1), :)/ &
1585 (field_results(rtc%print_pol_elements(k, 2), :) + &
1586 1.0e-10*field_results(rtc%print_pol_elements(k, 2), 2))
1587 END DO
1588 ! Print to the file
1589 ft_unit = cp_print_key_unit_nr(logger, pol_section, extension=file_extensions(i), &
1590 file_form="FORMATTED", file_position="REWIND")
1591 IF (ft_unit > 0) THEN
1592 ALLOCATE (headers(2*n_elems + 1))
1593 DO k = 1, n_elems
1594 WRITE (headers(2*k), "(A16,I2,I2)") "real pol. elem.", &
1595 rtc%print_pol_elements(k, 1), &
1596 rtc%print_pol_elements(k, 2)
1597 WRITE (headers(2*k + 1), "(A16,I2,I2)") "imag pol. elem.", &
1598 rtc%print_pol_elements(k, 1), &
1599 rtc%print_pol_elements(k, 2)
1600 END DO
1601 ! Write header
1602 IF (info_unit == ft_unit) THEN
1603 headers(1) = "# Energy [eV]"
1604 prefix = " POLARIZABILITY|"
1605 prefix_format = "(A16)"
1606 CALL print_rt_file(ft_unit, headers, omegas, pol_results, &
1607 prefix, prefix_format, evolt)
1608 ELSE
1609 headers(1) = "# omega [at.u.]"
1610 CALL print_rt_file(ft_unit, headers, omegas, pol_results)
1611 END IF
1612 DEALLOCATE (headers)
1613 END IF
1614 CALL cp_print_key_finished_output(ft_unit, logger, pol_section)
1615 ! Static polarizability alpha(0): pol_results at the FFT-grid omega
1616 ! closest to zero. Re is alpha(0); Im should be machine-zero (sanity).
1617 IF (info_unit > 0) THEN
1618 idx_omega_zero = minloc(abs(omegas), dim=1)
1619 IF (i == 1) THEN
1620 WRITE (info_unit, '(A,T22,A,T28,A,T36,A,T59,A)') &
1621 " STATIC_POL|", "spin", "element", "Re [a.u.]", "Im [a.u.]"
1622 END IF
1623 DO k_static = 1, n_elems
1624 WRITE (info_unit, '(A,T22,I4,T28,I3,",",I3,T36,ES22.10E3,T59,ES22.10E3)') &
1625 " STATIC_POL|", i, &
1626 rtc%print_pol_elements(k_static, 1), &
1627 rtc%print_pol_elements(k_static, 2), &
1628 REAL(pol_results(k_static, idx_omega_zero), kind=dp), &
1629 aimag(pol_results(k_static, idx_omega_zero))
1630 END DO
1631 END IF
1632 END DO
1633 ! Spin-summed total polarizability (open shell only; inert for nspin=1).
1634 ! Field is spin-independent, so (sum_s moments_s)/field == sum_s (moments_s/field).
1635 IF (nspin > 1) THEN
1636 ALLOCATE (pol_results_spin_total(n_elems, n))
1637 pol_results_spin_total(:, :) = (0.0_dp, 0.0_dp)
1638 DO k = 1, n_elems
1639 DO i = 1, nspin
1640 pol_results_spin_total(k, :) = pol_results_spin_total(k, :) + &
1641 results(3*(i - 1) + rtc%print_pol_elements(k, 1), :)
1642 END DO
1643 pol_results_spin_total(k, :) = pol_results_spin_total(k, :)/ &
1644 (field_results(rtc%print_pol_elements(k, 2), :) + &
1645 1.0e-10*field_results(rtc%print_pol_elements(k, 2), 2))
1646 END DO
1647 ft_unit = cp_print_key_unit_nr(logger, pol_section, extension="_SPIN_TOTAL.dat", &
1648 file_form="FORMATTED", file_position="REWIND")
1649 IF (ft_unit > 0) THEN
1650 ALLOCATE (headers(2*n_elems + 1))
1651 DO k = 1, n_elems
1652 WRITE (headers(2*k), "(A16,I2,I2)") "real pol. elem.", &
1653 rtc%print_pol_elements(k, 1), &
1654 rtc%print_pol_elements(k, 2)
1655 WRITE (headers(2*k + 1), "(A16,I2,I2)") "imag pol. elem.", &
1656 rtc%print_pol_elements(k, 1), &
1657 rtc%print_pol_elements(k, 2)
1658 END DO
1659 IF (info_unit == ft_unit) THEN
1660 headers(1) = "# Energy [eV]"
1661 prefix = " POLARIZABILITY|"
1662 prefix_format = "(A16)"
1663 CALL print_rt_file(ft_unit, headers, omegas, pol_results_spin_total, &
1664 prefix, prefix_format, evolt)
1665 ELSE
1666 headers(1) = "# omega [at.u.]"
1667 CALL print_rt_file(ft_unit, headers, omegas, pol_results_spin_total)
1668 END IF
1669 DEALLOCATE (headers)
1670 END IF
1671 CALL cp_print_key_finished_output(ft_unit, logger, pol_section)
1672 ! Static polarizability total (header row already emitted by the per-spin block)
1673 IF (info_unit > 0) THEN
1674 idx_omega_zero = minloc(abs(omegas), dim=1)
1675 DO k_static = 1, n_elems
1676 WRITE (info_unit, '(A,T22,A,T28,I3,",",I3,T36,ES22.10E3,T59,ES22.10E3)') &
1677 " STATIC_POL|", "TOT", &
1678 rtc%print_pol_elements(k_static, 1), &
1679 rtc%print_pol_elements(k_static, 2), &
1680 REAL(pol_results_spin_total(k_static, idx_omega_zero), kind=dp), &
1681 aimag(pol_results_spin_total(k_static, idx_omega_zero))
1682 END DO
1683 END IF
1684 DEALLOCATE (pol_results_spin_total)
1685 END IF
1686 END IF
1687
1688 ! Padé polarizability
1689 IF (rtc%pade_requested .AND. do_polarizability) THEN
1690 ! Start with the field pade
1691 ALLOCATE (field_results_pade(3, n_pade))
1692 IF (rtc%apply_delta_pulse) THEN
1693 DO k = 1, 3
1694 field_results_pade(k, :) = cmplx(delta_vec(k), 0.0, kind=dp)
1695 END DO
1696 ELSE
1697 DO k = 1, 3
1698 CALL greenx_refine_ft(rtc%pade_fit_e_min, rtc%pade_fit_e_max, &
1699 omegas_complex, field_results(k, :), &
1700 omegas_pade, field_results_pade(k, :))
1701 END DO
1702 END IF
1703 ! Allocate polarisation pade
1704 ALLOCATE (pol_results_pade(n_elems, n_pade))
1705 ! Refine
1706 DO i = 1, nspin
1707 DO k = 1, n_elems
1708 ! NOTE : Regularization to small value
1709 pol_results_pade(k, :) = results_pade(3*(i - 1) + rtc%print_pol_elements(k, 1), :)/( &
1710 field_results_pade(rtc%print_pol_elements(k, 2), :) + &
1711 field_results_pade(rtc%print_pol_elements(k, 2), 2)*1.0e-10_dp)
1712 END DO
1713 ! Print to the file
1714 ft_unit = cp_print_key_unit_nr(logger, pol_section, extension="_PADE"//file_extensions(i), &
1715 file_form="FORMATTED", file_position="REWIND")
1716 IF (ft_unit > 0) THEN
1717 ALLOCATE (headers(2*n_elems + 1))
1718 DO k = 1, n_elems
1719 WRITE (headers(2*k), "(A16,I2,I2)") "re,pade,pol.", &
1720 rtc%print_pol_elements(k, 1), &
1721 rtc%print_pol_elements(k, 2)
1722 WRITE (headers(2*k + 1), "(A16,I2,I2)") "im,pade,pol.", &
1723 rtc%print_pol_elements(k, 1), &
1724 rtc%print_pol_elements(k, 2)
1725 END DO
1726 ! Write header
1727 IF (info_unit == ft_unit) THEN
1728 headers(1) = "# Energy [eV]"
1729 prefix = " POLARIZABILITY_PADE|"
1730 prefix_format = "(A21)"
1731 CALL print_rt_file(ft_unit, headers, omegas_pade_real, pol_results_pade, &
1732 prefix, prefix_format, evolt)
1733 ELSE
1734 headers(1) = "# omega [at.u.]"
1735 CALL print_rt_file(ft_unit, headers, omegas_pade_real, pol_results_pade)
1736 END IF
1737 DEALLOCATE (headers)
1738 END IF
1739 CALL cp_print_key_finished_output(ft_unit, logger, pol_section)
1740 END DO
1741 ! Spin-summed total Padé polarizability (open shell only; inert for nspin=1)
1742 IF (nspin > 1) THEN
1743 ALLOCATE (pol_results_pade_spin_total(n_elems, n_pade))
1744 pol_results_pade_spin_total(:, :) = (0.0_dp, 0.0_dp)
1745 DO k = 1, n_elems
1746 DO i = 1, nspin
1747 pol_results_pade_spin_total(k, :) = pol_results_pade_spin_total(k, :) + &
1748 results_pade(3*(i - 1) + rtc%print_pol_elements(k, 1), :)
1749 END DO
1750 pol_results_pade_spin_total(k, :) = pol_results_pade_spin_total(k, :)/( &
1751 field_results_pade(rtc%print_pol_elements(k, 2), :) + &
1752 field_results_pade(rtc%print_pol_elements(k, 2), 2)*1.0e-10_dp)
1753 END DO
1754 ft_unit = cp_print_key_unit_nr(logger, pol_section, extension="_PADE_SPIN_TOTAL.dat", &
1755 file_form="FORMATTED", file_position="REWIND")
1756 IF (ft_unit > 0) THEN
1757 ALLOCATE (headers(2*n_elems + 1))
1758 DO k = 1, n_elems
1759 WRITE (headers(2*k), "(A16,I2,I2)") "re,pade,pol.", &
1760 rtc%print_pol_elements(k, 1), &
1761 rtc%print_pol_elements(k, 2)
1762 WRITE (headers(2*k + 1), "(A16,I2,I2)") "im,pade,pol.", &
1763 rtc%print_pol_elements(k, 1), &
1764 rtc%print_pol_elements(k, 2)
1765 END DO
1766 IF (info_unit == ft_unit) THEN
1767 headers(1) = "# Energy [eV]"
1768 prefix = " POLARIZABILITY_PADE|"
1769 prefix_format = "(A21)"
1770 CALL print_rt_file(ft_unit, headers, omegas_pade_real, pol_results_pade_spin_total, &
1771 prefix, prefix_format, evolt)
1772 ELSE
1773 headers(1) = "# omega [at.u.]"
1774 CALL print_rt_file(ft_unit, headers, omegas_pade_real, pol_results_pade_spin_total)
1775 END IF
1776 DEALLOCATE (headers)
1777 END IF
1778 CALL cp_print_key_finished_output(ft_unit, logger, pol_section)
1779 DEALLOCATE (pol_results_pade_spin_total)
1780 END IF
1781 DEALLOCATE (field_results_pade)
1782 DEALLOCATE (pol_results_pade)
1783 END IF
1784
1785 IF (rtc%pade_requested .AND. (do_moments_ft .OR. do_polarizability)) THEN
1786 DEALLOCATE (omegas_complex)
1787 DEALLOCATE (omegas_pade)
1788 DEALLOCATE (omegas_pade_real)
1789 DEALLOCATE (results_pade)
1790 END IF
1791
1792 IF (do_polarizability) THEN
1793 DEALLOCATE (pol_results)
1794 DEALLOCATE (field_results)
1795 END IF
1796
1797 IF (do_moments_ft .OR. do_polarizability) THEN
1798 DEALLOCATE (results)
1799 DEALLOCATE (omegas)
1800 END IF
1801
1802 CALL timestop(handle)
1803
1804 END SUBROUTINE print_ft
1805
1806! **************************************************************************************************
1807!> \brief ...
1808!> \param rt_unit ...
1809!> \param headers ...
1810!> \param xvals ...
1811!> \param yvals ...
1812!> \param prefix ...
1813!> \param prefix_format ...
1814!> \param xscale_opt ...
1815!> \param comp_opt ...
1816! **************************************************************************************************
1817 SUBROUTINE print_rt_file(rt_unit, headers, xvals, yvals, prefix, prefix_format, xscale_opt, comp_opt)
1818 INTEGER, INTENT(IN) :: rt_unit
1819 CHARACTER(len=20), DIMENSION(:), INTENT(IN), &
1820 OPTIONAL :: headers
1821 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: xvals
1822 COMPLEX(kind=dp), DIMENSION(:, :), INTENT(IN) :: yvals
1823 CHARACTER(len=21), INTENT(IN), OPTIONAL :: prefix
1824 CHARACTER(len=5), INTENT(IN), OPTIONAL :: prefix_format
1825 REAL(kind=dp), INTENT(IN), OPTIONAL :: xscale_opt
1826 INTEGER, OPTIONAL :: comp_opt
1827
1828 INTEGER :: do_comp, i, j, ncols, nrows
1829 LOGICAL :: do_headers, do_prefix
1830 REAL(kind=dp) :: xscale
1831
1832 do_prefix = .false.
1833 IF (PRESENT(prefix)) THEN
1834 IF (PRESENT(prefix_format)) THEN
1835 do_prefix = .true.
1836 ELSE
1837 cpabort("Printing of prefix with missing format!")
1838 END IF
1839 END IF
1840
1841 xscale = 1.0_dp
1842 IF (PRESENT(xscale_opt)) xscale = xscale_opt
1843
1844 ncols = SIZE(yvals, 1)
1845 nrows = SIZE(yvals, 2)
1846
1847 ! Check whether printing complex data (default) or just a component
1848 do_comp = rt_file_comp_both
1849 IF (PRESENT(comp_opt)) do_comp = comp_opt
1850
1851 do_headers = PRESENT(headers)
1852 ! Check whether enough headers for yvals and xvals is present
1853 IF (do_headers) THEN
1854 IF ((do_comp /= rt_file_comp_both .AND. SIZE(headers) < ncols + 1) &
1855 .OR. (do_comp == rt_file_comp_both .AND. SIZE(headers) < 2*ncols + 1)) THEN
1856 cpabort("Not enought headers to print the file!")
1857 END IF
1858 END IF
1859
1860 IF (SIZE(xvals) < nrows) THEN
1861 cpabort("Not enough xvals to print all yvals!")
1862 END IF
1863
1864 IF (rt_unit > 0) THEN
1865 ! Print headers
1866 IF (do_headers) THEN
1867 ! If prefix is present, write prefix
1868 IF (do_prefix) THEN
1869 WRITE (rt_unit, prefix_format, advance="no") prefix
1870 END IF
1871 WRITE (rt_unit, "(A20)", advance="no") headers(1)
1872 ! Print the rest of the headers
1873 SELECT CASE (do_comp)
1874 CASE (rt_file_comp_both)
1875 ! Complex case
1876 DO j = 1, 2*ncols - 1
1877 WRITE (rt_unit, "(A20)", advance="no") headers(j + 1)
1878 END DO
1879 WRITE (rt_unit, "(A20)") headers(2*ncols + 1)
1880 CASE DEFAULT
1881 ! For other cases, just one component is printed
1882 DO j = 1, ncols - 1
1883 WRITE (rt_unit, "(A20)", advance="no") headers(j + 1)
1884 END DO
1885 WRITE (rt_unit, "(A20)") headers(ncols + 1)
1886 END SELECT
1887 END IF
1888 ! Done with the headers, print actual data
1889 DO i = 1, nrows
1890 ! If prefix is present, write prefix
1891 IF (do_prefix) THEN
1892 WRITE (rt_unit, prefix_format, advance="no") prefix
1893 END IF
1894 WRITE (rt_unit, "(E20.8E3)", advance="no") xvals(i)*xscale
1895 DO j = 1, ncols - 1
1896 SELECT CASE (do_comp)
1897 CASE (rt_file_comp_real)
1898 WRITE (rt_unit, "(E20.8E3)", advance="no") &
1899 REAL(yvals(j, i))
1900 CASE (rt_file_comp_imag)
1901 WRITE (rt_unit, "(E20.8E3)", advance="no") &
1902 aimag(yvals(j, i))
1903 CASE DEFAULT
1904 ! Print both components
1905 WRITE (rt_unit, "(E20.8E3,E20.8E3)", advance="no") &
1906 REAL(yvals(j, i)), aimag(yvals(j, i))
1907 END SELECT
1908 END DO
1909 ! Print the final column(s)
1910 SELECT CASE (do_comp)
1911 CASE (rt_file_comp_real)
1912 WRITE (rt_unit, "(E20.8E3)") real(yvals(j, i))
1913 CASE (rt_file_comp_imag)
1914 WRITE (rt_unit, "(E20.8E3)") aimag(yvals(j, i))
1915 CASE DEFAULT
1916 ! Print both components
1917 WRITE (rt_unit, "(E20.8E3,E20.8E3)") &
1918 REAL(yvals(j, i)), aimag(yvals(j, i))
1919 END SELECT
1920 END DO
1921 END IF
1922 END SUBROUTINE print_rt_file
1923
1924END MODULE rt_propagation_output
cp_dbcsr_api::dbcsr_create
Definition cp_dbcsr_api.F:194
cp_dbcsr_operations::dbcsr_allocate_matrix_set
Definition cp_dbcsr_operations.F:77
cp_dbcsr_operations::dbcsr_deallocate_matrix_set
Definition cp_dbcsr_operations.F:85
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:88
cp_log_handling::cp_to_string
Definition cp_log_handling.F:90
kahan_sum::accurate_sum
Definition kahan_sum.F:41
parallel_gemm_api::parallel_gemm
Definition parallel_gemm_api.F:39
pw_methods::pw_zero
Definition pw_methods.F:83
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
cell_types
Handles all functions related to the CELL.
Definition cell_types.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_scale
subroutine, public dbcsr_scale(matrix, alpha_scalar)
...
Definition cp_dbcsr_api.F:1188
cp_dbcsr_api::dbcsr_deallocate_matrix
subroutine, public dbcsr_deallocate_matrix(matrix)
...
Definition cp_dbcsr_api.F:235
cp_dbcsr_api::dbcsr_iterator_blocks_left
logical function, public dbcsr_iterator_blocks_left(iterator)
...
Definition cp_dbcsr_api.F:944
cp_dbcsr_api::dbcsr_iterator_stop
subroutine, public dbcsr_iterator_stop(iterator)
...
Definition cp_dbcsr_api.F:1050
cp_dbcsr_api::dbcsr_desymmetrize
subroutine, public dbcsr_desymmetrize(matrix_a, matrix_b)
...
Definition cp_dbcsr_api.F:525
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_multiply
subroutine, public dbcsr_multiply(transa, transb, alpha, matrix_a, matrix_b, beta, matrix_c, first_row, last_row, first_column, last_column, first_k, last_k, retain_sparsity, filter_eps, flop)
...
Definition cp_dbcsr_api.F:1096
cp_dbcsr_api::dbcsr_get_info
subroutine, public dbcsr_get_info(matrix, nblkrows_total, nblkcols_total, nfullrows_total, nfullcols_total, nblkrows_local, nblkcols_local, nfullrows_local, nfullcols_local, my_prow, my_pcol, local_rows, local_cols, proc_row_dist, proc_col_dist, row_blk_size, col_blk_size, row_blk_offset, col_blk_offset, distribution, name, matrix_type, group)
...
Definition cp_dbcsr_api.F:808
cp_dbcsr_api::dbcsr_init_p
subroutine, public dbcsr_init_p(matrix)
...
Definition cp_dbcsr_api.F:207
cp_dbcsr_api::dbcsr_iterator_next_block
subroutine, public dbcsr_iterator_next_block(iterator, row, column, block, block_number_argument_has_been_removed, row_size, col_size, row_offset, col_offset, transposed)
...
Definition cp_dbcsr_api.F:971
cp_dbcsr_api::dbcsr_filter
subroutine, public dbcsr_filter(matrix, eps)
...
Definition cp_dbcsr_api.F:658
cp_dbcsr_api::dbcsr_binary_write
subroutine, public dbcsr_binary_write(matrix, filepath)
...
Definition cp_dbcsr_api.F:290
cp_dbcsr_api::dbcsr_get_occupation
real(kind=dp) function, public dbcsr_get_occupation(matrix)
...
Definition cp_dbcsr_api.F:893
cp_dbcsr_api::dbcsr_iterator_start
subroutine, public dbcsr_iterator_start(iterator, matrix, shared, dynamic, dynamic_byrows)
...
Definition cp_dbcsr_api.F:1012
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_checksum
real(kind=dp) function, public dbcsr_checksum(matrix, pos)
Calculates the checksum of a DBCSR matrix.
Definition cp_dbcsr_contrib.F:628
cp_dbcsr_contrib::dbcsr_trace
subroutine, public dbcsr_trace(matrix, trace)
Computes the trace of the given matrix, also known as the sum of its diagonal elements.
Definition cp_dbcsr_contrib.F:412
cp_dbcsr_cp2k_link
Routines that link DBCSR and CP2K concepts together.
Definition cp_dbcsr_cp2k_link.F:14
cp_dbcsr_cp2k_link::cp_dbcsr_alloc_block_from_nbl
subroutine, public cp_dbcsr_alloc_block_from_nbl(matrix, sab_orb, desymmetrize)
allocate the blocks of a dbcsr based on the neighbor list
Definition cp_dbcsr_cp2k_link.F:445
cp_dbcsr_operations
DBCSR operations in CP2K.
Definition cp_dbcsr_operations.F:18
cp_dbcsr_operations::cp_dbcsr_sm_fm_multiply
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
Definition cp_dbcsr_operations.F:548
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_and_add
subroutine, public cp_fm_scale_and_add(alpha, matrix_a, beta, matrix_b)
calc A <- alpha*A + beta*B optimized for alpha == 1.0 (just add beta*B) and beta == 0....
Definition cp_fm_basic_linalg.F:170
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_double
subroutine, public cp_fm_struct_double(fmstruct, struct, context, col, row)
creates a struct with twice the number of blocks on each core. If matrix A has to be multiplied with ...
Definition cp_fm_struct.F:527
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_set_all
subroutine, public cp_fm_set_all(matrix, alpha, beta)
set all elements of a matrix to the same value, and optionally the diagonal to a different one
Definition cp_fm_types.F:534
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_unit_nr
recursive integer function, public cp_logger_get_default_unit_nr(logger, local, skip_not_ionode)
asks the default unit number of the given logger. try to use cp_logger_get_unit_nr
Definition cp_log_handling.F:576
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
cp_output_handling
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
Definition cp_output_handling.F:25
cp_output_handling::cp_iter_string
character(len=default_string_length) function, public cp_iter_string(iter_info, print_key, for_file)
returns the iteration string, a string that is useful to create unique filenames (once you trim it)
Definition cp_output_handling.F:517
cp_output_handling::cp_print_key_unit_nr
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)
...
Definition cp_output_handling.F:853
cp_output_handling::cp_print_key_finished_output
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,...
Definition cp_output_handling.F:1073
cp_output_handling::cp_p_file
integer, parameter, public cp_p_file
Definition cp_output_handling.F:103
cp_output_handling::cp_printkey_is_on
logical function, public cp_printkey_is_on(iteration_info, print_key)
returns true if the printlevel activates this printkey does not look if this iteration it should be p...
Definition cp_output_handling.F:429
cp_output_handling::cp_print_key_should_output
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...
Definition cp_output_handling.F:345
cp_realspace_grid_cube
A wrapper around pw_to_cube() which accepts particle_list_type.
Definition cp_realspace_grid_cube.F:12
cp_realspace_grid_cube::cp_pw_to_cube
subroutine, public cp_pw_to_cube(pw, unit_nr, title, particles, zeff, stride, max_file_size_mb, zero_tails, silent, mpi_io)
...
Definition cp_realspace_grid_cube.F:48
efield_utils
all routins needed for a nonperiodic electric field
Definition efield_utils.F:12
efield_utils::make_field
subroutine, public make_field(dft_control, field, sim_step, sim_time)
computes the amplitude of the efield within a given envelop
Definition efield_utils.F:120
greenx_interface
Interface to the Greenx library.
Definition greenx_interface.F:13
greenx_interface::greenx_refine_ft
subroutine, public greenx_refine_ft(fit_e_min, fit_e_max, x_fit, y_fit, x_eval, y_eval, n_pade_opt)
Refines the FT grid using Padé approximants.
Definition greenx_interface.F:205
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:18
input_constants::ehrenfest
integer, parameter, public ehrenfest
Definition input_constants.F:131
input_constants::real_time_propagation
integer, parameter, public real_time_propagation
Definition input_constants.F:131
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
input_section_types::section_get_ivals
integer function, dimension(:), pointer, public section_get_ivals(section_vals, keyword_name)
...
Definition input_section_types.F:1001
input_section_types::section_vals_get_subs_vals
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
Definition input_section_types.F:735
kahan_sum
sums arrays of real/complex numbers with much reduced round-off as compared to a naive implementation...
Definition kahan_sum.F:29
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kinds::default_path_length
integer, parameter, public default_path_length
Definition kinds.F:58
machine
Machine interface based on Fortran 2003 and POSIX.
Definition machine.F:17
machine::m_flush
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Definition machine.F:124
mathconstants
Definition of mathematical constants and functions.
Definition mathconstants.F:16
mathconstants::one
real(kind=dp), parameter, public one
Definition mathconstants.F:123
mathconstants::twopi
real(kind=dp), parameter, public twopi
Definition mathconstants.F:112
mathconstants::zero
real(kind=dp), parameter, public zero
Definition mathconstants.F:123
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:27
parallel_gemm_api
basic linear algebra operations for full matrixes
Definition parallel_gemm_api.F:14
particle_list_types
represent a simple array based list of the given type
Definition particle_list_types.F:15
particle_methods
Define methods related to particle_type.
Definition particle_methods.F:14
particle_methods::get_particle_set
subroutine, public get_particle_set(particle_set, qs_kind_set, first_sgf, last_sgf, nsgf, nmao, basis, ncgf)
Get the components of a particle set.
Definition particle_methods.F:123
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
physcon
Definition of physical constants:
Definition physcon.F:68
physcon::femtoseconds
real(kind=dp), parameter, public femtoseconds
Definition physcon.F:153
physcon::evolt
real(kind=dp), parameter, public evolt
Definition physcon.F:183
pw_env_types
container for various plainwaves related things
Definition pw_env_types.F:14
pw_env_types::pw_env_get
subroutine, public pw_env_get(pw_env, pw_pools, cube_info, gridlevel_info, auxbas_pw_pool, auxbas_grid, auxbas_rs_desc, auxbas_rs_grid, rs_descs, rs_grids, xc_pw_pool, vdw_pw_pool, poisson_env, interp_section)
returns the various attributes of the pw env
Definition pw_env_types.F:113
pw_methods
Definition pw_methods.F:25
pw_pool_types
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Definition pw_pool_types.F:24
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:542
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_kind_types::get_qs_kind_set
subroutine, public get_qs_kind_set(qs_kind_set, all_potential_present, tnadd_potential_present, gth_potential_present, sgp_potential_present, paw_atom_present, dft_plus_u_atom_present, maxcgf, maxsgf, maxco, maxco_proj, maxgtops, maxlgto, maxlprj, maxnset, maxsgf_set, ncgf, npgf, nset, nsgf, nshell, maxpol, maxlppl, maxlppnl, maxppnl, nelectron, maxder, max_ngrid_rad, max_sph_harm, maxg_iso_not0, lmax_rho0, basis_rcut, do_mtlr_present, basis_type, total_zeff_corr, npgf_seg, cneo_potential_present, nkind_q, natom_q)
Get attributes of an atomic kind set.
Definition qs_kind_types.F:984
qs_linres_current
given the response wavefunctions obtained by the application of the (rxp), p, and ((dk-dl)xp) operato...
Definition qs_linres_current.F:17
qs_linres_current::calculate_jrho_resp
subroutine, public calculate_jrho_resp(mat_d0, mat_jp, mat_jp_rii, mat_jp_riii, ib, idir, current_rs, current_gs, qs_env, current_env, soft_valid, retain_rsgrid)
Calculation of the idir component of the response current density in the presence of a constant magne...
Definition qs_linres_current.F:578
qs_linres_types
Type definitiona for linear response calculations.
Definition qs_linres_types.F:12
qs_mo_io
Definition and initialisation of the mo data type.
Definition qs_mo_io.F:21
qs_mo_io::write_rt_mos_to_restart
subroutine, public write_rt_mos_to_restart(mo_array, rt_mos, particle_set, dft_section, qs_kind_set)
...
Definition qs_mo_io.F:244
qs_moments
Calculates the moment integrals <a|r^m|b> and <a|r x d/dr|b>.
Definition qs_moments.F:14
qs_moments::build_local_moment_matrix
subroutine, public build_local_moment_matrix(qs_env, moments, nmoments, ref_point, ref_points, basis_type, all_images, minimum_image, neighbor_image, first_component)
...
Definition qs_moments.F:160
qs_neighbor_list_types
Define the neighbor list data types and the corresponding functionality.
Definition qs_neighbor_list_types.F:17
qs_operators_ao
Definition qs_operators_ao.F:13
qs_operators_ao::build_lin_mom_matrix
subroutine, public build_lin_mom_matrix(qs_env, matrix, minimum_image)
Calculation of the linear momentum matrix <mu|∂|nu> over Cartesian Gaussian functions.
Definition qs_operators_ao.F:275
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
qs_scf_post_gpw
Does all kind of post scf calculations for GPW/GAPW.
Definition qs_scf_post_gpw.F:15
qs_scf_post_gpw::write_mo_free_results
subroutine, public write_mo_free_results(qs_env)
Write QS results always available (if switched on through the print_keys) Can be called from ls_scf.
Definition qs_scf_post_gpw.F:2599
qs_scf_post_gpw::qs_scf_post_moments
subroutine, public qs_scf_post_moments(input, logger, qs_env, output_unit)
Computes and prints electric moments.
Definition qs_scf_post_gpw.F:1671
qs_scf_post_gpw::write_mo_dependent_results
subroutine, public write_mo_dependent_results(qs_env, scf_env)
Write QS results available if MO's are present (if switched on through the print_keys) Writes only MO...
Definition qs_scf_post_gpw.F:2264
qs_scf_post_tb
Does all kind of post scf calculations for DFTB.
Definition qs_scf_post_tb.F:15
qs_scf_post_tb::scf_post_calculation_tb
subroutine, public scf_post_calculation_tb(qs_env, tb_type, no_mos)
collects possible post - scf calculations and prints info / computes properties.
Definition qs_scf_post_tb.F:172
qs_scf_types
module that contains the definitions of the scf types
Definition qs_scf_types.F:14
qs_subsys_types
types that represent a quickstep subsys
Definition qs_subsys_types.F:12
qs_subsys_types::qs_subsys_get
subroutine, public qs_subsys_get(subsys, atomic_kinds, atomic_kind_set, particles, particle_set, local_particles, molecules, molecule_set, molecule_kinds, molecule_kind_set, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, natom, nparticle, ncore, nshell, nkind, atprop, virial, results, cell, cell_ref, use_ref_cell, energy, force, qs_kind_set, cp_subsys, nelectron_total, nelectron_spin)
...
Definition qs_subsys_types.F:137
rt_projection_mo_utils
Function related to MO projection in RTP calculations.
Definition rt_projection_mo_utils.F:12
rt_projection_mo_utils::compute_and_write_proj_mo
subroutine, public compute_and_write_proj_mo(qs_env, mos_new, proj_mo, n_proj)
Compute the projection of the current MO coefficients on reference ones and write the results.
Definition rt_projection_mo_utils.F:352
rt_propagation_ft
Separation of Fourier transform utilities into separate file.
Definition rt_propagation_ft.F:12
rt_propagation_ft::multi_fft
subroutine, public multi_fft(time_series, value_series, result_series, omega_series, damping_opt, t0_opt, subtract_initial_opt)
Calculates the Fourier transform - couples to FFT libraries in CP2K, if available.
Definition rt_propagation_ft.F:99
rt_propagation_output
Routine for the real time propagation output.
Definition rt_propagation_output.F:13
rt_propagation_output::print_ft
subroutine, public print_ft(rtp_section, moments, times, fields, rtc, info_opt, cell)
Calculate and print the Fourier transforms + polarizabilites from moment trace.
Definition rt_propagation_output.F:1342
rt_propagation_output::report_density_occupation
subroutine, public report_density_occupation(filter_eps, rho)
Reports the sparsity pattern of the complex density matrix.
Definition rt_propagation_output.F:780
rt_propagation_output::print_rt_file
subroutine, public print_rt_file(rt_unit, headers, xvals, yvals, prefix, prefix_format, xscale_opt, comp_opt)
...
Definition rt_propagation_output.F:1818
rt_propagation_output::rt_prop_output
subroutine, public rt_prop_output(qs_env, run_type, delta_iter, used_time)
...
Definition rt_propagation_output.F:136
rt_propagation_output::print_moments
subroutine, public print_moments(moments_section, info_unit, moments, time, imag_opt, append_opt)
Print the dipole moments into a file.
Definition rt_propagation_output.F:1173
rt_propagation_output::rt_convergence
subroutine, public rt_convergence(rtp, matrix_s, delta_mos, delta_eps)
computes the convergence criterion for RTP and EMD
Definition rt_propagation_output.F:568
rt_propagation_output::rt_file_comp_both
integer, parameter, public rt_file_comp_both
Definition rt_propagation_output.F:122
rt_propagation_output::calc_local_moment
subroutine, public calc_local_moment(moment_matrices, density_matrices, work, moment, imag_opt)
Calculate the values of real/imaginary parts of moments in all directions.
Definition rt_propagation_output.F:1297
rt_propagation_output::rt_file_comp_imag
integer, parameter, public rt_file_comp_imag
Definition rt_propagation_output.F:122
rt_propagation_output::rt_file_comp_real
integer, parameter, public rt_file_comp_real
Definition rt_propagation_output.F:122
rt_propagation_output::rt_convergence_density
subroutine, public rt_convergence_density(rtp, delta_p, delta_eps)
computes the convergence criterion for RTP and EMD based on the density matrix
Definition rt_propagation_output.F:679
rt_propagation_types
Types and set_get for real time propagation depending on runtype and diagonalization method different...
Definition rt_propagation_types.F:21
rt_propagation_types::get_rtp
subroutine, public get_rtp(rtp, exp_h_old, exp_h_new, h_last_iter, rho_old, rho_next, rho_new, mos, mos_new, mos_old, mos_next, s_inv, s_half, s_minus_half, b_mat, c_mat, propagator_matrix, mixing, mixing_factor, s_der, dt, nsteps, sinvh, sinvh_imag, sinvb, admm_mos)
...
Definition rt_propagation_types.F:416
rt_propagation_utils
Routines needed for EMD.
Definition rt_propagation_utils.F:13
rt_propagation_utils::write_rtp_mos_to_output_unit
subroutine, public write_rtp_mos_to_output_unit(qs_env, rtp)
...
Definition rt_propagation_utils.F:635
rt_propagation_utils::write_rtp_mo_cubes
subroutine, public write_rtp_mo_cubes(qs_env, rtp)
Write the time dependent amplitude of the MOs in real grid. Very close to qs_scf_post_gpw/qs_scf_post...
Definition rt_propagation_utils.F:735
rt_propagation_velocity_gauge
Routines to perform the RTP in the velocity gauge.
Definition rt_propagation_velocity_gauge.F:12
rt_propagation_velocity_gauge::velocity_gauge_nl_commutator
subroutine, public velocity_gauge_nl_commutator(qs_env, matrix_cre, matrix_cim, vec_pot)
Calculate the commutator [Vnl~, r_d] (d=x,y,z) of the gauge-transformed nonlocal pseudopotential with...
Definition rt_propagation_velocity_gauge.F:540
atomic_kind_types::atomic_kind_type
Provides all information about an atomic kind.
Definition atomic_kind_types.F:45
cell_types::cell_type
Type defining parameters related to the simulation cell.
Definition cell_types.F:60
cp_control_types::dft_control_type
Definition cp_control_types.F:758
cp_control_types::rtp_control_type
Definition cp_control_types.F:94
cp_dbcsr_api::dbcsr_distribution_type
Definition cp_dbcsr_api.F:182
cp_dbcsr_api::dbcsr_iterator_type
Definition cp_dbcsr_api.F:188
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
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:127
message_passing::mp_comm_type
Definition message_passing.F:178
particle_list_types::particle_list_type
represent a list of objects
Definition particle_list_types.F:46
particle_types::particle_type
Definition particle_types.F:35
pw_env_types::pw_env_type
contained for different pw related things
Definition pw_env_types.F:70
pw_pool_types::pw_pool_type
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Definition pw_pool_types.F:83
pw_types::pw_c1d_gs_type
Definition pw_types.F:127
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:224
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:190
qs_linres_types::current_env_type
Definition qs_linres_types.F:94
qs_neighbor_list_types::neighbor_list_set_p_type
Definition qs_neighbor_list_types.F:67
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
qs_scf_types::qs_scf_env_type
Definition qs_scf_types.F:99
qs_subsys_types::qs_subsys_type
Definition qs_subsys_types.F:56
rt_propagation_types::rt_prop_type
Definition rt_propagation_types.F:76