PHYS 580 COMPUTATIONAL PHYSICS
Sample programs for numerical algorithms
Clink on the links below; you can save these files to disk and compile them.
Note: very important! To view properly, you must view as page source (open the file, and then ctrl+u) You can also right-click on the link and choose "Save Target As..."
Exercise: use the above programs (or write your own) and using xmgrace or some other graphing program and try to make graphs similar to those above. Try to write as either a postscript file or as a jpeg file (hint: you have to use “print” you cannot simply “save”).
Quadrature (integration of functions)
(to compile, do gfortran –o qtest.x quadmaster.f90 quadraturesubs.f90 or use your favorite compiler)
Exercise: Modify the above program but compare for a different function (you choose). Make graphs.
Exercise: Write a subroutine for Bode’s rule and add to above program. Make graphs.
Alternate version using standard fortran function random_number
For C/C++ use rand() in <stdlib.h>
Exercise: modify the example calling routine above but instead of computing the average and std deviation of a uniform distribution (ran3 or random_number) use the Gaussian distribution gaussvar . Do you get the average and (importantly) the std dev you expect? For experts: Compute the 4th moment, both analytically and numerically.
algorithm: here are 2 codes that compute
via the Metropolis
algorithm (see the recommended texts for details, also these notes
on the Metropolis algorithm).
Both codes use ranlib.f and metrolib.f90.
The subroutine metropolis (in metrolib.f90) writes the functions FUNC0 and FUNC1 to files, so that one can look at the correlation function C(k) = < f(xi)f(xi+k) > - <f>2
< x**2 exp(-x**2) > = < FUNC1 >
< exp(-x**2) > = < FUNC0>
For version 0, FUNC0 =1 and correlation function is trivial. Correlation function for FUNC1.
For version 1, correlation functions for both FUNC1 and FUNC0
(The program I used to compute the correlation function is here, but it is not well documented)
Another simple Monte Carlo example -- computing
volume of n-dimensional sphere
Matrix multiplication (for benchmarking speed and optimization) uses Fortran inbuilt function random_number()
Notes on Broyden's method by Dan Mickelsen; for an example application see Phys Rev. C 78, 014318 (2008).
For a derivative-free minimizer POUNDerS see section III in Phys. Rev C 82, 024313 (2010).
Ordinary differential equations
Solving a simple, first-order DE:
More 4th-order Runge-Kutta: a central force problem (must link to rungekuttalib.f) which calculates orbit (r vs. theta, or x vs y). Some results: (1) "gravitational" force (1/r2 force) with different angular momenta; (2) gravitational force with different step sizes in theta; (3) different central forces (1/rn) notice the perihelion shift.
Partial differential equations
Sample code for Neumann boundary condition du/dx = 0 at x = 0, L
Relaxation for simple inhomogenous differential equations:
Dirichlet boundary conditions with simple relaxation sweep
Combined Dirichlet and Neumann boundary conditions (with improved relaxation sweep) and output graph
Last modified Nov 15, 2017