El método de Brent.

Algoritmo de búsqueda de raíces

En 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 confiabilidad de la bisección pero puede ser tan rápido como algunos de los métodos menos confiables. El algoritmo intenta utilizar el método secante potencialmente de rápida convergencia o la interpolación cuadrática inversa si es posible, pero recurre al método de bisección 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 método de Brent-Dekker .

Las mejoras modernas al 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]} Método de Ridders , que realiza interpolaciones exponenciales en lugar de cuadráticas y proporciona 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 óptimas en el peor de los casos y asintóticas.

El método de Dekker.

La idea de combinar el método de la bisección con el método de la 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 ₀ ) tienen 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 para 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 ( segundo _k )| debe ser menor o igual que | f ( a _k )|, de modo que b _k es una mejor estimación para 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 viene dado por interpolación lineal, también conocido 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 está dado por el método de la 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 sigue siendo el mismo: 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 ( ak _{+1 )} | < | f ( b _{k +1} )|, entonces a _{k +1} es probablemente una mejor suposición para la solución que b _{k +1} y, por lo tanto, se intercambian los valores de a _{k +1} y b _{k +1 .}

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 secante, pero las iteraciones b _k convergen muy lentamente (en particular, | b _k − b _{k −1} | puede ser arbitrariamente pequeño). En este caso, el método de Dekker requiere muchas más iteraciones que el método de bisección.

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 se acepte el resultado del método secante como la siguiente iteración. Deben satisfacerse 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 el paso anterior realizó la interpolación, entonces la desigualdad se usa para realizar la siguiente acción (para elegir) la 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 realiza el método de bisección y su resultado se utiliza para la siguiente iteración. Si el paso anterior realizó una interpolación, entonces se utiliza la desigualdad. ${\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 paso de interpolación consecutivos se reduzcan a la mitad cada dos iteraciones, y después de como máximo las iteraciones, el tamaño del paso será menor que , 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á mediante interpolación cuadrática inversa o lineal, en cuyo caso convergerá de forma 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 de la secante). Si f ( b _k ), f ( a _k ) y f ( b _{k −1} ) son distintos, aumenta ligeramente la eficiencia. Como consecuencia, la condición para aceptar s (el valor propuesto por interpolación lineal o interpolación cuadrática inversa) debe cambiarse: s debe estar entre (3 a _k + b _k ) / 4 y b _k .

Algoritmo

ingrese  a , b y (un puntero a) una función para f
calcule f ( a )calcular f ( b ) si  f ( a ) f ( b ) ≥ 0 entonces  Salir de la función porque la raíz no está entre corchetes.terminar si si | f ( a )| < | f ( segundo )| luego intercambie ( a , b ) finalice si c  : = a set mflag repita hasta  f ( b o s ) = 0 o | segundo - un | es lo 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 ) termina 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á desactivado 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 ) establece mflag ; de lo contrario, borra mflag end 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 si f ( a ) f ( s ) < 0 entonces b  := s else a  := s final si si | 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 ( segundo )| luego intercambie ( a , b ) final si final repite la salida  b  o s (devuelve la raíz)

Ejemplo

Supongamos que buscamos 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 las condiciones f ( a ₀ ) f ( b ₀ ) < 0 y | f ( segundo ₀ )| ≤ | f ( un ₀ )| estan satisfechos.

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

1. En la primera iteración, usamos 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. Este 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, usamos interpolación cuadrática inversa entre ( a ₁ , f ( a ₁ )) = (−4, −25) y ( b ₀ , f ( b ₀ )) = (1.33333, 0.48148) y ( b ₁ , f ( segundo ₁ )) = (1,23256, 0,22891). Esto produce 1,14205, que se encuentra entre (3 a ₁ + b ₁ ) / 4 y b ₁ . Además, la desigualdad |1.14205 − b ₁ | ≤ | segundo ₀ - segundo _-1 | /2 se cumple, 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, usamos interpolación cuadrática inversa entre ( a ₂ , f ( a ₂ )) = (−4, −25) y ( b ₁ , f ( b ₁ )) = (1.23256, 0.22891) y ( b ₂ , f ( segundo ₂ )) = (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 ₂ | ≤ | segundo ₁ - segundo ₀ | /2 no se cumple, por lo que se rechaza este valor. En cambio, 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, usamos interpolación cuadrática inversa entre ( a ₃ , f ( a ₃ )) = (−4, −25) y ( b ₂ , f ( b ₂ )) = (1.14205, 0.083582) y ( b ₃ , f ( segundo ₃ )) = (−1,42897, 9,26891). Esto produce 1,15448, que no está en el intervalo entre (3 a ₃ + b ₃ ) / 4 y b ₃ ). Por 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 la desigualdad |−3.45500 − b ₄ | ≤ | segundo ₄ - segundo ₃ | / 2 necesitan estar satisfechos. Esta desigualdad es falsa, por lo que usamos 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 sigue siendo la misma ( b ₅ = b ₄ ).
6. En la sexta iteración, no podemos usar la interpolación cuadrática inversa porque b ₅ = b ₄ . Por lo tanto, usamos 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, entonces establecemos a ₇ = b ₆ y b ₇ = −3.00219 ( a ₇ y b ₇ se intercambian de modo que la condición | f ( b ₇ )| ≤ | f ( a ₇ ) | está satisfecho). ( Correcto  : interpolación lineal ) ${\displaystyle s=-2.99436,f(s)=0.089961}$
8. En la octava iteración, no podemos usar la interpolación cuadrática inversa porque a ₇ = b ₆ . La interpolación lineal produce s = −2,99994, lo cual se acepta. ( Correcto  : ) ${\displaystyle s=-2.9999,f(s)=0.0016}$
9. En las siguientes iteraciones, 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 en Fortran de esta implementación con ligeras modificaciones.
• El método PARI/GPsolve implementa el método.
• Se pueden encontrar otras implementaciones del algoritmo (en C++, C y Fortran) en los libros Numerical Recipes .
• 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 de matemáticas:
  1. Minimización de funciones en minima.hpp con un ejemplo de localización de mínimos de funciones.
  2. La búsqueda de raíces implementa el TOMS748 más nuevo, 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 informática 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 C# alojada en Code Project.

Referencias

1. Brent 1973
2. Dekker 1969
3. Chandrupatla, Tirupathi R. (1997). "Un nuevo algoritmo híbrido cuadrático/bisección para encontrar el cero de una función no lineal sin utilizar derivadas". Avances en software de ingeniería . 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 el cero de una función", Algoritmos de minimización sin derivadas , Englewood Cliffs, Nueva Jersey: 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

Otras lecturas

• Atkinson, Kendall E. (1989). "Sección 2.8". Introducción al análisis numérico (2ª ed.). John Wiley e hijos. ISBN 0-471-50023-2.
• Prensa, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP (2007). "Sección 9.3. Método Van Wijngaarden-Dekker-Brent". Recetas numéricas: el arte de la informática científica (3ª ed.). Nueva York: Cambridge University Press. ISBN 978-0-521-88068-8. Archivado desde el original el 11 de agosto de 2011 . Consultado el 28 de febrero de 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 GSL.
• Impulsar la implementación de C++.
• Implementación de Python (Scipy)