about
Jim: the dft/ and diffeq/ directories under this folder are from jesse's home directory, work done in the first month of the term. The mathematica files were done in the last month.
Jesse
Below you will see many files which have been uploaded. I wasnt sure how else to get them off of
my harddrive and onto the cs account. Anyway, youve seen the harmonic stuff. Its calculating
the harmonic oscillator states from the schrodinger equation. After i finished that, i worked on
calculating perturbed eigenstates of a perturbed SHO potential. After many tries to get it working i have worked out the solutions to several perturbed hamiltonians. Then there is also densityplots, which works out the hydrogen atom wavefunctions and plots the density in 2d. The actual commands may not work for you, as i downloaded an opengl 3d plot viewer module for my computer, and the commands are specific to that module. You can still see the 2d slices though. So thats what ive been doing.
In a little more detail:
Numerical Calculation of Harmonic Oscillator Eigenstates:
I started off with the file harmonic.nb. In this file, i started with the schrodinger equation for the simple harmonic oscillator. I solve this equation numerically using mathematicas built in ndsolve function over the range from -L to 0. using the boundary conditions psi(-L) = 0, and psi'(-L) = 1. This corresponds to an exponentially decaying solution to the left superimposed with one exponentially increasing to the left. The problem is that we only want the decreasing function as a solution because we need to be able to normalize psi, and an exponentially increasing function cant be normalized. If we take L large enough that V(-L)>>E, there will be negligible error because neither solution blows up as we integrate through our range of -L to 0. What we get out of ndsolve is an interpolating function which is actually a nested list. We use this as a replacement for psi. We also need to get an equation for psi' as we need that later. To find the even pairity solutions we need to solve for what energy makes the derivative of the wavefunction equal to zero at x=0. So by using the findroot function, and giving it a starting guess for the energy it finds the eigenvalue. Similarly for odd solutions, we need to find the energy that makes the value of the wavefunction equal to zero at x=0. Again we use the findroot function and give it a guess for the energy value. Then we need the eigenfunction which corresponds to the eigenvalue we found. So we again use the replacement rules, but instead of a variable energy value, we give it the eigenvalue we just found. This then gives half of the eigenfunction, for -L to 0. Now depending on whether we have an odd or even function, we define the second half differently. For even, we say that psi(x) = eigenfunction(x) for x<0 and psi(x) = eigenfunction(-x) for x > 0. And for odd functions we say that psi(x) = eigenfunction(x) for x<0 and psi(x) = -eigenfunction(-x) for x > 0. We then have an eigenfunction and an eigenvalue associated with it, the only problem is that its not normalized. So we take the integral over all space of the square of the eigenfunction. We take this value and divide the eigenfunction by the its squareroot, this is then the normalized eigenfunction. So this file you have to specifically go through and manually enter the value of the energy range your looking for, and the length of integration to calculate. So i went through 4 cases, the ground state and the first three excited states and plotted the normalized wavefunctions.
In my second iteration of this process harmonic2.nb, i attempted to automate the whole thing, make one big function which called the variable to plot, as well as the number of energy states you wanted plotted. I used if statements to determine if it was an even or odd state, and for some reason it never would plot out the values, it claims that my functions are non numerical, though its the same code as before, just set into a module. Everything below the main module section is the same as the last three sections of the first file, it just got left there.
In my third file harmonic3.nb, i tried to get a similar function to the previous file, however i made it so it just tried to calculate the ground state energy and eigenfunction, however this still was giving me the same error as before, nonnumerical. I also tried adding in more and less things to the module local definitions section, figuring that might have something to do with it. The first block on this file is just junk copied out to make sure i didnt just lose the code, sloppy i know, im sorry.
The fourth file, which is actually called harmonic5.nb was where i really started to make good progress. I set up the whole system as a table instead of as a function. This way, all you have to do is define the what the last step is in the table creation to decide how many solutions you want. The other thing i did differently was to just calculate one pairity of solution, in this case the even solutions. So the output is a table of the eigenstates plotted together for however many you want. The first piece on the top is something i got from the differential geometry book ive been using. It essentially blocks all graphical output unless you put S@ before the plot command. This way you can store graphics into memory without them being plotted on the screen, then plot them together later in a more cohesive way, though i dont really utilize it to its fullest extent here.
In the next file, harmonic6.nb, i tried to make a table of graphical output, by making two loops inside, one for the even and one for the odd. For some reason i started getting exactly the same errors again as when i was trying to calculate the two types of states together in previous notebooks. Even though there were two entirely seperate sets of variables, neither of them interacting, something was going on and i still am not sure why it didnt work.
In harmonic7.nb, i tried a whole new approach. In this one i made pretty much each line of the previous code into its own function of the value of the energy state. Each function calls upon the previous lines and functions, so when you get to the end the last function in essence calls all the other functions by use of the delayed equality statement. And im back to using if statements to gauge when the energy state is even or odd, but when you get to calling the final function, you still seem to get the same exact errors as ive been getting. Non numerical integrand when calculating the normalization constant.
In harmonic8.nb, i started to think about what each function was a function of, and started using dummy variables in the definitions of the functions, only to use a replacement at the end, but this still didnt seem to make much of a difference because it still gives me the non numerical integrand error when calculating the normalization constant.
In harmonic9.nb i have to admit right off that i made a stupid error that i never took out. The normalization constant is just stupid in this one, i dont know what i was thinking when i wrote that. Other than that, i made an adjustment to the evalue function which changes the value of L for different ranges of energy states, the higher the state, the larger L needs to be to get an accurate result. I also realized in this notebook that i could plot the non normalized psi functions and see their shape anyway. I wasnt sure why then, but my plots were wrong for all odd functions, it was because my normalization function was just stupid.
harmonicfinal.nb is essentially the same thing as harmonic9.nb, just cleaned up a little bit.
The real final work here was harmonicfinalb.nb. This one utilized all the methods of harmonic9.nb, while getting the normalization function correct. This time i was able to plot out any amount of solutions i wanted, as well as their squares. I did this, and kinda wrote in a little bit about the solution, though ive done i much better job here, as well as making the plots look nice. So this is it, solving for the SHO eigenvalues and eigenfunctions numerically with the goal of plotting them with their respective probabilities.
Calculating the eigenfunctions and eigenvalues of a pertubation to the SHO potential by use of the eigenstates of the SHO hardcoded in:
I started with notharmonic.nb. The basic stratagy was to code in the eigenstates of the SHO as a function of energy state. i then define what the hamiltonian used is. In this case i used 0.1Exp[-u^2/2] + Hsho. Then i wrote out two functions which calculated the dot product of two eigenfunctions around each part of the total hamiltonian. So theres one which calculates Integral[ym[u]*Hsho[u]*yn[u],{u}] and there one that calculates Integral[ym[u]*Hp[u]*yn[u],{u}]
Where H = Hsho + Hp. Then, by making a matrix of dotproducts for 0
In the second file, notharmonic2.nb, i combine the two dotproduct functions into one, so it calculates just one matrix to begin with. I do find the eigenvalues and plot them on a diagonal matrix, and i also find the first 15 eigenvectors, which i then use to plot the first 15 perturbed eigenstates of this new system.
In notharmonic3.nb, i added in a function which calculated the set of eigenvalues for a given number of used eigenstates to calculate, and one which calculated the set of eigenvectors as well. I then used the function eigenvectors as the coefficients when plotting the perturbed eigenstates.
notharmonic4.nb is essentially the same thing as notharmonic3, however in this savefile i have run the eigenvalue and eigenvector calculations for a 15x15 matrix. Then i plot out the first 15 perturbed eigenstates on top of the original SHO states. Below each of those graphs, i plot out the difference between them. This gives those 15 eigenstates at a pretty high accuracy, of the order of magnitude of 1.
notharmonic5.nb is where i made the most progress. In this i made several calculations like in notharmonic4, however i made calculations with different pertubations. starting with no pertubation, which is just trivial since the solutions are the eigenstates we want to write this out in terms of. Then 0.5*Exp[-u^2/.1], then 1*Exp[-u^2/.1], then 8*Exp[-u^2/.1], then 8*Exp[-u^2/.1]-5*Exp[-u^2/.01]. The first several are simple additions to the potential, a slight bump in the center of a quadratic well, first with small size, then with added strength. The last pertubation is the same bump as the previous run, but there is another small well built into the bump thats in the well. A bit confusing, i just wanted to see what would happen if we changed it that much.
notharmonic6.nb is the same data, just graphed in a more legible and aesthetic way. I added in legends as well as plots of the potential that is being invstigated in each case. The difference graph was also added in on top of psi and psi^2.
notharmonic7.nb is the final draft. Everything calculated and plotted out as in the previous two files, but i also added in information and documentation about what i did for each step into the file. That one should follow through pretty self explainitory.
