Справка по API

Основных функций две: find_zero для нахождения нуля функции при некотором начальном значении или ограничивающем интервале и find_zeros для эвристического нахождения всех нулей в заданном интервале.

# Roots.find_zero — Function

Интерфейс одного из нескольких методов нахождения нулей одномерной функции, например решение .

Аргументы

Позиционные аргументы

Именованные аргументы

Расширенная справка

Начальное значение

Для большинства методов x0 — это скалярное значение, указывающее начальное значение в итеративной процедуре. (Методы секущих могут иметь кортеж, задающий их начальные значения.) Значения должны быть подтипом Number и иметь методы, определенные для float, real и oneunit.

Для заключаемых в скобки интервалов x0 задается с помощью кортежа, вектора или любого итерируемого объекта с определенными экстремумами (extrema). Заключаемый в скобки интервал ] представляет собой интервал, где и имеют разные знаки.

Возвращаемое значение

Если алгоритм завершается успешно, возвращается найденный приблизительный корень. Если алгоритм завершается сбоем, возникает ошибка ConvergenceFailed . В случае сбоя альтернативная форма solve(ZeroProblem(f,x0), M) возвращает NaN .

Указание метода

Метод указывается для того, чтобы определить, какой алгоритм следует использовать.

Дополнительные сведения см. на странице справки по каждому методу (например, ?Order1). Неэкспортируемые методы должны быть снабжены именем модуля, как в ?Roots.Schroder.

Если метод не указан, используемый по умолчанию метод зависит от x0.

По умолчанию выбираются надежные методы. Они могут быть не столь эффективны, как некоторые другие.

Указание функции

Функции передаются в качестве первых аргументов.

Для тех немногих методов, которые используют одну или несколько производных (Newton, Halley, Schroder, LithBoonkkampIJzerman(S,D) и т. д.), применяется кортеж функций. Для классических алгоритмов можно использовать функцию, возвращающую (f(x), f(x)/f'(x), [f'(x)/f''(x)]).

Необязательные аргументы (допуски, определение ограничений, трассировка)

См. раздел справки по Roots.assess_convergence со сведениями о сходимости См. страницу справки по Roots.default_tolerances(method) со сведениями о заданных по умолчанию допусках.

В целом при работе с числами с плавающей запятой сходимость не следует понимать как безусловное утверждение. Даже если математически α является ответом, а xstar — реализацией плавающей запятой, может оказаться, что f(xstar) - f(α) ≈ xstar ⋅ f'(α) ⋅ eps(α), поэтому необходимо принимать во внимание допуски, а иногда и указывать их.

Для методов Bisection сходимость гарантирована для значений Float64, поэтому допуски по умолчанию заданы равными .

Если метод заключения в скобки передается после спецификации метода, то всякий раз, когда в ходе работы алгоритма будет обнаружена скобка, метод будет переключаться на метод заключения в скобки для нахождения нуля. (Методы заключения в скобки математически гарантированно сходятся, методы без заключения в скобки могут сходиться, а могут и нет.) Это то, что Order0 делает по умолчанию — первоначальный метод секущей переключается на метод AlefeldPotraShi при обнаружении скобки.

Примечание. На порядок выполнения метода указывает схема именования. Схема имеет порядок r, если при eᵢ = xᵢ - α, eᵢ₊₁ = C⋅eᵢʳ. Если погрешность eᵢ достаточно мала, то, по сути, на каждом шаге она будет набирать в r раз больше начальных нулей. Однако если погрешность велика, этого не произойдет. Без хороших начальных предположений метод высокого порядка может сходиться медленно, если вообще будет сходиться. Методы OrderN имеют ряд эвристических правил, используемых для получения более широкого диапазона сходимости ценой не совсем точной реализации метода, хотя они доступны через неэкспортированные методы.

Примеры:

Методы по умолчанию.

При указании метода

Изменение допусков.

Трассировка

При передаче verbose=true будет показана подробная информация о шагах алгоритма. Алгоритм tracks позволяет передать объект Roots.Tracks для записи значений x и f(x), используемых в алгоритме.

# Roots.find_zeros — Function

Выполняет поиск нулей f в интервале [a,b] с помощью эвристического алгоритма.

Возвращает вектор нулей в отсортированном порядке, возможно, пустой.

Расширенная справка

Примеры

Основной алгоритм проверяет наличие нулей в конечных точках, затем разбивает интервал (a,b) на подынтервалы no_pts-1, а затем переходит к поиску нулей с помощью метода деления пополам или метода без производной. Поскольку проверка на наличие заключаемого в скобки интервала относительно дешева, а метод деления пополам гарантированно сходится, в каждом интервале есть k пар промежуточных точек, проверяемых на наличие скобки.

При нахождении нулей алгоритм использует их для разбиения интервала (a,b) на подынтервалы. Каждый подынтервал уменьшается так, чтобы конечные точки не были нулями, и процесс повторяется в этом подынтервале. Если начальный интервал слишком велик, примитивная проверка на наличие нулей может оказаться безрезультатной, и нули не будут обнаружены. Если есть близлежащие нули, при уменьшении интервала они могут быть пропущены, хотя, как видно из примеров, близлежащие корни могут быть найдены правильно, хотя для действительно близких точек или очень плоских функций может помочь увеличение no_pts.

