Effacer les filtres
Effacer les filtres

how to get the correct plot of function sin(N+1)*p​i*x^2/sin(​pi*x)^2 with different N

7 vues (au cours des 30 derniers jours)
Daniel Niu
Daniel Niu le 10 Sep 2023
Commenté : Daniel Niu le 10 Sep 2023
% 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
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.
Daniel Niu
Daniel Niu le 10 Sep 2023
how do you perform that and get the right plot. thank you

Connectez-vous pour commenter.

Réponses (3)

Steven Lord
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])
ans = 1×3
NaN NaN NaN

Sam Chak
Sam Chak le 10 Sep 2023
I checked using the limit() function.
syms x N
y = limit((sin((N + 1)*pi*x)^2)/(sin(pi*x)^2), x, 0)
y = 
N = sym(9999);
y = limit((sin((N + 1)*pi*x)^2)/(sin(pi*x)^2), x, 0)
y = 
100000000
  1 commentaire
Daniel Niu
Daniel Niu le 10 Sep 2023
thank you,Sam,you are right I can not understand why the maximum is different at 0 and 1 when N=10000

Connectez-vous pour commenter.


Paul
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)
ans =
NaN
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)
i1 =
250001
and its corresponding value of y
y(i1)
ans =
5.237999158081999e+03
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)]
ans = 1×2
1.0e-28 * 0.785574047890805 0.000149975978266
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)
i1 =
250001
y(i1)
ans =
2.599231705387884e+08
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.

Tags

Community Treasure Hunt

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

Start Hunting!

Translated by