View source for Halley's method

{{Short description|A method of numerically finding roots of a function}}
{{Use dmy dates|date=October 2023}}
In [[numerical analysis]], '''Halley's method''' is a [[root-finding algorithm]] used for [[function (mathematics)|functions]] of one real variable with a continuous [[second derivative]]. [[Edmond Halley]] was an English mathematician and astronomer who introduced the method now called by his name. 

The algorithm is second in the class of [[Householder's method]]s, after [[Newton's method]]. Like the latter, it iteratively produces a sequence of approximations to the root; their [[rate of convergence]] to the root is cubic. Multidimensional versions of this method exist.<ref>{{cite journal |last1=Gundersen |first1=Geir |last2=Steihaug |first2=Trond |title=On large scale unconstrained optimization problems and higher order methods |journal=Optimization Methods & Software |date=2010 |volume=25 |issue=3 |pages=337–358 |doi=10.1080/10556780903239071 |url=https://optimization-online.org/wp-content/uploads/2007/03/1610.pdf |access-date=2 September 2024}}</ref>

Halley's method exactly finds the roots of a linear-over-linear [[Padé approximant|Padé approximation]] to the function, in contrast to [[Newton's method]] or the [[Secant method]] which approximate the function linearly, or [[Muller's method]] which approximates the function quadratically.<ref>{{Cite journal | last1 = Boyd | first1 = John P. | title = Finding the Zeros of a Univariate Equation: Proxy Rootfinders, Chebyshev Interpolation, and the Companion Matrix | doi = 10.1137/110838297 | journal = SIAM Review | volume = 55 | issue = 2 | pages = 375–396 | year = 2013 | url = http://epubs.siam.org/doi/pdf/10.1137/110838297 }}</ref>

==Method==
Halley's method is a numerical algorithm for solving the nonlinear equation {{math|1=''f''(''x'') = 0}}. In this case, the function ''f'' has to be a function of one real variable. The method consists of a sequence of iterations:

:<math>x_{n+1} = x_n - \frac {2 f(x_n) f'(x_n)} {2 {[f'(x_n)]}^2 - f(x_n) f''(x_n)} </math>

