scf_driver Subroutine

public subroutine scf_driver(basis, infos, molGrid)

Uses

  • proc~~scf_driver~~UsesGraph proc~scf_driver scf_driver module~basis_tools basis_tools proc~scf_driver->module~basis_tools module~constants constants proc~scf_driver->module~constants module~dft dft proc~scf_driver->module~dft module~guess guess proc~scf_driver->module~guess module~int2_compute int2_compute proc~scf_driver->module~int2_compute module~io_constants io_constants proc~scf_driver->module~io_constants module~mathlib mathlib proc~scf_driver->module~mathlib module~messages messages proc~scf_driver->module~messages module~mod_dft_molgrid mod_dft_molgrid proc~scf_driver->module~mod_dft_molgrid module~oqp_tagarray_driver oqp_tagarray_driver proc~scf_driver->module~oqp_tagarray_driver module~precision precision proc~scf_driver->module~precision module~printing printing proc~scf_driver->module~printing module~scf_converger scf_converger proc~scf_driver->module~scf_converger module~types types proc~scf_driver->module~types module~util util proc~scf_driver->module~util module~basis_tools->module~constants 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~parallel parallel module~basis_tools->module~parallel module~constants->module~precision module~dft->module~basis_tools module~dft->module~io_constants module~dft->module~messages module~dft->module~mod_dft_molgrid module~dft->module~precision module~guess->module~precision module~oqp_linalg oqp_linalg module~guess->module~oqp_linalg 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~precision module~mathlib->module~oqp_linalg 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~mod_dft_molgrid->module~precision module~bragg_slater_radii bragg_slater_radii module~mod_dft_molgrid->module~bragg_slater_radii module~lebedev lebedev module~mod_dft_molgrid->module~lebedev module~mod_grid_storage mod_grid_storage module~mod_dft_molgrid->module~mod_grid_storage 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~scf_converger->module~precision 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~bragg_slater_radii->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~constants module~int2e_libint->module~precision module~int2e_libint->iso_c_binding 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~constants module~int2e_rys->module~precision module~lebedev->module~precision module~mod_grid_storage->module~precision module~blas_wrap blas_wrap module~oqp_linalg->module~blas_wrap module~lapack_wrap lapack_wrap module~oqp_linalg->module~lapack_wrap module~parallel->module~precision module~parallel->iso_c_binding module~parallel->iso_fortran_env mpi mpi module~parallel->mpi module~blas_wrap->module~messages module~blas_wrap->module~precision module~mathlib_types mathlib_types module~blas_wrap->module~mathlib_types module~lapack_wrap->module~messages module~lapack_wrap->module~precision module~lapack_wrap->module~mathlib_types module~libint_f->iso_c_binding

Arguments

Type IntentOptional Attributes Name
type(basis_set), intent(in) :: basis
type(information), intent(inout), target :: infos
type(dft_grid_t), intent(in) :: molGrid

Calls

