If you haven't done so already, modify your program from HW2 to use
2nd order accurate surface flux terms (fit a quadratic curve to the
surface temperature and the two volume centers in from the
surface). Note that I'm requiring a second order
evaluation of the temperature derivative at the surface in your flux
term. I am not requiring that you get a conduction equation in
the surface volumes that is second order accurate. If your HW2
implementation was first order,
provide me with a revised set of equations for the 9 volume
types.

Run the HW2 problem with 3x3, 9x9, and 27x27
volumes. Apply Richardson Extrapolation analysis to temperatures
at
the
center point in the rod, and points located midway along one x-facing
surface and midway along one y-facing surface. Note that the last
two evaluations are surface not volume centered temperatures and you
will need to use an equation consistent with your second order surface
flux evaluation to get them. At each of these three points in
space give me the predicted
order of accuracy, the predicted error on the finest (27x27) mesh, and
predicted error on the 9x9 mesh.

You will probably need more decimal digits of precision for your
temperatures than currently available in the file
steady.out. Just after the line:

21 FORMAT(10F8.2)

in subroutine edit, add the following lines:

WRITE(11,*)In addition you will need to change the formulas used to generate txs and tys, to get values that are second order accurate.

WRITE(11,*) ' Middle of y facing edge T = ', tys(nx/2+1,1)

WRITE(11,*) ' Middle of x facing edge T = ', txs(ny/2+1,1)

WRITE(11,*) ' Center Point in the rod T = ', t(nx/2+1,ny/2+1)