(git:691081d)
Loading...
Searching...
No Matches
qs_tddfpt2_methods.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8MODULE qs_tddfpt2_methods
9 USE admm_methods, ONLY: admm_fit_mo_coeffs
10 USE admm_types, ONLY: admm_type,&
11 get_admm_env
12 USE atomic_kind_types, ONLY: atomic_kind_type
13 USE bibliography, ONLY: grimme2013,&
14 grimme2016,&
15 hernandez2025,&
16 iannuzzi2005,&
17 cite_reference
18 USE cell_types, ONLY: cell_type
19 USE cp_blacs_env, ONLY: cp_blacs_env_type
20 USE cp_control_types, ONLY: dft_control_type,&
21 rixs_control_type,&
22 tddfpt2_control_type
23 USE cp_dbcsr_api, ONLY: dbcsr_create,&
24 dbcsr_deallocate_matrix,&
25 dbcsr_p_type,&
26 dbcsr_set,&
27 dbcsr_type,&
28 dbcsr_type_antisymmetric,&
29 dbcsr_type_symmetric
30 USE cp_dbcsr_cp2k_link, ONLY: cp_dbcsr_alloc_block_from_nbl
31 USE cp_dbcsr_operations, ONLY: cp_dbcsr_sm_fm_multiply,&
32 dbcsr_deallocate_matrix_set
33 USE cp_fm_pool_types, ONLY: fm_pool_create_fm
34 USE cp_fm_struct, ONLY: cp_fm_struct_create,&
35 cp_fm_struct_release,&
36 cp_fm_struct_type
37 USE cp_fm_types, ONLY: cp_fm_copy_general,&
38 cp_fm_create,&
39 cp_fm_get_element,&
40 cp_fm_get_info,&
41 cp_fm_release,&
42 cp_fm_to_fm,&
43 cp_fm_to_fm_submat,&
44 cp_fm_type
45 USE cp_log_handling, ONLY: cp_get_default_logger,&
46 cp_logger_get_default_io_unit,&
47 cp_logger_type
48 USE cp_output_handling, ONLY: cp_add_iter_level,&
49 cp_iterate,&
50 cp_print_key_finished_output,&
51 cp_print_key_unit_nr,&
52 cp_rm_iter_level
53 USE exstates_types, ONLY: excited_energy_type
54 USE header, ONLY: tddfpt_header,&
55 tddfpt_soc_header
56 USE hfx_admm_utils, ONLY: aux_admm_init
57 USE hfx_types, ONLY: compare_hfx_sections,&
58 hfx_create
59 USE input_constants, ONLY: &
60 do_admm_aux_exch_func_none, do_admm_basis_projection, do_admm_exch_scaling_none, &
61 do_admm_purify_none, do_potential_truncated, no_sf_tddfpt, oe_none, &
62 tddfpt_dipole_scf_moment, tddfpt_dipole_velocity, tddfpt_kernel_full, tddfpt_kernel_none, &
63 tddfpt_kernel_stda, tddfpt_sf_col, tddfpt_sf_noncol
64 USE input_section_types, ONLY: section_vals_get,&
65 section_vals_get_subs_vals,&
66 section_vals_type,&
67 section_vals_val_get,&
68 section_vals_val_set
69 USE kinds, ONLY: dp
70 USE kpoint_methods, ONLY: rskp_transform
71 USE kpoint_types, ONLY: get_kpoint_info,&
72 kpoint_env_p_type,&
73 kpoint_env_type,&
74 kpoint_type
75 USE lri_environment_methods, ONLY: lri_print_stat
76 USE lri_environment_types, ONLY: lri_density_release,&
77 lri_env_release
78 USE machine, ONLY: m_flush
79 USE message_passing, ONLY: mp_para_env_type
80 USE min_basis_set, ONLY: create_minbas_set
81 USE parallel_gemm_api, ONLY: parallel_gemm
82 USE particle_types, ONLY: particle_type
83 USE physcon, ONLY: evolt
84 USE qs_environment_types, ONLY: get_qs_env,&
85 qs_environment_type
86 USE qs_kernel_methods, ONLY: create_fxc_kernel,&
87 create_kernel_env
88 USE qs_kernel_types, ONLY: full_kernel_env_type,&
89 kernel_env_type,&
90 release_kernel_env
91 USE qs_kind_types, ONLY: qs_kind_type
92 USE qs_ks_types, ONLY: qs_ks_env_type
93 USE qs_mo_types, ONLY: get_mo_set,&
94 mo_set_type
95 USE qs_moments, ONLY: qs_moment_kpoints_scf_mos
96 USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type
97 USE qs_overlap, ONLY: build_overlap_matrix
98 USE qs_rho_types, ONLY: qs_rho_type
99 USE qs_scf_methods, ONLY: eigensolver
100 USE qs_scf_types, ONLY: qs_scf_env_type
101 USE qs_tddfpt2_assign, ONLY: assign_state
102 USE qs_tddfpt2_densities, ONLY: tddfpt_construct_aux_fit_density,&
103 tddfpt_construct_ground_state_orb_density
104 USE qs_tddfpt2_eigensolver, ONLY: tddfpt_davidson_solver,&
105 tddfpt_orthogonalize_psi1_psi0,&
106 tddfpt_orthonormalize_psi1_psi1
107 USE qs_tddfpt2_forces, ONLY: tddfpt_forces_main
108 USE qs_tddfpt2_fprint, ONLY: tddfpt_print_forces
109 USE qs_tddfpt2_lri_utils, ONLY: tddfpt2_lri_init
110 USE qs_tddfpt2_properties, ONLY: tddfpt_dipole_operator,&
111 tddfpt_print_excitation_analysis,&
112 tddfpt_print_exciton_descriptors,&
113 tddfpt_print_nto_analysis,&
114 tddfpt_print_summary
115 USE qs_tddfpt2_restart, ONLY: tddfpt_read_restart,&
116 tddfpt_write_newtonx_output,&
117 tddfpt_write_restart
118 USE qs_tddfpt2_smearing_methods, ONLY: tddfpt_smeared_occupation
119 USE qs_tddfpt2_soc, ONLY: tddfpt_soc
120 USE qs_tddfpt2_stda_types, ONLY: allocate_stda_env,&
121 deallocate_stda_env,&
122 stda_env_type,&
123 stda_init_param
124 USE qs_tddfpt2_stda_utils, ONLY: get_lowdin_mo_coefficients,&
125 stda_init_matrices
126 USE qs_tddfpt2_subgroups, ONLY: tddfpt_sub_env_init,&
127 tddfpt_sub_env_release,&
128 tddfpt_subgroup_env_type
129 USE qs_tddfpt2_types, ONLY: hfxsr_create_work_matrices,&
130 stda_create_work_matrices,&
131 tddfpt_create_work_matrices,&
132 tddfpt_ground_state_mos,&
133 tddfpt_release_work_matrices,&
134 tddfpt_work_matrices
135 USE qs_tddfpt2_utils, ONLY: tddfpt_guess_vectors,&
136 tddfpt_init_mos,&
137 tddfpt_oecorr,&
138 tddfpt_release_ground_state_mos
139 USE rixs_types, ONLY: rixs_env_type,&
140 tddfpt2_valence_type
141 USE string_utilities, ONLY: integer_to_string
142 USE util, ONLY: sort
143 USE xc_write_output, ONLY: xc_write
144#include "./base/base_uses.f90"
145
146 IMPLICIT NONE
147
148 PRIVATE
149
150 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_tddfpt2_methods'
151
152 LOGICAL, PARAMETER, PRIVATE :: debug_this_module = .false.
153 ! number of first derivative components (3: d/dx, d/dy, d/dz)
154 INTEGER, PARAMETER, PRIVATE :: nderivs = 3
155 INTEGER, PARAMETER, PRIVATE :: maxspins = 2
156
157 PUBLIC :: tddfpt, tddfpt_energies, tddfpt_input
158
159! **************************************************************************************************
160
161CONTAINS
162
163! **************************************************************************************************
164!> \brief Perform TDDFPT calculation. If calc_forces then it also builds the response vector for the
165!> Z-vector method and calculates some contributions to the force
166!> \param qs_env Quickstep environment
167!> \param calc_forces ...
168!> \param rixs_env ...
169!> \par History
170!> * 05.2016 created [Sergey Chulkov]
171!> * 06.2016 refactored to be used with Davidson eigensolver [Sergey Chulkov]
172!> * 03.2017 cleaned and refactored [Sergey Chulkov]
173!> \note Based on the subroutines tddfpt_env_init(), and tddfpt_env_deallocate().
174! **************************************************************************************************
175 SUBROUTINE tddfpt(qs_env, calc_forces, rixs_env)
176 TYPE(qs_environment_type), POINTER :: qs_env
177 LOGICAL, INTENT(IN) :: calc_forces
178 TYPE(rixs_env_type), OPTIONAL, POINTER :: rixs_env
179
180 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt'
181
182 INTEGER :: handle, ispin, istate, log_unit, mult, &
183 my_state, nao, nocc, nspins, &
184 nstate_max, nstates, nvirt, old_state
185 INTEGER, DIMENSION(maxspins) :: nactive
186 LOGICAL :: do_admm, do_exck, do_hfx, do_hfxlr, &
187 do_hfxsr, do_kpoints, do_rixs, do_sf, &
188 do_soc, lmult_tmp, state_change
189 REAL(kind=dp) :: gsmin, gsval, xsval
190 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: evals, ostrength
191 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
192 TYPE(cell_type), POINTER :: cell
193 TYPE(cp_blacs_env_type), POINTER :: blacs_env
194 TYPE(cp_fm_struct_type), POINTER :: matrix_struct
195 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:) :: my_active, my_mos
196 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: dipole_op_mos_occ, evects, s_evects
197 TYPE(cp_logger_type), POINTER :: logger
198 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_ks_oep, matrix_s, &
199 matrix_s_aux_fit, &
200 matrix_s_aux_fit_vs_orb
201 TYPE(dft_control_type), POINTER :: dft_control
202 TYPE(excited_energy_type), POINTER :: ex_env
203 TYPE(full_kernel_env_type), TARGET :: full_kernel_env, kernel_env_admm_aux
204 TYPE(kernel_env_type) :: kernel_env
205 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos, mos_aux_fit
206 TYPE(mp_para_env_type), POINTER :: para_env
207 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
208 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
209 TYPE(qs_scf_env_type), POINTER :: scf_env
210 TYPE(rixs_control_type), POINTER :: rixs_control
211 TYPE(section_vals_type), POINTER :: hfxsr_section, kernel_section, &
212 lri_section, soc_section, &
213 tddfpt_print_section, tddfpt_section, &
214 xc_section
215 TYPE(stda_env_type), TARGET :: stda_kernel
216 TYPE(tddfpt2_control_type), POINTER :: tddfpt_control
217 TYPE(tddfpt2_valence_type), POINTER :: valence_state
218 TYPE(tddfpt_ground_state_mos), DIMENSION(:), &
219 POINTER :: gs_mos
220 TYPE(tddfpt_subgroup_env_type) :: sub_env
221 TYPE(tddfpt_work_matrices) :: work_matrices
222
223 CALL timeset(routinen, handle)
224
225 NULLIFY (logger)
226 logger => cp_get_default_logger()
227
228 NULLIFY (tddfpt_section, tddfpt_control)
229
230 CALL get_qs_env(qs_env, &
231 dft_control=dft_control, &
232 do_rixs=do_rixs)
233 do_kpoints = dft_control%nimages > 1
234
235 IF (do_rixs) THEN
236 tddfpt_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%RIXS%TDDFPT")
237 NULLIFY (rixs_control, valence_state)
238 rixs_control => dft_control%rixs_control
239 tddfpt_control => rixs_control%tddfpt2_control
240 valence_state => rixs_env%valence_state
241 ELSE
242 tddfpt_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%TDDFPT")
243 tddfpt_control => dft_control%tddfpt2_control
244 END IF
245
246 ! input section print/xc
247 CALL tddfpt_input(qs_env, tddfpt_section, tddfpt_control, do_hfx, do_admm, do_exck, &
248 do_hfxsr, do_hfxlr, xc_section, tddfpt_print_section, &
249 lri_section, hfxsr_section)
250
251 log_unit = cp_print_key_unit_nr(logger, tddfpt_print_section, "PROGRAM_BANNER", &
252 extension=".tddfptLog")
253
254 tddfpt_control%do_hfx = do_hfx
255 tddfpt_control%do_admm = do_admm
256 tddfpt_control%do_hfxsr = do_hfxsr
257 tddfpt_control%hfxsr_primbas = 0
258 tddfpt_control%hfxsr_re_int = .true.
259 tddfpt_control%do_hfxlr = do_hfxlr
260 tddfpt_control%do_exck = do_exck
261 do_sf = tddfpt_control%spinflip /= no_sf_tddfpt
262 IF (do_sf) CALL cite_reference(hernandez2025)
263 IF (tddfpt_control%do_hfxlr) THEN
264 kernel_section => section_vals_get_subs_vals(tddfpt_section, "XC%HFX_KERNEL%HFXLR")
265 CALL section_vals_val_get(kernel_section, "RCUT", r_val=tddfpt_control%hfxlr_rcut)
266 CALL section_vals_val_get(kernel_section, "SCALE", r_val=tddfpt_control%hfxlr_scale)
267 END IF
268
269 soc_section => section_vals_get_subs_vals(tddfpt_section, "SOC")
270 CALL section_vals_get(soc_section, explicit=do_soc)
271
272 IF (do_soc) THEN
273 ! start with multiplicity that is not specified in input
274 ! so that excited-state gradient is for multiplicity given in input
275 lmult_tmp = tddfpt_control%rks_triplets
276 tddfpt_control%rks_triplets = .NOT. (tddfpt_control%rks_triplets)
277 END IF
278
279 CALL cite_reference(iannuzzi2005)
280 IF (tddfpt_control%kernel == tddfpt_kernel_stda) THEN
281 CALL cite_reference(grimme2013)
282 CALL cite_reference(grimme2016)
283 END IF
284
285 CALL tddfpt_header(log_unit)
286 CALL kernel_info(log_unit, dft_control, tddfpt_control, xc_section)
287
288 IF (do_kpoints) THEN
289 IF (calc_forces) THEN
290 cpabort("TDDFPT forces are not implemented for k-points")
291 END IF
292 IF (do_rixs) THEN
293 cpabort("RIXS/TDDFPT is not implemented for k-points")
294 END IF
295 IF (do_soc) THEN
296 cpabort("TDDFPT-SOC is not implemented for k-points")
297 END IF
298 CALL tddfpt_kpoint_independent_particle(qs_env, logger, tddfpt_control)
299 CALL cp_print_key_finished_output(log_unit, &
300 logger, &
301 tddfpt_print_section, &
302 "PROGRAM_BANNER")
303 CALL timestop(handle)
304 RETURN
305 END IF
306
307 CALL get_qs_env(qs_env, &
308 blacs_env=blacs_env, &
309 cell=cell, &
310 matrix_ks=matrix_ks, &
311 matrix_s=matrix_s, &
312 mos=mos, &
313 scf_env=scf_env)
314
315 ! obtain occupied and virtual (unoccupied) ground-state Kohn-Sham orbitals
316 NULLIFY (gs_mos)
317 CALL tddfpt_init_mos(qs_env, gs_mos, log_unit)
318
319 ! obtain smeared occupation numbers
320 IF (tddfpt_control%do_smearing) THEN
321 CALL tddfpt_smeared_occupation(qs_env, gs_mos, log_unit)
322 END IF
323
324 ! obtain corrected KS-matrix
325 CALL tddfpt_oecorr(qs_env, gs_mos, matrix_ks_oep)
326
327 IF ((tddfpt_control%do_lrigpw) .AND. &
328 (tddfpt_control%kernel /= tddfpt_kernel_full)) THEN
329 CALL cp_abort(__location__, "LRI only implemented for full kernel")
330 END IF
331
332 IF (ASSOCIATED(matrix_ks_oep)) matrix_ks => matrix_ks_oep
333
334 ! determine active orbitals
335 ! default is all occupied MOs
336 CALL init_res_method(qs_env, gs_mos, tddfpt_control, tddfpt_section, log_unit)
337
338 ! components of the dipole operator
339 CALL tddfpt_dipole_operator(dipole_op_mos_occ, &
340 tddfpt_control, &
341 gs_mos, &
342 qs_env)
343
344 nspins = SIZE(gs_mos)
345 ! multiplicity of molecular system
346 IF (nspins > 1) THEN
347 mult = abs(SIZE(gs_mos(1)%evals_occ) - SIZE(gs_mos(2)%evals_occ)) + 1
348 IF (mult > 2) THEN
349 CALL cp_warn(__location__, "There is a convergence issue for multiplicity >= 3")
350 END IF
351 ELSE
352 IF (tddfpt_control%rks_triplets) THEN
353 mult = 3
354 ELSE
355 mult = 1
356 END IF
357 END IF
358
359 ! split mpi communicator
360 ALLOCATE (my_mos(nspins), my_active(nspins))
361 DO ispin = 1, nspins
362 my_mos(ispin) = gs_mos(ispin)%mos_occ
363 my_active(ispin) = gs_mos(ispin)%mos_active
364 END DO
365 CALL tddfpt_sub_env_init(sub_env, qs_env, &
366 mos_occ=my_mos(:), mos_active=my_active(:), &
367 kernel=tddfpt_control%kernel)
368 DEALLOCATE (my_mos, my_active)
369
370 IF (tddfpt_control%kernel == tddfpt_kernel_full) THEN
371 ! create environment for Full Kernel
372 IF (dft_control%qs_control%xtb) THEN
373 cpabort("TDDFPT: xTB only works with sTDA Kernel")
374 END IF
375
376 IF (tddfpt_control%do_hfxsr) THEN
377 kernel_section => section_vals_get_subs_vals(tddfpt_section, "XC%HFX_KERNEL")
378 CALL section_vals_val_get(kernel_section, "HFXSR_PRIMBAS", &
379 i_val=tddfpt_control%hfxsr_primbas)
380 ! basis set
381 CALL create_minbas_set(qs_env, log_unit, basis_type="TDA_HFX", &
382 primitive=tddfpt_control%hfxsr_primbas)
383 ! admm control
384 ALLOCATE (full_kernel_env%admm_control)
385 full_kernel_env%admm_control%purification_method = do_admm_purify_none
386 full_kernel_env%admm_control%method = do_admm_basis_projection
387 full_kernel_env%admm_control%scaling_model = do_admm_exch_scaling_none
388 full_kernel_env%admm_control%aux_exch_func = do_admm_aux_exch_func_none
389 ! hfx section
390 full_kernel_env%hfxsr_section => hfxsr_section
391 !
392 CALL aux_admm_init(qs_env, mos, full_kernel_env%admm_env, &
393 full_kernel_env%admm_control, "TDA_HFX")
394 CALL get_admm_env(full_kernel_env%admm_env, mos_aux_fit=mos_aux_fit, &
395 matrix_s_aux_fit=matrix_s_aux_fit, &
396 matrix_s_aux_fit_vs_orb=matrix_s_aux_fit_vs_orb)
397 CALL admm_fit_mo_coeffs(full_kernel_env%admm_env, matrix_s_aux_fit, &
398 matrix_s_aux_fit_vs_orb, mos, mos_aux_fit, .true.)
399 ! x_data
400 CALL get_qs_env(qs_env, cell=cell, atomic_kind_set=atomic_kind_set, &
401 qs_kind_set=qs_kind_set, particle_set=particle_set, &
402 para_env=para_env)
403 CALL hfx_create(full_kernel_env%x_data, para_env, hfxsr_section, atomic_kind_set, &
404 qs_kind_set, particle_set, dft_control, cell, orb_basis="TDA_HFX")
405 END IF
406
407 ! allocate pools and work matrices
408 nstates = tddfpt_control%nstates
409 !! Too many states can lead to Problems
410 !! You should be warned if there are more states
411 !! than occ-virt Combinations!!
412 CALL cp_fm_get_info(gs_mos(1)%mos_occ, ncol_global=nocc)
413 IF (tddfpt_control%spinflip == no_sf_tddfpt) THEN
414 CALL cp_fm_get_info(gs_mos(1)%mos_virt, ncol_global=nvirt)
415 ELSE
416 CALL cp_fm_get_info(gs_mos(2)%mos_virt, ncol_global=nvirt)
417 END IF
418 nstate_max = nocc*nvirt
419 IF (nstates > nstate_max) THEN
420 cpwarn("NUMBER OF EXCITED STATES COULD LEAD TO PROBLEMS!")
421 cpwarn("Experimental: CHANGED NSTATES TO ITS MAXIMUM VALUE!")
422 nstates = nstate_max
423 tddfpt_control%nstates = nstate_max
424 END IF
425 CALL tddfpt_create_work_matrices(work_matrices, gs_mos, nstates, &
426 do_hfx, do_admm, do_hfxlr, do_exck, do_sf, qs_env, sub_env)
427
428 ! create full_kernel and admm_kernel within tddfpt_energies
429 kernel_env%full_kernel => full_kernel_env
430 kernel_env%admm_kernel => kernel_env_admm_aux
431 NULLIFY (kernel_env%stda_kernel)
432 IF (do_hfxsr) THEN
433 ! work matrices for SR HFX
434 CALL hfxsr_create_work_matrices(work_matrices, qs_env, full_kernel_env%admm_env)
435 END IF
436 IF (do_hfxlr) THEN
437 ! calculate S_half and Lowdin MO coefficients
438 CALL get_lowdin_mo_coefficients(qs_env, sub_env, work_matrices)
439 END IF
440 ELSE IF (tddfpt_control%kernel == tddfpt_kernel_stda) THEN
441 ! setup for kernel_stda outside tddfpt_energies
442 CALL cp_fm_get_info(gs_mos(1)%mos_occ, nrow_global=nao)
443 nactive = tddfpt_control%nactive
444 CALL allocate_stda_env(qs_env, stda_kernel, nao, nactive)
445 ! sTDA parameters
446 CALL stda_init_param(qs_env, stda_kernel, tddfpt_control%stda_control)
447 ! allocate pools and work matrices
448 nstates = tddfpt_control%nstates
449 CALL stda_create_work_matrices(work_matrices, gs_mos, nstates, qs_env, sub_env)
450 !
451 CALL stda_init_matrices(qs_env, stda_kernel, sub_env, &
452 work_matrices, tddfpt_control)
453 !
454 kernel_env%stda_kernel => stda_kernel
455 NULLIFY (kernel_env%full_kernel)
456 NULLIFY (kernel_env%admm_kernel)
457 ELSE IF (tddfpt_control%kernel == tddfpt_kernel_none) THEN
458 ! allocate pools and work matrices
459 nstates = tddfpt_control%nstates
460 CALL stda_create_work_matrices(work_matrices, gs_mos, nstates, qs_env, sub_env)
461 NULLIFY (kernel_env%full_kernel)
462 NULLIFY (kernel_env%admm_kernel)
463 NULLIFY (kernel_env%stda_kernel)
464 END IF
465
466 IF (do_sf) THEN
467 ! only alpha -> beta excitations are considered in spin-flip TDDFT
468 ALLOCATE (evects(1, nstates))
469 ELSE
470 ALLOCATE (evects(nspins, nstates))
471 END IF
472 ALLOCATE (evals(nstates))
473 ALLOCATE (s_evects(SIZE(evects, 1), nstates))
474
475 DO istate = 1, nstates
476 DO ispin = 1, SIZE(evects, 1)
477 CALL fm_pool_create_fm( &
478 work_matrices%fm_pool_ao_mo_active(ispin)%pool, &
479 s_evects(ispin, istate))
480 END DO
481 END DO
482
483 IF (.NOT. do_soc) THEN
484 ! compute tddfpt excitation energies of multiplicity mult
485 CALL tddfpt_energies(qs_env, nstates, nspins, work_matrices, &
486 tddfpt_control, logger, tddfpt_print_section, evects, evals, &
487 gs_mos, tddfpt_section, s_evects, matrix_s, kernel_env, matrix_ks, &
488 sub_env, ostrength, dipole_op_mos_occ, mult, xc_section, full_kernel_env, &
489 kernel_env_admm_aux)
490 ELSE
491 CALL tddfpt_soc_energies(qs_env, nstates, work_matrices, &
492 tddfpt_control, logger, tddfpt_print_section, &
493 evects, evals, ostrength, &
494 gs_mos, tddfpt_section, s_evects, matrix_s, kernel_env, matrix_ks, &
495 sub_env, dipole_op_mos_occ, lmult_tmp, xc_section, full_kernel_env, &
496 kernel_env_admm_aux)
497 END IF
498
499 !print forces for selected states
500 IF (calc_forces) THEN
501 CALL tddfpt_print_forces(qs_env, evects, evals, ostrength, &
502 tddfpt_print_section, gs_mos, &
503 kernel_env, sub_env, work_matrices)
504 END IF
505
506 ! excited state potential energy surface
507 IF (qs_env%excited_state) THEN
508 IF (sub_env%is_split) THEN
509 CALL cp_abort(__location__, &
510 "Excited state forces not possible when states"// &
511 " are distributed to different CPU pools.")
512 END IF
513 ! for gradients unshifted KS matrix
514 IF (ASSOCIATED(matrix_ks_oep)) CALL get_qs_env(qs_env, matrix_ks=matrix_ks)
515 CALL get_qs_env(qs_env, exstate_env=ex_env)
516 state_change = .false.
517 IF (ex_env%state > 0) THEN
518 my_state = ex_env%state
519 ELSE IF (ex_env%state < 0) THEN
520 ! state following
521 ALLOCATE (my_mos(nspins))
522 DO ispin = 1, nspins
523 my_mos(ispin) = gs_mos(ispin)%mos_occ
524 END DO
525 my_state = abs(ex_env%state)
526 CALL assign_state(qs_env, matrix_s, evects, my_mos, ex_env%wfn_history, my_state)
527 DEALLOCATE (my_mos)
528 IF (my_state /= abs(ex_env%state)) THEN
529 state_change = .true.
530 old_state = abs(ex_env%state)
531 END IF
532 ex_env%state = -my_state
533 ELSE
534 CALL cp_warn(__location__, &
535 "Active excited state not assigned. Use the first state.")
536 my_state = 1
537 END IF
538 cpassert(my_state > 0)
539 IF (my_state > nstates) THEN
540 CALL cp_warn(__location__, &
541 "There were not enough excited states calculated.")
542 cpabort("excited state potential energy surface")
543 END IF
544 !
545 ! energy
546 ex_env%evalue = evals(my_state)
547 ! excitation vector
548 CALL cp_fm_release(ex_env%evect)
549 ALLOCATE (ex_env%evect(SIZE(evects, 1)))
550 DO ispin = 1, SIZE(evects, 1)
551 CALL cp_fm_get_info(matrix=evects(ispin, 1), &
552 matrix_struct=matrix_struct)
553 CALL cp_fm_create(ex_env%evect(ispin), matrix_struct)
554 CALL cp_fm_to_fm(evects(ispin, my_state), ex_env%evect(ispin))
555 END DO
556
557 IF (log_unit > 0) THEN
558 gsval = ex_env%wfn_history%gsval
559 gsmin = ex_env%wfn_history%gsmin
560 xsval = ex_env%wfn_history%xsval
561 WRITE (log_unit, "(1X,A,T40,F10.6,A,T62,F10.6,A)") "Ground state orbital alignment:", &
562 gsmin, "[MinVal]", gsval, "[Average]"
563 WRITE (log_unit, "(1X,A,T71,F10.6)") "Excitation vector alignment:", xsval
564 IF (state_change) THEN
565 WRITE (log_unit, "(1X,A,I5,T60,A14,T76,I5)") &
566 "Target state has been changed from state ", &
567 old_state, " to new state ", my_state
568 END IF
569 WRITE (log_unit, "(1X,A,I4,A,F12.5,A)") "Calculate properties for state:", &
570 my_state, " with excitation energy ", ex_env%evalue*evolt, " eV"
571 END IF
572
573 ! Calculate response vector
574 IF (calc_forces) THEN
575 CALL tddfpt_forces_main(qs_env, gs_mos, ex_env, kernel_env, &
576 sub_env, work_matrices)
577 END IF
578 END IF
579
580 ! share evals, evects and mo_coefs with rixs
581 IF (do_rixs) THEN
582 ! copy evals
583 valence_state%nstates = nstates
584 ALLOCATE (valence_state%evals(SIZE(evals)))
585 valence_state%evals(:) = evals(:)
586
587 ALLOCATE (valence_state%evects(nspins, nstates))
588 ALLOCATE (valence_state%mos_active(nspins))
589 DO ispin = 1, nspins
590 ! copy evects
591 DO istate = 1, nstates
592 CALL cp_fm_get_info(matrix=evects(ispin, istate), &
593 matrix_struct=matrix_struct)
594 CALL cp_fm_create(valence_state%evects(ispin, istate), matrix_struct)
595 CALL cp_fm_to_fm(evects(ispin, istate), valence_state%evects(ispin, istate))
596 END DO
597 ! copy mos_occ
598 CALL cp_fm_get_info(matrix=gs_mos(ispin)%mos_active, &
599 matrix_struct=matrix_struct)
600 CALL cp_fm_create(valence_state%mos_active(ispin), matrix_struct)
601 CALL cp_fm_to_fm(gs_mos(ispin)%mos_active, valence_state%mos_active(ispin))
602 END DO
603 END IF
604
605 ! clean up
606 CALL cp_fm_release(evects)
607 CALL cp_fm_release(s_evects)
608
609 CALL cp_print_key_finished_output(log_unit, &
610 logger, &
611 tddfpt_print_section, &
612 "PROGRAM_BANNER")
613
614 DEALLOCATE (evals, ostrength)
615
616 IF (tddfpt_control%kernel == tddfpt_kernel_full) THEN
617 IF (do_admm) CALL release_kernel_env(kernel_env%admm_kernel)
618 IF (tddfpt_control%do_lrigpw) THEN
619 CALL lri_env_release(kernel_env%full_kernel%lri_env)
620 DEALLOCATE (kernel_env%full_kernel%lri_env)
621 CALL lri_density_release(kernel_env%full_kernel%lri_density)
622 DEALLOCATE (kernel_env%full_kernel%lri_density)
623 END IF
624 CALL release_kernel_env(kernel_env%full_kernel)
625 ELSE IF (tddfpt_control%kernel == tddfpt_kernel_stda) THEN
626 CALL deallocate_stda_env(stda_kernel)
627 ELSE IF (tddfpt_control%kernel == tddfpt_kernel_none) THEN
628 !
629 ELSE
630 cpabort('Unknown kernel type')
631 END IF
632 CALL tddfpt_release_work_matrices(work_matrices, sub_env)
633 CALL tddfpt_sub_env_release(sub_env)
634
635 CALL cp_fm_release(dipole_op_mos_occ)
636
637 DO ispin = nspins, 1, -1
638 CALL tddfpt_release_ground_state_mos(gs_mos(ispin))
639 END DO
640 DEALLOCATE (gs_mos)
641
642 IF (ASSOCIATED(matrix_ks_oep)) THEN
643 CALL dbcsr_deallocate_matrix_set(matrix_ks_oep)
644 END IF
645
646 CALL timestop(handle)
647
648 END SUBROUTINE tddfpt
649
650! **************************************************************************************************
651!> \brief TDDFPT input
652!> \param qs_env Quickstep environment
653!> \param tddfpt_section ...
654!> \param tddfpt_control ...
655!> \param do_hfx ...
656!> \param do_admm ...
657!> \param do_exck ...
658!> \param do_hfxsr ...
659!> \param do_hfxlr ...
660!> \param xc_section ...
661!> \param tddfpt_print_section ...
662!> \param lri_section ...
663!> \param hfxsr_section ...
664! **************************************************************************************************
665 SUBROUTINE tddfpt_input(qs_env, tddfpt_section, tddfpt_control, do_hfx, do_admm, do_exck, &
666 do_hfxsr, do_hfxlr, xc_section, tddfpt_print_section, lri_section, &
667 hfxsr_section)
668 TYPE(qs_environment_type), POINTER :: qs_env
669 TYPE(section_vals_type), POINTER :: tddfpt_section
670 TYPE(tddfpt2_control_type), POINTER :: tddfpt_control
671 LOGICAL, INTENT(INOUT) :: do_hfx, do_admm, do_exck, do_hfxsr, &
672 do_hfxlr
673 TYPE(section_vals_type), POINTER :: xc_section, tddfpt_print_section, &
674 lri_section, hfxsr_section
675
676 CHARACTER(len=20) :: nstates_str
677 LOGICAL :: exar, exf, exgcp, exhf, exhfxk, exk, &
678 explicit, explicit_root, expot, exvdw, &
679 exwfn, found, same_hfx, use_real_wfn
680 REAL(kind=dp) :: c_hf
681 TYPE(dft_control_type), POINTER :: dft_control
682 TYPE(kpoint_type), POINTER :: kpoints
683 TYPE(section_vals_type), POINTER :: hfx_section, hfx_section_gs, input, &
684 print_sub, xc_root, xc_sub
685
686 NULLIFY (dft_control, input, kpoints)
687 CALL get_qs_env(qs_env, dft_control=dft_control, input=input, kpoints=kpoints)
688
689 IF (dft_control%nimages > 1) THEN
690 IF (tddfpt_control%kernel /= tddfpt_kernel_none) THEN
691 cpabort("TDDFPT with k-points currently supports only KERNEL NONE")
692 END IF
693 CALL get_kpoint_info(kpoints, use_real_wfn=use_real_wfn)
694 IF (use_real_wfn) THEN
695 cpabort("K-point TDDFPT requires complex wavefunctions")
696 END IF
697 IF (tddfpt_control%spinflip /= no_sf_tddfpt) THEN
698 cpabort("Spin-flip TDDFPT is not implemented for k-points")
699 END IF
700 IF (tddfpt_control%do_smearing) THEN
701 cpabort("Smeared-occupation TDDFPT is not implemented for k-points")
702 END IF
703 IF (tddfpt_control%oe_corr /= oe_none) THEN
704 cpabort("Orbital-energy-corrected TDDFPT is not implemented for k-points")
705 END IF
706 IF (tddfpt_control%dipole_form /= 0 .AND. &
707 tddfpt_control%dipole_form /= tddfpt_dipole_velocity .AND. &
708 tddfpt_control%dipole_form /= tddfpt_dipole_scf_moment) THEN
709 cpabort("K-point TDDFPT supports only velocity-form or SCF_MOMENT transition dipoles")
710 END IF
711 END IF
712
713 IF (tddfpt_control%nstates <= 0) THEN
714 CALL integer_to_string(tddfpt_control%nstates, nstates_str)
715 CALL cp_warn(__location__, "TDDFPT calculation was requested for "// &
716 trim(nstates_str)//" excited states: nothing to do.")
717 RETURN
718 END IF
719
720 NULLIFY (tddfpt_print_section)
721 tddfpt_print_section => section_vals_get_subs_vals(tddfpt_section, "PRINT")
722
723 IF (dft_control%nimages > 1) THEN
724 IF (tddfpt_control%do_exciton_descriptors .OR. &
725 tddfpt_control%do_directional_exciton_descriptors) THEN
726 cpabort("Exciton descriptors are not implemented for k-point TDDFPT")
727 END IF
728 print_sub => section_vals_get_subs_vals(tddfpt_print_section, "NTO_ANALYSIS")
729 CALL section_vals_get(print_sub, explicit=explicit)
730 IF (explicit) cpabort("NTO analysis is not implemented for k-point TDDFPT")
731 print_sub => section_vals_get_subs_vals(tddfpt_print_section, "NAMD_PRINT")
732 CALL section_vals_get(print_sub, explicit=explicit)
733 IF (explicit) cpabort("NAMD_PRINT is not implemented for k-point TDDFPT")
734 END IF
735
736 IF (tddfpt_control%kernel == tddfpt_kernel_full) THEN
737 NULLIFY (xc_root)
738 xc_root => section_vals_get_subs_vals(tddfpt_section, "XC")
739 CALL section_vals_get(xc_root, explicit=explicit_root)
740 NULLIFY (xc_section)
741 IF (explicit_root) THEN
742 ! No ADIABATIC_RESCALING option possible
743 NULLIFY (xc_sub)
744 xc_sub => section_vals_get_subs_vals(xc_root, "ADIABATIC_RESCALING")
745 CALL section_vals_get(xc_sub, explicit=exar)
746 IF (exar) THEN
747 CALL cp_warn(__location__, "TDDFPT Kernel with ADIABATIC_RESCALING not possible.")
748 cpabort("TDDFPT Input")
749 END IF
750 ! No GCP_POTENTIAL option possible
751 NULLIFY (xc_sub)
752 xc_sub => section_vals_get_subs_vals(xc_root, "GCP_POTENTIAL")
753 CALL section_vals_get(xc_sub, explicit=exgcp)
754 IF (exgcp) THEN
755 CALL cp_warn(__location__, "TDDFPT Kernel with GCP_POTENTIAL not possible.")
756 cpabort("TDDFPT Input")
757 END IF
758 ! No VDW_POTENTIAL option possible
759 NULLIFY (xc_sub)
760 xc_sub => section_vals_get_subs_vals(xc_root, "VDW_POTENTIAL")
761 CALL section_vals_get(xc_sub, explicit=exvdw)
762 IF (exvdw) THEN
763 CALL cp_warn(__location__, "TDDFPT Kernel with VDW_POTENTIAL not possible.")
764 cpabort("TDDFPT Input")
765 END IF
766 ! No WF_CORRELATION option possible
767 NULLIFY (xc_sub)
768 xc_sub => section_vals_get_subs_vals(xc_root, "WF_CORRELATION")
769 CALL section_vals_get(xc_sub, explicit=exwfn)
770 IF (exwfn) THEN
771 CALL cp_warn(__location__, "TDDFPT Kernel with WF_CORRELATION not possible.")
772 cpabort("TDDFPT Input")
773 END IF
774 ! No XC_POTENTIAL option possible
775 NULLIFY (xc_sub)
776 xc_sub => section_vals_get_subs_vals(xc_root, "XC_POTENTIAL")
777 CALL section_vals_get(xc_sub, explicit=expot)
778 IF (expot) THEN
779 CALL cp_warn(__location__, "TDDFPT Kernel with XC_POTENTIAL not possible.")
780 cpabort("TDDFPT Input")
781 END IF
782 !
783 NULLIFY (xc_sub)
784 xc_sub => section_vals_get_subs_vals(xc_root, "XC_FUNCTIONAL")
785 CALL section_vals_get(xc_sub, explicit=exf)
786 NULLIFY (xc_sub)
787 xc_sub => section_vals_get_subs_vals(xc_root, "XC_KERNEL")
788 CALL section_vals_get(xc_sub, explicit=exk)
789 IF ((exf .AND. exk) .OR. .NOT. (exf .OR. exk)) THEN
790 CALL cp_warn(__location__, "TDDFPT Kernel needs XC_FUNCTIONAL or XC_KERNEL section.")
791 cpabort("TDDFPT Input")
792 END IF
793 NULLIFY (xc_sub)
794 xc_sub => section_vals_get_subs_vals(xc_root, "HF")
795 CALL section_vals_get(xc_sub, explicit=exhf)
796 NULLIFY (xc_sub)
797 xc_sub => section_vals_get_subs_vals(xc_root, "HFX_KERNEL")
798 CALL section_vals_get(xc_sub, explicit=exhfxk)
799 !
800 xc_section => xc_root
801 hfx_section => section_vals_get_subs_vals(xc_section, "HF")
802 CALL section_vals_get(hfx_section, explicit=do_hfx)
803 IF (do_hfx) THEN
804 CALL section_vals_val_get(hfx_section, "FRACTION", r_val=c_hf)
805 do_hfx = (c_hf /= 0.0_dp)
806 END IF
807 !TDDFPT only works if the kernel has the same HF section as the DFT%XC one
808 IF (do_hfx) THEN
809 hfx_section_gs => section_vals_get_subs_vals(input, "DFT%XC%HF")
810 CALL compare_hfx_sections(hfx_section, hfx_section_gs, same_hfx)
811 IF (.NOT. same_hfx) THEN
812 cpabort("TDDFPT Kernel must use the same HF section as DFT%XC or no HF at all.")
813 END IF
814 END IF
815
816 do_admm = do_hfx .AND. dft_control%do_admm
817 IF (do_admm) THEN
818 ! 'admm_env%xc_section_primary' and 'admm_env%xc_section_aux' need to be redefined
819 CALL cp_abort(__location__, &
820 "ADMM is not implemented for a TDDFT kernel XC-functional which is different from "// &
821 "the one used for the ground-state calculation. A ground-state 'admm_env' cannot be reused.")
822 END IF
823 ! SET HFX_KERNEL and/or XC_KERNEL
824 IF (exk) THEN
825 do_exck = .true.
826 ELSE
827 do_exck = .false.
828 END IF
829 IF (exhfxk) THEN
830 xc_sub => section_vals_get_subs_vals(xc_root, "HFX_KERNEL")
831 CALL section_vals_val_get(xc_sub, "DO_HFXSR", l_val=do_hfxsr)
832 xc_sub => section_vals_get_subs_vals(xc_root, "HFX_KERNEL%HFXLR")
833 CALL section_vals_get(xc_sub, explicit=do_hfxlr)
834 ELSE
835 do_hfxsr = .false.
836 do_hfxlr = .false.
837 END IF
838 ELSE
839 xc_section => section_vals_get_subs_vals(input, "DFT%XC")
840 hfx_section => section_vals_get_subs_vals(xc_section, "HF")
841 CALL section_vals_get(hfx_section, explicit=do_hfx)
842 IF (do_hfx) THEN
843 CALL section_vals_val_get(hfx_section, "FRACTION", r_val=c_hf)
844 do_hfx = (c_hf /= 0.0_dp)
845 END IF
846 do_admm = do_hfx .AND. dft_control%do_admm
847 do_exck = .false.
848 do_hfxsr = .false.
849 do_hfxlr = .false.
850 END IF
851 ELSE
852 do_hfx = .false.
853 do_admm = .false.
854 do_exck = .false.
855 do_hfxsr = .false.
856 do_hfxlr = .false.
857 END IF
858
859 ! reset rks_triplets if UKS is in use
860 IF (tddfpt_control%rks_triplets .AND. dft_control%nspins > 1) THEN
861 tddfpt_control%rks_triplets = .false.
862 CALL cp_warn(__location__, "Keyword RKS_TRIPLETS has been ignored for spin-polarised calculations")
863 END IF
864
865 ! lri input
866 IF (tddfpt_control%do_lrigpw) THEN
867 lri_section => section_vals_get_subs_vals(tddfpt_section, "LRIGPW")
868 END IF
869
870 ! set defaults for short range HFX
871 NULLIFY (hfxsr_section)
872 IF (do_hfxsr) THEN
873 hfxsr_section => section_vals_get_subs_vals(tddfpt_section, "XC%HFX_KERNEL%HF")
874 CALL section_vals_get(hfxsr_section, explicit=found)
875 IF (.NOT. found) THEN
876 cpabort("HFXSR option needs &HF section defined")
877 END IF
878 CALL section_vals_val_get(hfxsr_section, "INTERACTION_POTENTIAL%POTENTIAL_TYPE", explicit=found)
879 IF (.NOT. found) THEN
880 CALL section_vals_val_set(hfxsr_section, "INTERACTION_POTENTIAL%POTENTIAL_TYPE", &
881 i_val=do_potential_truncated)
882 END IF
883 CALL section_vals_val_get(hfxsr_section, "INTERACTION_POTENTIAL%CUTOFF_RADIUS", explicit=found)
884 IF (.NOT. found) THEN
885 CALL section_vals_val_set(hfxsr_section, "INTERACTION_POTENTIAL%CUTOFF_RADIUS", r_val=7.5589_dp)
886 END IF
887 CALL section_vals_val_get(hfxsr_section, "RI%_SECTION_PARAMETERS_", l_val=found)
888 IF (found) THEN
889 CALL cp_abort(__location__, "Short range TDA kernel with RI not possible")
890 END IF
891 END IF
892
893 END SUBROUTINE tddfpt_input
894
895! **************************************************************************************************
896!> \brief ...
897!> \param log_unit ...
898!> \param dft_control ...
899!> \param tddfpt_control ...
900!> \param xc_section ...
901! **************************************************************************************************
902 SUBROUTINE kernel_info(log_unit, dft_control, tddfpt_control, xc_section)
903 INTEGER, INTENT(IN) :: log_unit
904 TYPE(dft_control_type), POINTER :: dft_control
905 TYPE(tddfpt2_control_type), POINTER :: tddfpt_control
906 TYPE(section_vals_type), POINTER :: xc_section
907
908 CHARACTER(LEN=4) :: ktype
909 LOGICAL :: lsd
910
911 lsd = (dft_control%nspins > 1)
912 IF (tddfpt_control%kernel == tddfpt_kernel_full) THEN
913 ktype = "FULL"
914 IF (log_unit > 0) THEN
915 WRITE (log_unit, "(T2,A,T77,A4)") "KERNEL|", trim(ktype)
916 CALL xc_write(log_unit, xc_section, lsd)
917 IF (tddfpt_control%do_hfx) THEN
918 IF (tddfpt_control%do_admm) THEN
919 WRITE (log_unit, "(T2,A,T62,A19)") "KERNEL|", "ADMM Exact Exchange"
920 IF (tddfpt_control%admm_xc_correction) THEN
921 WRITE (log_unit, "(T2,A,T60,A21)") "KERNEL|", "Apply ADMM Kernel XC Correction"
922 END IF
923 IF (tddfpt_control%admm_symm) THEN
924 WRITE (log_unit, "(T2,A,T60,A21)") "KERNEL|", "Symmetric ADMM Kernel"
925 END IF
926 ELSE
927 WRITE (log_unit, "(T2,A,T67,A14)") "KERNEL|", "Exact Exchange"
928 END IF
929 END IF
930 IF (tddfpt_control%do_hfxsr) THEN
931 WRITE (log_unit, "(T2,A,T43,A38)") "KERNEL|", "Short range HFX approximation"
932 END IF
933 IF (tddfpt_control%do_hfxlr) THEN
934 WRITE (log_unit, "(T2,A,T43,A38)") "KERNEL|", "Long range HFX approximation"
935 END IF
936 IF (tddfpt_control%do_lrigpw) THEN
937 WRITE (log_unit, "(T2,A,T42,A39)") "KERNEL|", "LRI approximation of transition density"
938 END IF
939 END IF
940 ELSE IF (tddfpt_control%kernel == tddfpt_kernel_stda) THEN
941 ktype = "sTDA"
942 IF (log_unit > 0) THEN
943 WRITE (log_unit, "(T2,A,T77,A4)") "KERNEL|", trim(ktype)
944 IF (tddfpt_control%stda_control%do_ewald) THEN
945 WRITE (log_unit, "(T2,A,T78,A3)") "KERNEL| Coulomb term uses Ewald summation"
946 ELSE
947 WRITE (log_unit, "(T2,A,T78,A3)") "KERNEL| Coulomb term uses direct summation (MIC)"
948 END IF
949 IF (tddfpt_control%stda_control%do_exchange) THEN
950 WRITE (log_unit, "(T2,A,T78,A3)") "KERNEL| Exact exchange term", "YES"
951 WRITE (log_unit, "(T2,A,T71,F10.3)") "KERNEL| Short range HFX fraction:", &
952 tddfpt_control%stda_control%hfx_fraction
953 ELSE
954 WRITE (log_unit, "(T2,A,T79,A2)") "KERNEL| Exact exchange term", "NO"
955 END IF
956 WRITE (log_unit, "(T2,A,T66,E15.3)") "KERNEL| Transition density filter", &
957 tddfpt_control%stda_control%eps_td_filter
958 END IF
959 ELSE IF (tddfpt_control%kernel == tddfpt_kernel_none) THEN
960 ktype = "NONE"
961 IF (log_unit > 0) THEN
962 WRITE (log_unit, "(T2,A,T77,A4)") "KERNEL|", trim(ktype)
963 END IF
964 ELSE
965 !CPABORT("Unknown kernel")
966 END IF
967 !
968 IF (log_unit > 0) THEN
969 IF (tddfpt_control%rks_triplets) THEN
970 WRITE (log_unit, "(T2,A,T74,A7)") "KERNEL| Spin symmetry of excitations", "Triplet"
971 ELSE IF (lsd) THEN
972 ! Spin-conserving excitations where requested
973 IF (tddfpt_control%spinflip == no_sf_tddfpt) THEN
974 WRITE (log_unit, "(T2,A,T69,A12)") "KERNEL| Spin symmetry of excitations", "Unrestricted"
975 ! Spin-flip excitations with collinear exchange-correlation kernel requested
976 ELSE IF (tddfpt_control%spinflip == tddfpt_sf_col) THEN
977 WRITE (log_unit, "(T2,A,T72,A9)") "KERNEL| Spin flip", "Collinear"
978 ! Spin-flip excitations with noncollinear exchange-correlation kernel requested
979 ELSE IF (tddfpt_control%spinflip == tddfpt_sf_noncol) THEN
980 WRITE (log_unit, "(T2,A,T69,A12)") "KERNEL| Spin flip", "Noncollinear"
981 END IF
982 ELSE
983 WRITE (log_unit, "(T2,A,T74,A7)") "KERNEL| Spin symmetry of excitations", "Singlet"
984 END IF
985 WRITE (log_unit, "(T2,A,T73,I8)") "TDDFPT| Number of states calculated", tddfpt_control%nstates
986 WRITE (log_unit, "(T2,A,T73,I8)") "TDDFPT| Number of Davidson iterations", tddfpt_control%niters
987 WRITE (log_unit, "(T2,A,T66,E15.3)") "TDDFPT| Davidson iteration convergence", tddfpt_control%conv
988 WRITE (log_unit, "(T2,A,T73,I8)") "TDDFPT| Max. number of Krylov space vectors", tddfpt_control%nkvs
989 END IF
990
991 END SUBROUTINE kernel_info
992
993! **************************************************************************************************
994!> \brief Print independent-particle vertical transitions for k-point calculations.
995!> \param qs_env ...
996!> \param logger ...
997!> \param tddfpt_control ...
998! **************************************************************************************************
999 SUBROUTINE tddfpt_kpoint_independent_particle(qs_env, logger, tddfpt_control)
1000 TYPE(qs_environment_type), POINTER :: qs_env
1001 TYPE(cp_logger_type), POINTER :: logger
1002 TYPE(tddfpt2_control_type), POINTER :: tddfpt_control
1003
1004 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt_kpoint_independent_particle'
1005
1006 COMPLEX(KIND=dp), ALLOCATABLE, &
1007 DIMENSION(:, :, :, :, :) :: kpoint_dipole
1008 INTEGER :: handle, ideriv, ikp, ikp_local, iocc, ispin, istate, itrans, ivirt, log_unit, &
1009 nao, nkp, nkp_local, nspins, nstates, ntrans_kpoint, ntrans_spin, ntrans_total, &
1010 spin_offset, trans_index
1011 INTEGER, ALLOCATABLE, DIMENSION(:) :: inds
1012 INTEGER, DIMENSION(2) :: kp_range
1013 INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
1014 INTEGER, DIMENSION(maxspins) :: homo_spin, nao_spin, nmo_spin, nvirt_spin
1015 LOGICAL :: my_kpgrp, use_scf_moment_dipoles
1016 REAL(kind=dp) :: checksum, dipole_im, dipole_re, fsum, &
1017 gap, oscillator_factor, spin_factor
1018 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: eigenvalues_kp, evals, &
1019 oscillator_strength, transition_energy
1020 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: transition_dipole_im, &
1021 transition_dipole_re
1022 REAL(kind=dp), DIMENSION(:), POINTER :: eigenvalues, wkp
1023 REAL(kind=dp), DIMENSION(nderivs) :: transition_dipole_abs
1024 TYPE(cp_blacs_env_type), POINTER :: blacs_env, blacs_env_all
1025 TYPE(cp_fm_struct_type), POINTER :: fm_struct, moment_struct
1026 TYPE(cp_fm_type) :: fm_dummy, fm_tmp, mo_coeff_im_global, &
1027 mo_coeff_re_global, moment_im, &
1028 moment_re
1029 TYPE(cp_fm_type), POINTER :: mo_coeff_im, mo_coeff_re
1030 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: overlap_deriv
1031 TYPE(dbcsr_type), POINTER :: cmatrix, rmatrix
1032 TYPE(dft_control_type), POINTER :: dft_control
1033 TYPE(kpoint_env_p_type), DIMENSION(:), POINTER :: kp_env
1034 TYPE(kpoint_env_type), POINTER :: kp
1035 TYPE(kpoint_type), POINTER :: kpoints
1036 TYPE(mo_set_type), DIMENSION(:, :), POINTER :: mos_kp
1037 TYPE(mp_para_env_type), POINTER :: para_env, para_env_inter_kp, para_env_kp
1038 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
1039 POINTER :: sab_kp, sab_orb
1040 TYPE(qs_ks_env_type), POINTER :: ks_env
1041
1042 CALL timeset(routinen, handle)
1043
1044 NULLIFY (blacs_env, blacs_env_all, cell_to_index, cmatrix, dft_control, eigenvalues, &
1045 fm_struct, kp, kp_env, kpoints, ks_env, mo_coeff_im, mo_coeff_re, &
1046 moment_struct, mos_kp, overlap_deriv, para_env, para_env_inter_kp, para_env_kp, &
1047 rmatrix, sab_kp, sab_orb, wkp)
1048 CALL get_qs_env(qs_env, dft_control=dft_control, kpoints=kpoints, ks_env=ks_env, &
1049 sab_orb=sab_orb)
1050 cpassert(ASSOCIATED(kpoints))
1051
1052 CALL get_kpoint_info(kpoints, nkp=nkp, kp_range=kp_range, kp_env=kp_env, &
1053 para_env=para_env, blacs_env_all=blacs_env_all, &
1054 para_env_inter_kp=para_env_inter_kp, para_env_kp=para_env_kp, &
1055 blacs_env=blacs_env, wkp=wkp, cell_to_index=cell_to_index, &
1056 sab_nl=sab_kp)
1057 cpassert(ASSOCIATED(para_env))
1058 cpassert(ASSOCIATED(para_env_inter_kp))
1059 cpassert(ASSOCIATED(para_env_kp))
1060 cpassert(ASSOCIATED(blacs_env_all))
1061 cpassert(ASSOCIATED(blacs_env))
1062 cpassert(ASSOCIATED(kp_env))
1063 cpassert(ASSOCIATED(ks_env))
1064 cpassert(ASSOCIATED(sab_orb))
1065 cpassert(ASSOCIATED(sab_kp))
1066 cpassert(ASSOCIATED(cell_to_index))
1067
1068 nspins = dft_control%nspins
1069 nmo_spin = 0
1070 homo_spin = 0
1071 nao_spin = 0
1072 nkp_local = max(0, kp_range(2) - kp_range(1) + 1)
1073 IF (nkp_local > 0) THEN
1074 kp => kp_env(1)%kpoint_env
1075 mos_kp => kp%mos
1076 cpassert(ASSOCIATED(mos_kp))
1077 cpassert(SIZE(mos_kp, 2) == nspins)
1078 DO ispin = 1, nspins
1079 CALL get_mo_set(mos_kp(1, ispin), nmo=nmo_spin(ispin), homo=homo_spin(ispin), &
1080 nao=nao_spin(ispin))
1081 END DO
1082 END IF
1083 CALL para_env%max(nmo_spin)
1084 CALL para_env%max(homo_spin)
1085 CALL para_env%max(nao_spin)
1086
1087 ntrans_kpoint = 0
1088 DO ispin = 1, nspins
1089 nvirt_spin(ispin) = nmo_spin(ispin) - homo_spin(ispin)
1090 IF (homo_spin(ispin) <= 0 .OR. nvirt_spin(ispin) <= 0) THEN
1091 cpabort("At least one occupied and one unoccupied MO are required for k-point TDDFPT")
1092 END IF
1093 ntrans_kpoint = ntrans_kpoint + homo_spin(ispin)*nvirt_spin(ispin)
1094 END DO
1095 ntrans_total = nkp*ntrans_kpoint
1096 IF (ntrans_total <= 0) THEN
1097 cpabort("No independent-particle k-point transitions available")
1098 END IF
1099
1100 ALLOCATE (transition_energy(ntrans_total), transition_dipole_re(ntrans_total, nderivs), &
1101 transition_dipole_im(ntrans_total, nderivs), oscillator_strength(ntrans_total), &
1102 inds(ntrans_total))
1103 transition_energy = 0.0_dp
1104 transition_dipole_re = 0.0_dp
1105 transition_dipole_im = 0.0_dp
1106 oscillator_strength = 0.0_dp
1107 use_scf_moment_dipoles = (tddfpt_control%dipole_form == tddfpt_dipole_scf_moment)
1108 IF (use_scf_moment_dipoles) THEN
1109 CALL cp_warn(__location__, "SCF_MOMENT k-point dipoles use direct SCF MO matrix "// &
1110 "elements; compare folded energy blocks, not individual degenerate states.")
1111 CALL qs_moment_kpoints_scf_mos(qs_env, kpoint_dipole)
1112 END IF
1113
1114 IF (.NOT. use_scf_moment_dipoles) THEN
1115 CALL build_overlap_matrix(ks_env, matrixkp_s=overlap_deriv, nderivative=1, &
1116 basis_type_a="ORB", basis_type_b="ORB", sab_nl=sab_orb, &
1117 ext_kpoints=kpoints)
1118
1119 ALLOCATE (rmatrix, cmatrix)
1120 CALL dbcsr_create(rmatrix, template=overlap_deriv(1, 1)%matrix, &
1121 matrix_type=dbcsr_type_symmetric)
1122 CALL dbcsr_create(cmatrix, template=overlap_deriv(1, 1)%matrix, &
1123 matrix_type=dbcsr_type_antisymmetric)
1124 CALL cp_dbcsr_alloc_block_from_nbl(rmatrix, sab_kp)
1125 CALL cp_dbcsr_alloc_block_from_nbl(cmatrix, sab_kp)
1126 END IF
1127
1128 DO ikp = 1, nkp
1129 my_kpgrp = (ikp >= kp_range(1) .AND. ikp <= kp_range(2))
1130 IF (my_kpgrp) THEN
1131 ikp_local = ikp - kp_range(1) + 1
1132 kp => kp_env(ikp_local)%kpoint_env
1133 mos_kp => kp%mos
1134 ELSE
1135 NULLIFY (kp, mos_kp)
1136 END IF
1137 spin_offset = 0
1138 DO ispin = 1, nspins
1139 nao = nao_spin(ispin)
1140 ALLOCATE (eigenvalues_kp(nmo_spin(ispin)))
1141 eigenvalues_kp = 0.0_dp
1142
1143 IF (.NOT. use_scf_moment_dipoles) THEN
1144 CALL cp_fm_struct_create(fm_struct, nrow_global=nao, ncol_global=nmo_spin(ispin), &
1145 para_env=para_env, context=blacs_env_all)
1146 CALL cp_fm_create(mo_coeff_re_global, fm_struct)
1147 CALL cp_fm_create(mo_coeff_im_global, fm_struct)
1148 CALL cp_fm_create(fm_tmp, fm_struct)
1149 CALL cp_fm_struct_release(fm_struct)
1150 CALL cp_fm_struct_create(moment_struct, nrow_global=nmo_spin(ispin), &
1151 ncol_global=nmo_spin(ispin), para_env=para_env, &
1152 context=blacs_env_all)
1153 CALL cp_fm_create(moment_re, moment_struct)
1154 CALL cp_fm_create(moment_im, moment_struct)
1155 CALL cp_fm_struct_release(moment_struct)
1156 END IF
1157
1158 IF (my_kpgrp) THEN
1159 CALL get_mo_set(mos_kp(1, ispin), eigenvalues=eigenvalues)
1160 cpassert(ASSOCIATED(eigenvalues))
1161 IF (para_env_kp%is_source()) THEN
1162 eigenvalues_kp(1:nmo_spin(ispin)) = eigenvalues(1:nmo_spin(ispin))
1163 END IF
1164 IF (.NOT. use_scf_moment_dipoles) THEN
1165 CALL get_mo_set(mos_kp(1, ispin), mo_coeff=mo_coeff_re)
1166 CALL get_mo_set(mos_kp(2, ispin), mo_coeff=mo_coeff_im)
1167 cpassert(ASSOCIATED(mo_coeff_re))
1168 cpassert(ASSOCIATED(mo_coeff_im))
1169 CALL cp_fm_copy_general(mo_coeff_re, mo_coeff_re_global, para_env)
1170 CALL cp_fm_copy_general(mo_coeff_im, mo_coeff_im_global, para_env)
1171 END IF
1172 ELSE IF (.NOT. use_scf_moment_dipoles) THEN
1173 CALL cp_fm_copy_general(fm_dummy, mo_coeff_re_global, para_env)
1174 CALL cp_fm_copy_general(fm_dummy, mo_coeff_im_global, para_env)
1175 END IF
1176 CALL para_env%sum(eigenvalues_kp)
1177
1178 spin_factor = 1.0_dp
1179 IF (nspins == 1) THEN
1180 IF (tddfpt_control%rks_triplets) THEN
1181 spin_factor = 0.0_dp
1182 ELSE
1183 spin_factor = 2.0_dp
1184 END IF
1185 END IF
1186
1187 DO ideriv = 1, nderivs
1188 IF (.NOT. use_scf_moment_dipoles) THEN
1189 CALL dbcsr_set(rmatrix, 0.0_dp)
1190 CALL dbcsr_set(cmatrix, 0.0_dp)
1191 CALL rskp_transform(rmatrix=rmatrix, cmatrix=cmatrix, rsmat=overlap_deriv, &
1192 ispin=ideriv + 1, xkp=kpoints%xkp(:, ikp), &
1193 cell_to_index=cell_to_index, sab_nl=sab_kp)
1194
1195 CALL cp_dbcsr_sm_fm_multiply(rmatrix, mo_coeff_re_global, fm_tmp, nmo_spin(ispin))
1196 CALL parallel_gemm("T", "N", nmo_spin(ispin), nmo_spin(ispin), nao, &
1197 1.0_dp, mo_coeff_re_global, fm_tmp, 0.0_dp, moment_re)
1198 CALL parallel_gemm("T", "N", nmo_spin(ispin), nmo_spin(ispin), nao, &
1199 1.0_dp, mo_coeff_im_global, fm_tmp, 0.0_dp, moment_im)
1200
1201 CALL cp_dbcsr_sm_fm_multiply(rmatrix, mo_coeff_im_global, fm_tmp, nmo_spin(ispin))
1202 CALL parallel_gemm("T", "N", nmo_spin(ispin), nmo_spin(ispin), nao, &
1203 1.0_dp, mo_coeff_re_global, fm_tmp, 1.0_dp, moment_im)
1204 CALL parallel_gemm("T", "N", nmo_spin(ispin), nmo_spin(ispin), nao, &
1205 -1.0_dp, mo_coeff_im_global, fm_tmp, 1.0_dp, moment_re)
1206
1207 CALL cp_dbcsr_sm_fm_multiply(cmatrix, mo_coeff_re_global, fm_tmp, nmo_spin(ispin))
1208 CALL parallel_gemm("T", "N", nmo_spin(ispin), nmo_spin(ispin), nao, &
1209 1.0_dp, mo_coeff_re_global, fm_tmp, 1.0_dp, moment_im)
1210 CALL parallel_gemm("T", "N", nmo_spin(ispin), nmo_spin(ispin), nao, &
1211 -1.0_dp, mo_coeff_im_global, fm_tmp, 1.0_dp, moment_re)
1212
1213 CALL cp_dbcsr_sm_fm_multiply(cmatrix, mo_coeff_im_global, fm_tmp, nmo_spin(ispin))
1214 CALL parallel_gemm("T", "N", nmo_spin(ispin), nmo_spin(ispin), nao, &
1215 -1.0_dp, mo_coeff_re_global, fm_tmp, 1.0_dp, moment_re)
1216 CALL parallel_gemm("T", "N", nmo_spin(ispin), nmo_spin(ispin), nao, &
1217 -1.0_dp, mo_coeff_im_global, fm_tmp, 1.0_dp, moment_im)
1218 END IF
1219
1220 DO iocc = 1, homo_spin(ispin)
1221 DO ivirt = homo_spin(ispin) + 1, nmo_spin(ispin)
1222 trans_index = (ikp - 1)*ntrans_kpoint + spin_offset + &
1223 (iocc - 1)*nvirt_spin(ispin) + ivirt - homo_spin(ispin)
1224 gap = eigenvalues_kp(ivirt) - eigenvalues_kp(iocc)
1225 IF (gap <= 0.0_dp) THEN
1226 cpabort("K-point TDDFPT requires positive occupied-virtual energy gaps")
1227 END IF
1228 IF (use_scf_moment_dipoles) THEN
1229 oscillator_factor = sqrt(spin_factor*wkp(ikp))
1230 dipole_re = real(kpoint_dipole(ispin, ikp, ideriv, iocc, ivirt), kind=dp)
1231 dipole_im = aimag(kpoint_dipole(ispin, ikp, ideriv, iocc, ivirt))
1232 ELSE
1233 oscillator_factor = sqrt(spin_factor*wkp(ikp))/gap
1234 CALL cp_fm_get_element(moment_re, ivirt, iocc, dipole_re)
1235 CALL cp_fm_get_element(moment_im, ivirt, iocc, dipole_im)
1236 END IF
1237 transition_dipole_re(trans_index, ideriv) = oscillator_factor*dipole_re
1238 transition_dipole_im(trans_index, ideriv) = oscillator_factor*dipole_im
1239 END DO
1240 END DO
1241 END DO
1242
1243 DO iocc = 1, homo_spin(ispin)
1244 DO ivirt = homo_spin(ispin) + 1, nmo_spin(ispin)
1245 trans_index = (ikp - 1)*ntrans_kpoint + spin_offset + &
1246 (iocc - 1)*nvirt_spin(ispin) + ivirt - homo_spin(ispin)
1247 transition_energy(trans_index) = eigenvalues_kp(ivirt) - eigenvalues_kp(iocc)
1248 oscillator_strength(trans_index) = 2.0_dp/3.0_dp*transition_energy(trans_index)* &
1249 sum(transition_dipole_re(trans_index, :)**2 + &
1250 transition_dipole_im(trans_index, :)**2)
1251 END DO
1252 END DO
1253 IF (.NOT. use_scf_moment_dipoles) THEN
1254 CALL cp_fm_release(moment_im)
1255 CALL cp_fm_release(moment_re)
1256 CALL cp_fm_release(fm_tmp)
1257 CALL cp_fm_release(mo_coeff_im_global)
1258 CALL cp_fm_release(mo_coeff_re_global)
1259 END IF
1260 DEALLOCATE (eigenvalues_kp)
1261 spin_offset = spin_offset + homo_spin(ispin)*nvirt_spin(ispin)
1262 END DO
1263 END DO
1264
1265 IF (any(transition_energy <= 0.0_dp)) THEN
1266 cpabort("K-point TDDFPT KERNEL NONE requires positive occupied-virtual energy gaps")
1267 END IF
1268
1269 CALL sort(transition_energy, ntrans_total, inds)
1270 nstates = min(tddfpt_control%nstates, ntrans_total)
1271 IF (tddfpt_control%nstates > ntrans_total) THEN
1272 cpwarn("Requested more TDDFPT states than independent-particle k-point transitions")
1273 END IF
1274
1275 ALLOCATE (evals(nstates))
1276 evals(1:nstates) = transition_energy(1:nstates)
1277 checksum = sqrt(sum(evals**2))
1278
1279 log_unit = cp_logger_get_default_io_unit(logger)
1280 IF (log_unit > 0) THEN
1281 WRITE (log_unit, "(1X,A)") "", &
1282 "-------------------------------------------------------------------------------", &
1283 "- TDDFPT K-point Independent-particle Transitions -", &
1284 "-------------------------------------------------------------------------------"
1285 WRITE (log_unit, "(1X,A)") &
1286 "Only KERNEL NONE is active for k-point TDDFPT; transition dipole magnitudes are shown."
1287 WRITE (log_unit, '(/,T10,A,T19,A,T37,A,T69,A)') "State", "Excitation", &
1288 "Transition dipole (a.u.)", "Oscillator"
1289 WRITE (log_unit, '(T10,A,T19,A,T37,A,T49,A,T61,A,T67,A)') "number", "energy (eV)", &
1290 "x", "y", "z", "strength (a.u.)"
1291 WRITE (log_unit, '(T10,72("-"))')
1292 END IF
1293
1294 fsum = 0.0_dp
1295 DO istate = 1, nstates
1296 itrans = inds(istate) - 1
1297 ikp = itrans/ntrans_kpoint + 1
1298 itrans = mod(itrans, ntrans_kpoint)
1299 spin_offset = 0
1300 DO ispin = 1, nspins
1301 ntrans_spin = homo_spin(ispin)*nvirt_spin(ispin)
1302 IF (itrans < spin_offset + ntrans_spin) THEN
1303 itrans = itrans - spin_offset
1304 iocc = itrans/nvirt_spin(ispin) + 1
1305 ivirt = mod(itrans, nvirt_spin(ispin)) + homo_spin(ispin) + 1
1306 EXIT
1307 END IF
1308 spin_offset = spin_offset + ntrans_spin
1309 END DO
1310
1311 IF (log_unit > 0) THEN
1312 transition_dipole_abs(1:nderivs) = &
1313 sqrt(transition_dipole_re(inds(istate), 1:nderivs)**2 + &
1314 transition_dipole_im(inds(istate), 1:nderivs)**2)
1315 WRITE (log_unit, '(1X,A,T9,I7,T19,F11.5,T31,3(1X,ES11.4E2),T69,ES12.5E2)') &
1316 "TDDFPT|", istate, evals(istate)*evolt, transition_dipole_abs, &
1317 oscillator_strength(inds(istate))
1318 fsum = fsum + oscillator_strength(inds(istate))**2
1319 WRITE (log_unit, '(1X,A,T18,I7,T28,I7,T38,I7,T50,I7,T62,I7,T74,F10.5)') &
1320 "TDDFPT_KPOINT|", istate, ikp, ispin, iocc, ivirt, wkp(ikp)
1321 END IF
1322 END DO
1323
1324 IF (log_unit > 0) THEN
1325 WRITE (log_unit, '(/,T2,A,E16.8)') 'TDDFPT : CheckSum E = ', checksum
1326 WRITE (log_unit, '(T2,A,E16.8)') 'TDDFPT : CheckSum F = ', sqrt(fsum)
1327 WRITE (log_unit, "(1X,A)") &
1328 "-------------------------------------------------------------------------------"
1329 END IF
1330
1331 IF (use_scf_moment_dipoles) THEN
1332 DEALLOCATE (kpoint_dipole)
1333 ELSE
1334 CALL dbcsr_deallocate_matrix(rmatrix)
1335 CALL dbcsr_deallocate_matrix(cmatrix)
1336 CALL dbcsr_deallocate_matrix_set(overlap_deriv)
1337 END IF
1338
1339 DEALLOCATE (evals, inds, oscillator_strength, transition_dipole_im, transition_dipole_re, &
1340 transition_energy)
1341
1342 CALL timestop(handle)
1343
1344 END SUBROUTINE tddfpt_kpoint_independent_particle
1345
1346! **************************************************************************************************
1347!> \brief The energy calculation has been moved to its own subroutine
1348!> \param qs_env ...
1349!> \param nstates ...
1350!> \param nspins ...
1351!> \param work_matrices ...
1352!> \param tddfpt_control ...
1353!> \param logger ...
1354!> \param tddfpt_print_section ...
1355!> \param evects ...
1356!> \param evals ...
1357!> \param gs_mos ...
1358!> \param tddfpt_section ...
1359!> \param S_evects ...
1360!> \param matrix_s ...
1361!> \param kernel_env ...
1362!> \param matrix_ks ...
1363!> \param sub_env ...
1364!> \param ostrength ...
1365!> \param dipole_op_mos_occ ...
1366!> \param mult ...
1367!> \param xc_section ...
1368!> \param full_kernel_env ...
1369!> \param kernel_env_admm_aux ...
1370! **************************************************************************************************
1371 SUBROUTINE tddfpt_energies(qs_env, nstates, nspins, work_matrices, &
1372 tddfpt_control, logger, tddfpt_print_section, evects, evals, &
1373 gs_mos, tddfpt_section, S_evects, matrix_s, kernel_env, matrix_ks, &
1374 sub_env, ostrength, dipole_op_mos_occ, mult, xc_section, full_kernel_env, &
1375 kernel_env_admm_aux)
1376
1377 TYPE(qs_environment_type), POINTER :: qs_env
1378 INTEGER :: nstates, nspins
1379 TYPE(tddfpt_work_matrices) :: work_matrices
1380 TYPE(tddfpt2_control_type), POINTER :: tddfpt_control
1381 TYPE(cp_logger_type), POINTER :: logger
1382 TYPE(section_vals_type), POINTER :: tddfpt_print_section
1383 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: evects
1384 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: evals
1385 TYPE(tddfpt_ground_state_mos), DIMENSION(:), &
1386 POINTER :: gs_mos
1387 TYPE(section_vals_type), POINTER :: tddfpt_section
1388 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: s_evects
1389 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
1390 TYPE(kernel_env_type) :: kernel_env
1391 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks
1392 TYPE(tddfpt_subgroup_env_type) :: sub_env
1393 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: ostrength
1394 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: dipole_op_mos_occ
1395 INTEGER :: mult
1396 TYPE(section_vals_type), POINTER :: xc_section
1397 TYPE(full_kernel_env_type), TARGET :: full_kernel_env, kernel_env_admm_aux
1398
1399 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt_energies'
1400
1401 CHARACTER(len=20) :: nstates_str
1402 INTEGER :: energy_unit, handle, iter, log_unit, &
1403 niters, nocc, nstate_max, &
1404 nstates_read, nvirt
1405 LOGICAL :: do_admm, do_exck, do_soc, explicit
1406 REAL(kind=dp) :: conv
1407 TYPE(admm_type), POINTER :: admm_env
1408 TYPE(cp_blacs_env_type), POINTER :: blacs_env
1409 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks_oep
1410 TYPE(dft_control_type), POINTER :: dft_control
1411 TYPE(qs_rho_type), POINTER :: rho0_struct
1412 TYPE(section_vals_type), POINTER :: lri_section, namd_print_section, &
1413 soc_section
1414
1415 CALL timeset(routinen, handle)
1416
1417 CALL get_qs_env(qs_env, dft_control=dft_control)
1418
1419 NULLIFY (admm_env, matrix_ks_oep)
1420 do_admm = tddfpt_control%do_admm
1421 IF (do_admm) CALL get_qs_env(qs_env, admm_env=admm_env)
1422
1423 ! setup for full_kernel and admm_kernel within tddfpt_energies due to dependence on multiplicity
1424 IF (tddfpt_control%kernel == tddfpt_kernel_full) THEN
1425
1426 CALL tddfpt_construct_ground_state_orb_density( &
1427 rho_orb_struct=work_matrices%rho_orb_struct_sub, &
1428 rho_xc_struct=work_matrices%rho_xc_struct_sub, &
1429 is_rks_triplets=tddfpt_control%rks_triplets, &
1430 qs_env=qs_env, sub_env=sub_env, &
1431 wfm_rho_orb=work_matrices%rho_ao_orb_fm_sub)
1432 IF (dft_control%qs_control%gapw_xc) THEN
1433 rho0_struct => work_matrices%rho_xc_struct_sub
1434 ELSE
1435 rho0_struct => work_matrices%rho_orb_struct_sub
1436 END IF
1437
1438 IF (do_admm) THEN
1439 ! Full kernel with ADMM
1440 IF (tddfpt_control%admm_xc_correction) THEN
1441 CALL create_kernel_env(kernel_env=full_kernel_env, &
1442 rho_struct_sub=rho0_struct, &
1443 xc_section=admm_env%xc_section_primary, &
1444 is_rks_triplets=tddfpt_control%rks_triplets, &
1445 sub_env=sub_env, qs_env=qs_env)
1446 ELSE
1447 CALL create_kernel_env(kernel_env=full_kernel_env, &
1448 rho_struct_sub=rho0_struct, &
1449 xc_section=xc_section, &
1450 is_rks_triplets=tddfpt_control%rks_triplets, &
1451 sub_env=sub_env, qs_env=qs_env)
1452 END IF
1453
1454 CALL tddfpt_construct_aux_fit_density( &
1455 rho_orb_struct=work_matrices%rho_orb_struct_sub, &
1456 rho_aux_fit_struct=work_matrices%rho_aux_fit_struct_sub, &
1457 local_rho_set=sub_env%local_rho_set_admm, &
1458 qs_env=qs_env, sub_env=sub_env, &
1459 wfm_rho_orb=work_matrices%rho_ao_orb_fm_sub, &
1460 wfm_rho_aux_fit=work_matrices%rho_ao_aux_fit_fm_sub, &
1461 wfm_aux_orb=work_matrices%wfm_aux_orb_sub)
1462
1463 CALL create_kernel_env(kernel_env=kernel_env_admm_aux, &
1464 rho_struct_sub=work_matrices%rho_aux_fit_struct_sub, &
1465 xc_section=admm_env%xc_section_aux, &
1466 is_rks_triplets=tddfpt_control%rks_triplets, &
1467 sub_env=sub_env, qs_env=qs_env)
1468 kernel_env%full_kernel => full_kernel_env
1469 kernel_env%admm_kernel => kernel_env_admm_aux
1470 ELSE
1471 ! Full kernel
1472 CALL create_kernel_env(kernel_env=full_kernel_env, &
1473 rho_struct_sub=rho0_struct, &
1474 xc_section=xc_section, &
1475 is_rks_triplets=tddfpt_control%rks_triplets, &
1476 sub_env=sub_env, qs_env=qs_env)
1477 kernel_env%full_kernel => full_kernel_env
1478 NULLIFY (kernel_env%admm_kernel)
1479 END IF
1480 ! Fxc from kernel definition
1481 do_exck = tddfpt_control%do_exck
1482 kernel_env%full_kernel%do_exck = do_exck
1483 ! initilize xc kernel
1484 IF (do_exck) THEN
1485 CALL create_fxc_kernel(work_matrices%rho_orb_struct_sub, work_matrices%fxc_rspace_sub, &
1486 xc_section, tddfpt_control%rks_triplets, sub_env, qs_env)
1487 END IF
1488 END IF
1489
1490 ! lri input
1491 IF (tddfpt_control%do_lrigpw) THEN
1492 lri_section => section_vals_get_subs_vals(tddfpt_section, "LRIGPW")
1493 CALL tddfpt2_lri_init(qs_env, kernel_env, lri_section, &
1494 tddfpt_print_section)
1495 END IF
1496
1497 !! Too many states can lead to Problems
1498 !! You should be warned if there are more states
1499 !! than occ-virt Combinations!!
1500 CALL cp_fm_get_info(gs_mos(1)%mos_occ, ncol_global=nocc)
1501 IF (tddfpt_control%spinflip == no_sf_tddfpt) THEN
1502 CALL cp_fm_get_info(gs_mos(1)%mos_virt, ncol_global=nvirt)
1503 ELSE
1504 CALL cp_fm_get_info(gs_mos(2)%mos_virt, ncol_global=nvirt)
1505 END IF
1506 nstate_max = nocc*nvirt
1507 IF ((SIZE(gs_mos, 1) == 2) .AND. (tddfpt_control%spinflip == no_sf_tddfpt)) THEN
1508 CALL cp_fm_get_info(gs_mos(2)%mos_occ, ncol_global=nocc)
1509 CALL cp_fm_get_info(gs_mos(2)%mos_virt, ncol_global=nvirt)
1510 nstate_max = nocc*nvirt + nstate_max
1511 END IF
1512 IF (nstates > nstate_max) THEN
1513 cpwarn("NUMBER OF EXCITED STATES COULD LEAD TO PROBLEMS!")
1514 cpwarn("Experimental: CHANGED NSTATES TO ITS MAXIMUM VALUE!")
1515 nstates = nstate_max
1516 END IF
1517
1518 soc_section => section_vals_get_subs_vals(tddfpt_section, "SOC")
1519 CALL section_vals_get(soc_section, explicit=do_soc)
1520
1521 ! reuse Ritz vectors from the previous calculation if available
1522 IF (tddfpt_control%is_restart .AND. .NOT. do_soc) THEN
1523 CALL get_qs_env(qs_env, blacs_env=blacs_env)
1524
1525 nstates_read = tddfpt_read_restart( &
1526 evects=evects, &
1527 evals=evals, &
1528 gs_mos=gs_mos, &
1529 logger=logger, &
1530 tddfpt_section=tddfpt_section, &
1531 tddfpt_print_section=tddfpt_print_section, &
1532 fm_pool_ao_mo_active=work_matrices%fm_pool_ao_mo_active, &
1533 blacs_env_global=blacs_env)
1534 ELSE
1535 nstates_read = 0
1536 END IF
1537
1538 ! build the list of missed singly excited states and sort them in ascending order
1539 ! according to their excitation energies
1540 log_unit = cp_print_key_unit_nr(logger, tddfpt_print_section, &
1541 "GUESS_VECTORS", extension=".tddfptLog")
1542 CALL tddfpt_guess_vectors(evects=evects, evals=evals, &
1543 gs_mos=gs_mos, log_unit=log_unit, tddfpt_control=tddfpt_control, &
1544 fm_pool_ao_mo_active=work_matrices%fm_pool_ao_mo_active, &
1545 qs_env=qs_env, nspins=nspins)
1546 CALL cp_print_key_finished_output(log_unit, logger, &
1547 tddfpt_print_section, "GUESS_VECTORS")
1548
1549 CALL tddfpt_orthogonalize_psi1_psi0(evects, work_matrices%S_C0_C0T, qs_env, &
1550 gs_mos, evals, tddfpt_control, work_matrices%S_C0)
1551 CALL tddfpt_orthonormalize_psi1_psi1(evects, nstates, s_evects, matrix_s(1)%matrix)
1552
1553 niters = tddfpt_control%niters
1554 IF (niters > 0) THEN
1555 log_unit = cp_print_key_unit_nr(logger, tddfpt_print_section, &
1556 "ITERATION_INFO", extension=".tddfptLog")
1557 energy_unit = cp_print_key_unit_nr(logger, &
1558 tddfpt_print_section, &
1559 "DETAILED_ENERGY", &
1560 extension=".tddfptLog")
1561
1562 IF (log_unit > 0) THEN
1563 WRITE (log_unit, "(1X,A)") "", &
1564 "-------------------------------------------------------------------------------", &
1565 "- TDDFPT WAVEFUNCTION OPTIMIZATION -", &
1566 "-------------------------------------------------------------------------------"
1567
1568 WRITE (log_unit, '(/,T11,A,T27,A,T40,A,T62,A)') "Step", "Time", "Convergence", "Conv. states"
1569 WRITE (log_unit, '(1X,79("-"))')
1570 END IF
1571
1572 CALL cp_add_iter_level(logger%iter_info, "TDDFT_SCF")
1573
1574 DO
1575 ! *** perform Davidson iterations ***
1576 conv = tddfpt_davidson_solver( &
1577 evects=evects, &
1578 evals=evals, &
1579 s_evects=s_evects, &
1580 gs_mos=gs_mos, &
1581 tddfpt_control=tddfpt_control, &
1582 matrix_ks=matrix_ks, &
1583 qs_env=qs_env, &
1584 kernel_env=kernel_env, &
1585 sub_env=sub_env, &
1586 logger=logger, &
1587 iter_unit=log_unit, &
1588 energy_unit=energy_unit, &
1589 tddfpt_print_section=tddfpt_print_section, &
1590 work_matrices=work_matrices)
1591
1592 ! at this point at least one of the following conditions are met:
1593 ! a) convergence criteria has been achieved;
1594 ! b) maximum number of iterations has been reached;
1595 ! c) Davidson iterations must be restarted due to lack of Krylov vectors
1596
1597 CALL cp_iterate(logger%iter_info, increment=0, iter_nr_out=iter)
1598 ! terminate the loop if either (a) or (b) is true ...
1599 IF ((conv <= tddfpt_control%conv) .OR. iter >= niters) EXIT
1600
1601 ! ... otherwise restart Davidson iterations
1602 evals = 0.0_dp
1603 IF (log_unit > 0) THEN
1604 WRITE (log_unit, '(1X,25("-"),1X,A,1X,25("-"))') "Restart Davidson iterations"
1605 CALL m_flush(log_unit)
1606 END IF
1607 END DO
1608
1609 ! write TDDFPT restart file at the last iteration if requested to do so
1610 CALL cp_iterate(logger%iter_info, increment=0, last=.true.)
1611 CALL tddfpt_write_restart(evects=evects, &
1612 evals=evals, &
1613 gs_mos=gs_mos, &
1614 logger=logger, &
1615 tddfpt_print_section=tddfpt_print_section)
1616
1617 CALL cp_rm_iter_level(logger%iter_info, "TDDFT_SCF")
1618
1619 ! print convergence summary
1620 IF (log_unit > 0) THEN
1621 CALL integer_to_string(iter, nstates_str)
1622 IF (conv <= tddfpt_control%conv) THEN
1623 WRITE (log_unit, "(1X,A)") "", &
1624 "-------------------------------------------------------------------------------", &
1625 "- TDDFPT run converged in "//trim(nstates_str)//" iteration(s) ", &
1626 "-------------------------------------------------------------------------------"
1627 ELSE
1628 WRITE (log_unit, "(1X,A)") "", &
1629 "-------------------------------------------------------------------------------", &
1630 "- TDDFPT run did NOT converge after "//trim(nstates_str)//" iteration(s) ", &
1631 "-------------------------------------------------------------------------------"
1632 END IF
1633 END IF
1634
1635 CALL cp_print_key_finished_output(energy_unit, logger, &
1636 tddfpt_print_section, "DETAILED_ENERGY")
1637 CALL cp_print_key_finished_output(log_unit, logger, &
1638 tddfpt_print_section, "ITERATION_INFO")
1639 ELSE
1640 CALL cp_warn(__location__, &
1641 "Skipping TDDFPT wavefunction optimization")
1642 END IF
1643
1644 IF (ASSOCIATED(matrix_ks_oep)) THEN
1645 IF (tddfpt_control%dipole_form == tddfpt_dipole_velocity) THEN
1646 CALL cp_warn(__location__, &
1647 "Transition dipole moments and oscillator strengths are likely to be incorrect "// &
1648 "when computed using an orbital energy correction XC-potential together with "// &
1649 "the velocity form of dipole transition integrals")
1650 END IF
1651 END IF
1652
1653 ! *** print summary information ***
1654 log_unit = cp_logger_get_default_io_unit(logger)
1655
1656 namd_print_section => section_vals_get_subs_vals( &
1657 tddfpt_print_section, &
1658 "NAMD_PRINT")
1659 CALL section_vals_get(namd_print_section, explicit=explicit)
1660 IF (explicit) THEN
1661 CALL tddfpt_write_newtonx_output(evects, &
1662 evals, &
1663 gs_mos, &
1664 logger, &
1665 tddfpt_print_section, &
1666 matrix_s(1)%matrix, &
1667 s_evects, &
1668 sub_env)
1669 END IF
1670 ALLOCATE (ostrength(nstates))
1671 ostrength = 0.0_dp
1672 CALL tddfpt_print_summary(log_unit, &
1673 evects, &
1674 evals, &
1675 gs_mos, &
1676 ostrength, &
1677 mult, &
1678 dipole_op_mos_occ, &
1679 tddfpt_control%dipole_form)
1680 CALL tddfpt_print_excitation_analysis( &
1681 log_unit, &
1682 evects, &
1683 evals, &
1684 gs_mos, &
1685 matrix_s(1)%matrix, &
1686 tddfpt_control%spinflip, &
1687 min_amplitude=tddfpt_control%min_excitation_amplitude)
1688 CALL tddfpt_print_nto_analysis(qs_env, &
1689 evects, evals, &
1690 ostrength, &
1691 gs_mos, &
1692 matrix_s(1)%matrix, &
1693 tddfpt_print_section)
1694 IF (tddfpt_control%do_exciton_descriptors) THEN
1695 CALL tddfpt_print_exciton_descriptors( &
1696 log_unit, &
1697 evects, &
1698 gs_mos, &
1699 matrix_s(1)%matrix, &
1700 tddfpt_control%do_directional_exciton_descriptors, &
1701 tddfpt_control%do_directional_exciton_crosscorrelation, &
1702 qs_env)
1703 END IF
1704
1705 IF (tddfpt_control%do_lrigpw) THEN
1706 CALL lri_print_stat(qs_env, &
1707 ltddfpt=.true., &
1708 tddfpt_lri_env=kernel_env%full_kernel%lri_env)
1709 END IF
1710
1711 CALL timestop(handle)
1712 END SUBROUTINE tddfpt_energies
1713
1714! **************************************************************************************************
1715!> \brief Perform singlet and triplet computations for subsequent TDDFPT-SOC calculation.
1716!> \param qs_env Quickstep environment
1717!> \param nstates number of requested exited states
1718!> \param work_matrices ...
1719!> \param tddfpt_control ...
1720!> \param logger ...
1721!> \param tddfpt_print_section ...
1722!> \param evects Eigenvector of the requested multiplicity
1723!> \param evals Eigenvalue of the requested multiplicity
1724!> \param ostrength Oscillatorstrength
1725!> \param gs_mos ...
1726!> \param tddfpt_section ...
1727!> \param S_evects ...
1728!> \param matrix_s ...
1729!> \param kernel_env ...
1730!> \param matrix_ks ...
1731!> \param sub_env ...
1732!> \param dipole_op_mos_occ ...
1733!> \param lmult_tmp ...
1734!> \param xc_section ...
1735!> \param full_kernel_env ...
1736!> \param kernel_env_admm_aux ...
1737!> \par History
1738!> * 02.2023 created [Jan-Robert Vogt]
1739!> \note Based on tddfpt2_methods and xas_tdp_utils.
1740!> \note only the values of one multiplicity will be passed back for force calculations!
1741! **************************************************************************************************
1742
1743 SUBROUTINE tddfpt_soc_energies(qs_env, nstates, work_matrices, &
1744 tddfpt_control, logger, tddfpt_print_section, &
1745 evects, evals, ostrength, &
1746 gs_mos, tddfpt_section, S_evects, matrix_s, kernel_env, matrix_ks, &
1747 sub_env, dipole_op_mos_occ, lmult_tmp, xc_section, full_kernel_env, &
1748 kernel_env_admm_aux)
1749
1750 TYPE(qs_environment_type), INTENT(IN), POINTER :: qs_env
1751 INTEGER, INTENT(in) :: nstates
1752 TYPE(tddfpt_work_matrices) :: work_matrices
1753 TYPE(tddfpt2_control_type), POINTER :: tddfpt_control
1754 TYPE(cp_logger_type), POINTER :: logger
1755 TYPE(section_vals_type), POINTER :: tddfpt_print_section
1756 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: evects
1757 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: evals, ostrength
1758 TYPE(tddfpt_ground_state_mos), DIMENSION(:), &
1759 POINTER :: gs_mos
1760 TYPE(section_vals_type), POINTER :: tddfpt_section
1761 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: s_evects
1762 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
1763 TYPE(kernel_env_type) :: kernel_env
1764 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks
1765 TYPE(tddfpt_subgroup_env_type) :: sub_env
1766 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: dipole_op_mos_occ
1767 LOGICAL, INTENT(in) :: lmult_tmp
1768 TYPE(section_vals_type), POINTER :: xc_section
1769 TYPE(full_kernel_env_type), TARGET :: full_kernel_env, kernel_env_admm_aux
1770
1771 CHARACTER(LEN=*), PARAMETER :: routinen = 'tddfpt_soc_energies'
1772
1773 INTEGER :: handle, ispin, istate, log_unit, mult, &
1774 nspins
1775 LOGICAL :: do_sf
1776 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: evals_mult, ostrength_mult
1777 TYPE(cp_fm_type), ALLOCATABLE, DIMENSION(:, :) :: evects_mult
1778
1779 CALL timeset(routinen, handle)
1780
1781 log_unit = cp_print_key_unit_nr(logger, tddfpt_print_section, &
1782 "PROGRAM_BANNER", &
1783 extension=".tddfptLog")
1784 CALL tddfpt_soc_header(log_unit)
1785
1786 nspins = SIZE(gs_mos)
1787 ALLOCATE (evects_mult(nspins, nstates))
1788 ALLOCATE (evals_mult(nstates))
1789 do_sf = tddfpt_control%spinflip /= no_sf_tddfpt
1790
1791 ! First multiplicity
1792 IF (lmult_tmp) THEN
1793 IF (log_unit > 0) THEN
1794 WRITE (log_unit, "(1X,A)") "", &
1795 "-------------------------------------------------------------------------------", &
1796 "- TDDFPT SINGLET ENERGIES -", &
1797 "-------------------------------------------------------------------------------"
1798 END IF
1799 mult = 1
1800 ELSE
1801 IF (log_unit > 0) THEN
1802 WRITE (log_unit, "(1X,A)") "", &
1803 "-------------------------------------------------------------------------------", &
1804 "- TDDFPT TRIPLET ENERGIES -", &
1805 "-------------------------------------------------------------------------------"
1806 END IF
1807 mult = 3
1808 END IF
1809
1810 CALL tddfpt_energies(qs_env, nstates, nspins, work_matrices, tddfpt_control, logger, &
1811 tddfpt_print_section, evects_mult, evals_mult, &
1812 gs_mos, tddfpt_section, s_evects, matrix_s, &
1813 kernel_env, matrix_ks, sub_env, ostrength_mult, &
1814 dipole_op_mos_occ, mult, xc_section, full_kernel_env, &
1815 kernel_env_admm_aux)
1816
1817 ! Clean up in between for full kernel
1818 IF (tddfpt_control%kernel == tddfpt_kernel_full) THEN
1819 IF (tddfpt_control%do_admm) CALL release_kernel_env(kernel_env%admm_kernel)
1820 CALL release_kernel_env(kernel_env%full_kernel)
1821 CALL tddfpt_release_work_matrices(work_matrices, sub_env)
1822 CALL tddfpt_create_work_matrices(work_matrices, gs_mos, nstates, &
1823 tddfpt_control%do_hfx, &
1824 tddfpt_control%do_admm, tddfpt_control%do_hfxlr, &
1825 tddfpt_control%do_exck, do_sf, qs_env, sub_env)
1826 END IF
1827
1828 DO istate = 1, nstates
1829 DO ispin = 1, nspins
1830 CALL cp_fm_release(s_evects(ispin, istate))
1831 END DO
1832 END DO
1833
1834 DO istate = 1, nstates
1835 DO ispin = 1, nspins
1836 CALL fm_pool_create_fm( &
1837 work_matrices%fm_pool_ao_mo_active(ispin)%pool, &
1838 s_evects(ispin, istate))
1839 END DO
1840 END DO
1841
1842 tddfpt_control%rks_triplets = lmult_tmp
1843
1844 ! Second multiplicity
1845 IF (lmult_tmp) THEN
1846 IF (log_unit > 0) THEN
1847 WRITE (log_unit, "(1X,A)") "", &
1848 " singlet excitations finished ", &
1849 " ", &
1850 "-------------------------------------------------------------------------------", &
1851 "- TDDFPT TRIPLET ENERGIES -", &
1852 "-------------------------------------------------------------------------------"
1853 END IF !log_unit
1854 mult = 3
1855 ELSE
1856 IF (log_unit > 0) THEN
1857 WRITE (log_unit, "(1X,A)") "", &
1858 " triplet excitations finished ", &
1859 " ", &
1860 "-------------------------------------------------------------------------------", &
1861 "- TDDFPT SINGLET ENERGIES -", &
1862 "-------------------------------------------------------------------------------"
1863 END IF !log_unit
1864 mult = 1
1865 END IF
1866
1867 CALL tddfpt_energies(qs_env, nstates, nspins, work_matrices, tddfpt_control, logger, &
1868 tddfpt_print_section, evects, evals, &
1869 gs_mos, tddfpt_section, s_evects, matrix_s, &
1870 kernel_env, matrix_ks, sub_env, ostrength, &
1871 dipole_op_mos_occ, mult, xc_section, full_kernel_env, &
1872 kernel_env_admm_aux)
1873
1874 ! Compute perturbative SOC correction
1875 ! Order should always be singlet triplet in tddfpt_soc
1876 IF (lmult_tmp) THEN
1877 CALL tddfpt_soc(qs_env, evals_mult, evals, evects_mult, evects, gs_mos) !mult=singlet
1878 ELSE
1879 CALL tddfpt_soc(qs_env, evals, evals_mult, evects, evects_mult, gs_mos) !mult=triplet
1880 END IF
1881
1882 ! deallocate the additional multiplicity
1883 DO ispin = 1, SIZE(evects_mult, 1)
1884 DO istate = 1, SIZE(evects_mult, 2)
1885 CALL cp_fm_release(evects_mult(ispin, istate))
1886 END DO
1887 END DO
1888 DEALLOCATE (evects_mult, evals_mult, ostrength_mult)
1889
1890 CALL timestop(handle)
1891
1892 END SUBROUTINE tddfpt_soc_energies
1893
1894! **************************************************************************************************
1895!> \brief ...
1896!> \param qs_env ...
1897!> \param gs_mos ...
1898!> \param tddfpt_control ...
1899!> \param tddfpt_section ...
1900!> \param iounit ...
1901! **************************************************************************************************
1902 SUBROUTINE init_res_method(qs_env, gs_mos, tddfpt_control, tddfpt_section, iounit)
1903
1904 TYPE(qs_environment_type), POINTER :: qs_env
1905 TYPE(tddfpt_ground_state_mos), DIMENSION(:), &
1906 POINTER :: gs_mos
1907 TYPE(tddfpt2_control_type), POINTER :: tddfpt_control
1908 TYPE(section_vals_type), POINTER :: tddfpt_section
1909 INTEGER, INTENT(IN) :: iounit
1910
1911 CHARACTER(LEN=*), PARAMETER :: routinen = 'init_res_method'
1912
1913 INTEGER :: handle, i, io, ispin, nao, nmo, nmol, &
1914 nspins
1915 INTEGER, ALLOCATABLE, DIMENSION(:, :) :: orblist
1916 INTEGER, DIMENSION(:), POINTER :: mollist
1917 LOGICAL :: do_res, do_sf, ew1, ew2, ew3, ewcut, lms
1918 REAL(kind=dp) :: eclow, ecup, eint, emo
1919 REAL(kind=dp), DIMENSION(:), POINTER :: rvint
1920 TYPE(cp_blacs_env_type), POINTER :: blacs_env
1921 TYPE(cp_fm_struct_type), POINTER :: fm_struct
1922 TYPE(section_vals_type), POINTER :: res_section
1923
1924 CALL timeset(routinen, handle)
1925
1926 res_section => section_vals_get_subs_vals(tddfpt_section, "REDUCED_EXCITATION_SPACE")
1927 CALL section_vals_val_get(res_section, "_SECTION_PARAMETERS_", l_val=do_res)
1928
1929 ! spin flip TDA
1930 IF (tddfpt_control%spinflip == no_sf_tddfpt) THEN
1931 do_sf = .false.
1932 ELSE
1933 do_sf = .true.
1934 END IF
1935
1936 nspins = SIZE(gs_mos)
1937 IF (.NOT. do_res) THEN
1938 DO ispin = 1, nspins
1939 nmo = gs_mos(ispin)%nmo_occ
1940 tddfpt_control%nactive(ispin) = nmo
1941 gs_mos(ispin)%nmo_active = nmo
1942 ALLOCATE (gs_mos(ispin)%index_active(nmo))
1943 DO i = 1, nmo
1944 gs_mos(ispin)%index_active(i) = i
1945 END DO
1946 END DO
1947 ELSE
1948 IF (iounit > 0) THEN
1949 WRITE (iounit, "(/,1X,27('='),A,26('='))") ' REDUCED EXCITATION SPACE '
1950 END IF
1951 CALL section_vals_val_get(res_section, "ENERGY_WINDOW", explicit=ew1)
1952 CALL section_vals_val_get(res_section, "UPPER_ENERGY_CUTOFF", explicit=ew2)
1953 CALL section_vals_val_get(res_section, "LOWER_ENERGY_CUTOFF", explicit=ew3)
1954 ewcut = (ew1 .OR. ew2 .OR. ew3)
1955 CALL section_vals_val_get(res_section, "MOLECULE_LIST", explicit=lms)
1956
1957 CALL section_vals_val_get(res_section, "ENERGY_WINDOW", r_vals=rvint)
1958 cpassert(SIZE(rvint) == 2)
1959 eclow = rvint(1)
1960 ecup = rvint(2)
1961 CALL section_vals_val_get(res_section, "UPPER_ENERGY_CUTOFF", r_val=eint)
1962 ecup = min(ecup, eint)
1963 CALL section_vals_val_get(res_section, "LOWER_ENERGY_CUTOFF", r_val=eint)
1964 eclow = max(eclow, eint)
1965 IF (ewcut .AND. (iounit > 0)) THEN
1966 IF (eclow < -1.e8_dp .AND. ecup > 1.e8_dp) THEN
1967 WRITE (iounit, "(1X,A,T51,A10,T71,A10)") &
1968 'Orbital Energy Window [eV]', " -Inf", " Inf"
1969 ELSE IF (eclow < -1.e8_dp) THEN
1970 WRITE (iounit, "(1X,A,T51,A10,T71,F10.4)") &
1971 'Orbital Energy Window [eV]', " -Inf", evolt*ecup
1972 ELSE IF (ecup > 1.e8_dp) THEN
1973 WRITE (iounit, "(1X,A,T51,F10.4,T71,A10)") &
1974 'Orbital Energy Window [eV]', evolt*eclow, " Inf"
1975 ELSE
1976 WRITE (iounit, "(1X,A,T51,F10.4,T71,F10.4)") &
1977 'Orbital Energy Window [eV]', evolt*eclow, evolt*ecup
1978 END IF
1979 END IF
1980
1981 nmol = 0
1982 IF (lms) THEN
1983 CALL section_vals_val_get(res_section, "MOLECULE_LIST", i_vals=mollist)
1984 nmol = SIZE(mollist)
1985 WRITE (iounit, "(1X,A)") 'List of Selected Molecules'
1986 WRITE (iounit, "(1X,15(I5))") mollist(1:nmol)
1987 END IF
1988
1989 DO ispin = 1, nspins
1990 tddfpt_control%nactive(ispin) = gs_mos(ispin)%nmo_occ
1991 END DO
1992 nmo = maxval(tddfpt_control%nactive)
1993 ALLOCATE (orblist(nmo, nspins))
1994 orblist = 0
1995
1996 IF (lms) THEN
1997 ! ignore for now
1998 orblist = 1
1999 DO ispin = 1, nspins
2000 cpassert(.NOT. ASSOCIATED(gs_mos(ispin)%evals_occ_matrix))
2001 END DO
2002 ELSE IF (ewcut) THEN
2003 ! Filter orbitals wrt energy window
2004 DO ispin = 1, nspins
2005 DO i = 1, gs_mos(ispin)%nmo_occ
2006 emo = gs_mos(ispin)%evals_occ(i)
2007 IF (emo > eclow .AND. emo < ecup) orblist(i, ispin) = 1
2008 END DO
2009 END DO
2010 ELSE
2011 orblist = 1
2012 END IF
2013
2014 ! count active orbitals
2015 nmo = sum(orblist)
2016 IF (nmo == 0) THEN
2017 cpabort("RSE TDA: no active orbitals selected.")
2018 END IF
2019 DO ispin = 1, nspins
2020 nmo = sum(orblist(:, ispin))
2021 tddfpt_control%nactive(ispin) = nmo
2022 gs_mos(ispin)%nmo_active = nmo
2023 ALLOCATE (gs_mos(ispin)%index_active(nmo))
2024 io = 0
2025 DO i = 1, SIZE(orblist, 1)
2026 IF (orblist(i, ispin) == 1) THEN
2027 io = io + 1
2028 gs_mos(ispin)%index_active(io) = i
2029 END IF
2030 END DO
2031 END DO
2032 DEALLOCATE (orblist)
2033
2034 IF (lms) THEN
2035 ! output information
2036 ELSE
2037 IF (iounit > 0) THEN
2038 WRITE (iounit, "(1X,A)") 'List of Selected States'
2039 IF (nspins == 1) THEN
2040 WRITE (iounit, "(A,T67,A)") ' Active State Orbital', 'Orbital Energy'
2041 DO i = 1, gs_mos(1)%nmo_active
2042 io = gs_mos(1)%index_active(i)
2043 WRITE (iounit, "(T8,I6,T21,I6,T61,F20.4)") i, io, gs_mos(1)%evals_occ(io)*evolt
2044 END DO
2045 ELSE
2046 DO ispin = 1, nspins
2047 WRITE (iounit, "(1X,A,I2)") 'Spin ', ispin
2048 WRITE (iounit, "(A,T67,A)") ' Active State Orbital', 'Orbital Energy'
2049 DO i = 1, gs_mos(ispin)%nmo_active
2050 io = gs_mos(ispin)%index_active(i)
2051 WRITE (iounit, "(T8,I6,T21,I6,T61,F20.4)") i, io, gs_mos(ispin)%evals_occ(io)*evolt
2052 END DO
2053 END DO
2054 END IF
2055 END IF
2056 END IF
2057
2058 IF (do_sf) THEN
2059 cpabort("Restricted Active Space with spin flip TDA NYA")
2060 END IF
2061
2062 IF (iounit > 0) THEN
2063 WRITE (iounit, "(1X,79('='))")
2064 END IF
2065 END IF
2066
2067 ! Allocate mos_active
2068 DO ispin = 1, nspins
2069 CALL get_qs_env(qs_env, blacs_env=blacs_env)
2070 CALL cp_fm_get_info(gs_mos(ispin)%mos_occ, nrow_global=nao)
2071 nmo = gs_mos(ispin)%nmo_active
2072 CALL cp_fm_struct_create(fm_struct, template_fmstruct=gs_mos(ispin)%mos_occ%matrix_struct, &
2073 ncol_global=nmo, context=blacs_env)
2074 NULLIFY (gs_mos(ispin)%mos_active)
2075 ALLOCATE (gs_mos(ispin)%mos_active)
2076 CALL cp_fm_create(gs_mos(ispin)%mos_active, fm_struct)
2077 CALL cp_fm_struct_release(fm_struct)
2078 ! copy the active orbitals
2079 IF (gs_mos(ispin)%nmo_active == gs_mos(ispin)%nmo_occ) THEN
2080 DO i = 1, gs_mos(ispin)%nmo_active
2081 cpassert(i == gs_mos(ispin)%index_active(i))
2082 END DO
2083 CALL cp_fm_to_fm_submat(gs_mos(ispin)%mos_occ, gs_mos(ispin)%mos_active, &
2084 nao, nmo, 1, 1, 1, 1)
2085 ELSE
2086 DO i = 1, gs_mos(ispin)%nmo_active
2087 io = gs_mos(ispin)%index_active(i)
2088 CALL cp_fm_to_fm_submat(gs_mos(ispin)%mos_occ, gs_mos(ispin)%mos_active, &
2089 nao, 1, 1, io, 1, i)
2090 END DO
2091 END IF
2092 END DO
2093
2094 CALL timestop(handle)
2095
2096 END SUBROUTINE init_res_method
2097
2098END MODULE qs_tddfpt2_methods
cp_dbcsr_api::dbcsr_create
Definition cp_dbcsr_api.F:194
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_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:83
parallel_gemm_api::parallel_gemm
Definition parallel_gemm_api.F:39
admm_methods
Contains ADMM methods which require molecular orbitals.
Definition admm_methods.F:15
admm_methods::admm_fit_mo_coeffs
subroutine, public admm_fit_mo_coeffs(admm_env, matrix_s_aux_fit, matrix_s_mixed, mos, mos_aux_fit, geometry_did_change)
...
Definition admm_methods.F:871
admm_types
Types and set/get functions for auxiliary density matrix methods.
Definition admm_types.F:15
admm_types::get_admm_env
subroutine, public get_admm_env(admm_env, mo_derivs_aux_fit, mos_aux_fit, sab_aux_fit, sab_aux_fit_asymm, sab_aux_fit_vs_orb, matrix_s_aux_fit, matrix_s_aux_fit_kp, matrix_s_aux_fit_vs_orb, matrix_s_aux_fit_vs_orb_kp, task_list_aux_fit, matrix_ks_aux_fit, matrix_ks_aux_fit_kp, matrix_ks_aux_fit_im, matrix_ks_aux_fit_dft, matrix_ks_aux_fit_hfx, matrix_ks_aux_fit_dft_kp, matrix_ks_aux_fit_hfx_kp, rho_aux_fit, rho_aux_fit_buffer, admm_dm)
Get routine for the ADMM env.
Definition admm_types.F:599
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
bibliography
collects all references to literature in CP2K as new algorithms / method are included from literature...
Definition bibliography.F:21
bibliography::grimme2016
integer, save, public grimme2016
Definition bibliography.F:36
bibliography::iannuzzi2005
integer, save, public iannuzzi2005
Definition bibliography.F:36
bibliography::hernandez2025
integer, save, public hernandez2025
Definition bibliography.F:36
bibliography::grimme2013
integer, save, public grimme2013
Definition bibliography.F:36
cell_types
Handles all functions related to the CELL.
Definition cell_types.F:15
cp_blacs_env
methods related to the blacs parallel environment
Definition cp_blacs_env.F:15
cp_control_types
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
Definition cp_control_types.F:12
cp_dbcsr_api
Definition cp_dbcsr_api.F:8
cp_dbcsr_api::dbcsr_deallocate_matrix
subroutine, public dbcsr_deallocate_matrix(matrix)
...
Definition cp_dbcsr_api.F:235
cp_dbcsr_api::dbcsr_set
subroutine, public dbcsr_set(matrix, alpha)
...
Definition cp_dbcsr_api.F:1204
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_pool_types
pool for for elements that are retained and released
Definition cp_fm_pool_types.F:14
cp_fm_pool_types::fm_pool_create_fm
subroutine, public fm_pool_create_fm(pool, element, name)
returns an element, allocating it if none is in the pool
Definition cp_fm_pool_types.F:185
cp_fm_struct
represent the structure of a full matrix
Definition cp_fm_struct.F:14
cp_fm_struct::cp_fm_struct_create
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
Definition cp_fm_struct.F:132
cp_fm_struct::cp_fm_struct_release
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition cp_fm_struct.F:355
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_copy_general
subroutine, public cp_fm_copy_general(source, destination, para_env)
General copy of a fm matrix to another fm matrix. Uses non-blocking MPI rather than ScaLAPACK.
Definition cp_fm_types.F:1524
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_get_element
subroutine, public cp_fm_get_element(matrix, irow_global, icol_global, alpha, local)
returns an element of a fm this value is valid on every cpu using this call is expensive
Definition cp_fm_types.F:627
cp_fm_types::cp_fm_to_fm_submat
subroutine, public cp_fm_to_fm_submat(msource, mtarget, nrow, ncol, s_firstrow, s_firstcol, t_firstrow, t_firstcol)
copy just a part ot the matrix
Definition cp_fm_types.F:1455
cp_fm_types::cp_fm_create
subroutine, public cp_fm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
creates a new full matrix with the given structure
Definition cp_fm_types.F:165
cp_log_handling
various routines to log and control the output. The idea is that decisions about where to log should ...
Definition cp_log_handling.F:41
cp_log_handling::cp_logger_get_default_io_unit
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
Definition cp_log_handling.F:523
cp_log_handling::cp_get_default_logger
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
Definition cp_log_handling.F:234
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_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_iterate
subroutine, public cp_iterate(iteration_info, last, iter_nr, increment, iter_nr_out)
adds one to the actual iteration
Definition cp_output_handling.F:574
cp_output_handling::cp_rm_iter_level
subroutine, public cp_rm_iter_level(iteration_info, level_name, n_rlevel_att)
Removes an iteration level.
Definition cp_output_handling.F:705
cp_output_handling::cp_add_iter_level
subroutine, public cp_add_iter_level(iteration_info, level_name, n_rlevel_new)
Adds an iteration level.
Definition cp_output_handling.F:662
exstates_types
Types for excited states potential energies.
Definition exstates_types.F:14
header::tddfpt_header
subroutine, public tddfpt_header(iw)
...
Definition header.F:122
header::tddfpt_soc_header
subroutine, public tddfpt_soc_header(iw)
...
Definition header.F:143
hfx_admm_utils
Utilities for hfx and admm methods.
Definition hfx_admm_utils.F:17
hfx_admm_utils::aux_admm_init
subroutine, public aux_admm_init(qs_env, mos, admm_env, admm_control, basis_type)
Minimal setup routine for admm_env No forces No k-points No DFT correction terms.
Definition hfx_admm_utils.F:405
hfx_types
Types and set/get functions for HFX.
Definition hfx_types.F:16
hfx_types::hfx_create
subroutine, public hfx_create(x_data, para_env, hfx_section, atomic_kind_set, qs_kind_set, particle_set, dft_control, cell, orb_basis, ri_basis, nelectron_total, nkp_grid)
This routine allocates and initializes all types in hfx_data
Definition hfx_types.F:607
hfx_types::compare_hfx_sections
subroutine, public compare_hfx_sections(hfx_section1, hfx_section2, is_identical, same_except_frac)
Compares the non-technical parts of two HFX input section and check whether they are the same Ignore ...
Definition hfx_types.F:3036
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:18
input_constants::do_admm_purify_none
integer, parameter, public do_admm_purify_none
Definition input_constants.F:906
input_constants::tddfpt_sf_col
integer, parameter, public tddfpt_sf_col
Definition input_constants.F:626
input_constants::tddfpt_kernel_none
integer, parameter, public tddfpt_kernel_none
Definition input_constants.F:623
input_constants::tddfpt_dipole_scf_moment
integer, parameter, public tddfpt_dipole_scf_moment
Definition input_constants.F:744
input_constants::oe_none
integer, parameter, public oe_none
Definition input_constants.F:629
input_constants::do_admm_basis_projection
integer, parameter, public do_admm_basis_projection
Definition input_constants.F:915
input_constants::do_admm_aux_exch_func_none
integer, parameter, public do_admm_aux_exch_func_none
Definition input_constants.F:921
input_constants::tddfpt_dipole_velocity
integer, parameter, public tddfpt_dipole_velocity
Definition input_constants.F:744
input_constants::do_potential_truncated
integer, parameter, public do_potential_truncated
Definition input_constants.F:855
input_constants::tddfpt_kernel_full
integer, parameter, public tddfpt_kernel_full
Definition input_constants.F:623
input_constants::tddfpt_sf_noncol
integer, parameter, public tddfpt_sf_noncol
Definition input_constants.F:626
input_constants::tddfpt_kernel_stda
integer, parameter, public tddfpt_kernel_stda
Definition input_constants.F:623
input_constants::no_sf_tddfpt
integer, parameter, public no_sf_tddfpt
Definition input_constants.F:626
input_constants::do_admm_exch_scaling_none
integer, parameter, public do_admm_exch_scaling_none
Definition input_constants.F:933
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_vals_val_set
subroutine, public section_vals_val_set(section_vals, keyword_name, i_rep_section, i_rep_val, val, l_val, i_val, r_val, c_val, l_vals_ptr, i_vals_ptr, r_vals_ptr, c_vals_ptr)
sets the requested value
Definition input_section_types.F:1233
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
input_section_types::section_vals_get
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
Definition input_section_types.F:708
input_section_types::section_vals_val_get
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Definition input_section_types.F:1052
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kpoint_methods
Routines needed for kpoint calculation.
Definition kpoint_methods.F:15
kpoint_methods::rskp_transform
subroutine, public rskp_transform(rmatrix, cmatrix, rsmat, ispin, xkp, cell_to_index, sab_nl, is_complex, rs_sign)
Transformation of real space matrices to a kpoint.
Definition kpoint_methods.F:1175
kpoint_types
Types and basic routines needed for a kpoint calculation.
Definition kpoint_types.F:15
kpoint_types::get_kpoint_info
subroutine, public get_kpoint_info(kpoint, kp_scheme, nkp_grid, kp_shift, symmetry, verbose, full_grid, use_real_wfn, eps_geo, parallel_group_size, kp_range, nkp, xkp, wkp, para_env, blacs_env_all, para_env_kp, para_env_inter_kp, blacs_env, kp_env, kp_aux_env, mpools, iogrp, nkp_groups, kp_dist, cell_to_index, index_to_cell, sab_nl, sab_nl_nosym, inversion_symmetry_only, symmetry_backend, symmetry_reduction_method, gamma_centered, lattice_fft)
Retrieve information from a kpoint environment.
Definition kpoint_types.F:544
lri_environment_methods
Calculates integral matrices for LRIGPW method lri : local resolution of the identity.
Definition lri_environment_methods.F:18
lri_environment_methods::lri_print_stat
subroutine, public lri_print_stat(qs_env, ltddfpt, tddfpt_lri_env)
...
Definition lri_environment_methods.F:1248
lri_environment_types
contains the types and subroutines for dealing with the lri_env lri : local resolution of the identit...
Definition lri_environment_types.F:18
lri_environment_types::lri_density_release
subroutine, public lri_density_release(lri_density)
releases the given lri_density
Definition lri_environment_types.F:554
lri_environment_types::lri_env_release
subroutine, public lri_env_release(lri_env)
releases the given lri_env
Definition lri_environment_types.F:433
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
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:27
min_basis_set
generate or use from input minimal basis set
Definition min_basis_set.F:14
min_basis_set::create_minbas_set
subroutine, public create_minbas_set(qs_env, unit_nr, basis_type, primitive)
...
Definition min_basis_set.F:51
parallel_gemm_api
basic linear algebra operations for full matrixes
Definition parallel_gemm_api.F:14
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::evolt
real(kind=dp), parameter, public evolt
Definition physcon.F:183
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_kernel_methods
Definition qs_kernel_methods.F:8
qs_kernel_methods::create_kernel_env
subroutine, public create_kernel_env(kernel_env, xc_section, is_rks_triplets, rho_struct_sub, sub_env, qs_env)
Create kernel environment.
Definition qs_kernel_methods.F:61
qs_kernel_methods::create_fxc_kernel
subroutine, public create_fxc_kernel(rho_struct, fxc_rspace, xc_section, is_rks_triplets, sub_env, qs_env)
Create the xc kernel potential for the approximate Fxc kernel model.
Definition qs_kernel_methods.F:143
qs_kernel_types
Definition qs_kernel_types.F:8
qs_kernel_types::release_kernel_env
subroutine, public release_kernel_env(kernel_env)
Release kernel environment.
Definition qs_kernel_types.F:93
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_ks_types
Definition qs_ks_types.F:15
qs_mo_types
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
qs_mo_types::get_mo_set
subroutine, public get_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, mo_coeff, mo_coeff_b, uniform_occupation, kts, mu, flexible_electron_count)
Get the components of a MO set data structure.
Definition qs_mo_types.F:401
qs_moments
Calculates the moment integrals <a|r^m|b> and <a|r x d/dr|b>.
Definition qs_moments.F:14
qs_moments::qs_moment_kpoints_scf_mos
subroutine, public qs_moment_kpoints_scf_mos(qs_env, dipole, rcc, nmo_spin_out)
Calculates interband k-point dipoles in the existing SCF MO basis.
Definition qs_moments.F:3183
qs_neighbor_list_types
Define the neighbor list data types and the corresponding functionality.
Definition qs_neighbor_list_types.F:17
qs_overlap
Calculation of overlap matrix, its derivatives and forces.
Definition qs_overlap.F:19
qs_overlap::build_overlap_matrix
subroutine, public build_overlap_matrix(ks_env, matrix_s, matrixkp_s, matrix_name, nderivative, basis_type_a, basis_type_b, sab_nl, calculate_forces, matrix_p, matrixkp_p, ext_kpoints)
Calculation of the overlap matrix over Cartesian Gaussian functions.
Definition qs_overlap.F:121
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_scf_methods
groups fairly general SCF methods, so that modules other than qs_scf can use them too split off from ...
Definition qs_scf_methods.F:19
qs_scf_methods::eigensolver
subroutine, public eigensolver(matrix_ks_fm, mo_set, ortho, work, cholesky_method, do_level_shift, level_shift, matrix_u_fm, use_jacobi)
Diagonalise the Kohn-Sham matrix to get a new set of MO eigen- vectors and MO eigenvalues....
Definition qs_scf_methods.F:167
qs_scf_types
module that contains the definitions of the scf types
Definition qs_scf_types.F:14
qs_tddfpt2_assign
Definition qs_tddfpt2_assign.F:8
qs_tddfpt2_assign::assign_state
subroutine, public assign_state(qs_env, matrix_s, evects, psi0, wfn_history, my_state)
...
Definition qs_tddfpt2_assign.F:56
qs_tddfpt2_densities
Definition qs_tddfpt2_densities.F:8
qs_tddfpt2_densities::tddfpt_construct_aux_fit_density
subroutine, public tddfpt_construct_aux_fit_density(rho_orb_struct, rho_aux_fit_struct, local_rho_set, qs_env, sub_env, wfm_rho_orb, wfm_rho_aux_fit, wfm_aux_orb)
Project a charge density expressed in primary basis set into the auxiliary basis set.
Definition qs_tddfpt2_densities.F:162
qs_tddfpt2_densities::tddfpt_construct_ground_state_orb_density
subroutine, public tddfpt_construct_ground_state_orb_density(rho_orb_struct, rho_xc_struct, is_rks_triplets, qs_env, sub_env, wfm_rho_orb)
Compute the ground-state charge density expressed in primary basis set.
Definition qs_tddfpt2_densities.F:72
qs_tddfpt2_eigensolver
Definition qs_tddfpt2_eigensolver.F:8
qs_tddfpt2_eigensolver::tddfpt_orthonormalize_psi1_psi1
subroutine, public tddfpt_orthonormalize_psi1_psi1(evects, nvects_new, s_evects, matrix_s)
Make new TDDFPT trial vectors orthonormal to all previous trial vectors.
Definition qs_tddfpt2_eigensolver.F:280
qs_tddfpt2_eigensolver::tddfpt_orthogonalize_psi1_psi0
subroutine, public tddfpt_orthogonalize_psi1_psi0(evects, s_c0_c0t, qs_env, gs_mos, evals, tddfpt_control, s_c0)
Make TDDFPT trial vectors orthogonal to all occupied molecular orbitals.
Definition qs_tddfpt2_eigensolver.F:107
qs_tddfpt2_eigensolver::tddfpt_davidson_solver
real(kind=dp) function, public tddfpt_davidson_solver(evects, evals, s_evects, gs_mos, tddfpt_control, matrix_ks, qs_env, kernel_env, sub_env, logger, iter_unit, energy_unit, tddfpt_print_section, work_matrices)
Perform Davidson iterations.
Definition qs_tddfpt2_eigensolver.F:830
qs_tddfpt2_forces
Definition qs_tddfpt2_forces.F:8
qs_tddfpt2_forces::tddfpt_forces_main
subroutine, public tddfpt_forces_main(qs_env, gs_mos, ex_env, kernel_env, sub_env, work_matrices)
Perform TDDFPT gradient calculation. This routine calculates the response vector R of Eq....
Definition qs_tddfpt2_forces.F:163
qs_tddfpt2_fprint
Definition qs_tddfpt2_fprint.F:8
qs_tddfpt2_fprint::tddfpt_print_forces
subroutine, public tddfpt_print_forces(qs_env, evects, evals, ostrength, print_section, gs_mos, kernel_env, sub_env, work_matrices)
Calculate and print forces of selected excited states.
Definition qs_tddfpt2_fprint.F:102
qs_tddfpt2_lri_utils
Definition qs_tddfpt2_lri_utils.F:8
qs_tddfpt2_lri_utils::tddfpt2_lri_init
subroutine, public tddfpt2_lri_init(qs_env, kernel_env, lri_section, tddfpt_print_section)
Initialize LRI environment, basis, neighborlists and matrices.
Definition qs_tddfpt2_lri_utils.F:79
qs_tddfpt2_methods
Definition qs_tddfpt2_methods.F:8
qs_tddfpt2_methods::tddfpt
subroutine, public tddfpt(qs_env, calc_forces, rixs_env)
Perform TDDFPT calculation. If calc_forces then it also builds the response vector for the Z-vector m...
Definition qs_tddfpt2_methods.F:176
qs_tddfpt2_methods::tddfpt_input
subroutine, public tddfpt_input(qs_env, tddfpt_section, tddfpt_control, do_hfx, do_admm, do_exck, do_hfxsr, do_hfxlr, xc_section, tddfpt_print_section, lri_section, hfxsr_section)
TDDFPT input.
Definition qs_tddfpt2_methods.F:668
qs_tddfpt2_methods::tddfpt_energies
subroutine, public tddfpt_energies(qs_env, nstates, nspins, work_matrices, tddfpt_control, logger, tddfpt_print_section, evects, evals, gs_mos, tddfpt_section, s_evects, matrix_s, kernel_env, matrix_ks, sub_env, ostrength, dipole_op_mos_occ, mult, xc_section, full_kernel_env, kernel_env_admm_aux)
The energy calculation has been moved to its own subroutine.
Definition qs_tddfpt2_methods.F:1376
qs_tddfpt2_properties
Definition qs_tddfpt2_properties.F:8
qs_tddfpt2_properties::tddfpt_dipole_operator
subroutine, public tddfpt_dipole_operator(dipole_op_mos_occ, tddfpt_control, gs_mos, qs_env)
Compute the action of the dipole operator on the ground state wave function.
Definition qs_tddfpt2_properties.F:154
qs_tddfpt2_properties::tddfpt_print_nto_analysis
subroutine, public tddfpt_print_nto_analysis(qs_env, evects, evals, ostrength, gs_mos, matrix_s, print_section)
Print natural transition orbital analysis.
Definition qs_tddfpt2_properties.F:970
qs_tddfpt2_properties::tddfpt_print_exciton_descriptors
subroutine, public tddfpt_print_exciton_descriptors(log_unit, evects, gs_mos, matrix_s, do_directional_exciton_descriptors, do_directional_exciton_crosscorrelation, qs_env)
Print exciton descriptors, cf. Mewes et al., JCTC 14, 710-725 (2018).
Definition qs_tddfpt2_properties.F:1236
qs_tddfpt2_properties::tddfpt_print_excitation_analysis
subroutine, public tddfpt_print_excitation_analysis(log_unit, evects, evals, gs_mos, matrix_s, spinflip, min_amplitude)
Print excitation analysis.
Definition qs_tddfpt2_properties.F:672
qs_tddfpt2_properties::tddfpt_print_summary
subroutine, public tddfpt_print_summary(log_unit, evects, evals, gs_mos, ostrength, mult, dipole_op_mos_occ, dipole_form)
Print final TDDFPT excitation energies and oscillator strengths.
Definition qs_tddfpt2_properties.F:528
qs_tddfpt2_restart
Definition qs_tddfpt2_restart.F:8
qs_tddfpt2_restart::tddfpt_write_newtonx_output
subroutine, public tddfpt_write_newtonx_output(evects, evals, gs_mos, logger, tddfpt_print_section, matrix_s, s_evects, sub_env)
Write Ritz vectors to a binary restart file.
Definition qs_tddfpt2_restart.F:318
qs_tddfpt2_restart::tddfpt_write_restart
subroutine, public tddfpt_write_restart(evects, evals, gs_mos, logger, tddfpt_print_section)
Write Ritz vectors to a binary restart file.
Definition qs_tddfpt2_restart.F:78
qs_tddfpt2_restart::tddfpt_read_restart
integer function, public tddfpt_read_restart(evects, evals, gs_mos, logger, tddfpt_section, tddfpt_print_section, fm_pool_ao_mo_active, blacs_env_global)
Initialise initial guess vectors by reading (un-normalised) Ritz vectors from a binary restart file.
Definition qs_tddfpt2_restart.F:162
qs_tddfpt2_smearing_methods
Definition qs_tddfpt2_smearing_methods.F:8
qs_tddfpt2_smearing_methods::tddfpt_smeared_occupation
subroutine, public tddfpt_smeared_occupation(qs_env, gs_mos, log_unit)
...
Definition qs_tddfpt2_smearing_methods.F:66
qs_tddfpt2_soc
Definition qs_tddfpt2_soc.F:7
qs_tddfpt2_soc::tddfpt_soc
subroutine, public tddfpt_soc(qs_env, evals_a, evals_b, evects_a, evects_b, gs_mos)
Perform TDDFPT-SOC calculation.
Definition qs_tddfpt2_soc.F:140
qs_tddfpt2_stda_types
Simplified Tamm Dancoff approach (sTDA).
Definition qs_tddfpt2_stda_types.F:10
qs_tddfpt2_stda_types::stda_init_param
subroutine, public stda_init_param(qs_env, stda_kernel, stda_control)
Get the parameters needed for an sTDA calculation.
Definition qs_tddfpt2_stda_types.F:119
qs_tddfpt2_stda_types::deallocate_stda_env
subroutine, public deallocate_stda_env(stda_kernel)
Deallocate the sTDA environment.
Definition qs_tddfpt2_stda_types.F:255
qs_tddfpt2_stda_types::allocate_stda_env
subroutine, public allocate_stda_env(qs_env, stda_kernel, n_ao, nactive)
Allocate the sTDA environment.
Definition qs_tddfpt2_stda_types.F:220
qs_tddfpt2_stda_utils
Simplified Tamm Dancoff approach (sTDA).
Definition qs_tddfpt2_stda_utils.F:10
qs_tddfpt2_stda_utils::get_lowdin_mo_coefficients
subroutine, public get_lowdin_mo_coefficients(qs_env, sub_env, work)
Calculate Lowdin MO coefficients.
Definition qs_tddfpt2_stda_utils.F:330
qs_tddfpt2_stda_utils::stda_init_matrices
subroutine, public stda_init_matrices(qs_env, stda_kernel, sub_env, work, tddfpt_control)
Calculate sTDA matrices.
Definition qs_tddfpt2_stda_utils.F:103
qs_tddfpt2_subgroups
Definition qs_tddfpt2_subgroups.F:8
qs_tddfpt2_subgroups::tddfpt_sub_env_init
subroutine, public tddfpt_sub_env_init(sub_env, qs_env, mos_occ, mos_active, kernel)
Split MPI communicator to create a set of parallel (sub)groups.
Definition qs_tddfpt2_subgroups.F:185
qs_tddfpt2_subgroups::tddfpt_sub_env_release
subroutine, public tddfpt_sub_env_release(sub_env)
Release parallel group environment.
Definition qs_tddfpt2_subgroups.F:444
qs_tddfpt2_types
Definition qs_tddfpt2_types.F:8
qs_tddfpt2_types::hfxsr_create_work_matrices
subroutine, public hfxsr_create_work_matrices(work_matrices, qs_env, admm_env)
Allocate work matrices for hfxsr.
Definition qs_tddfpt2_types.F:571
qs_tddfpt2_types::tddfpt_create_work_matrices
subroutine, public tddfpt_create_work_matrices(work_matrices, gs_mos, nstates, do_hfx, do_admm, do_hfxlr, do_exck, do_sf, qs_env, sub_env)
Allocate work matrices for full kernel.
Definition qs_tddfpt2_types.F:232
qs_tddfpt2_types::tddfpt_release_work_matrices
subroutine, public tddfpt_release_work_matrices(work_matrices, sub_env)
Release work matrices.
Definition qs_tddfpt2_types.F:783
qs_tddfpt2_types::stda_create_work_matrices
subroutine, public stda_create_work_matrices(work_matrices, gs_mos, nstates, qs_env, sub_env)
Allocate work matrices for sTDA kernel.
Definition qs_tddfpt2_types.F:650
qs_tddfpt2_utils
Definition qs_tddfpt2_utils.F:8
qs_tddfpt2_utils::tddfpt_release_ground_state_mos
subroutine, public tddfpt_release_ground_state_mos(gs_mos)
Release molecular orbitals.
Definition qs_tddfpt2_utils.F:571
qs_tddfpt2_utils::tddfpt_oecorr
subroutine, public tddfpt_oecorr(qs_env, gs_mos, matrix_ks_oep)
Callculate orbital corrected KS matrix for TDDFPT.
Definition qs_tddfpt2_utils.F:630
qs_tddfpt2_utils::tddfpt_guess_vectors
subroutine, public tddfpt_guess_vectors(evects, evals, gs_mos, log_unit, tddfpt_control, fm_pool_ao_mo_active, qs_env, nspins)
Generate missed guess vectors.
Definition qs_tddfpt2_utils.F:988
qs_tddfpt2_utils::tddfpt_init_mos
subroutine, public tddfpt_init_mos(qs_env, gs_mos, iounit)
Prepare MOs for TDDFPT Calculations.
Definition qs_tddfpt2_utils.F:105
rixs_types
Define Resonant Inelastic XRAY Scattering (RIXS) control type and associated create,...
Definition rixs_types.F:13
string_utilities
Utilities for string manipulations.
Definition string_utilities.F:16
string_utilities::integer_to_string
subroutine, public integer_to_string(inumber, string)
Converts an integer number to a string. The WRITE statement will return an error message,...
Definition string_utilities.F:3239
util
All kind of helpful little routines.
Definition util.F:14
xc_write_output
Writes information on XC functionals to output.
Definition xc_write_output.F:11
xc_write_output::xc_write
subroutine, public xc_write(iounit, xc_section, lsd)
...
Definition xc_write_output.F:44
admm_types::admm_type
stores some data used in wavefunction fitting
Definition admm_types.F:120
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_blacs_env::cp_blacs_env_type
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
Definition cp_blacs_env.F:53
cp_control_types::dft_control_type
Definition cp_control_types.F:761
cp_control_types::rixs_control_type
Definition cp_control_types.F:745
cp_control_types::tddfpt2_control_type
Definition cp_control_types.F:637
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
exstates_types::excited_energy_type
Contains information on the excited states energy.
Definition exstates_types.F:53
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:127
kpoint_types::kpoint_env_p_type
Definition kpoint_types.F:108
kpoint_types::kpoint_env_type
Keeps information about a specific k-point.
Definition kpoint_types.F:90
kpoint_types::kpoint_type
Contains information about kpoints.
Definition kpoint_types.F:187
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:768
particle_types::particle_type
Definition particle_types.F:35
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:224
qs_kernel_types::full_kernel_env_type
Collection of variables required to evaluate adiabatic TDDFPT kernel.
Definition qs_kernel_types.F:47
qs_kernel_types::kernel_env_type
Type to hold environments for the different kernels.
Definition qs_kernel_types.F:78
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:190
qs_ks_types::qs_ks_env_type
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
Definition qs_ks_types.F:137
qs_mo_types::mo_set_type
Definition qs_mo_types.F:47
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_tddfpt2_stda_types::stda_env_type
Definition qs_tddfpt2_stda_types.F:45
qs_tddfpt2_subgroups::tddfpt_subgroup_env_type
Parallel (sub)group environment.
Definition qs_tddfpt2_subgroups.F:103
qs_tddfpt2_types::tddfpt_ground_state_mos
Ground state molecular orbitals.
Definition qs_tddfpt2_types.F:85
qs_tddfpt2_types::tddfpt_work_matrices
Set of temporary ("work") matrices.
Definition qs_tddfpt2_types.F:117
rixs_types::rixs_env_type
Definition rixs_types.F:29
rixs_types::tddfpt2_valence_type
Valence state coming from the qs_tddfpt2 routine.
Definition rixs_types.F:40