#############################################################
# Building 'bvp1d' and sample output on the SGIs
#############################################################

einstein% pwd; ls
/usr2/people/phy329/linsys/ex2
Makefile   bvp1d.f

einstein% make
  f77 -g -c bvp1d.f
  f77 -g -L/usr/local/lib bvp1d.o \
         -lp329f -llapack -lblas -o bvp1d

einstein% bvp1d
 usage: bvp1d <level> [<option>]

        Specify option .ne. 0 for output
        of error instead of solution

einstein% bvp1d 4
   0.0000000000000000E+00  3.2196467714129540E-15
   6.2500000000000000E-02  0.3876394685723126    
   0.1250000000000000      0.7162643420150208    
   0.1875000000000000      0.9358444623383715    
   0.2500000000000000       1.012950746721882    
   0.3125000000000000      0.9358444623383710    
   0.3750000000000000      0.7162643420150199    
   0.4375000000000000      0.3876394685723112    
   0.5000000000000000      1.4432899320127035E-15
   0.5625000000000000     -0.3876394685723084    
   0.6250000000000000     -0.7162643420150170    
   0.6875000000000000     -0.9358444623383682    
   0.7500000000000000      -1.012950746721879    
   0.8125000000000000     -0.9358444623383686    
   0.8750000000000000     -0.7162643420150178    
   0.9375000000000000     -0.3876394685723095    
    1.000000000000000      0.0000000000000000E+00
 rmserr =   8.8841389573660889E-03

#############################################################
# Simple convergence test: solve BVP on a sequence of 
# levels (h, h/2, h/4, h/16, etc.), redirect stdout to 
# /dev/null so that only the overall RMS error appears on 
# terminal.  Note how RMS error goes down by very nearly 
# a factor of 4 at each successive level, indicating 
# O(h^2) convergence.
#############################################################
einstein% foreach level (4 5 6 7 8 9 10)
foreach? bvp1d $level > /dev/null
foreach? end
 rmserr =   8.8841389573660889E-03
 rmserr =   2.2413991373375825E-03
 rmserr =   5.6382739826356723E-04
 rmserr =   1.4145099550242024E-04
 rmserr =   3.5428279660453927E-05
 rmserr =   8.8654982385453620E-06
 rmserr =   2.2174426933897138E-06

#############################################################
# Making output files for subsequent plotting via gnuplot.
# See Class Notes for postscript.
#############################################################
einstein% bvp1d 8 > out8
 rmserr =   3.5428279660453927E-05

einstein% bvp1d 5 1 > err5
 rmserr =   2.2413991373375825E-03

einstein% bvp1d 6 1 > err6
 rmserr =   5.6382739826356723E-04

einstein% bvp1d 7 1 > err7
 rmserr =   1.4145099550242024E-04

#############################################################
# Gnuplot 'script' (gpsoln8) for making plot of level-8 
# solution
#############################################################
einstein% cat gpsoln8

set terminal postscript portrait
set size 0.760,1.0
set output "soln8.ps"
plot [0:1] [-1:1] "out8"
quit

#############################################################
# Make the plot
#############################################################
einstein% gnuplot < gpsoln8

#############################################################
# Gnuplot 'script' (gperr) for making plot of error from 
# level 5, 6 and 7 calculations
#############################################################
einstein% cat gperr

set terminal postscript portrait
set size 0.760,1.0
set output "err456.ps"
plot "err5", "err6", "err7"
quit

#############################################################
# Make the plot
#############################################################
einstein% gnuplot < gperr