Tell me more ×
Computational Science Stack Exchange is a question and answer site for scientists using computers to solve scientific problems. It's 100% free, no registration required.
up vote 0 down vote favorite

This post is the ccontinuation of this pose enter link description here

deDLMMddq1=2*infec1x1*q5tp^2*cos(q3t) - (L1^2*q3tp^2*roL*cos(q3t))/2 - ...  
infec1x1*q3tp^2*q5t^2*cos(q3t) +...
2*infec2q72x2*q7tp^2*cos(q3t)*cos(sita) - 2*infe2x2*q5tp^2*cos(q3t)*sin(sita) -...
2*infe2x2*q5tp^2*cos(sita)*sin(q3t) - 2*infec2q72x2*q7tp^2*sin(q3t)*sin(sita) -...
2*infe1x1*q3tp*q5tp*cos(q3t) + infe1x1*q3tp^2*q5t*sin(q3t) -...
4*infec1x1*q3tp*q5t*q5tp*sin(q3t) -...
L1*L2*q3tp^2*roL*cos(q3t) - infec2q54x2*q3tp^2*q5t^4*cos(q3t)*cos(sita) +...
12*infec2q54x2*q5t^2*q5tp^2*cos(q3t)*cos(sita) -...
infec2q72x2*q3tp^2*q7t^2*cos(q3t)*cos(sita) + ...
(L2^2*q3tp^2*roL*sin(q3t)*sin(sita))/2 + infe2x2*q3tp^2*q5t^2*cos(q3t)*sin(sita) +...
infe2x2*q3tp^2*q5t^2*cos(sita)*sin(q3t)+infec2q54x2*q3tp^2*q5t^4*sin(q3t)*sin(sita)-...
12*infec2q54x2*q5t^2*q5tp^2*sin(q3t)*sin(sita) +... 
infec2q72x2*q3tp^2*q7t^2*sin(q3t)*sin(sita) -...
2*infe3x2*q3tp*q7tp*cos(q3t)*cos(sita) + 2*infe3x2*q3tp*q7tp*sin(q3t)*sin(sita) +...
2*infec2q52q7x2*q5tp^2*q7t*cos(q3t)*cos(sita) +...
infe3x2*q3tp^2*q7t*cos(q3t)*sin(sita)+infe3x2*q3tp^2*q7t*cos(sita)*sin(q3t) -...
2*infec2q52q7x2*q5tp^2*q7t*sin(q3t)*sin(sita)-(L2^2*q3tp^2*roL*cos(q3t)*cos(sita))/2-...
8*infec2q54x2*q3tp*q5t^3*q5tp*cos(q3t)*sin(sita)-...
8*infec2q54x2*q3tp*q5t^3*q5tp*cos(sita)*sin(q3t) -...
2*infec2q52q7x2*q3tp*q5t^2*q7tp*cos(q3t)*sin(sita) -...
2*infec2q52q7x2*q3tp*q5t^2*q7tp*cos(sita)*sin(q3t) -... 
infec2q52q7x2*q3tp^2*q5t^2*q7t*cos(q3t)*cos(sita) + ...
infec2q52q7x2*q3tp^2*q5t^2*q7t*sin(q3t)*sin(sita) -... 
4*infe2x2*q3tp*q5t*q5tp*cos(q3t)*cos(sita) +...
4*infec2q52q7x2*q5t*q5tp*q7tp*cos(q3t)*cos(sita) - ...
4*infec2q72x2*q3tp*q7t*q7tp*cos(q3t)*sin(sita) - ...
4*infec2q72x2*q3tp*q7t*q7tp*cos(sita)*sin(q3t) + ...
4*infe2x2*q3tp*q5t*q5tp*sin(q3t)*sin(sita) - ...
4*infec2q52q7x2*q5t*q5tp*q7tp*sin(q3t)*sin(sita) -...
4*infec2q52q7x2*q3tp*q5t*q5tp*q7t*cos(q3t)*sin(sita) -...
4*infec2q52q7x2*q3tp*q5t*q5tp*q7t*cos(sita)*sin(q3t);


