780.20: 1094 Session 3
Handouts: Background notes for Session 3;
Landau/Paez (LP) chapter 4 excerpts;
Gnuplot fitting example;
C++ formatting; printouts of codes
Your goals for today and next time (in order of priority):
 Finish leftover tasks from Session 2
 Look at a quick demo about comparing floatingpoint numbers
 Use a makefile to compile a project with multiple .cpp files
and a header (.h) file
 Duplicate the results in LP Fig. 4.3 with a loglog plot,
then fit the slopes with Gnuplot
 Practice coding an algorithm:
Add a new subroutine with the 3/8 rule and analyze
 Add an appropriate integration routine from the GSL and analyze
the error
Please work in pairs.
The instructors will bounce around 1094 to ask and answer questions.
You should make or change to a 780 subdirectory,
download session03.zip from the 780 homepage and unzip it.
Leftover Tasks from Session 2
Here are the priority tasks to finish from Session 2. Don't worry about
other items; when you complete these, move on to Session 3.
 Summing in Different Orders. Make sure
to finish the
entire section. Be careful in part 5: There is one bug (a typo) that
cannot be caught by the compiler.
 Bessel 13 (you'll finish this in PS#1 if
necessary).
Comparing Floating Point Numbers
A common task in a computational problem is to check whether two
floating point numbers are the same. In C++, the comparison operator
is == (not just one =),
so it might seem that a simple if statement
would do the job. Here's an example of why this fails in general:
 Look at the file number_comparison.cpp in an editor. What does
it do?
 Compile, link, and run the
program (there is a makefile).
Why doesn't it give the answer you want? Add statements
to print out x1 and x2 to check your response.
What did you add?
 Suggest a better way to compare two numbers based
on the idea that the relative inaccuracy of any number can
be as large as the machine precision.
Numerical Derivatives: Pass 1
The Session 3 notes have a short introduction to numerical
differentiation (see the HjorthJensen notes online for more detail).
These are among the simplest algorithms for us to derive and to verify
the theoretical approximation and roundoff errors (the errors for
most realworld algorithms are noisy).
 Look at the file derivative_test_simple.cpp
while reviewing the discussion in the notes. (Also look at the
last section in the notes, which describes pointers to functions,
which are used here.)
What part(s) of the code do you not understand? What is being printed?
 Use the makefile to compile, link and run the code, generating the file
derivative_test_simple.dat.
There is a missing #include statement you need
to add (compare with previous codes). Add a statement to print
a header line to the output file with a label for each column
(and start the line with # so that it doesn't screw up gnuplot).
Show an instructor your output file.
 Make a graph with gnuplot
with two plots:
the logarithm of the relative error for forward difference vs. the
logarithm of h (which is Deltax) and the analogous plot for
central difference.
Print the graph and attach.
 Are the slopes in each region consistent with
the analysis in the text? Which is the better algorithm? Explain.
What value of h is optimal for each algorithm?
Makefiles for multiple project files (including header file)
Many of the example programs started as
"allinone" C programs from the Landau/Paez
text. One of these was integ.c, which we converted to C++ and
then split up into:
 integ_test.cpp, which has the main program and the function
to be integrated
 integ_routines.cpp, which has the integration functions themselves
 integ_routines.h, which has function prototypes
(a prototype tells the compiler the function return type and the
type of all its arguments)
There is also a function in gauss.cpp and there is make_integ_test to compile
it all.
In a subsequent step, we modified the codes so that the function is
passed as an argument to the integration functions.
The idea is that the integration routines should be isolated in a file
by themselves. The main program just invokes these routines.
The header file conveys the prototype information about the integration
routines to the main program and to any other functions
that might call the routines. (Note: Later we'll consider going
further and defining an integration class.)
 Look at the original file (integ.c) and see how it was converted
and split up
into integ_test.cpp, integ_routines.cpp, and integ_routines.h
(the cpp files are printed out for you).
Also look in the ORIGINALS subdirectory
to see the change to passing the function for
the integrand. Why do these changes improve the code?
 Create the executable
integ_test using the makefile make_integ_test and run it to generate
the integ.dat output file.

Change integ_test.cpp to output relative (rather
than absolute) errors and add column headings to the output file.
Note that only the files that have changed
since the last compilation are recompiled each time!
 Use Gnuplot to reproduce Figure 4.3 on page 59 of the LandauPaez
handout. Briefly explain what you
can learn from the plot (remember slopes!).
Consider all regions of the graphs.
 Bonus: Change the loop in integ_test.cpp so that the points on
the loglog plot are evenly spaced. What did you change?
Finding the Approximation Error From a LogLog Plot
From the plot in the previous section, we can estimate the approximation
errors by eye. Now we want to actually fit lines to find the slope
using Gnuplot. Use the handout as a guide.
 Modify the code so that it outputs the logarithm base 10 (log10)
of N and the relative errors.
 Find the slopes of the trapezoid and Simpson's rule plots in
the regions where they are linear.
 Are the results consistent with
the analysis in the text? Can you fit the roundoff error region?
Coding an Algorithm
This is primarily practice converting an algorithm to working code.
Your goal is to add a new function to integ_routines.cpp, called three_eighth,
which implements the 3/8 rule from LP Table 4.1.
 From the rule and the discussion in the text, write some
pseudocode to explain how you would integrate a function with this
rule.
 Implement your pseudocode in C++ by adding a new routine (called
three_eighth) to integ_routines.cpp (be sure to add a prototype to
integ_routines.h).
 Modify integ_test.cpp to output to a new file the results from
the new routine, and add them to the previous plot. Explain
the approximation error. [Warning: If you don't change the way that
integ_test.cpp loops through the number of intervals, you will likely
run into a subtle error that prevents you from getting the
approximation predicted theoretically!]
GSL Scientific Library Yet Again
 Go to the web page with the GSL manual (it is linked from the
780.20 web page). Find an appropriate integration routine for the
test integral we've been working on. Which one did you choose,
and why?
 Add another calculation to integ_test.cpp (with output to a file) using
this GSL routine.
 Compare the accuracy of the GSL routine to the others on an error plot.
780.20: 1094 Session 3.
Last modified: 07:47 am, January 12, 2009.
furnstahl.1@osu.edu