Sidi's generalized secant method

Sidi's generalized secant method is a root-finding algorithm, that is, a numerical method for solving equations of the form $f(x)=0$ . The method was published by Avram Sidi.^[1]

The method is a generalization of the secant method. Like the secant method, it is an iterative method which requires one evaluation of $f$ in each iteration and no derivatives of $f$. The method can converge much faster though, with an order which approaches 2 provided that $f$ satisfies the regularity conditions described below.

Algorithm

We call $\alpha$ the root of $f$, that is, $f(\alpha)=0$. Sidi's method is an iterative method which generates a sequence $\{ x_i \}$ of approximations of $\alpha$. Starting with k + 1 initial approximations $x_1 , \dots , x_{k+1}$, the approximation $x_{k+2}$ is calculated in the first iteration, the approximation $x_{k+3}$ is calculated in the second iteration, etc. Each iteration takes as input the last k + 1 approximations and the value of $f$ at those approximations. Hence the nth iteration takes as input the approximations $x_n , \dots , x_{n+k}$ and the values $f(x_n) , \dots , f(x_{n+k})$.

The number k must be 1 or larger: k = 1, 2, 3, .... It remains fixed during the execution of the algorithm. In order to obtain the starting approximations $x_1 , \dots , x_{k+1}$ one could carry out a few initializing iterations with a lower value of k.

The approximation $x_{n+k+1}$ is calculated as follows in the nth iteration. A polynomial of interpolation $p_{n,k} (x)$ of degree k is fitted to the k + 1 points $(x_n, f(x_n)), \dots (x_{n+k}, f(x_{n+k}))$. With this polynomial, the next approximation $x_{n+k+1}$ of $\alpha$ is calculated as

$x_{n+k+1} = x_{n+k} - \frac{f(x_{n+k})}{p_{n,k}'(x_{n+k})}$

(1)

with $p_{n,k}'(x_{n+k})$ the derivative of $p_{n,k}$ at $x_{n+k}$. Having calculated $x_{n+k+1}$ one calculates $f(x_{n+k+1})$ and the algorithm can continue with the (n + 1)th iteration. Clearly, this method requires the function $f$ to be evaluated only once per iteration; it requires no derivatives of $f$.

The iterative cycle is stopped if an appropriate stop-criterion is met. Typically the criterion is that the last calculated approximation is close enough to the sought-after root $\alpha$.

To execute the algorithm effectively, Sidi's method calculates the interpolating polynomial $p_{n,k} (x)$ in its Newton form.

Convergence

Sidi showed that if the function $f$ is (k + 1)-times continuously differentiable in an open interval $I$ containing $\alpha$ (that is, $f \in C^{k+1} (I)$), $\alpha$ is a simple root of $f$ (that is, $f'(\alpha) \neq 0$) and the initial approximations $x_1 , \dots , x_{k+1}$ are chosen close enough to $\alpha$, then the sequence $\{ x_i \}$ converges to $\alpha$, meaning that the following limit holds: $\lim\limits_{n \to \infty} x_n = \alpha$.

Sidi furthermore showed that

$\lim_{n\to\infty} \frac{x_{n +1}-\alpha}{\prod^k_{i=0}(x_{n-i}-\alpha)} = L = \frac{(-1)^{k+1}} {(k+1)!}\frac{f^{(k+1)}(\alpha)}{f'(\alpha)},$

and that the sequence converges to $\alpha$ of order $\psi_k$, i.e.

$\lim\limits_{n \to \infty} \frac{|x_{n+1}-\alpha|}{|x_n-\alpha|^{\psi_k}} = |L|^{(\psi_k-1)/k}$

The order of convergence $\psi_k$ is the only positive root of the polynomial

$s^{k+1} - s^k - s^{k-1} - \dots - s - 1$

We have e.g. $\psi_1 = (1+\sqrt{5})/2$ ≈ 1.6180, $\psi_2$ ≈ 1.8393 and $\psi_3$ ≈ 1.9276. The order approaches 2 from below if k becomes large: $\lim\limits_{k \to \infty} \psi_k = 2$ ^{[2] [3]}

Related algorithms

Sidi's method reduces to the secant method if we take k = 1. In this case the polynomial $p_{n,1} (x)$ is the linear approximation of $f$ around $\alpha$ which is used in the nth iteration of the secant method.

We can expect that the larger we choose k, the better $p_{n,k} (x)$ is an approximation of $f(x)$ around $x=\alpha$. Also, the better $p_{n,k}' (x)$ is an approximation of $f'(x)$ around $x=\alpha$. If we replace $p_{n,k}'$ with $f'$ in (1) we obtain that the next approximation in each iteration is calculated as

$x_{n+k+1} = x_{n+k} - \frac{f(x_{n+k})}{f'(x_{n+k})}$

(2)

This is the Newton–Raphson method. It starts off with a single approximation $x_1$ so we can take k = 0 in (2). It does not require an interpolating polynomial but instead one has to evaluate the derivative $f'$ in each iteration. Depending on the nature of $f$ this may not be possible or practical.

Once the interpolating polynomial $p_{n,k} (x)$ has been calculated, one can also calculate the next approximation $x_{n+k+1}$ as a solution of $p_{n,k} (x)=0$ instead of using (1). For k = 1 these two methods are identical: it is the secant method. For k = 2 this method is known as Muller's method.^[3] For k = 3 this approach involves finding the roots of a cubic function, which is unattractively complicated. This problem becomes worse for even larger values of k. An additional complication is that the equation $p_{n,k} (x)=0$ will in general have multiple solutions and a prescription has to be given which of these solutions is the next approximation $x_{n+k+1}$. Muller does this for the case k = 2 but no such prescriptions appear to exist for k > 2.

References

1. Sidi, Avram, "Generalization Of The Secant Method For Nonlinear Equations", Applied Mathematics E-notes 8 (2008), 115–123, http://www.math.nthu.edu.tw/~amen/2008/070227-1.pdf
2. Traub, J.F., "Iterative Methods for the Solution of Equations", Prentice Hall, Englewood Cliffs, N.J. (1964)
3. Muller, David E., "A Method for Solving Algebraic Equations Using an Automatic Computer", Mathematical Tables and Other Aids to Computation 10 (1956), 208–215