#############################################################
# Building 'integral' and sample output on the lnx machines
#############################################################
lnx1% pwd; ls
/home/phys410/ode/integral
Makefile  fcn.f  integral.f

lnx1% make
pgf77 -g -c integral.f
pgf77 -g -c fcn.f
pgf77 -g -L/usr/local/PGI/lib integral.o \
               fcn.o -lp410f -lodepack \
               -llinpack -lblas -o integral

#############################################################
# Usage
#############################################################
lnx1% integral
 usage: integral <xs> <xf> [<tol>]

#############################################################
# We can check the results using the following Maple 
# code (or similar)
# 
# > Digits := 25;
# > evalf(int(exp(-x^2),x=0.0..5.0));
#
#       .8862 2692 5451 3954 7538 24606
#############################################################
lnx1% integral 0.0 5.0
   0.8862 2692 5451 8431
                    ^

#############################################################
# > evalf(int(exp(-x^2),x=0.0..100.0));
#
#       .8862 2692 5452 7580 1364 90837
#############################################################
lnx1% integral 0.0 100
   0.8862 2692 5447 2388  
                 ^

#############################################################
# Repeat previous computation with less stringent tolerance,
# note that answer is (roughly) correspondingly less 
# accurate.
#############################################################
lnx1% integral 0.0 100.0 1.0d-6
   0.8862 2897 2928 0249
           ^