am trying to plot the Euler parameters constraint equation vs revolutions and plot nutation, precession, and spin angles vs time. I am unable to run this code with the [q1; 0; 0; q4;0;1;0] initial conditions. I cannot run the code with initial condition [q1; 0; 0; q4;0;1;0] or [q1; 0; 0; q4;0;-1;0]. But it works fine for any higher whole number such as [q1; 0; 0; q4;0;2;0]. Why is that? How can I fix this?     clc; clear all; q1 = sin((6/2)*pi/180); q4 = cos((6/2)*pi/180); t = 0:0.01:(2*pi*4); options=odeset('RelTol',1e-12,'AbsTol',1e-12); [t,y] = ode45(@euler_eqns, t, [q1; 0; 0; q4;0;1;0], options); q = y(:,1:4) w = y(:,5:7) K = q(:,1).^2 + q(:,2).^2 + q(:,3).^2 + q(:,4).^2; K0 = 1; c = K - K0; v = t / (2*pi); plot(v,c); xlabel('revs '); ylabel('K-K0'); nutation = acos(1 - 2*q(:,1).^2 - 2*q(:,3).^2); precession = asin((2*(q(:,1).*q(:,2) - q(:,3).*q(:,4)))./(sin(nutation))); spin = asin((2*(q(:,1).*q(:,2) + q(:,3).*q(:,4)))./(sin(nutation))); psi = rad2deg(precession); theta = rad2deg(nutation); phi = rad2deg(spin); % Plot the angles over time plot(v, psi, v, theta, v, phi) legend('Precession', 'Nutation', 'Spin') xlabel('revs') ylabel('Angles (deg)')   function dqwdt = euler_eqns(t,qwd, q)    q = qwd(1:4);    w = qwd(5:7);    I = 500;    J = 125;    K = (I-J)/I;    T = 35;    ohm = 1;    qstar = 1/ohm;        dqdt = zeros(4,1);    dqdt(1) = pi*(q(4)*w(1) + q(3)*(qstar - w(2) - ohm) + q(2)*w(3));    dqdt(2) = pi*(q(3)*w(1) + q(4)*(w(2) - ohm - qstar) - q(1)*w(3));    dqdt(3) = pi*(-q(2)*w(1) - q(1)*(qstar - w(2) - ohm) + q(4)*w(3));    dqdt(4) = pi*(-q(1)*w(1) - q(2)*(w(2) - ohm - qstar) - q(3)*w(3));   dwdt = zeros(3,1);    dwdt(1) = 12*pi*K*(q(2)*q(3) + q(1)*q(4))*(1 - 2*q(1)^2 - 2*q(2)^2) ...               - 2*pi*K*w(2)*w(3) + 2*pi*qstar*w(3);    dwdt(2) = 0;    dwdt(3) = -24*pi*K*(q(3)*q(1) - q(2)*q(4))*(q(2)*q(3) + q(1)*q(4)) ...               + 2*pi*K*w(1)*w(2) + 2*pi*qstar*w(1);      % Combine the time derivatives into a single vector    dqwdt = [dqdt; dwdt];     end

Computer Networking: A Top-Down Approach (7th Edition)
7th Edition
ISBN:9780133594140
Author:James Kurose, Keith Ross
Publisher:James Kurose, Keith Ross
Chapter1: Computer Networks And The Internet
Section: Chapter Questions
Problem R1RQ: What is the difference between a host and an end system? List several different types of end...
icon
Related questions
Question

I am trying to plot the Euler parameters constraint equation vs revolutions and plot nutation, precession, and spin angles vs time. I am unable to run this code with the [q1; 0; 0; q4;0;1;0] initial conditions. I cannot run the code with initial condition [q1; 0; 0; q4;0;1;0] or [q1; 0; 0; q4;0;-1;0]. But it works fine for any higher whole number such as [q1; 0; 0; q4;0;2;0]. Why is that? How can I fix this?

 

 

clc;
clear all;

