Roots.jl

Документация по Roots.jl

Общие сведения

Roots — это пакет Julia для нахождения нулей непрерывных скалярных функций одной вещественной переменной с использованием чисел с плавающей запятой. То есть он служит для решения уравнений вида для с учетом специфических особенностей чисел с плавающей запятой.

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

Методы наподобие метода деления пополам. Для функций с известным ограничивающим интервалом (таким, в котором и чередуются по знаку) есть несколько ограничивающих методов, включая Bisection. Для большинства типов чисел с плавающей запятой деление пополам происходит с использованием соглашений о хранении значений с плавающей запятой, в результате чего определяется точный нуль или ограничивающий интервал минимальной длины, которой можно достигнуть при вычислениях с плавающей запятой. К другим методам относятся A42, AlefeldPotraShi, Roots.Brent, Roots.Chandrapatlu, Roots.ITP, Roots.Ridders и разновидности метода FalsePosition. Ограничивающий метод по умолчанию для базовых типов с плавающей запятой — Bisection, так как он более устойчив в случае некоторых входных данных, однако A42 и AlefeldPotraShi обычно сходятся за меньшее число итераций и поэтому более эффективны.
Несколько методов без вычисления производных. Задаются посредством методов Order0, Order1 (метод секущей), Order2 (метод Стеффенсена), Order5, Order8 и Order16. Число приблизительно указывает порядок сходимости. Метод по умолчанию — Order0. Он же является наиболее надежным, так как сводится к ограничивающему методу при достижении ограничения. Методы более высокого порядка могут обеспечивать сходимость более высокого порядка (более быструю), однако для получения результата не всегда требуют меньше вызовов функций, чем Order1 или Order2. Методы Roots.Order1B и Roots.Order2B — это сверхлинейные и квадратически сходящиеся методы, не зависящие от кратности нуля.
Методы, требующие одной или нескольких производных: Классические варианты — Roots.Newton и Roots.Halley. Помимо них, есть методы Roots.QuadraticInverse, Roots.ChebyshevLike и Roots.SuperHalley. Roots.Schroder предоставляет квадратичный метод, подобный методу Ньютона, который не зависит от кратности нуля. Методы Roots.ThukralXB, где X=2, 3, 4 или 5, также не зависят от кратности. X означает количество определяемых производных. Методы Roots.LithBoonkkampIJzerman{S,D} запоминают S шагов и используют D производных.

Базовое использование

Рассмотрим полиномиальную функцию . Корни (или нули) этого многочлена можно было бы определить с помощью функции roots из пакета Polynomials. Однако для нахождения значений даже такой функции применяется численный метод, так как решение с радикалами отсутствует. То есть даже в случае с многочленами для решения требуются нелинейные методы нахождения корней. (Хотя полиномиальные методы нахождения корней могут использовать некоторые свойства, недоступные для универсальных нелинейных функций.)

Пакет Roots предоставляет различные алгоритмы для решения этой задачи. В этом кратком обзоре представлены лишь стандартные алгоритмы.

Для функции на простом графике видно, что один нуль находится где-то между и , а еще два — возле .

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

julia> using Roots

julia> f(x) =  x^5 - x + 1/2
f (generic function with 1 method)

julia> find_zero(f, (-1.2,  -1)) ≈ -1.0983313019186336
true

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

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

julia> find_zero(f,  0.6) ≈ 0.550606579334135
true

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

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

julia> find_zeros(f, -5,  5)
3-element Vector{Float64}:
 -1.0983313019186334
  0.550606579334135
  0.7690997031778959