int1e Subroutine

public subroutine int1e(infos)

Uses

  • proc~~int1e~~UsesGraph proc~int1e int1e module~basis_tools basis_tools proc~int1e->module~basis_tools module~constants constants proc~int1e->module~constants module~int1 int1 proc~int1e->module~int1 module~io_constants io_constants proc~int1e->module~io_constants module~messages messages proc~int1e->module~messages module~oqp_tagarray_driver oqp_tagarray_driver proc~int1e->module~oqp_tagarray_driver module~physical_constants physical_constants proc~int1e->module~physical_constants module~precision precision proc~int1e->module~precision module~printing printing proc~int1e->module~printing module~strings strings proc~int1e->module~strings module~types types proc~int1e->module~types 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~int1->module~basis_tools module~int1->module~messages module~int1->iso_fortran_env module~mod_1e_primitives mod_1e_primitives module~int1->module~mod_1e_primitives module~mod_shell_tools mod_shell_tools module~int1->module~mod_shell_tools 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 iso_c_binding iso_c_binding module~oqp_tagarray_driver->iso_c_binding tagarray tagarray module~oqp_tagarray_driver->tagarray module~physical_constants->iso_fortran_env module~precision->iso_fortran_env module~printing->module~precision module~strings->iso_c_binding 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~atomic_structure_m->iso_c_binding module~functionals->module~precision module~functionals->iso_c_binding xc_f03_lib_m xc_f03_lib_m module~functionals->xc_f03_lib_m module~mod_1e_primitives->module~constants module~mod_1e_primitives->iso_fortran_env module~mod_1e_primitives->module~mod_shell_tools module~mod_gauss_hermite mod_gauss_hermite module~mod_1e_primitives->module~mod_gauss_hermite module~rys rys module~mod_1e_primitives->module~rys module~xyz_order xyz_order module~mod_1e_primitives->module~xyz_order module~mod_shell_tools->module~basis_tools module~mod_shell_tools->module~precision module~parallel->module~precision module~parallel->iso_c_binding module~parallel->iso_fortran_env mpi mpi module~parallel->mpi module~mod_gauss_hermite->module~precision module~rys->module~constants module~rys->module~precision module~rys_lut rys_lut module~rys->module~rys_lut

@brief Calculate the basic H, S, and T 1e-integrals

Arguments

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

Calls

proc~~int1e~~CallsGraph proc~int1e int1e interface~data_has_tags data_has_tags proc~int1e->interface~data_has_tags interface~tagarray_get_data tagarray_get_data proc~int1e->interface~tagarray_get_data proc~omp_hst omp_hst proc~int1e->proc~omp_hst proc~print_module_info print_module_info proc~int1e->proc~print_module_info proc~print_sym_labeled print_sym_labeled proc~int1e->proc~print_sym_labeled remove_records remove_records proc~int1e->remove_records reserve_data reserve_data proc~int1e->reserve_data proc~omp_hst->proc~print_sym_labeled interface~bas_norm_matrix bas_norm_matrix proc~omp_hst->interface~bas_norm_matrix interface~int1_coul int1_coul proc~omp_hst->interface~int1_coul none~alloc~2 shpair_t%alloc proc~omp_hst->none~alloc~2 none~bcast par_env_t%bcast proc~omp_hst->none~bcast none~fetch_by_id shell_t%fetch_by_id proc~omp_hst->none~fetch_by_id none~init~14 par_env_t%init proc~omp_hst->none~init~14 none~shell_pair shpair_t%shell_pair proc~omp_hst->none~shell_pair proc~add_ecpint add_ecpint proc~omp_hst->proc~add_ecpint proc~comp_kin_ovl_int1_prim comp_kin_ovl_int1_prim proc~omp_hst->proc~comp_kin_ovl_int1_prim proc~comp_lz_int1_prim comp_lz_int1_prim proc~omp_hst->proc~comp_lz_int1_prim proc~update_triang_matrix update_triang_matrix proc~omp_hst->proc~update_triang_matrix none~bf_label basis_set%bf_label proc~print_sym_labeled->none~bf_label 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 mpi_comm_rank mpi_comm_rank none~init~14->mpi_comm_rank mpi_comm_size mpi_comm_size none~init~14->mpi_comm_size interface~compute_integrals compute_integrals proc~add_ecpint->interface~compute_integrals interface~free_integrator free_integrator proc~add_ecpint->interface~free_integrator interface~free_result free_result proc~add_ecpint->interface~free_result interface~init_integrator init_integrator proc~add_ecpint->interface~init_integrator interface~init_integrator_instance init_integrator_instance proc~add_ecpint->interface~init_integrator_instance interface~set_ecp_basis set_ecp_basis proc~add_ecpint->interface~set_ecp_basis none~bf_to_shell basis_set%bf_to_shell proc~add_ecpint->none~bf_to_shell proc~doquadgausshermite doQuadGaussHermite proc~comp_kin_ovl_int1_prim->proc~doquadgausshermite proc~comp_lz_int1_prim->proc~doquadgausshermite 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 abrt abrt proc~doquadgausshermite->abrt

