/* A demonstration program for the class "Scientific Computing" taught
   at the Courant Institute of Mathematical Sciences of NYU.  It is 
   supposed to demonstrate ideas about good programming practice for
   numerical computing.  If you want to debate what's here, send email
   to Jonathan Goodman at goodman@cims.nyu.edu.  
   
   This program was downloaded from the web page
   <a href="http://math.nyu.edu/faculty/goodman/teaching/Sci_Comp_99/"
            http://math.nyu.edu/faculty/goodman/teaching/Sci_Comp_99/</a>
            
------------------------------------------------------------------------*/  

#include <iostream.h>

//  A procedure to solve a recurrence relation, defined more fully in the
//  header file.

double RecurrenceRelationComputation( double a, int n, int*   errorFlag) {

   *errorFlag = 0;  // An optimistic default.
   if ( n < 0 ) {
      cout << "In routine RecurrenceRelationComputation, stopping because"
           << " n is negative.  Got n = " << n << endl;
      *errorFlag = 1;
      return 0.0; }
   
   if ( n == 0 ) return 1.0;  // Some dummy might do this, and it isn't 
   if ( n == 1 ) return 1.0;  // exactly wrong ... so be prepared.
   
   double x_k, x_km1, x_km2;  // x_k means "x sub k".  km1 means "k minus 1".
   x_km1 = x_km2 = 1.0;
   for ( int k = 2; k <= n; k++ ) {
      x_k   = a*x_km1 + x_km2;  // There is a more efficient way to program 
      x_km2 = x_km1;            // this.  In the future, compilers will figure
      x_km1 = x_k;              // it out.
    }
    
   return x_k;
 }