Integrate and NIntegrate

An Oscillatory Example

Consider the integral of a "bad" function

Plot[Sin[1/x], {x, 0, Pi/2}, AspectRatio -> 1/6];

Plot[Sin[1/x], {x, 0, 0.01}, AspectRatio -> 1/6];

Integrate[Sin[1/x], {x, 0, Pi/2}]

N[%]

NIntegrate[Sin[1/x], {x, 0, Pi/2}] // Timing

NIntegrate[Sin[1/x], {x, 0, Pi/2}, MaxRecursion -> 12] // 
      Timing

NIntegrate[Sin[1/x], {x, 0, Pi/2}, Method -> DoubleExponential, 
        MaxRecursion -> 12] // Timing

N[Integrate[Sin[1/x], {x, 0, 10^(10)}]]

N[Integrate[Sin[1/x], {x, 0, 10^(100)}]]

Integrate[Sin[1/x], {x, 0, Infinity}]

A Localized Example

Consider now the integral of a "good" function

Plot[Exp[-x^2], {x, -10, 10}, 
        PlotRange -> {0, 1}];

Integrate[Exp[-x^2], {x, -3, 3}]

N[%]

Integrate[Exp[-x^2], {x, -Infinity, Infinity}]

N[%]

Pretty close, isn't it?

But, suppose we were trying to integrate this numerically.

NIntegrate[Exp[-x^2], {x, -3, 3}]

Again, pretty close! So, how about improving the precision by extending the limits:

NIntegrate[Exp[-x^2], {x, -10000, 10000}]

What gives !?!

NIntegrate uses recusrion to "zero in" on regions of interest...

NIntegrate[Exp[-x^2], {x, -10000, 10000}, 
      MaxRecursion -> 1000]

Nope. That's not it.

NIntegrate[Exp[-x^2], {x, -10000, 10000}, 
      Method -> Oscillatory]

Since it is not sampling the function sufficiently close to zero, let's help it:

NIntegrate[
      Exp[-x^2], {x, -100000, -1, 0, 1, 100000}]

Quite a bit better...

NIntegrate[
      Exp[-x^2], {x, -100000, -3, 0, 3, 100000}]

Back to "pretty good"...

NIntegrate[
      Exp[-x^2], {x, -100000, -5, 0, 5, 100000}]

Finally.

Integration Path

Having subdivided the path of an integral,

Integrate[1, {x, -5, 5}]

can we integrate back-and-forth?

Integrate[1, {x, -5, 0, -1, 5}]

How about in the complex plane?

Integrate[1, {x, -5 + I, 5 - I}]

Indeed, the integration contour is a slanted straight line, extending from -5 to +5 in the real direction, and +i to -i in the imaginary direction. These two don't add, but are summed separately. This integral is not measuring the length of the contour in the usual, Euclidean metric, where the length would be \(TraditionalForm\`ds = \@\(\((d\ x)\)\^2 + \((d\ y)\)\^2\)\), but in the "taxicab metric" (also called the Manhattan metric), where you add up the number of blocks you need to drive North-South plus the number of blocks East-West---and backtracking doesn't count.

But, since we can do contour integrations, how about a closed-contour integral about the origin:

Integrate[1/z, {z, 1, I, -1, -I, 1}]

Neat! Mathematica knows residues!

Integrate[1/z^2, {z, 1, I, -1, -I, 1}]

Integrate[1/z^3, {z, 1, I, -1, -I, 1}]

That's because 1/z has a simple pole at z = 0, and hence a residue (= the coefficient of the (z^(-1))term in the Laurent expansion ( = power expansion over all integral powers, from -∞ to +∞ )) is 1; 1/(z^m), for m > 1 has a pole of order m, and so the coefficient of the  (z^(-1))term is, of course, zero.

How about multiply-winding contours?

Integrate[
      1/z, {z, 1, I, -1, -I, 2, 2  I, -2, -2 I, 
        1}]

Not bad for a 'puter!

So, we can define a function calculating residues, except Mathematica already has one:

Residue[1/z, {z, 0}]

And it's pretty good:

{Residue[Sin[1/z], {z, 1}], Residue[Exp[1/z], {z, 1}]}

But has problems with essential singularities:

{Residue[Sin[1/z], {z, 0}], Residue[Exp[1/z], {z, 0}]}

Speaking of poles and residues:

Integrate[1/(x^3 - 1), {x, -Infinity, Infinity}]

?? Integrate

Integrate[1/(x^3 - 1), {x, -Infinity, Infinity}, 
        PrincipalValue -> True] // Timing

Integrate[Exp[I x]/(x^3 - 1), {x, -Infinity, Infinity}, 
        PrincipalValue -> True] // Timing

Pretty good, even if you have to wait a bit over two minutes.

Finally, you can actually specify the contour as a parametric curve, specified in the complex plane as z = (e^(i t))(
          2  cos(t) + 1):

ParametricPlot[{Re[Exp[I t] (2 Cos[t] + 1)], 
          Im[Exp[I t] (2 Cos[t] + 1)]}, {t, 0, 2  Pi}];

curve[t_] := Exp[I t] (2 Cos[t] + 1)

Integrate[1/(curve[t] - 1/2), {t, 0, 2  Pi}]

Now, this is not right...

...since

Residue[1/(z - 1/2), {z, 1/2}]

and the the limaçon goes around the residue twice. A-ha.

Integrate[1/(curve[t] - 1/2), {t, 0, Pi}]

ParametricPlot[{Re[Exp[I t] (2 Cos[t] + 1)], 
          Im[Exp[I t] (2 Cos[t] + 1)]}, {t, 0, Pi}];

Integrate[1/(curve[t] - 1/2), {t, 0, Pi/2}]

ParametricPlot[{Re[Exp[I t] (2 Cos[t] + 1)], 
          Im[Exp[I t] (2 Cos[t] + 1)]}, {t, 0, Pi/2}];

Oh, well... forget about that. Let's do this numerically.

Numerical Contour Integration

ContourIntegral[f_, par : (z_ -> g_), {t_, a_, b_}] := 
      NIntegrate[Evaluate[(f /. par) D[g, t]], {t, a, b}]

ContourIntegral[1/(z - 1/2), 
      z -> curve[t], {t, 0, 2 Pi}]

Chop[%]

N[%/(4 Pi)]