Метод Бренца - Википедия - Brents method

В числовой анализ, Метод Брента это алгоритм поиска корней объединение метод деления пополам, то секущий метод и обратная квадратичная интерполяция. Он надежен на деление пополам, но может быть таким же быстрым, как и некоторые из менее надежных методов. Алгоритм пытается использовать потенциально быстро сходящийся метод секущей или обратную квадратичную интерполяцию, если это возможно, но при необходимости возвращается к более надежному методу деления пополам. Метод Брента обусловлен Ричард Брент^[1] и основан на более раннем алгоритме Теодорус Деккер.^[2] Следовательно, метод также известен как Метод Брента – Деккера.

Метод Чандрупатлы - это вариант, который проще и быстрее сходится для функций, которые имеют плоскую форму вокруг своих корней (что означает, что они имеют несколько корней или близко расположенные корни).^[3][4]

Метод Деккера

Идея объединить метод деления пополам с методом секущей восходит к Деккер (1969).

Предположим, что мы хотим решить уравнение ж(Икс) = 0. Как и в случае с методом деления пополам, нам нужно инициализировать метод Деккера двумя точками, скажем а₀ и б₀, так что ж(а₀) и ж(б₀) имеют противоположные знаки. Если ж непрерывно на [а₀, б₀], теорема о промежуточном значении гарантирует существование решения между а₀ и б₀.

В каждой итерации участвуют три момента:

• б_k это текущая итерация, то есть текущее предположение для корня ж.
• а_k является «контрапунктом», т. е. такой точкой, что ж(а_k) и ж(б_k) имеют противоположные знаки, поэтому интервал [а_k, б_k] содержит решение. Кроме того, |ж(б_k) | должно быть меньше или равно |ж(а_k) |, так что б_k это лучшее предположение для неизвестного решения, чем а_k.
• б_k−1 это предыдущая итерация (для первой итерации мы устанавливаем б_k−1 = а₀).

Вычисляются два предварительных значения для следующей итерации. Первый задается линейной интерполяцией, также известной как метод секущей:

${ displaystyle s = { begin {case} 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 {else}} end {случаи }}}$

а второй - методом деления пополам

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

Если результат метода секущих, s, лежит строго между б_k и м, то следующая итерация (б_k+1 = s), в противном случае используется средняя точка (б_k+1 = м).

Затем значение нового контрапункта выбирается так, чтобы ж(а_k+1) и ж(б_k+1) имеют противоположные знаки. Если ж(а_k) и ж(б_k+1) имеют противоположные знаки, то контрапункт остается прежним: а_k+1 = а_k. Иначе, ж(б_k+1) и ж(б_k) имеют противоположные знаки, поэтому новый контрапункт становится а_k+1 = б_k.

Наконец, если |ж(а_k+1)| < |ж(б_k+1) |, то а_k+1 вероятно лучшее предположение для решения, чем б_k+1, а значит, и значения а_k+1 и б_k+1 обмениваются.

На этом заканчивается описание единственной итерации метода Деккера.

Метод Деккера работает хорошо, если функция ж ведет себя достаточно хорошо. Однако есть обстоятельства, при которых каждая итерация использует метод секущей, но итерация повторяет б_k сходятся очень медленно (в частности, |б_k − б_k−1| может быть сколь угодно малым). В этом случае метод Деккера требует гораздо большего количества итераций, чем метод деления пополам.

Метод Брента

Брент (1973) предложила небольшую модификацию, чтобы избежать этой проблемы. Он вставил дополнительный тест, который должен быть удовлетворен, прежде чем результат метода секущей будет принят в качестве следующей итерации. Одновременно должны выполняться два неравенства: с учетом определенного числового допуска ${ displaystyle delta}$, если на предыдущем шаге использовался метод деления пополам, неравенство

${ displaystyle | delta | <| b_ {k} -b_ {k-1} |}$

должен удерживаться для выполнения интерполяции, в противном случае выполняется метод деления пополам и его результат используется для следующей итерации.

Если на предыдущем шаге выполнялась интерполяция, то выполняется неравенство

