How to plot a function which is defined on different subintervals

2021 3 7
1 个回答
回答已采纳
更新时间：2021 3 10
3 次查看（30 天）

0 个投票

I am trying to plot the function ,
But I don't know how to write the code for the definition of the function f which is given on different subintervals.

4 个评论
显示 2更早的评论 隐藏 2更早的评论

Jan
Jan 2021-3-7
At the points i, f has two different values. Shouldn't one of the intervals exclude i?
Cris19
Cris19 2021-3-7
编辑：Cris19 2021-3-7
@JanThank you for the comment. Indeed, f had two different values at i and I have corrected the function. Now it should be continuous.

请先登录，再进行评论。

请先登录，再回答此问题。

 采纳的回答

Walter Roberson
Walter Roberson 2021-3-7

1 个投票

There are several methods available. The most straight forward is to write a function that loops over the inputs, testing each one to decide what the result should be.
function y = f(X)
y = zeros(size(X));
for K = 1 : numel(X)
x = X(K);
if x < 2
y(K) = 0;
elseif x <= 5
y(K) = x.^2;
elseif x <= 8
y(K) = x-x.^3;
else
y(K) = 0;
end
end
With regards to those intervals you need, think about floor(x+1/2)

11 个评论
显示 9更早的评论 隐藏 9更早的评论

