# Split Boundary Value Problems with Algebraic Equations

Is it possible at all to solve with NDSolve (or other built in function) a split boundary value problem with algebraic equations?

Please look at the following example:

s = NDSolve[
{
x'[t] == y[t] + x[t] y[t] + z[t],
z'[t] == 2*y[t],
2 x[t] + y[t] - z[t] == 1,
x[0] == 0,
z[0] == 0}
, {x, y, z}, {t, 0, 1}];


This example solves well, and outputs Interpolation Functions.

I compute the final value of variable z:

In[50]:= z[1] /. s
Out[51]:= {0.694658}


When I try the same example with split conditions:

s1 = NDSolve[
{
x'[t] == y[t] + x[t] y[t] + z[t],
z'[t] == 2*y[t],
2 x[t] + y[t] - z[t] == 1,
x[0] == 0,
z[1] == 0.694658}
, {x, y, z}, {t, 0, 1}];


I get the following error:

NDSolve::bvdae: Differential-algebraic equations must be given as initial value problems.


I am trying to solve a complex problem with this structure, is there any way to tackle it in mathematica?

Edit: Unlike this Minimum Working Example, in the original problem I am trying to solve the variable y cannot be eliminated from the system

Indeed, NDSolve cannot solve this equation as written. However, it is easy enough to eliminate y from the system.

{x'[t] == y[t] + x[t] y[t] + z[t], z'[t] == 2*y[t]} /. y[t] -> 1 + z[t] - 2 x[t]


and then solve and plot

s1 = NDSolve[{Derivative[1][x][t] == 1 - 2 x[t] + 2 z[t] + x[t] (1 - 2 x[t] + z[t]),
Derivative[1][z][t] == 2 (1 - 2 x[t] + z[t]), x[0] == 0,
z[1] == 0.694658}, {x, z}, {t, 0, 1}]
Plot[Evaluate[{x[t], 1 + z[t] - 2 x[t], z[t]} /. s1], {t, 0, 1},
AxesLabel -> {t, "x, y, z"}]


A similar problem is 63824. Its approach can be applied as follows.

s3 = ParametricNDSolve[{x'[t] == y[t] + x[t] y[t] + z[t], z'[t] == 2*y[t],
2 x[t] + y[t] - z[t] == 1, x[0] == 0, z[0] == a}, {x, y, z}, {t, 0, 1}, {a}];
f[w_?NumericQ] := z[w][1] /. s3
w0 = w /. FindRoot[f[w] == 0.694658, {w, .1}]
(* -2.47361*10^-6 *)
Plot[Evaluate[{x[w0][t], y[w0][t], z[w0][t]} /. s3], {t, 0, 1},
AxesLabel -> {t, "x, y, z"}]


which gives the same plot shown above. Note that the value of w0 is zero to roundoff.

• Thank you, unfortunately I am trying to solve a more complex system where the equivalent of y cannot be eliminated. – antonio Feb 4 '16 at 17:23
• @Antonio You may find the discussion here useful. – bbgodfrey Feb 4 '16 at 17:35
• @Antonio In response to your comment, I was able to find an alternative approach that does not involve eliminating y from the equations. Basically, it is the standard shooting method. Why NDSolve does not use it in this case is unclear to me. – bbgodfrey Feb 4 '16 at 21:55