El método de Brent

Algoritmo de búsqueda de raíces

En el análisis numérico , el método de Brent es un algoritmo híbrido de búsqueda de raíces que combina el método de bisección , el método de la secante y la interpolación cuadrática inversa . Tiene la fiabilidad de la bisección, pero puede ser tan rápido como algunos de los métodos menos fiables. El algoritmo intenta utilizar el método de la secante, que puede converger rápidamente, o la interpolación cuadrática inversa si es posible, pero recurre al método de bisección, que es más robusto, si es necesario. El método de Brent se debe a Richard Brent ^[1] y se basa en un algoritmo anterior de Theodorus Dekker ^[2] . En consecuencia, el método también se conoce como el método de Brent-Dekker .

Las mejoras modernas del método de Brent incluyen el método de Chandrupatla, que es más simple y rápido para funciones que son planas alrededor de sus raíces; ^{[3] [4]} el método de Ridders , que realiza interpolaciones exponenciales en lugar de cuadráticas, proporcionando una fórmula cerrada más simple para las iteraciones; y el método ITP , que es un híbrido entre regula-falsi y bisección que logra garantías asintóticas y en el peor de los casos óptimos.

El método de Dekker

La idea de combinar el método de bisección con el método secante se remonta a Dekker (1969).

Supongamos que queremos resolver la ecuación f ( x ) = 0. Al igual que con el método de bisección, necesitamos inicializar el método de Dekker con dos puntos, digamos a ₀ y b ₀ , tales que f ( a ₀ ) y f ( b ₀ ) tengan signos opuestos. Si f es continua en [ a ₀ , b ₀ ], el teorema del valor intermedio garantiza la existencia de una solución entre a ₀ y b ₀ .

En cada iteración intervienen tres puntos:

• b _k es la iteración actual, es decir, la estimación actual de la raíz de f .
• a _k es el "contrapunto", es decir, un punto tal que f ( a _k ) y f ( b _k ) tienen signos opuestos, por lo que el intervalo [ a _k , b _k ] contiene la solución. Además, | f ( b _k )| debe ser menor o igual que | f ( a _k )|, de modo que b _k es una mejor estimación de la solución desconocida que a _k .
• b _{k −1} es la iteración anterior (para la primera iteración, establecemos b _{k −1} = a ₀ ).

Se calculan dos valores provisionales para la siguiente iteración. El primero se obtiene mediante interpolación lineal, también conocida como método de la secante:

${\displaystyle s={\begin{cases}b_{k}-{\frac {b_{k}-b_{k-1}}{f(b_{k})-f(b_{k-1})}}f(b_{k}),&{\mbox{if }}f(b_{k})\neq f(b_{k-1})\\m&{\mbox{otherwise }}\end{cases}}}$

y el segundo se da por el método de bisección

${\displaystyle m={\frac {a_{k}+b_{k}}{2}}.}$

Si el resultado del método secante, s , se encuentra estrictamente entre b _k y m , entonces se convierte en la siguiente iteración ( b _{k +1} = s ); de lo contrario, se utiliza el punto medio ( b _{k +1} = m ).

Luego, el valor del nuevo contrapunto se elige de manera que f ( a _{k +1} ) y f ( b _{k +1} ) tengan signos opuestos. Si f ( a _k ) y f ( b _{k +1} ) tienen signos opuestos, entonces el contrapunto permanece igual: a _{k +1} = a _k . De lo contrario, f ( b _{k +1} ) y f ( b _k ) tienen signos opuestos, por lo que el nuevo contrapunto se convierte en a _{k +1} = b _k .

Finalmente, si | f ( a _{k +1} )| < | f ( b _{k +1} )|, entonces a _{k +1} es probablemente una mejor estimación de la solución que b _{k +1} y, por lo tanto, los valores de a _{k +1} y b _{k +1} se intercambian.

Con esto finaliza la descripción de una única iteración del método de Dekker.

El método de Dekker funciona bien si la función f se comporta razonablemente bien. Sin embargo, hay circunstancias en las que cada iteración emplea el método de la secante, pero las iteraciones b _k convergen muy lentamente (en particular, | b _k − b _{k −1} | puede ser arbitrariamente pequeño). El método de Dekker requiere muchas más iteraciones que el método de bisección en este caso.

El método de Brent

Brent (1973) propuso una pequeña modificación para evitar el problema del método de Dekker. Inserta una prueba adicional que debe cumplirse antes de que el resultado del método de la secante sea aceptado como la siguiente iteración. Se deben cumplir simultáneamente dos desigualdades:

