72#include "./base/base_uses.f90"
78 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_tddfpt2_utils'
80 LOGICAL,
PARAMETER,
PRIVATE :: debug_this_module = .false.
82 INTEGER,
PARAMETER,
PRIVATE :: nderivs = 3
83 INTEGER,
PARAMETER,
PRIVATE :: maxspins = 2
103 INTEGER,
INTENT(IN) :: iounit
105 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_init_mos'
107 INTEGER :: handle, ispin, nmo_avail, nmo_occ, &
109 INTEGER,
DIMENSION(2, 2) :: moc, mvt
110 LOGICAL :: print_virtuals_newtonx
111 REAL(kind=
dp),
DIMENSION(:),
POINTER :: evals_virt_spin
115 TYPE(
cp_fm_type),
ALLOCATABLE,
DIMENSION(:), &
118 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_ks, matrix_s
125 CALL timeset(routinen, handle)
127 CALL get_qs_env(qs_env, blacs_env=blacs_env, cell=cell, dft_control=dft_control, &
128 matrix_ks=matrix_ks, matrix_s=matrix_s, mos=mos, scf_env=scf_env)
129 tddfpt_control => dft_control%tddfpt2_control
131 cpassert(.NOT.
ASSOCIATED(gs_mos))
133 nspins = dft_control%nspins
134 ALLOCATE (gs_mos(nspins))
138 CALL section_vals_val_get(print_section,
"NAMD_PRINT%PRINT_VIRTUALS", l_val=print_virtuals_newtonx)
143 NULLIFY (evals_virt, evals_virt_spin, mos_virt_spin)
144 IF (
ASSOCIATED(scf_env))
THEN
145 IF ((scf_env%method ==
ot_method_nr .AND. tddfpt_control%nlumo > 0) .OR. &
146 (scf_env%method ==
ot_method_nr .AND. print_virtuals_newtonx))
THEN
152 CALL get_mo_set(mos(ispin), nmo=nmo_avail, homo=nmo_occ)
153 nmo_virt = min(nmo_virt, nmo_avail - nmo_occ)
156 nmo_virt = tddfpt_control%nlumo - nmo_virt
157 IF (.NOT. print_virtuals_newtonx)
THEN
158 IF (nmo_virt > 0)
THEN
159 ALLOCATE (evals_virt(nspins), mos_virt(nspins))
161 CALL make_lumo_gpw(qs_env, scf_env, mos_virt, evals_virt, nmo_virt, nmo_avail)
168 IF (
ASSOCIATED(evals_virt))
THEN
169 evals_virt_spin => evals_virt(ispin)%array
171 NULLIFY (evals_virt_spin)
173 IF (
ALLOCATED(mos_virt))
THEN
174 mos_virt_spin => mos_virt(ispin)
176 NULLIFY (mos_virt_spin)
179 nlumo=tddfpt_control%nlumo, &
181 matrix_ks=matrix_ks(ispin)%matrix, matrix_s=matrix_s(1)%matrix, &
182 mos_virt=mos_virt_spin, evals_virt=evals_virt_spin, &
189 CALL cp_fm_get_info(gs_mos(ispin)%mos_occ, nrow_global=moc(1, ispin), ncol_global=moc(2, ispin))
190 CALL cp_fm_get_info(gs_mos(ispin)%mos_virt, nrow_global=mvt(1, ispin), ncol_global=mvt(2, ispin))
193 WRITE (iounit,
"(T2,A,T36,A)")
"TDDFPT| Molecular Orbitals:", &
194 " Spin AOs Occ Virt Total"
196 WRITE (iounit,
"(T2,A,T37,I4,4I10)")
"TDDFPT| ", ispin, moc(1, ispin), moc(2, ispin), &
197 mvt(2, ispin), moc(2, ispin) + mvt(2, ispin)
201 IF (
ASSOCIATED(evals_virt))
THEN
202 DO ispin = 1,
SIZE(evals_virt)
203 IF (
ASSOCIATED(evals_virt(ispin)%array))
DEALLOCATE (evals_virt(ispin)%array)
205 DEALLOCATE (evals_virt)
210 CALL timestop(handle)
236 mos_virt, evals_virt, qs_env)
239 INTEGER,
INTENT(in) :: nlumo
241 INTEGER,
INTENT(in) :: cholesky_method
242 TYPE(
dbcsr_type),
POINTER :: matrix_ks, matrix_s
243 TYPE(
cp_fm_type),
INTENT(IN),
POINTER :: mos_virt
244 REAL(kind=
dp),
DIMENSION(:),
POINTER :: evals_virt
247 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_init_ground_state_mos'
248 REAL(kind=
dp),
PARAMETER :: eps_dp = epsilon(0.0_dp)
250 INTEGER :: cholesky_method_inout, handle, icol_global, icol_local, imo, iounit, irow_global, &
251 irow_local, nao, ncol_local, nelectrons, nmo_occ, nmo_scf, nmo_virt, nrow_local, sign_int
252 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: minrow_neg_array, minrow_pos_array, &
254 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
255 LOGICAL :: do_eigen, print_phases
256 REAL(kind=
dp) :: element, maxocc
257 REAL(kind=
dp),
CONTIGUOUS,
DIMENSION(:, :), &
259 REAL(kind=
dp),
DIMENSION(:),
POINTER :: mo_evals_extended, mo_occ_extended, &
262 ao_mo_virt_fm_struct, wfn_fm_struct
263 TYPE(
cp_fm_type) :: matrix_ks_fm, ortho_fm, work_fm, &
265 TYPE(
cp_fm_type),
POINTER :: mo_coeff_extended
271 CALL timeset(routinen, handle)
277 CALL blacs_env%get(para_env=para_env)
279 CALL get_mo_set(mo_set, nao=nao, nmo=nmo_scf, homo=nmo_occ, maxocc=maxocc, &
280 nelectron=nelectrons, occupation_numbers=mo_occ_scf)
285 nmo_virt = nao - nmo_occ
287 nmo_virt = min(nmo_virt, nlumo)
290 CALL cp_abort(__location__, &
291 'At least one unoccupied molecular orbital is required to calculate excited states.')
295 IF (
ASSOCIATED(evals_virt))
THEN
296 cpassert(
ASSOCIATED(mos_virt))
297 IF (nmo_virt > nmo_scf - nmo_occ +
SIZE(evals_virt)) do_eigen = .true.
299 IF (nmo_virt > nmo_scf - nmo_occ) do_eigen = .true.
303 NULLIFY (ao_mo_occ_fm_struct, ao_mo_virt_fm_struct)
304 CALL cp_fm_struct_create(ao_mo_occ_fm_struct, nrow_global=nao, ncol_global=nmo_occ, context=blacs_env)
305 CALL cp_fm_struct_create(ao_mo_virt_fm_struct, nrow_global=nao, ncol_global=nmo_virt, context=blacs_env)
307 NULLIFY (gs_mos%mos_occ, gs_mos%mos_virt, gs_mos%evals_occ_matrix)
308 ALLOCATE (gs_mos%mos_occ, gs_mos%mos_virt)
310 CALL cp_fm_create(gs_mos%mos_virt, ao_mo_virt_fm_struct)
312 ALLOCATE (gs_mos%evals_occ(nmo_occ))
313 ALLOCATE (gs_mos%evals_virt(nmo_virt))
314 ALLOCATE (gs_mos%phases_occ(nmo_occ))
315 ALLOCATE (gs_mos%phases_virt(nmo_virt))
318 NULLIFY (ao_ao_fm_struct, wfn_fm_struct)
319 NULLIFY (mos_extended, mo_coeff_extended, mo_evals_extended, mo_occ_extended)
320 CALL cp_fm_struct_create(ao_ao_fm_struct, nrow_global=nao, ncol_global=nao, context=blacs_env)
324 CALL cp_fm_struct_create(wfn_fm_struct, nrow_global=nao, ncol_global=nmo_occ + nmo_virt, context=blacs_env)
325 ALLOCATE (mos_extended)
326 CALL allocate_mo_set(mos_extended, nao, nmo_occ + nmo_virt, nelectrons, &
327 REAL(nelectrons,
dp), maxocc, flexible_electron_count=0.0_dp)
328 CALL init_mo_set(mos_extended, fm_struct=wfn_fm_struct, name=
"mos-extended")
330 CALL get_mo_set(mos_extended, mo_coeff=mo_coeff_extended, &
331 eigenvalues=mo_evals_extended, occupation_numbers=mo_occ_extended)
339 mo_occ_extended(imo) = mo_occ_scf(imo)
341 mo_occ_extended(nmo_scf + 1:) = 0.0_dp
354 cpabort(
'CHOLESKY DBCSR_INVERSE is not implemented in TDDFT.')
356 cpabort(
'CHOLESKY OFF is not implemented in TDDFT.')
365 cholesky_method_inout = cholesky_method
366 CALL eigensolver(matrix_ks_fm=matrix_ks_fm, mo_set=mos_extended, ortho=ortho_fm, &
367 work=work_fm, cholesky_method=cholesky_method_inout, &
368 do_level_shift=.false., level_shift=0.0_dp, use_jacobi=.false.)
376 CALL get_mo_set(mo_set, mo_coeff=mo_coeff_extended, &
377 eigenvalues=mo_evals_extended, occupation_numbers=mo_occ_extended)
386 CALL cp_fm_to_fm(mo_coeff_extended, work_fm, ncol=nmo_occ, source_start=1, target_start=1)
387 CALL cp_fm_get_info(work_fm, nrow_local=nrow_local, ncol_local=ncol_local, &
388 row_indices=row_indices, col_indices=col_indices, local_data=my_block)
390 ALLOCATE (minrow_neg_array(nmo_occ), minrow_pos_array(nmo_occ), sum_sign_array(nmo_occ))
391 minrow_neg_array(:) = nao
392 minrow_pos_array(:) = nao
393 sum_sign_array(:) = 0
394 DO icol_local = 1, ncol_local
395 icol_global = col_indices(icol_local)
397 DO irow_local = 1, nrow_local
398 element = my_block(irow_local, icol_local)
401 IF (element >= eps_dp)
THEN
403 ELSE IF (element <= -eps_dp)
THEN
407 sum_sign_array(icol_global) = sum_sign_array(icol_global) + sign_int
409 irow_global = row_indices(irow_local)
410 IF (sign_int > 0)
THEN
411 IF (minrow_pos_array(icol_global) > irow_global) &
412 minrow_pos_array(icol_global) = irow_global
413 ELSE IF (sign_int < 0)
THEN
414 IF (minrow_neg_array(icol_global) > irow_global) &
415 minrow_neg_array(icol_global) = irow_global
420 CALL para_env%sum(sum_sign_array)
421 CALL para_env%min(minrow_neg_array)
422 CALL para_env%min(minrow_pos_array)
424 DO icol_local = 1, nmo_occ
425 IF (sum_sign_array(icol_local) > 0)
THEN
427 gs_mos%phases_occ(icol_local) = 1.0_dp
428 ELSE IF (sum_sign_array(icol_local) < 0)
THEN
430 gs_mos%phases_occ(icol_local) = -1.0_dp
433 IF (minrow_pos_array(icol_local) <= minrow_neg_array(icol_local))
THEN
436 gs_mos%phases_occ(icol_local) = 1.0_dp
439 gs_mos%phases_occ(icol_local) = -1.0_dp
444 DEALLOCATE (minrow_neg_array, minrow_pos_array, sum_sign_array)
447 CALL cp_fm_to_fm(mo_coeff_extended, gs_mos%mos_occ, ncol=nmo_occ, source_start=1, target_start=1)
448 gs_mos%evals_occ(1:nmo_occ) = mo_evals_extended(1:nmo_occ)
450 IF (
ASSOCIATED(evals_virt) .AND. (.NOT. do_eigen) .AND. nmo_virt > nmo_scf - nmo_occ)
THEN
451 CALL cp_fm_to_fm(mo_coeff_extended, gs_mos%mos_virt, ncol=nmo_scf - nmo_occ, &
452 source_start=nmo_occ + 1, target_start=1)
453 CALL cp_fm_to_fm(mos_virt, gs_mos%mos_virt, ncol=nmo_virt - (nmo_scf - nmo_occ), &
454 source_start=1, target_start=nmo_scf - nmo_occ + 1)
455 gs_mos%evals_virt(1:nmo_scf - nmo_occ) = evals_virt(nmo_occ + 1:nmo_occ + nmo_scf)
456 gs_mos%evals_virt(nmo_scf - nmo_occ + 1:nmo_virt) = evals_virt(1:nmo_virt - (nmo_scf - nmo_occ))
458 CALL cp_fm_to_fm(mo_coeff_extended, gs_mos%mos_virt, ncol=nmo_virt, source_start=nmo_occ + 1, target_start=1)
459 gs_mos%evals_virt(1:nmo_virt) = mo_evals_extended(nmo_occ + 1:nmo_occ + nmo_virt)
462 IF (print_phases)
THEN
468 CALL cp_fm_to_fm(gs_mos%mos_virt, work_fm_virt, ncol=nmo_virt, source_start=1, target_start=1)
469 CALL cp_fm_get_info(work_fm_virt, nrow_local=nrow_local, ncol_local=ncol_local, &
470 row_indices=row_indices, col_indices=col_indices, local_data=my_block)
472 ALLOCATE (minrow_neg_array(nmo_virt), minrow_pos_array(nmo_virt), sum_sign_array(nmo_virt))
473 minrow_neg_array(:) = nao
474 minrow_pos_array(:) = nao
475 sum_sign_array(:) = 0
476 DO icol_local = 1, ncol_local
477 icol_global = col_indices(icol_local)
479 DO irow_local = 1, nrow_local
480 element = my_block(irow_local, icol_local)
483 IF (element >= eps_dp)
THEN
485 ELSE IF (element <= -eps_dp)
THEN
489 sum_sign_array(icol_global) = sum_sign_array(icol_global) + sign_int
491 irow_global = row_indices(irow_local)
492 IF (sign_int > 0)
THEN
493 IF (minrow_pos_array(icol_global) > irow_global) &
494 minrow_pos_array(icol_global) = irow_global
495 ELSE IF (sign_int < 0)
THEN
496 IF (minrow_neg_array(icol_global) > irow_global) &
497 minrow_neg_array(icol_global) = irow_global
502 CALL para_env%sum(sum_sign_array)
503 CALL para_env%min(minrow_neg_array)
504 CALL para_env%min(minrow_pos_array)
505 DO icol_local = 1, nmo_virt
506 IF (sum_sign_array(icol_local) > 0)
THEN
508 gs_mos%phases_virt(icol_local) = 1.0_dp
509 ELSE IF (sum_sign_array(icol_local) < 0)
THEN
511 gs_mos%phases_virt(icol_local) = -1.0_dp
514 IF (minrow_pos_array(icol_local) <= minrow_neg_array(icol_local))
THEN
517 gs_mos%phases_virt(icol_local) = 1.0_dp
520 gs_mos%phases_virt(icol_local) = -1.0_dp
524 DEALLOCATE (minrow_neg_array, minrow_pos_array, sum_sign_array)
533 DEALLOCATE (mos_extended)
536 CALL timestop(handle)
549 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_release_ground_state_mos'
553 CALL timeset(routinen, handle)
555 IF (
ALLOCATED(gs_mos%phases_occ)) &
556 DEALLOCATE (gs_mos%phases_occ)
558 IF (
ALLOCATED(gs_mos%evals_virt)) &
559 DEALLOCATE (gs_mos%evals_virt)
561 IF (
ALLOCATED(gs_mos%evals_occ)) &
562 DEALLOCATE (gs_mos%evals_occ)
564 IF (
ALLOCATED(gs_mos%phases_virt)) &
565 DEALLOCATE (gs_mos%phases_virt)
567 IF (
ASSOCIATED(gs_mos%evals_occ_matrix))
THEN
569 DEALLOCATE (gs_mos%evals_occ_matrix)
572 IF (
ASSOCIATED(gs_mos%mos_virt))
THEN
574 DEALLOCATE (gs_mos%mos_virt)
577 IF (
ASSOCIATED(gs_mos%mos_occ))
THEN
579 DEALLOCATE (gs_mos%mos_occ)
582 CALL timestop(handle)
595 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_ks_oep
597 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_oecorr'
599 INTEGER :: handle, iounit, ispin, nao, nmo_occ, &
604 mo_occ_mo_occ_fm_struct
613 CALL timeset(routinen, handle)
619 CALL get_qs_env(qs_env, blacs_env=blacs_env, dft_control=dft_control, matrix_ks=matrix_ks)
620 tddfpt_control => dft_control%tddfpt2_control
624 nspins =
SIZE(gs_mos)
625 NULLIFY (matrix_ks_oep)
626 IF (tddfpt_control%oe_corr /=
oe_none)
THEN
628 WRITE (iounit,
"(1X,A)")
"", &
629 "-------------------------------------------------------------------------------", &
630 "- Orbital Eigenvalue Correction Started -", &
631 "-------------------------------------------------------------------------------"
634 CALL cp_warn(__location__, &
635 "Orbital energy correction potential is an experimental feature. "// &
636 "Use it with extreme care")
641 CALL cp_abort(__location__, &
642 "Implementation of orbital energy correction XC-potentials is "// &
643 "currently incompatible with exact-exchange functionals")
649 CALL dbcsr_copy(matrix_ks_oep(ispin)%matrix, matrix_ks(ispin)%matrix)
655 NULLIFY (xc_fun_empty)
660 IF (dft_control%qs_control%semi_empirical)
THEN
661 cpabort(
"TDDFPT with SE not possible")
662 ELSEIF (dft_control%qs_control%dftb)
THEN
663 cpabort(
"TDDFPT with DFTB not possible")
664 ELSEIF (dft_control%qs_control%xtb)
THEN
666 ext_ks_matrix=matrix_ks_oep)
669 ext_ks_matrix=matrix_ks_oep)
672 IF (tddfpt_control%oe_corr ==
oe_saop .OR. &
673 tddfpt_control%oe_corr ==
oe_lb .OR. &
674 tddfpt_control%oe_corr ==
oe_gllb)
THEN
676 WRITE (iounit,
"(T2,A)")
" Orbital energy correction of SAOP type "
678 CALL add_saop_pot(matrix_ks_oep, qs_env, tddfpt_control%oe_corr)
679 ELSE IF (tddfpt_control%oe_corr ==
oe_shift)
THEN
681 WRITE (iounit,
"(T2,A,T71,F10.3)") &
682 " Virtual Orbital Eigenvalue Shift [eV] ", tddfpt_control%ev_shift*
evolt
683 WRITE (iounit,
"(T2,A,T71,F10.3)") &
684 " Open Shell Orbital Eigenvalue Shift [eV] ", tddfpt_control%eos_shift*
evolt
686 CALL ev_shift_operator(qs_env, gs_mos, matrix_ks_oep, &
687 tddfpt_control%ev_shift, tddfpt_control%eos_shift)
689 CALL cp_abort(__location__, &
690 "Unimplemented orbital energy correction potential")
697 NULLIFY (mo_occ_mo_occ_fm_struct)
699 nmo_occ =
SIZE(gs_mos(ispin)%evals_occ)
700 CALL cp_fm_struct_create(mo_occ_mo_occ_fm_struct, nrow_global=nmo_occ, ncol_global=nmo_occ, &
702 ALLOCATE (gs_mos(ispin)%evals_occ_matrix)
703 CALL cp_fm_create(gs_mos(ispin)%evals_occ_matrix, mo_occ_mo_occ_fm_struct)
711 work_fm, ncol=nmo_occ, alpha=1.0_dp, beta=0.0_dp)
712 CALL parallel_gemm(
'T',
'N', nmo_occ, nmo_occ, nao, 1.0_dp, gs_mos(ispin)%mos_occ, work_fm, &
713 0.0_dp, gs_mos(ispin)%evals_occ_matrix)
717 WRITE (iounit,
"(1X,A)") &
718 "-------------------------------------------------------------------------------"
723 CALL timestop(handle)
737 INTEGER(kind=int_8) :: nstates_total
739 INTEGER :: ispin, nspins
742 nspins =
SIZE(gs_mos)
745 nstates_total = nstates_total + &
746 SIZE(gs_mos(ispin)%evals_occ, kind=
int_8)* &
747 SIZE(gs_mos(ispin)%evals_virt, kind=
int_8)
764 SUBROUTINE ev_shift_operator(qs_env, gs_mos, matrix_ks, ev_shift, eos_shift)
770 REAL(kind=
dp),
INTENT(IN) :: ev_shift, eos_shift
772 CHARACTER(len=*),
PARAMETER :: routinen =
'ev_shift_operator'
774 INTEGER :: handle, ispin, n_spins, na, nb, nhomo, &
775 nl, nos, nrow, nu, nvirt
783 CALL timeset(routinen, handle)
785 n_spins =
SIZE(gs_mos)
786 cpassert(n_spins ==
SIZE(matrix_ks))
788 IF (eos_shift /= 0.0_dp .AND. n_spins > 1)
THEN
789 cpabort(
"eos_shift not implemented")
790 CALL get_qs_env(qs_env, mos=mos, matrix_s=matrix_s)
791 smat => matrix_s(1)%matrix
797 DO ispin = 1, n_spins
798 coeff => gs_mos(ispin)%mos_occ
799 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nhomo)
800 IF (nhomo == nu)
THEN
809 alpha=-eos_shift, keep_sparsity=.true.)
814 coeff => gs_mos(ispin)%mos_virt
815 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nvirt)
817 cpassert(nvirt >= nos)
821 alpha=eos_shift, keep_sparsity=.true.)
826 IF (ev_shift /= 0.0_dp)
THEN
827 DO ispin = 1, n_spins
828 CALL dbcsr_add(matrix_ks(ispin)%matrix, smat, &
829 alpha_scalar=1.0_dp, beta_scalar=ev_shift)
830 coeff => gs_mos(ispin)%mos_occ
831 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nhomo)
835 alpha=-ev_shift, keep_sparsity=.true.)
839 coeff => gs_mos(ispin)%mos_virt
840 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nvirt)
841 cpassert(nvirt >= nos)
845 alpha=-ev_shift, keep_sparsity=.true.)
852 IF (ev_shift /= 0.0_dp)
THEN
853 CALL get_qs_env(qs_env, mos=mos, matrix_s=matrix_s)
854 smat => matrix_s(1)%matrix
855 DO ispin = 1, n_spins
856 CALL dbcsr_add(matrix_ks(ispin)%matrix, smat, &
857 alpha_scalar=1.0_dp, beta_scalar=ev_shift)
858 coeff => gs_mos(ispin)%mos_occ
859 CALL cp_fm_get_info(coeff, matrix_struct=fmstruct, ncol_global=nhomo)
863 alpha=-ev_shift, keep_sparsity=.true.)
869 IF (eos_shift == 0.0_dp .OR. n_spins == 1)
THEN
870 DO ispin = 1, n_spins
871 IF (
ALLOCATED(gs_mos(ispin)%evals_virt))
THEN
872 gs_mos(ispin)%evals_virt(:) = gs_mos(ispin)%evals_virt(:) + ev_shift
882 IF (
ALLOCATED(gs_mos(1)%evals_occ))
THEN
883 gs_mos(1)%evals_occ(nl + 1:nu) = gs_mos(1)%evals_occ(nl + 1:nu) - eos_shift
885 IF (
ALLOCATED(gs_mos(1)%evals_virt))
THEN
886 gs_mos(1)%evals_virt(:) = gs_mos(1)%evals_virt(:) + ev_shift
888 IF (
ALLOCATED(gs_mos(2)%evals_virt))
THEN
889 gs_mos(2)%evals_virt(1:nos) = gs_mos(2)%evals_virt(1:nos) + eos_shift
890 gs_mos(2)%evals_virt(nos + 1:) = gs_mos(2)%evals_virt(nos + 1:) + ev_shift
893 IF (
ALLOCATED(gs_mos(1)%evals_virt))
THEN
894 gs_mos(1)%evals_virt(1:nos) = gs_mos(1)%evals_virt(1:nos) + eos_shift
895 gs_mos(1)%evals_virt(nos + 1:) = gs_mos(1)%evals_virt(nos + 1:) + ev_shift
897 IF (
ALLOCATED(gs_mos(2)%evals_occ))
THEN
898 gs_mos(2)%evals_occ(nl + 1:nu) = gs_mos(2)%evals_occ(nl + 1:nu) - eos_shift
900 IF (
ALLOCATED(gs_mos(2)%evals_virt))
THEN
901 gs_mos(2)%evals_virt(:) = gs_mos(2)%evals_virt(:) + ev_shift
906 CALL timestop(handle)
908 END SUBROUTINE ev_shift_operator
931 TYPE(
cp_fm_type),
DIMENSION(:, :),
INTENT(inout) :: evects
932 REAL(kind=
dp),
DIMENSION(:),
INTENT(inout) :: evals
934 DIMENSION(SIZE(evects, 1)),
INTENT(in) :: gs_mos
935 INTEGER,
INTENT(in) :: log_unit
937 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tddfpt_guess_vectors'
939 CHARACTER(len=5) :: spin_label
940 INTEGER :: handle, imo_occ, imo_virt, ind, ispin, &
941 istate, jspin, nspins, nstates, &
942 nstates_occ_virt_alpha, &
944 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: inds
945 INTEGER,
DIMENSION(maxspins) :: nmo_occ_avail, nmo_occ_selected, &
947 REAL(kind=
dp) :: e_occ
948 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: e_virt_minus_occ
951 CALL timeset(routinen, handle)
953 nspins =
SIZE(evects, 1)
954 nstates =
SIZE(evects, 2)
956 IF (debug_this_module)
THEN
957 cpassert(nstates > 0)
958 cpassert(nspins == 1 .OR. nspins == 2)
963 nmo_occ_avail(ispin) =
SIZE(gs_mos(ispin)%evals_occ)
966 nmo_occ_selected(ispin) = min(nmo_occ_avail(ispin), nstates)
967 nmo_virt_selected(ispin) = min(
SIZE(gs_mos(ispin)%evals_virt), nstates)
969 CALL cp_fm_get_info(gs_mos(ispin)%mos_occ, matrix_struct=fm_struct_evects(ispin)%struct)
974 nstates_selected = dot_product(nmo_occ_selected(1:nspins), nmo_virt_selected(1:nspins))
976 ALLOCATE (inds(nstates_selected))
977 ALLOCATE (e_virt_minus_occ(nstates_selected))
981 DO imo_occ = 1, nmo_occ_selected(ispin)
983 e_occ = gs_mos(ispin)%evals_occ(nmo_occ_avail(ispin) - imo_occ + 1)
985 DO imo_virt = 1, nmo_virt_selected(ispin)
987 e_virt_minus_occ(istate) = gs_mos(ispin)%evals_virt(imo_virt) - e_occ
992 IF (debug_this_module)
THEN
993 cpassert(istate == nstates_selected)
996 CALL sort(e_virt_minus_occ, nstates_selected, inds)
998 IF (nspins == 1)
THEN
1003 nstates_occ_virt_alpha = nmo_occ_selected(1)*nmo_virt_selected(1)
1004 IF (log_unit > 0)
THEN
1005 WRITE (log_unit,
"(1X,A)")
"", &
1006 "-------------------------------------------------------------------------------", &
1007 "- TDDFPT Initial Guess -", &
1008 "-------------------------------------------------------------------------------"
1009 WRITE (log_unit,
'(T11,A)')
"State Occupied -> Virtual Excitation"
1010 WRITE (log_unit,
'(T11,A)')
"number orbital orbital energy (eV)"
1011 WRITE (log_unit,
'(1X,79("-"))')
1014 DO istate = 1, nstates
1015 IF (
ASSOCIATED(evects(1, istate)%matrix_struct))
THEN
1017 WRITE (log_unit,
'(T7,I8,T28,A19,T60,F14.5)') &
1018 istate,
"*** restarted ***", evals(istate)*
evolt
1020 ind = inds(istate) - 1
1021 IF (nspins > 1)
THEN
1022 IF (ind < nstates_occ_virt_alpha)
THEN
1024 spin_label =
'(alp)'
1027 ind = ind - nstates_occ_virt_alpha
1028 spin_label =
'(bet)'
1032 imo_occ = nmo_occ_avail(ispin) - ind/nmo_virt_selected(ispin)
1033 imo_virt = mod(ind, nmo_virt_selected(ispin)) + 1
1034 evals(istate) = e_virt_minus_occ(istate)
1037 WRITE (log_unit,
'(T7,I8,T24,I8,T37,A5,T45,I8,T54,A5,T60,F14.5)') &
1038 istate, imo_occ, spin_label, nmo_occ_avail(ispin) + imo_virt, spin_label, e_virt_minus_occ(istate)*
evolt
1040 DO jspin = 1, nspins
1042 CALL cp_fm_create(evects(jspin, istate), fm_struct_evects(jspin)%struct)
1045 IF (jspin == ispin) &
1046 CALL cp_fm_to_fm(gs_mos(ispin)%mos_virt, evects(ispin, istate), &
1047 ncol=1, source_start=imo_virt, target_start=imo_occ)
1052 IF (log_unit > 0)
THEN
1054 WRITE (log_unit,
"(1X,A)") &
1055 "-------------------------------------------------------------------------------"
1058 DEALLOCATE (e_virt_minus_occ)
1061 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,...
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_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
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
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, 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, 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, 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)
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)
routine where the real calculations are made: the KS matrix is calculated
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.
pure integer(kind=int_8) function, public tddfpt_total_number_of_states(gs_mos)
Compute the number of possible singly excited states (occ -> virt)
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_guess_vectors(evects, evals, gs_mos, log_unit)
Generate missed guess vectors.
subroutine, public tddfpt_oecorr(qs_env, gs_mos, matrix_ks_oep)
Callculate orbital corrected KS matrix for TDDFPT.
subroutine, public tddfpt_init_mos(qs_env, gs_mos, iounit)
Prepare MOs for TDDFPT Calculations.
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...
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...
stores all the informations relevant to an mpi environment
Ground state molecular orbitals.