Deriving the IDE we get

$$
-n b\psi(b) + c_0 \psi''(b) + (b^2+k)\psi'(b) + 2b\psi(b) = 0
$$

or

$$
c_0 \psi''(b) + (b^2+k)\psi'(b)+b(2-n)\psi(b) = 0
$$
 Now if $n < 2,\ \ c_0 > 0,\ \ k > 0$ it looks as an stable ODE.

Follows a possible script to solve the complex case ($\psi = \psi_r+i\psi_i$) according to the comment.

    parms = {c0 -> 10, k -> 10, n -> 1, psir0 -> 1, dpsir0 -> 0, psii0 -> -1, dpsii0 -> 0};
    bmax = 10;
    sol = NDSolve[{c0 psir''[b] + (b^2 + k) psii'[b] + b (n - 2) psii[b] == 0, -c0 psii''[b] + (b^2 + k) psir'[b] + b (n - 2) psir[b] == 0, psir[0] == psir0, psir'[0] == dpsir0, psii[0] == psii0, psii'[0] == dpsii0} /. parms, {psir, psii}, {b, 0, bmax}]
    ParametricPlot[Evaluate[{psir[b], psii[b]} /. sol], {b, 0, bmax}]