Dada una tolerancia numérica específica , si en el paso anterior se utilizó el método de bisección, la desigualdad debe cumplirse para realizar la interpolación; de lo contrario, se realiza el método de bisección y su resultado se utiliza para la siguiente iteración. ${\displaystyle \delta }$ ${\textstyle |\delta |<|b_{k}-b_{k-1}|}$

Si en el paso anterior se realizó una interpolación, entonces se utiliza la desigualdad para realizar la siguiente acción (a elegir): interpolación (cuando la desigualdad es verdadera) o el método de bisección (cuando la desigualdad no es verdadera). ${\textstyle |\delta |<|b_{k-1}-b_{k-2}|}$

Además, si en el paso anterior se utilizó el método de bisección, la desigualdad debe cumplirse; de ​​lo contrario, se aplica el método de bisección y su resultado se utiliza para la siguiente iteración. Si en el paso anterior se utilizó la interpolación, se utiliza la desigualdad en su lugar. ${\textstyle |s-b_{k}|<{\begin{matrix}{\frac {1}{2}}\end{matrix}}|b_{k}-b_{k-1}|}$ ${\textstyle |s-b_{k}|<{\begin{matrix}{\frac {1}{2}}\end{matrix}}|b_{k-1}-b_{k-2}|}$

Esta modificación garantiza que en la k-ésima iteración, se realizará un paso de bisección en, como máximo, iteraciones adicionales, porque las condiciones anteriores obligan a que los tamaños de los pasos de interpolación consecutivos se reduzcan a la mitad cada dos iteraciones y, después de, como máximo, iteraciones, el tamaño del paso será menor que , lo que invoca un paso de bisección. Brent demostró que su método requiere, como máximo, N ² iteraciones, donde N denota el número de iteraciones para el método de bisección. Si la función f se comporta bien, entonces el método de Brent generalmente procederá por interpolación cuadrática inversa o lineal, en cuyo caso convergerá de manera superlineal . ${\displaystyle 2\log _{2}(|b_{k-1}-b_{k-2}|/\delta )}$ ${\displaystyle 2\log _{2}(|b_{k-1}-b_{k-2}|/\delta )}$ ${\displaystyle \delta }$

Además, el método de Brent utiliza interpolación cuadrática inversa en lugar de interpolación lineal (como la que utiliza el método secante). Si f ( b _k ), f ( a _k ) y f ( b _{k −1} ) son distintas, aumenta ligeramente la eficiencia. Como consecuencia, la condición para aceptar s (el valor propuesto por la interpolación lineal o la interpolación cuadrática inversa) tiene que cambiarse: s tiene que estar entre (3 a _k + b _k ) / 4 y b _k .

Algoritmo

entrada  a , b y (un puntero a) una función para f
calcular f ( a )Calcular f ( b ) si  f ( a ) f ( b ) ≥ 0 entonces  función de salida porque la raíz no está entre corchetes.fin si si | f ( a )| < | f ( b )| entonces intercambia ( a , b ) fin si c  := a establecer mflag repetir hasta  f ( b o s ) = 0 o | b − a | es suficientemente pequeño (convergencia)  si  f ( a ) ≠ f ( c ) y  f ( b ) ≠ f ( c ) entonces  ( interpolación cuadrática inversa ) de lo contrario ( método secante ) fin si si (condición 1) s no está entre y b o (condición 2) (mflag está establecido y | s − b | ≥ | b − c |/2) o (condición 3) (mflag está borrado y | s − b | ≥ | c − d |/2) o (condición 4) (mflag está establecido y | b − c | < | δ |) o (condición 5) (mflag está borrado y | c − d | < | δ |) entonces ( método de bisección ) establecer mflag de lo contrario borrar mflag fin si calcular f ( s ) d  := c (d se asigna por primera vez aquí; no se usará arriba en la primera iteración porque mflag está configurado) c  := b if f ( a ) f ( s ) < 0 then b  := s else a  := s end if if | f ( a )| < | ${\textstyle s:={\frac {af(b)f(c)}{(f(a)-f(b))(f(a)-f(c))}}+{\frac {bf(a)f(c)}{(f(b)-f(a))(f(b)-f(c))}}+{\frac {cf(a)f(b)}{(f(c)-f(a))(f(c)-f(b))}}}$    ${\textstyle s:=b-f(b){\frac {b-a}{f(b)-f(a)}}}$        ${\displaystyle (3a+b)/4}$               ${\textstyle s:={\frac {a+b}{2}}}$            f ( b )| entonces intercambia ( a , b ) fin si fin repetir salida  b  o s (devuelve la raíz)

Ejemplo

Supongamos que estamos buscando un cero de la función definida por f ( x ) = ( x + 3)( x − 1) ² .