Cris19
Cris19 2021-3-8
@Walter Roberson Thank you, but I still don't know how to write the code for f(x) on [2,500]. After defining and plotting the function f(x) on the interval [2,500], I also want to plot on [2,500] the solution of the second order ODE y''+f(x)y'+y=0 with the initial data y(2)=0.1, y'(2)=0.1. This is why I think we should know the function f(x) on the whole interval [2,500].
Walter Roberson
Walter Roberson 2021-3-8
The function is defined in sections. You would have a difficult time writing it as a formula unless you use piecewise or an infinite sum of heaviside.
Cris19
Cris19 2021-3-8
How exactly can I use the piecewise syntax? Would you please write the code? I am not too familiarized with Matlab.
Cris19
Cris19 2021-3-8
编辑：Cris19 2021-3-8
I tried to use this code:
for i= 2:500
Bounds_1 = [i - (1/i^2),i];
Bounds_2 = [i,i + (1/i^2)];
Bounds = [Bounds_1 Bounds_2];
Min = min(Bounds);
Max = max(Bounds);
f1 = @(x) (i^2.*x - i^3 + 1).*(Bounds_1(1) < x & x <= Bounds_1(2));
f2 = @(x) (-i^2.*x + i^3 + 1).*(Bounds_2(1) < x & x <= Bounds_2(2));
f = @(x) f1(x) + f2(x);
end
fplot(f);
xlim([Min-2 Max+2]);
ylim([0 1.1]);
But the only plotted 'spike' is the last one, for i=500. How can I keep all spikes in the expression of f?
Walter Roberson
Walter Roberson 2021-3-8
for i= 2:500
i is a positive integer
Bounds_1 = [i - (1/i^2),i];
With i being real valued, 1/i^2 is real valued and positive (or infinite and positive), and subtracting a positive amount from a real value always gives you a result that is less than the original amount.
Bounds_2 = [i,i + (1/i^2)];
With i being real-valued, 1/i^2 is real valued and positive (or infinite and positive), and adding a positive amount to a real value always gives you a result that is greater than the original amount
Bounds = [Bounds_1 Bounds_2];
That would be a list [i - (1/i^2),i,i,i + (1/i^2)]
Min = min(Bounds);
As demonstrated above, i-1/i^2 is always less than i, and i+1/i^2 is always greater than i, so we can quickly see that the minimum value must be i-1/i^2
Max = max(Bounds);
and the maximum value must be 1+1/i^2 . So there is no point in creating the Bounds list and taking the min and max: we already know which will be the min and which will be the max
end
You do not use Min and Max inside the loop, so you overwrite the value each time. Therefore at the end of the loop, the value will be whatever was calculated on the last iteration, the one for i=500
fplot(f);
xlim([Min-2 Max+2]);
and you create your xlim according to only that iteration.
f1 = @(x) (i^2.*x - i^3 + 1).*(Bounds_1(1) < x & x <= Bounds_1(2));
f2 = @(x) (-i^2.*x + i^3 + 1).*(Bounds_2(1) < x & x <= Bounds_2(2));
f = @(x) f1(x) + f2(x);
end
You overwrite all of f1, f2, and f each iteration, so in that code there is no point doing any iteration except the last one.
Walter Roberson
Walter Roberson 2021-3-8
Ran in:
Ran in:
You could try this:
f = @(x) zeros(size(x));
for i= 2:500
Bounds_1 = [i - (1/i^2),i];
Bounds_2 = [i,i + (1/i^2)];
f1 = @(x) (i^2.*x - i^3 + 1).*(Bounds_1(1) < x & x <= Bounds_1(2));
f2 = @(x) (-i^2.*x + i^3 + 1).*(Bounds_2(1) < x & x <= Bounds_2(2));
f = @(x) f(x) + f1(x) + f2(x);
end
fplot(f, [1 6.5]); ylim([-0.1 1.1])
It is difficult to plot much further than this, as the spikes get so narrow that fplot misses them. By 500 the spikes are only 1/2500 wide, and assuming you want at least 5 points around 500, that implies you would have to use a sampling frequency no less than 125000 Hz, times your range of at least 501 (you need to include 1 and 501 exactly in the sample so that the fourier transform does not accidentally think that the pattern repeats indefinitely when instead it stops at 1 and 501.)
Cris19
Cris19 2021-3-9
Thank you so very much for your thorough comments! The code works very well!
Cris19
Cris19 2021-3-9
编辑：Cris19 2021-3-9
One more question, if possible... I try to solve to second order ODE y''+f(x)y'+y=0, whre f(x) is the function above. I am interested to plot on [0,10] the solution of this ODE. for the initial data y(0)=y'(0)=0.1. I used the function (where f(x) is defined with the code you wrote)
function dwdt=susp(t,w)
f = @(x) zeros(size(x));
for i= 2:10
Bounds_1 = [i - (1/i^2),i];
Bounds_2 = [i,i + (1/i^2)];
f1 = @(x) (i^2.*x - i^3 + 1).*(Bounds_1(1) < x & x <= Bounds_1(2));
f2 = @(x) (-i^2.*x + i^3 + 1).*(Bounds_2(1) < x & x <= Bounds_2(2));
f = @(x) f(x) + f1(x) + f2(x);
end
dwdt=zeros(2,1);
dwdt(1)=w(2);
dwdt(2)=-f*w(2)-w(1);
end
and the commands
tspan = [0 10];
z0=[0.1 0.1];
[t,z] = ode45(@(t,z) susp(t,z), tspan, z0);
figure
plot(t,z(:,1),'r');
legend('x(t)');
But I get several errors...
Unary operator '-' is not supported for operand of type 'function_handle'.
Error in susp (line 12)
dwdt(2)=-f*w(2)-w(1);
Error in @(t,z)susp(t,z)
Error in odearguments (line 90)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 115)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
How could I fix them?
Walter Roberson
Walter Roberson 2021-3-9
It is not possible to use the ode*() functions with discontinuous functions (or rather, if it works at all it will give the wrong answer.)
The exception to that is if the ode*() call has continuous second derivatives for the entire call then it is allowed. So you would have to break up the time into segments: [2, 2+1/4], [2+1/4, 3-1/9], [3-1/9, 3], [3,3+1/9], [3+1/9, 4-1/16], [4-1/16, 4], [4,4+1/6] and so on.
Cris19
Cris19 2021-3-10
Thank you. I wonder if there is an alternative solution. I am interested to find the asymptotic behavior at +infinity of the solution of that ODE. So is this why I think it is interesting to see the plotting of the solution on intervals larger and larger, such as [0,10], [0,500] etc.
I really don't know how to code this...
Walter Roberson
Walter Roberson 2021-3-10
At asymptopic behaviour is
syms n f(x)
D = symsum(dirac(x-n), n, 1, 200)
df = diff(f)
d2f = diff(df)
eqn = d2f + D*df + f == 0
dsolve([eqn, f(0)==0])
but in practice you will not get any solution. Even if you reduce it down to dirac at one particular integer you are not going to get a solution.

请先登录，再进行评论。

更多回答（0 个）

类别

在 帮助中心 和 File Exchange 中查找有关 Programming 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by