Ebook Chapter 5 Problem Solving and Program Design in C

• 5.9 Figure 5.19 Finding a Function Root Using the Bisection Method

#include <cstdlib>
#include <iostream>
#include <math.h>

#define FALSE 0
#define TRUE 1
#define NO_ROOT -99999.0

double bisect(double x_left,double x_right, double epsilon, double f( double farg));
double g(double x);
double h(double x);

using namespace std;

int main(void)
{
    double x_left,x_right,epsilon,root;
   
    cout<<"Enter interval endpoints> ";
    cin>>x_left>>x_right;
    cout<<"Enter tolerance> ";
    cin>>epsilon;
   
    cout<<"Function g"<<endl;
    root = bisect(x_left,x_right,epsilon,g);
    if(root != NO_ROOT)
            cout<<"g       ="<<root<<g(root)<<endl<<endl;
            cout<<"Function h"<<endl;
            root= bisect(x_left,x_right,epsilon,h);
            if(root != NO_ROOT)
                    cout<<"h       =n"<<root<<h(root);
                   
    return(0);
}

double bisect(double x_left, double x_right, double epsilon, double f(double farg))
{
       double x_mid, f_left, f_mid, f_right;
      
       int root_found;
       f_left = f(x_left);
       f_right = f(x_right);
      
       if(f_left * f_right>0){
                 cout<<"May be no root in       "<<x_left<<x_right<<endl;
                 return NO_ROOT;
                 }
                
       root_found = FALSE;
       while(fabs(x_right - x_left) > epsilon && !root_found)
       {
           x_mid = (x_left + x_right)/2.0;
           f_mid = f(x_mid);
          
       if(f_mid == 0.0){
                root_found = TRUE;
                }
       else if (f_left * f_mid < 0.0) {
            x_right = x_mid;
            }
       else {
            x_left + x_mid;
            }
           
       if(root_found)
                     cout<<"Root found at x =      ,midpoint of       "<<x_mid<<x_right<<endl;
       else
           cout<<"New interval is       "<<x_left<<x_right;
           }
      
       return ((x_left + x_right)/2.0);
       }
      
       double g(double x) {
              return (5 * pow(x, 3.0) -2 *pow(x,2.0)+3);
              }
             
       double h(double x){
              return(pow(x, 4.0) - 3 * pow(x,2.0) - 8);
              };