4.3 Локальная интерполяция
В алгебраическом интерполировании при росте числа узлов увеличивается, как правило, степень интерполяционного многочлена. Кроме того, когда интерполирование выполняется для функций, не являющихся достаточно гладкими, интерполирование высокого порядка нецелесообразно. В этих случаях лучший результат может дать локальная (кусочная) интерполяция.
В отличие от глобальной интерполяции, когда интерполирующая функция g(x) строится сразу для всего рассматриваемого интервала изменения x, при локальной интерполяции такая функция ищется отдельно для разных частей этого интервала, а затем как бы сшивается из соответствующих отдельных кусков. При локальной интерполяции должно сохраняться основное условие интерполирования: в узлах xi значения интерполирующей функции должны быть равны значениям данной функции yi = f(xi).
Простейшим и часто используемым видом локальной интерполяции является линейная (кусочно-линейная) интерполяция. Она состоит в том, что заданные точки (xk, yk) (k = 1, 2 ,…, n) соединяются прямолинейными отрезками, и функция f(x) приближается ломаной с вершинами в этих точках.
Уравнения каждого отрезка ломаной в общем случае разные. Поскольку имеется n – 1 интервалов (xk, xk+1), то для каждого из них в качестве уравнения интерполяционного многочлена используется уравнение прямой, проходящей через две заданные точки. В частности для k-го интервала можно написать уравнение прямой, проходящей через две точки (xk, yk) и (xk+1, yk+1) в виде
,
xkxxk+1.
(4.5)
Следовательно, при использовании линейной интерполяции сначала нужно определить интервал, в который попадает значение аргумента, а затем подставить его в формулу (4.5) и найти приближенное значение в данной точке.
Приведем реализацию этого алгоритма в среде Maple для таблично заданной функции из предыдущего примера.
> restart;
> X:=[0.5,2,3,4,5];Y:=[0.3,1,2,2,1];n:=5: # исходные данные
Программируем формулу (4.5), причем величину y из этой формулы определяем как функцию (с именем h) не только аргумента x, но и номера интервала k:
> h:=(x,k)->(Y[k+1]-Y[k])*(x-X[k])/(X[k+1]-X[k])+Y[k];
Интерполирующую функцию на всем числовом промежутке задаем кусочным образом:
> w:=x->piecewise(x<X[2],h(x,1),x<X[3],h(x,2), x<X[4],h(x,3),h(x,4));
> simplify(w(x));
Итак, искомая интерполирующая функция есть w(x). Далее – графический вывод результата интерполирования:
> g3:=plot(w(x),x=X[1]..X[5]):
> g2:=plot([[X[i],Y[i]]$i=1..5],x=0..5,style=POINT,
symbol=CIRCLE,color=black);
> plots[display](g3,g2);
Рис. 4.2.
Рассмотрим теперь случай квадратичной локальной интерполяции. В качестве интерполяционной функции на отрезке [xk, xk+2] принимается квадратичный трехчлен. Такую интерполяцию называют также параболической. Если отвлечься от того, что отрезок [xk, xk+2] является подмножеством отрезка [x1,xn], и рассмотреть его отдельно, то квадратичную интерполяцию можно трактовать как глобальную с тремя узлами. Интерполирующую функцию на этом отрезке можно выразить явно как многочлен Лагранжа L2(x):
Ниже приводится фрагмент программы на языке Maple, реализующий кусочно-квадратичную интерполяцию (подразумевается, что переменные X, Y, n, g2, введенные ранее, остались в оперативной памяти).
> h:=(x,k)->Y[k]*(x-X[k+1])*(x-X[k+2])/(X[k]-X[k+1])/(X[k]-X[k+2])+Y[k+1]*(x-X[k])*(x-X[k+2])/(X[k+1]-X[k])/(X[k+1]-X[k+2])+Y[k+2]*(x-X[k])*(x-X[k+1])/(X[k+2]-X[k])/(X[k+2]-X[k+1]);
> w:=x->piecewise( x<X[3],h(x,1),x<=X[5],h(x,3));
> simplify(w(x));
> g3:=plot(w(x),x=X[1]..X[5]):
> plots[display](g3,g2);
Рис. 4.3.
Если присмотреться, из последнего графика можно увидеть главный недостаток данного способа кусочной интерполяции: в граничных точках отдельных участков, где осуществляется интерполирование, производная аппроксимирующей функции терпит разрыв. От этого недостатка свободны другие виды локальной интерполяции, в частности, сплайн-интерполяция и кусочная интерполяция многочленами Эрмита.