q1 = sin((6/2)*pi/180);
q4 = cos((6/2)*pi/180);
t = 0:0.01:(2*pi*4);
options=odeset('RelTol',1e-12,'AbsTol',1e-12);
[t,y] = ode45(@euler_eqns, t, [q1; 0; 0; q4;0;1;0], options);
q = y(:,1:4)
w = y(:,5:7)

K = q(:,1).^2 + q(:,2).^2 + q(:,3).^2 + q(:,4).^2;
K0 = 1;
c = K - K0;
v = t / (2*pi);

plot(v,c);
xlabel('revs ');
ylabel('K-K0');


nutation = acos(1 - 2*q(:,1).^2 - 2*q(:,3).^2);
precession = asin((2*(q(:,1).*q(:,2) - q(:,3).*q(:,4)))./(sin(nutation)));
spin = asin((2*(q(:,1).*q(:,2) + q(:,3).*q(:,4)))./(sin(nutation)));

psi = rad2deg(precession);
theta = rad2deg(nutation);
phi = rad2deg(spin);

% Plot the angles over time
plot(v, psi, v, theta, v, phi)
legend('Precession', 'Nutation', 'Spin')
xlabel('revs')
ylabel('Angles (deg)')

 

function dqwdt = euler_eqns(t,qwd, q)
   q = qwd(1:4);
   w = qwd(5:7);
   I = 500;
   J = 125;
   K = (I-J)/I;
   T = 35;
   ohm = 1;
   qstar = 1/ohm;
   
   dqdt = zeros(4,1);
   dqdt(1) = pi*(q(4)*w(1) + q(3)*(qstar - w(2) - ohm) + q(2)*w(3));
   dqdt(2) = pi*(q(3)*w(1) + q(4)*(w(2) - ohm - qstar) - q(1)*w(3));
   dqdt(3) = pi*(-q(2)*w(1) - q(1)*(qstar - w(2) - ohm) + q(4)*w(3));
   dqdt(4) = pi*(-q(1)*w(1) - q(2)*(w(2) - ohm - qstar) - q(3)*w(3));

  dwdt = zeros(3,1);
   dwdt(1) = 12*pi*K*(q(2)*q(3) + q(1)*q(4))*(1 - 2*q(1)^2 - 2*q(2)^2) ...
              - 2*pi*K*w(2)*w(3) + 2*pi*qstar*w(3);
   dwdt(2) = 0;
   dwdt(3) = -24*pi*K*(q(3)*q(1) - q(2)*q(4))*(q(2)*q(3) + q(1)*q(4)) ...
              + 2*pi*K*w(1)*w(2) + 2*pi*qstar*w(1);

 

   % Combine the time derivatives into a single vector
   dqwdt = [dqdt; dwdt];
   
end

 

Expert Solution
trending now

Trending now

This is a popular solution!

steps

Step by step

Solved in 3 steps

Blurred answer
Recommended textbooks for you
Computer Networking: A Top-Down Approach (7th Edi…
Computer Networking: A Top-Down Approach (7th Edi…
Computer Engineering
ISBN:
9780133594140
Author:
James Kurose, Keith Ross
Publisher:
PEARSON
Computer Organization and Design MIPS Edition, Fi…
Computer Organization and Design MIPS Edition, Fi…
Computer Engineering
ISBN:
9780124077263
Author:
David A. Patterson, John L. Hennessy
Publisher:
Elsevier Science
Network+ Guide to Networks (MindTap Course List)
Network+ Guide to Networks (MindTap Course List)
Computer Engineering
ISBN:
9781337569330
Author:
Jill West, Tamara Dean, Jean Andrews
Publisher:
Cengage Learning
Concepts of Database Management
Concepts of Database Management
Computer Engineering
ISBN:
9781337093422
Author:
Joy L. Starks, Philip J. Pratt, Mary Z. Last
Publisher:
Cengage Learning
Prelude to Programming
Prelude to Programming
Computer Engineering
ISBN:
9780133750423
Author:
VENIT, Stewart
Publisher:
Pearson Education
Sc Business Data Communications and Networking, T…
Sc Business Data Communications and Networking, T…
Computer Engineering
ISBN:
9781119368830
Author:
FITZGERALD
Publisher:
WILEY