I'm glad I was of help
Actually, after doing a quick search on the web for Newton's Method, I changed my solve() and findRoot() functions. It's much simpler, and (I think) will yield better results for the same number of iterations. The code for solveNewton() and findRootNewton() is:
Code:
double findRootNewton(const Polynomial& p, const Polynomial& dp, double x1, double x2, int prec = 100)
{
double x = (x1 + x2) / 2;
for(int i = 0; i < prec; ++i)
x -= (p.evaluate(x) / dp.evaluate(x));
return x;
}
std::deque<double> solveNewton(Polynomial& p)
{
std::deque<double> ret;
if(p.order() < 1)
return ret;
if(p.order() == 1)
{
ret.push_back(-(p.coefficient(0) / p.coefficient(1)));
return ret;
}
Polynomial der = p.derivative();
std::deque<double> extrema = solveNewton(der);
if(extrema.size() > 0)
{
extrema.push_front(extrema.front() - 1.);
extrema.push_back(extrema.back() + 1.);
}
else
{
extrema.push_front(-1.);
extrema.push_back(1.);
}
std::deque<double>::size_type i, j;
j = extrema.size() - 1;
for(i = 0; i < j; ++i)
{
double root = findRootNewton(p, der, extrema[i], extrema[i + 1]);
if(root != std::numeric_limits<double>::min())
ret.push_back(root);
}
return ret;
}
The only potential problem I see, is that I have no idea how Newton's Method really works (just grabbed the formula "x_(n+1) = x - f(x)/f'(x)" from somewhere) - and so I don't know for sure if this method will only converge on a root between the extrema on the left and right of the current x value, or if it will potentially find the same root several times while omitting others. So far as I have tested, however, it seems to hold up OK.