I'd like to extend my comments to an answer. First of all, for those in _v12_ or higher, `FiniteElement` is the best choice for this problem, as shown in [user21's answer](https://mathematica.stackexchange.com/a/219608/1871). But, if you're in a version lower than _v12_ but higher than _v9_, it becomes a little troublesome, because 1. nonlinear `FiniteElement` isn't implemented yet. 2. `Shooting` method can't handle the problem well, which is [an arguable backslide](https://mathematica.stackexchange.com/a/130373/1871). [![enter image description here][1]][1] _As we can see, though `ndsz` warning is generated, `NDSolve` manages to find the desired result in v9._ OK, so what to do? Well, to be honest I don't know if the following solution will cause other problem in _v10.0_, because _v10.0_ is a quite unstable version, but it does work in _v9_ and _v12.1_: rin = 10^-30; psol = ParametricNDSolveValue[{D[χ[x], {x, 2}] == (χ[x]^3/2)/Sqrt[x], χ'[ 10] == d, χ[10] == 0}, χ, {x, rin, 10}, d] drule = FindRoot[psol[d][rin] == 1, {d, 0}] // Quiet Plot[psol[d /. drule][r] // Evaluate, {r, rin, 10}] Alternatively, we can turn to finite difference method (FDM). I'll use [`pdetoae`](http://mathematica.stackexchange.com/a/127997/1871) for the generation of finite difference equations: rin = 0; eq = D[χ[x], {x, 2}] Sqrt[x] == (χ[x]^3/2); bc = {χ[rin] == 1, χ[10] == 0}; points = 25; domain = {rin, 10}; grid = Array[# &, points, domain]; difforder = 2; del = #[[2 ;; -2]] &; (* Definition of pdetoae isn't included in this post, please find it in the link above. *) ptoafunc = pdetoae[χ[x], grid, difforder]; ae = ptoafunc@eq // del; initialguess[x_] = 0; solrule = FindRoot[{ae, bc}, Table[{χ[x], initialguess[x]}, {x, grid}]] sol = ListInterpolation[solrule[[All, -1]], grid] [1]: https://i.sstatic.net/eKZrT.png