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?
bind = 0
beforexdot
,xdot
is defined to be0
. Ifbind
is defined afterxdot
, thenxdot
is defined in terms ofxdot
, which causes an infinite recursion. $\endgroup$ – Michael E2 Apr 26 '16 at 0:46xdot
now (and definebind=0
at the appropriate time). However, the central issue remains unresolved. $\endgroup$ – clr66 Apr 26 '16 at 17:45xdot
whenbind==1
in the Piecewise function, NDsolve doesn't generate any dynamics inx[t]
after the WhenEvent triggers. $\endgroup$ – clr66 Apr 26 '16 at 17:54