/* --------------------------------------------------------------------- 
// Figure 8.11:  Using Newton's method to find a root.
// Figure 8.12:  Sliding downhill.
// --------------------------------------------------------------------- 
*/

#include "tools.h"
void newton( void );          

#define EPSILON    .000001         /* Convergence tolerance.         */
#define TRIES       100            /* Give up after 100 iterations.  */

/* ------------------------------------------- function to be solved */    
double f (double x)  { return (x*x + 3*x - 3); }
/* ------------------------------------------------- derivative of f */
double fprime (double x)  {  return (2*x + 3); }
/* ----------------------------------------------------------------- */      

void main( void )                                                                                                              
{  
    int do_it_again;                /* repeat-or-stop switch */ 
 
    banner();                                                        
    puts("\n Finding root of function f(x) = x*x + 3*x - 3 by Newton's method.");
    
    do { newton();
         printf( " \n\n Enter '1' to continue or '0' to quit: " ); 
         scanf( "%i", &do_it_again );           
    } while (do_it_again != 0);
    bye();            
} 

/* -----------------------------------------------------------------       
// Find one root of function f using Newton's method.                
*/  
void  
newton( void )          
{  
    int k;              /* Iteration counter. */
    double x;           /* Current guess for root. */                  
    double delta_x;     /* Difference between successive guesses. */
    
    printf( " Enter x0, the initial guess: " );   
    scanf( "%lf", &x );                                        

    for (k = 0; k < TRIES; ++k) {    /* Give up after TRIES iterations. */
        delta_x = -f(x) / fprime(x);
        if (!finite( delta_x )) break;            /* Overflow or underflow. */
        x += delta_x; 
        if (fabs( delta_x ) < EPSILON) break; /* Process has converged. */
    }

    /* Analyze reason for leaving loop and print results. ------------- */
    if (k == TRIES)            
        printf( "\t Process did not converge in %i iterations.\n", TRIES );

    else if (!finite( delta_x )) 
        printf( "\t Overflow on delta_x.\n" );

    else 
        printf( "\t Process converged.\n" );

    printf( "\t x = %g,  delta_x = %g. \n", x,  delta_x );                  
    printf( "\t f(x) = %g,  f'(x) = %g. \n", f(x), fprime(x) );                  
}