proc~~scf_driver~~CallsGraph proc~scf_driver scf_driver interface~data_has_tags data_has_tags proc~scf_driver->interface~data_has_tags interface~show_message show_message proc~scf_driver->interface~show_message interface~tagarray_get_data tagarray_get_data proc~scf_driver->interface~tagarray_get_data interface~unpack_matrix unpack_matrix proc~scf_driver->interface~unpack_matrix none~add_data scf_conv%add_data proc~scf_driver->none~add_data none~bcast par_env_t%bcast proc~scf_driver->none~bcast none~clean~18 int2_compute_t%clean proc~scf_driver->none~clean~18 none~get_fock scf_conv_result%get_fock proc~scf_driver->none~get_fock none~init~16 int2_compute_t%init proc~scf_driver->none~init~16 none~run~6 int2_compute_t%run proc~scf_driver->none~run~6 none~set_screening~2 int2_compute_t%set_screening proc~scf_driver->none~set_screening~2 proc~build_pfon_density build_pfon_density proc~scf_driver->proc~build_pfon_density proc~dftexcor dftexcor proc~scf_driver->proc~dftexcor proc~e_charge_repulsion e_charge_repulsion proc~scf_driver->proc~e_charge_repulsion proc~form_rohf_fock form_rohf_fock proc~scf_driver->proc~form_rohf_fock proc~get_ab_initio_density get_ab_initio_density proc~scf_driver->proc~get_ab_initio_density proc~get_ab_initio_orbital get_ab_initio_orbital proc~scf_driver->proc~get_ab_initio_orbital proc~matrix_invsqrt matrix_invsqrt proc~scf_driver->proc~matrix_invsqrt proc~measure_time measure_time proc~scf_driver->proc~measure_time proc~mo_reorder mo_reorder proc~scf_driver->proc~mo_reorder proc~mo_to_ao mo_to_ao proc~scf_driver->proc~mo_to_ao proc~pfon_occupations pfon_occupations proc~scf_driver->proc~pfon_occupations proc~print_mo_range print_mo_range proc~scf_driver->proc~print_mo_range proc~print_scf_energy print_scf_energy proc~scf_driver->proc~print_scf_energy proc~traceprod_sym_packed traceprod_sym_packed proc~scf_driver->proc~traceprod_sym_packed proc~unpack_f90 UNPACK_F90 interface~unpack_matrix->proc~unpack_f90 none~add_data->interface~unpack_matrix interface~pack_matrix pack_matrix none~add_data->interface~pack_matrix proc~antisymmetrize_matrix antisymmetrize_matrix none~add_data->proc~antisymmetrize_matrix proc~oqp_dgemm_i64 oqp_dgemm_i64 none~add_data->proc~oqp_dgemm_i64 proc~oqp_dsymm_i64 oqp_dsymm_i64 none~add_data->proc~oqp_dsymm_i64 none~par_env_t_bcast_byte par_env_t%par_env_t_bcast_byte none~bcast->none~par_env_t_bcast_byte none~par_env_t_bcast_c_bool par_env_t%par_env_t_bcast_c_bool none~bcast->none~par_env_t_bcast_c_bool none~par_env_t_bcast_dp_1d par_env_t%par_env_t_bcast_dp_1d none~bcast->none~par_env_t_bcast_dp_1d none~par_env_t_bcast_dp_2d par_env_t%par_env_t_bcast_dp_2d none~bcast->none~par_env_t_bcast_dp_2d none~par_env_t_bcast_dp_3d par_env_t%par_env_t_bcast_dp_3d none~bcast->none~par_env_t_bcast_dp_3d none~par_env_t_bcast_dp_4d par_env_t%par_env_t_bcast_dp_4d none~bcast->none~par_env_t_bcast_dp_4d none~par_env_t_bcast_dp_scalar par_env_t%par_env_t_bcast_dp_scalar none~bcast->none~par_env_t_bcast_dp_scalar none~par_env_t_bcast_int32_1d par_env_t%par_env_t_bcast_int32_1d none~bcast->none~par_env_t_bcast_int32_1d none~par_env_t_bcast_int32_scalar par_env_t%par_env_t_bcast_int32_scalar none~bcast->none~par_env_t_bcast_int32_scalar none~par_env_t_bcast_int64_1d par_env_t%par_env_t_bcast_int64_1d none~bcast->none~par_env_t_bcast_int64_1d none~par_env_t_bcast_int64_scalar par_env_t%par_env_t_bcast_int64_scalar none~bcast->none~par_env_t_bcast_int64_scalar 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~build_pfon_density->interface~pack_matrix proc~dmatd_blk dmatd_blk proc~dftexcor->proc~dmatd_blk proc~form_rohf_fock->interface~unpack_matrix dsymm dsymm proc~form_rohf_fock->dsymm proc~form_rohf_fock->interface~pack_matrix proc~orthogonal_transform2 orthogonal_transform2 proc~form_rohf_fock->proc~orthogonal_transform2 proc~orthogonal_transform_sym orthogonal_transform_sym proc~form_rohf_fock->proc~orthogonal_transform_sym proc~triangular_to_full triangular_to_full proc~form_rohf_fock->proc~triangular_to_full proc~get_ab_initio_density->interface~show_message proc~orb_to_dens orb_to_dens proc~get_ab_initio_density->proc~orb_to_dens proc~get_ab_initio_orbital->interface~show_message proc~get_ab_initio_orbital->interface~unpack_matrix proc~diag_symm_full diag_symm_full proc~get_ab_initio_orbital->proc~diag_symm_full proc~get_ab_initio_orbital->proc~oqp_dgemm_i64 proc~orthogonal_transform orthogonal_transform proc~get_ab_initio_orbital->proc~orthogonal_transform proc~matrix_invsqrt->interface~show_message proc~diag_symm_packed diag_symm_packed proc~matrix_invsqrt->proc~diag_symm_packed dgemm dgemm proc~mo_reorder->dgemm proc~reordermos reordermos proc~mo_reorder->proc~reordermos proc~mo_to_ao->interface~unpack_matrix proc~mo_to_ao->interface~pack_matrix proc~mo_to_ao->proc~oqp_dgemm_i64 proc~mo_to_ao->proc~oqp_dsymm_i64 proc~print_eigvec_vals_labeled print_eigvec_vals_labeled proc~print_mo_range->proc~print_eigvec_vals_labeled 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 mpi_bcast mpi_bcast none~par_env_t_bcast_byte->mpi_bcast none~par_env_t_bcast_c_bool->mpi_bcast none~par_env_t_bcast_dp_1d->mpi_bcast none~par_env_t_bcast_dp_2d->mpi_bcast none~par_env_t_bcast_dp_3d->mpi_bcast none~par_env_t_bcast_dp_4d->mpi_bcast none~par_env_t_bcast_dp_scalar->mpi_bcast none~par_env_t_bcast_int32_1d->mpi_bcast none~par_env_t_bcast_int32_scalar->mpi_bcast none~par_env_t_bcast_int64_1d->mpi_bcast none~par_env_t_bcast_int64_scalar->mpi_bcast 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~diag_symm_full->interface~show_message dsyev dsyev proc~diag_symm_full->dsyev proc~diag_symm_packed->interface~show_message dspev dspev proc~diag_symm_packed->dspev dspevx dspevx proc~diag_symm_packed->dspevx proc~dmatd_blk->none~init~14 none~clean~14 xc_consumer_ks_t%clean proc~dmatd_blk->none~clean~14 proc~run_xc run_xc proc~dmatd_blk->proc~run_xc 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~oqp_dgemm_i64->interface~show_message proc~oqp_dgemm_i64->dgemm proc~oqp_dsymm_i64->interface~show_message proc~oqp_dsymm_i64->dsymm proc~orb_to_dens->interface~pack_matrix proc~oqp_dsyr2k_i64 oqp_dsyr2k_i64 proc~orb_to_dens->proc~oqp_dsyr2k_i64 proc~orthogonal_transform->interface~show_message proc~orthogonal_transform->proc~oqp_dgemm_i64 proc~orthogonal_transform2->interface~show_message proc~orthogonal_transform2->proc~oqp_dgemm_i64 proc~orthogonal_transform_sym->interface~show_message proc~orthogonal_transform_sym->proc~oqp_dgemm_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~print_eigvec_vals_labeled->interface~data_has_tags proc~print_eigvec_vals_labeled->interface~tagarray_get_data none~bf_label basis_set%bf_label proc~print_eigvec_vals_labeled->none~bf_label dswap dswap proc~reordermos->dswap proc~triangular_to_full->interface~show_message 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 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~oqp_dsyr2k_i64->interface~show_message dsyr2k dsyr2k proc~oqp_dsyr2k_i64->dsyr2k 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~pack_f90->interface~show_message proc~pack_f90->proc~oqp_dtrttp_i64 proc~run_xc->parallel_start proc~run_xc->parallel_stop proc~run_xc->update none~allreduce par_env_t%allreduce proc~run_xc->none~allreduce none~compaos xc_engine_t%compAOs proc~run_xc->none~compaos none~compxc xc_engine_t%compXC proc~run_xc->none~compxc none~getslicenonzero dft_grid_t%getSliceNonZero proc~run_xc->none~getslicenonzero none~getstats xc_engine_t%getStats proc~run_xc->none~getstats none~init~18 xc_engine_t%init proc~run_xc->none~init~18 none~pruneaos xc_engine_t%pruneAOs proc~run_xc->none~pruneaos none~resetpointers xc_engine_t%resetPointers proc~run_xc->none~resetpointers postUpdate postUpdate proc~run_xc->postUpdate none~par_env_t_allreduce_byte par_env_t%par_env_t_allreduce_byte none~allreduce->none~par_env_t_allreduce_byte none~par_env_t_allreduce_c_bool par_env_t%par_env_t_allreduce_c_bool none~allreduce->none~par_env_t_allreduce_c_bool none~par_env_t_allreduce_dp_1d par_env_t%par_env_t_allreduce_dp_1d none~allreduce->none~par_env_t_allreduce_dp_1d none~par_env_t_allreduce_dp_2d par_env_t%par_env_t_allreduce_dp_2d none~allreduce->none~par_env_t_allreduce_dp_2d none~par_env_t_allreduce_dp_3d par_env_t%par_env_t_allreduce_dp_3d none~allreduce->none~par_env_t_allreduce_dp_3d none~par_env_t_allreduce_dp_4d par_env_t%par_env_t_allreduce_dp_4d none~allreduce->none~par_env_t_allreduce_dp_4d none~par_env_t_allreduce_dp_scalar par_env_t%par_env_t_allreduce_dp_scalar none~allreduce->none~par_env_t_allreduce_dp_scalar none~par_env_t_allreduce_int32_1d par_env_t%par_env_t_allreduce_int32_1d none~allreduce->none~par_env_t_allreduce_int32_1d none~par_env_t_allreduce_int32_scalar par_env_t%par_env_t_allreduce_int32_scalar none~allreduce->none~par_env_t_allreduce_int32_scalar none~par_env_t_allreduce_int64_1d par_env_t%par_env_t_allreduce_int64_1d none~allreduce->none~par_env_t_allreduce_int64_1d none~par_env_t_allreduce_int64_scalar par_env_t%par_env_t_allreduce_int64_scalar none~allreduce->none~par_env_t_allreduce_int64_scalar none~aoval basis_set%aoval none~compaos->none~aoval none~compmos xc_engine_t%compMOs none~compxc->none~compmos none~comprhoall xc_engine_t%compRhoAll none~compxc->none~comprhoall proc~compute~2 xc_libxc_t%compute none~compxc->proc~compute~2 none~getbyid list_grid_3d_t%getByID none~getslicenonzero->none~getbyid none~getenergy xc_lib_t%getEnergy none~getstats->none~getenergy proc~init xc_libxc_t%init none~init~18->proc~init none~resetprunedpointers xc_engine_t%resetPrunedPointers none~pruneaos->none~resetprunedpointers none~resetorbpointers~2 xc_engine_t%resetOrbPointers none~resetpointers->none~resetorbpointers~2 none~resetxcpointers xc_engine_t%resetXCPointers none~resetpointers->none~resetxcpointers proc~setpts xc_libxc_t%setPts none~resetpointers->proc~setpts none~compaov basis_set%compAOv none~aoval->none~compaov none~compaovg basis_set%compAOvg none~aoval->none~compaovg none~compaovgg basis_set%compAOvgg none~aoval->none~compaovgg proc~mo_tran_gemm_ mo_tran_gemm_ none~compmos->proc~mo_tran_gemm_ proc~mo_tran_symm_ mo_tran_symm_ none~compmos->proc~mo_tran_symm_ mpi_allreduce mpi_allreduce none~par_env_t_allreduce_byte->mpi_allreduce none~par_env_t_allreduce_c_bool->mpi_allreduce none~par_env_t_allreduce_dp_1d->mpi_allreduce none~par_env_t_allreduce_dp_2d->mpi_allreduce none~par_env_t_allreduce_dp_3d->mpi_allreduce none~par_env_t_allreduce_dp_4d->mpi_allreduce none~par_env_t_allreduce_dp_scalar->mpi_allreduce none~par_env_t_allreduce_int32_1d->mpi_allreduce none~par_env_t_allreduce_int32_scalar->mpi_allreduce none~par_env_t_allreduce_int64_1d->mpi_allreduce none~par_env_t_allreduce_int64_scalar->mpi_allreduce none~calc_evfxc functional_t%calc_evfxc proc~compute~2->none~calc_evfxc none~calc_evxc functional_t%calc_evxc proc~compute~2->none~calc_evxc none~calc_xc functional_t%calc_xc proc~compute~2->none~calc_xc none~scalexc xc_lib_t%scalexc proc~compute~2->none~scalexc none~clean~13 xc_lib_t%clean proc~init->none~clean~13 none~resetenergy xc_lib_t%resetEnergy proc~init->none~resetenergy none~calc_evfxc->interface~show_message xc_f03_func_info_get_family xc_f03_func_info_get_family none~calc_evfxc->xc_f03_func_info_get_family xc_f03_func_info_get_kind xc_f03_func_info_get_kind none~calc_evfxc->xc_f03_func_info_get_kind xc_f03_func_info_get_name xc_f03_func_info_get_name none~calc_evfxc->xc_f03_func_info_get_name xc_f03_gga_exc_vxc_fxc xc_f03_gga_exc_vxc_fxc none~calc_evfxc->xc_f03_gga_exc_vxc_fxc xc_f03_lda_exc_vxc_fxc xc_f03_lda_exc_vxc_fxc none~calc_evfxc->xc_f03_lda_exc_vxc_fxc xc_f03_mgga_exc_vxc_fxc xc_f03_mgga_exc_vxc_fxc none~calc_evfxc->xc_f03_mgga_exc_vxc_fxc none~calc_evxc->interface~show_message none~calc_evxc->xc_f03_func_info_get_family none~calc_evxc->xc_f03_func_info_get_kind none~calc_evxc->xc_f03_func_info_get_name xc_f03_gga_exc_vxc xc_f03_gga_exc_vxc none~calc_evxc->xc_f03_gga_exc_vxc xc_f03_lda_exc_vxc xc_f03_lda_exc_vxc none~calc_evxc->xc_f03_lda_exc_vxc xc_f03_mgga_exc_vxc xc_f03_mgga_exc_vxc none~calc_evxc->xc_f03_mgga_exc_vxc none~calc_xc->interface~show_message none~calc_xc->xc_f03_func_info_get_family none~calc_xc->xc_f03_func_info_get_kind none~calc_xc->xc_f03_func_info_get_name xc_f03_gga_exc_vxc_fxc_kxc xc_f03_gga_exc_vxc_fxc_kxc none~calc_xc->xc_f03_gga_exc_vxc_fxc_kxc xc_f03_lda_exc_vxc_fxc_kxc xc_f03_lda_exc_vxc_fxc_kxc none~calc_xc->xc_f03_lda_exc_vxc_fxc_kxc xc_f03_mgga_exc_vxc_fxc_kxc xc_f03_mgga_exc_vxc_fxc_kxc none~calc_xc->xc_f03_mgga_exc_vxc_fxc_kxc proc~mo_tran_gemm_->proc~oqp_dgemm_i64 proc~mo_tran_symm_->proc~oqp_dsymm_i64

