Not enough input arguments ode45

function dy = rigid(t,y)
% a column vector
dy = zeros(3,1);    
w1=sin(t)*exp(-t);
w2=sin(2*t)*exp(-2*t);
w3=sin(3*t)*exp(-3*t);
dy(1) =(1/cos(y(2)))* (cos(y(3))*w1-sin(y(3))*w2);
dy(2) =(sin(y(3))*w1+cos(y(3))*w2);
dy(3) =w3+(1/cos(y(2)))*(-w1*cos(y(3))*sin(y(2))+w2*sin(y(2))*sin(y(3)));

syms a b c ad bd cd om1 om2 om3 t real
%rotations for 1(a)-2(b)-3(c)
X=[1 0 0;0 cos(a) sin(a);0 -sin(a) cos(a)];
Y=[cos(b) 0 -sin(b);0 1 0;sin(b) 0 cos(b)];
Z=[cos(c) sin(c) 0;-sin(c) cos(c) 0;0 0 1];
%Last rotation 3(c)
L=Z
%Last 3(c) and middle rotation 2(b)
LM=Z*Y
%Last rotation angle derivative
dv1=[0; 0; cd]
%Middle rotation angle derivative
dv2=[0; bd; 0]
%First rotation angle derivative
dv3=[ad; 0; 0]
w11=dv1;
w22=L*dv2;
w33=LM*dv3;
%summation
w=w11+w22+w33
[S] = equationsToMatrix([w], [ad, bd, cd])

%Convert to differential equation form
%kinematic equation
IS=simplify (inv(S))

%Calculate differential equations
%angular velocity vectors
om1=sin(t)*exp(-t);
om2=sin(2*t)*exp(-2*t);
om3=sin(3*t)*exp(-3*t);
%Differential equation
angd=IS*[om1 om2 om3]'

options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-4 1e-5]);
[T,Y] = ode45(@rigid,[0 2],[0 0 0],options);
plot(T,Y(:,1),'*',T,Y(:,2),'.',T,Y(:,3),'+')
legend('a','b','c')

[Aleyna]

Question: Calculate Euler angles differential equations for 1(a)-2(b)-3(c) rotation sequence

Answer 1

[Dyuman Joshi]

I am not sure how the different symbolic variables and expressions come into play here for solving the ODE.

Nonetheless, Any functions defined in a script must be at the end of the script.

Also, the call to the ODE solver (ode45 in this case) should be outside the ODE funciton.

syms a b c ad bd cd om1 om2 om3 t real
%rotations for 1(a)-2(b)-3(c)
X=[1 0 0;0 cos(a) sin(a);0 -sin(a) cos(a)];
Y=[cos(b) 0 -sin(b);0 1 0;sin(b) 0 cos(b)];
Z=[cos(c) sin(c) 0;-sin(c) cos(c) 0;0 0 1];
%Last rotation 3(c)
L=Z;
%Last 3(c) and middle rotation 2(b)
LM=Z*Y;
%Last rotation angle derivative
dv1=[0; 0; cd];
%Middle rotation angle derivative
dv2=[0; bd; 0];
%First rotation angle derivative
dv3=[ad; 0; 0];
w11=dv1;
w22=L*dv2;
w33=LM*dv3;
%summation
w=w11+w22+w33;
[S] = equationsToMatrix([w], [ad, bd, cd]);

%Convert to differential equation form
%kinematic equation
IS=simplify (inv(S));

%Calculate differential equations
%angular velocity vectors
om1=sin(t)*exp(-t);
om2=sin(2*t)*exp(-2*t);
om3=sin(3*t)*exp(-3*t);
%Differential equation
angd=IS*[om1 om2 om3]';


%% Call to ODE
options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-4 1e-5]);
[T,Y] = ode45(@rigid,[0 2],[0 0 0],options);

%% Plotting the ODE
plot(T,Y(:,1),'*',T,Y(:,2),'.',T,Y(:,3),'+')
legend('a','b','c')

%% Move the ODE function to the bottom/end of the script

function dy = rigid(t,y)
% a column vector
dy = zeros(3,1);    
w1=sin(t)*exp(-t);
w2=sin(2*t)*exp(-2*t);
w3=sin(3*t)*exp(-3*t);
dy(1) =(1/cos(y(2)))* (cos(y(3))*w1-sin(y(3))*w2);
dy(2) =(sin(y(3))*w1+cos(y(3))*w2);
dy(3) =w3+(1/cos(y(2)))*(-w1*cos(y(3))*sin(y(2))+w2*sin(y(2))*sin(y(3)));
end

Answer 2

[Sulaymon Eshkabilov]

Don't use the MATLAB's built in function names as a variable name - cd. Why to display all iterations that can be also skipped with ";". Here is the corrected code:

options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-4 1e-5]);
[T,Y] = ode45(@rigid,[0 2],[0 0 0],options);

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

angd = 

plot(T,Y(:,1),'*',T,Y(:,2),'.',T,Y(:,3),'+')
legend('a','b','c')

function dy = rigid(t,y)
% a column vector
dy = zeros(3,1);    
w1=sin(t)*exp(-t);
w2=sin(2*t)*exp(-2*t);
w3=sin(3*t)*exp(-3*t);
dy(1) =(1/cos(y(2)))* (cos(y(3))*w1-sin(y(3))*w2);
dy(2) =(sin(y(3))*w1+cos(y(3))*w2);
dy(3) =w3+(1/cos(y(2)))*(-w1*cos(y(3))*sin(y(2))+w2*sin(y(2))*sin(y(3)));

syms a b c ad bd cdd  t real  % om1 om2 om3 should be skipped
%rotations for 1(a)-2(b)-3(c)
X=[1 0 0;0 cos(a) sin(a);0 -sin(a) cos(a)];
Y=[cos(b) 0 -sin(b);0 1 0;sin(b) 0 cos(b)];
Z=[cos(c) sin(c) 0;-sin(c) cos(c) 0;0 0 1];
%Last rotation 3(c)
L=Z;
%Last 3(c) and middle rotation 2(b)
LM=Z*Y;
%Last rotation angle derivative
dv1=[0; 0; cdd];
%Middle rotation angle derivative
dv2=[0; bd; 0];
%First rotation angle derivative
dv3=[ad; 0; 0];
w11=dv1;
w22=L*dv2;
w33=LM*dv3;
%summation
w=w11+w22+w33;
[S] = equationsToMatrix([w], [ad, bd, cdd]);

%Convert to differential equation form
%kinematic equation
IS=simplify (inv(S));

%Calculate differential equations
%angular velocity vectors
om1=sin(t)*exp(-t);
om2=sin(2*t)*exp(-2*t);
om3=sin(3*t)*exp(-3*t);
%Differential equation
angd=IS*[om1 om2 om3]'
end