> restart:#"m06_p56"

Un depósito de 30 m3 aislado térmicamente contiene inicialmente amoníaco a 20 ºC en estado bifásico al 50% en volumen. A partir de un cierto instante se extrae vapor hasta que la temperatura interior baja a 10 ºC. Se pide:
a) Presión y masa iniciales en el interior, y situación del estado en el diagrama T-s.
b) Estado del vapor que se extrae, entalpía, y representación en los diagramas T-s y p-h.
c) Plantear el balance energético del depósito, indicando los datos y las incógnitas.
d) Presión y masa finales en el interior

Datos:

> read"../therm_eq.m":read"../therm_proc.m":with(therm_proc):

> su1:="NH3":dat:=[V=30*m_^3,Vv=15*m_^3,T1=(20+273)*K_,T2=(10+273)*K_];

`:=`(dat, [V = `+`(`*`(30, `*`(`^`(m_, 3)))), Vv = `+`(`*`(15, `*`(`^`(m_, 3)))), T1 = `+`(`*`(293, `*`(K_))), T2 = `+`(`*`(283, `*`(K_)))])

Esquema:

> `:=`(Sistemas, [elemento])

Image

> `:=`(Estados, [1 = inicial, 2 = final])

Eqs. const.:

> gdat:=get_gas_data(su1):ldat:=get_liq_data(su1):dat:=op(dat),gdat,ldat,Const,SI2,SI1:get_pv_data(su1):

a) Presión y masa iniciales en el interior, y situación del estado en el diagrama T-s.

> p1=p[v](T1);p1_:=subs(dat,pv(T1)):p1=subs(dat,evalf(%,3));ml:=rho*(V-Vv);ml_:=subs(dat,%):'ml'=evalf(%,4);rho=subs(dat,rho);mv:=p1*Vv/(R*T1);mv_:=subs(p1=p1_,dat,mv):'mv'=evalf(%,2);x1:='mv/(ml+mv)';x1_:=mv_/(ml_+mv_):'x1'=evalf(%,2);m1='ml+mv';m1_:=ml_+mv_:'m1'=evalf(%,3);

p1 = p[v](T1)

p1 = `+`(`*`(0.854e6, `*`(Pa_)))

`:=`(ml, `*`(rho, `*`(`+`(V, `-`(Vv)))))

ml = `+`(`*`(0.1046e5, `*`(kg_)))

rho = `+`(`/`(`*`(697., `*`(kg_)), `*`(`^`(m_, 3))))

`:=`(mv, `/`(`*`(p1, `*`(Vv)), `*`(R, `*`(T1))))

mv = `+`(`*`(89., `*`(kg_)))

`:=`(x1, `/`(`*`(mv), `*`(`+`(ml, mv))))

x1 = 0.85e-2

m1 = `+`(ml, mv)

m1 = `+`(`*`(0.105e5, `*`(kg_)))

Conviene comprobar que el modelo de gas ideal es adecuado.

> pR_:=evalf(subs(dat,p1_/p[cr]),3);TR_:=evalf(subs(dat,T1/T[cr]),3);

`:=`(pR_, 0.756e-1)

`:=`(TR_, .721)

y mirando en el diagrama generalizado Z(pR,TR) se obtiene Z1v=0,92, i.e.la masa de vapor sería 97 kg en vez de 89 kg, pero la incertidumbre es menor que en la modelización de la fase líquida.

b) Estado del vapor que se extrae, entalpía, y representación en los diagramas T-s y p-h.

Elijo la entalpía del líquido en el estado triple como referencia.

El estado irá variando desde el 1V al 2V al ir extrayendo. Los valores inicial y final serán:

> Tref_:=subs(dat,T[f]);hv1=href+hlvref+c[p]*(T1-Tref);T[b]=subs(dat,T[b]);h[lv0]=subs(dat,h[lv0]);hlvT:=hlv(T);hlvref_:=subs(dat,T=Tref_,dat,hlvT):'hlref'=evalf(%,3);hv1_:=subs(dat,hlvref_+c[p]*(T1-Tref_)):'hv1'=evalf(%,3);hv2_:=subs(dat,hlvref_+c[p]*(T2-Tref_)):'hv2'=evalf(%,3);p2=p[v](T2);p2_:=subs(dat,pv(T2)):p2=subs(dat,evalf(%,3));

