# ODE45 producing NaN values for a set of differential equations

3 views (last 30 days)
Mishal Mohanlal on 26 Sep 2021
Commented: Mishal Mohanlal on 28 Sep 2021
Hi Guys
Please refer to the below code.
When I try and excute the function using ODE45 solver all the values produced are NaN.
What is the source of the issue? I have checked the code numerous times and there are no mistakes in equations.
clc;
clear;
r0 = 1;
L = 2.5;
k = 353160*2;
m1 = 3500/2;
m2 = 3500/2;
SolverOptions = odeset('RelTol',1e-5,'AbsTol',1e-5);
[t,y] = ode45(@Sspace,[0:0.001:10],[0 0 0 0 0 0],SolverOptions,r0,L,k,m1,m2);
y
function xp = Sspace(t,x,r0,L,k,m1,m2)
% r2 = x(5);
% r2_dot = x(6);
% theta1 = x(1);
% theta1_dot = x(2);
% theta2 = x(3);
% theta2_dot = x(4);
one = sqrt(L+x(5)*cos(x(1))-x(5)*cos(x(3)));
two = sqrt(L-x(5)*cos(x(1))+x(5)*cos(x(3)));
r1 = (L+one*two+x(5)*sin(x(3)))/sin(x(1));
one_dot = x(6)*cos(x(1))-x(5)*x(2)*sin(x(1))-x(6)*cos(x(3))+x(5)*x(4)*sin(x(3));
two_dot = -x(6)*cos(x(1))+x(6)*x(2)*sin(x(1))+x(6)*cos(x(3))-x(5)*x(4)*sin(x(3));
R1_dot = one_dot*two/(2*one)+two_dot*one/(2*two);
r1_dot = ((R1_dot+x(6)*sin(x(3))+x(5)*x(4)*cos(x(3)))*sin(x(1))-(x(2)*cos(x(1))*(L+one*two+2*x(5)*sin(x(3)))))/sin(x(1))^2;
theta1_dotdot = (-2*m1*r1*r1_dot*x(2)-m1*9.81*r1*sin(x(1)))/(m1*r1^2);
theta2_dotdot = (-2*m2*x(5)*x(6)*x(4)-m2*9.81*x(5)*sin(x(3)))/(m2*x(5)^2);
r2_dotdot = (m2*x(5)*x(4)^2-k*(x(5)-r0)+m2*9.81*cos(x(3)))/m2;
xp = [x(2);
theta2_dotdot;
x(4);
theta2_dotdot;
x(6);
r2_dotdot];
end

Cris LaPierre on 26 Sep 2021
Your numerator and denominators are both 0. And 0/0 = NaN
x=[0 0 0 0 0 0];
r0 = 1;
L = 2.5;
k = 353160*2;
m1 = 3500/2;
m2 = 3500/2;
one = sqrt(L+x(5)*cos(x(1))-x(5)*cos(x(3)));
two = sqrt(L-x(5)*cos(x(1))+x(5)*cos(x(3)));
r1 = (L+one*two+x(5)*sin(x(3)))/sin(x(1));
one_dot = x(6)*cos(x(1))-x(5)*x(2)*sin(x(1))-x(6)*cos(x(3))+x(5)*x(4)*sin(x(3));
two_dot = -x(6)*cos(x(1))+x(6)*x(2)*sin(x(1))+x(6)*cos(x(3))-x(5)*x(4)*sin(x(3));
R1_dot = one_dot*two/(2*one)+two_dot*one/(2*two);
r1_dot = ((R1_dot+x(6)*sin(x(3))+x(5)*x(4)*cos(x(3)))*sin(x(1))-(x(2)*cos(x(1))*(L+one*two+2*x(5)*sin(x(3)))))/sin(x(1))^2
r1_dot = NaN
% numerator
((R1_dot+x(6)*sin(x(3))+x(5)*x(4)*cos(x(3)))*sin(x(1))-(x(2)*cos(x(1))*(L+one*two+2*x(5)*sin(x(3)))))
ans = 0
% denominator
((R1_dot+x(6)*sin(x(3))+x(5)*x(4)*cos(x(3)))*sin(x(1))-(x(2)*cos(x(1))*(L+one*two+2*x(5)*sin(x(3)))))
ans = 0
% result of 0/0
0/0
ans = NaN
Mishal Mohanlal on 28 Sep 2021
Thank you kindly

Paul on 26 Sep 2021
One problem appears to be that the initial condition of x(5) is 0, which results in a divide by zero in the computation of theta2_dotdot.

R2021a

### Community Treasure Hunt

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

Start Hunting!