${ displaystyle | delta | <| b_ {k-1} -b_ {k-2} |}$

вместо этого используется для выполнения следующего действия (для выбора) интерполяции (когда неравенство истинно) или метода деления пополам (когда неравенство не верно).

Кроме того, если на предыдущем шаге использовался метод деления пополам, неравенство

${ displaystyle | s-b_ {k} | <{ begin {matrix} { frac {1} {2}} end {matrix}} | b_ {k} -b_ {k-1} |}$

должен выполняться, иначе будет выполнен метод деления пополам, и его результат будет использоваться для следующей итерации. Если на предыдущем шаге выполнялась интерполяция, то выполняется неравенство

${ displaystyle | s-b_ {k} | <{ begin {matrix} { frac {1} {2}} end {matrix}} | b_ {k-1} -b_ {k-2} |}$

вместо этого используется.

Эта модификация гарантирует, что на k-й итерации шаг деления пополам будет выполнен не более чем за ${ Displaystyle 2 журнал _ {2} (| b_ {k-1} -b_ {k-2} | / delta)}$ дополнительных итераций, потому что вышеуказанные условия вынуждают размеры последовательных шагов интерполяции уменьшаться вдвое каждые две итерации, и не более чем после ${ Displaystyle 2 журнал _ {2} (| b_ {k-1} -b_ {k-2} | / delta)}$ итераций размер шага будет меньше чем ${ displaystyle delta}$, который вызывает шаг пополам. Брент доказал, что его метод требует не более N² итераций, где N обозначает количество итераций для метода деления пополам. Если функция ж хорошо себя ведет, то метод Брента обычно использует либо обратную квадратичную, либо линейную интерполяцию, и в этом случае он сходится сверхлинейно.

Кроме того, метод Брента использует обратная квадратичная интерполяция вместо линейная интерполяция (как используется секущим методом). Если ж(б_k), ж(а_k) и ж(б_k−1) отчетливы, это немного увеличивает эффективность. Как следствие, условие принятия s (значение, предложенное либо линейной интерполяцией, либо обратной квадратичной интерполяцией) необходимо изменить: s должен находиться между (3а_k + б_k) / 4 и б_k.

Алгоритм

Эта секция возможно содержит оригинальные исследования. Пожалуйста Улучши это к проверка заявленные претензии и добавление встроенные цитаты. Заявления, содержащие только оригинальные исследования, следует удалить. (Февраль 2019 г.) (Узнайте, как и когда удалить этот шаблон сообщения)

Вход а, б, и (указатель на) функцию для жвычислить ж(а) вычислить ж(б)если ж(а)ж(б) ≥ 0 тогда     функция выхода, потому что корень не заключен в квадратные скобки.конец, еслиесли |ж(а)| < |ж(б)| тогда    замена (а,б)конец, еслиc := анабор mflagповторять, пока ж(б или же s) = 0 или же |б − а| является достаточно мал (конвергенция)    если ж(а) ≠ ж(c) и ж(б) ≠ ж(c) тогда        ${ displaystyle 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))}}}$ (обратная квадратичная интерполяция )    еще        ${ Displaystyle s: = b-f (b) { frac {b-a} {f (b) -f (a)}}}$ (секущий метод )    конец, если    если (условие 1) s не является между ${ displaystyle (3a + b) / 4}$ и б или же       (условие 2) (mflag является набор и |s−б| ≥ |б−c|/2) или же       (условие 3) (mflag является очищен и |s−б| ≥ |c−d|/2) или же       (условие 4) (mflag является набор и |б−c| < |δ|) или же       (условие 5) (mflag является очищен и |c−d| < |δ|) тогда        ${ displaystyle s: = { frac {a + b} {2}}}$ (метод деления пополам )        набор mflag еще        Чисто mflag конец, если    вычислить ж(s)    d := c  (d здесь назначается впервые; он не будет использоваться выше на первой итерации, поскольку установлен mflag)    c := б    если ж(а)ж(s) < 0 тогда        б := s     еще        а := s     конец, если    если |ж(а)| < |ж(б)| тогда        замена (а,б)     конец, есликонец повторениявыход б или s (вернуть корень)

