Now that we have the functionality to solve the 1 dimensional Schrödinger equation numerically, all we need to do is to give it something to solve by specifying the potential energy Vfun! For the first example we will explore the quantum simple harmonic oscillator. # obtain neigs solutions from the sparse matrix # create the Hamiltonian Operator matrix: X = np.linspace(xmin, xmax, Nx) # x axis grid Note that the function schrodinger1D has the following arguments: xmin = minimum x grid value, xmax = maximum x grid value, Nx = number of grid points, Vfun = potential energy function object, neigs = number of eigenvalues to solve for, and findpsi = flag to determine what the function returns.
#Schrödinger code#
The code to construct H and calculate its eigenvectors and eigenvalues is presented below. We then normalise the eigenvectors at the end. We also choose to obtain the smallest magnitude solutions with which='SM'. We choose the number of solutions to obtain by setting the value of k=neigs. Also using the sparse matrix eigenvalue calculators allow us to control the number of solutions to find, rather than calculating every single possible solution which will take a lot of time and resources.Īfter constructing H we then calculate its eigenvectors and eigenvalues using = sla.eigs(H, k=neigs, which='SM'). the matrix is sparse), and using the sparse matrix libraries will help to speed up computation for extremely large systems. This is because the Hamiltonian matrix H is composed mainly of 0s (i.e. Instead of using standard linear algebra libraries, we use the sparse matrix linear algebra library. This should work for most situations in general, but might break down for situations which require other forms of boundary conditions such as periodic boundary conditions. Note that setting up the Hamiltonian matrix in such a manner results in Dirichlet boundary conditions being intrinsically applied to the numerical system. We then add to this the potential energy V to complete the Hamiltonian matrix H. Therefore if we scale the values of V and E by a factor of 2, the Hamiltonian matrix H takes the form of a tridiagonal matrix with 2 /dx 2 in the main diagonal, and -1 /dx 2 in the first diagonals above and below the main diagonal. Note that we have set ħ = m = 1 for simplicity’s sake, and that dx is the step size of the spatial grid used.
By using second order central finite differences for 1 dimension, the wave function is differentiated numerically as: -1/(2 dx 2) as can be seen in the discretised equation below.
#Schrödinger how to#
The first step is to determine how to implement the derivative numerically. Higher dimensions will be treated next time! Also for this post we will solve the Schrödinger equation numerically for 1 dimension only. For simplicity’s sake, we will assume that the Planck constant ħ and particle mass m are both equal to 1. Therefore, if we can determine the form of the matrix H, we can then numerically solve the system to get the both quantum mechanical wave function Ψ, as well as the corresponding total particle energy E.Īs can be seen in the equation at the top of this post, the matrix H is the sum of a second order derivative (physically this is the kinetic energy) and a potential energy function V. Physically this means that operating the Hamiltonian which is the sum of kinetic and potential energies returns the total energy of the particle! For the less mathematically inclined, this relation simply means that multiplying matrix H on vector Ψ gives the same result as multiplying the scalar value E on vector Ψ. This means that it can be cast in matrix form as: HΨ = E Ψ, where H is the Hamiltonian matrix (the Hamiltonian is essentially the sum of a particle’s kinetic and potential energies), Ψ is the wave function vector and E is the energy eigenvalue. The time independent Schrödinger equation is an eigenvalue problem. I will put some links below at the bottom which you can refer to if you would like to know more! The main aim of this post is not to teach you quantum physics, but how to use Python programming to numerically solve complex physics problems.
So why Python you might ask? Well, because we will be using Python in this post to numerically solve the time independent Schrödinger equation shown above! Yes I know, snakes are usually not associated with quantum physics. Unfortunately there is no programming language I know of named “Cat” (there should be!), so today’s post will have to make do with the title “Schrödinger’s Python”.
When one talks about quantum mechanics, the animal one usually associates with the subject is the cat. Please check here for the latest version of this post. Update: due to various difficulties encountered in writing Python code and mathematical equations in WordPress, I have decided to start migrating most of my content to Github.