44 USE dbcsr_api,
ONLY: dbcsr_p_type
77 #include "./base/base_uses.f90"
83 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_pdos'
91 INTEGER :: maxl, nlist
92 LOGICAL :: separate_components
93 INTEGER,
DIMENSION(:),
POINTER :: list_index
94 REAL(KIND=
dp),
DIMENSION(:, :), &
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
109 TYPE(ldos_type),
POINTER :: ldos
113 TYPE(r_ldos_type),
POINTER :: ldos
114 END TYPE r_ldos_p_type
138 dft_section, ispin, xas_mittle, external_matrix_shalf)
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
151 CHARACTER(len=*),
PARAMETER :: routinen =
'calculate_projected_dos'
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, &
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, &
170 LOGICAL,
DIMENSION(:, :),
POINTER :: read_r
171 REAL(kind=
dp) :: dh(3, 3), dvol, e_fermi, r(3), r_vec(3), &
173 REAL(kind=
dp),
DIMENSION(:),
POINTER :: eigenvalues, evals_virt, &
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
198 ionode = logger%para_env%is_source()
205 IF ((.NOT. should_output))
RETURN
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)
212 CALL timeset(routinen, handle)
213 iterstep = logger%iter_info%iteration(logger%iter_info%n_rlevel)
215 IF (output_unit > 0)
WRITE (unit=output_unit, fmt=
'(/,(T3,A,T61,I10))') &
216 " Calculate PDOS at iteration step ", iterstep
222 nkind =
SIZE(atomic_kind_set)
224 CALL get_mo_set(mo_set=mo_set, mo_coeff=mo_coeff, homo=homo, nao=nao, nmo=nmo, &
227 context=context, para_env=para_env, &
228 nrow_global=nrow_global, &
229 ncol_global=ncol_global)
232 IF (out_each == -1) out_each = nao + 1
234 IF (nlumo == -1) nlumo = nao - homo
235 do_virt = (nlumo > (nmo - homo))
236 nvirt = nlumo - (nmo - homo)
239 IF (
PRESENT(ispin))
THEN
245 CALL generate_virtual_mo(qs_env, mo_set, evals_virt, mo_virt, nvirt, ispin=my_spin)
252 IF (
PRESENT(external_matrix_shalf)) calc_matsh = .false.
256 NULLIFY (matrix_shalf)
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")
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)
267 matrix_shalf => external_matrix_shalf
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)
279 IF (output_unit > 0)
WRITE (unit=output_unit, fmt=
'(/,(T3,A,T14,I10,T27,A))') &
280 " Compute ", nvirt,
" additional unoccupied KS orbitals"
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)
290 CALL cp_fm_release(matrix_shalf)
291 DEALLOCATE (matrix_shalf)
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
302 ALLOCATE (ldos_p(nldos))
303 ALLOCATE (ldos_index(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)
312 ldos_p(ildos)%ldos%nlist = 0
317 IF (
ASSOCIATED(
list))
THEN
318 CALL reallocate(ldos_p(ildos)%ldos%list_index, 1, ldos_p(ildos)%ldos%nlist +
SIZE(
list))
320 ldos_p(ildos)%ldos%list_index(i + ldos_p(ildos)%ldos%nlist) =
list(i)
322 ldos_p(ildos)%ldos%nlist = ldos_p(ildos)%ldos%nlist +
SIZE(
list)
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"
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))
336 ALLOCATE (ldos_p(ildos)%ldos%pdos_array(0:maxlgto, nmo + nvirt))
338 ldos_p(ildos)%ldos%pdos_array = 0.0_dp
339 ldos_p(ildos)%ldos%maxl = -1
347 IF (n_r_ldos > 0)
THEN
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))
355 dft_control=dft_control, &
357 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool, &
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)
371 r_ldos_p(ildos)%ldos%nlist = 0
376 IF (
ASSOCIATED(
list))
THEN
377 CALL reallocate(r_ldos_p(ildos)%ldos%list_index, 1, r_ldos_p(ildos)%ldos%nlist +
SIZE(
list))
379 r_ldos_p(ildos)%ldos%list_index(i + r_ldos_p(ildos)%ldos%nlist) =
list(i)
381 r_ldos_p(ildos)%ldos%nlist = r_ldos_p(ildos)%ldos%nlist +
SIZE(
list)
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
394 r_ldos_p(ildos)%ldos%x_range)
396 ALLOCATE (r_ldos_p(ildos)%ldos%x_range(2))
397 r_ldos_p(ildos)%ldos%x_range = 0.0_dp
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
402 r_ldos_p(ildos)%ldos%y_range)
404 ALLOCATE (r_ldos_p(ildos)%ldos%y_range(2))
405 r_ldos_p(ildos)%ldos%y_range = 0.0_dp
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
410 r_ldos_p(ildos)%ldos%z_range)
412 ALLOCATE (r_ldos_p(ildos)%ldos%z_range(2))
413 r_ldos_p(ildos)%ldos%z_range = 0.0_dp
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
419 r_vals=r_ldos_p(ildos)%ldos%eval_range)
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)
426 bo => wf_r%pw_grid%bounds_local
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))
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)
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)
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))
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
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
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
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
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"
501 IF (separate_components)
THEN
502 ALLOCATE (pdos_array(
nsoset(maxlgto), nkind, nmo + nvirt))
504 ALLOCATE (pdos_array(0:maxlgto, nkind, nmo + nvirt))
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)
514 eigenvalues => mo_set%eigenvalues
515 occupation_numbers => mo_set%occupation_numbers
520 ALLOCATE (vecbuffer(1, nao))
522 ALLOCATE (firstrow(natom))
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)
533 IF (output_unit > 0)
WRITE (unit=output_unit, fmt=
'(/,(T15,A))') &
534 "---- PDOS: start iteration on the KS states --- "
536 DO imo = 1, nmo + nvirt
538 IF (output_unit > 0 .AND. mod(imo, out_each) == 0)
WRITE (unit=output_unit, fmt=
'((T20,A,I10))') &
539 " KS state index : ", imo
543 nao, 1, transpose=.true.)
546 nao, 1, transpose=.true.)
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)
560 IF (separate_components)
THEN
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)
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)
592 DO il = 1, ldos_p(ildos)%ldos%nlist
593 iatom = ldos_p(ildos)%ldos%list_index(il)
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)
604 ldos_p(ildos)%ldos%maxl = max(ldos_p(ildos)%ldos%maxl, maxl)
605 IF (ldos_p(ildos)%ldos%separate_components)
THEN
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)
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)
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
643 wf_r, wf_g, atomic_kind_set, qs_kind_set, cell, dft_control, particle_set, &
647 wf_r, wf_g, atomic_kind_set, qs_kind_set, cell, dft_control, particle_set, &
650 r_ldos_p(ildos)%ldos%pdos_array(imo) = 0.0_dp
651 DO il = 1, r_ldos_p(ildos)%ldos%npoints
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)
659 r_ldos_p(ildos)%ldos%pdos_array(imo) = r_ldos_p(ildos)%ldos%pdos_array(imo)*dvol
664 CALL cp_fm_release(matrix_shalfc)
665 DEALLOCATE (vecbuffer)
668 IF (append .AND. iterstep > 1)
THEN
676 NULLIFY (orb_basis_set)
678 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
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)))
688 my_mittle = trim(spin(ispin))//
"_k"//trim(adjustl(cp_to_string(ikind)))
692 my_mittle =
"k"//trim(adjustl(cp_to_string(ikind)))
697 extension=
".pdos", file_position=my_pos, file_action=my_act, &
698 file_form=
"FORMATTED", middle_name=trim(my_mittle))
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)
707 WRITE (unit=fmtstr1(15:16), fmt=
"(I2)") maxl + 1
708 WRITE (unit=fmtstr2(6:7), fmt=
"(I2)") maxl + 1
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))
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))
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)
750 IF (ldos_p(ildos)%ldos%maxl > 0)
THEN
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))
756 my_mittle = trim(spin(ispin))//
"_list"//trim(ldos_index(ildos))
760 my_mittle =
"list"//trim(ldos_index(ildos))
765 extension=
".pdos", file_position=my_pos, file_action=my_act, &
766 file_form=
"FORMATTED", middle_name=trim(my_mittle))
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)
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
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))
786 DO j = 0, ldos_p(ildos)%ldos%maxl
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))
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)
817 DO ildos = 1, n_r_ldos
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))
827 my_mittle = trim(spin(ispin))//
"_r_list"//trim(r_ldos_index(ildos))
831 my_mittle =
"r_list"//trim(r_ldos_index(ildos))
836 extension=
".pdos", file_position=my_pos, file_action=my_act, &
837 file_form=
"FORMATTED", middle_name=trim(my_mittle))
839 fmtstr1 =
"(I8,2X,2F16.6, (2X,F16.8))"
840 fmtstr2 =
"(A42, (10X,A8))"
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
862 DEALLOCATE (pdos_array)
863 DEALLOCATE (firstrow)
866 DEALLOCATE (ldos_p(ildos)%ldos%pdos_array)
867 DEALLOCATE (ldos_p(ildos)%ldos%list_index)
868 DEALLOCATE (ldos_p(ildos)%ldos)
871 DEALLOCATE (ldos_index)
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)
881 IF (.NOT. read_r(2, ildos))
THEN
882 DEALLOCATE (r_ldos_p(ildos)%ldos%y_range)
884 IF (.NOT. read_r(3, ildos))
THEN
885 DEALLOCATE (r_ldos_p(ildos)%ldos%z_range)
887 IF (.NOT. read_r(4, ildos))
THEN
888 DEALLOCATE (r_ldos_p(ildos)%ldos%eval_range)
890 DEALLOCATE (r_ldos_p(ildos)%ldos)
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)
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)
906 CALL timestop(handle)
928 SUBROUTINE generate_virtual_mo(qs_env, mo_set, evals_virt, mo_virt, &
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
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
947 ALLOCATE (evals_virt(nvirt))
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)
956 nrow_global=nrow_global, ncol_global=nvirt)
957 CALL cp_fm_create(mo_virt, fm_struct_tmp, name=
"virtual")
961 NULLIFY (local_preconditioner)
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, &
967 iter_max=scf_control%max_iter_lumos, &
968 size_ortho_space=nmo)
970 CALL calculate_subspace_eigenvalues(mo_virt, ks_matrix(ispin)%matrix, &
973 END SUBROUTINE generate_virtual_mo
subroutine pbc(r, r_pbc, s, s_pbc, a, b, c, alpha, beta, gamma, debug, info, pbc0, h, hinv)
...
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
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.
methods related to the blacs parallel environment
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...
subroutine, public cp_fm_power(matrix, work, exponent, threshold, n_dependent, verbose, eigvals)
...
represent the structure of a full matrix
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
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
represent a full matrix distributed on many processors
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
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,...
subroutine, public cp_fm_init_random(matrix, ncol, start_col)
fills a matrix with random numbers
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
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...
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
An array-based list which grows on demand. When the internal array is full, a new array of twice the ...
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
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.
computes preconditioners, and implements methods to apply them currently used in qs_ot
container for various plainwaves related things
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
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
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.
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 and initialisation of the mo data type.
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.
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...
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.
parameters that control an scf iteration