> restart:#"m07_p42"

Se quiere estudiar la influencia del gas de trabajo en el empuje de un pequeño cohete de control de actitud. Para el gas de trabajo se van a considerar mezclas de hidrógeno y argón. La tobera (convergente-divergente) se va a considerar isoentrópica, con condiciones de entrada fijas, p1=1 MPa y T1=300 K y baja velocidad (correspondientes al estado inicial del depósito de propulsante), y salida a través de un área de salida de A2=1 mm2, adaptada a una presión ambiente de 100 kPa. Se pide, para una fracción molar de hidrógeno genérica, y para los valores límite de este parámetro:

a) Masa molar y capacidad térmica específica a presión constante.

b) Relación de capacidades térmicas y temperatura de salida.

c) Densidades de entrada y salida.

d) Velocidad de salida y número de Mach correspondiente.

e) Gasto másico y empuje.

Datos:

> read`../therm_eq.m`:read`../therm_const.m`:read`../therm_proc.m`:with(therm_proc):Digits:=5:

> su1:="H2":su2:="Ar":dat:=[p1=1e6*Pa_,T1=300*K_,p2=100e3*Pa_,A2=1e-6*m_^2];

[p1 = `+`(`*`(0.1e7, `*`(Pa_))), T1 = `+`(`*`(300, `*`(K_))), p2 = `+`(`*`(0.100e6, `*`(Pa_))), A2 = `+`(`*`(0.1e-5, `*`(`^`(m_, 2))))]

Image

> g1dat:=get_gas_data(su1):g2dat:=get_gas_data(su2):dat:=op(dat),Const,SI2,SI1:Molar_mass_H2_Ar:=[subs(g1dat,M),subs(g2dat,M)];

[`+`(`/`(`*`(0.2e-2, `*`(kg_)), `*`(mol_))), `+`(`/`(`*`(0.40e-1, `*`(kg_)), `*`(mol_)))]

a) Masa molar y capacidad térmica específica a presión constante.

Sea x1 la fracción molar de H2.

> Mm:=x1*M1+x2*M2;cpm:=x1*cp1+x2*cp2;gammam:='cpm/(cpm-R[u])';Mm:=subs(x2=1-x1,x1*M1+x2*M2);Mm_:=subs(M1=M,g1dat,M2=M,g2dat,Mm);plot(subs(SI0,Mm_),x1=0..1,'Mm'=0..0.05);cpm_:=subs(x2=1-x1,cp1=c[p]*M,g1dat,cp2=c[p]*M,g2dat,cpm);plot(subs(SI0,cpm_),x1=0..1,'c[p,molar]'=0..30);plot(subs(SI0,cpm_/Mm_),x1=0..1,'c[p,mass]'=0..15000);cp_mass_H2_Ar:=[subs(x1=1,cpm_/Mm_),subs(x1=0,cpm_/Mm_)];

`+`(`*`(x1, `*`(M1)), `*`(x2, `*`(M2)))
`+`(`*`(x1, `*`(cp1)), `*`(x2, `*`(cp2)))
`/`(`*`(cpm), `*`(`+`(cpm, `-`(R[u]))))
`+`(`*`(x1, `*`(M1)), `*`(`+`(1, `-`(x1)), `*`(M2)))
`+`(`/`(`*`(0.2e-2, `*`(x1, `*`(kg_))), `*`(mol_)), `/`(`*`(0.40e-1, `*`(`+`(1, `-`(x1)), `*`(kg_))), `*`(mol_)))
Plot_2d
`+`(`/`(`*`(28.400, `*`(x1, `*`(J_))), `*`(K_, `*`(mol_))), `/`(`*`(20.920, `*`(`+`(1, `-`(x1)), `*`(J_))), `*`(K_, `*`(mol_))))
Plot_2d
Plot_2d
[`+`(`/`(`*`(14200., `*`(J_)), `*`(kg_, `*`(K_)))), `+`(`/`(`*`(523.00, `*`(J_)), `*`(kg_, `*`(K_))))]

La masa molar de la mezcla varía linealmente con x1 desde 40 g/mol (todo Ar) hasta 2 g/mol (todo H2).

