mulQuadGaussHermite Subroutine

public subroutine mulQuadGaussHermite(tint, t, rij, ri, rj, r, ni, nj, m)

@brief Gauss-Hermite quadrature using minimum point formula @details Compute: xint = sum( w(1:npts,npts) * (h(1:npts,npts)t+dxi)(ni-1) * (h(1:npts,npts)t+dxj)(nj-1) ) yint = sum( w(1:npts,npts) * (h(1:npts,npts)*t+dyi)(ni-1) * (h(1:npts,npts)t+dyj)(nj-1) ) zint = sum( w(1:npts,npts) * (h(1:npts,npts)t+dzi)(ni-1) * (h(1:npts,npts)*t+dzj)(nj-1) )

Note

Use of I functions will let NI run up to 7 (S=1, P=2, ..., I=7) Use of I functions will let NJ run up to 9 (for K.E. ints) Use of I functions requires NPTS=8 to do kinetic energy integrals. @param[out] xint x-component of the integral @param[out] yint y-component of the integral @param[out] zint z-component of the integral @param[out] t inverse square root of total exponent @param[out] x0 x-coord. of the center of primitive pair @param[out] y0 y-coord. of the center of primitive pair @param[out] z0 z-coord. of the center of primitive pair @param[out] xi x-coord. of the center of 1st primitive @param[out] yi y-coord. of the center of 1st primitive @param[out] zi z-coord. of the center of 1st primitive @param[out] xj x-coord. of the center of 2nd primitive @param[out] yj y-coord. of the center of 2nd primitive @param[out] zj z-coord. of the center of 2nd primitive @param[out] ni current angular momentum on center i @param[out] nj current angular momentum on center j

Note

based on STVINT from int1.src @author Vladimir Mironov @date Sep, 2018 Initial release

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(out) :: tint(3,0:*)
real(kind=dp), intent(in) :: t
real(kind=dp), intent(in) :: rij(3)
real(kind=dp), intent(in) :: ri(3)
real(kind=dp), intent(in) :: rj(3)
real(kind=dp), intent(in) :: r(3)
integer, intent(in) :: ni
integer, intent(in) :: nj
integer, intent(in) :: m

Calls

proc~~mulquadgausshermite~~CallsGraph proc~mulquadgausshermite mulQuadGaussHermite abrt abrt proc~mulquadgausshermite->abrt

Called by

proc~~mulquadgausshermite~~CalledByGraph proc~mulquadgausshermite mulQuadGaussHermite proc~comp_allmult_int1_prim comp_allmult_int1_prim proc~comp_allmult_int1_prim->proc~mulquadgausshermite proc~comp_mult_int1_prim comp_mult_int1_prim proc~comp_mult_int1_prim->proc~mulquadgausshermite proc~multipole_integrals multipole_integrals proc~multipole_integrals->proc~comp_allmult_int1_prim proc~electric_moments electric_moments proc~electric_moments->proc~multipole_integrals proc~get_td_transition_dipole get_td_transition_dipole proc~get_td_transition_dipole->proc~multipole_integrals proc~get_transition_dipole get_transition_dipole proc~get_transition_dipole->proc~multipole_integrals proc~tdhf_energy tdhf_energy proc~tdhf_energy->proc~get_td_transition_dipole proc~tdhf_mrsf_energy tdhf_mrsf_energy proc~tdhf_mrsf_energy->proc~get_transition_dipole proc~tdhf_sf_energy tdhf_sf_energy proc~tdhf_sf_energy->proc~get_transition_dipole proc~tdhf_energy_c tdhf_energy_C proc~tdhf_energy_c->proc~tdhf_energy proc~tdhf_mrsf_energy_c tdhf_mrsf_energy_C proc~tdhf_mrsf_energy_c->proc~tdhf_mrsf_energy proc~tdhf_sf_energy_c tdhf_sf_energy_C proc~tdhf_sf_energy_c->proc~tdhf_sf_energy