bessel_i (2, %i*x) or bessel_i(2, x) with integrate --> noun form



The conditional return is cute, but it will flummox diff, solve, and ...

Examples:

(%i62) integrate(bessel_i(a,x),x);
(%o62) if equal(a,-2) then (hypergeometric([3/2],[5/2,3],x^2/4)*x^3)/24 else (hypergeometric([(a+1)/2],[a+1,(a+3)/2],x^2/4)*x^(a+1))/(2^a*gamma(a+2))

(%i63) integrate(bessel_i(5,x),x);
(%o63) (hypergeometric([3],[4,6],x^2/4)*x^6)/23040

(%i64) integrate(bessel_i(-2,x),x);
(%o64) (hypergeometric([3/2],[5/2,3],x^2/4)*x^3)/24

What's the story with all this rat stuff....

(%i68) integrate(bessel_i(2,x),x,0,10.0);
rat: replaced 10.0 by 10/1 = 10.0
rat: replaced 10.0 by 10/1 = 10.0
rat: replaced 10.0 by 10/1 = 10.0
rat: replaced 2348.932064591061 by 241940/103 = 2348.932038834952
(%o68) 241940/103

The code passes a few tests such as

(%i94) integrate(bessel_i(-2,x),x);
(%o94) (hypergeometric([3/2],[5/2,3],x^2/4)*x^3)/24
(%i95) bessel_i(-2,x) - diff(%,x);
(%o95) -(hypergeometric([5/2],[7/2,4],x^2/4)*x^4)/240-(hypergeometric([3/2],[5/2,3],x^2/4)*x^2)/8+bessel_i(2,x)
(%i96) taylor(%,x,0,13);
(%o96)/T/ 0+...

;;----start---------------------------
(defun bessel-i-integral (a z)
  (mfuncall 'mcond
	    (meqp a -2)
	    (mul (opcons 'mexpt z 3)
		 (div 1 24)
		 (opcons '$hypergeometric
			 (opcons 'mlist (div 3 2))
			 (opcons 'mlist  (div 5 2) 3)
			 (div (mul z z) 4)))
	    t 
	    (mul
	     (opcons 'mexpt z (add a 1))
	     (opcons '$hypergeometric
		     (opcons 'mlist (div (add a 1) 2))
		     (opcons 'mlist (div (add a 3) 2) (add a 1))
		     (div (mul z z) 4))
	     (div 
	      1
	      (mul 
	       (opcons 'mexpt 2 a)
	       (opcons '%gamma (add a 2)))))))
	       		      
(putprop '%bessel_i `((a z) nil ,#'bessel-i-integral) 'integral) 
 
;;---end---------------------------------------------------
--Barton