@brief this will calculatte the density matrix
Type | Intent | Optional | 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 |
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