# NIntegrate of a highly oscillatory integral double exponential oscillatory

I checked NIntegrate integration strategies and found that my integral is a double exponential oscillatory function. I tried increasing the precision and recursions and specifying the method but I still get the following error

NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small. >>


Here is the code that i used.

g[σ_, center_, q_] := q^5 + q^4 + q;
f[σ_, center_, energy_, mass_, q_] :=
E^(I (g[σ, center, q])) E^(-((q - center)^2/(4 σ^2))) E^(I energy q)
val[σ_, center_, energy_, mass_, q_] :=
NIntegrate[
f[σ, center, energy, mass, q], {q, -∞, ∞}, MaxRecursion -> 20,
WorkingPrecision -> 50, Method -> "DoubleExponentialOscillatory"];
val[3, 10, 500, 1, q]


I'm basically having a problem on how to deal with the highly oscillatory behavior of the integrand such that I get a nice numerical result.

• Isn't the integral of both the real and imaginary components 0 by symmetry? – Feyre Aug 27 '16 at 11:26
• i dont think its zero – PhilCsar Aug 27 '16 at 11:30
• Any intuition on how large it should be? – Michael E2 Aug 27 '16 at 12:21
• i really have no idea on what the value should be – PhilCsar Aug 27 '16 at 12:23
• How about an order of magnitude? – Feyre Aug 27 '16 at 12:54

I see the problem from another point of view. Basically, the integral is a Gaussian multiplying a highly oscillatory function (two functions indeed, a Sin and a Cos). First, I made a variable change to center the function at $q=0$, thus, for the example:

f[3, 10, 500, 1, q]
(*Out[1]= *E^(-(1/36) (-10 + q)^2 + 500 I q + I (q + q^4 + q^5))*)
E^(-(1/36) (-10 + q)^2)*ExpToTrig[E^(+500 I q + I (q + q^4 + q^5))] /. {(q - 10) -> x, q -> (10 + x)} // ExpandAll


If we plot the real part, these are two functions, the exponential and the cosine:

GraphicsRow[Plot[#[[1]], {x, -15, 15},PlotRange -> #[[2]]] & /@ {{E^(-(x^2/36)) , {-15, 15}}, {Cos[115010 + 54501 x + 10600 x^2 + 1040 x^3 + 51 x^4 + x^5], {-0.001, 0.001}}}]


We can see that the integrand is not odd. The same happens to the imaginary part.

Thus, I look for another strategy. As we have a moderate rapidly decreasing high oscillatory function, I have calculated the integral value for the sequence in which the integration interval increases (it takes sometime):

Chop@Table[NIntegrate[E^(-(x^2/
36)) (Cos[115010 + 54501 x + 10600 x^2 + 1040 x^3 + 51 x^4 + x^5] +
I Sin[115010 + 54501 x + 10600 x^2 + 1040 x^3 + 51 x^4 +
x^5]) , {x, -a, a}, Method -> {"GlobalAdaptive", "MaxErrorIncreases" -> 1000}, MaxRecursion -> 20, WorkingPrecision -> 25], {a, 1, 40, 4}]
(*{0.00003599894018448947571489283, 0.00001013416214467363360152173,-0.00007020239541982327974633850, -5.291033137885836075741994*10^-7, -1.046793294458294734393493*10^-8, 0, 0, 0, 0, 0}*)


and plot the points:

Therefore, as the sequence converges, I conclude that the integral is zero.

I realise that my post is basic and lacks formalism as well as a deep knowledge of the MMA numerical integration techniques. However, I hope this helps.

I thought that you might consider it as a Fourier-Transform (from q- to energy-space) and try something along the lines of

I1[σ_, center_, q_] :=
Exp[-((-center + q)^2/(4 σ^2))] Cos[(q + q^4 + q^5)]
I2[σ_, center_, q_] :=
I Exp[-((-center + q)^2/(4 σ^2))] Sin[(q + q^4 + q^5)]
NFourierTransform[I1[3, 10, q], q, 500] + NFourierTransform[I2[3, 10, q], q, 500]
(* Out -3.15636*10^-13 + 2.38794*10^-13 I *)


(Be careful, I have not fixed any conventions and don't know how well NFourierTransform works)

• Which means 0, using NIntegrate[..,MinRecursion -> 20, MaxRecursion -> 200, WorkingPrecision -> 50] also yields 0 (3.5354832 ...*^-23 - 2.3618059 ...*^-24 I). The Integral is 0. – Feyre Aug 27 '16 at 11:36
• think so as well.... but not 100% sure – Armin Aug 27 '16 at 11:38
• @PhilCsar actually... if you consider a smaller energy, the integral is non-zero (numerically). Just fishing in the dust here: are you trying to do something with electron wave packages in natural units? if so, you may have to check your unit conversion... in other words: what is the origin of your computation? – Armin Aug 27 '16 at 12:48
• im pretty sure that theres nothing wrong in the unit conversion. also, what I posted is something that just looks similar to my original problem and I cant post my original problem here for particular reasons. – PhilCsar Aug 27 '16 at 12:57
• @PhilCsar In that case, does MinRecursion -> 20, MaxRecursion -> 200, WorkingPrecision -> 50 not help? – Feyre Aug 27 '16 at 13:11