Code covered by the BSD License

# Clacar

by

### Lacika (view profile)

descompunere_CaCO3

[sys x0]=descompunere_CaCO3(t,x,u,flag)
```function [sys x0]=descompunere_CaCO3(t,x,u,flag)
if abs(flag)==1

% Intrari
Qg=u(1);
Qs=u(2);
Qc=u(3);
Ts0=u(4);
TgL=u(5);
xCaCO3_0=u(6);
xCaO_0=0.1;
xCO2_L=0.01;

% Iesiri (necunoscute)
for i=1:30
xCaCO3(i)=x(i);
xCaO(i)=x(30+i);
xCO2(i)=x(60+i);
Ts(i)=x(90+i);
Tg(i)=x(120+i);
end;

% Constantele proceseului
M_CaCO3=100;                    %[kg/kmol]
M_CaO=56;
M_CO2=44;
R=2*4.185;                      % [J/(kmolK)]
Ea=0.74*10^6;                   % [kJ/kmol]
A_CaCO3=1.25*3.5*1.26*10^26;    % factor preexponential, legea Arrhenius [1/s]
dH_CaCO3=420*4.185;             % [kJ/kg]
Cpg=0.25*4.185;
Cps=0.35*4.185;                 % [J/kg]
Kt=1/3600*440*4.185;            % [kJ/(sm2K)]
At=300;                         % [m2]
ro=2830;                        % [kg/m3]
rog=1.2;
raport_ro=1.75;

% Parametrii cuptorului
hr=6;                           % inaltimea zonei reactie [m]
d=4.5;                          % diametrul cuptorului [m]
vg=1.5;                         % viteza gazului [m/s]
vs=0.4062/3600;                 % viteza solidului [m/s]

% Ecuatii algebrice
dz=hr/30;                       % grosimea elementului dz
A=pi*d^2/4;

% Ecuatiile diferentiale
k=A_CaCO3*exp(-Ea/(R*Ts(1)));
Hgs=At*Kt*(Tg(1)-Ts(1)) /A;
Hr=k*ro*xCaCO3(1)/M_CaCO3*dH_CaCO3;
Hc=Qc*5600*4.185/(hr*A);

dxCaCO3(1)=-vs*(xCaCO3(1)-xCaCO3_0)/dz-xCaCO3(1)*k;
dxCaO(1)=-vs*(xCaO(1)-xCaO_0)/dz + xCaCO3(1)*k*(M_CaO/M_CaCO3)*raport_ro;
dxCO2(1)=-vg*(xCO2(1)-xCO2(2))/dz +xCaCO3(1)*k*(M_CO2/M_CaCO3)*ro/rog;
dTs(1)=-vs*(Ts(1)-Ts0)/dz + vs/(Qs*Cps)*(Hgs-Hr);
dTg(1)=-vg*(Tg(1)-Tg(2))/dz-vg/(Qs*Cps)*(Hgs-Hc);

for i=2:29
k=A_CaCO3*exp(-Ea/(R*Ts(i)));
Hgs=At*Kt*(Tg(i)-Ts(i))/A;
Hr=k*ro*xCaCO3(i)/M_CaCO3*dH_CaCO3;
Hc=Qc*5600*4.185/(hr*A);

dxCaCO3(i)=-vs*(xCaCO3(i+1)-xCaCO3(i-1))/(2*dz)-xCaCO3(i)*k;
dxCaO(i)=-vs*(xCaO(i+1)-xCaO(i-1))/(2*dz)+xCaCO3(i)*k*(M_CaO/M_CaCO3)*raport_ro;
dxCO2(i)=-vg*(xCO2(i-1)-xCO2(i+1))/(2*dz)+xCaCO3(i)*k*(M_CO2/M_CaCO3)*ro/rog;
dTs(i)=-vs*(Ts(i+1)-Ts(i-1))/(2*dz) + vs/(Qs*Cps)*(Hgs-Hr);
dTg(i)=-vg*(Tg(i-1)-Tg(i+1))/(2*dz)-vg/(Qs*Cps)*(Hgs-Hc);
end;

k=A_CaCO3*exp(-Ea/(R*Ts(30)));
Hgs=At*Kt*(Tg(30)-Ts(30))/A;
Hr=k*ro*xCaCO3(30)/M_CaCO3*dH_CaCO3;
Hc=Qc*5600*4.185/(hr*A);

dxCaCO3(30)=-vs*(xCaCO3(30)-xCaCO3(29))/dz-xCaCO3(30)*k;
dxCaO(30)=-vs*(xCaO(30)-xCaO(29))/dz+xCaCO3(30)*k*(M_CaO/M_CaCO3)*raport_ro;
dxCO2(30)=-vg*(xCO2(30)-xCO2_L)/dz+ xCaCO3(30)*k*(M_CO2/M_CaCO3)*ro/rog;
dTs(30)=-vs*(Ts(30)-Ts(29))/dz+vs/(Qs*Cps)*(Hgs-Hr);
dTg(30)=-vg*(Tg(30)-TgL)/dz-vg/(Qs*Cps)*(Hgs-Hc);

sys=[dxCaCO3 dxCaO dxCO2 dTs dTg];
elseif flag==3
sys=x;
elseif flag==0
sys=[150 0 150 6 0 0];
%for i=1:30
%    x1(i)=0.9-(0.9-0.07)*i/30;
%    x2(i)=0.01+(0.89)*i/30;
%    x3(i)=0.38-(0.37)*i/30;
%    x4(i)=1033+100*i/30;
%    x5(i)=1083+100*i/30;
% end;