I have the following ODE for a function $F(x)$: $F''-\frac{1}{x}F'-aF=0$ with the following boundary conditions: $F(x\to0) = 1$, $F(x\to\infty)=0$. It can be solved analytically: $F = \sqrt{a}xK_1(\sqrt{a}x)$, where $K_1$ is a modified Bessel function. However, I want to solve it numerically, because my original equation is more difficult. The question is: how to specify boundary conditions?

I tried to implement this solution. However, it seems that it won't work, because it gives Power::infy: Infinite expression 1/0. encountered error, which, I suppose, is connected with the fact that $F'(x) \sim \frac{1}{x}$ when $x\to 0$ (this limit can be obtained from the exact solution).

So, how to specify boundary conditions in this case?

  • $\begingroup$ What about changing variable x -> 1/(x+1), solving on [0,1] and changing variables back? $\endgroup$
    – Ihor
    Commented Nov 9, 2018 at 6:13

1 Answer 1


You can find an approximative solution for your problem as follows. It's some kind of shootingmethod (see NDSolve[...,Method-> "Shooting"]:

The initial slope Fs is adapted to force F[T]==0

T = 10;                  (* "infinity" *)
\[CurlyEpsilon] = .0001; (*avoid singularity x=0 *)

fa = ParametricNDSolveValue[{  F''[x] - F'[x]/x - a F[x] == 0,F[\[CurlyEpsilon]] == 1, F'[\[CurlyEpsilon]] == Fs},F, {x, \[CurlyEpsilon], T}, {Fs, a  },Method -> "StiffnessSwitching"]   

In the next step the optimal slope for different ais calculated

fit = Table[{fs, a} /. FindRoot[fa[fs, a][T] ==  0, {fs, -\[CurlyEpsilon]}], {a, 1/2, 5,1/2}]

and shown

Plot[Map[ fa[#[[1]], #[[2]]] [x] &, fit], {x, \[CurlyEpsilon], T},PlotRange -> {0, 1}]

enter image description here


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.