error code drying of solid
Show older comments
when running the program it reflects the error to me
code
%% PARAMETROS DE ENTRADA AL MODELO
a=0.2755; b=-0.34; D0=0.008 ;r0=D0/2; as=3/r0;
%Cálculo del flujo másico de aire
l1=0.53;
l2=0.51;
A=l1*l2;
L1=100/100;
V=A*L1;
v=4; %Velocidad del aire
ma=0.028;
Ga=(ma/A)*3600;
N=4;
L=0.5;
x=L/N;
h=0.5;
tfinal=24*4.5;
Mi=0.2409
v=1;
Sp=0.685*0.51;
%Coeficientes de transferencia de calor por conducción
lamb_madera=0.13;
lamb_vidrio=0.05;
lamb_yeso=0.35;
em=2.5/100;
ev=6/100;
ey=1.5/100;
hcp=(lamb_madera/em)+(lamb_vidrio/ev)+(lamb_yeso/ey);
hvv=(5.67+3.86*v);
hvfpi=4.1;
%Coeficiente de Fourier
a0=29.86
a11=-2.997
b1 =-6.595
a22 =0.2317
b2 =1.949
a33 =-0.2656
b3 =-0.4193
a4 =-0.3443
b4 =-0.3929
a5 =0.2462
b5 =0.3422
a6 =-0.1787
b6 =-0.4957
a7 =-0.3067
b7 =0.2473
w =0.1488
%Temperatura de ambiente
a0a=25.91
a1a=-2.422
b1a =-3.467
a2a =0.4955
b2a =1.07
a3a =-0.4853
b3a =0.4878
a4a =0.1001
b4a =0.07658
a5a =0.2295
b5a =-0.04174
a6a =-0.0225
b6a =0.02728
a7a =-0.07581
b7a =-0.1319
a8a =-0.1274
b8a =0.1078
wa =0.1298
%% DETERMINACION DEL NUMERO DE NODOS PARA LA SIMULACION
N=L/x;
N=ceil(N);
X=tfinal/h;
X=ceil(X);
%% CONDICIONES INICIALES
M0=56;
M0=100*M0/(100-M0);
for n=2:N
M(n,13)=M0;
Tg(n,13)=20;
end
% load('Tf0.mat');
for i=14:X
H(1,i)=0.01;
end
for ia=14:(X)
T(1,ia)=a0+a11*cos(ia*w)+b1*sin(ia*w)+a22*cos(2*ia*w)+b2*sin(2*ia*w)+a33*cos(3*ia*w)+b3*sin(3*ia*w)+a4*cos(4*ia*w)+b4*sin(4*ia*w)+a5*cos(5*ia*w)+b5*sin(5*ia*w)+a6*cos(6*ia*w)+b6*sin(6*ia*w)+a7*cos(7*ia*w)+b7*sin(7*ia*w);
Ta(1,ia)=a0a+a1a*cos(ia*wa)+b1a*sin(ia*wa)+a2a*cos(2*ia*wa)+b2a*sin(2*ia*wa)+a3a*cos(3*ia*wa)+b3a*sin(3*ia*wa)+a4a*cos(4*ia*wa)+b4a*sin(4*ia*wa)+a5a*cos(5*ia*wa)+b5a*sin(5*ia*wa)+a6a*cos(6*ia*wa)+b6a*sin(6*ia*wa)+a7a*cos(7*ia*wa)+b7a*sin(7*ia*wa)+a8a*cos(8*ia*wa)+b8a*sin(8*ia*wa);
Tsky(1,ia)=Ta(1,ia)-6;
end
%% CALCULO DE LOS PARAMETROS INTERMEDIOS y H, M,T,Tg
for j=14:X
for n=2:N
end
%Calculo de M(n,j)
RH(n-1,j)=1.626*H(n-1,j)*(10^5)/((1+1.608*H(n-1,j))*6892.4*exp(54.633-12301.69/(1.8*T(n-1,j)+492)-5.169*log(1.8*T(n-1,j)+492)));
RH(n-1,j)=abs(RH(n-1,j));
Me(n-1,j)=100*(0.01087324+0.23758842*RH(n-1,j)-0.63848052*(RH(n-1,j)^2)+0.56197989*(RH(n-1,j)^3))*exp((0.02965424-0.33506858*RH(n-1,j)+0.95779665*(RH(n-1,j)^2)-1.27182954*(RH(n-1,j)^3)+0.573953345*(RH(n-1,j)^4))*(T(n-1,j)-57.2913401));
Me(n-1,j)=abs(Me(n-1,j));
K(n-1,j)=(pi*pi/(9))*as*as*3600*8.64*(10^-4)*exp(-39.94*1000/(8.314*(Tg(n,j-1)+273)));
K(n-1,j)=abs(K(n-1,j));
k1(n-1,j)=-K(n-1,j)*(M(n,j-1)-Me(n-1,j));
k2(n-1,j)=-K(n-1,j)*((M(n,j-1)+0.05*k1(n-1,j))-Me(n-1,j));
k3(n-1,j)=-K(n-1,j)*((M(n,j-1)+0.05*k2(n-1,j))-Me(n-1,j));
k4(n-1,j)=-K(n-1,j)*((M(n,j-1)+0.1*k3(n-1,j))-Me(n-1,j));
M(n,j)=M(n,j-1)+0.0167*(k1(n-1,j)+2*k2(n-1,j)+2*k3(n-1,j)+k4(n-1,j));
dM(n,j)=0.0167*(k1(n-1,j)+2*k2(n-1,j)+2*k3(n-1,j)+k4(n-1,j));
%Calculo de H(n,j)
rhop(n,j-1)=726.4-4.135*-4.135*(M(n,j-1)/(100+M(n,j-1)));
rhop(n,j-1)=abs(rhop(n,j-1));
H(n,j)=H(n-1,j)-0.01*(rhop(n,j-1)*x/(Ga*2*h))*(M(n,j)-M(n,j-1));
%Calculo de T(n,j)
Ca(n-1,j)=1003.4+0.178*T(n-1,j);
Ca(n-1,j)=abs(Ca(n-1,j));
Cv(n-1,j)=1859+0.236*T(n-1,j);
Cv(n-1,j)=abs(Cv(n-1,j));
ua(n-1,j)=(1.691*(10^-5)+4.987*(10^-8)*T(n-1,j)-3.187*(10^-11)*(T(n-1,j)^2)+1.319*(10^-14)*(T(n-1,j)^3));
ua(n-1,j)=abs(ua(n-1,j));
ha(n-1,j)=172.2*(Ga^0.5)*as;
ha(n-1,j)=abs(ha(n-1,j));
hrc(1,j)=0.93*(5.670373*(10^-8))*(T(n-1,j)+Tsky(1,j))*((T(n-1,j)^2)+(Tsky(1,j)^2));
Tap(1,j)=(Ta(1,j)*hvv+Tsky(1,j)*hrc)/(hvv+hrc);
Tpi(1,j)=((hvfpi*(hcp+hrc(1,j)+hvv)*T(n-1,j))+(hcp*(hrc(1,j)+hvv)*Tap(1,j)))/((hcp*(hrc(1,j)+hvv))+(hvfpi*(hcp+hrc(1,j)+hvv)));
per=x*3600*(hvfpi*Sp*(T(n-1,j)-Tpi(1,j)))/(Ga*Ca(n-1,j)+Ga*H(n,j)*Cv(n-1,j));
T(n,j)=T(n-1,j)-(ha(n-1,j)*x/(Ga*Ca(n-1,j)+Ga*H(n,j)*Cv(n-1,j)))*(T(n-1,j)-Tg(n,j-1));
%Calculo de Tg(n,j)
Cp(n,j-1)=1755+23.45*0.01*M(n,j-1);
Cp(n,j-1)=abs(Cp(n,j-1));
hfg(n,j-1)=2501.6+(1840-4180)*(Tg(n,j-1)+273);
hfg(n,j-1)=abs(hfg(n,j-1));
Cw=4200;
Tg(n,j)=Tg(n,j-1)+(10/6)*(h*ha(n-1,j)/(rhop(n,j-1)*Cp(n,j-1)+rhop(n,j-1)*Cw*0.01*M(n,j-1)))*(T(n-1,j)-Tg(n,j-1))-(10/6)*Ga*h*(hfg(n,j-1)+(T(n,j)-Tg(n,j-1))*Cv(n-1,j))*(H(n,j)-H(n-1,j))/(x*(rhop(n,j-1)*Cp(n,j-1)+rhop(n,j-1)*Cw*0.01*M(n,j-1)));
LV(n,j)=hfg(n,j-1)+(T(n,j)-Tg(n,j-1))*Cv(n-1,j);
end
t=6:h:(tfinal-0.5);
t1=6:h:(tfinal-1);
M4=M(4,:);
M4=M4(13:1:X);
M3=M(3,:);
M3=M3(13:1:X);
M2=M(2,:);
M2=M2(13:1:X);
plot(t,M4,'b','linewidth',3)
hold on
plot(t,M3,'r','linewidth',3)
plot(t,M2,'k','linewidth',3)
legend('Bandeja 3','Bandeja 2','Bandeja 1')
xlabel('tiempo(hrs del dia)')
ylabel('Humedad del cacao(%db)')
grid on
figure(2)
dM1=dM(2,:);
dM1=dM1(2:1:X);
divM=M/M0;
divM1=divM(2,:);
divM1=divM1(1:1:X-1);
plot(divM1(13:1:X-1),-dM1(13:1:X-1)/0.5,'k','linewidth',3)
hold on
dM2=dM(3,:);
dM2=dM2(2:1:X);
divM2=M/M0;
divM2=divM2(3,:);
divM2=divM2(1:1:X-1);
plot(divM2(13:1:X-1),-dM2(13:1:X-1)/0.5,'r','linewidth',3)
dM3=dM(4,:);
dM3=dM3(2:1:X);
divM3=M/M0;
divM3=divM3(4,:);
divM3=divM3(1:1:X-1);
plot(divM3(13:1:X-1),-dM3(13:1:X-1)/0.5,'b','linewidth',3)
legend('Bandeja 1','Bandeja 2','Bandeja 3')
xlabel('M/M0 (db)')
ylabel('-dM/dt (kg de agua/kg de cacao *h^-1)')
grid on
figure(3)
dM11=dM(2,:);
dM11=dM11(14:1:X);
dM22=dM(3,:);
dM22=dM22(14:1:X);
dM33=dM(4,:);
dM33=dM33(14:1:X);
plot(t1,-dM11/0.5,'k','linewidth',3)
hold on
plot(t1,-dM22/0.5,'r','linewidth',3)
plot(t1,-dM33/0.5,'b','linewidth',3)
legend('Bandeja 1','Bandeja 2','Bandeja 3')
xlabel('tiempo(hrs del dia)')
ylabel('-dM/dt (kg de agua/kg de cacao *h^-1)')
grid on
figure(4)
plot(H(2,:),'.k')
hold on
plot(H(3,:),'.r')
plot(H(4,:),'.b')
legend('layer 1','layer 2','layer 3','layer 4')
xlabel('tiempo(hrs)')
ylabel('Humedad especifica del aire(kg/kg)')
grid on
figure(5)
Tent=T(1,:);
Tent=Tent(14:1:X)
T1=T(2,:);
T1=T1(14:1:X);
T2=T(3,:);
T2=T2(14:1:X);
T3=T(4,:);
T3=T3(14:1:X);
plot(t1,Tent,'g','linewidth',2)
hold on
plot(t1,T1,'k','linewidth',2)
plot(t1,T2,'r','linewidth',2)
plot(t1,T3,'b','linewidth',2)
legend('Entrada a la camara','Bandeja 1','Bandeja 2','Bandeja 3')
xlabel('tiempo(hrs)')
ylabel('Temperatura del aire desecante (C)')
grid on
figure(6)
Tg1=Tg(2,:);
Tg1=Tg1(13:1:X)
Tg2=Tg(3,:);
Tg2=Tg2(13:1:X)
Tg3=Tg(4,:);
Tg3=Tg3(13:1:X)
plot(t,Tg1,'k','linewidth',3)
hold on
plot(t,Tg2,'r','linewidth',3)
plot(t,Tg3,'b','linewidth',3)
legend('Bandeja 1','Bandeja 2','Bandeja 3')
xlabel('tiempo(hrs)')
ylabel('Temperatura del cacao (C)')
grid on
figure(7)
H1=H(2,:);
H1=H1(14:1:X);
H2=H(3,:);
H2=H2(14:1:X);
H3=H(4,:);
H3=H3(14:1:X);
plot(t1,H1,'k','linewidth',3)
hold on
plot(t1,H2,'r','linewidth',3)
plot(t1,H3,'b','linewidth',3)
legend('Bandeja 1','Bandeja 2','Bandeja 3')
xlabel('tiempo(hrs)')
ylabel('Humedad especifica del aire desecante')
grid on
axis([6 75 0.005 0.035])
error:
Index in position 1 exceeds array bounds (must not exceed 1).
Error in secado (line 94)
RH(n-1,j)=1.626*H(n-1,j)*(10^5)/((1+1.608*H(n-1,j))*6892.4*exp(54.633-12301.69/(1.8*T(n-1,j)+492)-5.169*log(1.8*T(n-1,j)+492)));
solutions please :(
Accepted Answer
More Answers (0)
Categories
Find more on Programming in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!