eq1 = ((A-B+C-D*x6d)/(-F))-x2d;
eq2 = ((G+H*I+L*x6d)/(-K))-x4d;
eq3 = ((M-N*O-P*x4d)/(-R))-x6d;
[sol1, sol2, sol3] = solve([eq1, eq2, eq3],[x2d,x4d,x6d]);
X0 = [1 0 1 0 1 0];
t = 0:60:5676.98*1.6;
[time_lin, solution_lin] = ode45(@(t,x) Lin_func(t,x), t, X0);
[time_nonlin, solution_nonlin] = ode45(@(t,x) nonLin_func(t,x), t, X0);
function d_dt = Lin_func(t,X)
Ix = 6;
Iy = 8;
Iz = 4;
mu = (3.986*10^14);
Rad = 6878140;
w = sqrt(mu/Rad^3);
d_dt = X*NaN;
d_dt(1) = X(2);
d_dt(2) = (w*(Ix+Iz-Iy)*X(6)-4*(w^2)*(Iy-Iz)*X(1))/Ix;
d_dt(3) = X(4);
d_dt(4) = 3*(w^2)*(Iz-Ix)*X(3)/Iy;
d_dt(5) = X(6);
d_dt(6) = ((w^2)*(Ix-Iy)*X(5)-w*(Iz+Ix-Iy)*X(2))/Iz;
end
function d_dt = nonLin_func(t,X)
Ix = 6;
Iy = 8;
Iz = 4;
mu = (3.986*10^14);
Rad = 6878140;
w = sqrt(mu/Rad^3);
A = (Iy-Iz)*(w*cosd(X(5))*sind(X(1))-w*cosd(X(1))*sind(X(3))*sind(X(5))-sind(X(1))*X(4)+cosd(X(3))*cosd(X(1))*X(6))*(w*cosd(X(1))*cosd(X(5))+w*sind(X(3))*sind(X(1))*sind(X(5))-cosd(X(1))*X(4)-cosd(X(3))*sind(X(1))*X(6));
B = Ix*w*cosd(X(3))*cosd(X(5))*X(6);
C = Ix*X(4)*(w*sind(X(3))*sind(X(5))-cosd(X(3))*X(6));
D = Ix*sind(X(3));
E = (3*mu/(2*Rad^3))*(sind(2*X(1))*(cosd(X(3))^2)*(Iy-Iz));
F = Ix;
G = (Iz-Ix)*(w*cosd(X(5))*sind(X(1))-w*cosd(X(1))*sind(X(3))*sind(X(5))-sind(X(1))*X(4)+cosd(X(3))*cosd(X(1))*X(6))*(w*cosd(X(3))*sind(X(5))-X(2)+sind(X(3))*X(6));
H = Iy;
I = (-w*cosd(X(5))*sind(X(3))*sind(X(1))*X(6)+w*cosd(X(1))*sind(X(5))*X(6)+X(2)*(w*cosd(X(5))*sind(X(1))-w*cosd(X(1))*sind(X(3))*sind(X(5))+cosd(X(3))*cosd(X(1))*X(6))-sind(X(1))*X(4)*(w*cosd(X(3))*sind(X(5))+X(2)+sind(X(3))*X(6)));
J = (3*mu/(2*Rad^3))*(sind(2*X(3))*cosd(X(1))*(Ix-Iz));
K = cosd(X(1))*Iy;
L = Iy*cosd(X(3))*sind(X(1));
M = (Iy-Ix)*(-w*cosd(X(3))*sind(X(5))+X(2)-sind(X(3))*X(6))*(-w*(cosd(X(1))*cosd(X(5))+sind(X(3))*sind(X(1))*sind(X(5)))+cosd(X(1))*X(4)+cosd(X(3))*sind(X(1))*X(6));
N = Iz;
O = (w*cosd(X(1))*cosd(X(5))*sind(X(3))*X(6)+w*sind(X(1))*sind(X(5))*X(6)+cosd(X(1))*X(4)*(w*cosd(X(3))*sind(X(5))+X(2)+sind(X(3))*X(6))-X(2)*(w*cosd(X(1))*cosd(X(5))+w*sind(X(3))*sind(X(1))*sind(X(5))-cosd(X(3))*sind(X(1))*X(6)));
P = sind(X(1))*Iz;
Q = (3*mu/(2*Rad^3))*(sind(2*X(3))*sind(X(1))*(Iy-Ix));
R = cosd(X(3))*cosd(X(1))*Iz;
d_dt = X*NaN;
d_dt(1) = X(2);
d_dt(2) = -(D*G*P + D*K*M + A*L*P + A*K*R - B*L*P - D*J*P - B*K*R + C*L*P + C*K*R - D*K*Q - E*L*P - E*K*R + D*H*I*P - D*K*N*O)/(F*(L*P + K*R));
d_dt(3) = X(4);
d_dt(4) = -(G*R - L*M - J*R + L*Q + H*I*R + L*N*O)/(L*P + K*R);
d_dt(5) = X(6);
d_dt(6) = -(G*P + K*M - J*P - K*Q + H*I*P - K*N*O)/(L*P + K*R);
end
2 Comments
Direct link to this comment
https://fr.mathworks.com/matlabcentral/answers/495624-solving-nonlinear-equation-and-linearized-eqaution-using-ode45#comment_776060
Direct link to this comment
https://fr.mathworks.com/matlabcentral/answers/495624-solving-nonlinear-equation-and-linearized-eqaution-using-ode45#comment_776060
Direct link to this comment
https://fr.mathworks.com/matlabcentral/answers/495624-solving-nonlinear-equation-and-linearized-eqaution-using-ode45#comment_776183
Direct link to this comment
https://fr.mathworks.com/matlabcentral/answers/495624-solving-nonlinear-equation-and-linearized-eqaution-using-ode45#comment_776183
Sign in to comment.