Re: NIntegrate and speed

*To*: mathgroup at smc.vnet.net*Subject*: [mg116785] Re: NIntegrate and speed*From*: Andrew Moylan <amoylan at wolfram.com>*Date*: Mon, 28 Feb 2011 04:59:32 -0500 (EST)

You're right, the inner integral is a Bessel function: In[1]:= Integrate[Cos[k R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}, Assumptions -> {R > 0, k > 0}] Out[1]= \[Pi] BesselJ[0, k R] So, you can use a one-dimensional NIntegrate: In[2]:= R = 8000; Z = 1; rd = 3500; In[4]:= NIntegrate[ Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5 (\[Pi] BesselJ[0, k R]), {k, 0, \[Infinity]}] // Timing Out[4]= {0.474709, 0.000424068} You can try to speed up the result by choosing a different method. In this case, NIntegrate selects "ExtrapolatingOscillatory" as the method: In[5]:= NIntegrate[ Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5 (\[Pi] BesselJ[0, k R]), {k, 0, \[Infinity]}, Method -> "ExtrapolatingOscillatory"] // Timing Out[5]= {0.434378, 0.000424068} The new "LevinRule" general oscillatory method is actually faster for this integral: In[6]:= NIntegrate[ Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5 (\[Pi] BesselJ[0, k R]), {k, 0, \[Infinity]}, Method -> "LevinRule"] // Timing Out[6]= {0.135486, 0.000424068} To do many NIntegrates in a loop, Parallelize would be the first thing to try: In[7]:= Clear[R]; Parallelize[ Table[NIntegrate[ Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5 (\[Pi] BesselJ[0, k R]), {k, 0, \[Infinity]}, Method -> "LevinRule"], {R, 6000, 9000, 1000}]] Out[8]= {0.00052313, 0.000470404, 0.000424068, 0.000383576} I hope this helps, Andrew Moylan Wolfram Research On Feb 27, 2011, at 8:35 PM, Marco Masi wrote: > I have the following problems with NIntegrate. > > 1) I would like to make the following double numerical integral converge without errors > > R = 8000; Z = 1; rd = 3500; > NIntegrate[Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5 (NIntegrate[Cos[k R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}]), {k, 0, \[Infinity]}] > > It tells non numerical values present and I don't understand why, since it evaluates finally a numerical value? 0.000424067 > > 2) Isn't the second integrand a cylindrical Bessel function of order 0? So, I expected that > NIntegrate[Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5 BesselJZero[0, k R], {k, 0, \[Infinity]}] doing the same job. But it fails to converge and gives 0.00185584- i4.96939*10^-18. Trying with WorkingPrecision didn't make things better. How can this be fixed? > > 3) The above Nintegrals will go into a loop and should be evaluated as fast as possible. How? With Compile, CompilationTarget -> "C", Paralleization, etc.? > > Any suggestions? > > Marco. >