Error/convergence analysis using Taylor expansion for a root-finding method

Question

I want to prove that this root-finding scheme has a quartic convergence rate $$\left\{ \begin{array}{c} y_{n}=x_{n}-\frac{f(x_{n})}{f^{\prime }(x_{n})}, n=0,1,2,...., \\ v_{n}=x_{n}-\frac{f(x_{n})-\frac{9}{30}f(y_{n})}{f(x_{n})-\frac{39}{30}% f(y_{n})}\left( \frac{f(x_{n})}{f^{\prime }(x_{n})}\right) , \\ x_{n+1}=v_{n}-\frac{f(v_{n})}{\frac{f(v_{n})-f(x_{n})}{v_{n}-x_{n}}}.% \end{array}% \right. $$ and the error below: $$\mu _{n+1}=\frac{7b_{2}^{3}}{10}\mu ^{4}+O\left( \mu ^{5}\right), $$ where $\mu =x_{n}-\alpha $ and $b_{k}=\frac{1}{k!}\frac{ f^{(k)}(\alpha )}{f^{\prime }(\alpha )}$, $k=2,3,4,\ldots $.

For proof of that, we will use Taylor expansion first for $y_{n}$, then for $v_{n}$, and finally for $x_{n+1}$. I got it for $y_{n}$ successfully but for $v_{n}$, it was written that

$$v(x_{n})=\alpha +% \frac{7}{10}b_{2}^{2}\mu ^{3}+\left( \frac{9b_{2}b_{3}}{5}-\frac{% 159b_{2}^{3}}{100}\right) \mu ^{4}+O\left( \mu ^{5}\right).$$

Now, how to get this relation for $v(x_{n})$?

My try for this purpose, I used Taylor expansion as follows

$$f(x_{n})=f^{\prime }(\alpha )\left( \mu +b_{2}\mu^{2}+b_{3}\mu ^{3}+b_{4}\mu ^{4}+b_{5}\mu ^{5}\right)+O( \mu ^{6})$$ and $$f^{\prime }(x_{n})=f^{\prime }(\alpha )\left( 1 +2b_{2}\mu+3b_{3}\mu ^{2}+4b_{4}\mu ^{3}+5b_{5}\mu ^{5}\right)+O( \mu ^{5})$$ then I applied the binomial theorem with negative power to get that $$\frac{f(x_{n})}{f^{\prime }(x_{n})}=\mu - b_{2}\mu ^{2}+\left( 2b_{2}^{2}-2b_{3}\right) \mu ^{3}+\mu ^{4}\left( -4b_{2}^{3}+7b_{3}b_{2}-3b_{4}\right) +O( \mu ^{5})$$ So $$f(y_{n})=f\left( x_{n}-\frac{f(x_{n})}{f^{\prime }(x_{n})}\right) =f^{\prime }(\alpha )\left( b_{2}\mu ^{2}+\left( -2b_{2}^{2}+2b_{3}\right) \mu ^{3}+\left( 5b_{2}^{3}-7b_{3}b_{2}+3b_{4}\right) \mu ^{4}\right)+O( \mu ^{5}) $$

Then I got stuck here and I could not get that for $v(x_{n})$.

Answer 1

rough outline of the method

The next value $x_{+1}$ is the secant root computed for the points $(x,f(x))$ and $(v,f(v))$. It is well-known that the error propagation for this is $$ x_{+1}-\alpha = C\,(x-\alpha)(v-\alpha)~~\text{ or }~~ f(x_{+1})=C'\, f(x)f(v) $$ To get an order $4$ error reduction, that is, $f(x_{+1})=O(f(x)^4)$ we need that $v$ is an order 3 approximation, $f(v)=O(f(x)^3)$. The Newton value $y=x+s$ has order 2. Its value $f(y)$ can be used to remove derivatives from third order methods like Halley's method $$ v=x-\frac{f(x)f'(x)}{f'(x)^2-\frac12f''(x)f(x)}=x+\frac{f(x)}{f(x)-f(y)}s $$ or the Euler-Chebyshev method $$ v=x+s-\frac12f'(x)^{-1}f''(x)s^2=x+\frac{f(x)+f(y)}{f(x)}s $$ or any variant in-between like in the proposed method.

Euler-Chebyshev method and variants in detail

Start with the Newton update $s=y-x=-\frac{f(x)}{f'(x)}$. Then any correction to it computes as $v=x+s+\gamma s$ with $\gamma =O( f(x))$ and the function value $$ f(y+\gamma s)=f(x)+f'(x)(s+\gamma s)+\frac12f''(x)(1+\gamma)^2s^2+... \\ =-\gamma f(x)+\frac12f''(x)(1+2\gamma+\gamma^2)s^2+.. $$ This allows in a first step to replace $\frac12f''(x)s^2+...$ with $f(y)$. The remaining terms on the lowest order $O(f(x)^2)$ are $-\gamma f(x)+f(y)$, so set $\gamma =\frac{f(y)}{f(x)}$, with modifications of order $O(f(x)^2)$ possible, like $\gamma=\frac{f(y)}{f(x)+\kappa f(y)}$.

The first variant inserted gives the new leading error terms, now in order higher, as \begin{align} \frac{f(y)}{f(x)^2}&=\frac12\frac{f''(x)}{f'(x)^2}+\frac16\frac{f'''(x)f(x)}{f'(x)^3}+\dots \\ f(v)&=-f(y)+\frac12f''(x)\frac{(f(x)+f(y))^2}{f'(x)^2}+\frac16f'''(x)\frac{(f(x)+f(y))^3}{f'(x)^3}+\dots \\ &=-f(y)+\left(\frac{f(y)}{f(x)^2}-\frac16\frac{f'''(x)f(x)}{f'(x)^3}\right)(f(x)+f(y))^2+\frac16f'''(x)\frac{(f(x)+f(y))^3}{f'(x)^3}+\dots \\ &=2\frac{f(y)^2}{f(x)}\;+\;\frac{f(y)^3}{f(x)^2}+\frac16\frac{f'''(x)f(x)^2f(y)}{f'(x)^3}+\dots \end{align} Only the first term is of order 3.

Secant root

Translating the secant root error formula into function values gives $$ f(x_{+1})=\frac{f''(x)}{2f'(x)^2}f(x)f(v)+\dots =2\frac{f(y)^3}{f(x)^2}+\dots $$

Now you are free to replace $f(x)=\mu+b_2\mu^2+b_3\mu^3+...$, so that $f(y)=b_2\mu^2$ and $$f(x_{+1})=-2b_2^3\mu^4+...$$

I do not see where factors $5$ and $7$ or the value $\kappa=-\frac{39}{30}$ as an optimal parameter choice would originate.

Using CAS Magma for Taylor series computations

Using computer algebra for the Taylor series operations confirms the claim against above insights/objections.

CF<b2,b3,b4,b5>:=FunctionField(Rationals(),4);
A<k>:=FunctionField(CF);
PS<u>:=PowerSeriesRing(A);

f := u+b2*u^2+b3*u^3+b4*u^4+b5*u^5;
df := 1+2*b2*u+3*b3*u^2+4*b4*u^3+5*b5*u^4;

x := u+O(u^6);
fx := Evaluate(f,x);
dfx := Evaluate(df,x);

"Newton step";
s := PS!(-fx/dfx); "s = ", s;
y := x+s; "y = ", y;
fy := Evaluate(fx,y); "f(y) = ", fy;

"Chebyshev step";
gamma := 1+PS!(fy/(fx+k*fy)); "1+gamma = ", gamma;

v := x+gamma*s; "v = ", v;
fv := Evaluate(fx,v); "f(v) = ", fv;

"Secant root";
sr := PS!((u*fv-v*fx)/(fv-fx)); "secant root = ", sr;
"with value = ", Evaluate(fx,sr);

with output

Newton step

s =  -u + b2*u^2 + (-2*b2^2 + 2*b3)*u^3 + (4*b2^3 - 7*b2*b3 + 3*b4)*u^4 + 
    (-8*b2^4 + 20*b2^2*b3 - 10*b2*b4 - 6*b3^2 + 4*b5)*u^5 + O(u^6)
y =  b2*u^2 + (-2*b2^2 + 2*b3)*u^3 + (4*b2^3 - 7*b2*b3 + 3*b4)*u^4 + (-8*b2^4 + 
    20*b2^2*b3 - 10*b2*b4 - 6*b3^2 + 4*b5)*u^5 + O(u^6)
f(y) =  b2*u^2 + (-2*b2^2 + 2*b3)*u^3 + (5*b2^3 - 7*b2*b3 + 3*b4)*u^4 + 
    (-12*b2^4 + 24*b2^2*b3 - 10*b2*b4 - 6*b3^2 + 4*b5)*u^5 + O(u^6)

Chebyshev step

1+gamma =  1 + b2*u + (-b2^2*k - 3*b2^2 + 2*b3)*u^2 + (b2^3*k^2 + (6*b2^3 - 
    4*b2*b3)*k + (8*b2^3 - 10*b2*b3 + 3*b4))*u^3 + (-b2^4*k^3 + (-9*b2^4 + 
    6*b2^2*b3)*k^2 + (-25*b2^4 + 32*b2^2*b3 - 6*b2*b4 - 4*b3^2)*k - 20*b2^4 + 
    37*b2^2*b3 - 14*b2*b4 - 8*b3^2 + 4*b5)*u^4 + O(u^5)
v =  (b2^2*k + 2*b2^2)*u^3 + (-b2^3*k^2 + (-7*b2^3 + 4*b2*b3)*k - 9*b2^3 + 
    7*b2*b3)*u^4 + (b2^4*k^3 + (10*b2^4 - 6*b2^2*b3)*k^2 + (33*b2^4 - 38*b2^2*b3
    + 6*b2*b4 + 4*b3^2)*k + (30*b2^4 - 44*b2^2*b3 + 10*b2*b4 + 6*b3^2))*u^5 + 
    O(u^6)
f(v) =  (b2^2*k + 2*b2^2)*u^3 + (-b2^3*k^2 + (-7*b2^3 + 4*b2*b3)*k - 9*b2^3 + 
    7*b2*b3)*u^4 + (b2^4*k^3 + (10*b2^4 - 6*b2^2*b3)*k^2 + (33*b2^4 - 38*b2^2*b3
    + 6*b2*b4 + 4*b3^2)*k + (30*b2^4 - 44*b2^2*b3 + 10*b2*b4 + 6*b3^2))*u^5 + 
    O(u^6)

Secant root

secant root =  (b2^3*k + 2*b2^3)*u^4 + (-b2^4*k^2 + (-8*b2^4 + 5*b2^2*b3)*k - 
    11*b2^4 + 9*b2^2*b3)*u^5 + O(u^6)
with value =  (b2^3*k + 2*b2^3)*u^4 + (-b2^4*k^2 + (-8*b2^4 + 5*b2^2*b3)*k - 
    11*b2^4 + 9*b2^2*b3)*u^5 + O(u^6)

Higher order?

Looking at the leading coefficient, the variant with $k=\kappa=-2$ increases the order in $v$. Then computing the secant from the points $y$ and $v$ gives a method of order $6$, which is better than 2 Newton steps with less evaluation effort for it. 2 Newton steps followed by a secant step gives the same order, but with more evaluations.