Приложение А

Программа описания динамики ОУ

dt=10;

tk=10000;

t=[0:dt:tk];

n=length(t);

% uravnenie materialnogo balansa dlya teplonositelya

kv1=1; % propusknaya sposobnost klapana goryachego teplonositelya

kv2=0; % propusknaya sposobnost klapana holodnogo teplonositelya

Pt=230000; % davlenie teplonositelya pered klapanom

Pa=101325; % davlenie teplonositelya v rubashke

ro2=1000; % plotnost teplonositelya

Gt=0.1*(kv1+kv2)*((Pt-Pa)/ro2)^0.5/3600; % rashod teplonositelya

N=zeros(1,n);

A=30; % koefficient

E=40000; % energia aktivacii

R=8.314; % gasovaya const

T=zeros(1,n); % temperatura smesi v reaktore

Tvh=15; % temperatura smesi na vhode v reaktor

ca=zeros(1,n);

cavh=20; % konzentraziya komponenta "A" na vxode v reaktor

cb=zeros(1,n);

cbvh=100; % konzentraziya komponenta "B" na vxode v reaktor

cc=zeros(1,n);

K=1000; % koefficient teploprovodnosti

d=1.0; % diametr reaktora

H=1.515; % visota reaktora

s=3.14*d*H; % ploshad bokovoj poverhnosti reaktora

Tt=zeros(1,n); % temperatura smesi teplonositela

Tt1=90; % temperatura teplonositelya 1 na vhode v rubashku

deltaH=10000; % teplovoi effekt reakcii

v1=3.14*H*d^2/4; % ob`em reaktora

ro1=800; % plotnost reakcionnoi massi

cp1=3500; % teploemkost smesi v reaktore

cp2=4100; % teploemkost teplonositelya

a1=1;

D=1.1;

v2=3.14*H*D^2/4-v1; % ob'em rubashki

% nachalnie uslovia

T(1,1)=Tvh; % temperatura smesi v reaktore

ca(1,1)=cavh; % konzentraziya komponenta "A" v reaktore

cb(1,1)=cbvh; % konzentraziya komponenta "B" v reaktore

cc(1,1)=0; % konzentraziya komponenta "C" v reaktore

Tt(1,1)=Tt1; % temperatura teplonositelya v rubashke

for i=2:n

% skorost reakcii

N(1,i)=A*exp(-E/(R*T(1,i-1)))*ca(1,i-1)*cb(1,i-1);

% pokomponentniy materialniy balans dlya reaktora

ca(1,i)=-N(1,i)*dt+ca(1,i-1); % konzentraziya komponenta "A" v reaktore

cb(1,i)=-N(1,i)*dt+cb(1,i-1); % konzentraziya komponenta "B" v reaktore

cc(1,i)=N(1,i)*dt+cc(1,i-1); % konzentraziya komponenta "C" v reaktore

% teplovoy balans (T-reaktor, Tt-rubashka)

T(1,i)=(K*s*(Tt(1,i-1)-T(1,i-1))+deltaH*N(1,i)*v1)*dt/(ro1*v1*cp1)+T(1,i-1);

Tt(1,i)=Tt(1,i-1)+dt*(Gt*ro2*cp2*a1*(Tt1-Tt(1,i-1))-K*s*(Tt(1,i-1)-T(1,i-1)))/(ro2*v2*cp2);

end

% postroenie grafikov

figure(1)

plot(t,T,'b');

hold on;

plot(t,Tt,'r');

legend('T','Tt');

hold on;

grid on;

title('GRAFIKI TEMPERATUR');

figure(2)

plot(t,ca,'b');

hold on;

plot(t,cb,'r');

hold on;

plot(t,cc,'g');

legend('ca','cb','cc');

hold on;

grid on;

title('GRAFIKI KONCENTRACIY');

Приложение Б

Программа для построения переходных характеристик при кусочно-линейном воздействии на объект в динамическом режиме

dt=1;

tk=10000;

t=[0:dt:tk];

n=length(t);

% uravnenie materialnogo balansa dlya teplonositelya

kv1=1; % propusknaya sposobnost klapana goryachego teplonositelya

kv2=0; % propusknaya sposobnost klapana holodnogo teplonositelya

Pt=230000; % davlenie teplonositelya pered klapanom

Pa=101325; % davlenie teplonositelya v rubashke

ro2=1000; % plotnost teplonositelya

Gt=0.1*(kv1+kv2)*((Pt-Pa)/ro2)^0.5/3600; % rashod teplonositelya

N=zeros(1,n);

A=30; % koefficient

E=4000; % energia aktivacii

R=8.314; % gasovaya const

T=zeros(1,n); % temperatura smesi v reaktore

Tvh=15; % temperatura smesi na vhode v reaktor

ca=zeros(1,n);

cavh=20; % konzentraziya komponenta "A" na vxode v reaktor

cb=zeros(1,n);

cbvh=100; % konzentraziya komponenta "B" na vxode v reaktor

cc=zeros(1,n);

K=1000; % koefficient teploprovodnosti

d=1.0; % diametr reaktora

H=1.515; % visota reaktora

s=3.14*d*H; % ploshad bokovoj poverhnosti reaktora

Tt=zeros(1,n); % temperatura smesi teplonositela

Tt1=90; % temperatura teplonositelya 1 na vhode v rubashku

Tt2=15; % temperatura teplonositelya 2 na vhode v rubashku

deltaH=10000; % teplovoi effekt reakcii

v1=3.14*H*d^2/4; % ob`em reaktora