Called by

proc~~scf_driver~~CalledByGraph proc~scf_driver scf_driver proc~hf_energy hf_energy proc~hf_energy->proc~scf_driver

Source Code

  subroutine scf_driver(basis, infos, molGrid)

  ! Main drirver for HF and DFT with RHF, ROHF and UHF
  ! DIIS is the main algorithm for SCF

     USE precision, only: dp
     use oqp_tagarray_driver
     use types, only: information
     use int2_compute, only: int2_compute_t, int2_fock_data_t, int2_rhf_data_t, int2_urohf_data_t
     use dft, only: dftexcor
     use mod_dft_molgrid, only: dft_grid_t
     use messages,  only: show_message, WITH_ABORT
     use guess, only: get_ab_initio_density, get_ab_initio_orbital
     use util, only: measure_time, e_charge_repulsion
     use printing, only: print_mo_range
     use mathlib, only: traceprod_sym_packed, matrix_invsqrt
     use mathlib, only: unpack_matrix
     use io_constants, only: IW
     use constants, only: kB_HaK 
     use basis_tools, only: basis_set
     use scf_converger, only: scf_conv_result, scf_conv, &
             conv_cdiis, conv_ediis

     implicit none

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

     type(basis_set), intent(in) :: basis
     type(information), target, intent(inout) :: infos
     type(dft_grid_t), intent(in) :: molGrid
     integer :: i, ii, iter, nschwz, nbf, nbf_tri, nbf2, ok, maxit
     real(kind=dp) :: ehf, ehf1, nenergy, etot, diffe, e_old, psinrm, &
                 scalefactor,vne, vnn, vee, vtot, virial, tkin
     real(kind=dp), allocatable :: tempvec(:), lwrk(:), lwrk2(:)
     real(kind=dp), allocatable, target :: smat_full(:,:), pdmat(:,:), pfock(:,:), rohf_bak(:)
     real(kind=dp), allocatable, target :: dold(:,:), fold(:,:)
     integer :: nfocks, diis_reset
