I get this error:"This system does not seem to be linear."
1 view (last 30 days)
Show older comments
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:
0 Comments
Accepted Answer
Torsten
on 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 Comments
Torsten
on 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
More Answers (1)
madhan ravi
on 3 Jan 2019
Just follow the same way showed in your previous question?
3 Comments
madhan ravi
on 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.
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!