TISE Solution: Shooting Method to Solve TISE - Computational Physics Formulation & Simulation

Program 1

Shooting Method to Solve TISE

Python 3.11 Ready

'''Solve 1-dimensional time-independent Schroedinger Equation (TISE) for
particle in 1D box'''

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

def V(x):
    return 0

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

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

#Identify eigenvalues (where psi(L) crosses zero)
Y = np.array([psiboundval(E) for E in Erange])
plt.plot(Erange, Y, label="$\Psi(L)$")
plt.axhline(y=0, color='red', linestyle='--', label="Boundary Condition $\Psi(L)$=0")
plt.ylabel(r"$\Psi(L)$")
plt.xlabel(r"$E$")
plt.title(r"$\Psi(L) vs. E$")
plt.legend()
plt.savefig("plot.png")
plt.show()

# Determine eigenvalues by Bisection method
En = []
for i in range(len(Y) - 1):
    if np.sign(Y[i]) != np.sign(Y[i + 1]):  # Check for sign change
        E_eigen = bisect(psiboundval, Erange[i], Erange[i + 1])
        En.append(E_eigen)
print("Energy Eigen Values:", np.round(En,4))

Click "Run Code" to execute script and render figures.

Mathematical Problem Formulation

Shooting Method to solve TISE

Consider a particle of mass $ m $ confined in a one-dimensional box of length $ L $. The potential energy $ V(x) $ inside the box is zero, while it is infinite outside the box. The goal is to solve the time-independent SchrÃ¶dinger equation (TISE) for this system, which is given by:

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

where:

$ \hbar $ is the reduced Planck's constant,
$ E $ is the energy of the particle,
$ \Psi(x) $ is the wave function.

Implement a numerical solution using the shooting method to find the energy eigenvalues of the particle in the box.

Theoretical Background & Computational Mechanics

Shooting Method

The shooting method is a numerical technique to solve boundary value problems by converting them into initial value problems. In this code:

The TISE is transformed into a system of first-order ordinary differential equations (ODEs).
The boundary conditions are specified (i.e., $ \Psi(0) = 0 $ and $ \Psi(L) = 0 $).

The differential equations derived from the TISE are:

Let $ \Psi(x) = y $ (the wave function).

The first derivative $ \frac{d \Psi}{dx} = z $.

Thus, the system becomes:

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

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

This can be represented in the code as:

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

u[0] = 0: This sets the initial value of the wave function $ \Psi(x) $ at the starting position (which is usually $ x = 0 $ in a 1D box). Setting $ \Psi(0) = 0 $ means that the wave function starts at zero at this boundary.

u[1] = 0.01: This sets the initial value of the first derivative of the wave function $ \frac{d \Psi}{dx} $ at the starting position $ x = 0 $. A small positive value (0.01) is chosen here, which indicates that the wave function has a slight slope at the beginning. This is necessary to ensure that the wave function can evolve from this point.

The odeint function from SciPy is used to solve the ODEs numerically:

sol = odeint(f, u, x, args=(E, ))

Boundary Value Calculation

The function psiboundval(E) calculates the wave function at the right boundary $ x = L $:

def psiboundval(E):
    sol = odeint(f, u, x, args=(E, ))
    return sol[:, 0][-1]  # Return the wave function at x = L

Visualization

Finally, the code plots $ \Psi(L) $ against the energy values $ E $:

plt.plot(Erange, Y, label="$\Psi(L)$")
plt.axhline(y=0, color='red', linestyle='--', label="Boundary Condition $\Psi(L)=0$")

Eigenvalue Determination

Identifying Zeros

The eigenvalues are identified by locating the energy values $ E $ for which the wave function $ \Psi(L) $ crosses zero. The crossing points correspond to valid eigenvalues of the system.

Bisection Method

The code uses the bisection method to find these eigenvalues, where it checks for changes in the sign of $ \Psi(L) $:

for i in range(len(Y) - 1):
    if np.sign(Y[i]) != np.sign(Y[i + 1]):
        E_eigen = bisect(psiboundval, Erange[i], Erange[i + 1])
        En.append(E_eigen)

This method effectively narrows down the interval in which an eigenvalue exists by looking for points where the function changes from positive to negative (or vice versa).

Theoretical Principles & Computational Physics Guide

Shooting Method to Solve TISE: Quantum Energy Eigenvalue Search

Dr. Alorika Chatterjee & Dr. Anirban Shaw Science & Physics Learners Python 3 / NumPy / SciPy

1. Core Principles & Theoretical Framework

Eigenvalue Quantization: Transforms the boundary value problem into an initial value problem parameterized by trial energy $E$.

Target Matching: Integrates from $x = 0$ to boundary $x = L$ and searches for trial energies that satisfy $\psi(L; E) = 0$.

Node Counting: The $n$-th excited bound state contains exactly $n$ internal nodes in its wavefunction.

2. Mathematical Formulation & Algorithmic Pipeline

1D Time-Independent Schrödinger Equation (TISE):

$$ -\frac{\hbar^2}{2m} \frac{d^2\psi}{dx^2} + V(x)\psi = E\psi \implies \frac{d^2\psi}{dx^2} = -\frac{2m}{\hbar^2}\big(E - V(x)\big)\psi(x) $$

3. Numerical Mechanics & Implementation Logic

Trial Bracket: Choose energy bracket $[E_1, E_2]$ such that boundary values $\psi(L; E_1)$ and $\psi(L; E_2)$ have opposite signs.
ODE Integration: Integrate the TISE using RK4 or Numerov's method from $x=0$ to $x=L$ for test energy $E_{\text{mid}}$.
Bisection Convergence: Halve the energy bracket until $|E_2 - E_1| < 10^{-6}\text{ eV}$.
Normalization: Scale $\psi(x)$ such that $\int_0^L |\psi(x)|^2 dx = 1$.

4. Relevant Computational Practice Exercises

Exercise 1

Compute the ground state and first two excited state energies for an electron in an infinite square well of width $L = 1\text{ nm}$.

Exercise 2

Verify that the computed eigenvalues match the analytical formula $E_n = \frac{n^2 \pi^2 \hbar^2}{2m L^2}$ within $0.01\%$.

Exercise 3

Plot the normalized probability density $|\psi_n(x)|^2$ for $n = 1, 2, 3$.

Standard Literature & References:

David J. Griffiths, "Introduction to Quantum Mechanics", Cambridge University Press (3rd ed., 2018).
Mark Newman, "Computational Physics", University of Michigan (2013).

Ch. 19.1 Shooting Method to Solve TISE