Интерполирование сплайнами

Интерполяционные формулы Лагранжа, Ньютона и Стирлинга и др. при использовании большого числа узлов интерполяции на всем отрезке [a, b] часто приводят к плохому приближению из-за накопления погрешностей в процессе вычислений [2]. Кроме того, из-за расходимости процесса интерполяции увеличение числа узлов не обязательно приводит к повышению точности. Для снижения погрешностей весь отрезок [a, b] разбивается на частичные отрезки и на каждом из них функцию f(x) заменяют приближенно полиномом невысокой степени. Это называется кусочно-полиномиальной интерполяцией.

Один из способов интерполирования на всем отрезке [a, b] является интерполирование сплайнами.

Сплайном называется кусочно-полиномиальная функция, определенная на отрезке [a, b] и имеющая на этом отрезке некоторое количество непрерывных производных. Преимущества интерполяции сплайнами по сравнению с обычными методами интерполяции – в сходимости и устойчивости вычислительного процесса.

На практике предпочитают для получения сходящегося интерполяционного процесса использовать кусочно-полиноминальную интерполяцию. Для этого отрезок [а,b] разделяют на части точками и интерполируют не при помощи одного многочлена на всем отрезке [а,b] а для каждого частичного отрезка избирают свои узлы и строят свой интерполирующий многочлен невысокой степени. Однако этот способ приближения имеет недостаток: в точках «стыка» двух соседних многочленов производная, как правило, имеет разрыв. Часто это обстоятельство не играет большой роли. Вместе с тем нередко требуется, чтобы аппроксимирующая функция была гладкой и тогда простейшая кусочно-гладкая интерполяция становиться неприемлемой. Естественна потребность в наличии аппроксимирующих функций, которые сочетали бы в себе локальную простоту многочлена невысокой степени и глобальную на всем отрезке [а,b] гладкость, привела к появлению так называемых гладких сплайн-функций или сплайнов – специальным образом построенных гладких кусочно-многочленных функций.

Дадим строгое определение сплайна. Пусть отрезок [а,b] разбит точками a=x₀<x₁<x_n=b на n частичных отрезков [x_i-1,x_i]. Сплайном S_m(x) степени m называется определенная на [a,b] функция, принадлежащая классу C^p[a,b] p

раз непрерывно дифференцируемых функций, такая, что на каждом частичном промежутке [x_i-1,x_i] i=1,…, n – это многочлен m-ой степени.

Разность d=m-p между степенью сплайна и наивысшим порядком непрерывной на отрезке [а,b] производной называется дефектом сплайна.

Пусть функция y=f(x) задана таблицей своих значений (1). Сплайн S_m(x) называется интерполяционным, если S_m(x_i)= y(x_i) для всех i= 0, 1,…, n. Значение s_i= S¢_m(x_i) называется наклоном сплайна в точке x_i.

Простейший пример сплайна дает непрерывная кусочно-линейная функция, являющаяся сплайном первой степени (линейным сплайном) с дефектом, равным единице. Действительно на отрезке [а,b] сама функция (нулевая производная) непрерывна. В то же время на каждом частичном отрезке S₁(x) совпадает с некоторым многочленом первой степени.

Наиболее широкое распространение на практике получили сплайны S₃(x) третьей степени (кубические сплайны) с дефектом, равным 1 или 2. Такие сплайны на каждом из частичных отрезков [x_i-1,x_i] совпадают с кубическим многочленом. Дадим определение кубического сплайна.

Кубическим сплайном дефекта 1, интерполирующим на отрезке [а,b] данную функцию f(x), называется функция

g(x)={g_i(x)=a_i+b_i(x-x_i-1)+c_i(x-x_i-1)²+d _i(x-x_i-1)³ при xÎ[x_i-1,x_i] и i=1, …, n},

удовлетворяющая следующим условиям:

1. условиям интерполяции - g(x_i)=f_i, i=0, …, n;