Tomamos [ a ₀ , b ₀ ] = [−4, 4/3] como nuestro intervalo inicial.

Tenemos f ( a ₀ ) = −25 y f ( b ₀ ) = 0,48148 (todos los números en esta sección están redondeados), por lo que se cumplen las condiciones f ( a ₀ ) f ( b ₀ ) < 0 y | f ( b ₀ )| ≤ | f ( a ₀ )|.

Gráfica de f ( x ) = ( x + 3)( x − 1) ²

1. En la primera iteración, utilizamos interpolación lineal entre ( b ₋₁ , f ( b ₋₁ )) = ( a ₀ , f ( a ₀ )) = (−4, −25) y ( b ₀ , f ( b ₀ )) = (1.33333, 0.48148), lo que produce s = 1.23256. Esto se encuentra entre (3 a ₀ + b ₀ ) / 4 y b ₀ , por lo que se acepta este valor. Además, f (1.23256) = 0.22891, por lo que establecemos a ₁ = a ₀ y b ₁ = s = 1.23256.
2. En la segunda iteración, utilizamos la interpolación cuadrática inversa entre ( a ₁ , f ( a ₁ )) = (−4, −25) y ( b ₀ , f ( b ₀ )) = (1.33333, 0.48148) y ( b ₁ , f ( b ₁ )) = (1.23256, 0.22891). Esto produce 1.14205, que se encuentra entre (3 a ₁ + b ₁ ) / 4 y b ₁ . Además, se satisface la desigualdad |1.14205 − b ₁ | ≤ | b ₀ − b ₋₁ | / 2, por lo que se acepta este valor. Además, f (1,14205) = 0,083582, por lo que establecemos a ₂ = a ₁ y b ₂ = 1,14205.
3. En la tercera iteración, utilizamos la interpolación cuadrática inversa entre ( a ₂ , f ( a ₂ )) = (−4, −25) y ( b ₁ , f ( b ₁ )) = (1,23256, 0,22891) y ( b ₂ , f ( b ₂ )) = (1,14205, 0,083582). Esto produce 1,09032, que se encuentra entre (3 a ₂ + b ₂ ) / 4 y b ₂ . Pero aquí entra en juego la condición adicional de Brent: la desigualdad |1,09032 − b ₂ | ≤ | b ₁ − b ₀ | / 2 no se satisface, por lo que este valor se rechaza. En su lugar, se calcula el punto medio m = −1,42897 del intervalo [ a ₂ , b ₂ ]. Tenemos f ( m ) = 9,26891, por lo que establecemos a ₃ = a ₂ y b ₃ = −1,42897.
4. En la cuarta iteración, utilizamos la interpolación cuadrática inversa entre ( a ₃ , f ( a ₃ )) = (−4, −25) y ( b ₂ , f ( b ₂ )) = (1,14205, 0,083582) y ( b ₃ , f ( b ₃ )) = (−1,42897, 9,26891). Esto produce 1,15448, que no está en el intervalo entre (3 a ₃ + b ₃ ) / 4 y b ₃ ). Por lo tanto, se reemplaza por el punto medio m = −2,71449. Tenemos f ( m ) = 3,93934, por lo que establecemos a ₄ = a ₃ y b ₄ = −2,71449.
5. En la quinta iteración, la interpolación cuadrática inversa produce −3,45500, que se encuentra en el intervalo requerido. Sin embargo, la iteración anterior fue un paso de bisección, por lo que se debe satisfacer la desigualdad |−3,45500 − b ₄ | ≤ | b ₄ − b ₃ | / 2 . Esta desigualdad es falsa, por lo que utilizamos el punto medio m = −3,35724. Tenemos f ( m ) = −6,78239, por lo que m se convierte en el nuevo contrapunto ( a ₅ = −3,35724) y la iteración permanece igual ( b ₅ = b ₄ ).
6. En la sexta iteración, no podemos usar la interpolación cuadrática inversa porque b ₅ = b ₄ . Por lo tanto, usamos la interpolación lineal entre ( a ₅ , f ( a ₅ )) = (−3.35724, −6.78239) y ( b ₅ , f ( b ₅ )) = (−2.71449, 3.93934). El resultado es s = −2.95064, que satisface todas las condiciones. Pero como la iteración no cambió en el paso anterior, rechazamos este resultado y volvemos a la bisección. Actualizamos s = -3.03587 y f ( s ) = -0.58418.
7. En la séptima iteración, podemos volver a utilizar la interpolación cuadrática inversa. El resultado es s = −3,00219, que satisface todas las condiciones. Ahora, f ( s ) = −0,03515, por lo que establecemos a ₇ = b ₆ y b ₇ = −3,00219 ( a ₇ y b ₇ se intercambian de modo que se satisface la condición | f ( b ₇ )| ≤ | f ( a ₇ )|). ( Correcto  : interpolación lineal ⁠ ⁠ ${\displaystyle s=-2.99436,f(s)=0.089961}$ )
8. En la octava iteración, no podemos utilizar la interpolación cuadrática inversa porque a ₇ = b ₆ . La interpolación lineal da como resultado s = −2,99994, que es aceptado. ( Correcto  : ⁠ ⁠ ${\displaystyle s=-2.9999,f(s)=0.0016}$ )
9. En las iteraciones siguientes, la raíz x = −3 se aproxima rápidamente: b ₉ = −3 + 6·10 ⁻⁸ y b ₁₀ = −3 − 3·10 ⁻¹⁵ . ( Correcto  : Iter 9 : f ( s ) = −1,4 × 10 ⁻⁷ , Iter 10 : f ( s ) = 6,96 × 10 ⁻¹² )