Source Code

  subroutine int1e(infos)

    use types, only: information
    use oqp_tagarray_driver
    use precision, only: dp
    use io_constants, only: iw
    use constants, only: tol_int
    use int1, only: omp_hst
    use basis_tools, only: basis_set
    use printing, only: print_sym_labeled
    use messages, only: show_message, WITH_ABORT
    use strings, only: Cstring, fstring
    use physical_constants, only: BOHR_TO_ANGSTROM
    use printing, only: print_module_info

    implicit none

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

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

    real(kind=dp) :: tol
    integer :: i, nbf, nat, nbf2

    ! tagarray
    real(kind=dp), contiguous, pointer :: &
      hcore(:), tmat(:), smat(:)
    character(len=*), parameter :: tags_general(3) = (/ character(len=80) :: &
      OQP_SM, OQP_TM, OQP_Hcore /)

    logical dbg
    dbg = .false.

!   Files open:
!   LOG: Read and Write: Main output file
    open(unit=iw,  file=infos%log_filename, position="append")

!   Load basis set
    basis => infos%basis
    basis%atoms => infos%atoms
!
    call print_module_info('int1e','Computing H, S and T Matrices')

!   Print out the Cartesian Coordinates.
    write(iw, fmt="(&
              &/21X,32('=')&
              &/21X,a&
              &/21X,32('=')&
              &/8X,'ATOM     ZNUC',11X,'X',14X,'Y',14X,'Z'&
              &/6X,62('-'))") "Cartesian Coordinate in Angstrom"

    do i = 1, size(basis%atoms%zn(:))
       write(iw,'(7x,i4,5x,f4.1,3(x,f15.9))') &
               i, basis%atoms%zn(i), basis%atoms%xyz(1:3,i)*BOHR_TO_ANGSTROM
    end do

!   Allocate H, S and T matrices
    nbf2 = basis%nbf*(basis%nbf+1)/2

    call infos%dat%remove_records(tags_general)

    call infos%dat%reserve_data(OQP_SM, TA_TYPE_REAL64, nbf2, comment=OQP_SM_comment)
    call infos%dat%reserve_data(OQP_TM, TA_TYPE_REAL64, nbf2, comment=OQP_TM_comment)
    call infos%dat%reserve_data(OQP_Hcore, TA_TYPE_REAL64, nbf2, comment=OQP_Hcore_comment)

    call data_has_tags(infos%dat, tags_general, module_name, subroutine_name, WITH_ABORT)
    call tagarray_get_data(infos%dat, OQP_SM, smat)
    call tagarray_get_data(infos%dat, OQP_TM, tmat)
    call tagarray_get_data(infos%dat, OQP_Hcore, Hcore)

!   Create arrays of atomic coordinates and charges for one-electron code
    nbf = basis%nbf
    nat = ubound(infos%atoms%zn,1)

!   Compute conventional H, S, and T integrals
    tol = log(10.0d0)*tol_int
    call omp_hst(basis, infos%atoms%xyz, infos%atoms%zn - infos%basis%ecp_zn_num, hcore, smat, tmat,&
            logtol=tol, comm=infos%mpiinfo%comm, usempi=infos%mpiinfo%usempi)

    if (dbg) then
        write(iw,'(/"BARE NUCLEUS HAMILTONIAN INTEGRALS (H=T+V)")')
        call print_sym_labeled(hcore,nbf, basis)

        write(iw,'(/"OVERLAP MATRIX")')
        call print_sym_labeled(Smat,nbf, basis)

        write(iw,'(/"KINETIC ENERGY INTEGRALS")')
        call print_sym_labeled(tmat,nbf, basis)
    end if

    write(iw,"(/1X,'...... End Of One Electron Integrals ......'/)")

    close(iw)

  end subroutine int1e