I get this error:"This system does not seem to be linear."

4 Ansichten (letzte 30 Tage)
Arda Nova
Arda Nova am 3 Jan. 2019
Kommentiert: Arda Nova am 4 Jan. 2019
When I try this code for simple pendulum, I get results:
syms pt(t) th(t)
m=1; l=0.5; g=9.81;
e1= diff(th)*(m*l^2)==pt;
e2= diff(pt)==-m*g*l*sin(th);
vars = [pt(t); th(t)];
V = odeToVectorField([e1,e2]);
M = matlabFunction(V, 'vars', {'t','Y'});
interval = [0 5];
y0 = [0; pi/4];
ySol = ode45(M,interval,y0);
tValues = linspace(interval(1),interval(2),1000);
a= deval(ySol,tValues,1)/(m*l^2);
plot(tValues,a)
But when I use it for triple pendulum, it gives error. Couldn't solve it. Sorry if it's simple to figure out. Really new here.
syms theta1(t) theta2(t) theta3(t) p1(t) p2(t) p3(t)
m1=1; m2=1; m3=1; l1=1; l2=1;
l3=1; g=9.81; tau1=0; tau2=0; tau3=0;
I1=0; I2=0; I3=0;
e1= diff(theta1)*(I1+(m1+m2+m3)*l1^2)...
+diff(theta2)*(m2+m3)*l1*l2*cos(theta1-theta2)+...
diff(theta3)*m3*l1*l3*cos(theta1-theta3)==p1;
e2= diff(theta1)*(m2+m3)*l1*l2*cos(theta1-theta2)+...
diff(theta2)*(I2+(m2+m3)*l2^2)+...
diff(theta3)*m3*l2*l3*cos(theta2-theta3)==p2;
e3= diff(theta1)*m3*l1*l3*cos(theta1-theta3)+...
diff(theta2)*m3*l2*l3*cos(theta2-theta3)+...
diff(theta3)*(I3+m3*l3^2)==p3;
e4= diff(p1)== tau1-tau2-(m2+m3)*diff(theta1)*diff(theta2)*sin(theta1-theta2)...
-m3*diff(theta1)*diff(theta3)*l1*l3*sin(theta1-theta3)...
-(m1+m2+m3)*g*l1*cos(theta1);
e5= diff(p2)==tau2-tau3+(m2+m3)*diff(theta1)*theta2*l1*l2*sin(theta1-theta2)...
-m3*diff(theta2)*diff(theta3)*l2*l3*sin(theta2-theta3)...
-(m2+m3)*g*l2*cos(theta2);
e6= diff(p3)==tau3+(m3)*diff(theta1)*diff(theta3)*l1*l3*sin(theta1-theta3)...
+m3*diff(theta2)*diff(theta3)*l2*l3*sin(theta2-theta3)...
-m3*g*l3*cos(theta3);
vars= [theta1(t);theta2(t);theta3(t);p1(t);p2(t);p3(t)];
V = odeToVectorField([e1,e2,e3,e4,e5,e6]);
M = matlabFunction(V,'vars', {'t','Y'});
I have get that error in simple pendulum too, it was pt==diff(th)*(m*l^2), then I put the pt to the end, and it's solved. In triple pendulum I tried leaving diff(theta1) alone didn't work, tried to this code too, but nothing changed. Original equations are:
Adsız1.png
Adsız2.png

Akzeptierte Antwort

Torsten
Torsten am 4 Jan. 2019
The product of differentials in your equations (diff(theta1)*diff(theta3), e.g.) makes it impossible to use ODE45.
I don't know if it can be applied directly, but ODE15I is the correct solver to use in this case.
Best wishes
Torsten.
  4 Kommentare
Torsten
Torsten am 4 Jan. 2019
function main
y0 = [pi/2; pi/2; pi/2; 0; 0; 0];
yp0=[0; 0; 0; 0; 0; 0;];
[y0,yp0] = decic(@odennotfunatall,0,y0,[pi/2 pi/2 pi/2 0 0 0],yp0,[]);
[t,y] = ode15i(@odennotfunatall,[0 5],y0,yp0);
plot(t,y)
end
function hell2 = odennotfunatall(~,y,yp)
m1=1; m2=1; m3=1; l1=1; l2=1;
l3=1; g=9.81; tau1=0; tau2=0; tau3=0; I1=0; I2=0; I3=0;
hell2=zeros(6,1);
hell2(1)=yp(1)*(I1+(m1+m2+m3)*l1^2)...
+yp(2)*(m2+m3)*l1*l2*cos(y(1)-y(2))+...
yp(3)*m3*l1*l3*cos(y(1)-y(3))-y(4);
hell2(2)=yp(1)*(m2+m3)*l1*l2*cos(y(1)-y(2))+...
yp(2)*(I2+(m2+m3)*l2^2)+...
yp(3)*m3*l2*l3*cos(y(2)-y(3))-y(5);
hell2(3)=yp(1)*m3*l1*l3*cos(y(1)-y(3))+...
yp(2)*m3*l2*l3*cos(y(2)-y(3))+...
yp(3)*(I3+m3*l3^2)-y(6);
hell2(4)=-yp(4)+tau1-tau2-(m2+m3)*yp(1)*yp(2)*sin(y(1)-y(2))...
-m3*yp(1)*yp(3)*l1*l3*sin(y(1)-y(3))...
-(m1+m2+m3)*g*l1*cos(y(1));
hell2(5)=-yp(5)+tau2-tau3+(m2+m3)*yp(1)*yp(2)*l1*l2*sin(y(1)-y(2))...
-m3*yp(2)*yp(3)*l2*l3*sin(y(2)-y(3))...
-(m2+m3)*g*l2*cos(y(2));
hell2(6)=-yp(6)+tau3+(m3)*yp(1)*yp(3)*l1*l3*sin(y(1)-y(3))...
+m3*yp(2)*yp(3)*l2*l3*sin(y(2)-y(3))...
-m3*g*l3*cos(y(3));
end
Arda Nova
Arda Nova am 4 Jan. 2019
Yeah the zeros, couldn't find a way to put it, thank you very much.

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (1)

madhan ravi
madhan ravi am 3 Jan. 2019
Just follow the same way showed in your previous question?
  3 Kommentare
madhan ravi
madhan ravi am 3 Jan. 2019
No problem , please recheck your equations take your time there are only four equations whereas you have 6 equations in your code.
Arda Nova
Arda Nova am 3 Jan. 2019
Oh sorry, upper picture is streched out. At the right of "(60)" (on the picture), there is two more equations.

Melden Sie sich an, um zu kommentieren.

Kategorien

Mehr zu Numerical Integration and Differential Equations finden Sie in Help Center und File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by