/* 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>
            
------------------------------------------------------------------------  

A program to solve a simple two term recurrence relation for a sequence of
parameter values and write the output to a file that can be read by a 
plotting program.

                              */          
            

#include <math.h>  /* Declarations of the routines in the math library.
                      I always put this in, even though some numerical 
                      programs don't need it; it's a habbit.    */

#include <iostream.h> /* The C++ io library.  This makes "cin" and "cout"
                         possible.  */
                         
#include <fstream.h>  /* The C++ routines that communicate with files.  
                         These make "ifstream" and "ofstream" possible. */
                         
#include <iomanip.h>  /* So that we can format the output.              */
                      
#include "recurrence.h"  /* This has the declaration of the procedures,
                            so I don't have to clutter up the main program
                            with them.            */
                            
#define N_MAX 30   /* The farthest in the recurrence relation ever to go. 
                      The first thing the compiler does is substitute "30"
                      for every "N_MAX" in the code below. */
                      
#define N_VALS 20  /* The number of points to plot.  */
                      
#define   A0      1   /* Often add white space to make things line up and */
#define DELTA_A  .01  /* be easier to read and check for correctness.     */ 

#define OUTPUT_FILE "recurrence_run_output"    
                          
int main() {

   double x0, x1; x0 = 0; x1 = 0; // It is possible to put more than one
                                  // statement on a line.  Some people think
                                  // it's bad programming practice, however.
                                  
   double a,       // It's also possible to put one statement on more than
          b,       // one line, particularly if that makes the program more
          c;       // readable.    
                           
   double xn[N_VALS];  // An array of N_MAX doubles, for holding the n-th
                       // term in the series for the i-th parameter value.
                       // Making N_MAX all capital letters indicates that it
                       // is not a variable, but is defined by a #define above.
                      
                      
   int errorFlag;      // Another style for multi-word variable names.
   for ( int i = 0; i < N_VALS; i++ ) {
   
      a = A0 + i*DELTA_A ;  // 
      xn[i] = RecurrenceRelationComputation( a, N_MAX, &errorFlag);  // Is this name too long?
      
         if (errorFlag) {  // Always test the error flag, even if you are sure
                           // nothing can go wrong.  In Java, this test is less
                           // clumsy.
            cout << "In main program, " // Error messages must always say where 
                                        // they were generated ... 
                 << "got error flag " << errorFlag   // ... and what it was.
                 << " with i = " << i << ", and a = " << a 
                 << ", and return value = " << xn[i] << endl;
            cout << " Stopping here and returning." << endl;
           return 1; }
     }  // End of "for ( int i = 0; i < N_VALS; i++ )" loop.
     
   
//  Write the data to a file.  This uses the C++ file io mechanism.  It could
//  also have been done in plain C.

   ofstream outFile;  
   outFile.open(OUTPUT_FILE, ios::out ); // Overwrite previous data, if any.
      if ( ! outFile ) { //  Check that the file opened properly.
         cout << "In main program, outfile.open could not open file "
              << OUTPUT_FILE << ".  Stopping here and returning." << endl;
         return 1;
       }
       
/*  Set up the format so that floating point numbers are printed in a
    uniform way.  */ 
      
   outFile << setprecision( 7 );  
           // Print six digits after the decimal point.
           
   outFile.setf( ios::scientific, ios::floatfield );  
           //  Scientific notation with all zeros printed.
           
   for ( i = 0; i < N_VALS; i++ ) {
   
           /*  These if statements are there to make the output line up.
               I wish I knew a better way to do this in C++.  If someone 
               knows a way that does not use a C io routine, please let me
               know.          */
               
      if ( i < 10    ) 
          outFile << "    " << i << "  " << xn[i] << endl;
      else if 
         ( i < 100   ) outFile << "   "  << i << "  " << xn[i] << endl;
      else if  
         ( i < 1000  ) outFile << "  "   << i << "  " << xn[i] << endl;
      else if  
         ( i < 10000 ) outFile << "  "   << i << "  " << xn[i] << endl;
    }
     
     
    return 0;
    
 }