!    For DIIS
     real(kind=dp) :: diis_error
     real(kind=dp), parameter :: ethr_cdiis_big = 2.0_dp
     real(kind=dp), parameter :: ethr_ediis = 1.0_dp
     integer :: diis_nfocks, diis_stat, maxdiis
     character(len=6), dimension(5) :: diis_name
!    Vshift
     real(kind=dp) :: vshift, H_U_gap, H_U_gap_crit
     logical :: vshift_last_iter

     type(scf_conv) :: conv
     class(scf_conv_result), allocatable :: conv_res
     character(16) :: scf_name = ""
     real(kind=dp) :: eexc, totele, totkin
     logical :: is_dft, diis_reset_condition
     integer :: scf_type
!    MOM
     logical :: do_mom, step_0_mom, do_mom_flag
     real(kind=dp), allocatable, dimension(:,:) :: mo_a_for_mom, mo_b_for_mom, work

     integer :: nelec, nelec_a, nelec_b
     integer, parameter :: scf_rhf  = 1, &
                           scf_uhf  = 2, &
                           scf_rohf = 3
     real(kind=dp), allocatable :: pfxc(:,:), qmat(:,:)

     type(int2_compute_t) :: int2_driver
     class(int2_fock_data_t), allocatable :: int2_data
!    pFON 
     logical :: do_pfon
     real(kind=dp) :: beta_pfon, start_temp, end_temp, temp_pfon
     real(kind=dp) :: electron_sum_a, electron_sum_b, pfon_start_temp 
     real(kind=dp), allocatable :: occ_a(:), occ_b(:)
     real(kind=dp) :: sum_occ_alpha, sum_occ_beta, cooling_rate
     real(kind=dp) :: pfon_cooling_rate, pfon_nsmear
     real(kind=dp) :: last_cooled_temp = 0.0_dp 
     integer :: nsmear 
  ! tagarray
     real(kind=dp), contiguous, pointer :: &
       dmat_a(:), dmat_b(:), fock_a(:), fock_b(:), hcore(:), mo_b(:,:), &
       smat(:), tmat(:), mo_a(:,:), &
       mo_energy_b(:), mo_energy_a(:), mo_energy_a_for_mom(:)
     character(len=*), parameter :: tags_general(3) = (/ character(len=80) :: &
       OQP_SM, OQP_TM, OQP_Hcore /)
     character(len=*), parameter :: tags_alpha(4) = (/ character(len=80) :: &
       OQP_FOCK_A, OQP_DM_A, OQP_E_MO_A, OQP_VEC_MO_A /)
     character(len=*), parameter :: tags_beta(4) = (/ character(len=80) :: &
       OQP_FOCK_B, OQP_DM_B, OQP_E_MO_B, OQP_VEC_MO_B /)

  !  Default values

  !  MOM settings
  !  Current MOM option works for both RHF and ROHF
     do_mom = infos%control%mom
  !  Vshift settings
  !  Current VSHIFT only works for ROHF
     vshift = infos%control%vshift
     vshift_last_iter=.false.
     H_U_gap_crit=0.02_dp

  !  pFON settings
     do_pfon = .false. 
     do_pfon = infos%control%pfon 
     start_temp = infos%control%pfon_start_temp
     cooling_rate = infos%control%pfon_cooling_rate
     if (start_temp <= 0.0_dp) then 
         start_temp = 2000.0_dp 
     end if 
     temp_pfon  = start_temp
     if (temp_pfon < 1.0_dp) temp_pfon = 1.0_dp

     beta_pfon = 1.0_dp / (kB_HaK * temp_pfon)

  !  DIIS options
  !  none IS NOT recommended!
  !  c-DIIS: Default commutator DIIS
  !  e-DIIS: Energy-Weighted DIIS
  !  a-DIIS: Augmented DIIS
  !  v-DIIS: Vshift with DIIS, which will be eventually become c-DIIS
  !  MOM: Maximum Overlap Method for better convergency

     diis_error = 2.0_dp
     diis_name = [character(len=6) :: "none", "c-DIIS","e-DIIS","a-DIIS","v-DIIS"]
  !  Read calculation metadata from infos
     is_dft = infos%control%hamilton >= 20
     select case (infos%control%scftype)
     case (1)
       scf_type = scf_rhf
     case (2)
       scf_type = scf_uhf
     case (3)
       scf_type = scf_rohf
     end select
  !  Total number of electrons
     nelec = infos%mol_prop%nelec
     nelec_a = infos%mol_prop%nelec_a
     nelec_b = infos%mol_prop%nelec_b
  !  Matrix size
     nbf = basis%nbf
     nbf_tri = nbf*(nbf+1)/2
     nbf2 = nbf*nbf
     maxit = infos%control%maxit
     maxdiis = infos%control%maxdiis
     diis_reset = infos%control%diis_reset_mod
  !  DFT
     if (.not.is_dft) then
         scalefactor = 1.0_dp
     else
         scalefactor = infos%dft%HFscale
     end if
  !  SCF Type
     select case (scf_type)
     case (scf_rhf)
       scf_name = "RHF"
       nfocks = 1
       diis_nfocks = 1
     case (scf_uhf)
       scf_name = "UHF"
       nfocks = 2
       diis_nfocks = 2
     case (scf_rohf)
       scf_name = "ROHF"
       nfocks = 2
       diis_nfocks = 1
     end select

  !  Now we are allocating dynamic memories of tag arrays
  !  Tag arrays
     call data_has_tags(infos%dat, tags_general, module_name, subroutine_name, WITH_ABORT)
     call tagarray_get_data(infos%dat, OQP_Hcore, hcore)
     call tagarray_get_data(infos%dat, OQP_SM, smat)
     call tagarray_get_data(infos%dat, OQP_TM, tmat)

     call data_has_tags(infos%dat, tags_alpha, module_name, subroutine_name, WITH_ABORT)
     call tagarray_get_data(infos%dat, OQP_FOCK_A, fock_a)
     call tagarray_get_data(infos%dat, OQP_DM_A, dmat_a)
     call tagarray_get_data(infos%dat, OQP_E_MO_A, mo_energy_a)
     call tagarray_get_data(infos%dat, OQP_VEC_MO_A, mo_a)

     if (nfocks > 1) then
       call data_has_tags(infos%dat, tags_beta, module_name, subroutine_name, WITH_ABORT)
       call tagarray_get_data(infos%dat, OQP_FOCK_B, fock_b)
       call tagarray_get_data(infos%dat, OQP_DM_B, dmat_b)
       call tagarray_get_data(infos%dat, OQP_E_MO_B, mo_energy_b)
       call tagarray_get_data(infos%dat, OQP_VEC_MO_B, mo_b)
     endif

     allocate(smat_full(nbf,nbf), pdmat(nbf_tri,nfocks), pfock(nbf_tri,nfocks), &
              qmat(nbf,nbf), &
              stat=ok, &
              source=0.0_dp)
  !  Alloc work matrices for MOM
     if (do_mom) then
        step_0_mom = .true.
        do_mom_flag=.false.
        allocate(mo_a_for_mom(nbf,nbf), &
                 mo_b_for_mom(nbf,nbf), &
                 mo_energy_a_for_mom(nbf), &
                 work(nbf,nbf), &
                 source=0.0_dp)
     end if
     if(ok/=0) call show_message('Cannot allocate memory for SCF',WITH_ABORT)
  !  Incremental Fock building for Direct SCF
     if (infos%control%scf_incremental /= 0) then
       allocate( dold(nbf_tri,nfocks), fold(nbf_tri,nfocks), &
                stat=ok, &
                source=0.0_dp)
       if(ok/=0) call show_message('Cannot allocate memory for SCF',WITH_ABORT)
     end if
  !  DFT
     if (is_dft) then
       allocate(pfxc(nbf_tri,nfocks), &
                stat=ok, &
                source=0.0_dp)
       if(ok/=0) call show_message('Cannot allocate memory for temporary vectors',WITH_ABORT)
     end if
  !  ROHF temporary arrays
     if (scf_type == scf_rohf .or. is_dft) then
       allocate(tempvec(nbf2), &
                stat=ok, &
                source=0.0_dp)
       if(ok/=0) call show_message('Cannot allocate memory for temporary vectors',WITH_ABORT)
     end if
     if (scf_type == scf_rohf) then
       allocate(lwrk(nbf), lwrk2(nbf2), rohf_bak(nbf_tri), &
                stat=ok, &
                source=0.0_dp)
       if(ok/=0) call show_message('Cannot allocate memory for ROHF temporaries',WITH_ABORT)
     end if

  !  Allocating dynamic memories done

     call measure_time(print_total=1, log_unit=iw)

     call matrix_invsqrt(smat, qmat, nbf)

  !  First Compute Nuclear-Nuclear energy
     nenergy = e_charge_repulsion(infos%atoms%xyz, infos%atoms%zn - infos%basis%ecp_zn_num)

  !  During guess, the Hcore, Q nd Overlap matrices were formed.
  !  Using these, the initial orbitals (VEC) and density (Dmat) were subsequently computed.
  !  Now we are going to calculate ERI(electron repulsion integrals) to form a new FOCK
  !  matrix.

  !  Vshift settings
  !   vshift = infos%control%vshift
  !   vshift_last_iter=.false.
  !   H_U_gap_crit=0.02_dp

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

  !  The main Do loop of SCF
     select case (scf_type)
     case (scf_rhf)
         pdmat(:,1) = dmat_a
         allocate(int2_rhf_data_t :: int2_data)
         int2_data = int2_rhf_data_t(nfocks=1, d=pdmat, scale_exchange=scalefactor)
     case (scf_uhf, scf_rohf)
         pdmat(:,1) = dmat_a
         pdmat(:,2) = dmat_b
         allocate(int2_urohf_data_t :: int2_data)
         int2_data = int2_urohf_data_t(nfocks=2, d=pdmat, scale_exchange=scalefactor)
     end select
     call unpack_matrix(smat,smat_full,nbf,'U')

