Finite Square Well: Program 1

Particle in a Finite Potential Well

Consider a particle of mass $ m $ in a finite potential well $ V(x) $ defined by:

\[ V(x) = \begin{cases} -V_0 & \text{for } -\frac{L}{2} < x < \frac{L}{2} \\ 0 & \text{otherwise} \end{cases} \]

where $ V_0 $ is the depth of the potential well and $ L $ is the width of the well.

The Schrödinger equation governing the particle's wavefunction $ \Psi(x) $ is given by:

\[ \frac{d^2 \Psi}{dx^2} = \frac{2m}{\hbar^2} \left( V(x) - E \right) \Psi \]

where $ E $ is the energy eigenvalue of the particle.

Solve for the Energy Eigenvalues:
Define a function that uses the shooting method to determine the energy eigenvalues of the particle in the potential well. The shooting method iteratively finds energy values $ E $ such that the boundary conditions at $ x = \pm \frac{L}{2} $ are satisfied.

The eigenvalues are obtained by searching within a specified energy range $ E_0 \leq E \leq E_n $.
For each energy eigenvalue $ E $, compute the corresponding wavefunction $ \Psi(x) $ using the differential equation.

import numpy as np
from scipy.integrate import odeint, simps
from scipy.optimize import bisect
import matplotlib.pyplot as plt

def V(x):
    if(x > -L/2 and x < L/2):
        V = -V0
    else:
        V = 0
    return V

V0 = 10    
L, m, hbar, w, E0, En = 4, 1, 1, 1, -V0, 0
Erange = np.arange(E0, En, 0.1)
x = np.linspace(-L, L, 1000)
u = [0, 0.1]
    
def f(u, x, E):
    y, z = u
    f1, f2 = z, ((2*m)/(hbar)**2)*(V(x) -E)*y
    return (f1, f2)

def psiboundval(E):
    sol = odeint(f, u, x, args = (E, ))
    return sol[:, 0][-1]

def shoot(Erange): 
    global E
    Y = np.array([psiboundval(E) for E in Erange])
    eigval = np.array([bisect(psiboundval, Erange[i], Erange[i + 1])
                       for i in np.where(np.diff(np.signbit(Y)))[0]])
    return eigval

def psinorm(E):
    sol = odeint(f, u, x, args = (eigen, ))
    psi = sol[:, 0]
    return psi/np.sqrt(simps(psi*psi, x))

En = shoot(Erange)
print("Energy Eigen Values:", En)

P = np.array([V(i)for i in x])
plt.figure(figsize=(10,20))
st = 1
for eigen in En:
    psi = psinorm(eigen)
    normpsi = eigen + psi
    plt.subplot(len(En),2,st)
    plt.ylabel(f"$\Psi${st}(x)", fontsize=16)
    plt.xlabel(r'$x$', fontsize=16)
    plt.title(f'Energy level: n = {st}, Energy eigen value ={eigen: .4f}')
    plt.plot(x, normpsi)
    plt.plot(x, P)
    plt.axhline(y = eigen, color = 'g')
    plt.axvline(x = 0)
    st = st + 1
plt.savefig('psi.png')

TISE Solution

Finite Square Well: Program 1

Particle in a Finite Potential Well

Output 1