Міністерство освіти та науки України
Національний університет “Львівська політехніка”
Контрольна робота
З курсу
«Алгоритмізація і програмування»
на тему:
«Залежність густини ртуті від температури»
Мета: побудувати графік залежності густини ртуті від температури; отримати розв’язок різними методами розв’язування нелінійних рівнянь.
Завдання:
0 c ≤ t ≤ 30 c
Знайти, в який момент часу t після зміни тиску відносна зміна вихідної величини дорівнюватиме (t)=0,6.
Залежність густини ртуті від температури
Табулювання
Блок-схема алгоритму
Програма мовою С
//tabyljyvannja
#include<stdio.h>
#include<conio.h>
#include<math.h>
float fi(float t)
{return 1/25*t*exp(-(1+0.04)/2/25*t);}
main()
{clrscr();
for(float t=0;t<=30;t+=3)
printf("\n t=%.f fi=%6.4f",t,fi(t));
getchar();
}
Результат виконання
t=0 fi=0.0 t=12 fi=0.3740 t=24 fi=0.5827
t=3 fi=0.1127 t=15 fi=0.4392 t=27 fi=0.6159
t=6 fi=0.2118 t=18 fi=0.4951 t=30 fi=0.6430
t=9 fi=0.2985 t=21 fi=0.5427
Програма в середовищі MATLAB
%rozraxa2
clear,clc
t=[0:3:30];
fi=kti2(t);
disp([t(:) fi(:)]);
figure(1);
plot(t,fi,'k');
xlabel('t');ylabel('fi(t)');
grid;
Функція-процедура
%kti2
function y=kti2(t)
T0=25; Kp=0.04; V0=1;
y=V0/T0.*t.*exp(-(1+Kp)/2./T0.*t);
Результат виконання програми
0 0 12.0000 0.3740 24.0000 0.5827
3.0000 0.1127 15.0000 0.4392 27.0000 0.6159
6.0000 0.2118 18.0000 0.4951 30.0000 0.6430
9.0000 0.2985 21.0000 0.5427
Графік залежності густини ртуті від температури
Виділення кореня
fi(t)=1/25∙t∙exp(-(1+0.04)/2/25*t)-0,6
t([0;30]. fi(0)=0; fi(30)= 0.6430
Блок-схема алгоритму
Програма мовою С
//vudilennja_korenja
#include<stdio.h>
#include<conio.h>
#include<math.h>
float fi(float t)
{return 1/25*t*exp(-(1+0.04)/2/25*t)-0,6;}
main()
{clrscr();
for(float t=0;t<=30;t+=3)
printf("\n t=%.f fi=%6.4f",t,fi(t));
getchar();
}
Результат виконання
t=0 f=-0.6000
t=3 f=-0.4872
t=6 f=-0.3881
t=9 f=-0.3014
t=12 f=-0.2260
t=15 f=-0.1608
t=18 f=-0.1048
t=21 f=-0.0572
t=24 f=-0.0172
t=27 f= 0.0159
t=30 f= 0.0429
Програма в середовищі MATLAB
Функція-процедура
%vud_k
function y=vud_k(t)
T0=25; Kp=0.04; V0=1;
y=V0/T0.*t.*exp(-(1+Kp)/2./T0.*t)-0.6;
%vud_kor
clear,clc
t=[0:3:30];
fi=vud_k(t);
disp([t(:) fi(:)]);
format long;
z=fzero('vud_k',0)
figure(1);
plot(t,fi,'k',z,0,'*k');
xlabel('t');ylabel('fi(t)');
grid;
Результат виконання програми
0 -0.6000 12.0000 -0.2260 24.0000 -0.0172
3.0000 -0.4872 15.0000 -0.1608 27.0000 0.0159
6.0000 -0.3881 18.0000 -0.1048 30.0000 0.0429
9.0000 -0.3014 21.0000 -0.0572
z =25.48761202714193
Метод половинного ділення
Блок-схема алгоритму
Програма мовою С
//pol_dil
#include<stdio.h>
#include<conio.h>
#include<math.h>
float fi(float t)
{return 1/25*t*exp(-(1+0,04)/2/25*t);}
main()
{clrscr();
float a=0,b=30,t0,eps=1e-4;
int k=0;
A:t0=(a+b)/2;
if(fabs(a-b)>eps)
{if(fi(a)*fi(t0)>0) a=t0; else b=t0;
k++; goto A;
}
printf("\n t0=%.7f fi(t0)=%.7f",t0,fi(t0));
printf(" k=%d eps=%.4f a-b=%.7f",k,eps,fabs(a-b));
getchar();
}
Результат виконання програми
t0=0.00002548 fi(t0)=0.0 k=19 eps=0.0001 a-b=0.0000572
Програма в середовищі MATLAB
%pol_dil
clear;clc
a=0; b=30;
k=0; eps=1e-4;
t0=(a+b)/2
while abs(a-b)>eps
if vud_k(a)*vud_k(t0)>0 a=t0;
else b=t0;
end;
t0=(a+b)/2;
k=k+1;
end;
format short;
disp(' t0 ')
disp([ t0 ])
disp(' k ')
disp([ k ])
disp(' abs(a-b) ')
disp([abs(a-b)])
Результат виконання програми
t0=25.4876 k=19 abs(a-b)=5.7220e-005
Метод простих ітерацій
fi(t)=1/25∙t∙exp(-(1+0.04)/2/25*t)-0.6;
fi´(t)=1/25*exp(-(1+0.04)/2/25*t)*(1-(1+0.04)/2/25*t);
φ(t)=t-f(t)∕k , де |k|≥Q∕2
Q=max|f´(t)|
[0;30]
fi´(0)=4∙10-2=Q
fi´(30)=0.81∙10-2
fi´(t)>0 => k>0
k=2*10-2
φ(t)=t-fi(t)/k;
q=max|φ´(t)|
[0;30]
φ´(t)=1- fi´(t)/k;
φ´(0)=0.6777; φ´(30)=-0.6 => q=0.6777<1
Умова виходу:|t1-t|<∆(1-q)∕q; ε=∆∙(1-q)∕q
Блок-схема алгоритму
Програма мовою С
//prostux_iteracij
#include<stdio.h>
#include<conio.h>
#include<math.h>
float f(float t)
{return 1/25*t*exp(-(1+0.04)/2/25*t)-0.6;}
main()
{clrscr();
float a=0,b=30,t,t1,det=1e-4,eps=det*(1-0.6777)/0.6777;
int k=0;
t=(a+b)/2;
A:t1= t-(1/25*t*exp(-(1+0.04)/2/25*t)-0.6)/2/1e-2;
if(fabs(t1-t)>eps)
{t=t1; k++; goto A;}
else
{printf("\n t1=%.4f f(t1)=%.7f",t1,f(t1));
printf(" k=%d eps=%.7f ",k,eps);}
getchar();
}
Результат виконання програми
t1=25.4876 f(t1)=0.0000000 k=14 eps=0.0001
Програма в середовищі MATLAB
%pro_itr
clear;clc
a=0; b=30;
k=0; eps=1e-4;
t=(a+b)/2;
t1=t-(1/25*t*exp(-(1+0.04)/2/25*t)-0.6)/2/1e-2;
while abs(t1-t)>eps*(1-0.6777)/0.6777
t=t1;
t1=t-(1./25.*t.*exp(-(1+0.04)./2./25*t)-0.6)/2/1e-2;
k=k+1;
end;
disp(' a ') disp([ a ])
disp(' t1 ') disp([ t1 ])
disp(' b ') disp([ b ])
disp(' k ') disp([ k ])
disp(' eps ') disp([ eps ])
Результат виконання програми
a=0 t1=25.4876 b=30 k=14 eps=1.0000e-004
Метод хорд
fi(a)<0 fi(b)>0 fi´(t)>0 fi´´(t)>0
t=a-fi(a)∙(b-a)/(fi(b)-fi(a)) ; γ=fi(t0)∙(t0-tfi)/(fi(t0)-fi(tfi)) ;
t=t0-γ;
Умова виходу: |γ|≤ε; ε=∆∙m/(M-m)
m=min |f´(t)|=1.83675∙10-4; M=max |f´(t)|=1.94287∙10-4;
Програма мовою С
//metod_xord
#include<stdio.h>
#include<conio.h>
#include<math.h>
float fi(float t)
{return 1/25∙t∙exp(-(1+0.04)/2/25*t)-0.6;}
main()
{clrscr();
float a=0,b=30,t0,t,tf,ha,det=1e-4,
m=81*1e-4,M=400*1e-4,eps=m/(M-m)*det;
int k=0;
t=a-fi(a)/(fi(b)-fi(a))*(b-a);
if(fi(a)*fi(t)>0) {t0=a;tf=b;} else {t0=b;tf=a;}
A:ha=fi(t0)*(t0-tf)/(fi(t0)-fi(tf));
t=t0-ha;
if(fabs(ha)>eps)
{t0=t; k++; goto A;}
else
{printf("\n t=%.4f fi(t)=%.7f",t,fi(t));
printf(" k=%d eps=%.7f fabs(ha)=%.7f",k,eps,fabs(ha));}
getchar();
}
Результат виконання програми
t=25.4876 f(t)=-0.0000 k=8 eps=0.00002539 abs(ha)=0.0000089559
Блок-схема алгоритму
Програма в середовищі MATLAB
%Metod hord
clear;clc
k=0; t0=0; tf=30;
eps=1e-4; M=400*1e-4; m=81*1e-4;
det=eps*m/(M-m);
ha=vud_k(t0)*(t0-tf)/(vud_k(t0)-vud_k(tf));
t=t0-ha;
while abs(ha)>eps;
t0=t;
ha=vud_k(t0)*(t0-tf)/(vud_k(t0)-vud_k(tf));
t=t0-ha;
k=k+1;
end;
disp(' t '); disp([ t ]);
disp(' k '); disp([ k ]);
format long;
disp(' eps '); disp([ eps ]);
disp(' abs(ha) '); disp([ abs(ha) ]);
Результат виконання програми
t=25.4876 k=8 eps=2.53918e-005 abs(ha)=8.955989e-006
Метод дотичних
fi´(t)=1/25*exp(-(1+0.04)/2/25*t)*(1-(1+0.04)/2/25*t);
t=b-fi(b)/ fi´(b)
γ=fi(t0)/ fi´(t0)
t=t0-γ;
ε=√2∙m1∙∆/M2
Умова виходу: γ≤√2∙m1∙∆/M2
m1=min|fi´(t)| => m1=fi´(30)=81∙10-4
M2=max|fi´´(t)|
fi´´(t)=-(1+0.04)/2/25/25*exp(-(1+0.04)/2/25*t)*(2-(1+0.04)/2/25*t);
fi´´(0)= -16.64∙10-4
fi´´(30)= -6.13∙10-4
M2=16.64∙10-4
Програма мовою С
//metog_dotu4nux
#include<stdio.h>
#include<conio.h>
#include<math.h>
float fi(float t)
{return 1/25*t*exp(-(1+0.04)/2/25*t)-0.6;}
main()
{clrscr();
float a=0,b=30,t0,t,tf,ha,det=1e-4,m1=81*1e-4,
M2=16.64*1e-4,eps=sqrt(2*m1*det/M2);
int k=0;
t=b-fi(b)/(1/25*exp(-(1+0.04)/2/25*t)*(1-(1+0.04)/2/25*t));
if(fi(a)*fi(t)>0) {t0=b;tf=a;} else {t0=a;tf=b;}
A:ha=fi(t0)/(1/25*exp(-(1+0.04)/2/25*t)*(1-(1+0.04)/2/25*t));
t=t0-ha;
if(fabs(ha)>eps)
{t0=t; k++; goto A;}
else
{printf("\n t=%.7f fi(t)=%.7f",t,fi(t));
printf(" k=%d eps=%.7f ha=%.7f ",k,eps,ha);}
getchar();
}
Результат виконання програми
t=25.4876 k=2 eps=0.03120005 ha=-0.02104984
Блок-схема алгоритму
Програма в середовищі MATLAB
Функція-процедура
%poxidna
function y=poxidna(t)
y=1./25.*exp(-(1+0.04)./2./25.*t).*(1-(1+0.04)./2./25.*t);
%metod_dotuchnux
clear;clc
k=0; tf=0; t0=30;
M2=16.64*1e-4; m1=81*1e-4;
det=1e-4;
eps=sqrt(2*m1*det/M2);
ha=vud_k(t0)/poxidna(t0);
t=t0-ha;
while abs(ha)>eps;
t0=t;
k=k+1;
ha=vud_k(t0)/poxidna(t0);
t=t0-ha;
end;
format long;
disp(' t '); disp([ t ]);
disp(' k '); disp([ k ]);
disp(' eps '); disp([ eps ]);
disp(' ha '); disp([ ha ]);
Результат виконання програми
t=25.4876 k=2 eps=0.0312018 ha=-0.02104984
Комбінований метод
Наближення зліва (метод хорд):
γ=fi(a)∙(a-b)/(fi(a)-fi(b)) ; t1=a-γ.
Наближення справа (метод дотичних):
γ=fi(b)/ fi´(b) ; t2=b-γ.
Блок-схема алгоритму
Програма мовою С
//kombinovanuj_metod
#include<stdio.h>
#include<conio.h>
#include<math.h>
float fi(float t)
{return 1/25*t*exp(-(1+0.04)/2/25*t)-0.6;}
main()
{clrscr();
float a=0,b=30,t2,t1,t,ha,det=1e-4;
int k=0;
A:ha=fi(a)*(a-b)/(fi(a)-fi(b));
t1=a-ha;
ha=fi(b)/(1/25*exp(-(1+0.04)/2/25*t)*(1-(1+0.04)/2/25*t));
t2=b-ha;
if(fabs(b-a)>det)
{a=t1;b=t2; k++; goto A;}
else
{t=(t1+t2)/2;
printf("\n t=%.7f",t);
printf(" k=%d eps=%.7f ",k,fabs(b-a));}
getchar();
}
Результат виконання роботи
t=25.4876 k=3 eps= 0.0000000
Програма в середовищі MATLAB
%Kombinovanui metod
clear;clc
a=0; b=30; eps=1e-4; k=0;
b1=b-vud_k(b)/poxidna(b);
a1=a-(vud_k(a)*(a-b1))/(vud_k(a)-vud_k(b1));
while (abs(b-a)>eps)
b1=b-vud_k(b)/poxidna(b);
a1=a-(vud_k(a)*(a-b1))/(vud_k(a)-vud_k(b1));
b=b1;
a=a1;
k=k+1;
end;
t=(a+b)/2;
format long;
disp(' t '); disp([ t ]);
disp(' k '); disp([ k ]);
disp('abs(b-a)'); disp([ abs(b-a) ]);
Результат виконання роботи
t=25.4876 k=3
abs(b-a)=1.442278010799214e-005
Порівняльна таблиця результатів обчислень
Метод уточнення кореня
Задана похибка
Значення кореня
Кількість ітерацій
Половинного ділення
10-4
25.4876
19
Простих ітерацій
10-4
25.4876
14
Хорд
10-4
25.4876
8
Дотичних
10-4
25.4876
2
Комбінований
10-4
25.4876
3