I need to solve a differential equation of the type

$\qquad \partial_{x_1}y(x_1,x_2)= y(x_1,x_1)\,y(x_1,x_2)$

with initial condition

$\qquad y(0,x_2)=x_2$.

Now if I try to solve this with NDSolve I get an error that tells me "the arguments should be ordered consistently" (NDSolve::conarg), here is my code:

s = 
    {D[y[x1, x2], x1] == y[x1, x1] y[x1, x2], y[0, x2] == x2},
    y, {x1, 0, 1}, {x2, 0, 1}]

On my daylong search I haven't found any helpful solution for this seemingly easy problem.

  • 1
    $\begingroup$ Since it's not mentioned in the question, I'll point out that the cause of the error is that the arguments of y[x1, x1] and y[x1, x2] are not the same. $\endgroup$
    – Michael E2
    Sep 10 '18 at 16:35
  • $\begingroup$ The y[x1,x1] is the point of the problem, otherwise it would be trivial. Sorry, I should have clarified that. I want it to work this way. $\endgroup$
    – A G P
    Sep 10 '18 at 17:16
  • $\begingroup$ Just wanted to make sure the meaning of the error was clear. I wonder, what makes this problem seem easy? Is it a standard type? Are there standard methods for numerically integrating it? $\endgroup$
    – Michael E2
    Sep 10 '18 at 22:09
  • 1
    $\begingroup$ Welcome to Mathematica.SE! I hope you will become a regular contributor. To get started, 1) take the introductory tour now, 2) 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, 3) remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign, and 4) give help too, by answering questions in your areas of expertise. $\endgroup$
    – bbgodfrey
    Sep 12 '18 at 3:55

This PDE can be solved symbolically. To begin, replace y[x1, x1] by f[x1], which DSolve can handle without difficulty.

DSolve[{D[y[x1, x2], x1] == f[x1] y[x1, x2], y[0, x2] == x2}, 
    y[x1, x2], {x1, x2}] [[1, 1]]
(* y[x1, x2] -> E^(Integrate[-f[K[1]], {K[1], 1, 0}] 
   - Integrate[-f[K[1]], {K[1], 1, x1}])*x2

From this solution, it is clear that y[x1, x2]/x2 is independent of x2. Consequently, y[x1, x1]/x1 == y[x1, x2]/x2. Substituting this identity into the PDE then yields,

s = DSolve[{D[y[x1, x2], x1] == (x1/x2) y[x1, x2]^2, y[0, x2] == x2}, 
    y[x1, x2], {x1, x2}] [[1, 1]]
(* y[x1, x2] -> -((2 x2)/(-2 + x1^2)) *)

Plot3D[s//Last, {x1, 0, 1}, {x2, 0, 1}, AxesLabel -> {x1, x2, y}, 
    ImageSize -> Large, LabelStyle -> {Bold, Black, Medium}]

enter image description here

Not surprisingly, the plot closely resembles the numerical result of Michael E2. Note that this method works with this PDE for any initial condition.


This solution can, if desired, be validated by back-substitution:

{D[y[x1, x2], x1] == y[x1, x1] y[x1, x2], y[0, x2] == x2} /. 
    y -> Function[{x1, x2}, -((2 x2)/(-2 + x1^2))]
(* {True, True} *)
  • $\begingroup$ Writing down f[x1], you implicitly claim it a known function, but it is not so. Therefore, your solution is built on the sand. Am I not right? $\endgroup$
    – user64494
    Sep 12 '18 at 4:57
  • 1
    $\begingroup$ @user64494 No. I do not claim that f[x1] is known, only that y[x1, x1] is a function only of x1. $\endgroup$
    – bbgodfrey
    Sep 12 '18 at 5:02
  • $\begingroup$ @user64494 Please see the addendum just added. $\endgroup$
    – bbgodfrey
    Sep 12 '18 at 5:07
  • $\begingroup$ OK, Function[{x1, x2}, -((2 x2)/(-2 + x1^2))] is a solution. I still have doubts concerning your DSolving. $\endgroup$
    – user64494
    Sep 12 '18 at 5:20

One could try a iterative approach approximating y[x1, x1] with the previous, if the iterations converge to solution:

foo::div = "Error increased. Divergent?";
foo::slwcon = "Slow convergence. Error `` after `` steps.";
tol = 10^-6;
 Function[{s1, s2, err1, err2, iter},
    With[{s = 
       NDSolveValue[{D[y[x1, x2], x1] == s1[x1, x1] y[x1, x2], 
         y[0, x2] == x2}, y, {x1, 0, 1}, {x2, 0, 1}]},
     {s, s1,
      Flatten[s1["ValuesOnGrid"], 1] - s @@@ Flatten[s1["Grid"], 1] //
         Abs // Max,
      iter + 1}
     ]] @@ # &,
 {NDSolveValue[{D[y[x1, x2], x1] == 1*y[x1, x2], y[0, x2] == x2}, 
   y, {x1, 0, 1}, {x2, 0, 1}],
  1 &,
  10., 100., 0},
 Function[{s1, s2, err1, err2, iter},
    If[err1 > err2,
     Message[foo::div]; False,
     If[err1 > tol,
      If[iter == 10, Message[foo::slwcon, err1, iter]];
    ] @@ # &,
 1, 100
s = {y -> First[%]};

Mathematica graphics

Plot3D[y[x1, x2] /. s, {x1, 0, 1}, {x2, 0, 1}, AxesLabel -> Automatic]

Mathematica graphics

The error is so-so:

Plot3D[D[y[x1, x2], x1] - y[x1, x1] y[x1, x2] /. s // Evaluate,
 {x1, 0, 1}, {x2, 0, 1}, AxesLabel -> Automatic]

Mathematica graphics

  • $\begingroup$ Thanks for your effort! $\endgroup$
    – A G P
    Sep 11 '18 at 14:38

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.