Алгоритм Де Бура

В математической области анализа численного алгоритм де Бура ^[1] представляет собой полиномиальный и численно стабильный алгоритм для оценки сплайновых кривых в форме B-сплайна . Это обобщение алгоритма де Кастельжо для кривых Безье . Алгоритм был разработан немецко-американским математиком Карлом Р. де Буром . Были созданы упрощенные, потенциально более быстрые варианты алгоритма де Бура, но они имеют сравнительно меньшую стабильность. ^[2]^[3]

Введение

Общее введение в B-сплайны дано в основной статье . Здесь мы обсуждаем алгоритм де Бура, эффективную и численно стабильную схему для оценки сплайновой кривой. $\mathbf {S} (x)$ на позиции $x$ . Кривая строится из суммы B-сплайн-функций. $B_{i,p}(x)$ умноженный на потенциально векторные константы $\mathbf {c} _{i}$ , называемые контрольными точками, $\mathbf {S} (x)=\sum _{i}\mathbf {c} _{i}B_{i,p}(x).$ B-сплайны порядка $p+1$ — связные кусочно-полиномиальные функции степени $p$ определяется по сетке узлов ${t_{0},\dots ,t_{i},\dots ,t_{m}}$ (далее мы всегда используем индексы, отсчитываемые от нуля). Алгоритм Де Бура использует O (p ²) + O (p) операций для оценки сплайновой кривой. Примечание: основная статья о B-сплайнах и классические публикации. ^[1] используйте другое обозначение: B-сплайн индексируется как $B_{i,n}(x)$ с $n=p+1$ .

Местная поддержка

B-сплайны имеют локальную поддержку, что означает, что полиномы положительны только в конечной области и равны нулю в других местах. Формула рекурсии Кокса-де Бура ^[4] показывает это: $B_{i,0}(x):={\begin{cases}1&{\text{if }}\quad t_{i}\leq x<t_{i+1}\\0&{\text{otherwise}}\end{cases}}$ $B_{i,p}(x):={\frac {x-t_{i}}{t_{i+p}-t_{i}}}B_{i,p-1}(x)+{\frac {t_{i+p+1}-x}{t_{i+p+1}-t_{i+1}}}B_{i+1,p-1}(x).$

Пусть индекс $k$ определить интервал узла, который содержит позицию, $x\in [t_{k},t_{k+1})$ . В формуле рекурсии мы видим, что только B-сплайны с $i=k-p,\dots ,k$ отличны от нуля для этого узла. Таким образом, сумма уменьшается до: $\mathbf {S} (x)=\sum _{i=k-p}^{k}\mathbf {c} _{i}B_{i,p}(x).$

Это следует из $i\geq 0$ что $k\geq p$ . Точно так же мы видим в рекурсии, что самое высокое запрошенное местоположение узла находится по индексу $k+1+p$ . Это означает, что любой интервал узла $[t_{k},t_{k+1})$ который фактически используется, должен иметь как минимум $p$ дополнительные узлы до и после. В компьютерной программе это обычно достигается путем повторения первого и последнего использованного местоположения узла. $p$ раз. Например, для $p=3$ и реальные места узлов $(0,1,2)$ , можно было бы дополнить вектор узла до $(0,0,0,0,1,2,2,2,2)$ .

Алгоритм

Имея эти определения, мы теперь можем описать алгоритм де Бура. Алгоритм не вычисляет функции B-сплайна. $B_{i,p}(x)$ напрямую. Вместо этого он оценивает $\mathbf {S} (x)$ через эквивалентную рекурсивную формулу.

Позволять $\mathbf {d} _{i,r}$ быть новыми контрольными точками с $\mathbf {d} _{i,0}:=\mathbf {c} _{i}$ для $i=k-p,\dots ,k$ . Для $r=1,\dots ,p$ применяется следующая рекурсия: $\mathbf {d} _{i,r}=(1-\alpha _{i,r})\mathbf {d} _{i-1,r-1}+\alpha _{i,r}\mathbf {d} _{i,r-1};\quad i=k-p+r,\dots ,k$ $\alpha _{i,r}={\frac {x-t_{i}}{t_{i+1+p-r}-t_{i}}}.$

После завершения итераций мы имеем $\mathbf {S} (x)=\mathbf {d} _{k,p}$ , это означает, что $\mathbf {d} _{k,p}$ является желаемым результатом.

Алгоритм Де Бура более эффективен, чем явный расчет B-сплайнов. $B_{i,p}(x)$ с рекурсивной формулой Кокса-де Бура, поскольку она не вычисляет члены, которые гарантированно умножаются на ноль.

Оптимизации

