%Codigo para del análisis termico de CCP con flujo constante de calor
%Especificaciones del tubo absorbedor:
%Flujo masico (kg/s)
m=5.85e-3;
%Diametro externo (m)
Do=0.0254;
%Diametro interno (m)
Di=0.0234;
%Longitud (m)
L=2.44;
%Absortancia
abs=0.9;
%Emitancia
em=0.9;
%Temperatura de entrada (K)
Ti=295.15;
%Temperatura de salida (K)
%To= Variables a calcular
%Area externa del receptor
Aro=pi*Do*L;
%Area interna del receptor
Ari=pi*Di*L;
%Especificaciones del concentrador,misma longitud del receptor:
%Foco de la parabola (m)
f=0.266;
%ángulo de borde (grados)
rad=360/(2*pi);
ab=90/rad;
%Reflectancia
Refle=0.92;
%Apertura de la parabola (m)
W=4*f*tan(ab/2);
%Area de captación del concentrador (m^2)
Ac=W*L;
%Concentración geometrica
C=Ac/Aro;
%Radio maximo de concentrador (m)
RMC=(2*f)/(1+cos(ab));
%Otros parametros opticos
%Errores ópticos (rad)
sig=0.011;
%Error de desalineación reflector (grados)
bet0=0.5;
bet=bet0/rad;
%Error de enfoque (m)
dr=0.0075;
%Factor de intercepción
sigma=sig*C;
beta=bet*C;
d=dr/Do;
syms psii
M=((sin(ab)*(1+ cos(psii))*(1-(2*d)*sin(psii)))-((pi*beta)*(1+cos(ab))))/((sqrt(2)*pi*sigma)*(1+cos(ab)));
N=((sin(ab)*(1+ cos(psii))*(1+(2*d)*sin(psii)))+((pi*beta)*(1+cos(ab))))/((sqrt(2)*pi*sigma)*(1+cos(ab)));
q=double(int(((erf(M)-erf(-N))/(1+cos(psii))),0,ab));
gama=((1 + cos(ab))/(2*sin(ab)))*q
%Condiciones ambientales
%Temparatura ambiente (K)
Ta=27+273.15;
%Irradiación directa (W/m^2)
Gb=950;
%Velocidad del viento (m/s)
Vv=1.25;
%Constantes
%Stefan-Bolzmant (W/m^2K^4)
SB=5.67e-8;
%Balance de energía en el receptor
% Qu=mCp(To-Ti); Considerando que Cp es constante, el la ecuacion de anterior
% tenemos que hay cuatro posibles variables, si fijamos Qu,m y Ti podemos
% calcular To.
%Propiedad del fluido de trabajo agua
%Capacidad calorifica (J/kg-K)
cpf=4209;
%Densidad (kg/m3)
denf=983.2;
%Viscosidad (kg/m-s)
visf=0.000306;
%Conductividad termica del fluido (W/m.K)
kf=0.677;
%Numero de Prandt
Prf=(cpf*visf)/kf;
%Propiedad del aire
%Capacidad calorifica (J/kg-K)
cpa=1;
%Densidad (kg/m3)
dena=1;
%Viscosidad (kg/m-s)
visa=1;
%Conductividad termica del fluido (W/m.K)
ka=26.3e-3;
%Numero de Prandt
Pra=(15.89e-6)/(22.5e-6);
syms T
Ts=200
%Coeficiente radiación
hr=4*SB*em*(Ts^3);
%Calculamos el coeficiente de tranferencia del fluido
%Numero Reynolds
Ref=(4*m)/(pi*visf*Di);
%Numero de Nusselt
if Ref>4000
Nu=0.23*(Ref^.8)*(Prf^.4);
else
Nu=4.36;
end
hf=(kf/Di)*Nu;
%Coeficiente convectivo del viento
Rev=(Vv*Do)/15.89e-6;
if Rev>=0.4 & Rev<4
B=0.989;
n=0.33;
elseif Rev>=4 & Rev<40
B=0.911;
n=0.385;
elseif Rev>=40 & Rev<4000
B=0.683;
n=0.466;
elseif Rev>=4000 & Rev<40000
B=0.193;
n=0.618;
elseif Rev>=40000 & Rev<400000
B=0.027;
n=0.805;
end
Nuv=(B)*(Rev^n)*(Pra^(1/3))*(0.95);
hv=Nuv*(ka/Do);
%Calor util
Qu=(gama*abs*Ac*Refle*Gb)-(Aro*(hr+hv)*(Ts-Ta));
%Calculo de Ts
qs=Qu/Ari;
To=Ti+(Qu/(m*cpf));
%Se estima la Tsii y Tsio
Tsii=Ti+(qs/hf);
Tsio=To+(qs/hf);
T=(Tsii+Tsio)/2;
while T~=Ts
Ts=T;
%Coeficiente radiación
hr=4*SB*em*(Ts^3);
%Calculamos el coeficiente de tranferencia del fluido
%Numero Reynolds
Ref=(4*m)/(pi*visf*Di);
%Numero de Nusselt
if Ref>4000
Nu=0.23*(Ref^.8)*(Prf^.4);
else
Nu=4.36;
end
hf=(kf/Di)*Nu;
%Coeficiente convectivo del viento
Rev=(Vv*Do)/15.89e-6;
if Rev>=0.4 & Rev<4
B=0.989;
n=0.33;
elseif Rev>=4 & Rev<40
B=0.911;
n=0.385;
elseif Rev>=40 & Rev<4000
B=0.683;
n=0.466;
elseif Rev>=4000 & Rev<40000
B=0.193;
n=0.618;
elseif Rev>=40000 & Rev<400000
B=0.027;
n=0.805;
end
Nuv=(B)*(Rev^n)*(Pra^(1/3))*(0.95);
hv=Nuv*(ka/Do);
%Calor util
Qu=(gama*abs*Refle*Ac*Gb)-(Aro*(hr+hv)*(Ts-Ta))
%Calculo de Ts
qs=Qu/Ari;
To=Ti+(Qu/(m*cpf));
%Se estima la Tsii y Tsio
Tsii=Ti+(qs/hf);
Tsio=To+(qs/hf);
T=(Tsii+Tsio)/2;
end
disp(T)
