# How can I code a contour integral representation of cosine?

7 visualizzazioni (ultimi 30 giorni)
Richard Mc Enteggart il 23 Ott 2022
Modificato: David Goodmanson il 25 Ott 2022
% "Contour Integral Representation of Cosine"
%==========================%
% making 'z' a possible number such as pi
N = 10000000; % example number
z_lower = 0;
z_upper = 6*pi;
%==========================%
z1 = linspace(z_lower,z_upper,N);
y = 1;
Sl = y - 10*1i;
Su = y + 10*1i;
ds = (Su - Sl)/N; % size of each step
sum = 0.0;
%==========================%
for z = linspace(z_lower,z_upper,N)
for Sn = linspace(Sl,Su,N)
sum = sum + ((exp(Sn) - (z.^2/4*Sn))/sqrt(Sn))*ds;
end
end
sum = sum*(sqrt(pi)/2*pi*1i);
plot(Sn,sum)
The above is my poor attempt! I am trying to write a code to numerically approximate the integral along a suitable contour and plot the approximation against the exact value of cos z. I would like to do this twice - once for z ∈ [0, 6π] and once for complex valued z in the range z ∈ [0 + i, 6π + i]. Then I intend to plot both the real and imaginary part of the computed cos z.
I have a following lost of steps that I am trying to implement!
• Choose γ, Sl, Su , N.
• Step through z from z lower to z upper (using a different number of steps (other than N) for this).
• For each value of z compute cos z by stepping along the contour in N discrete steps from Sl to Su .
• For each value of sn along the contour compute the integrand "((exp(Sn) - (z.^2/4*Sn))/sqrt(Sn))" [I have attached an image to make this integrand clearer to read] and add it to the rolling sum.
• Plot the result.
• Repeat using a better integration rule
Image of the integrand!
I am pretty new to general forms of computation and would really appreciate any help! Thanks!
##### 8 CommentiMostra 6 commenti meno recentiNascondi 6 commenti meno recenti
Torsten il 23 Ott 2022
I don't know. As you can see, it runs under MATLAB R2022b. Which MATLAB version do you have in use ?
Richard Mc Enteggart il 24 Ott 2022
Modificato: Richard Mc Enteggart il 24 Ott 2022
I use MATLAB R2022b as well! I ran the code suggested by you and it works! My machine just did not compute it I guess haha
Thank you very much :D

Accedi per commentare.

### Risposte (1)

David Goodmanson il 24 Ott 2022
Modificato: David Goodmanson il 24 Ott 2022
Hi Richard,
this integral converges very slowly along the vertical path gamma-i*inf to gamma+i*inf.
Convergence depends on what happens for large values of s. In that case, you can drop the z^2/s term in the exponential and look at exp(s). When the path of integration has large imaginary part, exp(s) is oscillatory and only 1/sqrt(s) is making the integrand smaller. Convergence is slow.
Rather than stick to the vertical path, which will prove the point philosophically but is difficult to achieve numerically, it helps to use a different one. The path is moved over from the original path in a continuous fashion. Here the new path is a rectangle:
-inf+ia a+ib
------<-------
o ^ o = origin
------>------^
-inf-ia a-ib
The idea is that at the left hand parts of the path, s has a large negative real part and exp(s) dies off quickly rather than oscillating.
In the code below it would have been nice to use the set
I1 = integral(@(t) fun(t,z), -inf-i*a, b-i*a);
I2 = integral(@(t) fun(t,z), b-i*a, b+i*a);
I3 = integral(@(t) fun(t,z), b+i*a, -inf+i*a);
but the integral function can't handle complex end points at infinity, so the vertical offset is inserted manually in I1 and I3. With the choices for a and b the code does well up to z = 7*pi. You can mess around with a and b to see what works.
z = pi/3;
a = 3
b = 1
I1 = integral(@(t) fun(t,z,-i*a), -inf, b);
I2 = integral(@(t) fun(t,z,0 ), b-i*a, b+i*a);
I3 = integral(@(t) fun(t,z, i*a), b, -inf);
I = sqrt(pi)/(2*pi*i)*(I1+I2+I3)
delta = I - cos(z) % check
function y = fun(t,z,u)
s = t+u;
y = (exp(s-z.^2./(4*s))./sqrt(s));
end
I = 0.5000 + 0.0000i
delta = -2.7756e-16 + 6.2638e-17i
##### 2 CommentiMostra NessunoNascondi Nessuno
Richard Mc Enteggart il 24 Ott 2022
Thanks for getting back to me! That is incredibly interesting! Your expalantions are making sense. I did not know my integrand was getting significantly smaller.
Is there any way I could plot this bloc of code?
David Goodmanson il 24 Ott 2022
Modificato: David Goodmanson il 25 Ott 2022
Hi Richard, I expanded the original answer to hopefully provide some clarification. As to plotting, if you use the 'arrayvalued' option in the integral function then z can be a vector of values and you can plot the resulting output vs z.

Accedi per commentare.

### Categorie

Scopri di più su Numerical Integration and Differentiation in Help Center e File Exchange

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by