Приведенный выше алгоритм не оптимизирован для реализации на компьютере. Для этого требуется память $(p+1)+p+\dots +1=(p+1)(p+2)/2$ временные контрольные точки $\mathbf {d} _{i,r}$ . Каждая временная контрольная точка записывается ровно один раз и читается дважды. Обратив итерацию $i$ (считая вниз, а не вверх), мы можем запустить алгоритм с памятью всего за $p+1$ временные контрольные точки, позволяя $\mathbf {d} _{i,r}$ повторно использовать память для $\mathbf {d} _{i,r-1}$ . Аналогично, существует только одно значение $\alpha$ используется на каждом этапе, поэтому мы также можем повторно использовать память.

Кроме того, удобнее использовать индекс, отсчитываемый от нуля. $j=0,\dots ,p$ для временных контрольных пунктов. Связь с предыдущим индексом $i=j+k-p$ . Таким образом мы получаем улучшенный алгоритм:

Позволять $\mathbf {d} _{j}:=\mathbf {c} _{j+k-p}$ для $j=0,\dots ,p$ . Итерация для $r=1,\dots ,p$ : $\mathbf {d} _{j}:=(1-\alpha _{j})\mathbf {d} _{j-1}+\alpha _{j}\mathbf {d} _{j};\quad j=p,\dots ,r\quad$ $\alpha _{j}:={\frac {x-t_{j+k-p}}{t_{j+1+k-r}-t_{j+k-p}}}.$ Обратите внимание, что $j$ необходимо отсчитывать. После завершения итераций результат: $\mathbf {S} (x)=\mathbf {d} _{p}$ .

Пример реализации

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

def deBoor(k: int, x: int, t, c, p: int):
    """Evaluates S(x).

    Arguments
    ---------
    k: Index of knot interval that contains x.
    x: Position.
    t: Array of knot positions, needs to be padded as described above.
    c: Array of control points.
    p: Degree of B-spline.
    """
    d = [c[j + k - p] for j in range(0, p + 1)] 

    for r in range(1, p + 1):
        for j in range(p, r - 1, -1):
            alpha = (x - t[j + k - p]) / (t[j + 1 + k - r] - t[j + k - p]) 
            d[j] = (1.0 - alpha) * d[j - 1] + alpha * d[j]

    return d[p]

См. также

Внешние ссылки

Компьютерный код

PPPACK : содержит множество сплайн-алгоритмов на Фортране.
Научная библиотека GNU : C-библиотека, содержит подбиблиотеку для сплайнов, перенесенных из PPPACK.
SciPy : библиотека Python, содержит подбиблиотеку scipy.interpolate со сплайн-функциями на основе FITPACK.
TinySpline : C-библиотека для сплайнов с оболочкой C++ и привязками для C#, Java, Lua, PHP, Python и Ruby.
Einspline : C-библиотека для сплайнов в 1, 2 и 3 измерениях с оболочками Фортрана.

Ссылки

^ Jump up to: ^а ^б К. де Бур [1971], «Пакет подпрограмм для вычислений с помощью B-сплайнов», Techn.Rep. LA-4728-MS, Лос-Аламосская научная лаборатория, Лос-Аламос, Нью-Мексико; п. 109, 121.
^ Ли, ETY (декабрь 1982 г.). «Упрощенная процедура расчета B-сплайна». Вычисление . 29 (4). Спрингер-Верлаг: 365–371. дои : 10.1007/BF02246763 . S2CID 2407104 .
^ Ли, ETY (1986). «Комментарии к некоторым алгоритмам B-сплайна». Вычисление . 36 (3). Спрингер-Верлаг: 229–238. дои : 10.1007/BF02240069 . S2CID 7003455 .
^ К. де Бур, с. 90

Цитируемые работы

Карл де Бур (2003). Практическое руководство по сплайнам, переработанное издание . Спрингер-Верлаг. ISBN 0-387-95366-3 .

[de_boor_paper-1] Jump up to: ^а ^б К. де Бур [1971], «Пакет подпрограмм для вычислений с помощью B-сплайнов», Techn.Rep. LA-4728-MS, Лос-Аламосская научная лаборатория, Лос-Аламос, Нью-Мексико; п. 109, 121.

[2] Ли, ETY (декабрь 1982 г.). «Упрощенная процедура расчета B-сплайна». Вычисление . 29 (4). Спрингер-Верлаг: 365–371. дои : 10.1007/BF02246763 . S2CID 2407104 .

[3] Ли, ETY (1986). «Комментарии к некоторым алгоритмам B-сплайна». Вычисление . 36 (3). Спрингер-Верлаг: 229–238. дои : 10.1007/BF02240069 . S2CID 7003455 .

[4] К. де Бур, с. 90

[1]

[2]

[3]

[4]