program bestime
use fvn_fnlib
implicit none

real(kind=dp_kind) :: x,xmin,xmax,xstep,t1,t2,bes,b(51)
integer(kind=ip_kind) :: npoints,iter,n,i,it,nn


iter=10000
n=10
npoints=200
xmin=-50.
xmax=50.  
xstep=(xmax-xmin)/dble(npoints)

write(*,*) "Computation time comparison between bsjn and besrj"
write(*,*) "bsjn is faster when computing J_n for only one value of n"
write(*,*) "besrj is faster when computing J_n for n=0 to n"

call cpu_time(t1)
do it=1,iter 
  do i=1,npoints
      x=xmin+i*xstep
      do nn=0,n
        bes=bsjn(nn,x)
      end do
  end do
end do
call cpu_time(t2)
write (*,'("Calculating "(I10)" values of J_n for n=0 to "(I3)" with bsjn :")') iter*npoints,n
write(*,*) t2-t1

call cpu_time(t1)
do it=1,iter 
  do i=1,npoints
      x=xmin+i*xstep
        bes=bsjn(n,x)
  end do
end do
call cpu_time(t2)
write (*,'("Calculating "(I10)" values of J_n for n="(I3)" with bsjn :")') iter*npoints,n
write(*,*) t2-t1


! Calcultate J_n(x) on 200 points for n=0 to 50 with besrj
! do it iter times to have a sufficient time value
call cpu_time(t1)
do it=1,iter
  do i=1,npoints
      x=xmin+i*xstep
        call besrj(x,n+1,b)
  end do
end do
call cpu_time(t2)
write (*,'("Calculating "(I10)" values of J_n for n=0 to "(I3)" with besrj :")') iter*npoints,n
write(*,*) t2-t1




end program bestime