Лабораторные работы / Работа 2 - Вариант 4 - 2003 / dyna
.docРХТУ им. Д.И. Менделеева
Кафедра информатики и компьютерного моделирования
Построение динамических моделей простых гидравлических систем
Вариант 4.
Выполнил:
гр. МЦ-43
Москва 2003
Схема гидравлической системы:
Известные величины:
Давления на входе и выходе системы p1 - p6;
Коэффициенты сужающих устройств k1 - k7;
Геометрические высоты емкостей HG1, HG2.
Неизвестные величины:
Высоты столбов жидкости в емкостях h1, h2;
Давления жидкости на дно сосуда p7, p8;
Давления газа над жидкостью p9, p10;
Расходы жидкости v1 - v7.
Основные допущения, накладываемые на модель:
1) Система работает в стационарном режиме;
2) Трубопроводы находятся на одном уровне;
3) Газ - идеальный pV=nRT;
4) Режим - изотермический T=const;
6) Изменение температуры газа в емкости мало;
7) Жидкость несжимаема;
8) Течение жидкости - слева направо;
9) Цилиндрическая форма закрытой емкости с площадью поперечного сечения
S и геометрической высотой HG;
10) Одинаковое давление газа PN в емкостях, не заполненных жидкостью;
Система уравнений математического описания:
A) Уравнения расходов жидкости через клапаны:
1. v1=k1(p1-p7) 1/2
2. v2=k2(p2-p7) 1/2
3. v3=k3(p3-p7) 1/2
4. v4=k4(p8-p4) 1/2
5. v5=k5(p7-p5) 1/2
6. v6=k6(p7-p5) 1/2
7. v7=k7(p7-p8) 1/2
B) Определение давлений жидкости и газа в емкостях:
8. p7=p9+ρgh1
9. p9=PN HG1/(HG1-h1)
10. p8=p10+ρgh2
11. p10=PNHG2/(HG2-h2)
C) Дифференциальные уравнения заполнения емкостей:
12. dh1/dt=(v1+ v2+ v3 –v7– v5–v6)/s1
13. dh2/dt=(v7– v4)/s2
D) Начальные условия:
14. h1(t=0)=h10
15. h2(t=0)=h20
Информационная матрица системы уравнений:
N п/п |
v1 |
v2 |
v3 |
v4 |
v5 |
v6 |
v7 |
p7 |
p8 |
p9 |
p10 |
h1 |
h01 |
h2 |
h02 |
|
1 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
2 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
3 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
4 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
5 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
6 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
7 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
8 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
9 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
10 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
11 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
12 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
13 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
14 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
15 |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
Блок-схема алгоритма расчета:
Текст программы:
program dina;
const m=2 (*количество решаемых ОДУ*);np=10 ;nk=7 ;nv=15 ;
type areal=array [1..m] of real;
var del (*шаг*),k,x0(*начальное значение аргумента*),x1: real;
y0(*начальные значения функций*),y1:areal;
n(*число шагов*),kv(*кратность вывода*),i,kl,kl1,kl2:integer;
vm,v,ak:array[1..nk] of real; p:array[1..np] of real;
hg,h,s,pr:areal; ro,pn,g,x:real;
label 10,20;
FUNCTION SGN (ab:real):integer;
begin
if ab>0 then sgn:=1;
if ab<0 then sgn:=-1;
if ab=0 then sgn:=0;
end;
PROCEDURE dydx (x:real;y:areal;var pr:areal);
(*расчет производных (pr) для текущих x,y*)
label 300,400;
var vm:array[1..7] of real; i:integer;
begin
g:=9.81;
h:=y;
p[9]:=pn*hg[1]/(hg[1]-h[1]);
p[7]:=p[9]+ro*g*h[1]*0.000001;
p[10]:=pn*hg[2]/(hg[2]-h[2]);
p[8]:=p[10]+ro*g*h[2]*0.000001;
v[1]:=ak[1]*sgn(p[1]-p[7])*sqrt(abs(p[1]-p[7]));
v[2]:=ak[2]*sgn(p[2]-p[7])*sqrt(abs(p[2]-p[7]));
v[3]:=ak[3]*sgn(p[3]-p[7])*sqrt(abs(p[3]-p[7]));
v[4]:=ak[4]*sgn(p[8]-p[4])*sqrt(abs(p[8]-p[4]));
v[5]:=ak[5]*sgn(p[7]-p[5])*sqrt(abs(p[7]-p[5]));
v[6]:=ak[6]*sgn(p[7]-p[6])*sqrt(abs(p[7]-p[6]));
v[7]:=ak[7]*sgn(p[7]-p[8])*sqrt(abs(p[7]-p[8]));
pr[1]:=(v[1]+v[2]+v[3]-v[7]-v[5]-v[6])/s[1];
pr[2]:=(v[7]-v[4])/s[2];
for i:=1 to 7 do vm[i]:=v[i]*ro;
if kl=0 then goto 400;
if kl=1 then goto 300;
writeln; writeln(' p(7-10) vm ');
writeln(' ',p[7]:nv:5, vm[1]:nv:5);
writeln(' ',p[8]:nv:5, vm[2]:nv:5);
writeln(' ',p[9]:nv:5, vm[3]:nv:5);
writeln(' ',p[10]:nv:5 ,vm[4]:nv:5);
writeln(' ',' ',vm[5]:nv:5);
writeln(' ',' ',vm[6]:nv:5);
writeln(' ',' ',vm[7]:nv:5);
300:writeln('x= ',x:nv:2, ' y[1]= ',y[1]:nv:5, ' y[2]= ',y[2]:nv:5);
writeln('pr[1]= ',pr[1]:nv:5, ' pr[2]= ',pr[2]:nv:5);
400:
end;
PROCEDURE STEP (x:real;y:areal;del:real;var x1:real;var y1:areal);
(*шаг интегрирования:по X и Y расчет x1 и y1*)
var pr,y12:areal;i:integer;
begin
dydx(x,y,pr);
for i:=1 to m do y12[i]:=y[i]+del*pr[i]/2;
dydx(x+del/2,y12,pr);
x1:=x+del;
for i:=1 to m do y1[i]:=y[i]+del*pr[i];
end;
begin
writeln('ВЫСОТА ЕМКОСТИ (1,2) /М/ = '); read ( hg[1], hg[2]);
writeln('ПЛОЩАДИ ОСНОВАНИЙ (1,2) /M^2/ = '); read (s[1],s[2]);
writeln('ПЛОТНОСТЬ /кг/м^3/ = '); read (ro);
writeln('НАЧ.ДАВЛЕНИЕ /МПа/ = '); read (pn);
writeln('ДАВЛЕНИЕ (1-6) /МПа/ = ');
for i:=1 to 6 do read(p[i]);
writeln('КОЭФ. ПРОПУСКНОЙ СПОСОБНОСТИ (1-7) = ');
for i:=1 to 7 do read (ak[i]);
10: writeln('н.у.:x0,y0[1],y0[2]');read(x0,y0[1],y0[2]);
20: writeln('число шагов,кратность вывода,del'); read(n,kv,del);
writeln('вывод промежуточ.: нет - 0 частич. - 1 полн. - 2 ' );
read(kl1);
writeln('вывод результатов: частич. - 1 полн. - 2 ' );
read(kl2);
for i:=1 to n do
begin
kl:=kl1; step(x0,y0,del,x1,y1);
x0:=x1; y0:=y1;
if (i div kv)=(i/kv) then
begin
write('ПРИ i= ',i:3);
if kl2=1 then writeln(' x= ',x0:10:5,' y[1]= ',y0[1]:nv:5, 'y[2]= ',y0[2]:nv:5);
if kl=2 then begin kl:=kl2;dydx(x0,y0,pr);end
end
end;
writeln('прекращение счета - 0; продолжение - 1; счет с новымми н.у. - 2 ');
read(kl);
if kl=1 then goto 20; if kl=2 then goto 10;
end.
Результаты моделирования:
Параметры гидравлической системы:
HG1=10 м; HG2=10 м;
S1=1 м2; S2=1 м2;
Давления(МПа): p1=1.0; p2=1.0 ; p3=0.7 ; p4=0.7; p5=0.1; p6=0.1; p7=0.1; p8=0.1; pN=0.1
Коэффициенты сужающих устройств:
для номинального режима k[1-9]=0.010
для закрытого вентиля k[i]=0.000
Номинальный режим |
Заполнение емкостей (k4, k5, k6 =0) |
Опорожнение (k1, k2, k3 =0) |
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|
|
|
Графики: