I have this code, which I think should work fine

Param = {Mo -> 112, mu1 -> 275.52,k -> 4.144,m2 -> 14,m3 -> 24, m4 -> 12, 
         R3 -> 0.16, r3 -> 0.12, rho3 -> 0.14, r1 -> 0.06};
eq = (m3*r1^2 + m3 *(r1*(r3+R3)/(2*r3+R3)) ^2 + m3*rho3*(r1/(2*r3+R3))^2 + 
      m4*(r1*(r3+R3)/(2*r3+R3))^2)*(phi1''[t]*Sin[phi1[t]]^2 - 
Q = Mo - k*phi1'[t] + mu1*phi1'[t]*4*(r1*(r3+R3)/(2*r3+R3));
sol = NDSolve [{eq==Q, phi1[0] == 0, phi1'[0] == 23.375} /.Param,
               {phi1},{t, 0,0.2688}] 

However, it gives me warning

NDSolve::ndnum: Encountered non-numerical value for a derivative at t == 0.`.

There are a lot of similar questions, but none of them were helpful since I don't understand what that error means and what causes it. I believe it has something to do with phi1[t] being inside a sine function.


2 Answers 2


Let see what happens to the ode at $t=0$

In your case

Param = {Mo -> 112, mu1 -> 275.52, k -> 4.144, m2 -> 14, m3 -> 24, 
  m4 -> 12, R3 -> 0.16, r3 -> 0.12, rho3 -> 0.14, r1 -> 0.06}
eq = (m3*r1^2 + m3*(r1*(r3 + R3)/(2*r3 + R3))^2 + 
    m3*rho3*(r1/(2*r3 + R3))^2 + 
    m4*(r1*(r3 + R3)/(2*r3 + R3))^2)*(phi1''[t]*Sin[phi1[t]]^2 - 
Q = Mo - k*phi1'[t] + mu1*phi1'[t]*4*(r1*(r3 + R3)/(2*r3 + R3));
ode = eq == Q 
ode = ode /. Param // Rationalize
ode = First@ode - Last@ode == 0

Mathematica graphics

You can see that at $t=0$ the above reduces to

Mathematica graphics

Because phi1[t]=0 at $t=0$, and hence $\sin(0)=0$ so all of the last term goes away.

Now lets apply the other initial condition which is phi1'[0]== 187/8 and see if it satisfies the above at $t=0$

ode /. phi1[t] -> 0

Mathematica graphics

% /. phi1'[t] -> Rationalize[23.375]

Mathematica graphics

So at $t=0$, the initial conditions are not consistent with the ode.

So I think the problem is that your initial conditions are not consistent. Normally NDSolve will say that as another warning message ibcinc , but may be in this case it did not?

  • $\begingroup$ ibcinc is only for PDE. If we force NDSolve to use "Residual" method for equation simplification, we'll see ivres warning. (See my answer for more info. ) $\endgroup$
    – xzczd
    Jan 12, 2023 at 3:52

The problem is strongly related to e.g.

Numerical solution of Bessel-like equation using NDSolve

NDSolve divide-by-zero trouble

Trouble in the second order ordinary differential equations with second oder coupled iterm

But it may be hard for beginners to notice the underlying relationship, so let me elaborate a bit. As discussed in this post, NDSolve by default transforms the ODE to a first order system in standard form. This standard form can be obtained with to1storder:

internalsys = 
 to1storder[eq == Q /. Param, phi1, t, "standardize" -> True, "form" -> Subscript]

enter image description here

As shown above, the standard form involves a Csc[Subscript[phi1, 0][t]] term. (Subscript[phi1, 0] stands for phi1 in original equation. ) So the initial condition phi1[0] == 0 will lead to non-numerical ComplexInfinity:

internalsys[[1]] /. t -> 0 /. Subscript[phi1, 0][0] -> 0
(* {Subscript[phi1, 1]'[0] == ComplexInfinity} *)

This is the reason why ndnum warning pops up.

The ndnum warning can be circumvented by, as discussed in linked posts above, SolveDelayed -> True:

NDSolve[{eq == Q, phi1[0] == 0, phi1'[0] == 23.375} /. Param, {phi1}, {t, 0, 0.2688},
  SolveDelayed -> True]

But this leads to the warning

NDSolve::ivres: NDSolve has computed initial values that give a zero residual for the differential-algebraic system, but some components are different from those specified. If you need them to be satisfied, giving initial conditions for all dependent variables and their derivatives is recommended.

The reason is already elaborated in Nasser's answer, so I'd like not to repeat. You should double check the underlying model of your equation.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.