4
$\begingroup$

If c1 cnd c2 are zero, the code responds very well; but if one of c1 or c2 is nonzero, the code doesn't work correctly. I get the error

NDSolve::ndinid: Initial condition {0} is not in the range specified by the discrete variable NDSolve`s$147246. >>

I think it's because of the Sign function. Without the Sign function, the code works properly. I would like it to work correctly, and also I would like to know why this code didn't work properly.

ClearAll[θ1, θ2];
m1 = 1;  m2 = 1; lc1 = 1;  lc2 = 1; l1 = 2; l2 = 2; Ix2 = 1; Iy2 = 1;  Iz1 = 1;
Iz2 = 1;  g = 9.81;  Im1 = 1; Im2 = 1; η1 = 28;  η2 = 28; b1 = 1; b2 = 2;
c1 = 0; c2 = 1;

τ1[t_] := 10;
τ2[t_] := 10;

u1[t_] := τ1[t] - b1*θ1'[t] - c1*Sign[θ1'[t]];
u2[t_] := τ2[t] - b2*θ2'[t] - c2*Sign[θ2'[t]];


m11[t_] := m1*(lc1)^2 + m2*(l1 + lc2*Cos[θ2[t]])^2 + 
           Ix2*(Sin[θ2[t]])^2 + Iy2*(Cos[θ2[t]])^2 + Iz1;
m22[t_] := m2*(lc2)^2 + Iz2;
V11[t_] := 0;
V12[t_] := (-2*m2*lc2*(l1 + lc2*Cos[θ2[t]])* Sin[θ2[t]] + 
           2 (Ix2 - Iy2)*Sin[θ2[t]]*Cos[θ2[t]])*θ1'[t];
V21[t_] := (m2*lc2*(l1 + lc2*Cos[θ2[t]])* Sin[θ2[t]] - 
           (Ix2 - Iy2)*Sin[θ2[t]]* Cos[θ2[t]])*θ1'[t];
V22[t_] := 0;
M = {{m11[t], 0}, {0, m22[t]}} + {{Im1*(η1)^2, 0}, {0, Im2*(η2)^2}};
V = {{V11[t], V12[t]}, {V21[t], V22[t]}};
G = {{0}, {m2*lc2*Cos[θ2[t]]*g}};
initu = q[0] == {{0}, {.1}};
initv = q'[0] == {{0}, {0}};
q[t_] := {θ1[t], θ2[t]};
sol = First@  NDSolve[{Thread[M.q''[t] + V.q'[t] + G == {u1[t], u2[t]}], initu, 
     initv}, q[t], {t, 0, 100}];

Plot[{θ1[t] /. sol, θ2[t] /. sol}, {t, 0, 100}, PlotTheme -> "Detailed"]
$\endgroup$
1
  • $\begingroup$ This problem is quite tricky, it seems to me. $\endgroup$
    – Michael E2
    Jul 15, 2015 at 0:39

1 Answer 1

5
$\begingroup$

It's a good idea to include the error you got in your question. It gives a clue and might prompt someone to investigate:

NDSolve::ndinid: Initial condition {0} is not in the range specified by the discrete variable NDSolve`s$147246. >>

Now Sign is discontinuous and will cause NDSolve to invoke special processing of the ODE. I suspect that the strange internal variable NDSolve`s$147246 is added to the system along with an event to implement the change in Sign. Note that the message is helpful in calling this variable a discrete variable, which is just what one would need to implement the discontinuous changes in Sign. (Normally one wouldn't expect to see an internal local variable, and I'm still undecided whether you really should. It would have been nicer if the basename were not the single letter s but something like Sign, say.)

In any case, the problem is with the initial condition for q, that is to say, for {θ1, θ2}. It defines θ1 and θ2 to be 1D vectors instead of scalars, so that Sign[θ1'[t]] will be a vector. But apparently NDSolve thought Sign was going to be a scalar. The mismatch causes the error.

If it is acceptable for θ1 and θ2 to be scalars, then here is a fix:

initu = q[0] == Flatten@{{0}, {.1}};  (* <--- N.B. *)
initv = q'[0] == Flatten@{{0}, {0}}; (* <--- N.B. *)
sol = First@NDSolve[{Thread[M.q''[t] + V.q'[t] + G == {u1[t], u2[t]}], 
    initu, initv}, q[t], {t, 0, 100}];

Plot[{θ1[t] /. sol, θ2[t] /. sol}, {t, 0, 100}, PlotTheme -> "Detailed"]

Mathematica graphics

Note: It is not essential, but what I often do, especially when debugging code on this site, is to examine the pieces. The error message complains about the initial conditions. Since I didn't write the code, one of the first things I'm curious about in NDSolve problems is what the ODE system actually looks like. So I evaluate it:

 {Thread[M.q''[t] + V.q'[t] + G == {u1[t], u2[t]}], initu, initv}

Mathematica graphics

The first thing that struck me were those strange (to me) initial conditions for θ1 and θ2. I often do this sort of examination when I've written the code for my own system and it does not work as expected.

$\endgroup$
1
  • $\begingroup$ +1 For this nice way to examine what the ODE actually looks like! $\endgroup$
    – Janosh
    Nov 11, 2017 at 10:32

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.