(git:26ffdda)
Loading...
Searching...
No Matches
bse_print.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Routines for printing information in context of the BSE calculation
10!> \par History
11!> 10.2024 created [Maximilian Graml]
12! **************************************************************************************************
13MODULE bse_print
14
15 USE atomic_kind_types, ONLY: get_atomic_kind
16 USE bibliography, ONLY: graml2026,&
17 cite_reference
18 USE bse_properties, ONLY: compute_and_print_absorption_spectrum,&
19 exciton_descr_type
20 USE bse_util, ONLY: filter_eigvec_contrib,&
21 get_bse_spin_block_layout
22 USE cp_fm_types, ONLY: cp_fm_get_info,&
23 cp_fm_type
24 USE input_constants, ONLY: bse_screening_alpha,&
25 bse_screening_rpa,&
26 bse_screening_tdhf,&
27 bse_screening_w0
28 USE kinds, ONLY: dp
29 USE mp2_types, ONLY: mp2_type
30 USE particle_types, ONLY: particle_type
31 USE physcon, ONLY: angstrom,&
32 evolt
33 USE qs_environment_types, ONLY: get_qs_env,&
34 qs_environment_type
35#include "./base/base_uses.f90"
36
37 IMPLICIT NONE
38
39 PRIVATE
40
41 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'bse_print'
42
43 PUBLIC :: print_bse_start_flag, fm_write_thresh, print_excitation_energies, &
44 print_output_header, print_transition_amplitudes, print_optical_properties, &
45 print_exciton_descriptors
46
47CONTAINS
48
49! **************************************************************************************************
50!> \brief ...
51!> \param bse_tda ...
52!> \param bse_abba ...
53!> \param unit_nr ...
54! **************************************************************************************************
55 SUBROUTINE print_bse_start_flag(bse_tda, bse_abba, unit_nr)
56
57 LOGICAL, INTENT(IN) :: bse_tda, bse_abba
58 INTEGER, INTENT(IN) :: unit_nr
59
60 CHARACTER(LEN=*), PARAMETER :: routinen = 'print_BSE_start_flag'
61
62 INTEGER :: handle
63
64 CALL timeset(routinen, handle)
65 CALL cite_reference(graml2026)
66
67 IF (unit_nr > 0) THEN
68 WRITE (unit_nr, *) ' '
69 WRITE (unit_nr, '(T2,A79)') '*******************************************************************************'
70 WRITE (unit_nr, '(T2,A79)') '** **'
71 WRITE (unit_nr, '(T2,A79)') '** Bethe Salpeter equation (BSE) for excitation energies **'
72 IF (bse_tda .AND. bse_abba) THEN
73 WRITE (unit_nr, '(T2,A79)') '** solved with and without Tamm-Dancoff approximation (TDA) **'
74 ELSE IF (bse_tda) THEN
75 WRITE (unit_nr, '(T2,A79)') '** solved with Tamm-Dancoff approximation (TDA) **'
76 ELSE
77 WRITE (unit_nr, '(T2,A79)') '** solved without Tamm-Dancoff approximation (TDA) **'
78 END IF
79
80 WRITE (unit_nr, '(T2,A79)') '** **'
81 WRITE (unit_nr, '(T2,A79)') '*******************************************************************************'
82 WRITE (unit_nr, *) ' '
83 END IF
84
85 CALL timestop(handle)
86
87 END SUBROUTINE print_bse_start_flag
88
89! **************************************************************************************************
90!> \brief ...
91!> \param homo ...
92!> \param virtual ...
93!> \param homo_irred ...
94!> \param flag_TDA ...
95!> \param multiplet ...
96!> \param alpha ...
97!> \param mp2_env ...
98!> \param unit_nr ...
99! **************************************************************************************************
100 SUBROUTINE print_output_header(homo, virtual, homo_irred, flag_TDA, &
101 multiplet, alpha, mp2_env, unit_nr)
102
103 INTEGER, INTENT(IN) :: homo, virtual, homo_irred
104 LOGICAL, INTENT(IN) :: flag_tda
105 CHARACTER(LEN=10), INTENT(IN) :: multiplet
106 REAL(kind=dp), INTENT(IN) :: alpha
107 TYPE(mp2_type), INTENT(INOUT) :: mp2_env
108 INTEGER, INTENT(IN) :: unit_nr
109
110 CHARACTER(LEN=*), PARAMETER :: routinen = 'print_output_header'
111
112 INTEGER :: handle
113
114 CALL timeset(routinen, handle)
115
116 IF (unit_nr > 0) THEN
117 WRITE (unit_nr, '(T2,A4)') 'BSE|'
118 WRITE (unit_nr, '(T2,A4)') 'BSE|'
119 IF (flag_tda) THEN
120 WRITE (unit_nr, '(T2,A4,T7,A74)') 'BSE|', '**************************************************************************'
121 WRITE (unit_nr, '(T2,A4,T7,A74)') 'BSE|', '* Bethe Salpeter equation (BSE) with Tamm Dancoff approximation (TDA) *'
122 WRITE (unit_nr, '(T2,A4,T7,A74)') 'BSE|', '**************************************************************************'
123 WRITE (unit_nr, '(T2,A4)') 'BSE|'
124 WRITE (unit_nr, '(T2,A4,T7,A48,A23)') 'BSE|', 'The excitations are calculated by diagonalizing ', &
125 'the BSE within the TDA:'
126 WRITE (unit_nr, '(T2,A4)') 'BSE|'
127 WRITE (unit_nr, '(T2,A4,T29,A16)') 'BSE|', 'A X^n = Ω^n X^n'
128 WRITE (unit_nr, '(T2,A4)') 'BSE|'
129 WRITE (unit_nr, '(T2,A4,T7,A23)') 'BSE|', 'i.e. in index notation:'
130 WRITE (unit_nr, '(T2,A4)') 'BSE|'
131 WRITE (unit_nr, '(T2,A4,T7,A41)') 'BSE|', 'sum_jb ( A_ia,jb X_jb^n ) = Ω^n X_ia^n'
132 WRITE (unit_nr, '(T2,A4)') 'BSE|'
133 WRITE (unit_nr, '(T2,A4,T7,A14)') 'BSE|', 'Ref.: Eq. (18)'
134 WRITE (unit_nr, '(T2,A4,T7,A61)') 'BSE|', 'in PRB 113,205152 (2026); https://doi.org/10.1103/38k2-d55h .'
135 ELSE
136 WRITE (unit_nr, '(T2,A4,T7,A74)') 'BSE|', '**************************************************************************'
137 WRITE (unit_nr, '(T2,A4,T7,A74)') 'BSE|', '* Full ("ABBA") Bethe Salpeter equation (BSE) (i.e. without TDA) *'
138 WRITE (unit_nr, '(T2,A4,T7,A74)') 'BSE|', '**************************************************************************'
139 WRITE (unit_nr, '(T2,A4)') 'BSE|'
140 WRITE (unit_nr, '(T2,A4,T7,A48,A24)') 'BSE|', 'The excitations are calculated by diagonalizing ', &
141 'the BSE without the TDA:'
142 WRITE (unit_nr, '(T2,A4)') 'BSE|'
143 WRITE (unit_nr, '(T2,A4,T22,A30)') 'BSE|', '|A B| |X^n| |1 0| |X^n|'
144 WRITE (unit_nr, '(T2,A4,T22,A31)') 'BSE|', '|B A| |Y^n| = Ω^n |0 -1| |Y^n|'
145 WRITE (unit_nr, '(T2,A4)') 'BSE|'
146 WRITE (unit_nr, '(T2,A4,T7,A23)') 'BSE|', 'i.e. in index notation:'
147 WRITE (unit_nr, '(T2,A4)') 'BSE|'
148 WRITE (unit_nr, '(T2,A4,T7,A62)') 'BSE|', ' sum_jb ( A_ia,jb X_jb^n + B_ia,jb Y_jb^n ) = Ω^n X_ia^n'
149 WRITE (unit_nr, '(T2,A4,T7,A62)') 'BSE|', '- sum_jb ( B_ia,jb X_jb^n + A_ia,jb Y_jb^n ) = Ω^n Y_ia^n'
150 WRITE (unit_nr, '(T2,A4)') 'BSE|'
151 WRITE (unit_nr, '(T2,A4,T7,A13)') 'BSE|', 'Ref.: Eq. (8)'
152 WRITE (unit_nr, '(T2,A4,T7,A61)') 'BSE|', 'in PRB 113,205152 (2026); https://doi.org/10.1103/38k2-d55h .'
153 END IF
154 WRITE (unit_nr, '(T2,A4)') 'BSE|'
155 WRITE (unit_nr, '(T2,A4)') 'BSE|'
156 WRITE (unit_nr, '(T2,A4,T7,A4,T18,A42,T70,A1,I4,A1,I4,A1)') 'BSE|', 'i,j:', &
157 'occupied molecular orbitals, i.e. state in', '[', homo_irred - homo + 1, ',', homo_irred, ']'
158 WRITE (unit_nr, '(T2,A4,T7,A4,T18,A44,T70,A1,I4,A1,I4,A1)') 'BSE|', 'a,b:', &
159 'unoccupied molecular orbitals, i.e. state in', '[', homo_irred + 1, ',', homo_irred + virtual, ']'
160 WRITE (unit_nr, '(T2,A4,T7,A2,T18,A16)') 'BSE|', 'n:', 'Excitation index'
161 WRITE (unit_nr, '(T2,A4)') 'BSE|'
162 IF (mp2_env%bse%screening_method == bse_screening_w0) THEN
163 WRITE (unit_nr, '(T2,A4,T7,A)') 'BSE|', 'A_ia,jb = (ε_a-ε_i) δ_ij δ_ab + α * v_ia,jb - W_ij,ab'
164 ELSE IF (mp2_env%bse%screening_method == bse_screening_rpa) THEN
165 WRITE (unit_nr, '(T2,A4,T7,A)') 'BSE|', 'A_ia,jb = (ε_a-ε_i) δ_ij δ_ab + α * v_ia,jb'
166 END IF
167 IF (.NOT. flag_tda) THEN
168 IF (mp2_env%bse%screening_method == bse_screening_w0) THEN
169 WRITE (unit_nr, '(T2,A4,T7,A)') 'BSE|', 'B_ia,jb = α * v_ia,jb - W_ib,aj'
170 ELSE IF (mp2_env%bse%screening_method == bse_screening_rpa) THEN
171 WRITE (unit_nr, '(T2,A4,T7,A)') 'BSE|', 'B_ia,jb = α * v_ia,jb'
172 END IF
173 WRITE (unit_nr, '(T2,A4)') 'BSE|'
174 WRITE (unit_nr, '(T2,A4,T7,A14)') 'BSE|', 'Ref.: Eq. (10)'
175 WRITE (unit_nr, '(T2,A4,T7,A61)') 'BSE|', 'in PRB 113,205152 (2026); https://doi.org/10.1103/38k2-d55h .'
176 END IF
177 IF (.NOT. flag_tda) THEN
178 WRITE (unit_nr, '(T2,A4)') 'BSE|'
179 WRITE (unit_nr, '(T2,A4,T7,A74)') 'BSE|', 'The BSE is solved for Ω^n and X_ia^n as a hermitian problem, e.g. Eq.(42)'
180 WRITE (unit_nr, '(T2,A4,T7,A71)') 'BSE|', 'in PRB 92,045209 (2015); http://dx.doi.org/10.1103/PhysRevB.92.045209 .'
181 END IF
182 WRITE (unit_nr, '(T2,A4)') 'BSE|'
183 WRITE (unit_nr, '(T2,A4,T7,A7,T31,A23)') 'BSE|', 'ε_...:', 'GW quasiparticle energy'
184 WRITE (unit_nr, '(T2,A4,T7,A7,T31,A15)') 'BSE|', 'δ_...:', 'Kronecker delta'
185 WRITE (unit_nr, '(T2,A4,T7,A3,T31,A21)') 'BSE|', 'α:', 'spin-dependent factor (Singlet/Triplet)'
186 WRITE (unit_nr, '(T2,A4,T7,A6,T30,A34)') 'BSE|', 'v_...:', 'Electron-hole exchange interaction'
187 IF (mp2_env%bse%screening_method == bse_screening_w0) THEN
188 WRITE (unit_nr, '(T2,A4,T7,A,T31,A)') 'BSE|', 'W_... = 1/ϵ v_...:', &
189 'Direct interaction screened by '
190 WRITE (unit_nr, '(T2,A4,T30,A)') 'BSE|', &
191 'dielectric function ϵ(ω=0)'
192 ELSE IF (mp2_env%bse%screening_method == bse_screening_tdhf) THEN
193 WRITE (unit_nr, '(T2,A4,T7,A,T30,A)') 'BSE|', 'W_... = v_...:', 'Direct interaction without screening'
194 ELSE IF (mp2_env%bse%screening_method == bse_screening_alpha) THEN
195 WRITE (unit_nr, '(T2,A4,T7,A,T31,A,F5.2)') 'BSE|', 'W_... = γ v_...:', &
196 'Direct interaction with artificial screening γ=', mp2_env%bse%screening_factor
197 END IF
198 WRITE (unit_nr, '(T2,A4)') 'BSE|'
199 WRITE (unit_nr, '(T2,A4)') 'BSE|'
200 WRITE (unit_nr, '(T2,A4,T7,A47,A7,A9,F3.1)') 'BSE|', &
201 'The spin-dependent factor is for the requested ', multiplet, " is α = ", alpha
202 WRITE (unit_nr, '(T2,A4)') 'BSE|'
203 END IF
204
205 CALL timestop(handle)
206
207 END SUBROUTINE print_output_header
208
209! **************************************************************************************************
210!> \brief ...
211!> \param Exc_ens ...
212!> \param homo ...
213!> \param virtual ...
214!> \param flag_TDA ...
215!> \param multiplet ...
216!> \param info_approximation ...
217!> \param mp2_env ...
218!> \param unit_nr ...
219! **************************************************************************************************
220 SUBROUTINE print_excitation_energies(Exc_ens, homo, virtual, flag_TDA, multiplet, &
221 info_approximation, mp2_env, unit_nr)
222
223 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: exc_ens
224 INTEGER, INTENT(IN) :: homo, virtual
225 LOGICAL, INTENT(IN) :: flag_tda
226 CHARACTER(LEN=10), INTENT(IN) :: multiplet, info_approximation
227 TYPE(mp2_type), INTENT(INOUT) :: mp2_env
228 INTEGER, INTENT(IN) :: unit_nr
229
230 CHARACTER(LEN=*), PARAMETER :: routinen = 'print_excitation_energies'
231
232 INTEGER :: handle, i_exc
233
234 CALL timeset(routinen, handle)
235
236 IF (unit_nr > 0) THEN
237 IF (flag_tda) THEN
238 WRITE (unit_nr, '(T2,A4,T7,A56)') 'BSE|', 'Excitation energies from solving the BSE within the TDA:'
239 ELSE
240 WRITE (unit_nr, '(T2,A4,T7,A57)') 'BSE|', 'Excitation energies from solving the BSE without the TDA:'
241 END IF
242 WRITE (unit_nr, '(T2,A4)') 'BSE|'
243 WRITE (unit_nr, '(T2,A4,T11,A12,T30,A7,T44,A8,T55,A27)') 'BSE|', &
244 'Excitation n', multiplet, 'TDA/ABBA', 'Excitation energy Ω^n (eV)'
245 END IF
246 !prints actual energies values
247 IF (unit_nr > 0) THEN
248 DO i_exc = 1, min(homo*virtual, mp2_env%bse%num_print_exc)
249 WRITE (unit_nr, '(T2,A4,T7,I16,T30,A7,T46,A6,T59,F22.4)') &
250 'BSE|', i_exc, multiplet, info_approximation, exc_ens(i_exc)*evolt
251 END DO
252 END IF
253
254 CALL timestop(handle)
255
256 END SUBROUTINE print_excitation_energies
257
258! **************************************************************************************************
259!> \brief ...
260!> \param fm_eigvec_X ...
261!> \param homo ...
262!> \param virtual ...
263!> \param homo_irred ...
264!> \param info_approximation ...
265!> \param mp2_env ...
266!> \param unit_nr ...
267!> \param fm_eigvec_Y ...
268! **************************************************************************************************
269 SUBROUTINE print_transition_amplitudes(fm_eigvec_X, homo, virtual, homo_irred, &
270 info_approximation, mp2_env, unit_nr, fm_eigvec_Y)
271
272 TYPE(cp_fm_type), INTENT(IN) :: fm_eigvec_x
273 INTEGER, DIMENSION(:), INTENT(IN) :: homo, virtual, homo_irred
274 CHARACTER(LEN=10), INTENT(IN) :: info_approximation
275 TYPE(mp2_type), INTENT(INOUT) :: mp2_env
276 INTEGER, INTENT(IN) :: unit_nr
277 TYPE(cp_fm_type), INTENT(IN), OPTIONAL :: fm_eigvec_y
278
279 CHARACTER(LEN=*), PARAMETER :: routinen = 'print_transition_amplitudes'
280
281 INTEGER :: handle, i_exc, isp, n_ov_joint, nspins
282 INTEGER, ALLOCATABLE, DIMENSION(:) :: n_ov, offsets
283
284 CALL timeset(routinen, handle)
285
286 nspins = SIZE(homo)
287 ALLOCATE (n_ov(nspins), offsets(nspins))
288 CALL get_bse_spin_block_layout(homo, virtual, n_ov, offsets, n_ov_joint)
289
290 IF (unit_nr > 0) THEN
291 WRITE (unit_nr, '(T2,A4)') 'BSE|'
292 WRITE (unit_nr, '(T2,A4,T7,A61)') &
293 'BSE|', "Single-particle transitions are built up by (de-)excitations,"
294 WRITE (unit_nr, '(T2,A4,T7,A18)') &
295 'BSE|', "which we denote by"
296 WRITE (unit_nr, '(T2,A4)') 'BSE|'
297 WRITE (unit_nr, '(T2,A4,T20,A2,T30,A40)') &
298 'BSE|', "=>", "for excitations, i.e. entries of X_ia^n,"
299 WRITE (unit_nr, '(T2,A4,T20,A2,T30,A42)') &
300 'BSE|', "<=", "for deexcitations, i.e. entries of Y_ia^n."
301 WRITE (unit_nr, '(T2,A4)') &
302 'BSE|'
303 WRITE (unit_nr, '(T2,A4,T7,A73)') &
304 'BSE|', "The following single-particle transitions have significant contributions,"
305 WRITE (unit_nr, '(T2,A4,T7,A16,F5.3,A15,F5.3,A16)') &
306 'BSE|', "i.e. |X_ia^n| > ", mp2_env%bse%eps_x, " or |Y_ia^n| > ", &
307 mp2_env%bse%eps_x, ", respectively :"
308
309 IF (nspins == 1) THEN
310 WRITE (unit_nr, '(T2,A4,T15,A27,I5,A13,I5,A3)') 'BSE|', '-- Quick reminder: HOMO i =', &
311 homo_irred(1), ' and LUMO a =', homo_irred(1) + 1, " --"
312 WRITE (unit_nr, '(T2,A4)') 'BSE|'
313 WRITE (unit_nr, '(T2,A4,T7,A12,T30,A1,T32,A5,T42,A1,T49,A8,T64,A17)') &
314 "BSE|", "Excitation n", "i", "=>/<=", "a", 'TDA/ABBA', "|X_ia^n|/|Y_ia^n|"
315 ELSE
316 ! bare A (no width) for sigma-bearing literals: explicit widths count bytes, and the
317 ! 2-byte UTF-8 sigma would otherwise truncate.
318 DO isp = 1, nspins
319 WRITE (unit_nr, '(T2,A4,T15,A,I2,A,I5,A,I5,A)') 'BSE|', &
320 '-- Quick reminder: σ =', isp, ', HOMO i =', homo_irred(isp), &
321 ' and LUMO a =', homo_irred(isp) + 1, " --"
322 END DO
323 WRITE (unit_nr, '(T2,A4)') 'BSE|'
324 WRITE (unit_nr, '(T2,A4,T7,A12,T22,A,T30,A1,T32,A5,T42,A1,T49,A8,T64,A)') &
325 "BSE|", "Excitation n", "σ", "i", "=>/<=", "a", 'TDA/ABBA', "|X_iaσ^n|/|Y_iaσ^n|"
326 END IF
327 END IF
328 DO i_exc = 1, min(n_ov_joint, mp2_env%bse%num_print_exc)
329 IF (unit_nr > 0) THEN
330 WRITE (unit_nr, '(T2,A4)') 'BSE|'
331 END IF
332 !Iterate through eigenvector and print values above threshold
333 CALL print_transition_amplitudes_core(fm_eigvec_x, "=>", info_approximation, &
334 i_exc, virtual, homo, homo_irred, &
335 unit_nr, mp2_env, offsets)
336 IF (PRESENT(fm_eigvec_y)) THEN
337 CALL print_transition_amplitudes_core(fm_eigvec_y, "<=", info_approximation, &
338 i_exc, virtual, homo, homo_irred, &
339 unit_nr, mp2_env, offsets)
340 END IF
341 END DO
342
343 DEALLOCATE (n_ov, offsets)
344 CALL timestop(handle)
345
346 END SUBROUTINE print_transition_amplitudes
347
348! **************************************************************************************************
349!> \brief ...
350!> \param Exc_ens ...
351!> \param oscill_str ...
352!> \param trans_mom_bse ...
353!> \param polarizability_residues ...
354!> \param homo ...
355!> \param virtual ...
356!> \param homo_irred ...
357!> \param flag_TDA ...
358!> \param info_approximation ...
359!> \param mp2_env ...
360!> \param unit_nr ...
361!> \param open_shell if .TRUE., print spin-summed (UKS) dipole formula instead of the sqrt(2) one
362! **************************************************************************************************
363 SUBROUTINE print_optical_properties(Exc_ens, oscill_str, trans_mom_bse, polarizability_residues, &
364 homo, virtual, homo_irred, flag_TDA, &
365 info_approximation, mp2_env, unit_nr, open_shell)
366
367 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: exc_ens, oscill_str
368 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :) :: trans_mom_bse, polarizability_residues
369 INTEGER, INTENT(IN) :: homo, virtual, homo_irred
370 LOGICAL, INTENT(IN) :: flag_tda
371 CHARACTER(LEN=10), INTENT(IN) :: info_approximation
372 TYPE(mp2_type), INTENT(INOUT) :: mp2_env
373 INTEGER, INTENT(IN) :: unit_nr
374 LOGICAL, INTENT(IN), OPTIONAL :: open_shell
375
376 CHARACTER(LEN=*), PARAMETER :: routinen = 'print_optical_properties'
377
378 INTEGER :: handle, i_exc
379 LOGICAL :: my_open_shell
380
381 CALL timeset(routinen, handle)
382
383 my_open_shell = .false.
384 IF (PRESENT(open_shell)) my_open_shell = open_shell
385
386 ! Discriminate between singlet and triplet, since triplet state can't couple to light
387 ! and therefore calculations of dipoles etc are not necessary
388 IF (mp2_env%bse%bse_spin_config == 0) THEN
389 IF (unit_nr > 0) THEN
390 WRITE (unit_nr, '(T2,A4)') 'BSE|'
391 WRITE (unit_nr, '(T2,A4,T7,A60)') &
392 'BSE|', "Transition moments d_r^n (with r∈(x,y,z), in atomic units)"
393 WRITE (unit_nr, '(T2,A4,T7,A67)') &
394 'BSE|', "and oscillator strength f^n of excitation level n are obtained from"
395 WRITE (unit_nr, '(T2,A4)') 'BSE|'
396 IF (my_open_shell) THEN
397 IF (flag_tda) THEN
398 WRITE (unit_nr, '(T2,A4,T10,A)') &
399 'BSE|', "d_r^n = sum_σ sum_ia < ψ_iσ | r | ψ_aσ > X_iaσ^n"
400 ELSE
401 WRITE (unit_nr, '(T2,A4,T10,A)') &
402 'BSE|', "d_r^n = sum_σ sum_ia < ψ_iσ | r | ψ_aσ > ( X_iaσ^n + Y_iaσ^n )"
403 END IF
404 ELSE
405 IF (flag_tda) THEN
406 WRITE (unit_nr, '(T2,A4,T10,A)') &
407 'BSE|', "d_r^n = sqrt(2) sum_ia < ψ_i | r | ψ_a > X_ia^n"
408 ELSE
409 WRITE (unit_nr, '(T2,A4,T10,A)') &
410 'BSE|', "d_r^n = sum_ia sqrt(2) < ψ_i | r | ψ_a > ( X_ia^n + Y_ia^n )"
411 END IF
412 END IF
413 WRITE (unit_nr, '(T2,A4)') 'BSE|'
414 WRITE (unit_nr, '(T2,A4,T14,A)') &
415 'BSE|', "f^n = 2/3 * Ω^n sum_r∈(x,y,z) ( d_r^n )^2"
416 WRITE (unit_nr, '(T2,A4)') 'BSE|'
417 WRITE (unit_nr, '(T2,A4,T7,A19)') &
418 'BSE|', "where we introduced"
419 WRITE (unit_nr, '(T2,A4)') 'BSE|'
420 WRITE (unit_nr, '(T2,A4,T7,A5,T15,A28)') &
421 'BSE|', "ψ_i:", "occupied molecular orbitals,"
422 WRITE (unit_nr, '(T2,A4,T7,A5,T15,A28)') &
423 'BSE|', "ψ_a:", "empty molecular orbitals and"
424 WRITE (unit_nr, '(T2,A4,T9,A2,T14,A18)') &
425 'BSE|', "r:", "position operator."
426 WRITE (unit_nr, '(T2,A4)') 'BSE|'
427 WRITE (unit_nr, '(T2,A4,T7,A21)') &
428 'BSE|', "Ref.: Eqs. (21), (D3)"
429 WRITE (unit_nr, '(T2,A4,T7,A69)') &
430 'BSE|', "in Phys. Rev. B 113, 205152 (2026); https://doi.org/10.1103/38k2-d55h"
431 WRITE (unit_nr, '(T2,A4)') 'BSE|'
432 IF (flag_tda) THEN
433 WRITE (unit_nr, '(T2,A4,T7,A55)') 'BSE|', &
434 'Optical properties from solving the BSE within the TDA:'
435 ELSE
436 WRITE (unit_nr, '(T2,A4,T7,A56)') 'BSE|', &
437 'Optical properties from solving the BSE without the TDA:'
438 END IF
439 WRITE (unit_nr, '(T2,A4)') 'BSE|'
440 WRITE (unit_nr, '(T2,A4,T8,A12,T22,A8,T38,A5,T48,A5,T58,A5,T64,A17)') 'BSE|', &
441 'Excitation n', "TDA/ABBA", "d_x^n", "d_y^n", "d_z^n", 'Osc. strength f^n'
442 DO i_exc = 1, min(homo*virtual, mp2_env%bse%num_print_exc)
443 WRITE (unit_nr, '(T2,A4,T8,I12,T24,A6,T35,F8.3,T45,F8.3,T55,F8.3,T65,F16.3)') &
444 'BSE|', i_exc, info_approximation, trans_mom_bse(1, 1, i_exc), trans_mom_bse(2, 1, i_exc), &
445 trans_mom_bse(3, 1, i_exc), oscill_str(i_exc)
446 END DO
447 WRITE (unit_nr, '(T2,A4)') 'BSE|'
448 WRITE (unit_nr, '(T2,A4,T7,A)') 'BSE|', &
449 'Check for Thomas-Reiche-Kuhn sum rule'
450 WRITE (unit_nr, '(T2,A4)') 'BSE|'
451 WRITE (unit_nr, '(T2,A4,T35,A15)') 'BSE|', &
452 'N_e = Σ_n f^n'
453 WRITE (unit_nr, '(T2,A4)') 'BSE|'
454 ! Open shell: caller passes homo_irred = n_alpha + n_beta (total electrons).
455 ! Closed shell: homo_irred = n_occ, i.e. 2 electrons per occupied orbital.
456 WRITE (unit_nr, '(T2,A4,T7,A24,T65,I16)') 'BSE|', &
457 'Number of electrons N_e:', merge(homo_irred, homo_irred*2, my_open_shell)
458 WRITE (unit_nr, '(T2,A4,T7,A40,T66,F16.3)') 'BSE|', &
459 'Sum over oscillator strengths Σ_n f^n :', sum(oscill_str)
460 WRITE (unit_nr, '(T2,A4)') 'BSE|'
461 IF (mp2_env%bse%bse_cutoff_occ > 0 .OR. mp2_env%bse%bse_cutoff_empty > 0) THEN
462 CALL cp_warn(__location__, &
463 "Accuracy of TRK sum rule might suffer from cutoffs.")
464 END IF
465 END IF
466
467 ! Compute and print the absorption spectrum to external file
468 IF (mp2_env%bse%bse_print_spectrum) THEN
469 CALL compute_and_print_absorption_spectrum(oscill_str, polarizability_residues, exc_ens, &
470 info_approximation, unit_nr, mp2_env)
471 END IF
472
473 ELSE
474 IF (unit_nr > 0) THEN
475 WRITE (unit_nr, '(T2,A4)') 'BSE|'
476 WRITE (unit_nr, '(T2,A4)') 'BSE|'
477 CALL cp_warn(__location__, &
478 "Requested triplet excitation cannot couple to light. "// &
479 "Skipping calculation of transition moments, "// &
480 "oscillator strengths, and spectrum.")
481 END IF
482 END IF
483
484 CALL timestop(handle)
485
486 END SUBROUTINE print_optical_properties
487
488! **************************************************************************************************
489!> \brief ...
490!> \param fm_eigvec ...
491!> \param direction_excitation ...
492!> \param info_approximation ...
493!> \param i_exc ...
494!> \param virtual ...
495!> \param homo ...
496!> \param homo_irred ...
497!> \param unit_nr ...
498!> \param mp2_env ...
499!> \param offsets ...
500! **************************************************************************************************
501 SUBROUTINE print_transition_amplitudes_core(fm_eigvec, direction_excitation, info_approximation, &
502 i_exc, virtual, homo, homo_irred, &
503 unit_nr, mp2_env, offsets)
504
505 TYPE(cp_fm_type), INTENT(IN) :: fm_eigvec
506 CHARACTER(LEN=2), INTENT(IN) :: direction_excitation
507 CHARACTER(LEN=10), INTENT(IN) :: info_approximation
508 INTEGER, INTENT(IN) :: i_exc
509 INTEGER, DIMENSION(:), INTENT(IN) :: virtual, homo, homo_irred
510 INTEGER, INTENT(IN) :: unit_nr
511 TYPE(mp2_type), INTENT(INOUT) :: mp2_env
512 INTEGER, DIMENSION(:), INTENT(IN) :: offsets
513
514 CHARACTER(LEN=*), PARAMETER :: routinen = 'print_transition_amplitudes_core'
515 CHARACTER(LEN=2), DIMENSION(2), PARAMETER :: spin_label = ["α", "β"]
516
517 INTEGER :: handle, isp, k, num_entries
518 INTEGER, ALLOCATABLE, DIMENSION(:) :: idx_homo, idx_spin, idx_virt
519 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: eigvec_entries
520
521! 2-byte UTF-8 glyphs (LEN=1 would truncate both alpha/beta to the shared 0xCE byte)
522
523 CALL timeset(routinen, handle)
524
525 ! direction_excitation can be either => (means excitation; from fm_eigvec_X)
526 ! or <= (means deexcitation; from fm_eigvec_Y)
527 IF (SIZE(homo) == 1) THEN
528 CALL filter_eigvec_contrib(fm_eigvec, idx_homo, idx_virt, eigvec_entries, &
529 i_exc, virtual(1), num_entries, mp2_env)
530 IF (unit_nr > 0) THEN
531 DO k = 1, num_entries
532 WRITE (unit_nr, '(T2,A4,T14,I5,T26,I5,T35,A2,T38,I5,T51,A6,T65,F16.4)') &
533 "BSE|", i_exc, homo_irred(1) - homo(1) + idx_homo(k), direction_excitation, &
534 homo_irred(1) + idx_virt(k), info_approximation, abs(eigvec_entries(k))
535 END DO
536 END IF
537 ELSE
538 CALL filter_eigvec_contrib(fm_eigvec, idx_homo, idx_virt, eigvec_entries, &
539 i_exc, virtual(1), num_entries, mp2_env, &
540 offsets=offsets, virtual_per_spin=virtual, idx_spin=idx_spin)
541 IF (unit_nr > 0) THEN
542 DO k = 1, num_entries
543 isp = idx_spin(k)
544 WRITE (unit_nr, '(T2,A4,T14,I5,T22,A2,T26,I5,T35,A2,T38,I5,T51,A6,T65,F16.4)') &
545 "BSE|", i_exc, spin_label(isp), &
546 homo_irred(isp) - homo(isp) + idx_homo(k), direction_excitation, &
547 homo_irred(isp) + idx_virt(k), info_approximation, abs(eigvec_entries(k))
548 END DO
549 END IF
550 DEALLOCATE (idx_spin)
551 END IF
552 DEALLOCATE (idx_homo)
553 DEALLOCATE (idx_virt)
554 DEALLOCATE (eigvec_entries)
555 CALL timestop(handle)
556
557 END SUBROUTINE print_transition_amplitudes_core
558
559! **************************************************************************************************
560!> \brief Prints exciton descriptors, cf. Mewes et al., JCTC 14, 710-725 (2018)
561!> \param exc_descr Exciton descriptors with size of num_print_exc_descr
562!> \param ref_point_multipole Reference point for computation of multipole moments, e.g. center of mass
563!> \param unit_nr ...
564!> \param num_print_exc_descr Number of excitation levels for which descriptors are printed
565!> \param print_checkvalue Flag, which determines if values for regtests should be printed
566!> \param print_directional_exc_descr Flag, which activates printing of directional descriptors
567!> \param print_directional_crosscorrelation Flag, which adds the crosscorrelation matrix to the directional descriptors
568!> \param prefix_output String, which is put in front of prints, i.e. "BSE|" or "" for TDDFPT
569!> \param qs_env ...
570! **************************************************************************************************
571 SUBROUTINE print_exciton_descriptors(exc_descr, ref_point_multipole, unit_nr, &
572 num_print_exc_descr, print_checkvalue, print_directional_exc_descr, &
573 print_directional_crosscorrelation, prefix_output, qs_env)
574
575 TYPE(exciton_descr_type), ALLOCATABLE, &
576 DIMENSION(:) :: exc_descr
577 REAL(kind=dp), ALLOCATABLE, DIMENSION(:), &
578 INTENT(IN) :: ref_point_multipole
579 INTEGER, INTENT(IN) :: unit_nr, num_print_exc_descr
580 LOGICAL, INTENT(IN) :: print_checkvalue, print_directional_exc_descr, &
581 print_directional_crosscorrelation
582 CHARACTER(LEN=4), INTENT(IN) :: prefix_output
583 TYPE(qs_environment_type), POINTER :: qs_env
584
585 CHARACTER(LEN=*), PARAMETER :: routinen = 'print_exciton_descriptors'
586
587 CHARACTER(LEN=1), DIMENSION(3) :: array_direction_str
588 CHARACTER(LEN=5) :: method_name
589 INTEGER :: handle, i_dir, i_exc
590 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
591
592 IF (prefix_output == 'BSE|') THEN
593 method_name = 'BSE'
594 ELSE
595 method_name = 'TDDFT'
596 END IF
597
598 CALL timeset(routinen, handle)
599 CALL get_qs_env(qs_env, particle_set=particle_set)
600 IF (unit_nr > 0) THEN
601 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
602 'Exciton descriptors for excitation level n are given by'
603 WRITE (unit_nr, '(T2,A4)') prefix_output
604 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
605 'd_eh = | <r_h - r_e>_exc |'
606 WRITE (unit_nr, '(T2,A4)') prefix_output
607 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
608 'σ_e = sqrt( <r_e^2>_exc - <r_e>_exc^2 )'
609 WRITE (unit_nr, '(T2,A4)') prefix_output
610 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
611 'σ_h = sqrt( <r_h^2>_exc - <r_h>_exc^2 )'
612 WRITE (unit_nr, '(T2,A4)') prefix_output
613 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
614 'COV_eh = <r_e r_h>_exc - <r_e>_exc <r_h>_exc'
615 WRITE (unit_nr, '(T2,A4)') prefix_output
616 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
617 'd_exc = sqrt( | < |r_h - r_e|^2 >_exc )'
618 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
619 ' = sqrt( d_eh^2 + σ_e^2 + σ_h^2 - 2 * COV_eh )'
620 WRITE (unit_nr, '(T2,A4)') prefix_output
621 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
622 'R_eh = COV_eh / (σ_e * σ_h)'
623 WRITE (unit_nr, '(T2,A4)') prefix_output
624 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
625 'where the expectation values <.>_exc are taken with respect to the '
626 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
627 'exciton wavefunction of excitation n:'
628 WRITE (unit_nr, '(T2,A4)') prefix_output
629
630 IF (exc_descr(1)%flag_TDA) THEN
631 WRITE (unit_nr, '(T2,A4,T20,A)') prefix_output, &
632 '𝚿_n(r_e,r_h) = Σ_{i,a} X_ia^n ψ_i(r_h) ψ_a(r_e) ,'
633 ELSE
634 WRITE (unit_nr, '(T2,A4,T20,A)') prefix_output, &
635 '𝚿_n(r_e,r_h) = Σ_{i,a} X_ia^n ψ_i(r_h) ψ_a(r_e)'
636 WRITE (unit_nr, '(T2,A4,T40,A)') prefix_output, &
637 '+ Y_ia^n ψ_a(r_h) ψ_i(r_e) ,'
638 END IF
639 WRITE (unit_nr, '(T2,A4)') prefix_output
640 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
641 'i.e.'
642 WRITE (unit_nr, '(T2,A4)') prefix_output
643 WRITE (unit_nr, '(T2,A4,T20,A)') prefix_output, &
644 '< O >_exc = < 𝚿_n | O | 𝚿_n > / < 𝚿_n | 𝚿_n > ,'
645 WRITE (unit_nr, '(T2,A4)') prefix_output
646 IF (exc_descr(1)%flag_TDA) THEN
647 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
648 'where c_n = < 𝚿_n | 𝚿_n > = 1 within TDA.'
649 ELSE
650 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
651 'where c_n = < 𝚿_n | 𝚿_n > ≥ 1 without TDA.'
652 END IF
653 WRITE (unit_nr, '(T2,A4)') prefix_output
654 WRITE (unit_nr, '(T2,A4)') prefix_output
655 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
656 'Here, we introduced'
657 WRITE (unit_nr, '(T2,A4)') prefix_output
658 WRITE (unit_nr, '(T2,A4,T7,A5,T15,A)') &
659 prefix_output, "ψ_i:", "occupied molecular orbitals,"
660 WRITE (unit_nr, '(T2,A4,T7,A5,T15,A)') &
661 prefix_output, "ψ_a:", "empty molecular orbitals and"
662 WRITE (unit_nr, '(T2,A4,T9,A2,T14,A)') &
663 prefix_output, "r:", "position operator."
664 WRITE (unit_nr, '(T2,A4)') prefix_output
665 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
666 'Ref.: Eqs. (28)-(32)'
667 WRITE (unit_nr, '(T2,A4,T7,A,A)') prefix_output, &
668 'Phys. Rev. B 113, 205152 (2026); ', &
669 'https://doi.org/10.1103/38k2-d55h'
670 WRITE (unit_nr, '(T2,A4)') prefix_output
671 WRITE (unit_nr, '(T2,A4)') prefix_output
672 IF (exc_descr(1)%flag_TDA) THEN
673 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
674 'Exciton descriptors from solving the ', method_name, ' within the TDA:'
675 ELSE
676 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
677 'Exciton descriptors from solving the ', method_name, ' without the TDA:'
678 END IF
679 WRITE (unit_nr, '(T2,A4)') prefix_output
680 WRITE (unit_nr, '(T2,A4,T10,A1,6X,A3,1X,4X,A10,5X,A10,5X,A10,3X,A11,8X,A4)') prefix_output, &
681 'n', 'c_n', 'd_eh [Å]', 'σ_e [Å]', 'σ_h [Å]', 'd_exc [Å]', 'R_eh'
682 DO i_exc = 1, num_print_exc_descr
683 WRITE (unit_nr, '(T2,A4,T7,I4,4X,F5.3,1X,5(2X,F10.4))') &
684 prefix_output, i_exc, exc_descr(i_exc)%norm_XpY, &
685 exc_descr(i_exc)%diff_r_abs*angstrom, &
686 exc_descr(i_exc)%sigma_e*angstrom, exc_descr(i_exc)%sigma_h*angstrom, &
687 exc_descr(i_exc)%diff_r_sqr*angstrom, exc_descr(i_exc)%corr_e_h
688 END DO
689 WRITE (unit_nr, '(T2,A4)') prefix_output
690 ! For debug runs, print first d_exc separately to allow the regtests to read in
691 IF (print_checkvalue) THEN
692 WRITE (unit_nr, '(T2)')
693 WRITE (unit_nr, '(T2,A28,T65,F16.4)') 'Checksum exciton descriptors', &
694 exc_descr(1)%diff_r_sqr*angstrom
695 WRITE (unit_nr, '(T2)')
696 END IF
697 WRITE (unit_nr, '(T2,A4)') prefix_output
698 ! Print exciton descriptor resolved per direction
699 IF (print_directional_exc_descr) THEN
700 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
701 'We can restrict the exciton descriptors to a specific direction,'
702 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
703 'e.g. the x-components are:'
704 WRITE (unit_nr, '(T2,A4)') prefix_output
705 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
706 'd_eh^x = | <x_h - x_e>_exc |'
707 WRITE (unit_nr, '(T2,A4)') prefix_output
708 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
709 'σ_e^x = sqrt( <x_e^2>_exc - <x_e>_exc^2 )'
710 WRITE (unit_nr, '(T2,A4)') prefix_output
711 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
712 'σ_h^x = sqrt( <x_h^2>_exc - <x_h>_exc^2 )'
713 WRITE (unit_nr, '(T2,A4)') prefix_output
714 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
715 "COV_eh^{μμ'} = <r^μ_e r^μ'_h>_exc - <r^μ_e>_exc <r^μ'_h>_exc"
716 WRITE (unit_nr, '(T2,A4)') prefix_output
717 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
718 'd_exc^x = sqrt( | < |x_h - x_e|^2 >_exc )'
719 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
720 ' = sqrt( (d_eh^x)^2 + (σ_e^x)^2'
721 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
722 " + (σ_h^x)^2 - 2 * (COV_eh^{xx}) )"
723 WRITE (unit_nr, '(T2,A4)') prefix_output
724 IF (print_directional_crosscorrelation) THEN
725 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
726 "Subsequently, the cross-correlation matrix R_eh^{μμ'} is printed"
727 WRITE (unit_nr, '(T2,A4)') prefix_output
728 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
729 "R_eh^{μμ'} = COV_eh^{μμ'}/(σ^μ_e σ^μ'_h)"
730 WRITE (unit_nr, '(T2,A4)') prefix_output
731 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
732 "Rows: electron direction μ, columns: hole direction μ',"
733 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
734 "so in general R_eh^{μμ'} ≠ R_eh^{μ'μ}."
735 WRITE (unit_nr, '(T2,A4)') prefix_output
736 END IF
737 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
738 'Ref.: Eqs. (33)-(36)'
739 WRITE (unit_nr, '(T2,A4,T7,A,A)') prefix_output, &
740 'Phys. Rev. B 113, 205152 (2026); ', &
741 'https://doi.org/10.1103/38k2-d55h'
742 WRITE (unit_nr, '(T2,A4)') prefix_output
743 WRITE (unit_nr, '(T2,A4)') prefix_output
744 IF (exc_descr(1)%flag_TDA) THEN
745 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
746 'Exciton descriptors per direction from solving the ', method_name, ' within the TDA:'
747 ELSE
748 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
749 'Exciton descriptors per direction from solving the ', method_name, ' without the TDA:'
750 END IF
751 WRITE (unit_nr, '(T2,A4)') prefix_output
752 WRITE (unit_nr, '(T2,A4,T12,A1,2X,A9,5X,A12,5X,A12,5X,A12,3X,A13)') prefix_output, &
753 'n', 'r = x/y/z', 'd_eh^r [Å]', 'σ_e^r [Å]', 'σ_h^r [Å]', 'd_exc^r [Å]'
754 DO i_exc = 1, num_print_exc_descr
755 DO i_dir = 1, 3
756 array_direction_str = ["x", "y", "z"]
757 WRITE (unit_nr, '(T2,A4,T9,I4,10X,A1,1X,4(4X,F10.4))') &
758 prefix_output, i_exc, array_direction_str(i_dir), &
759 exc_descr(i_exc)%d_eh_dir(i_dir)*angstrom, &
760 exc_descr(i_exc)%sigma_e_dir(i_dir)*angstrom, &
761 exc_descr(i_exc)%sigma_h_dir(i_dir)*angstrom, &
762 exc_descr(i_exc)%d_exc_dir(i_dir)*angstrom
763 END DO
764 WRITE (unit_nr, '(T2,A4)') prefix_output
765 END DO
766 IF (print_directional_crosscorrelation) THEN
767 WRITE (unit_nr, '(T2,A4)') prefix_output
768 WRITE (unit_nr, '(T2,A4)') prefix_output
769 IF (exc_descr(1)%flag_TDA) THEN
770 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
771 'Crosscorrelation matrix from solving the ', method_name, ' within the TDA:'
772 ELSE
773 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
774 'Crosscorrelation matrix from solving the ', method_name, ' without the TDA:'
775 END IF
776 WRITE (unit_nr, '(T2,A4)') prefix_output
777 WRITE (unit_nr, '(T2,A4,T12,A1,T23,A3,T44,A1,T62,A1,T80,A1)') prefix_output, &
778 'n', 'e\h', 'x', 'y', 'z'
779 ! rows: electron direction, columns: hole direction
780 DO i_exc = 1, num_print_exc_descr
781 WRITE (unit_nr, '(T2,A4,T9,I4,T25,A1,1X,3(11X,F7.4))') prefix_output, i_exc, &
782 array_direction_str(1), exc_descr(i_exc)%corr_e_h_matrix(1, :)
783 DO i_dir = 2, 3
784 WRITE (unit_nr, '(T2,A4,T25,A1,1X,3(11X,F7.4))') prefix_output, &
785 array_direction_str(i_dir), exc_descr(i_exc)%corr_e_h_matrix(i_dir, :)
786 END DO
787 WRITE (unit_nr, '(T2,A4)') prefix_output
788 END DO
789 END IF
790 END IF
791 ! Print the reference atomic geometry for the exciton descriptors
792 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
793 'With the center of charge as reference point r_0,'
794 WRITE (unit_nr, '(T2,A4,T15,A7,F10.4,A2,F10.4,A2,F10.4,A1)') prefix_output, &
795 'r_0 = (', ref_point_multipole(1)*angstrom, ', ', ref_point_multipole(2)*angstrom, ', ', &
796 ref_point_multipole(3)*angstrom, ')'
797 IF (exc_descr(1)%flag_TDA) THEN
798 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
799 'we further obtain r_e and r_h from solving the ', method_name, ' within the TDA'
800 ELSE
801 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
802 'we further obtain r_e and r_h from solving the ', method_name, ' without the TDA'
803 END IF
804 WRITE (unit_nr, '(T2,A4)') prefix_output
805 WRITE (unit_nr, '(T2,A4,T8,A12,1X,13X,A9,13X,A9,13X,A9)') prefix_output, &
806 'Excitation n', 'x_e [Å]', 'y_e [Å]', 'z_e [Å]'
807 DO i_exc = 1, num_print_exc_descr
808 WRITE (unit_nr, '(T2,A4,T8,I12,1X,3(5X,F15.4))') &
809 prefix_output, i_exc, &
810 exc_descr(i_exc)%r_e_shift(:)*angstrom
811 END DO
812 WRITE (unit_nr, '(T2,A4)') prefix_output
813 WRITE (unit_nr, '(T2,A4,T8,A12,1X,13X,A9,13X,A9,13X,A9)') prefix_output, &
814 'Excitation n', 'x_h [Å]', 'y_h [Å]', 'z_h [Å]'
815 DO i_exc = 1, num_print_exc_descr
816 WRITE (unit_nr, '(T2,A4,T8,I12,1X,3(5X,F15.4))') &
817 prefix_output, i_exc, &
818 exc_descr(i_exc)%r_h_shift(:)*angstrom
819 END DO
820 WRITE (unit_nr, '(T2,A4)') prefix_output
821 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
822 'The reference atomic geometry for these values is given by'
823 END IF
824 CALL write_qs_particle_coordinates_bse(particle_set, unit_nr, prefix_output)
825 IF (unit_nr > 0) THEN
826 WRITE (unit_nr, '(T2,A4)') prefix_output
827 END IF
828 CALL timestop(handle)
829
830 END SUBROUTINE print_exciton_descriptors
831
832! **************************************************************************************************
833!> \brief Debug function to write elements of a full matrix to file, if they are larger than a given threshold
834!> \param fm ...
835!> \param thresh ...
836!> \param header ...
837!> \param unit_nr ...
838!> \param abs_vals ...
839! **************************************************************************************************
840 SUBROUTINE fm_write_thresh(fm, thresh, header, unit_nr, abs_vals)
841
842 TYPE(cp_fm_type), INTENT(IN) :: fm
843 REAL(kind=dp), INTENT(IN) :: thresh
844 CHARACTER(LEN=*), INTENT(IN) :: header
845 INTEGER, INTENT(IN) :: unit_nr
846 LOGICAL, OPTIONAL :: abs_vals
847
848 CHARACTER(LEN=*), PARAMETER :: my_footer = " | ENDING WRITING OF MATRIX", &
849 routinen = 'fm_write_thresh'
850
851 INTEGER :: handle, i, j, ncol_local, nrow_local
852 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
853 LOGICAL :: my_abs_vals
854
855 CALL timeset(routinen, handle)
856
857 IF (PRESENT(abs_vals)) THEN
858 my_abs_vals = abs_vals
859 ELSE
860 my_abs_vals = .false.
861 END IF
862
863 CALL cp_fm_get_info(matrix=fm, &
864 nrow_local=nrow_local, &
865 ncol_local=ncol_local, &
866 row_indices=row_indices, &
867 col_indices=col_indices)
868
869 IF (unit_nr > 0) THEN
870 WRITE (unit_nr, *) header
871 END IF
872 IF (my_abs_vals) THEN
873 DO i = 1, nrow_local
874 DO j = 1, ncol_local
875 IF (abs(fm%local_data(i, j)) > thresh) THEN
876 IF (unit_nr > 0) THEN
877 WRITE (unit_nr, "(A7,T10,I5,T20,I5,T30,F13.5)") header, row_indices(i), col_indices(j), &
878 abs(fm%local_data(i, j))
879 END IF
880 END IF
881 END DO
882 END DO
883 ELSE
884 DO i = 1, nrow_local
885 DO j = 1, ncol_local
886 IF (abs(fm%local_data(i, j)) > thresh) THEN
887 IF (unit_nr > 0) THEN
888 WRITE (unit_nr, "(A7,T10,I5,T20,I5,T30,F13.5)") header, row_indices(i), col_indices(j), &
889 fm%local_data(i, j)
890 END IF
891 END IF
892 END DO
893 END DO
894 END IF
895 CALL fm%matrix_struct%para_env%sync()
896 IF (unit_nr > 0) THEN
897 WRITE (unit_nr, *) my_footer
898 END IF
899
900 CALL timestop(handle)
901
902 END SUBROUTINE fm_write_thresh
903
904! **************************************************************************************************
905!> \brief Write the atomic coordinates to the output unit.
906!> \param particle_set ...
907!> \note Adapted from particle_methods.F [MG]
908!> \param unit_nr ...
909!> \param prefix_output ...
910! **************************************************************************************************
911 SUBROUTINE write_qs_particle_coordinates_bse(particle_set, unit_nr, prefix_output)
912
913 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
914 INTEGER, INTENT(IN) :: unit_nr
915 CHARACTER(LEN=4), INTENT(IN) :: prefix_output
916
917 CHARACTER(len=*), PARAMETER :: routinen = 'write_qs_particle_coordinates_bse'
918
919 CHARACTER(LEN=2) :: element_symbol
920 INTEGER :: handle, iatom, natom
921
922 CALL timeset(routinen, handle)
923
924 IF (unit_nr > 0) THEN
925 WRITE (unit_nr, '(T2,A4)') prefix_output
926 WRITE (unit_nr, '(T2,A4,T13,A7,16X,A7,15X,A7,15X,A7)') prefix_output, &
927 'Element', 'x [Å]', 'y [Å]', 'z [Å]'
928 natom = SIZE(particle_set)
929 DO iatom = 1, natom
930 CALL get_atomic_kind(atomic_kind=particle_set(iatom)%atomic_kind, &
931 element_symbol=element_symbol)
932 WRITE (unit_nr, '(T2,A4,T8,A12,1X,3(5X,F15.4))') &
933 prefix_output, element_symbol, particle_set(iatom)%r(1:3)*angstrom
934 END DO
935 END IF
936
937 CALL timestop(handle)
938
939 END SUBROUTINE write_qs_particle_coordinates_bse
940
941END MODULE bse_print
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
atomic_kind_types::get_atomic_kind
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
Definition atomic_kind_types.F:138
bibliography
collects all references to literature in CP2K as new algorithms / method are included from literature...
Definition bibliography.F:21
bibliography::graml2026
integer, save, public graml2026
Definition bibliography.F:36
bse_print
Routines for printing information in context of the BSE calculation.
Definition bse_print.F:13
bse_print::print_output_header
subroutine, public print_output_header(homo, virtual, homo_irred, flag_tda, multiplet, alpha, mp2_env, unit_nr)
...
Definition bse_print.F:102
bse_print::print_bse_start_flag
subroutine, public print_bse_start_flag(bse_tda, bse_abba, unit_nr)
...
Definition bse_print.F:56
bse_print::print_optical_properties
subroutine, public print_optical_properties(exc_ens, oscill_str, trans_mom_bse, polarizability_residues, homo, virtual, homo_irred, flag_tda, info_approximation, mp2_env, unit_nr, open_shell)
...
Definition bse_print.F:366
bse_print::print_transition_amplitudes
subroutine, public print_transition_amplitudes(fm_eigvec_x, homo, virtual, homo_irred, info_approximation, mp2_env, unit_nr, fm_eigvec_y)
...
Definition bse_print.F:271
bse_print::fm_write_thresh
subroutine, public fm_write_thresh(fm, thresh, header, unit_nr, abs_vals)
Debug function to write elements of a full matrix to file, if they are larger than a given threshold.
Definition bse_print.F:841
bse_print::print_excitation_energies
subroutine, public print_excitation_energies(exc_ens, homo, virtual, flag_tda, multiplet, info_approximation, mp2_env, unit_nr)
...
Definition bse_print.F:222
bse_print::print_exciton_descriptors
subroutine, public print_exciton_descriptors(exc_descr, ref_point_multipole, unit_nr, num_print_exc_descr, print_checkvalue, print_directional_exc_descr, print_directional_crosscorrelation, prefix_output, qs_env)
Prints exciton descriptors, cf. Mewes et al., JCTC 14, 710-725 (2018).
Definition bse_print.F:574
bse_properties
Routines for computing excitonic properties, e.g. exciton diameter, from the BSE.
Definition bse_properties.F:13
bse_properties::compute_and_print_absorption_spectrum
subroutine, public compute_and_print_absorption_spectrum(oscill_str, polarizability_residues, exc_ens, info_approximation, unit_nr, mp2_env)
Computes and returns absorption spectrum for the frequency range and broadening provided by the user....
Definition bse_properties.F:228
bse_util
Auxiliary routines for GW + Bethe-Salpeter for computing electronic excitations.
Definition bse_util.F:13
bse_util::filter_eigvec_contrib
subroutine, public filter_eigvec_contrib(fm_eigvec, idx_homo, idx_virt, eigvec_entries, i_exc, virtual, num_entries, mp2_env, offsets, virtual_per_spin, idx_spin)
Filters eigenvector entries above a given threshold to describe excitations in the singleparticle bas...
Definition bse_util.F:1027
bse_util::get_bse_spin_block_layout
subroutine, public get_bse_spin_block_layout(homo_red, virt_red, n_ov, offsets, n_ov_joint)
Spin-block layout for the open-shell (joint) BSE matrix: per-spin OV-pair counts and the block offset...
Definition bse_util.F:1180
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_get_info
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
Definition cp_fm_types.F:1038
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:18
input_constants::bse_screening_tdhf
integer, parameter, public bse_screening_tdhf
Definition input_constants.F:1398
input_constants::bse_screening_w0
integer, parameter, public bse_screening_w0
Definition input_constants.F:1398
input_constants::bse_tda
integer, parameter, public bse_tda
Definition input_constants.F:1393
input_constants::bse_screening_alpha
integer, parameter, public bse_screening_alpha
Definition input_constants.F:1398
input_constants::bse_screening_rpa
integer, parameter, public bse_screening_rpa
Definition input_constants.F:1398
input_constants::bse_abba
integer, parameter, public bse_abba
Definition input_constants.F:1393
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
mp2_types
Types needed for MP2 calculations.
Definition mp2_types.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
physcon::angstrom
real(kind=dp), parameter, public angstrom
Definition physcon.F:144
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
bse_properties::exciton_descr_type
Definition bse_properties.F:67
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:114
mp2_types::mp2_type
Definition mp2_types.F:343
particle_types::particle_type
Definition particle_types.F:35
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:224