Применение численных методов для решения уравнений с частными производными

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

САНКТ-ПЕТЕРБУРГСКИЙ УНИВЕРСИТЕТ ПУТЕЙ СООБЩЕНИЯ

Кафедра «Прикладная математика»

ОТЧЕТ

ПО ВЫПОЛНЕННОЙ КУРСОВОЙ РАБОТЕ

Предмет «Численные методы»

«Применение численных методов для решения Уравнений с частными производными»

Санкт-Петербург 2008г.

Лабораторная работа N1 "Интерполирование алгебраическими многочленами"

Для решения задачи локального интерполирования алгебраическими многочленами в системе MATLAB предназначены функции polyfit (POLYnomial FITting - аппроксимация многочленом) и polyval (POLYnomial VALue - значение многочлена).

Функция polyfit (X,Y,n) находит коэффициенты многочлена степени n , построенного по данным вектора Х, который аппроксимирует данные вектора Y в смысле наименьшего квадрата отклонения. Если число элементов векторов X и Y равно n+1, то функция polyfit (X,Y,n) решает задачу интерполирования многочленом степени n.

Функция polyval (P,z) вычисляет значения полинома, коэффициенты которого являются элементами вектора P, от аргумента z . Если z - вектор или матрица, то полином вычисляется во всех точках z.

Воспользуемся указанными функциями системы MATLAB для решения задачи локального интерполирования алгебраическими многочленами функции, заданной таблицей своих значений

X	0.0	1.0	2.0	3.0	4.0
Y	1.0	1.8	2.2	1.4	1.0

и вычисления ее приближенного значения в точке x* = 2.2 .

Задача 1 (задача локального интерполирования многочленами)

Построить интерполяционные многочлены 1-ой, 2-ой и 3-ей степени.

Вычислить их значения при x=x*.

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

Решение задачи средствами системы MATLAB:

X=[0.0000 0.5000 1.0000 1.5000 2.0000 2.5000];

Y=[0.0378 0.0653 0.3789 1.0353 0.5172 0.9765];

xzv=1.61;

P1=polyfit(X(4:5),Y(4:5),1) Коэффициенты многочлена P1

P2=polyfit(X(3:5),Y(3:5),2) Коэффициенты многочлена P2

P3=polyfit(X(3:6),Y(3:6),3) Коэффициенты многочлена P3

Полученные таким образом коэффициенты интерполяционных многочленов и значения этих многочленов при x=x* :

P1 = -1.0362 2.5896

P2 = -2.3490 7.1853 -4.4574

P3 = 2.8692 -15.2604 25.8351 -13.0650

z1 = 0.9213

z2 = 1.0221

z3 = 0.9470

многочлены P1, P2, P3

P₁ = -1.0362*X+2.5896

P₂ = -2.3490*X²+7.1853*X+-4.4574

P₃ = 2.8692*X³ -15.2604*X² + 25.8351 + -13.0650

Для построения графиков интерполяционных многочленов следует создать векторы xi1, xi2, xi3, моделирующие интервалы (X(3):X(4)), (X(2):X(4)),(X(2):X(5)), соответственно, и вычислить значения многочленов P1, P2, P3 для элементов векторов xi1, xi2, xi3, соответственно:

xi1=X(4):0.05:X(5);

xi2=X(3):0.05:X(5);

xi3=X(3):0.05:X(6);

y1=polyval(P1,xi1);

y2=polyval(P2,xi2);

y3=polyval(P3,xi3);

plot(X,Y,'*k',xi1,y1,xi2,y2,xi3,y3);grid

Интерполирование нелинейной функцией Y=A*exp(-B*X)

y_l=log(Y)

Pu=polyfit(X(4:5),y_l(4:5),1)

z_l=(exp(Pu(2))*exp(Pu(1)*xzv))

Y= 8.3040*exp(-1.3880*X)

Функция plot с указанными аргументами строит табличные значения функции черными звездочками('*k'), а также графики многочленов P1 (по векторам xi1 и y1), P2 (по векторам xi2 и y2) и P3 (по векторам xi3 и y3), и функцией Y=A*exp(-B*X), соответственно синей, красной и зеленой кривыми.

plot(X,Y,'*k',xi1,y1,xi2,y2,xi3,y3,xi1,exp(Pu(2))*exp(Pu(1)*xi1));grid

Оценка погрешности интерполирования

При оценке погрешности решения задачи интерполирования в точке x* за погрешность epsk интерполяционного многочлена степени k принимается модуль разности значений этого многочлена и многочлена степени k+1 в точке x*.

С помощью уже полученных значений мы можем оценить погрешности интерполяционных многочленов P1 и P2 в точке x* , используя функцию abs системы MATLAB для вычисления модуля:

eps1 = abs(z1-z2)

eps1 = 0.1008

eps2 = abs(z2-z3)

eps2 = 0.0751

Для оценки погрешности многочлена P3 необходимо предварительно вычислить

значение z4=P4(x*), а затем - eps3.

P4=polyfit(X,Y,4);z4=polyval(P4,xzv);

eps3=abs(z4-z3)

eps3 = 0.1450

«Построение сплайна»

X=[0.0000 0.5000 1.0000 1.5000 2.0000 2.5000];

Y=[0.0378 0.0653 0.3789 1.0353 0.5172 0.9765];

cs = spline(X,[0 Y 0]);

xx = linspace(0,2.5);

plot(X,Y,'*m',xx,ppval(cs,xx),'-k');

h=0.5

esstestvennii spline

A=[4 2 0 0 0 0

1 4 1 0 0 0

0 1 4 1 0 0

0 0 1 4 1 0

0 0 0 1 4 1

0 0 0 0 2 4]

B=[6*(Y(2)-Y(1))/h 0 0 0 0 6*(Y(length(Y))-Y(length(Y)-1))/h]

for i = 2:(length(Y)-1)

B(i)=(3/h)*(Y(i+1)-Y(i-1))

end

S=inv(A)*B'

otsutstvie uzla

A1=[1 0 -1 0 0 0

1 4 1 0 0 0

0 1 4 1 0 0

0 0 1 4 1 0

0 0 0 1 4 1

0 0 0 1 0 -1]

B1=[2*(2*Y(2)-Y(1)-Y(3))/h 0 0 0 0 2*(2*Y(length(Y)-1)-Y(length(Y))-Y(length(Y)-2))/h]

for i = 2:(length(Y)-1)

B1(i)=(3/h)*(Y(i+1)-Y(i-1))

end

S1=inv(A1)*B1'

c1 = spline(X,[S(2) Y S(5)]);

x1 = linspace(0,2.5,101);

c2 = spline(X,[S1(2) Y S1(5)]);

x2 = linspace(0,2.5,101);

plot(X,Y,'ob',xx,ppval(cs,xx),'-',x1,ppval(c1,x1),'*',x2,ppval(c2,x2),'^',xx,spline(X,Y,xx));

h = 0.5000

A =

4 2 0 0 0 0

1 4 1 0 0 0

0 1 4 1 0 0

0 0 1 4 1 0

0 0 0 1 4 1

0 0 0 0 2 4

B = 0.3300 0 0 0 0 5.5116

B = 0.3300 2.0466 0 0 0 5.5116

B = 0.3300 2.0466 5.8200 0 0 5.5116

B = 0.3300 2.0466 5.8200...