#include<stdio.h>
#include<iostream>
using namespace std;

double f (double x) {
    return (x * (1-x));
}

int ternary_search(double lo, double hi) {
    /* We know that f() is a unimodal function in the range [lo,hi] --
       Specifically there are points x0=lo <= x1 <= x2 <= x3=hi where the
       region [x0,x1] is increasing, the region [x1,x2] is constant
       and the region [x2,x3] is decreasing.  Your algorithm knows x0
       and x3. The goal is to find a point of maximum value in
       [x0,x3].  Each of these three regions is possibly of zero
       length.
       
       It works both when the range is continuous or discrete, and the
       function value is continuous or discrete.

       In this example we use the continuous case for both.
    */
    
    double x0,x1,x2,x3,f1,f2;
    x0 = lo;
    x3 = hi;

    /* Induction hypothesis at start of loop:
       (1) x0 < x3
       (1) the function is unimodal in the range [x0,x3]
       (2) the range [x0,x3] contains a maximum of the function
           over the range [lo,hi] 
    */
    for (int count = 0; count < 100; count++) {
	x1 = (2*x0 + x3)/3;
	x2 = (x0 + 2*x3)/3;
	f1 = f(x1);
	f2 = f(x2);
	if (f1 < f2) x0 = x1;
	else x3 = x2;
    }
    return (x0 + x3)/2;
}

int main() {
    printf ("searching in [%.2f, %.2f] minimum at %f\n", 0.0, 1.0, ternary_search(0.0, 1.0));
    printf ("searching in [%.2f, %.2f] minimum at %f\n", 1.0, 10.0, ternary_search(1.0, 10.0));
    printf ("searching in [%.2f, %.2f] minimum at %f\n", -10.0, 10.0, ternary_search(-10.0, 10.0));	    
    return 0;
}