huckel_guess Subroutine

public subroutine huckel_guess(ovl, orbitals, infos, basis, huckel_basis)

Uses

  • proc~~huckel_guess~~UsesGraph proc~huckel_guess huckel_guess module~basis_tools basis_tools proc~huckel_guess->module~basis_tools module~constants constants proc~huckel_guess->module~constants module~guess guess proc~huckel_guess->module~guess module~int1 int1 proc~huckel_guess->module~int1 module~mathlib mathlib proc~huckel_guess->module~mathlib module~messages messages proc~huckel_guess->module~messages module~types types proc~huckel_guess->module~types module~basis_tools->module~constants 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~io_constants io_constants module~basis_tools->module~io_constants module~parallel parallel module~basis_tools->module~parallel module~precision precision module~basis_tools->module~precision module~constants->module~precision module~oqp_linalg oqp_linalg module~guess->module~oqp_linalg module~guess->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~mathlib->module~oqp_linalg module~mathlib->module~precision comm_IOFILE comm_IOFILE module~messages->comm_IOFILE comm_PAR comm_PAR module~messages->comm_PAR module~messages->module~io_constants module~messages->module~precision module~types->module~basis_tools iso_c_binding iso_c_binding module~types->iso_c_binding module~types->module~atomic_structure_m module~functionals functionals module~types->module~functionals module~types->module~parallel module~types->module~precision tagarray tagarray module~types->tagarray module~atomic_structure_m->iso_c_binding module~functionals->iso_c_binding module~functionals->module~precision 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~blas_wrap blas_wrap module~oqp_linalg->module~blas_wrap module~lapack_wrap lapack_wrap module~oqp_linalg->module~lapack_wrap module~parallel->iso_c_binding module~parallel->iso_fortran_env module~parallel->module~precision mpi mpi module~parallel->mpi module~precision->iso_fortran_env 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~mod_gauss_hermite->module~precision module~rys->module~constants module~rys->module~precision module~rys_lut rys_lut module~rys->module~rys_lut

Arguments

Type IntentOptional Attributes Name
real(kind=dp) :: ovl(*)
real(kind=dp) :: orbitals(*)
type(information), intent(in) :: infos
type(basis_set), intent(in) :: basis
type(basis_set), intent(in) :: huckel_basis

Calls

proc~~huckel_guess~~CallsGraph proc~huckel_guess huckel_guess interface~show_message show_message proc~huckel_guess->interface~show_message interface~unpack_matrix unpack_matrix proc~huckel_guess->interface~unpack_matrix proc~basis_overlap basis_overlap proc~huckel_guess->proc~basis_overlap proc~corresponding_orbital_projection corresponding_orbital_projection proc~huckel_guess->proc~corresponding_orbital_projection proc~diag_symm_full diag_symm_full proc~huckel_guess->proc~diag_symm_full proc~matrix_invsqrt matrix_invsqrt proc~huckel_guess->proc~matrix_invsqrt proc~mksphar mksphar proc~huckel_guess->proc~mksphar proc~oqp_dgemm_i64 oqp_dgemm_i64 proc~huckel_guess->proc~oqp_dgemm_i64 proc~oqp_dgeqrf_i64 oqp_dgeqrf_i64 proc~huckel_guess->proc~oqp_dgeqrf_i64 proc~oqp_dormqr_i64 oqp_dormqr_i64 proc~huckel_guess->proc~oqp_dormqr_i64 proc~oqp_dsymm_i64 oqp_dsymm_i64 proc~huckel_guess->proc~oqp_dsymm_i64 proc~orthogonal_transform_sym orthogonal_transform_sym proc~huckel_guess->proc~orthogonal_transform_sym proc~overlap overlap proc~huckel_guess->proc~overlap proc~unpack_f90 UNPACK_F90 interface~unpack_matrix->proc~unpack_f90 none~alloc2 shpair_t%alloc2 proc~basis_overlap->none~alloc2 none~fetch_by_id shell_t%fetch_by_id proc~basis_overlap->none~fetch_by_id none~shell_pair2 shpair_t%shell_pair2 proc~basis_overlap->none~shell_pair2 proc~comp_kin_ovl_int1_prim comp_kin_ovl_int1_prim proc~basis_overlap->proc~comp_kin_ovl_int1_prim proc~update_rectangular_matrix update_rectangular_matrix proc~basis_overlap->proc~update_rectangular_matrix proc~corresponding_orbital_projection->interface~show_message proc~corresponding_orbital_projection->proc~diag_symm_full proc~corresponding_orbital_projection->proc~oqp_dgemm_i64 proc~oqp_dswap_i64 oqp_dswap_i64 proc~corresponding_orbital_projection->proc~oqp_dswap_i64 proc~oqp_dsyrk_i64 oqp_dsyrk_i64 proc~corresponding_orbital_projection->proc~oqp_dsyrk_i64 proc~schmd schmd proc~corresponding_orbital_projection->proc~schmd proc~diag_symm_full->interface~show_message dsyev dsyev proc~diag_symm_full->dsyev proc~matrix_invsqrt->interface~show_message proc~diag_symm_packed diag_symm_packed proc~matrix_invsqrt->proc~diag_symm_packed proc~mksphar->interface~show_message proc~oqp_dgemm_i64->interface~show_message dgemm dgemm proc~oqp_dgemm_i64->dgemm proc~oqp_dgeqrf_i64->interface~show_message dgeqrf dgeqrf proc~oqp_dgeqrf_i64->dgeqrf proc~oqp_dormqr_i64->interface~show_message dormqr dormqr proc~oqp_dormqr_i64->dormqr proc~oqp_dsymm_i64->interface~show_message dsymm dsymm proc~oqp_dsymm_i64->dsymm 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 none~alloc~2 shpair_t%alloc proc~overlap->none~alloc~2 proc~overlap->none~fetch_by_id none~shell_pair shpair_t%shell_pair proc~overlap->none~shell_pair proc~overlap->proc~comp_kin_ovl_int1_prim proc~update_triang_matrix update_triang_matrix proc~overlap->proc~update_triang_matrix proc~doquadgausshermite doQuadGaussHermite proc~comp_kin_ovl_int1_prim->proc~doquadgausshermite proc~diag_symm_packed->interface~show_message dspev dspev proc~diag_symm_packed->dspev dspevx dspevx proc~diag_symm_packed->dspevx proc~oqp_dswap_i64->interface~show_message dswap dswap proc~oqp_dswap_i64->dswap proc~oqp_dsyrk_i64->interface~show_message dsyrk dsyrk proc~oqp_dsyrk_i64->dsyrk 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~schmd->interface~show_message proc~schmd->proc~oqp_dgeqrf_i64 proc~oqp_dorgqr_i64 oqp_dorgqr_i64 proc~schmd->proc~oqp_dorgqr_i64 proc~unpack_f90->interface~show_message proc~unpack_f90->proc~oqp_dtpttr_i64 abrt abrt proc~doquadgausshermite->abrt proc~oqp_dorgqr_i64->interface~show_message dorgqr dorgqr proc~oqp_dorgqr_i64->dorgqr

