/* Find an approximate solution of the equation F(x) = 0 in a given
   interval, using the bisection method.

   Input:   Desired accuracy of approximation, endpoints of an
            interval containing a solution
   Output:  Prompts to the user and the approximate solution

Notes:
   1. The solution will be within DesiredAccuracy of the exact
      solution.
   2. To find the solution for a different F(x), redefine F() below.
------------------------------------------------------------------*/

double F(double x)
{
   return x*x*x + x - 5;       // the function F for which the
}                              // equation F(x) = 0 is being solved

//-----------------------------------------------------------------------

#include <iostream.h>

double Bisector(double, double, double);

int main(void)
{
   cout << "\nThis program uses the bisection method to find an\n"
        << "\tapproximate solution to the equation F(x) = 0.\n";

   double
      DesiredAccuracy;      // the accuracy desired

   cout << "\nEnter the accuracy desired (e.g. .001): ";
   cin >> DesiredAccuracy;

   double
      Left,
      Right;

   do                       // get the interval containing a solution
   {
      cout << "Enter the x-values of interval containing solution: ";
      cin >> Left >> Right;
   }
   while (F(Left) * F(Right) >= 0.0);


   double                   // solve the equation
      Solution = Bisector(Left, Right, DesiredAccuracy);

   cout << "\n--> "<< Solution
        << " is an approximate solution of F(x) = 0, to within "
        << DesiredAccuracy << ".\n\n";

   return 0;
}

/*------------------------------------------------------------------
Bisector performs the bisection algorithm.

   Receive:  Left, the left endpoint of the original interval
             Right, the right endpoint of the original interval
             Accuracy, the desired accuracy of the approximation
   Return:   MidPt, the middle of the final interval
-------------------------------------------------------------------*/

double Bisector(double Left, double Right, double Accuracy)
{
   double
      Width = Right - Left,              // the interval width
      MidPt,                             // the midpoint of the interval
      fMid;                              // value of F at Midpoint

   while (Width/2.0 > Accuracy)          // While loop
   {
      MidPt = (Left + Right) / 2.0;      //   compute midpoint
      fMid = F(MidPt);                   //   compute function at midpoint

      if (F(Left) * fMid < 0.0 )         //   solution is in left half
         Right = MidPt;
      else                               //   solution is in right half
         Left = MidPt;

      Width /= 2.0;                      //   split the interval
   }                                     // End loop

   return MidPt;
}