Why does DSolve fail in this first order differential equation?

Question

I am trying to solve the following first order differential equation: $$ g'(R)=-2 \sqrt{R^2-g(R)}+2 R$$

By direct substitution it can be verified that an obvious solution is $$g(R)=R^2$$

ode=2R-2Sqrt[R^2-g[R]]-Derivative[1][g][R]
g[R_]:=R^2
ode==0

Giving True as expected.

However, using DSolve (Mathematica 11.2 for Windows 10):

Clear[g]
DSolve[ode==0, g[R], R]//FullSimplify//Expand

we find

{{g[R]->-1+E^(2 C[1])-2 R-2 Sqrt[-E^(2 C[1]) (1+R)^2]},
{g[R]->-1+E^(2 C[1])-2 R+2 Sqrt[-E^(2 C[1]) (1+R)^2]}}

As seen, Mathematica completely missed the correct solution!

Edit 1: For more details on the mathematics aspects of the solution (and the two branches of the solution), see my post at the Mathematics stack.

Edit 2: Adding initial conditions like $g(0)=0$ or $g(1)=1$, as suggested in the comments, unfortunately doesn't help.

Edit 3: Just for future reference I add the way I found to solve the ode and recover the $R^2$ solution. To do that define a new function $y(R)$ that should be zero when $g(R)=R^2$, in the following manner

g[R_] := R^2 (1 - y[R]^2)
ode1 = ode // FullSimplify // PowerExpand // FullSimplify;
DSolve[ode1 == 0, y[R], R] // Expand

which gives

{{y[R] -> 0}, {y[R] -> 1 + C[1]/R}}

The solution $y[R] \rightarrow 0$ implies $g(R)=R^2$ as expected.

After massaging the ODE as shown above, I can get Mathematica to give the correct solution, however I still don't understand why it failed in the first place.

So, my question is, why did DSolve fail in this simple example?

Answer 1

I can no longer test V11.2, but in V13, DSolve has been improved. It still needs some help. The option IncludeSingularSolutions yields more solutions, but we need to rationalize the ode to get all of them. Of course, the solutions to the rationalized ode need to be checked on the original ode.

The constant of integration shows up in several places as Exp[2 C[1]]. One can apply constSimplify or just do the substitution C[1] -> Log[C[1]]/2 directly to get a nicer looking result. We check the solutions and find conditions for which they are valid.

dsol = DSolve[ode, g[R], R, IncludeSingularSolutions -> True] /. 
   C[1] -> Log[C[1]]/2 // Simplify
(*
{{g[R] -> -1 - 2 R + C[1] - 2 Sqrt[-(1 + R)^2 C[1]]},
 {g[R] -> -1 - 2 R + C[1] + 2 Sqrt[-(1 + R)^2 C[1]]},
 {g[R] -> 0},
 {g[R] -> -1 - 2 R}}
*)

cond = ode /. DSolve`DSolveToPureFunction@dsol // Map[Reduce[#, Reals] &]
(*
{R > -1 && C[1] <= 0,
 (R < -1 && C[1] <= -1 - 2 R - R^2) ||
  (R > -1 && -1 - 2 R - R^2 <= C[1] <= 0), 
 R >= 0,
 R >= -1}
*)

Likewise, for the rationalization oderat of the original ode, we solve the differential equation and find conditions for which they are valid.

oderat = Eliminate[{ode /. Sqrt[R^2 - g[R]] -> dummy, dummy^2 == R^2 - g[R]}, 
  dummy]

(*  (4 R - g'[R]) g'[R] == 4 g[R]  *)

dsolrat = DSolve[oderat, g[R], R, IncludeSingularSolutions -> True]

(*  {{g[R] -> R C[1] - C[1]^2/4}, {g[R] -> R^2}}  *)

condrat = 
 ode /. DSolve`DSolveToPureFunction@dsolrat // Map[Reduce[#, Reals] &]

(*  {R >= C[1]/2, True}  *)

We can then assemble the final solutions, using ConditionExpression to restrict each solution (optional). All seem to be valid for some values of C[1] and interval for R.

allsols = {g[R] -> #} & /@ MapThread[ConditionalExpression,
   {g[R] /. Join[dsol, dsolrat], Join[cond, condrat]}]
(*
{{g[R] -> ConditionalExpression[
     -1 - 2*R + C[1] - 2*Sqrt[(-(1 + R)^2)*C[1]], 
     R > -1 && C[1] <= 0]}, 
 {g[R] -> ConditionalExpression[
     -1 - 2*R + C[1] + 2*Sqrt[(-(1 + R)^2)*C[1]], 
     (R > -1 && 2*R + R^2 + C[1] >= -1 && C[1] <= 0) || 
      (R < -1 && 2*R + R^2 + C[1] <= -1)]}, 
 {g[R] -> ConditionalExpression[ 0, R >= 0]}, 
 {g[R] -> ConditionalExpression[ -1 - 2*R, R >= -1]}, 
 {g[R] -> ConditionalExpression[ R*C[1] - C[1]^2/4, 
     R >= C[1]/2]},
 {g[R] -> R^2}}
*)

Re-check solutions, assuming the conditions are satisfied for each one:

Map[
 With[{conditions = 
     Cases[#, ConditionalExpression[_, assum_] :> assum, Infinity]}, 
   Simplify[Reduce[Flatten@{ode /. #, conditions}], conditions]] &,
 DSolve`DSolveToPureFunction@allsols]

(*  {True, True, True, True, True, True}  *)

Answer 2

Usually a solution to an ODE needs as many constants as their degree. In this case, the solution

$$ g(R) = R^2 $$

doesn't requires such a constant. Analyzing the phase plane for this ODE we can observe that the initial conditions for a real solution should obey

$$ R_0 \ge 0,\ \ g(R_0) \lt R_0^2 $$

and for $g(R)= R^2$ we have a particular solution. Follows a plot showing the stream plot (phase plane) for this ODE, in red the particular solution and in green two solutions: for $g(0) = 0$ and for $g(2+\epsilon^2) = 4$

gr1 = StreamPlot[{1, -2 (Sqrt[x^2 - y] - x)}, {x, 0, 10}, {y, -1, 9}];
gr2 = ContourPlot[y == x^2, {x, 0, 10}, {y, -1, 9}, ContourStyle -> Red];
sol3 = NDSolve[{g'[r] == -2 Sqrt[-g[r] + r^2] + 2 r, g[0] == 0}, g, {r, 0, 10}][[1]];
gr3 = Plot[Evaluate[g[r] /. sol3], {r, 0, 10}, PlotStyle -> Green];
sol4 = NDSolve[{g'[r] == -2 Sqrt[-g[r] + r^2] + 2 r, g[2.01] == 4}, g, {r, 0, 10}][[1]];
gr4 = Plot[Evaluate[g[r] /. sol4], {r, 0, 10}, PlotStyle -> Green];
Show[gr1, gr4, gr3, gr2]

Answer 3

Not an answer, but might help getting help...

I evaluated the integral using the package NCAlgebra, which considers g as a noncommutative variable. You can see (it you click on the image...) that it transforms the DSolve into a Solve with two integration constants.

<< NC`;
<< NCAlgebra`;
DSolve[2 R - 2 Sqrt[R^2 - g[R]] - Derivative[1][g][R] == 0, g[R], R]

Maybe some gurus can infer what is happening from there (probably MMA transforms the ODE).