Explanation of Error Noted in Question
The computation in the question fails, because there is no real value of F'[0.01]
for F[0.01] == 0.01
and gamma == -0.9
. This can be seen by solving the ODE for F[0.01]
as a function of F'[0.01]
.
sf = f /. Flatten@Solve[((1 + F'[x])^2*(1 - (1 + F'[x])^gamma)^2/gamma^2 ==
x*F'[x] - F[x]) /. gamma -> -0.9 /. F'[x] -> fp /. F[x] -> f /. x -> 0.01, f]
(* 0.01 fp - 1.23457 (1. + fp)^2 (1. - 1./(1. + fp)^0.9)^2 *)
FindMaximum[sf, fp]
(* {0.0000249876, {fp -> 0.00499627}} *)
So, choose F[0.01]
less than that maximum, for instance,
s = ParametricNDSolveValue[{(1 + F'[x])^2*(1 - (1 + F'[x])^gamma)^2/
gamma^2 == x*F'[x] - F[x], F[0.01] == 10^-5}, F, {x, 0.01, 100}, {gamma},
Method -> {"EquationSimplification" -> "Residual"}]
Plot[s[-.9][x], {x, 0.01, 100}]

(It is convenient to use ParametricNDSolveValue
for problems like this.)
General Linear Solution
Differentiating the ODE reveals that there are two categories of solutions.
der2 = Collect[Subtract @@ D[(1 + F'[x])^2*(1 - (1 + F'[x])^gamma)^2/gamma^2 ==
x*F'[x] - F[x], x], F''[x], Simplify]
(* (-x + (2 (1 + F'[x])^(1 + gamma) (-1 + (1 + F'[x])^gamma))/gamma +
(2 (1 + F'[x]) (-1 + (1 + F'[x])^gamma)^2)/gamma^2) F''[x] *)
Solutions given by F''[x] == 0
are, of course, linear, F[x] == a x + b
. Insert this into the ODE to obtain b
as a function of a
.
sb = -First[((1 + F'[x])^2*(1 - (1 + F'[x])^gamma)^2/gamma^2 ==
x*F'[x] - F[x]) /. F[x] -> a x + b /. F'[x] -> a]
(* ((1 + a)^2 (1 - (1 + a)^gamma)^2)/gamma^2 *)
ParametricPlot[{sb /. gamma -> -0.9, a}, {a, -1, 1}, AxesLabel -> {"F[0]", "F'[0]"},
ImageSize -> Large, LabelStyle -> Directive[12, Black, Bold],
AspectRatio -> 1/GoldenRatio]

The linear solutions are a one-dimensional family, parameterized by a
or b
. The solution obtained in the previous section is an example.
Nonlinear Solutions
The other factor of der2
is an algebraic expression for F'[x]
.
sx = First[der2]
(* -x + (2 (1 + F'[x])^(1 + gamma) (-1 + (1 + F'[x])^gamma))/gamma +
(2 (1 + F'[x]) (-1 + (1 + F'[x])^gamma)^2)/gamma^2 *)
ParametricPlot[{x + sx /. F'[x] -> fp /. gamma -> -0.9, fp}, {fp, -1, 1},
AxesLabel -> {x, "F'[x]"}, ImageSize -> Large,
LabelStyle -> Directive[12, Black, Bold], AspectRatio -> 1/GoldenRatio]

For x > -1
there are two solutions. Visibly, F'[0] == 0
for one of them, and the corresponding value of F[0]
is given by
((1 + F'[x])^2*(1 - (1 + F'[x])^gamma)^2/gamma^2 ==
x*F'[x] - F[x]) /. x -> 0 /. gamma -> -0.9 /. F'[0] -> 0
(* 0. == -F[0] *)
With these values, the solution is given by
fprime[gg_?NumericQ, xx_?NumericQ] :=
fp /. FindRoot[sx /. {F'[x] -> fp, gamma -> gg, x -> xx}, {fp, xx/2}]
Flatten@NDSolve[{F'[x] == fprime[-.9, x], F[10^-4] == 10^-8/4}, F, {x, 10^-4, 100}];
Plot[F[x] /. %, {x, 10^-4, 100}, AxesLabel -> {x, F}, ImageSize -> Large,
LabelStyle -> Directive[12, Black, Bold]]

F[0]
is determined for the second solution by
((1 + F'[x])^2*(1 - (1 + F'[x])^gamma)^2/gamma^2 ==
x*F'[x] - F[x]) /. x -> 0 /. gamma -> -0.9
(* 1.23457 (1 + F'[0])^2 (1 - 1/(1 + F'[0])^0.9)^2 == -F[0] *)
F'[x] /. FindRoot[0 == sx + x /. gamma -> -0.9, {F'[x], -.8, -.9}];
f01 = -First[%%] /. F'[0] -> %
(* -0.599484 *)
fprime1[gg_?NumericQ, xx_?NumericQ] :=
fp /. FindRoot[sx /. {F'[x] -> fp, gamma -> gg, x -> xx}, {fp, -.8, -.9}]
Flatten@NDSolve[{F'[x] == fprime1[-.9, x], F[10^-4] == f01}, F, {x, 10^-4, 100}];
Plot[F[x] /. %, {x, 10^-4, 100}, AxesLabel -> {x, F}, ImageSize -> Large,
LabelStyle -> Directive[12, Black, Bold]]

It is nearly linear, because fprime1[-0.9, x]
is nearly constant.
NDSolve::ivcon: The given initial conditions were not consistent with the differential-algebraic equations. NDSolve will attempt to correct the values.
$\endgroup$