0
$\begingroup$

I want to integrate the function

Log [ (x^2 + t^2 + m^2 + 50^2)/( x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))]

for x,t and y between -50 and 50. This function jumps suddenly at x=-y, so the code

NIntegrate[SetPrecision[Log [ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-(x + y)^2/(2 0.00001)], 10], {x, -50, 50}, {t, -50, 50}, {y, -50, 50}, Method -> {"GlobalAdaptive", "SingularityDepth" -> 8, "SingularityHandler" -> "IMT"}, MinRecursion -> 30, MaxRecursion -> 40, WorkingPrecision -> 10 ]

gives the following output:

NIntegrate::errprec: Catastrophic loss of precision in the global error estimate due to insufficient WorkingPrecision or divergent integral.

Then, i thought about partitioning the region, which is a square, in 6 sub-regions, as follows:

NIntegrate[ SetPrecision[(50^2)/Log[ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, -50, 0}, {t, -50, 50}, {y, -50, 0}, Method -> {"GlobalAdaptive", "SingularityDepth" -> 8, "SingularityHandler" -> "IMT"}, MinRecursion -> 30, MaxRecursion -> 40, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log[ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, -50, 0}, {t, -50, 50}, {y, 0, -x}, Method -> {"GlobalAdaptive", "SingularityDepth" -> 8, "SingularityHandler" -> "IMT"}, MinRecursion -> 30, MaxRecursion -> 40, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log[ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, -50, 0}, {t, -50, 50}, {y, -x, 50}, Method -> {"GlobalAdaptive", "SingularityDepth" -> 8, "SingularityHandler" -> "IMT"}, MinRecursion -> 30, MaxRecursion -> 40, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log[ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, 0, 50}, {t, -50, 50}, {y, 0, 50}, Method -> {"GlobalAdaptive", "SingularityDepth" -> 8, "SingularityHandler" -> "IMT"}, MinRecursion -> 30, MaxRecursion -> 40, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log[ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, 0, 50}, {t, -50, 50}, {y, -x, 0}, Method -> {"GlobalAdaptive", "SingularityDepth" -> 8, "SingularityHandler" -> "IMT"}, MinRecursion -> 30, MaxRecursion -> 40, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log[ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, 0, 50}, {t, -50, 50}, {y, -50, -x}, Method -> {"GlobalAdaptive", "SingularityDepth" -> 8, "SingularityHandler" -> "IMT"}, MinRecursion -> 30, MaxRecursion -> 40, AccuracyGoal -> 10, WorkingPrecision -> 10]

Mathematica is able to compute this integral and gives the answer 129869.8932. But, I decided to compute the same thing using other code:

NIntegrate[SetPrecision[(50^2)/Log[(x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, -50, 0}, {t, -50, 50}, {y, -50, 0}, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log [ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, -50, 0}, {t, -50, 50}, {y, 0, -x}, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log[(x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, -50, 0}, {t, -50, 50}, {y, -x, 50}, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log[ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, 0, 50}, {t, -50, 50}, {y, 0, 50}, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log [(x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, 0, 50}, {t, -50, 50}, {y, -x, 0}, AccuracyGoal -> 10, WorkingPrecision -> 10] + NIntegrate[SetPrecision[(50^2)/Log[ (x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-((x + y)^2/(2 0.00001))], 10], {x, 0, 50}, {t, -50, 50}, {y, -50, -x}, AccuracyGoal -> 10, WorkingPrecision -> 10]

and then the answer is 2.738034446*10^-18! Mathematica does not warn be about any issues in the computation, the calculations are done in seconds, so it seems like everything is ok, but the results differ by more than 20 orders of magnitude! Anyone knows what is happening? And how should I deal with the singularity of my integral?

Edit: I must be able to compute the integral of the function Log [(y^2 + t^2 + m^2 + 50^2)/(y^2 + t^2 + m^2)]/ Log [(x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-(x + y)^2/( 2 0.00001)] as well, and the trick suggested in the answers does not work in this case.

$\endgroup$
0
2
$\begingroup$

Your function doesn't "jump" at x=-y, there is only a "continous peak"

Looking at the integrand I observe that you could separate the integration concerning y and t:

eps=1/100000;
inty = Integrate[Exp[-((x + y)^2/(2 eps))], {y, -50, 50}];
(* function of x and m *)

intt = Simplify[Integrate[Log[(x^2 + t^2 + m^2 + 50^2)] , {t, -50, 50}] - 
Integrate[Log[(x^2 + t^2 + m^2)] , {t,  -50, 50}], 
{Element[{x, m},Reals], x^2 + m^2 > 0, -50 <= x <= 50}]  
(* function of x and m*)

Now you can evaluate the remaining integral

int[m_?NumericQ] := NIntegrate[intt inty , {x, -50, 50}]

which depends on m:

Plot[int[m], {m, -10, 10}, PlotRange -> {0, Automatic}]

enter image description here

$\endgroup$
7
  • $\begingroup$ Thanks for your answer, it was very helpful! Could you please explain to me what "m_?NumericQ" means? Say I wanted to integrate from -l to l, instead of -50 to 50, so that the final result would depend on l as well. Then intt and inty would depend on m and l, and int[m_?NumericQ, l_?NumericQ] := NIntegrate[intt inty , {x, -l, l}] would give me a function of m and l? $\endgroup$
    – Fisher
    Jun 15 '18 at 13:31
  • $\begingroup$ Another question: I know that inty is a function of m and x. Say I want to convert it to a function of m and y, so that I can use it later without having to evaluate it again. Is that possible? $\endgroup$
    – Fisher
    Jun 15 '18 at 13:44
  • $\begingroup$ m_?NumericQ restricts the function argument to be numerical.That's why NIntegrate is only called, if this condition is fullfilled. I'm not sure, wether I got your second question. x,y , are independent integration variables. I think it's not possible to convert x to y . $\endgroup$ Jun 15 '18 at 14:07
  • 1
    $\begingroup$ Btw, I need to compute a similar integral, but the function now is Log [(y^2 + t^2 + m^2 + 50^2)/(y^2 + t^2 + m^2)]/ Log [(x^2 + t^2 + m^2 + 50^2)/(x^2 + t^2 + m^2)] Exp[-(x + y)^2/( 2 0.00001)], also for all variables ranging from -50 to 50. But I dont think that I can use your trick this time... Should I ask another question or edit this one? $\endgroup$
    – Fisher
    Jun 15 '18 at 14:47
  • $\begingroup$ In this case only the integral over t is separable...I'm not sure wether the x,y integral converges. $\endgroup$ Jun 15 '18 at 15:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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