More integrals for special functions



I have implemented a small extension to get more integrals for special
functions. The idea is to support a special routine which is looked up
in the routine INTFORM, if Maxima does not get a result for the desired
special function.

This is the implementation in INTFORM:

...

((intform (cadr expres)))
;; Maxima had no sucess to integrate an integrand with a
;; special function. Try to look up a special integration routine.
;; At this point of the algorithm this works only for functions 
;; with an integral on the property list.
((setq op (safe-get (mop expres) 'integrate-special))
 (funcall op exp var))
(t nil)))

...

And this is a small extension to get the integral for the Incomplete
Gamma function and a power function:

(defprop %gamma_incomplete integrate-gamma-incomplete integrate-special)

(defun integrate-gamma-incomplete (expr var)
  (let (part fact rest gamma l)
    ;; Decompose the integrand into [fact,gamma_incomplete(a,z)]
    (setq part (if (mtimesp expr)
                   ($partition expr (setq gamma 
                                          (isinop expr '%
gamma_incomplete)))
                   (list '(mlist simp) 1 expr)))
    (setq fact (cadr part)
          rest (caddr part))
    (cond ((and (alike1 gamma rest) ; no power of Incomplete gamma
                (setq l (m2-c*z^a fact))) ; match the factor
           ;; Integrate z^(v-1)*gamma_incomplete(a,z)
           (let ((v (add 1 (cdras 'a l)))
                 (a (cadr gamma)))
             (if (not (equal v 0))
                 (mul (inv v)
                      (sub (mul (power var v)
                                (take '(%gamma_incomplete) a var))
                           (take '(%gamma_incomplete) (add a v) var)))
                 nil)))
          (t nil))))

The implementation of the special integration routine is done in a way
which is independent from the global special variables of the
integrator. (With the exception of the global special variable VAR.) The
only assumption is that the integrand only has factors. This is true
when we are in the routine INTFORM.

This is an example for a new result:

(%i4) integrate(x*gamma_incomplete(a,x),x);
(%o4) (gamma_incomplete(a,x)*x^2-gamma_incomplete(a+2,x))/2


Furthermore, I have implemented more integrals for the Bessel J
function. The following types are new:

integrate(z^a*bessel_j(v,b*z^r),z)

integrate(z^a*exp(%i*b*z^r)*bessel_j(v,b*z^r),z)
integrate(z^a*exp(-%i*b*z^r)*bessel_j(v,b*z^r),z)

integrate(z^a*sin(b*z^r+d)*bessel_j(v,b*z^r),z)
integrate(z^a*cos(b*z^r+d)*bessel_j(v,b*z^r),z)

integrate(z^a*bessel_j(v,b*z^r)^2,z)
integrate(z^a*bessel_j(u,b*z^r)*bessel_j(v,b*z^r),z)

The general results are in terms of hypergeometric functions. To get
some more simple results special cases are implemented too. The special
cases might be not necessary, if we improve the code to simplify
hypergeometric functions in the future.

These are only some arbitrary examples:

(%i5) integrate(x^3*bessel_j(2,x),x);
(%o5) bessel_j(3,x)*x^3

(%i6) integrate(x^3*bessel_j(2,x)^2,x);
(%o6) hypergeometric([5/2,4],[3,5,5],-x^2)*x^8/512

(%i7) integrate(x^2*bessel_j(2,x)^2,x);
(%o7) hypergeometric([5/2,7/2],[3,9/2,5],-x^2)*x^7/448

(%i8) integrate(x*bessel_j(2,x)^2,x);
(%o8) (bessel_j(2,x)^2-bessel_j(1,x)*bessel_j(3,x))*x^2/2

(%i9) integrate(sqrt(x)*bessel_j(1,x)*bessel_j(2,x),x);
(%o9) hypergeometric([9/4,5/2],[3,13/4,4],-x^2)*x^(9/2)/72

(%i27) integrate(sin(x)*bessel_j(0,x),x);
(%o27) hypergeometric([3/4,5/4],[3/2,3/2,2],-x^2)*x^2/2

(%i28) integrate(sin(x+%pi/3)*bessel_j(0,x),x);
(%o28) x*(hypergeometric([3/4,5/4],[3/2,3/2,2],-x^2)*x/2
         +sqrt(3)*hypergeometric([1/4,3/4],[1/2,1,3/2],-x^2))
        /2

The testsuite has no problems. We get no slow down, because the new code
is only executed, if a special function is part of the integrand. We
only get new integrals.

If the implementation is of interest I would like to commit the code.

Dieter Kaiser