Hi SB
You can use the principal value for the integration,which is defined as
(1/2) * [ (integral with path on one side of the pole) + (integral with path on other side of the pole)]
You have
pole 3 pole 2 pole 1
-4 -1 +1
In this case there are two poles that get crossed, so there are 4 paths to consider
pole 2 pole 1 case for example:
path passes to L L a path of radius 1 shifted to left by 1/2
L R b path of radius 2
R L c path of radius .8
R R d path of radius 1 shifted to right by 1/2
and the answer is 1/2 x 1/2 = 1/4 of the sum of the integrals along these paths. By explicit integration on those paths:
f = @(z) (z.^2+3*z+2i)./((z+4).*(z-1).*(z+1));
% case a
Ra = 1; shift = -1/2;
z = @(theta) Ra*exp(i*theta) + shift;
integranda = @(theta) f(z(theta)).*(Ra*exp(i*theta))*i;
Ia = integral(integranda,0,2*pi);
% case b
Rb = 2;
z = @(theta) Rb*exp(i*theta);
integrandb = @(theta) f(z(theta)).*(Rb*exp(i*theta))*i;
Ib = integral(integrandb,0,2*pi);
% case c
Rc = .8; % this integral will be zero since no poles are enclosed
z = @(theta) Rc*exp(i*theta);
integrandc = @(theta) f(z(theta)).*(Rc*exp(i*theta))*i;
Ic = integral(integrandc,0,2*pi);
% case d
Rd = 1; shift = 1/2
z = @(theta) Rd*exp(i*theta) + shift;
integrandd = @(theta) f(z(theta)).*(Rd*exp(i*theta))*i;
Id = integral(integrandd,0,2*pi);
I = (Ia+Ib+Ic+Id)/4
I = 0.4189 + 2.3038i
which agrees with the result from residue theory, no explicit integrartion.
A simpler way, still with numerical integration, is to draw two closely spaced vertical lines up the y axis, which splits the unit circle into a D and a backwards D. The ccw path around the unit circle is the same as the sum of ccw paths around each of the two D's since the integrals along the straight sides of the D's are in opposite directions and cancel each other. For the right hand D, the path going to just the left of the pole results in 0 since the new D encloses no poles. The path going just to the right of the pole does enclose the pole. You can alter the path to be a small circle enclosing the pole and keep the factor of (1/2) from the principal value definition. The backwards D works similarly which results in
% case a
Ra = 1/3; % for example
z = @(theta) Ra*exp(i*theta) + 1;
integranda = @(theta) f(z(theta)).*(Ra*exp(i*theta))*i;
Ia = (1/2)*integral(integranda,0,2*pi);
% case b
Rb = 1/3; % for example
z = @(theta) Rb*exp(i*theta) - 1;
integrandb = @(theta) f(z(theta)).*(Rb*exp(i*theta))*i;
Ib = (1/2)*integral(integrandb,0,2*pi);
I = Ia+Ib
and gives the same result.