Пример

Предположим, что мы ищем нуль функции, определенной формулой ж(Икс) = (Икс + 3)(Икс − 1)².

Мы принимаем [а₀, б₀] = [−4, 4/3] как наш начальный интервал.

У нас есть ж(а₀) = −25 и ж(б₀) = 0,48148 (все числа в этом разделе округлены), поэтому условия ж(а₀) ж(б₀) <0 и |ж(б₀)| ≤ |ж(а₀) | довольны.

График ж(Икс) = (Икс + 3)(Икс − 1)²

1. В первой итерации мы используем линейную интерполяцию между (б₋₁, ж(б₋₁)) = (а₀, ж(а₀)) = (−4, −25) и (б₀, ж(б₀)) = (1,33333, 0,48148), что дает s = 1,23256. Это находится между (3а₀ + б₀) / 4 и б₀, поэтому это значение принято. Более того, ж(1.23256) = 0.22891, поэтому мы полагаем а₁ = а₀ и б₁ = s = 1.23256.
2. Во второй итерации мы используем обратную квадратичную интерполяцию между (а₁, ж(а₁)) = (−4, −25) и (б₀, ж(б₀)) = (1,33333, 0,48148) и (б₁, ж(б₁)) = (1,23256, 0,22891). Это дает 1,14205, что находится между (3а₁ + б₁) / 4 и б₁. Кроме того, неравенство | 1.14205 - б₁| ≤ |б₀ − б₋₁| / 2 удовлетворяется, поэтому это значение принято. Более того, ж(1,14205) = 0,083582, поэтому мы полагаем а₂ = а₁ и б₂ = 1.14205.
3. В третьей итерации мы используем обратную квадратичную интерполяцию между (а₂, ж(а₂)) = (−4, −25) и (б₁, ж(б₁)) = (1,23256, 0,22891) и (б₂, ж(б₂)) = (1,14205, 0,083582). Это дает 1,09032, что находится между (3а₂ + б₂) / 4 и б₂. Но здесь вступает в действие дополнительное условие Брента: неравенство | 1.09032 - б₂| ≤ |б₁ − б₀| / 2 не выполняется, поэтому это значение отклоняется. Вместо этого середина м = −1,42897 интервала [а₂, б₂] вычисляется. У нас есть ж(м) = 9.26891, поэтому положим а₃ = а₂ и б₃ = −1.42897.
4. В четвертой итерации мы используем обратную квадратичную интерполяцию между (а₃, ж(а₃)) = (−4, −25) и (б₂, ж(б₂)) = (1,14205, 0,083582) и (б₃, ж(б₃)) = (-1,42897, 9,26891). Это дает 1,15448, что не находится в интервале между (3а₃ + б₃) / 4 и б₃). Следовательно, он заменяется средней точкой м = −2,71449. У нас есть ж(м) = 3.93934, поэтому положим а₄ = а₃ и б₄ = −2.71449.
5. На пятой итерации обратная квадратичная интерполяция дает −3,45500, что лежит в требуемом интервале. Однако предыдущая итерация была делением пополам, поэтому неравенство | −3.45500 - б₄| ≤ |б₄ − б₃| / 2 должны быть удовлетворены. Это неравенство неверно, поэтому мы используем среднюю точку м = -3,35724. У нас есть ж(м) = −6.78239, поэтому м становится новым контрапунктом (а₅ = −3,35724) и итерация остается прежней (б₅ = б₄).
6. На шестой итерации мы не можем использовать обратную квадратичную интерполяцию, потому что б₅ = б₄. Следовательно, мы используем линейную интерполяцию между (а₅, ж(а₅)) = (−3,35724, −6,78239) и (б₅, ж(б₅)) = (-2,71449, 3,93934). Результат s = −2.95064, что удовлетворяет всем условиям. Но поскольку итерация не изменилась на предыдущем шаге, мы отклоняем этот результат и возвращаемся к делению пополам. Мы обновляем s = -3,03587, и ж(s) = -0.58418.
7. В седьмой итерации мы снова можем использовать обратную квадратичную интерполяцию. Результат s = −3.00219, что удовлетворяет всем условиям. Сейчас же, ж(s) = −0,03515, поэтому положим а₇ = б₆ и б₇ = −3.00219 (а₇ и б₇ меняются местами, так что условие |ж(б₇)| ≤ |ж(а₇) | доволен). (Правильный : линейная интерполяция ${ displaystyle s = -2,99436, f (s) = 0,089961}$)
8. На восьмой итерации мы не можем использовать обратную квадратичную интерполяцию, потому что а₇ = б₆. Линейная интерполяция дает s = −2,99994, что принято. (Правильный : ${ displaystyle s = -2,9999, f (s) = 0,0016}$)
9. В следующих итерациях корень Икс = −3 приближается быстро: б₉ = −3 + 6·10⁻⁸ и б₁₀ = −3 − 3·10⁻¹⁵. (Правильный : Iter 9: ж(s) = −1.4 × 10⁻⁷, Iter 10: ж(s) = 6.96 × 10⁻¹²)