Called by

proc~~huckel_guess~~CalledByGraph proc~huckel_guess huckel_guess proc~guess_huckel guess_huckel proc~guess_huckel->proc~huckel_guess proc~guess_huckel_c guess_huckel_C proc~guess_huckel_c->proc~guess_huckel

Source Code

 subroutine huckel_guess(ovl, orbitals, infos, basis, huckel_basis)

   use constants, only: tol_int
   use types,     only: information
   use messages,  only: show_message, WITH_ABORT
   use mathlib, only: matrix_invsqrt
   use basis_tools, only: basis_set
   use int1, only: basis_overlap
   use guess, only: corresponding_orbital_projection

   implicit none

   type(information), intent(in) :: infos
   type(basis_set), intent(in) :: basis, huckel_basis

   real(kind=dp) :: ovl(*), orbitals(*)
   integer :: nat, i, ok, l0, l0co, nbf, nbf2, nbf_co, nact, ndoc, nproj

   real(kind=dp), allocatable :: scr(:)
   real(kind=dp), allocatable :: q(:)
   real(kind=dp), allocatable :: vec(:,:)
   real(kind=dp), allocatable :: sco(:,:)

   nbf = basis%nbf
   nbf2 = nbf*(nbf+1)/2
   nat = infos%mol_prop%natom

!  Number of orbitals in MINI basis used in Huckel
   nbf_co = huckel_basis%nbf

   allocate(scr(nbf), &
            q(nbf*nbf), &
            vec(nbf_co,nbf_co),     &
            sco(nbf_co,nbf),     &
            stat=ok)
   if (ok/=0) call show_message('Cannot allocate memory', WITH_ABORT)

!  Get overlap between the minimal basis set and the input basis
   call basis_overlap(sco, basis, huckel_basis, tol=log(10.0d0)*tol_int)

   do i = 1, nbf
     sco(:,i) = sco(:,i)*basis%bfnrm(i) * huckel_basis%bfnrm
   end do

!  Determine which orbitals should be projected
   if (infos%control%scftype == 1) then
     ndoc = infos%mol_prop%nelec/2
     nact = 0
   else if (infos%control%scftype >= 2) then
     ndoc = infos%mol_prop%nelec_b
     nact = infos%mol_prop%nelec_a-infos%mol_prop%nelec_b
   end if

!  Extended Huckel calculation in mini basis set
   call huckel_calc(huckel_basis, vec, l0co, nat, infos%atoms%zn, tol_int)

!  Do at most 5 virtuals from the huckel
   nproj = min(l0co,ndoc+nact+5)

!  Get canonical orbitals in input basis space
   call matrix_invsqrt(ovl, q, nbf, qrnk=l0)
   orbitals(1:nbf*nbf) = q(1:nbf*nbf)

!  Project minimal basis set guess onto the input canonical orbitals,
   call corresponding_orbital_projection(vec, sco, orbitals, ndoc, nact, nproj, nbf, nbf_co, l0)

   deallocate(vec, sco)

   call orthogonalize_orbitals(q, ovl, orbitals, nproj, l0, nbf, nbf)

 end subroutine huckel_guess