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) .OR. (tddfpt_control%do_bse_w_only) .OR. &
137 (tddfpt_control%do_bse_gw_only))
THEN
138 NULLIFY (ks_env, ex_env)
139 CALL get_qs_env(qs_env, exstate_env=ex_env, ks_env=ks_env)
140 CALL dbcsr_copy(matrix_ks(1)%matrix, ex_env%matrix_ks(1)%matrix)
144 cpassert(.NOT.
ASSOCIATED(gs_mos))
146 nspins = dft_control%nspins
147 ALLOCATE (gs_mos(nspins))
151 CALL section_vals_val_get(print_section,
"NAMD_PRINT%PRINT_VIRTUALS", l_val=print_virtuals_newtonx)
156 NULLIFY (evals_virt, evals_virt_spin, mos_virt_spin)
157 IF (
ASSOCIATED(scf_env))
THEN
158 IF ((scf_env%method ==
ot_method_nr .AND. tddfpt_control%nlumo > 0) .OR. &
159 (scf_env%method ==
ot_method_nr .AND. print_virtuals_newtonx))
THEN
165 CALL get_mo_set(mos(ispin), nmo=nmo_avail, homo=nmo_occ)
166 nmo_virt = min(nmo_virt, nmo_avail - nmo_occ)
169 nmo_virt = tddfpt_control%nlumo - nmo_virt
170 IF (.NOT. print_virtuals_newtonx)
THEN
171 IF (nmo_virt > 0)
THEN
172 ALLOCATE (evals_virt(nspins), mos_virt(nspins))
174 CALL make_lumo_gpw(qs_env, scf_env, mos_virt, evals_virt, nmo_virt, nmo_avail)
181 IF (
ASSOCIATED(evals_virt))
THEN
182 evals_virt_spin => evals_virt(ispin)%array
184 NULLIFY (evals_virt_spin)
186 IF (
ALLOCATED(mos_virt))
THEN
187 mos_virt_spin => mos_virt(ispin)
189 NULLIFY (mos_virt_spin)
192 nlumo=tddfpt_control%nlumo, &
194 matrix_ks=matrix_ks(ispin)%matrix, matrix_s=matrix_s(1)%matrix, &
195 mos_virt=mos_virt_spin, evals_virt=evals_virt_spin, &
202 CALL cp_fm_get_info(gs_mos(ispin)%mos_occ, nrow_global=moc(1, ispin), ncol_global=moc(2, ispin))
203 CALL cp_fm_get_info(gs_mos(ispin)%mos_virt, nrow_global=mvt(1, ispin), ncol_global=mvt(2, ispin))
206 WRITE (iounit,
"(T2,A,T36,A)")
"TDDFPT| Molecular Orbitals:", &
207 " Spin AOs Occ Virt Total"
209 WRITE (iounit,
"(T2,A,T37,I4,4I10)")
"TDDFPT| ", ispin, moc(1, ispin), moc(2, ispin), &
210 mvt(2, ispin), moc(2, ispin) + mvt(2, ispin)
214 IF (
ASSOCIATED(evals_virt))
THEN
215 DO ispin = 1,
SIZE(evals_virt)
216 IF (
ASSOCIATED(evals_virt(ispin)%array))
DEALLOCATE (evals_virt(ispin)%array)
218 DEALLOCATE (evals_virt)
223 CALL timestop(handle)
249 mos_virt, evals_virt, qs_env)
252 INTEGER,
INTENT(in) :: nlumo
254 INTEGER,
INTENT(in) :: cholesky_method
255 TYPE(
dbcsr_type),
POINTER :: matrix_ks, matrix_s
256 TYPE(
cp_fm_type),
INTENT(IN),
POINTER :: mos_virt
257 REAL(kind=
dp),
DIMENSION(:),
POINTER :: evals_virt
260 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_init_ground_state_mos'
261 REAL(kind=
dp),
PARAMETER :: eps_dp = epsilon(0.0_dp)
263 INTEGER :: cholesky_method_inout, handle, icol_global, icol_local, imo, iounit, irow_global, &
264 irow_local, nao, ncol_local, nelectrons, nmo_occ, nmo_scf, nmo_virt, nrow_local, sign_int
265 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: minrow_neg_array, minrow_pos_array, &
267 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
268 LOGICAL :: do_eigen, print_phases
269 REAL(kind=
dp) :: element, maxocc
270 REAL(kind=
dp),
CONTIGUOUS,
DIMENSION(:, :), &
272 REAL(kind=
dp),
DIMENSION(:),
POINTER :: mo_evals_extended, mo_occ_extended, &
275 ao_mo_virt_fm_struct, wfn_fm_struct
276 TYPE(
cp_fm_type) :: matrix_ks_fm, ortho_fm, work_fm, &
278 TYPE(
cp_fm_type),
POINTER :: mo_coeff_extended
284 CALL timeset(routinen, handle)
290 CALL blacs_env%get(para_env=para_env)
292 CALL get_mo_set(mo_set, nao=nao, nmo=nmo_scf, homo=nmo_occ, maxocc=maxocc, &
293 nelectron=nelectrons, occupation_numbers=mo_occ_scf)
298 nmo_virt = nao - nmo_occ
300 nmo_virt = min(nmo_virt, nlumo)
303 CALL cp_abort(__location__, &
304 'At least one unoccupied molecular orbital is required to calculate excited states.')
308 IF (
ASSOCIATED(evals_virt))
THEN
309 cpassert(
ASSOCIATED(mos_virt))
310 IF (nmo_virt > nmo_scf - nmo_occ +
SIZE(evals_virt)) do_eigen = .true.
312 IF (nmo_virt > nmo_scf - nmo_occ) do_eigen = .true.
316 NULLIFY (ao_mo_occ_fm_struct, ao_mo_virt_fm_struct)
318 CALL cp_fm_struct_create(ao_mo_occ_fm_struct, template_fmstruct=mo_set%mo_coeff%matrix_struct, &
319 ncol_global=nmo_occ, context=blacs_env)
320 CALL cp_fm_struct_create(ao_mo_virt_fm_struct, template_fmstruct=mo_set%mo_coeff%matrix_struct, &
321 ncol_global=nmo_virt, context=blacs_env)
323 NULLIFY (gs_mos%mos_occ, gs_mos%mos_virt, gs_mos%evals_occ_matrix)
324 ALLOCATE (gs_mos%mos_occ, gs_mos%mos_virt)
326 CALL cp_fm_create(gs_mos%mos_virt, ao_mo_virt_fm_struct)
327 gs_mos%nmo_occ = nmo_occ
329 ALLOCATE (gs_mos%evals_occ(nmo_occ))
330 ALLOCATE (gs_mos%evals_virt(nmo_virt))
331 ALLOCATE (gs_mos%phases_occ(nmo_occ))
332 ALLOCATE (gs_mos%phases_virt(nmo_virt))
335 NULLIFY (ao_ao_fm_struct, wfn_fm_struct)
336 NULLIFY (mos_extended, mo_coeff_extended, mo_evals_extended, mo_occ_extended)
340 CALL cp_fm_struct_create(ao_ao_fm_struct, nrow_global=nao, ncol_global=nao, context=blacs_env)
341 CALL cp_fm_struct_create(wfn_fm_struct, nrow_global=nao, ncol_global=nmo_occ + nmo_virt, context=blacs_env)
342 ALLOCATE (mos_extended)
343 CALL allocate_mo_set(mos_extended, nao, nmo_occ + nmo_virt, nelectrons, &
344 REAL(nelectrons,
dp), maxocc, flexible_electron_count=0.0_dp)
345 CALL init_mo_set(mos_extended, fm_struct=wfn_fm_struct, name=
"mos-extended")
347 CALL get_mo_set(mos_extended, mo_coeff=mo_coeff_extended, &
348 eigenvalues=mo_evals_extended, occupation_numbers=mo_occ_extended)
356 mo_occ_extended(imo) = mo_occ_scf(imo)
358 mo_occ_extended(nmo_scf + 1:) = 0.0_dp
371 cpabort(
'CHOLESKY DBCSR_INVERSE is not implemented in TDDFT.')
373 cpabort(
'CHOLESKY OFF is not implemented in TDDFT.')
382 cholesky_method_inout = cholesky_method
383 CALL eigensolver(matrix_ks_fm=matrix_ks_fm, mo_set=mos_extended, ortho=ortho_fm, &
384 work=work_fm, cholesky_method=cholesky_method_inout, &
385 do_level_shift=.false., level_shift=0.0_dp, use_jacobi=.false.)
393 CALL get_mo_set(mo_set, mo_coeff=mo_coeff_extended, &
394 eigenvalues=mo_evals_extended, occupation_numbers=mo_occ_extended)
403 CALL cp_fm_to_fm(mo_coeff_extended, work_fm, ncol=nmo_occ, source_start=1, target_start=1)
404 CALL cp_fm_get_info(work_fm, nrow_local=nrow_local, ncol_local=ncol_local, &
405 row_indices=row_indices, col_indices=col_indices, local_data=my_block)
407 ALLOCATE (minrow_neg_array(nmo_occ), minrow_pos_array(nmo_occ), sum_sign_array(nmo_occ))
408 minrow_neg_array(:) = nao
409 minrow_pos_array(:) = nao
410 sum_sign_array(:) = 0
411 DO icol_local = 1, ncol_local
412 icol_global = col_indices(icol_local)
414 DO irow_local = 1, nrow_local
415 element = my_block(irow_local, icol_local)
418 IF (element >= eps_dp)
THEN
420 ELSE IF (element <= -eps_dp)
THEN
424 sum_sign_array(icol_global) = sum_sign_array(icol_global) + sign_int
426 irow_global = row_indices(irow_local)
427 IF (sign_int > 0)
THEN
428 IF (minrow_pos_array(icol_global) > irow_global) &
429 minrow_pos_array(icol_global) = irow_global
430 ELSE IF (sign_int < 0)
THEN
431 IF (minrow_neg_array(icol_global) > irow_global) &
432 minrow_neg_array(icol_global) = irow_global
437 CALL para_env%sum(sum_sign_array)
438 CALL para_env%min(minrow_neg_array)
439 CALL para_env%min(minrow_pos_array)
441 DO icol_local = 1, nmo_occ
442 IF (sum_sign_array(icol_local) > 0)
THEN
444 gs_mos%phases_occ(icol_local) = 1.0_dp
445 ELSE IF (sum_sign_array(icol_local) < 0)
THEN
447 gs_mos%phases_occ(icol_local) = -1.0_dp
450 IF (minrow_pos_array(icol_local) <= minrow_neg_array(icol_local))
THEN
453 gs_mos%phases_occ(icol_local) = 1.0_dp
456 gs_mos%phases_occ(icol_local) = -1.0_dp
461 DEALLOCATE (minrow_neg_array, minrow_pos_array, sum_sign_array)
464 CALL cp_fm_to_fm(mo_coeff_extended, gs_mos%mos_occ, ncol=nmo_occ, source_start=1, target_start=1)
465 gs_mos%evals_occ(1:nmo_occ) = mo_evals_extended(1:nmo_occ)
467 IF (
ASSOCIATED(evals_virt) .AND. (.NOT. do_eigen) .AND. nmo_virt > nmo_scf - nmo_occ)
THEN
468 CALL cp_fm_to_fm(mo_coeff_extended, gs_mos%mos_virt, ncol=nmo_scf - nmo_occ, &
469 source_start=nmo_occ + 1, target_start=1)
470 CALL cp_fm_to_fm(mos_virt, gs_mos%mos_virt, ncol=nmo_virt - (nmo_scf - nmo_occ), &
471 source_start=1, target_start=nmo_scf - nmo_occ + 1)
472 gs_mos%evals_virt(1:nmo_scf - nmo_occ) = evals_virt(nmo_occ + 1:nmo_occ + nmo_scf)
473 gs_mos%evals_virt(nmo_scf - nmo_occ + 1:nmo_virt) = evals_virt(1:nmo_virt - (nmo_scf - nmo_occ))
475 CALL cp_fm_to_fm(mo_coeff_extended, gs_mos%mos_virt, ncol=nmo_virt, source_start=nmo_occ + 1, target_start=1)
476 gs_mos%evals_virt(1:nmo_virt) = mo_evals_extended(nmo_occ + 1:nmo_occ + nmo_virt)
479 IF (print_phases)
THEN
485 CALL cp_fm_to_fm(gs_mos%mos_virt, work_fm_virt, ncol=nmo_virt, source_start=1, target_start=1)
486 CALL cp_fm_get_info(work_fm_virt, nrow_local=nrow_local, ncol_local=ncol_local, &
487 row_indices=row_indices, col_indices=col_indices, local_data=my_block)
489 ALLOCATE (minrow_neg_array(nmo_virt), minrow_pos_array(nmo_virt), sum_sign_array(nmo_virt))
490 minrow_neg_array(:) = nao
491 minrow_pos_array(:) = nao
492 sum_sign_array(:) = 0
493 DO icol_local = 1, ncol_local
494 icol_global = col_indices(icol_local)
496 DO irow_local = 1, nrow_local
497 element = my_block(irow_local, icol_local)
500 IF (element >= eps_dp)
THEN
502 ELSE IF (element <= -eps_dp)
THEN
506 sum_sign_array(icol_global) = sum_sign_array(icol_global) + sign_int
508 irow_global = row_indices(irow_local)
509 IF (sign_int > 0)
THEN
510 IF (minrow_pos_array(icol_global) > irow_global) &
511 minrow_pos_array(icol_global) = irow_global
512 ELSE IF (sign_int < 0)
THEN
513 IF (minrow_neg_array(icol_global) > irow_global) &
514 minrow_neg_array(icol_global) = irow_global
519 CALL para_env%sum(sum_sign_array)
520 CALL para_env%min(minrow_neg_array)
521 CALL para_env%min(minrow_pos_array)
522 DO icol_local = 1, nmo_virt
523 IF (sum_sign_array(icol_local) > 0)
THEN
525 gs_mos%phases_virt(icol_local) = 1.0_dp
526 ELSE IF (sum_sign_array(icol_local) < 0)
THEN
528 gs_mos%phases_virt(icol_local) = -1.0_dp
531 IF (minrow_pos_array(icol_local) <= minrow_neg_array(icol_local))
THEN
534 gs_mos%phases_virt(icol_local) = 1.0_dp
537 gs_mos%phases_virt(icol_local) = -1.0_dp
541 DEALLOCATE (minrow_neg_array, minrow_pos_array, sum_sign_array)
550 DEALLOCATE (mos_extended)
553 CALL timestop(handle)
566 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_release_ground_state_mos'
570 CALL timeset(routinen, handle)
572 IF (
ALLOCATED(gs_mos%phases_occ)) &
573 DEALLOCATE (gs_mos%phases_occ)
575 IF (
ALLOCATED(gs_mos%evals_virt)) &
576 DEALLOCATE (gs_mos%evals_virt)
578 IF (
ALLOCATED(gs_mos%evals_occ)) &
579 DEALLOCATE (gs_mos%evals_occ)
581 IF (
ALLOCATED(gs_mos%phases_virt)) &
582 DEALLOCATE (gs_mos%phases_virt)
584 IF (
ALLOCATED(gs_mos%index_active)) &
585 DEALLOCATE (gs_mos%index_active)
587 IF (
ASSOCIATED(gs_mos%evals_occ_matrix))
THEN
589 DEALLOCATE (gs_mos%evals_occ_matrix)
592 IF (
ASSOCIATED(gs_mos%mos_virt))
THEN
594 DEALLOCATE (gs_mos%mos_virt)
597 IF (
ASSOCIATED(gs_mos%mos_occ))
THEN
599 DEALLOCATE (gs_mos%mos_occ)
602 IF (
ASSOCIATED(gs_mos%mos_active))
THEN
604 DEALLOCATE (gs_mos%mos_active)
607 CALL timestop(handle)
621 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_ks_oep
623 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_oecorr'
625 INTEGER :: handle, iounit, ispin, nao, nmo_occ, &
630 mo_occ_mo_occ_fm_struct
639 CALL timeset(routinen, handle)
645 CALL get_qs_env(qs_env, blacs_env=blacs_env, dft_control=dft_control, matrix_ks=matrix_ks)
646 tddfpt_control => dft_control%tddfpt2_control
650 nspins =
SIZE(gs_mos)
651 NULLIFY (matrix_ks_oep)
652 IF (tddfpt_control%oe_corr /=
oe_none)
THEN
654 WRITE (iounit,
"(1X,A)")
"", &
655 "-------------------------------------------------------------------------------", &
656 "- Orbital Eigenvalue Correction Started -", &
657 "-------------------------------------------------------------------------------"
660 CALL cp_warn(__location__, &
661 "Orbital energy correction potential is an experimental feature. "// &
662 "Use it with extreme care")
667 CALL cp_abort(__location__, &
668 "Implementation of orbital energy correction XC-potentials is "// &
669 "currently incompatible with exact-exchange functionals")
675 CALL dbcsr_copy(matrix_ks_oep(ispin)%matrix, matrix_ks(ispin)%matrix)
681 NULLIFY (xc_fun_empty)
686 IF (dft_control%qs_control%semi_empirical)
THEN
687 cpabort(
"TDDFPT with SE not possible")
688 ELSEIF (dft_control%qs_control%dftb)
THEN
689 cpabort(
"TDDFPT with DFTB not possible")
690 ELSEIF (dft_control%qs_control%xtb)
THEN
692 ext_ks_matrix=matrix_ks_oep)
695 ext_ks_matrix=matrix_ks_oep)
698 IF (tddfpt_control%oe_corr ==
oe_saop .OR. &
699 tddfpt_control%oe_corr ==
oe_lb .OR. &
700 tddfpt_control%oe_corr ==
oe_gllb)
THEN
702 WRITE (iounit,
"(T2,A)")
" Orbital energy correction of SAOP type "
704 CALL add_saop_pot(matrix_ks_oep, qs_env, tddfpt_control%oe_corr)
705 ELSE IF (tddfpt_control%oe_corr ==
oe_shift)
THEN
707 WRITE (iounit,
"(T2,A,T71,F10.3)") &
708 " Virtual Orbital Eigenvalue Shift [eV] ", tddfpt_control%ev_shift*
evolt
709 WRITE (iounit,
"(T2,A,T71,F10.3)") &
710 " Open Shell Orbital Eigenvalue Shift [eV] ", tddfpt_control%eos_shift*
evolt
712 CALL ev_shift_operator(qs_env, gs_mos, matrix_ks_oep, &
713 tddfpt_control%ev_shift, tddfpt_control%eos_shift)
715 CALL cp_abort(__location__, &
716 "Unimplemented orbital energy correction potential")
723 NULLIFY (mo_occ_mo_occ_fm_struct)
725 nmo_occ =
SIZE(gs_mos(ispin)%evals_occ)
726 CALL cp_fm_struct_create(mo_occ_mo_occ_fm_struct, nrow_global=nmo_occ, ncol_global=nmo_occ, &
728 ALLOCATE (gs_mos(ispin)%evals_occ_matrix)
729 CALL cp_fm_create(gs_mos(ispin)%evals_occ_matrix, mo_occ_mo_occ_fm_struct)
737 work_fm, ncol=nmo_occ, alpha=1.0_dp, beta=0.0_dp)
738 CALL parallel_gemm(
'T',
'N', nmo_occ, nmo_occ, nao, 1.0_dp, gs_mos(ispin)%mos_occ, work_fm, &
739 0.0_dp, gs_mos(ispin)%evals_occ_matrix)
743 WRITE (iounit,
"(1X,A)") &
744 "-------------------------------------------------------------------------------"
749 CALL timestop(handle)
765 INTEGER(kind=int_8) :: nstates_total
767 INTEGER :: ispin, nspins
770 nspins =
SIZE(gs_mos)
775 nstates_total = nstates_total + &
776 gs_mos(ispin)%nmo_active* &
777 SIZE(gs_mos(ispin)%evals_virt, kind=
int_8)
781 nstates_total = gs_mos(1)%nmo_active* &
782 SIZE(gs_mos(2)%evals_virt, kind=
int_8)
799 SUBROUTINE ev_shift_operator(qs_env, gs_mos, matrix_ks, ev_shift, eos_shift)
805 REAL(kind=
dp),
INTENT(IN) :: ev_shift, eos_shift
807 CHARACTER(len=*),
PARAMETER :: routinen =
'ev_shift_operator'
809 INTEGER :: handle, ispin, n_spins, na, nb, nhomo, &
810 nl, nos, nrow, nu, nvirt
818 CALL timeset(routinen, handle)
820 n_spins =
SIZE(gs_mos)
821 cpassert(n_spins ==
SIZE(matrix_ks))
823 IF (eos_shift /= 0.0_dp .AND. n_spins > 1)
THEN
824 cpabort(
"eos_shift not implemented")
825 CALL get_qs_env(qs_env, mos=mos, matrix_s=matrix_s)
826 smat => matrix_s(1)%matrix
832 DO ispin = 1, n_spins
833 coeff => gs_mos(ispin)%mos_occ
834 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nhomo)
835 IF (nhomo == nu)
THEN
844 alpha=-eos_shift, keep_sparsity=.true.)
849 coeff => gs_mos(ispin)%mos_virt
850 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nvirt)
852 cpassert(nvirt >= nos)
856 alpha=eos_shift, keep_sparsity=.true.)
861 IF (ev_shift /= 0.0_dp)
THEN
862 DO ispin = 1, n_spins
863 CALL dbcsr_add(matrix_ks(ispin)%matrix, smat, &
864 alpha_scalar=1.0_dp, beta_scalar=ev_shift)
865 coeff => gs_mos(ispin)%mos_occ
866 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nhomo)
870 alpha=-ev_shift, keep_sparsity=.true.)
874 coeff => gs_mos(ispin)%mos_virt
875 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nvirt)
876 cpassert(nvirt >= nos)
880 alpha=-ev_shift, keep_sparsity=.true.)
887 IF (ev_shift /= 0.0_dp)
THEN
888 CALL get_qs_env(qs_env, mos=mos, matrix_s=matrix_s)
889 smat => matrix_s(1)%matrix
890 DO ispin = 1, n_spins
891 CALL dbcsr_add(matrix_ks(ispin)%matrix, smat, &
892 alpha_scalar=1.0_dp, beta_scalar=ev_shift)
893 coeff => gs_mos(ispin)%mos_occ
894 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nhomo)
898 alpha=-ev_shift, keep_sparsity=.true.)
904 IF (eos_shift == 0.0_dp .OR. n_spins == 1)
THEN
905 DO ispin = 1, n_spins
906 IF (
ALLOCATED(gs_mos(ispin)%evals_virt))
THEN
907 gs_mos(ispin)%evals_virt(:) = gs_mos(ispin)%evals_virt(:) + ev_shift
917 IF (
ALLOCATED(gs_mos(1)%evals_occ))
THEN
918 gs_mos(1)%evals_occ(nl + 1:nu) = gs_mos(1)%evals_occ(nl + 1:nu) - eos_shift
920 IF (
ALLOCATED(gs_mos(1)%evals_virt))
THEN
921 gs_mos(1)%evals_virt(:) = gs_mos(1)%evals_virt(:) + ev_shift
923 IF (
ALLOCATED(gs_mos(2)%evals_virt))
THEN
924 gs_mos(2)%evals_virt(1:nos) = gs_mos(2)%evals_virt(1:nos) + eos_shift
925 gs_mos(2)%evals_virt(nos + 1:) = gs_mos(2)%evals_virt(nos + 1:) + ev_shift
928 IF (
ALLOCATED(gs_mos(1)%evals_virt))
THEN
929 gs_mos(1)%evals_virt(1:nos) = gs_mos(1)%evals_virt(1:nos) + eos_shift
930 gs_mos(1)%evals_virt(nos + 1:) = gs_mos(1)%evals_virt(nos + 1:) + ev_shift
932 IF (
ALLOCATED(gs_mos(2)%evals_occ))
THEN
933 gs_mos(2)%evals_occ(nl + 1:nu) = gs_mos(2)%evals_occ(nl + 1:nu) - eos_shift
935 IF (
ALLOCATED(gs_mos(2)%evals_virt))
THEN
936 gs_mos(2)%evals_virt(:) = gs_mos(2)%evals_virt(:) + ev_shift
941 CALL timestop(handle)
943 END SUBROUTINE ev_shift_operator
970 fm_pool_ao_mo_active, qs_env, nspins)
971 TYPE(
cp_fm_type),
DIMENSION(:, :),
INTENT(inout) :: evects
972 REAL(kind=
dp),
DIMENSION(:),
INTENT(inout) :: evals
973 INTEGER,
INTENT(in) :: nspins
977 INTEGER,
INTENT(in) :: log_unit
981 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_guess_vectors'
983 CHARACTER(len=5) :: spin_label1, spin_label2
984 INTEGER :: handle, i, imo_occ, imo_virt, ind, ispin, istate, j, jspin, k, no, nstates, &
985 nstates_occ_virt_alpha, nstates_selected, nv, spin1, spin2
986 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: inds
987 INTEGER,
ALLOCATABLE,
DIMENSION(:, :) :: reverse_index
988 INTEGER,
DIMENSION(maxspins) :: nmo, nmo_occ_avail, nmo_occ_selected, &
990 REAL(kind=
dp) :: e_occ
991 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: e_virt_minus_occ, ev_occ, ev_virt
994 CALL timeset(routinen, handle)
996 nstates =
SIZE(evects, 2)
998 IF (debug_this_module)
THEN
999 cpassert(nstates > 0)
1000 cpassert(nspins == 1 .OR. nspins == 2)
1006 DO ispin = 1, nspins
1008 nmo_occ_avail(ispin) = gs_mos(ispin)%nmo_active
1009 nmo(ispin) = gs_mos(ispin)%nmo_occ
1012 nmo_occ_selected(ispin) = min(nmo_occ_avail(ispin), nstates)
1013 nmo_virt_selected(ispin) = min(
SIZE(gs_mos(ispin)%evals_virt), nstates)
1019 nstates_selected = dot_product(nmo_occ_selected(1:nspins), nmo_virt_selected(1:nspins))
1021 nstates_selected = nmo_occ_selected(1)*nmo_virt_selected(2)
1024 ALLOCATE (inds(nstates_selected))
1025 ALLOCATE (e_virt_minus_occ(nstates_selected))
1030 DO ispin = 1, nspins
1031 no = nmo_occ_selected(ispin)
1032 nv = nmo_virt_selected(ispin)
1033 ALLOCATE (ev_virt(nv), ev_occ(no))
1035 IF ((tddfpt_control%do_bse) .OR. (tddfpt_control%do_bse_w_only) .OR. &
1036 (tddfpt_control%do_bse_gw_only))
THEN
1037 ev_virt(1:nv) = ex_env%gw_eigen(nmo(ispin) + 1:nmo(ispin) + nv)
1039 j = nmo_occ_avail(ispin) - i + 1
1040 k = gs_mos(ispin)%index_active(j)
1041 ev_occ(i) = ex_env%gw_eigen(k)
1044 ev_virt(1:nv) = gs_mos(ispin)%evals_virt(1:nv)
1046 j = nmo_occ_avail(ispin) - i + 1
1047 k = gs_mos(ispin)%index_active(j)
1048 ev_occ(i) = gs_mos(ispin)%evals_occ(k)
1052 DO imo_occ = 1, nmo_occ_selected(ispin)
1054 e_occ = ev_occ(imo_occ)
1056 DO imo_virt = 1, nmo_virt_selected(ispin)
1058 e_virt_minus_occ(istate) = ev_virt(imo_virt) - e_occ
1062 DEALLOCATE (ev_virt, ev_occ)
1066 DO imo_occ = 1, nmo_occ_selected(1)
1068 i = gs_mos(1)%nmo_active - imo_occ + 1
1069 k = gs_mos(1)%index_active(i)
1070 e_occ = gs_mos(1)%evals_occ(k)
1072 DO imo_virt = 1, nmo_virt_selected(2)
1074 e_virt_minus_occ(istate) = gs_mos(2)%evals_virt(imo_virt) - e_occ
1079 IF (debug_this_module)
THEN
1080 cpassert(istate == nstates_selected)
1083 CALL sort(e_virt_minus_occ, nstates_selected, inds)
1086 IF (nspins == 1)
THEN
1090 spin_label2 = spin_label1
1095 spin_label1 =
'(alp)'
1096 spin_label2 =
'(bet)'
1101 nstates_occ_virt_alpha = nmo_occ_selected(1)*nmo_virt_selected(1)
1104 nstates_occ_virt_alpha = nmo_occ_selected(1)*nmo_virt_selected(2)
1106 IF (log_unit > 0)
THEN
1107 WRITE (log_unit,
"(1X,A)")
"", &
1108 "-------------------------------------------------------------------------------", &
1109 "- TDDFPT Initial Guess -", &
1110 "-------------------------------------------------------------------------------"
1111 WRITE (log_unit,
'(T11,A)')
"State Occupied -> Virtual Excitation"
1112 WRITE (log_unit,
'(T11,A)')
"number orbital orbital energy (eV)"
1113 WRITE (log_unit,
'(1X,79("-"))')
1117 ALLOCATE (reverse_index(i, nspins))
1119 DO ispin = 1, nspins
1120 DO i = 1,
SIZE(gs_mos(ispin)%index_active)
1121 j = gs_mos(ispin)%index_active(i)
1122 reverse_index(j, ispin) = i
1126 DO istate = 1, nstates
1127 IF (
ASSOCIATED(evects(1, istate)%matrix_struct))
THEN
1130 WRITE (log_unit,
'(T7,I8,T28,A19,T60,F14.5)') &
1131 istate,
"*** restarted ***", evals(istate)*
evolt
1136 ind = inds(istate) - 1
1139 IF ((nspins > 1) .AND. (tddfpt_control%spinflip ==
no_sf_tddfpt))
THEN
1140 IF (ind < nstates_occ_virt_alpha)
THEN
1143 spin_label1 =
'(alp)'
1144 spin_label2 =
'(alp)'
1146 ind = ind - nstates_occ_virt_alpha
1149 spin_label1 =
'(bet)'
1150 spin_label2 =
'(bet)'
1156 i = ind/nmo_virt_selected(spin2) + 1
1157 j = nmo_occ_avail(spin1) - i + 1
1158 imo_occ = gs_mos(spin1)%index_active(j)
1159 imo_virt = mod(ind, nmo_virt_selected(spin2)) + 1
1161 evals(istate) = e_virt_minus_occ(istate)
1164 WRITE (log_unit,
'(T7,I8,T24,I8,T37,A5,T45,I8,T54,A5,T60,F14.5)') &
1165 istate, imo_occ, spin_label1, nmo(spin2) + imo_virt, spin_label2, e_virt_minus_occ(istate)*
evolt
1167 DO jspin = 1,
SIZE(evects, 1)
1172 IF (jspin == spin1)
THEN
1175 i = reverse_index(imo_occ, spin1)
1176 CALL cp_fm_to_fm(gs_mos(spin2)%mos_virt, evects(spin1, istate), &
1177 ncol=1, source_start=imo_virt, target_start=i)
1183 DEALLOCATE (reverse_index)
1185 IF (log_unit > 0)
THEN
1186 WRITE (log_unit,
'(/,T7,A,T50,I24)')
'Number of active states:', &
1188 WRITE (log_unit,
"(1X,A)") &
1189 "-------------------------------------------------------------------------------"
1192 DEALLOCATE (e_virt_minus_occ)
1195 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, xcint_weights, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
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, exc_accint, 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, xcint_weights, 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.