Numerically integrate function with sudden jumps

 Numerically integrate function with sudden jumps

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?

Let’s block ads! (Why?)

Recent Questions – Mathematica Stack Exchange