(git:9cbee8b)
Loading...
Searching...
No Matches
qs_pdos.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 Calculation and writing of projected density of states
10!> The DOS is computed per angular momentum and per kind
11!> \par History
12!> -
13!> \author Marcella (29.02.2008,MK)
14! **************************************************************************************************
15MODULE qs_pdos
16 USE atomic_kind_types, ONLY: atomic_kind_type,&
17 get_atomic_kind,&
18 get_atomic_kind_set
19 USE basis_set_types, ONLY: get_gto_basis_set,&
20 gto_basis_set_type
21 USE cell_types, ONLY: cell_type,&
22 pbc
23 USE cp_array_utils, ONLY: cp_1d_r_p_type
24 USE cp_blacs_env, ONLY: cp_blacs_env_type
25 USE cp_cfm_types, ONLY: cp_cfm_create,&
26 cp_cfm_get_submatrix,&
27 cp_cfm_release,&
28 cp_cfm_type
29 USE cp_control_types, ONLY: dft_control_type
30 USE cp_dbcsr_api, ONLY: dbcsr_p_type
31 USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm
32 USE cp_fm_diag, ONLY: cp_fm_power
33 USE cp_fm_struct, ONLY: cp_fm_struct_create,&
34 cp_fm_struct_release,&
35 cp_fm_struct_type
36 USE cp_fm_types, ONLY: cp_fm_create,&
37 cp_fm_get_info,&
38 cp_fm_get_submatrix,&
39 cp_fm_release,&
40 cp_fm_type
41 USE cp_log_handling, ONLY: cp_get_default_logger,&
42 cp_logger_get_default_io_unit,&
43 cp_logger_type,&
44 cp_to_string
45 USE cp_output_handling, ONLY: cp_p_file,&
46 cp_print_key_finished_output,&
47 cp_print_key_should_output,&
48 cp_print_key_unit_nr
49 USE input_section_types, ONLY: section_vals_get,&
50 section_vals_get_subs_vals,&
51 section_vals_type,&
52 section_vals_val_get
53 USE kinds, ONLY: default_string_length,&
54 dp
55 USE kpoint_methods, ONLY: lowdin_kp_mo_coeff
56 USE kpoint_types, ONLY: kpoint_env_type,&
57 kpoint_type
58 USE memory_utilities, ONLY: reallocate
59 USE message_passing, ONLY: mp_para_env_type
60 USE orbital_pointers, ONLY: nso,&
61 nsoset
62 USE orbital_symbols, ONLY: l_sym,&
63 sgf_symbol
64 USE parallel_gemm_api, ONLY: parallel_gemm
65 USE particle_types, ONLY: particle_type
66 USE pw_env_types, ONLY: pw_env_get,&
67 pw_env_type
68 USE pw_pool_types, ONLY: pw_pool_p_type,&
69 pw_pool_type
70 USE pw_types, ONLY: pw_c1d_gs_type,&
71 pw_r3d_rs_type
72 USE qs_collocate_density, ONLY: calculate_wavefunction
73 USE qs_dos_utils, ONLY: &
74 add_broadened_value, broadening_cutoff, broadening_function, dos_density_scale, &
75 dos_energy_label, dos_energy_scale, dos_energy_unit_ev, dos_energy_zero_absolute, &
76 dos_energy_zero_auto, dos_energy_zero_hoco, dos_energy_zero_label, &
77 dos_resolve_energy_zero, write_broadening_info
78 USE qs_environment_types, ONLY: get_qs_env,&
79 qs_environment_type
80 USE qs_kind_types, ONLY: get_qs_kind,&
81 get_qs_kind_set,&
82 qs_kind_type
83 USE qs_mo_types, ONLY: get_mo_set,&
84 mo_set_type
85 USE qs_scf_diagonalization, ONLY: diag_kp_smat
86 USE qs_scf_types, ONLY: qs_scf_env_type
87#include "./base/base_uses.f90"
88
89 IMPLICIT NONE
90
91 PRIVATE
92
93 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_pdos'
94
95! **************************************************************************************************
96 ! *** Public subroutines ***
97
98 PUBLIC :: calculate_projected_dos, calculate_projected_dos_kp
99
100 TYPE ldos_type
101 INTEGER :: maxl = -1, nlist = -1
102 LOGICAL :: separate_components = .false.
103 INTEGER, DIMENSION(:), POINTER :: list_index => null()
104 REAL(KIND=dp), DIMENSION(:, :), &
105 POINTER :: pdos_array => null()
106 END TYPE ldos_type
107
108 TYPE r_ldos_type
109 INTEGER :: nlist = -1, npoints = -1
110 INTEGER, DIMENSION(:, :), POINTER :: index_grid_local => null()
111 INTEGER, DIMENSION(:), POINTER :: list_index => null()
112 REAL(KIND=dp), DIMENSION(:), POINTER :: x_range => null(), y_range => null(), z_range => null()
113 REAL(KIND=dp), DIMENSION(:), POINTER :: eval_range => null()
114 REAL(KIND=dp), DIMENSION(:), &
115 POINTER :: pdos_array => null()
116 END TYPE r_ldos_type
117
118 TYPE ldos_p_type
119 TYPE(ldos_type), POINTER :: ldos => null()
120 END TYPE ldos_p_type
121
122 TYPE r_ldos_p_type
123 TYPE(r_ldos_type), POINTER :: ldos => null()
124 END TYPE r_ldos_p_type
125CONTAINS
126
127! **************************************************************************************************
128!> \brief Compute and write projected density of states
129!> \param mo_set ...
130!> \param atomic_kind_set ...
131!> \param qs_kind_set ...
132!> \param particle_set ...
133!> \param qs_env ...
134!> \param dft_section ...
135!> \param ispin ...
136!> \param xas_mittle ...
137!> \param external_matrix_shalf ...
138!> \param unoccupied_orbs ...
139!> \param unoccupied_evals ...
140!> \param pdos_print_key ...
141!> \param write_pdos ...
142!> \param write_pdos_curve ...
143!> \date 26.02.2008
144!> \par History:
145!> - Added optional external matrix_shalf to avoid recomputing it (A. Bussy, 09.2019)
146!> \par Variables
147!> -
148!> -
149!> \author MI
150!> \version 1.0
151! **************************************************************************************************
152 SUBROUTINE calculate_projected_dos(mo_set, atomic_kind_set, qs_kind_set, particle_set, qs_env, &
153 dft_section, ispin, xas_mittle, external_matrix_shalf, &
154 unoccupied_orbs, unoccupied_evals, pdos_print_key, write_pdos, write_pdos_curve)
155
156 TYPE(mo_set_type), INTENT(IN) :: mo_set
157 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
158 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
159 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
160 TYPE(qs_environment_type), POINTER :: qs_env
161 TYPE(section_vals_type), POINTER :: dft_section
162 INTEGER, INTENT(IN), OPTIONAL :: ispin
163 CHARACTER(LEN=default_string_length), INTENT(IN), &
164 OPTIONAL :: xas_mittle
165 TYPE(cp_fm_type), INTENT(IN), OPTIONAL, TARGET :: external_matrix_shalf, unoccupied_orbs
166 TYPE(cp_1d_r_p_type), INTENT(IN), OPTIONAL, TARGET :: unoccupied_evals
167 CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: pdos_print_key
168 LOGICAL, INTENT(IN), OPTIONAL :: write_pdos, write_pdos_curve
169
170 CHARACTER(len=*), PARAMETER :: routinen = 'calculate_projected_dos'
171
172 CHARACTER(LEN=16) :: energy_label, fmtstr2
173 CHARACTER(LEN=27) :: fmtstr1
174 CHARACTER(LEN=32) :: zero_label
175 CHARACTER(LEN=6), ALLOCATABLE, DIMENSION(:, :, :) :: tmp_str
176 CHARACTER(LEN=default_string_length) :: kind_name, my_act, my_mittle, my_pos, &
177 my_print_key, spin(2)
178 CHARACTER(LEN=default_string_length), &
179 ALLOCATABLE, DIMENSION(:) :: ldos_index, r_ldos_index
180 INTEGER :: broaden_type, energy_unit, energy_zero, handle, homo, i, iatom, ikind, il, ildos, &
181 im, imo, imo_ref, in_x, in_y, in_z, ir, irow, iset, isgf, ishell, iso, ispin_ref, &
182 iterstep, iw, j, jx, jy, jz, k, lcomponent, lshell, maxl, maxlgto, my_spin, n_dependent, &
183 n_r_ldos, n_rep, nao, natom, ncol_global, ndigits, nkind, nldos, nmo, nmo_ref, np_tot, &
184 npoints, nrow_global, nset, nsgf, nvirt, out_each, output_unit, resolved_energy_zero
185 INTEGER, ALLOCATABLE, DIMENSION(:) :: firstrow
186 INTEGER, DIMENSION(:), POINTER :: list, nshell
187 INTEGER, DIMENSION(:, :), POINTER :: bo, l
188 LOGICAL :: append, calc_matsh, do_curve, do_ldos, do_r_ldos, do_virt, fractional_occupation, &
189 ionode, separate_components, should_output, write_curve, write_pdos_file
190 LOGICAL, DIMENSION(:, :), POINTER :: read_r
191 REAL(kind=dp) :: broaden_width, de, dh(3, 3), dvol, e_fermi, e_fermi_ref(2), energy_factor, &
192 energy_ref, ev_factor, hoco, hoco_ref(2), r(3), r_vec(3), ratom(3), voigt_mixing
193 REAL(kind=dp), DIMENSION(:), POINTER :: eigenvalues, eval_ref, evals_virt, &
194 occ_ref, occupation_numbers
195 REAL(kind=dp), DIMENSION(:, :), POINTER :: vecbuffer
196 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: pdos_array
197 TYPE(cell_type), POINTER :: cell
198 TYPE(cp_blacs_env_type), POINTER :: context
199 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
200 TYPE(cp_fm_type) :: matrix_shalfc, matrix_work
201 TYPE(cp_fm_type), POINTER :: matrix_shalf, mo_coeff, mo_virt
202 TYPE(cp_logger_type), POINTER :: logger
203 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: s_matrix
204 TYPE(dft_control_type), POINTER :: dft_control
205 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
206 TYPE(ldos_p_type), DIMENSION(:), POINTER :: ldos_p
207 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos_ref
208 TYPE(mp_para_env_type), POINTER :: para_env
209 TYPE(pw_c1d_gs_type) :: wf_g
210 TYPE(pw_env_type), POINTER :: pw_env
211 TYPE(pw_pool_p_type), DIMENSION(:), POINTER :: pw_pools
212 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
213 TYPE(pw_r3d_rs_type) :: wf_r
214 TYPE(r_ldos_p_type), DIMENSION(:), POINTER :: r_ldos_p
215 TYPE(section_vals_type), POINTER :: curve_section, ldos_section
216
217 NULLIFY (logger, mos_ref, eval_ref, occ_ref)
218 logger => cp_get_default_logger()
219 ionode = logger%para_env%is_source()
220 my_print_key = "PRINT%PDOS"
221 IF (PRESENT(pdos_print_key)) my_print_key = trim(pdos_print_key)
222 write_pdos_file = .true.
223 IF (PRESENT(write_pdos)) write_pdos_file = write_pdos
224 curve_section => section_vals_get_subs_vals(dft_section, trim(my_print_key)//"%CURVE")
225 CALL section_vals_get(curve_section, explicit=write_curve)
226 IF (PRESENT(write_pdos_curve)) write_curve = write_pdos_curve
227 should_output = btest(cp_print_key_should_output(logger%iter_info, dft_section, &
228 trim(my_print_key)), cp_p_file)
229 output_unit = cp_logger_get_default_io_unit(logger)
230
231 spin(1) = "ALPHA"
232 spin(2) = "BETA"
233 IF ((.NOT. should_output)) RETURN
234
235 NULLIFY (context, s_matrix, orb_basis_set, para_env, pdos_array)
236 NULLIFY (eigenvalues, fm_struct_tmp, mo_coeff, vecbuffer, mo_virt)
237 NULLIFY (curve_section, ldos_section, list, cell, pw_env, auxbas_pw_pool, evals_virt)
238 NULLIFY (occupation_numbers, ldos_p, r_ldos_p, dft_control, occupation_numbers)
239
240 CALL timeset(routinen, handle)
241 iterstep = logger%iter_info%iteration(logger%iter_info%n_rlevel)
242
243 IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T3,A,T61,I10))') &
244 " Calculate PDOS at iteration step ", iterstep
245 CALL get_qs_env(qs_env=qs_env, &
246 matrix_s=s_matrix, &
247 dft_control=dft_control)
248
249 CALL get_atomic_kind_set(atomic_kind_set, natom=natom)
250 CALL get_qs_kind_set(qs_kind_set, nsgf=nsgf, maxlgto=maxlgto)
251 nkind = SIZE(atomic_kind_set)
252
253 CALL get_mo_set(mo_set=mo_set, mo_coeff=mo_coeff, homo=homo, nao=nao, nmo=nmo, &
254 mu=e_fermi)
255 CALL cp_fm_get_info(mo_coeff, &
256 context=context, para_env=para_env, &
257 nrow_global=nrow_global, &
258 ncol_global=ncol_global)
259
260 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%OUT_EACH_MO", i_val=out_each)
261 IF (out_each == -1) out_each = nao + 1
262 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%DELTA_E", r_val=de)
263 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%BROADEN%TYPE", i_val=broaden_type)
264 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%BROADEN%WIDTH", r_val=broaden_width)
265 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%BROADEN%VOIGT_MIXING", r_val=voigt_mixing)
266 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%NDIGITS", i_val=ndigits)
267 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%ENERGY_UNIT", i_val=energy_unit)
268 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%ENERGY_ZERO", i_val=energy_zero)
269 ndigits = min(max(ndigits, 1), 10)
270 IF (write_curve .AND. de <= 0.0_dp) THEN
271 cpwarn("Broadened PDOS output requires DELTA_E > 0 and will be skipped")
272 write_curve = .false.
273 END IF
274 IF (write_curve .AND. broaden_width <= 0.0_dp) THEN
275 cpwarn("Broadened PDOS output requires a finite WIDTH and will be skipped")
276 write_curve = .false.
277 END IF
278 do_curve = write_curve .AND. (broaden_width > 0.0_dp)
279 IF (do_curve) de = max(de, 0.00001_dp)
280 nvirt = 0
281 NULLIFY (evals_virt)
282 IF (PRESENT(unoccupied_orbs) .AND. PRESENT(unoccupied_evals)) THEN
283 IF (ASSOCIATED(unoccupied_evals%array)) THEN
284 nvirt = SIZE(unoccupied_evals%array)
285 IF (nvirt > 0) THEN
286 mo_virt => unoccupied_orbs
287 evals_virt => unoccupied_evals%array
288 END IF
289 END IF
290 END IF
291 do_virt = (nvirt > 0)
292
293 calc_matsh = .true.
294 IF (PRESENT(external_matrix_shalf)) calc_matsh = .false.
295
296 ! Create S^1/2 : from sparse to full matrix, if no external available
297 IF (calc_matsh) THEN
298 NULLIFY (matrix_shalf)
299 CALL cp_fm_struct_create(fm_struct_tmp, para_env=para_env, context=context, &
300 nrow_global=nrow_global, ncol_global=nrow_global)
301 ALLOCATE (matrix_shalf)
302 CALL cp_fm_create(matrix_shalf, fm_struct_tmp, name="matrix_shalf")
303 CALL cp_fm_create(matrix_work, fm_struct_tmp, name="matrix_work")
304 CALL cp_fm_struct_release(fm_struct_tmp)
305 CALL copy_dbcsr_to_fm(s_matrix(1)%matrix, matrix_shalf)
306 CALL cp_fm_power(matrix_shalf, matrix_work, 0.5_dp, epsilon(0.0_dp), n_dependent)
307 CALL cp_fm_release(matrix_work)
308 ELSE
309 matrix_shalf => external_matrix_shalf
310 END IF
311
312 ! Multiply S^(1/2) time the mOS coefficients to get orthonormalized MOS
313 CALL cp_fm_struct_create(fm_struct_tmp, para_env=para_env, context=context, &
314 nrow_global=nrow_global, ncol_global=ncol_global)
315 CALL cp_fm_create(matrix_shalfc, fm_struct_tmp, name="matrix_shalfc")
316 CALL parallel_gemm("N", "N", nrow_global, ncol_global, nrow_global, &
317 1.0_dp, matrix_shalf, mo_coeff, 0.0_dp, matrix_shalfc)
318 CALL cp_fm_struct_release(fm_struct_tmp)
319
320 IF (do_virt) THEN
321 IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T3,A,T14,I10,T27,A))') &
322 " Use ", nvirt, " additional unoccupied KS orbitals"
323 CALL cp_fm_struct_create(fm_struct_tmp, para_env=para_env, context=context, &
324 nrow_global=nrow_global, ncol_global=nvirt)
325 CALL cp_fm_create(matrix_work, fm_struct_tmp, name="matrix_shalfc")
326 CALL parallel_gemm("N", "N", nrow_global, nvirt, nrow_global, &
327 1.0_dp, matrix_shalf, mo_virt, 0.0_dp, matrix_work)
328 CALL cp_fm_struct_release(fm_struct_tmp)
329 END IF
330
331 IF (calc_matsh) THEN
332 CALL cp_fm_release(matrix_shalf)
333 DEALLOCATE (matrix_shalf)
334 END IF
335 ! Array to store the PDOS per kind and angular momentum
336 do_ldos = .false.
337 ldos_section => section_vals_get_subs_vals(dft_section, trim(my_print_key)//"%LDOS")
338
339 CALL section_vals_get(ldos_section, n_repetition=nldos)
340 IF (nldos > 0) THEN
341 IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T3,A,T61,I10))') &
342 " Prepare the list of atoms for LDOS. Number of lists: ", nldos
343 do_ldos = .true.
344 ALLOCATE (ldos_p(nldos))
345 ALLOCATE (ldos_index(nldos))
346 DO ildos = 1, nldos
347 WRITE (ldos_index(ildos), '(I0)') ildos
348 ALLOCATE (ldos_p(ildos)%ldos)
349 NULLIFY (ldos_p(ildos)%ldos%pdos_array)
350 NULLIFY (ldos_p(ildos)%ldos%list_index)
351
352 CALL section_vals_val_get(ldos_section, "LIST", i_rep_section=ildos, n_rep_val=n_rep)
353 IF (n_rep > 0) THEN
354 ldos_p(ildos)%ldos%nlist = 0
355 DO ir = 1, n_rep
356 NULLIFY (list)
357 CALL section_vals_val_get(ldos_section, "LIST", i_rep_section=ildos, i_rep_val=ir, &
358 i_vals=list)
359 IF (ASSOCIATED(list)) THEN
360 CALL reallocate(ldos_p(ildos)%ldos%list_index, 1, ldos_p(ildos)%ldos%nlist + SIZE(list))
361 DO i = 1, SIZE(list)
362 ldos_p(ildos)%ldos%list_index(i + ldos_p(ildos)%ldos%nlist) = list(i)
363 END DO
364 ldos_p(ildos)%ldos%nlist = ldos_p(ildos)%ldos%nlist + SIZE(list)
365 END IF
366 END DO
367 ELSE
368 ! stop, LDOS without list of atoms is not implemented
369 END IF
370
371 IF (output_unit > 0) WRITE (unit=output_unit, fmt='((T10,A,T18,I6,T25,A,T36,I10,A))') &
372 " List ", ildos, " contains ", ldos_p(ildos)%ldos%nlist, " atoms"
373 CALL section_vals_val_get(ldos_section, "COMPONENTS", i_rep_section=ildos, &
374 l_val=ldos_p(ildos)%ldos%separate_components)
375 IF (ldos_p(ildos)%ldos%separate_components) THEN
376 ALLOCATE (ldos_p(ildos)%ldos%pdos_array(nsoset(maxlgto), nmo + nvirt))
377 ELSE
378 ALLOCATE (ldos_p(ildos)%ldos%pdos_array(0:maxlgto, nmo + nvirt))
379 END IF
380 ldos_p(ildos)%ldos%pdos_array = 0.0_dp
381 ldos_p(ildos)%ldos%maxl = -1
382
383 END DO
384 END IF
385
386 do_r_ldos = .false.
387 ldos_section => section_vals_get_subs_vals(dft_section, trim(my_print_key)//"%R_LDOS")
388 CALL section_vals_get(ldos_section, n_repetition=n_r_ldos)
389 IF (n_r_ldos > 0) THEN
390 do_r_ldos = .true.
391 IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T3,A,T61,I10))') &
392 " Prepare the list of points for R_LDOS. Number of lists: ", n_r_ldos
393 ALLOCATE (r_ldos_p(n_r_ldos))
394 ALLOCATE (r_ldos_index(n_r_ldos))
395 CALL get_qs_env(qs_env=qs_env, &
396 cell=cell, &
397 dft_control=dft_control, &
398 pw_env=pw_env)
399 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool, &
400 pw_pools=pw_pools)
401
402 CALL auxbas_pw_pool%create_pw(wf_r)
403 CALL auxbas_pw_pool%create_pw(wf_g)
404 ALLOCATE (read_r(4, n_r_ldos))
405 DO ildos = 1, n_r_ldos
406 WRITE (r_ldos_index(ildos), '(I0)') ildos
407 ALLOCATE (r_ldos_p(ildos)%ldos)
408 NULLIFY (r_ldos_p(ildos)%ldos%pdos_array)
409 NULLIFY (r_ldos_p(ildos)%ldos%list_index)
410
411 CALL section_vals_val_get(ldos_section, "LIST", i_rep_section=ildos, n_rep_val=n_rep)
412 IF (n_rep > 0) THEN
413 r_ldos_p(ildos)%ldos%nlist = 0
414 DO ir = 1, n_rep
415 NULLIFY (list)
416 CALL section_vals_val_get(ldos_section, "LIST", i_rep_section=ildos, i_rep_val=ir, &
417 i_vals=list)
418 IF (ASSOCIATED(list)) THEN
419 CALL reallocate(r_ldos_p(ildos)%ldos%list_index, 1, r_ldos_p(ildos)%ldos%nlist + SIZE(list))
420 DO i = 1, SIZE(list)
421 r_ldos_p(ildos)%ldos%list_index(i + r_ldos_p(ildos)%ldos%nlist) = list(i)
422 END DO
423 r_ldos_p(ildos)%ldos%nlist = r_ldos_p(ildos)%ldos%nlist + SIZE(list)
424 END IF
425 END DO
426 ELSE
427 ! stop, LDOS without list of atoms is not implemented
428 END IF
429
430 ALLOCATE (r_ldos_p(ildos)%ldos%pdos_array(nmo + nvirt))
431 r_ldos_p(ildos)%ldos%pdos_array = 0.0_dp
432 read_r(1:3, ildos) = .false.
433 CALL section_vals_val_get(ldos_section, "XRANGE", i_rep_section=ildos, explicit=read_r(1, ildos))
434 IF (read_r(1, ildos)) THEN
435 CALL section_vals_val_get(ldos_section, "XRANGE", i_rep_section=ildos, r_vals= &
436 r_ldos_p(ildos)%ldos%x_range)
437 ELSE
438 ALLOCATE (r_ldos_p(ildos)%ldos%x_range(2))
439 r_ldos_p(ildos)%ldos%x_range = 0.0_dp
440 END IF
441 CALL section_vals_val_get(ldos_section, "YRANGE", i_rep_section=ildos, explicit=read_r(2, ildos))
442 IF (read_r(2, ildos)) THEN
443 CALL section_vals_val_get(ldos_section, "YRANGE", i_rep_section=ildos, r_vals= &
444 r_ldos_p(ildos)%ldos%y_range)
445 ELSE
446 ALLOCATE (r_ldos_p(ildos)%ldos%y_range(2))
447 r_ldos_p(ildos)%ldos%y_range = 0.0_dp
448 END IF
449 CALL section_vals_val_get(ldos_section, "ZRANGE", i_rep_section=ildos, explicit=read_r(3, ildos))
450 IF (read_r(3, ildos)) THEN
451 CALL section_vals_val_get(ldos_section, "ZRANGE", i_rep_section=ildos, r_vals= &
452 r_ldos_p(ildos)%ldos%z_range)
453 ELSE
454 ALLOCATE (r_ldos_p(ildos)%ldos%z_range(2))
455 r_ldos_p(ildos)%ldos%z_range = 0.0_dp
456 END IF
457
458 CALL section_vals_val_get(ldos_section, "ERANGE", i_rep_section=ildos, explicit=read_r(4, ildos))
459 IF (read_r(4, ildos)) THEN
460 CALL section_vals_val_get(ldos_section, "ERANGE", i_rep_section=ildos, &
461 r_vals=r_ldos_p(ildos)%ldos%eval_range)
462 ELSE
463 ALLOCATE (r_ldos_p(ildos)%ldos%eval_range(2))
464 r_ldos_p(ildos)%ldos%eval_range(1) = -huge(0.0_dp)
465 r_ldos_p(ildos)%ldos%eval_range(2) = +huge(0.0_dp)
466 END IF
467
468 bo => wf_r%pw_grid%bounds_local
469 dh = wf_r%pw_grid%dh
470 dvol = wf_r%pw_grid%dvol
471 np_tot = wf_r%pw_grid%npts(1)*wf_r%pw_grid%npts(2)*wf_r%pw_grid%npts(3)
472 ALLOCATE (r_ldos_p(ildos)%ldos%index_grid_local(3, np_tot))
473
474 r_ldos_p(ildos)%ldos%npoints = 0
475 DO jz = bo(1, 3), bo(2, 3)
476 DO jy = bo(1, 2), bo(2, 2)
477 DO jx = bo(1, 1), bo(2, 1)
478 !compute the position of the grid point
479 i = jx - wf_r%pw_grid%bounds(1, 1)
480 j = jy - wf_r%pw_grid%bounds(1, 2)
481 k = jz - wf_r%pw_grid%bounds(1, 3)
482 r(3) = k*dh(3, 3) + j*dh(3, 2) + i*dh(3, 1)
483 r(2) = k*dh(2, 3) + j*dh(2, 2) + i*dh(2, 1)
484 r(1) = k*dh(1, 3) + j*dh(1, 2) + i*dh(1, 1)
485
486 DO il = 1, r_ldos_p(ildos)%ldos%nlist
487 iatom = r_ldos_p(ildos)%ldos%list_index(il)
488 ratom = particle_set(iatom)%r
489 r_vec = pbc(ratom, r, cell)
490 IF (cell%orthorhombic) THEN
491 IF (cell%perd(1) == 0) r_vec(1) = modulo(r_vec(1), cell%hmat(1, 1))
492 IF (cell%perd(2) == 0) r_vec(2) = modulo(r_vec(2), cell%hmat(2, 2))
493 IF (cell%perd(3) == 0) r_vec(3) = modulo(r_vec(3), cell%hmat(3, 3))
494 END IF
495
496 in_x = 0
497 in_y = 0
498 in_z = 0
499 IF (r_ldos_p(ildos)%ldos%x_range(1) /= 0.0_dp) THEN
500 IF (r_vec(1) > r_ldos_p(ildos)%ldos%x_range(1) .AND. &
501 r_vec(1) < r_ldos_p(ildos)%ldos%x_range(2)) THEN
502 in_x = 1
503 END IF
504 ELSE
505 in_x = 1
506 END IF
507 IF (r_ldos_p(ildos)%ldos%y_range(1) /= 0.0_dp) THEN
508 IF (r_vec(2) > r_ldos_p(ildos)%ldos%y_range(1) .AND. &
509 r_vec(2) < r_ldos_p(ildos)%ldos%y_range(2)) THEN
510 in_y = 1
511 END IF
512 ELSE
513 in_y = 1
514 END IF
515 IF (r_ldos_p(ildos)%ldos%z_range(1) /= 0.0_dp) THEN
516 IF (r_vec(3) > r_ldos_p(ildos)%ldos%z_range(1) .AND. &
517 r_vec(3) < r_ldos_p(ildos)%ldos%z_range(2)) THEN
518 in_z = 1
519 END IF
520 ELSE
521 in_z = 1
522 END IF
523 IF (in_x*in_y*in_z > 0) THEN
524 r_ldos_p(ildos)%ldos%npoints = r_ldos_p(ildos)%ldos%npoints + 1
525 r_ldos_p(ildos)%ldos%index_grid_local(1, r_ldos_p(ildos)%ldos%npoints) = jx
526 r_ldos_p(ildos)%ldos%index_grid_local(2, r_ldos_p(ildos)%ldos%npoints) = jy
527 r_ldos_p(ildos)%ldos%index_grid_local(3, r_ldos_p(ildos)%ldos%npoints) = jz
528 EXIT
529 END IF
530 END DO
531 END DO
532 END DO
533 END DO
534 CALL reallocate(r_ldos_p(ildos)%ldos%index_grid_local, 1, 3, 1, r_ldos_p(ildos)%ldos%npoints)
535 npoints = r_ldos_p(ildos)%ldos%npoints
536 CALL para_env%sum(npoints)
537 IF (output_unit > 0) WRITE (unit=output_unit, fmt='((T10,A,T18,I6,T25,A,T36,I10,A))') &
538 " List ", ildos, " contains ", npoints, " grid points"
539 END DO
540 END IF
541
542 IF (trim(my_print_key) == "PRINT%DOS") THEN
543 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%PDOS%COMPONENTS", l_val=separate_components)
544 ELSE
545 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%COMPONENTS", l_val=separate_components)
546 END IF
547 IF (separate_components) THEN
548 ALLOCATE (pdos_array(nsoset(maxlgto), nkind, nmo + nvirt))
549 ELSE
550 ALLOCATE (pdos_array(0:maxlgto, nkind, nmo + nvirt))
551 END IF
552 IF (do_virt) THEN
553 ALLOCATE (eigenvalues(nmo + nvirt))
554 eigenvalues(1:nmo) = mo_set%eigenvalues(1:nmo)
555 eigenvalues(nmo + 1:nmo + nvirt) = evals_virt(1:nvirt)
556 ALLOCATE (occupation_numbers(nmo + nvirt))
557 occupation_numbers(:) = 0.0_dp
558 occupation_numbers(1:nmo) = mo_set%occupation_numbers(1:nmo)
559 ELSE
560 eigenvalues => mo_set%eigenvalues
561 occupation_numbers => mo_set%occupation_numbers
562 END IF
563
564 hoco = -huge(0.0_dp)
565 fractional_occupation = .false.
566 DO imo = 1, nmo + nvirt
567 IF (occupation_numbers(imo) > 1.0e-10_dp) hoco = max(hoco, eigenvalues(imo))
568 IF (abs(occupation_numbers(imo) - real(nint(occupation_numbers(imo)), kind=dp)) > &
569 1.0e-8_dp) fractional_occupation = .true.
570 END DO
571 IF (hoco < -0.5_dp*huge(0.0_dp)) hoco = e_fermi
572 IF (PRESENT(ispin) .AND. dft_control%nspins == 2) THEN
573 e_fermi_ref(:) = 0.0_dp
574 hoco_ref(:) = -huge(0.0_dp)
575 CALL get_qs_env(qs_env=qs_env, mos=mos_ref)
576 IF (ASSOCIATED(mos_ref)) THEN
577 DO ispin_ref = 1, dft_control%nspins
578 CALL get_mo_set(mo_set=mos_ref(ispin_ref), nmo=nmo_ref, mu=e_fermi_ref(ispin_ref))
579 eval_ref => mos_ref(ispin_ref)%eigenvalues
580 occ_ref => mos_ref(ispin_ref)%occupation_numbers
581 DO imo_ref = 1, nmo_ref
582 IF (occ_ref(imo_ref) > 1.0e-10_dp) THEN
583 hoco_ref(ispin_ref) = max(hoco_ref(ispin_ref), eval_ref(imo_ref))
584 END IF
585 IF (abs(occ_ref(imo_ref) - real(nint(occ_ref(imo_ref)), kind=dp)) > &
586 1.0e-8_dp) fractional_occupation = .true.
587 END DO
588 IF (hoco_ref(ispin_ref) < -0.5_dp*huge(0.0_dp)) hoco_ref(ispin_ref) = e_fermi_ref(ispin_ref)
589 END DO
590 ELSE
591 e_fermi_ref(:) = e_fermi
592 hoco_ref(:) = hoco
593 END IF
594 ELSE
595 e_fermi_ref(:) = e_fermi
596 hoco_ref(:) = hoco
597 END IF
598 resolved_energy_zero = dos_resolve_energy_zero(energy_zero, dft_control%smear, fractional_occupation)
599 SELECT CASE (resolved_energy_zero)
600 CASE (dos_energy_zero_absolute)
601 energy_ref = 0.0_dp
602 CASE (dos_energy_zero_hoco)
603 energy_ref = maxval(hoco_ref(1:dft_control%nspins))
604 CASE DEFAULT
605 energy_ref = maxval(e_fermi_ref(1:dft_control%nspins))
606 END SELECT
607 energy_factor = dos_energy_scale(energy_unit)
608 energy_label = dos_energy_label(energy_unit)
609 zero_label = dos_energy_zero_label(resolved_energy_zero)
610 IF (energy_zero == dos_energy_zero_auto) zero_label = "AUTO -> "//trim(zero_label)
611 ev_factor = dos_energy_scale(dos_energy_unit_ev)
612
613 pdos_array = 0.0_dp
614 nao = mo_set%nao
615 ALLOCATE (vecbuffer(1, nao))
616 vecbuffer = 0.0_dp
617 ALLOCATE (firstrow(natom))
618 firstrow = 0
619
620 !Adjust energy range for r_ldos
621 DO ildos = 1, n_r_ldos
622 IF (eigenvalues(1) > r_ldos_p(ildos)%ldos%eval_range(1)) THEN
623 r_ldos_p(ildos)%ldos%eval_range(1) = eigenvalues(1)
624 END IF
625 IF (eigenvalues(nmo + nvirt) < r_ldos_p(ildos)%ldos%eval_range(2)) THEN
626 r_ldos_p(ildos)%ldos%eval_range(2) = eigenvalues(nmo + nvirt)
627 END IF
628 END DO
629
630 IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T15,A))') &
631 "---- PDOS: start iteration on the KS states --- "
632
633 DO imo = 1, nmo + nvirt
634
635 IF (output_unit > 0 .AND. mod(imo, out_each) == 0) WRITE (unit=output_unit, fmt='((T20,A,I10))') &
636 " KS state index : ", imo
637 ! Extract the eigenvector from the distributed full matrix
638 IF (imo > nmo) THEN
639 CALL cp_fm_get_submatrix(matrix_work, vecbuffer, 1, imo - nmo, &
640 nao, 1, transpose=.true.)
641 ELSE
642 CALL cp_fm_get_submatrix(matrix_shalfc, vecbuffer, 1, imo, &
643 nao, 1, transpose=.true.)
644 END IF
645
646 ! Calculate the pdos for all the kinds
647 irow = 1
648 DO iatom = 1, natom
649 firstrow(iatom) = irow
650 NULLIFY (orb_basis_set)
651 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
652 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
653 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
654 nset=nset, &
655 nshell=nshell, &
656 l=l, maxl=maxl)
657 IF (separate_components) THEN
658 isgf = 1
659 DO iset = 1, nset
660 DO ishell = 1, nshell(iset)
661 lshell = l(ishell, iset)
662 DO iso = 1, nso(lshell)
663 lcomponent = nsoset(lshell - 1) + iso
664 pdos_array(lcomponent, ikind, imo) = &
665 pdos_array(lcomponent, ikind, imo) + &
666 vecbuffer(1, irow)*vecbuffer(1, irow)
667 irow = irow + 1
668 END DO ! iso
669 END DO ! ishell
670 END DO ! iset
671 ELSE
672 isgf = 1
673 DO iset = 1, nset
674 DO ishell = 1, nshell(iset)
675 lshell = l(ishell, iset)
676 DO iso = 1, nso(lshell)
677 pdos_array(lshell, ikind, imo) = &
678 pdos_array(lshell, ikind, imo) + &
679 vecbuffer(1, irow)*vecbuffer(1, irow)
680 irow = irow + 1
681 END DO ! iso
682 END DO ! ishell
683 END DO ! iset
684 END IF
685 END DO ! iatom
686
687 ! Calculate the pdos for all the lists
688 DO ildos = 1, nldos
689 DO il = 1, ldos_p(ildos)%ldos%nlist
690 iatom = ldos_p(ildos)%ldos%list_index(il)
691
692 irow = firstrow(iatom)
693 NULLIFY (orb_basis_set)
694 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
695 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
696
697 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
698 nset=nset, &
699 nshell=nshell, &
700 l=l, maxl=maxl)
701 ldos_p(ildos)%ldos%maxl = max(ldos_p(ildos)%ldos%maxl, maxl)
702 IF (ldos_p(ildos)%ldos%separate_components) THEN
703 isgf = 1
704 DO iset = 1, nset
705 DO ishell = 1, nshell(iset)
706 lshell = l(ishell, iset)
707 DO iso = 1, nso(lshell)
708 lcomponent = nsoset(lshell - 1) + iso
709 ldos_p(ildos)%ldos%pdos_array(lcomponent, imo) = &
710 ldos_p(ildos)%ldos%pdos_array(lcomponent, imo) + &
711 vecbuffer(1, irow)*vecbuffer(1, irow)
712 irow = irow + 1
713 END DO ! iso
714 END DO ! ishell
715 END DO ! iset
716 ELSE
717 isgf = 1
718 DO iset = 1, nset
719 DO ishell = 1, nshell(iset)
720 lshell = l(ishell, iset)
721 DO iso = 1, nso(lshell)
722 ldos_p(ildos)%ldos%pdos_array(lshell, imo) = &
723 ldos_p(ildos)%ldos%pdos_array(lshell, imo) + &
724 vecbuffer(1, irow)*vecbuffer(1, irow)
725 irow = irow + 1
726 END DO ! iso
727 END DO ! ishell
728 END DO ! iset
729 END IF
730 END DO !il
731 END DO !ildos
732
733 ! Calculate the DOS projected in a given volume in real space
734 DO ildos = 1, n_r_ldos
735 IF (r_ldos_p(ildos)%ldos%eval_range(1) <= eigenvalues(imo) .AND. &
736 r_ldos_p(ildos)%ldos%eval_range(2) >= eigenvalues(imo)) THEN
737
738 IF (imo > nmo) THEN
739 CALL calculate_wavefunction(mo_virt, imo - nmo, &
740 wf_r, wf_g, atomic_kind_set, qs_kind_set, cell, dft_control, particle_set, &
741 pw_env)
742 ELSE
743 CALL calculate_wavefunction(mo_coeff, imo, &
744 wf_r, wf_g, atomic_kind_set, qs_kind_set, cell, dft_control, particle_set, &
745 pw_env)
746 END IF
747 r_ldos_p(ildos)%ldos%pdos_array(imo) = 0.0_dp
748 DO il = 1, r_ldos_p(ildos)%ldos%npoints
749 j = j + 1
750 jx = r_ldos_p(ildos)%ldos%index_grid_local(1, il)
751 jy = r_ldos_p(ildos)%ldos%index_grid_local(2, il)
752 jz = r_ldos_p(ildos)%ldos%index_grid_local(3, il)
753 r_ldos_p(ildos)%ldos%pdos_array(imo) = r_ldos_p(ildos)%ldos%pdos_array(imo) + &
754 wf_r%array(jx, jy, jz)*wf_r%array(jx, jy, jz)
755 END DO
756 r_ldos_p(ildos)%ldos%pdos_array(imo) = r_ldos_p(ildos)%ldos%pdos_array(imo)*dvol
757 END IF
758 END DO
759 END DO ! imo
760
761 CALL cp_fm_release(matrix_shalfc)
762 DEALLOCATE (vecbuffer)
763
764 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%APPEND", l_val=append)
765 IF (append .AND. iterstep > 1) THEN
766 my_pos = "APPEND"
767 ELSE
768 my_pos = "REWIND"
769 END IF
770 my_act = "WRITE"
771 IF (write_pdos_file) THEN
772 DO ikind = 1, nkind
773
774 NULLIFY (orb_basis_set)
775 CALL get_atomic_kind(atomic_kind_set(ikind), name=kind_name)
776 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
777 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, maxl=maxl)
778
779 ! basis none has no associated maxl, and no pdos
780 IF (maxl < 0) cycle
781
782 IF (PRESENT(ispin)) THEN
783 IF (PRESENT(xas_mittle)) THEN
784 my_mittle = trim(xas_mittle)//trim(spin(ispin))//"_k"//trim(adjustl(cp_to_string(ikind)))
785 ELSE
786 my_mittle = trim(spin(ispin))//"_k"//trim(adjustl(cp_to_string(ikind)))
787 END IF
788 my_spin = ispin
789 ELSE
790 my_mittle = "k"//trim(adjustl(cp_to_string(ikind)))
791 my_spin = 1
792 END IF
793
794 IF (write_pdos_file .AND. do_curve) THEN
795 IF (separate_components) THEN
796 CALL write_broadened_pdos(logger, dft_section, trim(my_mittle)//"_curve", my_pos, my_act, iterstep, &
797 e_fermi, hoco, energy_ref, trim(zero_label), &
798 "Projected DOS for atomic kind "//trim(kind_name), &
799 maxl, .true., pdos_array(1:nsoset(maxl), ikind, :), &
800 eigenvalues, nmo + nvirt, de, broaden_type, broaden_width, &
801 voigt_mixing, ndigits, pdos_print_key=trim(my_print_key))
802 ELSE
803 CALL write_broadened_pdos(logger, dft_section, trim(my_mittle)//"_curve", my_pos, my_act, iterstep, &
804 e_fermi, hoco, energy_ref, trim(zero_label), &
805 "Projected DOS for atomic kind "//trim(kind_name), &
806 maxl, .false., pdos_array(0:maxl, ikind, :), &
807 eigenvalues, nmo + nvirt, de, broaden_type, broaden_width, &
808 voigt_mixing, ndigits, pdos_print_key=trim(my_print_key))
809 END IF
810 END IF
811
812 IF (write_pdos_file) THEN
813 iw = cp_print_key_unit_nr(logger, dft_section, trim(my_print_key), &
814 extension=".pdos", file_position=my_pos, file_action=my_act, &
815 file_form="FORMATTED", middle_name=trim(my_mittle))
816 IF (iw > 0) THEN
817
818 fmtstr1 = "(I8,2X,2F16.6, (2X,F16.8))"
819 fmtstr2 = "(A42, (10X,A8))"
820 IF (separate_components) THEN
821 WRITE (unit=fmtstr1(15:16), fmt="(I2)") nsoset(maxl) + 1
822 WRITE (unit=fmtstr2(6:7), fmt="(I2)") nsoset(maxl) + 1
823 ELSE
824 WRITE (unit=fmtstr1(15:16), fmt="(I2)") maxl + 2
825 WRITE (unit=fmtstr2(6:7), fmt="(I2)") maxl + 2
826 END IF
827
828 WRITE (unit=iw, fmt="(A,I0)") &
829 "# Projected DOS for atomic kind "//trim(kind_name)//" at iteration step i = ", iterstep
830 WRITE (unit=iw, fmt="(A,F12.6,A,F12.6,A)") &
831 "# E(Fermi) = ", e_fermi, " a.u. = ", e_fermi*ev_factor, " eV"
832 WRITE (unit=iw, fmt="(A,F12.6,A,F12.6,A)") &
833 "# E(HOCO) = ", hoco, " a.u. = ", hoco*ev_factor, " eV"
834 IF (separate_components) THEN
835 ALLOCATE (tmp_str(0:0, 0:maxl, -maxl:maxl))
836 tmp_str = ""
837 DO j = 0, maxl
838 DO i = -j, j
839 tmp_str(0, j, i) = sgf_symbol(0, j, i)
840 END DO
841 END DO
842
843 WRITE (unit=iw, fmt=fmtstr2) &
844 "# MO Energy[a.u.] Occupation", "Total", &
845 ((trim(tmp_str(0, il, im)), im=-il, il), il=0, maxl)
846 DO imo = 1, nmo + nvirt
847 WRITE (unit=iw, fmt=fmtstr1) imo, eigenvalues(imo), &
848 occupation_numbers(imo), sum(pdos_array(1:nsoset(maxl), ikind, imo)), &
849 (pdos_array(lshell, ikind, imo), lshell=1, nsoset(maxl))
850 END DO
851 DEALLOCATE (tmp_str)
852 ELSE
853 WRITE (unit=iw, fmt=fmtstr2) &
854 "# MO Energy[a.u.] Occupation", "Total", &
855 (trim(l_sym(il)), il=0, maxl)
856 DO imo = 1, nmo + nvirt
857 WRITE (unit=iw, fmt=fmtstr1) imo, eigenvalues(imo), &
858 occupation_numbers(imo), sum(pdos_array(0:maxl, ikind, imo)), &
859 (pdos_array(lshell, ikind, imo), lshell=0, maxl)
860 END DO
861 END IF
862 END IF
863 CALL cp_print_key_finished_output(iw, logger, dft_section, &
864 trim(my_print_key))
865 END IF
866
867 END DO ! ikind
868 END IF
869
870 ! write the pdos for the lists, each ona different file,
871 ! the filenames are indexed with the list number
872 DO ildos = 1, nldos
873 ! basis none has no associated maxl, and no pdos
874 IF (ldos_p(ildos)%ldos%maxl > 0) THEN
875
876 IF (PRESENT(ispin)) THEN
877 IF (PRESENT(xas_mittle)) THEN
878 my_mittle = trim(xas_mittle)//trim(spin(ispin))//"_list"//trim(ldos_index(ildos))
879 ELSE
880 my_mittle = trim(spin(ispin))//"_list"//trim(ldos_index(ildos))
881 END IF
882 my_spin = ispin
883 ELSE
884 my_mittle = "list"//trim(ldos_index(ildos))
885 my_spin = 1
886 END IF
887
888 IF (do_curve) THEN
889 IF (ldos_p(ildos)%ldos%separate_components) THEN
890 CALL write_broadened_pdos(logger, dft_section, trim(my_mittle)//"_curve", my_pos, my_act, iterstep, &
891 e_fermi, hoco, energy_ref, trim(zero_label), &
892 "Projected DOS for atom list "//trim(ldos_index(ildos)), &
893 ldos_p(ildos)%ldos%maxl, .true., &
894 ldos_p(ildos)%ldos%pdos_array(1:nsoset(ldos_p(ildos)%ldos%maxl), :), &
895 eigenvalues, nmo + nvirt, de, broaden_type, broaden_width, &
896 voigt_mixing, ndigits, pdos_print_key=trim(my_print_key))
897 ELSE
898 CALL write_broadened_pdos(logger, dft_section, trim(my_mittle)//"_curve", my_pos, my_act, iterstep, &
899 e_fermi, hoco, energy_ref, trim(zero_label), &
900 "Projected DOS for atom list "//trim(ldos_index(ildos)), &
901 ldos_p(ildos)%ldos%maxl, .false., &
902 ldos_p(ildos)%ldos%pdos_array(0:ldos_p(ildos)%ldos%maxl, :), &
903 eigenvalues, nmo + nvirt, de, broaden_type, broaden_width, &
904 voigt_mixing, ndigits, pdos_print_key=trim(my_print_key))
905 END IF
906 END IF
907
908 iw = cp_print_key_unit_nr(logger, dft_section, trim(my_print_key), &
909 extension=".pdos", file_position=my_pos, file_action=my_act, &
910 file_form="FORMATTED", middle_name=trim(my_mittle))
911 IF (iw > 0) THEN
912
913 fmtstr1 = "(I8,2X,2F16.6, (2X,F16.8))"
914 fmtstr2 = "(A42, (10X,A8))"
915 IF (ldos_p(ildos)%ldos%separate_components) THEN
916 WRITE (unit=fmtstr1(15:16), fmt="(I2)") nsoset(ldos_p(ildos)%ldos%maxl) + 1
917 WRITE (unit=fmtstr2(6:7), fmt="(I2)") nsoset(ldos_p(ildos)%ldos%maxl) + 1
918 ELSE
919 WRITE (unit=fmtstr1(15:16), fmt="(I2)") ldos_p(ildos)%ldos%maxl + 2
920 WRITE (unit=fmtstr2(6:7), fmt="(I2)") ldos_p(ildos)%ldos%maxl + 2
921 END IF
922
923 WRITE (unit=iw, fmt="(A,I0,A,I0,A,I0)") &
924 "# Projected DOS for list ", ildos, " of ", ldos_p(ildos)%ldos%nlist, &
925 " atoms, at iteration step i = ", iterstep
926 WRITE (unit=iw, fmt="(A,F12.6,A,F12.6,A)") &
927 "# E(Fermi) = ", e_fermi, " a.u. = ", e_fermi*ev_factor, " eV"
928 WRITE (unit=iw, fmt="(A,F12.6,A,F12.6,A)") &
929 "# E(HOCO) = ", hoco, " a.u. = ", hoco*ev_factor, " eV"
930 IF (ldos_p(ildos)%ldos%separate_components) THEN
931 ALLOCATE (tmp_str(0:0, 0:ldos_p(ildos)%ldos%maxl, -ldos_p(ildos)%ldos%maxl:ldos_p(ildos)%ldos%maxl))
932 tmp_str = ""
933 DO j = 0, ldos_p(ildos)%ldos%maxl
934 DO i = -j, j
935 tmp_str(0, j, i) = sgf_symbol(0, j, i)
936 END DO
937 END DO
938
939 WRITE (unit=iw, fmt=fmtstr2) &
940 "# MO Energy[a.u.] Occupation", "Total", &
941 ((trim(tmp_str(0, il, im)), im=-il, il), il=0, ldos_p(ildos)%ldos%maxl)
942 DO imo = 1, nmo + nvirt
943 WRITE (unit=iw, fmt=fmtstr1) imo, eigenvalues(imo), &
944 occupation_numbers(imo), &
945 sum(ldos_p(ildos)%ldos%pdos_array(1:nsoset(ldos_p(ildos)%ldos%maxl), imo)), &
946 (ldos_p(ildos)%ldos%pdos_array(lshell, imo), &
947 lshell=1, nsoset(ldos_p(ildos)%ldos%maxl))
948 END DO
949 DEALLOCATE (tmp_str)
950 ELSE
951 WRITE (unit=iw, fmt=fmtstr2) &
952 "# MO Energy[a.u.] Occupation", "Total", &
953 (trim(l_sym(il)), il=0, ldos_p(ildos)%ldos%maxl)
954 DO imo = 1, nmo + nvirt
955 WRITE (unit=iw, fmt=fmtstr1) imo, eigenvalues(imo), &
956 occupation_numbers(imo), &
957 sum(ldos_p(ildos)%ldos%pdos_array(0:ldos_p(ildos)%ldos%maxl, imo)), &
958 (ldos_p(ildos)%ldos%pdos_array(lshell, imo), lshell=0, ldos_p(ildos)%ldos%maxl)
959 END DO
960 END IF
961 END IF
962 CALL cp_print_key_finished_output(iw, logger, dft_section, &
963 trim(my_print_key))
964 END IF ! maxl>0
965 END DO ! ildos
966
967 ! write the pdos for the lists, each ona different file,
968 ! the filenames are indexed with the list number
969 DO ildos = 1, n_r_ldos
970
971 npoints = r_ldos_p(ildos)%ldos%npoints
972 CALL para_env%sum(npoints)
973 CALL para_env%sum(np_tot)
974 CALL para_env%sum(r_ldos_p(ildos)%ldos%pdos_array)
975 IF (PRESENT(ispin)) THEN
976 IF (PRESENT(xas_mittle)) THEN
977 my_mittle = trim(xas_mittle)//trim(spin(ispin))//"_r_list"//trim(r_ldos_index(ildos))
978 ELSE
979 my_mittle = trim(spin(ispin))//"_r_list"//trim(r_ldos_index(ildos))
980 END IF
981 my_spin = ispin
982 ELSE
983 my_mittle = "r_list"//trim(r_ldos_index(ildos))
984 my_spin = 1
985 END IF
986
987 iw = cp_print_key_unit_nr(logger, dft_section, trim(my_print_key), &
988 extension=".pdos", file_position=my_pos, file_action=my_act, &
989 file_form="FORMATTED", middle_name=trim(my_mittle))
990 IF (iw > 0) THEN
991 fmtstr1 = "(I8,2X,2F16.6, (2X,F16.8))"
992 fmtstr2 = "(A42, (10X,A8))"
993
994 WRITE (unit=iw, fmt="(A,I0,A,F12.6,F12.6,A)") &
995 "# Projected DOS in real space, using ", npoints, &
996 " points of the grid, and eval in the range", r_ldos_p(ildos)%ldos%eval_range(1:2), " Hartree"
997 WRITE (unit=iw, fmt="(A,F12.6,A,F12.6,A)") &
998 "# E(Fermi) = ", e_fermi, " a.u. = ", e_fermi*ev_factor, " eV"
999 WRITE (unit=iw, fmt="(A,F12.6,A,F12.6,A)") &
1000 "# E(HOCO) = ", hoco, " a.u. = ", hoco*ev_factor, " eV"
1001 WRITE (unit=iw, fmt="(A)") &
1002 "# MO Energy[a.u.] Occupation LDOS"
1003 DO imo = 1, nmo + nvirt
1004 IF (r_ldos_p(ildos)%ldos%eval_range(1) <= eigenvalues(imo) .AND. &
1005 r_ldos_p(ildos)%ldos%eval_range(2) >= eigenvalues(imo)) THEN
1006 WRITE (unit=iw, fmt="(I8,2X,2F16.6,E20.10,E20.10)") imo, &
1007 eigenvalues(imo), occupation_numbers(imo), &
1008 r_ldos_p(ildos)%ldos%pdos_array(imo), r_ldos_p(ildos)%ldos%pdos_array(imo)*np_tot
1009 END IF
1010 END DO
1011
1012 END IF
1013 CALL cp_print_key_finished_output(iw, logger, dft_section, &
1014 trim(my_print_key))
1015 END DO
1016
1017 ! deallocate local variables
1018 DEALLOCATE (pdos_array)
1019 DEALLOCATE (firstrow)
1020 IF (do_ldos) THEN
1021 DO ildos = 1, nldos
1022 DEALLOCATE (ldos_p(ildos)%ldos%pdos_array)
1023 DEALLOCATE (ldos_p(ildos)%ldos%list_index)
1024 DEALLOCATE (ldos_p(ildos)%ldos)
1025 END DO
1026 DEALLOCATE (ldos_p)
1027 DEALLOCATE (ldos_index)
1028 END IF
1029 IF (do_r_ldos) THEN
1030 DO ildos = 1, n_r_ldos
1031 DEALLOCATE (r_ldos_p(ildos)%ldos%index_grid_local)
1032 DEALLOCATE (r_ldos_p(ildos)%ldos%pdos_array)
1033 DEALLOCATE (r_ldos_p(ildos)%ldos%list_index)
1034 IF (.NOT. read_r(1, ildos)) THEN
1035 DEALLOCATE (r_ldos_p(ildos)%ldos%x_range)
1036 END IF
1037 IF (.NOT. read_r(2, ildos)) THEN
1038 DEALLOCATE (r_ldos_p(ildos)%ldos%y_range)
1039 END IF
1040 IF (.NOT. read_r(3, ildos)) THEN
1041 DEALLOCATE (r_ldos_p(ildos)%ldos%z_range)
1042 END IF
1043 IF (.NOT. read_r(4, ildos)) THEN
1044 DEALLOCATE (r_ldos_p(ildos)%ldos%eval_range)
1045 END IF
1046 DEALLOCATE (r_ldos_p(ildos)%ldos)
1047 END DO
1048 DEALLOCATE (read_r)
1049 DEALLOCATE (r_ldos_p)
1050 DEALLOCATE (r_ldos_index)
1051 CALL auxbas_pw_pool%give_back_pw(wf_r)
1052 CALL auxbas_pw_pool%give_back_pw(wf_g)
1053 END IF
1054 IF (do_virt) THEN
1055 CALL cp_fm_release(matrix_work)
1056 DEALLOCATE (eigenvalues)
1057 DEALLOCATE (occupation_numbers)
1058 END IF
1059
1060 CALL timestop(handle)
1061
1062 END SUBROUTINE calculate_projected_dos
1063
1064! **************************************************************************************************
1065!> \brief Compute and write broadened projected density of states for k-point calculations.
1066!> \param qs_env ...
1067!> \param dft_section ...
1068!> \param pdos_print_key ...
1069!> \param write_pdos ...
1070!> \param write_pdos_curve ...
1071! **************************************************************************************************
1072 SUBROUTINE calculate_projected_dos_kp(qs_env, dft_section, pdos_print_key, write_pdos, write_pdos_curve)
1073
1074 TYPE(qs_environment_type), POINTER :: qs_env
1075 TYPE(section_vals_type), POINTER :: dft_section
1076 CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: pdos_print_key
1077 LOGICAL, INTENT(IN), OPTIONAL :: write_pdos, write_pdos_curve
1078
1079 CHARACTER(len=*), PARAMETER :: routinen = 'calculate_projected_dos_kp'
1080
1081 CHARACTER(LEN=32) :: zero_label
1082 CHARACTER(LEN=default_string_length) :: kind_name, my_act, my_mittle, my_pos, &
1083 my_print_key, spin(2)
1084 COMPLEX(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: zvecbuffer
1085 INTEGER :: broaden_type, energy_zero, fractional_occupation_int, handle, icomp, ik, ikind, &
1086 imo, ispin, iterstep, maxl, maxlgto, n_r_ldos, nao, ncomp, ndigits, nhist, nkind, nldos, &
1087 nmo_kp, nspins, output_unit, resolved_energy_zero
1088 INTEGER, ALLOCATABLE, DIMENSION(:) :: ao_comp, ao_kind, ao_l, kind_maxl
1089 LOGICAL :: append, fractional_occupation, &
1090 separate_components, should_output, &
1091 write_curve, write_pdos_file
1092 REAL(kind=dp) :: broaden_cutoff, broaden_width, de, e1, &
1093 e2, e_fermi(2), emax, emin, &
1094 energy_ref(2), hoco(2), voigt_mixing, &
1095 wkp
1096 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: ao_weight
1097 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: proj_weight, vecbuffer
1098 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :, :) :: pdos_curve
1099 REAL(kind=dp), DIMENSION(:), POINTER :: eigenvalues, occupation_numbers
1100 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
1101 TYPE(cp_cfm_type) :: cshalfc
1102 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
1103 TYPE(cp_fm_type) :: shalfc
1104 TYPE(cp_logger_type), POINTER :: logger
1105 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_s_kp
1106 TYPE(dft_control_type), POINTER :: dft_control
1107 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
1108 TYPE(kpoint_env_type), POINTER :: kp
1109 TYPE(kpoint_type), POINTER :: kpoints
1110 TYPE(mo_set_type), POINTER :: mo_set
1111 TYPE(mp_para_env_type), POINTER :: para_env
1112 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1113 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1114 TYPE(qs_scf_env_type), POINTER :: scf_env
1115 TYPE(section_vals_type), POINTER :: curve_section, ldos_section
1116
1117 NULLIFY (logger, kpoints, dft_control, para_env, atomic_kind_set, qs_kind_set, particle_set)
1118 NULLIFY (matrix_s_kp, scf_env)
1119 NULLIFY (kp, mo_set, eigenvalues, fm_struct_tmp, orb_basis_set, curve_section, ldos_section)
1120 logger => cp_get_default_logger()
1121 my_print_key = "PRINT%PDOS"
1122 IF (PRESENT(pdos_print_key)) my_print_key = trim(pdos_print_key)
1123 write_pdos_file = .true.
1124 IF (PRESENT(write_pdos)) write_pdos_file = write_pdos
1125 curve_section => section_vals_get_subs_vals(dft_section, trim(my_print_key)//"%CURVE")
1126 CALL section_vals_get(curve_section, explicit=write_curve)
1127 IF (PRESENT(write_pdos_curve)) write_curve = write_pdos_curve
1128 should_output = btest(cp_print_key_should_output(logger%iter_info, dft_section, &
1129 trim(my_print_key)), cp_p_file)
1130 output_unit = cp_logger_get_default_io_unit(logger)
1131 IF ((.NOT. should_output)) RETURN
1132
1133 CALL timeset(routinen, handle)
1134 iterstep = logger%iter_info%iteration(logger%iter_info%n_rlevel)
1135
1136 IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T3,A,T61,I10))') &
1137 " Calculate k-point PDOS at iteration step ", iterstep
1138
1139 CALL get_qs_env(qs_env=qs_env, &
1140 kpoints=kpoints, &
1141 dft_control=dft_control, &
1142 matrix_s_kp=matrix_s_kp, &
1143 scf_env=scf_env, &
1144 atomic_kind_set=atomic_kind_set, &
1145 qs_kind_set=qs_kind_set, &
1146 particle_set=particle_set)
1147 para_env => kpoints%para_env_inter_kp
1148 nspins = dft_control%nspins
1149 nkind = SIZE(atomic_kind_set)
1150 CALL get_qs_kind_set(qs_kind_set, maxlgto=maxlgto)
1151 IF (.NOT. ASSOCIATED(kpoints%kp_env)) THEN
1152 cpwarn("No local k points available for k-point PDOS")
1153 CALL timestop(handle)
1154 RETURN
1155 END IF
1156 IF (SIZE(kpoints%kp_env) == 0) THEN
1157 cpwarn("No local k points available for k-point PDOS")
1158 CALL timestop(handle)
1159 RETURN
1160 END IF
1161
1162 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%DELTA_E", r_val=de)
1163 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%APPEND", l_val=append)
1164 IF (trim(my_print_key) == "PRINT%DOS") THEN
1165 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%PDOS%COMPONENTS", l_val=separate_components)
1166 ELSE
1167 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%COMPONENTS", l_val=separate_components)
1168 END IF
1169 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%NDIGITS", i_val=ndigits)
1170 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%BROADEN%TYPE", i_val=broaden_type)
1171 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%BROADEN%WIDTH", r_val=broaden_width)
1172 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%BROADEN%VOIGT_MIXING", r_val=voigt_mixing)
1173 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%ENERGY_ZERO", i_val=energy_zero)
1174 ndigits = min(max(ndigits, 1), 10)
1175 IF (write_curve .AND. de <= 0.0_dp) THEN
1176 cpwarn("Broadened k-point PDOS output requires DELTA_E > 0 and will be skipped")
1177 CALL timestop(handle)
1178 RETURN
1179 END IF
1180 IF (write_curve) de = max(de, 0.00001_dp)
1181
1182 ldos_section => section_vals_get_subs_vals(dft_section, trim(my_print_key)//"%LDOS")
1183 CALL section_vals_get(ldos_section, n_repetition=nldos)
1184 ldos_section => section_vals_get_subs_vals(dft_section, trim(my_print_key)//"%R_LDOS")
1185 CALL section_vals_get(ldos_section, n_repetition=n_r_ldos)
1186 IF (nldos > 0 .OR. n_r_ldos > 0) THEN
1187 cpwarn("LDOS/R_LDOS are not implemented for k-point PDOS and will be ignored")
1188 END IF
1189 IF (write_pdos_file) THEN
1190 cpwarn("State-resolved k-point PDOS output is not implemented yet")
1191 END IF
1192 IF (.NOT. write_curve) THEN
1193 CALL timestop(handle)
1194 RETURN
1195 END IF
1196 IF (broaden_width <= 0.0_dp) THEN
1197 cpwarn("Broadened k-point PDOS output requires a finite WIDTH and will be skipped")
1198 CALL timestop(handle)
1199 RETURN
1200 END IF
1201
1202 IF (separate_components) THEN
1203 ncomp = nsoset(maxlgto)
1204 ELSE
1205 ncomp = maxlgto + 1
1206 END IF
1207 ALLOCATE (kind_maxl(nkind))
1208 kind_maxl = -1
1209 DO ikind = 1, nkind
1210 NULLIFY (orb_basis_set)
1211 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
1212 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, maxl=maxl)
1213 kind_maxl(ikind) = maxl
1214 END DO
1215
1216 emin = huge(0.0_dp)
1217 emax = -huge(0.0_dp)
1218 e_fermi(:) = 0.0_dp
1219 hoco(:) = -huge(0.0_dp)
1220 fractional_occupation = .false.
1221 IF (kpoints%nkp /= 0) THEN
1222 DO ik = 1, SIZE(kpoints%kp_env)
1223 kp => kpoints%kp_env(ik)%kpoint_env
1224 DO ispin = 1, nspins
1225 mo_set => kp%mos(1, ispin)
1226 CALL get_mo_set(mo_set=mo_set, nmo=nmo_kp, mu=e_fermi(ispin))
1227 eigenvalues => mo_set%eigenvalues
1228 occupation_numbers => mo_set%occupation_numbers
1229 DO imo = 1, nmo_kp
1230 IF (occupation_numbers(imo) > 1.0e-10_dp) hoco(ispin) = max(hoco(ispin), eigenvalues(imo))
1231 IF (abs(occupation_numbers(imo) - real(nint(occupation_numbers(imo)), kind=dp)) > &
1232 1.0e-8_dp) fractional_occupation = .true.
1233 END DO
1234 e1 = minval(eigenvalues(1:nmo_kp))
1235 e2 = maxval(eigenvalues(1:nmo_kp))
1236 emin = min(emin, e1)
1237 emax = max(emax, e2)
1238 END DO
1239 END DO
1240 END IF
1241 CALL para_env%min(emin)
1242 CALL para_env%max(emax)
1243 CALL para_env%max(e_fermi)
1244 CALL para_env%max(hoco)
1245 fractional_occupation_int = merge(1, 0, fractional_occupation)
1246 CALL para_env%max(fractional_occupation_int)
1247 fractional_occupation = (fractional_occupation_int /= 0)
1248 DO ispin = 1, nspins
1249 IF (hoco(ispin) < -0.5_dp*huge(0.0_dp)) hoco(ispin) = e_fermi(ispin)
1250 END DO
1251 resolved_energy_zero = dos_resolve_energy_zero(energy_zero, dft_control%smear, fractional_occupation)
1252 SELECT CASE (resolved_energy_zero)
1253 CASE (dos_energy_zero_absolute)
1254 energy_ref(:) = 0.0_dp
1255 CASE (dos_energy_zero_hoco)
1256 energy_ref(:) = maxval(hoco(1:nspins))
1257 CASE DEFAULT
1258 energy_ref(:) = maxval(e_fermi(1:nspins))
1259 END SELECT
1260 zero_label = dos_energy_zero_label(resolved_energy_zero)
1261 IF (energy_zero == dos_energy_zero_auto) zero_label = "AUTO -> "//trim(zero_label)
1262 broaden_cutoff = broadening_cutoff(broaden_type, broaden_width)
1263 emin = emin - broaden_cutoff
1264 emax = emax + broaden_cutoff
1265 nhist = nint((emax - emin)/de) + 1
1266 ALLOCATE (pdos_curve(nhist, ncomp, nkind, nspins))
1267 pdos_curve = 0.0_dp
1268
1269 ! Ensure that S(k)^1/2 is available for the Lowdin projection.
1270 ! This is normally only constructed for Lowdin population/DFT+U paths.
1271 CALL diag_kp_smat(matrix_s_kp, kpoints, scf_env%scf_work1)
1272
1273 ! Use the first local k point to construct the AO -> kind/l/component map.
1274 kp => kpoints%kp_env(1)%kpoint_env
1275 mo_set => kp%mos(1, 1)
1276 CALL get_mo_set(mo_set=mo_set, nao=nao)
1277 CALL build_pdos_ao_map(qs_kind_set, particle_set, nao, ao_kind, ao_l, ao_comp)
1278 ALLOCATE (ao_weight(nao), proj_weight(ncomp, nkind))
1279 ALLOCATE (vecbuffer(1, nao), zvecbuffer(1, nao))
1280
1281 IF (kpoints%nkp /= 0) THEN
1282 DO ik = 1, SIZE(kpoints%kp_env)
1283 kp => kpoints%kp_env(ik)%kpoint_env
1284 wkp = kp%wkp
1285 DO ispin = 1, nspins
1286 mo_set => kp%mos(1, ispin)
1287 CALL get_mo_set(mo_set=mo_set, nao=nao, nmo=nmo_kp)
1288 eigenvalues => mo_set%eigenvalues
1289 CALL cp_fm_get_info(mo_set%mo_coeff, matrix_struct=fm_struct_tmp)
1290 IF (kpoints%use_real_wfn) THEN
1291 CALL cp_fm_create(shalfc, fm_struct_tmp, name="shalfc")
1292 CALL lowdin_kp_mo_coeff(kp, ispin, kpoints%use_real_wfn, shalfc=shalfc)
1293 ELSE
1294 CALL cp_cfm_create(cshalfc, fm_struct_tmp, name="cshalfc")
1295 CALL lowdin_kp_mo_coeff(kp, ispin, kpoints%use_real_wfn, cshalfc=cshalfc)
1296 END IF
1297
1298 DO imo = 1, nmo_kp
1299 IF (kpoints%use_real_wfn) THEN
1300 CALL cp_fm_get_submatrix(shalfc, vecbuffer, 1, imo, nao, 1, transpose=.true.)
1301 ao_weight(:) = vecbuffer(1, 1:nao)**2
1302 ELSE
1303 CALL cp_cfm_get_submatrix(cshalfc, zvecbuffer, 1, imo, nao, 1, transpose=.true.)
1304 ao_weight(:) = real(conjg(zvecbuffer(1, 1:nao))*zvecbuffer(1, 1:nao), kind=dp)
1305 END IF
1306 proj_weight = 0.0_dp
1307 CALL accumulate_pdos_weights(ao_weight, ao_kind, ao_l, ao_comp, &
1308 separate_components, proj_weight)
1309 DO ikind = 1, nkind
1310 IF (kind_maxl(ikind) < 0) cycle
1311 IF (separate_components) THEN
1312 DO icomp = 1, nsoset(kind_maxl(ikind))
1313 CALL add_broadened_value(pdos_curve(:, icomp, ikind, ispin), &
1314 emin, de, eigenvalues(imo), &
1315 wkp*proj_weight(icomp, ikind), &
1316 broaden_type, broaden_width, voigt_mixing)
1317 END DO
1318 ELSE
1319 DO icomp = 1, kind_maxl(ikind) + 1
1320 CALL add_broadened_value(pdos_curve(:, icomp, ikind, ispin), &
1321 emin, de, eigenvalues(imo), &
1322 wkp*proj_weight(icomp, ikind), &
1323 broaden_type, broaden_width, voigt_mixing)
1324 END DO
1325 END IF
1326 END DO
1327 END DO
1328
1329 IF (kpoints%use_real_wfn) THEN
1330 CALL cp_fm_release(shalfc)
1331 ELSE
1332 CALL cp_cfm_release(cshalfc)
1333 END IF
1334 END DO
1335 END DO
1336 END IF
1337 CALL para_env%sum(pdos_curve)
1338
1339 IF (append .AND. iterstep > 1) THEN
1340 my_pos = "APPEND"
1341 ELSE
1342 my_pos = "REWIND"
1343 END IF
1344 my_act = "WRITE"
1345 spin(1) = "ALPHA"
1346 spin(2) = "BETA"
1347 IF (write_pdos_file) THEN
1348 DO ikind = 1, nkind
1349 IF (kind_maxl(ikind) < 0) cycle
1350 CALL get_atomic_kind(atomic_kind_set(ikind), name=kind_name)
1351 DO ispin = 1, nspins
1352 IF (nspins == 2) THEN
1353 my_mittle = trim(spin(ispin))//"_k"//trim(adjustl(cp_to_string(ikind)))
1354 ELSE
1355 my_mittle = "k"//trim(adjustl(cp_to_string(ikind)))
1356 END IF
1357 CALL write_broadened_pdos_curve(logger, dft_section, trim(my_mittle)//"_curve", my_pos, my_act, &
1358 iterstep, e_fermi(ispin), hoco(ispin), energy_ref(ispin), &
1359 trim(zero_label), &
1360 "K-point projected DOS for atomic kind "//trim(kind_name), &
1361 kind_maxl(ikind), separate_components, &
1362 pdos_curve(:, :, ikind, ispin), emin, de, &
1363 broaden_type, broaden_width, voigt_mixing, ndigits, pdos_print_key=trim(my_print_key))
1364 END DO
1365 END DO
1366 END IF
1367
1368 DEALLOCATE (ao_comp, ao_kind, ao_l, ao_weight, kind_maxl, pdos_curve, proj_weight, &
1369 vecbuffer, zvecbuffer)
1370
1371 CALL timestop(handle)
1372
1373 END SUBROUTINE calculate_projected_dos_kp
1374
1375! **************************************************************************************************
1376!> \brief Build AO mapping arrays for PDOS accumulation.
1377!> \param qs_kind_set ...
1378!> \param particle_set ...
1379!> \param nao ...
1380!> \param ao_kind ...
1381!> \param ao_l ...
1382!> \param ao_comp ...
1383! **************************************************************************************************
1384 SUBROUTINE build_pdos_ao_map(qs_kind_set, particle_set, nao, ao_kind, ao_l, ao_comp)
1385
1386 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1387 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1388 INTEGER, INTENT(IN) :: nao
1389 INTEGER, ALLOCATABLE, DIMENSION(:), INTENT(OUT) :: ao_kind, ao_l, ao_comp
1390
1391 INTEGER :: iatom, ikind, irow, iset, ishell, iso, &
1392 lshell, maxl, nset
1393 INTEGER, DIMENSION(:), POINTER :: nshell
1394 INTEGER, DIMENSION(:, :), POINTER :: l
1395 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
1396
1397 ALLOCATE (ao_kind(nao), ao_l(nao), ao_comp(nao))
1398 irow = 0
1399 DO iatom = 1, SIZE(particle_set)
1400 NULLIFY (orb_basis_set)
1401 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
1402 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
1403 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1404 nset=nset, &
1405 nshell=nshell, &
1406 l=l, maxl=maxl)
1407 DO iset = 1, nset
1408 DO ishell = 1, nshell(iset)
1409 lshell = l(ishell, iset)
1410 DO iso = 1, nso(lshell)
1411 irow = irow + 1
1412 cpassert(irow <= nao)
1413 ao_kind(irow) = ikind
1414 ao_l(irow) = lshell
1415 ao_comp(irow) = nsoset(lshell - 1) + iso
1416 END DO
1417 END DO
1418 END DO
1419 END DO
1420 cpassert(irow == nao)
1421
1422 END SUBROUTINE build_pdos_ao_map
1423
1424! **************************************************************************************************
1425!> \brief Accumulate AO weights into kind/l or kind/component projected weights.
1426!> \param ao_weight ...
1427!> \param ao_kind ...
1428!> \param ao_l ...
1429!> \param ao_comp ...
1430!> \param separate_components ...
1431!> \param proj_weight ...
1432! **************************************************************************************************
1433 SUBROUTINE accumulate_pdos_weights(ao_weight, ao_kind, ao_l, ao_comp, &
1434 separate_components, proj_weight)
1435
1436 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: ao_weight
1437 INTEGER, DIMENSION(:), INTENT(IN) :: ao_kind, ao_l, ao_comp
1438 LOGICAL, INTENT(IN) :: separate_components
1439 REAL(kind=dp), DIMENSION(:, :), INTENT(INOUT) :: proj_weight
1440
1441 INTEGER :: iao, icomp, ikind
1442
1443 DO iao = 1, SIZE(ao_weight)
1444 ikind = ao_kind(iao)
1445 IF (separate_components) THEN
1446 icomp = ao_comp(iao)
1447 ELSE
1448 icomp = ao_l(iao) + 1
1449 END IF
1450 proj_weight(icomp, ikind) = proj_weight(icomp, ikind) + ao_weight(iao)
1451 END DO
1452
1453 END SUBROUTINE accumulate_pdos_weights
1454
1455! **************************************************************************************************
1456!> \brief Write a broadened k-point PDOS curve.
1457!> \param logger ...
1458!> \param dft_section ...
1459!> \param middle_name ...
1460!> \param file_position ...
1461!> \param file_action ...
1462!> \param iterstep ...
1463!> \param e_fermi ...
1464!> \param hoco ...
1465!> \param energy_ref ...
1466!> \param zero_label ...
1467!> \param title ...
1468!> \param maxl ...
1469!> \param separate_components ...
1470!> \param pdos_curve ...
1471!> \param emin ...
1472!> \param de ...
1473!> \param broaden_type ...
1474!> \param broaden_width ...
1475!> \param voigt_mixing ...
1476!> \param ndigits ...
1477!> \param pdos_print_key ...
1478! **************************************************************************************************
1479 SUBROUTINE write_broadened_pdos_curve(logger, dft_section, middle_name, file_position, &
1480 file_action, iterstep, e_fermi, hoco, energy_ref, zero_label, &
1481 title, maxl, &
1482 separate_components, pdos_curve, emin, de, &
1483 broaden_type, broaden_width, voigt_mixing, ndigits, pdos_print_key)
1484
1485 TYPE(cp_logger_type), POINTER :: logger
1486 TYPE(section_vals_type), POINTER :: dft_section
1487 CHARACTER(LEN=*), INTENT(IN) :: middle_name, file_position, file_action
1488 INTEGER, INTENT(IN) :: iterstep
1489 REAL(kind=dp), INTENT(IN) :: e_fermi, hoco, energy_ref
1490 CHARACTER(LEN=*), INTENT(IN) :: zero_label, title
1491 INTEGER, INTENT(IN) :: maxl
1492 LOGICAL, INTENT(IN) :: separate_components
1493 REAL(kind=dp), DIMENSION(:, :), INTENT(IN) :: pdos_curve
1494 REAL(kind=dp), INTENT(IN) :: emin, de
1495 INTEGER, INTENT(IN) :: broaden_type
1496 REAL(kind=dp), INTENT(IN) :: broaden_width, voigt_mixing
1497 INTEGER, INTENT(IN) :: ndigits
1498 CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: pdos_print_key
1499
1500 CHARACTER(LEN=16) :: energy_label
1501 CHARACTER(LEN=20) :: fmtstr_data
1502 CHARACTER(LEN=6), ALLOCATABLE, DIMENSION(:, :, :) :: tmp_str
1503 CHARACTER(LEN=default_string_length) :: my_print_key
1504 INTEGER :: energy_unit, i, icomp, il, im, iw, &
1505 ncomponents, nhist
1506 REAL(kind=dp) :: density_factor, energy_factor, &
1507 ev_factor, eval
1508
1509 my_print_key = "PRINT%PDOS"
1510 IF (PRESENT(pdos_print_key)) my_print_key = trim(pdos_print_key)
1511
1512 nhist = SIZE(pdos_curve, 1)
1513 IF (separate_components) THEN
1514 ncomponents = nsoset(maxl)
1515 ELSE
1516 ncomponents = maxl + 1
1517 END IF
1518 CALL section_vals_val_get(dft_section, trim(my_print_key)//"%CURVE%ENERGY_UNIT", i_val=energy_unit)
1519 energy_factor = dos_energy_scale(energy_unit)
1520 density_factor = dos_density_scale(energy_unit)
1521 energy_label = dos_energy_label(energy_unit)
1522 ev_factor = dos_energy_scale(dos_energy_unit_ev)
1523 iw = cp_print_key_unit_nr(logger, dft_section, trim(my_print_key), &
1524 extension=".pdos", file_position=file_position, file_action=file_action, &
1525 file_form="FORMATTED", middle_name=trim(middle_name))
1526 IF (iw > 0) THEN
1527 WRITE (unit=iw, fmt="(A,I0)") "# "//trim(title)//" at iteration step i = ", iterstep
1528 WRITE (unit=iw, fmt="(A,F12.6,A,F12.6,A)") &
1529 "# E(Fermi) = ", e_fermi, " a.u. = ", e_fermi*ev_factor, " eV"
1530 WRITE (unit=iw, fmt="(A,F12.6,A,F12.6,A)") &
1531 "# E(HOCO) = ", hoco, " a.u. = ", hoco*ev_factor, " eV"
1532 WRITE (unit=iw, fmt="(A,A)") "# Energy zero: ", trim(zero_label)
1533 CALL write_broadening_info(iw, broaden_type, broaden_width, voigt_mixing)
1534 WRITE (unit=iw, fmt="(A)", advance="NO") "# "//trim(energy_label)
1535 WRITE (unit=iw, fmt="(2X,A)", advance="NO") "total"
1536 IF (separate_components) THEN
1537 ALLOCATE (tmp_str(0:0, 0:maxl, -maxl:maxl))
1538 tmp_str = ""
1539 DO il = 0, maxl
1540 DO im = -il, il
1541 tmp_str(0, il, im) = sgf_symbol(0, il, im)
1542 WRITE (unit=iw, fmt="(2X,A)", advance="NO") trim(tmp_str(0, il, im))
1543 END DO
1544 END DO
1545 DEALLOCATE (tmp_str)
1546 ELSE
1547 DO il = 0, maxl
1548 WRITE (unit=iw, fmt="(2X,A)", advance="NO") trim(l_sym(il))
1549 END DO
1550 END IF
1551 WRITE (unit=iw, fmt="()")
1552 WRITE (unit=fmtstr_data, fmt="(A,I0,A)") "(2X,F20.", ndigits, ")"
1553 DO i = 1, nhist
1554 eval = (emin + (i - 1)*de - energy_ref)*energy_factor
1555 WRITE (unit=iw, fmt="(F15.8)", advance="NO") eval
1556 WRITE (unit=iw, fmt=fmtstr_data, advance="NO") sum(pdos_curve(i, 1:ncomponents))*density_factor
1557 DO icomp = 1, ncomponents
1558 WRITE (unit=iw, fmt=fmtstr_data, advance="NO") pdos_curve(i, icomp)*density_factor
1559 END DO
1560 WRITE (unit=iw, fmt="()")
1561 END DO
1562 END IF
1563 CALL cp_print_key_finished_output(iw, logger, dft_section, trim(my_print_key))
1564
1565 END SUBROUTINE write_broadened_pdos_curve
1566
1567! **************************************************************************************************
1568!> \brief Write a broadened PDOS curve for a projected weight matrix.
1569!> \param logger ...
1570!> \param dft_section ...
1571!> \param middle_name ...
1572!> \param file_position ...
1573!> \param file_action ...
1574!> \param iterstep ...
1575!> \param e_fermi ...
1576!> \param hoco ...
1577!> \param energy_ref ...
1578!> \param zero_label ...
1579!> \param title ...
1580!> \param maxl ...
1581!> \param separate_components ...
1582!> \param weights ...
1583!> \param eigenvalues ...
1584!> \param nstates ...
1585!> \param de ...
1586!> \param broaden_type ...
1587!> \param broaden_width ...
1588!> \param voigt_mixing ...
1589!> \param ndigits ...
1590!> \param pdos_print_key ...
1591! **************************************************************************************************
1592 SUBROUTINE write_broadened_pdos(logger, dft_section, middle_name, file_position, file_action, &
1593 iterstep, e_fermi, hoco, energy_ref, zero_label, title, maxl, separate_components, weights, &
1594 eigenvalues, nstates, de, broaden_type, broaden_width, &
1595 voigt_mixing, ndigits, pdos_print_key)
1596
1597 TYPE(cp_logger_type), POINTER :: logger
1598 TYPE(section_vals_type), POINTER :: dft_section
1599 CHARACTER(LEN=*), INTENT(IN) :: middle_name, file_position, file_action
1600 INTEGER, INTENT(IN) :: iterstep
1601 REAL(kind=dp), INTENT(IN) :: e_fermi, hoco, energy_ref
1602 CHARACTER(LEN=*), INTENT(IN) :: zero_label, title
1603 INTEGER, INTENT(IN) :: maxl
1604 LOGICAL, INTENT(IN) :: separate_components
1605 REAL(kind=dp), DIMENSION(:, :), INTENT(IN) :: weights
1606 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: eigenvalues
1607 INTEGER, INTENT(IN) :: nstates
1608 REAL(kind=dp), INTENT(IN) :: de
1609 INTEGER, INTENT(IN) :: broaden_type
1610 REAL(kind=dp), INTENT(IN) :: broaden_width, voigt_mixing
1611 INTEGER, INTENT(IN) :: ndigits
1612 CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: pdos_print_key
1613
1614 INTEGER :: i, icomp, imo, ncomponents, nhist
1615 REAL(kind=dp) :: cutoff, emax, emin, eval, line_shape
1616 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: pdos_curve
1617
1618 IF (broaden_width <= 0.0_dp) RETURN
1619
1620 ncomponents = SIZE(weights, 1)
1621 cutoff = broadening_cutoff(broaden_type, broaden_width)
1622 emin = minval(eigenvalues(1:nstates)) - cutoff
1623 emax = maxval(eigenvalues(1:nstates)) + cutoff
1624 nhist = nint((emax - emin)/de) + 1
1625 ALLOCATE (pdos_curve(nhist, ncomponents))
1626 pdos_curve = 0.0_dp
1627
1628 DO imo = 1, nstates
1629 DO i = max(1, floor((eigenvalues(imo) - cutoff - emin)/de) + 1), &
1630 min(nhist, ceiling((eigenvalues(imo) + cutoff - emin)/de) + 1)
1631 eval = emin + (i - 1)*de
1632 line_shape = broadening_function(eval - eigenvalues(imo), broaden_type, broaden_width, &
1633 voigt_mixing)
1634 DO icomp = 1, ncomponents
1635 pdos_curve(i, icomp) = pdos_curve(i, icomp) + weights(icomp, imo)*line_shape
1636 END DO
1637 END DO
1638 END DO
1639
1640 CALL write_broadened_pdos_curve(logger, dft_section, middle_name, file_position, file_action, &
1641 iterstep, e_fermi, hoco, energy_ref, zero_label, title, maxl, separate_components, &
1642 pdos_curve, emin, de, broaden_type, broaden_width, &
1643 voigt_mixing, ndigits, pdos_print_key=pdos_print_key)
1644
1645 DEALLOCATE (pdos_curve)
1646
1647 END SUBROUTINE write_broadened_pdos
1648
1649END MODULE qs_pdos
modulo
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
Definition grid_common.h:125
cell_types::pbc
Definition cell_types.F:103
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:88
cp_log_handling::cp_to_string
Definition cp_log_handling.F:90
memory_utilities::reallocate
Definition memory_utilities.F:27
parallel_gemm_api::parallel_gemm
Definition parallel_gemm_api.F:39
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
atomic_kind_types::get_atomic_kind_set
subroutine, public get_atomic_kind_set(atomic_kind_set, atom_of_kind, kind_of, natom_of_kind, maxatom, natom, nshell, fist_potential_present, shell_present, shell_adiabatic, shell_check_distance, damping_present)
Get attributes of an atomic kind set.
Definition atomic_kind_types.F:223
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
basis_set_types
Definition basis_set_types.F:15
basis_set_types::get_gto_basis_set
subroutine, public get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, npgf_sum, nshell_sum, maxder, short_kind_radius, npgf_seg_sum, ccon)
...
Definition basis_set_types.F:643
cell_types
Handles all functions related to the CELL.
Definition cell_types.F:15
cp_array_utils
various utilities that regard array of different kinds: output, allocation,... maybe it is not a good...
Definition cp_array_utils.F:18
cp_blacs_env
methods related to the blacs parallel environment
Definition cp_blacs_env.F:15
cp_cfm_types
Represents a complex full matrix distributed on many processors.
Definition cp_cfm_types.F:12
cp_cfm_types::cp_cfm_get_submatrix
subroutine, public cp_cfm_get_submatrix(fm, target_m, start_row, start_col, n_rows, n_cols, transpose)
Extract a sub-matrix from the full matrix: op(target_m)(1:n_rows,1:n_cols) = fm(start_row:start_row+n...
Definition cp_cfm_types.F:356
cp_cfm_types::cp_cfm_release
subroutine, public cp_cfm_release(matrix)
Releases a full matrix.
Definition cp_cfm_types.F:182
cp_cfm_types::cp_cfm_create
subroutine, public cp_cfm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
Creates a new full matrix with the given structure.
Definition cp_cfm_types.F:125
cp_control_types
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
Definition cp_control_types.F:12
cp_dbcsr_api
Definition cp_dbcsr_api.F:8
cp_dbcsr_operations
DBCSR operations in CP2K.
Definition cp_dbcsr_operations.F:18
cp_dbcsr_operations::copy_dbcsr_to_fm
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
Definition cp_dbcsr_operations.F:212
cp_fm_diag
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
Definition cp_fm_diag.F:17
cp_fm_diag::cp_fm_power
subroutine, public cp_fm_power(matrix, work, exponent, threshold, n_dependent, verbose, eigvals)
...
Definition cp_fm_diag.F:1127
cp_fm_struct
represent the structure of a full matrix
Definition cp_fm_struct.F:14
cp_fm_struct::cp_fm_struct_create
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
Definition cp_fm_struct.F:132
cp_fm_struct::cp_fm_struct_release
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition cp_fm_struct.F:355
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_get_info
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
Definition cp_fm_types.F:1038
cp_fm_types::cp_fm_create
subroutine, public cp_fm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
creates a new full matrix with the given structure
Definition cp_fm_types.F:165
cp_fm_types::cp_fm_get_submatrix
subroutine, public cp_fm_get_submatrix(fm, target_m, start_row, start_col, n_rows, n_cols, transpose)
gets a submatrix of a full matrix op(target_m)(1:n_rows,1:n_cols) =fm(start_row:start_row+n_rows,...
Definition cp_fm_types.F:925
cp_log_handling
various routines to log and control the output. The idea is that decisions about where to log should ...
Definition cp_log_handling.F:41
cp_log_handling::cp_logger_get_default_io_unit
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
Definition cp_log_handling.F:523
cp_log_handling::cp_get_default_logger
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
Definition cp_log_handling.F:234
cp_output_handling
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
Definition cp_output_handling.F:25
cp_output_handling::cp_print_key_unit_nr
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
Definition cp_output_handling.F:853
cp_output_handling::cp_print_key_finished_output
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
Definition cp_output_handling.F:1073
cp_output_handling::cp_p_file
integer, parameter, public cp_p_file
Definition cp_output_handling.F:103
cp_output_handling::cp_print_key_should_output
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
Definition cp_output_handling.F:345
input_section_types
objects that represent the structure of input sections and the data contained in an input section
Definition input_section_types.F:15
input_section_types::section_vals_get_subs_vals
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
Definition input_section_types.F:735
input_section_types::section_vals_get
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
Definition input_section_types.F:708
input_section_types::section_vals_val_get
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Definition input_section_types.F:1052
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kinds::default_string_length
integer, parameter, public default_string_length
Definition kinds.F:57
kpoint_methods
Routines needed for kpoint calculation.
Definition kpoint_methods.F:15
kpoint_methods::lowdin_kp_mo_coeff
subroutine, public lowdin_kp_mo_coeff(kp, ispin, use_real_wfn, shalfc, cshalfc)
Calculate S(k)^1/2 C(k) for real or complex k-point wavefunctions.
Definition kpoint_methods.F:1697
kpoint_types
Types and basic routines needed for a kpoint calculation.
Definition kpoint_types.F:15
list
An array-based list which grows on demand. When the internal array is full, a new array of twice the ...
Definition list.F:24
memory_utilities
Utility routines for the memory handling.
Definition memory_utilities.F:14
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:23
orbital_pointers
Provides Cartesian and spherical orbital pointers and indices.
Definition orbital_pointers.F:15
orbital_pointers::nsoset
integer, dimension(:), allocatable, public nsoset
Definition orbital_pointers.F:42
orbital_pointers::nso
integer, dimension(:), allocatable, public nso
Definition orbital_pointers.F:42
orbital_symbols
orbital_symbols
Definition orbital_symbols.F:14
orbital_symbols::sgf_symbol
character(len=6) function, public sgf_symbol(n, l, m)
Build a spherical orbital symbol (orbital labels for printing).
Definition orbital_symbols.F:107
orbital_symbols::l_sym
character(len=1), dimension(0:11), parameter, public l_sym
Definition orbital_symbols.F:24
parallel_gemm_api
basic linear algebra operations for full matrixes
Definition parallel_gemm_api.F:14
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
pw_env_types
container for various plainwaves related things
Definition pw_env_types.F:14
pw_env_types::pw_env_get
subroutine, public pw_env_get(pw_env, pw_pools, cube_info, gridlevel_info, auxbas_pw_pool, auxbas_grid, auxbas_rs_desc, auxbas_rs_grid, rs_descs, rs_grids, xc_pw_pool, vdw_pw_pool, poisson_env, interp_section)
returns the various attributes of the pw env
Definition pw_env_types.F:113
pw_pool_types
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Definition pw_pool_types.F:24
pw_types
Definition pw_types.F:24
qs_collocate_density
Calculate the plane wave density by collocating the primitive Gaussian functions (pgf).
Definition qs_collocate_density.F:38
qs_collocate_density::calculate_wavefunction
subroutine, public calculate_wavefunction(mo_vectors, ivector, rho, rho_gspace, atomic_kind_set, qs_kind_set, cell, dft_control, particle_set, pw_env, basis_type)
maps a given wavefunction on the grid
Definition qs_collocate_density.F:2867
qs_dos_utils
Utilities for broadened DOS and PDOS output.
Definition qs_dos_utils.F:11
qs_dos_utils::broadening_function
pure real(kind=dp) function, public broadening_function(delta_e, broaden_type, broaden_width, voigt_mixing)
Normalized broadening function. BROADEN_WIDTH is FWHM.
Definition qs_dos_utils.F:281
qs_dos_utils::write_broadening_info
subroutine, public write_broadening_info(iw, broaden_type, broaden_width, voigt_mixing)
Write broadening metadata.
Definition qs_dos_utils.F:320
qs_dos_utils::add_broadened_value
subroutine, public add_broadened_value(curve, emin, de, eig, weight, broaden_type, broaden_width, voigt_mixing)
Add a broadened spectral line with a scalar weight to a curve.
Definition qs_dos_utils.F:232
qs_dos_utils::dos_energy_zero_hoco
integer, parameter, public dos_energy_zero_hoco
Definition qs_dos_utils.F:31
qs_dos_utils::dos_energy_zero_auto
integer, parameter, public dos_energy_zero_auto
Definition qs_dos_utils.F:31
qs_dos_utils::dos_energy_unit_ev
integer, parameter, public dos_energy_unit_ev
Definition qs_dos_utils.F:29
qs_dos_utils::dos_energy_label
character(len=16) function, public dos_energy_label(energy_unit)
Return the energy-column label for DOS-like output.
Definition qs_dos_utils.F:81
qs_dos_utils::dos_resolve_energy_zero
integer function, public dos_resolve_energy_zero(energy_zero, smearing_enabled, fractional_occupation)
Resolve AUTO energy-zero selection for DOS-like output.
Definition qs_dos_utils.F:125
qs_dos_utils::dos_energy_scale
real(kind=dp) function, public dos_energy_scale(energy_unit)
Return the conversion factor from internal energy units to the selected DOS energy unit.
Definition qs_dos_utils.F:48
qs_dos_utils::dos_energy_zero_absolute
integer, parameter, public dos_energy_zero_absolute
Definition qs_dos_utils.F:31
qs_dos_utils::broadening_cutoff
pure real(kind=dp) function, public broadening_cutoff(broaden_type, broaden_width)
Broadening cutoff used for numerical accumulation.
Definition qs_dos_utils.F:259
qs_dos_utils::dos_density_scale
real(kind=dp) function, public dos_density_scale(energy_unit)
Return the DOS-density conversion factor for the selected energy unit.
Definition qs_dos_utils.F:67
qs_dos_utils::dos_energy_zero_label
character(len=16) function, public dos_energy_zero_label(energy_zero)
Return the label for the selected DOS energy zero.
Definition qs_dos_utils.F:100
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, 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:533
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_kind_types::get_qs_kind
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, cneo_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zatom, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, u_of_dft_plus_u, j_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, monovalent, floating, name, element_symbol, pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
Definition qs_kind_types.F:468
qs_kind_types::get_qs_kind_set
subroutine, public get_qs_kind_set(qs_kind_set, all_potential_present, tnadd_potential_present, gth_potential_present, sgp_potential_present, paw_atom_present, dft_plus_u_atom_present, maxcgf, maxsgf, maxco, maxco_proj, maxgtops, maxlgto, maxlprj, maxnset, maxsgf_set, ncgf, npgf, nset, nsgf, nshell, maxpol, maxlppl, maxlppnl, maxppnl, nelectron, maxder, max_ngrid_rad, max_sph_harm, maxg_iso_not0, lmax_rho0, basis_rcut, basis_type, total_zeff_corr, npgf_seg, cneo_potential_present, nkind_q, natom_q)
Get attributes of an atomic kind set.
Definition qs_kind_types.F:913
qs_mo_types
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
qs_mo_types::get_mo_set
subroutine, public get_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, mo_coeff, mo_coeff_b, uniform_occupation, kts, mu, flexible_electron_count)
Get the components of a MO set data structure.
Definition qs_mo_types.F:401
qs_pdos
Calculation and writing of projected density of states The DOS is computed per angular momentum and p...
Definition qs_pdos.F:15
qs_pdos::calculate_projected_dos_kp
subroutine, public calculate_projected_dos_kp(qs_env, dft_section, pdos_print_key, write_pdos, write_pdos_curve)
Compute and write broadened projected density of states for k-point calculations.
Definition qs_pdos.F:1073
qs_pdos::calculate_projected_dos
subroutine, public calculate_projected_dos(mo_set, atomic_kind_set, qs_kind_set, particle_set, qs_env, dft_section, ispin, xas_mittle, external_matrix_shalf, unoccupied_orbs, unoccupied_evals, pdos_print_key, write_pdos, write_pdos_curve)
Compute and write projected density of states.
Definition qs_pdos.F:155
qs_scf_diagonalization
Different diagonalization schemes that can be used for the iterative solution of the eigenvalue probl...
Definition qs_scf_diagonalization.F:15
qs_scf_diagonalization::diag_kp_smat
subroutine, public diag_kp_smat(matrix_s, kpoints, fmwork)
Kpoint diagonalization routine Transforms matrices to kpoint, distributes kpoint groups,...
Definition qs_scf_diagonalization.F:2135
qs_scf_types
module that contains the definitions of the scf types
Definition qs_scf_types.F:14
atomic_kind_types::atomic_kind_type
Provides all information about an atomic kind.
Definition atomic_kind_types.F:45
basis_set_types::gto_basis_set_type
Definition basis_set_types.F:72
cell_types::cell_type
Type defining parameters related to the simulation cell.
Definition cell_types.F:60
cp_array_utils::cp_1d_r_p_type
represent a pointer to a 1d array
Definition cp_array_utils.F:99
cp_blacs_env::cp_blacs_env_type
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
Definition cp_blacs_env.F:53
cp_cfm_types::cp_cfm_type
Represent a complex full matrix.
Definition cp_cfm_types.F:69
cp_control_types::dft_control_type
Definition cp_control_types.F:757
cp_dbcsr_api::dbcsr_p_type
Definition cp_dbcsr_api.F:172
cp_fm_struct::cp_fm_struct_type
keeps the information about the structure of a full matrix
Definition cp_fm_struct.F:82
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:114
cp_log_handling::cp_logger_type
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Definition cp_log_handling.F:140
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:127
kpoint_types::kpoint_env_type
Keeps information about a specific k-point.
Definition kpoint_types.F:83
kpoint_types::kpoint_type
Contains information about kpoints.
Definition kpoint_types.F:175
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:746
particle_types::particle_type
Definition particle_types.F:35
pw_env_types::pw_env_type
contained for different pw related things
Definition pw_env_types.F:70
pw_pool_types::pw_pool_p_type
to create arrays of pools
Definition pw_pool_types.F:135
pw_pool_types::pw_pool_type
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Definition pw_pool_types.F:83
pw_types::pw_c1d_gs_type
Definition pw_types.F:127
pw_types::pw_r3d_rs_type
Definition pw_types.F:62
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:222
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:177
qs_mo_types::mo_set_type
Definition qs_mo_types.F:47
qs_scf_types::qs_scf_env_type
Definition qs_scf_types.F:97