2. двойной непрерывной дифференцируемостью - g(x)ÎC²[a,b];

3. краевым условиям.

Краевые условия могут быть весьма разнообразными. Можно задавать в граничных точках конкретные наклоны графика сплайна:

g^¢(x₀) = f ^¢(a), g^¢(x_n) = f ¢(b)

или конкретную кривизну:

g^¢¢(x₀) = f ^¢¢(a), g^¢¢(x_n) = f ¢¢(b).

Часто задают условия “свободного провисания”:

g^¢¢(x₀) = 0, g^¢¢(x_n) = 0.

Эти условия определяют так называемый естественный кубический сплайн (или чертежный сплайн).

Рассмотрим построение естественного кубического сплайна. Для этого надо найти 4n его коэффициентов { a_i, b_i, c_i, d_i; i=1,…, n}. Для нахождения этих коэффициентов имеем следующие соотношения:

1. из условий интерполяции g₁(x₀)=f₀, g_i(x_i)=f_i, i=1,2,…, n; (уравнений (n+1))

2. непрерывность сплайна g_i-1(x_i-1)=g_i(x_i-1), i= 2, 3, …,n; (уравнений (n-1))

3. непрерывность производной g^¢_i-1(x_i-1)=g^¢_i(x_i-1), i= 2, …,n; (уравнений (n-1))

4. непрерывность второй производной g^¢¢_i-1(x_i-1)=g^¢¢_i(x_i-1), i= 2,…,n; (n-1)

5. краевые условия g^¢¢₁(x₀)=0, g^¢¢_n(x_n)=0. (два уравнения)

Для i-ой составляющей сплайна имеем:

g_i(x)=a_i+b_i(x-x_i)+c_i(x-x_i)²+d _i(x-x_i)³,

g^¢_i(x)=b_i+2c_i(x-x_i)+3d₁(x-x_i)²,

g^¢¢_i(x)=2c_i+6d_i(x-x_i).

Полагая для краткости узлы равноотстоящими: h=x_i-x_i-1, i=1,2,…, n, получим детализированную систему уравнений:

a₁-b₁h+c₁h²-d₁h³=f₀ , a_i=f_i, i=1,2,…, n;

a_i-1= a_i-b_ih+c_ih²-d₁h³, i= 2, 3, …,n;

b_i-1=b_i-2c_ih+3d_ih², i= 2, 3, …,n;

c_i-1=c_i-3d_ih, i= 2, 3, …,n;

c₁-3d₁h=0, с_n=0.

Всего получается 4n уравнений т.е. столько же сколько и неизвестных.

Теорема о степени приближения кубическим сплайном:

Если функция f(x) принадлежит классу четырежды дифференцируемых на отрезке [а,b] функций: f(x)ÎC⁴[a,b], то кубический сплайн g(x) удовлетворяет неравенствам:

max½f^(m)(x)-g^(m)(x)½£ Kh^4-m max½f⁽⁴⁾(x)½ при a£x£b и m=0, 1, 2;

т.е. ½f(x)-g(x)½= O(h⁴), ½f ¢(x)-g¢(x)½= O(h³), ½f ¢¢(x)-g¢¢(x)½= O(h²).

Замечание: Благодаря большей простоте записи естественные сплайны получили значительное распространение. Однако искусственное наложение условий f ^¢¢(a)=0 и f ^¢¢(b)=0 при интерполяции функций, которые этим условиям не удовлетворяют, приводит к значительной потере точности. Вместо четвертого порядка естественный сплайн обладает лишь вторым порядком точности. Если использование естественного сплайна не вызвано какими-либо специальными причинами, то следует, по-видимому, отказаться от него в пользу кубического сплайна с граничными условиями другого типа.

45.Метод наименьших квадратов.