La capacidad térmica específica a presión constante no varía linealmente con x1 (la capacidad térmica molar sí lo haría). Los límites son cp=14000 J/(kg·K) para todo H2 y 523 J/(kg·K) para todo Ar.

b) Relación de capacidades térmicas y temperatura de salida.

Nótese que como la entrada es a baja velocidad, las condiciones totales serán prácticamente las estáticas: T1t?T1 y p1t=p1, pero no así a la salida.

> gammam_:=subs(dat,x2=1-x1,cp1=c[p]*M,g1dat,cp2=c[p]*M,g2dat,gammam);plot(subs(SI0,gammam_),x1=0..1,'gamma'=1..1.7);gamma_H2_Ar:=[subs(x1=1,gammam_),subs(x1=0,gammam_)];eqT2:=T2=T1*(p2/p1)^((gamma-1)/gamma);T2_:=subs(dat,simplify(subs(gamma=gammam_,dat,rhs(eqT2))));plot(subs(SI0,T2_),x1=0..1,T2=0..300);T2_H2_Ar:=evalf([subs(x1=1,T2_),subs(x1=0,T2_)]);

`/`(`*`(`+`(`/`(`*`(28.400, `*`(x1, `*`(J_))), `*`(K_, `*`(mol_))), `/`(`*`(20.920, `*`(`+`(1, `-`(x1)), `*`(J_))), `*`(K_, `*`(mol_))))), `*`(`+`(`/`(`*`(28.400, `*`(x1, `*`(J_))), `*`(K_, `*`(mol_))...
Plot_2d
[1.4139, 1.6595]
T2 = `*`(T1, `*`(`^`(`/`(`*`(p2), `*`(p1)), `/`(`*`(`+`(gamma, `-`(1))), `*`(gamma)))))
`+`(`*`(300, `*`(K_, `*`(exp(`+`(`-`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.))))))))))
Plot_2d
[`+`(`*`(152.89, `*`(K_))), `+`(`*`(120.14, `*`(K_)))]

La relación de capacidades térmicas no es lineal en x1 (es un cociente de funciones lineales, i.e. varía hiperbólicamente con x1), yendo de 1,40 para H2 puro hasta 1,67 para Ar puro.

La temperatura de salida es muy baja debido al gran enfriamiento que conlleva la expansión isentrópica.

c) Densidades de entrada y salida.

> eqrho:=rho='p*M/(R[u]*T)';eqrho1:=expand(subs(p=p1,T=T1,M=Mm_,dat,%));plot(subs(SI0,rhs(eqrho1)),x1=0..1,rho[1]=0..20);rho1_H2_Ar:=[subs(x1=1,rhs(eqrho1)),subs(x1=0,rhs(eqrho1))];eqrho2:=expand(subs(p=p2,T=T2_,M=Mm_,dat,eqrho));plot(subs(T2=T2_,SI0,rhs(eqrho2)),x1=0..1,rho[2]=0..4);rho2_H2_Ar:=evalf([subs(x1=1,rhs(eqrho2)),subs(x1=0,rhs(eqrho2))]);

