Computational Physics II
• Assignment 2 - Schrödinger's Equation Solutions via Conjugate Gradient

Due: February 18

Code

• I used the conjugate gradient solver provided by Professor Dawson to solve the 2D Schrodinger Equation {-∇² + V(x,y)} Ψ(x,y) = E Ψ(x,y) for a two diimensional cylindrical square well with a potential of radius 1 and energy -20. Here is a picture of the potential:

• For each iteration of the solver, the conjugate gradient marches along a steepest descent line while preserving the global gradient. In order to do this, the system must be initialized with an arbitrary search vector. For this assignment, the vector was set to the 0 vector. In addition, the matrix equation that results from the Schrodinger Equation must be initialized to some source vector b in order to solve a more general set of equations. As is normal with most matrix equations, the form of the Schrodinger equation becomes Ax=b, where A = (E-H) from above, x = Ψ, and b is a randomly initiatlized source vector.
• The system is discretized onto a 50x50 lattice, with a lattice spacing of 0.08. This ensures that any wave functions resulting on the lattice are well represented. The vector x is then solved for using the conjugate gradient algorithm.
• In practice, a range of E values is allowed to be input to the system. The program can either scan this in from command line parameters or prompt for them during execution. They three paramters are E_min, E_max, and ΔE. The program will sweep from E_min - E_max in steps of ΔE and solve the system for x at each iteration.
• The goal of this problem is to discover the eigenvalues of the matrix equation Ax=b, from the Schrodinger Equation. This will tell us the energy levels associated with this potential. As a small aside, since the potential is a well, it will have bound eigenstates.
• We can use something called the Rayleigh quotient to use the approximate solution x_approx to obtain an approximate eigenvalue λ. The reason this is useful is because there is a bound on the nearness of λ and its closest actual eigenvalue. A very obvious use of the Squeeze Theorem allows us to plot the lower bound and upper bound to find where they are effectively equal. This will tell us where the eigenvalues are.
• For the potential given at the beginning of this section, this approached was used to find all of the bound states (aka states which have negative energy). Here is a plot of the lower and upper bounds as a function of E:

• Notice the function periodically pinches at the points where we expect eigenvalues to occur. These pinches occur at roughly, -16.2, -10.5, -3.8, and -2. Also notice that at some pinches (-10.5 and -3.8) it actually appears there are multiple pinches contained in some small range around those points. We will next investigate this phenomenon.
• For the first three eigenvalues, the eigenstates will be found. In terms of the Schrodinger Equation, this means that the ground state and first two excited state energies will be found. The eigenstates that correspond to these eigenvalues will be the Schrodinger wave function at the particular energy.
• The ground state was relatively easy. The program was run many times, and at each interval I manually resized the search space until I zeroed in on the exact energy where the bounds pinched in. This value was -16.1908 and the wave functions that resulted are pictured below. The interesting note here is that it seems there are actually two wavefunctions present, but this is only because we can take any linear combination of the actual wave function and it will be a solution to the Schrodinger Equation because it is in fact a linear equation. The two wave functions you see here are just a biproduct of this property:

• The next two energy levels have an added intricacy. These are the levels that I pointed out earlier as having possibly two eigenvalues where it looks as though there should be only one. But this is not abnormal in the world of Qunatum Mechanics! Take the hyperfine splitting, and the Zeeman effect. These are instances where there are broken degeneracies. In this case there is not actual broken degeneracy. There are actual degenerate eigenvalues of the system, but what causes the splitting. Remember in writing the program that we were forced to discretize a lattice of some spacing a. Because of this discretization we cannot hope to represent an actual continuous potential. Therefore we inadvertantly break some of the degeneracies of the system because we have perturbed the actual system by discretization. This is the same principle behind the hyperfine splitting. There is some small correction term that is added to the potential of a hydrogen atom that causes its energy levels to split apart. The only difference is we didn't actually mean to do it.
• There actually is no splitting around the -10.5 eigenvalue, but there is still degeneracy resulting from the linearity effect we saw in the ground state. The difference is that there is less symmetry in the first excited state (i.e. not rotationally symmetric for all rotations), so there will be more distinct wave functions:

• The second excited state was labeled as having energy -3.8. When inspected closer, we can see the splitting that I was just referring to earlier:

• The splitting produces two eigenvalues -3.76077 and -3.46268. Although these are distinct eigenvalues, we we look at the eigenstates they correspond to, we see that these are actually the same wave function and the splitting is just a biproduct of the discretization of space:

Wave Function with Energy -3.46268	Wave Function with Energy -3.46268
Wave Function with Energy -3.76077	Wave Function with Energy -3.76077

• The only question remains as to why some degeneracies are not split on the lattice. This is again related to the splitting we normally see in QM. If the lattice spacing were decreased to 0, we would have a continous representation of the potential. In this case, there would be no splitting of the eigenvalues, because a degeneracy would be a degeneracy. That means that the amount of splitting that occurs must be proportional to a, and for the particular a chosen in this example those energy levels were not degenerate.

2. Rotationally Symmetric Linear Potential

• Using the same processes and program as in problem 1, we can calculate the eigenvalues (energies) of a new potential. For this problem the potential becomes V(x,y) = √(x²+y²). Here is a picture of what it looks like:

• We do the same process as before wherein we scan a wide range of energies to pick out the first 6 in the spectrum. Since this potential is positive everywhere, there will be no bound states and we need to consider negative energy values. Here is a scan from E=0-22:

• Notice there is splitting in this case too as you can see for the first and second excited states. The first 6 eigenvalues can be identified through the finer resolutions as 1.40, (3.89, 3.92), (6.37, 6.50), 11.30, (13.72, 13.77, 13.80, 13.85, 13.91, 13.97), and 21.18. There is a lot of splitting around the 13.7-14.0 region as you can see in this picture:

• Finally, here is a plot of the ground state wave function (and remember this is just one of the possible linear combinations that still solve the problem):