Brief explanation of the problem: I use the Newton-Raphson algorithm to search for roots in polynomials and does not work in some cases. why?
I took the Newton Raphson hybrid algorithm from "numeric recipes in C ++", which bisects if New-Raph does not converge properly (with a low derivative value or convergence rate is not fast).
I checked the algorithm with several polynomials and it worked. Now I am testing inside the software that I have and I always had a bug with a specific polynomial. My problem is that I do not know why this polynomial simply does not achieve a result when many others do it. Since I want to improve the algorithm for any polynomial y, I need to know which one is the reason for the lack of convergence, so I can handle it correctly.
The following I will post all the information I can provide about the algorithm and polynomial in which I have an error.
Polynomial:
f(t)= t^4 + 0,557257315256597*t^3 - 3,68254086033178*t^2 + + 0,139389107255627*t + 1,75823776590795
First derivative:
f'(t)= 4*t^3 + 1.671771945769790*t^2 - 7.365081720663563*t + 0.139389107255627
Plot: 
Roots (by Matlab):
-2.133112008595826 1.371976341295347 0.883715461977390 -0.679837109933505
Algorithm:
double rtsafe(double* coeffs, int degree, double x1, double x2,double xacc,double xacc2) { int j; double df,dx,dxold,f,fh,fl; double temp,xh,xl,rts; double* dcoeffs=dvector(0,degree); for(int i=0;i<=degree;i++) dcoeffs[i]=0.0; PolyDeriv(coeffs,dcoeffs,degree); evalPoly(x1,coeffs,degree,&fl); evalPoly(x2,coeffs,degree,&fh); evalPoly(x2,dcoeffs,degree-1,&df); if ((fl > 0.0 && fh > 0.0) || (fl < 0.0 && fh < 0.0)) nrerror("Root must be bracketed in rtsafe"); if (fl == 0.0) return x1; if (fh == 0.0) return x2; if (fl < 0.0) { // Orient the search so that f(xl) < 0. xl=x1; xh=x2; } else { xh=x1; xl=x2; } rts=0.5*(x1+x2); //Initialize the guess for root, dxold=fabs(x2-x1); //the "stepsize before last," dx=dxold; //and the last step evalPoly(rts,coeffs,degree,&f); evalPoly(rts,dcoeffs,degree-1,&dx); for (j=1;j<=MAXIT;j++) { //Loop over allowed iterations if ((((rts-xh)*df-f)*((rts-xl)*df-f) > 0.0) //Bisect if Newton out of range, || (fabs(2.0*f) > fabs(dxold*df))) { //or not decreasing fast enough. dxold=dx; dx=0.5*(xh-xl); rts=xl+dx; if (xl == rts) return rts; //Change in root is negligible. } else {// Newton step acceptable. Take it. dxold=dx; dx=f/df; temp=rts; rts -= dx; if (temp == rts) return rts; } if (fabs(dx) < xacc) return rts;// Convergence criterion evalPoly(rts,coeffs,degree,&f); evalPoly(rts,dcoeffs,degree-1,&dx); //The one new function evaluation per iteration. if (f < 0.0) //Maintain the bracket on the root. xl=rts; else xh=rts; } //As the Accuracy asked to the algorithm is really high (but usually easily reached) //the results precission is checked again, but with a less exigent result dx=f/df; if(fabs(dx)<xacc2) return rts; nrerror("Maximum number of iterations exceeded in rtsafe"); return 0.0;// Never get here. }
The algorithm is called with the following variables:
x1=0.019 x2=1.05 xacc=1e-10 xacc2=0.1 degree=4 MAXIT=1000 coeffs[0]=1.75823776590795; coeffs[1]=0.139389107255627; coeffs[2]=-3.68254086033178; coeffs[3]=0.557257315256597; coeffs[4]=1.0;
The problem is that the algorithm is superior to the iterations of amximum and there is an upwind to the root aproximatedly 0.15 .
So, my direct and abbreviated question is: why does this polynomial fail to achieve an exact error when many (1000, at least) other very similar polynomials (up to 1e-10 accuracy and a few iterations!)
I know that the question is complex, and it may not have a really direct answer, but I was stuck with this for several days, and I do not know how to solve it. Thank you very much for taking the time to read my question.