how to get the correct plot of function sin(N+1)*pi*x^2/sin(pi*x)^2 with different N
4 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
% Define the value of N
format long
N = 1000;
% Define a range of x values
%x = linspace(-0.01*pi,0.01*pi, 1000); % You can adjust the range and number of points as needed
x = -1.5:0.00001:1.5; % You can adjust the range and number of points as needed
% Calculate the function values
y = (sin((N+1) *pi* x).^2) ./ (sin(pi*x).^2);
% Plot the function
plot(x, y)
title('Plot of (sin((N+1) *pi* x).^2) ./ (sin(pi*x).^2)')
grid on
when I change N=10000;
% Define the value of N
format long
N = 10000;
% Define a range of x values
%x = linspace(-0.01*pi,0.01*pi, 1000); % You can adjust the range and number of points as needed
x = -1.5:0.00001:1.5; % You can adjust the range and number of points as needed
% Calculate the function values
y = (sin((N+1) *pi* x).^2) ./ (sin(pi*x).^2);
% Plot the function
plot(x, y)
title('Plot of (sin((N+1) *pi* x).^2) ./ (sin(pi*x).^2)')
grid on
the limit of (sin((N+1) *pi* x).^2) ./ (sin(pi*x).^2) is (N+1)^2, which means for N=1000, the maximum is 1000000, which is consistent with the plot
but for N=10000, the maximum is 1e8, which has something wrong.
4 commentaires
Torsten
le 10 Sep 2023
Modifié(e) : Torsten
le 10 Sep 2023
The precision of usual floating point arithmetics reaches its limits. If you want a more precise plot, use symbolic computation.
And don't evaluate your y exactly at integer points. It will give a division by zero for which the result you get is unpredictable.
Réponses (3)
Steven Lord
le 10 Sep 2023
To eliminate the effect of the floating point approximation to π, try using sinpi instead of sin.
% Define the value of N
format long
N = 1000;
% Define a range of x values
%x = linspace(-0.01*pi,0.01*pi, 1000); % You can adjust the range and number of points as needed
x = -1.5:0.00001:1.5; % You can adjust the range and number of points as needed
% Calculate the function values
f = @(x) (sinpi((N+1)*x).^2) ./ (sinpi(x).^2);
y = f(x);
% Plot the function
plot(x, y)
title('Plot of (sin((N+1) *pi* x).^2) ./ (sin(pi*x).^2)')
grid on
Why the spikes at the integer values? Well, those peaks aren't exactly at the integer values but very close to them, since the function is undefined (0/0) at exact integer values of x.
f([-1 0 1])
1 commentaire
Bruno Luong
le 10 Sep 2023
For older version one can use sind
f = @(x) (sind(180*(N+1)*x).^2) ./ (sind(180*x).^2)
Sam Chak
le 10 Sep 2023
Hi @Daniel Niu
syms x N
y = limit((sin((N + 1)*pi*x)^2)/(sin(pi*x)^2), x, 0)
N = sym(9999);
y = limit((sin((N + 1)*pi*x)^2)/(sin(pi*x)^2), x, 0)
Paul
le 10 Sep 2023
Modifié(e) : Paul
le 10 Sep 2023
Hi Daniel,
I think you are correct that the limit of y as x approaches any integer is (N+1)^2. However, all the compuations are numeric, and numeric computations don't know about limits, they just compute what's given with the finite precisicion of the data types.
Here's the first case:
format long
N = 1000;
% Define a range of x values
x = -1.5:0.00001:1.5; % You can adjust the range and number of points as needed
% Calculate the function values
y = (sin((N+1) *pi* x).^2) ./ (sin(pi*x).^2);
There is a value of x that is exactly equal to 0. For x == 0, the value of y is
y(x==0)
and if you zoom in on the plot around x = 0 you'll see a gap.
Plot what's happening around x = 1, and add a point for the expected value.
plot(x,y,'-o',1,(N+1)^2,'x'),grid,xlim([0.9999 1.0001])
Clearly, we see something fishy at x = 1.
There is a value in x that is exactly equal to 1
i1 = find(x==1)
and its corresponding value of y
y(i1)
The problem is that, in double precision, the Matlab function pi cannot exactly represent mathematical pi. Hence, the numerator of y, which should ideally be zero, isn't zero. Same for the denominator
[(sin((N+1) *pi* x(i1)).^2) , (sin(pi*x(i1)).^2)]
Those values get divided and we get what we get, which is basically junk. If you zoom in further, you'll see that the plot at x = 1 is y(i1).
Here's the same thing for the second case, and you'll see that
N = 10000;
% Define a range of x values
%x = linspace(-0.01*pi,0.01*pi, 1000); % You can adjust the range and number of points as needed
x = -1.5:0.00001:1.5; % You can adjust the range and number of points as needed
% Calculate the function values
y = (sin((N+1) *pi* x).^2) ./ (sin(pi*x).^2);
i1 = find(x == 1)
y(i1)
figure
plot(x,y,'-o',1,(N+1)^2,'x'),grid,xlim([0.9999 1.0001])
Basically, we're suffering from finite precision and not getting exact values of 0 in the numerator and denominator of y for nonzero, integer values of x.
0 commentaires
Voir également
Catégories
En savoir plus sur Number Theory dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!