a problem in numerical integration in Mathematica

Question

I am trying to numerically compute an integral of the following type in Mathematica. I can show that this integral can be computed analytically. Still, I want to make sure Mathematica can also give me the same result using NIntegrate since I will compute some complicated integral of the same kind numerically later. Here is the integral,

$$\int_0^{3} d\omega e^{-i\omega t}\omega\Big(\int_0 ^{\infty}d\omega' \frac{e^{i\omega' t}}{\omega'-\omega-i\epsilon}-\int_0 ^{\infty}d\omega' \frac{e^{-i\omega' t}}{\omega'+\omega+i\epsilon}\Big),$$

where $0\leq\epsilon\ll 1$ after the computation of the integral. One can change the variable of integration in the second integral (inside the parenthesis) such that $\omega' \to -\omega'$ and then combine the two integrals with the result $$ \Big(\int_0 ^{\infty}d\omega' \frac{e^{i\omega' t}}{\omega'-\omega-i\epsilon}-\int_0 ^{\infty}d\omega' \frac{e^{-i\omega' t}}{\omega'+\omega+i\epsilon}\Big)=\int_{-\infty} ^{\infty}d\omega' \frac{e^{i\omega' t}}{\omega'-\omega-i\epsilon}.$$

Then, we can use complex integration to solve this integral. The answer is $2\pi i e^{i\omega t}$. Plugging this result into the original double integral, we find $$\int_0^{3} d\omega e^{-i\omega t}\omega (2\pi i e^{-i\omega t})=2\pi i\int_0^3 d\omega \omega .$$

So, we can see that this integral is $t$ independent. However, if we use the NIntegrate and compute the original integral numerically, we get a t-dependent result. I am going to make this clear with an example. We know that $$2\pi i\int_0^3 d\omega \omega =0. + 28.2743 I$$. Next, we try numerical integration. The original integral can be written as $$\int_0^{3} d\omega e^{-i\omega t}\omega\Big(\int_0 ^{\infty}d\omega' \frac{e^{i\omega' t }}{\omega'-\omega-i\epsilon}-\int_0 ^{\infty}d\omega' \frac{e^{-i\omega' t}}{\omega'+\omega+i\epsilon}\Big)$$. The first integral can be written as a principal value integral plus a Dirac delta part: $$\int_0^{3} d\omega e^{-i\omega t}\omega \Bigg(PV(\int_0 ^{\infty}d\omega' \frac{e^{i\omega' t }}{\omega'-\omega})+I\pi e^{i\omega t}-\int_0 ^{\infty}d\omega' \frac{e^{-i\omega' t}}{\omega'+\omega+i\epsilon}\Bigg)$$

We can simply set $\epsilon=0$ in the last integral since it doesn't have any singularity. I computed these integrals using NIntegrate and for t=1, the result was 9.22151244253655*^-13 + 28.274333882309378 I. For t=10, the answer is 1.83864*10^-7 + 28.2743 I. The answers are t dependent and also the real part is not equal to the correct answer. As you can see, the correct answer is a pure imaginary number, but the numerical answer has a small real part. I have included my code here. I will appreciate any advice or suggestion. Thanks!

I1[t] is the last integral (in the last equation). I have performed the integral over y first. That gave me an exponential integral.

I1[t_]:=NIntegrate[Exp[-I*x*t](E^(I t x) ExpIntegralEi[-I t(x+Infinity)]-E^(I t x) ExpIntegralEi[-I t (x + 0)])*(x), {x, 0, 3}] 

The contribution of the second term ($i\pi e^{i\omega t}$) to the integral has been computed by

I2 = NIntegrate[I*Pi*x, {x, 0, 3}]

This part is t independent and is equal to 0. + 14.1372 I. The principal value integral is computed by

I3[ts_]:=NIntegrate[NIntegrate[Exp[-I*x*ts]*(Exp[I*y*ts]/(y - x))*(x), {y, 0, x, 3},PrincipalValue -> True], {x, 0, 3}]

Since we $3<y<\infty$, there is no singularly in the integral so there is no need to use the principal value. For that piece, I used

I4[ts_]:=NIntegrate[Exp[-i*x*t] (x)*(E^(I t x) ExpIntegralEi[I t (-x + Infinity)] -E^(I t x) ExpIntegralEi[I t (-x + 3)]), {x, 0, 3}].

Answer 1

Here's a straightforward way:

i10[t_, w_?NumericQ] := NIntegrate[
   Exp[I v t]/(v - w), {v, -Infinity, w, Infinity},
   Method -> "PrincipalValue",
   WorkingPrecision -> MachinePrecision];
i1[t_] := NIntegrate[E^(-I t w) w*i10[t, w], {w, 0, 3},
  WorkingPrecision -> MachinePrecision]

i1[1/2]

(*  1.69022*10^-8 + 14.1372 I  *)

Notes: (1) I think your 28.27 is off by a factor of 2. (2) I don't think ExpIntegralEi represents your integral.

---

You can check the interior integral this way:

Assuming[t > 0 && w > 0,
 Sqrt[2 Pi] FourierTransform[1/(v - w), v, t]
 ]

(*  I E^(I t w) \[Pi]  *)

i2[t_] := 
 NIntegrate[E^(-I t w) w*(I E^(I t w) \[Pi]), {w, 0, 3}, 
  WorkingPrecision -> MachinePrecision]

i2[1/2]

(*  0. + 14.1372 I  *)

And i2[t] may be computed exactly, of course.