deDLMMddq2=g*mB + 2*infec1x1*q5tp^2*sin(q3t) + L1*g*roL + L2*g*roL -...
(L1^2*q3tp^2*roL*sin(q3t))/2 +...
2*infe2x2*q5tp^2*cos(q3t)*cos(sita) - infec1x1*q3tp^2*q5t^2*sin(q3t) +...
2*infec2q72x2*q7tp^2*cos(q3t)*sin(sita) + 2*infec2q72x2*q7tp^2*cos(sita)*sin(q3t) - ...
2*infe2x2*q5tp^2*sin(q3t)*sin(sita) - 2*infe1x1*q3tp*q5tp*sin(q3t) -...
infe1x1*q3tp^2*q5t*cos(q3t) -...
(L2^2*q3tp^2*roL*cos(q3t)*sin(sita))/2 - (L2^2*q3tp^2*roL*cos(sita)*sin(q3t))/2 - ...
infe2x2*q3tp^2*q5t^2*cos(q3t)*cos(sita) - L1*L2*q3tp^2*roL*sin(q3t) - ...
infec2q54x2*q3tp^2*q5t^4*cos(q3t)*sin(sita) -...
infec2q54x2*q3tp^2*q5t^4*cos(sita)*sin(q3t) +...
12*infec2q54x2*q5t^2*q5tp^2*cos(q3t)*sin(sita) + ...
12*infec2q54x2*q5t^2*q5tp^2*cos(sita)*sin(q3t) -...
infec2q72x2*q3tp^2*q7t^2*cos(q3t)*sin(sita) -...
infec2q72x2*q3tp^2*q7t^2*cos(sita)*sin(q3t) + ...
infe2x2*q3tp^2*q5t^2*sin(q3t)*sin(sita) - 2*infe3x2*q3tp*q7tp*cos(q3t)*sin(sita) - ...
2*infe3x2*q3tp*q7tp*cos(sita)*sin(q3t) - infe3x2*q3tp^2*q7t*cos(q3t)*cos(sita) +...
2*infec2q52q7x2*q5tp^2*q7t*cos(q3t)*sin(sita) +...
2*infec2q52q7x2*q5tp^2*q7t*cos(sita)*sin(q3t) +...
infe3x2*q3tp^2*q7t*sin(q3t)*sin(sita) + 4*infec1x1*q3tp*q5t*q5tp*cos(q3t) + ...
8*infec2q54x2*q3tp*q5t^3*q5tp*cos(q3t)*cos(sita) +...
2*infec2q52q7x2*q3tp*q5t^2*q7tp*cos(q3t)*cos(sita) -...
8*infec2q54x2*q3tp*q5t^3*q5tp*sin(q3t)*sin(sita) -...
2*infec2q52q7x2*q3tp*q5t^2*q7tp*sin(q3t)*sin(sita) -...
infec2q52q7x2*q3tp^2*q5t^2*q7t*cos(q3t)*sin(sita) - ...
infec2q52q7x2*q3tp^2*q5t^2*q7t*cos(sita)*sin(q3t) +...
4*infec2q72x2*q3tp*q7t*q7tp*cos(q3t)*cos(sita) -...
4*infe2x2*q3tp*q5t*q5tp*cos(q3t)*sin(sita) - ...
4*infe2x2*q3tp*q5t*q5tp*cos(sita)*sin(q3t) + ...
4*infec2q52q7x2*q5t*q5tp*q7tp*cos(q3t)*sin(sita) +...
4*infec2q52q7x2*q5t*q5tp*q7tp*cos(sita)*sin(q3t) -...
4*infec2q72x2*q3tp*q7t*q7tp*sin(q3t)*sin(sita) +...
4*infec2q52q7x2*q3tp*q5t*q5tp*q7t*cos(q3t)*cos(sita) - ...
4*infec2q52q7x2*q3tp*q5t*q5tp*q7t*sin(q3t)*sin(sita);

