/*
* FIXED-POINT ALGORITHM 2.2
*
* To find a solution to p = g(p) given an
* initial approximation p0
* This code is to solve Illustration problem on Page 61 of the textbook. 
* Formulation (d) is used.
*
* INPUT:  initial approximation; tolerance TOL;
*         maximum number of iterations NO.
*
* OUTPUT: approximate solution p or 
*         a message that the method fails.
*/

#include <stdio.h>
#include <math.h> 
#include <iostream>
#include <iomanip>

#define true 1
#define false 0

static   double    G(double );             /* function g */

using namespace std;

int main()
{
      double TOL, P0, P;
      int I,   NO,  FLAG = 2, OK;
      char AA;
      FILE *OUP[1];

      // P0 = 1.5;              // init. guess 
      P0 = 2.5;              // init. guess 
      TOL = 1.0e-7;         // tolerance for |p-p0|
      NO =  200;             // max. number of iterations

      /* STEP 1 */ 
      I = 1;  
      /* STEP 2 */ 
      OK = true;
      while(I<=NO && OK) 
      {
         /* STEP 3 */ 
         /* compute P(I) */  
         P = G(P0);
         if (FLAG == 2) 
            cout<<setw(3)<<I <<"  " << setiosflags(ios::fixed) << setprecision(8) << P <<endl;
         /* STEP 4 */ 
         if (fabs((P-P0)/P) < TOL) { 
            /* procedure completed successfully */
            cout<<"\nApproximate solution P = " << setiosflags(ios::fixed) << setprecision(8) << P<<endl;
            cout<<"Number of iterations = "<< I <<endl; 
            cout<<"(    Tolerance = "<<setiosflags(ios::fixed) << setprecision(8) << TOL <<endl;
            cout<<"|p-p0| = "<<setiosflags(ios::fixed) << setprecision(8) << fabs(P-P0) <<endl; 
            OK = false;
         }
         else {
            /* STEP 5 */ 
            I++;
            /* STEP 6 */ 
            /* update P0 */ 
            P0 = P; 
         }
      }
      if (true == OK) {
         /* STEP 7 */  
         /* procedure completed unsuccessfully */ 
         cout<<"\nIteration number "<< NO <<endl;
         cout<<" gave approximation "<< setiosflags(ios::fixed) << setprecision(8) << P <<endl; 
         cout<<"not within tolerance "<<setiosflags(ios::fixed) << setprecision(8) << TOL <<endl; 
      }

   return 1;
}
/* Change function G for a new problem */
static double G(double X)
{
   double g; 

   // g = sqrt(10.0 / (4.0 + X)); 
   // g = X*X*X + 4*X*X + X - 10; 
   // g = sqrt(10/X - 4*X); 
   // g = 0.5*sqrt(10 - X*X*X); 
   g = 5.0/(X*X) + 2.0; 
   return g;
}