mathematica: evaluation order in numerical integration

0
votes

Hi i have a problem concerning a numerical integration in mathematica. Here is my test function

Table[NIntegrate[
Boole[Sqrt[1 - cosk^2]*Sqrt[1 - cosk2^2] > Abs[a - cosk*cosk2]]/
Sqrt[(1 - cosk^2)*(1 - cosk2^2) - (a - cosk*cosk2)^2],
{cosk, -1,1}, {cosk2, -1, 1}, Method -> "GlobalAdaptive"], {a, -.9, .9, .1}]

The integration yields complex values although due to the boolean function in the integrand the argument in the sqrt should be always positive and hence only result in real values. Is it possible to first evaluate the boole function and only if it is true then start to numerically integrate?

If i calculate the same integral using a monte carlo integration strategy

Table[NIntegrate[
Boole[Sqrt[1 - cosk^2]*Sqrt[1 - cosk2^2] > Abs[a - cosk*cosk2]]/
Sqrt[(1 - cosk^2)*(1 - cosk2^2) - (a - cosk*cosk2)^2], {cosk, -1, 
1}, {cosk2, -1, 1}, Method -> {"MonteCarlo", "MaxPoints" -> 10^8, 
"SymbolicProcessing" -> None}], {a, -.9, .9, .1}]

how can I find out if it sums up a lot of zeroes due to the boolean function? I think the evaluation can save a lot of computation time if it first evaluates the boolean function for each sample point of the monte carlo grid. If i replace "MonteCarlo" by "AdaptiveMonteCarlo" the result goes totally wrong.

2
mathwolfram-mathematicanumerical-integration

2 Answers

0
votes

What happens if you convert to a true guard? (The Boole form may be too clever by half (Mathematica might not treat false --> 0 as a real guard).)

Table[NIntegrate[
If[Sqrt[1 - cosk^2]*Sqrt[1 - cosk2^2] > Abs[a - cosk*cosk2],
1/Sqrt[(1 - cosk^2)*(1 - cosk2^2) - (a - cosk*cosk2)^2], 0],
{cosk, -1,1}, {cosk2, -1, 1}], {a, -.9, .9, .1}]
0
votes

here is how to use EvaluationMonitor to learn what the sampling methods do:

out = Reap[With[{a = -.9},
   NIntegrate[
    v = Boole[
       Sqrt[1 - cosk^2]*Sqrt[1 - cosk2^2] > Abs[a - cosk*cosk2]]/
      Sqrt[(1 - cosk^2)*(1 - cosk2^2) - (a - 
           cosk*cosk2)^2], {cosk, -1, 1}, {cosk2, -1, 1}, 
    Method -> "LocalAdaptive", 
    EvaluationMonitor :> Sow[{cosk, cosk2, v}]]]]

Note LocalAdaptive is very slow, but is likely the most accurate purely numerical result: 

     Graphics@Point[out[[2, 1]][[All, {1, 2}]]]

enter image description here

here is global adaptive.

enter image description here

MonteCarlo just samples everywhere without regard to the integrand value.

Your best bet probably is to use global adaptive and set MinRecursion. Satisfy yourself the imaginary parts are negligible and take the real part of the result. Of course it would be far better to analytically solve for your integration limits and eliminate the Boole all together.