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
Categories
Find more on Assembly in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!