# NDSolve with embedded piecewise function for state-variable

Given ODE:

$$\dot{x}=-2 e^{-x^2} x + U$$

where:

$$U=-u_{as}$$

$$u_{as} = \begin{cases} \frac{e}{\psi_+-e}, & \mbox{if } e\mbox{>=0} \\ \frac{e}{e-\psi_-}, & \mbox{if } e\mbox{<0} \end{cases}$$

$$e=x$$

$$\psi_+=\frac{1-\delta}{T \cdot t+1}+\delta$$

$$\psi_-=-\frac{0.5+\delta}{3T \cdot t+1}-\delta$$

There is my code:

(***)

Clear["Derivative"]

ClearAll["Global*"]

pars = {xs = -1/4, xe = 1, T = 5, \[Delta] = 0.1}

e = x[t]

\[Psi]plus = (1 - \[Delta])/(T t + 1) + \[Delta]

\[Psi]minus = -(1/2 + \[Delta])/(3 T t + 1) - \[Delta]

uas = Piecewise[{{e/(\[Psi]plus - e), e >= 0}, {e/(e - \[Psi]minus),
e < 0}}]

U = -uas

sys = NDSolve[{x'[t] == -2 E^-x[t]^2 x[t] + U, x[0] == xs}, {x}, {t,
0, 3}]

Plot[{Evaluate[x[t] /. sys], \[Psi]plus, \[Psi]minus}, {t, 0, 3},
PlotRange -> Full, PlotPoints -> 100]

Plot[{Evaluate[U /. sys]}, {t, 0, 1}, PlotRange -> Full,
PlotPoints -> 100]


When I run NDSolve, I get a error on the 3rd second of the calculation.

NDSolve::smpf: Failure to project onto the discontinuity surface when computing Filippov continuation at time 2.359060463990643.


I want to understand what the reason for this error is and how to avoid it in the numerical calculations.

I would be grateful for the help.

• It looks like a numerical problem at t=2.35. Decrease the range {t, 0, 2.35} and NDSolve evaluates without error message May 11, 2021 at 8:33
• @UlrichNeumann Yes, if we reduce the time interval of the calculation. But I need calculation at large intervals {t, 0,100} and above. It is necessary to correct the error and avoid it in the future.
– dtn
May 11, 2021 at 8:36
• @UlrichNeumann The problem disappeared when I added to the code Method -> {"DiscontinuityProcessing" -> False}. But I still did not understand what was the reason for the error.
– dtn
May 11, 2021 at 8:42
• Try AccuracyGoal -> 32 or even AccuracyGoal -> Infinity. May 11, 2021 at 14:47
• +1. Consider reporting it to WRI support. I feel there is room for improvement here. There's no real discontinuity in the vector field, just a suspected one arising from Piecewise. But the vector field changes continuously as crosses the equilibrium $x=0$. I suspect it's this singularity in the vector field at the boundary between the pieces that leads to the failed projection. (Raising AccuracyGoal keeps the numerical solution from crossing the boundary, and the true solution does not cross it either.) May 11, 2021 at 20:58

Perhaps StreamPlot shows you the reason:

StreamPlot[{1, -2 E^-x ^2 x  -
Which[x >= 0, x /(0.1 + (0.9/(1 + 5 t)) - x ),x < 0, x /(0.1 + (0.6/(1 + 15 t)) + x ), True, 0]}, {t, 0,5}, {x, -.5, .5}]


The trajectory x->0 seems to be unstable!

• The problem disappeared when I added to the code Method -> {"DiscontinuityProcessing" -> False}.
– dtn
May 11, 2021 at 8:56
• Yes that's already clear. Trajectory x~0 isn't attractive, because of the discontinuity of the ode near x~0 May 11, 2021 at 9:04
• Ok, so the reason is internal, not algorithmic.
– dtn
May 11, 2021 at 9:12
• The reason is your discontinuous ode near x==0. Method -> {"DiscontinuityProcessing" -> False}seems to ignore this fact! May 11, 2021 at 9:16
• That's what I'm talking about.
– dtn
May 11, 2021 at 9:18