# NDSolve with Piecewise solving BVP gives warning bvdisc [duplicate]

Backslide introduced in 9.0, persisting through 11.0.1

I try to solve a boundary value problem of coupled ODEs. For reference, when I solve this set of equations it works without any problems:

eq1 = {y'[x] == Piecewise[{{z[x], x < 5}, {y[x], x > 5}}]};
eq2 = {z'[x] == Piecewise[{{2 z[x] + 1, x < 5}, {2 y[x] + 1, x > 5}}]};
eqs = {eq1, eq2, y == 1, z == 0};
sol1 = NDSolve[eqs, {y, z}, {x, 0, 10},
Method -> {"Shooting",
"StartingInitialConditions" -> {y == 1, z == 1}}]


When I try to use more than 2 sections in the piecewise function I receive the following error:

eq1 = {y'[x] ==
Piecewise[{{z[x], x < 5}, {y[x], 8 > x > 5}, {2 z[x], x > 8}}]};
eq2 = {z'[x] == Piecewise[{{2 z[x] + 1, x < 5}, {2 y[x] + 1, x > 5}}]};
eqs = {eq1, eq2, y == 1, z == 0};
sol1 = NDSolve[eqs, {y, z}, {x, 0, 10},
Method -> {"Shooting",
"StartingInitialConditions" -> {y == 1, z == 1}}]

NDSolve::bvdisc: NDSolve is not currently able to solve boundary value problems with discrete variables. >>


I searched for a solution to this error online but non of the offered solutions helped. I don't understand why mathematica reads at as discrete variables when more than 2 sections are included. Same error appears whenever a boundary of the form a<x<b is set in the piecewise function. It is important to note that for initial boundary problems, there are no issues. For example:

eq1 = {y'[x] ==
Piecewise[{{z[x], x < 5}, {y[x], 8 > x > 5}, {2 z[x], x > 8}}]};
eq2 = {z'[x] == Piecewise[{{2 z[x] + 1, x < 5}, {2 y[x] + 1, x > 5}}]};
eqs = {eq1, eq2, y == 1, z == 0};
sol1 = NDSolve[eqs, {y, z}, {x, 0, 10}]


Yields

{{y -> InterpolatingFunction[{{0., 10.}}, <>],
z -> InterpolatingFunction[{{0., 10.}}, <>]}}


Does anyone have a solution to this?

It's a backslide introduced since v9:  and another issue of "DiscontinuityProcessing" so we at least have 2 workarounds.

One is to turn off the "DiscontinuityProcessing":

eq1 = {y'[x] == Piecewise[{{z[x], x < 5}, {y[x], 8 > x > 5}, {2 z[x], x > 8}}]};
eq2 = {z'[x] == Piecewise[{{2 z[x] + 1, x < 5}, {2 y[x] + 1, x > 5}}]};
eqs = {eq1, eq2, y == 1, z == 0};
sol1 = NDSolve[eqs, {y, z}, {x, 0, 10}, Method -> {"DiscontinuityProcessing" -> False}]


But as pointed out by Michael E2 in his answer, this is dangerous.

Another workaround is to transform the Piecewise into UnitStep, with the help of the undocumented SimplifyPWToUnitStep:

eq1 = {y'[x] == Piecewise[{{z[x], x < 5}, {y[x], 8 > x > 5}, {2 z[x], x > 8}}]};
eq2 = {z'[x] == Piecewise[{{2 z[x] + 1, x < 5}, {2 y[x] + 1, x > 5}}]};
eqs = {eq1, eq2, y == 1, z == 0};
sol1 = NDSolve[eqs // SimplifyPWToUnitStep, {y, z}, {x, 0, 10}]

• NDSolve implements discontinuities with discrete variables for Piecewise with more than two pieces, not counting the default. This approach may solve some problems, but causes one here. Why do you consider this a bug, since there is a proper warning, instead of it being a limitation? Dec 14, 2016 at 13:49
• @MichaelE2 After a second thought I modified the tag to "backslide", given NDSolve doesn't output a wrong answer in this case. BTW you mean NDSolve is trying to handle the discontinuity with DiscreteVariables in this case? (At least in v9 the warning bvdisc isn't documented at all. ) Maybe you can elaborate a little with an answer? Dec 14, 2016 at 14:31
• I won't have time, at least for several days. Here's a closely related answer: mathematica.stackexchange.com/a/106193/4999 -- I get the same behavior in V9 as V11 (and the message is not documented in V11, either). Dec 14, 2016 at 14:31
• You can see a little of what's going on with {state} = NDSolveProcessEquations[<ndsolve code>] and with state["NumericalFunction"] and state["Variables"]. Image: i.stack.imgur.com/OWNZ8.png Dec 14, 2016 at 14:37

Interesting: Changing the conditions in the Piecewise function to simple one-sided inequalities allows NDSolve to work.

eq1 = {y'[x] == Piecewise[{{z[x], x < 5}, {2 z[x], x > 8}, {y[x], x > 5}}]};
eq2 = {z'[x] == Piecewise[{{2 z[x] + 1, x < 5}, {2 y[x] + 1, x > 5}}]};
eqs = {eq1, eq2, y == 1, z == 0};
sol1 = NDSolve[eqs, {y, z}, {x, 0, 10},
Method -> {"Shooting",
"StartingInitialConditions" -> {y == 1, z == 1}}]


You still get a NDSolve::berr message (scaled boundary error), but that it not surprising in a discontinuous BVP.

• Actually the Piecewise should be modified to Piecewise[{{z[x], x < 5}, {2 z[x], x > 8}, {y[x], x > 5}}]`. Your solution still works of course :) . Dec 14, 2016 at 14:44
• @xzczd Thanks. Much too busy today to be careful. :) Dec 14, 2016 at 14:46
• Thank you very much for the answers @xzczd and Michael E2, definitely solved my problem. Dec 15, 2016 at 7:50
• @Eng.Technion Now you can consider accepting the answer you like by clicking the checkmark sign :) Dec 15, 2016 at 8:05