tdhf_mrsf_energy Subroutine

public subroutine tdhf_mrsf_energy(infos)

Uses

  • proc~~tdhf_mrsf_energy~~UsesGraph proc~tdhf_mrsf_energy tdhf_mrsf_energy module~basis_tools basis_tools proc~tdhf_mrsf_energy->module~basis_tools module~int2_compute int2_compute proc~tdhf_mrsf_energy->module~int2_compute module~io_constants io_constants proc~tdhf_mrsf_energy->module~io_constants module~mathlib mathlib proc~tdhf_mrsf_energy->module~mathlib module~messages messages proc~tdhf_mrsf_energy->module~messages module~oqp_linalg oqp_linalg proc~tdhf_mrsf_energy->module~oqp_linalg module~oqp_tagarray_driver oqp_tagarray_driver proc~tdhf_mrsf_energy->module~oqp_tagarray_driver module~precision precision proc~tdhf_mrsf_energy->module~precision module~printing printing proc~tdhf_mrsf_energy->module~printing module~strings strings proc~tdhf_mrsf_energy->module~strings module~tdhf_lib tdhf_lib proc~tdhf_mrsf_energy->module~tdhf_lib module~tdhf_mrsf_lib tdhf_mrsf_lib proc~tdhf_mrsf_energy->module~tdhf_mrsf_lib module~tdhf_sf_lib tdhf_sf_lib proc~tdhf_mrsf_energy->module~tdhf_sf_lib module~types types proc~tdhf_mrsf_energy->module~types module~util util proc~tdhf_mrsf_energy->module~util module~basis_tools->module~io_constants module~basis_tools->module~precision iso_fortran_env iso_fortran_env module~basis_tools->iso_fortran_env module~atomic_structure_m atomic_structure_m module~basis_tools->module~atomic_structure_m module~constants constants module~basis_tools->module~constants module~parallel parallel module~basis_tools->module~parallel module~int2_compute->module~basis_tools module~int2_compute->module~messages module~int2_compute->module~precision module~int2_compute->module~atomic_structure_m module~int2_pairs int2_pairs module~int2_compute->module~int2_pairs module~int2e_libint int2e_libint module~int2_compute->module~int2e_libint module~int2e_rys int2e_rys module~int2_compute->module~int2e_rys module~int2_compute->module~parallel module~mathlib->module~oqp_linalg module~mathlib->module~precision module~messages->module~io_constants module~messages->module~precision comm_IOFILE comm_IOFILE module~messages->comm_IOFILE comm_PAR comm_PAR module~messages->comm_PAR module~blas_wrap blas_wrap module~oqp_linalg->module~blas_wrap module~lapack_wrap lapack_wrap module~oqp_linalg->module~lapack_wrap iso_c_binding iso_c_binding module~oqp_tagarray_driver->iso_c_binding tagarray tagarray module~oqp_tagarray_driver->tagarray module~precision->iso_fortran_env module~printing->module~precision module~strings->iso_c_binding module~tdhf_lib->module~basis_tools module~tdhf_lib->module~int2_compute module~tdhf_lib->module~oqp_linalg module~tdhf_lib->module~precision module~tdhf_mrsf_lib->module~basis_tools module~tdhf_mrsf_lib->module~int2_compute module~tdhf_mrsf_lib->module~oqp_linalg module~tdhf_mrsf_lib->module~precision module~tdhf_sf_lib->module~oqp_linalg module~types->module~basis_tools module~types->module~precision module~types->iso_c_binding module~types->module~atomic_structure_m module~functionals functionals module~types->module~functionals module~types->module~parallel module~types->tagarray module~util->module~precision module~atomic_structure_m->iso_c_binding module~blas_wrap->module~messages module~blas_wrap->module~precision module~mathlib_types mathlib_types module~blas_wrap->module~mathlib_types module~constants->module~precision module~functionals->module~precision module~functionals->iso_c_binding xc_f03_lib_m xc_f03_lib_m module~functionals->xc_f03_lib_m module~int2_pairs->module~precision module~int2e_libint->module~precision module~int2e_libint->iso_c_binding module~int2e_libint->module~constants module~int2e_libint->module~int2_pairs module~libint_f libint_f module~int2e_libint->module~libint_f module~int2e_rys->module~basis_tools module~int2e_rys->module~precision module~int2e_rys->module~constants module~lapack_wrap->module~messages module~lapack_wrap->module~precision module~lapack_wrap->module~mathlib_types module~parallel->module~precision module~parallel->iso_c_binding module~parallel->iso_fortran_env mpi mpi module~parallel->mpi module~libint_f->iso_c_binding

Arguments

Type IntentOptional Attributes Name
type(information), intent(inout), target :: infos

Calls