rho = `/`(`*`(p, `*`(M)), `*`(R[u], `*`(T)))
rho = `+`(`-`(`/`(`*`(15.235, `*`(x1, `*`(kg_))), `*`(`^`(m_, 3)))), `/`(`*`(16.037, `*`(kg_)), `*`(`^`(m_, 3))))
Plot_2d
[`+`(`/`(`*`(.802, `*`(kg_)), `*`(`^`(m_, 3)))), `+`(`/`(`*`(16.037, `*`(kg_)), `*`(`^`(m_, 3))))]
rho = `+`(`-`(`/`(`*`(1.5235, `*`(x1, `*`(kg_))), `*`(`^`(m_, 3), `*`(exp(`+`(`-`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.)))))))))), `/`(`*`(1.6037, `*`(kg_)), `*`(`^`(m_, 3), `*`(exp(`+`(`-...
Plot_2d
[`+`(`/`(`*`(.15737, `*`(kg_)), `*`(`^`(m_, 3)))), `+`(`/`(`*`(4.0044, `*`(kg_)), `*`(`^`(m_, 3))))]

La densidad de entrada es lineal en x1 porque lo es en la masa molar, yendo de 0,8 kg/m3 con H2 puro a 16 kg/m3 con Ar puro.

La densidad de salida varía casi linealmente desde 0,16 kg/m3 con H2 puro a 4 kg/m3 con Ar puro.

d) Velocidad de salida y número de Mach correspondiente.

> eqBE:=h1=h2;eqBE:=cp*T1=cp*T2+v2^2/2;v2_:=sqrt('2*(cpm_/Mm_)*(T1-T2_)');v2__:=simplify(subs(T2=T2_,dat,%));plot(subs(SI0,v2__),x1=0..1,v=0..2000);v2_H2_Ar:=evalf([subs(x1=1,dat,v2_),subs(x1=0,dat,v2_)]);

h1 = h2
`*`(cp, `*`(T1)) = `+`(`*`(cp, `*`(T2)), `*`(`/`(1, 2), `*`(`^`(v2, 2))))
`*`(`^`(2, `/`(1, 2)), `*`(`^`(`/`(`*`(cpm_, `*`(`+`(T1, `-`(T2_)))), `*`(Mm_)), `/`(1, 2))))
`+`(`*`(109.54, `*`(`^`(`/`(`*`(`+`(`-`(1.), exp(`+`(`-`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.))))))), `*`(`^`(m_, 2), `*`(`+`(`*`(187., `*`(x1)), 523.)))), `*`(`^`(s_, 2), `*`(`+`(`*`(19....
Plot_2d
[`+`(`*`(2044.0, `*`(`^`(`/`(`*`(`^`(m_, 2)), `*`(`^`(s_, 2))), `/`(1, 2))))), `+`(`*`(433.73, `*`(`^`(`/`(`*`(`^`(m_, 2)), `*`(`^`(s_, 2))), `/`(1, 2)))))]

> eqM:=Mach=v/c;eqc:=c=sqrt(gamma*R*T);eqc:=c=sqrt(gamma*R[u]*T/M);eqc_:=simplify(subs(gamma=gammam_,T=T2_,M=Mm_,dat,eqc));plot(subs(SI0,rhs(eqc_)),x1=0..1,c=0..1000);c_H2_Ar:=evalf([subs(x1=1,dat,rhs(eqc_)),subs(x1=0,dat,rhs(eqc_))]);eqM_:=simplify(subs(eqc_,v=v2__,dat,eqM));plot(subs(SI0,rhs(eqM_)),x1=0..1,Mach=0..3);Mach_H2_Ar:=evalf([subs(x1=1,dat,rhs(eqM_)),subs(x1=0,dat,rhs(eqM_))]);

Mach = `/`(`*`(v), `*`(c))
c = `*`(`^`(`*`(gamma, `*`(R, `*`(T))), `/`(1, 2)))
c = `*`(`^`(`/`(`*`(gamma, `*`(R[u], `*`(T))), `*`(M)), `/`(1, 2)))
c = `+`(`*`(4994.3, `*`(`^`(`/`(`*`(kg_, `*`(`^`(m_, 4), `*`(`+`(`*`(187., `*`(x1)), 523.), `*`(exp(`+`(`-`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.)))))))))), `*`(`^`(s_, 2), `*`(`+`(`-`(`*`...
c = `+`(`*`(4994.3, `*`(`^`(`/`(`*`(kg_, `*`(`^`(m_, 4), `*`(`+`(`*`(187., `*`(x1)), 523.), `*`(exp(`+`(`-`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.)))))))))), `*`(`^`(s_, 2), `*`(`+`(`-`(`*`...
c = `+`(`*`(4994.3, `*`(`^`(`/`(`*`(kg_, `*`(`^`(m_, 4), `*`(`+`(`*`(187., `*`(x1)), 523.), `*`(exp(`+`(`-`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.)))))))))), `*`(`^`(s_, 2), `*`(`+`(`-`(`*`...
Plot_2d
[`+`(`*`(947.97, `*`(`^`(`/`(`*`(`^`(m_, 2)), `*`(`^`(s_, 2))), `/`(1, 2))))), `+`(`*`(203.57, `*`(`^`(`/`(`*`(`^`(m_, 2)), `*`(`^`(s_, 2))), `/`(1, 2)))))]
Mach = `+`(`/`(`*`(0.72743e-1, `*`(`^`(`/`(`*`(`+`(`-`(1.), exp(`+`(`-`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.))))))), `*`(`^`(m_, 2), `*`(`+`(`*`(187., `*`(x1)), 523.)))), `*`(`^`(s_, 2), ...
Plot_2d
[2.1561, 2.1304]

La velocidad de salida varía mucho con x1, desde 434 m/s con Ar puro a 2044 m/s con H2 puro, pero la velocidad local del sonido en el gas varía de modo muy parecido, por lo que el número de Mach apenas varía con x1, desde 2,13 en Ar a 2,16 en H2.

e) Gasto másico y empuje.

> eqm:=m=rho*v*A;eqm:=m=rho2*v2*A2;eqm_:=subs(rho2=rho,eqrho2,v2=v2__,dat,eqm);plot(subs(SI0,rhs(eqm_)),x1=0..1,m=0..0.002);m_H2_Ar:=evalf([subs(x1=1,rhs(eqm_)),subs(x1=0,rhs(eqm_))]);

m = `*`(rho, `*`(v, `*`(A)))
m = `*`(rho2, `*`(v2, `*`(A2)))
m = `+`(`*`(0.10954e-3, `*`(`+`(`-`(`/`(`*`(1.5235, `*`(x1, `*`(kg_))), `*`(`^`(m_, 3), `*`(exp(`+`(`-`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.)))))))))), `/`(`*`(1.6037, `*`(kg_)), `*`(`^`(...
m = `+`(`*`(0.10954e-3, `*`(`+`(`-`(`/`(`*`(1.5235, `*`(x1, `*`(kg_))), `*`(`^`(m_, 3), `*`(exp(`+`(`-`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.)))))))))), `/`(`*`(1.6037, `*`(kg_)), `*`(`^`(...
Plot_2d
[`+`(`/`(`*`(0.32165e-3, `*`(kg_, `*`(`^`(`/`(`*`(`^`(m_, 2)), `*`(`^`(s_, 2))), `/`(1, 2))))), `*`(m_))), `+`(`/`(`*`(0.17368e-2, `*`(kg_, `*`(`^`(`/`(`*`(`^`(m_, 2)), `*`(`^`(s_, 2))), `/`(1, 2)))))...

el gasto másico sería de 0,3 g/s para H2 puro y de 1,7 g/s para Ar puro.

Por estar la tobera adaptada, el empuje es el gasto másico por la volocidad de salida (en general el empuje de un cohete sería F=m*v2+A2*(p2-p0), siendo p0 la presión ambiente.

> F:=m*v2;F:=rho2*v2^2*A2;F_:=simplify(subs(rho2=rho,eqrho2,v2=v2__,dat,F));plot(subs(SI0,F_),x1=0..1,Thrust=0..1);F_H2_Ar:=subs(dat,evalf([subs(x1=1,dat,F_),subs(x1=0,dat,F_)]));

`*`(m, `*`(v2))
`*`(rho2, `*`(`^`(v2, 2), `*`(A2)))
`+`(`-`(`/`(`*`(0.11999e-5, `*`(exp(`+`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.))))), `*`(kg_, `*`(`+`(`*`(15235., `*`(x1)), `-`(16037.)), `*`(`+`(`-`(1.), exp(`+`(`-`(`/`(`*`(478.60), `*`(`...
`+`(`-`(`/`(`*`(0.11999e-5, `*`(exp(`+`(`/`(`*`(478.60), `*`(`+`(`*`(187., `*`(x1)), 523.))))), `*`(kg_, `*`(`+`(`*`(15235., `*`(x1)), `-`(16037.)), `*`(`+`(`-`(1.), exp(`+`(`-`(`/`(`*`(478.60), `*`(`...
Plot_2d
[`+`(`*`(.65741, `*`(N_))), `+`(`*`(.75329, `*`(N_)))]

i.e. el empuje apenas depende de x1: el H2 puro da 0,66 N y el Ar puro 0,75.

Se concluye que los motores cohete de gas frío (i.e. sin reacciones químicas) conviene usar gases ligeros, puesto que dan casi el mismo empuje con un gasto másico mucho menor (si bien las toberas para estos últimos son más delicadas por la menor garganta requerida para conseguir las altas velocidades implicadas).

>