Changing the equation using WhenEvent in NDSolve

Question

I want to change part of a set of equations by using WhenEvent in NDSolve. I have made a much less complicated minimum working example as follows:

sol = x /. First@NDSolve[{
      x''[t] +
        If[x[t] < 0,
         50 x[t] ,
         16 50 x[t] ] == 0,
      x'[0] == 0,
      x[0] == 1
      },
     {x}, {t, 0, 2}];
Plot[sol[t], {t, 0, 1}]

This does not use WhenEvent but does switch a term in the equation. Now I try and introduce a WhenEvent. I also introduce a module that contains the two versions of the term in the equation.

ClearAll[ff];
ff[status_, t_?NumericQ, x_?NumericQ] := 
  If[x < 0,
   Sow[{t, x, -1, status}, tag1]; 50 x ,
   Sow[{t, x,   1, status}, tag1]; 16 50 x ];

status = "More";
a = Reap[
   res = Reap[
      sol = x /. First@NDSolve[{
            x''[t] +
              ff[status, t, x[t]] == 0,
            x'[0] == 0, x[0] == 1,

            WhenEvent[x[t] == 0, 
             If[x[t] < 0, status = "Less", status = "More"]; 
             Sow[{t, x[t], status, "WhenEvent"}, tag2]]},
           {x}, {t, 0, 2}];,
      tag1];,
   tag2];
Plot[sol[t], {t, 0, 1}, 
 Epilog -> {Line[res[[2, 1, All, {1, 3}]]], Red, PointSize[0.01], 
   Point[a[[2, 1, All, {1, 2}]]]}]

This has also worked and the two Sow functions collect various values which shows everything is working. I have introduced a status value which shows which of the two states I am in. I would like to use the status value for controlling the switching of the term in the differential equation. So I modify the function ff so that the status is examined rather than the value of the variable. This does not work.

ClearAll[ff];
ff[status_, t_?NumericQ, x_?NumericQ] := 
  If[status == "Less",
   Sow[{t, x, -1, status}, tag1]; 50 x ,
   Sow[{t, x,   1, status}, tag1]; 16 50 x ];

status = "More";
a = Reap[
   res = Reap[
      sol = x /. First@NDSolve[{
            x''[t] +
              ff[status, t, x[t]] == 0,
            x'[0] == 0, x[0] == 1,

            WhenEvent[x[t] == 0, 
             If[x[t] < 0, status = "Less", status = "More"]; 
             Sow[{t, x[t], status, "WhenEvent"}, tag2]]},
           {x}, {t, 0, 2}];,
      tag1];,
   tag2];
Plot[sol[t], {t, 0, 1}, 
 Epilog -> {Line[res[[2, 1, All, {1, 3}]]], Red, PointSize[0.01], 
   Point[a[[2, 1, All, {1, 2}]]]}]

Why does this not work? How can I fix it?

Answer 1

I think you need to make status a DiscreteVariable in order to change it within NDSolve. In that case, you'll also need to use numeric values instead of "More" and "Less".

Does this work for you?

ClearAll[ff];
ff[status_, t_?NumericQ, x_?NumericQ] := 
  If[status == -1, Sow[{t, x, -1, status}, tag1]; 50 x, 
   Sow[{t, x, 1, status}, tag1]; 16 50 x];

Clear[status];
a = Reap[res = 
     Reap[sol = 
        x /. First@
          NDSolve[{x''[t] + ff[status[t], t, x[t]] == 0, x'[0] == 0, 
            x[0] == 1, status[0] == 1,            
            WhenEvent[
             x[t] < 0, {status[t] -> -1, 
              Sow[{t, x[t], status[t], "WhenEvent"}, tag2]}],
            WhenEvent[
             x[t] > 0, {status[t] -> 1, 
              Sow[{t, x[t], status[t], "WhenEvent"}, tag2]}]}, {x, 
            status}, {t, 0, 2}, DiscreteVariables -> {status}];, 
      tag1];, tag2];

Plot[sol[t], {t, 0, 1}, 
 Epilog -> {Line[res[[2, 1, All, {1, 3}]]], Red, PointSize[0.01], 
   Point[a[[2, 1, All, {1, 2}]]]}]