`:=`(Tref_, `+`(`*`(195., `*`(K_))))

hv1 = `+`(href, hlvref, `*`(c[p], `*`(`+`(T1, `-`(Tref)))))

T[b] = `+`(`*`(239., `*`(K_)))

h[lv0] = `+`(`/`(`*`(0.1357e7, `*`(J_)), `*`(kg_)))

`:=`(hlvT, `+`(h[lv0], `*`(`+`(c[p], `-`(c)), `*`(`+`(T, `-`(T[b]))))))

hlref = `+`(`/`(`*`(0.146e7, `*`(J_)), `*`(kg_)))

hv1 = `+`(`/`(`*`(0.168e7, `*`(J_)), `*`(kg_)))

hv2 = `+`(`/`(`*`(0.166e7, `*`(J_)), `*`(kg_)))

p2 = p[v](T2)

p2 = `+`(`*`(0.608e6, `*`(Pa_)))

c) Plantear el balance energético del depósito, indicando los datos y las incógnitas.

A la vista de lo anterior, la entalpía del vapor puede considerarse constante (i.e. sus variaciones despreciables frente a los otros términos).

> eqBE:=d(m*u)=dQ+dW+hv*dm;eqBE:=m2*u2-m1*u1=hv*(m2-m1);m*u=m*h-p*V;eqBE:=m2*h2-p2*V-m1*h1+p1*V=hv*(m2-m1);h=href+c*(T-Tref)+x*'hlvT';hlvT1_:=subs(dat,T=T1,dat,hlv(T)):hlvT1=evalf(%,3);h1_:=subs(dat,c*(T1-Tref_)+x1_*hlvT1_):'h1'=evalf(%,2);hlvT2_:=subs(dat,T=T2,dat,hlv(T)):hlvT2=evalf(%,3);h2_:=subs(dat,c*(T2-Tref_)+x2*hlvT2_):'h2'=evalf(%,3);

`:=`(eqBE, d(`*`(m, `*`(u))) = `+`(dQ, dW, `*`(hv, `*`(dm))))

`:=`(eqBE, `+`(`*`(m2, `*`(u2)), `-`(`*`(m1, `*`(u1)))) = `*`(hv, `*`(`+`(m2, `-`(m1)))))

`*`(m, `*`(u)) = `+`(`*`(m, `*`(h)), `-`(`*`(p, `*`(V))))

`:=`(eqBE, `+`(`*`(m2, `*`(h2)), `-`(`*`(p2, `*`(V))), `-`(`*`(m1, `*`(h1))), `*`(p1, `*`(V))) = `*`(hv, `*`(`+`(m2, `-`(m1)))))

h = `+`(href, `*`(c, `*`(`+`(T, `-`(Tref)))), `*`(x, `*`(hlvT)))

hlvT1 = `+`(`/`(`*`(0.123e7, `*`(J_)), `*`(kg_)))

h1 = `+`(`/`(`*`(0.46e6, `*`(J_)), `*`(kg_)))

hlvT2 = `+`(`/`(`*`(0.125e7, `*`(J_)), `*`(kg_)))

h2 = `+`(`/`(`*`(0.405e6, `*`(J_)), `*`(kg_)), `/`(`*`(0.125e7, `*`(x2, `*`(`^`(m_, 2)))), `*`(`^`(s_, 2))))

i.e., hay 2 incógnitas en el balance energético, la masa total final, m2, y la fracción másica de vapor final, x2. Para resolver necesitamos otra ecuación: la del volumen.

d) Presión y masa finales en el interior

