# Numerical Solution to Schrodinger Equation w/ Coulomb Potential

• Tertius
Tertius
Homework Statement
I am trying to numerically solve the time-independent Schrödinger equation for the tunneling problem encountered in fusion processes (with a Coulomb barrier), specifically avoiding the use of the WKB approximation.
Relevant Equations
Schrodinger Equation
Coulomb Potential
I am doing this to have my own solution for customization and understanding. I also want to manually check the WKB approximation accuracy at various energies against this static solution.

I've split the problem into 3 regions and am solving it in 1D, but am having problems with how to define the initial boundary of Region 1.
$$r_{ctp}$$ is the classical turning point where incoming kinetic energy E equals coulomb potential V.
$$r_{nuc}$$ is the point where the nuclear strong force takes over at approximately the nuclear radius.

Region 1:
$$-\infty < x <= r_{ctp}$$
Region 1 has some of V(x) before reaching ##r_{ctp}##.
Equation:
$$-\frac{\hbar^{2}}{2m} \frac{d^{2}}{dx^{2}} \psi_I(x) +V(x)\psi_I(x)=E\psi_I(x)$$
BCs:
$$\psi_I(-100) = 0$$
$$\psi_I'(-100) = 0.0001$$ (this is nonsense, and it propagates through the solutions).
If it is set to 0, then the solution is trivially 0. But the wave function should approach 0 going to ##-\infty##.

Region 2:
$$r_{ctp} < x <= r_{nuc}$$
Region 2 has negative kinetic energy because E < V.
Equation:
$$-\frac{\hbar^{2}}{2m} \frac{d^{2}}{dx^{2}} \psi_{II}(x) +V(x)\psi_{II}(x)=E\psi_{II}(x)$$
BCs:
$$\psi_I(r_{ctp}) = \psi_{II}(r_{ctp})$$
$$\psi_I'(r_{ctp}) = \psi_{II}'(r_{ctp})$$

Region 3:
$$r_{nuc} < x <= 0$$
Region 3 is 'freely propagating' where ##V(x)=0##
Equation:
$$-\frac{\hbar^{2}}{2m} \frac{d^{2}}{dx^{2}} \psi_{III}(x) =E\psi_{III}(x)$$
BCs:
$$\psi_{II}(r_{nuc}) = \psi_{III}(r_{nuc})$$

$$\psi_{II}'(r_{nuc}) = \psi_{III}'(r_{nuc})$$

This setup is producing the attached graph using NDSolve in Mathematica.

The two questions i have are:
1) what are the appropriate BCs for Region 1?
2) why is the Region 2 solution oscillatory when it should be exponential? Usually that's a minus sign problem but wondering if there is something else going on.

Thanks

#### Attachments

• tunneling v1.jpeg
21.9 KB · Views: 6
Last edited:
I figured out where my approach was wrong. The boundary conditions for this are a little tricky because it is not a bound system.
So I discretized it manually and just solved $$H\psi = E\psi$$, where E is the kinetic energy of the incoming particle. E in this case is just a single eigenvalue because it is not a bound system, so it is possible to solve it as a simplified linear system as $$(H-EI)\psi = 0$$. This made it fast and give the same answer. But I did need to give the LinearSolver a small non-zero initial wave function value so it didn't give the trivial answer. Hope this helps someone.

jonjacson

Replies
5
Views
146
Replies
4
Views
554
Replies
6
Views
1K
Replies
5
Views
1K
Replies
1
Views
1K
Replies
11
Views
1K
Replies
2
Views
2K
Replies
1
Views
857