deDLMMddq3=2*inx2fe2x2*q5tp^2 + g*inx1*cos(q3t) - 2*infe1x1fec1x1*q5t*q5tp^2 - ...
2*infe3x2fec2q72x2*q7t*q7tp^2 - 10*infe2x2fec2q54x2*q5t^4*q5tp^2 - ...
2*infe2x2fec2q72x2*q5t^2*q7tp^2 +...
2*infe2x2fec2q72x2*q5tp^2*q7t^2 - 2*infe3x2fec2q52q7x2*q5tp^2*q7t^2 + ...
2*infe1x1pf*q3tp*q5t*q5tp +...
4*inx1fec1x1*q3tp*q5t*q5tp + 2*infe3x2pf*q3tp*q7t*q7tp +...
4*inx2fec2q72x2*q3tp*q7t*q7tp  +...
4*infe2x2pf*q3tp*q5t^3*q5tp + 4*infec1x1pf*q3tp*q5t^3*q5tp + ...
8*infec2q54x2pf*q3tp*q5t^7*q5tp +...
8*inx2fec2q54x2*q3tp*q5t^3*q5tp + 2*infe2x2fe3x2*q3tp*q5t^2*q7tp +...
2*infec2q52q7x2fec2q54x2*q3tp*q5t^6*q7tp + 2*inx2fec2q52q7x2*q3tp*q5t^2*q7tp +...
4*infec2q72x2pf*q3tp*q7t^3*q7tp - 4*infe2x2fec2q52q7x2*q5t^3*q5tp*q7tp - ...
g*infe1x1*q5t*sin(q3t) + ...
2*L1*infe2x2*q5tp^2*cos(sita) + 2*L1*infec2q72x2*q7tp^2*sin(sita) +...
g*infec1x1*q5t^2*cos(q3t) +...
g*inx2*cos(q3t)*cos(sita) - 12*infe3x2fec2q54x2*q5t^2*q5tp^2*q7t -... 
g*inx2*sin(q3t)*sin(sita) +...
2*L1*infec2q52q7x2*q5tp^2*q7t*sin(sita) + 4*infec2q52q7x2pf*q3tp*q5t^3*q5tp*q7t^2 + ...
8*infec2q54x2fec2q72x2*q3tp*q5t^3*q5tp*q7t^2 + ...
6*infec2q52q7x2fec2q72x2*q3tp*q5t^2*q7t^2*q7tp -...
g*infe3x2*q7t*cos(q3t)*sin(sita) - g*infe3x2*q7t*cos(sita)*sin(q3t) +... 
4*infe2x2fe3x2*q3tp*q5t*q5tp*q7t +...
4*inx2fec2q52q7x2*q3tp*q5t*q5tp*q7t - 4*infe3x2fec2q52q7x2*q5t*q5tp*q7t*q7tp +...
12*L1*infec2q54x2*q5t^2*q5tp^2*sin(sita) + L1*L2*g*roL*cos(q3t) + ...
g*infec2q54x2*q5t^4*cos(q3t)*cos(sita) + g*infec2q72x2*q7t^2*cos(q3t)*cos(sita) -...
g*infe2x2*q5t^2*cos(q3t)*sin(sita) - g*infe2x2*q5t^2*cos(sita)*sin(q3t) - ...
g*infec2q54x2*q5t^4*sin(q3t)*sin(sita) - g*infec2q72x2*q7t^2*sin(q3t)*sin(sita) -...
2*L1*infe3x2*q3tp*q7tp*sin(sita) + 4*infec2q52q7x2fec2q72x2*q3tp*q5t*q5tp*q7t^3 +...
12*infec2q52q7x2fec2q54x2*q3tp*q5t^5*q5tp*q7t +...
2*infec2q52q7x2pf*q3tp*q5t^4*q7t*q7tp +...
4*infec2q54x2fec2q72x2*q3tp*q5t^4*q7t*q7tp +...
g*infec2q52q7x2*q5t^2*q7t*cos(q3t)*cos(sita)  - ...
g*infec2q52q7x2*q5t^2*q7t*sin(q3t)*sin(sita) + ...
4*L1*infec2q72x2*q3tp*q7t*q7tp*cos(sita) - ...
4*L1*infe2x2*q3tp*q5t*q5tp*sin(sita) + 4*L1*infec2q52q7x2*q5t*q5tp*q7tp*sin(sita) + ...
8*L1*infec2q54x2*q3tp*q5t^3*q5tp*cos(sita) + ...
2*L1*infec2q52q7x2*q3tp*q5t^2*q7tp*cos(sita) +...
4*L1*infec2q52q7x2*q3tp*q5t*q5tp*q7t*cos(sita);

