Глобальная и кусочно-полиномиальная интерполяция
Пусть функция f(x) интерполируется на отрезке [a, b]. Метод решения этой задачи с помощью единого многочлена для всего отрезка называют глобальной полиномиальной интерполяцией. В вычислительной практике такой поход применяется редко в силу различных причин. Одна из причин v необходимо задать стратегию выбора узлов при интерполяции функции f многочленами все возрастающей степени n.
Теорема Фабера. Какова бы ни была стратегия выбора узлов интерполяции, найдется непрерывная на отрезке [a, b] функция f, для которой при . Таким образом, теорема Фабера отрицает существование единой для всех непрерывных функций стратегии выбора узлов. Проиллюстрируем сказанное примером.
Предположим, что выбираем равноотстоящие узлы, то есть , i = 0, 1,...n, где . Покажем, что для функции Рунге такая стратегия является неудачной.
ПРИМЕР 1. Глобальная интерполяция функции Рунге (рис. 10.1).
% Интерполяция функции Рунге
% Введём функцию Рунге
f = inline('1./(1+25*x.^2)');
% Вычислим таблицу значений
x = linspace(-1, 1, 10);
y = f(x);
% Проинтерполируем функцию Рунге многочленами Лагранжа
p = polyfit(x, y, 10);
xx = linspace(-1, 1, 100);
yy = polyval(p, xx);
axes('NextPlot', 'Add');
% Покажем, что глобальная аппроксимация плохо работает для функции Рунге
plot(x, y);
plot(xx, yy, 'Color', 'r');
Рис. 10.1 - Глобальная интерполяция функции Рунге
% С увеличением узлов сетки, ситуация только ухудщается
% Вычислим таблицу значений. 20 узлов сетки
x = linspace(-1, 1, 20);
y = f(x);
% Проинтерполируем функцию Рунге многочленами Лагранжа
p = polyfit(x, y, 20);
xx = linspace(-1, 1, 100);
yy = polyval(p, xx);
figure
axes('NextPlot', 'Add');
% Покажем, что глобальная аппроксимация плохо работает для функции Рунге
plot(x, y);
plot(xx, yy, 'Color', 'r');
Рис. 10.2 - интерполяция функции Рунге многочленами Лагранжа
На практике чаще используют кусочно-полиномиальную интерполяцию: исходный отрезок разбивается на части и на каждом отрезке малой длины исходная функция заменяется многочленом невысокой степени. Система Mathcad предоставляет возможность аппроксимации двумя важными классами функций: кусочно-линейной и сплайнами.
При кусочно-линейной интерполяции узловые точки соединяются отрезками прямых, то есть через каждые две точки , проводится полином первой степени.
ПРИМЕР 2. Кусочно-линейная интерполяция функции Рунге.
% Кусочно-линейная интерполяция функции Рунге
% Введём функцию Рунге
f = inline('1./(1+25*x.^2)');
% Вычислим таблицу значений
x = linspace(-1, 1, 10);
y = f(x);
% Начертим график кусочно-линейной аппроксимации (рис. 10.3)
plot(x, y);
Рис. 10.3 - график кусочно-линейной аппроксимации
Как видно из приведенного примера этот способ приближения также имеет недостаток: в точках «стыка» двух соседних многочленов производная, как правило, имеет разрыв.
Если исходная функция была гладкой и требуется, чтобы и аппроксимирующая функция была гладкой, то кусочно-полиномиальная интерполяция неприемлема. В этом случае применяют сплайны v специальным образом построенные гладкие кусочно-многочленные функции.
Интерполяция сплайнами
Пусть отрезок [a,b] разбит точками на n частичных отрезков . Сплайном степени m называется функция , обладающая следующими свойствами:
1) функция непрерывна на отрезке [a,b] вместе со своими производными до некоторого порядка p;
2) на каждом частичном отрезке функция совпадает с некоторым алгебраическим многочленом степени m.
Разность m-p между степенью сплайна и наивысшим порядком непрерывной на отрезке [a ,b] производной называют дефектом сплайна. Кусочно-линейная функция является сплайном первой степени с дефектом, равным единице. Действительно, на отрезке [a ,b] сама функция (нулевая производная) непрерывна. В то же время на каждом частичном отрезке совпадает с некоторым многочленом первой степени.
ПРИМЕР 3. Построение параболического сплайна.
Пусть дан фрагмент таблицы значений функции:
x | -1 | ||
y | 1.5 | 0.5 | 2.5 |
Требуется построить параболический сплайн дефекта 1.
Так как строится сплайн , то он будет представлен двумя полиномами 2-ой степени:
.
Функция должна удовлетворять условиям:
- это есть условие интерполяции;
- это есть условие непрерывности первой производной.
Таким образом, получили 5 условий для нахождения 6-сти неизвестных. Два условия дополнительно накладывают на сплайн в граничных точках.
Возьмем, например дополнительное граничное условие следующего вида .
Тогда получим систему уравнений относительно неизвестных коэффициентов :
Эта система легко решается:
, , , , , .
Таким образом:
.
Наиболее широкое распространение получили сплайны 3 степени (кубические сплайны) с дефектом равным 1 или 2. Система для осуществления сплайн-интерполяции кубическими полиномами предусматривает несколько встроенных функций. Одна из них рассмотрена в примере.
ПРИМЕР 4 . Построение сплайн-интерполяции (рис. 10.4).
% Построить интерполяцию сплайнами функции Рунге
% Введём функцию Рунге
f = inline('1./(1+25*x.^2)');
% Вычислим таблицу значений
x = linspace(-1, 1, 10);
y = f(x);
% Вычислим сплайн-интерполяцию
xx = linspace(-1, 1, 100);
yy = spline(x, y, xx);
% Начертим графики
axes('NextPlot', 'Add');
plot(x, y, 'LineWidth', 2);
% Красным на графике - аппроксимация, жирным - исходная функция.
plot(xx, yy, 'Color', 'r');
Рис. 10.4 - построение сплайн-интерполяции
Погрешность приближения кубическими сплайнами.
Пусть функция f имеет на отрезке [a,b] непрерывную производную четвертого порядка и .
Тогда для интерполяционного кубического сплайна справедлива оценка погрешности: .
Задание для самостоятельной работы
1. Функция y = f(x) задана таблицей своих значений.
X | 0.2 | 0.4 | 0.6 | 0.8 | ||
Y | 0.75 | 1.1 | 1.35 | 1.25 | 1.05 | 0.8 |
Предложить способы интерполирования для нахождения значений функции в точках 0.24, 0.5, 0.96.
2. Функция y = f(x) задана таблицей своих значений.
X | |||
Y |
Построить интерполяционный кубический сплайн с граничными условиями , .
3. Проинтерполировать функцию задачи 2 методом кусочно-линейной интерполяции и построить график исходной функции и найденных многочленов.
Вопросы
1. Объясните разницу между глобальной и кусочно-полиномиальной интерполяцией. Почему на практике чаще используется кусочно-полиномиальная интерполяция.
2. Дайте определение интерполяционного сплайна m-ой степени.
3. Что такое дефект сплайна.
4. Запишите формулу сплайна первой степени с дефектом единица.
Численное дифференцирование
Производная функции есть предел отношения приращения функции к приращению независимой переменной при стремлении к нулю приращения независимой переменной:
.
При численном нахождении производной заменим отношение бесконечно малых приращений функций и аргумента отношением конечных разностей. Очевидно, что чем меньше будет приращение аргумента, тем точнее численное значение производной.