oqp_apply_basis Subroutine

public subroutine oqp_apply_basis(infos)

Uses

  • proc~~oqp_apply_basis~~UsesGraph proc~oqp_apply_basis oqp_apply_basis basis_api basis_api proc~oqp_apply_basis->basis_api iso_c_binding iso_c_binding proc~oqp_apply_basis->iso_c_binding module~atomic_structure_m atomic_structure_m proc~oqp_apply_basis->module~atomic_structure_m module~messages messages proc~oqp_apply_basis->module~messages module~oqp_tagarray_driver oqp_tagarray_driver proc~oqp_apply_basis->module~oqp_tagarray_driver module~parallel parallel proc~oqp_apply_basis->module~parallel module~strings strings proc~oqp_apply_basis->module~strings module~types types proc~oqp_apply_basis->module~types module~atomic_structure_m->iso_c_binding comm_IOFILE comm_IOFILE module~messages->comm_IOFILE comm_PAR comm_PAR module~messages->comm_PAR module~io_constants io_constants module~messages->module~io_constants module~precision precision module~messages->module~precision module~oqp_tagarray_driver->iso_c_binding tagarray tagarray module~oqp_tagarray_driver->tagarray module~parallel->iso_c_binding iso_fortran_env iso_fortran_env module~parallel->iso_fortran_env module~parallel->module~precision mpi mpi module~parallel->mpi module~strings->iso_c_binding module~types->iso_c_binding module~types->module~atomic_structure_m module~types->module~parallel module~basis_tools basis_tools module~types->module~basis_tools module~functionals functionals module~types->module~functionals module~types->module~precision module~types->tagarray module~basis_tools->module~atomic_structure_m module~basis_tools->module~parallel module~basis_tools->iso_fortran_env module~basis_tools->module~io_constants module~basis_tools->module~precision module~constants constants module~basis_tools->module~constants module~functionals->iso_c_binding module~functionals->module~precision xc_f03_lib_m xc_f03_lib_m module~functionals->xc_f03_lib_m module~precision->iso_fortran_env module~constants->module~precision

Arguments

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

Calls

proc~~oqp_apply_basis~~CallsGraph proc~oqp_apply_basis oqp_apply_basis interface~data_has_tags data_has_tags proc~oqp_apply_basis->interface~data_has_tags interface~tagarray_get_data tagarray_get_data proc~oqp_apply_basis->interface~tagarray_get_data map_shell2basis_set map_shell2basis_set proc~oqp_apply_basis->map_shell2basis_set none~bcast par_env_t%bcast proc~oqp_apply_basis->none~bcast none~init~14 par_env_t%init proc~oqp_apply_basis->none~init~14 print_basis print_basis proc~oqp_apply_basis->print_basis 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 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

Called by

proc~~oqp_apply_basis~~CalledByGraph proc~oqp_apply_basis oqp_apply_basis proc~apply_basis_c apply_basis_C proc~apply_basis_c->proc~oqp_apply_basis

Source Code

  subroutine oqp_apply_basis(infos)
    use messages, only: with_abort
    use types, only: information
    use atomic_structure_m, only: atomic_structure
    use strings, only: fstring
    use oqp_tagarray_driver
    use iso_c_binding, only: c_char
    use parallel, only: par_env_t
    use basis_api, only: map_shell2basis_set, print_basis

    implicit none
    type(information), intent(inout) :: infos
    type(par_env_t) :: pe
    character(len=:), allocatable :: basis_file
    integer :: iw, i
    logical :: err
  !
  ! Section of Tagarray for the basis filename
  ! We are getting basis file name from Python via tagarray
  !
    character(len=1,kind=c_char), contiguous, pointer :: basis_filename(:)
    character(len=*), parameter :: subroutine_name = "oqp_apply_basis"
    character(len=*), parameter :: tags_general(1) = (/ character(len=80) :: &
          OQP_basis_filename /)
    call data_has_tags(infos%dat, tags_general, module_name, subroutine_name, with_abort)
    call tagarray_get_data(infos%dat, OQP_basis_filename, basis_filename)
    allocate(character(ubound(basis_filename,1)) :: basis_file)
    do i = 1, ubound(basis_filename,1)
       basis_file(i:i) = basis_filename(i)
    end do
!
!  ! Files open
!  ! 3. LOG: Write: Main output file
!  ! 5. BAS: read: Basis set library (internally)
!
    open (newunit=iw, file=infos%log_filename, position="append")
   !
    write(iw,'(/,20x,"++++++++++++++++++++++++++++++++++++++++")')
    write(iw,'(  22X,"MODULE: apply_basis ")')
    write(iw,'(  22X,"Setting up basis set information")')
    write(iw,'(20x,"++++++++++++++++++++++++++++++++++++++++")')
    call pe%init(infos%mpiinfo%comm, infos%mpiinfo%usempi)
    call map_shell2basis_set(infos, infos%basis)
!    if (pe%rank == 0) then
!      call infos%basis%from_file(basis_file, infos%atoms, err)
!      infos%control%basis_set_issue = err
!    endif

!    call infos%basis%basis_broadcast(infos%mpiinfo%comm, infos%mpiinfo%usempi)

    if (sum(infos%basis%ecp_zn_num)>0) then
      call pe%bcast(infos%mol_prop%nelec, 1)
      call pe%bcast(infos%mol_prop%nelec_A, 1)
      call pe%bcast(infos%mol_prop%nelec_B, 1)
      call pe%bcast(infos%mol_prop%nocc, 1)
    end if

! Checking error of basis set reading..
!    call pe%bcast(infos%control%basis_set_issue, 1)

    write(iw,'(/5X,"Basis Sets options"/&
                  &5X,18("-")/&
                  &5X,"Basis Sets: ",A/&
                  &5X,"Number of Shells  =",I8,5X,"Number of Primitives  =",I8/&
                  &5X,"Number of Basis Set functions  =",I8/&
                  &5X,"Maximum Angluar Momentum =",I8/)') &
                    trim(basis_file), &
                    infos%basis%nshell, infos%basis%nprim, &
                    infos%basis%nbf, infos%basis%mxam
    close (iw)

    call print_basis(infos)


  end subroutine oqp_apply_basis