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 Nov 9 '18 at 6:13

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, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.