Q:

Switching Differential Equation in NDSolve

I am trying to solve the following system of differential equations using NDSolve:
$\dot{z}_t=.5(1-z_t)$
$\dot{y}_t=.05y_t+z_t-x_t$
subject to the following constraint:
$-y_t-z_t\le0$
where $z_t$ essentially follows a predetermined path, and $x_t=0.75$ until the constraint becomes binding, then dynamically adjusts to keep $y_t$ at the level constrained by $z_t$.
I am using the following code:
zdot=.5*(1-z[t]);
ydot=.05*y[t]+z[t]-x[t];
xdotbind=D[Solve[-ydot-zdot==0,x[t]][[1,1,2]],t];
xdot=Piecewise[{{0,bind==0},{xdotbind,bind==1}}];
bind=0;
sol=NDSolve[{x'[t]==xdot,y'[t]==ydot,z'[t]==zdot,x[0]==.75,y[0]==0,z[0]==0.1,WhenEvent[-y[t]-z[t]==0,{x[t]->.05*y[t]+z[t],bind=1,"RemoveEvent"}]},{x,y,z},{t,0,25}];
Grid[{{Plot[Evaluate[x[t]/.sol],{t,0,2.5},PlotRange->All,AxesLabel->{"t","x"}],Plot[{Evaluate[y[t]/.sol],-Evaluate[z[t]/.sol]},{t,0,2.5},PlotRange->All,AxesLabel->{"t","y"},PlotStyle->{Automatic,{Gray,Dashed}}],Plot[Evaluate[z[t]/.sol],{t,0,25},PlotRange->All,AxesLabel->{"t","z"}]}}]

which generates the following:

Basically, I am using a WhenEvent trigger to change the value of $x_t$ once $y_t$ hits the constraint (represented by the dotted line).  This works fine, but then I try to force NDSolve to switch the expression used for $\dot{x}_t$ to the one that would make $y_t$ follow the constraint path going forward by using a Piecewise function, but this doesn't seem to work for some reason.
Alternatively, I've been able to achieve what I want using "StopIntegration" in a WhenEvent trigger to shut down NDSolve when it hits the constraint, then using the values of the state variables at this point as the initial values in a second NDSolve with the appropriate constraint expression for $\dot{x}_t$.  However, I then have to stitch together the two sets of results from the two NDSolve commands, which is a bit cumbersome, so I'm hoping there's a way to do this within a single NDSolve.
Any ideas?

A:

Changing parameter values during integration works better with DiscreteVariables.
But I think the problem with OP's code, in the question and the OP's answer, has more to do with Mathematica numerics.
My solution
Clear[bind];
zdot = 1/2 (1 - z[t]);
ydot = 1/20*y[t] + z[t] - x[t];
xdotbind = 
  D[Solve[-ydot - zdot == 0, x[t]][[1, 1, 2]], t] /. {y'[t] -> ydot, 
    z'[t] -> zdot};
xdot = Piecewise[{{0, bind[t] == 0}, {xdotbind, bind[t] == 1}}];
{sol, {data}} = Reap@NDSolve[{
     x'[t] == xdot, y'[t] == ydot, z'[t] == zdot,
     bind[0] == 0, x[0] == 3/4, y[0] == 0, z[0] == 1/10,
     WhenEvent[-y[t] - z[t] == 0,
      {x[t] -> 1/20*y[t] + z[t] + 1/2*(1 - z[t]), bind[t] -> 1, 
       "RemoveEvent"}]},
    {x, y, z}, {t, 0, 5},
    DiscreteVariables -> {bind \[Element] {0, 1}},
    StepMonitor :> 
     Sow[{bind[t], t, {x[t], y[t], z[t]}, {xdot, ydot, -zdot}}]];

Grid[{{Plot[Evaluate[x[t] /. sol], {t, 0, 5}, PlotRange -> All, 
    AxesLabel -> {"t", "x"}], 
   Plot[{Evaluate[y[t] /. sol], Evaluate[-z[t] /. sol]}, {t, 0, 5}, 
    PlotRange -> All, AxesLabel -> {"t", "y"}, 
    PlotStyle -> {Automatic, {Gray, Dashed}}], 
   Plot[Evaluate[z[t] /. sol], {t, 0, 5}, PlotRange -> All, 
    AxesLabel -> {"t", "z"}]}}]

Numerical issues
The OP puts a finger on the problem with the observation in the OP's answer:

if I change the region condition in the definition of xdot from y[t] > -z[t] to y[t]+z[t] > 0 (which seems equivalent to me...)

The inequalities are mathematically equivalent but not numerically equivalent.  This can be demonstrated by the following:
1. + $MachineEpsilon > 1.
(*  False  *)

(1. + $MachineEpsilon) - 1. > 0
(*  True  *)

This is because there is a tolerance in comparing numbers.  If they are approximately equal, as defined by Internal`$EqualTolerance, then strict comparison will return False.
Block[{Internal`$EqualTolerance = 0.},
 1. + $MachineEpsilon > 1.
 ]
(*  True  *)

The documentation page for Less points out another problem,
rounding error:
0.00001 < 2.00006 - 2.00005
(*  True  *)

0.00001 + 2.00005 == 2.00006
(*  True  *)

The problem is almost bound to arise in the OP's answer with the variant
xdot = Piecewise[{{0, y[t] + z[t] > 0}}, xdotbind]

because NDSolve is trying to tiptoe along the path where y[t] == z[t].
If we look at the two forms of the condition in xdot for the OP's first solution (the correct one), we can see that the comparison is an issue:
grid = x["Grid"] /. sol // Flatten; (* the times of the steps *)

x'["ValuesOnGrid"] /. sol // First; (* values of the derivative   x'  *)
Split[Flatten@Position[%, 0.],      (* positions where  x'[t] == 0  *)
  Subtract[##] == -1 &] /. {a_Integer, ___, b_Integer} :> a ;; b
(*  {1 ;; 38}  *)

y[grid[[38 ;; -1]]] > -z[grid[[38 ;; -1]]] /. sol // First // Thread
y[grid[[38 ;; -1]]] + z[grid[[38 ;; -1]]] > 0 /. sol // First // Thread

(*  y[t] > -z[t]  from the correct solution
{True, False, False, False, False, False, False, False, False, False, 
... 30 False's ...
False, False, False, False, False, False, False}

    y[t] + z[t] > 0  from the alternative
{True, False, False, False, False, False, False, False, False, True, 
... 33 True's ...
True, True, True, True}
*)

The first True to False change is where the jump in x[t] occurs around t == 0.73.  In the correct solution, the condition is False at all subsequent steps, which means the x'[t] == xdotbind and the solution has y track z.  The alternative condition is False for a few steps and then is True.  This would set x'[t] == 0 and y would stop tracking z.  This is what happens when NDSolve is run with the alternative condition defining xdot.

