Integrating products of Bessel functions with an additional exponential or rational factor

We provide two Matlab programs to compute integrals of the form∫0∞e−cxxm∏i=1kJνi(aix)dxand∫0∞xmr2+x2∏i=1kJνi(aix)dx with Jνi(x)Jνi(x) the Bessel function of the first kind and (real) order νiνi. The parameter m   is a real number such that ∑iνi+m>−1∑iνi+m>−1 (to assure integrability near zero), r is real and the numbers c   and aiai are all strictly positive. The program can deliver accurate error estimates.Program summaryProgram title: BESSELINTR, BESSELINTCCatalogue identifier: AEAH_v1_0Program summary URL:http://cpc.cs.qub.ac.uk/summaries/AEAH_v1_0.htmlProgram obtainable from: CPC Program Library, Queen's University, Belfast, N. IrelandLicensing provisions: Standard CPC licence, http://cpc.cs.qub.ac.uk/licence/licence.htmlNo. of lines in distributed program, including test data, etc.: 1601No. of bytes in distributed program, including test data, etc.: 13 161Distribution format: tar.gzProgramming language:   Matlab (version ⩾6.5), Octave (version ⩾2.1.69⩾2.1.69)Computer: All supporting Matlab or OctaveOperating system: All supporting Matlab or OctaveRAM: For k   Bessel functions our program needs approximately (500+140k500+140k) double precision variablesClassification: 4.11Nature of problem: The problem consists in integrating an arbitrary product of Bessel functions with an additional rational or exponential factor over a semi-infinite interval. Difficulties arise from the irregular oscillatory behaviour and the possible slow decay of the integrand, which prevents truncation at a finite point.Solution method: The interval of integration is split into a finite and infinite part. The integral over the finite part is computed using Gauss–Legendre quadrature. The integrand on the infinite part is approximated using asymptotic expansions and this approximation is integrated exactly with the aid of the upper incomplete gamma function. In the case where a rational factor is present, this factor is first expanded in a Taylor series around infinity.Restrictions:   Some (and eventually all) numerical accuracy is lost when one or more of the parameters r,c,air,c,ai or vivi grow very large, or when r becomes small.Running time: Less than 0.02 s for a simple problem (two Bessel functions, small parameters), a few seconds for a more complex problem (more than six Bessel functions, large parameters), in Matlab 7.4 (R2007a) on a 2.4 GHz AMD Opteron Processor 250.References:[1]J. Van Deun, R. Cools, Algorithm 858: Computing infinite range integrals of an arbitrary product of Bessel functions, ACM Trans. Math. Software 32 (4) (2006) 580–596.