qstar Subroutine

public subroutine qstar(m, n, c, Ck, Ck1, Qs, Qt)

Compute Q*mn(-ic) for oblate radial functions with a small argument

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: m
integer, intent(in) :: n
real(kind=wp), intent(in) :: c
real(kind=wp), intent(in) :: Ck(200)
real(kind=wp), intent(in) :: Ck1
real(kind=wp), intent(out) :: Qs
real(kind=wp), intent(out) :: Qt

Called by

proc~~qstar~~CalledByGraph proc~qstar specfun_module::qstar proc~rmn2so specfun_module::rmn2so proc~rmn2so->proc~qstar proc~rswfo specfun_module::rswfo proc~rswfo->proc~rmn2so

Source Code

   subroutine qstar(m,n,c,Ck,Ck1,Qs,Qt)

      integer,intent(in) :: m
      integer,intent(in) :: n
      real(wp),intent(in) :: c
      real(wp),intent(in) :: Ck(200)
      real(wp),intent(in) :: Ck1
      real(wp),intent(out) :: Qs
      real(wp),intent(out) :: Qt

      real(wp) :: ap(200) , qs0 , r , s , sk
      integer :: i , ip , k , l

      ip = 1
      if ( n-m==2*int((n-m)/2) ) ip = 0
      r = 1.0_wp/Ck(1)**2
      ap(1) = r
      do i = 1 , m
         s = 0.0_wp
         do l = 1 , i
            sk = 0.0_wp
            do k = 0 , l
               sk = sk + Ck(k+1)*Ck(l-k+1)
            enddo
            s = s + sk*ap(i-l+1)
         enddo
         ap(i+1) = -r*s
      enddo
      qs0 = ap(m+1)
      do l = 1 , m
         r = 1.0_wp
         do k = 1 , l
            r = r*(2.0_wp*k+ip)*(2.0_wp*k-1.0_wp+ip)/(2.0_wp*k)**2
         enddo
         qs0 = qs0 + ap(m-l+1)*r
      enddo
      Qs = (-1)**ip*Ck1*(Ck1*qs0)/c
      Qt = -2.0_wp/Ck1*Qs

   end subroutine qstar