> eqV:=V/m2=vl2+x2*vlv2;x2:=(V/m2-vl2)/(vv2-vl2);vl2=1/rho;vl2:=subs(dat,1/rho):'vl2'=evalf(%,2);vv2=R*T2/p2;vv2:=subs(dat,R*T2/p2_):'vv2'=evalf(%,2);'eqBE'=eqBE;eqBE_:=subs(dat,m2*h2_-p2_*V-m1_*h1_+p1_*V=hv1_*(m2-m1_)):'eqBE_'=evalf(%,2);m2_:=subs(dat,solve(eqBE_,m2)):'m2'=evalf(%,4);Dm='m2-m1';Dm_:=m2_-m1_:Dm=evalf(%,2);x2_:=subs(m2=m2_,dat,x2):'x2'=evalf(%,2);Vl2='m2/rho';Vl2_:=subs(dat,m2_/rho):'Vl2'=evalf(%,4);m2l=m2_*(1-x2_):evalf(%,4);

`:=`(eqV, `/`(`*`(V), `*`(m2)) = `+`(vl2, `*`(x2, `*`(vlv2))))

`:=`(x2, `/`(`*`(`+`(`/`(`*`(V), `*`(m2)), `-`(vl2))), `*`(`+`(vv2, `-`(vl2)))))

vl2 = `/`(1, `*`(rho))

vl2 = `+`(`/`(`*`(0.14e-2, `*`(`^`(m_, 3))), `*`(kg_)))

vv2 = `/`(`*`(R, `*`(T2)), `*`(p2))

vv2 = `+`(`/`(`*`(.23, `*`(`^`(m_, 3))), `*`(kg_)))

eqBE = (`+`(`*`(m2, `*`(h2)), `-`(`*`(p2, `*`(V))), `-`(`*`(m1, `*`(h1))), `*`(p1, `*`(V))) = `*`(hv, `*`(`+`(m2, `-`(m1)))))

