# Problem solving a nonlinear differential equation with specific boundary conditions

I am trying to solve the following differential equation:

DSolve[{y''[x] == A*B (y[x])^2, y' == -B*C, y'[-D] == 0}, y[x], x]


For my system, this differential equation would be valid in the range of -D < x < 0. Mathematica says:

For some branches of the general solution, unable to solve the conditions.

Without boundary conditions, I get a solution consisting of a Weierstrass elliptic function, but Mathematica is not able to solve for the boundary conditions, even if I use a simpler series expansion of the solution.

In general, I am not directly interested in y[x], but in y/C. I have experimental data for sets of y, B, C and D.

Ideally, I would like to have a general solution for y/C to be able to analyze it and to compare the behavior with my experimental data, but I am not able to solve the boundary value problem. Maybe a part of the problem is that the Weierstrass elliptic functions are periodic and I only need one branch of that ?

If a general solution is not possible, a numerical solution could also be helpful. How could I calculate the value of A for the data sets of y, B, C and D ? I appreciate any help.

• Welcome to Mathematica.SE! I suggest the following: 1) As you receive help, try to give it too, by answering questions in your area of expertise. 2) Read the faq! 3) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! – Michael E2 Dec 3 '14 at 17:51
• 1. It's better to avoid using capitals as variable names as they might conflict with builtins (C and D are builtins) 2. For a numerical solution, ParametricNDSolve might be helpful. – Szabolcs Dec 3 '14 at 18:24
• The Weierstrass Elliptic Function contains 2nd order poles, and DSolve has no way of knowing whether your two boundary points lie on opposite sides of a pole. Perhaps, specifying the allowed range of values for D would help. – bbgodfrey Dec 5 '14 at 2:38

## 1 Answer

For a numerical solution, we can use FindRoot to find the a for given b, c, d, and y == yinit such that the solution to the IVP y == yinit, y' == -b*c, satisfies the boundary condition y'[-d] == 0.

SeedRandom;          (* random data *)
b0 = RandomReal[1, 3];
c0 = RandomReal[1, 3];
d0 = RandomReal[1, 3];
y0 = RandomReal[1, 3];

{solp} = ParametricNDSolve[
{y''[x] == a*b (y[x])^2, y' == -b*c, y == yinit},
y, {x, -d, 0}, {a, b, c, d, yinit}];

FindRoot[
y[a, Sequence @@ #]'[           (* y' for the parameters *)
-#[]     (* -d *)
] /. solp,
{a, 0.01}
] & /@ Transpose[{b0, c0, d0, y0}]
(*
{{a -> -5.20439}, {a -> -14.6164}, {a -> -13.7974}}
*)


Remark: (1) To find a single a for a given b = b1, c = c1, d = d1, y = y1, the call to FindRoot has the form

FindRoot[y[a, b1, c1, d1, y1]'[-d] /. solp, {a, 0.01}]


(2) For some functions the starting value in FindRoot for a is important. For the random data above, it does not happen to matter. (I chose 0.01 to avoid the value 0 in case there was a problem.)