Ode45 I get wrong results.

function dydt = idn (t, y)
dydt = zeros(2,1);
dydt(1)=y(2);
dydt(2)=-((9.81*sin(y(1)))/0.5);


dydt = [0.78 0];
[t1,y1] = ode45(@idn, [0 20], dydt);
plot(t1,y1(:,2));
grid on;

This is a code for simple pendulum. When I use this code group, acceleration comes out wrong(it goes between max 3.. when it should be 13.). What am I doing wrong? When I use sind and put 45 degrees instead of sin and 0.78 radia, this time I get 30 ish something as a max result when it shouldn't change.

Original Eqn: Theta''=-g*sin(Theta)/ l

[madhan ravi]

why 13 any reason?

[Arda Nova]

dydt(2)=-((9.81*sin(y(1)))/0.5);

At the bottom, acceleration is maximum. If you put 270 degrees acceleration comes 19.62 rad/sec^2 (Sorry I had made a mistake at taking the degree. Tried it too, wrong result). One think is for sure, acceleration should not be affected from me using radian or degrees. I don't even know where does the program gets the 3.something. If I knew then I would maybe solve the issue but, couldn't figure it out.

Answer 1

[Stephan]

Hi,

ode45 integrates your second order ode two times. So the result of ode45 will be angle and velocity over time. To get acceleration you have to calculate it seperatly:

% initial conditions of angle and velocity
dydt = [pi/2; 0];
% integrat ode to get angle an velocity
[t1,y1] = ode45(@idn, [0 20], dydt);
% calculate acceleration
y1(:,3) = -9.81.*sin(y1(:,1))./0.5;
% plot results
plot(t1,y1(:,1),t1,y1(:,2),t1,y1(:,3));
grid on;

function dydt = idn(t, y)
dydt = zeros(2,1);
dydt(1)=y(2);
dydt(2)=-((9.81*sin(y(1)))/0.5);
end

Best regards

Stephan

[Arda Nova]

I am very thankful to you, good sir. You have no idea how much you have helped.