Задание:
- Аппроксимировать модельную функцию f(x) из таб. 1.1 двумя методами согласно таб. 8.1.
- Определить дифференциальную и интегральную погрешности аппроксимации для всех вариантов.
- Построить графики аппроксимируемой и аппроксимирующей функций, а также графики погрешностей от числа узлов аппроксимации.
Функция:
3x4+4x3-12x*x+1=0
Методы согласно таблице 8.1:
Метод 1: Метод наименьших квадратов.
Метод 2: Метод Стирлинга.
Метод решения СЛАУ: метод с обратной матрицей.
Метод интегрирования: Метод Ньютона-Котеса (n=6).
Число узлов сетки: 12, 16, 25.
Для приближенного поиска первого положительного корня был построен график функции:
Как видно из графика корень лежит в интервале от 0.2 до 0.4. Для уточнения корня используем метод деления пополам. С помощью этого метода был найден корень:
уточнённый корень: 0.30893109560193.
Отрезок аппроксимации: [0.30893109560193-2, 0.30893109560193+3]
Разработка алгоритмов для –заданных видов аппроксимации
Алгоритм аппроксимации методом наименьших квадратов:
Алгоритм аппроксимации методом Стирлинга:
Также для решения СЛАУ был разработан метод с обратной матрицей:
Алгоритм Ньютона-Котеса (n=6)
Один шаг:
Полный алгоритм:
Тексты программ:
Файл «DiffPogr.m»:
1 2 3 4 5 6 7 8 9 10 11 |
function result=DiffPogr(F,FA,xl,xr,N) h=(xr-xl)/N; x=xl;i=0;result=0; while (x<=xr) i=i+1; result=result+abs(feval(F,x)-feval(FA,x)); x=xl+h*i end; result result=result/N; |
Файл «DN.m»:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 |
function result=DN(X,Y,N,x); %eta procedure vichislyaet znacheniya DY (5-i shtuk) DY=DNone(Y,N); DY(N)=0; DY'; D2Y=DNone(DY,N-1);D2Y(N)=0; D2Y(N-1)=0;D2Y'; D3Y=DNone(D2Y,N-2);D3Y(N)=0; D3Y(N-1)=0;D3Y(N-2)=0;D3Y'; D4Y=DNone(D3Y,N-3);D4Y(N)=0; D4Y(N-1)=0;D4Y(N-2)=0;D4Y(N-3)=0;D4Y'; D5Y=DNone(D4Y,N-4);D5Y(N)=0; D5Y(N-1)=0;D5Y(N-2)=0;D5Y(N-3)=0;D5Y(N-4);D5Y'; i=1; while (X(i)<x) index=i; i=i+1; end; res=zeros(N,5); res(:,1)=DY';%(index);poldecrement(index); res(:,2)=D2Y';%(index);poldecrement(index); res(:,3)=D3Y';%(index);poldecrement(index); res(:,4)=D4Y';%(index);poldecrement(index); res(:,5)=D5Y';%(index);poldecrement(index); result=res; |
Файл «DNone.m»:
1 2 3 4 5 6 |
function result=DNone(Y,N); DY=0; for i=1:(N-1) DY(i)=-Y(i+1)+Y(i); end; result=DY; |
Файл «func.m»:
1 2 |
function result=func(x) result=3*x^4+4*x^3-12*x*x+1; |
Файл «funcappms.m»:
1 2 3 4 5 6 7 |
function y=funcappms(x,AI) AI=[ 1.0000 -0.0000 -12.0000 4.0000 3.0000]; y=AI(1)+AI(2)*x+AI(3)*x^2+AI(4)*x^3+AI(5)*x^4; |
Файл «graphic.m»:
1 2 3 4 5 6 7 8 |
i=0;Y=0;x=0.2;X=0; while (x<2) i=i+1; Y(i)=func(x); x=x+0.01; X(i)=x; end; plot(X,Y) |
Файл «lab5.m»:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 |
F=@func; clc xl=0.30893109560193-2; xr=0.30893109560193+3; % kolichestvo uzlov 12, 16, 25 N=16; % minsquare(F,N,xl,xr,0.1) Stirling(F,N,xl,xr,0.000001); FA=@funcappms; % I1=NyutonKotes6(xl,xr,FA,1/N) % I2=NyutonKotes6(xl,xr,F,0.1/N) % inteps=I1-I2 % difeps=DiffPogr(F,FA,xl,xr,N) |
Файл «minsquare.m»:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 |
function result=minsquare(F,n,xl,xr,h2) h=(xr-xl)/n; x=xl;i=0; while(x<=xr) i=i+1; X(i)=x; Y(i)=feval(F,x); x=xl+h*i; end; % X=0; % Y=0; % X=[0.75,1.5,2.25,3,3.75] % Y=[2.5,1.2,1.12,2.25,4.28] % m=2;n=4; mc=5; A=zeros(mc,mc) B=zeros(mc,1); nc=n+1; for k=1:mc for l=1:mc for i=1:nc i A(k,l)=A(k,l)+(X(i))^(k-1+l-1); end; end; for i=1:nc B(k)=B(k)+X(i)^(k-1)*Y(i); end; end; A B %bilo polucheno SLAU teper reshim SLAY i poluchim aproksimiruyaushuyu %funkciyu AI=SLAU(A,B) Y1=0; x=xl; ysize=(xr-xl)/h2+1; i=0; while (x<xr) i=i+1; X1(i)=x; Y1(i)=AI(1)+AI(2)*X1(i)+AI(3)*X1(i)^2+AI(4)*X1(i)^3+AI(5)*X1(i)^4; x=xl+h2*i; end; X1; Y1; hold on; plot(X,Y,'g.-') plot(X1,Y1,'b.-') result=AI; |
Файл «NyutonKotes6.m»:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 |
function RES=NyutonKotes6(A,B,F,e) % function RES=NyutonKotes6(A,B,F,X,Y,h,middle,DY,e) number_of_iteration=1; ce=e*2; I=0;dx=B-A; x=A; I=NyutonKotes6step(x,x+dx,F); % I=NyutonKotes6step(x,x+dx,F,X,Y,h,middle,DY); while (ce>e) x=A; dx=dx*0.5; II=0; while (x<B) II=II+NyutonKotes6step(x,x+dx,F); % II=II+NyutonKotes6step(x,x+dx,F,X,Y,h,middle,DY); x=x+dx; end; ce= abs(I-II) I=II; % number_of_iteration=number_of_iteration+1; end; RES=I; |
Файл «NyutonKotes6step.m»:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 |
function RES=NyutonKotes6step(A,B,F) % function RES=NyutonKotes6step(A,B,F,X,Y,h,middle,DY) H=[41, 216, 27, 272, 27, 216, 41]; N=840;m=1; x=A; n=6;dx=(B-A)/n; I=0; for i=1:n+1 I=I+H(i)*feval(F,x); % q=(x-X(middle))/h; % % yy=Y(middle)+q*(DY(middle-1,1)+DY(middle,1))/2+q^2/2*DY(middle-1,2)+q*(q^2-1)/factorial(3)*(DY(middle-2,3)+DY(middle-1,3))/2+q^2*(q^2-1)/factorial(4)*DY(middle-2,4)+q*(q^2-1)*(q^2-2^2)/factorial(5)*(DY(middle-3,5)+DY(middle-2,5))/2; % I=I+H(i)*yy; x=A+dx*i; end; I; I=(I/N)*n*dx; RES=I; |
Файл «poldel.m»:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 |
function poldel; x=0.3; xl=0.2; xr=0.4; eps=0.000000000000000001;y=func(x); while (abs(y)>eps) y=func(x); if (y>0) if (func(xl)<0) xr=x; x=(x+xl)/2; else xl=x; x=(x+xr)/2; end; else if (func(xl)>0) xr=x; x=(x+xl)/2; else xl=x; x=(x+xr)/2; end; end; xl x xr func(xl) func(x) func(xr) end; y x |
Файл «SLAU.m»:
1 2 3 4 5 |
function X=SLAU(A,B); A B inv(A) X=inv(A)*B |
Файл «Stirling.m»:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 |
function result=Stirling(F,N,xl,xr,h2); if (mod(N,2)==1) N=N+1;end; h=(xr-xl)/N; x=xl;i=0; while(x<=xr) i=i+1; X(i)=x; Y(i)=feval(F,x); x=xl+h*i; end; % X(i)=x; x=xl+h*(i-1); DY=DN(X,Y,N+1,x); Y1=zeros(N+1); middle=round((N+1)/2); Y1=0 Y i=0; xl xr x=xl X1=0; Y1=0; ysize=(xr-xl)/h2+1; difpog=0; while (x<xr) i=i+1; q=(x-X(middle))/h; X1(i)=x; Y1(ysize-i)=Y(middle)+q*(DY(middle-1,1)+DY(middle,1))/2+q^2/2*DY(middle-1,2)+q*(q^2-1)/factorial(3)*(DY(middle-2,3)+DY(middle-1,3))/2+q^2*(q^2-1)/factorial(4)*DY(middle-2,4)+q*(q^2-1)*(q^2-2^2)/factorial(5)*(DY(middle-3,5)+DY(middle-2,5))/2; difpog=difpog+abs(feval(F,xr-(x-xl))-Y1(ysize-i)) x=xl+h2*i; end; % % Y1(1)=Y(1); % X1(ysize-2)=X(N+1); % Y1(ysize-2)=Y(N+1); X1; Y1; hold on; plot(X,Y,'g.-'); plot(X1,Y1,'b'); % LLX=[xl,xl]; % YY=[Y(1),Y(N)]; % RRX=[xr,xr]; % plot(LLX,YY,'y'); % plot(RRX,YY,'y'); % difpog=difpog/i % intpog=NyutonKotes6(xl,xr,F,X,Y,h,middle,DY,0.1) intpog=NyutonKotes6(xl,xr,F,0.0001) result=Y1; |
Результаты аппроксимации функции:
(на графиках зелёным цветом обозначен график реальной функции, а синим – график аппроксимирующей функции)
Для метода наименьших квадратов:
Число узлов 12:
Дифференциальная погрешность = 4.699944137579830e-015
Интегральная погрешность = 0.0767
Число узлов 16:
Дифференциальная погрешность = 5.141720382795256e-015
Интегральная погрешность = 0.0358
Число узлов 25:
Дифференциальная погрешность = 1.596639487289053e-015
Интегральная погрешность = 0.0179
Для метода Стирлинга:
Число узлов 12:
Дифференциальная погрешность = 3.287053575018462e-014
Интегральная погрешность = 0,0000000000002
Для числа узлов 16:
Дифференциальная погрешность = 2.252785354639464e-013
Интегральная погрешность = 0,0000000000006
Для числа узлов 25:
Дифференциальная погрешность =5.491511713503260e-012
Интегральная погрешность = 0,0000000000027
График зависимости интегральной погрешности от числа узлов сетки для метода наименьших квадратов:
График зависимости дифференциальной погрешности от числа узлов сетки для метода наименьших квадратов:
График зависимости дифференциальной погрешности для метода Стирлинга от числа узлов сетки:
График зависимости интегральной погрешности от числа узлов сетки для метода Стирлинга:
Анализ полученных результатов: проведенное исследование метода наименьших квадратов и метода Стирлинга показало, что при помощи этих методов можно аппроксимировать функцию.
Анализ зависимости погрешности от количества узлов показал, что для метода Стирлинга дифференциальная и интегральная погрешность растут с увеличением числа узлов, а для метода наименьших квадратов наблюдается обратная зависимость, что можно объяснить следующим образом: при увеличении числа узлов для метода Стирлинга происходит рост относительной погрешности определения значений конечных разностей, в то время как в случае со методом наименьших квадратов происходит увеличение числа многочленов, график каждого из которых в итоге всё более приближается к графику того участка исходной кривой, за аппроксимацию которого он отвечает.
Просмотр графиков аппроксимирующей и аппроксимируемой функции показал, что «на глаз» они практически неотличимы на заданном отрезке аппроксимации.
Выводы: используя 3x4+4x3-12x*x+1=0 в качестве аппроксимируемой функции функцию, мы изучили два метода аппроксимации: метод Стирлинга и метод наименьших квадратов, построили для обоих методов графики аппроксимирующей и аппроксимируемой функций для различного количества узлов. Были оценены зависимости погрешности от числа узлов для обоих методов.