# Indeterminate solution of 2nd order nonlinear differential equation

I seek to prove a solution in power series in y or by numerical methods for the second order nonlinear differential equation:

nlde = {Sqrt[1 - f[y]^2] f''[y] + y == 0, f'[0] == 0, f[0] == 1};


It seems that there is an analytical answer in power series for this ODE (see appendix in this paper):

My failed attempts to prove the analytical solution with Mathematica:

Power series:

AsymptoticDSolveValue[nlde, f, {y, 0, 4}]
(* Indeterminate expression 0 ComplexInfinity encountered*)

Series[DifferentialRoot[Function @@ {{f, y}, nlde}][y], {y, 0, 4}]
(* supplied equation is not a
linear differential equation with polynomial coefficient *)


Numerical:

NDSolveValue[nlde, f, {y, 0, 4}]
(* Infinite expression encountered*)


Anyone has any suggestions on how to solve the equation? Thanks!

• The ode is singular for t=0! – Ulrich Neumann Feb 1 at 13:58
• If you look at the result of SolveAlways[With[{f = C[0] + Sum[C[k] y^k, {k, 1, 6}] + O[y]^7}, Sqrt[1 - f^2] D[f, {y, 2}] + y == 0], y], you might get an idea on what constraints your series coefficients have to follow for a power series to satisfy your nonlinear ODE. – J. M.'s ennui Feb 1 at 14:03
• Yes, it is singular. Somehow there appears to be an analytic solution in power series (see edited question). The approach of using SolveAlways doesn't seem to give the same coefficients as the expected results – Ferca Feb 1 at 14:15
• Using the excerpt you have shown, you could use SolveAlways[] with formula A6: 1 + Sum[C[k] y^k, {k, 2, 10}] /. With[{g = -Sum[C[k] y^(k - 2), {k, 2, 10}] + O[y]^(11 - 2)}, Last[SolveAlways[g (2 - y^2 g) (y^2 D[g, {y, 2}] + 4 y D[g, y] + 2 g)^2 == 1, y]]] – J. M.'s ennui Feb 2 at 2:19
• Yes, your approach is correct. It seems that SolveAlways is picking even the other solutions for the imaginary coefficients! – Ferca Feb 2 at 17:53

If we rationalize the DE, we can get an asymptotic solution:

nlde2 = {(1 - f[y]^2) f''[y]^2 - y^2 == 0, f'[0] == 0, f[0] == 1};
asol = AsymptoticDSolveValue[nlde2, f, {y, 0, 8}]


AsymptoticDSolveValue::asdb: There are multiple solution branches for the equations, but AsymptoticDSolveValue will return only one.

1 - y^2/2 - y^4/104 - (67 y^6)/83824 - (49463 y^8)/496908672


It returns a solution that satisfies the original nlde for y >= 0:

sol = {f -> Function[y, Evaluate[asol + O[y]^(Exponent[asol, y] + 1)]]}

First@nlde /. Equal -> Subtract /. sol
Simplify[%, y >= 0]


### Alternative

Implementing the method in the paper (sans rationalization):

nlde3 = nlde /. f -> Function[y, 1 - y^2 g[y]] //
FullSimplify[#, y > 0] &

1 - y^2 AsymptoticDSolveValue[{First@nlde3, Reduce[nlde3 /. y -> 0]},
g, {y, 0, 6}] // Expand

• I think the point is that Sqrt[] is not analytic at 0, which is where Sqrt[1 - f[y]^2] starts. – Michael E2 Feb 1 at 14:31
• Thanks! It seems the trick is to rationalise the ODE with the square root term. The power series solution works now. The numerical approach with NDSolveValue[nlde2, f, {y, 0, 4}] still complains about indeterminate expressions encountered. – Ferca Feb 1 at 15:01

A direct way to get the solution, known from perturbation theory, follows with assumptions f=1+y+...and y>0

cl = CoefficientList[Simplify[Normal[Series[Sqrt[1 - f[y]^2] f''[y] + y
/.f -> Function[y, 1 + Sum[C[k] y^k, {k, 2, 6}]], {y, 0, 5}]],y > 0], y]

sol=Solve[0 == cl , {C[2], C[3], C[4], C[5], C[6]}]
(*{{C[2] -> -(1/2), C[3] -> 0, C[4] -> -(1/104), C[5] -> 0,C[6] -> -(67/83824)}}*)

Function[y, 1 + Sum[C[k] y^k, {k, 2, 6}]][y]/.sol
(*1 - y^2/2 - y^4/104 - (67 y^6)/83824*)

• Thanks! Your approach works. It seems to be critical to set the two first terms 1+0 and then add the coefficients for this work. – Ferca Feb 1 at 15:15
• @Ferca The method works for the general case C[0],C[1] too. Analyzing the coefficientlist let you derive conditions for the existence of a power serie. – Ulrich Neumann Feb 1 at 18:38
• cl = CoefficientList[ Simplify[ Normal[Series[ Sqrt[1 - f[y]^2] f''[y] + y /. f -> Function[y, C[0] + Sum[C[k] y^k, {k, 1, 6}]], {y, 0, 5}]], y > 0], y]; sol = Solve[0 == cl, {C[0], C[1], C[2], C[3], C[4], C[5], C[6]}]; Function[y, C[0] + Sum[C[k] y^k, {k, 1, 6}]][y] /. sol; % /. C[0] -> 1 leads to Indeterminate – Ferca Feb 1 at 19:40
• @Ferca No , you have to solve it step by step: First look at cl[[1]] which vanishs for C[0]->1 (or C[2]->0) , Second set C[0]=1 and look at cl[[1]] which only vanishs for C[1]->1. And so on... – Ulrich Neumann Feb 2 at 6:58
• @Ferca In this way you might find a second solution for arbitrary C[0],C[1] : Solve[0 == cl[[1 ;; 3]], {C[2], C[3], C[4]}] // Simplify (*{{C[2] -> 0, C[3] -> -(1/(6 Sqrt[1 - C[0]^2])), C[4] -> -((C[0] C[1])/(12 (1 - C[0]^2)^(3/2)))}}*) – Ulrich Neumann Feb 2 at 7:31