Министерство образования и науки Украины
Национальный аэрокосмический университет им. Н.Е. Жуковского
«Харьковский авиационный институт»
Кафедра 304
Лабораторная работа №5
По курсу «Численные методы»
По теме: «Интерполирование для таблиц с постоянным шагом.
Численное дифференцирование. Обратное интерполирование.»
Выполнил:
студент 325 группы
Меняйлов Евгений
Проверила:
Яовая О.В.
___________________
Харьков 2009
Цель. Для табличной функции f(x), построить аналитическую функцию в виде интерполяционного многочлена Ньютона, Гаусса.
x |
0.0 |
0.1 |
0.2 |
0.3 |
0.4 |
0.5 |
0.6 |
0.7 |
0.8 |
0.9 |
f(x) |
0.96285 |
1.04482 |
1.12351 |
1.19871 |
1.27023 |
1.33790 |
1.40159 |
1.46117 |
1.51655 |
1.56760 |
Реализация в MathCad
1. Функция построения матрицы конечных разностей. |
2. Первая и вторая формулы Ньютона. |
3. Первая и вторая формулы Гаусса |
4. Проверка основного условия интерполирования |
5. Вычисление значений функции в заданных точках.
6. График интерполяционного многочлена.
__________________________________________________
Реализация в MatLab
%Функция EvalNuton.m вычисляет значения указанного полинома в заданном множестве точек:
function yy = EvalNuton(X,Y,xx,formNumber)
n=length(X)-1;
if (formNumber==1)
for i=1:length(xx);
yy(i)=Nuton1(xx(i),X,Y,n);
end;
else
for i=1:length(xx);
yy(i)=Nuton2(xx(i),X,Y,n+1);
end;
end;
return
------------------------------------------------------------------
%X,Y – табличные значения; x – точка интерполяции; n – количество узловых точек.
function res = Nuton1(x,X,Y,n);
res=Y(1);
% определение шага таблицы
h=X(2)-X(1);
for i=1:n
% построение матрицы конечных разностей
konrazn=diff(Y,i);
konrazn=konrazn(1);
% вычисление факториала
fact=prod(1:i);
% вычисление (x-x0)(x-x1)...(x-xn-1)
Mult=1;
for j=1:i
Mult=Mult*(x-X(j));
end;
res=res+(konrazn/(fact * h^i))*Mult;
end;
return
------------------------------------------------------------------
function res = Nuton2(x,X,Y,n);
res=Y(n);
% определение шага таблицы
h=X(2)-X(1);
%определение q
q=(x-X(n))/h;
for i=1:(n-1)
% построение матрицы конечных разностей
konrazn=diff(Y,i);
konrazn=konrazn(n-i);
% вычисление факториала
fact=prod(1:i);
% вычисление q (q-1)(q-2)...(q-(i-1))
Mult=1;
for j=0:(i-1)
Mult=Mult*(q+j);
end;
res=res+(konrazn/fact)*Mult;
end;
return
-------------------------------------------------
Построение графика интерполяционного многочлена.
xx=linspace(0,1,1000);
yy=EvalNuton(x,y,xx,1);
figure('Color','w');
hold on;
plot(xx,yy,'-b');
plot(x,y,'-.b*');
-------------------------------------------------------------------
Результаты.
>> x = [0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9];
>> y = [0.96285 1.04482 1.12351 1.19871 1.27023 1.33790 1.40159 1.46117 1.51655 1.5676];
>> EvalNuton(x,y,0.02,1)
ans =
0.97948373715354
>> EvalNuton(x,y,0.89,2)
ans =
1.56270742560460
function yy = EvalGauss(X,Y,xx,formNumber)
n=length(X)-1;
if (formNumber==1)
for i=1:length(xx);
yy(i)=Gauss1(xx(i),X,Y,n);
end;
else
for i=1:length(xx);
yy(i)=Gauss2(xx(i),X,Y,n+1);
end;
end;
return
---------------------------------------------------
function res = Gauss1(x, X, Y, n)
k = round(n/2);
res = Y(k);
% определение шага таблицы
h=X(2)-X(1);
for i=1:n
% построение матрицы конечных разностей
konrazn=diff(Y,i);
konrazn=konrazn(k-floor(i/2));
% вычисление факториала
fact=prod(1:i);
Mult=1;
for j=1:i
Mult=Mult*(x-X(k+(-1)^j*floor(j/2)));
end;
res=res+(konrazn/(fact * h^i))*Mult;
end;
return
---------------------------------------------------
function res = Gauss2(x, X, Y, n)
k = round(n/2);
res = Y(k);
% определение шага таблицы
h=X(2)-X(1);
for i=1:(n-1)
% построение матрицы конечных разностей
konrazn=diff(Y,i);
konrazn=konrazn(k-floor(i/2));
% вычисление факториала
fact=prod(1:i);
Mult=1;
for j=0:(i-1)
Mult=Mult*(x-X(k+(-1)^j*floor((i-j)/2)));
end;
res=res+(konrazn/(fact * h^i))*Mult;
end;
return
------------------------------------------------------------
Результаты.
>> EvalGauss(x,y,0.32,1)
ans =
1.21331463672730
>> EvalGauss(x,y,0.32,2)
ans =
1.21331463672730
Вывод. Для большого количества точек интерполирования можно прерывать вычисления конечных разностей с учетом заданной точности. При этом выбор формулы для интерполирования зависит от расположения заданной точки по отношению к точкам интерполирования : если точка находится вначале таблицы или в конце, вычисления проводят по первой и второй формуле Ньютона соответственно, если точка лежит в середине таблицы и ниже диагональных элементов, то работают с первой формулой Гаусса, если выше диагональных элементов – со второй формулой Гаусса. Это позволяем уменьшить погрешность вычислений.
Уважаемый посетитель!
Чтобы распечатать файл, скачайте его (в формате Word).
Ссылка на скачивание - внизу страницы.