An implementation of an improved & simplified Brent's Method. Calculates root of a function f(x) in the interval [a,b]. Zhang, Z. (2011). An Improvement to the Brent’s Method. International Journal of Experimental Algorithms (IJEA), (2), 21–26. Retrieved from http://www.cscjournals.org/csc/manuscript/Journals/IJEA/volume2/Issue1/IJEA-7.pdf

gistfile1.cpp

template <class T>
inline void __swap(T *a, T *b)
{
	T temp = *a;
	*a = *b;
	*b = temp;
}

/*
An implementation of an improved & simplified Brent's Method.

Calculates root of a function f(x) in the interval [a,b].

Zhang, Z. (2011). An Improvement to the Brent’s Method. International Journal of Experimental Algorithms (IJEA), (2), 21–26. 
Retrieved from http://www.cscjournals.org/csc/manuscript/Journals/IJEA/volume2/Issue1/IJEA-7.pdf

f -> Function pointer to the function to be evaluated
a -> Starting point of interval
b -> Ending point of interval
root -> This will contain the root when the functiuon is complete
tol -> tolerance . Set to a low value like 1e-6

Returns zero on success
Returns -1 on failure
*/
template <class F>
int cpu_brent_improved(double (*)(double), double a, double b, double& root, double tol)
{
	double s, c, f_s, f_a, f_b, f_c;
	f_a = f(a);
	f_b = f(b);
	
	if(a>b || f_a * f_b >= 0)	// The algorithm assumes a<b
	{
		return -1; // Root not bound by the given guesses
	}

	do{
		c = (a+b)/2;
		f_c = f(c);
		if(fabs(f_a-f_c) > tol && fabs(f_b-f_c) > tol) 	// f(a)!=f(c) and f(b)!=f(c)
		{
			// Inverse quadratic interpolation
			s = a*f_b*f_c/((f_a-f_b)*(f_a-f_c)) + b*f_a*f_c/((f_b-f_a)*(f_b-f_c)) + c*f_a*f_b/((f_c-f_a)*(f_c-f_b));
		}else{
			// Secant rule
			s = b - f_b*(b-a)/(f_b-f_a);
		}
		f_s = f(s);
		if(c>s)
		{
			__swap(&s, &c);
		}
		if(f_c*f_s < 0)
		{
			a = s;
			b = c;
		}else if(f_s*f_b < 0){
			a = c;
		}else{
			b = s;
		}
	}while(f_s < tol || fabs(b-a)<tol);	// Convergence conditions
	
	root = s;
	return 0;
}

There's another article by Steven Stage that points out a couple of problems with the Zhang algorithm. It also profiles its performance vs. Brent v. Regula Falsi. Check out: https://www.cscjournals.org/manuscript/Journals/IJEA/Volume4/Issue1/IJEA-33.pdf

I copied this and tweaked it slightly for my environment, fixed the templated arg issue, then modified per Stage. However, I didn't get proper convergence until I fixed a further problem with the above code, namely that after a and b are calculated f_a and/or f_b need to be recalculated. If you do this at the top of the do loop, there may be redundant calculation, whereas if you do it after they are modified, the recalculation might be superfluous.