(git:58e3e09)
qs_pdos.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 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 ! **************************************************************************************************
15 MODULE qs_pdos
16  USE atomic_kind_types, ONLY: atomic_kind_type,&
20  gto_basis_set_type
21  USE cell_types, ONLY: cell_type,&
22  pbc
23  USE cp_blacs_env, ONLY: cp_blacs_env_type
24  USE cp_control_types, ONLY: dft_control_type
26  USE cp_fm_diag, ONLY: cp_fm_power
29  cp_fm_struct_type
30  USE cp_fm_types, ONLY: cp_fm_create,&
34  cp_fm_release,&
35  cp_fm_type
38  cp_logger_type,&
39  cp_to_string
40  USE cp_output_handling, ONLY: cp_p_file,&
44  USE dbcsr_api, ONLY: dbcsr_p_type
47  section_vals_type,&
49  USE kinds, ONLY: default_string_length,&
50  dp
51  USE memory_utilities, ONLY: reallocate
52  USE message_passing, ONLY: mp_para_env_type
53  USE orbital_pointers, ONLY: nso,&
54  nsoset
55  USE orbital_symbols, ONLY: l_sym,&
57  USE parallel_gemm_api, ONLY: parallel_gemm
58  USE particle_types, ONLY: particle_type
59  USE preconditioner_types, ONLY: preconditioner_type
60  USE pw_env_types, ONLY: pw_env_get,&
61  pw_env_type
62  USE pw_pool_types, ONLY: pw_pool_p_type,&
63  pw_pool_type
64  USE pw_types, ONLY: pw_c1d_gs_type,&
65  pw_r3d_rs_type
67  USE qs_environment_types, ONLY: get_qs_env,&
68  qs_environment_type
69  USE qs_kind_types, ONLY: get_qs_kind,&
71  qs_kind_type
72  USE qs_mo_methods, ONLY: calculate_subspace_eigenvalues
73  USE qs_mo_types, ONLY: get_mo_set,&
74  mo_set_type
76  USE scf_control_types, ONLY: scf_control_type
77 #include "./base/base_uses.f90"
78 
79  IMPLICIT NONE
80 
81  PRIVATE
82 
83  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_pdos'
84 
85 ! **************************************************************************************************
86  ! *** Public subroutines ***
87 
88  PUBLIC :: calculate_projected_dos
89 
90  TYPE ldos_type
91  INTEGER :: maxl, nlist
92  LOGICAL :: separate_components
93  INTEGER, DIMENSION(:), POINTER :: list_index
94  REAL(KIND=dp), DIMENSION(:, :), &
95  POINTER :: pdos_array
96  END TYPE ldos_type
97 
98  TYPE r_ldos_type
99  INTEGER :: nlist, npoints
100  INTEGER, DIMENSION(:, :), POINTER :: index_grid_local
101  INTEGER, DIMENSION(:), POINTER :: list_index
102  REAL(KIND=dp), DIMENSION(:), POINTER :: x_range, y_range, z_range
103  REAL(KIND=dp), DIMENSION(:), POINTER :: eval_range
104  REAL(KIND=dp), DIMENSION(:), &
105  POINTER :: pdos_array
106  END TYPE r_ldos_type
107 
108  TYPE ldos_p_type
109  TYPE(ldos_type), POINTER :: ldos
110  END TYPE ldos_p_type
111 
112  TYPE r_ldos_p_type
113  TYPE(r_ldos_type), POINTER :: ldos
114  END TYPE r_ldos_p_type
115 CONTAINS
116 
117 ! **************************************************************************************************
118 !> \brief Compute and write projected density of states
119 !> \param mo_set ...
120 !> \param atomic_kind_set ...
121 !> \param qs_kind_set ...
122 !> \param particle_set ...
123 !> \param qs_env ...
124 !> \param dft_section ...
125 !> \param ispin ...
126 !> \param xas_mittle ...
127 !> \param external_matrix_shalf ...
128 !> \date 26.02.2008
129 !> \par History:
130 !> - Added optional external matrix_shalf to avoid recomputing it (A. Bussy, 09.2019)
131 !> \par Variables
132 !> -
133 !> -
134 !> \author MI
135 !> \version 1.0
136 ! **************************************************************************************************
137  SUBROUTINE calculate_projected_dos(mo_set, atomic_kind_set, qs_kind_set, particle_set, qs_env, &
138  dft_section, ispin, xas_mittle, external_matrix_shalf)
139 
140  TYPE(mo_set_type), INTENT(IN) :: mo_set
141  TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
142  TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
143  TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
144  TYPE(qs_environment_type), POINTER :: qs_env
145  TYPE(section_vals_type), POINTER :: dft_section
146  INTEGER, INTENT(IN), OPTIONAL :: ispin
147  CHARACTER(LEN=default_string_length), INTENT(IN), &
148  OPTIONAL :: xas_mittle
149  TYPE(cp_fm_type), INTENT(IN), OPTIONAL, TARGET :: external_matrix_shalf
150 
151  CHARACTER(len=*), PARAMETER :: routinen = 'calculate_projected_dos'
152 
153  CHARACTER(LEN=16) :: fmtstr2
154  CHARACTER(LEN=27) :: fmtstr1
155  CHARACTER(LEN=6), ALLOCATABLE, DIMENSION(:, :, :) :: tmp_str
156  CHARACTER(LEN=default_string_length) :: kind_name, my_act, my_mittle, my_pos, &
157  spin(2)
158  CHARACTER(LEN=default_string_length), &
159  ALLOCATABLE, DIMENSION(:) :: ldos_index, r_ldos_index
160  INTEGER :: handle, homo, i, iatom, ikind, il, ildos, im, imo, in_x, in_y, in_z, ir, irow, &
161  iset, isgf, ishell, iso, iterstep, iw, j, jx, jy, jz, k, lcomponent, lshell, maxl, &
162  maxlgto, my_spin, n_dependent, n_r_ldos, n_rep, nao, natom, ncol_global, nkind, nldos, &
163  nlumo, nmo, np_tot, npoints, nrow_global, nset, nsgf, nvirt, out_each, output_unit
164  INTEGER, ALLOCATABLE, DIMENSION(:) :: firstrow
165  INTEGER, DIMENSION(:), POINTER :: list, nshell
166  INTEGER, DIMENSION(:, :), POINTER :: bo, l
167  LOGICAL :: append, calc_matsh, do_ldos, do_r_ldos, &
168  do_virt, ionode, separate_components, &
169  should_output
170  LOGICAL, DIMENSION(:, :), POINTER :: read_r
171  REAL(kind=dp) :: dh(3, 3), dvol, e_fermi, r(3), r_vec(3), &
172  ratom(3)
173  REAL(kind=dp), DIMENSION(:), POINTER :: eigenvalues, evals_virt, &
174  occupation_numbers
175  REAL(kind=dp), DIMENSION(:, :), POINTER :: vecbuffer
176  REAL(kind=dp), DIMENSION(:, :, :), POINTER :: pdos_array
177  TYPE(cell_type), POINTER :: cell
178  TYPE(cp_blacs_env_type), POINTER :: context
179  TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
180  TYPE(cp_fm_type) :: matrix_shalfc, matrix_work, mo_virt
181  TYPE(cp_fm_type), POINTER :: matrix_shalf, mo_coeff
182  TYPE(cp_logger_type), POINTER :: logger
183  TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: s_matrix
184  TYPE(dft_control_type), POINTER :: dft_control
185  TYPE(gto_basis_set_type), POINTER :: orb_basis_set
186  TYPE(ldos_p_type), DIMENSION(:), POINTER :: ldos_p
187  TYPE(mp_para_env_type), POINTER :: para_env
188  TYPE(pw_c1d_gs_type) :: wf_g
189  TYPE(pw_env_type), POINTER :: pw_env
190  TYPE(pw_pool_p_type), DIMENSION(:), POINTER :: pw_pools
191  TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
192  TYPE(pw_r3d_rs_type) :: wf_r
193  TYPE(r_ldos_p_type), DIMENSION(:), POINTER :: r_ldos_p
194  TYPE(section_vals_type), POINTER :: ldos_section
195 
196  NULLIFY (logger)
197  logger => cp_get_default_logger()
198  ionode = logger%para_env%is_source()
199  should_output = btest(cp_print_key_should_output(logger%iter_info, dft_section, &
200  "PRINT%PDOS"), cp_p_file)
201  output_unit = cp_logger_get_default_io_unit(logger)
202 
203  spin(1) = "ALPHA"
204  spin(2) = "BETA"
205  IF ((.NOT. should_output)) RETURN
206 
207  NULLIFY (context, s_matrix, orb_basis_set, para_env, pdos_array)
208  NULLIFY (eigenvalues, fm_struct_tmp, mo_coeff, vecbuffer)
209  NULLIFY (ldos_section, list, cell, pw_env, auxbas_pw_pool, evals_virt)
210  NULLIFY (occupation_numbers, ldos_p, r_ldos_p, dft_control, occupation_numbers)
211 
212  CALL timeset(routinen, handle)
213  iterstep = logger%iter_info%iteration(logger%iter_info%n_rlevel)
214 
215  IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T3,A,T61,I10))') &
216  " Calculate PDOS at iteration step ", iterstep
217  CALL get_qs_env(qs_env=qs_env, &
218  matrix_s=s_matrix)
219 
220  CALL get_atomic_kind_set(atomic_kind_set, natom=natom)
221  CALL get_qs_kind_set(qs_kind_set, nsgf=nsgf, maxlgto=maxlgto)
222  nkind = SIZE(atomic_kind_set)
223 
224  CALL get_mo_set(mo_set=mo_set, mo_coeff=mo_coeff, homo=homo, nao=nao, nmo=nmo, &
225  mu=e_fermi)
226  CALL cp_fm_get_info(mo_coeff, &
227  context=context, para_env=para_env, &
228  nrow_global=nrow_global, &
229  ncol_global=ncol_global)
230 
231  CALL section_vals_val_get(dft_section, "PRINT%PDOS%OUT_EACH_MO", i_val=out_each)
232  IF (out_each == -1) out_each = nao + 1
233  CALL section_vals_val_get(dft_section, "PRINT%PDOS%NLUMO", i_val=nlumo)
234  IF (nlumo == -1) nlumo = nao - homo
235  do_virt = (nlumo > (nmo - homo))
236  nvirt = nlumo - (nmo - homo)
237  ! Generate virtual orbitals
238  IF (do_virt) THEN
239  IF (PRESENT(ispin)) THEN
240  my_spin = ispin
241  ELSE
242  my_spin = 1
243  END IF
244 
245  CALL generate_virtual_mo(qs_env, mo_set, evals_virt, mo_virt, nvirt, ispin=my_spin)
246  ELSE
247  NULLIFY (evals_virt)
248  nvirt = 0
249  END IF
250 
251  calc_matsh = .true.
252  IF (PRESENT(external_matrix_shalf)) calc_matsh = .false.
253 
254  ! Create S^1/2 : from sparse to full matrix, if no external available
255  IF (calc_matsh) THEN
256  NULLIFY (matrix_shalf)
257  CALL cp_fm_struct_create(fm_struct_tmp, para_env=para_env, context=context, &
258  nrow_global=nrow_global, ncol_global=nrow_global)
259  ALLOCATE (matrix_shalf)
260  CALL cp_fm_create(matrix_shalf, fm_struct_tmp, name="matrix_shalf")
261  CALL cp_fm_create(matrix_work, fm_struct_tmp, name="matrix_work")
262  CALL cp_fm_struct_release(fm_struct_tmp)
263  CALL copy_dbcsr_to_fm(s_matrix(1)%matrix, matrix_shalf)
264  CALL cp_fm_power(matrix_shalf, matrix_work, 0.5_dp, epsilon(0.0_dp), n_dependent)
265  CALL cp_fm_release(matrix_work)
266  ELSE
267  matrix_shalf => external_matrix_shalf
268  END IF
269 
270  ! Multiply S^(1/2) time the mOS coefficients to get orthonormalized MOS
271  CALL cp_fm_struct_create(fm_struct_tmp, para_env=para_env, context=context, &
272  nrow_global=nrow_global, ncol_global=ncol_global)
273  CALL cp_fm_create(matrix_shalfc, fm_struct_tmp, name="matrix_shalfc")
274  CALL parallel_gemm("N", "N", nrow_global, ncol_global, nrow_global, &
275  1.0_dp, matrix_shalf, mo_coeff, 0.0_dp, matrix_shalfc)
276  CALL cp_fm_struct_release(fm_struct_tmp)
277 
278  IF (do_virt) THEN
279  IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T3,A,T14,I10,T27,A))') &
280  " Compute ", nvirt, " additional unoccupied KS orbitals"
281  CALL cp_fm_struct_create(fm_struct_tmp, para_env=para_env, context=context, &
282  nrow_global=nrow_global, ncol_global=nvirt)
283  CALL cp_fm_create(matrix_work, fm_struct_tmp, name="matrix_shalfc")
284  CALL parallel_gemm("N", "N", nrow_global, nvirt, nrow_global, &
285  1.0_dp, matrix_shalf, mo_virt, 0.0_dp, matrix_work)
286  CALL cp_fm_struct_release(fm_struct_tmp)
287  END IF
288 
289  IF (calc_matsh) THEN
290  CALL cp_fm_release(matrix_shalf)
291  DEALLOCATE (matrix_shalf)
292  END IF
293  ! Array to store the PDOS per kind and angular momentum
294  do_ldos = .false.
295  ldos_section => section_vals_get_subs_vals(dft_section, "PRINT%PDOS%LDOS")
296 
297  CALL section_vals_get(ldos_section, n_repetition=nldos)
298  IF (nldos > 0) THEN
299  IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T3,A,T61,I10))') &
300  " Prepare the list of atoms for LDOS. Number of lists: ", nldos
301  do_ldos = .true.
302  ALLOCATE (ldos_p(nldos))
303  ALLOCATE (ldos_index(nldos))
304  DO ildos = 1, nldos
305  WRITE (ldos_index(ildos), '(I0)') ildos
306  ALLOCATE (ldos_p(ildos)%ldos)
307  NULLIFY (ldos_p(ildos)%ldos%pdos_array)
308  NULLIFY (ldos_p(ildos)%ldos%list_index)
309 
310  CALL section_vals_val_get(ldos_section, "LIST", i_rep_section=ildos, n_rep_val=n_rep)
311  IF (n_rep > 0) THEN
312  ldos_p(ildos)%ldos%nlist = 0
313  DO ir = 1, n_rep
314  NULLIFY (list)
315  CALL section_vals_val_get(ldos_section, "LIST", i_rep_section=ildos, i_rep_val=ir, &
316  i_vals=list)
317  IF (ASSOCIATED(list)) THEN
318  CALL reallocate(ldos_p(ildos)%ldos%list_index, 1, ldos_p(ildos)%ldos%nlist + SIZE(list))
319  DO i = 1, SIZE(list)
320  ldos_p(ildos)%ldos%list_index(i + ldos_p(ildos)%ldos%nlist) = list(i)
321  END DO
322  ldos_p(ildos)%ldos%nlist = ldos_p(ildos)%ldos%nlist + SIZE(list)
323  END IF
324  END DO
325  ELSE
326  ! stop, LDOS without list of atoms is not implemented
327  END IF
328 
329  IF (output_unit > 0) WRITE (unit=output_unit, fmt='((T10,A,T18,I6,T25,A,T36,I10,A))') &
330  " List ", ildos, " contains ", ldos_p(ildos)%ldos%nlist, " atoms"
331  CALL section_vals_val_get(ldos_section, "COMPONENTS", i_rep_section=ildos, &
332  l_val=ldos_p(ildos)%ldos%separate_components)
333  IF (ldos_p(ildos)%ldos%separate_components) THEN
334  ALLOCATE (ldos_p(ildos)%ldos%pdos_array(nsoset(maxlgto), nmo + nvirt))
335  ELSE
336  ALLOCATE (ldos_p(ildos)%ldos%pdos_array(0:maxlgto, nmo + nvirt))
337  END IF
338  ldos_p(ildos)%ldos%pdos_array = 0.0_dp
339  ldos_p(ildos)%ldos%maxl = -1
340 
341  END DO
342  END IF
343 
344  do_r_ldos = .false.
345  ldos_section => section_vals_get_subs_vals(dft_section, "PRINT%PDOS%R_LDOS")
346  CALL section_vals_get(ldos_section, n_repetition=n_r_ldos)
347  IF (n_r_ldos > 0) THEN
348  do_r_ldos = .true.
349  IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T3,A,T61,I10))') &
350  " Prepare the list of points for R_LDOS. Number of lists: ", n_r_ldos
351  ALLOCATE (r_ldos_p(n_r_ldos))
352  ALLOCATE (r_ldos_index(n_r_ldos))
353  CALL get_qs_env(qs_env=qs_env, &
354  cell=cell, &
355  dft_control=dft_control, &
356  pw_env=pw_env)
357  CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool, &
358  pw_pools=pw_pools)
359 
360  CALL auxbas_pw_pool%create_pw(wf_r)
361  CALL auxbas_pw_pool%create_pw(wf_g)
362  ALLOCATE (read_r(4, n_r_ldos))
363  DO ildos = 1, n_r_ldos
364  WRITE (r_ldos_index(ildos), '(I0)') ildos
365  ALLOCATE (r_ldos_p(ildos)%ldos)
366  NULLIFY (r_ldos_p(ildos)%ldos%pdos_array)
367  NULLIFY (r_ldos_p(ildos)%ldos%list_index)
368 
369  CALL section_vals_val_get(ldos_section, "LIST", i_rep_section=ildos, n_rep_val=n_rep)
370  IF (n_rep > 0) THEN
371  r_ldos_p(ildos)%ldos%nlist = 0
372  DO ir = 1, n_rep
373  NULLIFY (list)
374  CALL section_vals_val_get(ldos_section, "LIST", i_rep_section=ildos, i_rep_val=ir, &
375  i_vals=list)
376  IF (ASSOCIATED(list)) THEN
377  CALL reallocate(r_ldos_p(ildos)%ldos%list_index, 1, r_ldos_p(ildos)%ldos%nlist + SIZE(list))
378  DO i = 1, SIZE(list)
379  r_ldos_p(ildos)%ldos%list_index(i + r_ldos_p(ildos)%ldos%nlist) = list(i)
380  END DO
381  r_ldos_p(ildos)%ldos%nlist = r_ldos_p(ildos)%ldos%nlist + SIZE(list)
382  END IF
383  END DO
384  ELSE
385  ! stop, LDOS without list of atoms is not implemented
386  END IF
387 
388  ALLOCATE (r_ldos_p(ildos)%ldos%pdos_array(nmo + nvirt))
389  r_ldos_p(ildos)%ldos%pdos_array = 0.0_dp
390  read_r(1:3, ildos) = .false.
391  CALL section_vals_val_get(ldos_section, "XRANGE", i_rep_section=ildos, explicit=read_r(1, ildos))
392  IF (read_r(1, ildos)) THEN
393  CALL section_vals_val_get(ldos_section, "XRANGE", i_rep_section=ildos, r_vals= &
394  r_ldos_p(ildos)%ldos%x_range)
395  ELSE
396  ALLOCATE (r_ldos_p(ildos)%ldos%x_range(2))
397  r_ldos_p(ildos)%ldos%x_range = 0.0_dp
398  END IF
399  CALL section_vals_val_get(ldos_section, "YRANGE", i_rep_section=ildos, explicit=read_r(2, ildos))
400  IF (read_r(2, ildos)) THEN
401  CALL section_vals_val_get(ldos_section, "YRANGE", i_rep_section=ildos, r_vals= &
402  r_ldos_p(ildos)%ldos%y_range)
403  ELSE
404  ALLOCATE (r_ldos_p(ildos)%ldos%y_range(2))
405  r_ldos_p(ildos)%ldos%y_range = 0.0_dp
406  END IF
407  CALL section_vals_val_get(ldos_section, "ZRANGE", i_rep_section=ildos, explicit=read_r(3, ildos))
408  IF (read_r(3, ildos)) THEN
409  CALL section_vals_val_get(ldos_section, "ZRANGE", i_rep_section=ildos, r_vals= &
410  r_ldos_p(ildos)%ldos%z_range)
411  ELSE
412  ALLOCATE (r_ldos_p(ildos)%ldos%z_range(2))
413  r_ldos_p(ildos)%ldos%z_range = 0.0_dp
414  END IF
415 
416  CALL section_vals_val_get(ldos_section, "ERANGE", i_rep_section=ildos, explicit=read_r(4, ildos))
417  IF (read_r(4, ildos)) THEN
418  CALL section_vals_val_get(ldos_section, "ERANGE", i_rep_section=ildos, &
419  r_vals=r_ldos_p(ildos)%ldos%eval_range)
420  ELSE
421  ALLOCATE (r_ldos_p(ildos)%ldos%eval_range(2))
422  r_ldos_p(ildos)%ldos%eval_range(1) = -huge(0.0_dp)
423  r_ldos_p(ildos)%ldos%eval_range(2) = +huge(0.0_dp)
424  END IF
425 
426  bo => wf_r%pw_grid%bounds_local
427  dh = wf_r%pw_grid%dh
428  dvol = wf_r%pw_grid%dvol
429  np_tot = wf_r%pw_grid%npts(1)*wf_r%pw_grid%npts(2)*wf_r%pw_grid%npts(3)
430  ALLOCATE (r_ldos_p(ildos)%ldos%index_grid_local(3, np_tot))
431 
432  r_ldos_p(ildos)%ldos%npoints = 0
433  DO jz = bo(1, 3), bo(2, 3)
434  DO jy = bo(1, 2), bo(2, 2)
435  DO jx = bo(1, 1), bo(2, 1)
436  !compute the position of the grid point
437  i = jx - wf_r%pw_grid%bounds(1, 1)
438  j = jy - wf_r%pw_grid%bounds(1, 2)
439  k = jz - wf_r%pw_grid%bounds(1, 3)
440  r(3) = k*dh(3, 3) + j*dh(3, 2) + i*dh(3, 1)
441  r(2) = k*dh(2, 3) + j*dh(2, 2) + i*dh(2, 1)
442  r(1) = k*dh(1, 3) + j*dh(1, 2) + i*dh(1, 1)
443 
444  DO il = 1, r_ldos_p(ildos)%ldos%nlist
445  iatom = r_ldos_p(ildos)%ldos%list_index(il)
446  ratom = particle_set(iatom)%r
447  r_vec = pbc(ratom, r, cell)
448  IF (cell%orthorhombic) THEN
449  IF (cell%perd(1) == 0) r_vec(1) = modulo(r_vec(1), cell%hmat(1, 1))
450  IF (cell%perd(2) == 0) r_vec(2) = modulo(r_vec(2), cell%hmat(2, 2))
451  IF (cell%perd(3) == 0) r_vec(3) = modulo(r_vec(3), cell%hmat(3, 3))
452  END IF
453 
454  in_x = 0
455  in_y = 0
456  in_z = 0
457  IF (r_ldos_p(ildos)%ldos%x_range(1) /= 0.0_dp) THEN
458  IF (r_vec(1) > r_ldos_p(ildos)%ldos%x_range(1) .AND. &
459  r_vec(1) < r_ldos_p(ildos)%ldos%x_range(2)) THEN
460  in_x = 1
461  END IF
462  ELSE
463  in_x = 1
464  END IF
465  IF (r_ldos_p(ildos)%ldos%y_range(1) /= 0.0_dp) THEN
466  IF (r_vec(2) > r_ldos_p(ildos)%ldos%y_range(1) .AND. &
467  r_vec(2) < r_ldos_p(ildos)%ldos%y_range(2)) THEN
468  in_y = 1
469  END IF
470  ELSE
471  in_y = 1
472  END IF
473  IF (r_ldos_p(ildos)%ldos%z_range(1) /= 0.0_dp) THEN
474  IF (r_vec(3) > r_ldos_p(ildos)%ldos%z_range(1) .AND. &
475  r_vec(3) < r_ldos_p(ildos)%ldos%z_range(2)) THEN
476  in_z = 1
477  END IF
478  ELSE
479  in_z = 1
480  END IF
481  IF (in_x*in_y*in_z > 0) THEN
482  r_ldos_p(ildos)%ldos%npoints = r_ldos_p(ildos)%ldos%npoints + 1
483  r_ldos_p(ildos)%ldos%index_grid_local(1, r_ldos_p(ildos)%ldos%npoints) = jx
484  r_ldos_p(ildos)%ldos%index_grid_local(2, r_ldos_p(ildos)%ldos%npoints) = jy
485  r_ldos_p(ildos)%ldos%index_grid_local(3, r_ldos_p(ildos)%ldos%npoints) = jz
486  EXIT
487  END IF
488  END DO
489  END DO
490  END DO
491  END DO
492  CALL reallocate(r_ldos_p(ildos)%ldos%index_grid_local, 1, 3, 1, r_ldos_p(ildos)%ldos%npoints)
493  npoints = r_ldos_p(ildos)%ldos%npoints
494  CALL para_env%sum(npoints)
495  IF (output_unit > 0) WRITE (unit=output_unit, fmt='((T10,A,T18,I6,T25,A,T36,I10,A))') &
496  " List ", ildos, " contains ", npoints, " grid points"
497  END DO
498  END IF
499 
500  CALL section_vals_val_get(dft_section, "PRINT%PDOS%COMPONENTS", l_val=separate_components)
501  IF (separate_components) THEN
502  ALLOCATE (pdos_array(nsoset(maxlgto), nkind, nmo + nvirt))
503  ELSE
504  ALLOCATE (pdos_array(0:maxlgto, nkind, nmo + nvirt))
505  END IF
506  IF (do_virt) THEN
507  ALLOCATE (eigenvalues(nmo + nvirt))
508  eigenvalues(1:nmo) = mo_set%eigenvalues(1:nmo)
509  eigenvalues(nmo + 1:nmo + nvirt) = evals_virt(1:nvirt)
510  ALLOCATE (occupation_numbers(nmo + nvirt))
511  occupation_numbers(:) = 0.0_dp
512  occupation_numbers(1:nmo) = mo_set%occupation_numbers(1:nmo)
513  ELSE
514  eigenvalues => mo_set%eigenvalues
515  occupation_numbers => mo_set%occupation_numbers
516  END IF
517 
518  pdos_array = 0.0_dp
519  nao = mo_set%nao
520  ALLOCATE (vecbuffer(1, nao))
521  vecbuffer = 0.0_dp
522  ALLOCATE (firstrow(natom))
523  firstrow = 0
524 
525  !Adjust energy range for r_ldos
526  DO ildos = 1, n_r_ldos
527  IF (eigenvalues(1) > r_ldos_p(ildos)%ldos%eval_range(1)) &
528  r_ldos_p(ildos)%ldos%eval_range(1) = eigenvalues(1)
529  IF (eigenvalues(nmo + nvirt) < r_ldos_p(ildos)%ldos%eval_range(2)) &
530  r_ldos_p(ildos)%ldos%eval_range(2) = eigenvalues(nmo + nvirt)
531  END DO
532 
533  IF (output_unit > 0) WRITE (unit=output_unit, fmt='(/,(T15,A))') &
534  "---- PDOS: start iteration on the KS states --- "
535 
536  DO imo = 1, nmo + nvirt
537 
538  IF (output_unit > 0 .AND. mod(imo, out_each) == 0) WRITE (unit=output_unit, fmt='((T20,A,I10))') &
539  " KS state index : ", imo
540  ! Extract the eigenvector from the distributed full matrix
541  IF (imo > nmo) THEN
542  CALL cp_fm_get_submatrix(matrix_work, vecbuffer, 1, imo - nmo, &
543  nao, 1, transpose=.true.)
544  ELSE
545  CALL cp_fm_get_submatrix(matrix_shalfc, vecbuffer, 1, imo, &
546  nao, 1, transpose=.true.)
547  END IF
548 
549  ! Calculate the pdos for all the kinds
550  irow = 1
551  DO iatom = 1, natom
552  firstrow(iatom) = irow
553  NULLIFY (orb_basis_set)
554  CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
555  CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
556  CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
557  nset=nset, &
558  nshell=nshell, &
559  l=l, maxl=maxl)
560  IF (separate_components) THEN
561  isgf = 1
562  DO iset = 1, nset
563  DO ishell = 1, nshell(iset)
564  lshell = l(ishell, iset)
565  DO iso = 1, nso(lshell)
566  lcomponent = nsoset(lshell - 1) + iso
567  pdos_array(lcomponent, ikind, imo) = &
568  pdos_array(lcomponent, ikind, imo) + &
569  vecbuffer(1, irow)*vecbuffer(1, irow)
570  irow = irow + 1
571  END DO ! iso
572  END DO ! ishell
573  END DO ! iset
574  ELSE
575  isgf = 1
576  DO iset = 1, nset
577  DO ishell = 1, nshell(iset)
578  lshell = l(ishell, iset)
579  DO iso = 1, nso(lshell)
580  pdos_array(lshell, ikind, imo) = &
581  pdos_array(lshell, ikind, imo) + &
582  vecbuffer(1, irow)*vecbuffer(1, irow)
583  irow = irow + 1
584  END DO ! iso
585  END DO ! ishell
586  END DO ! iset
587  END IF
588  END DO ! iatom
589 
590  ! Calculate the pdos for all the lists
591  DO ildos = 1, nldos
592  DO il = 1, ldos_p(ildos)%ldos%nlist
593  iatom = ldos_p(ildos)%ldos%list_index(il)
594 
595  irow = firstrow(iatom)
596  NULLIFY (orb_basis_set)
597  CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
598  CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
599 
600  CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
601  nset=nset, &
602  nshell=nshell, &
603  l=l, maxl=maxl)
604  ldos_p(ildos)%ldos%maxl = max(ldos_p(ildos)%ldos%maxl, maxl)
605  IF (ldos_p(ildos)%ldos%separate_components) THEN
606  isgf = 1
607  DO iset = 1, nset
608  DO ishell = 1, nshell(iset)
609  lshell = l(ishell, iset)
610  DO iso = 1, nso(lshell)
611  lcomponent = nsoset(lshell - 1) + iso
612  ldos_p(ildos)%ldos%pdos_array(lcomponent, imo) = &
613  ldos_p(ildos)%ldos%pdos_array(lcomponent, imo) + &
614  vecbuffer(1, irow)*vecbuffer(1, irow)
615  irow = irow + 1
616  END DO ! iso
617  END DO ! ishell
618  END DO ! iset
619  ELSE
620  isgf = 1
621  DO iset = 1, nset
622  DO ishell = 1, nshell(iset)
623  lshell = l(ishell, iset)
624  DO iso = 1, nso(lshell)
625  ldos_p(ildos)%ldos%pdos_array(lshell, imo) = &
626  ldos_p(ildos)%ldos%pdos_array(lshell, imo) + &
627  vecbuffer(1, irow)*vecbuffer(1, irow)
628  irow = irow + 1
629  END DO ! iso
630  END DO ! ishell
631  END DO ! iset
632  END IF
633  END DO !il
634  END DO !ildos
635 
636  ! Calculate the DOS projected in a given volume in real space
637  DO ildos = 1, n_r_ldos
638  IF (r_ldos_p(ildos)%ldos%eval_range(1) <= eigenvalues(imo) .AND. &
639  r_ldos_p(ildos)%ldos%eval_range(2) >= eigenvalues(imo)) THEN
640 
641  IF (imo > nmo) THEN
642  CALL calculate_wavefunction(mo_virt, imo - nmo, &
643  wf_r, wf_g, atomic_kind_set, qs_kind_set, cell, dft_control, particle_set, &
644  pw_env)
645  ELSE
646  CALL calculate_wavefunction(mo_coeff, imo, &
647  wf_r, wf_g, atomic_kind_set, qs_kind_set, cell, dft_control, particle_set, &
648  pw_env)
649  END IF
650  r_ldos_p(ildos)%ldos%pdos_array(imo) = 0.0_dp
651  DO il = 1, r_ldos_p(ildos)%ldos%npoints
652  j = j + 1
653  jx = r_ldos_p(ildos)%ldos%index_grid_local(1, il)
654  jy = r_ldos_p(ildos)%ldos%index_grid_local(2, il)
655  jz = r_ldos_p(ildos)%ldos%index_grid_local(3, il)
656  r_ldos_p(ildos)%ldos%pdos_array(imo) = r_ldos_p(ildos)%ldos%pdos_array(imo) + &
657  wf_r%array(jx, jy, jz)*wf_r%array(jx, jy, jz)
658  END DO
659  r_ldos_p(ildos)%ldos%pdos_array(imo) = r_ldos_p(ildos)%ldos%pdos_array(imo)*dvol
660  END IF
661  END DO
662  END DO ! imo
663 
664  CALL cp_fm_release(matrix_shalfc)
665  DEALLOCATE (vecbuffer)
666 
667  CALL section_vals_val_get(dft_section, "PRINT%PDOS%APPEND", l_val=append)
668  IF (append .AND. iterstep > 1) THEN
669  my_pos = "APPEND"
670  ELSE
671  my_pos = "REWIND"
672  END IF
673  my_act = "WRITE"
674  DO ikind = 1, nkind
675 
676  NULLIFY (orb_basis_set)
677  CALL get_atomic_kind(atomic_kind_set(ikind), name=kind_name)
678  CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
679  CALL get_gto_basis_set(gto_basis_set=orb_basis_set, maxl=maxl)
680 
681  ! basis none has no associated maxl, and no pdos
682  IF (maxl < 0) cycle
683 
684  IF (PRESENT(ispin)) THEN
685  IF (PRESENT(xas_mittle)) THEN
686  my_mittle = trim(xas_mittle)//trim(spin(ispin))//"_k"//trim(adjustl(cp_to_string(ikind)))
687  ELSE
688  my_mittle = trim(spin(ispin))//"_k"//trim(adjustl(cp_to_string(ikind)))
689  END IF
690  my_spin = ispin
691  ELSE
692  my_mittle = "k"//trim(adjustl(cp_to_string(ikind)))
693  my_spin = 1
694  END IF
695 
696  iw = cp_print_key_unit_nr(logger, dft_section, "PRINT%PDOS", &
697  extension=".pdos", file_position=my_pos, file_action=my_act, &
698  file_form="FORMATTED", middle_name=trim(my_mittle))
699  IF (iw > 0) THEN
700 
701  fmtstr1 = "(I8,2X,2F16.6, (2X,F16.8))"
702  fmtstr2 = "(A42, (10X,A8))"
703  IF (separate_components) THEN
704  WRITE (unit=fmtstr1(15:16), fmt="(I2)") nsoset(maxl)
705  WRITE (unit=fmtstr2(6:7), fmt="(I2)") nsoset(maxl)
706  ELSE
707  WRITE (unit=fmtstr1(15:16), fmt="(I2)") maxl + 1
708  WRITE (unit=fmtstr2(6:7), fmt="(I2)") maxl + 1
709  END IF
710 
711  WRITE (unit=iw, fmt="(A,I0,A,F12.6,A)") &
712  "# Projected DOS for atomic kind "//trim(kind_name)//" at iteration step i = ", &
713  iterstep, ", E(Fermi) = ", e_fermi, " a.u."
714  IF (separate_components) THEN
715  ALLOCATE (tmp_str(0:0, 0:maxl, -maxl:maxl))
716  tmp_str = ""
717  DO j = 0, maxl
718  DO i = -j, j
719  tmp_str(0, j, i) = sgf_symbol(0, j, i)
720  END DO
721  END DO
722 
723  WRITE (unit=iw, fmt=fmtstr2) &
724  "# MO Eigenvalue [a.u.] Occupation", &
725  ((trim(tmp_str(0, il, im)), im=-il, il), il=0, maxl)
726  DO imo = 1, nmo + nvirt
727  WRITE (unit=iw, fmt=fmtstr1) imo, eigenvalues(imo), occupation_numbers(imo), &
728  (pdos_array(lshell, ikind, imo), lshell=1, nsoset(maxl))
729  END DO
730  DEALLOCATE (tmp_str)
731  ELSE
732  WRITE (unit=iw, fmt=fmtstr2) &
733  "# MO Eigenvalue [a.u.] Occupation", &
734  (trim(l_sym(il)), il=0, maxl)
735  DO imo = 1, nmo + nvirt
736  WRITE (unit=iw, fmt=fmtstr1) imo, eigenvalues(imo), occupation_numbers(imo), &
737  (pdos_array(lshell, ikind, imo), lshell=0, maxl)
738  END DO
739  END IF
740  END IF
741  CALL cp_print_key_finished_output(iw, logger, dft_section, &
742  "PRINT%PDOS")
743 
744  END DO ! ikind
745 
746  ! write the pdos for the lists, each ona different file,
747  ! the filenames are indexed with the list number
748  DO ildos = 1, nldos
749  ! basis none has no associated maxl, and no pdos
750  IF (ldos_p(ildos)%ldos%maxl > 0) THEN
751 
752  IF (PRESENT(ispin)) THEN
753  IF (PRESENT(xas_mittle)) THEN
754  my_mittle = trim(xas_mittle)//trim(spin(ispin))//"_list"//trim(ldos_index(ildos))
755  ELSE
756  my_mittle = trim(spin(ispin))//"_list"//trim(ldos_index(ildos))
757  END IF
758  my_spin = ispin
759  ELSE
760  my_mittle = "list"//trim(ldos_index(ildos))
761  my_spin = 1
762  END IF
763 
764  iw = cp_print_key_unit_nr(logger, dft_section, "PRINT%PDOS", &
765  extension=".pdos", file_position=my_pos, file_action=my_act, &
766  file_form="FORMATTED", middle_name=trim(my_mittle))
767  IF (iw > 0) THEN
768 
769  fmtstr1 = "(I8,2X,2F16.6, (2X,F16.8))"
770  fmtstr2 = "(A42, (10X,A8))"
771  IF (ldos_p(ildos)%ldos%separate_components) THEN
772  WRITE (unit=fmtstr1(15:16), fmt="(I2)") nsoset(ldos_p(ildos)%ldos%maxl)
773  WRITE (unit=fmtstr2(6:7), fmt="(I2)") nsoset(ldos_p(ildos)%ldos%maxl)
774  ELSE
775  WRITE (unit=fmtstr1(15:16), fmt="(I2)") ldos_p(ildos)%ldos%maxl + 1
776  WRITE (unit=fmtstr2(6:7), fmt="(I2)") ldos_p(ildos)%ldos%maxl + 1
777  END IF
778 
779  WRITE (unit=iw, fmt="(A,I0,A,I0,A,I0,A,F12.6,A)") &
780  "# Projected DOS for list ", ildos, " of ", ldos_p(ildos)%ldos%nlist, &
781  " atoms, at iteration step i = ", iterstep, &
782  ", E(Fermi) = ", e_fermi, " a.u."
783  IF (ldos_p(ildos)%ldos%separate_components) THEN
784  ALLOCATE (tmp_str(0:0, 0:ldos_p(ildos)%ldos%maxl, -ldos_p(ildos)%ldos%maxl:ldos_p(ildos)%ldos%maxl))
785  tmp_str = ""
786  DO j = 0, ldos_p(ildos)%ldos%maxl
787  DO i = -j, j
788  tmp_str(0, j, i) = sgf_symbol(0, j, i)
789  END DO
790  END DO
791 
792  WRITE (unit=iw, fmt=fmtstr2) &
793  "# MO Eigenvalue [a.u.] Occupation", &
794  ((trim(tmp_str(0, il, im)), im=-il, il), il=0, ldos_p(ildos)%ldos%maxl)
795  DO imo = 1, nmo + nvirt
796  WRITE (unit=iw, fmt=fmtstr1) imo, eigenvalues(imo), occupation_numbers(imo), &
797  (ldos_p(ildos)%ldos%pdos_array(lshell, imo), lshell=1, nsoset(ldos_p(ildos)%ldos%maxl))
798  END DO
799  DEALLOCATE (tmp_str)
800  ELSE
801  WRITE (unit=iw, fmt=fmtstr2) &
802  "# MO Eigenvalue [a.u.] Occupation", &
803  (trim(l_sym(il)), il=0, ldos_p(ildos)%ldos%maxl)
804  DO imo = 1, nmo + nvirt
805  WRITE (unit=iw, fmt=fmtstr1) imo, eigenvalues(imo), occupation_numbers(imo), &
806  (ldos_p(ildos)%ldos%pdos_array(lshell, imo), lshell=0, ldos_p(ildos)%ldos%maxl)
807  END DO
808  END IF
809  END IF
810  CALL cp_print_key_finished_output(iw, logger, dft_section, &
811  "PRINT%PDOS")
812  END IF ! maxl>0
813  END DO ! ildos
814 
815  ! write the pdos for the lists, each ona different file,
816  ! the filenames are indexed with the list number
817  DO ildos = 1, n_r_ldos
818 
819  npoints = r_ldos_p(ildos)%ldos%npoints
820  CALL para_env%sum(npoints)
821  CALL para_env%sum(np_tot)
822  CALL para_env%sum(r_ldos_p(ildos)%ldos%pdos_array)
823  IF (PRESENT(ispin)) THEN
824  IF (PRESENT(xas_mittle)) THEN
825  my_mittle = trim(xas_mittle)//trim(spin(ispin))//"_r_list"//trim(r_ldos_index(ildos))
826  ELSE
827  my_mittle = trim(spin(ispin))//"_r_list"//trim(r_ldos_index(ildos))
828  END IF
829  my_spin = ispin
830  ELSE
831  my_mittle = "r_list"//trim(r_ldos_index(ildos))
832  my_spin = 1
833  END IF
834 
835  iw = cp_print_key_unit_nr(logger, dft_section, "PRINT%PDOS", &
836  extension=".pdos", file_position=my_pos, file_action=my_act, &
837  file_form="FORMATTED", middle_name=trim(my_mittle))
838  IF (iw > 0) THEN
839  fmtstr1 = "(I8,2X,2F16.6, (2X,F16.8))"
840  fmtstr2 = "(A42, (10X,A8))"
841 
842  WRITE (unit=iw, fmt="(A,I0,A,F12.6,F12.6,A,F12.6,A)") &
843  "# Projected DOS in real space, using ", npoints, &
844  " points of the grid, and eval in the range", r_ldos_p(ildos)%ldos%eval_range(1:2), &
845  " Hartree, E(Fermi) = ", e_fermi, " a.u."
846  WRITE (unit=iw, fmt="(A)") &
847  "# MO Eigenvalue [a.u.] Occupation LDOS"
848  DO imo = 1, nmo + nvirt
849  IF (r_ldos_p(ildos)%ldos%eval_range(1) <= eigenvalues(imo) .AND. &
850  r_ldos_p(ildos)%ldos%eval_range(2) >= eigenvalues(imo)) THEN
851  WRITE (unit=iw, fmt="(I8,2X,2F16.6,E20.10,E20.10)") imo, eigenvalues(imo), occupation_numbers(imo), &
852  r_ldos_p(ildos)%ldos%pdos_array(imo), r_ldos_p(ildos)%ldos%pdos_array(imo)*np_tot
853  END IF
854  END DO
855 
856  END IF
857  CALL cp_print_key_finished_output(iw, logger, dft_section, &
858  "PRINT%PDOS")
859  END DO
860 
861  ! deallocate local variables
862  DEALLOCATE (pdos_array)
863  DEALLOCATE (firstrow)
864  IF (do_ldos) THEN
865  DO ildos = 1, nldos
866  DEALLOCATE (ldos_p(ildos)%ldos%pdos_array)
867  DEALLOCATE (ldos_p(ildos)%ldos%list_index)
868  DEALLOCATE (ldos_p(ildos)%ldos)
869  END DO
870  DEALLOCATE (ldos_p)
871  DEALLOCATE (ldos_index)
872  END IF
873  IF (do_r_ldos) THEN
874  DO ildos = 1, n_r_ldos
875  DEALLOCATE (r_ldos_p(ildos)%ldos%index_grid_local)
876  DEALLOCATE (r_ldos_p(ildos)%ldos%pdos_array)
877  DEALLOCATE (r_ldos_p(ildos)%ldos%list_index)
878  IF (.NOT. read_r(1, ildos)) THEN
879  DEALLOCATE (r_ldos_p(ildos)%ldos%x_range)
880  END IF
881  IF (.NOT. read_r(2, ildos)) THEN
882  DEALLOCATE (r_ldos_p(ildos)%ldos%y_range)
883  END IF
884  IF (.NOT. read_r(3, ildos)) THEN
885  DEALLOCATE (r_ldos_p(ildos)%ldos%z_range)
886  END IF
887  IF (.NOT. read_r(4, ildos)) THEN
888  DEALLOCATE (r_ldos_p(ildos)%ldos%eval_range)
889  END IF
890  DEALLOCATE (r_ldos_p(ildos)%ldos)
891  END DO
892  DEALLOCATE (read_r)
893  DEALLOCATE (r_ldos_p)
894  DEALLOCATE (r_ldos_index)
895  CALL auxbas_pw_pool%give_back_pw(wf_r)
896  CALL auxbas_pw_pool%give_back_pw(wf_g)
897  END IF
898  IF (do_virt) THEN
899  DEALLOCATE (evals_virt)
900  CALL cp_fm_release(mo_virt)
901  CALL cp_fm_release(matrix_work)
902  DEALLOCATE (eigenvalues)
903  DEALLOCATE (occupation_numbers)
904  END IF
905 
906  CALL timestop(handle)
907 
908  END SUBROUTINE calculate_projected_dos
909 
910 ! **************************************************************************************************
911 !> \brief Compute additional virtual states starting from the available MOS
912 !> \param qs_env ...
913 !> \param mo_set ...
914 !> \param evals_virt ...
915 !> \param mo_virt ...
916 !> \param nvirt ...
917 !> \param ispin ...
918 !> \date 08.03.2008
919 !> \par History:
920 !> -
921 !> \par Variables
922 !> -
923 !> -
924 !> \author MI
925 !> \version 1.0
926 ! **************************************************************************************************
927 
928  SUBROUTINE generate_virtual_mo(qs_env, mo_set, evals_virt, mo_virt, &
929  nvirt, ispin)
930 
931  TYPE(qs_environment_type), POINTER :: qs_env
932  TYPE(mo_set_type), INTENT(IN) :: mo_set
933  REAL(kind=dp), DIMENSION(:), POINTER :: evals_virt
934  TYPE(cp_fm_type), INTENT(OUT) :: mo_virt
935  INTEGER, INTENT(IN) :: nvirt, ispin
936 
937  INTEGER :: nmo, nrow_global
938  TYPE(cp_blacs_env_type), POINTER :: context
939  TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
940  TYPE(cp_fm_type), POINTER :: mo_coeff
941  TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: ks_matrix, s_matrix
942  TYPE(mp_para_env_type), POINTER :: para_env
943  TYPE(preconditioner_type), POINTER :: local_preconditioner
944  TYPE(scf_control_type), POINTER :: scf_control
945 
946  NULLIFY (evals_virt)
947  ALLOCATE (evals_virt(nvirt))
948 
949  CALL get_qs_env(qs_env, matrix_ks=ks_matrix, matrix_s=s_matrix, &
950  scf_control=scf_control)
951  CALL get_mo_set(mo_set=mo_set, mo_coeff=mo_coeff, nmo=nmo)
952  CALL cp_fm_get_info(mo_coeff, context=context, para_env=para_env, &
953  nrow_global=nrow_global)
954 
955  CALL cp_fm_struct_create(fm_struct_tmp, para_env=para_env, context=context, &
956  nrow_global=nrow_global, ncol_global=nvirt)
957  CALL cp_fm_create(mo_virt, fm_struct_tmp, name="virtual")
958  CALL cp_fm_struct_release(fm_struct_tmp)
959  CALL cp_fm_init_random(mo_virt, nvirt)
960 
961  NULLIFY (local_preconditioner)
962 
963  CALL ot_eigensolver(matrix_h=ks_matrix(ispin)%matrix, matrix_s=s_matrix(1)%matrix, &
964  matrix_c_fm=mo_virt, matrix_orthogonal_space_fm=mo_coeff, &
965  eps_gradient=scf_control%eps_lumos, &
966  preconditioner=local_preconditioner, &
967  iter_max=scf_control%max_iter_lumos, &
968  size_ortho_space=nmo)
969 
970  CALL calculate_subspace_eigenvalues(mo_virt, ks_matrix(ispin)%matrix, &
971  evals_virt)
972 
973  END SUBROUTINE generate_virtual_mo
974 
975 END MODULE qs_pdos
976 
subroutine pbc(r, r_pbc, s, s_pbc, a, b, c, alpha, beta, gamma, debug, info, pbc0, h, hinv)
...
Definition: dumpdcd.F:1203
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:117
Define the atomic kind types and their sub types.
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.
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.
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)
...
Handles all functions related to the CELL.
Definition: cell_types.F:15
methods related to the blacs parallel environment
Definition: cp_blacs_env.F:15
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
DBCSR operations in CP2K.
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
Definition: cp_fm_diag.F:17
subroutine, public cp_fm_power(matrix, work, exponent, threshold, n_dependent, verbose, eigvals)
...
Definition: cp_fm_diag.F:896
represent the structure of a full matrix
Definition: cp_fm_struct.F:14
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:125
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition: cp_fm_struct.F:320
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
Definition: cp_fm_types.F:1016
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:901
subroutine, public cp_fm_init_random(matrix, ncol, start_col)
fills a matrix with random numbers
Definition: cp_fm_types.F:452
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
Definition: cp_fm_types.F:167
various routines to log and control the output. The idea is that decisions about where to log should ...
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...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
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)
...
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,...
integer, parameter, public cp_p_file
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...
objects that represent the structure of input sections and the data contained in an input section
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
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
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
Defines the basic variable types.
Definition: kinds.F:23
integer, parameter, public dp
Definition: kinds.F:34
integer, parameter, public default_string_length
Definition: kinds.F:57
An array-based list which grows on demand. When the internal array is full, a new array of twice the ...
Definition: list.F:24
Utility routines for the memory handling.
Interface to the message passing library MPI.
Provides Cartesian and spherical orbital pointers and indices.
integer, dimension(:), allocatable, public nsoset
integer, dimension(:), allocatable, public nso
orbital_symbols
character(len=6) function, public sgf_symbol(n, l, m)
Build a spherical orbital symbol (orbital labels for printing).
character(len=1), dimension(0:11), parameter, public l_sym
basic linear algebra operations for full matrixes
Define the data structure for the particle information.
types of preconditioners
computes preconditioners, and implements methods to apply them currently used in qs_ot
container for various plainwaves related things
Definition: pw_env_types.F:14
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
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
Calculate the plane wave density by collocating the primitive Gaussian functions (pgf).
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, external_vector)
maps a given wavefunction on the grid
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, 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_nonbond, sab_almo, sab_kp, sab_kp_nosym, 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, 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, ecoul_1c, rho0_s_rs, rho0_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, 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, rhs)
Get the QUICKSTEP environment.
Define the quickstep kind type and their sub types.
Definition: qs_kind_types.F:23
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, 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_r3d_rs_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, floating, name, element_symbol, pao_basis_size, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
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)
Get attributes of an atomic kind set.
collects routines that perform operations directly related to MOs
Definition: qs_mo_methods.F:14
Definition and initialisation of the mo data type.
Definition: qs_mo_types.F:22
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:397
an eigen-space solver for the generalised symmetric eigenvalue problem for sparse matrices,...
subroutine, public ot_eigensolver(matrix_h, matrix_s, matrix_orthogonal_space_fm, matrix_c_fm, preconditioner, eps_gradient, iter_max, size_ortho_space, silent, ot_settings)
...
Calculation and writing of projected density of states The DOS is computed per angular momentum and p...
Definition: qs_pdos.F:15
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)
Compute and write projected density of states.
Definition: qs_pdos.F:139
parameters that control an scf iteration