!    Variable DIIS options:
!      a) If we use VSHIFT option, the combination of e-DIIS and c-DIIS is used,
!      since e-DIIS can better reduce total energy than c-DIIS.
!      Once sufficiently converged, c-DIIS is used to finalize the SCF.
!      b) If VSHIFT is not set, c-DIIS is default.
!      c) If vdiis (diis_type=5) is chosen, the VSHIFT is initally turned on.
!      d) if MOM=.true., MOM turns on if DIIS error < mom_switch
     if (infos%control%diis_type == 5) then
        call conv%init(ldim=nbf, maxvec=maxdiis, &
             subconvergers=[conv_cdiis, conv_ediis, conv_cdiis], &
             thresholds=[ethr_cdiis_big, ethr_ediis, infos%control%vdiis_cdiis_switch], &
             overlap=smat_full, overlap_sqrt=qmat, &
             num_focks=diis_nfocks, verbose=1)
        if (infos%control%vshift == 0.0_dp) then
           infos%control%vshift=0.1_dp
           vshift=0.1_dp
           call show_message("")
           call show_message('Setting VSHIFT=0.1, since VDIIS is chosen without VSHIFT value.')
        endif
     elseif (infos%control%vshift /= 0.0_dp) then
        call conv%init(ldim=nbf, maxvec=maxdiis, &
             subconvergers=[conv_cdiis, conv_ediis, conv_cdiis], &
             thresholds=[ethr_cdiis_big, ethr_ediis, infos%control%vshift_cdiis_switch], &
             overlap=smat_full, overlap_sqrt=qmat, &
             num_focks=diis_nfocks, verbose=1)
     else