Реализации

• Брент (1973) опубликовал Алгол 60 выполнение.
• Netlib содержит перевод этой реализации на Фортран с небольшими изменениями.
• В PARI / GP метод решать реализует метод.
• Другие реализации алгоритма (на C ++, C и Fortran) можно найти в Числовые рецепты книги.
• В Apache Commons Математическая библиотека реализует алгоритм в Ява.
• В SciPy Модуль optimize реализует алгоритм в Python (язык программирования)
• Стандартная библиотека Modelica реализует алгоритм в Modelica.
• В униор функция реализует алгоритм в R (программное обеспечение).
• В fzero функция реализует алгоритм в MATLAB.
• В Boost (библиотеки C ++) реализует два алгоритма на основе метода Брента в C ++ в инструментарии Math:

1. Минимизация функции при minima.hpp с примером минимумы функции определения местоположения.
2. В поиске корней реализован новый TOMS748, более современный и эффективный алгоритм, чем оригинал Брента, на TOMS748, и Поиск корней Boost.Math который использует TOMS748 внутри с Примеры.

• В Optim.jl пакет реализует алгоритм в Юля (язык программирования)

Рекомендации

• Брент, Р. П. (1973), "Глава 4: Алгоритм с гарантированной сходимостью для нахождения нуля функции", Алгоритмы минимизации без производных, Энглвуд Клиффс, Нью-Джерси: Прентис-Холл, ISBN  0-13-022335-2
• Деккер, Т. Дж. (1969), «Нахождение нуля с помощью последовательной линейной интерполяции», у Dejon, B .; Хенрици, П. (ред.), Конструктивные аспекты основной теоремы алгебры, Лондон: Wiley-Interscience, ISBN  978-0-471-20300-1

дальнейшее чтение

• Аткинсон, Кендалл Э. (1989). «Раздел 2.8.». Введение в численный анализ (2-е изд.). Джон Уайли и сыновья. ISBN  0-471-50023-2.
• Press, W. H .; Теукольский, С. А .; Vetterling, W. T .; Фланнери, Б. П. (2007). «Раздел 9.3. Метод Ван Вейнгаардена – Деккера – Брента». Числовые рецепты: искусство научных вычислений (3-е изд.). Нью-Йорк: Издательство Кембриджского университета. ISBN  978-0-521-88068-8.
• Алефельд, Г. Э .; Potra, F.A .; Ши Исюнь (сентябрь 1995 г.). «Алгоритм 748: заключительные нули непрерывных функций». Транзакции ACM на математическом ПО. 21 (3): 327–344. Дои:10.1145/210089.210111.

внешняя ссылка

• zeroin.f в Netlib.
• модуль brent в C ++ (также C, Fortran, Matlab) Джон Буркардт
• GSL выполнение.
• Boost C ++ выполнение.
• Python (Scipy) выполнение