15
$\begingroup$

I am trying to evaluate a highly oscillatory integral using NIntegrate. I fear that due to limited resources (time and/or memory), I will not be able to evaluate the integral to the desired precision. Thus, I would like to programmatically access the error estimates that are e.g. reported by the messages NIntegrate::maxp, NIntegrate::ncvb, or NIntegrate::eincr. I could not find an option of NIntegrate that would directly make these error estimates available. However, given that I have to evaluate a multitude of integrals, it is impractical to obtain the errors from the warnings by hand.

The following example generates the NIntegrate::maxp message (obviously this very integral has an analytical solution):

NIntegrate[Sin[x]/Sqrt[x], {x, 0, 100}, Method -> "MonteCarlo",PrecisionGoal -> 6]

NIntegrate::maxp: The integral failed to converge after 50100 integrand evaluations. NIntegrate obtained 1.1787733508261242and 0.07678430788995934 for the integral and error estimates.

How to get (if necessary, extract) the error estimate (0.07678430788995934`)?

Remark: The example from the help of NIntegrate::eincr, i.e. ref/message/NIntegrate/eincr, does not produce the expected message in version 8.0; unfortunate my integrals still do.

Improve this question
$\endgroup$
2
  • 2
    $\begingroup$ "...evaluate a highly oscillatory integral using NIntegrate[]." - then, why "MonteCarlo"? There's "DoubleExponential" or "ClenshawCurtisOscillatoryRule" which you could have used... unless your actual integrals are in fact multidimensional, and you've just grossly oversimplified. $\endgroup$
    – J. M.'s eventual burnout ♦
    Nov 13, 2012 at 11:39
  • 2
    $\begingroup$ Yes. I am grossly simplifying and the actual integral is multidimensional (4D). In particular, I have just chosen the method since it generates one of the messages in question. It turned out that generating the NIntegrate::eincr message with a simple 1D integral is unexpectedly (given my mathematical naivety) difficult, i.e. NIntegrate is very robust (see my remark). For the actual integral, Method->”MonteCarlo” is in fact my best bet. $\endgroup$
    – dan
    Nov 13, 2012 at 12:54

2 Answers 2

Reset to default
14
$\begingroup$

It seems to me that there's a better approach, but one way is to define your own DownValue for this particular message. For example:

Unprotect[Message];
Message[NIntegrate::maxp, its_, int_, err_] := Sow[err]

Then

NIntegrate[Sin[x]/Sqrt[x], {x, 0, 100}, 
  Method -> "MonteCarlo", PrecisionGoal -> 6] // Reap

(* Out: {1.07721, {{0.0761274}}} *)
Improve this answer
$\endgroup$
2
  • 1
    $\begingroup$ Excellent. For NIntegrate::eincr the same approach applies. For NIntegrate::ncvb, Message[NIntegrate::ncvb, nr_, var_, varlist_, at_, int_, err_] := Sow[err] seems to do the job. In general, Message[NIntegrate::ncvb | NIntegrate::eincr | NIntegrate::maxp, l___] := Sow[Last@{l}] could probably be used. $\endgroup$
    – dan
    Nov 13, 2012 at 13:25
  • 2
    $\begingroup$ In order to retain the original message, one could probably use something along the lines of: Block[{original = False}, Unprotect[Message]; Message[NIntegrate::maxp, l___] /; Not[original] := (Sow[Last@{l}]; original = True; Message[NIntegrate::maxp, l]); NIntegrate[Sin[x]/Sqrt[x], {x, 0, 100}, Method -> "MonteCarlo", PrecisionGoal -> 6] // Reap ] $\endgroup$
    – dan
    Nov 13, 2012 at 13:33
6
$\begingroup$

From the description of the question it seems to me that using the (undocumented) option IntegrationMonitor to obtain integration intervals and estimates might be very useful.

Here is an example:

t = Reap[NIntegrate[Sin[x]/Sqrt[x], {x, 0, 100}, PrecisionGoal -> 6, 
    Method -> "MonteCarlo", 
    IntegrationMonitor -> (Sow[
        Map[{#1@"Boundaries", #1@"Integral", #1@"Error"} &, #1]] &)]];
res = t[[1]];
t = t[[2, 1]];
Take[t, -4]

enter image description here

More examples and explanations about the use of IntegrationMonitor can be found in the notebook "Finding the applied NIntegrate methods.nb" attached to the community.wolfram.com discussion "Integration method used in NIntegrate".

Improve this answer
$\endgroup$
5
  • $\begingroup$ How could one out of this example get the total error of the numerical integral (and not for each subdivision)? $\endgroup$
    – hal
    Dec 23, 2020 at 11:52
  • 1
    $\begingroup$ @hal Just replace Map with Total@Map. Please post a new MSE question. $\endgroup$
    – Anton Antonov
    Dec 23, 2020 at 13:54
  • $\begingroup$ I will. Thank you! $\endgroup$
    – hal
    Dec 23, 2020 at 22:42
  • $\begingroup$ @AntonAntonov What is the meaning of "Boundaries" in this example? Thanks! $\endgroup$
    – Ulrich Neumann
    Apr 30, 2021 at 6:51
  • $\begingroup$ @UlrichNeumann "Boundaries" is for the boundaries of the integration region of the integration region object. Each integration region object has its own integration region, integration function, singularity handler, and integral and error estimates. $\endgroup$
    – Anton Antonov
    Apr 30, 2021 at 12:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.