Responsive Navbar with Google Search
TISE Solution: Particle in 1D Potential Box

TISE Solution

Particle in 1D Potential Box: Program 1

Particle in a One-Dimensional Box

Consider a particle confined in a one-dimensional box with length \( L \). The 1-dimensional time-independent Schrödinger equation (TISE) for a particle in a potential \( V(x) \) is given by:

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

For a particle in an infinite potential box, the potential \( V(x) \) inside the box is zero:

\[ V(x) = 0 \quad \text{for } 0 \leq x \leq L \]

The boundary conditions are:

\[ \Psi(0) = 0 \quad \text{and} \quad \Psi(L) = 0 \]

  1. Determine the Energy Eigenvalues:

    Use the shooting method to find the eigenvalues \( E \) of the particle in the box. Search for eigenvalues in the range \( E_0 \leq E \leq E_n \).

  2. Normalize the Wavefunctions:

    For each eigenvalue \( E \), compute the corresponding wavefunction \( \Psi(x) \) by solving the TISE.

    Normalize the wavefunction \( \Psi(x) \) such that:

    \[ \int_0^L |\Psi(x)|^2 \, dx = 1 \]

Output 1

Quantum Mechanics: Particle in a Box

1. Defining the Potential Function

Python code:

def V(x):
    return 0

Since the potential \( V(x) = 0 \) inside the box, this function simply returns 0 for any input \( x \).

2. Setting Parameters and Initial Setup

Python code:

L, m, hbar, w, E0, En = 10, 1, 1, 1, 0.001, 5
plt.figure(figsize=(10, 20))
Erange = np.arange(E0, En, 0.1)
x = np.linspace(0, L, 100)
u = [0, 0.01]

Here we define:

  • \( L \): Length of the box.
  • \( m \): Mass of the particle (set to 1 for simplicity).
  • \( \hbar \): Reduced Planck's constant (set to 1 for simplicity).
  • \( E_0 \) and \( E_n \): Initial and final values for the energy range.
  • \( x \): Array of positions within the box, with 100 points from \( 0 \) to \( L \).
  • \( u \): Initial conditions: \( y(0) = 0 \) and \( y'(0) = 0.01 \). The small initial slope ensures we start solving the differential equation from a non-zero slope.

3. Defining the System of ODEs

Python code:

def f(u, x, E):
    y, z = u
    f1, f2 = z, ((2 * m) / (hbar ** 2)) * (V(x) - E) * y
    return (f1, f2)

This function defines the system of first-order ODEs by breaking the TISE into two equations:

\[ f_1 = \frac{dy}{dx} = z \]

\[ f_2 = \frac{dz}{dx} = \frac{2m}{\hbar^2} (V(x) - E) y \]

Here, \( y \) represents the wavefunction \( \Psi(x) \), and \( z \) represents \( \frac{d \Psi}{dx} \).

4. Defining the Boundary Condition Checker Function

Python code:

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

This function uses odeint to integrate the system from \( x = 0 \) to \( x = L \) for a given energy \( E \). It returns the value of \( \Psi(L) \) i.e. the right most value of \( \Psi(x) \). If \( \Psi(L) \approx 0 \), \( E \) is likely an eigenvalue. Therefore, we have to determine the enrgy values for which \( \Psi(L) \approx 0 \) and those energies will be the eigen values.

5. Shooting Method for Finding Eigenvalues

Python code:

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

This function applies the shooting method to find energy eigenvalues:

  • \( Y \) stores the values of \( \Psi(L) \) for each energy \( E \) in Erange.

  • The code np.where(np.diff(np.signbit(Y)))[0] is used to identify indices in the array Y where a sign change occurs between consecutive values.

    This is important because a sign change in the boundary value of the wavefunction \( \Psi(L) \) (stored in Y) suggests the presence of an eigenvalue between those points in the energy range.

  • bisect then finds the precise eigenvalue between each detected sign change by finding the zero crossing of psiboundval.

6. Normalizing and Plotting the Eigenfunctions

Python code:

En = np.round(shoot(Erange), 4)
print("Energy Eigen Values:", En)

st = 1
for eigen in En:
    sol = odeint(f, u, x, args=(eigen,))
    psi = sol[:, 0]
    normpsi = psi / np.sqrt(simps(psi * psi, x))
    plt.subplot(len(En), 2, st)
    plt.ylabel(r"$\Psi$"+str(st)+"(x)", fontsize=16)
    plt.xlabel(r'$x$', fontsize=16)
    plt.title('Energy level: n= '+str(st)+', Energy eigen value ='+str(eigen))
    plt.plot(x, normpsi)
    plt.grid()
    st += 1
plt.tight_layout()
plt.savefig('psi.png')
plt.show()

Here:

  • \( E_n \) stores the energy eigenvalues found by shoot.
  • Each eigenvalue corresponds to a particular eigenstate, or wavefunction.
  • odeint(f, u, x, args=(eigen, )) solves the TISE for each eigenvalue.
  • Normalization: \( \text{normpsi} = \psi / \sqrt{\int \psi^2 dx} \) normalizes each wavefunction by dividing it by the square root of its integral (using simps for integration).
  • The eigenfunctions are then plotted, labeled by energy level and eigenvalue.