bisection.h
Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00034 #ifndef MAT_BISECTION
00035 #define MAT_BISECTION
00036 #include <cmath>
00037 namespace mat {
00045 template<typename Treal>
00046 inline int sign(Treal value) {
00047 if (value > 0)
00048 return 1;
00049 else if (value < 0)
00050 return -1;
00051 else
00052 return 0;
00053 }
00054
00055
00067 template<typename Treal, typename Tfun>
00068 Treal bisection(Tfun const & fun, Treal min, Treal max, Treal const tol) {
00069 int sign_min = sign(fun.eval(min));
00070 int sign_max = sign(fun.eval(max));
00071 if (sign_min == sign_max)
00072 throw Failure("bisection(Tfun&, Treal, Treal, Treal): interval "
00073 "incorrect");
00074 Treal middle = (max + min) / 2;
00075 int sign_middle = sign(fun.eval(middle));
00076 while (template_blas_fabs(max - min) > tol * 2 && sign_middle != 0) {
00077 if (sign_middle == sign_min) {
00078 min = middle;
00079 sign_min = sign_middle;
00080 }
00081 else {
00082 max = middle;
00083 sign_max = sign_middle;
00084 }
00085 middle = (max + min) / 2;
00086 sign_middle = sign(fun.eval(middle));
00087 }
00088 return middle;
00089 }
00090
00091 }
00092 #endif