I wonder why NaN appears at Y. How can I fix it?
1 view (last 30 days)
Show older comments
Student
on 1 Oct 2023
Commented: Walter Roberson
on 1 Oct 2023
I was solving this differential equation, but NaN came out. :(
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1498929/image.png)
syms Z(R)
g = 9.80665;
%b = ;
r = 0.072;
la = 997;
lb = 1.29;
%k = (la - lb) * g * b^2 / r;
k = 20;
b = sqrt(k * r / (g * (la - lb)));
ode = diff(Z, 2) / (1 + (diff(Z)^2))^1.5 + diff(Z) / (R * (1 + (diff(Z)^2))^2) == 2 + k * Z;
[V] = odeToVectorField(ode);
M = matlabFunction(V, 'vars', {'R', 'Y'});
[R, Y]= ode45(M,[0, 3],[0, 0]);
plot(R, Y);
0 Comments
Accepted Answer
Walter Roberson
on 1 Oct 2023
Your ode is a function of R, and divides by R, and you started R at 0.
1 Comment
Walter Roberson
on 1 Oct 2023
Look at the second term of the expression. You have defined boundary values at R = 0 with Z'(0)=0. The numerator of the second term is Z' and if we evaluate at R=0 then by assumption in the initial conditions, that is 0 so the numerator would be 0.
Meanwhile the denominator of the second term involves R. At R=0 that would be 0. The other part of the denominator is sqrt(1+Z'^2) and at 0 by way of the initial conditions Z' is 0 there so sqrt(1+0^2) = 1. So the denominator is 0*1
Therefore at R=0 the second term is 0/0 which gives problems. Numeric solvers cannot solve that problem.
Can theory solve the problem? Well, if we understand the equation in terms of limits then if the Z' of the numerator approaches 0 "faster" than the linear R approaches 0 then in the limit the second term is potentially "fast 0 divided by slow 0" which in calculus is 0.
If we temporarily substitute 0 for the second term and look at the first term at 0 with Z'(0)=0 then the first term simplifies to Z''(0) and the right hand side simplifies to 2+k*0 so Z''(0) would have to equal 2. So start exploring equations in which Z''(0)=2, Z' (0)=0, Z(0)=0. Can you make the simplest such equation work with the ode?
More Answers (0)
See Also
Categories
Find more on Ordinary Differential Equations in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!