eqBE_ = (`+`(`*`(m2, `*`(`+`(`/`(`*`(0.40e6, `*`(J_)), `*`(kg_)), `/`(`*`(0.56e7, `*`(`+`(`/`(`*`(30., `*`(`^`(m_, 3))), `*`(m2)), `-`(`/`(`*`(0.14e-2, `*`(`^`(m_, 3))), `*`(kg_)))), `*`(kg_))), `*`(m...
eqBE_ = (`+`(`*`(m2, `*`(`+`(`/`(`*`(0.40e6, `*`(J_)), `*`(kg_)), `/`(`*`(0.56e7, `*`(`+`(`/`(`*`(30., `*`(`^`(m_, 3))), `*`(m2)), `-`(`/`(`*`(0.14e-2, `*`(`^`(m_, 3))), `*`(kg_)))), `*`(kg_))), `*`(m...

m2 = `+`(`*`(0.1015e5, `*`(kg_)))

Dm = `+`(m2, `-`(m1))

Dm = `+`(`-`(`*`(0.39e3, `*`(kg_))))

x2 = 0.68e-2

Vl2 = `/`(`*`(m2), `*`(rho))

Vl2 = `+`(`*`(14.56, `*`(`^`(m_, 3))))

m2l = `+`(`*`(0.1008e5, `*`(kg_)))

i.e. se han extraido 390 kg de vapor de amoníaco y el nivel de líquido bajará un poco.

Otro procedimiento de resolver este problema es asimilarlo a una expansión isoentrópica de la masa de control inicial, y luego descartar la parte que en realidad ha salido (ver Esquema equivalente).

Por cierto, este modelo enseña que el punto 2 está a la izquierda del punto 1 en el diagrama T-s, ya que la parte de vapor que se elimina tiene mayor entropía específica que el resto.

Tomando la misma referencia para las entropías que para las entalpías:

> eqBM:=m1=m2;eqBM:=m1l+m1v=m2l+m2v;eqBS:=S1=S2;eqBS:=m1l*s1l+m1v*s1v=m2l*s2l+m2v*s2v;sl:=sref+c*ln(T/Tref);s1l_:=subs(dat,evalf(subs(dat,c*ln(T1/Tref_)))):'s1l'=evalf(%,2);s2l_:=subs(dat,evalf(subs(dat,c*ln(T2/Tref_)))):'s2l'=evalf(%,2);sv:=sref+hlvref/Tref+c[p]*ln(T/Tref)-R*ln(p/pref);pref=p[v](Tref);pref_:=pv(Tref_):'pref'=evalf(%,2);s1v_:=subs(dat,evalf(subs(dat,sref=0,Tref=Tref_,hlvref=hlvref_,T=T1,p=p1_,pref=pref_,dat,sv))):'s1v'=evalf(%,2);s2v_:=subs(dat,evalf(subs(dat,sref=0,Tref=Tref_,hlvref=hlvref_,T=T2,p=p2_,pref=pref_,dat,sv))):'s2v'=evalf(%,2);eqBM_:=ml_+mv_=m2l+m2v:'eqBM'=evalf(%,2);eqBS_:=ml_*s1l_+mv_*s1v_=m2l*s2l_+m2v*s2v_:'eqBS'=evalf(%,2);sol_:=solve({eqBM_,eqBS_},{m2l,m2v}):m2l_:=evalf(subs(sol_,m2l),4);m2v_:=evalf(subs(sol_,m2v),3);V2l='m2l/rho';V2l_:=subs(sol_,dat,m2l/rho):'V2l'=evalf(%,4);

`:=`(eqBM, m1 = m2)

`:=`(eqBM, `+`(m1l, m1v) = `+`(m2l, m2v))

`:=`(eqBS, S1 = S2)

`:=`(eqBS, `+`(`*`(m1l, `*`(s1l)), `*`(m1v, `*`(s1v))) = `+`(`*`(m2l, `*`(s2l)), `*`(m2v, `*`(s2v))))

`:=`(sl, `+`(sref, `*`(c, `*`(ln(`/`(`*`(T), `*`(Tref)))))))

s1l = `+`(`/`(`*`(0.19e4, `*`(J_)), `*`(kg_, `*`(K_))))

s2l = `+`(`/`(`*`(0.17e4, `*`(J_)), `*`(kg_, `*`(K_))))

`:=`(sv, `+`(sref, `/`(`*`(hlvref), `*`(Tref)), `*`(c[p], `*`(ln(`/`(`*`(T), `*`(Tref))))), `-`(`*`(R, `*`(ln(`/`(`*`(p), `*`(pref))))))))

pref = p[v](Tref)

pref = `+`(`*`(0.59e4, `*`(Pa_)))

s1v = `+`(`/`(`*`(0.60e4, `*`(J_)), `*`(kg_, `*`(K_))))

s2v = `+`(`/`(`*`(0.61e4, `*`(J_)), `*`(kg_, `*`(K_))))

eqBM = (`+`(`*`(0.11e5, `*`(kg_))) = `+`(m2l, m2v))

eqBS = (`+`(`/`(`*`(0.20e8, `*`(J_)), `*`(K_))) = `+`(`/`(`*`(0.17e4, `*`(m2l, `*`(J_))), `*`(kg_, `*`(K_))), `/`(`*`(0.61e4, `*`(m2v, `*`(J_))), `*`(kg_, `*`(K_)))))

`:=`(m2l_, `+`(`*`(0.1007e5, `*`(kg_))))

`:=`(m2v_, `+`(`*`(472., `*`(kg_))))

V2l = `/`(`*`(m2l), `*`(rho))

V2l = `+`(`*`(14.45, `*`(`^`(m_, 3))))

pero ahora el vapor que realmente quedaría dentro sería hasta completar el volumen:

> m2v_dentro:=p2*(V-V2l)/(R*T2);m2v_dentro_:=subs(p2=p2_,V2l=V2l_,dat,%):'m2v_dentro'=evalf(%,2);;mv_sale:=m2v-m2v_dentro;mv_sale_:=m2v_-m2v_dentro_:'mv_sale'=evalf(%,2);

`:=`(m2v_dentro, `/`(`*`(p2, `*`(`+`(V, `-`(V2l)))), `*`(R, `*`(T2))))

m2v_dentro = `+`(`*`(69., `*`(kg_)))

`:=`(mv_sale, `+`(m2v, `-`(`/`(`*`(p2, `*`(`+`(V, `-`(V2l)))), `*`(R, `*`(T2))))))

mv_sale = `+`(`*`(0.40e3, `*`(kg_)))

que es aproximadamente igual al del otro modelo.

>