Suppose we want to solve a simple schrodinger equation under columb potential,namely,
\Psi''=(V-E)\Psi and V(r)=-\frac{1}{r},
The wavefunction u(r) of ground state is: \sim a\times r\cdot e^{-br} which should be the exact solution no matter how far the interval.But the numerical solution turns to be weired. For very large r_{max}, the wavefunction goes to large again for a certain point.
The codes goes as following:
using DifferentialEquations,Plots,Roots
function columbode!(du,u,p,r)
v(r)=-1/r
du[1]=u[2]
du[2]=(v(r)-p)*u[1]
end
function plotf(rmax)
u0=[0.0,1.0]
rspan=(1e-10,rmax)
prob(e)=ODEProblem(columbode!,u0,rspan,e)
solf(e)=solve(prob(e), Tsit5(), reltol=1e-8, abstol=1e-8)
e_eigen=find_zero(e->solf(e).u[end][1],(-1.0,-0.1))
plot(solf(e_eigen),vars=[1])
end
plot(plotf(50.0),plotf(60.0),plotf(70.0),plotf(80.0),layout=(2,2))
The result is as following:
Any idea on how to solve this problem?By the way, I don’t want to use callback to set the wavefunction to artificial zero after I need to find the root of wavefunction with respect to trial energy. Thanks!