deDLMMddq4=infe1x1pppf*q5t + 2*infe2x2pppf*q5t^3 + g*infe1x1*cos(q3t) -... 
infe1x1pf*q3tp^2*q5t - 2*inx1fec1x1*q3tp^2*q5t + 4*infe2x2pf*q5t*q5tp^2 +...
4*infec1x1pf*q5t*q5tp^2 - 2*infe2x2pf*q3tp^2*q5t^3 - 2*infec1x1pf*q3tp^2*q5t^3 - ...
4*infec2q54x2pf*q3tp^2*q5t^7 -  4*inx2fec2q54x2*q3tp^2*q5t^3 +...
48*infec2q54x2pf*q5t^5*q5tp^2 + 8*infec2q54x2fec2q72x2*q5t^3*q7tp^2 +...
2*infe2x2ppfe3x2pp*q5t*q7t - 2*infec2q52q7x2pf*q3tp^2*q5t^3*q7t^2 - ...
4*infec2q54x2fec2q72x2*q3tp^2*q5t^3*q7t^2 - 2*infe2x2fe3x2*q3tp^2*q5t*q7t -...
8*infe3x2fec2q54x2*q3tp*q5t^3*q7tp + 4*infe2x2fec2q52q7x2*q3tp*q5t^3*q7tp - ...
2*inx2fec2q52q7x2*q3tp^2*q5t*q7t + 16*infec2q52q7x2fec2q54x2*q5t^4*q5tp*q7tp + ...
4*infec2q52q7x2fec2q72x2*q5t*q7t*q7tp^2 + 2*g*infec1x1*q5t*sin(q3t) -...
2*infec2q52q7x2fec2q72x2*q3tp^2*q5t*q7t^3 - ...
6*infec2q52q7x2fec2q54x2*q3tp^2*q5t^5*q7t +...
4*infec2q52q7x2pf*q5t*q5tp^2*q7t^2 + 32*infec2q52q7x2fec2q54x2*q5t^3*q5tp^2*q7t +...
2*L1*infe2x2*q3tp^2*q5t*sin(sita) + 2*g*infe2x2*q5t*cos(q3t)*cos(sita) - ...
2*g*infe2x2*q5t*sin(q3t)*sin(sita) - 4*L1*infec2q54x2*q3tp^2*q5t^3*cos(sita) + ...
8*infe2x2fec2q72x2*q3tp*q5t*q7t*q7tp - 4*infe3x2fec2q52q7x2*q3tp*q5t*q7t*q7tp + ...
4*g*infec2q54x2*q5t^3*cos(q3t)*sin(sita) + 4*g*infec2q54x2*q5t^3*cos(sita)*sin(q3t) + ...
8*infec2q52q7x2pf*q5t^2*q5tp*q7t*q7tp - 2*L1*infec2q52q7x2*q3tp^2*q5t*q7t*cos(sita) +...
2*g*infec2q52q7x2*q5t*q7t*cos(q3t)*sin(sita) +... 
2*g*infec2q52q7x2*q5t*q7t*cos(sita)*sin(q3t);

