76#include "./base/base_uses.f90"
82 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_tddfpt2_utils'
84 LOGICAL,
PARAMETER,
PRIVATE :: debug_this_module = .false.
86 INTEGER,
PARAMETER,
PRIVATE :: nderivs = 3
87 INTEGER,
PARAMETER,
PRIVATE :: maxspins = 2
107 INTEGER,
INTENT(IN) :: iounit
109 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_init_mos'
111 INTEGER :: handle, ispin, nmo_avail, nmo_occ, &
113 INTEGER,
DIMENSION(2, 2) :: moc, mvt
114 LOGICAL :: print_virtuals_newtonx
115 REAL(kind=
dp),
DIMENSION(:),
POINTER :: evals_virt_spin
119 TYPE(
cp_fm_type),
ALLOCATABLE,
DIMENSION(:), &
122 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_ks, matrix_s
131 CALL timeset(routinen, handle)
133 CALL get_qs_env(qs_env, blacs_env=blacs_env, cell=cell, dft_control=dft_control, &
134 matrix_ks=matrix_ks, matrix_s=matrix_s, mos=mos, scf_env=scf_env)
135 tddfpt_control => dft_control%tddfpt2_control
136 IF (tddfpt_control%do_bse)
THEN
137 NULLIFY (ks_env, ex_env)
138 CALL get_qs_env(qs_env, exstate_env=ex_env, ks_env=ks_env)
139 CALL dbcsr_copy(matrix_ks(1)%matrix, ex_env%matrix_ks(1)%matrix)
143 cpassert(.NOT.
ASSOCIATED(gs_mos))
145 nspins = dft_control%nspins
146 ALLOCATE (gs_mos(nspins))
150 CALL section_vals_val_get(print_section,
"NAMD_PRINT%PRINT_VIRTUALS", l_val=print_virtuals_newtonx)
155 NULLIFY (evals_virt, evals_virt_spin, mos_virt_spin)
156 IF (
ASSOCIATED(scf_env))
THEN
157 IF ((scf_env%method ==
ot_method_nr .AND. tddfpt_control%nlumo > 0) .OR. &
158 (scf_env%method ==
ot_method_nr .AND. print_virtuals_newtonx))
THEN
164 CALL get_mo_set(mos(ispin), nmo=nmo_avail, homo=nmo_occ)
165 nmo_virt = min(nmo_virt, nmo_avail - nmo_occ)
168 nmo_virt = tddfpt_control%nlumo - nmo_virt
169 IF (.NOT. print_virtuals_newtonx)
THEN
170 IF (nmo_virt > 0)
THEN
171 ALLOCATE (evals_virt(nspins), mos_virt(nspins))
173 CALL make_lumo_gpw(qs_env, scf_env, mos_virt, evals_virt, nmo_virt, nmo_avail)
180 IF (
ASSOCIATED(evals_virt))
THEN
181 evals_virt_spin => evals_virt(ispin)%array
183 NULLIFY (evals_virt_spin)
185 IF (
ALLOCATED(mos_virt))
THEN
186 mos_virt_spin => mos_virt(ispin)
188 NULLIFY (mos_virt_spin)
191 nlumo=tddfpt_control%nlumo, &
193 matrix_ks=matrix_ks(ispin)%matrix, matrix_s=matrix_s(1)%matrix, &
194 mos_virt=mos_virt_spin, evals_virt=evals_virt_spin, &
201 CALL cp_fm_get_info(gs_mos(ispin)%mos_occ, nrow_global=moc(1, ispin), ncol_global=moc(2, ispin))
202 CALL cp_fm_get_info(gs_mos(ispin)%mos_virt, nrow_global=mvt(1, ispin), ncol_global=mvt(2, ispin))
205 WRITE (iounit,
"(T2,A,T36,A)")
"TDDFPT| Molecular Orbitals:", &
206 " Spin AOs Occ Virt Total"
208 WRITE (iounit,
"(T2,A,T37,I4,4I10)")
"TDDFPT| ", ispin, moc(1, ispin), moc(2, ispin), &
209 mvt(2, ispin), moc(2, ispin) + mvt(2, ispin)
213 IF (
ASSOCIATED(evals_virt))
THEN
214 DO ispin = 1,
SIZE(evals_virt)
215 IF (
ASSOCIATED(evals_virt(ispin)%array))
DEALLOCATE (evals_virt(ispin)%array)
217 DEALLOCATE (evals_virt)
222 CALL timestop(handle)
248 mos_virt, evals_virt, qs_env)
251 INTEGER,
INTENT(in) :: nlumo
253 INTEGER,
INTENT(in) :: cholesky_method
254 TYPE(
dbcsr_type),
POINTER :: matrix_ks, matrix_s
255 TYPE(
cp_fm_type),
INTENT(IN),
POINTER :: mos_virt
256 REAL(kind=
dp),
DIMENSION(:),
POINTER :: evals_virt
259 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_init_ground_state_mos'
260 REAL(kind=
dp),
PARAMETER :: eps_dp = epsilon(0.0_dp)
262 INTEGER :: cholesky_method_inout, handle, icol_global, icol_local, imo, iounit, irow_global, &
263 irow_local, nao, ncol_local, nelectrons, nmo_occ, nmo_scf, nmo_virt, nrow_local, sign_int
264 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: minrow_neg_array, minrow_pos_array, &
266 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
267 LOGICAL :: do_eigen, print_phases
268 REAL(kind=
dp) :: element, maxocc
269 REAL(kind=
dp),
CONTIGUOUS,
DIMENSION(:, :), &
271 REAL(kind=
dp),
DIMENSION(:),
POINTER :: mo_evals_extended, mo_occ_extended, &
274 ao_mo_virt_fm_struct, wfn_fm_struct
275 TYPE(
cp_fm_type) :: matrix_ks_fm, ortho_fm, work_fm, &
277 TYPE(
cp_fm_type),
POINTER :: mo_coeff_extended
283 CALL timeset(routinen, handle)
289 CALL blacs_env%get(para_env=para_env)
291 CALL get_mo_set(mo_set, nao=nao, nmo=nmo_scf, homo=nmo_occ, maxocc=maxocc, &
292 nelectron=nelectrons, occupation_numbers=mo_occ_scf)
297 nmo_virt = nao - nmo_occ
299 nmo_virt = min(nmo_virt, nlumo)
302 CALL cp_abort(__location__, &
303 'At least one unoccupied molecular orbital is required to calculate excited states.')
307 IF (
ASSOCIATED(evals_virt))
THEN
308 cpassert(
ASSOCIATED(mos_virt))
309 IF (nmo_virt > nmo_scf - nmo_occ +
SIZE(evals_virt)) do_eigen = .true.
311 IF (nmo_virt > nmo_scf - nmo_occ) do_eigen = .true.
315 NULLIFY (ao_mo_occ_fm_struct, ao_mo_virt_fm_struct)
317 CALL cp_fm_struct_create(ao_mo_occ_fm_struct, template_fmstruct=mo_set%mo_coeff%matrix_struct, &
318 ncol_global=nmo_occ, context=blacs_env)
319 CALL cp_fm_struct_create(ao_mo_virt_fm_struct, template_fmstruct=mo_set%mo_coeff%matrix_struct, &
320 ncol_global=nmo_virt, context=blacs_env)
322 NULLIFY (gs_mos%mos_occ, gs_mos%mos_virt, gs_mos%evals_occ_matrix)
323 ALLOCATE (gs_mos%mos_occ, gs_mos%mos_virt)
325 CALL cp_fm_create(gs_mos%mos_virt, ao_mo_virt_fm_struct)
326 gs_mos%nmo_occ = nmo_occ
328 ALLOCATE (gs_mos%evals_occ(nmo_occ))
329 ALLOCATE (gs_mos%evals_virt(nmo_virt))
330 ALLOCATE (gs_mos%phases_occ(nmo_occ))
331 ALLOCATE (gs_mos%phases_virt(nmo_virt))
334 NULLIFY (ao_ao_fm_struct, wfn_fm_struct)
335 NULLIFY (mos_extended, mo_coeff_extended, mo_evals_extended, mo_occ_extended)
339 CALL cp_fm_struct_create(ao_ao_fm_struct, nrow_global=nao, ncol_global=nao, context=blacs_env)
340 CALL cp_fm_struct_create(wfn_fm_struct, nrow_global=nao, ncol_global=nmo_occ + nmo_virt, context=blacs_env)
341 ALLOCATE (mos_extended)
342 CALL allocate_mo_set(mos_extended, nao, nmo_occ + nmo_virt, nelectrons, &
343 REAL(nelectrons,
dp), maxocc, flexible_electron_count=0.0_dp)
344 CALL init_mo_set(mos_extended, fm_struct=wfn_fm_struct, name=
"mos-extended")
346 CALL get_mo_set(mos_extended, mo_coeff=mo_coeff_extended, &
347 eigenvalues=mo_evals_extended, occupation_numbers=mo_occ_extended)
355 mo_occ_extended(imo) = mo_occ_scf(imo)
357 mo_occ_extended(nmo_scf + 1:) = 0.0_dp
370 cpabort(
'CHOLESKY DBCSR_INVERSE is not implemented in TDDFT.')
372 cpabort(
'CHOLESKY OFF is not implemented in TDDFT.')
381 cholesky_method_inout = cholesky_method
382 CALL eigensolver(matrix_ks_fm=matrix_ks_fm, mo_set=mos_extended, ortho=ortho_fm, &
383 work=work_fm, cholesky_method=cholesky_method_inout, &
384 do_level_shift=.false., level_shift=0.0_dp, use_jacobi=.false.)
392 CALL get_mo_set(mo_set, mo_coeff=mo_coeff_extended, &
393 eigenvalues=mo_evals_extended, occupation_numbers=mo_occ_extended)
402 CALL cp_fm_to_fm(mo_coeff_extended, work_fm, ncol=nmo_occ, source_start=1, target_start=1)
403 CALL cp_fm_get_info(work_fm, nrow_local=nrow_local, ncol_local=ncol_local, &
404 row_indices=row_indices, col_indices=col_indices, local_data=my_block)
406 ALLOCATE (minrow_neg_array(nmo_occ), minrow_pos_array(nmo_occ), sum_sign_array(nmo_occ))
407 minrow_neg_array(:) = nao
408 minrow_pos_array(:) = nao
409 sum_sign_array(:) = 0
410 DO icol_local = 1, ncol_local
411 icol_global = col_indices(icol_local)
413 DO irow_local = 1, nrow_local
414 element = my_block(irow_local, icol_local)
417 IF (element >= eps_dp)
THEN
419 ELSE IF (element <= -eps_dp)
THEN
423 sum_sign_array(icol_global) = sum_sign_array(icol_global) + sign_int
425 irow_global = row_indices(irow_local)
426 IF (sign_int > 0)
THEN
427 IF (minrow_pos_array(icol_global) > irow_global) &
428 minrow_pos_array(icol_global) = irow_global
429 ELSE IF (sign_int < 0)
THEN
430 IF (minrow_neg_array(icol_global) > irow_global) &
431 minrow_neg_array(icol_global) = irow_global
436 CALL para_env%sum(sum_sign_array)
437 CALL para_env%min(minrow_neg_array)
438 CALL para_env%min(minrow_pos_array)
440 DO icol_local = 1, nmo_occ
441 IF (sum_sign_array(icol_local) > 0)
THEN
443 gs_mos%phases_occ(icol_local) = 1.0_dp
444 ELSE IF (sum_sign_array(icol_local) < 0)
THEN
446 gs_mos%phases_occ(icol_local) = -1.0_dp
449 IF (minrow_pos_array(icol_local) <= minrow_neg_array(icol_local))
THEN
452 gs_mos%phases_occ(icol_local) = 1.0_dp
455 gs_mos%phases_occ(icol_local) = -1.0_dp
460 DEALLOCATE (minrow_neg_array, minrow_pos_array, sum_sign_array)
463 CALL cp_fm_to_fm(mo_coeff_extended, gs_mos%mos_occ, ncol=nmo_occ, source_start=1, target_start=1)
464 gs_mos%evals_occ(1:nmo_occ) = mo_evals_extended(1:nmo_occ)
466 IF (
ASSOCIATED(evals_virt) .AND. (.NOT. do_eigen) .AND. nmo_virt > nmo_scf - nmo_occ)
THEN
467 CALL cp_fm_to_fm(mo_coeff_extended, gs_mos%mos_virt, ncol=nmo_scf - nmo_occ, &
468 source_start=nmo_occ + 1, target_start=1)
469 CALL cp_fm_to_fm(mos_virt, gs_mos%mos_virt, ncol=nmo_virt - (nmo_scf - nmo_occ), &
470 source_start=1, target_start=nmo_scf - nmo_occ + 1)
471 gs_mos%evals_virt(1:nmo_scf - nmo_occ) = evals_virt(nmo_occ + 1:nmo_occ + nmo_scf)
472 gs_mos%evals_virt(nmo_scf - nmo_occ + 1:nmo_virt) = evals_virt(1:nmo_virt - (nmo_scf - nmo_occ))
474 CALL cp_fm_to_fm(mo_coeff_extended, gs_mos%mos_virt, ncol=nmo_virt, source_start=nmo_occ + 1, target_start=1)
475 gs_mos%evals_virt(1:nmo_virt) = mo_evals_extended(nmo_occ + 1:nmo_occ + nmo_virt)
478 IF (print_phases)
THEN
484 CALL cp_fm_to_fm(gs_mos%mos_virt, work_fm_virt, ncol=nmo_virt, source_start=1, target_start=1)
485 CALL cp_fm_get_info(work_fm_virt, nrow_local=nrow_local, ncol_local=ncol_local, &
486 row_indices=row_indices, col_indices=col_indices, local_data=my_block)
488 ALLOCATE (minrow_neg_array(nmo_virt), minrow_pos_array(nmo_virt), sum_sign_array(nmo_virt))
489 minrow_neg_array(:) = nao
490 minrow_pos_array(:) = nao
491 sum_sign_array(:) = 0
492 DO icol_local = 1, ncol_local
493 icol_global = col_indices(icol_local)
495 DO irow_local = 1, nrow_local
496 element = my_block(irow_local, icol_local)
499 IF (element >= eps_dp)
THEN
501 ELSE IF (element <= -eps_dp)
THEN
505 sum_sign_array(icol_global) = sum_sign_array(icol_global) + sign_int
507 irow_global = row_indices(irow_local)
508 IF (sign_int > 0)
THEN
509 IF (minrow_pos_array(icol_global) > irow_global) &
510 minrow_pos_array(icol_global) = irow_global
511 ELSE IF (sign_int < 0)
THEN
512 IF (minrow_neg_array(icol_global) > irow_global) &
513 minrow_neg_array(icol_global) = irow_global
518 CALL para_env%sum(sum_sign_array)
519 CALL para_env%min(minrow_neg_array)
520 CALL para_env%min(minrow_pos_array)
521 DO icol_local = 1, nmo_virt
522 IF (sum_sign_array(icol_local) > 0)
THEN
524 gs_mos%phases_virt(icol_local) = 1.0_dp
525 ELSE IF (sum_sign_array(icol_local) < 0)
THEN
527 gs_mos%phases_virt(icol_local) = -1.0_dp
530 IF (minrow_pos_array(icol_local) <= minrow_neg_array(icol_local))
THEN
533 gs_mos%phases_virt(icol_local) = 1.0_dp
536 gs_mos%phases_virt(icol_local) = -1.0_dp
540 DEALLOCATE (minrow_neg_array, minrow_pos_array, sum_sign_array)
549 DEALLOCATE (mos_extended)
552 CALL timestop(handle)
565 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_release_ground_state_mos'
569 CALL timeset(routinen, handle)
571 IF (
ALLOCATED(gs_mos%phases_occ)) &
572 DEALLOCATE (gs_mos%phases_occ)
574 IF (
ALLOCATED(gs_mos%evals_virt)) &
575 DEALLOCATE (gs_mos%evals_virt)
577 IF (
ALLOCATED(gs_mos%evals_occ)) &
578 DEALLOCATE (gs_mos%evals_occ)
580 IF (
ALLOCATED(gs_mos%phases_virt)) &
581 DEALLOCATE (gs_mos%phases_virt)
583 IF (
ALLOCATED(gs_mos%index_active)) &
584 DEALLOCATE (gs_mos%index_active)
586 IF (
ASSOCIATED(gs_mos%evals_occ_matrix))
THEN
588 DEALLOCATE (gs_mos%evals_occ_matrix)
591 IF (
ASSOCIATED(gs_mos%mos_virt))
THEN
593 DEALLOCATE (gs_mos%mos_virt)
596 IF (
ASSOCIATED(gs_mos%mos_occ))
THEN
598 DEALLOCATE (gs_mos%mos_occ)
601 IF (
ASSOCIATED(gs_mos%mos_active))
THEN
603 DEALLOCATE (gs_mos%mos_active)
606 CALL timestop(handle)
620 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_ks_oep
622 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_oecorr'
624 INTEGER :: handle, iounit, ispin, nao, nmo_occ, &
629 mo_occ_mo_occ_fm_struct
638 CALL timeset(routinen, handle)
644 CALL get_qs_env(qs_env, blacs_env=blacs_env, dft_control=dft_control, matrix_ks=matrix_ks)
645 tddfpt_control => dft_control%tddfpt2_control
649 nspins =
SIZE(gs_mos)
650 NULLIFY (matrix_ks_oep)
651 IF (tddfpt_control%oe_corr /=
oe_none)
THEN
653 WRITE (iounit,
"(1X,A)")
"", &
654 "-------------------------------------------------------------------------------", &
655 "- Orbital Eigenvalue Correction Started -", &
656 "-------------------------------------------------------------------------------"
659 CALL cp_warn(__location__, &
660 "Orbital energy correction potential is an experimental feature. "// &
661 "Use it with extreme care")
666 CALL cp_abort(__location__, &
667 "Implementation of orbital energy correction XC-potentials is "// &
668 "currently incompatible with exact-exchange functionals")
674 CALL dbcsr_copy(matrix_ks_oep(ispin)%matrix, matrix_ks(ispin)%matrix)
680 NULLIFY (xc_fun_empty)
685 IF (dft_control%qs_control%semi_empirical)
THEN
686 cpabort(
"TDDFPT with SE not possible")
687 ELSEIF (dft_control%qs_control%dftb)
THEN
688 cpabort(
"TDDFPT with DFTB not possible")
689 ELSEIF (dft_control%qs_control%xtb)
THEN
691 ext_ks_matrix=matrix_ks_oep)
694 ext_ks_matrix=matrix_ks_oep)
697 IF (tddfpt_control%oe_corr ==
oe_saop .OR. &
698 tddfpt_control%oe_corr ==
oe_lb .OR. &
699 tddfpt_control%oe_corr ==
oe_gllb)
THEN
701 WRITE (iounit,
"(T2,A)")
" Orbital energy correction of SAOP type "
703 CALL add_saop_pot(matrix_ks_oep, qs_env, tddfpt_control%oe_corr)
704 ELSE IF (tddfpt_control%oe_corr ==
oe_shift)
THEN
706 WRITE (iounit,
"(T2,A,T71,F10.3)") &
707 " Virtual Orbital Eigenvalue Shift [eV] ", tddfpt_control%ev_shift*
evolt
708 WRITE (iounit,
"(T2,A,T71,F10.3)") &
709 " Open Shell Orbital Eigenvalue Shift [eV] ", tddfpt_control%eos_shift*
evolt
711 CALL ev_shift_operator(qs_env, gs_mos, matrix_ks_oep, &
712 tddfpt_control%ev_shift, tddfpt_control%eos_shift)
714 CALL cp_abort(__location__, &
715 "Unimplemented orbital energy correction potential")
722 NULLIFY (mo_occ_mo_occ_fm_struct)
724 nmo_occ =
SIZE(gs_mos(ispin)%evals_occ)
725 CALL cp_fm_struct_create(mo_occ_mo_occ_fm_struct, nrow_global=nmo_occ, ncol_global=nmo_occ, &
727 ALLOCATE (gs_mos(ispin)%evals_occ_matrix)
728 CALL cp_fm_create(gs_mos(ispin)%evals_occ_matrix, mo_occ_mo_occ_fm_struct)
736 work_fm, ncol=nmo_occ, alpha=1.0_dp, beta=0.0_dp)
737 CALL parallel_gemm(
'T',
'N', nmo_occ, nmo_occ, nao, 1.0_dp, gs_mos(ispin)%mos_occ, work_fm, &
738 0.0_dp, gs_mos(ispin)%evals_occ_matrix)
742 WRITE (iounit,
"(1X,A)") &
743 "-------------------------------------------------------------------------------"
748 CALL timestop(handle)
764 INTEGER(kind=int_8) :: nstates_total
766 INTEGER :: ispin, nspins
769 nspins =
SIZE(gs_mos)
774 nstates_total = nstates_total + &
775 gs_mos(ispin)%nmo_active* &
776 SIZE(gs_mos(ispin)%evals_virt, kind=
int_8)
780 nstates_total = gs_mos(1)%nmo_active* &
781 SIZE(gs_mos(2)%evals_virt, kind=
int_8)
798 SUBROUTINE ev_shift_operator(qs_env, gs_mos, matrix_ks, ev_shift, eos_shift)
804 REAL(kind=
dp),
INTENT(IN) :: ev_shift, eos_shift
806 CHARACTER(len=*),
PARAMETER :: routinen =
'ev_shift_operator'
808 INTEGER :: handle, ispin, n_spins, na, nb, nhomo, &
809 nl, nos, nrow, nu, nvirt
817 CALL timeset(routinen, handle)
819 n_spins =
SIZE(gs_mos)
820 cpassert(n_spins ==
SIZE(matrix_ks))
822 IF (eos_shift /= 0.0_dp .AND. n_spins > 1)
THEN
823 cpabort(
"eos_shift not implemented")
824 CALL get_qs_env(qs_env, mos=mos, matrix_s=matrix_s)
825 smat => matrix_s(1)%matrix
831 DO ispin = 1, n_spins
832 coeff => gs_mos(ispin)%mos_occ
833 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nhomo)
834 IF (nhomo == nu)
THEN
843 alpha=-eos_shift, keep_sparsity=.true.)
848 coeff => gs_mos(ispin)%mos_virt
849 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nvirt)
851 cpassert(nvirt >= nos)
855 alpha=eos_shift, keep_sparsity=.true.)
860 IF (ev_shift /= 0.0_dp)
THEN
861 DO ispin = 1, n_spins
862 CALL dbcsr_add(matrix_ks(ispin)%matrix, smat, &
863 alpha_scalar=1.0_dp, beta_scalar=ev_shift)
864 coeff => gs_mos(ispin)%mos_occ
865 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nhomo)
869 alpha=-ev_shift, keep_sparsity=.true.)
873 coeff => gs_mos(ispin)%mos_virt
874 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nvirt)
875 cpassert(nvirt >= nos)
879 alpha=-ev_shift, keep_sparsity=.true.)
886 IF (ev_shift /= 0.0_dp)
THEN
887 CALL get_qs_env(qs_env, mos=mos, matrix_s=matrix_s)
888 smat => matrix_s(1)%matrix
889 DO ispin = 1, n_spins
890 CALL dbcsr_add(matrix_ks(ispin)%matrix, smat, &
891 alpha_scalar=1.0_dp, beta_scalar=ev_shift)
892 coeff => gs_mos(ispin)%mos_occ
893 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nhomo)
897 alpha=-ev_shift, keep_sparsity=.true.)
903 IF (eos_shift == 0.0_dp .OR. n_spins == 1)
THEN
904 DO ispin = 1, n_spins
905 IF (
ALLOCATED(gs_mos(ispin)%evals_virt))
THEN
906 gs_mos(ispin)%evals_virt(:) = gs_mos(ispin)%evals_virt(:) + ev_shift
916 IF (
ALLOCATED(gs_mos(1)%evals_occ))
THEN
917 gs_mos(1)%evals_occ(nl + 1:nu) = gs_mos(1)%evals_occ(nl + 1:nu) - eos_shift
919 IF (
ALLOCATED(gs_mos(1)%evals_virt))
THEN
920 gs_mos(1)%evals_virt(:) = gs_mos(1)%evals_virt(:) + ev_shift
922 IF (
ALLOCATED(gs_mos(2)%evals_virt))
THEN
923 gs_mos(2)%evals_virt(1:nos) = gs_mos(2)%evals_virt(1:nos) + eos_shift
924 gs_mos(2)%evals_virt(nos + 1:) = gs_mos(2)%evals_virt(nos + 1:) + ev_shift
927 IF (
ALLOCATED(gs_mos(1)%evals_virt))
THEN
928 gs_mos(1)%evals_virt(1:nos) = gs_mos(1)%evals_virt(1:nos) + eos_shift
929 gs_mos(1)%evals_virt(nos + 1:) = gs_mos(1)%evals_virt(nos + 1:) + ev_shift
931 IF (
ALLOCATED(gs_mos(2)%evals_occ))
THEN
932 gs_mos(2)%evals_occ(nl + 1:nu) = gs_mos(2)%evals_occ(nl + 1:nu) - eos_shift
934 IF (
ALLOCATED(gs_mos(2)%evals_virt))
THEN
935 gs_mos(2)%evals_virt(:) = gs_mos(2)%evals_virt(:) + ev_shift
940 CALL timestop(handle)
942 END SUBROUTINE ev_shift_operator
969 fm_pool_ao_mo_active, qs_env, nspins)
970 TYPE(
cp_fm_type),
DIMENSION(:, :),
INTENT(inout) :: evects
971 REAL(kind=
dp),
DIMENSION(:),
INTENT(inout) :: evals
972 INTEGER,
INTENT(in) :: nspins
976 INTEGER,
INTENT(in) :: log_unit
980 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_guess_vectors'
982 CHARACTER(len=5) :: spin_label1, spin_label2
983 INTEGER :: handle, i, imo_occ, imo_virt, ind, ispin, istate, j, jspin, k, no, nstates, &
984 nstates_occ_virt_alpha, nstates_selected, nv, spin1, spin2
985 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: inds
986 INTEGER,
ALLOCATABLE,
DIMENSION(:, :) :: reverse_index
987 INTEGER,
DIMENSION(maxspins) :: nmo, nmo_occ_avail, nmo_occ_selected, &
989 REAL(kind=
dp) :: e_occ
990 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: e_virt_minus_occ, ev_occ, ev_virt
993 CALL timeset(routinen, handle)
995 nstates =
SIZE(evects, 2)
997 IF (debug_this_module)
THEN
998 cpassert(nstates > 0)
999 cpassert(nspins == 1 .OR. nspins == 2)
1005 DO ispin = 1, nspins
1007 nmo_occ_avail(ispin) = gs_mos(ispin)%nmo_active
1008 nmo(ispin) = gs_mos(ispin)%nmo_occ
1011 nmo_occ_selected(ispin) = min(nmo_occ_avail(ispin), nstates)
1012 nmo_virt_selected(ispin) = min(
SIZE(gs_mos(ispin)%evals_virt), nstates)
1018 nstates_selected = dot_product(nmo_occ_selected(1:nspins), nmo_virt_selected(1:nspins))
1020 nstates_selected = nmo_occ_selected(1)*nmo_virt_selected(2)
1023 ALLOCATE (inds(nstates_selected))
1024 ALLOCATE (e_virt_minus_occ(nstates_selected))
1029 DO ispin = 1, nspins
1030 no = nmo_occ_selected(ispin)
1031 nv = nmo_virt_selected(ispin)
1032 ALLOCATE (ev_virt(nv), ev_occ(no))
1034 IF (tddfpt_control%do_bse)
THEN
1035 ev_virt(1:nv) = ex_env%gw_eigen(nmo(ispin) + 1:nmo(ispin) + nv)
1037 j = nmo_occ_avail(ispin) - i + 1
1038 k = gs_mos(ispin)%index_active(j)
1039 ev_occ(i) = ex_env%gw_eigen(k)
1042 ev_virt(1:nv) = gs_mos(ispin)%evals_virt(1:nv)
1044 j = nmo_occ_avail(ispin) - i + 1
1045 k = gs_mos(ispin)%index_active(j)
1046 ev_occ(i) = gs_mos(ispin)%evals_occ(k)
1050 DO imo_occ = 1, nmo_occ_selected(ispin)
1052 e_occ = ev_occ(imo_occ)
1054 DO imo_virt = 1, nmo_virt_selected(ispin)
1056 e_virt_minus_occ(istate) = ev_virt(imo_virt) - e_occ
1060 DEALLOCATE (ev_virt, ev_occ)
1064 DO imo_occ = 1, nmo_occ_selected(1)
1066 i = gs_mos(1)%nmo_active - imo_occ + 1
1067 k = gs_mos(1)%index_active(i)
1068 e_occ = gs_mos(1)%evals_occ(k)
1070 DO imo_virt = 1, nmo_virt_selected(2)
1072 e_virt_minus_occ(istate) = gs_mos(2)%evals_virt(imo_virt) - e_occ
1077 IF (debug_this_module)
THEN
1078 cpassert(istate == nstates_selected)
1081 CALL sort(e_virt_minus_occ, nstates_selected, inds)
1084 IF (nspins == 1)
THEN
1088 spin_label2 = spin_label1
1093 spin_label1 =
'(alp)'
1094 spin_label2 =
'(bet)'
1099 nstates_occ_virt_alpha = nmo_occ_selected(1)*nmo_virt_selected(1)
1102 nstates_occ_virt_alpha = nmo_occ_selected(1)*nmo_virt_selected(2)
1104 IF (log_unit > 0)
THEN
1105 WRITE (log_unit,
"(1X,A)")
"", &
1106 "-------------------------------------------------------------------------------", &
1107 "- TDDFPT Initial Guess -", &
1108 "-------------------------------------------------------------------------------"
1109 WRITE (log_unit,
'(T11,A)')
"State Occupied -> Virtual Excitation"
1110 WRITE (log_unit,
'(T11,A)')
"number orbital orbital energy (eV)"
1111 WRITE (log_unit,
'(1X,79("-"))')
1115 ALLOCATE (reverse_index(i, nspins))
1117 DO ispin = 1, nspins
1118 DO i = 1,
SIZE(gs_mos(ispin)%index_active)
1119 j = gs_mos(ispin)%index_active(i)
1120 reverse_index(j, ispin) = i
1124 DO istate = 1, nstates
1125 IF (
ASSOCIATED(evects(1, istate)%matrix_struct))
THEN
1128 WRITE (log_unit,
'(T7,I8,T28,A19,T60,F14.5)') &
1129 istate,
"*** restarted ***", evals(istate)*
evolt
1134 ind = inds(istate) - 1
1137 IF ((nspins > 1) .AND. (tddfpt_control%spinflip ==
no_sf_tddfpt))
THEN
1138 IF (ind < nstates_occ_virt_alpha)
THEN
1141 spin_label1 =
'(alp)'
1142 spin_label2 =
'(alp)'
1144 ind = ind - nstates_occ_virt_alpha
1147 spin_label1 =
'(bet)'
1148 spin_label2 =
'(bet)'
1154 i = ind/nmo_virt_selected(spin2) + 1
1155 j = nmo_occ_avail(spin1) - i + 1
1156 imo_occ = gs_mos(spin1)%index_active(j)
1157 imo_virt = mod(ind, nmo_virt_selected(spin2)) + 1
1159 evals(istate) = e_virt_minus_occ(istate)
1162 WRITE (log_unit,
'(T7,I8,T24,I8,T37,A5,T45,I8,T54,A5,T60,F14.5)') &
1163 istate, imo_occ, spin_label1, nmo(spin2) + imo_virt, spin_label2, e_virt_minus_occ(istate)*
evolt
1165 DO jspin = 1,
SIZE(evects, 1)
1170 IF (jspin == spin1)
THEN
1173 i = reverse_index(imo_occ, spin1)
1174 CALL cp_fm_to_fm(gs_mos(spin2)%mos_virt, evects(spin1, istate), &
1175 ncol=1, source_start=imo_virt, target_start=i)
1181 DEALLOCATE (reverse_index)
1183 IF (log_unit > 0)
THEN
1184 WRITE (log_unit,
'(/,T7,A,T50,I24)')
'Number of active states:', &
1186 WRITE (log_unit,
"(1X,A)") &
1187 "-------------------------------------------------------------------------------"
1190 DEALLOCATE (e_virt_minus_occ)
1193 CALL timestop(handle)
Handles all functions related to the CELL.
various utilities that regard array of different kinds: output, allocation,... maybe it is not a good...
methods related to the blacs parallel environment
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
subroutine, public dbcsr_get_info(matrix, nblkrows_total, nblkcols_total, nfullrows_total, nfullcols_total, nblkrows_local, nblkcols_local, nfullrows_local, nfullcols_local, my_prow, my_pcol, local_rows, local_cols, proc_row_dist, proc_col_dist, row_blk_size, col_blk_size, row_blk_offset, col_blk_offset, distribution, name, matrix_type, group)
...
subroutine, public dbcsr_init_p(matrix)
...
subroutine, public dbcsr_add(matrix_a, matrix_b, alpha_scalar, beta_scalar)
...
DBCSR operations in CP2K.
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
subroutine, public cp_dbcsr_plus_fm_fm_t(sparse_matrix, matrix_v, matrix_g, ncol, alpha, keep_sparsity, symmetry_mode)
performs the multiplication sparse_matrix+dense_mat*dens_mat^T if matrix_g is not explicitly given,...
Basic linear algebra operations for full matrices.
subroutine, public cp_fm_triangular_invert(matrix_a, uplo_tr)
inverts a triangular matrix
various cholesky decomposition related routines
subroutine, public cp_fm_cholesky_decompose(matrix, n, info_out)
used to replace a symmetric positive def. matrix M with its cholesky decomposition U: M = U^T * U,...
pool for for elements that are retained and released
subroutine, public fm_pool_create_fm(pool, element, name)
returns an element, allocating it if none is in the pool
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_create(matrix, matrix_struct, name, use_sp, nrow, ncol, set_zero)
creates a new full matrix with the given structure
subroutine, public cp_fm_to_fm_submat(msource, mtarget, nrow, ncol, s_firstrow, s_firstcol, t_firstrow, t_firstcol)
copy just a part ot the matrix
subroutine, public cp_fm_set_all(matrix, alpha, beta)
set all elements of a matrix to the same value, and optionally the diagonal to a different one
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
Types for excited states potential energies.
Defines the basic variable types.
integer, parameter, public int_8
integer, parameter, public dp
Interface to the message passing library MPI.
basic linear algebra operations for full matrixes
Definition of physical constants:
real(kind=dp), parameter, public evolt
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, 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.
routines that build the Kohn-Sham matrix (i.e calculate the coulomb and xc parts
subroutine, public qs_ks_build_kohn_sham_matrix(qs_env, calculate_forces, just_energy, print_active, ext_ks_matrix, ext_xc_section)
routine where the real calculations are made: the KS matrix is calculated
subroutine, public set_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, vppl, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, kpoints, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, task_list, task_list_soft, subsys, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env)
...
Definition and initialisation of the mo data type.
subroutine, public allocate_mo_set(mo_set, nao, nmo, nelectron, n_el_f, maxocc, flexible_electron_count)
Allocates a mo set and partially initializes it (nao,nmo,nelectron, and flexible_electron_count are v...
subroutine, public deallocate_mo_set(mo_set)
Deallocate a wavefunction data structure.
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.
subroutine, public init_mo_set(mo_set, fm_pool, fm_ref, fm_struct, name)
initializes an allocated mo_set. eigenvalues, mo_coeff, occupation_numbers are valid only after this ...
groups fairly general SCF methods, so that modules other than qs_scf can use them too split off from ...
subroutine, public eigensolver(matrix_ks_fm, mo_set, ortho, work, cholesky_method, do_level_shift, level_shift, matrix_u_fm, use_jacobi)
Diagonalise the Kohn-Sham matrix to get a new set of MO eigen- vectors and MO eigenvalues....
Does all kind of post scf calculations for GPW/GAPW.
subroutine, public make_lumo_gpw(qs_env, scf_env, unoccupied_orbs, unoccupied_evals, nlumo, nlumos)
Gets the lumos, and eigenvalues for the lumos.
module that contains the definitions of the scf types
integer, parameter, public ot_method_nr
subroutine, public tddfpt_release_ground_state_mos(gs_mos)
Release molecular orbitals.
subroutine, public tddfpt_init_ground_state_mos(gs_mos, mo_set, nlumo, blacs_env, cholesky_method, matrix_ks, matrix_s, mos_virt, evals_virt, qs_env)
Generate all virtual molecular orbitals for a given spin by diagonalising the corresponding Kohn-Sham...
subroutine, public tddfpt_oecorr(qs_env, gs_mos, matrix_ks_oep)
Callculate orbital corrected KS matrix for TDDFPT.
subroutine, public tddfpt_guess_vectors(evects, evals, gs_mos, log_unit, tddfpt_control, fm_pool_ao_mo_active, qs_env, nspins)
Generate missed guess vectors.
subroutine, public tddfpt_init_mos(qs_env, gs_mos, iounit)
Prepare MOs for TDDFPT Calculations.
pure integer(kind=int_8) function, public tddfpt_total_number_of_states(tddfpt_control, gs_mos)
Compute the number of possible singly excited states (occ -> virt)
All kind of helpful little routines.
Calculate the saop potential.
subroutine, public add_saop_pot(ks_matrix, qs_env, oe_corr)
...
Calculation of KS matrix in xTB Reference: Stefan Grimme, Christoph Bannwarth, Philip Shushkov JCTC 1...
subroutine, public build_xtb_ks_matrix(qs_env, calculate_forces, just_energy, ext_ks_matrix)
...
Type defining parameters related to the simulation cell.
represent a pointer to a 1d array
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
to create arrays of pools
keeps the information about the structure of a full matrix
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Contains information on the excited states energy.
stores all the informations relevant to an mpi environment
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
Ground state molecular orbitals.