Finding the root of a function by Bisection Method The program has to look for a root in an interval [a,b]. The root should be declared with a certain accuracy eps. I.e it should look for a part-interval in [a,b], which has the length of eps.

The bisection algorithm should be:


Save the interval boundaries
Look if [a,b] has a root. (original given interval)
look if a-b &lt; eps. If yes, part-interval found. 
If no, divide [a,b] in half and continue with point 2. etc.


(We can assume that there is already a root in the given original interval [a,b])

I've tested it with following functions:

log x + 2

x*x*x

-3.*x*x + x + 2.

Are there any improvements possible or did I make a big mistake?

Bisection.cpp

#include "Bisection.hpp"
#include &lt;cmath&gt;

using namespace std;


void bisection(std::vector&lt;double&gt;&amp; interval, double eps) {
  double a;
  double b;
  double c;

  a = interval[0];
  b = interval[1];

  if (b - a &lt; eps) // interval already 'found'
  {
    return;
  }

  do
  {
    c = (a + b) / 2.0;

    if (f(c) == 0)
    {
      a = c;
      b = c; 
    }
    else 
    {
      if (f(a) * f(c) &gt; 0)
      {
        a = c;
      }
      else
      {
        b = c;
      }                  
    }


  } while (f(c) != 0 || b - a &gt; eps);

  /*found part-interval with the length of &lt; eps*/
  interval[0] = a;
  interval[1] = b;
}
A std::vector&lt;T&gt; vec is an arbitrary collection of Ts. vec might contain none, one or many elements. However, you want to model an interval, which consists of exactly two elements. A std::vector&lt;T&gt; therefore isn't the correct container for an interval.

You could use std::pair&lt;double, double&gt; instead, or just write your own variant:

struct interval_type {
    double lower_bound;
    double upper_bound;
};


Next, try to initialize the variables instead of just declaring them. While we're at it, try to come up with better names, e.g. instead of

  double a;
  double b;
  double c;

  a = interval[0];
  b = interval[1];


write

  double lower = interval.lower_bound; // if you use the struct from above
  double upper = interval.upper_bound;


c's scope is too large, by the way. It's just necessary inside the loop, so we can limit its life time, but we need to break in this case:

  do
  {
    const double midpoint = (lower + upper) / 2.0;

    if (f(midpoint) == 0)
    {
      a = midpoint;
      b = midpoint; 
      break;
    }
    else 
    {
      if (f(lower) * f(midpoint) &gt; 0)
      {
        lower = midpoint;
      }
      else
      {
        upper = midpoint;
      }                  
    }
  } while (upper - lower &gt; eps);


Other than that, the test f(value) == 0 is usually too strict, since floating point arithmetic can be a finicky beast. std::abs(f(value)) &lt;= some_eps works well enough most of the times.

Other than that, I suggest you to allow other functions than a fixed f. Something along

template &lt;typename F&gt;
interval_type bisection(F &amp;&amp; f, interval_type interval, double eps = 1e-8){

}


comes to mind. As you can see, there is no reason to change the interval in place. Just return a new one.