Математика / Интерполяция

Периодический кубический сплайн по циклической системе

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

Опубликовано: Обновлено:

Формула

$$M_{i-1}h_{i-1}+2M_i(h_{i-1}+h_i)+M_{i+1}h_i=6\left(\frac{y_{i+1}-y_i}{h_i}-\frac{y_i-y_{i-1}}{h_{i-1}}\right),\quad M_0=M_n$$
Схема Один период кубического сплайна
00.511.522.53-0.500.51yxпериодический сплайн

Узлы периода: (0;0), (1;1), (2;0), (3;0); крайние значения совпадают, и кривая замыкается.

Обозначения

$x_i$
узел сетки, единица аргумента
$y_i$
значение функции в узле, единица функции
$h_i=x_{i+1}-x_i$
длина i-го интервала, единица аргумента
$M_i=S″(x_i)$
вторая производная сплайна в узле, единица функции / аргумент²
$S_i(x)$
кубический участок сплайна на [x_i,x_{i+1}], единица функции
$n$
число интервалов одного периода, —

Условия применения

Таблица должна содержать один полный период с одинаковыми значениями на его концах: x_n=x_0+T и y_n=y_0. Сетка может быть неравномерной, но узлы должны строго возрастать. Периодическое продолжение требует согласованного поведения функции на стыке; для гладкого периодического сигнала совпадение значений, первой и второй производных на концах естественно. При равномерном шаге система для M упрощается до циклической трёхдиагональной матрицы.

Ограничения

Замыкание только по значениям без периодических производных создаёт излом или скачок кривизны на стыке. Если измеренный последний отсчёт отличается от первого из-за шума, его нельзя молча считать точным периодическим концом: нужно согласовать данные или применить сглаживание. Сплайн интерполирует каждую точку и поэтому может повторять шум. Выбранная периодичность также не определяет корректность периода физической модели; ошибочный T даст гладкую, но неверно фазированную кривую.

Подробное объяснение

Для каждого промежутка кубический участок определяется концами и значениями M_i и M_{i+1}. После дифференцирования и требования непрерывности первой производной в узлах получается линейное уравнение: соседние кривизны умножаются на длины смежных интервалов, а правая часть равна шестиразности наклонов соседних хорд. В периодическом случае соседями первого узла становятся узлы последнего интервала, поэтому матрица имеет циклические угловые элементы, а M_n отождествляется с M_0. Решив систему, участок на [x_i,x_{i+1}] можно вычислить по формуле через две линейные базисные части и два кубических поправочных слагаемых.

Для ручного счёта достаточно сначала выписать наклоны секущих d_i=(y_{i+1}−y_i)/h_i, затем правые части 6(d_i−d_{i−1}) и циклические уравнения. В равномерной сетке коэффициенты становятся 1,4,1, если уравнение разделить на h. В примере с четырьмя узлами 0,1,2,3 и значениями 0,1,0,0 правая часть равна (6,−12,6); условие M_3=M_0 даёт решение (2,−4,2,2). Подстановка этих значений в кубический участок между 1 и 2 даёт середину 0,625.

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

Как пользоваться формулой

  1. Проверьте, что данные покрывают ровно один заявленный период и первый и последний отсчёты равны.
  2. Вычислите длины интервалов h_i и наклоны секущих d_i между соседними узлами.
  3. Составьте циклические уравнения для M_i, учитывая M_n=M_0 и соседство последнего интервала с первым.
  4. Решите систему и подставьте две соседние кривизны в формулу участка, содержащего целевой x.
  5. Проверьте значение в узле и согласованность стыка по наклону и второй производной.

Историческая справка

Кубические сплайны выросли из более ранней теории кусочно-полиномиальной аппроксимации. В 1946 году Исаак Якоб Шёнберг ввёл термин spline function в работах об аппроксимации равномерных данных аналитическими функциями; позднее теория получила систематическое развитие в численном анализе. Периодические сплайны стали отдельной конструкцией для циклических данных: статья П. В. Галкина «О возможности периодической сплайн-интерполяции» 1970 года прямо рассматривает разрешимость классической задачи периодической интерполяции.

Современный вычислительный способ задаёт периодические условия для производных на обоих концах и решает систему для узловых кривизн. Учебник Парвиза Моина *Fundamentals of Engineering Numerical Analysis* отдельно описывает периодическое условие и объясняет, что оно предназначено для одного периода известного сигнала. Поэтому конкретная формула циклической системы — часть общей теории кубических сплайнов, а не открытие, которое корректно приписать одному автору. В статье рассматривается именно интерполяционная версия: она проходит через данные, а не сглаживает их.