Implementaciones

• Brent (1973) publicó una implementación de Algol 60 .
• Netlib contiene una traducción Fortran de esta implementación con ligeras modificaciones.
• El método PARI/GPsolve implementa el método.
• Otras implementaciones del algoritmo (en C++, C y Fortran) se pueden encontrar en los libros de Recetas numéricas .
• La biblioteca Apache Commons Math implementa el algoritmo en Java .
• El módulo de optimización de SciPy implementa el algoritmo en Python (lenguaje de programación)
• La biblioteca estándar de Modelica implementa el algoritmo en Modelica .
• La unirootfunción implementa el algoritmo en R (software) .
• La fzerofunción implementa el algoritmo en MATLAB .
• Boost (bibliotecas de C++) implementa dos algoritmos basados ​​en el método de Brent en C++ en el kit de herramientas Math:
  1. Minimización de funciones en minima.hpp con un ejemplo que localiza los mínimos de la función.
  2. La búsqueda de raíces implementa el nuevo TOMS748, un algoritmo más moderno y eficiente que el original de Brent, en TOMS748, y la búsqueda de raíces Boost.Math que utiliza TOMS748 internamente con ejemplos.
• El paquete Optim.jl implementa el algoritmo en Julia (lenguaje de programación)
• El sistema de álgebra computacional Emmy (escrito en Clojure (lenguaje de programación) ) implementa una variante del algoritmo diseñado para la minimización de funciones univariadas.
• Búsqueda de raíces en la biblioteca de C# alojada en Code Project.

Referencias

1. Brent 1973
2. Dekker 1969
3. Chandrupatla, Tirupathi R. (1997). "Un nuevo algoritmo híbrido cuadrático/de bisección para hallar el cero de una función no lineal sin utilizar derivadas". Advances in Engineering Software . 28 (3): 145–149. doi :10.1016/S0965-9978(96)00051-8.
4. "Diez pequeños algoritmos, parte 5: interpolación de extremos cuadráticos y método de Chandrupatla - Jason Sachs".

• Brent, RP (1973), "Capítulo 4: Un algoritmo con convergencia garantizada para encontrar un cero de una función", Algoritmos para la minimización sin derivadas , Englewood Cliffs, NJ: Prentice-Hall, ISBN 0-13-022335-2
• Dekker, TJ (1969), "Encontrar un cero mediante interpolación lineal sucesiva", en Dejon, B.; Henrici, P. (eds.), Aspectos constructivos del teorema fundamental del álgebra , Londres: Wiley-Interscience, ISBN 978-0-471-20300-1

Lectura adicional

• Atkinson, Kendall E. (1989). "Sección 2.8". Introducción al análisis numérico (2.ª ed.). John Wiley and Sons. ISBN 0-471-50023-2.
• Press, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP (2007). "Sección 9.3. Método de Van Wijngaarden–Dekker–Brent". Recetas numéricas: el arte de la computación científica (3.ª ed.). Nueva York: Cambridge University Press. ISBN 978-0-521-88068-8Archivado desde el original el 11-08-2011 . Consultado el 28-02-2012 .
• Alefeld, GE; Potra, FA; Shi, Yixun (septiembre de 1995). "Algoritmo 748: encerrar ceros de funciones continuas". Transacciones ACM sobre software matemático . 21 (3): 327–344. doi : 10.1145/210089.210111 . S2CID  207192624.

Enlaces externos

• zeroin.f en Netlib .
• módulo brent en C++ (también C, Fortran, Matlab) Archivado el 5 de abril de 2018 en Wayback Machine por John Burkardt
• Implementación de GSL.
• Impulsar la implementación de C++.
• Implementación de Python (Scipy)