Cauchy integral theorem vs. integral formula

Let's have the following integral:

Here y and η are real and positive number. We can solve this integral using Cachy integral theorem by converting this integeral to a complex integral. Let's assume a contour over upper half of the complex plane. Then, we will have:

this integral have two simple poles at and . The residues at both poles are and , respectively. Hence, the total integral will be:

Now, I want to check that if this integral is actually zero, using numerical integration. For that I write the following code:

clear; clc;

a = -1000000; %close to -infinity
b = 1000000; %close to +infinity

n = rand(1); %\eta: some randome number
y = rand(1)*10; %some randome number

f = @(x, n, y) 1 ./ ((x - 1i*n) .* (x - y - 1i*n));


x = linspace(a,b,1e5); dx =x(2)-x(1);

yy = f(x,n,y);

%plot(x,yy)

sum = sum(yy)*dx

this gives me number for one set of random n and y.

1.3303 + 0.2697i

which is not near to zero. what could be the reason?

Answer 1

[Paul]

Hi Luqman,

The numerical integration is using a step size that is proably too large to capture the rapid variation in the integrand near the poles. And the limits of the the numerical integration are going way out into the tails of the integrand that are not really contributing mugh to the integral.

Original code

clear; clc;
a = -1000000; %close to -infinity
b = 1000000; %close to +infinity
rng(101)
n = rand(1); %\eta: some randome number
y = rand(1)*10; %some randome number
f = @(x, n, y) 1 ./ ((x - 1i*n) .* (x - y - 1i*n));
x = linspace(a,b,1e5); dx =x(2)-x(1);
dx  % large relative to 1/x

dx = 20.0002

yy = f(x,n,y);
yy([1 end])

ans = 
   1.0e-11 *

   0.1000 - 0.0000i   0.1000 + 0.0000i

%plot(x,yy)
% don't override names of built-in functions
thesum = sum(yy)*dx

thesum = 0.6772 + 0.0686i

Tighten the limits, which also makes the dx smaller

x = linspace(-30,60,1e5); dx =x(2)-x(1);
dx

dx = 9.0001e-04

yy = f(x,n,y);
yy([1 end])

ans = 
   1.0e-03 *

   0.9328 - 0.0296i   0.3069 + 0.0056i

%plot(x,yy)
thesum = sum(yy)*dx

thesum = -0.0480 + 0.0003i

Also, we can use integral

g = @(x) f(x,n,y);
integral(g,-inf,inf)

ans = -1.6653e-16 + 3.3307e-16i

[Luqman Saleem]

Thank you.