Effective implementation of algorithms (Master Thesis)
Effective and error-free implementation of algorithms
|
00001 #ifndef H_MATH_BINSEARCH_FUNCTION_BINSEARCH 00002 #define H_MATH_BINSEARCH_FUNCTION_BINSEARCH 00003 00009 #include "utils/preconditions/preconditions.h" 00010 #include <cmath> 00011 00012 namespace math { 00013 namespace binsearch { 00014 00018 template<typename T> class Function 00019 { 00020 public: 00021 virtual T operator()(T x) = 0; 00022 }; 00023 00031 template<typename T> 00032 class ConvexFunction: public Function<T> { 00033 }; 00034 00035 template<typename T> 00036 class FunctionBinsearch { 00037 private: 00049 int number_of_iterations(T left, T right, T precision, double 00050 division) { 00051 Preconditions::check(right > left, "invalid range"); 00052 Preconditions::check(precision > 0, "precision must be positive"); 00053 00054 double tmp = log(right - left) - log(precision); 00055 return 2 + (int) (tmp / log(division)); 00056 } 00057 00058 public: 00072 T root(Function<T>* f, T left, T right, T precision) 00073 { 00074 Preconditions::check(right >= left, "invalid range"); 00075 Preconditions::check(precision > 0, "precision must be positive"); 00076 Preconditions::check((*f)(left) * (*f)(right) <= 0, 00077 "endpoint should have different sign"); 00078 00079 if (right - left <= precision) { 00080 return left; 00081 } 00082 00083 assert(left < right); 00084 00085 int sgn = (*f)(left) > 0 ? 1 : -1; 00086 00087 int iterations = number_of_iterations(left, right, precision, 2.0); 00088 while (iterations--) { 00089 T middle = (left + right) / 2; 00090 if (sgn * (*f)(middle) > 0) { 00091 left = middle; 00092 } else { 00093 right = middle; 00094 } 00095 } 00096 return left; 00097 } 00098 00114 T convex_min(ConvexFunction<T>* f, 00115 T left, T right, T precision) 00116 { 00117 Preconditions::check(right >= left); 00118 Preconditions::check(precision > 0); 00119 00120 if (right - left <= precision) { 00121 return left; 00122 } 00123 00124 assert(left < right); 00125 00126 int iterations = number_of_iterations(left, right, precision, 1.5); 00127 while (iterations--) { 00128 T mid1 = (left * 2 + right ) / 3; 00129 T mid2 = (left + right * 2) / 3; 00130 00131 if ((*f)(mid1) < (*f)(mid2)) { 00132 right = mid2; 00133 } else { 00134 left = mid1; 00135 } 00136 } 00137 00138 return left; 00139 } 00140 }; 00141 00142 } // namespace binsearch 00143 } // namespace math 00144 #endif