/*  Bisection Method (15 points)

MATH BACKGROUND:
================
	
The bisection method is one of the simplest algorithms for finding the roots 
of a function f(x) (i.e., for finding when f(x)=0 ). To use this method one 
must provide an initial interval [a,b] that contains or brackets the root, 
such that f(a) and f(b) have opposite signs.

The algorithm repeatedly divides the interval [a,b] into two equal halves 
and focuses its attention on the half that continues to bracket the root. 
This condition is met if the function continues to have opposite signs at the 
beginning and at the end of the subinterval.

More formally the iterative algorithm consists of the following steps:
	
  Step 1: calculate the mid-point c = (a+b)/2 of the interval [a,b]

  Step 2: if f(c) = 0, then c is the root and we are done
          else if f(a)*f(c) < 0, then make  b = c 
          else if f(c)*f(b) < 0, then make  a = c

  Step 3: if |a-b| < EPSILON, then stop

  Step 4: Otherwise go to Step 1 using the new values for and and b


PROBLEM DESCRIPTION:
====================

Write a complete C program that implements the bisection method and uses it
to calculate the root of a given function in a user specified interval 
(e.g., f(x) = x*x -2 in the interval [1,2]). You can hardcode several functions
in your code and uncomment one of them for testing. Something like this will do:

 double f(double x)
 {
	//return x*x -2;  // the root is sqrt(2)
	return x*x -4;  // the root is 2
 }

Format the output similar to the sample runs provided below. 
These values were obtained with EPSILON=0.00001 and using %.7f to print them.


SAMPLE RUN #1:  ( f(x) = x*x -2; the root should be sqrt(2) )
	
Please enter a: 1
Please enter b: 2
      a           b           b-a         f(a)       f(b)
-----------------------------------------------------------
[ 1.0000000,  2.0000000]   1.0000000  -1.0000000  2.0000000
[ 1.0000000,  1.5000000]   0.5000000  -1.0000000  0.2500000
[ 1.2500000,  1.5000000]   0.2500000  -0.4375000  0.2500000
[ 1.3750000,  1.5000000]   0.1250000  -0.1093750  0.2500000
[ 1.3750000,  1.4375000]   0.0625000  -0.1093750  0.0664063
[ 1.4062500,  1.4375000]   0.0312500  -0.0224609  0.0664063
[ 1.4062500,  1.4218750]   0.0156250  -0.0224609  0.0217285
[ 1.4140625,  1.4218750]   0.0078125  -0.0004272  0.0217285
[ 1.4140625,  1.4179688]   0.0039063  -0.0004272  0.0106354
[ 1.4140625,  1.4160156]   0.0019531  -0.0004272  0.0051003
[ 1.4140625,  1.4150391]   0.0009766  -0.0004272  0.0023355
[ 1.4140625,  1.4145508]   0.0004883  -0.0004272  0.0009539
[ 1.4140625,  1.4143066]   0.0002441  -0.0004272  0.0002633
[ 1.4141846,  1.4143066]   0.0001221  -0.0000820  0.0002633
[ 1.4141846,  1.4142456]   0.0000610  -0.0000820  0.0000906
[ 1.4141846,  1.4142151]   0.0000305  -0.0000820  0.0000043
[ 1.4141998,  1.4142151]   0.0000153  -0.0000388  0.0000043

The root of the function in this interval is x=1.41421


SAMPLE RUN #2:  ( f(x) = x*x -4; the root should be 2 )

Please enter a: -1
Please enter b: 3
      a           b           b-a         f(a)       f(b)
-----------------------------------------------------------
[-1.0000000,  3.0000000]   4.0000000  -3.0000000  5.0000000
[ 1.0000000,  3.0000000]   2.0000000  -3.0000000  5.0000000
Found it at c=2.0000000!

The root of the function in this interval is x=2

*/ 

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

double f(double x)
{
	return x*x -4;  // the root is 2
	//return x*x -2;  // the root is sqrt(2)
}

#define EPSILON 0.00001

int main(void)
{
	double a, b,c;
	
	do
	{
		printf("Please enter a: "); 	scanf("%lf", &a);
		printf("Please enter b: "); 	scanf("%lf", &b);
		if(a>b)
			printf("  Error: a must be smaller than b.\n\n");
		if(f(a) * f(b) > 0)
			printf("  Error: f(a)=%g and f(b)=%g must have different signs.\n\n", f(a), f(b));
	} while((a>b) || (f(a) * f(b) > 0));
	
	
	printf("      a           b           b-a         f(a)       f(b)\n");
	printf("-----------------------------------------------------------\n");
	
	do
 	{
		printf("[% .7f, % .7f]  % .7f  % .7f % .7f\n",a, b, b-a, f(a), f(b));
		
		c = (a+b)/2.0;

		if( f(a) * f(c) <0 )
			b=c;
		else if( f(c)*f(b)<0 )
			a=c;
		else
		{
			printf("Found it at c=%.7f!\n", c);
			break;   // We are at the root so we can stop
		}
		
	} while( fabs(a - b) > EPSILON);
	
	printf("\nThe root of the function in this interval is x=%g\n\n", c);
	
    system("pause");
}