$
\def\figA#1{\def\nextfig##1{\figB##1}\def#1{1}#1}
\def\figB#1{\def\nextfig##1{\figC##1}\def#1{2}#1}
\def\figC#1{\def\nextfig##1{\figD##1}\def#1{3}#1}
\def\figD#1{\def\nextfig##1{\figE##1}\def#1{4}#1}
\def\figE#1{\def\nextfig##1{\figF##1}\def#1{5}#1}
\def\figF#1{\def\nextfig##1{\figG##1}\def#1{6}#1}
\def\figG#1{\def\nextfig##1{\figH##1}\def#1{7}#1}
\def\figH#1{\def\nextfig##1{\figI##1}\def#1{8}#1}
\def\figI#1{\def\nextfig##1{\figJ##1}\def#1{9}#1}
\def\figJ#1{\def\nextfig##1{\figK##1}\def#1{10}#1}
\def\figK#1{\def\nextfig##1{\figL##1}\def#1{11}#1}
\def\nextfig#1{\figA#1}
$
Introduction
In this answer we will develop a geometric picture of the OP's differential equation. It will illustrate the behavior of DSolve
on such equations. We will show that the singular solution arises as the discriminant of the ODE viewed as an equation in three variables, and that the discriminant curve necessarily is the same as the envelope of the general solution shown in the answer by Michael Seifert. The approach will reveal two more solutions, one of which might be considered spurious according to one's point of view. We will see that the general solution is valid only over a restricted domain. Finally, our analysis will suggest reasons why DSolve
does not, and probably does not attempt to, find the singular solution.
We will give a general development of the ideas and return from time to time to the OP's ODE as a particular example.
Set-up
A solution $y=y(x)$ of a first-order ODE $F(x,y(x),y'(x))=0$ parametrizes a curve on the surface
$$F(x,y,p)=0\,,$$
where we think of the coordinate $p$ as $p = dy/dx$. The projection $(x,y,p)\mapsto(x,y)$ maps the curve on the surface back to the solution curve in the $xy$ plane. In what follows, we will assume $F$ is sufficiently smooth and, for the code, that $F$ is a polynomial (or rational). That is not required by mathematics, but if the ODE is a polynomial at least in $p$, then there is a simple function in Mathematica that can yield the answer.
In the OP's case, we need to rationalize the ODE to eliminate the radical. There are many ways to do this (by hand is probably faster in this case than trying to figure out how to do it in Mathematica), such as GroebnerBasis
.
Clear[x, y, p, F];
ode = y'[x] == y[x]/x - 2/x^2 + (2*Sqrt[1 - x y[x]])/x^2;
odeToFn = {y[x] -> y, y'[x] -> p, Equal -> Subtract};
fn = ode /. odeToFn;
F[x_, y_, p_] = Module[{u},
With[{rads = DeleteDuplicates@
Cases[fn, Power[a_, b_Rational] :> u[a, b], Infinity]},
First@GroebnerBasis[
Flatten@{fn /. Power[a_, b_Rational] :> u[a, b],
Map[
Numerator@ Together[ (* gets rid of denominators in neg. powers *)
#^Denominator[Last@#] - First[#]^Numerator[Last@#]] &,
rads]},
rads]
]]
F[x_, y_, p_] = Numerator@ Together@ F[x, y, p]
(* p^2 + (4 p)/x^2 - (2 p y)/x + y^2/x^2 *)
(* 4 p + p^2 x^2 - 2 p x y + y^2 *)
Clearing the denominators with Numerator@Together@F[x, y, p]
can help with numerics.
Rationalizing usually introduces extra solutions arising from conjugate equations. We will return to this issue after we have developed a clearer picture of the ODE.
Terminology
For convenience we will call the solution returned by DSolve
the general solution. A particular solution of the general solution we will call a regular solution. The differential equation ode
we will call the OP's ODE. The rationalized equation $F(x, y(x), y'(x))=0$ we will call the differential equation $F = 0$; the surface $F(x,y,p)=0$ we will call the surface $F=0$. We will use these terms generally; the function F[x, y, p]
defined above is a particular example.
Simplest answer
When there is a singular solution, it is simply the discriminant of the equation:
Discriminant[F[x, y, p], p] == 0 // Simplify
(* x y == 1 *)
Now, every surface that folds over itself will have such a discriminant, but it will not turn always out to be a solution to the corresponding differential equation $F(x,y,y')=0$. We will discuss why and when the discriminant yields a singular solution below.
What's going on?
The singular points of the projection from the surface $F=0$ to the $xy$ plane form the so-called criminant. It is the set of points on the surface for which $\partial F/\partial p = 0$. The discriminant curve is the projection of the criminant onto the $xy$ plane. (See, for instance, Arnold, Geometric Methods in the Theory of Differential Equations.)
The discriminant may be obtained by eliminating $p$ from the system
$$(1)\ \ F(x,y,p)=0,\quad\text{and}\quad(2)\ \ F_p(x,y,p)=0\,.$$
Eliminate[{F[x, y, p] == 0, D[F[x, y, p], p] == 0}, p]
(* x y == 1 *)
Geometrically, this corresponds to where the tangent plane is vertical (parallel to the $p$ axis):
Fig. $\nextfig\criminant$. The criminant (left) and its projection onto the discriminant curve (right).
GraphicsRow[{#, Show[#, ViewPoint -> {0, 0, Infinity}]}] &[
ContourPlot3D[F[x, y, p] == 0,
{x, 0.02, 2}, {y, 0.02, 2}, {p, -4, 0},
ContourStyle -> Opacity[0.5],
MeshFunctions -> {Function[{x, y, p}, x^2 p + 1]},
Mesh -> {{0}}, MeshStyle -> {Red, Thick]},
AxesLabel -> Automatic]
]
(* Because x y <= 1 always, we use a different form x^2 p == -1 in
the MeshFunctions to draw the criminant. *)
If the discriminant defines $y$ as a smooth function of $x$, then we can deduce a relation for its derivative:
$$
0={d \over dx}\,F{(x,y,y')}=F_x+F_y\,y'+F_p\,y''=F_x+y'F_y \implies
y'=-F_x/F_y\,.
\tag{1}
$$
Before discussing when the discriminant is a solution to the ODE, let us first look at its connection to the envelope of the general solution.
Why is the discriminant the same as the envelope of the general solution?
Suppose we have a smooth parametrization of the general solution in the form $y = y(x,c)$. Since
$$
F(x,\, y(x,c),\, y'(x,c))=0
$$
(where $y'(x,c)$ is the derivative with respect to $x$), then we have at the (dis)criminant where $F_p=0$
$$
{\partial \over \partial c}\,F(x,\, y(x,c),\, y'(x,c)) =
F_y\,{\partial\over\partial c}\,y(x,c) = 0
$$
Thus, provided $F_y \ne 0$, which is generally true†, the discriminant satisfies the conditions $y = y(x,c),\ 0 = {\partial y / \partial c}$, which are the conditions for the envelope.
[†Otherwise $F_y=F_p=0$ and $dF=0$ would imply $F_x=0$ and the surface would have a singular point.]
Fig. $\nextfig\family$. The family of solutions for different values of the constant, plotted on the surface $F=0$ (left) and showing the envelope (right).
{sol} = DSolve[ode, y, x];
GraphicsRow[{#, Show[#, ViewPoint -> {0, 0, Infinity}]}] &[
ParametricPlot3D[{x, y[x], y'[x]} /. sol /. C[1] -> 1/c // Simplify // Evaluate,
{c, -1, -0.01}, {x, Max[0.02, -((1 + 2 c)/(2 c^2))], Min[2, (-0.005 - c)/c^2]},
PlotPoints -> 35, PlotStyle -> Opacity[0.5],
MeshFunctions -> {Function[{x, y, p, c, x2}, c]},
MeshStyle -> {Thick}, BoxRatios -> {1, 1, 1},
AxesLabel -> {x, y, p}] /. {foo = 0; mesh_Line :> {ColorData[97, ++foo], mesh}}
]
What is more, if any curve on the surface $F=0$ passes through the criminant, its projection onto the $xy$ plane intersects the discriminant in one of two ways. Either it is tangent or it has a singularity at the point of intersection. For at a point on the criminant we have exactly as in equation (1)
$$0={dF}=F_x\;dx+F_y\;dy+F_p\;dp \implies dy/dx=-F_x/F_y\,,$$
and the solution will be tangent to the discriminant, provided it is nonsingular.
What's really going on?
We can construct a direction field on the surface $F=0$ by a system of planes, called a contact structure, as follows. At a point $(x,y,p)$ on the surface, the contact plane is given by $dx = p\; dy$, i.e., the plane through $(x,y,p)$ with normal $(1,-p,0)$. It is vertical and it intersects the $xy$ plane in a line with a slope $dy/dx$ that is given by the height $p$. A curve on the surface $F=0$ that is tangent to the contact plane at each point it passes through is called an integral curve of the contact structure. Its projection onto the $xy$ plane is an integral curve of the ODE $F=0$.
Fig. $\nextfig\contact$. Contact structure on the surface: $F=0$ (top), $p^2=y$ (left), and $p^2=x$ (right). The integral curves are shown on each surface and their projections, which are integral curves of the ODE are shown below.
At a point on the criminant there are two possibilities. The contact plane intersects the tangent plane in a vertical line, or the contact plane is identical to the tangent plane.
The second case implies that the slope of the discriminant is given by $p = dy/dx$ and satisfies the ODE.
If each contact plane of the criminant is tangent to the criminant, the discriminant curve will be an integral curve of the ODE.
It should be clear that the contact planes coinciding with the tangent planes is an exceptional situation and to have all of them line up along the criminant is really unlikely. In that case, the contact planes do not define a direction on the surface, and technically, the criminant is not an integral curve, but a singular locus of the direction field.
A generic, small, smooth perturbation of the surface $F=0$ will transform the contact system to the first case.
For instance, if the surface is translated up, the slope $p$ increases and the contact planes rotate but the tangent planes stay the same (Fig. $\nextfig\deformation$). If the surface is rotated, then the tangent planes rotate but not the contact planes.
Fig. $\deformation$. Deformation of the contact structure (via translation). A slight raising or lowering changes the direction field, which is given by the height $p$ at each point. In turn the integral curves are altered, most dramatically near the discriminant curve, which itself is destroyed as a solution.
After such a deformation, the discriminant curve will no longer be solution of the differential equation. The contact plane of the transformed surface will define a vertical direction at the criminant. The integral curves of the contact structure will pass through the criminant, and their projections onto the $xy$ plane will have a cusp on the discriminant.
Indeed this is the norm. Generically at the criminant an ODE can locally be put in the normal form (Fig. $\contact$, right)
$$P^2=X \quad (P = dY/dX)$$
by suitable change of coordinates. Only exceptionally can it be put in the form (Fig. $\contact$, left)
$$P^2=Y$$
for which the $X$ axis is the criminant and the exceptional/singular solution. Indeed the OP's ODE may be put in this form by the change of coordinates
$$y= {1 \over x}-x\,Y,\quad x = 2/X\,.$$
(* change of variables *)
alg2Fn = HoldPattern[y_[x_] -> fx_] :> (y -> Function[{x}, fx]);
x2X = First@Solve[x == 2/X, x]; (* change x to X *)
X2x = First@Solve[x == 2/X, X]; (* change X to x *)
y2Y = First@Solve[y[x] == 1/x - x Y[2/x], y[x]] /. alg2Fn; (* change y to Y *)
Y2y = First@Solve[y[x] == 1/x - x Y[2/x] /. x2X, Y[X]] /. alg2Fn; (* change Y to y *)
(* mapping ode to P^2 == Y *)
odeY = F[x, y[x], y'[x]] == 0 /. y2Y /. x2X // Simplify
(* Y[X] == Y'[X]^2 *)
Note that DSolve
does not find the singular solution Y -> Function[{X}, 0]
to this equation either, even though it is "obvious." However, the general and singular solution transform under the change of coordinates back to the general and singular solution of of the ODE $F=0$, after adjusting the constants of integration:
Ysol = DSolve[odeY, Y, X]
y[x] /. y2Y /. x2X (* change y[x] to terms of Y, X *)
% /. Ysol /. X2x /. C[1] -> 2 C[1] // Expand
y[x] /. DSolve[F[x, y[x], y'[x]] == 0, y[x], x] /. C[1] -> 2/C[1] // Expand
y[x] /. DSolve[ode, y[x], x] /. C[1] -> 2/C[1] // Expand
(*
{{Y -> Function[{X}, 1/4 (X^2 - 2 X C[1] + C[1]^2)]},
{Y -> Function[{X}, 1/4 (X^2 + 2 X C[1] + C[1]^2)]}}
X/2 - (2 Y[X])/X <-- y[x] in terms of Y, X
{2 C[1] - x C[1]^2, -2 C[1] - x C[1]^2} <-- Ysol in terms of y, x
{2 C[1] - x C[1]^2, -2 C[1] - x C[1]^2} <-- DSolve @ F == 0
{-2 C[1] - x C[1]^2} <-- DSolve @ ode (OP's)
*)
There are a couple of things worth observing about these results.
First, the DSolve
solution of the ODE $F = 0$ returns two formulas, while the solution to ode
has only one. This is due to the rationalization of ode
that produced F
. The two formulas are equivalent and may be reconciled by replacing C[1] -> - C[1]
.
Second, the missing singular solution of odeY
is mapped to the missing singular solution of the OP's ode
:
Y[X] == 0 /. Y2y /. X2x // Simplify
(* y[x] == 1/x *)
On solving for the singular solutions
The singular solution begins to seem special. Let us consider how one might solve for it. The simplest way to find the singular solution is to check whether the discriminant curve satisfies the differential equation. Another approach is to integrate the contact structure. At the end we will review the desirability of automatically solving for the singular solution.
Checking the discriminant
Instead of Eliminate
, we can use three-argument form of Solve
in which the third argument is a List
of variables to be eliminated. This will return the discriminant as a list of formulas for y
. We then can pick the ones that satisfy the ODE.
ClearAll[singularSolve];
singularSolve[de_Equal, form : y_[x_] | y_, x_] :=
Module[{p, q, u, fn, rads, a, b, ode, discriminant},
fn = de /. {y[x] -> y, y'[x] -> p, Equal -> Subtract};
rads = DeleteDuplicates@ Cases[fn, Power[a_, b_Rational], Infinity];
ode = First@ GroebnerBasis[
Flatten@ {
fn /. Power[a_, b_Rational] :> u[a, b],
rads /. Power[a_, b_Rational] :> u[a, b]^Denominator[b] - a^Numerator[b]},
rads /. Power[a_, b_Rational] :> u[a, b]];
discriminant = Solve[{ode == 0, D[ode, p] == 0}, {y}, {p}];
Pick[
Which[
MatchQ[form, y[x]], {y[x] -> (y /. #)},
MatchQ[form, y], {y -> Function @@ {{x}, y /. #}},
True, # /. p -> y'[x]
] & /@ discriminant,
TrueQ@Simplify[ode == 0 /. Flatten@{#, p -> D[y /. #, x]}] & /@ discriminant
]
]
Examples:
singularSolve[ode, y, x] (* OP's ode *)
(* {{y -> Function[{x}, 1/x]}} *)
singularSolve[y'[x]^2 == y[x], y[x], x] (* generic singular solution *)
(* {{y[x] -> 0}} *)
singularSolve[y'[x]^2 == x, y, x] (* generic normal form (no singular curve) *)
(* {} *)
Integrating the contact structure
By solving the contact structure we would find four solutions to the ODE $F=0$:
$$
\begin{align}
y &= 1/x & &\text{singular solution} \\
y &= -4 (c + x)/c^2 & &\text{general solution} \\
y &= 0 & &\text{limit of general solution} \\
x &= 0 & &\text{projective closure of general solution} \\
\end{align}
$$
We will explain how these arise below. The fourth one may or may not be admitted depending on one's viewpoint.
Let's illustrate how the singular curve arises in the simple case of $(dY/dX)^2=Y$. Practically by inspection one can see that the singular solution is $Y=0$, which is not represented by the general solution $Y={1\over4}(X-c)^2$ returned by DSolve
nor a limit of it. If we solve the contact system, we will see how the singular solution falls out. The corresponding equations are
$$
\begin{align}
P^2 &= Y & &\text{contact manifold} \\
2P \; dP &= dY & &\text{tangent field} \\
dY &= P \; dX & &\text{contact field} \\
\end{align}
$$
whence $P(2\;dP-dX)=0$. One factor gives $P=0$, from which follows the singular curve $Y=0$. Integrating the other factor yields $2P=X-c$ or after squaring, $Y={1\over4}(X-c)^2$.
In general, the contact system to be solved is
$$
\begin{align}
F(x,y,p) &= 0 & &\text{contact manifold} \\
dF(x,y,p) &= 0 & &\text{tangent field} \\
dy &= p \; dx & &\text{contact field} \\
\end{align}
$$
The general approach is to reduce the direction field equations, that is, the tangent and contact fields, to disjoint systems of algebraic and differential equations (via Reduce
). The contact manifold $F(x,y,p)=0$ is used to eliminate variables where necessary. The differential equations are integrated. Then we can get $y$ in terms of $x$ and a constant of integration, eliminating $p$ as needed. It turns out the Reduce
may produce a component of the system with more than one differential equation. For instance, in the OP's case I obtained four disjoint systems, the first two of which have two differential equations:
(p == 0 && dy == 0 && -2 + x y != 0 && dp == 0 && dx != 0) || (*=> y == 0 *)
(p != 0 && dx == dy/p && 2 + p x^2 - x y != 0 && dp == 0 && dy != 0) || (*=> y == (4(c-x))/c^2 *)
(x != 0 && p == (-2 + x y)/x^2 && p != 0 && dx == dy/p && dy != 0) || (*=> y == 1/x *)
(y != 0 && x == 2/y && p == 0 && dy == 0 && dx != 0) (* => not a solution *)
Their integrations produce independent constants that need to be reconciled.
It is worth emphasizing that solving this system yields both the general and the singular solutions. In the OP's case, it turned up two solutions not yet mentioned. One is
y -> Function[{x}, 0]
This is the limit of the general solution returned by DSolve
as C[1] -> Infinity
. Another is x == 0
: It can be considered an integral curve of the projective closure of the contact structure at p == ∞
. It makes no sense as a solution for y
, but if we think of the differential equation in terms of (1) x
as a function of y
and (2) p == 1/x'[y]
, then it is a solution of the differential equation
F[x[y], y, 1/x'[y]] == 0
(after clearing denominators). In fact, F
is symmetric in this way and this differential equation is essentially the same as the original with x
and y
switched.
DSolve[Numerator@Together@F[x[y], y, 1/x'[y]] == 0, x, y]
(*
{{x -> Function[{y}, (4 (-y + C[1]))/C[1]^2]},
{x -> Function[{y}, -((4 (y + C[1]))/C[1]^2)]}}
*)
One can also view the $y$ axis $x == 0$ as a limit of the solution as C[1] -> Infinity
, or in terms of the original general solution as C[1] -> 0
, which is indeed how it appears in Michael Seifert's animation.
Note: The results obtained depended on how the Reduce
command was set up. I have had trouble developing an approach that did not need tweaking on a case by case basis. In the big picture it does not seem that important to do so.
The OP's ODE revisited
We promised to say something about the OP's original ode
and the rationalized F[x, y[x], y'[x]] == 0
. Two observations raise questions about the DSolve
general solution and its relationship to the contact structure. One is that the direction field does not exist along the criminant. The other is that the surface $F=0$ folds over itself along the criminant, and the surface that lies above criminant corresponds to the OP's ode
while the surface below corresponds to the conjugate ODE; and yet the traces of the general solution on $F=0$ cross the criminant.
The general solution is valid only over a restricted domain
In fact the general solution of the OP's ode
is valid only up to the point it contacts the discriminant.
{sol} = DSolve[ode, y, x];
Reduce[ode && $y == y[x] /. sol, {y, x}, {C[1]}, Reals] (* compute domain *)
(* ($y < 0 && (1/$y <= x < 0 || x > 0)) || ($y > 0 && (x < 0 || 0 < x <= 1/$y)) *)
If we check the residual error of ode
for the solution that passes through $(1,1)$, we see it is zero up until x == 1
at which point y == 1
and the solution contacts the discriminant.
Fig. $\nextfig\residual$. After the solution touches the discriminant at x == 1
, it no longer satisfies the ODE ode
(left), but it does satisfy the ODE F[x, y[x], y'[x]] == 0
(right).
Block[{C},
C[1] = C[1] /. First@Solve[y[1] == 1 /. sol, C[1]];
GraphicsRow[{
Plot[ode /. Equal -> Subtract /. sol, {x, 0.5, 2},
PlotRange -> All, Frame -> True],
Plot[F[x, y[x], y'[x]] /. sol, {x, 0.5, 2}, PlotRange -> 0.8*^-15,
Frame -> True, WorkingPrecision -> 16]
}, ImageSize -> 533]
]
Using Reduce
instead of Solve
within DSolve
One can get Mathematica to produce this restriction by setting the option Method -> Reduce
to Solve
as I did here.
With[{opts = Options[Solve]},
Internal`WithLocalSettings[
SetOptions[Solve, Method -> Reduce],
{redsol} = DSolve[ode, y, x],
SetOptions[Solve, opts]
]]
Solve::useq: The answer found by Solve contains equational condition(s)....A likely reason for this is that the solution set depends on branch cuts of Wolfram Language functions. >>
(*
{{y -> Function[{x},
ConditionalExpression[-((4 (x + C[1]))/C[1]^2),
-((2 x)/C[1]) + Sqrt[(2 x + C[1])^2/C[1]^2] == 1]]}}
*)
We simplify the condition further.
redsol /. Cases[
redsol,
ConditionalExpression[val_, cond_] :>
(HoldPattern[ConditionalExpression[val, cond]] ->
With[{c = Reduce[cond, x, Reals]}, ConditionalExpression[val, c]]),
Infinity]
(*
{y -> Function[{x},
ConditionalExpression[-((4 (x + C[1]))/C[1]^2),
(2 x + C[1] >= 0 && C[1] > 0) || (2 x + C[1] <= 0 && C[1] < 0)]]}
*)
The condition 2 x + C[1] == 0
is equivalent to the solution touching the discriminant y = 1/x
:
Reduce[y[x] == 1/x /. redsol]
(* C[1] != 0 && x == -(C[1]/2) *)
Fig. $\nextfig\domain$ shows a plot of the valid solutions to ode
. (One may observe the uniqueness of a solution with a given initial condition.)
Fig. $\nextfig\domain$. The general solution from DSolve
is valid only up until it contacts the singular solution.
Show[
Plot[1/x, {x, 0, 2}, PlotStyle -> Gray, PlotRange -> {0, 2}, AspectRatio -> Automatic],
Table[Plot[
Evaluate@
Block[{C},
C[1] = C[1] /. First@Solve[y[x0] == 1/x0 /. redsol, C[1]];
y[x] /. redsol],
{x, 0, 2}, PlotStyle -> ColorData[97, Round[(x0 - 0.4)/0.2]],
MeshFunctions -> {Function[{x, y}, x y]}, Mesh -> {{0.99999}},
MeshStyle -> {PointSize[Medium], ColorData[97, Round[(x0 - 0.4)/0.2]]}],
{x0, 0.6, 1.8, 0.2}]
]
Connection to the perturbed ODE
One can see in the animation above (Fig. $\deformation$) how the general solution is "cut" at the discriminant and reflected back along the discriminant as the translated ODE $F(x,y,y'-a)=0$ moves through $a=0$. For a small perturbation, two branches of a regular solution meet in a cusp at the discriminant, one of them closely tracking a regular solution (line) of the orginal ODE ($a=0$) and the other closely tracking the discriminant.
Fig. $\nextfig\perturbed$. The perturbed ODE $F(x,y,y'-a)=0$ for $a = \pm0.01$. Each solution has a cusp on the discriminant. One branch follows closely a line of the original general solution (shown lightly shaded); the other follows closely the discriminant.
It is interesting that one can see the singular solution as a limit, not of the general solution to the OP's ODE, but of the ODEs in a neighborhood of $F(x,y,y')=0$. When the criminant aligns with the contact field, the singular solution appears as an independent solution and the solutions pass by the discriminant (or through the criminant).
Conclusion
The upshot is that the rationalized ODE admits as solutions the complete lines in the general solution to the OP's ode
, which only admits rays starting at the singular solution.
Aside: Connection to numerical integration
The picture of the ODE as a contact system also illustrates the numerical integration of the ODE.
- When the integration of a regular solution by
NDSolve
reaches the discriminant/singular solution, the system becomes stiff and integration can continue only along the singular solution. This is true for both ode
and $F=0$; in the case of F[x, y[x], y'[x]] == 0
, NDSolve
first solves for y'[x]
and effectively creates two ODEs, equivalent to ode
and its conjugate.
- The discriminant appears as boundary in
ode
in the form Sqrt[1 - x y[x]]
, which causes complex numbers to enter the calculations as round-off error and approximate steps lead to the boundary being violated.
- Whether
NDSolve
tracks the singular solution on the boundary or leaves it to follow a regular solution depends on which branch of the ODE $F=0$ is being integrated and the direction of integration.
- Rectifying (straightening) the singular solution (by transforming the ODE to the form
Y'[X]^2 == Y[X]
) aligns the steps taken by NDSolve
with the direction of the singular solution Y[X] == 0
. In this case, NDSolve
naturally follows the singular solution and does not leave it, provided the solution lies exactly on Y == 0
; otherwise, the effects of error will be similar to those in integrating ode
.
- If we integrate the contact structure to construct parametrizations of the integral curves,
NDSolve
seems to always integrate across the singularity at the criminant.
Fig. $\nextfig\ndsolve$. Integrating the ODE $F=0$ with an initial condition on the singular solution. NSolve
computes two solutions (dark green, dark red). The integral curve (light yellow) is traced until it meets the singular solution (light blue). Then it traces the discriminant, with increasing error the further from the initial condition. The values along the singular solution have a small imaginary component; only the real part is plotted.
solF = NDSolve[
{F[x, y[x], y'[x]] == 0, y[1] == 1}, y, {x, 0.5, 2},
"ExtrapolationHandler" -> {Indeterminate &, "WarningMessage" -> False},
Method -> "StiffnessSwitching" (* aids integrations along singular solution *)
];
Plot @@ {
Flatten[{1/x, 2 - x, Re[y[x] /. solF]}], Flatten[{x, 0.5, 2}],
PlotStyle -> {Directive[Thickness[0.015], Opacity[0.5]],
Directive[Thickness[0.015], Opacity[0.5]],
Directive[Darker[Red, 0.5], Opacity[0.8]],
Directive[Darker[Green, 0.66], Opacity[0.8]]},
GridLines -> {{1}, {1}}
}
The numerical integration of a related problem is discussed here: NDSolve solves this ordinary differential equation only "half-way"
DSolve
and the complete solution
Why is it that DSolve
misses a solution of this differential equation?
A little reflection will suggest why DSolve
does not solve for the singular solution.
Why singular solutions are missing in the first place
- Singular solutions are not a member or a limit of the general solution.
- A singular solution is a rare coincidence.
- The singular solution in the example fell out from algebraic equations obtained from differentiation and not from integrating differential equations.
- Using
Reduce
does not guarantee all solutions to the contact system will be found. (This is an observed fact. I do not understand why it happened or whether another approach might always be successful.)
Why they are not automatically sought by DSolve
- Singular solutions are rare, so one has to know when to look for them, or it will be mainly a waste of time.
- Reducing systems can take a lot of time. Is it really worth it if most of the time the result will be the null set?
How to be sure that DSolve
returned all possible solutions?
This is the OP's primary question. Leaving the difficulties of symbolic integration aside, there seem to be three steps to getting the complete answer (or as complete as possible in Mathematica.)
- Use
DSolve
to get the general solution.
- Check the limits of the general solution as the constants of integration approach the boundaries of their natural domains. (This is the issue with the ODE in DSolve not finding solution I expected.)
- Check the discriminant or the envelope of the general solution.
DSolve can find general solutions for linear and weakly nonlinear partial differential equations. Truly nonlinear partial differential equations usually admit no general solutions.
Since your equation is not linear there should be no expectation of finding all possible solutions. $\endgroup$