Clear Filters
Clear Filters

Solving system of nonlinear differential equations using ode45

2 views (last 30 days)
I have a question for system of ordinary differential equations, because Matlab gives some strange solution as output. There is the code:
Jo1=1;
Jo2=2;
Jo3=3;
Mo1=1;
Mo2=1;
Mo3=1;
f=@(t,x)[x(4).*sin(x(3))./sin(x(2))+x(5).*cos(x(3))./sin(x(2));
cos(x(3)).*x(4)-sin(x(3)).*x(5);
-x(4).*sin(x(3)).*cos(x(2))./sin(x(2))-x(6).*cos(x(3)).*cos(x(2))./sin(x(2));
(Mo1-(Jo3-Jo2).*x(6).*x(5))./Jo1;
(Mo2-(Jo1-Jo3).*x(4).*x(6))./Jo2;
(Mo3-(Jo2-Jo1).*x(5).*x(4))./Jo3];
[t, x]= ode45(f, [0,1],[0,0,0,0,0,0]);
I get solution for x4, x5 and x6, but for x1, x2 and x3 the solution is NaN.
So if someone have any advice it would help.

Answers (3)

Francisco J. Triveno Vargas
Edited: Torsten on 12 Jul 2024
A simples tests is change the initial values like this:
clear
close all
Jo1=1;
Jo2=2;
Jo3=3;
Mo1=1;
Mo2=1;
Mo3=1;
f=@(t,x)[x(4).*sin(x(3))./sin(x(2))+x(5).*cos(x(3))./sin(x(2));
cos(x(3)).*x(4)-sin(x(3)).*x(5);
-x(4).*sin(x(3)).*cos(x(2))./sin(x(2))-x(6).*cos(x(3)).*cos(x(2))./sin(x(2));
(Mo1-(Jo3-Jo2).*x(6).*x(5))./Jo1;
(Mo2-(Jo1-Jo3).*x(4).*x(6))./Jo2;
(Mo3-(Jo2-Jo1).*x(5).*x(4))./Jo3];
[t, x]= ode45(f, [0,1],[0.01,0.01,0.01,0,0,0]);
figure(10)
plot(t,x)

Sam Chak
Sam Chak on 11 Jul 2024
Edited: Sam Chak on 11 Jul 2024
It appears that your current formulations may be erroneous. If you are not comfortable memorizing the required formulas, I would suggest referring to established textbooks to copy the most commonly used and accepted expressions.
Additionally, since you are applying constant torques of Mo1 = 1, Mo2 = 1, and Mo3 = 1 to the system, the system is not going to remain in the desired equilibrium point of [0, 0, 0, 0, 0, 0] (also the initial values). This may indicate the need to further review your system dynamics and control strategies.
Rotational kinematics in terms of a 3-2-1 Euler rotation sequence:
Here is a simple demo:
%% The System
function dx = ode(t, x)
% the parameters
Jo1 = 1;
Jo2 = 2;
Jo3 = 3;
% the inputs (need to be designed)
Mo1 = - 1*x(1) - 2*x(4);
Mo2 = - 1*x(2) - 2*x(5);
Mo3 = - 1*x(3) - 2*x(6);
% the differential equations
dx = [x(4) + sin(x(1))*tan(x(2))*x(5) + cos(x(1))*tan(x(2))*x(6)
cos(x(1))* x(5) - sin(x(1))* x(6)
sin(x(1))*sec(x(2))*x(5) + cos(x(1))*sec(x(2))*x(6)
(Mo1 + (Jo2 - Jo3)*x(2)*x(3))/Jo1
(Mo2 + (Jo3 - Jo1)*x(3)*x(1))/Jo2
(Mo3 + (Jo1 - Jo2)*x(1)*x(2))/Jo3];
end
%% Run the simulation
tspan = [0, 15];
x0 = deg2rad([3, 6, 9, 0, 0, 0]);
[t, x] = ode45(@ode, tspan, x0);
plot(t, rad2deg(x(:,1:3))), grid on, xlabel('t'), ylabel('Angle (degree)')
title('Time responses of Euler angles')
legend('\theta_{1}', '\theta_{2}', '\theta_{3}', 'fontsize', 14)
  2 Comments
Miroslav Jovic
Miroslav Jovic on 17 Jul 2024
I have tried for constant values at first, and now i've calculate real torques and put it in code. Thank you for your help!
Sam Chak
Sam Chak on 17 Jul 2024
Hi @Miroslav Jovic, I'm glad to hear that. If you found the example and code helpful, please consider clicking 'Accept' ✔️ on the answer. Additionally, you can show your appreciation by voting 👍 other answers as a token of support for knowledge sharing. Your support is greatly appreciated!

Sign in to comment.


Muhammed abdulmalek
Muhammed abdulmalek on 11 Jul 2024
Jo1 = 1;
Jo2 = 2;
Jo3 = 3;
Mo1 = 1;
Mo2 = 1;
Mo3 = 1;
f = @(t,x) [x(4) + sin(x(1)*tan(x(2))*x(5) + cos(x(1))*tan(x(2))*x(6));
cos(x(1)*x(5) - sin(x(1))*x(6));
sin(x(1)*sn(x(2))*x(5) + cos(x(1))*sn(x(2))*x(6));
(Mo1 - (Jo3 - Jo2)*x(3)*x(2))/Jo1;
(Mo2 - (Jo1 - Jo3)*x(1)*x(3))/Jo2;
(Mo3 - (Jo2 - Jo1)*x(2)*x(1))/Jo3];
[t, x]= ode45(f, [0, 1], [0, 0, 0, 0, 0, 0]);
plot(t, x), ızgara açık, xlabel('t')

Categories

Find more on Particle & Nuclear Physics in Help Center and File Exchange

Products


Release

R2015a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!