deDLMMddq5=infe3x2pppf*q7t + 2*infe2x2fe3x2*q5tp^2 + infe2x2ppfe3x2pp*q5t^2 -...
infe3x2pf*q3tp^2*q7t -...
2*inx2fec2q72x2*q3tp^2*q7t + 4*infec2q72x2pf*q7t*q7tp^2 - infe2x2fe3x2*q3tp^2*q5t^2 - ...
infec2q52q7x2fec2q54x2*q3tp^2*q5t^6 - inx2fec2q52q7x2*q3tp^2*q5t^2 -...
2*infec2q72x2pf*q3tp^2*q7t^3 +...
12*infec2q52q7x2fec2q54x2*q5t^4*q5tp^2 + 2*infec2q52q7x2fec2q72x2*q5t^2*q7tp^2 +...
4*infec2q52q7x2fec2q72x2*q5tp^2*q7t^2 - 3*infec2q52q7x2fec2q72x2*q3tp^2*q5t^2*q7t^2 +...
8*infe3x2fec2q54x2*q3tp*q5t^3*q5tp - 4*infe2x2fec2q52q7x2*q3tp*q5t^3*q5tp + ...
4*infec2q52q7x2pf*q5t^3*q5tp*q7tp + L1*infe3x2*q3tp^2*sin(sita) + ...
g*infe3x2*cos(q3t)*cos(sita) -...
infec2q52q7x2pf*q3tp^2*q5t^4*q7t - 2*infec2q54x2fec2q72x2*q3tp^2*q5t^4*q7t +...
2*infec2q52q7x2pf*q5t^2*q5tp^2*q7t + 24*infec2q54x2fec2q72x2*q5t^2*q5tp^2*q7t -...
g*infe3x2*sin(q3t)*sin(sita) - 2*L1*infec2q72x2*q3tp^2*q7t*cos(sita) + ...
2*g*infec2q72x2*q7t*cos(q3t)*sin(sita) + 2*g*infec2q72x2*q7t*cos(sita)*sin(q3t) -...
L1*infec2q52q7x2*q3tp^2*q5t^2*cos(sita) - 8*infe2x2fec2q72x2*q3tp*q5t*q5tp*q7t +...
4*infe3x2fec2q52q7x2*q3tp*q5t*q5tp*q7t + 8*infec2q52q7x2fec2q72x2*q5t*q5tp*q7t*q7tp +...
g*infec2q52q7x2*q5t^2*cos(q3t)*sin(sita) + g*infec2q52q7x2*q5t^2*cos(sita)*sin(q3t);



DLsubMMddq=[deDLMMddq1;deDLMMddq2;deDLMMddq3;deDLMMddq4;deDLMMddq5];

fty1=Q-DLsubMMddq;

fty=InvMty22*fty1;

odefun=@(t,y)[y(6); y(7); y(8);y(9);y(10);subs(fty,...
{q1t,q2t,q3t,q5t,q7t,q1tp,q2tp,q3tp,q5tp,q7tp,sigmu},...
{y(1),y(2),y(3),y(4),y(5),y(6),y(7),y(8),y(9),y(10),...
-sign(y(6) - L1*y(8)*sin(y(3)) - fe3L2*y(10)*cos(sita)*sin(y(3)) -...
fe3L2*y(10)*cos(y(3))*sin(sita) -...
L2*y(8)*cos(sita)*sin(y(3)) - L2*y(8)*cos(y(3))*sin(sita) +... 
2*fec2q72L2*y(10)*y(5)*cos(sita)*cos(y(3)) -...
fe3L2*y(5)*y(8)*cos(sita)*cos(y(3)) - 2*fe2L2*y(4)*y(9)*cos(sita)*sin(y(3)) -...
2*fe2L2*y(4)*y(9)*cos(y(3))*sin(sita) -...
2*fec2q72L2*y(10)*y(5)*sin(sita)*sin(y(3)) + fe3L2*y(5)*y(8)*sin(sita)*sin(y(3)) + ...
fec2q52q7L2*y(10)*y(4)^2*cos(sita)*cos(y(3)) -...
fe2L2*y(4)^2*y(8)*cos(sita)*cos(y(3)) + 4*fec2q54L2*y(4)^3*y(9)*cos(sita)*cos(y(3)) -...
fec2q54L2*y(4)^4*y(8)*cos(sita)*sin(y(3))-fec2q54L2*y(4)^4*y(8)*cos(y(3))*sin(sita) -...
fec2q72L2*y(5)^2*y(8)*cos(sita)*sin(y(3)) -fec2q72L2*y(5)^2*y(8)*cos(y(3))*sin(sita) -...  
fec2q52q7L2*y(10)*y(4)^2*sin(sita)*sin(y(3))    +...
fe2L2*y(4)^2*y(8)*sin(sita)*sin(y(3)) - 4*fec2q54L2*y(4)^3*y(9)*sin(sita)*sin(y(3)) -...
2*fec2q52q7L2*y(4)*y(5)*y(9)*sin(sita)*sin(y(3)) -... 
fec2q52q7L2*y(4)^2*y(5)*y(8)*cos(sita)*sin(y(3)) -...
fec2q52q7L2*y(4)^2*y(5)*y(8)*cos(y(3))*sin(sita) + ...
2*fec2q52q7L2*y(4)*y(5)*y(9)*cos(sita)*cos(y(3)))})];


