ODE45 what to put into the function argument

16 views (last 30 days)
Rick on 14 Aug 2014
Commented: Star Strider on 14 Aug 2014
Hello, In a file a have an attachment to a problem of modeling a pendulum. It has a format in the code that I have to follow exactly as is, and I was wondering if I could get some help with what to put into the ODE 45 part. Here is the code
function ydot = InvPend1(y,g,L,k1,k2)
% initialize ydot array
ydot = zeros(2,1);
u = k1*y(1)+k2*y(2);
ydot(1) = y(2);
ydot(2) = (g/L)*sin(y(1))-(1/L)*cos(y(1))*u;
g = 9.8; L=1;k1=12;k2=2;
% Set initial conditions
y0 =[10*pi/180;0];
% Set Simulation time span
tspan = [0 8];
% Define function a handle
[T,Y] = ode45(@(y,g,L,u) InvPend1,tspan,y0)
But I have a lot of difficulty figuring out what to put in the InvPed1 argument, and in general I am struggling on how to use the ode45 function because I am having a difficult time figuring out what to put as its first argument.
I could also use some help with the subplots. How do I plot theta vs. time? plot(T,Y)?? But what is U array? I think my problem here is that I don't fully comprehend what the output [T, Y] will be. Will Y be the original matrix [Y(1);Y(2)], meaning this ode45 will output the theta and derivative of theta for all times in the time span?
% Use subplot to create top plot
subplot(2,1,1);
% plot theta (degrees) vs. time
plot(T,Y);
% Set xlabel
xlabel('time-seconds');
% Set yLabel
ylabel('theta-deg');
% Compute U array
U = k1*Y(1) + k2*Y(2);
% plot u vs. time
plot(T,U);
% Set xlabel
xlabel('time-seconds');
% Set ylabel

Star Strider on 14 Aug 2014
The ODE function arguments have a specific format, and have to include your independent variable (usually t) as the first argument (even if you don’t use it):
function ydot = InvPend1(t,y,g,L,k1,k2)
[T,Y] = ode45(@(t,y) InvPend1(t,y,g,L,k1,k2),tspan,y0)
I don’t know what you’re doing with this line, so I suggest you delete it:
I got your code to work, with ‘InvPend1’ as an anonymous function. Your function file is otherwise fine as it is (once you change the first line to include t), so my code will work as an illustration. (Anonymous functions are easier for me to work with.)
You calculate ‘u’ within your ‘InvPend1’, so you don’t need to pass it as an argument. If you want to use it in your plot, you simply need to calculate it:
u = k1*Y(:,1)+k2*Y(:,2);
since your calculated values for Y is a (81x2) matrix (when I ran it), you need to refer to your values by column references.
I don’t know what you want to do w.r.t. your ‘subplots’ because I don’t know your requirements. I plotted everything in one plot here.
InvPend1 = @(t,y,g,L,k1,k2) [y(2); (g/L)*sin(y(1))-(1/L)*cos(y(1))*(k1*y(1)+k2*y(2))];
g = 9.8; L=1;k1=12;k2=2;
% Set initial conditions
y0 =[10*pi/180;0];
% Set Simulation time span
tspan = [0 8];
[T,Y] = ode45(@(t,y) InvPend1(t,y,g,L,k1,k2),tspan,y0);
u = k1*Y(:,1)+k2*Y(:,2);
figure(1)
plot(T, Y, T, u)
grid
legend('Y(t)', 'Y''(t)', 'u(t)', 'Location','NE')
Star Strider on 14 Aug 2014
Great!
Correct. Only the call to ode45 needs to.