Update refrigerants images for mEoS documentation including PIP plot #1

.script/plotMEoS.py

@@ -117,15 +117,15 @@
     extended = np.logspace(np.log10(Tc), np.log10(Tmax), 10)
 else:
     extended = np.linspace(Tc, Tmax, 5)
-#     extended = np.r_[Tc, Tmax]
+    # extended = np.r_[Tc, Tmax]
 isoT = np.concatenate([
     np.linspace(Tmin, Tc, 5, endpoint=False), extended])
 isoT_kw = {"ls": ":", "color": "red", "lw": 0.8}
 
 # Isobar lines to plot
 isoP = np.logspace(np.log10(Pt.P), np.log10(Pmax), 10)
-isoP_PIP = np.r_[0.1, 0.2, 0.5, 1, 2, 5, fluid.Pc.MPa, 15, 20, 25, 30, 35, 50,
-             100, 200, 500, 1000]*1e6
+isoP_PIP = np.r_[0.1, 0.2, 0.5, 1, 2, 5, fluid.Pc.MPa, 15, 20, 25, 30, 35,
+                 50, 100, 200, 500, 1000]*1e6
 isoP_kw = {"ls": ":", "color": "blue", "lw": 0.8}
 
 # # Isoenthalpic lines to plot
@@ -315,7 +315,7 @@
 smax = max(svap)
 smax += abs(smax * 0.1)
 
-# ax_Ideal.set_ylim(bottom=0.01, top=)
+# ax_Ideal.set_ylim(bottom=1e-6)
 ax_PIP.set_ylim(-10, 15)
 ax_PIP.set_xlim(Tmin, min(2*fluid.Tc, Tmax))
 ax_Ph.set_xlim(Pt.h.kJkg, hmax)
@@ -445,9 +445,15 @@
     if P < fluid.Pc:
         sat_pnt = []
         for x in np.linspace(1, 0.1, 10):
-            sat_pnt.append(fluid(P=P, x=x))
+            if sat_pnt and sat_pnt[-1].status == 1:
+                sat_pnt.append(fluid(T=sat_pnt[-1].T, x=x))
+            else:
+                sat_pnt.append(fluid(P=P, x=x))
         for x in np.linspace(0.1, 0, 11):
-            sat_pnt.append(fluid(P=P, x=x))
+            if sat_pnt and sat_pnt[-1].status == 1:
+                sat_pnt.append(fluid(T=sat_pnt[-1].T, x=x))
+            else:
+                sat_pnt.append(fluid(P=P, x=x))
         sat = True
     else:
         sat = False
@@ -457,10 +463,9 @@
         point = fluid(P=P, T=t)
 
         if point.status == 1:
-#             print(t, point.T, point.v, point.u)
 
-#             if point.rhoM > point._constants["rhomax"]:
-#                 continue
+            # if point.rhoM > point._constants["rhomax"]:
+            #     continue
 
             # Discard point below the melting line
             if fluid._melting:
@@ -479,17 +484,18 @@ def f_p(t):
                 for p in sat_pnt[::-1]:
                     if p.status == 1:
                         pts.append(p)
-#                         print("SAT", t, p.T, p.v, p.u)
+                        # print("SAT", t, p.T, p.v, p.u, p.x)
                 sat = False
+            # print(t, point.T, point.v, point.u, point.x, point.status)
             pts.append(point)
 
     # Clean triple point when create a loop in isotherm
-#     if P == isoP[0]:
-#         for p in pts:
-#             print(p.u, p.v)
-#         for p in sat_pnt:
-#             print("sat", p.uv)
-#         del pts[0]
+    # if P == isoP[0]:
+    #     for p in pts:
+    #         print(p.u, p.v)
+    #     for p in sat_pnt:
+    #         print("sat", p.uv)
+    #     del pts[0]
 
     h = []
     T = []
@@ -501,8 +507,9 @@ def f_p(t):
         if p.status == 1:
             if p.x in (0, 1):
                 if not PIP and (p.PIP < 0 or p.PIP > 10):
-                    continue
-                PIP.append(p.PIP)
+                    PIP.append(None)
+                else:
+                    PIP.append(p.PIP)
             else:
                 PIP.append(None)
             h.append(p.h.kJkg)
@@ -760,10 +767,10 @@ def f_v(ti):
     for t in T:
         # Optional limits for ideal curves to avoid false values out of range
 #         if curva in ["ideal", "boyle"]:
-#             if t > 5*fluid.Tc:
+#             if t > 3*fluid.Tc:
 #                 continue
 #         if curva == "joule-thomson":
-#             if t > 8*fluid.Tc:
+#             if t > 4.5*fluid.Tc:
 #                 continue
 
         P = fluid()._IdealCurve(curva, t)

lib/mEoS/C5F12.py

@@ -61,7 +61,7 @@ class C5F12(MEoS):
         "cp": CP1,
         "ref": "NBP",
 
-        "Tmin": Tt, "Tmax": 500, "Pmax": 10000,
+        "Tmin": Tt, "Tmax": 500, "Pmax": 100000,
 
         "nr1": [0.036359237, 1.1222855, -0.63721964, -1.3599906, 0.24995429],
         "d1": [4, 1, 1, 2, 3],

lib/mEoS/D2O.py

@@ -223,8 +223,8 @@ def _visco0(self, rho, T, fase):
                "method": "_thermo0",
                "__name__": "IAPWS (2021)",
                "__doi__": {
-                   "autor": "Kestin, J., Sengers, J.V., Kamgar-Parsi, B., "
-                            "Levelt Sengers, J.M.H.",
+                   "autor": "Huber, M.L., Perkins, R.A., Assael, M.J., "
+                            "Monogenidou, S.A., Hellmann, R., Sengers, J.V.",
                    "title": "New International Formulation for the Thermal "
                             "Conductivity of Heavy Water",
                    "ref": "J. Phys. Chem. Ref. Data 51(1) (2022) 013102",

lib/mEoS/R12.py

@@ -135,7 +135,7 @@ class R12(MEoS):
         "t": [0.57, 0.72, 0.89, 1.07, 1.25]}
     _vapor_Density = {
         "eq": 2,
-        "n": [-3.153, -6.4734, -17.346, -15.918e2, -32.492, -120.72],
+        "n": [-3.153, -6.4734, -17.346, -15.918, -32.492, -120.72],
         "t": [0.418, 1.32, 3.3, 6.6, 7.0, 15.0]}
 
     thermo0 = {"__name__": "Krauss (1989)",

lib/mEoS/R13.py

@@ -153,3 +153,10 @@ class R13(MEoS):
 
     _viscosity = (trnECS, )
     _thermal = (trnECS, )
+
+
+if __name__ == "__main__":
+    st = R13(P=256, T=92)
+    print(st.status, st.msg)
+    print(st.P)
+    print(st.v, st.u)

lib/mEoS/R13I1.py

@@ -22,7 +22,7 @@
 
 from lib import unidades
 from lib.meos import MEoS
-from lib.mEoS import R134a
+from lib.mEoS.R134a import R134a
 
 
 class R13I1(MEoS):