NIntegrate and speed
Marco Masi
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.
Albert Retey
>
> 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
many numeric functions in Mathematica try to evaluate their "body"
symbolically. In your case that will call NIntegrate with a symbolic k
which results in the messages you see. Actually these can be ignored in
that case, but of course it is nicer to avoid them, e.g. like this:
R = 8000; Z = 1; rd = 3500;
f[k_?NumericQ] :=
NIntegrate[Cos[k R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}]
NIntegrate[
Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5 (f[k]), {k, 0, \[Infinity]}]
> 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?
If I use the result Mathematica gives me for:
Integrate[Cos[k R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}]
the result is correct:
NIntegrate[
Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5 (\[Pi] BesselJ[0, R k]), {k,
0, \[Infinity]}]
so I guess you just didn't choose the correct Bessel function. Of course
that approach is much faster than the nested NIntegrate...
> 3) The above Nintegrals will go into a loop and should be evaluated
> as fast as possible. How? With Compile, CompilationTarget -> "C",
> Paralleization, etc.?
NIntegrate itself can't be compiled but the integrand will usually be
compiled automatically (see the Compiled option for NIntegrate). I think
ParallelDo will be your friend, but you have of course to take care that
there are no dependencies in the evaluations of the body for different
values of the iteration variable...
hth,
albert
Andrew Moylan
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.
>
Daniel Lichtblau
----- Original Message -----
> From: "Marco Masi" <marco...@ymail.com>
> To: math...@smc.vnet.net
> Sent: Sunday, February 27, 2011 3:35:46 AM
> Subject: NIntegrate and speed
> 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
You presented it as an iterated integral. Mathematically that is fine but from a language semantics viewpoint you now have a small problem. It is that the outer integral cannot correctly do symbolic analysis of its integrand but it may try to do so anyway. In essence, the outer integrand "looks" to be nonnumerical until actual values areplugged in for the outer variable of integration.
There are (at least) two ways to work around this. One is to recast as a double (as opposed to iterated) integral.
Timing[i1 =
NIntegrate[
Exp[-k*Abs[Z]]/(1 + (k*rd)^2)^(3/2)*
Cos[k*R*Sin[\[Theta]]], {\[Theta], 0, \[Pi]}, {k,
0, \[Infinity]}]]
{39.733, 0.0004240679194556497}
An alternative is to define the inner function as a black box that only evaluates for numeric input. In that situation the outer NIntegrate will not attempt to get cute with its integrand.
ff[t_?NumericQ] :=
NIntegrate[Cos[t* R*Sin[\[Theta]]], {\[Theta], 0, \[Pi]}]
In[90]:= Timing[
i2 = NIntegrate[
Exp[-k* Abs[Z]]/(1 + (k* rd)^2)^1.5 *ff[k], {k, 0, \[Infinity]}]]
Out[90]= {26.63, 0.0004240673399701612}
> 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?
Use the correct symbolic form of the inner integral. It involves BesselJ rather than BesselJZero.
In[91]:= ff2[t_] =
Integrate[Cos[t* Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
Assumptions -> Element[t, Reals]]
Out[91]= \[Pi] BesselJ[0, Abs[t]]
In[92]:= Timing[
i3 = NIntegrate[
Exp[-k Abs[Z]]/(1 + (k *rd)^2)^(3/2)* ff2[k*R], {k,
0, \[Infinity]}]]
Out[92]= {0.7019999999999982, 0.0004240679192434893}
Not surprisingly this is much faster, and will help to get you past the speed bumps you allude to below.
>
> 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.
Compile will not help because most of the time will be spent in NIntegrate code called from the virtual machine of the run time library (that latter if you compile to C). Evaluating in parallel should help. Also there might be option settings that allow NIntegrate to handle this faster than by default but without significant degradation in quality of results. Here is a set of timings using a few different methods, and have PrecisionGoal set fairly low (three digits).
In[109]:= Table[
Timing[NIntegrate[
Exp[-k Abs[Z]]/(1 + (k *rd)^2)^(3/2)* \[Pi] BesselJ[0,
Abs[k*R]], {k, 0, \[Infinity]}, PrecisionGoal -> 3,
Method -> meth]], {meth, {Automatic, "DoubleExponential",
"Trapezoidal", "RiemannRule"}}]
During evaluation of In[109]:= NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in k near {k} = {0.0002724458978988764}. NIntegrate obtained 0.00042483953211734914` and 0.000012161444876769691` for the integral and error estimates. >>
Out[109]= {{0.6709999999999923,
0.0004240678889181539}, {0.0150000000000432,
0.0004240644189596502}, {0.03100000000000591,
0.0004240644189596502}, {0.04699999999996862,
0.0004248395321173491}}
I rather suspect there are more option tweaks that could make this faster still without appreciable degradation in quality of results.
Daniel Lichtblau
Wolfram Research
Bill Rowe
>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
You are getting an error since when the inner integral is
initially sampled, k hasn't been assigned a value. You can do
the double integral without the error by doing:
In[3]:= NIntegrate[
Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5 Cos[k R Sin[\[Theta]]], {\[Theta],
0, \[Pi]}, {k, 0, \[Infinity]}]
Out[3]= 0.000424068
That is, to do a double integral you do not need to type
NIntegrate twice. And by using this syntax, the algorithm being
used can immediately see k is to be assigned numerical values
during the integration.
Bob Hanlon
$Version
"8.0 for Mac OS X x86 (64-bit) (November 6, 2010)"
Clear[R, k]
Assuming[Element[{k, R}, Reals],
Integrate[Cos[k R Sin[t]], {t, 0, Pi}]]
Pi*BesselJ[0, Abs[k]*Abs[R]]
R = 8000; Z = 1; rd = 3500;
NIntegrate[Exp[-k Abs[Z]]/(1 + (k rd)^2)^1.5
Pi* BesselJ[0, Abs[k] Abs[R]], {k, 0, Infinity}]
0.000424068
I get the same results with version 7.01.0
Bob Hanlon
---- Marco Masi <marco...@ymail.com> 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.
DrMajorBob
(Bob Hanlon?), but didn't notice it in yours.
It's Short Attention Span Theater, here.
Bobby
On Mon, 28 Feb 2011 11:43:52 -0600, Daniel Lichtblau <da...@wolfram.com>
wrote:
> DrMajorBob wrote:
>> The inner integration is much faster done symbolically:
>> Timing[Clear[ff];
>> ff[t_] =
>> Integrate[Cos[t*R*Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
>> Assumptions -> {t >= 0}];
>> {ff[t], i2 =
>> NIntegrate[
>> Exp[-k*Abs[Z]]/(1 + (k*rd)^2)^1.5*ff[k], {k, 0, \[Infinity]}]}
>> ]
>> {0.725109, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
>> Bobby
>
> Didn't I do that? Under item (2) somewhere?
>
> Daniel
>
>
>> On Mon, 28 Feb 2011 03:59:54 -0600, Daniel Lichtblau <da...@wolfram.com>
>> wrote:
>>
>>>
>>>
>>> ----- Original Message -----
>>>> From: "Marco Masi" <marco...@ymail.com>
>>>> To: math...@smc.vnet.net
>>>> Sent: Sunday, February 27, 2011 3:35:46 AM
>>>> Subject: NIntegrate and speed
>>>> 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
>>>
>>> You presented it as an iterated integral. Mathematically that is fine
>>> but from a language semantics viewpoint you now have a small problem.
>>> It is that the outer integral cannot correctly do symbolic analysis of
>>> its integrand but it may try to do so anyway. In essence, the outer
>>> integrand "looks" to be nonnumerical until actual values areplugged in
>>> for the outer variable of integration.
>>>
>>> There are (at least) two ways to work around this. One is to recast as
>>> a double (as opposed to iterated) integral.
>>>
>>> Timing[i1 =
>>> NIntegrate[
>>> Exp[-k*Abs[Z]]/(1 + (k*rd)^2)^(3/2)*
>>> Cos[k*R*Sin[\[Theta]]], {\[Theta], 0, \[Pi]}, {k,
>>> 0, \[Infinity]}]]
>>> {39.733, 0.0004240679194556497}
>>>
>>> An alternative is to define the inner function as a black box that
>>> only evaluates for numeric input. In that situation the outer
>>> NIntegrate will not attempt to get cute with its integrand.
>>>
>>> ff[t_?NumericQ] :=
>>> NIntegrate[Cos[t* R*Sin[\[Theta]]], {\[Theta], 0, \[Pi]}]
>>>
>>> In[90]:= Timing[
>>> i2 = NIntegrate[
>>> Exp[-k* Abs[Z]]/(1 + (k* rd)^2)^1.5 *ff[k], {k, 0, \[Infinity]}]]
>>> Out[90]= {26.63, 0.0004240673399701612}
>>>
>>>
>>>> 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?
>>>
>>> Use the correct symbolic form of the inner integral. It involves
>>> BesselJ rather than BesselJZero.
>>>
>>> In[91]:= ff2[t_] =
>>> Integrate[Cos[t* Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
>>> Assumptions -> Element[t, Reals]]
>>> Out[91]= \[Pi] BesselJ[0, Abs[t]]
>>>
>>> In[92]:= Timing[
>>> i3 = NIntegrate[
>>> Exp[-k Abs[Z]]/(1 + (k *rd)^2)^(3/2)* ff2[k*R], {k,
>>> 0, \[Infinity]}]]
>>> Out[92]= {0.7019999999999982, 0.0004240679192434893}
>>>
>>> Not surprisingly this is much faster, and will help to get you past
>>> the speed bumps you allude to below.
>>>
>>>>
>>>> 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.
>>>
Daniel Lichtblau
> Mea culpa... yes, you did. I saw it in another post after writing mine
> (Bob Hanlon?), but didn't notice it in yours.
>
> It's Short Attention Span Theater, here.
>
> Bobby
Yez could'a written "Mea maxima culpa". Then it would show up as "Mea
culpa in another program".
Daniel
DrMajorBob
Timing[Clear[ff];
ff[t_] =
Integrate[Cos[t*R*Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
Assumptions -> {t >= 0}];
{ff[t], i2 =
NIntegrate[
Exp[-k*Abs[Z]]/(1 + (k*rd)^2)^1.5*ff[k], {k, 0, \[Infinity]}]}
]
{0.725109, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
Bobby
On Mon, 28 Feb 2011 03:59:54 -0600, Daniel Lichtblau <da...@wolfram.com>
DrMajorBob
rather than 1.5. I see no reason that should make a huge difference, but
it obviously does!
The same change is no help in the double integral.
These are timings on my machine (3-yr-old iMac 24/2.8/4GB/1TB/SD).
R = 8000; Z = 1; rd = 3500;
Timing[Clear[ff];
ff[t_] =
Integrate[Cos[t*R*Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
Assumptions -> {t >= 0}];
{ff[t], i2 =
NIntegrate[
Exp[-k*Abs[Z]]/(1 + (k*rd)^2)^(3/2)*ff[k], {k, 0, \[Infinity]}]}]
{0.712826, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
(*Symbolic evaluation of ff ala Bobby*)
Timing[Clear[ff];
ff[t_] =
Integrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
Assumptions -> {t >= 0}]; {ff[t],
i2 = NIntegrate[(Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k,
0, \[Infinity]}]}]
{20.7246, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
(*Numerical evaluation of ff ala Kevin*)
R = 8000; Z = 1; rd = 3500;
Clear[ff];
ff[t_?NumericQ] :=
NIntegrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}];
Timing[{ff[t],
i2 = NIntegrate[(Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^(3/2), {k,
0, \[Infinity]}]}]
{34.8137, {ff[t], 0.000424067}}
Bobby
On Tue, 01 Mar 2011 05:10:53 -0600, Kevin J. McCann
<Kevin....@umbc.edu> wrote:
> A couple of comments. (1) Bobby, your times are amazingly fast. I am
> running a 3.6 GHz quad-core i5 processor with 4G of memory and my times
> are quite a bit longer than yours; (2) my timings with ff done
> numerically are about twice as as fast as the symbolic:
>
> (* Symbolic evaluation of ff ala Bobby *)
> R = 8000; Z = 1; rd = 3500;
> Timing[Clear[ff];
> ff[t_] = Integrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
> Assumptions -> {t >= 0}]; {ff[t],
> i2 = NIntegrate[(
> Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k, 0, \[Infinity]}]}]
>
> {22.907, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
>
> (* Numerical evaluation of ff ala Kevin *)
> R = 8000; Z = 1; rd = 3500;
> Clear[ff];
> ff[t_?NumericQ] :=
> NIntegrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}];
> Timing[{ff[t],
> i2 = NIntegrate[(
> Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k, 0, \[Infinity]}]}]
>
> {13.484, {ff[t], 0.000424067}}
>
> I guess it depends a lot on the computer. Incidentally, I am running
> 64-bit mma8 under 64-bit Vista.
>
> Kevin
DrMajorBob
whether Mathematica has cached the Integrate result.
NIntegrate is fast, but Integrate can take 20 seconds.
Bobby
On Tue, 01 Mar 2011 14:10:19 -0600, DrMajorBob <btr...@austin.rr.com>
wrote:
> I may have posted the wrong code, but the only difference is using 3/2
> rather than 1.5. I see no reason that should make a huge difference, but
> it obviously does!
>
> The same change is no help in the double integral.
>
> These are timings on my machine (3-yr-old iMac 24/2.8/4GB/1TB/SD).
>
> R = 8000; Z = 1; rd = 3500;
> Timing[Clear[ff];
> ff[t_] =
> Integrate[Cos[t*R*Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
> Assumptions -> {t >= 0}];
> {ff[t], i2 =
> NIntegrate[
> Exp[-k*Abs[Z]]/(1 + (k*rd)^2)^(3/2)*ff[k], {k, 0, \[Infinity]}]}]
>
> {0.712826, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
>
> (*Symbolic evaluation of ff ala Bobby*)
> Timing[Clear[ff];
> ff[t_] =
> Integrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
> Assumptions -> {t >= 0}]; {ff[t],
> i2 = NIntegrate[(Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k,
> 0, \[Infinity]}]}]
>
> {20.7246, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
>
> (*Numerical evaluation of ff ala Kevin*)
> R = 8000; Z = 1; rd = 3500;
> Clear[ff];
> ff[t_?NumericQ] :=
> NIntegrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}];
> Timing[{ff[t],
> i2 = NIntegrate[(Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^(3/2), {k,
> 0, \[Infinity]}]}]
>
> {34.8137, {ff[t], 0.000424067}}
>
> Bobby
>
> On Tue, 01 Mar 2011 05:10:53 -0600, Kevin J. McCann
> <Kevin....@umbc.edu> wrote:
>
>> A couple of comments. (1) Bobby, your times are amazingly fast. I am
>> running a 3.6 GHz quad-core i5 processor with 4G of memory and my times
>> are quite a bit longer than yours; (2) my timings with ff done
>> numerically are about twice as as fast as the symbolic:
>>
>> (* Symbolic evaluation of ff ala Bobby *)
>> R = 8000; Z = 1; rd = 3500;
>> Timing[Clear[ff];
>> ff[t_] = Integrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
>> Assumptions -> {t >= 0}]; {ff[t],
>> i2 = NIntegrate[(
>> Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k, 0, \[Infinity]}]}]
>>
>> {22.907, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
>>
>> (* Numerical evaluation of ff ala Kevin *)
>> R = 8000; Z = 1; rd = 3500;
>> Clear[ff];
>> ff[t_?NumericQ] :=
>> NIntegrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}];
>> Timing[{ff[t],
>> i2 = NIntegrate[(
>> Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k, 0, \[Infinity]}]}]
>>
>> {13.484, {ff[t], 0.000424067}}
>>
>> I guess it depends a lot on the computer. Incidentally, I am running
>> 64-bit mma8 under 64-bit Vista.
>>
>> Kevin
>>
>> On 3/1/2011 5:29 AM, DrMajorBob wrote:
Kevin J. McCann
running a 3.6 GHz quad-core i5 processor with 4G of memory and my times
are quite a bit longer than yours; (2) my timings with ff done
numerically are about twice as as fast as the symbolic:
(* Symbolic evaluation of ff ala Bobby *)
R = 8000; Z = 1; rd = 3500;
Timing[Clear[ff];
ff[t_] = Integrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
Assumptions -> {t >= 0}]; {ff[t],
i2 = NIntegrate[(
Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k, 0, \[Infinity]}]}]
{22.907, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
(* Numerical evaluation of ff ala Kevin *)
R = 8000; Z = 1; rd = 3500;
Clear[ff];
ff[t_?NumericQ] :=
NIntegrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}];
Timing[{ff[t],
i2 = NIntegrate[(
Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k, 0, \[Infinity]}]}]
{13.484, {ff[t], 0.000424067}}
I guess it depends a lot on the computer. Incidentally, I am running
64-bit Mathematica 8 under 64-bit Vista.
Kevin
On 3/1/2011 5:29 AM, DrMajorBob wrote:
Fahim Chandurwala
> A couple of comments. (1) Bobby, your times are amazingly fast. I am
> running a 3.6 GHz quad-core i5 processor with 4G of memory and my times
> are quite a bit longer than yours; (2) my timings with ff done
> numerically are about twice as as fast as the symbolic:
>
> (* Symbolic evaluation of ff ala Bobby *)
> R = 8000; Z = 1; rd = 3500;
> Timing[Clear[ff];
> ff[t_] = Integrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
> Assumptions -> {t >= 0}]; {ff[t],
> Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k, 0, \[Infinity]}]}]
>
> {22.907, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
>
> (* Numerical evaluation of ff ala Kevin *)
> R = 8000; Z = 1; rd = 3500;
> Clear[ff];
> ff[t_?NumericQ] :=
> NIntegrate[Cos[t R Sin[\[Theta]]], {\[Theta], 0, \[Pi]}];
> Timing[{ff[t],
> Exp[-k Abs[Z]] ff[k])/(1 + (k rd)^2)^1.5, {k, 0, \[Infinity]}]}]
>
> {13.484, {ff[t], 0.000424067}}
>
> I guess it depends a lot on the computer. Incidentally, I am running
> 64-bit Mathematica 8 under 64-bit Vista.
>
> Kevin
>
> On 3/1/2011 5:29 AM, DrMajorBob wrote:
>
>
>
>
>
>
>
> > The inner integration is much faster done symbolically:
>
> > Timing[Clear[ff];
> > ff[t_] =
> > Integrate[Cos[t*R*Sin[\[Theta]]], {\[Theta], 0, \[Pi]},
> > Assumptions -> {t>= 0}];
> > {ff[t], i2 =
> > NIntegrate[
> > Exp[-k*Abs[Z]]/(1 + (k*rd)^2)^1.5*ff[k], {k, 0, \[Infinity]}]}
> > ]
>
> > {0.725109, {\[Pi] BesselJ[0, 8000 t], 0.000424068}}
>
> > Bobby
>
> >> (that latter if you compile to C). Evaluating in parallel should help.
> >> faster than by default but without significant degradation in quality of
> >> results. Here is a set of timings using a few different methods, and
> >> have PrecisionGoal set fairly low (three digits).
>
> >> In[109]:= Table[
> >> Timing[NIntegrate[
> >> Exp[-k Abs[Z]]/(1 + (k *rd)^2)^(3/2)* \[Pi] BesselJ[0,
> >> Abs[k*R]], {k, 0, \[Infinity]}, PrecisionGoal -> 3,
> >> Method -> meth]], {meth, {Automatic, "DoubleExponential",
> >> "Trapezoidal", "RiemannRule"}}]
>
> >> converge to prescribed accuracy after 9 recursive bisections in k near
> >> {k} = {0.0002724458978988764}.NIntegrateobtained
> >> 0.00042483953211734914` and 0.000012161444876769691` for the integral
> >> and error estimates.>>
>
> >> Out[109]= {{0.6709999999999923,
> >> 0.0004240678889181539}, {0.0150000000000432,
> >> 0.0004240644189596502}, {0.03100000000000591,
> >> 0.0004240644189596502}, {0.04699999999996862,
> >> 0.0004248395321173491}}
>
> >> I rather suspect there are more option tweaks that could make this
> >> faster still without appreciable degradation in quality of results.
>
> >> Daniel Lichtblau
> >> Wolfram Research
I am trying to numerically integrate, in Mathematica 8, the following:
NIntegrate[y Exp[-t],{y,-1,1},{t,-1,1}]. It gives me slow convergence
errors and is slow. Of course, it works fine if I just symbolically,
integrate this function. But my concern is timing. Is there a way to
have Mathematica give me the best answer it can reach in a specified
time? Setting PrecisionGoal and AccuracyGoals arent vvery helpful in
speeding up the process. Technically I am trying to do quadruple
integral about 11^4 times. Therefore speed is a huge concern, and
accuracy isn't. Especially if something is close to Zero. Thank you.
DrMajorBob
> NIntegrate[y Exp[-t],{y,-1,1},{t,-1,1}].
These are all pretty fast:
Clear[i]
Timing[
(i[y_, t_] = Integrate[y Exp[-t], y, t]);
list = i @@@ Tuples@{{-1, 1}, {-1, 1}};
list.{1, -1, -1, 1}
]
{0.001188, 0}
Timing@Integrate[y Exp[-t], {y, -1, 1}, {t, -1, 1}]
{0.042455, 0}
Timing@NIntegrate[y Exp[-t], {y, -1, 1}, {t, -1, 1}]
{0.005376, 0.}
... and here's a sanity check with different limits:
Clear[i]
yLimits = {-1/2, 1};
tLimits = {-1, 2};
Timing[
(i[y_, t_] = Integrate[y Exp[-t], y, t]);
list = i @@@ Tuples@{yLimits, tLimits};
list.{1, -1, -1, 1} // N
]
{0.001236, 0.968605}
Timing@N@Integrate[
y Exp[-t], {y, Sequence @@ yLimits}, {t, Sequence @@ tLimits}]
{0.043329, 0.968605}
Timing@NIntegrate[y Exp[-t], Evaluate@{y, Sequence @@ yLimits},
Evaluate@{t, Sequence @@ tLimits}]
{0.005648, 0.968605}
Bobby
On Sun, 01 May 2011 05:22:14 -0500, Fahim Chandurwala <fcha...@gmail.com>
wrote: