get_ab_initio_density Subroutine

public subroutine get_ab_initio_density(alpha_density, alpha_orbital, beta_density, beta_orbital, infos, basis)

Uses

  • proc~~get_ab_initio_density~~UsesGraph proc~get_ab_initio_density get_ab_initio_density module~basis_tools basis_tools proc~get_ab_initio_density->module~basis_tools module~mathlib mathlib proc~get_ab_initio_density->module~mathlib module~messages messages proc~get_ab_initio_density->module~messages module~types types proc~get_ab_initio_density->module~types 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~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~oqp_linalg oqp_linalg 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~constants->module~precision module~functionals->iso_c_binding module~functionals->module~precision xc_f03_lib_m xc_f03_lib_m module~functionals->xc_f03_lib_m 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

@brief this will calculatte the density matrix

Arguments

Type IntentOptional Attributes Name
real(kind=dp) :: alpha_density(:)
real(kind=dp) :: alpha_orbital(:,:)
real(kind=dp), optional :: beta_density(:)
real(kind=dp), optional :: beta_orbital(:,:)
type(information), intent(in) :: infos
type(basis_set), intent(in) :: basis

Calls

proc~~get_ab_initio_density~~CallsGraph proc~get_ab_initio_density get_ab_initio_density interface~show_message show_message proc~get_ab_initio_density->interface~show_message proc~orb_to_dens orb_to_dens proc~get_ab_initio_density->proc~orb_to_dens interface~pack_matrix pack_matrix proc~orb_to_dens->interface~pack_matrix proc~oqp_dsyr2k_i64 oqp_dsyr2k_i64 proc~orb_to_dens->proc~oqp_dsyr2k_i64 proc~pack_f90 PACK_F90 interface~pack_matrix->proc~pack_f90 proc~oqp_dsyr2k_i64->interface~show_message dsyr2k dsyr2k proc~oqp_dsyr2k_i64->dsyr2k proc~pack_f90->interface~show_message proc~oqp_dtrttp_i64 oqp_dtrttp_i64 proc~pack_f90->proc~oqp_dtrttp_i64 proc~oqp_dtrttp_i64->interface~show_message dtrttp dtrttp proc~oqp_dtrttp_i64->dtrttp

Called by

proc~~get_ab_initio_density~~CalledByGraph proc~get_ab_initio_density get_ab_initio_density proc~guess_hcore guess_hcore proc~guess_hcore->proc~get_ab_initio_density proc~guess_huckel guess_huckel proc~guess_huckel->proc~get_ab_initio_density proc~guess_json guess_json proc~guess_json->proc~get_ab_initio_density proc~scf_driver scf_driver proc~scf_driver->proc~get_ab_initio_density proc~guess_hcore_c guess_hcore_C proc~guess_hcore_c->proc~guess_hcore proc~guess_huckel_c guess_huckel_C proc~guess_huckel_c->proc~guess_huckel proc~guess_json_c guess_json_C proc~guess_json_c->proc~guess_json proc~hf_energy hf_energy proc~hf_energy->proc~scf_driver

Source Code

subroutine get_ab_initio_density(alpha_density,alpha_orbital,beta_density,beta_orbital, infos, basis)
   use mathlib, only: orb_to_dens
   use messages, only: show_message, with_abort
   use types, only: information
   use basis_tools, only: basis_set
   implicit none
   type(information), intent(in) :: infos
   type(basis_set), intent(in) :: basis
   integer    :: ok
   integer    :: scftype, nocc, na, nb, nbasis, na2
   real(dp)   :: alpha_density(:), &
                 alpha_orbital(:,:)
   real(dp), optional :: beta_density(:), &
                         beta_orbital(:,:)
   real(dp), dimension(:), allocatable :: occno
!
   scftype = infos%control%scftype
   na   = infos%mol_prop%nelec_a
   nb   = infos%mol_prop%nelec_b
   nbasis = basis%nbf
   nocc = infos%mol_prop%nocc

   allocate(occno(nocc), stat=ok)
   if(ok/=0) call show_message('allocation of occno fails', with_abort)

   occno = 2
   if (scftype/=1) occno = 1

   na2 = nocc
   if (scftype/=1) na2 = na

!  alpha density is always calculated for RHF/ROHF/UHF
   call orb_to_dens(alpha_density,alpha_orbital,occno,na2,nbasis,nbasis)

   if (nb == 0) then
     beta_density = 0
     if (scftype == 2) beta_orbital = 0

   else
     select case (scftype)
     case (1)

     case (2)
       call orb_to_dens(beta_density,beta_orbital,occno,nb,nbasis,nbasis)

     case (3)
       call orb_to_dens(beta_density,alpha_orbital,occno,nb,nbasis,nbasis)
     end select
   end if

end subroutine get_ab_initio_density