vvsa Subroutine

public subroutine vvsa(Va, x, Pv)

Compute parabolic cylinder function Vv(x) for small argument

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: Va

Order

real(kind=wp), intent(in) :: x

Argument

real(kind=wp), intent(out) :: Pv

Vv(x)


Calls

proc~~vvsa~~CallsGraph proc~vvsa specfun_module::vvsa proc~gamma2 specfun_module::gamma2 proc~vvsa->proc~gamma2

Called by

proc~~vvsa~~CalledByGraph proc~vvsa specfun_module::vvsa proc~pbvv specfun_module::pbvv proc~pbvv->proc~vvsa

Source Code

   subroutine vvsa(Va,x,Pv)

      real(wp),intent(in) :: Va !! Order
      real(wp),intent(in) :: x !! Argument
      real(wp),intent(out) :: Pv !! Vv(x)

      real(wp) :: a0 , ep , fac , g1 , ga0 , gm , gw ,  &
                  r , r1 , sv , sv0 , v1 , va0 ,   &
                  vb0 , vm
      integer :: m

      real(wp),parameter :: eps = 1.0e-15_wp

      ep = exp(-0.25_wp*x*x)
      va0 = 1.0_wp + 0.5_wp*Va
      if ( x/=0.0_wp ) then
         a0 = 2.0_wp**(-0.5_wp*Va)*ep/(twopi)
         sv = sin(-(Va+0.5_wp)*pi)
         v1 = -0.5_wp*Va
         call gamma2(v1,g1)
         Pv = (sv+1.0_wp)*g1
         r = 1.0_wp
         fac = 1.0_wp
         do m = 1 , 250
            vm = 0.5_wp*(m-Va)
            call gamma2(vm,gm)
            r = r*sq2*x/m
            fac = -fac
            gw = fac*sv + 1.0_wp
            r1 = gw*r*gm
            Pv = Pv + r1
            if ( abs(r1/Pv)<eps .and. gw/=0.0_wp ) exit
         enddo
         Pv = a0*Pv
      elseif ( va0<=0.0_wp .and. va0==int(va0) .or. Va==0.0_wp ) then
         Pv = 0.0_wp
      else
         vb0 = -0.5_wp*Va
         sv0 = sin(va0*pi)
         call gamma2(va0,ga0)
         Pv = 2.0_wp**vb0*sv0/ga0
      endif
   end subroutine vvsa