Problem with ode45 code ?
显示 更早的评论
Hi, i'm really new to Matlab. I am using ode45 to solve this differential equation:
y1prime = y2
y2prime = -y1 - (1/t)*(y2)
I am plotting it over the closed interval 1, 50.
It has the intial conditions:
y1(0)=0
y2(0)=0
y2prime(0)=0.
However i keep getting error messages about my m-script. I am using a code from quite an old book, so i'm not sure if that has anything to do with it. I would appreciate if anybody could have a quick look to see where i've gone wrong. Thank you
function yprime = pend(t,y)
yprime = [y(2); -y(1)-((y(2))/t)];
tspan = [1 50];
yazero = [0;0]; ybzero = [0;0]; yczero = [0;0];
[ta, ya] = ode45(@pend, tspan, yazero);
[tb, yb] = ode45(@pend, tspan, ybzero);
[tc, yc] = ode45(@pend, tspan, yczero);
[y1,y2] = meshgrid(-5:.5:5,-3:.5:3);
Dy1Dt = y2; Dy2Dt = -y(1)-((y(2))/t);
quiver(y1,y2,Dy1Dt,Dy2Dt)
hold on
plot(ya(:,1),ya(:,2),yb(:,1),yb(:,2),yc(:,1),yc(:,2))
axis equal, axis([-5 5 -3 3])
xlabel y_1(t), ylabel y_2(t), hold off
end
5 个评论
John D'Errico
2016-2-12
编辑：John D'Errico
2016-2-12
I seriously doubt that ANY old book has MATLAB code in it where you create a function that calls ODE45, inside passing ODE the name of the function that you just created, to create a recursive nightmare of a loop. So you have seriously hacked around with whatever code you did find.
As well, you state:
"I am plotting it over the closed interval 1,50.
It has the intial conditions: y1(0)=0 y2(0)=0 y2prime(0)=0."
But then specify initial conditions at 0, not 1.
I also have no idea why you think you need to solve it three times, all with the same initial conditions. Would you get different answers each time? No.
Finally, you have only two variables, not three. So specifying three initial conditions will be difficult. Anyway, y2prime appears to be trivially zero IF we assume that you really mean to start at 1 for t, NOT 0.
Jenny Briggs
2016-2-12
Brian Russell
2020-3-25
I see your problem. The old book you refer to is "MATLB Guide" by Higham and Higham, which I think is in general a brilliant book. However, they often leave out statements in their code which they assume their reader will know to put in. I also struggled with this example until I realized that their function is very short. In fact, it is only two lines:
function yprime = pend(t,y)
yprime = [y(2); -y(1)-((y(2))/t)];
There should be an end after these two lines. This short function should be saved in a file called pend.m. The rest of your code should be put in another program, say pend_run.m. This will fix your problem.
Brian Russell
2020-3-25
I meant "MATLAB Guide", of course!
回答（1 个）
John D'Errico
2016-2-12
So, just suppose I decide to try to solve your problem, assuming that you really intended to solve it over the interval [1,50]. That would take no more complicated code than this:
yp = @(y,t) [y(2);-y(1) - y(2)./t);
[tout,yout] = ode45(yp,[1,50],[0 0]);
The result would be rather boring though, since you have postulated a differential equation system of the form:
y1' = 0
y2' = 0
Note that since y1(1) = y2(1) = 0 at the beginning, they will STAY at zero forever. And ever. And ever.
Kind of boring.
类别
在 帮助中心 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息
2016-2-12
2020-3-25
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
