CoolProp
FluidLibrary.h
1 
2 #ifndef FLUIDLIBRARY_H
3 #define FLUIDLIBRARY_H
4 
5 #include "CoolPropFluid.h"
6 
7 #include "rapidjson_include.h"
8 #include "crossplatform_shared_ptr.h"
9 
10 #include <map>
11 #include <algorithm>
12 #include "Configuration.h"
13 #include "Backends/Cubics/CubicsLibrary.h"
14 #include "Helmholtz.h"
15 
16 namespace CoolProp{
17 
18 // Forward declaration of the necessary debug function to avoid including the whole header
19 extern int get_debug_level();
20 
22 
28 {
30  std::map<std::size_t, CoolPropFluid> fluid_map;
32  std::map<std::size_t, std::string> JSONstring_map;
33  std::vector<std::string> name_vector;
34  std::map<std::string, std::size_t> string_to_index_map;
35  bool _is_empty;
36 public:
37 
39  static ResidualHelmholtzContainer parse_alphar(rapidjson::Value &jsonalphar)
40  {
42 
43  for (rapidjson::Value::ValueIterator itr = jsonalphar.Begin(); itr != jsonalphar.End(); ++itr)
44  {
45  // A reference for code cleanness
46  rapidjson::Value &contribution = *itr;
47 
48  // Get the type (required!)
49  std::string type = contribution["type"].GetString();
50 
51  if (!type.compare("ResidualHelmholtzPower"))
52  {
53  std::vector<CoolPropDbl> n = cpjson::get_long_double_array(contribution["n"]);
54  std::vector<CoolPropDbl> d = cpjson::get_long_double_array(contribution["d"]);
55  std::vector<CoolPropDbl> t = cpjson::get_long_double_array(contribution["t"]);
56  std::vector<CoolPropDbl> l = cpjson::get_long_double_array(contribution["l"]);
57  assert(n.size() == d.size());
58  assert(n.size() == t.size());
59  assert(n.size() == l.size());
60 
61  alphar.GenExp.add_Power(n,d,t,l);
62  }
63  else if (!type.compare("ResidualHelmholtzGaussian"))
64  {
65  std::vector<CoolPropDbl> n = cpjson::get_long_double_array(contribution["n"]);
66  std::vector<CoolPropDbl> d = cpjson::get_long_double_array(contribution["d"]);
67  std::vector<CoolPropDbl> t = cpjson::get_long_double_array(contribution["t"]);
68  std::vector<CoolPropDbl> eta = cpjson::get_long_double_array(contribution["eta"]);
69  std::vector<CoolPropDbl> epsilon = cpjson::get_long_double_array(contribution["epsilon"]);
70  std::vector<CoolPropDbl> beta = cpjson::get_long_double_array(contribution["beta"]);
71  std::vector<CoolPropDbl> gamma = cpjson::get_long_double_array(contribution["gamma"]);
72  assert(n.size() == d.size());
73  assert(n.size() == t.size());
74  assert(n.size() == eta.size());
75  assert(n.size() == epsilon.size());
76  assert(n.size() == beta.size());
77  assert(n.size() == gamma.size());
78  alphar.GenExp.add_Gaussian(n,d,t,eta,epsilon,beta,gamma);
79  }
80  else if (!type.compare("ResidualHelmholtzGaoB"))
81  {
82  std::vector<CoolPropDbl> n = cpjson::get_long_double_array(contribution["n"]);
83  std::vector<CoolPropDbl> t = cpjson::get_long_double_array(contribution["t"]);
84  std::vector<CoolPropDbl> d = cpjson::get_long_double_array(contribution["d"]);
85  std::vector<CoolPropDbl> eta = cpjson::get_long_double_array(contribution["eta"]);
86  std::vector<CoolPropDbl> beta = cpjson::get_long_double_array(contribution["beta"]);
87  std::vector<CoolPropDbl> gamma = cpjson::get_long_double_array(contribution["gamma"]);
88  std::vector<CoolPropDbl> epsilon = cpjson::get_long_double_array(contribution["epsilon"]);
89  std::vector<CoolPropDbl> b = cpjson::get_long_double_array(contribution["b"]);
90  assert(n.size() == t.size());
91  assert(n.size() == d.size());
92  assert(n.size() == eta.size());
93  assert(n.size() == epsilon.size());
94  assert(n.size() == beta.size());
95  assert(n.size() == gamma.size());
96  assert(n.size() == b.size());
97  alphar.GaoB = ResidualHelmholtzGaoB(n,t,d,eta,beta,gamma,epsilon,b);
98  }
99  else if (!type.compare("ResidualHelmholtzNonAnalytic"))
100  {
101  if (alphar.NonAnalytic.N > 0){throw ValueError("Cannot add ");}
102  std::vector<CoolPropDbl> n = cpjson::get_long_double_array(contribution["n"]);
103  std::vector<CoolPropDbl> a = cpjson::get_long_double_array(contribution["a"]);
104  std::vector<CoolPropDbl> b = cpjson::get_long_double_array(contribution["b"]);
105  std::vector<CoolPropDbl> beta = cpjson::get_long_double_array(contribution["beta"]);
106  std::vector<CoolPropDbl> A = cpjson::get_long_double_array(contribution["A"]);
107  std::vector<CoolPropDbl> B = cpjson::get_long_double_array(contribution["B"]);
108  std::vector<CoolPropDbl> C = cpjson::get_long_double_array(contribution["C"]);
109  std::vector<CoolPropDbl> D = cpjson::get_long_double_array(contribution["D"]);
110  assert(n.size() == a.size());
111  assert(n.size() == b.size());
112  assert(n.size() == beta.size());
113  assert(n.size() == A.size());
114  assert(n.size() == B.size());
115  assert(n.size() == C.size());
116  assert(n.size() == D.size());
117  alphar.NonAnalytic = ResidualHelmholtzNonAnalytic(n,a,b,beta,A,B,C,D);
118  }
119  else if (!type.compare("ResidualHelmholtzLemmon2005"))
120  {
121  std::vector<CoolPropDbl> n = cpjson::get_long_double_array(contribution["n"]);
122  std::vector<CoolPropDbl> d = cpjson::get_long_double_array(contribution["d"]);
123  std::vector<CoolPropDbl> t = cpjson::get_long_double_array(contribution["t"]);
124  std::vector<CoolPropDbl> l = cpjson::get_long_double_array(contribution["l"]);
125  std::vector<CoolPropDbl> m = cpjson::get_long_double_array(contribution["m"]);
126  assert(n.size() == d.size());
127  assert(n.size() == t.size());
128  assert(n.size() == l.size());
129  assert(n.size() == m.size());
130  alphar.GenExp.add_Lemmon2005(n,d,t,l,m);
131  }
132  else if (!type.compare("ResidualHelmholtzExponential"))
133  {
134  std::vector<CoolPropDbl> n = cpjson::get_long_double_array(contribution["n"]);
135  std::vector<CoolPropDbl> d = cpjson::get_long_double_array(contribution["d"]);
136  std::vector<CoolPropDbl> t = cpjson::get_long_double_array(contribution["t"]);
137  std::vector<CoolPropDbl> g = cpjson::get_long_double_array(contribution["g"]);
138  std::vector<CoolPropDbl> l = cpjson::get_long_double_array(contribution["l"]);
139  assert(n.size() == d.size());
140  assert(n.size() == t.size());
141  assert(n.size() == g.size());
142  assert(n.size() == l.size());
143  alphar.GenExp.add_Exponential(n,d,t,g,l);
144  }
145  else if (!type.compare("ResidualHelmholtzAssociating"))
146  {
147  if (alphar.SAFT.disabled == false){throw ValueError("Cannot add ");}
148  CoolPropDbl a = cpjson::get_double(contribution,"a");
149  CoolPropDbl m = cpjson::get_double(contribution,"m");
150  CoolPropDbl epsilonbar = cpjson::get_double(contribution,"epsilonbar");
151  CoolPropDbl vbarn = cpjson::get_double(contribution,"vbarn");
152  CoolPropDbl kappabar = cpjson::get_double(contribution,"kappabar");
153  alphar.SAFT = ResidualHelmholtzSAFTAssociating(a,m,epsilonbar,vbarn,kappabar);
154  }
155  else
156  {
157  throw ValueError(format("Unsupported Residual helmholtz type: %s",type.c_str()));
158  }
159  }
160 
161  // Finish adding parts to the Generalized Exponential term, build other vectors
162  alphar.GenExp.finish();
163 
164  return alphar;
165  };
166 
168  static IdealHelmholtzContainer parse_alpha0(rapidjson::Value &jsonalpha0)
169  {
170  if (!jsonalpha0.IsArray()){throw ValueError();}
171 
173 
174  for (rapidjson::Value::ConstValueIterator itr = jsonalpha0.Begin(); itr != jsonalpha0.End(); ++itr)
175  {
176  // A reference for code cleanness
177  const rapidjson::Value &contribution = *itr;
178 
179  // Get the type (required!)
180  std::string type = contribution["type"].GetString();
181 
182  if (!type.compare("IdealGasHelmholtzLead"))
183  {
184  if (alpha0.Lead.is_enabled() == true){throw ValueError("Cannot add ");}
185  CoolPropDbl a1 = cpjson::get_double(contribution,"a1");
186  CoolPropDbl a2 = cpjson::get_double(contribution,"a2");
187 
188  alpha0.Lead = IdealHelmholtzLead(a1, a2);
189  }
190  else if (!type.compare("IdealGasHelmholtzPower"))
191  {
192  if (alpha0.Power.is_enabled() == true){throw ValueError("Cannot add ");}
193  std::vector<CoolPropDbl> n = cpjson::get_long_double_array(contribution["n"]);
194  std::vector<CoolPropDbl> t = cpjson::get_long_double_array(contribution["t"]);
195 
196  alpha0.Power = IdealHelmholtzPower(n, t);
197  }
198  else if (!type.compare("IdealGasHelmholtzLogTau"))
199  {
200  if (alpha0.LogTau.is_enabled() == true){throw ValueError("Cannot add ");}
201  CoolPropDbl a = cpjson::get_double(contribution,"a");
202 
203  alpha0.LogTau = IdealHelmholtzLogTau(a);
204  }
205  else if (!type.compare("IdealGasHelmholtzPlanckEinsteinGeneralized"))
206  {
207  // Retrieve the values
208  std::vector<CoolPropDbl> n = cpjson::get_long_double_array(contribution["n"]);
209  std::vector<CoolPropDbl> t = cpjson::get_long_double_array(contribution["t"]);
210 
211  std::vector<CoolPropDbl> c = cpjson::get_long_double_array(contribution["c"]);
212  std::vector<CoolPropDbl> d = cpjson::get_long_double_array(contribution["d"]);
213 
214  if (alpha0.PlanckEinstein.is_enabled() == true){
215  alpha0.PlanckEinstein.extend(n, t, c, d);
216  }
217  else{
218  alpha0.PlanckEinstein = IdealHelmholtzPlanckEinsteinGeneralized(n, t, c, d);
219  }
220  }
221  else if (!type.compare("IdealGasHelmholtzPlanckEinstein"))
222  {
223  // Retrieve the values
224  std::vector<CoolPropDbl> n = cpjson::get_long_double_array(contribution["n"]);
225  std::vector<CoolPropDbl> t = cpjson::get_long_double_array(contribution["t"]);
226  // Flip the sign of theta
227  for (std::size_t i = 0; i < t.size(); ++i){ t[i] *= -1;}
228  std::vector<CoolPropDbl> c(n.size(), 1);
229  std::vector<CoolPropDbl> d(c.size(), -1);
230 
231  if (alpha0.PlanckEinstein.is_enabled() == true){
232  alpha0.PlanckEinstein.extend(n, t, c, d);
233  }
234  else{
235  alpha0.PlanckEinstein = IdealHelmholtzPlanckEinsteinGeneralized(n, t, c, d);
236  }
237  }
238  else if (!type.compare("IdealGasHelmholtzCP0Constant"))
239  {
240  if (alpha0.CP0Constant.is_enabled() == true){throw ValueError("Cannot add another IdealGasHelmholtzCP0Constant term; join them together");}
241  CoolPropDbl cp_over_R = cpjson::get_double(contribution, "cp_over_R");
242  CoolPropDbl Tc = cpjson::get_double(contribution, "Tc");
243  CoolPropDbl T0 = cpjson::get_double(contribution, "T0");
244  alpha0.CP0Constant = IdealHelmholtzCP0Constant(cp_over_R, Tc, T0);
245  }
246  else if (!type.compare("IdealGasHelmholtzCP0PolyT"))
247  {
248  if (alpha0.CP0PolyT.is_enabled() == true){throw ValueError("Cannot add another CP0PolyT term; join them together");}
249  std::vector<CoolPropDbl> c = cpjson::get_long_double_array(contribution["c"]);
250  std::vector<CoolPropDbl> t = cpjson::get_long_double_array(contribution["t"]);
251  CoolPropDbl Tc = cpjson::get_double(contribution, "Tc");
252  CoolPropDbl T0 = cpjson::get_double(contribution, "T0");
253  alpha0.CP0PolyT = IdealHelmholtzCP0PolyT(c, t, Tc, T0);
254  }
255  else if (!type.compare("IdealGasHelmholtzCP0AlyLee"))
256  {
257 
258  std::vector<CoolPropDbl> constants = cpjson::get_long_double_array(contribution["c"]);
259  CoolPropDbl Tc = cpjson::get_double(contribution, "Tc");
260  CoolPropDbl T0 = cpjson::get_double(contribution, "T0");
261 
262  // Take the constant term if nonzero and set it as a polyT term
263  if (std::abs(constants[0]) > 1e-14){
264  std::vector<CoolPropDbl> c(1,constants[0]), t(1,0);
265  if (alpha0.CP0PolyT.is_enabled() == true){
266  alpha0.CP0PolyT.extend(c,t);
267  }
268  else{
269  alpha0.CP0PolyT = IdealHelmholtzCP0PolyT(c, t, Tc, T0);
270  }
271  }
272  std::vector<CoolPropDbl> n, c, d, t;
273  if (std::abs(constants[1]) > 1e-14){
274  // sinh term can be converted by setting a_k = C, b_k = 2*D, c_k = -1, d_k = 1
275  n.push_back(constants[1]);
276  t.push_back(-2*constants[2]/Tc);
277  c.push_back(1);
278  d.push_back(-1);
279  }
280  if (std::abs(constants[3]) > 1e-14){
281  // cosh term can be converted by setting a_k = C, b_k = 2*D, c_k = 1, d_k = 1
282  n.push_back(-constants[3]);
283  t.push_back(-2*constants[4]/Tc);
284  c.push_back(1);
285  d.push_back(1);
286  }
287  if (alpha0.PlanckEinstein.is_enabled() == true){
288  alpha0.PlanckEinstein.extend(n, t, c, d);
289  }
290  else{
291  alpha0.PlanckEinstein = IdealHelmholtzPlanckEinsteinGeneralized(n, t, c, d);
292  }
293  }
294  else if (!type.compare("IdealGasHelmholtzEnthalpyEntropyOffset"))
295  {
296  CoolPropDbl a1 = cpjson::get_double(contribution, "a1");
297  CoolPropDbl a2 = cpjson::get_double(contribution, "a2");
298  std::string reference = cpjson::get_string(contribution, "reference");
299  alpha0.EnthalpyEntropyOffsetCore = IdealHelmholtzEnthalpyEntropyOffset(a1, a2, reference);
300  }
301  else
302  {
303  std::cout << format("Unsupported ideal-gas Helmholtz type: %s\n",type.c_str());
304  //throw ValueError(format("Unsupported ideal-gas Helmholtz type: %s",type.c_str()));
305  }
306  }
307  return alpha0;
308  };
309 protected:
311  void parse_environmental(rapidjson::Value &json, CoolPropFluid &fluid)
312  {
313  fluid.environment.ASHRAE34 = cpjson::get_string(json,"ASHRAE34");
314  fluid.environment.GWP20 = cpjson::get_double(json,"GWP20");
315  fluid.environment.GWP100 = cpjson::get_double(json,"GWP100");
316  fluid.environment.GWP500 = cpjson::get_double(json,"GWP500");
317  fluid.environment.HH = cpjson::get_double(json,"HH");
318  fluid.environment.FH = cpjson::get_double(json,"FH");
319  fluid.environment.PH = cpjson::get_double(json,"PH");
320  fluid.environment.ODP = cpjson::get_double(json,"ODP");
321  }
322 
324  void parse_EOS(rapidjson::Value &EOS_json, CoolPropFluid &fluid)
325  {
326  EquationOfState E;
327  fluid.EOSVector.push_back(E);
328 
329  EquationOfState &EOS = fluid.EOSVector.at(fluid.EOSVector.size()-1);
330 
331  // Universal gas constant [J/mol/K]
332  EOS.R_u = cpjson::get_double(EOS_json,"gas_constant");
333  EOS.molar_mass = cpjson::get_double(EOS_json,"molar_mass");
334  EOS.acentric = cpjson::get_double(EOS_json,"acentric");
335 
336  EOS.pseudo_pure = cpjson::get_bool(EOS_json, "pseudo_pure");
337  EOS.limits.Tmax = cpjson::get_double(EOS_json, "T_max");
338  EOS.limits.pmax = cpjson::get_double(EOS_json, "p_max");
339 
340  rapidjson::Value &reducing_state = EOS_json["STATES"]["reducing"];
341  rapidjson::Value &satminL_state = EOS_json["STATES"]["sat_min_liquid"];
342  rapidjson::Value &satminV_state = EOS_json["STATES"]["sat_min_vapor"];
343 
344  // Reducing state
345  EOS.reduce.T = cpjson::get_double(reducing_state,"T");
346  EOS.reduce.rhomolar = cpjson::get_double(reducing_state,"rhomolar");
347  EOS.reduce.p = cpjson::get_double(reducing_state,"p");
348  EOS.reduce.hmolar = cpjson::get_double(reducing_state,"hmolar");
349  EOS.reduce.smolar = cpjson::get_double(reducing_state,"smolar");
350 
351  EOS.sat_min_liquid.T = cpjson::get_double(satminL_state, "T");
352  EOS.sat_min_liquid.p = cpjson::get_double(satminL_state, "p");
353  EOS.sat_min_liquid.rhomolar = cpjson::get_double(satminL_state, "rhomolar");
354  EOS.sat_min_vapor.T = cpjson::get_double(satminV_state, "T");
355  EOS.sat_min_vapor.p = cpjson::get_double(satminV_state, "p");
356  EOS.sat_min_vapor.rhomolar = cpjson::get_double(satminV_state, "rhomolar");
357 
359  EOS.limits.Tmin = cpjson::get_double(satminL_state, "T");
360  EOS.ptriple = cpjson::get_double(satminL_state, "p");
361  EOS.Ttriple = EOS.limits.Tmin;
362 
363  // BibTex keys
364  EOS.BibTeX_EOS = cpjson::get_string(EOS_json,"BibTeX_EOS");
365  EOS.BibTeX_CP0 = cpjson::get_string(EOS_json,"BibTeX_CP0");
366 
367  EOS.alphar = parse_alphar(EOS_json["alphar"]);
368  EOS.alpha0 = parse_alpha0(EOS_json["alpha0"]);
369 
370  if (EOS_json["STATES"].HasMember("hs_anchor")){
371  rapidjson::Value &hs_anchor = EOS_json["STATES"]["hs_anchor"];
372  EOS.hs_anchor.T = cpjson::get_double(hs_anchor, "T");
373  EOS.hs_anchor.p = cpjson::get_double(hs_anchor, "p");
374  EOS.hs_anchor.rhomolar = cpjson::get_double(hs_anchor, "rhomolar");
375  EOS.hs_anchor.hmolar = cpjson::get_double(hs_anchor, "hmolar");
376  EOS.hs_anchor.smolar = cpjson::get_double(hs_anchor, "smolar");
377  }
378 
379  if (EOS_json["STATES"].HasMember("pressure_max_sat")){
380  rapidjson::Value &s = EOS_json["STATES"]["pressure_max_sat"];
381  EOS.max_sat_p.T = cpjson::get_double(s, "T");
382  EOS.max_sat_p.p = cpjson::get_double(s, "p");
383  EOS.max_sat_p.rhomolar = cpjson::get_double(s, "rhomolar");
384  if (s.HasMember("hmolar")){
385  EOS.max_sat_p.hmolar = cpjson::get_double(s, "hmolar");
386  EOS.max_sat_p.smolar = cpjson::get_double(s, "smolar");
387  }
388  }
389 
390  if (EOS_json["STATES"].HasMember("temperature_max_sat")){
391  rapidjson::Value &s = EOS_json["STATES"]["temperature_max_sat"];
392  EOS.max_sat_T.T = cpjson::get_double(s, "T");
393  EOS.max_sat_T.p = cpjson::get_double(s, "p");
394  EOS.max_sat_T.rhomolar = cpjson::get_double(s, "rhomolar");
395  if (s.HasMember("hmolar")){
396  EOS.max_sat_T.hmolar = cpjson::get_double(s, "hmolar");
397  EOS.max_sat_T.smolar = cpjson::get_double(s, "smolar");
398  }
399  }
400 
401  if (EOS_json.HasMember("critical_region_splines")){
402  rapidjson::Value &spline = EOS_json["critical_region_splines"];
403  EOS.critical_region_splines.T_min = cpjson::get_double(spline, "T_min");
404  EOS.critical_region_splines.T_max = cpjson::get_double(spline, "T_max");
405  EOS.critical_region_splines.rhomolar_min = cpjson::get_double(spline, "rhomolar_min");
406  EOS.critical_region_splines.rhomolar_max = cpjson::get_double(spline, "rhomolar_max");
407  EOS.critical_region_splines.cL = cpjson::get_double_array(spline["cL"]);
408  EOS.critical_region_splines.cV = cpjson::get_double_array(spline["cV"]);
409  EOS.critical_region_splines.enabled = true;
410  }
411 
412  // Validate the equation of state that was just created
413  EOS.validate();
414  }
415 
417  void parse_EOS_listing(rapidjson::Value &EOS_array, CoolPropFluid & fluid)
418  {
419  for (rapidjson::Value::ValueIterator itr = EOS_array.Begin(); itr != EOS_array.End(); ++itr)
420  {
421  parse_EOS(*itr,fluid);
422  }
423  };
424 
426  void parse_dilute_viscosity(rapidjson::Value &dilute, CoolPropFluid & fluid)
427  {
428  if (dilute.HasMember("hardcoded")){
429  std::string target = cpjson::get_string(dilute, "hardcoded");
430  if (!target.compare("Ethane")){
431  fluid.transport.viscosity_dilute.type = CoolProp::ViscosityDiluteVariables::VISCOSITY_DILUTE_ETHANE; return;
432  }
433  else if (!target.compare("Cyclohexane")){
434  fluid.transport.viscosity_dilute.type = CoolProp::ViscosityDiluteVariables::VISCOSITY_DILUTE_CYCLOHEXANE; return;
435  }
436  else{
437  throw ValueError(format("hardcoded dilute viscosity [%s] is not understood for fluid %s",target.c_str(),fluid.name.c_str()));
438  }
439  }
440  std::string type = cpjson::get_string(dilute, "type");
441  if (!type.compare("collision_integral")){
442  // Get a reference to the entry in the fluid instance
443  CoolProp::ViscosityDiluteGasCollisionIntegralData &CI = fluid.transport.viscosity_dilute.collision_integral;
444 
445  // Set the type flag
446  fluid.transport.viscosity_dilute.type = CoolProp::ViscosityDiluteVariables::VISCOSITY_DILUTE_COLLISION_INTEGRAL;
447 
448  // Load up the values
449  CI.a = cpjson::get_long_double_array(dilute["a"]);
450  CI.t = cpjson::get_long_double_array(dilute["t"]);
451  CI.molar_mass = cpjson::get_double(dilute, "molar_mass");
452  CI.C = cpjson::get_double(dilute, "C");
453  }
454  else if (!type.compare("kinetic_theory")){
455  fluid.transport.viscosity_dilute.type = CoolProp::ViscosityDiluteVariables::VISCOSITY_DILUTE_KINETIC_THEORY;
456  }
457  else if (!type.compare("powers_of_T")){
458  // Get a reference to the entry in the fluid instance
459  CoolProp::ViscosityDiluteGasPowersOfT &CI = fluid.transport.viscosity_dilute.powers_of_T;
460 
461  // Load up the values
462  CI.a = cpjson::get_long_double_array(dilute["a"]);
463  CI.t = cpjson::get_long_double_array(dilute["t"]);
464 
465  fluid.transport.viscosity_dilute.type = CoolProp::ViscosityDiluteVariables::VISCOSITY_DILUTE_POWERS_OF_T;
466  }
467  else if (!type.compare("powers_of_Tr")){
468  // Get a reference to the entry in the fluid instance
469  CoolProp::ViscosityDiluteGasPowersOfTr &CI = fluid.transport.viscosity_dilute.powers_of_Tr;
470  // Load up the values
471  CI.a = cpjson::get_long_double_array(dilute["a"]);
472  CI.t = cpjson::get_long_double_array(dilute["t"]);
473  CI.T_reducing = cpjson::get_double(dilute, "T_reducing");
474  fluid.transport.viscosity_dilute.type = CoolProp::ViscosityDiluteVariables::VISCOSITY_DILUTE_POWERS_OF_TR;
475  }
476  else if (!type.compare("collision_integral_powers_of_Tstar")){
477  // Get a reference to the entry in the fluid instance
479 
480  // Load up the values
481  CI.a = cpjson::get_long_double_array(dilute["a"]);
482  CI.t = cpjson::get_long_double_array(dilute["t"]);
483  CI.T_reducing = cpjson::get_double(dilute,"T_reducing");
484  CI.C = cpjson::get_double(dilute,"C");
485 
487  }
488  else{
489  throw ValueError(format("type [%s] is not understood for fluid %s",type.c_str(),fluid.name.c_str()));
490  }
491  };
492 
494  void parse_initial_density_viscosity(rapidjson::Value &initial_density, CoolPropFluid & fluid)
495  {
496  std::string type = cpjson::get_string(initial_density, "type");
497  if (!type.compare("Rainwater-Friend")){
498  // Get a reference to the entry in the fluid instance
499  CoolProp::ViscosityRainWaterFriendData &RF = fluid.transport.viscosity_initial.rainwater_friend;
500 
501  // Load up the values
502  RF.b = cpjson::get_long_double_array(initial_density["b"]);
503  RF.t = cpjson::get_long_double_array(initial_density["t"]);
504 
505  // Set the type flag
507  }
508  else if (!type.compare("empirical")){
509  // Get a reference to the entry in the fluid instance
510  CoolProp::ViscosityInitialDensityEmpiricalData &EM = fluid.transport.viscosity_initial.empirical;
511 
512  // Load up the values
513  EM.n = cpjson::get_long_double_array(initial_density["n"]);
514  EM.d = cpjson::get_long_double_array(initial_density["d"]);
515  EM.t = cpjson::get_long_double_array(initial_density["t"]);
516  EM.T_reducing = cpjson::get_double(initial_density,"T_reducing");
517  EM.rhomolar_reducing = cpjson::get_double(initial_density,"rhomolar_reducing");
518 
519  // Set the type flag
521  }
522  else{
523  throw ValueError(format("type [%s] is not understood for fluid %s",type.c_str(),fluid.name.c_str()));
524  }
525  };
526 
528  void parse_higher_order_viscosity(rapidjson::Value &higher, CoolPropFluid & fluid)
529  {
530  // First check for hardcoded higher-order term
531  if (higher.HasMember("hardcoded")){
532  std::string target = cpjson::get_string(higher,"hardcoded");
533  if (!target.compare("Hydrogen")){
534  fluid.transport.viscosity_higher_order.type = CoolProp::ViscosityHigherOrderVariables::VISCOSITY_HIGHER_ORDER_HYDROGEN; return;
535  }
536  else if (!target.compare("n-Hexane")){
537  fluid.transport.viscosity_higher_order.type = CoolProp::ViscosityHigherOrderVariables::VISCOSITY_HIGHER_ORDER_HEXANE; return;
538  }
539  else if (!target.compare("n-Heptane")){
540  fluid.transport.viscosity_higher_order.type = CoolProp::ViscosityHigherOrderVariables::VISCOSITY_HIGHER_ORDER_HEPTANE; return;
541  }
542  else if (!target.compare("Toluene")){
543  fluid.transport.viscosity_higher_order.type = CoolProp::ViscosityHigherOrderVariables::VISCOSITY_HIGHER_ORDER_TOLUENE; return;
544  }
545  else if (!target.compare("Ethane")){
546  fluid.transport.viscosity_higher_order.type = CoolProp::ViscosityHigherOrderVariables::VISCOSITY_HIGHER_ORDER_ETHANE; return;
547  }
548  else if (!target.compare("Benzene")){
549  fluid.transport.viscosity_higher_order.type = CoolProp::ViscosityHigherOrderVariables::VISCOSITY_HIGHER_ORDER_BENZENE; return;
550  }
551  else{
552  throw ValueError(format("hardcoded higher order viscosity term [%s] is not understood for fluid %s",target.c_str(), fluid.name.c_str()));
553  }
554  }
555 
556  std::string type = cpjson::get_string(higher, "type");
557  if (!type.compare("modified_Batschinski_Hildebrand")){
558  // Get a reference to the entry in the fluid instance to simplify the code that follows
560 
561  // Set the flag for the type of this model
563 
564  BH.T_reduce = cpjson::get_double(higher, "T_reduce");
565  BH.rhomolar_reduce = cpjson::get_double(higher, "rhomolar_reduce");
566  // Load up the values
567  BH.a = cpjson::get_long_double_array(higher["a"]);
568  BH.t1 = cpjson::get_long_double_array(higher["t1"]);
569  BH.d1 = cpjson::get_long_double_array(higher["d1"]);
570  BH.gamma = cpjson::get_long_double_array(higher["gamma"]);
571  BH.l = cpjson::get_long_double_array(higher["l"]);
572  assert(BH.a.size() == BH.t1.size());
573  assert(BH.a.size() == BH.d1.size());
574  assert(BH.a.size() == BH.gamma.size());
575  assert(BH.a.size() == BH.l.size());
576  BH.f = cpjson::get_long_double_array(higher["f"]);
577  BH.t2 = cpjson::get_long_double_array(higher["t2"]);
578  BH.d2 = cpjson::get_long_double_array(higher["d2"]);
579  assert(BH.f.size() == BH.t2.size());
580  assert(BH.f.size() == BH.d2.size());
581  BH.g = cpjson::get_long_double_array(higher["g"]);
582  BH.h = cpjson::get_long_double_array(higher["h"]);
583  assert(BH.g.size() == BH.h.size());
584  BH.p = cpjson::get_long_double_array(higher["p"]);
585  BH.q = cpjson::get_long_double_array(higher["q"]);
586  assert(BH.p.size() == BH.q.size());
587  }
588  else if (!type.compare("friction_theory")){
589  // Get a reference to the entry in the fluid instance to simplify the code that follows
590  CoolProp::ViscosityFrictionTheoryData &F = fluid.transport.viscosity_higher_order.friction_theory;
591 
592  // Set the flag for the type of this model
593  fluid.transport.viscosity_higher_order.type = CoolProp::ViscosityHigherOrderVariables::VISCOSITY_HIGHER_ORDER_FRICTION_THEORY;
594 
595  // Always need these terms
596  F.Ai = cpjson::get_long_double_array(higher["Ai"]);
597  F.Aa = cpjson::get_long_double_array(higher["Aa"]);
598  F.Aaa = cpjson::get_long_double_array(higher["Aaa"]);
599  F.Ar = cpjson::get_long_double_array(higher["Ar"]);
600 
601 
602  F.Na = cpjson::get_integer(higher,"Na");
603  F.Naa = cpjson::get_integer(higher,"Naa");
604  F.Nr = cpjson::get_integer(higher,"Nr");
605  F.Nrr = cpjson::get_integer(higher,"Nrr");
606  F.c1 = cpjson::get_double(higher,"c1");
607  F.c2 = cpjson::get_double(higher,"c2");
608  assert(F.Aa.size() == 3);
609  assert(F.Aaa.size() == 3);
610  assert(F.Ar.size() == 3);
611 
612  F.T_reduce = cpjson::get_double(higher,"T_reduce");
613 
614  if (higher.HasMember("Arr") && !higher.HasMember("Adrdr")){
615  F.Arr = cpjson::get_long_double_array(higher["Arr"]);
616  assert(F.Arr.size() == 3);
617  }
618  else if (higher.HasMember("Adrdr") && !higher.HasMember("Arr")){
619  F.Adrdr = cpjson::get_long_double_array(higher["Adrdr"]);
620  assert(F.Adrdr.size() == 3);
621  }
622  else{
623  throw ValueError(format("can only provide one of Arr or Adrdr for fluid %s",fluid.name.c_str()));
624  }
625  if (higher.HasMember("Aii")){
626  F.Aii = cpjson::get_long_double_array(higher["Aii"]);
627  F.Nii = cpjson::get_integer(higher,"Nii");
628  }
629  if (higher.HasMember("Aaaa") && higher.HasMember("Arrr")){
630  F.Aaaa = cpjson::get_long_double_array(higher["Aaaa"]);
631  F.Arrr = cpjson::get_long_double_array(higher["Arrr"]);
632  F.Naaa = cpjson::get_integer(higher,"Naaa");
633  F.Nrrr = cpjson::get_integer(higher,"Nrrr");
634  }
635 
636  }
637  else{
638  throw ValueError(format("type [%s] is not understood for fluid %s",type.c_str(), fluid.name.c_str()));
639  }
640  };
641 
642  void parse_ECS_conductivity(rapidjson::Value &conductivity, CoolPropFluid & fluid)
643  {
644  fluid.transport.conductivity_ecs.reference_fluid = cpjson::get_string(conductivity,"reference_fluid");
645 
646  // Parameters for correction polynomials
647  fluid.transport.conductivity_ecs.psi_a = cpjson::get_long_double_array(conductivity["psi"]["a"]);
648  fluid.transport.conductivity_ecs.psi_t = cpjson::get_long_double_array(conductivity["psi"]["t"]);
649  fluid.transport.conductivity_ecs.psi_rhomolar_reducing = cpjson::get_double(conductivity["psi"],"rhomolar_reducing");
650  fluid.transport.conductivity_ecs.f_int_a = cpjson::get_long_double_array(conductivity["f_int"]["a"]);
651  fluid.transport.conductivity_ecs.f_int_t = cpjson::get_long_double_array(conductivity["f_int"]["t"]);
652  fluid.transport.conductivity_ecs.f_int_T_reducing = cpjson::get_double(conductivity["f_int"],"T_reducing");
653 
654  fluid.transport.conductivity_using_ECS = true;
655  }
656 
657  void parse_ECS_viscosity(rapidjson::Value &viscosity, CoolPropFluid & fluid)
658  {
659  fluid.transport.viscosity_ecs.reference_fluid = cpjson::get_string(viscosity,"reference_fluid");
660 
661  // Parameters for correction polynomial
662  fluid.transport.viscosity_ecs.psi_a = cpjson::get_long_double_array(viscosity["psi"]["a"]);
663  fluid.transport.viscosity_ecs.psi_t = cpjson::get_long_double_array(viscosity["psi"]["t"]);
664  fluid.transport.viscosity_ecs.psi_rhomolar_reducing = cpjson::get_double(viscosity["psi"],"rhomolar_reducing");
665 
666  fluid.transport.viscosity_using_ECS = true;
667  }
668 
669  void parse_Chung_viscosity(rapidjson::Value &viscosity, CoolPropFluid &fluid)
670  {
671  // These in base SI units
672  fluid.transport.viscosity_Chung.rhomolar_critical = cpjson::get_double(viscosity, "rhomolar_critical");
673  fluid.transport.viscosity_Chung.T_critical = cpjson::get_double(viscosity, "T_critical");
674  fluid.transport.viscosity_Chung.molar_mass = cpjson::get_double(viscosity, "molar_mass");
675  fluid.transport.viscosity_Chung.dipole_moment_D = cpjson::get_double(viscosity, "dipole_moment_D");
676  fluid.transport.viscosity_Chung.acentric = cpjson::get_double(viscosity, "acentric");
677  fluid.transport.viscosity_using_Chung = true;
678  }
679 
680  void parse_rhosr_viscosity(rapidjson::Value &viscosity, CoolPropFluid &fluid)
681  {
682  fluid.transport.viscosity_rhosr.C = cpjson::get_double(viscosity, "C");
683  fluid.transport.viscosity_rhosr.c_liq = cpjson::get_double_array(viscosity, "c_liq");
684  fluid.transport.viscosity_rhosr.c_vap = cpjson::get_double_array(viscosity, "c_vap");
685  fluid.transport.viscosity_rhosr.rhosr_critical = cpjson::get_double(viscosity, "rhosr_critical");
686  fluid.transport.viscosity_rhosr.x_crossover = cpjson::get_double(viscosity, "x_crossover");
687  fluid.transport.viscosity_using_rhosr = true;
688  }
689 
690 
692  void parse_viscosity(rapidjson::Value &viscosity, CoolPropFluid & fluid)
693  {
694  // If an array, use the first one, and then stop;
695  if (viscosity.IsArray()){
696  rapidjson::Value::ValueIterator itr = viscosity.Begin();
697  parse_viscosity(*itr, fluid);
698  return;
699  }
700 
701  // Load the BibTeX key
702  fluid.transport.BibTeX_viscosity = cpjson::get_string(viscosity,"BibTeX");
703 
704  // Set the Lennard-Jones 12-6 potential variables, or approximate them from method of Chung
705  if (!viscosity.HasMember("sigma_eta")|| !viscosity.HasMember("epsilon_over_k")){
706  default_transport(fluid);
707  }
708  else{
709  fluid.transport.sigma_eta = cpjson::get_double(viscosity, "sigma_eta");
710  fluid.transport.epsilon_over_k = cpjson::get_double(viscosity, "epsilon_over_k");
711  }
712 
713  // If it is using ECS, set ECS parameters and quit
714  if (viscosity.HasMember("type") && !cpjson::get_string(viscosity, "type").compare("ECS")){
715  parse_ECS_viscosity(viscosity, fluid);
716  return;
717  }
718 
719  // If it is using rho*sr CS, set parameters and quit
720  if (viscosity.HasMember("type") && !cpjson::get_string(viscosity, "type").compare("rhosr-CS")){
721  parse_rhosr_viscosity(viscosity, fluid);
722  return;
723  }
724 
725  // Use the method of Chung
726  // If it is using ECS, set ECS parameters and quit
727  if (viscosity.HasMember("type") && !cpjson::get_string(viscosity, "type").compare("Chung")){
728  parse_Chung_viscosity(viscosity, fluid);
729  return;
730  }
731 
732 
733  if (viscosity.HasMember("hardcoded")){
734  std::string target = cpjson::get_string(viscosity,"hardcoded");
735  if (!target.compare("Water")){
737  }
738  else if (!target.compare("HeavyWater")){
740  }
741  else if (!target.compare("Helium")){
743  }
744  else if (!target.compare("R23")){
746  }
747  else if (!target.compare("Methanol")){
749  }
750  else if (!target.compare("m-Xylene")) {
752  }
753  else if (!target.compare("o-Xylene")) {
755  }
756  else if (!target.compare("p-Xylene")) {
758  }
759  else{
760  throw ValueError(format("hardcoded viscosity [%s] is not understood for fluid %s",target.c_str(), fluid.name.c_str()));
761  }
762  }
763 
764  // Load dilute viscosity term
765  if (viscosity.HasMember("dilute")){
766  parse_dilute_viscosity(viscosity["dilute"], fluid);
767  }
768  // Load initial density term
769  if (viscosity.HasMember("initial_density")){
770  parse_initial_density_viscosity(viscosity["initial_density"], fluid);
771  }
772  // Load higher_order term
773  if (viscosity.HasMember("higher_order")){
774  parse_higher_order_viscosity(viscosity["higher_order"], fluid);
775  }
776  };
777 
779  void parse_dilute_conductivity(rapidjson::Value &dilute, CoolPropFluid & fluid)
780  {
781  if (dilute.HasMember("hardcoded")){
782  std::string target = cpjson::get_string(dilute, "hardcoded");
783  if (!target.compare("CO2")){
784  fluid.transport.conductivity_dilute.type = CoolProp::ConductivityDiluteVariables::CONDUCTIVITY_DILUTE_CO2; return;
785  }
786  else if (!target.compare("Ethane")){
787  fluid.transport.conductivity_dilute.type = CoolProp::ConductivityDiluteVariables::CONDUCTIVITY_DILUTE_ETHANE; return;
788  }
789  else if (!target.compare("none")){
790  fluid.transport.conductivity_dilute.type = CoolProp::ConductivityDiluteVariables::CONDUCTIVITY_DILUTE_NONE; return;
791  }
792  else{
793  throw ValueError(format("hardcoded dilute conductivity term [%s] is not understood for fluid %s",target.c_str(), fluid.name.c_str()));
794  }
795  }
796  std::string type = cpjson::get_string(dilute, "type");
797  if (!type.compare("ratio_of_polynomials")){
798  // Get a reference to the entry in the fluid instance
799  CoolProp::ConductivityDiluteRatioPolynomialsData &data = fluid.transport.conductivity_dilute.ratio_polynomials;
800 
801  // Set the type flag
802  fluid.transport.conductivity_dilute.type = CoolProp::ConductivityDiluteVariables::CONDUCTIVITY_DILUTE_RATIO_POLYNOMIALS;
803 
804  // Load up the values
805  data.A = cpjson::get_long_double_array(dilute["A"]);
806  data.B = cpjson::get_long_double_array(dilute["B"]);
807  data.n = cpjson::get_long_double_array(dilute["n"]);
808  data.m = cpjson::get_long_double_array(dilute["m"]);
809  data.T_reducing = cpjson::get_double(dilute, "T_reducing");
810  }
811  else if (!type.compare("eta0_and_poly")){
812  // Get a reference to the entry in the fluid instance
813  CoolProp::ConductivityDiluteEta0AndPolyData &data = fluid.transport.conductivity_dilute.eta0_and_poly;
814 
815  // Set the type flag
816  fluid.transport.conductivity_dilute.type = CoolProp::ConductivityDiluteVariables::CONDUCTIVITY_DILUTE_ETA0_AND_POLY;
817 
818  // Load up the values
819  data.A = cpjson::get_long_double_array(dilute["A"]);
820  data.t = cpjson::get_long_double_array(dilute["t"]);
821  }
822  else{
823  throw ValueError(format("type [%s] is not understood for fluid %s",type.c_str(),fluid.name.c_str()));
824  }
825  };
826 
828  void parse_residual_conductivity(rapidjson::Value &dilute, CoolPropFluid & fluid)
829  {
830  if (dilute.HasMember("hardcoded")){
831  std::string target = cpjson::get_string(dilute, "hardcoded");
832  if (!target.compare("CO2")){
833  fluid.transport.conductivity_residual.type = CoolProp::ConductivityResidualVariables::CONDUCTIVITY_RESIDUAL_CO2; return;
834  }
835  else{
836  throw ValueError(format("hardcoded residual conductivity term [%s] is not understood for fluid %s",target.c_str(), fluid.name.c_str()));
837  }
838  }
839  std::string type = cpjson::get_string(dilute, "type");
840  if (!type.compare("polynomial")){
841  // Get a reference to the entry in the fluid instance
842  CoolProp::ConductivityResidualPolynomialData &data = fluid.transport.conductivity_residual.polynomials;
843 
844  // Set the type flag
845  fluid.transport.conductivity_residual.type = CoolProp::ConductivityResidualVariables::CONDUCTIVITY_RESIDUAL_POLYNOMIAL;
846 
847  // Load up the values
848  data.B = cpjson::get_long_double_array(dilute["B"]);
849  data.d = cpjson::get_long_double_array(dilute["d"]);
850  data.t = cpjson::get_long_double_array(dilute["t"]);
851  data.T_reducing = cpjson::get_double(dilute, "T_reducing");
852  data.rhomass_reducing = cpjson::get_double(dilute, "rhomass_reducing");
853  }
854  else if (!type.compare("polynomial_and_exponential")){
855  // Get a reference to the entry in the fluid instance
856  CoolProp::ConductivityResidualPolynomialAndExponentialData &data = fluid.transport.conductivity_residual.polynomial_and_exponential;
857 
858  // Set the type flag
859  fluid.transport.conductivity_residual.type = CoolProp::ConductivityResidualVariables::CONDUCTIVITY_RESIDUAL_POLYNOMIAL_AND_EXPONENTIAL;
860 
861  // Load up the values
862  data.A = cpjson::get_long_double_array(dilute["A"]);
863  data.d = cpjson::get_long_double_array(dilute["d"]);
864  data.t = cpjson::get_long_double_array(dilute["t"]);
865  data.gamma = cpjson::get_long_double_array(dilute["gamma"]);
866  data.l = cpjson::get_long_double_array(dilute["l"]);
867  }
868  else{
869  throw ValueError(format("type [%s] is not understood for fluid %s",type.c_str(),fluid.name.c_str()));
870  }
871  };
872 
873  void parse_critical_conductivity(rapidjson::Value &critical, CoolPropFluid & fluid)
874  {
875  if (critical.HasMember("hardcoded")){
876  std::string target = cpjson::get_string(critical, "hardcoded");
877  if (!target.compare("R123")){
878  fluid.transport.conductivity_critical.type = CoolProp::ConductivityCriticalVariables::CONDUCTIVITY_CRITICAL_R123; return;
879  }
880  else if (!target.compare("Ammonia")){
881  fluid.transport.conductivity_critical.type = CoolProp::ConductivityCriticalVariables::CONDUCTIVITY_CRITICAL_AMMONIA; return;
882  }
883  else if (!target.compare("CarbonDioxideScalabrinJPCRD2006")){
884  fluid.transport.conductivity_critical.type = CoolProp::ConductivityCriticalVariables::CONDUCTIVITY_CRITICAL_CARBONDIOXIDE_SCALABRIN_JPCRD_2006; return;
885  }
886  else if (!target.compare("None")){
887  fluid.transport.conductivity_critical.type = CoolProp::ConductivityCriticalVariables::CONDUCTIVITY_CRITICAL_NONE; return;
888  }
889  else{
890  throw ValueError(format("critical conductivity term [%s] is not understood for fluid %s",target.c_str(), fluid.name.c_str()));
891  }
892  }
893  std::string type = cpjson::get_string(critical, "type");
894  if (!type.compare("simplified_Olchowy_Sengers")){
896  CoolProp::ConductivityCriticalSimplifiedOlchowySengersData &data = fluid.transport.conductivity_critical.Olchowy_Sengers;
897 
898  // Set the type flag
899  fluid.transport.conductivity_critical.type = CoolProp::ConductivityCriticalVariables::CONDUCTIVITY_CRITICAL_SIMPLIFIED_OLCHOWY_SENGERS;
900 
901  // Set values if they are found - otherwise fall back to default values
902  if (critical.HasMember("qD")){ data.qD = cpjson::get_double(critical,"qD"); }
903  if (critical.HasMember("zeta0")){ data.zeta0 = cpjson::get_double(critical,"zeta0"); }
904  if (critical.HasMember("GAMMA")){ data.GAMMA = cpjson::get_double(critical,"GAMMA"); }
905  if (critical.HasMember("gamma")){ data.gamma = cpjson::get_double(critical,"gamma"); }
906  if (critical.HasMember("R0")){ data.R0 = cpjson::get_double(critical,"R0"); }
907  if (critical.HasMember("T_ref")){ data.T_ref = cpjson::get_double(critical,"T_ref"); }
908  }
909  else{
910  throw ValueError(format("type [%s] is not understood for fluid %s",type.c_str(),fluid.name.c_str()));
911  }
912  };
913 
915  void parse_thermal_conductivity(rapidjson::Value &conductivity, CoolPropFluid & fluid)
916  {
917  // Load the BibTeX key
918  fluid.transport.BibTeX_conductivity = cpjson::get_string(conductivity,"BibTeX");
919 
920  // If it is using ECS, set ECS parameters and quit
921  if (conductivity.HasMember("type") && !cpjson::get_string(conductivity, "type").compare("ECS")){
922  parse_ECS_conductivity(conductivity, fluid);
923  return;
924  }
925 
926  if (conductivity.HasMember("hardcoded")){
927  std::string target = cpjson::get_string(conductivity, "hardcoded");
928  if (!target.compare("Water")){
930  }
931  else if (!target.compare("HeavyWater")){
933  }
934  else if (!target.compare("Methane")){
936  }
937  else if (!target.compare("R23")){
939  }
940  else if (!target.compare("Helium")){
942  }
943  else{
944  throw ValueError(format("hardcoded residual conductivity term [%s] is not understood for fluid %s",target.c_str(), fluid.name.c_str()));
945  }
946  }
947 
948  // Load dilute conductivity term
949  if (conductivity.HasMember("dilute")){
950  parse_dilute_conductivity(conductivity["dilute"], fluid);
951  }
952  // Load residual conductivity term
953  if (conductivity.HasMember("residual")){
954  parse_residual_conductivity(conductivity["residual"], fluid);
955  }
956  // Load critical conductivity term
957  if (conductivity.HasMember("critical")){
958  parse_critical_conductivity(conductivity["critical"], fluid);
959  }
960  };
961 
963  void parse_transport(rapidjson::Value &transport, CoolPropFluid & fluid)
964  {
965 
966  // Parse viscosity
967  if (transport.HasMember("viscosity")){
968  parse_viscosity(transport["viscosity"],fluid);
969  fluid.transport.viscosity_model_provided = true;
970  }
971 
972  // Parse thermal conductivity
973  if (transport.HasMember("conductivity")){
974  parse_thermal_conductivity(transport["conductivity"],fluid);
975  fluid.transport.conductivity_model_provided = true;
976  }
977  };
978 
979  void default_transport(CoolPropFluid & fluid)
980  {
981  // Use the method of Chung to approximate the values for epsilon_over_k and sigma_eta
982  // Chung, T.-H.; Ajlan, M.; Lee, L. L.; Starling, K. E. Generalized Multiparameter Correlation for Nonpolar and Polar Fluid Transport Properties. Ind. Eng. Chem. Res. 1988, 27, 671-679.
983  // rhoc needs to be in mol/L to yield a sigma in nm,
984  CoolPropDbl rho_crit_molar = fluid.EOS().reduce.rhomolar/1000.0;// [mol/m3 to mol/L]
985  CoolPropDbl Tc = fluid.EOS().reduce.T;
986  fluid.transport.sigma_eta = 0.809/pow(rho_crit_molar, static_cast<CoolPropDbl>(1.0/3.0))/1e9; // 1e9 is to convert from nm to m
987  fluid.transport.epsilon_over_k = Tc/1.2593; // [K]
988  }
989 
990  void parse_melting_line(rapidjson::Value &melting_line, CoolPropFluid & fluid)
991  {
992  fluid.ancillaries.melting_line.T_m = cpjson::get_double(melting_line, "T_m");
993  fluid.ancillaries.melting_line.BibTeX = cpjson::get_string(melting_line, "BibTeX");
994 
995  if (melting_line.HasMember("type"))
996  {
997  std::string type = cpjson::get_string(melting_line, "type");
998  if (!type.compare("Simon"))
999  {
1000  rapidjson::Value &parts = melting_line["parts"];
1002  for (rapidjson::Value::ValueIterator itr = parts.Begin(); itr != parts.End(); ++itr)
1003  {
1005  data.a = cpjson::get_double((*itr),"a");
1006  data.c = cpjson::get_double((*itr),"c");
1007  data.T_min = cpjson::get_double((*itr),"T_min");
1008  data.T_max = cpjson::get_double((*itr),"T_max");
1009  data.T_0 = cpjson::get_double((*itr),"T_0");
1010  data.p_0 = cpjson::get_double((*itr),"p_0");
1011  fluid.ancillaries.melting_line.simon.parts.push_back(data);
1012  }
1013  }
1014  else if (!type.compare("polynomial_in_Tr"))
1015  {
1016  rapidjson::Value &parts = melting_line["parts"];
1018  for (rapidjson::Value::ValueIterator itr = parts.Begin(); itr != parts.End(); ++itr)
1019  {
1021  data.a = cpjson::get_long_double_array((*itr),"a");
1022  data.t = cpjson::get_long_double_array((*itr),"t");
1023  data.T_min = cpjson::get_double((*itr),"T_min");
1024  data.T_max = cpjson::get_double((*itr),"T_max");
1025  data.T_0 = cpjson::get_double((*itr),"T_0");
1026  data.p_0 = cpjson::get_double((*itr),"p_0");
1027  fluid.ancillaries.melting_line.polynomial_in_Tr.parts.push_back(data);
1028  }
1029  }
1030  else if (!type.compare("polynomial_in_Theta"))
1031  {
1032  rapidjson::Value &parts = melting_line["parts"];
1034  for (rapidjson::Value::ValueIterator itr = parts.Begin(); itr != parts.End(); ++itr)
1035  {
1037  data.a = cpjson::get_long_double_array((*itr),"a");
1038  data.t = cpjson::get_long_double_array((*itr),"t");
1039  data.T_min = cpjson::get_double((*itr),"T_min");
1040  data.T_max = cpjson::get_double((*itr),"T_max");
1041  data.T_0 = cpjson::get_double((*itr),"T_0");
1042  data.p_0 = cpjson::get_double((*itr),"p_0");
1043  fluid.ancillaries.melting_line.polynomial_in_Theta.parts.push_back(data);
1044  }
1045  }
1046  else{
1047  throw ValueError(format("melting line type [%s] is not understood for fluid %s", type.c_str(), fluid.name.c_str()));
1048  }
1049  // Set the limits for the melting line curve
1050  fluid.ancillaries.melting_line.set_limits();
1051  }
1052  else{
1053  throw ValueError(format("melting line does not have \"type\" for fluid %s", fluid.name.c_str()));
1054  }
1055  };
1056 
1058  void parse_states(rapidjson::Value &states, CoolPropFluid & fluid)
1059  {
1060  if (!states.HasMember("critical")){ throw ValueError(format("fluid[\"STATES\"] [%s] does not have \"critical\" member",fluid.name.c_str())); }
1061  rapidjson::Value &crit = states["critical"];
1062  fluid.crit.T = cpjson::get_double(crit, "T");
1063  fluid.crit.p = cpjson::get_double(crit, "p");
1064  fluid.crit.rhomolar = cpjson::get_double(crit, "rhomolar");
1065  fluid.crit.hmolar = cpjson::get_double(crit, "hmolar");
1066  fluid.crit.smolar = cpjson::get_double(crit, "smolar");
1067 
1068  if (!states.HasMember("triple_liquid")){ throw ValueError(format("fluid[\"STATES\"] [%s] does not have \"triple_liquid\" member",fluid.name.c_str())); }
1069  rapidjson::Value &triple_liquid = states["triple_liquid"];
1070  if (triple_liquid.ObjectEmpty()){
1071  // State is empty - probably because the triple point temperature is below the minimum saturation temperature
1072  fluid.triple_liquid.T = -1;
1073  fluid.triple_liquid.p = -1;
1074  fluid.triple_liquid.rhomolar = -1;
1075  fluid.triple_liquid.hmolar = _HUGE;
1076  fluid.triple_liquid.smolar = _HUGE;
1077  }
1078  else{
1079  fluid.triple_liquid.T = cpjson::get_double(triple_liquid, "T");
1080  fluid.triple_liquid.p = cpjson::get_double(triple_liquid, "p");
1081  fluid.triple_liquid.rhomolar = cpjson::get_double(triple_liquid, "rhomolar");
1082  fluid.triple_liquid.hmolar = cpjson::get_double(triple_liquid, "hmolar");
1083  fluid.triple_liquid.smolar = cpjson::get_double(triple_liquid, "smolar");
1084  }
1085 
1086  if (!states.HasMember("triple_vapor")){ throw ValueError(format("fluid[\"STATES\"] [%s] does not have \"triple_vapor\" member",fluid.name.c_str())); }
1087  rapidjson::Value &triple_vapor = states["triple_vapor"];
1088  if (triple_vapor.ObjectEmpty()){
1089  // State is empty - probably because the triple point temperature is below the minimum saturation temperature
1090  fluid.triple_vapor.T = -1;
1091  fluid.triple_vapor.p = -1;
1092  fluid.triple_vapor.rhomolar = -1;
1093  fluid.triple_vapor.hmolar = _HUGE;
1094  fluid.triple_vapor.smolar = _HUGE;
1095  }
1096  else{
1097  fluid.triple_vapor.T = cpjson::get_double(triple_vapor, "T");
1098  fluid.triple_vapor.p = cpjson::get_double(triple_vapor, "p");
1099  fluid.triple_vapor.rhomolar = cpjson::get_double(triple_vapor, "rhomolar");
1100  fluid.triple_vapor.hmolar = cpjson::get_double(triple_vapor, "hmolar");
1101  fluid.triple_vapor.smolar = cpjson::get_double(triple_vapor, "smolar");
1102  }
1103  };
1104 
1106  void parse_ancillaries(rapidjson::Value &ancillaries, CoolPropFluid & fluid)
1107  {
1108  if (!ancillaries.HasMember("rhoL") || !ancillaries.HasMember("rhoV")){
1109  throw ValueError("Ancillary curves for either rhoL or rhoV are missing");
1110  }
1111  fluid.ancillaries.rhoL = SaturationAncillaryFunction(ancillaries["rhoL"]);
1112  fluid.ancillaries.rhoV = SaturationAncillaryFunction(ancillaries["rhoV"]);
1113 
1114  // If a pseudo-pure fluid, has pL and pV curves
1115  if (ancillaries.HasMember("pL") && ancillaries.HasMember("pV")){
1116  fluid.ancillaries.pL = SaturationAncillaryFunction(ancillaries["pL"]);
1117  fluid.ancillaries.pV = SaturationAncillaryFunction(ancillaries["pV"]);
1118  }
1119  // Otherwise has a single pS curve and not pL and not pV
1120  else if (!ancillaries.HasMember("pL") && !ancillaries.HasMember("pV") && ancillaries.HasMember("pS")){
1121  fluid.ancillaries.pL = SaturationAncillaryFunction(ancillaries["pS"]);
1122  fluid.ancillaries.pV = SaturationAncillaryFunction(ancillaries["pS"]);
1123  }
1124  else{
1125  throw ValueError("Pressure ancillary curves are missing or invalid");
1126  }
1127 
1128  if (ancillaries.HasMember("hL")){
1129  fluid.ancillaries.hL = SaturationAncillaryFunction(ancillaries["hL"]);
1130  }
1131  else{
1132  if (get_debug_level() > 0){ std::cout << "Missing hL ancillary for fluid " << fluid.name; }
1133  }
1134  if (ancillaries.HasMember("hLV")){
1135  fluid.ancillaries.hLV = SaturationAncillaryFunction(ancillaries["hLV"]);
1136  }
1137  else{
1138  if (get_debug_level() > 0){ std::cout << "Missing hLV ancillary for fluid " << fluid.name; }
1139  }
1140 
1141  if (ancillaries.HasMember("sL")){
1142  fluid.ancillaries.sL = SaturationAncillaryFunction(ancillaries["sL"]);
1143  }
1144  else{
1145  if (get_debug_level() > 0){ std::cout << "Missing sL ancillary for fluid " << fluid.name; }
1146  }
1147  if (ancillaries.HasMember("sLV")){
1148  fluid.ancillaries.sLV = SaturationAncillaryFunction(ancillaries["sLV"]);
1149  }
1150  else{
1151  if (get_debug_level() > 0){ std::cout << "Missing sLV ancillary for fluid " << fluid.name; }
1152  }
1153  if (!ValidNumber(fluid.ancillaries.sL.get_Tmin()) && get_debug_level()>0){
1154  std::cout << "Tmin invalid for sL for " << fluid.name << std::endl;
1155  }
1156 
1157  };
1158 
1160  void parse_surface_tension(rapidjson::Value &surface_tension, CoolPropFluid & fluid)
1161  {
1162  fluid.ancillaries.surface_tension = SurfaceTensionCorrelation(surface_tension);
1163  };
1164 
1166  void validate(CoolPropFluid & fluid)
1167  {
1168  assert(fluid.EOSVector.size() > 0);
1169  assert(fluid.CAS.length() > 0);
1170  assert(fluid.name.length() > 0);
1171  }
1172 public:
1173 
1174  // Default constructor;
1175  JSONFluidLibrary(){
1176  _is_empty = true;
1177  };
1178  bool is_empty(void){ return _is_empty;};
1179 
1181  static void add_many(const std::string &JSON_string);
1182 
1184  void add_many(rapidjson::Value &listing);
1185 
1186  void add_one(rapidjson::Value &fluid_json);
1187 
1188  std::string get_JSONstring(const std::string &key)
1189  {
1190  // Try to find it
1191  std::map<std::string, std::size_t>::const_iterator it = string_to_index_map.find(key);
1192  if (it != string_to_index_map.end()){
1193 
1194  std::map<std::size_t, std::string>::const_iterator it2 = JSONstring_map.find(it->second);
1195  if (it2 != JSONstring_map.end()){
1196  // Then, load the fluids we would like to add
1197  rapidjson::Document doc;
1198  cpjson::JSON_string_to_rapidjson(it2->second, doc);
1199  rapidjson::Document doc2; doc2.SetArray();
1200  doc2.PushBack(doc, doc.GetAllocator());
1201  return cpjson::json2string(doc2);
1202  }
1203  else{
1204  throw ValueError(format("Unable to obtain JSON string for this identifier [%d]", it->second));
1205  }
1206  }
1207  else{
1208  throw ValueError(format("Unable to obtain index for this identifier [%s]", key.c_str()));
1209  }
1210  }
1211 
1213 
1216  CoolPropFluid get(const std::string &key)
1217  {
1218  // Try to find it
1219  std::map<std::string, std::size_t>::const_iterator it = string_to_index_map.find(key);
1220  // If it is found
1221  if (it != string_to_index_map.end()){
1222  return get(it->second);
1223  }
1224  else{
1225  // Here we check for the use of a cubic Helmholtz energy transformation for a multi-fluid model
1226  std::vector<std::string> endings; endings.push_back("-SRK"); endings.push_back("-PengRobinson");
1227  for (std::vector<std::string>::const_iterator end = endings.begin(); end != endings.end(); ++end){
1228  if (endswith(key, *end)){
1229  std::string used_name = key.substr(0, key.size()-(*end).size());
1230  it = string_to_index_map.find(used_name);
1231  if (it != string_to_index_map.end()){
1232  // We found the name of the fluid within the library of multiparameter
1233  // Helmholtz-explicit models. We will load its parameters from the
1234  // multiparameter EOS
1235  //
1236  CoolPropFluid fluid = get(it->second);
1237  // Remove all the residual contributions to the Helmholtz energy
1238  fluid.EOSVector[0].alphar.empty_the_EOS();
1239  // Get the parameters for the cubic EOS
1240  CoolPropDbl Tc = fluid.EOSVector[0].reduce.T;
1241  CoolPropDbl pc = fluid.EOSVector[0].reduce.p;
1242  CoolPropDbl rhomolarc = fluid.EOSVector[0].reduce.rhomolar;
1243  CoolPropDbl acentric = fluid.EOSVector[0].acentric;
1244  CoolPropDbl R = 8.3144598; // fluid.EOSVector[0].R_u;
1245  // Set the cubic contribution to the residual Helmholtz energy
1246  shared_ptr<AbstractCubic> ac;
1247  if (*end == "-SRK"){
1248  ac.reset(new SRK(Tc, pc, acentric, R));
1249  }
1250  else if (*end == "-PengRobinson"){
1251  ac.reset(new PengRobinson(Tc, pc, acentric, R));
1252  }
1253  else {
1254  throw CoolProp::ValueError(format("Unable to match this ending [%s]", (*end).c_str()));
1255  }
1256  ac->set_Tr(Tc);
1257  ac->set_rhor(rhomolarc);
1258  fluid.EOSVector[0].alphar.cubic = ResidualHelmholtzGeneralizedCubic(ac);
1259  return fluid;
1260  }
1261  else{
1262  // Let's look in the library of cubic EOS
1263  CubicLibrary::CubicsValues vals = CubicLibrary::get_cubic_values(used_name);
1264  // Set the cubic contribution to the residual Helmholtz energy
1265  shared_ptr<AbstractCubic> ac;
1266  if (*end == "-SRK") {
1267  ac.reset(new SRK(vals.Tc, vals.pc, vals.acentric, get_config_double(R_U_CODATA)));
1268  }
1269  else if (*end == "-PengRobinson") {
1270  ac.reset(new PengRobinson(vals.Tc, vals.pc, vals.acentric, get_config_double(R_U_CODATA)));
1271  }
1272  else{
1273  throw CoolProp::ValueError(format("Unable to match this ending [%s]",(*end).c_str()));
1274  }
1275  ac->set_Tr(vals.Tc);
1276  if (vals.rhomolarc > 0){
1277  ac->set_rhor(vals.rhomolarc);
1278  }
1279  else{
1280  // Curve fit from all the pure fluids in CoolProp (thanks to recommendation of A. Kazakov)
1281  double v_c_Lmol = 2.14107171795*(vals.Tc/vals.pc*1000)+0.00773144012514; // [L/mol]
1282  ac->set_rhor(1/(v_c_Lmol/1000.0));
1283  }
1284  if (vals.alpha_type == "Twu") {
1285  std::vector<double> &c = vals.alpha_coeffs;
1286  ac->set_C_Twu(0, c[0], c[1], c[2]);
1287  }
1288  CoolPropFluid fluid;
1289  fluid.CAS = vals.CAS;
1290  EquationOfState E;
1291  E.acentric = vals.acentric;
1292  E.sat_min_liquid.T = _HUGE;
1293  E.sat_min_liquid.p = _HUGE;
1294  E.reduce.T = vals.Tc;
1295  E.reduce.p = vals.pc;
1296  E.reduce.rhomolar = ac->get_rhor();
1297  fluid.EOSVector.push_back(E);
1298  fluid.EOS().alphar.cubic = ResidualHelmholtzGeneralizedCubic(ac);
1299  fluid.EOS().alpha0 = vals.alpha0;
1300  fluid.crit.T = vals.Tc;
1301  fluid.crit.p = vals.pc;
1302  fluid.crit.rhomolar = ac->get_rhor();
1303 
1304  return fluid;
1305  }
1306  }
1307  }
1308  throw ValueError(format("key [%s] was not found in string_to_index_map in JSONFluidLibrary", key.c_str()));
1309  }
1310  };
1311 
1313 
1316  CoolPropFluid get(std::size_t key)
1317  {
1318  // Try to find it
1319  std::map<std::size_t, CoolPropFluid>::iterator it = fluid_map.find(key);
1320  // If it is found
1321  if (it != fluid_map.end()){
1322  return it->second;
1323  }
1324  else{
1325  throw ValueError(format("key [%d] was not found in JSONFluidLibrary",key));
1326  }
1327  };
1328  void set_fluid_enthalpy_entropy_offset(const std::string &fluid, double delta_a1, double delta_a2, const std::string &ref);
1330  std::string get_fluid_list(void)
1331  {
1332  return strjoin(name_vector, get_config_string(LIST_STRING_DELIMITER));
1333  };
1334 };
1335 
1338 
1340 std::string get_fluid_list(void);
1341 
1343 CoolPropFluid get_fluid(const std::string &fluid_string);
1344 
1346 std::string get_fluid_as_JSONstring(const std::string &indentifier);
1347 
1349 void set_fluid_enthalpy_entropy_offset(const std::string &fluid, double delta_a1, double delta_a2, const std::string &ref);
1350 
1351 } /* namespace CoolProp */
1352 #endif
Use TransportRoutines::viscosity_o_xylene_hardcoded.
Definition: CoolPropFluid.h:292
The core class for an equation of state.
Definition: CoolPropFluid.h:354
Definition: Ancillaries.h:143
Definition: CoolPropFluid.h:195
bool viscosity_using_rhosr
A flag for whether to use rho*sr CS model of Bell. False for no.
Definition: CoolPropFluid.h:321
ViscosityInitialDensityEmpiricalData empirical
Data for TransportRoutines::viscosity_initial_density_dependence_empirical.
Definition: CoolPropFluid.h:237
SimpleState sat_min_vapor
The saturated vapor state at the minimum saturation temperature.
Definition: CoolPropFluid.h:357
Use TransportRoutines::viscosity_hexane_higher_order_hardcoded.
Definition: CoolPropFluid.h:256
Use TransportRoutines::viscosity_p_xylene_hardcoded.
Definition: CoolPropFluid.h:293
ConductivityHardcodedEnum hardcoded_conductivity
Hardcoded flags for the conductivity.
Definition: CoolPropFluid.h:327
void parse_EOS_listing(rapidjson::Value &EOS_array, CoolPropFluid &fluid)
Parse the list of possible equations of state.
Definition: FluidLibrary.h:417
Use TransportRoutines::viscosity_heptane_higher_order_hardcoded.
Definition: CoolPropFluid.h:257
Ancillaries ancillaries
The set of ancillary equations for dewpoint, bubblepoint, surface tension, etc.
Definition: CoolPropFluid.h:501
static ResidualHelmholtzContainer parse_alphar(rapidjson::Value &jsonalphar)
Parse the contributions to the residual Helmholtz energy.
Definition: FluidLibrary.h:39
void parse_states(rapidjson::Value &states, CoolPropFluid &fluid)
Parse the critical state for the given EOS.
Definition: FluidLibrary.h:1058
bool pseudo_pure
Is a pseudo-pure fluid (true) or pure fluid (false)
Definition: CoolPropFluid.h:370
Variables for the dilute gas part.
Definition: CoolPropFluid.h:180
Use TransportRoutines::viscosity_m_xylene_hardcoded.
Definition: CoolPropFluid.h:291
CoolPropDbl epsilon_over_k
The Lennard-Jones 12-6 parameter.
Definition: CoolPropFluid.h:324
a polynomial in is in use
Definition: Ancillaries.h:212
Use TransportRoutines::viscosity_water_hardcoded.
Definition: CoolPropFluid.h:286
void parse_ancillaries(rapidjson::Value &ancillaries, CoolPropFluid &fluid)
Parse the critical state for the given EOS.
Definition: FluidLibrary.h:1106
Use TransportRoutines::conductivity_hardcoded_heavywater.
Definition: CoolPropFluid.h:298
void parse_dilute_viscosity(rapidjson::Value &dilute, CoolPropFluid &fluid)
Parse the transport properties.
Definition: FluidLibrary.h:426
double get_config_double(configuration_keys key)
Return the value of a double configuration key.
Definition: Configuration.cpp:92
Use TransportRoutines::viscosity_benzene_higher_order_hardcoded.
Definition: CoolPropFluid.h:259
void parse_dilute_conductivity(rapidjson::Value &dilute, CoolPropFluid &fluid)
Parse the transport properties.
Definition: FluidLibrary.h:779
void add_Gaussian(const std::vector< CoolPropDbl > &n, const std::vector< CoolPropDbl > &d, const std::vector< CoolPropDbl > &t, const std::vector< CoolPropDbl > &eta, const std::vector< CoolPropDbl > &epsilon, const std::vector< CoolPropDbl > &beta, const std::vector< CoolPropDbl > &gamma)
Add and convert an old-style Gaussian term to generalized form.
Definition: Helmholtz.h:305
std::vector< double > alpha_coeffs
The vector of coefficients for the alpha function.
Definition: CubicsLibrary.h:25
Definition: CoolPropFluid.h:100
This is generalized class that can be used to manage an ancillary curve, here they are ancillary curv...
Definition: Ancillaries.h:83
ViscosityDiluteCollisionIntegralPowersOfTstarData collision_integral_powers_of_Tstar
Data for TransportRoutines::viscosity_dilute_collision_integral_powers_of_T.
Definition: CoolPropFluid.h:213
JSONFluidLibrary & get_library(void)
Get a reference to the library instance.
Definition: FluidLibrary.cpp:316
ViscosityModifiedBatschinskiHildebrandData modified_Batschinski_Hildebrand
Data for TransportRoutines::viscosity_higher_order_modified_Batschinski_Hildebrand.
Definition: CoolPropFluid.h:265
int get_debug_level()
Get the debug level.
Definition: CoolProp.cpp:63
EOSLimits limits
Limits on the EOS.
Definition: CoolPropFluid.h:364
Definition: CoolPropFluid.h:223
ViscosityRainWaterFriendData rainwater_friend
Data for TransportRoutines::viscosity_initial_density_dependence_Rainwater_Friend.
Definition: CoolPropFluid.h:236
std::string get_config_string(configuration_keys key)
Return the value of a string configuration key.
Definition: Configuration.cpp:95
The surface tension correlation class uses correlations for the surface tension that are all of the f...
Definition: Ancillaries.h:24
void parse_EOS(rapidjson::Value &EOS_json, CoolPropFluid &fluid)
Parse the Equation of state JSON entry.
Definition: FluidLibrary.h:324
Use TransportRoutines::viscosity_initial_density_dependence_empirical.
Definition: CoolPropFluid.h:232
std::vector< EquationOfState > EOSVector
The equations of state that could be used for this fluid.
Definition: CoolPropFluid.h:486
std::string BibTeX_CP0
The bibtex key for the ideal gas specific heat correlation.
Definition: CoolPropFluid.h:373
double acentric
Acentric factor (-)
Definition: CubicsLibrary.h:15
double ptriple
Triple point pressure (Pa)
Definition: CoolPropFluid.h:365
std::string BibTeX_conductivity
The BibTeX key for the conductivity model.
Definition: CoolPropFluid.h:316
int type
The data needed for a melting curve formed of segments that are polynomials in .
Definition: Ancillaries.h:224
SimpleState sat_min_liquid
The saturated liquid state at the minimum saturation temperature.
Definition: CoolPropFluid.h:357
SimpleState max_sat_T
The state at the maximum saturation temperature for pseudo-pure.
Definition: CoolPropFluid.h:357
double R_u
The universal gas constant used for this EOS (usually, but not always, 8.314472 J/mol/K) ...
Definition: CoolPropFluid.h:365
void parse_viscosity(rapidjson::Value &viscosity, CoolPropFluid &fluid)
Parse the transport properties.
Definition: FluidLibrary.h:692
Use TransportRoutines::conductivity_hardcoded_methane.
Definition: CoolPropFluid.h:301
bool conductivity_model_provided
A flag for whether thermal conductivity model is provided. False for no.
Definition: CoolPropFluid.h:323
Use TransportRoutines::viscosity_dilute_powers_of_T.
Definition: CoolPropFluid.h:207
void parse_environmental(rapidjson::Value &json, CoolPropFluid &fluid)
Parse the environmental parameters (ODP, GWP, etc.)
Definition: FluidLibrary.h:311
A thermophysical property provider for critical and reducing values as well as derivatives of Helmhol...
Definition: CoolPropFluid.h:476
void add_Power(const std::vector< CoolPropDbl > &n, const std::vector< CoolPropDbl > &d, const std::vector< CoolPropDbl > &t, const std::vector< CoolPropDbl > &l)
Add and convert an old-style power (polynomial) term to generalized form.
Definition: Helmholtz.h:259
std::string BibTeX_EOS
The bibtex key for the equation of state.
Definition: CoolPropFluid.h:373
Term in the ideal-gas specific heat equation that is based on Aly-Lee formulation ** Specific heat is...
Definition: Helmholtz.h:1149
CriticalRegionSplines critical_region_splines
A cubic spline in the form T = f(rho) for saturated liquid and saturated vapor curves in the near-cri...
Definition: CoolPropFluid.h:375
Definition: Helmholtz.h:988
Definition: Helmholtz.h:826
Definition: GeneralizedCubic.h:541
double Ttriple
Triple point temperature (K)
Definition: CoolPropFluid.h:365
Definition: CubicsLibrary.h:14
Definition: Helmholtz.h:495
static IdealHelmholtzContainer parse_alpha0(rapidjson::Value &jsonalpha0)
Parse the contributions to the ideal-gas Helmholtz energy.
Definition: FluidLibrary.h:168
Definition: CoolPropFluid.h:191
Definition: CoolPropFluid.h:129
double get_Tmin(void)
Get the minimum temperature in K.
Definition: Ancillaries.h:131
Definition: Helmholtz.h:434
The leading term in the EOS used to set the desired reference state.
Definition: Helmholtz.h:750
ViscosityFrictionTheoryData friction_theory
Data for TransportRoutines::viscosity_higher_order_friction_theory.
Definition: CoolPropFluid.h:266
Use TransportRoutines::viscosity_hydrogen_higher_order_hardcoded.
Definition: CoolPropFluid.h:255
const EquationOfState & EOS() const
Get a reference to the equation of state.
Definition: CoolPropFluid.h:484
Use TransportRoutines::viscosity_dilute_kinetic_theory.
Definition: CoolPropFluid.h:204
Definition: Exceptions.h:26
Definition: Helmholtz.h:853
Definition: CoolPropFluid.h:246
CoolPropDbl sigma_eta
The Lennard-Jones 12-6 parameter.
Definition: CoolPropFluid.h:324
double molar_mass
The molar mass in kg/mol (note NOT kg/kmol)
Definition: CoolPropFluid.h:365
Definition: CoolPropFluid.h:219
std::string BibTeX
BibTeX key for the melting curve in use.
Definition: Ancillaries.h:219
void add_Lemmon2005(const std::vector< CoolPropDbl > &n, const std::vector< CoolPropDbl > &d, const std::vector< CoolPropDbl > &t, const std::vector< CoolPropDbl > &l, const std::vector< CoolPropDbl > &m)
Add and convert a term from Lemmon and Jacobsen (2005) used for R125.
Definition: Helmholtz.h:362
Use TransportRoutines::viscosity_higher_order_friction_theory.
Definition: CoolPropFluid.h:261
std::string get_fluid_list(void)
Return a comma-separated list of fluid names.
Definition: FluidLibrary.h:1330
SimpleState reduce
Reducing state used for the EOS (usually, but not always, the critical point)
Definition: CoolPropFluid.h:357
Use TransportRoutines::viscosity_dilute_powers_of_Tr.
Definition: CoolPropFluid.h:208
The evaluator class for a melting curve formed of segments in the form.
Definition: Ancillaries.h:159
CoolPropDbl T_m
Melting temperature at 1 atmosphere.
Definition: Ancillaries.h:220
double acentric
The acentric factor .
Definition: CoolPropFluid.h:365
static void add_many(const std::string &JSON_string)
Add all the fluid entries in the JSON-encoded string passed in.
Definition: FluidLibrary.cpp:82
std::string CAS
The CAS number of the fluid.
Definition: CoolPropFluid.h:490
A simon-type curve is in use.
Definition: Ancillaries.h:210
std::string name
The name of the fluid.
Definition: CoolPropFluid.h:488
void parse_transport(rapidjson::Value &transport, CoolPropFluid &fluid)
Parse the transport properties.
Definition: FluidLibrary.h:963
MeltingLinePiecewisePolynomialInTrData polynomial_in_Tr
The data used for a Simon-style curve.
Definition: Ancillaries.h:222
SimpleState hs_anchor
A fixed anchor state at Tc*1.1 and rhoc*0.9 used as a reference state for enthalpy and entropy ancill...
Definition: CoolPropFluid.h:357
double pc
Critical pressure (Pa)
Definition: CubicsLibrary.h:15
std::string get_fluid_as_JSONstring(const std::string &identifier)
Get the fluid as a JSON string, suitable for modification and reloading.
Definition: FluidLibrary.cpp:326
SimpleState triple_vapor
The saturated vapor state at the triple point temperature.
Definition: CoolPropFluid.h:503
Use TransportRoutines::viscosity_ethane_higher_order_hardcoded.
Definition: CoolPropFluid.h:258
double rhomolarc
Critical density (mol/m3) (initialized to an invalid negative number)
Definition: CubicsLibrary.h:15
The evaluator class for a melting curve formed of segments in the form.
Definition: Ancillaries.h:185
void parse_surface_tension(rapidjson::Value &surface_tension, CoolPropFluid &fluid)
Parse the surface_tension.
Definition: FluidLibrary.h:1160
double Tc
Critical temperature (K)
Definition: CubicsLibrary.h:15
Definition: Helmholtz.h:692
Use TransportRoutines::viscosity_higher_order_modified_Batschinski_Hildebrand.
Definition: CoolPropFluid.h:254
SimpleState triple_liquid
The saturated liquid state at the triple point temperature.
Definition: CoolPropFluid.h:503
Use TransportRoutines::viscosity_dilute_ethane.
Definition: CoolPropFluid.h:205
Use TransportRoutines::conductivity_hardcoded_helium.
Definition: CoolPropFluid.h:300
void parse_higher_order_viscosity(rapidjson::Value &higher, CoolPropFluid &fluid)
Parse the transport properties.
Definition: FluidLibrary.h:528
Definition: GeneralizedCubic.h:569
Use TransportRoutines::viscosity_toluene_higher_order_hardcoded.
Definition: CoolPropFluid.h:260
std::string BibTeX_viscosity
The BibTeX key for the viscosity model.
Definition: CoolPropFluid.h:316
MeltingLinePiecewisePolynomialInThetaData polynomial_in_Theta
The data needed for a melting curve formed of segments that are polynomials in .
Definition: Ancillaries.h:223
CoolPropFluid get_fluid(const std::string &fluid_string)
Get the fluid structure.
Definition: FluidLibrary.cpp:321
Use TransportRoutines::viscosity_dilute_collision_integral.
Definition: CoolPropFluid.h:202
void set_limits()
Evaluate the melting line to calculate the limits of the curve (Tmin/Tmax and pmin/pmax) ...
Definition: Ancillaries.cpp:119
bool viscosity_model_provided
A flag for whether viscosity model is provided. False for no.
Definition: CoolPropFluid.h:322
SimpleState crit
The state at the critical point.
Definition: CoolPropFluid.h:503
ViscosityDiluteGasPowersOfTr powers_of_Tr
Data for TransportRoutines::viscosity_dilute_powers_of_Tr.
Definition: CoolPropFluid.h:215
The term in the EOS used to shift the reference state of the fluid.
Definition: Helmholtz.h:781
Use TransportRoutines::viscosity_R23_hardcoded.
Definition: CoolPropFluid.h:289
void validate(CoolPropFluid &fluid)
Validate the fluid file that was just constructed.
Definition: FluidLibrary.h:1166
bool viscosity_using_ECS
A flag for whether to use extended corresponding states for viscosity. False for no.
Definition: CoolPropFluid.h:318
EnvironmentalFactorsStruct environment
The environmental variables for global warming potential, ODP, etc.
Definition: CoolPropFluid.h:500
Use TransportRoutines::viscosity_helium_hardcoded.
Definition: CoolPropFluid.h:288
void parse_thermal_conductivity(rapidjson::Value &conductivity, CoolPropFluid &fluid)
Parse the thermal conductivity data.
Definition: FluidLibrary.h:915
void parse_initial_density_viscosity(rapidjson::Value &initial_density, CoolPropFluid &fluid)
Parse the transport properties.
Definition: FluidLibrary.h:494
CoolPropDbl T_reducing
Reducing temperature [K[.
Definition: CoolPropFluid.h:187
Use TransportRoutines::viscosity_initial_density_dependence_Rainwater_Friend.
Definition: CoolPropFluid.h:231
ViscosityDiluteGasPowersOfT powers_of_T
Data for TransportRoutines::viscosity_dilute_powers_of_T.
Definition: CoolPropFluid.h:214
IdealHelmholtzContainer alpha0
The ideal Helmholtz energy.
Definition: CoolPropFluid.h:372
Use TransportRoutines::conductivity_hardcoded_water.
Definition: CoolPropFluid.h:297
Use TransportRoutines::viscosity_heavywater_hardcoded.
Definition: CoolPropFluid.h:287
Use TransportRoutines::conductivity_hardcoded_R23.
Definition: CoolPropFluid.h:299
This file contains flash routines in which the state is unknown, and a solver of some kind must be us...
Definition: AbstractState.h:19
CoolPropDbl C
Leading constant.
Definition: CoolPropFluid.h:187
Use TransportRoutines::viscosity_dilute_collision_integral_powers_of_T.
Definition: CoolPropFluid.h:203
bool conductivity_using_ECS
A flag for whether to use extended corresponding states for conductivity. False for no...
Definition: CoolPropFluid.h:319
bool viscosity_using_Chung
A flag for whether to use Chung model. False for no.
Definition: CoolPropFluid.h:320
ViscosityHardcodedEnum hardcoded_viscosity
Hardcoded flags for the viscosity.
Definition: CoolPropFluid.h:326
SimpleState max_sat_p
The state at the maximum saturation pressure for pseudo-pure.
Definition: CoolPropFluid.h:357
A container for the fluid parameters for the CoolProp fluids.
Definition: FluidLibrary.h:27
void parse_residual_conductivity(rapidjson::Value &dilute, CoolPropFluid &fluid)
Parse the transport properties.
Definition: FluidLibrary.h:828
ResidualHelmholtzContainer alphar
The residual Helmholtz energy.
Definition: CoolPropFluid.h:371
a polynomial in is in use
Definition: Ancillaries.h:211
std::string alpha_type
The type of alpha function.
Definition: CubicsLibrary.h:24
Definition: Helmholtz.h:955
void add_Exponential(const std::vector< CoolPropDbl > &n, const std::vector< CoolPropDbl > &d, const std::vector< CoolPropDbl > &t, const std::vector< CoolPropDbl > &g, const std::vector< CoolPropDbl > &l)
Add and convert an old-style exponential term to generalized form.
Definition: Helmholtz.h:283
Use TransportRoutines::viscosity_dilute_cyclohexane.
Definition: CoolPropFluid.h:206
IdealHelmholtzContainer alpha0
The ideal Helmholtz energy.
Definition: CubicsLibrary.h:26
void validate()
Validate the EOS that was just constructed.
Definition: CoolPropFluid.h:378
Use TransportRoutines::viscosity_methanol_hardcoded.
Definition: CoolPropFluid.h:290
ViscosityDiluteGasCollisionIntegralData collision_integral
Data for TransportRoutines::viscosity_dilute_collision_integral.
Definition: CoolPropFluid.h:212