Using ode45 to solve a Second Order Non-linear ODE

Question

Nivedita Tanksali 2021-4-4

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/792067-using-ode45-to-solve-a-second-order-non-linear-ode

评论： Alan Stevens 2021-4-4

在 MATLAB Online 中打开

I'm trying to solve a part of the following equation :

The part i'm trying to solve is

or

I've rewritten itt in the form of two first order equations, like so:

In vector form,

And my code is as follows:

tspan = 0:0.0033:100;
a=10*(pi/180);
b=0;
s0 = [a; b];
[t,s] = ode45(@pend_P,tspan,s0);
S1 = s(:,1);
S2 = s(:,2);
plot(t,S1*180/pi)
function sdot = pend_P(t,s)
a=10*(pi/180);
b=0;
% Expansion of s(1) or p :
p1 =0.000000001906*a^3 + (-0.0000007948)*a^2 + 0.00009188*a + (-0.003481);
p2 =0.00000915*a^2 + (-0.0009381)*a + 0.05331;
p3 =1*((-0.0001542)*a^2 + (-0.006078)*a + (-2.089));
p4 =a; 
s(1) =( (p1*t.^3) + (p2*t.^2) + (p3*t) + (p4) );
sdot = [s(2);  (2*s(2).^2/s(1))];     
end 

The curve I'm getting is a horixontal line though, as the value of S1 appears to be constant. Any idea where I've gone wrong?

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Answer 1

Alan Stevens 2021-4-4

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/792067-using-ode45-to-solve-a-second-order-non-linear-ode#answer_666527

在 MATLAB Online 中打开

Within your function pend you have

function sdot = pend_P(t,s)
a=10*(pi/180);
b=0;
...etc

This resets a and b to the initial conditions every time pend is called. You should have

function sdot = pend_P(t,s)
a=s(1);
b=s(2);
...etc

Also, if you start b (s(2)) at zero, then there is no way for anything to change given that

sdot = [s(2); (2*s(2).^2/s(1))];

i.e. both terms in sdot will stay at zero.

2 个评论
显示无隐藏无

Nivedita Tanksali 2021-4-4

I tried what you suggested, but it made no difference.

Also, 'a' is the initial displacement, it is an initial condition. It is not equal to s(1). And b is not equal to s(2) either, as it is the inital velocity.

Alan Stevens 2021-4-4

You pass a and b to the function as initial conditions, so, because b is zero, this means s(2) is intially zero, which means sdot returns 0, which means nothing changes! Because of the structure of sdot, the initial value of s(2) must be non-zero if you want anything to change.

If you intend that the values of p are calculated by the initial value of a only, and are unchanged after that, they are better calcuated outside of pend, and then passed in, otherwise they are recaculated on every call to pend, which is inefficient.

请先登录，再进行评论。