beginning with an initial guess {{math|''x''<sub>0</sub>}}.<ref>{{Cite journal|last1=Scavo|first1=T. R.|last2=Thoo|first2=J. B.|year=1995|title=On the geometry of Halley's method|journal=American Mathematical Monthly|volume=102|issue=5|pages=417–426|doi=10.2307/2975033|jstor=2975033}}</ref>

If ''f'' is a three times continuously differentiable function and ''a'' is a zero of ''f'' but not of its derivative, then, in a neighborhood of ''a'', the iterates {{math|''x''<sub>''n''</sub>}} satisfy:

:<math>| x_{n+1} - a | \le K \cdot {| x_n - a |}^3,\text{ for some }K > 0.</math>

This means that the iterates converge to the zero if the initial guess is sufficiently close, and that the convergence is cubic.<ref>{{Cite journal|last=Alefeld|first=G.|year=1981|title=On the convergence of Halley's method|jstor=2321760|journal=American Mathematical Monthly|volume=88|issue=7|pages=530–536|doi=10.2307/2321760}}</ref>

The following alternative formulation shows the similarity between Halley's method and Newton's method. The expression <math>f(x_n)/f'(x_n)</math> is computed only once, and it is particularly useful when <math>f''(x_n)/f'(x_n)</math> can be simplified:

:<math>x_{n+1} = x_n -\frac{f(x_n)}{f'(x_n) - \frac{f(x_n)}{f'(x_n)}\frac{f''(x_n)}{2}}= x_n - \frac {f(x_n)} {f'(x_n)} \left[1 - \frac {f(x_n)}{f'(x_n)} \cdot \frac {f''(x_n)} {2 f'(x_n)} \right]^{-1}.</math>

When the [[second derivative]] is very close to zero, the Halley's method iteration is almost the same as the Newton's method iteration.

==Derivation==

Consider the function

:<math>g(x) = \frac {f(x)} {\sqrt{|f'(x)|}}.</math>

Any root ''r'' of ''f'' that is ''not'' a root of its derivative is a root of ''g'' (i.e., <math>g(r)=0</math> when <math>f(r)=0 \ne {\sqrt{|f'(r)|}}</math>), and any root ''r'' of ''g'' must be a root of ''f'' provided the derivative of ''f'' at ''r'' is not infinite. Applying [[Newton's method]] to ''g'' gives

:<math>x_{n+1} = x_n - \frac{g(x_n)}{g'(x_n)}</math>

with 

:<math>g'(x) = \frac{2[f'(x)]^2 - f(x) f''(x)}{2 f'(x) \sqrt{|f'(x)|}}, </math>

and the result follows. Notice that if {{math|1=''f'' ′(''c'') = 0}}, then one cannot apply this at ''c'' because {{math|''g''(''c'')}} would be undefined. {{Explain|date=December 2023}}

==Cubic convergence==
Suppose ''a'' is a root of ''f'' but not of its derivative. And suppose that the third derivative of ''f'' exists and is continuous in a neighborhood of ''a'' and {{math|''x''<sub>''n''</sub>}} is in that neighborhood. Then [[Taylor's theorem]] implies:

:<math>0 = f(a) = f(x_n) + f'(x_n) (a - x_n) + \frac{f''(x_n)}{2} (a - x_n)^2 + \frac{f'''(\xi)}{6} (a - x_n)^3</math>

and also

:<math>0 = f(a) = f(x_n) + f'(x_n) (a - x_n) + \frac{f''(\eta)}{2} (a - x_n)^2,</math>

where ξ and η are numbers lying between ''a'' and {{math|''x''<sub>''n''</sub>}}. Multiply the first equation by <math>2f'(x_n)</math> and subtract from it the second equation times <math>f''(x_n)(a - x_n)</math> to give:

:<math>\begin{align}
0 &= 2 f(x_n) f'(x_n) + 2 [f'(x_n)]^2 (a - x_n) + f'(x_n) f''(x_n) (a - x_n)^2 + \frac{f'(x_n) f'''(\xi)}{3} (a - x_n)^3 \\
&\qquad- f(x_n) f''(x_n) (a - x_n) - f'(x_n) f''(x_n) (a - x_n)^2 - \frac{f''(x_n) f''(\eta)}{2} (a - x_n)^3.
\end{align}</math>

Canceling <math>f'(x_n) f''(x_n) (a - x_n)^2</math> and re-organizing terms yields:

:<math>0 = 2 f(x_n) f'(x_n) + \left(2 [f'(x_n)]^2 - f(x_n) f''(x_n) \right) (a - x_n) + \left(\frac{f'(x_n) f'''(\xi)}{3} - \frac{f''(x_n) f''(\eta)}{2} \right) (a - x_n)^3.</math>

Put the second term on the left side and divide through by

:<math> 2 [f'(x_n)]^2 - f(x_n) f''(x_n) </math>

to get:

:<math>a - x_n = \frac{-2f(x_n) f'(x_n)}{2[f'(x_n)]^2 - f(x_n) f''(x_n)} - \frac{2f'(x_n) f'''(\xi) - 3 f''(x_n) f''(\eta)}{6(2 [f'(x_n)]^2 - f(x_n) f''(x_n))} (a - x_n)^3.</math>

Thus:

:<math>a - x_{n+1} = - \frac{2 f'(x_n) f'''(\xi) - 3 f''(x_n) f''(\eta)}{12[f'(x_n)]^2 - 6 f(x_n) f''(x_n)} (a - x_n)^3.</math>

The limit of the coefficient on the right side as {{math|''x''<sub>''n''</sub> → ''a''}} is:

:<math>-\frac{2 f'(a) f'''(a) - 3 f''(a) f''(a)}{12 [f'(a)]^2 - 6 f(a) f''(a)}.</math>

If we take ''K'' to be a little larger than the absolute value of this, we can take absolute values of both sides of the formula and replace the absolute value of coefficient by its upper bound near ''a'' to get:

:<math>|a - x_{n+1}| \leq K |a - x_n|^3</math>

which is what was to be proved.

To summarize,

:<math>\Delta x_{i+1} =\frac{3(f'')^2 - 2f' f'''}{12(f')^2} (\Delta x_i)^3 + O[\Delta x_i]^4, \qquad \Delta x_i \triangleq x_i - a.</math><ref>{{Cite journal|last1=Proinov|first1=Petko D.|last2=Ivanov|first2=Stoil I.|year=2015|title=On the convergence of Halley's method for simultaneous computation of polynomial zeros|journal=J. Numer. Math.|volume=23|issue=4|pages=379–394|doi=10.1515/jnma-2015-0026|s2cid=10356202 }}</ref>

==Halley's irrational method==
{{Expand section|date=August 2024}}
Halley actually developed ''two'' third-order root-finding methods.  The above, using only a division, is referred to as ''Halley's rational method''.  A second, "irrational" method uses a square root as well:<ref>{{cite journal
 |first=Harry |last=Bateman |author-link=Harry Bateman
 |title=Halley's methods for solving equations
 |journal=The American Mathematical Monthly
 |volume=45 |issue=1 |date=January 1938 |pages=11–17
 |doi=10.2307/2303467
|jstor=2303467 }}</ref><ref name=Halley1694>{{cite journal
 |title=Methodus nova accurata & facilis inveniendi radices æqnationum quarumcumque generaliter, sine praviæ reductione
 |first=Edmond |last=Halley |author-link=Edmond Halley
 |journal=[[Philosophical Transactions of the Royal Society]]
 |volume=18 |issue=210 |date=May 1694 |pages=136–148
 |lang=la
 |doi=10.1098/rstl.1694.0029 |doi-access=free
}}  An English translation was published as {{cite book
 |chapter=A new, exact, and easy Method of finding the Roots of any Equations generally, and that without any previous Reduction
 |first=Edmond |last=Halley |author-link=Edmond Halley
 |title=The Philosophical Transactions of the Royal Society of London, from their commencement, in 1665, to the year 1800
 |editor1=C. Hutton |editor2=G. Shaw |editor3=R. Pearson
 |volume=III from 1683 to 1694
 |year=1809 |orig-date=May 1694 |pages=640–649
 |url=https://archive.org/details/s3id13654280
}}</ref><ref name=eggrobin>{{cite tech report
 |title=A correctly rounded binary64 cube root
 |first=Robin |last=Leroy
 |date=2021-06-21
 |url=https://raw.githubusercontent.com/mockingbirdnest/Principia/master/documentation/cbrt.pdf?raw=true
}}</ref>
:<math>x_{n+1} = x_n - \frac{f'(x_n) - \sqrt{[f'(x_n)]^2 - 2f(x_n)f''(x_n)}}{f''(x_n)}</math>
This iteration was "deservedly preferred" to the rational method by Halley{{r|Halley1694}} on the grounds that the denominator is smaller, making the division easier.  A second advantage is that it tends to have about half of the error of the rational method, a benefit which multiplies as it is iterated.  On a computer, it would appear to be slower as it has two slow operations (division and square root) instead of one, but on modern computers the reciprocal of the denominator can be computed at the same time as the square root via [[instruction pipelining]], so the [[Latency (engineering)#Computer hardware and software systems|latency]] of each iteration differs very little.{{r|eggrobin|p=24}}

==References==
<references/>

==External links==
* {{MathWorld|urlname=HalleysMethod|title=Halley’s method}}
* ''[http://numbers.computation.free.fr/Constants/Algorithms/newton.html Newton's method and high order iterations]'', Pascal Sebah and Xavier Gourdon, 2001 (the site has a link to a Postscript version for better formula display)

{{Root-finding algorithms}}
{{DEFAULTSORT:Halley's Method}}
[[Category:Root-finding algorithms]]