Историческая линия формулы

Периодические сплайны развивались внутри теории сплайнов, термин для которой Исаак Шёнберг ввёл в 1946 году. Галкину принадлежит опубликованная в 1970 году работа о разрешимости периодической сплайн-интерполяции. Циклическое условие в современном вычислительном алгоритме связывает конец таблицы с её началом.

Пример

Дан период [0;3] с узлами 0,1,2,3 и значениями 0,1,0,0; шаг равен 1, значения на концах совпадают. Наклоны секущих равны 1, −1 и 0. Циклическая система для уникальных кривизн M₀,M₁,M₂: 4M₀+M₁+M₂=6; M₀+4M₁+M₂=−12; M₁+4M₂+M₀=6, причём M₃=M₀. Решение: (M₀,M₁,M₂,M₃)=(2,−4,2,2). На участке [1;2] при x=1,5 применяем S(x)=M₁(2−x)³/6+M₂(x−1)³/6+(y₁−M₁/6)(2−x)+(y₂−M₂/6)(x−1). Получаем −4·0,5³/6+2·0,5³/6+(1+4/6)·0,5+(0−2/6)·0,5=0,625. Ответ: 0,625 единицы функции. Проверка: участок принимает y₁=1 и y₂=0 на концах, а производная и кривизна стыкуются с соседними периодическими сегментами.

Частая ошибка

Не заменяйте циклическое условие M_n=M_0 двумя нулевыми кривизнами: так получится натуральная, а не периодическая граница. Следите за циклическим правым соседом в первом и последнем уравнениях. Ещё одна ошибка — забыть, что периодическая таблица повторяет конечное значение, но конечный узел всё равно задаёт отдельный интервал для замыкания. Если первый и последний y различны, сначала выясните, являются ли они шумными измерениями одного значения или таблица не покрывает ровно период.

Практика

Задачи с решением

Производная на периодическом стыке

Условие. Для периода [0;3] с шагом 1 даны значения 0,1,0,0 и кривизны 2,−4,2,2. Проверьте равенство S′(0) и S′(3).

Решение. В начале первого участка производная равна секущей минус (2M₀+M₁)/6. Это 1−(4−4)/6=1. В конце последнего участка производная равна секущей плюс (M₂+2M₃)/6. Это 0+(2+4)/6=1.

Ответ. S′(0)=S′(3)=1 в единицах функции на единицу аргумента.

Замкнутый четырёхшаговый профиль

Условие. На узлах 0,1,2,3,4 заданы значения 0,1,0,−1,0. Найдите периодическое значение при x=0,5.

Решение. Секущие имеют наклоны 1,−1,−1,1. Циклическая система для M₀…M₃ даёт (0,−3,0,3), а M₄=M₀. На первом участке S(0,5)=−3·0,5³/6+(1−(−3)/6)·0,5=−0,0625+0,75=0,6875=11/16.

Ответ. 0,6875

Дополнительные источники

  • https://assets.cambridge.org/97805218/84327/excerpt/9780521884327_excerpt.pdf
  • https://www.mathnet.ru/php/archive.phtml?jrnid=mzm&option_lang=eng&paperid=7004&wshow=paper
  • https://pzs.dstu.dp.ua/DataMining/spline/bibl/Spline_Functions_Basic_Theory.pdf

Связанные формулы

Математика

Натуральный кубический сплайн: уравнения для вторых производных

$h_{i-1}M_{i-1}+2(h_{i-1}+h_i)M_i+h_iM_{i+1}=6(\delta_i-\delta_{i-1})$

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

Математика

Кубический сплайн с заданными наклонами на концах

$\begin{aligned}2h_0M_0+h_0M_1&=6(\delta_0-m_L),\\h_{n-1}M_{n-1}+2h_{n-1}M_n&=6(m_R-\delta_{n-1}).\end{aligned}$

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

Математика

Кубическая интерполяция Эрмита по значениям и наклонам

$H(x)=h_{00}(s)y_0+h_{10}(s)\ell m_0+h_{01}(s)y_1+h_{11}(s)\ell m_1,\quad s=\frac{x-x_0}{\ell},\quad \ell=x_1-x_0$

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

Математика

Линейная интерполяция: значение между двумя узлами

$P(x)=y_0+\frac{y_1-y_0}{x_1-x_0}(x-x_0)=(1-t)y_0+ty_1,\quad t=\frac{x-x_0}{x_1-x_0}$

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