+- HP Forums (https://www.hpmuseum.org/forum)
+-- Forum: HP Software Libraries (/forum-10.html)
+--- Forum: HP-41C Software Library (/forum-11.html)
+--- Thread: Thermodynamic Properties of Refrigerants (/thread-7525.html)
Thermodynamic Properties of Refrigerants - Ángel Martin - 01-07-2017 08:25 AM
Thermodynamic Properties of Refrigerants.
From the Author's Engineering Collection, included in the ETSII3 Module
Building on the "Principle of Corresponding States" program, this addition calculates the Thermodynamic properties Enthalpy and Entropy in the valid ranges for both saturated liquid & vapor dome and superheated vapor. The validity range is intended for pressures lower than the critical pressure PC.
These programs provide a replacement for the different refrigerant chart sheets – using a semi-empirical approach that combines the Martin Equation of State, the Antoine’s Equation, and polynomial expressions (in inverse T) for the isobaric specific heat capacity of the liquid and the gas - valid within the application range, defined as P < Pc.
The Martin Equation of State uses reduced magnitudes, and it’s cubic in the Volume:
Pr = Tr / [Zc.Vr – B ] – A / { Tr^n [ Zc.Vr – B + 1/8 ]^2 }
The data entry process involves all the parameters required, as obtained by the semi-empirical method used. The table below lists the parameters for the main refrigerant gases. Note that the parameters A, B, and N are dimensionless (“n” is not the number of moles). However all coefficients for the specific heats have dimensions - as required by the inversed polynomial expression: Cp = SUM ak/T^k ; with Cp [Energy/Mass*Temperature]
Equations and memory register map:-
Let Pr = P/Pc, Tr = T/Tc, and Vr = V/Vc the reduced pressure, temperature and volume over the corresponding critical values. Let To = 273 K, and Po = Pv(273).
The pressure is directly calculated by the EOS shown above, as:
Pr = Tr / [Zc.Vr – B ] – A / { Tr^n [ Zc.Vr – B + 1/8 ]^2 }
The specific volume is obtained as the maximum real root of the third degree equation given by:
a3.V^3 + a2.V^2 + a1.V + a0 = 0; where the coefficients are defined as:
a0 = A.B + (1/8 - B)^2 [ (Tr)^n+1 + B.Pr.(Tr)^n ]
a1 = (Zc/Vc). { A + Pr.(Tr)^n [3B^2 - B/2 + 1/64] - (1/8 - B)*2(Tr)^n+1 }
a2 = (Zc/Vc)^2 . [ (Tr)^n Pr.(1/4 - 3B) - (Tr)^n+1 ]
a3 = (Zc/Vc)^3 Pr (Tr)^n
The density of the saturated liquid is calculated by the formula:
Vs,l = Vsc. Vg (1-w.G) ; with the following auxiliary definitions:
Vsc = (R.Tc/PM.Pc).(0.292 – 0.0967 w)
G = 0.29607 – 0.09045 Tr – 0.048442 Tr^2
And Vg depends on the actual value of Tr, as follows:
• if: Tr < 0.2 or Tr > 1 => not in liquid phase
• if: 0.2 <= Tr < 0.8; then: Vg = 0.33593(1-Tr) + 1.51941 Tr^2 – 2.02512 Tr4^3 + 1.11422 Tr^4
• if: 0.8 <= Tr <= 1 ; then: Vg = 1 + 1.3 sqrt(1-Tr).log(1-Tr) – 0.50879(1-Tr) – 0.91534(1-Tr)^2
Enthalpy of Vapor: Defined as: H(T,P) = H*(T,P) + DisH(T,P)
- where H* is the value for an ideal gas and "Dis" indicates discrepancy from the perfect gas.
From the base expression for Enthalpy:
H(T,P) = U(T,P) – RT (Z-1) = U(T,P) – (P.V – R.T)
This is calculated using the Helmholtz Free Energy as auxiliary magnitude instead of the internal energy,
defined as: F = U – T.S ; thus: => U = F + TS
Hence the base equation results:
H(T,P)= [F(T,P) + T.S(T,P)] - Q(T,P); with: Q = P.V – R.T/PM
Moving along an isobaric line, we can relate these factors to the values in the saturated liquid - by adding the enthalpy of vaporization and the integral of the specific heat capacity (at constant pressure) between the boiling temperature and the final temperature. Doing that we get:
H(T,P) = [H*(Tb,P)-DisH(To,Po)] + Hvap + ITG {[Cpv]dT} + DisH(T,P)
and replacing the enthalpy by the free energy equivalents:
Hv(T,P) = {H*(Tb,P) + Hvap(Tb) + ITG[Cpv]dT } + [ Dis(F-Fo) - Dis(Q-Qo) + T.Dis(S-So) ]
The integral of the specific heat is easily obtained by direct integration:
ITG{Cpv} = ao(T-To) + a1.Ln(T/To) - a2(1/T - 1/To) - a3(1/T2 - 1/To2) /2 - a4(1/T3 - 1/To3) /3
The vaporization enthalpy is obtained using the Pitzer correlation shown below (per kg), where w is the acentric factor of the substance:
Hvap(To) = (R.Tc/PM) {10.95w.[(1-Tb/Tc)^0.456] + 7.08*(1-Tb/Tc)^0.354 }
and finally the discrepancies of free energy and entropy are calculated as follows, where Zo = P.V/R.T ; and: Zco = Pc.Vc/R.Tc ; in [mol]
DisF(T,P) = R.T.Ln(Z) + ITG{(P – R.T/V).dV} = (R.T/PM).Ln(PM.Zo) + Pc.Vc {A / [Zc.Trn.(Zc.Vr+1/8-B)] – Tr.[Ln(Vr)/PM.Zco - Ln(ZcVr –B)/Zc] }
DisS(T,P) = -∂[DisF]/∂T = (R/PM){ (n+1) + n.A.(PM.Zco/Zc) / [Zc.Vr + 1/8 - B].Tr^(n+1) – Ln(PM.Zo) + Ln[Vr / (Vr - B/Zc)^(PM.Zco/Zc)] }
Entropy of Vapor: Defined as: S(T,P) = S*(T,P) + DisS(T,P) ;
We'll follow the same approach as before, moving along the isobaric line at T=Tb (boiling temperature)l starting at the saturated liquid state.
Sv(T,P) = [S*(Tb,P)- DisS(Tb,P)] – R.Ln(P/Po) + ITG[(Cpv/T).dT] + DisS(T,P)
The integral of the specific heat over the temperature is easily obtained by direct integration of the inverted-polynomial expression:
ITG{Cpv/T} = a0.Ln(T/To) - a1(1/T - 1/To) - a2(1/T2 - 1/To2) / 2 - a3(1/T3 - 1/To3) /3 – a4(1/T4 - 1/To4) / 4
Using the expressions listed above, the final formula for the entropy is as follows:
Sv(T,P) = {S*(Tb,P) + Hvo(To)/To+ ITG[CPv/T]dT}–(R/PM).Ln(P/Po) + [DisS(T,P) – DisS(To,Po)]
Examples.
With the data from the table below you can test the program for a few commonly used refrigerants. Several refrigerant names ("R11", "R12", "R13", "R22", "NH3") can be entered at the initial prompt for an automated entry of all parameters - or alternatively each one of them is requested individually.
Code:
R-11 R-12 R-13 R-22 R-113 R-717
Parameter CFCl3 CF2Cl2 CF3Cl CHClF2 CCl2FCClF2 NH3
=====================================================================================
Zc 0.2766 0.2790 0.27739 0.267 0.2560 0.242004
A 0.421875 0.421875 0.421875 0.421875 0.421875 0.421875
B 0.05598 0.0578 0.0566 0.0488 0,0323 0.03
n -1.1 0.900 0.60 1.00 0.75 -1.1
Tc (C) 197.99 112.00 28.7708 96.00 214.09 132.2808
Pc (kp/cm2) 44.6003 41.9613 39.460 50.300 34.80 115.0036
Vc (cm3/g) 1.8018 1.79168 1.7212 1.9041 1.7301 4.247
PM (g/mol) 137.38 120.90 104.47 86.50 187.39 17.032
w acentric 0.188 0.1760 0.180 0.215 0.252 0.253
-------------------------------------------------------------------------------------
Ln(Pv) = x1 + x2/T ; with T in K and Pv in kp/cm2
x1 10.84436 10.1760 9.96485 10.6316 11.13964 11.71516
x2 3209.743 2469.5656 1901.3441 2460.1029 3568.1032 2803.591
-------------------------------------------------------------------------------------
Cp,v = SUM {ak / T^k} ; k= 0,1,..4 ; with T in K and Cp in kcal/kg.K
a0,v 0.1886 0.25659 0.30905 0.3485 0.25732 0.95239
a1,v -8.0924 -49.1350 -70.3763 -102.14 -15.946 -2.92521E2
a2,v -3.8118E3 5405.60 1.0052E4 1.6546E4 -1.6728E4 6.2791E4
a3,v 4.7327E5 -239860.0 -5.6219E5 -1.003E6 4.9196E6 -4.725E6
a4,v 0.0000 0.0000 0.0000 0.0000 -4.1739E8 0.0000
--------------------------------------------------------------------------------------
Cp,l = SUM {ak / T^k} ; k= 0,1,..4
a0,l 0.56642 9.9361 11.190 14.4110 4.5917 0.24694
a1,l -2.3715E2 -9.445E3 -8.4336E3 -1.3625E4 -5.058E3 -2.3079E4
a2,l 5.2605E4 3.4338E6 2.4332E6 4.9211E6 2.2046E6 8.4508E6
a3,l -4.0107E6 -5.5272E8 -3.1189E8 -7.8839E8 -4.2602E8 -1.3705E9
a4,l 0.0000 3.3182E10 1.4964E10 4.7189E10 3.0626E10 8.2854E10The programs work internally with the units reflected in the table but you can input and output the values in any other units. It comes without saying that these programs make extensive use of the Unit Management System (UMS), therefore the unit Conversion Module needs to be plugged into the calculator.
Example1.- Superheated Vapor.
For Freon-22, T= 25 deg C and P =0.2 atm
Vesp= 1.398 m3/kg
Hv = 324.646 kJ/kg
Sv = 1.440 kJ/kg.K