Overflow with NDSolve

Update:

Alex noted a typo in my answer. Sorry about this. You can still solve this equation with the nonlinear FEM solver, thought it's not a as straight forward as for the miss typed equation.

Using

NDSolveValue[{D[\[Chi][x], {x, 2}] == \[Chi][x]^(3/2)/Sqrt[x], \[Chi][
    rin] == 1, \[Chi][10] == 0}, \[Chi], {x, 0, 10}, 
 Method -> "FiniteElement"]

gives an error message:

This is a fairly general failure message. The reason it is so general is that for the code it's impossible to say why it failed. One cause can the that there is a transition from the real to the complex plain. An easy way to try this is to either give an complex valued initial seed or, alternatively, to add a 0. I complex component to the equation:

sol = NDSolveValue[{D[\[Chi][x], {x, 2}] == \[Chi][x]^(3/2)/
      Sqrt[x], \[Chi][rin] == 1, \[Chi][10] == 0}, \[Chi], {x, 0, 10},
    Method -> "FiniteElement", InitialSeeding -> \[Chi][x] == 0. I];
Plot[Re[sol[x]], {x, 0, 10}, PlotRange -> All]

Note the Re in the plot. Now, there is also a small complex component in the solution:

Plot[Im[sol[x]], {x, 0, 10}, PlotRange -> All]

But the solution compares favorably to other solutions presented here:

Plot[Evaluate[psol[d /. drule][r] - Re[sol[r]]], {r, rin, 10}, 
 PlotRange -> All]

Old answer:

How about:

rin = 0;
sol = NDSolve[{D[\[Chi][x], {x, 2}] == (\[Chi][x]^3/2)/
      Sqrt[x] , \[Chi][rin] == 1, \[Chi][10] == 0}, \[Chi], {x, rin, 
    10}, Method -> "FiniteElement"];
Chi[x_] := Evaluate[\[Chi][x] /. sol[[1]]];
Plot[Chi[x], {x, 0, 10}]

If you are interested in why I reformulated the equations then this section from the documentation is a good starting point.

I'd like to extend my comments to an answer. For those in v12 or higher, FiniteElement is a possible choice for this problem, as shown in user21's answer. But, if you're in a version lower than v12 but higher than v9, it becomes a bit more troublesome, because

1. nonlinear FiniteElement isn't implemented yet.

2. Shooting method can't handle the problem well, which is an arguable backslide.

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 (* choose -1/10 if in v9 *)}] // Quiet
(* {d -> -0.0116574} *)

Plot[psol[d /. drule][r] // Evaluate, {r, rin, 10}, PlotRange -> All]

Alternatively, we can turn to finite difference method (FDM). I'll use pdetoae 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]

There is also wavelet method for BVP. It is an example with Haar wavelets. It takes 0.36 s to solve this problem with 64 colocation points:

ClearAll["Global`*"]
L = 10; A = 0; B = 1; J = 5; M = 
 2^J; dx = (B - A)/(2 M); 
h1[x_] := Piecewise[{{1, A <= x <= B}, {0, True}}]; 
p1[x_, n_] := (1/n!)*(x - A)^n;
h[x_, k_, m_] := 
  Piecewise[{{1, 
     Inequality[k/m, LessEqual, x, Less, (1 + 2*k)/(2*m)]}, {-1, 
     Inequality[(1 + 2*k)/(2*m), LessEqual, x, Less, (1 + k)/m]}}, 0];
p[x_, k_, m_, n_] := 
  Piecewise[{{0, x < k/m}, {(-(k/m) + x)^n/n!, 
     Inequality[k/m, LessEqual, x, 
      Less, (1 + 2*k)/(2*m)]}, {((-(k/m) + x)^n - 
        2*(-((1 + 2*k)/(2*m)) + x)^n)/n!, (1 + 2*k)/(2*m) <= 
      x <= (1 + k)/
       m}, {((-(k/m) + x)^n + (-((1 + k)/m) + x)^n - 
        2*(-((1 + 2*k)/(2*m)) + x)^n)/n!, x > (1 + k)/m}}, 0];
xl = Table[A + l dx, {l, 0, 2 M}]; xcol = 
 Table[(xl[[l - 1]] + xl[[l]])/2, {l, 2, 2 M + 1}];
f2[x_] := 
 Sum[af[i, j] h[x, i, 2^j], {j, 0, J, 1}, {i, 0, 2^j - 1, 1}] + 
  a0 h1[x]; 
f1[x_] := 
 Sum[af[i, j] p[x, i, 2^j, 1], {j, 0, J, 1}, {i, 0, 2^j - 1, 1}] + 
  a0 p1[x, 1] + f10; 
f0[x_] := 
 Sum[af[i, j] p[x, i, 2^j, 2], {j, 0, J, 1}, {i, 0, 2^j - 1, 1}] + 
  a0 p1[x, 2] + f10 x + f00;
bc1 = {f0[0] == 1};
bc2 = {f0[1] == 0};
var = Flatten[Table[{af[i, j]}, {j, 0, J, 1}, {i, 0, 2^j - 1, 1}]];
varM = Join[{a0, f10, f00}, var];
eqq[x_] := Sqrt[x] f2[x]/L^(3/2) - f0[x]^(3/2);
eq = Flatten[Table[{eqq[x] == 0}, {x, xcol}]];
eqM = Join[eq, bc1, bc2];
sol = FindRoot[eqM, Table[{varM[[i]], 0.1}, {i, Length[varM]}], 
  MaxIterations -> 1000]; lst = 
 Table[{L x, Evaluate[f0[x] /. sol]}, {x, 0, 1, .01}];

ListLinePlot[lst, PlotRange -> All]

Now we can compare with solution by xzczd psol[d /. drule][r], here we show difference f0[x]-psol[d /. drule][x] And as we can see, combination of ParametricNDSolveValue[] and FindRoot[] is still good