/* Bisection Method for Finding the Square Root of a Positive Number */
	  
#include <stdio.h>
#include <assert.h>

const double ACCURACY = 0.0001;  /* desired accuracy of result */

/* square_root - Calculate the square root of a positive number
 *
 *  pre-conditions:  value>0
 *  post-conditions:  returns m such that |m*m - t| <= ACCURACY
 */
double
square_root (double value)
{
  double left, right, midpoint;  /* root will be in interval [left,right] */

  /* for f(x) = x^2 - t, differences represent values 
     f(left), f(right), and f(midpoint) */
  double diffLeft, diffRight, diffMid;  

  /* Set up initial interval for the bisection method */
  left = 0;
  right = (value > 2.0 ? 2.0 : value );
   
  diffLeft = left*left - value;
  diffRight = right*right - value;


  /* Iterate interval shrinking until sufficiently small */
  while ( (right - left) > ACCURACY)
  {
    /* f(x) = x^2 - t must have opposite signs at f(left) and f(right) */
    assert (diffLeft * diffRight < 0);  
      
    midpoint = (left + right) / 2.0;  /* midpoint of range [left,right] */
    diffMid = midpoint*midpoint - value;
    
    if (diffMid == 0.0) break;  /* stop loop if we have the exact root */
    
    if ((diffLeft * diffMid) < 0.0) {
      /* f(left) and f(mid) have opposite signs */
      right = midpoint;
      diffRight = diffMid;
    } else {
      left = midpoint;
      diffLeft = diffMid;
    }
  } // while

  return value;
  
} // square_root


/* Prompt the user for a number, calculate and print its square root */
int
main (void)
{
  double square;        /* we approximate the square root of this number */
  double root;
  int numAssigned;      /* return value for user input from scanf */

  /* Getting started */
  printf ("Program to compute a square root\n");
  printf ("Enter positive number: ");
  numAssigned = scanf ("%lf", &square);

  /* Verify user input */
  if (numAssigned == 0) {
    printf ("Error: Unable to read number\n");
    return 1;
  } else if (square <= 0) {
    printf ("Error: A positive number is required; %lf was entered.\n", square);
    return 1;
  }

  /* Calculate result */
  root = square_root (square);
  printf ("The square root of %lf is approximately %lf\n", square, root);

  return 0;
} // main