Міністерство Освіти та Науки України
Національний авіаційний університет
Кафедра інформаційних технологій
з дисципліни "Математичні моделі Інформаційних процесів та управління динамічних об'єктів"
Тема: Цифрова математична модель динаміки польоту літака Ту-І54Б в режимі автоматичного управління приладовою швидкістю через руль висоти
Виконавець: студент 3 курсу 303 групи ФІТ
Керівник: к.т.н., доцент Полухін Анатолій Васильович
Київ 2001 р.
Завдання на виконання курсового проекту
Дата видачи завдання на проектування: 01. 10. 2001 р.
Дата захисту проекту за навчальним планом: 01. 12. 2001 р.
Перелік питань, що підлягають розробці:
розробити математичну модель літака, та дослідити заданий закон автоматичного управління (стабілізація приладової швидкості через руль вистоти).
Початкові значення параметрів польоту літака:
Збурення, що діють на літак: Відмова центрального двигуна.
Метод чисельного інтегрування системи диференціальних рівнянь:
Рунги-Кутти 2 порядка з корекцією по середній похідній.
Дослідження динамічної подібності розробленої моделі та реального "вільного" об'єкту, використовуючи метод Рунге-Кутти другого порядку з корекцією у середній похідній.
Дослідження якості стабілізації нового заданого значення приладової швидкості при відсутності збурення
Дослідження якості стабілізації приладової швидкості при відмові центрального двигуна методом Рунге-Кутти другого порядку з корекцією у середній похідній.
uses crt;
type mass=array[1..7] of real;
matt=array[1..11] of real;
const
MzDZ = -0.0091; RAD= 180/PI; Ph = -0.287; CyDV = 0.005;
FIdv = 0; Pdg = 132.89; S = 201.45; Cy0 = -0.255;
_g=9.81; mzo=0.2; dT=0.05; CyFIst = 0.0132;
Po = 3500; KoefDZ = -0.14706; Pv = -9; cy_alfa=0.1008;
{____________________________________________________________________________}
{ Даннi для інтерполяції }
mzffist:mass=( -0.0465, -0.0475, -0.0480, -0.0498, -0.0520, -0.0565, -0.0580);
mzwzTab:mass=(-12.850, -12.910, -13.000, -13.900, -14.000, -14.900, -15.800);
mzalfa:mass=(-0.0318, -0.0320, -0.0325, -0.0340, -0.0360, -0.0440, -0.0510);
ChisloMahaTab:mass=(0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9);
CxTab:mass=(0.020, 0.022, 0.030, 0.040, 0.054, 0.082, 0.135);
Cy_Cx:mass=(0.000, 0.200, 0.400, 0.600, 0.800, 1.000, 1.200);
{Таблицi на 11 значень, для більшої точності }
time:matt=(0, 0.5, 1, 1.5, 2, 2.5, 2.8 , 3, 3.5, 4, 5);
pft:matt=(1686,1641.04,1475.25,1253.26, 899.2, 351.25, 0, -97, -240, -288, -299);
VzvukaTab:matt= (340.29, 339.91, 339.53, 339.14, 338.76, 337.98, 337.60, 337.21, 336.82, 336.43, 320.54);
HTab:matt= (0.0, 100.0, 200.0, 300.0, 400.0, 600.0, 700.0, 800.0, 900.0, 1000.0, 5000.0);
alfa_mzalt :matt =(0,3.600,5,6.2500,7.500,10,12.800,14.2,15,18,19.7);
Mzaltfal :matt =(-2,-2.2,-2.6,-3.1,-3.7,-5,-6.9,-8.8,-10.8,-13.6,-17);
CyTab__:matt=( 0.,-0.28, -0.07, 0.11, 0.30, 0.47, 0.63, 0.80, 0.93, 1.08, 1.18);
Cyal__:matt=( -1., 0., 2., 4., 6., 8., 10., 12., 14., 16., 18);
mzdvtab1:mass=(-0.01450,-0.014525,-0.01455,-0.014585, -0.01460,-0.01462, -0.01465);
mzdval1:mass=(0.000000, 2.5000, 5.000000, 8.5000, 10.00000, 12.000, 15.00000);
var
Buf,Cy ,DZ ,Cx ,CyAL ,AL ,MzFIst ,dVpr ,MzAL ,MzAL1 ,MzWz ,M ,MzDV ,Vprz
,dMzdv, _m ,G ,q ,Ro ,MzStat ,DV ,Xt ,FIst ,Pbal , MzVra ,Ba ,Mz ,Ye ,Iz ,P
,DGBal ,Xe ,Vpro ,Fun ,flag ,CyBal ,POtk ,TOtk ,ALBal ,DVBal ,Roh ,DG ,a
,Vpr ,Vprzad,Vpr0,Ny ,Tetab ,dTeta ,
delta1 ,delta ,Sigma , dVprzad ,CyAL1,T,
CxBal,XBal,YBal,Qdv,dvSau,wv,tv__:real;
pbaldv,dny:real;
kv,kv_,tv,tv_,kidv,qq,sss,P1:real;
mza,g0,dmz, deltavpr,flagg,flagbal,pdv,mz_cy,delta_ny,dv_ny:real;
CyB,cxb,gtop,Kh1,K_v,K_v1,T_v,T_v1,Kteta,Kh,Kwz,Ktet,Kint:real;
F :text;
x ,y ,d:array[1..20] of real;
vib,met,ssau,i ,j :integer;
{____________________________________________________________________________}
procedure menu;
begin
clrscr;
textcolor(3);
write('Курсова робота з дисциплiни: ');
writeln ('Математичнi моделi iнформацiйних процесiв та управлiння динамiчними обьєктами');
writeln('Виконав студент групи 303 ФiТ Єрофеєв Юрiй');
writeln;
writeln(' Варианти польоту: ');
writeln(' 1: Вiльний лiтак');
writeln(' 2: Вiдмова центрального двигуна');
writeln(' 3: Вiдмова центрального двигуна i стабiлiзацiя заданої прил. швидкостi через руль висоти');
writeln;
write(' Введiть варiант польоту: ');
readln(vib);
writeln;writeln;
textcolor(7);
{writeln('Пiдключити САУ ? 1-Так/2-Нi');
readln(ssau);
writeln; }
writeln('Метод чисельного iнтегрування');
writeln(' 1: Ейлера');
writeln(' 2: Рунге-Кутти з корекцiєю по середнiй похiднiй');
writeln;
write(' Ваш вибiр ');
readln(met);
end;
{____________________________________________________________________________}
{Интерполяція функциї 1 змінної методом Лагранжу}
function InterpolFunction1(x1:real; var xm,fm:mass ):real;
const n=6;
var i,j:integer;
p1,p2,p,r,r1,r2,q,s:real;
begin
p:=0.0; p1:=p; p2:=p;
for i:=1 to n do
begin
r:=1.0; r1:=0.0; r2:=r1;
for j:=1 to n do if i<>j then
begin
q:=xm[i]-xm[j]; s:=x1-xm[j]; r2:=(r2*s+2*r1)/q;r1:=(r1*s+r)/q;r:=r*s/q;
end;
p2:=p2+r2*fm[i];p1:=p1+r1*fm[i];p:=p+r*fm[i];
end;
InterpolFunction1:=p;
end;{InterpolFunction1}
{____________________________________________________________________________}
{Интерполяція функциї змінної методом Лагранжу таблиць на 11 значень}
function InterpolFunction11(x1:real; var xm,fm:matt ):real;
const n=10;
var i,j:integer;
p1,p2,p,r,r1,r2,q,s:real;
begin
p:=0.0; p1:=p; p2:=p;
for i:=1 to n do
begin
r:=1.0; r1:=0.0; r2:=r1;
for j:=1 to n do
if i<>j then
begin
q:=xm[i]-xm[j]; s:=x1-xm[j]; r2:=(r2*s+2*r1)/q;r1:=(r1*s+r)/q;r:=r*s/q;
end;
p2:=p2+r2*fm[i];p1:=p1+r1*fm[i];p:=p+r*fm[i];
end;
InterpolFunction11:=p;
end;{InterpolFunction1}
{____________________________________________________________________________}
{Процедура попередження користувача}
procedure attention;
begin
sound(100);
delay(2);
nosound;
end;
{____________________________________________________________________________}
{Процедура котролю за кутом атаки, закидом вертикального перевантаження, кількістю палива}
procedure control;
begin
dny:=ny-1;
if y[8]<0 then attention;
if (dny>0.2) or (dny<-0.2) then attention;
if y[5]>16 then attention;
end;
{____________________________________________________________________________}
{Вивід результатів на екран і в файл}
procedure PRINT;
begin
write('|',t:6:1); write(y[3]:9:2,' '); write(vpr:6:3,'');
write(y[6]:9:2,''); write(y[5]:7:3,' '); write(y[1]:7:3,' ');
write(P1:8:3, ''); write(dV:8:3,' '); writeln(ny:6:4,' |');
write(f,'|',t:6:1); write(f,y[3]:9:2,' '); write(f,vpr:7:3,' ');
write(f,y[6]:9:2,' '); write(f,y[5]:7:3,' '); write(f,y[1]:7:3,' ');
write(f,P1:8:3, ''); write(f,dV:8:3,' ');
write(f,ny:6:4,' ');
write(f,(2*Pdv+P1):6:4,' ');
writeln(f,y[4]:6:4,' |');
end;
{____________________________________________________________________________}
procedure prisvoenie;
begin
if flagg=1 then
begin
for i := 1 to 20 do X[i] := 0;
flagg := 0;
end;
Y[3] := 110; Y[6] := 1000; Gtop:=20000; G0:=53000; y[8]:=Gtop;
G :=Gtop+G0; Iz := 660000; DZ := 0; Ba := 5.285; FIst :=0;
Xt := 0.24; Ro := 0.1249; DG := 0; DV := 0; T := 0;
flagg:=1; flag:=1; flagbal:=1;
{SAU}
Kh1:= 0.4; K_v:= 1.8; K_v1:= 2.16; T_v:= 1.; T_v1:=1.; Kteta:= 2.0;
Kh:= 0.1; Kwz:= 3.0; Ktet:= 0.5; Kint:= 0.002; kidv:=0.04;
end;
{____________________________________________________________________________}
{ Чисельне інтегрування диференціальних рівнянь методом Эйлера }
procedure Eiler;
var i :integer;
begin
for i:=1 to 20 do y[i]:=y[i]+x[i]*dt;
end;{Integrirovanie}
{____________________________________________________________________________}
{Вивід на екран початкових значень польоту та показників стійкості}
procedure print_bal;
begin
clrscr;
writeln('Балансованi (початковi) значення польоту');
writeln;
writeln('Швидкiсть ',y[3]:5:3);
writeln('Висота ',y[6]:5:3);
writeln('Швидкiсть звуку ',a:6:3);
writeln('Густина повiтря ',roh:6:3);
writeln('Число маха ',M:6:3);
writeln('Cy ',cyb:6:3);
writeln('Cx ',cxb:6:3);
writeln('Кут атаки ',albal:6:3);
writeln('Руль висоти ',dvbal:6:3);
writeln('Дельта газу ',dgbal:6:3);
writeln('Ybal ',ybal:6:3);
writeln('Xbal ',xbal:6:3);
writeln('Початкова тяга силової установки ',pbal:6:3);
writeln('Початкова тяга двигуна ',p1:6:3);
writeln('Початкова приборна швидкiсть ',vpr0:6:3);
writeln;
textcolor(3);
writeln('Ступiнь повздовжноъ статичноъ стiйкостi ',mz_cy:6:3);
writeln('Запас стiйкостi по перевантаженню ',delta_ny:6:3);
writeln('Витрата руля висоти на одиницю перевантаження ',dv_ny:6:3);
writeln('Перiод коливань (по швидкостi) ',tv__:6:3);
textcolor(7);
readkey;
clrscr;
end;
{____________________________________________________________________________}
{ Балансуваня літака }
procedure Balansirovka;
begin
CyB:=(2*G)/(S*Roh*y[3]*y[3]);
albal:=InterpolFunction11(cyb,cytab__,cyal__);
CxB:=InterpolFunction1(CyB,Cy_Cx,CxTab);
Ybal := CyB*S*q;
Xbal := CxB*S*q;
Mzdv:=InterpolFunction1(albal,mzdval1,mzdvTab1);
DVBal := -(Mzo +MzAL*ALBal + CyB*(Xt-0.24) + MzDZ*DZ + MzFIst*FIst)/MzDV;
Pbal := Xbal/cos((ALBal+FIdv)/RAD);
DGBal:= (Pbal/3-Po-Pv*Y[3]-Ph*Y[6])/Pdg+66;
P1 := Pbal/3;
Y[5] := ALBal; Y[1] := Y[5] ;Tetab := ALBal;
Mzwz:=InterpolFunction1(M,ChisloMahaTab,mzwzTab);
mz_cy:=mzal/cy_alfa;
delta_ny:=mz_cy+(mzwz*ba*S*roh)/(2*_m);
mzdv:=mzdv*rad;
dv_ny:=-rad*(delta_ny*Cyb)/mzdv;
wv:=(_g/y[3])*sqrt((2*mz_cy)/delta_ny);
tv__:=(2*3.14)/wv;
vpr:=y[3]*sqrt(roh/ro);
vpr0:=vpr;
vprzad:=vpr0+10;
if vib=3 then vprzad:=vpr0;
print_bal;
end;
{____________________________________________________________________________}
{Закон стабілізаціїї приладової швидкості через руль висоти}
procedure SAU;
begin
deltavpr:=vprzad-vpr;
x[15]:=(1-y[15])/t_v;
x[16]:=(1-y[16])/t_v1;
x[17]:=kidv;
fun:=y[17]*deltavpr;
if(Fun>10) then Fun:=10;
if(Fun<-10) then Fun:=-10;
delta1:=(k_v*x[15]*deltavpr+x[16]*k_v1*100*deltavpr+fun);
if delta1 < -10 then delta1:=-10;
if delta1 > 10 then delta1:=10;
dTeta:=Y[1]-Tetab;
delta:=delta1+Kteta*dTeta;
if delta <-3.5 then delta:=-3.5;
if delta >3.5 then delta:=3.5;
Sigma:=delta+Kwz*Y[2];
if Sigma<-10 then Sigma:=-10;
if Sigma>10 then Sigma:=10;
X[9]:=0;
if delta<-2 then X[9]:=-0.6;
if delta>2 then X[9]:=0.6;
DVsau:=Sigma+Y[9];
end;
{____________________________________________________________________________}
{Система диференціальних рівнянь}
procedure dinamika;
begin
Roh := Ro - 0.0117*Y[6]/1000 + 0.000343*sqr(Y[6]/1000);
G:=G0+y[8];
_m := G/_g;
q := Roh*sqr(Y[3])/2;
a:=InterpolFunction11(y[6],HTab,VzvukaTab);
Vpr := y[3]*sqrt(Roh/Ro);
dVpr := Vpr - Vprzad;
M := Y[3]/a;
Qdv:=0.2586*M+0.4428;
mzfist:=InterpolFunction1(M,ChisloMahaTab,mzffist);
mzal:=InterpolFunction1(M,ChisloMahaTab,mzalfa);
if flagbal=1 then
begin
Balansirovka;
flagbal:=0;
end;
Mzal1:=InterpolFunction11(Y[5],alfa_mzalt , mzaltfal);
Mzdv:=InterpolFunction1(y[5],mzdval1,mzdvTab1);
Mzwz:=InterpolFunction1(M,ChisloMahaTab,mzwzTab);
Cy:=InterpolFunction11(y[5],cyal__,cytab__);
Cx:=InterpolFunction1(Cy,Cy_Cx,CxTab);
Ye := Cy*S*q;
Xe := Cx*S*q;
{ Відмова двигуна }
if (vib=2) or (vib=3) then begin
if ((T>=0) and (T<=4)) then
P1:=interpolfunction11(T,time,pft);
end;
if ((vib=3) or (ssau=1)) and (t>1.0) then SAU;
dv:=dvbal+dvsau;
MzStat:=mzo+mzal*Y[5]+mzdv*dv+mzdz*dz+mzfist*fist;
MzVra := (MzWz*Y[2] + MzAL1*X[5])*Ba/(RAD*Y[3]);
Mz := (MzStat + MzVra);
DG := 0;
Pdv := Po + Pv*Y[3] + Ph*Y[6] + Pdg*(DGBal+DG-66);
dMzdv := (Pbal/3 - Potk) * 0.5;
{____________________________________________________________________________}
x[1]:=y[2];
x[2]:=(Mz*S*Ba*q+dMzdv)*RAD/Iz;
x[3]:=(-G*sin(Y[4]/RAD)+(2*Pdv+P1)*cos((Y[5]+FIdv)/RAD)-Xe)/_m;
x[4]:=(-G*cos(Y[4]/RAD)+(2*Pdv+P1)*sin((Y[5]+FIdv)/RAD)+Ye)*RAD/(_m*Y[3]);
x[5]:=x[1]-x[4];
x[6]:=y[3]*sin(y[4]/RAD);
x[7]:=y[3]*cos(y[4]/RAD);
x[8]:=-3*qdv;
{____________________________________________________________________________}
{Добалансування літака}
if flag =1 then begin
d[1]:=-x[1];
d[2]:=-x[2];
d[3]:=-x[3];
d[4]:=-x[4];
d[5]:=-x[5];
d[6]:=-x[6];
flag:=0;
end;
x[1]:=x[1]+d[1];
x[2]:=x[2]+d[2];
x[3]:=x[3]+d[3];
x[4]:=x[4]+d[4];
x[5]:=x[5]+d[5];
x[6]:=x[6]+d[6];
{___________________________________________________________________}
{Вертикальне перевантаження }