proc~~tdhf_mrsf_energy~~CallsGraph proc~tdhf_mrsf_energy tdhf_mrsf_energy interface~data_has_tags data_has_tags proc~tdhf_mrsf_energy->interface~data_has_tags interface~show_message show_message proc~tdhf_mrsf_energy->interface~show_message interface~tagarray_get_data tagarray_get_data proc~tdhf_mrsf_energy->interface~tagarray_get_data interface~unpack_matrix unpack_matrix proc~tdhf_mrsf_energy->interface~unpack_matrix none~clean~18 int2_compute_t%clean proc~tdhf_mrsf_energy->none~clean~18 none~init~16 int2_compute_t%init proc~tdhf_mrsf_energy->none~init~16 none~run~6 int2_compute_t%run proc~tdhf_mrsf_energy->none~run~6 none~set_screening~2 int2_compute_t%set_screening proc~tdhf_mrsf_energy->none~set_screening~2 proc~get_mrsf_transition_density get_mrsf_transition_density proc~tdhf_mrsf_energy->proc~get_mrsf_transition_density proc~get_mrsf_transitions get_mrsf_transitions proc~tdhf_mrsf_energy->proc~get_mrsf_transitions proc~get_transition_density get_transition_density proc~tdhf_mrsf_energy->proc~get_transition_density proc~get_transition_dipole get_transition_dipole proc~tdhf_mrsf_energy->proc~get_transition_dipole proc~get_transitions get_transitions proc~tdhf_mrsf_energy->proc~get_transitions proc~iatogen iatogen proc~tdhf_mrsf_energy->proc~iatogen proc~inivec inivec proc~tdhf_mrsf_energy->proc~inivec proc~measure_time measure_time proc~tdhf_mrsf_energy->proc~measure_time proc~mntoia mntoia proc~tdhf_mrsf_energy->proc~mntoia proc~mrinivec mrinivec proc~tdhf_mrsf_energy->proc~mrinivec proc~mrsfcbc mrsfcbc proc~tdhf_mrsf_energy->proc~mrsfcbc proc~mrsfesum mrsfesum proc~tdhf_mrsf_energy->proc~mrsfesum proc~mrsfmntoia mrsfmntoia proc~tdhf_mrsf_energy->proc~mrsfmntoia proc~mrsfqroesum mrsfqroesum proc~tdhf_mrsf_energy->proc~mrsfqroesum proc~oqp_dgemm_i64 oqp_dgemm_i64 proc~tdhf_mrsf_energy->proc~oqp_dgemm_i64 proc~orthogonal_transform orthogonal_transform proc~tdhf_mrsf_energy->proc~orthogonal_transform proc~orthogonal_transform_sym orthogonal_transform_sym proc~tdhf_mrsf_energy->proc~orthogonal_transform_sym proc~print_module_info print_module_info proc~tdhf_mrsf_energy->proc~print_module_info proc~print_results~2 print_results proc~tdhf_mrsf_energy->proc~print_results~2 proc~rpaechk rpaechk proc~tdhf_mrsf_energy->proc~rpaechk proc~rpaeig rpaeig proc~tdhf_mrsf_energy->proc~rpaeig proc~rpanewb rpanewb proc~tdhf_mrsf_energy->proc~rpanewb proc~rpaprint rpaprint proc~tdhf_mrsf_energy->proc~rpaprint proc~rparedms rparedms proc~tdhf_mrsf_energy->proc~rparedms proc~rpavnorm rpavnorm proc~tdhf_mrsf_energy->proc~rpavnorm proc~sfqvec sfqvec proc~tdhf_mrsf_energy->proc~sfqvec proc~sfresvec sfresvec proc~tdhf_mrsf_energy->proc~sfresvec proc~trfrmb trfrmb proc~tdhf_mrsf_energy->proc~trfrmb remove_records remove_records proc~tdhf_mrsf_energy->remove_records reserve_data reserve_data proc~tdhf_mrsf_energy->reserve_data proc~unpack_f90 UNPACK_F90 interface~unpack_matrix->proc~unpack_f90 none~clean~5 int2_pair_storage%clean none~clean~18->none~clean~5 none~alloc int2_pair_storage%alloc none~init~16->none~alloc none~compute int2_pair_storage%compute none~init~16->none~compute none~init_shell_centers basis_set%init_shell_centers none~init~16->none~init_shell_centers none~init~14 par_env_t%init none~init~16->none~init~14 none~set int2_cutoffs_t%set none~init~16->none~set proc~libint_static_init libint_static_init none~init~16->proc~libint_static_init none~run~6->interface~show_message none~run_cam int2_compute_t%run_cam none~run~6->none~run_cam none~run_generic int2_compute_t%run_generic none~run~6->none~run_generic proc~ints_exchange ints_exchange none~set_screening~2->proc~ints_exchange proc~get_mrsf_transition_density->interface~show_message proc~get_trans_den get_trans_den proc~get_mrsf_transition_density->proc~get_trans_den proc~get_transition_density->interface~show_message proc~get_transition_density->proc~iatogen proc~get_transition_dipole->interface~show_message proc~get_transition_dipole->proc~orthogonal_transform interface~pack_matrix pack_matrix proc~get_transition_dipole->interface~pack_matrix proc~atomic_structure_center atomic_structure%atomic_structure_center proc~get_transition_dipole->proc~atomic_structure_center proc~multipole_integrals multipole_integrals proc~get_transition_dipole->proc~multipole_integrals proc~symmetrize_matrix symmetrize_matrix proc~get_transition_dipole->proc~symmetrize_matrix proc~traceprod_sym_packed traceprod_sym_packed proc~get_transition_dipole->proc~traceprod_sym_packed proc~mntoia->proc~oqp_dgemm_i64 proc~mrsfcbc->interface~show_message proc~mrsfcbc->proc~oqp_dgemm_i64 proc~mrsfesum->interface~show_message proc~mrsfesum->proc~oqp_dgemm_i64 proc~mrsfmntoia->interface~show_message proc~mrsfmntoia->proc~oqp_dgemm_i64 proc~oqp_dgemv_i64 oqp_dgemv_i64 proc~mrsfmntoia->proc~oqp_dgemv_i64 proc~oqp_dgemm_i64->interface~show_message dgemm dgemm proc~oqp_dgemm_i64->dgemm proc~orthogonal_transform->interface~show_message proc~orthogonal_transform->proc~oqp_dgemm_i64 proc~orthogonal_transform_sym->interface~show_message proc~orthogonal_transform_sym->proc~oqp_dgemm_i64 proc~oqp_dsymm_i64 oqp_dsymm_i64 proc~orthogonal_transform_sym->proc~oqp_dsymm_i64 proc~oqp_dtpttr_i64 oqp_dtpttr_i64 proc~orthogonal_transform_sym->proc~oqp_dtpttr_i64 proc~oqp_dtrttp_i64 oqp_dtrttp_i64 proc~orthogonal_transform_sym->proc~oqp_dtrttp_i64 proc~rpaeig->proc~oqp_dgemm_i64 proc~diag_symm_packed diag_symm_packed proc~rpaeig->proc~diag_symm_packed proc~rpaeig->proc~oqp_dtrttp_i64 proc~rparedms->proc~oqp_dgemm_i64 proc~sfresvec->proc~oqp_dgemm_i64 proc~trfrmb->proc~oqp_dgemm_i64 proc~pack_f90 PACK_F90 interface~pack_matrix->proc~pack_f90 mpi_comm_rank mpi_comm_rank none~init~14->mpi_comm_rank mpi_comm_size mpi_comm_size none~init~14->mpi_comm_size none~run_cam->none~run_generic none~run_generic->interface~show_message none~run_generic->proc~ints_exchange libint2_cleanup_eri libint2_cleanup_eri none~run_generic->libint2_cleanup_eri libint2_init_eri libint2_init_eri none~run_generic->libint2_init_eri none~clean~2 int2_rys_data_t%clean none~run_generic->none~clean~2 none~init~2 int2_rys_data_t%init none~run_generic->none~init~2 none~screen_ij int2_compute_data_t%screen_ij none~run_generic->none~screen_ij none~screen_ijkl int2_compute_data_t%screen_ijkl none~run_generic->none~screen_ijkl none~set_ids int2_rys_data_t%set_ids none~run_generic->none~set_ids parallel_start parallel_start none~run_generic->parallel_start parallel_stop parallel_stop none~run_generic->parallel_stop proc~genr22 genr22 none~run_generic->proc~genr22 proc~int2_rys_compute int2_rys_compute none~run_generic->proc~int2_rys_compute proc~libint_compute_eri libint_compute_eri none~run_generic->proc~libint_compute_eri proc~libint_print_eri libint_print_eri none~run_generic->proc~libint_print_eri proc~rys_print_eri rys_print_eri none~run_generic->proc~rys_print_eri update update none~run_generic->update proc~to_upper to_upper proc~atomic_structure_center->proc~to_upper proc~diag_symm_packed->interface~show_message dspev dspev proc~diag_symm_packed->dspev dspevx dspevx proc~diag_symm_packed->dspevx proc~ints_exchange->interface~show_message proc~ints_exchange->none~alloc proc~ints_exchange->none~compute proc~ints_exchange->none~set proc~ints_exchange->libint2_cleanup_eri proc~ints_exchange->libint2_init_eri proc~ints_exchange->none~clean~2 proc~ints_exchange->none~init~2 proc~ints_exchange->none~set_ids proc~ints_exchange->proc~genr22 proc~ints_exchange->proc~int2_rys_compute proc~ints_exchange->proc~libint_compute_eri libint2_static_init libint2_static_init proc~libint_static_init->libint2_static_init proc~multipole_integrals->interface~show_message interface~bas_norm_matrix bas_norm_matrix proc~multipole_integrals->interface~bas_norm_matrix none~alloc~2 shpair_t%alloc proc~multipole_integrals->none~alloc~2 none~fetch_by_id shell_t%fetch_by_id proc~multipole_integrals->none~fetch_by_id none~shell_pair shpair_t%shell_pair proc~multipole_integrals->none~shell_pair proc~comp_allmult_int1_prim comp_allmult_int1_prim proc~multipole_integrals->proc~comp_allmult_int1_prim proc~print_sym_labeled print_sym_labeled proc~multipole_integrals->proc~print_sym_labeled proc~update_triang_matrix update_triang_matrix proc~multipole_integrals->proc~update_triang_matrix proc~oqp_dgemv_i64->interface~show_message dgemv dgemv proc~oqp_dgemv_i64->dgemv proc~oqp_dsymm_i64->interface~show_message dsymm dsymm proc~oqp_dsymm_i64->dsymm proc~oqp_dtpttr_i64->interface~show_message dtpttr dtpttr proc~oqp_dtpttr_i64->dtpttr proc~oqp_dtrttp_i64->interface~show_message dtrttp dtrttp proc~oqp_dtrttp_i64->dtrttp proc~unpack_f90->interface~show_message proc~unpack_f90->proc~oqp_dtpttr_i64 unused_dummy unused_dummy none~screen_ij->unused_dummy none~screen_ijkl->unused_dummy proc~mulquadgausshermite mulQuadGaussHermite proc~comp_allmult_int1_prim->proc~mulquadgausshermite fgrid fgrid proc~genr22->fgrid rfinc rfinc proc~genr22->rfinc rmr rmr proc~genr22->rmr xgrid xgrid proc~genr22->xgrid none~evaluate rys_root_t%evaluate proc~int2_rys_compute->none~evaluate ab_x ab_x proc~libint_compute_eri->ab_x ab_y ab_y proc~libint_compute_eri->ab_y ab_z ab_z proc~libint_compute_eri->ab_z alpha1over_zetapluseta alpha1over_zetapluseta proc~libint_compute_eri->alpha1over_zetapluseta alpha1rho_over_zeta2 alpha1rho_over_zeta2 proc~libint_compute_eri->alpha1rho_over_zeta2 alpha2over_zetapluseta alpha2over_zetapluseta proc~libint_compute_eri->alpha2over_zetapluseta alpha2rho_over_zeta2 alpha2rho_over_zeta2 proc~libint_compute_eri->alpha2rho_over_zeta2 alpha3over_zetapluseta alpha3over_zetapluseta proc~libint_compute_eri->alpha3over_zetapluseta alpha3rho_over_eta2 alpha3rho_over_eta2 proc~libint_compute_eri->alpha3rho_over_eta2 alpha4over_zetapluseta alpha4over_zetapluseta proc~libint_compute_eri->alpha4over_zetapluseta alpha4rho_over_eta2 alpha4rho_over_eta2 proc~libint_compute_eri->alpha4rho_over_eta2 cd_x cd_x proc~libint_compute_eri->cd_x cd_y cd_y proc~libint_compute_eri->cd_y cd_z cd_z proc~libint_compute_eri->cd_z f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_0 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_0 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_0 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_1 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_1 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_1 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_10 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_10 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_10 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_11 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_11 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_11 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_12 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_12 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_12 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_13 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_13 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_13 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_14 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_14 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_14 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_15 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_15 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_15 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_16 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_16 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_16 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_17 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_17 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_17 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_18 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_18 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_18 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_19 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_19 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_19 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_2 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_2 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_2 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_20 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_20 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_20 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_3 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_3 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_3 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_4 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_4 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_4 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_5 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_5 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_5 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_6 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_6 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_6 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_7 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_7 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_7 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_8 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_8 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_8 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_9 f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_9 proc~libint_compute_eri->f_ab_s___0__s___1___twoprep_s___0__s___1___ab__up_9 proc~libint_compute_eri->fgrid igrid igrid proc~libint_compute_eri->igrid irgrd irgrd proc~libint_compute_eri->irgrd oo2e oo2e proc~libint_compute_eri->oo2e oo2z oo2z proc~libint_compute_eri->oo2z oo2ze oo2ze proc~libint_compute_eri->oo2ze pa_x pa_x proc~libint_compute_eri->pa_x pa_y pa_y proc~libint_compute_eri->pa_y pa_z pa_z proc~libint_compute_eri->pa_z qc_x qc_x proc~libint_compute_eri->qc_x qc_y qc_y proc~libint_compute_eri->qc_y qc_z qc_z proc~libint_compute_eri->qc_z proc~libint_compute_eri->rfinc rho12_over_alpha1 rho12_over_alpha1 proc~libint_compute_eri->rho12_over_alpha1 rho12_over_alpha2 rho12_over_alpha2 proc~libint_compute_eri->rho12_over_alpha2 rho34_over_alpha3 rho34_over_alpha3 proc~libint_compute_eri->rho34_over_alpha3 rho34_over_alpha4 rho34_over_alpha4 proc~libint_compute_eri->rho34_over_alpha4 proc~libint_compute_eri->rmr roe roe proc~libint_compute_eri->roe roz roz proc~libint_compute_eri->roz two_alpha0_bra two_alpha0_bra proc~libint_compute_eri->two_alpha0_bra two_alpha0_ket two_alpha0_ket proc~libint_compute_eri->two_alpha0_ket two_alpha1bra two_alpha1bra proc~libint_compute_eri->two_alpha1bra two_alpha1ket two_alpha1ket proc~libint_compute_eri->two_alpha1ket wp_x wp_x proc~libint_compute_eri->wp_x wp_y wp_y proc~libint_compute_eri->wp_y wp_z wp_z proc~libint_compute_eri->wp_z wq_x wq_x proc~libint_compute_eri->wq_x wq_y wq_y proc~libint_compute_eri->wq_y wq_z wq_z proc~libint_compute_eri->wq_z proc~pack_f90->interface~show_message proc~pack_f90->proc~oqp_dtrttp_i64 none~bf_label basis_set%bf_label proc~print_sym_labeled->none~bf_label abrt abrt proc~mulquadgausshermite->abrt

