(git:98357aa)
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! **************************************************************************************************
14
16 USE bibliography, ONLY: graml2026,&
17 cite_reference
22 USE cp_fm_types, ONLY: cp_fm_get_info,&
28 USE kinds, ONLY: dp
29 USE mp2_types, ONLY: mp2_type
31 USE physcon, ONLY: angstrom,&
32 evolt
35#include "./base/base_uses.f90"
36
37 IMPLICIT NONE
38
39 PRIVATE
40
41 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'bse_print'
42
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 prefix_output String, which is put in front of prints, i.e. "BSE|" or "" for TDDFPT
568!> \param qs_env ...
569! **************************************************************************************************
570 SUBROUTINE print_exciton_descriptors(exc_descr, ref_point_multipole, unit_nr, &
571 num_print_exc_descr, print_checkvalue, print_directional_exc_descr, &
572 prefix_output, qs_env)
573
574 TYPE(exciton_descr_type), ALLOCATABLE, &
575 DIMENSION(:) :: exc_descr
576 REAL(kind=dp), ALLOCATABLE, DIMENSION(:), &
577 INTENT(IN) :: ref_point_multipole
578 INTEGER, INTENT(IN) :: unit_nr, num_print_exc_descr
579 LOGICAL, INTENT(IN) :: print_checkvalue, &
580 print_directional_exc_descr
581 CHARACTER(LEN=4), INTENT(IN) :: prefix_output
582 TYPE(qs_environment_type), POINTER :: qs_env
583
584 CHARACTER(LEN=*), PARAMETER :: routinen = 'print_exciton_descriptors'
585
586 CHARACTER(LEN=1), DIMENSION(3) :: array_direction_str
587 CHARACTER(LEN=5) :: method_name
588 INTEGER :: handle, i_dir, i_exc
589 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
590
591 IF (prefix_output == 'BSE|') THEN
592 method_name = 'BSE'
593 ELSE
594 method_name = 'TDDFT'
595 END IF
596
597 CALL timeset(routinen, handle)
598 CALL get_qs_env(qs_env, particle_set=particle_set)
599 IF (unit_nr > 0) THEN
600 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
601 'Exciton descriptors for excitation level n are given by'
602 WRITE (unit_nr, '(T2,A4)') prefix_output
603 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
604 'd_eh = | <r_h - r_e>_exc |'
605 WRITE (unit_nr, '(T2,A4)') prefix_output
606 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
607 σ'_e = sqrt( <r_e^2>_exc - <r_e>_exc^2 )'
608 WRITE (unit_nr, '(T2,A4)') prefix_output
609 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
610 σ'_h = sqrt( <r_h^2>_exc - <r_h>_exc^2 )'
611 WRITE (unit_nr, '(T2,A4)') prefix_output
612 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
613 'COV_eh = <r_e r_h>_exc - <r_e>_exc <r_h>_exc'
614 WRITE (unit_nr, '(T2,A4)') prefix_output
615 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
616 'd_exc = sqrt( | < |r_h - r_e|^2 >_exc )'
617 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
618 σσ' = sqrt( d_eh^2 + _e^2 + _h^2 - 2 * COV_eh )'
619 WRITE (unit_nr, '(T2,A4)') prefix_output
620 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
621 σσ'R_eh = COV_eh / (_e * _h)'
622 WRITE (unit_nr, '(T2,A4)') prefix_output
623 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
624 'where the expectation values <.>_exc are taken with respect to the '
625 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
626 'exciton wavefunction of excitation n:'
627 WRITE (unit_nr, '(T2,A4)') prefix_output
628
629 IF (exc_descr(1)%flag_TDA) THEN
630 WRITE (unit_nr, '(T2,A4,T20,A)') prefix_output, &
631 𝚿Σψψ'_n(r_e,r_h) = _{i,a} X_ia^n _i(r_h) _a(r_e) ,'
632 ELSE
633 WRITE (unit_nr, '(T2,A4,T20,A)') prefix_output, &
634 𝚿Σψψ'_n(r_e,r_h) = _{i,a} X_ia^n _i(r_h) _a(r_e)'
635 WRITE (unit_nr, '(T2,A4,T40,A)') prefix_output, &
636 ψψ'+ Y_ia^n _a(r_h) _i(r_e) ,'
637 END IF
638 WRITE (unit_nr, '(T2,A4)') prefix_output
639 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
640 'i.e.'
641 WRITE (unit_nr, '(T2,A4)') prefix_output
642 WRITE (unit_nr, '(T2,A4,T20,A)') prefix_output, &
643 𝚿𝚿𝚿𝚿'< O >_exc = < _n | O | _n > / < _n | _n > ,'
644 WRITE (unit_nr, '(T2,A4)') prefix_output
645 IF (exc_descr(1)%flag_TDA) THEN
646 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
647 𝚿𝚿'where c_n = < _n | _n > = 1 within TDA.'
648 ELSE
649 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
650 𝚿𝚿≥'where c_n = < _n | _n > 1 without TDA.'
651 END IF
652 WRITE (unit_nr, '(T2,A4)') prefix_output
653 WRITE (unit_nr, '(T2,A4)') prefix_output
654 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
655 'Here, we introduced'
656 WRITE (unit_nr, '(T2,A4)') prefix_output
657 WRITE (unit_nr, '(T2,A4,T7,A5,T15,A)') &
658 prefix_output, ψ"_i:", "occupied molecular orbitals,"
659 WRITE (unit_nr, '(T2,A4,T7,A5,T15,A)') &
660 prefix_output, ψ"_a:", "empty molecular orbitals and"
661 WRITE (unit_nr, '(T2,A4,T9,A2,T14,A)') &
662 prefix_output, "r:", "position operator."
663 WRITE (unit_nr, '(T2,A4)') prefix_output
664 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
665 'Ref.: Eqs. (28)-(32)'
666 WRITE (unit_nr, '(T2,A4,T7,A,A)') prefix_output, &
667 'Phys. Rev. B 113, 205152 (2026); ', &
668 'https://doi.org/10.1103/38k2-d55h'
669 WRITE (unit_nr, '(T2,A4)') prefix_output
670 WRITE (unit_nr, '(T2,A4)') prefix_output
671 IF (exc_descr(1)%flag_TDA) THEN
672 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
673 'Exciton descriptors from solving the ', method_name, ' within the TDA:'
674 ELSE
675 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
676 'Exciton descriptors from solving the ', method_name, ' without the TDA:'
677 END IF
678 WRITE (unit_nr, '(T2,A4)') prefix_output
679 WRITE (unit_nr, '(T2,A4,T10,A1,6X,A3,1X,4X,A10,5X,A10,5X,A10,3X,A11,8X,A4)') prefix_output, &
680 'n', 'c_n', Å'd_eh []', σÅ'_e []', σÅ'_h []', Å'd_exc []', 'R_eh'
681 DO i_exc = 1, num_print_exc_descr
682 WRITE (unit_nr, '(T2,A4,T7,I4,4X,F5.3,1X,5(2X,F10.4))') &
683 prefix_output, i_exc, exc_descr(i_exc)%norm_XpY, &
684 exc_descr(i_exc)%diff_r_abs*angstrom, &
685 exc_descr(i_exc)%sigma_e*angstrom, exc_descr(i_exc)%sigma_h*angstrom, &
686 exc_descr(i_exc)%diff_r_sqr*angstrom, exc_descr(i_exc)%corr_e_h
687 END DO
688 WRITE (unit_nr, '(T2,A4)') prefix_output
689 ! For debug runs, print first d_exc separately to allow the regtests to read in
690 IF (print_checkvalue) THEN
691 WRITE (unit_nr, '(T2)')
692 WRITE (unit_nr, '(T2,A28,T65,F16.4)') 'Checksum exciton descriptors', &
693 exc_descr(1)%diff_r_sqr*angstrom
694 WRITE (unit_nr, '(T2)')
695 END IF
696 WRITE (unit_nr, '(T2,A4)') prefix_output
697 ! Print exciton descriptor resolved per direction
698 IF (print_directional_exc_descr) THEN
699 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
700 'We can restrict the exciton descriptors to a specific direction,'
701 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
702 'e.g. the x-components are:'
703 WRITE (unit_nr, '(T2,A4)') prefix_output
704 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
705 'd_eh^x = | <x_h - x_e>_exc |'
706 WRITE (unit_nr, '(T2,A4)') prefix_output
707 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
708 σ'_e^x = sqrt( <x_e^2>_exc - <x_e>_exc^2 )'
709 WRITE (unit_nr, '(T2,A4)') prefix_output
710 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
711 σ'_h^x = sqrt( <x_h^2>_exc - <x_h>_exc^2 )'
712 WRITE (unit_nr, '(T2,A4)') prefix_output
713 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
714 μμμμμμ"COV_eh^{'} = <r^_e r^'_h>_exc - <r^_e>_exc <r^'_h>_exc"
715 WRITE (unit_nr, '(T2,A4)') prefix_output
716 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
717 'd_exc^x = sqrt( | < |x_h - x_e|^2 >_exc )'
718 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
719 σ' = sqrt( (d_eh^x)^2 + (_e^x)^2'
720 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
721 σ" + (_h^x)^2 - 2 * (COV_eh^{xx}) )"
722 WRITE (unit_nr, '(T2,A4)') prefix_output
723 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
724 μμ"Subsequently, the cross-correlation matrix R_eh^{'} is printed"
725 WRITE (unit_nr, '(T2,A4)') prefix_output
726 WRITE (unit_nr, '(T2,A4,T15,A)') prefix_output, &
727 μμμμσμσμ"R_eh^{'} = COV_eh^{'}/(^_e ^_h) "
728 WRITE (unit_nr, '(T2,A4)') prefix_output
729 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
730 'Ref.: Eqs. (33)-(36)'
731 WRITE (unit_nr, '(T2,A4,T7,A,A)') prefix_output, &
732 'Phys. Rev. B 113, 205152 (2026); ', &
733 'https://doi.org/10.1103/38k2-d55h'
734 WRITE (unit_nr, '(T2,A4)') prefix_output
735 WRITE (unit_nr, '(T2,A4)') prefix_output
736 IF (exc_descr(1)%flag_TDA) THEN
737 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
738 'Exciton descriptors per direction from solving the ', method_name, ' within the TDA:'
739 ELSE
740 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
741 'Exciton descriptors per direction from solving the ', method_name, ' without the TDA:'
742 END IF
743 WRITE (unit_nr, '(T2,A4)') prefix_output
744 WRITE (unit_nr, '(T2,A4,T12,A1,2X,A9,5X,A12,5X,A12,5X,A12,3X,A13)') prefix_output, &
745 'n', 'r = x/y/z', Å'd_eh^r []', σÅ'_e^r []', σÅ'_h^r []', Å'd_exc^r []'
746 DO i_exc = 1, num_print_exc_descr
747 DO i_dir = 1, 3
748 array_direction_str = ["x", "y", "z"]
749 WRITE (unit_nr, '(T2,A4,T9,I4,10X,A1,1X,4(4X,F10.4))') &
750 prefix_output, i_exc, array_direction_str(i_dir), &
751 exc_descr(i_exc)%d_eh_dir(i_dir)*angstrom, &
752 exc_descr(i_exc)%sigma_e_dir(i_dir)*angstrom, &
753 exc_descr(i_exc)%sigma_h_dir(i_dir)*angstrom, &
754 exc_descr(i_exc)%d_exc_dir(i_dir)*angstrom
755 END DO
756 WRITE (unit_nr, '(T2,A4)') prefix_output
757 END DO
758 WRITE (unit_nr, '(T2,A4)') prefix_output
759 WRITE (unit_nr, '(T2,A4)') prefix_output
760 IF (exc_descr(1)%flag_TDA) THEN
761 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
762 'Crosscorrelation matrix from solving the ', method_name, ' within the TDA:'
763 ELSE
764 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
765 'Crosscorrelation matrix from solving the ', method_name, ' without the TDA:'
766 END IF
767 WRITE (unit_nr, '(T2,A4)') prefix_output
768 WRITE (unit_nr, '(T2,A4,T12,A1,8X,6(8X,A2))') prefix_output, &
769 'n', 'xx', 'yy', 'zz', 'xy', 'xz', 'yz'
770 DO i_exc = 1, num_print_exc_descr
771 WRITE (unit_nr, '(T2,A4,T9,I4,8X,6(3X,F7.4),3X,F7.4)') &
772 prefix_output, i_exc, &
773 exc_descr(i_exc)%corr_e_h_matrix(1, 1), &
774 exc_descr(i_exc)%corr_e_h_matrix(2, 2), &
775 exc_descr(i_exc)%corr_e_h_matrix(3, 3), &
776 exc_descr(i_exc)%corr_e_h_matrix(1, 2), &
777 exc_descr(i_exc)%corr_e_h_matrix(1, 3), &
778 exc_descr(i_exc)%corr_e_h_matrix(2, 3)
779 END DO
780 WRITE (unit_nr, '(T2,A4)') prefix_output
781 END IF
782 ! Print the reference atomic geometry for the exciton descriptors
783 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
784 'With the center of charge as reference point r_0,'
785 WRITE (unit_nr, '(T2,A4,T15,A7,F10.4,A2,F10.4,A2,F10.4,A1)') prefix_output, &
786 'r_0 = (', ref_point_multipole(1)*angstrom, ', ', ref_point_multipole(2)*angstrom, ', ', &
787 ref_point_multipole(3)*angstrom, ')'
788 IF (exc_descr(1)%flag_TDA) THEN
789 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
790 'we further obtain r_e and r_h from solving the ', method_name, ' within the TDA'
791 ELSE
792 WRITE (unit_nr, '(T2,A4,T7,A,A,A)') prefix_output, &
793 'we further obtain r_e and r_h from solving the ', method_name, ' without the TDA'
794 END IF
795 WRITE (unit_nr, '(T2,A4)') prefix_output
796 WRITE (unit_nr, '(T2,A4,T8,A12,1X,13X,A9,13X,A9,13X,A9)') prefix_output, &
797 'Excitation n', Å'x_e []', Å'y_e []', Å'z_e []'
798 DO i_exc = 1, num_print_exc_descr
799 WRITE (unit_nr, '(T2,A4,T8,I12,1X,3(5X,F15.4))') &
800 prefix_output, i_exc, &
801 exc_descr(i_exc)%r_e_shift(:)*angstrom
802 END DO
803 WRITE (unit_nr, '(T2,A4)') prefix_output
804 WRITE (unit_nr, '(T2,A4,T8,A12,1X,13X,A9,13X,A9,13X,A9)') prefix_output, &
805 'Excitation n', Å'x_h []', Å'y_h []', Å'z_h []'
806 DO i_exc = 1, num_print_exc_descr
807 WRITE (unit_nr, '(T2,A4,T8,I12,1X,3(5X,F15.4))') &
808 prefix_output, i_exc, &
809 exc_descr(i_exc)%r_h_shift(:)*angstrom
810 END DO
811 WRITE (unit_nr, '(T2,A4)') prefix_output
812 WRITE (unit_nr, '(T2,A4,T7,A)') prefix_output, &
813 'The reference atomic geometry for these values is given by'
814 END IF
815 CALL write_qs_particle_coordinates_bse(particle_set, unit_nr, prefix_output)
816 IF (unit_nr > 0) THEN
817 WRITE (unit_nr, '(T2,A4)') prefix_output
818 END IF
819 CALL timestop(handle)
820
821 END SUBROUTINE print_exciton_descriptors
822
823! **************************************************************************************************
824!> \brief Debug function to write elements of a full matrix to file, if they are larger than a given threshold
825!> \param fm ...
826!> \param thresh ...
827!> \param header ...
828!> \param unit_nr ...
829!> \param abs_vals ...
830! **************************************************************************************************
831 SUBROUTINE fm_write_thresh(fm, thresh, header, unit_nr, abs_vals)
832
833 TYPE(cp_fm_type), INTENT(IN) :: fm
834 REAL(kind=dp), INTENT(IN) :: thresh
835 CHARACTER(LEN=*), INTENT(IN) :: header
836 INTEGER, INTENT(IN) :: unit_nr
837 LOGICAL, OPTIONAL :: abs_vals
838
839 CHARACTER(LEN=*), PARAMETER :: my_footer = " | ENDING WRITING OF MATRIX", &
840 routinen = 'fm_write_thresh'
841
842 INTEGER :: handle, i, j, ncol_local, nrow_local
843 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
844 LOGICAL :: my_abs_vals
845
846 CALL timeset(routinen, handle)
847
848 IF (PRESENT(abs_vals)) THEN
849 my_abs_vals = abs_vals
850 ELSE
851 my_abs_vals = .false.
852 END IF
853
854 CALL cp_fm_get_info(matrix=fm, &
855 nrow_local=nrow_local, &
856 ncol_local=ncol_local, &
857 row_indices=row_indices, &
858 col_indices=col_indices)
859
860 IF (unit_nr > 0) THEN
861 WRITE (unit_nr, *) header
862 END IF
863 IF (my_abs_vals) THEN
864 DO i = 1, nrow_local
865 DO j = 1, ncol_local
866 IF (abs(fm%local_data(i, j)) > thresh) THEN
867 IF (unit_nr > 0) THEN
868 WRITE (unit_nr, "(A7,T10,I5,T20,I5,T30,F13.5)") header, row_indices(i), col_indices(j), &
869 abs(fm%local_data(i, j))
870 END IF
871 END IF
872 END DO
873 END DO
874 ELSE
875 DO i = 1, nrow_local
876 DO j = 1, ncol_local
877 IF (abs(fm%local_data(i, j)) > thresh) THEN
878 IF (unit_nr > 0) THEN
879 WRITE (unit_nr, "(A7,T10,I5,T20,I5,T30,F13.5)") header, row_indices(i), col_indices(j), &
880 fm%local_data(i, j)
881 END IF
882 END IF
883 END DO
884 END DO
885 END IF
886 CALL fm%matrix_struct%para_env%sync()
887 IF (unit_nr > 0) THEN
888 WRITE (unit_nr, *) my_footer
889 END IF
890
891 CALL timestop(handle)
892
893 END SUBROUTINE fm_write_thresh
894
895! **************************************************************************************************
896!> \brief Write the atomic coordinates to the output unit.
897!> \param particle_set ...
898!> \note Adapted from particle_methods.F [MG]
899!> \param unit_nr ...
900!> \param prefix_output ...
901! **************************************************************************************************
902 SUBROUTINE write_qs_particle_coordinates_bse(particle_set, unit_nr, prefix_output)
903
904 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
905 INTEGER, INTENT(IN) :: unit_nr
906 CHARACTER(LEN=4), INTENT(IN) :: prefix_output
907
908 CHARACTER(len=*), PARAMETER :: routinen = 'write_qs_particle_coordinates_bse'
909
910 CHARACTER(LEN=2) :: element_symbol
911 INTEGER :: handle, iatom, natom
912
913 CALL timeset(routinen, handle)
914
915 IF (unit_nr > 0) THEN
916 WRITE (unit_nr, '(T2,A4)') prefix_output
917 WRITE (unit_nr, '(T2,A4,T13,A7,16X,A7,15X,A7,15X,A7)') prefix_output, &
918 'Element', Å'x []', Å'y []', Å'z []'
919 natom = SIZE(particle_set)
920 DO iatom = 1, natom
921 CALL get_atomic_kind(atomic_kind=particle_set(iatom)%atomic_kind, &
922 element_symbol=element_symbol)
923 WRITE (unit_nr, '(T2,A4,T8,A12,1X,3(5X,F15.4))') &
924 prefix_output, element_symbol, particle_set(iatom)%r(1:3)*angstrom
925 END DO
926 END IF
927
928 CALL timestop(handle)
929
930 END SUBROUTINE write_qs_particle_coordinates_bse
931
932END MODULE bse_print
Define the atomic kind types and their sub types.
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.
collects all references to literature in CP2K as new algorithms / method are included from literature...
integer, save, public graml2026
Routines for printing information in context of the BSE calculation.
Definition bse_print.F:13
subroutine, public print_output_header(homo, virtual, homo_irred, flag_tda, multiplet, alpha, mp2_env, unit_nr)
...
Definition bse_print.F:102
subroutine, public print_bse_start_flag(bse_tda, bse_abba, unit_nr)
...
Definition bse_print.F:56
subroutine, public print_exciton_descriptors(exc_descr, ref_point_multipole, unit_nr, num_print_exc_descr, print_checkvalue, print_directional_exc_descr, prefix_output, qs_env)
Prints exciton descriptors, cf. Mewes et al., JCTC 14, 710-725 (2018)
Definition bse_print.F:573
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
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
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:832
subroutine, public print_excitation_energies(exc_ens, homo, virtual, flag_tda, multiplet, info_approximation, mp2_env, unit_nr)
...
Definition bse_print.F:222
Routines for computing excitonic properties, e.g. exciton diameter, from the BSE.
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....
Auxiliary routines for GW + Bethe-Salpeter for computing electronic excitations.
Definition bse_util.F:13
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:1026
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:1262
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
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
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public bse_screening_tdhf
integer, parameter, public bse_screening_w0
integer, parameter, public bse_tda
integer, parameter, public bse_screening_alpha
integer, parameter, public bse_screening_rpa
integer, parameter, public bse_abba
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Types needed for MP2 calculations.
Definition mp2_types.F:14
Define the data structure for the particle information.
Definition of physical constants:
Definition physcon.F:68
real(kind=dp), parameter, public evolt
Definition physcon.F:183
real(kind=dp), parameter, public angstrom
Definition physcon.F:144
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.
represent a full matrix