При обработке результатов измерений часто возникает необходимость построить эмпирическую формулу, дающую аналитическое выражение функциональной зависимости, заданной таблицей. Использование для этого интерполяционных полиномов не всегда целесообразно, так как зачастую результаты измерений помимо основной зависимости отражают разного рода малые случайные ошибки, и точное совпадение в узлах интерполяции построенной функции с табличными значениями может даже исказить основную зависимость. В других случаях речь может идти о построении наиболее простой формулы, достаточно хорошо отражающей указанную зависимость. Пусть результаты измерений представлены совокупностью измеренных значений функции в соответствующих точках f(x_i), i=0,1,…,n. Требуется найти функцию j(x) в виде некоторой эмпирической формулы, которая аппроксимирует функцию f(x), так чтобы отклонения v_i=j(x_i) - f(x_i), i=0,1,…, nоказались бы в каком-то смысле наименьшими. Наиболее распространенным критерием малости отклонений является следующий критерий, положенный в основу способа наименьших квадратов. Требуют, чтобы оказалась минимальной сумма квадратов отклонений Будем искать аппроксимирующую функцию в виде полинома степени m: j(x) = b₀+b₁x+…+b_mx^m (b_m¹0).Задача ставится следующим образом: подобрать коэффициенты полинома: b₀, b₁,…, b_m, так чтобы сумма квадратов отклонений S достигала минимума.Заметим, что величина S является неотрицательной функцией переменных b₀, b₁,…, b_m и, следовательно, всегда имеет минимум. Если m>n, то существует бесконечное множество полиномов степени m, обеспечивающих S=0. Для m=n равенство S=0 обеспечивается единственным полиномом, являющимся интерполяционным полиномом для заданной функции. Обычное соотношение между степенью полинома и числом узлов: m<<n, в этом случае при любых значениях коэффициентов полинома S>0. Таким образом, задача сводится к исследованию на минимум функции многих переменных S(b₀, b₁,…, b_m). Необходимым условием минимальности S является равенство нулю производных по всем параметрам: ¶S/¶b₀=0, ¶S/¶b₁=0, …, ¶S/¶b_m=0.

¶S/¶b₀=2åⁿ_i=0(b₀+b₁x_i+…+b_mx_i^m –f(x_i)),

¶S/¶b₁=2åⁿ_i=0(b₀+b₁x_i+…+b_mx_i^m- f(x_i))x_i,

… … …

¶S/¶b_m=2åⁿ_i=0(b₀+b₁x _i …+b_mx_i^m-f(x_i))x_i^m.

Преобразование уравнений дает систему, назывемую системой нормальных уравнений:

(n+1)b₀+b₁åⁿ_i=0x_i+ b₂åⁿ_i=0x_i²+…+b_måⁿ_i=0x_i^m = åⁿ_i=0f(x_i)

b₀åⁿ_i=0x_i +b₁åⁿ_i=0(x_i)²+…+b_måⁿ_i=0(x_i)^m+1 = åⁿ_i=0x_i×f(x_i)

… … …

b₀åⁿ_i=0x_i^m+b₁åⁿ_i=0x_i^m+¹+…+b_måⁿ_i=0x_i²^m = åⁿ_i=0x_i^m×f(x_i)

Решая эту систему уравнений с симметричной матрицей, найдем коэффициенты полинома: b₀,b₁,…,b_m. Если среди узлов нет совпадающих, то определитель системы отличен от нуля, и так как значения функции не могут быть все одновременно нулями, то система имеет единственное решение. Если известна погрешность измерения единичного значения функции - e, то величина максимального отклонения d=maxïv_iï (i=0,1,…, n) может служить обоснованием правильности выбора степени аппроксимирующего полинома m. Если максимальное отклонение полинома от функции много больше погрешности измерения: d>>e, то это означает слишком грубую аппроксимацию, и степень m необходимо увеличить. В противном случае - d<<e, степень m завышена и ее следует уменьшить. Только когда d»e, можно считать, что степень аппроксимирующего полинома выбрана правильно.