!    Normally, c-DIIS works best. But one can choose others (e-DIIS and a-DIIS).
        call conv%init(ldim=nbf, maxvec=maxdiis, &
             subconvergers=[infos%control%diis_type], &
             thresholds=[infos%control%diis_method_threshold], &
             overlap=smat_full, overlap_sqrt=qmat, &
             num_focks=diis_nfocks, verbose=1)
     endif

     eexc = 0.0_dp
     e_old = 0.0_dp

  !  SCF Options
     write(iw,'(/5X,"SCF options"/ &
                &5X,18("-")/ &
                &5X,"SCF type = ",A,5x,"MaxIT = ",I5/, &
                &5X,"MaxDIIS = ",I5,17x,"Conv = ",F14.10/, &
                &5X,"DIIS Type = ",A/, &
                &5X,"vDIIS_cDIIS_Switch = ",F8.5,3x,"vDIIS_vshift_Switch = ",F8.5/, &
                &5X,"DIIS Reset Mod = ",I5,10x,"DIIS Reset Conv = ",F12.8/, &
                &5X,"VShift = ",F8.5,15X,"VShift_cDIIS_Switch = ",F8.5)') &
                & scf_name, infos%control%maxit, &
                & infos%control%maxdiis, infos%control%conv, &
                & diis_name(infos%control%diis_type), &
                & infos%control%vdiis_cdiis_switch, infos%control%vdiis_vshift_switch, &
                & infos%control%diis_reset_mod, infos%control%diis_reset_conv, &
                & infos%control%vshift, infos%control%vshift_cdiis_switch
     write(iw,'(5X,"MOM = ",L5,21X,"MOM_Switch = ",F8.5)') &
                & infos%control%mom, infos%control%mom_switch 
     write(iw,'(5X,"pFON = ",L5,20X,"pFON Start Temp. = ",F9.2,/, &
                5X, "pFON Cooling Rate = ", F9.2,2X," pFON Num. Smearing = ",F8.5)') &
                infos%control%pfon, infos%control%pfon_start_temp, &
                infos%control%pfon_cooling_rate, infos%control%pfon_nsmear
     !  Initial message
     if (infos%control%pfon) then
           write(IW,fmt="&
                 &(/3x,'Direct SCF iterations begin.'/, &
                 &  3x,113('='),/ &
                 &  4x,'Iter',9x,'Energy',12x,'Delta E',9x,'Int Skip',5x,'DIIS Error',5x,'Shift',5x,'Method',5x,'pFON'/ &
                 &  3x,113('='))")
     else
          write(IW,fmt="&
                 &(/3x,'Direct SCF iterations begin.'/, &
                 &  3x,93('='),/ &
                 &  4x,'Iter',9x,'Energy',12x,'Delta E',9x,'Int Skip',5x,'DIIS Error',5x,'Shift',5x,'Method'/ &
                 &  3x,93('='))")
     endif
     call flush(iw)
     do iter = 1, maxit

  !     The main SCF iteration loop




!     pFON Cooling
      if (cooling_rate <= 0.0_dp) then
          cooling_rate = 50_dp 
      end if 

      if (do_pfon) then
          if (iter == maxit) then 
              temp_pfon = 0.0_dp 
          else if (abs(diis_error) < 10.0_dp * infos%control%conv) then 
              if (temp_pfon > 1.0_dp) then
                  last_cooled_temp = temp_pfon
              end if
              temp_pfon = 1.0_dp 
          else 
              if (temp_pfon == 1.0_dp .and. last_cooled_temp > 1.0_dp) then
                  temp_pfon = last_cooled_temp
              end if
              temp_pfon = temp_pfon - cooling_rate 
              if (temp_pfon < 1.0_dp) then
                  temp_pfon = 1.0_dp 
              end if
              last_cooled_temp = temp_pfon
          end if 

          if (temp_pfon > 1.0e-12_dp) then 
              beta_pfon = 1.0_dp / (kB_HaK * temp_pfon)
          else 
              beta_pfon = 1.0e20_dp 
          end if
      end if 

      pfock = 0.0_dp 

  !     Compute difference density matrix for incremental Fock build,
  !     which is the main advantage of direct SCF.
  !     It will provide much better ERI screening.

        if (infos%control%scf_incremental /= 0) then
          pdmat = pdmat - dold
        end if

        call int2_driver%run(int2_data, &
                cam=is_dft.and.infos%dft%cam_flag, &
                alpha=infos%dft%cam_alpha, &
                beta=infos%dft%cam_beta,&
                mu=infos%dft%cam_mu)
        nschwz = int2_driver%skipped

  !     Recover full Fock and density from difference matrices
        if (infos%control%scf_incremental /= 0) then
          pdmat = pdmat + dold
          int2_data%f(:,:,1) = int2_data%f(:,:,1) + fold
          fold = int2_data%f(:,:,1)
          dold = pdmat
        end if
  !     Scaling
        pfock(:,:) = 0.5_dp * int2_data%f(:,:,1)
        ii=0
        do i = 1, nbf
           ii = ii + i
           pfock(ii,1:nfocks) = 2.0_dp*pfock(ii,1:nfocks)
        end do

  !     Adding the skeleton H core to Fock for getting the new orbitals
        do i = 1, nfocks
          pfock(:,i)  = pfock(:,i) + hcore
        end do

  !     After this, we compute Energy.
        ehf = 0.0_dp
        ehf1 = 0.0_dp
        do i = 1, nfocks
          ehf1 = ehf1 + traceprod_sym_packed(pdmat(:,i),hcore,nbf)
          ehf = ehf + traceprod_sym_packed(pdmat(:,i),pfock(:,i),nbf)
        end do
        ehf = 0.5_dp*(ehf+ehf1)
        etot = ehf + nenergy
  !     DFT contribution
        if (is_dft) then
          if (scf_type == scf_rhf) then
            call dftexcor(basis,molgrid,1,pfxc,pfxc,mo_a,mo_a,nbf,nbf_tri,eexc,totele,totkin,infos)
          else if (scf_type == scf_uhf) then
             call dftexcor(basis,molgrid,2,pfxc(:,1),pfxc(:,2),mo_a,mo_b,nbf,nbf_tri,eexc,totele,totkin,infos)
          else if (scf_type == scf_rohf) then
!            ROHF does not have MO_B. So we copy MO_A to MO_B.
             mo_b = mo_a
             call dftexcor(basis,molgrid,2,pfxc(:,1),pfxc(:,2),mo_a,mo_b,nbf,nbf_tri,eexc,totele,totkin,infos)
          end if
          pfock = pfock + pfxc
          etot = etot + eexc
        end if

  !     Forming ROHF Fock by combing Alpha and Beta Focks.
        if (scf_type == scf_rohf) then
           rohf_bak = pfock(:,1)
           if (vshift_last_iter.eqv..true.) vshift = 0.0_dp
           call form_rohf_fock(pfock(:,1),pfock(:,2),tempvec, &
             mo_a,smat,lwrk2,nelec_a,nelec_b,nbf,vshift)
           pdmat(:,1) = pdmat(:,1) + pdmat(:,2)
        end if

  !     SCF Converger to get refined Fock Matrix
        call conv%add_data(f=pfock(:,1:diis_nfocks), &
                dens=pdmat(:,1:diis_nfocks), e=Etot)

  !     Run DIIS calculation, get DIIS error
        call conv%run(conv_res)
        diis_error = conv_res%error
  !     Checking the convergency
        diffe = etot-e_old
        if (infos%control%pfon) then
           write(IW,fmt="(4x,i4.1,2x,f17.10,1x,f17.10,1x,i13,1x,f14.8,5x,f5.3,5x,a,5x,a,f9.2)") &
                 iter, etot, diffe, nschwz, diis_error, vshift, &
                 trim(conv_res%active_converger_name), "Temp:", temp_pfon
           write(IW,fmt="(100x,a,f9.2)") "Beta:", beta_pfon
        else
           write(iw,'(4x,i4.1,2x,f17.10,1x,f17.10,1x,i13,1x,f14.8,5x,f5.3,5x,a)') &
                 iter, etot, diffe, nschwz, diis_error, vshift, &
                 trim(conv_res%active_converger_name)
        endif
        call flush(iw)
  !     VDIIS option
        if ((infos%control%diis_type.eq.5) &
           .and.(diis_error < infos%control%vdiis_vshift_switch)) then
           vshift=0.0_dp
        elseif ((infos%control%diis_type.eq.5) &
           .and.(diis_error.ge.infos%control%vdiis_vshift_switch)) then
           vshift=infos%control%vshift
        endif

        e_old = etot

  !     Exit if convergence criteria achieved
        if ((abs(diis_error)<infos%control%conv).and.(vshift==0.0_dp)) then
           exit
        elseif ((abs(diis_error)<infos%control%conv).and.(vshift/=0.0_dp)) then
           write(iw,"(3x,64('-')/10x,'Performing a last SCF with zero VSHIFT.')")
           vshift_last_iter=.true.
        elseif (vshift_last_iter.eqv..true.) then
  !        Only for ROHF case
           call get_ab_initio_orbital(pfock(:,1),mo_a,mo_energy_a,qmat)
           exit
        endif

  !     Reset DIIS for difficult case

  !     VSHIFT=0 and slow cases
        diis_reset_condition=(((iter/diis_reset).ge.1) &
           .and.(modulo(iter,diis_reset).eq.0) &
           .and.(diis_error.gt.infos%control%diis_reset_conv) &
           .and.(infos%control%vshift==0.0_dp))
        if (diis_reset_condition) then
!          Resetting DIIS
           write(iw,"(3x,64('-')/10x,'Resetting DIIS.')")
           call conv_res%get_fock(matrix=pfock(:,1:diis_nfocks), istat=diis_stat)
           call conv%init(ldim=nbf, maxvec=maxdiis, &
             subconvergers=[conv_cdiis], &
             thresholds=[ethr_cdiis_big], &
             overlap=smat_full, overlap_sqrt=qmat, &
             num_focks=diis_nfocks, verbose=1)
!          After resetting DIIS, we need to skip SD
           call conv%add_data(f=pfock(:,1:diis_nfocks), &
                dens=pdmat(:,1:diis_nfocks), e=Etot)
           call conv%run(conv_res)
        else
  !        Form the interpolated the Fock/density matrix
           call conv_res%get_fock(matrix=pfock(:,1:diis_nfocks), istat=diis_stat)
        endif

  !     Calculate new orbitals and density.
  !
        if (int2_driver%pe%rank == 0) then
           call get_ab_initio_orbital(pfock(:,1),mo_a,mo_energy_a,qmat)

           if (scf_type == scf_uhf .and. nelec_b /= 0) then
     !        Only UHF has beta orbitals.
               call get_ab_initio_orbital(pfock(:,2),mo_b,mo_energy_b,qmat)
           end if
        end if
        if (scf_type == scf_uhf .and. nelec_b /= 0) then
            call int2_driver%pe%bcast(mo_b, size(mo_b))
            call int2_driver%pe%bcast(mo_energy_b, size(mo_energy_b))
        end if
        call int2_driver%pe%bcast(mo_a, size(mo_a))
        call int2_driver%pe%bcast(mo_energy_a, size(mo_energy_a))

        ! pFON section
        do_pfon = infos%control%pfon
        if (do_pfon) then
            if (.not. allocated(occ_a)) allocate(occ_a(nbf))
            if (.not. allocated(occ_b)) allocate(occ_b(nbf))
            nsmear = infos%control%pfon_nsmear

            select case (scf_type)
            case (scf_rhf)
                call pfon_occupations(mo_energy_a, nbf, nelec, occ_a, beta_pfon, scf_type, nsmear)

            case (scf_uhf)
                call pfon_occupations(mo_energy_a, nbf, nelec, occ_a, beta_pfon, scf_type, nsmear, &
                                    is_beta=.false., nelec_a=nelec_a, nelec_b=nelec_b)
                if (nelec_b > 0) then
                    call pfon_occupations(mo_energy_b, nbf, nelec, occ_b, beta_pfon, scf_type, nsmear, &
                                        is_beta=.true., nelec_a=nelec_a, nelec_b=nelec_b)
                end if

            case (scf_rohf)
                call pfon_occupations(mo_energy_a, nbf, nelec, occ_a, beta_pfon, scf_type, nsmear, &
                                    is_beta=.false., nelec_a=nelec_a, nelec_b=nelec_b)
                if (nelec_b > 0) then
                    call pfon_occupations(mo_energy_a, nbf, nelec, occ_b, beta_pfon, scf_type, nsmear, &
                                        is_beta=.true., nelec_a=nelec_a, nelec_b=nelec_b)
                end if
            end select

            ! Alpha occupations
            sum_occ_alpha = sum(occ_a(1:nbf))

            ! Beta occupations
            sum_occ_beta = 0.0_dp
            if (scf_type == scf_uhf .or. scf_type == scf_rohf) then
                sum_occ_beta = sum(occ_b(1:nbf))
            end if

        end if


  !     MOM option works for RHF and ROHF
        if (do_mom .and. diis_error.lt.infos%control%mom_switch) do_mom_flag=.true.
        if (do_mom .and. do_mom_flag .and. .not. step_0_mom) then
           call mo_reorder(infos, mo_a_for_mom, mo_energy_a_for_mom, &
                           mo_a, mo_energy_a, smat_full)
        end if
        if (do_mom) then
           mo_a_for_mom = mo_a
           mo_energy_a_for_mom = mo_energy_a
        end if
        step_0_mom = .false.

  !     New density matrix in AO basis using MO.
        if (int2_driver%pe%rank == 0) then
            if (.not. do_pfon) then 
                call get_ab_initio_density(pdmat(:,1),mo_a,pdmat(:,2),mo_b,infos,basis)
            else 
                call build_pfon_density(pdmat, mo_a, mo_b, occ_a, occ_b, scf_type, nbf, nelec_a, nelec_b)
            end if
        end if 
        call int2_driver%pe%bcast(pdmat, size(pdmat))
 !     Checking the HOMO-LUMO gaps for predicting SCF convergency
        if ((iter > 10).and.(vshift==0.0_dp)) then
           select case (scf_type)
           case (scf_rhf)
              H_U_gap=mo_energy_a(nelec/2)-mo_energy_a(nelec/2-1)
           case (scf_uhf)
              H_U_gap=mo_energy_a(nelec/2)-mo_energy_a(nelec/2-1)
           case (scf_rohf)
              H_U_gap=mo_energy_a(nelec_a+1)-mo_energy_a(nelec_a)
           end select
        endif

  !  End of Iteration
     end do

     if (iter>maxit) then
       write(iw,"(3x,64('-')/10x,'SCF is not converged ....')")
       infos%mol_energy%SCF_converged=.false.
     else
       write(iw,"(3x,64('-')/10x,'SCF convergence achieved ....')")
       infos%mol_energy%SCF_converged=.true.
     end if

     write(iw,"(/' Final ',A,' energy is',F20.10,' after',I4,' iterations'/)") trim(scf_name), etot, iter

     if (is_dft) then
       write(iw,*)
       write(iw,"(' DFT: XC energy              = ',F20.10)") eexc
       write(iw,"(' DFT: total electron density = ',F20.10)") totele
       write(iw,"(' DFT: number of electrons    = ',I9,/)") nelec
     end if
  !
     if (scf_type == scf_uhf .and. nelec_b /= 0) then
         call int2_driver%pe%bcast(mo_b, size(mo_b))
         call int2_driver%pe%bcast(mo_energy_b, size(mo_energy_b))
     end if

     call int2_driver%pe%bcast(mo_a, size(mo_a))
     call int2_driver%pe%bcast(pdmat, size(pdmat))
     call int2_driver%pe%bcast(mo_energy_a, size(mo_energy_a))

     select case (scf_type)
     case (scf_rhf)
       fock_a = pfock(:,1)
       dmat_a = pdmat(:,1)
     case (scf_uhf)
       fock_a = pfock(:,1)
       fock_b = pfock(:,2)
       dmat_a = pdmat(:,1)
       dmat_b = pdmat(:,2)
     case (scf_rohf)
       fock_a = rohf_bak
       call mo_to_ao(fock_b, pfock(:,2), smat, mo_a, nbf, nbf)
       dmat_a = pdmat(:,1) - pdmat(:,2)
       dmat_b = pdmat(:,2)
       mo_b = mo_a
       mo_energy_b = mo_energy_a
     end select

  !  Print out the molecular orbitals
     call print_mo_range(basis, infos, mostart=1, moend=nbf)

  !  Print out the results
     psinrm = 0.0_dp
     tkin = 0.0_dp
     do i = 1, diis_nfocks
       psinrm = psinrm + traceprod_sym_packed(pdmat(:,i),smat,nbf)/nelec
       tkin = tkin + traceprod_sym_packed(pdmat(:,i),tmat,nbf)
     end do

  !  Writing out final results
     vne = ehf1 - tkin
     vee = etot - ehf1 - nenergy
     vnn = nenergy
     vtot = vne + vnn + vee
     virial = - vtot/tkin
     call print_scf_energy(psinrm, ehf1, nenergy, etot, vee, vne, vnn, vtot, tkin, virial)

  !  Save results to infos.
     infos%mol_energy%energy = etot
     infos%mol_energy%psinrm = psinrm
     infos%mol_energy%ehf1 = ehf1
     infos%mol_energy%vee = vee
     infos%mol_energy%nenergy = nenergy
     infos%mol_energy%vne = vne
     infos%mol_energy%vnn = vnn
     infos%mol_energy%vtot = vtot
     infos%mol_energy%tkin = tkin
     infos%mol_energy%virial = virial
     infos%mol_energy%energy = etot

  !  Clean up
     call int2_driver%clean()

     call measure_time(print_total=1, log_unit=iw)

  end subroutine scf_driver