Called by

proc~~tdhf_mrsf_energy~~CalledByGraph proc~tdhf_mrsf_energy tdhf_mrsf_energy proc~tdhf_mrsf_energy_c tdhf_mrsf_energy_C proc~tdhf_mrsf_energy_c->proc~tdhf_mrsf_energy

Source Code

  subroutine tdhf_mrsf_energy(infos)
    use io_constants, only: iw
    use oqp_tagarray_driver

    use types, only: information
    use strings, only: cstring, fstring
    use basis_tools, only: basis_set
    use messages, only: show_message, with_abort
    use util, only: measure_time

    use precision, only: dp
    use int2_compute, only: int2_compute_t
    use tdhf_mrsf_lib, only: int2_mrsf_data_t
    use tdhf_lib, only: int2_td_data_t
    use tdhf_lib, only: &
      iatogen, mntoia, rparedms, rpaeig, rpavnorm, &
      rpaechk, rpanewb, &
      rpaprint, inivec
    use tdhf_sf_lib, only: sfresvec, sfqvec, sfdmat, trfrmb, &
      get_transition_density, get_transitions, &
      get_transition_dipole, print_results
    use tdhf_mrsf_lib, only: &
      mrinivec, mrsfcbc, mrsfmntoia, mrsfesum, &
      mrsfqroesum, get_mrsf_transitions, &
      get_mrsf_transition_density
    use mathlib, only: orthogonal_transform, orthogonal_transform_sym, &
      unpack_matrix
    use oqp_linalg
    use printing, only: print_module_info

    implicit none

    character(len=*), parameter :: subroutine_name = "tdhf_mrsf_energy"

    type(basis_set), pointer :: basis
    type(information), target, intent(inout) :: infos

    integer :: s_size, ok

    real(kind=dp), allocatable :: scr2(:)
    real(kind=dp), allocatable :: wrk1(:,:), qvec(:,:)
    real(kind=dp), allocatable :: amo(:,:), wrk2(:,:)
    real(kind=dp), allocatable :: squared_S(:)
    real(kind=dp), allocatable :: amb(:,:), apb(:,:)
    real(kind=dp), allocatable, target :: vl(:), vr(:)
    real(kind=dp), pointer :: vl_p(:,:), vr_p(:,:)
    real(kind=dp), allocatable :: xm(:), scr(:)
    real(kind=dp), allocatable :: bvec_mo(:,:), for_trnsf_b_vec(:,:)
    real(kind=dp), allocatable, dimension(:,:) :: fa, fb
    real(kind=dp), allocatable, dimension(:) :: rnorm
    real(kind=dp), allocatable, dimension(:,:,:,:) :: trden
    integer, allocatable, dimension(:,:) :: trans
    real(kind=dp), allocatable, target :: mrsf_density(:,:,:,:)
    real(kind=dp), pointer :: fmrst2(:,:,:,:)
    real(kind=dp), allocatable, target :: fmrq1(:,:,:)
    real(kind=dp), allocatable :: dip(:,:,:), bvec_mo_tmp(:), eex(:)

    integer :: nocca, nvira, noccb, nvirb
    integer :: nbf, nbf2, xvec_dim
    integer :: mxvec, ist, jst, iend, nvec, novec
    integer :: iter, nv, iv, ivec
    integer :: mxiter
    logical :: tamm_dancoff
    integer :: imax
    integer :: ierr
    logical :: converged
    real(kind=dp) :: mxerr, cnvtol, scale_exch
    integer :: maxvec, mrst, nstates, target_state
    logical :: roref = .false.
    logical :: debug_mode

    type(int2_compute_t) :: int2_driver
    type(int2_mrsf_data_t), target :: int2_data_st
    type(int2_td_data_t), target :: int2_data_q

    logical :: dft = .false.
    integer :: scf_type, mol_mult

    ! tagarray
    real(kind=dp), contiguous, pointer :: &
      fock_a(:), dmat_a(:), mo_A(:,:), mo_energy_a(:), &
      fock_b(:), dmat_b(:), mo_b(:,:), mo_energy_b(:), &
      smat(:), ta(:), tb(:), td_t(:,:), bvec_mo_out(:,:), &
      mrsf_energies(:)
    character(len=*), parameter :: tags_alloc(3) = (/ character(len=80) :: &
      OQP_td_bvec_mo, OQP_td_t, OQP_td_energies /)
    character(len=*), parameter :: tags_required(9) = (/ character(len=80) :: &
      OQP_FOCK_A, OQP_DM_A, OQP_E_MO_A, OQP_VEC_MO_A, OQP_FOCK_B, OQP_DM_B, OQP_E_MO_B, OQP_VEC_MO_B, OQP_SM /)

  ! Readings
  ! Files open
  ! 3. LOG: Write: Main output file
    open (unit=iw, file=infos%log_filename, position="append")
  !
    call print_module_info('MRSF_TDHF_Energy','Computing Energy of MRSF-TDDFT')

  ! Load basis set
    basis => infos%basis
    basis%atoms => infos%atoms

   ! Input parameters
    dft = infos%control%hamilton == 20 ! dft or hf
    mrst = infos%tddft%mult
    nstates = infos%tddft%nstate
    target_state = infos%tddft%target_state
    maxvec = infos%tddft%maxvec
    cnvtol = infos%tddft%cnvtol
