numerov.procedure: Numerov's procedure
In npetraco/che302r: Handy functions for CHEM 302, Physical Chemistry II

numerov.procedure

R Documentation

Numerov's procedure

Description

Numerov's procedure to solve simple 1D Schodinger equations

Usage

numerov.procedure(
  xaxis,
  dxaxis,
  PE.function.name,
  nodes.in.state,
  E.guess,
  num.iterations,
  delay.time = 0,
  ...
)

Arguments

`xaxis`	The x-axis over which to find a solution to the Schrodinger equation with the requested potential
`dxaxis`	The "differential" increment for the x-axis
`PE.function.name`	Potential energy function name. Choices are: "box", "harmonic", "anharmonic", "radial"
`nodes.in.state`	Number of nodes that should be in the wave function.
`E.guess`	Guess for the energy of the wave function
`num.iterations`	Number of iterations to run the procedure
`delay.time`	Set this if you want to slow down the procedure and watch the solution evolve over the iterations.

Details

Numerov's procedure to solve simple 1D Schrodinger equations. Works for any of the above stated potential energy functions. See Ira Levine's Quantum Chemistry for a very nice explanation of the algorithm.

Value

A two column matrix. The first column is the x-axis over which the wave function was obtained. The second column is the un-normalized wave function solution.

References

Quantum Chemistry by Ira Levine. Any edition.

Examples

# Solve the Schrodinger equation for particle in a box:
#Domain (x-axis):
dx    <- 0.01                           # Increment on x-axis
x.min <- (0)                            # Start of x-axis. Well into (left) classically forbidden region
x.max <- 1.5                            # End of x-axis. Well into the (right) classically forbidden region
x     <- seq(from=x.min,to=x.max,by=dx) # x-axis

#Potential energy function
potential<-"box"
plot(x,V(x,potential))

#Numerov's procedure:
state        <- 0                # State you want, starting from 0. It is an integer: 0,1,2,3,...etc. I.E., number of nodes
Guess.Energy <- (0)              # Initial guess for energy of the state
max.iter     <- 40
psi.info     <- numerov.procedure(x,dx,potential,state,Guess.Energy,max.iter,delay.time=0.00)

#Approximately Normalize Psi:
Npsi.info <- approx.normalize(psi.info)
plot(Npsi.info[,1], Npsi.info[,2], xlab="x", ylab="N*psi(x)", typ="l")


# Solve the Schrodinger equation for radial wave functions (Basically Laguerre functions):
#Domain (x-axis):
dr    <- 0.1                            # Increment on r-axis
r.min <- 1e-15                          # Start of r-axis.
r.max <- 180                            # End of r-axis.
r     <- seq(from=r.min,to=r.max,by=dr) # r-axis

#Numerov's procedure to solve the (almost) Radial Schrodinger Equation:
n            <- 6      # n = 1,2,3,...
l            <- 4      # l = 0, ... , n-1 (e.g. s,p,d,f,g,.....)
Guess.Energy <- (-0.5) # Initial guess for energy of the state
max.iter     <- 50
state        <- (n-l-1)
Fr.info      <- numerov.procedure(r,dr,"radial",state,Guess.Energy,max.iter,0.00)

#Transform F(r) back into R(r), the proper radial wave function
psi.info     <- Fr.info
psi.info[,2] <- (r^-1)*Fr.info[,2]
plot(r, r^2 * psi.info[,2],typ="l")

#Approximately Normalize Psi:
Npsi.info <- approx.normalize(psi.info, include.jacobian = TRUE)

#Plot the normalized wave function:
plot(r, r^2*Npsi.info[,2], xlim=c(min(r),max(r)), ylab="N psi(x)", typ="l")

npetraco/che302r documentation built on April 17, 2025, 10:34 p.m.