محاسبات تقطیر در سیستم های بسته دو مخزنی در MATLAB
کد:
% Fil dist3.m
% Case 1: Konstant holdup
% zf=0.5;[t,x]=ode15s('dist3',[0 10],[zf*ones(1,22) 5 5 0])
% Case 2: PI-regulering topp
% Case 3: PI-regulering midt
% Case 4: P-regulering midt
% zf=0.5;[t,x]=ode15s('dist3',[0 10],[zf*ones(1,22) 0.01 9.99 0])
% Tid for aa naa gitt sammensetning i topp OG bunn, sammenligne anal.
% se litt paa effekt av PI-tuninger (KC, taui)
% se paa ulike verdier av zf og alpha
function xdot = dist3(t,x)
% Konstant relativ flyktighet
alfa = 2;
% Antall trinn
NT=22;
% Foedetrinn
NF=11;
% Oevrige data
F=0;
qF= 0;
zF=0.5;
V=10;
L0=10;
D=0;
B=0;
MD=x(NT+1);
MB=x(NT+2);
MK=0.0005;
KC=20;
%KC=0;
XDs=0.95;
XFs=0.5;
% Regulator
%err = x(NT)-XDs; % Case 2
err=x(NF)-XFs; % Case 3,4
LT=L0-KC*err; % Case 4
L0I=0; taui=0.05;
interr = x(NT+3);
%LT = L0I - KC*err - KC/taui*interr; % Case 2,3
%LT=10; % Case 1
% Damp-vaeske likevekt
i=1:NT-1;
y(i)=alfa*x(i)./(1+(alfa-1)*x(i));
% Dampstroem
i=1:NT-1;
V(i)=V*ones(1,NT-1);
i=NF:NT-1;
V(i)=V(i)+(1-qF)*F;
% Vaeskestroem
i=2:NT;
L(i)=LT*ones(1,NT-1);
i=2:NF;
L(i)=L(i)+qF*F;
% Total kondensator
y(NT)=x(NT);
% Holdup i kolonnen
i=2:NT-1;
M(i)=MK*ones(1,NT-2);
% Holdup i kokeren
MBdot=L(2)-V(1)-B;
% Holdup i kondensatoren
MDdot=V(NT-1)-L(NT)-D;
% Dynamiske ligninger (komponentbalanser)
% Kokeren
Mxdot(1) = L(2)*x(2)-V(1)*y(1)-B*x(1);
xdot(1)=(Mxdot(1)-x(1)*MBdot)/MB;
% Kolonnen
i=2:NF-1;
xdot(i)=(L(i+1).*x(i+1)' - L(i).*x(i)' + V(i-1).*y(i-1) - V(i).*y(i))./M(i);
i=NF+1:NT-1;
xdot(i)=(L(i+1).*x(i+1)' - L(i).*x(i)' + V(i-1).*y(i-1) - V(i).*y(i))./M(i);
% Foedetrinn
xdot(NF) =(L(NF+1)*x(NF+1)+V(NF-1)*y(NF-1)-L(NF)*x(NF)-V(NF)*y(NF)+F*zF)/MK;
% Kondensatoren
Mxdot(NT) = V(NT-1)*y(NT-1)-L(NT)*x(NT)-D*x(NT);
xdot(NT)=(Mxdot(NT)-x(NT)*MDdot)/MD;
% Funksjonen f(t,x) maa returnere en kolonnevektor
xdot(NT+1)=MDdot;
xdot(NT+2)=MBdot;
%Integrator in controller
xdot(NT+3)=err;
xdot=xdot(:);