Допуски используются для уменьшения интервалов, но не для поиска нулей в диапазоне поиска. При поисках метод деления пополам гарантированно сходится без указанного допуска. Для поиска без производной используется модификация метода Order0, который в худшем случае сравнивает |f(x)| <= 8*eps(x) для нахождения нуля. Алгоритм может определить несколько значений для нуля из-за аппроксимаций с плавающей запятой. Если потенциальная пара нулей удовлетворяет isapprox(a,b,atol=sqrt(xatol), rtol=sqrt(xrtol)), они объединяются.

Алгоритм может выполнять множество вызовов функции. При обнаружении нулей в интервале примитивный поиск выполняется в каждом подынтервале. Чтобы сократить количество вызовов функции, но с повышенной вероятностью пропуска нескольких нулей, адаптивный характер можно отменить, задав аргумент naive=true или уменьшив количество точек.

Алгоритм заимствован из PR у @djsegal..

Например, эта функция (из-за @truculentmath) особенно сложна, поскольку она положительна при каждом числе с плавающей запятой, но имеет два нуля (вертикальная асимптота при 15//11 отрицательна только в смежных значениях с плавающей запятой):

Менее экстремальным вариантом может быть следующий, где unique указывает на то, что метод деления пополам может быть полезен, и действительно find_zeros определяет эти значения:

Интерфейс «задача-алгоритм-решение» — это популярный шаблон в Julia, реализуемый набором пакетов DifferentialEquations.jl. Его можно использовать в качестве альтернативы функции find_zero. В отличие от find_zero, solve при несходимости возвращает NaN.

# CommonSolve.solve! — Function

Solve for the zero of a scalar-valued univariate function specified through ZeroProblem or ZeroProblemIterator using the CommonSolve interface.

The defaults for M and N depend on the ZeroProblem: if x0 is a number, then M=Secant() and N=AlefeldPotraShi() is used (Order0); if x0 has 2 (or more values) then it is assumed to be a bracketing interval and M=AlefeldPotraShi() is used.

The methods involved with this interface are:

The latter calls the following, which can be useful independently:

Returns NaN, not an error like find_zero, when the problem can not be solved. Tested for zero allocations.

Examples:

Or, if the iterable is required

keyword arguments can be used to adjust the default tolerances.

The above is equivalent to:

The keyword argument p may be used if the function(s) to be solved depend on a parameter in their second positional argument (e.g., f(x, p)). For example

This would be recommended, as there is no recompilation due to the function changing. For use with broadcasting, p may also be the last positional argument.

The argument verbose=true for init instructs that steps to be logged;

The iterator interface allows for the creation of hybrid solutions, such as is used when two methods are passed to solve. For example, this is essentially how the hybrid default is constructed:

# CommonSolve.solve — Function

Выполняет диспетчеризацию в solve!(init(fx, args...; kwargs...)). Дополнительные сведения см. в описании solve!.

# Roots.ZeroProblem — Type

Контейнер для функции и начального предположения, используемых с solve.

Начнем с описания классических методов, хотя они необязательно входят в число рекомендуемых, поскольку требуют больше усилий от пользователя. Однако это позволит получить представление о том, почему существует столько разнообразных методов.

Классические методы Ньютона и Галлея используют информацию о функции и ее производных для итеративного приближения к нулю при заданном начальном значении.

Метод Ньютона описывается легко:

Следующая за начальной точка в итерационном алгоритме находится как точка пересечения оси с касательной к в начальной точке. Этот процесс повторяется до сходимости или пока не станет понятно, что сходимость для начальной точки невозможна. Математически это можно записать так:

Для понимания методов, доступных в Roots, полезно знать ряд фактов:

# Roots.Newton — Type

Реализует метод Ньютона: xᵢ₊₁ = xᵢ - f(xᵢ)/f'(xᵢ). Это квадратично сходящийся метод, требующий одну производную и два вызова функции для каждого шага.

Примеры

Если вычисления функций являются дорогостоящими, можно передать функцию, которая возвращает (f, f/f') следующим образом:

Это может быть выгодно, если производная легко вычисляется из значения f, но в противном случае вычисление будет дорогостоящим.

Погрешность eᵢ = xᵢ - α может быть выражена в виде eᵢ₊₁ = f[xᵢ,xᵢ,α]/(2f[xᵢ,xᵢ])eᵢ² (методы Сиди, унифицированной обработки ложного положения, Ньютона-Рафсона, секущей и Стеффенсена для нелинейных уравнений).

# Roots.Halley — Type

Реализует метод Галлея `xᵢ₊₁ = xᵢ

сходится кубически, но требует вызова функции на шаг. Метод Галлея находит xₙ₊₁ как нуль гиперболы в точке (xₙ, f(xₙ)), совпадающей с первой и второй производными f.

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

Примеры

Это может быть выгодно, если производные легко вычисляются из вычислений для f, но в противном случае вычислять их отдельно будет дорого.

Погрешность eᵢ = xᵢ - α удовлетворяет условию eᵢ₊₁ ≈ -(2f'⋅f''' -3⋅(f'')²)/(12⋅(f'')²) ⋅ eᵢ³ (все вычисляется по α).

# Roots.QuadraticInverse — Type

Реализует квадратичный обратный метод, также известный как метод Чебышева, xᵢ₊₁ = xᵢ - (f/f')(xᵢ) * (1 + (f/f')(xᵢ) * (f''/f')(xᵢ) * 1/2). Этот метод сходится кубически и требует вызова функции на каждый шаг.

Пример

Если вычисления функций являются дорогостоящими, можно передать функцию, которая возвращает (f, f/f',f'/f''), следующим образом:

Это может быть выгодно, если производные легко вычисляются из вычислений для f, но в противном случае вычислять их отдельно будет дорого.

Погрешность eᵢ = xᵢ - α удовлетворяет условию eᵢ₊₁ ≈ (1/2⋅(f''/f')² - 1/6⋅f'''/f')) ⋅ eᵢ³ (все вычисляется по α).

# Roots.ChebyshevLike — Type

CHEBYSHEV-LIKE METHODS AND QUADRATIC EQUATIONS (J. A. EZQUERRO, J. M. GUTIÉRREZ, M. A. HERNÁNDEZ and M. A. SALANOVA)

# Roots.SuperHalley — Type

An acceleration of Newton’s method: Super-Halley method (J.M. Gutierrez, M.A. Hernandez)

Методы Ньютона и Галлея входят в это семейство методов:

# Roots.LithBoonkkampIJzerman — Type

Семейство различных методов, включающее метод секущей и метод Ньютона.

Определяет линейный многошаговый решатель с S шагами и D производными, как описано в статье Лит (Lith), Бункампа (Boonkkamp) и Айзермана (IJzerman) (https://doi.org/10.1016/j.amc.2017.09.003).

Расширенная справка

Примеры

Метод может быть более устойчивым к начальным условиям. Это пример из статьи (стр. 13). Метод Ньютона (случай S=1, D=1) завершается сбоем, если |x₀| ≥ 1.089, но методы с большим объемом памяти работают успешно.

Множество производных можно создать автоматически с помощью автоматического дифференцирования. Например,

Для случая D=1 метод заключения в скобки, основанный на этом подходе, реализован в LithBoonkkampIJzermanBracket.

Справка

В статье Лита, Бункампа и Айзермана (https://doi.org/10.1016/j.amc.2017.09.003) приведен анализ скорости сходимости при использовании линейных многошаговых методов для решения 0=f(x) с помощью f⁻¹(0) = x, когда f является достаточно сглаженной линейной функцией. Переформулировка, приписываемая Грау-Санчесу, находит дифференциальное уравнение для f⁻¹: dx/dy = [f⁻¹]′(y) = 1/f′(x) = F в виде x(0) = x₀ + ∫⁰_y₀ F(x(y)) dy.

Для численного решения этого уравнения используется линейный многошаговый метод. Пусть S — количество шагов памяти (S= 1,2,…​), а D — количество используемых производных, тогда при F(x) = dx/dy x_{n+S} = ∑_{k=0}^{S-1} aₖ x_{n+k} +∑d=1^D ∑_{k=1}^{S-1} aᵈ_{n+k}F⁽ᵈ⁾(x_{n+k}). Значения aₖ и aᵈₖ вычисляются на каждом шаге.

Эта таблица взята из таблиц 1 и 3 статьи, и в ней показана скорость сходимости для простых корней, указанных в ней:

То есть больший объем памяти приводит к более высокой скорости сходимости. Большее количество производных приводит к более высокой скорости сходимости. Однако интервал около α, нуль, где скорость сходимости гарантирована, может уменьшиться.

В методе секущей вместо члена производной, применяемого в методе Ньютона, используется угловой коэффициент секущей с двумя предыдущими значениями:

Хотя метод секущей имеет скорость сходимости порядка , то есть не является квадратичным, он требует лишь одного нового вызова функции за шаг и поэтому может быть очень эффективным. Кроме того, зачастую вычисление функций является самой трудоемкой частью, а в данном случае производная не требуется. Из-за своей потенциальной эффективности метод секущей применяется по умолчанию при вызове функции find_zero с одной начальной точкой.

Метод Стеффенсена — это метод с квадратичной сходимостью без вычисления производных, который использует секущую через и . Хотя он имеет более высокий порядок, он требует дополнительных вызовов функции за шаг и правильно выбранной начальной точки. Есть и другие методы без вычисления производных, в которых за счет увеличения вызовов функции достигается сходимость более высокого порядка. Они могут представлять интерес, когда нужна произвольная точность. Мерой эффективности является , где  — это порядок сходимости, а  — количество вызовов функции за шаг. Для метода секущей эта мера равна , а для метода Стеффенсена она меньше ( ).

# Roots.Secant — Type

Метод Order1() является псевдонимом для Secant. Он задает метод секущей. Этот метод хранит в своем состоянии два значения — xₙ и xₙ₋₁. Обновленная точка является точкой пересечения оси с секущей, сформированной из двух точек. Метод секущей использует вычисление функции на каждом шаге и имеет порядок φ≈ (1+sqrt(5))/2.

Погрешность eᵢ = xᵢ - α удовлетворяет условию eᵢ₊₂ = f[xᵢ₊₁,xᵢ,α] / f[xᵢ₊₁,xᵢ] * (xᵢ₊₁-α) * (xᵢ - α).

# Roots.Order1 — Type

Метод Order1() является псевдонимом для Secant. Он задает метод секущей. Этот метод хранит в своем состоянии два значения — xₙ и xₙ₋₁. Обновленная точка является точкой пересечения оси с секущей, сформированной из двух точек. Метод секущей использует вычисление функции на каждом шаге и имеет порядок φ≈ (1+sqrt(5))/2.

Погрешность eᵢ = xᵢ - α удовлетворяет условию eᵢ₊₂ = f[xᵢ₊₁,xᵢ,α] / f[xᵢ₊₁,xᵢ] * (xᵢ₊₁-α) * (xᵢ - α).

# Roots.Steffensen — Type

Квадратично сходящийся метод Стеффенсена используется для алгоритма Order2() без производных. В отличие от квадратично сходящегося метода Ньютона, производные не требуются, но, как и в методе Ньютона, на каждый шаг приходится два вызова функции. Алгоритм Стеффенсена более чувствителен, чем метод Ньютона, к плохим начальным предположениям, когда f(x) велика из-за того, как f'(x) аппроксимируется. Метод Order2 заменяет шаг Стеффенсена шагом секущей, когда значение f(x) велико.

Погрешность eᵢ - α удовлетворяет условию eᵢ₊₁ = f[xᵢ, xᵢ+fᵢ, α] / f[xᵢ,xᵢ+fᵢ] ⋅ (1 - f[xᵢ,α] ⋅ eᵢ²

# Roots.Order2 — Type

Steffensen с защитой на шаге секущей.

# Roots.Order5 — Type

Реализует алгоритм 5-го порядка из статьи A New Fifth Order Derivative Free Newton-Type Method for Solving Nonlinear Equations от Manoj Kumar, Akhilesh Kumar Singh и Akanksha, Appl. Math. Inf. Sci. 9, No. 3, 1507-1513 (2015), DOI: 10.12785/amis/090346. Для каждого шага требуется пять вызовов функции. Метод Order5 заменяет шаг Стеффенсена шагом секущей, когда значение f(x) велико.

Погрешность eᵢ = xᵢ - α удовлетворяет условию eᵢ₊₁ = K₁ ⋅ K₅ ⋅ M ⋅ eᵢ⁵ + O(eᵢ⁶)

# Roots.Order8 — Type

Реализует алгоритм 8-го порядка из статьи New Eighth-Order Derivative-Free Methods for Solving Nonlinear Equations от Rajinder Thukral, International Journal of Mathematics and Mathematical Sciences Volume 2012 (2012), Article ID 493456, 12 pages DOI: 10.1155/2012/493456. Для каждого шага требуется четыре вызова функции. Метод Order8 заменяет шаг Стеффенсена шагом секущей, когда значение f(x) велико.

Погрешность eᵢ = xᵢ - α выражается как eᵢ₊₁ = K ⋅ eᵢ⁸ в (2.25) статьи для явного K.

# Roots.Order16 — Type

Реализует алгоритм 16-го порядка из статьи New Sixteenth-Order Derivative-Free Methods for Solving Nonlinear Equations от R. Thukral, American Journal of Computational and Applied Mathematics p-ISSN: 2165-8935; e-ISSN: 2165-8943; 2012; 2(3): 112-118 DOI: 10.5923/j.ajcam.20120203.08.

Для каждого шага требуется пять вызовов функции. Хотя этот метод быстро сходится, в целом он не быстрее (меньше вызовов функций/шагов) других методов при использовании значений Float64, но может быть полезен для решения задач с BigFloat. Метод Order16 заменяет шаг Стеффенсена шагом секущей, когда значение f(x) велико.

Погрешность eᵢ = xᵢ - α выражается как eᵢ₊₁ = K⋅eᵢ¹⁶ для явного K в уравнении (50) из статьи.

Методом деления пополам находится нуль непрерывной функции между и , когда и имеют разные знаки. (Интервал ] называется ограничивающим, когда .) Базовый алгоритм очень простой: интервал ] делится по . Либо , либо один из интервалов ] или ] является ограничивающим. Он обозначается ]. Из этого описания видно, что ] имеет длину в раз большую, чем длина ], поэтому в конечном итоге будет найден нуль или либо интервалы сойдутся к нулю. Данная сходимость медленная (эффективность составляет всего ), но гарантированная. Для 16-, 32- и 64-битных значений с плавающей запятой переинтерпретация способа нахождения средней точки ( ) обеспечивает сходимость не более чем за итераций, в отличие от описанного выше способа, который в некоторых случаях может требовать гораздо большего числа шагов для сходимости.

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

При описанном выше способе выбора средней точки по умолчанию учитывается только знак функции . Алгоритмы, учитывающие форму функции, могут быть значительно эффективнее. Например, ограничивающий метод Roots.AlefeldPotraShi по Алефельду, Потра и Ши имеет эффективность . Этот метод также применяется по умолчанию в функции find_zero, когда дана одна начальная точка и определен ограничивающий интервал.

# Roots.Bisection — Type

Если возможно, используйте метод деления пополам для значений Float64. Метод заключения в скобки начинается с заключаемого в скобки интервала [a,b] и разбивает его на два интервала [a,c] и [c,b]. Если c не является нулем, то один из этих двух интервалов будет заключен в скобки и процесс продолжится. c вычисляется с помощью _middle, которая по-иному интерпретирует значения с плавающей запятой в целые числа без знака и выполняет разбиение. Его автором является Джейсон Меррилл (Jason Merrill). Этот метод позволяет избежать проблем с плавающей запятой и при нулевом значении допусков (по умолчанию) гарантирует «лучшее» решение (такое, в котором найден нуль или интервал скобок имеет тип [a, nextfloat(a)]).

Если заданы допуски, алгоритм завершается, когда длина интервала меньше или равна допуску max(δₐ, 2abs(u)δᵣ), причем u в {a,b} выбирается меньшим из |f(a)| и |f(b)|, или значение функции меньше max(tol, min(abs(a), abs(b)) * rtol). Последний вариант используется только в том случае, если настроены допуски по умолчанию (atol или rtol).

При решении для с помощью Bisection нельзя брать производную непосредственно через автоматическое дифференцирование, так как алгоритм не является дифференцируемым. Альтернативные варианты см. в разделе Чувствительность документации.

# Roots.A42 — Type

Метод заключения в скобки, который позволяет найти корень непрерывной функции в заданном заключаемом в скобки интервале [a, b], не требуя производных. Он основан на алгоритме 4.2, описанном в документе: G. E. Alefeld, F. A. Potra, and Y. Shi, "Algorithm 748: enclosing zeros of continuous functions," ACM Trans. Math. Softw. 21, 327—​344 (1995), DOI: 10.1145/210089.210111. Индекс асимптотической эффективности, , равен .

Автор: Джон Треверс (John Travers).

# Roots.AlefeldPotraShi — Type

Следует алгоритму 4.1 в статье ON ENCLOSING SIMPLE ROOTS OF NONLINEAR EQUATIONS от Alefeld, Potra, Shi; DOI: 10.1090/S0025-5718-1993-1192965-2.

Порядок сходимости — 2 + √5. Асимптотически на каждом шаге выполняется 3 вычисления функции. Индекс асимптотической эффективности равен . Менее эффективен, но может работать быстрее, чем связанный метод A42.

Автор: Джон Треверс (John Travers).

# Roots.Brent — Type

Реализация метода Брента (или Брента-Деккера). Этот метод использует обратную квадратичную интерполяцию или шаг секущей, при необходимости возвращаясь к шагу поиска делением пополам.

# Roots.Chandrapatla — Type

Использует алгоритм Чандрапатлы (ср. метод Шеррера) для решения .

Алгоритм Чандрапатлы выбирает обратный квадратичный шаг или шаг поиска делением пополам на основе вычисленного неравенства.

# Roots.Ridders — Type

Реализует метод Риддерса. Этот метод заключения в скобки находит среднюю точку x₁, затем интерполирует экспоненту, а затем использует ложное положение с интерполированным значением для нахождения c. Используется, если c и x₁ образуют скобку, в противном случае применяется подынтервал [a,c] или [c,b].

Пример:

Риддерс (Ridders) показал, что погрешность eₙ₊₁ ≈ 1/2 eₙeₙ₋₁eₙ₋₂ ⋅ (g^2-2fh)/f удовлетворяет условию для f=F', g=F''/2, h=F'''/6, что предполагает сходимость со скоростью ≈ 1.839.... Использует четыре вычисления функции на каждом шаге, поэтому порядок сходимости составляет ≈ 1.225....

# Roots.ITP — Type

Использует метод заключения в скобки ITP . Этот метод утверждает, что он «является первым алгоритмом поиска корней, который достигает суперлинейной сходимости метода секущей, сохраняя при этом оптимальную производительность метода деления пополам в наихудшем случае».

Значения κ1, κ₂ и n₀ являются параметрами настройки.

Примечание.

Предложено на форуме discourse пользователем @TheLateKronos, который предоставил исходную версию кода.

# Roots.FalsePosition — Type

Использует метод ложного положения для нахождения нуля для функции f в заключаемом в скобки интервале [a,b].

Метод ложного положения представляет собой модифицированный метод деления пополам, в котором средняя точка между [aₖ, bₖ] выбирается в качестве точки пересечения секущей линии с осью , а не среднее между двумя значениями.

Чтобы ускорить сходимость для вогнутых функций, этот алгоритм реализует понижающих коэффициентов Галдино (A family of regula falsi root-finding methods). Они задаются числом, как в FalsePosition(2), или одним из трех имен FalsePosition(:pegasus), FalsePosition(:illinois) или FalsePosition(:anderson_bjork) (по умолчанию). Вариант по умолчанию, как правило, имеет лучшую производительность, чем остальные, хотя бывают и исключения.

Для некоторых задач количество вызовов функции может быть больше, чем для метода Bisection, но в целом этот алгоритм выполняет меньше вызовов функции.

Примеры

# Roots.LithBoonkkampIJzermanBracket — Type

Метод заключения в скобки, являющийся модификацией метода Брента, согласно статье Лита, Бункампа и Айзермана (IJzerman). Наилучшая возможная скорость сходимости составляет 2,91.

Необходимо указать функцию, ее производную и заключаемый в скобки интервал.

Состояние включает 3 точки — скобку [a,b] (в b=xₙ f(b) ближе всего к 0) и c=xₙ₋₁ — и соответствующие значения функции и ее производной в этих трех точках.

Следующим предлагаемым шагом является выбор S=2 или S=3 для методов LithBoonkkampIJzerman с информацией о производной, включаемой только в тех случаях, если она может быть полезной. Предлагаемый шаг изменяется, если это не работает. Предлагаемый вариант сравнивается с шагом поиска делением пополам. В качестве следующего шага выбирается тот, который находится в скобках и имеет меньшее значение функции.

# Roots.BracketedHalley — Type

Для скобки [a,b] использует метод Roots.Halley, начиная со значения x, для которого fa или fb ближе всего к 0. Использует платформу Roots.AbstractAlefeldPotraShi для применения метода заключения в скобки, каждый раз выполняя дополнительный шаг двойной секущей.

# Roots.BracketedChebyshev — Type

Для скобки [a,b] использует метод Roots.QuadraticInverse, начиная со значения x, для которого fa или fb ближе всего к 0. Использует платформу Roots.AbstractAlefeldPotraShi для применения метода заключения в скобки, каждый раз выполняя дополнительный шаг двойной секущей.

# Roots.BracketedSchroder — Type

Для скобки [a,b] использует метод Roots.Schroder, начиная со значения x, для которого fa или fb ближе всего к 0. Использует платформу Roots.AbstractAlefeldPotraShi для применения метода заключения в скобки, каждый раз выполняя дополнительный шаг двойной секущей.

Порядок сходимости для большинства методов соответствует простым нулям, то есть значениям где , причем отлично от нуля. Если методы имеют порядок для непростых нулей, обычно требуется дополнительный вызов функции за шаг. Например, это можно заметить, если сравнить Roots.Newton и Roots.Schroder.

Для методов без вычисления производных для непростых нулей реализованы следующие типы:

# Roots.King — Type

Суперлинейная (порядка 1.6...) модификация метода секущей для нескольких корней. Представлено в статье A SECANT METHOD FOR MULTIPLE ROOTS, автор RICHARD F. KING, BIT 17 (1977), 321-328

Основная идея аналогична методу Шредера: применение метода секущей к f/f'. Однако здесь используется f' ~ fp = (fx - f(x-fx))/fx (шаг Стеффенсена). В этой реализации Order1B, когда значение fx слишком велико, используется шаг одной секущей f .

Асимптотическая погрешность eᵢ = xᵢ - α задается следующим образом: eᵢ₊₂ = 1/2⋅G''/G'⋅ eᵢ⋅eᵢ₊₁ + (1/6⋅G'''/G' - (1/2⋅G''/G'))^2⋅eᵢ⋅eᵢ₊₁⋅(eᵢ+eᵢ₊₁).

# Roots.Order1B — Type

Метод King с защищенным шагом секущей.

# Roots.Esser — Type

Метод Эссера. Это квадратично сходящийся метод, который, как и метод Шредера, не зависит от кратности нуля. Метод Шердера имеет обновленный шаг x - r2/(r2-r1) * r1, где ri = fⁱ⁻¹/fⁱ. Esser approximates Эссер считает, что f' ~ f[x-h, x+h], f'' ~ f[x-h,x,x+h], where , где h = fx, как и в методе Стеффенсена, требуется 3 вызова функции на каждый шаг. Реализация Order2B использует шаг секущей, когда значение |fx| считается слишком большим.

Esser, H. Computing (1975) 14: 367. DOI: 10.1007/BF02253547 Eine stets quadratisch konvergente Modification des Steffensen-Verfahrens

Примеры

# Roots.Order2B — Type

Метод Esser с защищенным шагом секущей.

Как показал Шредер, в случае с непростыми нулями можно использовать дополнительную производную, чтобы получить квадратичную сходимость на основе метода Ньютона:

# Roots.Schroder — Type

Метод Шредера, как и метод Галлея, использует f, f' и f''. В отличие от метода Галлея, он сходится квадратично, но не зависит от кратности нуля (ср. Schröder, E. Über unendlich viele Algorithmen zur Auflösung der Gleichungen. Math. Ann. 2, 317-365, 1870; mathworld).

Метод Шредера применяет метод Ньютона к f/f' — функции со всеми простыми нулями.

Пример

(В то время как при m=1 в методе Галлея выполняется 2 шага, а в методе Шредера — 3.)

Если вычисления функций являются дорогостоящими, можно передать функцию, которая возвращает (f, f/f',f'/f''), следующим образом:

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

Погрешность eᵢ = xᵢ - α такая же, как в методе Newton, с f замененной f/f'.

Семейство методов для непростых нулей, которые требуют производных порядка и сводятся к методу Шредера при , реализовано в следующем типе:

# Roots.AbstractThukralBMethod — Type

Абстрактный тип для методов ThukralXB для X, являющихся 2, 3, 4 или 5.

Это семейство методов, которые

функция, возвращающая коэффициенты: x -> (f(x), f(x)/f'(x), f'(x)/f''(x), …).

Примеры

Справка

Introduction to a family of Thukral -order method for finding multiple zeros of nonlinear equations, R. Thukral, JOURNAL OF ADVANCES IN MATHEMATICS 13(3):7230-7237, DOI: 10.24297/jam.v13i3.6146.

Полезной может оказаться такая стратегия: начинаем с неограничивающего метода, а затем, если встретится ограничение, переходим на ограничивающий метод. Это позволяет находить нули, не заключенные в ограничивающий интервал, и получать гарантированную сходимость при достижении ограничения. Такая стратегия используется по умолчанию функцией find_zero(f,a).

# Roots.Order0 — Type

Метод Order0 разработан как более надежная, хотя, возможно, и более медленная альтернатива другим методам нахождения корней без производных. Реализация примерно следует алгоритму, описанному в статье Personal Calculator Has Key to Solve Any Equation , кнопка SOLVE из HP-34C. Основная идея заключается в использовании шага секущей. Если по пути будет найдена скобка, следует переключиться на алгоритм заключения в скобки, используя AlefeldPotraShi. Если шаг секущей не уменьшает значение функции, используется шаг квадратичного метода до раз.

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

Порядок метода —  , где e_{i+1} \approx e_i^q. Известно, что метод Ньютона является квадратичным для простых корней; метод секущей порядка . Однако для метода Ньютона требуется вызовов, а для метода секущей — только . Асимптотическая эффективность равна ценой увеличения числа вызовов функций. Есть и другие методы порядка , требующие вызовов функции за шаг, например метод Галлея. Другие же требуют меньше шагов, как показано ниже. Многие методы используют обратные квадратичные шаги, а некоторые — обратные кубические. Они имеют порядок при решении уравнения ( для квадратичного случая). В случае с надежными методами для достижения скорости сходимости, как правило, требуется дополнительный вызов функции, что хорошо видно на примере метода Schroder.

Для определения того, сходится ли алгоритм или расходится, требуется указание допусков и условий сходимости.

В случае с точным делением пополам сходимость гарантируется математически. Для чисел с плавающей запятой либо находится точный нуль, либо ограничивающий интервал может быть разделен на ], где и  — соседние значения с плавающей запятой. То есть в случае с числами с плавающей запятой имеет минимальный возможный размер. Это можно рассматривать как условие остановки в . Для преждевременного завершения (обеспечивающего меньшую точность, но и меньшее число вызовов функции) можно задать допуск, чтобы при достаточно небольшом алгоритм успешно завершал выполнение. В случае со значениями с плавающей запятой для оценки равенства требуются два допуска: относительный допуск, так как минимальная разница значений с плавающей запятой зависит от размера и , и абсолютный допуск для значений в окрестности . Для определения близости в функцию Base.isapprox передаются значения xrtol и xatol.

Не во всех задачах можно полагаться на близость двух значений . В некоторых примерах разность может быть весьма малой и даже равна , хотя значение не близко к нулю . В связи с этим для неограничивающих методов также применяется проверка размера . При поиске аппроксимаций с плавающей запятой к , то есть нуля, необходимо учитывать значения, малые при малом . При использовании аппроксимации Тейлора можно ожидать, что они будут в районе . Таким образом, малая величина зависит от размера и производной в точке . Первая зависимость учитывается посредством относительной и абсолютной погрешностей (rtol и atol). Размер зависит от задачи, и для его учета можно увеличивать относительную или абсолютную погрешность.

Когда алгоритм возвращает значение NaN, его выполнение завершается. Это может происходить при приближении к сходимости либо указывать на наличие проблем. При досрочном завершении сходимость проверяется в пределах размера с увеличенной погрешностью, если задан аргумент strict=false (по умолчанию).

Для определения того, будет ли алгоритм завершать выполнение, подсчитывается количество итераций; значение по умолчанию настраивается с помощью аргумента maxiters. Так как большинство алгоритмов являются сверхлинейными, в окрестности решения сходимость достигается быстро, однако на приближение к решению любому алгоритму может потребоваться время, даже если его выполнение продвигается. В связи с этим максимум должен быть достаточно большим для учета линейных случаев, но не слишком во избежание слишком большого количества шагов, когда алгоритм не сходится.

Условия сходимости зависят от метода и оцениваются методами Roots.assess_convergence.

# Roots.assess_convergence — Function

Определяет сходимость алгоритма.

Возвращает флаг сходимости и логическое значение, указывающее, завершился ли алгоритм (со схождением или без схождения).

Если алгоритм не сошелся, возвращается (:not_converged, false).

Если алгоритм остановился или сошелся, возвращается флаг и true. Существуют следующие флаги:

Не проверяет ни количество выполненных шагов, ни количество вычислений функции.

В decide_convergence значения остановки (и :x_converged при strict=false) проверяются на сходимость с нестрогим допуском.

Погрешности по умолчанию задаются посредством методов Roots.default_tolerances.

# Roots.default_tolerances — Function

По умолчанию для большинства методов используются допуски xatol=eps(T), xrtol=eps(T), atol=4eps(S) и rtol=4eps(S) с соответствующими единицами (абсолютные допуски имеют единицы x и f(x); относительные допуски единицы не используют). Для значений Complex{T} используется T.

Количество итераций ограничивается maxiters=40.

Для метода Bisection, когда значения x имеют тип Float64, Float32 или Float16, допуски по умолчанию равны нулю, а количество итераций не ограничено. В этом случае алгоритм гарантированно сходится к точному нулю или точке, в которой функция изменяет знак одного из соседних значений ответа с плавающей запятой.

Для других типов по умолчанию заданы ненулевые допуски для xatol и xrtol.

Абстракции и многие проверки сходимости, используемые функцией find_zero, негативно сказываются на производительности. Если это критически важно, есть несколько «простых» методов, обеспечивающих более высокую производительность.

# Roots.secant_method — Function

Выполняет метод секущей для решения f(x) = 0.

Метод секущей — это итеративный метод с шагом обновления, заданным с помощью b - fb/m, где m — наклон секущей линии между (a,fa) и (b,fb).

Начальные значения могут быть указаны как пара двух значений, как в (x₀, x₁) или [x₀, x₁], или как одно значение x₁, в этом случае выбирается значение x₀.

Алгоритм возвращает m, когда abs(fm) <= max(atol, abs(m) * rtol). Если это не произойдет до шагов maxevals, или алгоритм столкнется с проблемой, будет возвращено значение NaN. Если выполнено слишком много шагов, текущее значение проверяется на предмет изменения знака для соседних значений с плавающей запятой.

Метод Order1 для find_zero также реализует метод секущей. Этот вариант должен быть немного быстрее, поскольку требуется меньше затрат на настройку.

Примеры:

# Roots.bisection — Function

Выполняет метод деления пополам для нахождения нуля непрерывной функции.

Предполагается, что (a,b) является скобкой, то есть функция имеет разные знаки в точках a и b. Интервал (a,b) преобразуется в значение с плавающей запятой и уменьшается, когда a или b являются бесконечными. В типичном случае функция f может быть бесконечной. Если f не является непрерывной, алгоритм может найти не только нули, но и точки перехода через ось x.

Если заданы нетривиальные допуски, процесс завершится, когда скобка (a,b) удовлетворяет условию isapprox(a, b, atol=xatol, rtol=xrtol) Для нулевых допусков, по умолчанию, для значений. For zero tolerances, the default, for Float64, Float32 или значений Float16 процесс завершится при значении x, при котором f(x)=0 или f(x)*f(prevfloat(x)) < 0 или f(x) * f(nextfloat(x)) < 0. For other number types, the Для других числовых типов используется метод Roots.A42.

# Roots.muller — Function

Отрывок и алгоритм взяты из статьи

Сходимость здесь определяется тем, что xᵢ₊₁ ≈ xᵢ с использованием указанных допусков, которые оба по умолчанию имеют значение eps(one(typeof(abs(xᵢ))))^4/5 в соответствующих единицах. Каждая итерация выполняет три вычисления f. Первый метод случайным образом выбирает две оставшиеся точки в относительной близости от xᵢ.

Обратите внимание, что метод может вернуть сложный результат даже для реальных исходных значений, так как это зависит от функции.

Примеры:

# Roots.newton — Function

Реализация метода Ньютона: xᵢ₊₁ = xᵢ - f(xᵢ)/f'(xᵢ).

Аргументы:

Пакет ForwardDiff позволяет вычислять производные автоматически. Например, при определении D(f) = x -> ForwardDiff.derivative(f, float(x)) можно использовать D(f) для первой производной.

Именованные аргументы передаются функции find_zero с помощью метода Roots.Newton().

Более простые реализации см. также в описании Roots.newton((f,fp), x0) и Roots.newton(fΔf, x0).

# Roots.dfree — Function

Более надежная реализация метода секущей

Решение для f(x) = 0 с использованием алгоритма из статьи Personal Calculator Has Key to Solve Any Equation f(x) = 0, кнопка SOLVE из HP-34C.

Также реализовано как метод Order0 для find_zero.

Начальные значения могут быть указаны как пара двух значений, как в (a,b) или [a,b], или как одно значение, в этом случае вычисляется значение b, возможно, из fb. Основная идея заключается в том, чтобы следовать методу секущих, за исключением случаев, когда

Сходимость наступает, когда f(m) == 0, происходит смена знака между m и соседним значением с плавающей запятой или f(m) <= 2^3*eps(m).

Значение NaN возвращается, если алгоритм делает слишком много шагов перед нахождением нуля.

Примеры

Согласно исходной схеме именования функция называлась fzero, а не fzeros, то есть так же, как аналогичная функция в MATLAB — fzero. Этот интерфейс не рекомендуется к использованию, но все еще поддерживается.

# Roots.fzero — Function

Находит нуль функции, используя один из нескольких итеративных алгоритмов.

Примеры:

# Roots.fzeros — Function

Выполняет поиск всех нулей f в интервале (a,b). Предполагается, что ни a, ни b не равны нулю.

Интерфейс совместимости для find_zeros.

При вызове функции find_zero можно добавить именованный аргумент verbose=true для получения подробных сведений о решении и данных по каждой итерации. Для сохранения этих данных можно передать объект Tracks в аргументе tracks.

# Roots.Tracks — Type

Экземпляр Tracks используется для записи хода выполнения алгоритма. T — это тип входных данных функции, а S — тип выходных данных функции. Оба по умолчанию имеют тип Float64. Обратите внимание, что поскольку этот тип не экспортируется, вам придется написать Roots.Tracks(), чтобы создать объект Tracks.

По умолчанию при нахождении корня отслеживание не ведется. Для изменения этого поведения создайте объект Tracks и передайте его именованному аргументу tracks. При этом изменится объект Tracks, хранящий значения входных данных и функции при каждой итерации, а также дополнительную информацию о процессе поиска корней.

Объекты Tracks отображаются в удобном для чтения формате. На внутреннем уровне хранится либо кортеж пар (x,f(x)), либо пары (aₙ, bₙ), последние — для методов заключения в скобки. (Эти особенности реализации могут быть изменены без уведомления.) Интерес могут представлять методы empty! для сброса объекта Tracks, get для получения отслеживаний, last для получения значения преобразования.

Если вы хотите только вывести только информацию, но в дальнейшем она вам не понадобится, можно передать verbose=true функции нахождения корней. Это не повлияет на возвращаемое значение, которое по-прежнему будет корнем функции.

Примеры