780.20: 1094 Session 5

Handouts: New: eigen_tridiagonal.cpp printout, eigen_basis.cpp printout, harmonic_oscillator.cpp printout, nan_test.cpp printout, and square_well.nb Mathematica notebook.
From session 4: eigen_test.cpp printout and GSL Eigensystems documentation.

Your goals for this session:

• Finish some leftover Session 4 tasks.
• See some examples of nan's and inf's.
• Use the eigen_tridiagonal program to look at eigenvalues of the harmonic oscillator potential.
• Find the lowest bound-state eigenvalues of two familiar potentials using the eigen_basis program.
• Examine (and try to understand) how the accuracy of your results depends on the size of the harmonic oscillator basis and the choice of the basis parameter b.

Please work in pairs (more or less). The instructors will bounce around 1094 and answer questions.

Leftover Tasks from Session 4

Work on these at most for the first 30 minutes of the period.

1. Pointer Games. The important part of this exercise is not understanding pointers in detail but to successfully modify a code based on an example. For now, just look at the derivative_test_new.cpp code, which is a solution to this exercise, and see if you have questions.
2. Linear Algebra with GSL Routines. You'll see the same things in Session 5. If you are running short of time, just go through the discussion in the Session 5 notes.

Nan's and Inf's

Just a quickie: Take a look at nan_test.cpp, use Dev-C++ to create nan_test and then run it. [Note: Dev-C++ uses .#INF for inf and .#IND for nan.]

• What do you think the result of 0.*(1./0.) will be? Predict and then modify the code to check it out. [Note: Use the variables "numerator" and "denominator" to calculate (1./0.).]
  
  
  
  

Bound States by Matrix Diagonalization in Coordinate Representation

The program eigen_tridiagonal.cpp uses the GSL library routines explored in eigen_test.cpp (session 4) to find the eigenvalues and lowest eigenvector of the harmonic oscillator using a method described in the Session 5 notes.

1. Create a project called eigen_tridiagonal, add eigen_tridiagonal.cpp, and compile and run a few times to see how it works. Try various values of Rmax and N. How can you check that the code is working?
   
   
   
   
   
2. Change the code so that only the lowest few eigenvalues are printed out. Look at the output file eigen_tridiagonal.dat and plot it with gnuplot. What is this function? If you want it to start at r=0, what is the appropriate value of the function there?
   
   
   
3. You need to pick a reasonable value of Rmax. Justify your choice:
   
   
   
   
4. For your choice of Rmax, try N = 4,8,16,32,64,128,256,512,1024. How does the relative error for the lowest eigenvalue scale with N? Attach an appropriate plot to validate your answer.
   
   
   
   

Bound States from Diagonalizing the Hamiltonian in a Basis

The program in eigen_basis.cpp uses the GSL library routines to diagonalize a Hamiltonian matrix in a basis of harmonic oscillator wave functions. You may want to refer to the GSL handout on eigensystems (there is also a printout of eigen_basis.cpp).

The eigen_basis program uses units in which the particle mass is 1 and hbar=1. The program asks you to choose

• a potential (Coulomb or square well);
• the parameter b for the harmonic oscillator basis (see harmonic_oscillator.cpp for the definition);
• the number of basis states to use.

The Coulomb potential is defined with Ze²=1, which means that the Bohr radius is also unity. This means that the bound energy levels are given by E_n = -1/2n², with n=1,2,...
The square well potential is defined with radius R=1 and depth V₀ = 50. You should find that there are three bound states.

Here are some subgoals. There won't be time for everything (as usual) but you'll have a chance to finish them as part of a future assignment.

1. Run the Mathematica notebook square_well.nb (make sure you understand what it is doing; e.g, look up FindRoot in the Help Browser). Find the bound-state energies for the square well parameters used here (you need to change the notebook parameters!).
   
   
   
   
2. Create a project called eigen_basis and compile and link the code (you'll need eigen_basis.cpp and also harmonic_oscillator.cpp). Run it a few times with each of the potentials to get familiar with it. If you try too large a basis size, the run time may be too long (so start small!). Look through the printout to see the basic idea of how the code works and find where the equation for the matrix element is implemented.
3. Based on the "exact" results from Mathematica, which predicted eigenvalue(s) do you think are most reliable? Which eigenvalues are calculated most effectively, those of the Coulomb potential or the square well potential? Why do you think this is?
   
   
   
   
   
   
4. You have under your control the size of the basis (i.e., the dimension of the matrix) and the harmonic oscillator parameter b (see harmonic_oscillator.cpp for the definition). For a fixed basis size (pick one that reproduces the ground state reasonably), how do you find the optimum b? (Hint: think gnuplot!) Can you qualitatively (or semi-quantitatively) account for your result? (Think about the potentials and guess what the lowest wave functions will look like and what changes about the basis when the harmonic oscillator parameter b is changed.)
   
   
   
   
   
   
5. If you now fix b (if you have time you can consider two or three different values in turn), how can you find how the accuracy of the ground state energy scales with the basis size? Make an appropriate plot.
   
   
   
   
   
   
6. Look at the code. How could you make it more efficient? (What do you think is the limiting factor based on the scaling of the time with the size of the basis?) For example, could you speed it up by almost a factor of two? (Hint, hint!)
   
   
   
   
   
   
7. How would you find the wave function that corresponds to a given state (e.g., the ground state)? Add code to generate the lowest wave function for the lowest bound state (hint: it involves the eigenvector).