html('<h2>Solución</h2>')
t_L_entra = 55
t_L_sale = 25
t_G_entra = 23
Y_R_entra = 0.2500
q_G_entra = 2000
L_entra = 950
K_Y_a_div_rho_G = 0.2*3600
p_t = 97.0
t, L = var ('t, L')
M_A = 18.01
M_B = 28.84
t_0 = 0
R = 8.314
p_A_sat(t) = 10^(10.79586*(1-273.16/(t+273.15)) + 5.02808*log(273.16/(t+273.15),10) + 1.50474*10^-4 * (1-10^(-8.29692*((t+273.15)/273.16 - 1))) + 4.2873*10^-4 * (10^(4.76955*(1-273.16/(t+273.15))) - 1) + 2.786118312) / 1000
C_pBG(t) = (1045.356 - 0.3161783*(t+273.15) + 7.083814*10^-4*(t+273.15)^2 - 2.705209*10^-7*(t+273.15)^3) / 1000
C_pAG(t) = (1360.5 + 2.31334*(t+273.15) - 2.46784*10^-10*(t+273.15)^5 + 5.91332 * 10^-13*(t+273.15)^6) / 1000
C_pAL(t) = (8155.99 - 28.0627*(t+273.15) + 0.0511283*(t+273.15)^2 - 2.17582*10^-13 * (t+273.15)^6)/1000
lambda_A(t) = (3483181.4 - 5862.7703*(t+273.15) + 12.139568*(t+273.15)^2 - 0.0140290431*(t+273.15)^3) / 1000
lambda_A0 = lambda_A(t_0)
C_s(t, Y) = C_pBG(t) + C_pAG(t)*Y
Y_s(t, p) = M_A*p_A_sat(t)/(M_B*(p-p_A_sat(t)))
H(t,Y) = C_s(t, Y) * (t-t_0) + lambda_A(t_0) * Y
H_s(t,p) = C_s(t, Y_s(t,p)) * (t-t_0) + lambda_A(t_0) * Y_s(t, p)
v_H(t,Y,p) = R*(t+273.15)/p * (1/M_B+Y/M_A)
rho_G(T, Y, p) = p*(M_B+M_A*Y)/(R*(T+273.15))
p_A_entra = Y_R_entra*p_A_sat(t_G_entra)
Y_entra = M_A/M_B * p_A_entra/(p_t-p_A_entra)
C_s_entra = C_s(t_G_entra, Y_entra)
H_entra = H(t_G_entra,Y_entra)
v_H_entra = v_H(t_G_entra, Y_entra, p_t)
rho_G_entra = rho_G(t_G_entra, Y_entra, p_t)
G_s = q_G_entra/v_H_entra
print "pA_entra = %8.3g kPa" %p_A_entra
print "Y'entra = %8.3g kg/kg" %Y_entra
print "Cs_entra = %8.4g kJ/(kg K)" %C_s_entra
print "H'entra = %8.3g kJ/kg" %H_entra
print "vH_entra = %8.3g m³/kg" %v_H_entra
print "G's = %8.4g kg/(m² h)" %G_s
H_G(t,L,Gs) = L*C_pAL(t)/Gs * (t-t_L_sale) + H_entra
H_sale = H_G(t_L_entra, L_entra, G_s)
LdivG_lim = var('LdivG_lim')
Relac_lim(t) = solve(H_G(t,LdivG_lim,1) == H_s(t,p_t), LdivG_lim)[0].rhs()
sol_Rel_lim = find_local_minimum(Relac_lim(t), t_L_sale, t_L_entra)
G_s_min = L_entra/sol_Rel_lim[0]
t_eq = sol_Rel_lim[1]
N = 30
H_eq_menos_H(t) = H_s(t,p_t)-H_G(t,L_entra,G_s)
Delta_H = (H_sale-H_entra)/N
lista_H_G = [H_entra + Delta_H*i for i in range(N+1)]
lista_t_I = [find_root(H_G(t,L_entra,G_s) == H_g, t_L_sale-0.1, t_L_entra+0.1) for H_g in lista_H_G]
lista_t_eq = [find_root(H_s(t,p_t) == H_g, 0, t_L_entra+0.1) for H_g in lista_H_G]
lista_H_I = [H_s(t_I,p_t) for t_I in lista_t_I]
lista_t_G = vector(RR, N+1)
lista_t_G[0] = t_G_entra
if lista_t_I[0] == lista_t_G[0]:
DH_DtG = 10^100
else:
DH_DtG = (lista_H_I[0]-lista_H_G[0])/(lista_t_I[0]-lista_t_G[0])
for i in range(1,N+1):
lista_t_G[i] = lista_t_I[i-1] - (lista_H_I[i-1] - lista_H_G[i])/DH_DtG
if lista_t_G[i] < lista_t_eq[i]:
lista_t_G[i] = lista_t_eq[i]
DH_DtG = (lista_H_I[i]-lista_H_G[i])/(lista_t_I[i]-lista_t_G[i])
t_G_sale = lista_t_G[N]
Y_sale = numerical_approx(solve(H(t_G_sale,Y) == lista_H_G[N], Y)[0].rhs())
Yn_sale = Y_sale * M_B/M_A
y_sale = Yn_sale/(1+Yn_sale)
p_A_sale = p_t*y_sale
Y_R = p_A_sale/p_A_sat(t_G_sale)
g(t) = derivative(H_G(t,L_entra, G_s)) / (H_s(t,p_t)-H_G(t,L_entra, G_s))
integral_g = numerical_integral(g, t_L_sale, t_L_entra)
graf_NTU = plot(g, xmin = t_L_sale, xmax = t_L_entra, axes_labels = ['$t$/($\degree$C)', r'$\frac{d H\'(t)}{H^{*}-H}$'], fill = 'axis', color = 'red')
N_tOG = integral_g[0]
H_tOG_entra = G_s*(1+Y_entra)/(K_Y_a_div_rho_G*rho_G_entra)
rho_G_sale = rho_G(t_G_sale, Y_sale, p_t)
H_tOG_sale = G_s*(1+Y_sale)/(K_Y_a_div_rho_G*rho_G_sale)
H_tOG_prom = mean([H_tOG_entra, H_tOG_sale])
graf_Hs = plot(H_s(t,p_t), xmin = min(t_L_sale,t_G_entra) , xmax = t_L_entra+5, legend_label = "$H\'_s(t)$", axes_labels=['$t$','$H\'_G$'])
graf_Gmin = line([(t_L_sale, H_entra),(t_L_entra, H_G(t_L_entra,Relac_lim(t_eq),1))], color = 'red', legend_label = "$\\frac{L C_{p,A,L}}{G\'_{s,min}}$")
graf_lop = line([(t_L_sale, H_entra),(t_L_entra, H_G(t_L_entra,L_entra,G_s))], color = 'purple', legend_label = "$H\'_G (t_L)$") \
graf_HG_tG = line(zip(lista_t_G,lista_H_G), color = 'green', legend_label = "$H\'_G (t_G)$")
graf_1 = graf_Hs + graf_lop + graf_HG_tG
graf_2 = graf_Hs + graf_lop + graf_Gmin
print "\n","a) Condiciones del aire de salida"
graf_1.show(figsize=(6,4))
graf_1.show(xmin = min(t_L_sale,t_G_entra), xmax = t_G_sale+0, ymax = H_s(t_G_sale+2,p_t), figsize=(6,4))
print "\n"
print "Y'sale = %8.3g kg/kg" %Y_sale
print "YR = %8.3g %%" %(Y_R*100)
print "H'sale = %8.3g kJ/kg" %H_sale
print "tG_sale = %8.3g °C (con %s deltas, ΔH = %.4g kJ/kg)" %(t_G_sale, N, Delta_H)
print "\n","b) Altura de relleno"
graf_NTU.show(figsize=(6,4))
print "\n","N_tOG = %s" %numerical_approx(N_tOG, digits=4)
print "H_tOG,entra = %s m" %numerical_approx(H_tOG_entra, digits=4)
print "H_c = %s m" %numerical_approx(N_tOG * H_tOG_entra, digits=4), "<--- Tomando solo H_tOG a la entrada"
print "H_tOG,sale = %s m" %numerical_approx(H_tOG_sale, digits=4)
print "H_tOG,prom = %s m" %numerical_approx(H_tOG_prom, digits=4)
print "H_c = %s m" %numerical_approx(N_tOG * H_tOG_prom, digits=4), "<--- Considerando H_tOG promedio"
print "\n","c) Relación Gs/Gsmin"
graf_2.show()
print "\n","G's/G's,mín = %s" %numerical_approx(G_s/G_s_min, digits=4), "<--- Respuesta"
print "G's,mín = %s kg/(m² h)" %numerical_approx(G_s_min, digits=4)