ro1=800; % plotnost reakcionnoi massi

cp1=3500; % teploemkost smesi v reaktore

cp2=4100; % teploemkost teplonositelya

a1=0;

a2=1;

D=1.1;

v2=3.14*H*D^2/4-v1; % ob'em rubashki

% nachalnie uslovia

T(1,1)=Tvh; % temperatura smesi v reaktore

ca(1,1)=cavh; % konzentraziya komponenta "A" v reaktore

cb(1,1)=cbvh; % konzentraziya komponenta "B" v reaktore

cc(1,1)=0; % konzentraziya komponenta "C" v reaktore

Tt(1,1)=Tt1; % temperatura teplonositelya v rubashke

% zadanie vozmusheniya

Tv=zeros(1,n);

tau1=4000;

tau2=7500;

tau3=10000;

t1=[0:dt:tau1];

t2=[0:dt:tau2];

t3=[0:dt:tau3];

n1=length(t1);

n2=length(t2);

n3=length(t3);

for i=1:n1

Tv(i)=(45/4000)*i*dt+15;

end

for i=(n1+1):n2

Tv(i)=60;

end

for i=(n2+1):n3

Tv(i)=-0.008*i*dt+120;

end

for i=2:n

% skorost reakcii

N(1,i)=A*exp(-E/(R*T(1,i-1)))*ca(1,i-1)*cb(1,i-1);

% pokomponentniy materialniy balans dlya reaktora

ca(1,i)=-N(1,i)*dt+ca(1,i-1); % konzentraziya komponenta "A" v reaktore

cb(1,i)=-N(1,i)*dt+cb(1,i-1); % konzentraziya komponenta "B" v reaktore

cc(1,i)=N(1,i)*dt+cc(1,i-1); % konzentraziya komponenta "C" v reaktore

% teplovoy balans (T-reaktor, Tt-rubashka)

T(1,i)=(K*s*(Tt(1,i-1)-Tv(1,i-1))+deltaH*N(1,i)*v1)*dt/(ro1*v1*cp1)+T(1,i-1);

Tt(1,i)=Tt(1,i-1)+dt*(Gt*ro2*cp2*(a1*(Tt1-Tt(1,i-1))+a2*(Tt2-Tt(1,i-1)))-K*s*(Tt(1,i-1)-Tv(1,i-1)))/(ro2*v2*cp2);

end

% postroenie grafikov

figure(1)

plot(t,T,'b');

hold on;

plot(t,Tt,'r');

hold on;

plot(t,Tv,'g');

legend('T','Tt','Tv');

hold on;

grid on;

title('GRAFIKI TEMPERATUR');

figure(2)

plot(t,ca,'b');

hold on;

plot(t,cb,'r');

hold on;

plot(t,cc,'g');

legend('ca','cb','cc');

hold on;

grid on;

title('GRAFIKI KONCENTRACIY');