% initial conditions
q1t0=0;
q2t0=sqrt(2)*L1-6*10^(-16);
q3t0=-pi/4;
q5t0=0.0001;
q7t0=-0.0001;

q1tp0=0;
q2tp0=-2;
q3tp0=0;
q5tp0=0.0001;
q7tp0=-0.0001;

y0=[q1t0 q2t0 q3t0 q5t0 q7t0 q1tp0 q2tp0 q3tp0 q5tp0 q7tp0];

tspan=[0 0.01];

options = odeset('Events', @event);

[t,y,te,ye,ie] = ode45(odefun, tspan, y0);

and the event functions as follows:

function [value,isterminal,direction] = event(t,y)

syms x20
sita=3*pi/2;
L1 = 0.2;
L2 = 0.2;

fe20=-1/(15*L1) + x20/(15*L1*L2);

fe30=(-3*x20)/L2 + (6*x20^2)/L2^2 - (4*x20^3)/L2^3;

fec2q54x20=-x20/(450*L1^2*L2^2);

fec2q52q7x20=x20/(5*L1*L2^2) - (2*x20^2)/(5*L1*L2^3) + (4*x20^3)/(15*L1*L2^4);

fec2q72x20=(-9*x20)/(2*L2^2) + (18*x20^2)/L2^3 - (36*x20^3)/L2^4;


fec2q54L200=subs(fec2q54x20,x20,L2);
fec2q52q7L200=subs(fec2q52q7x20,x20,L2);
fec2q72L200=subs(fec2q72x20,x20,L2);
fe2L200=subs(fe20,x20,L2);
fe3L200=subs(fe30,x20,L2);

R0Ey=y(2) + cos (sita)*cos (y(3))*fe2L200*y(4)^2 +... 
 cos (sita)*cos (y(3))*fe3L200*y(5) + L2*cos (y(3))*sin (sita) +...  
 cos (y(3))*fec2q54L200*y(4)^4*sin (sita) +...  
 cos (y(3))*fec2q52q7L200*y(4)^2*y(5)*sin (sita) + ... 
 cos (y(3))*fec2q72L200*y(5)^2*sin (sita) + L1*sin (y(3)) +...  
 L2*cos (sita)*sin (y(3)) + cos (sita)*fec2q54L200*y(4)^4*sin (y(3)) + ... 
 cos (sita)*fec2q52q7L200*y(4)^2*y(5)*sin (y(3)) + ... 
 cos (sita)*fec2q72L200*y(5)^2*sin (y(3)) - ... 
 fe2L200*y(4)^2*sin (sita)*sin (y(3)) - fe3L200*y(5)*sin (sita)*sin (y(3));

value = R0Ey;
isterminal = 1;
direction =1;

It should be noted that the event function is used to terminate the numerical integration.

share|improve this question
When your question is so long that SE forces you to break it up, that's a sign that you're doing it wrong. – David Ketcheson Sep 22 '12 at 6:21
2  
It is unlikely that anyone will run your code as is, especially when it is riddled with code smells. I suggest deleting the code, condensing your two questions into one, and sticking to the mathematics and numerical methods. If you can't describe your problem in a page or less, it's too long. A data dump asking us to ferret out what's wrong won't convince volunteers to answer your question; an adequately-posed question is the cost of convincing people to volunteer their time to help you. If it is difficult to do so, ask for advice on posing questions, and we'll do our best to assist. – Geoff Oxberry♦ Sep 22 '12 at 10:18

closed as not a real question by David Ketcheson, Geoff Oxberry♦ Sep 22 '12 at 10:19

It's difficult to tell what is being asked here. This question is ambiguous, vague, incomplete, overly broad, or rhetorical and cannot be reasonably answered in its current form. For help clarifying this question so that it can be reopened, see the FAQ.

Browse other questions tagged or ask your own question.