!   infos%tddft%debug_mode = .True.
    debug_mode = infos%tddft%debug_mode

    mol_mult = infos%mol_prop%mult
    if (mol_mult/=3) call show_message('MRSF-TDDFT are available for ROHF ref.&
        &with ONLY triplet multiplicity(mult=3)',with_abort)

    scf_type = infos%control%scftype
    if (scf_type==3) roref = .true.

    nbf = basis%nbf
    nbf2 = nbf*(nbf+1)/2
    s_size = (basis%nshell**2+basis%nshell)/2

    tamm_dancoff = .true.  ! tamm_dancoff: 0/1 means not doing/doing Tamm/Dancoff run

    nocca = infos%mol_prop%nelec_a
    nvira = nbf-nocca
    noccb = infos%mol_prop%nelec_b
    nvirb = nbf-noccb


    if (mrst==1 .or. mrst==3 ) then
      xvec_dim = nocca*nvirb
    else if (mrst==5) then
      xvec_dim = noccb*nvira
    end if


    if (mrst==1 ) then
      nstates = min(nstates, xvec_dim-1)
      mxvec = min(maxvec*nstates, xvec_dim-1, infos%control%maxit_dav*nstates)
    else if (mrst==3) then
      nstates = min(nstates, xvec_dim-3)
      mxvec = min(maxvec*nstates, xvec_dim-3, infos%control%maxit_dav*nstates)
    else if (mrst==5) then
      nstates = min(nstates, xvec_dim)
      mxvec = min(maxvec*nstates, xvec_dim, infos%control%maxit_dav*nstates)
    end if

    infos%tddft%nstate = nstates
    nvec = min(max(nstates,6), mxvec)

    call infos%dat%remove_records(tags_alloc)

    call infos%dat%reserve_data(OQP_td_bvec_mo, TA_TYPE_REAL64, &
        xvec_dim*nstates, (/xvec_dim, nstates/), comment=OQP_td_bvec_mo_comment)
    call infos%dat%reserve_data(OQP_td_t, TA_TYPE_REAL64, nbf2*2, (/ nbf2, 2 /), comment=OQP_td_t_comment)
    call infos%dat%reserve_data(OQP_td_energies, TA_TYPE_REAL64, nstates, comment=OQP_td_energies_comment)

    call data_has_tags(infos%dat, tags_alloc, module_name, subroutine_name, WITH_ABORT)
    call tagarray_get_data(infos%dat, OQP_td_bvec_mo, bvec_mo_out)

    call tagarray_get_data(infos%dat, OQP_td_t, td_t)
    call tagarray_get_data(infos%dat, OQP_td_energies, mrsf_energies)

    call data_has_tags(infos%dat, tags_required, module_name, subroutine_name, WITH_ABORT)
    call tagarray_get_data(infos%dat, OQP_SM, smat)
    call tagarray_get_data(infos%dat, OQP_FOCK_A, fock_a)
    call tagarray_get_data(infos%dat, OQP_FOCK_B, fock_b)
    call tagarray_get_data(infos%dat, OQP_DM_A, dmat_a)
    call tagarray_get_data(infos%dat, OQP_DM_B, dmat_b)
    call tagarray_get_data(infos%dat, OQP_E_MO_A, mo_energy_a)
    call tagarray_get_data(infos%dat, OQP_E_MO_B, mo_energy_b)
    call tagarray_get_data(infos%dat, OQP_VEC_MO_A, mo_a)
    call tagarray_get_data(infos%dat, OQP_VEC_MO_B, mo_b)

  ! Allocate temporary matrices for diagonalization
    allocate (fa(nbf,nbf), &
              fb(nbf,nbf), &
              stat=ok)
    if( ok/=0 ) call show_message('Cannot allocate memory',with_abort)

  ! Allocate TDDFT variables
    allocate(xm(xvec_dim),  &
             bvec_mo(xvec_dim,mxvec), &
             trden(nbf,nbf,nstates,nstates), &
             wrk1(nbf,nbf), &
             wrk2(nbf,nbf), &
             amo(xvec_dim,mxvec), &
             EEX(mxvec), &
             squared_S(nstates), &
             APB(mxvec,mxvec), &
             AMB(mxvec,mxvec), &
             VR(mxvec*mxvec), &
             VL(mxvec*mxvec), &
             for_trnsf_b_vec(mxvec,mxvec), & !
             dip(3,nstates,nstates), &
             scr(nbf2), &
             bvec_mo_tmp(xvec_dim), &
             scr2(mxvec*mxvec), &
             qvec(xvec_dim,nstates), &
             RNORM(nstates), &
             source=0.0_dp,stat=ok)
    allocate(trans(xvec_dim,2), &
             source=0, stat=ok)
    if( ok/=0 ) call show_message('Cannot allocate memory', with_abort)

    ta => td_t(:,1)
    tb => td_t(:,2)

    if (mrst==1 .or. mrst==3 ) then

      allocate(mrsf_density(nvec,7,nbf,nbf), &
               source=0.0_dp, &
               stat=ok)

    else if( mrst==5  )then

      allocate(fmrq1(nbf,nbf,nvec), &
               source=0.0_dp, &
               stat=ok)

    end if

    if( ok/=0 ) call show_message('Cannot allocate memory', with_abort)

    scale_exch = 1.0_dp
    if (infos%tddft%HFscale == -1.0_dp) &
          infos%tddft%HFscale = infos%dft%HFscale

    if (infos%dft%cam_flag) then
      if (infos%tddft%cam_alpha == -1.0_dp) &
            infos%tddft%cam_alpha = infos%dft%cam_alpha
      infos%tddft%HFscale = infos%tddft%cam_alpha
      if (infos%tddft%cam_beta == -1.0_dp) &
            infos%tddft%cam_beta = infos%dft%cam_beta
      if (infos%tddft%cam_mu == -1.0_dp) &
            infos%tddft%cam_mu = infos%dft%cam_mu
    end if
    if (dft) scale_exch = infos%tddft%HFscale
    ! set spin-pair coupling
    if (infos%tddft%spc_coco==-1.0_dp) &
          infos%tddft%spc_coco = infos%tddft%HFscale
    if (infos%tddft%spc_ovov==-1.0_dp) &
          infos%tddft%spc_ovov = infos%tddft%HFscale
    if (infos%tddft%spc_coov==-1.0_dp) &
          infos%tddft%spc_coov = infos%tddft%HFscale

    if(debug_mode.or..true.)then
      write(*,'(/,5x,"Input parameters:")')
      write(*,'(5x,"Number of states:                 ",1x,I0)') nstates
      write(*,'(5x,"Number of single excitations:     ",1x,I0)') xvec_dim
      write(*,'(5x,"Number of atomic orbitals:        ",1x,I0)') nbf
      write(*,'(5x,"Number of electrons:              ",1x,I0)') nocca+noccb
      write(*,'(5x,"Number of occupied alpha orbitals:",1x,I0)') nocca
      write(*,'(5x,"Number of occupied beta orbitals: ",1x,I0)') noccb
      write(*,'(5x,"Number of virtual alpha orbitals: ",1x,I0)') nvira
      write(*,'(5x,"Number of virtual beta orbitals:  ",1x,I0)') nvirb
      write(*,'(5x,"Maximum vectors:                  ",1x,I0)') mxvec
      write(*,'(5x,"Initial vectors:                  ",1x,I0)') nvec
      write(*, '(/7x,"Fitting parameters for MRSF-TDDFT")')
      if (.not.infos%dft%cam_flag) then
        write(*, '(10x,"Exact HF exchange:")')
        write(*, '(5x,"Reference: |", t20, f6.3, t29, "|")') infos%dft%HFscale
        write(*, '(5x,"Response:  |", t20, f6.3, t29, "|")') infos%tddft%HFscale
      else
        write(*, '(10x,"CAM parametres:")')
        write(*, '(16x,"|   alpha   |    beta   |     mu    |")')
        write(*, '(5x,"Reference: |", t20, f6.3, t29, "|", t32, f6.3, t41, "|", t44, f6.3, t53, "|")') &
           infos%dft%cam_alpha, infos%dft%cam_beta, infos%dft%cam_mu
        write(*, '(5x,"Response:  |", t20, f6.3, t29, "|", t32, f6.3, t41, "|", t44, f6.3, t53, "|")') &
           infos%tddft%cam_alpha, infos%tddft%cam_beta, infos%tddft%cam_mu
      end if
      write(*, '(10x,"Spin-pair coupling parametres:")')
      write(*, '(16x,"|   CO-CO   |   OV-OV   |   CO-OV   |")')
      write(*, '(16x,"|", t20, f6.3, t29, "|", t32, f6.3, t41, "|", t44, f6.3, t53, "|")') &
         infos%tddft%spc_coco, infos%tddft%spc_ovov, infos%tddft%spc_coov
    end if

    write(*,'(/,5x,46("="))')
    if (mrst==1) write(*,'(  5X,"Davidson algorithm for Singlet response states")')
    if (mrst==3) write(*,'(  5x,"Davidson algorithm for Triplet response states")')
    if (mrst==5) write(*,'(  5x,"Davidson algorithm for Quintet response states")')
    write(*,'(5x,46("="))')

    ! Initialize ERI (Electron Repulsion Integrals) calculations
    call int2_driver%init(basis, infos)
    call int2_driver%set_screening()
    call flush(iw)

  ! Prepare for ROHF
    if (roref) then
  !   Alapha
      call orthogonal_transform_sym(nbf, nbf, fock_a, mo_a, nbf, scr)
      call unpack_matrix(scr,fa)

  !   Beta
      call orthogonal_transform_sym(nbf, nbf, fock_b, mo_b, nbf, scr)
      call unpack_matrix(scr,fb)
    end if

  ! Construct TD trial vector
    if (mrst==1 .or. mrst==3) then

      call mrinivec(infos, mo_energy_a, mo_energy_a, bvec_mo, xm, nvec)

    else if (mrst==5) then

      call inivec(mo_energy_a,mo_energy_a,bvec_mo,xm,noccb,nocca,nvec)

    end if

    ist = 1
    iend = nvec
    iter = 0
    mxiter = infos%control%maxit_dav
    ierr = 0

    do iter = 1, mxiter
      nv = iend-ist+1

      if( mrst==1 .or. mrst==3 ) then

        mrsf_density = 0.0_dp   ! bo2v, bo1v, bco1, bco2, o21v, co12, ball

      else if( mrst==5 ) then

        fmrq1 = 0.0_dp

      end if

      do ivec = ist, iend

        iv = ivec-ist+1

        if (mrst==1 .or. mrst==3) then

          call iatogen(bvec_mo(:,ivec), wrk1, nocca, noccb)
          call mrsfcbc(infos, mo_a, mo_b, wrk1, mrsf_density(iv,:,:,:))

        else if (mrst==5) then

          call iatogen(bvec_mo(:,ivec), wrk1, noccb, nocca)
          call orthogonal_transform('t', nbf, mo_a, wrk1, fmrq1(:,:,iv), wrk2)

        end if

      end do

      if (mrst==1 .or. mrst==3) then

        int2_data_st = int2_mrsf_data_t( &
          d3 = mrsf_density(:iv,:,:,:), &
          tamm_dancoff = tamm_dancoff, &
          scale_exchange = scale_exch, &
          scale_coulomb = scale_exch)

        call int2_driver%run( &
          int2_data_st, &
          cam = dft.and.infos%dft%cam_flag, &
          alpha = infos%tddft%cam_alpha, &
          alpha_coulomb = infos%tddft%cam_alpha, &
          beta = infos%tddft%cam_beta, &
          beta_coulomb = infos%tddft%cam_beta, &
          mu = infos%tddft%cam_mu)

        fmrst2 => int2_data_st%f3(:,:,:,:,1) ! ado2v, ado1v, adco1, adco2, ao21v, aco12, agdlr

        ! Scaling factor if triplet
        if (mrst==3) fmrst2(:,1:6,:,:) = -fmrst2(:,1:6,:,:)

        ! Spin pair coupling
        if (infos%tddft%spc_coco /= infos%tddft%hfscale) &
           fmrst2(:,6,:,:) = fmrst2(:,6,:,:) * infos%tddft%spc_coco / infos%tddft%hfscale
        if (infos%tddft%spc_ovov /= infos%tddft%hfscale) &
           fmrst2(:,5,:,:) = fmrst2(:,5,:,:) * infos%tddft%spc_ovov / infos%tddft%hfscale
        if (infos%tddft%spc_coov /= infos%tddft%hfscale) &
           fmrst2(:,1:4,:,:) = fmrst2(:,1:4,:,:) * infos%tddft%spc_coov / infos%tddft%hfscale

      else if (mrst==5) then

        int2_data_q = int2_td_data_t( &
          d2=fmrq1(:,:,:iv), &
          int_apb = .false., &
          int_amb = .false., &
          tamm_dancoff = tamm_dancoff, &
          scale_exchange = scale_exch)
        call int2_driver%run( &
          int2_data_q, &
          cam = dft.and.infos%dft%cam_flag, &
          alpha = infos%tddft%cam_alpha, &
          beta = infos%tddft%cam_beta,&
          mu = infos%tddft%cam_mu)

      end if

      do ivec = ist, iend

        iv = ivec-ist+1

        if (mrst==1 .or. mrst==3) then

          ! Product (A-B)*X
          call mrsfmntoia(infos, fmrst2(iv,:,:,:), amo, mo_a, mo_b, ivec)

          call iatogen(bvec_mo(:,ivec), wrk1, nocca, noccb)

          call mrsfesum(infos, wrk1, fa, fb, amo, ivec)

        else if( mrst==5 )then

          call mntoia(int2_data_q%amb(:,:,iv,1), amo(:,ivec), mo_a, mo_b, noccb, nocca)

          ! Z(I+,A-)
          call iatogen(bvec_mo(:,ivec),wrk1,noccb,nocca)

          ! FB(I+,J+)*Z(J+,A-)
          call dgemm('n','n',noccb,nbf,noccb, &
                     1.0_dp,fb,nbf, &
                            wrk1,nbf, &
                     0.0_dp,wrk2,noccb)

          ! Z(I+,B-)*FA(B-,A-)
          call dgemm('n','n',noccb,nbf,nbf, &
                     1.0_dp,wrk1,nbf, &
                            fa,nbf, &
                    -1.0_dp,wrk2,noccb)

          call mrsfqroesum(wrk2,amo, &
                           nocca,noccb,nbf,ivec)

        end if
      end do

      vl_p(1:nvec, 1:nvec) => vl(1:nvec*nvec)
      vr_p(1:nvec, 1:nvec) => vr(1:nvec*nvec)
      call rparedms(bvec_mo,amo,amo,apb,amb,nvec,tamm_dancoff=.true.)
      call rpaeig(eex,vl_p,vr_p,apb,amb,scr2,tamm_dancoff=.true.)
      call rpavnorm(vr_p,vl_p,tamm_dancoff=.true.)
      call rpaechk(eex,nvec,nstates,imax,tamm_dancoff=.true.)

      for_trnsf_b_vec = vr_p
      call sfresvec(qvec,bvec_mo,amo,vr_p,eex,nvec,rnorm,nstates)
      call sfqvec(qvec,xm,eex,nstates)
      call rpaprint(eex, rnorm, cnvtol, iter, imax, nstates, do_neg=.true.)

      mxerr = maxval(rnorm)

!     Check convergence
      converged = mxerr<=cnvtol
      if (converged) exit

!     No space left for new vectors, exit
      if (nvec==mxvec) ierr = 1
      if (ierr/=0) exit

      call rpanewb(nstates,bvec_mo,qvec,novec,nvec,ierr,tamm_dancoff=.true.)

  !   ierr=1 nvec over mxvec: not converged case
      if (ierr/=0) exit

      ist = novec+1
      iend = nvec

    end do

    if (iter >= mxiter .and. .not. converged) ierr = -1

    select case (ierr)
    case (-1)
      write(*,'(/,2X,"MRSF-TD-DFT energies NOT CONVERGED after ",I4," iterations"/)') mxiter
      infos%mol_energy%Davidson_converged=.false.
    case (0)
      write(*,'(/,2X,"MRSF-TD-DFT energies converged in ",I4," iterations"/)') iter
      infos%mol_energy%Davidson_converged=.true.
    case (1)
      write(*,'(/,2X,"..something is wrong.. nvec = mxvec")')
      infos%mol_energy%Davidson_converged=.false.
    case (2)
      write(*,'(/,2x,"..something is wrong..  nvec > mxvec")')
      write(*,'(3x,"nvec/mxvec =",I4,"/",I4)') nvec, mxvec
      infos%mol_energy%Davidson_converged=.false.
    case (3)
      write(*,'(/,2x,"..something is wrong.. No vectors were added")')
      infos%mol_energy%Davidson_converged=.false.
    end select
    call flush(iw)

    call trfrmb(bvec_mo, for_trnsf_b_vec, nvec, nstates)

    select case (mrst)
      case(1)
        do ist = 1, nstates
          do jst = ist, nstates
            call get_mrsf_transition_density(infos,trden(:,:,ist,jst), bvec_mo, ist, jst)
          end do
        end do
        squared_S(:) = 0.0_dp
        call get_mrsf_transitions(trans, nocca, noccb, nbf)
        write(*,'(/,2x,35("="),/,2x,&
            &"Spin-adapted spin-flip excitations",/,2x,35("="))')
      case(3)
        do ist = 1, nstates
          do jst = ist, nstates
            call get_mrsf_transition_density(infos, trden(:,:,ist,jst), bvec_mo, ist, jst)
          end do
        end do
        squared_S(:) = 2.0_dp
        call get_mrsf_transitions(trans, nocca, noccb, nbf)
        write(*,'(/,2x,35("="),/,2x,&
            &"Spin-adapted spin-flip excitations",/,2x,35("="))')
      case(5)
        call get_transition_density(trden, bvec_mo, nbf, noccb, nocca, nstates)
        squared_S(:) = 6.0_dp
        call get_transitions(trans, noccb, nocca, nbf)
        write(*,'(/,2x,35("="),/,2x,&
            &"Beta -> Alpha spin-flip excitations",/,2x,35("="))')
      case default
        error stop "Unknown mrst value"
    end select

    call get_transition_dipole(basis, dip, mo_a, trden, nstates)

    mrsf_energies = eex(1:nstates)
    bvec_mo_out = bvec_mo(:,1:nstates)
    infos%mol_energy%excited_energy = mrsf_energies(infos%tddft%target_state)
    call print_results(infos, bvec_mo, eex, trans, dip, squared_S, nstates)
    call flush(iw)

    call int2_driver%clean()

    call measure_time(print_total=1, log_unit=iw)
    close(iw)

  end subroutine tdhf_mrsf_energy