I want to set up a model/simulation that involves numerically solving a set of coupled differential equations, but at some point in time I want some of the variables to quickly jump to a constant value for the remainder of the simulation.

More explicitly, the ideal scenario would be that each selected differential equation would look like $$\frac{\text dc}{\text dt} = \begin{cases} f(c,\dots), & t<t_0 \\ 0, & t\geq t_0 \end{cases}$$

and the function $c(t)$ would look like $$c(t)= \begin{cases} \text{solution to }\frac{\text dc}{\text dt}=f(c,\dots), & t<t_0 \\ c_0, & t\geq t_0 \end{cases} $$

where $f(c,\dots)$ is a function of the various dependent variables in the model, $t_0$ is the time at which the "switch" needs to occur, and $c_0$ is the constant value the variable needs to be set to after the "switch".

Of course, I could just make the differential equation look like $$\frac{\text d c}{\text dt}=H(t_0-t)f(c,\dots)+AH(t-t_0)(c_0-c(t))$$ where $H$ is the Heaviside step function and $A$ is very large so that $c$ moves quickly to $c_0$ and the derivative quickly decays to $0$ when $t\geq t_0$. But I didn't know if there was any good machinery in mathematica that could do this better. I know that step functions aren't always the best when numerically solving differential equations, so I feel like there is a faster and more eloquent way to do this.

As a simple example of a differential equation one could use $$\frac{\text dc}{\text dt} = \begin{cases} -c, & t<t_0 \\ 0, & t\geq t_0 \end{cases}$$

and the function $c(t)$ would look like $$c(t)= \begin{cases} \text{solution to }\frac{\text dc}{\text dt}=-c, & t<t_0 \\ 10, & t\geq t_0 \end{cases} $$ where the initial condition of $c(t)$ is not important. Of course this example doesn't require a numeric solution, but an example was requested in the comments. The reason this all needs to be done in the differential equation solver is because in my actual work other equations will depend on versions of this $c(t)$ example.

  • $\begingroup$ You might read the documentation for WhenEvent and see if you could model your problem to be similar to some of the working examples they show. $\endgroup$ – Bill May 2 at 18:32
  • $\begingroup$ @Bill I'll check it out. Thank you $\endgroup$ – Aaron Stevens May 2 at 18:33
  • $\begingroup$ @Bill It looks like WhenEvent is focused on triggering events based on when the dependent variable does something. I need to look at the independent variable, and it needs to be changed for all time, not a one time event after the independent variable passes a threshold. $\endgroup$ – Aaron Stevens May 2 at 22:06
  • $\begingroup$ @AaronStevens Use c'[t]=If[t<=t0,f[c[t]],0] $\endgroup$ – Alex Trounev May 3 at 12:22
  • 1
    $\begingroup$ @AaronStevens Then publish your equations. $\endgroup$ – Alex Trounev May 3 at 13:25

Here is a method for your simple example:

sol = NDSolveValue[
    u'[x]==Piecewise[{{-u[x], x<2}}], u[0]==10, 


Plot[sol[x], {x, 0, 4}]

enter image description here

  • $\begingroup$ Carl, I have had a busy summer and have just now gotten back to applying this to my own project. Could there possibly be some artifacts in numerically solving a set of coupled differential equations where some of the equations depend on the behavior of this function? Sometimes having large steps in differential equations messes up the numeric solving. $\endgroup$ – Aaron Stevens Jun 13 at 14:44

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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