Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic

Base functions as described in IAWPS/IF97

Information

Package description

Package BaseIF97/Basic computes the the fundamental functions for the 5 regions of the steam tables as described in the standards document IF97.pdf. The code of these functions has been generated using Mathematica and the add-on packages "Format" and "Optimize" to generate highly efficient, expression-optimized C-code from a symbolic representation of the thermodynamic functions. The C-code has than been transformed into Modelica code. An important feature of this optimization was to simultaneously optimize the functions and the directional derivatives because they share many common subexpressions.

Package contents

Version Info and Revision history

Author: Hubertus Tummescheit,
Modelon AB
Ideon Science Park
SE-22370 Lund, Sweden
email: hubertus@modelon.se
          

Extends from Modelica.Icons.Package (Icon for standard packages).

Package Content

NameDescription
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g1 g1 Gibbs function for region 1: g(p,T)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g2 g2 Gibbs function for region 2: g(p,T)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g2metastable g2metastable Gibbs function for metastable part of region 2: g(p,T)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.f3 f3 Helmholtz function for region 3: f(d,T)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g5 g5 base function for region 5: g(p,T)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.gibbs gibbs Gibbs function for region 1, 2 or 5: g(p,T,region)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g1pitau g1pitau derivative of g w.r.t. pi and tau
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g2pitau g2pitau derivative of g w.r.t. pi and tau
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g5pitau g5pitau derivative of g w.r.t. pi and tau
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.f3deltatau f3deltatau 1st derivatives of f w.r.t. delta and tau
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tph1 tph1 inverse function for region 1: T(p,h)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps1 tps1 inverse function for region 1: T(p,s)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tph2 tph2 reverse function for region 2: T(p,h)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2a tps2a reverse function for region 2a: T(p,s)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2b tps2b reverse function for region 2b: T(p,s)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2c tps2c reverse function for region 2c: T(p,s)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2 tps2 reverse function for region 2: T(p,s)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tsat tsat region 4 saturation temperature as a function of pressure
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.dtsatofp dtsatofp derivative of saturation temperature w.r.t. pressure
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tsat_der tsat_der derivative function for tsat
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.psat psat region 4 saturation pressure as a functionx of temperature
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.dptofT dptofT derivative of pressure w.r.t. temperature along the saturation pressure curve
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.psat_der psat_der derivative function for psat
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p1_hs p1_hs pressure as a function of ehtnalpy and entropy in region 1
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.h2ab_s h2ab_s boundary between regions 2a and 2b
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p2a_hs p2a_hs pressure as a function of enthalpy and entropy in subregion 2a
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p2b_hs p2b_hs pressure as a function of enthalpy and entropy in subregion 2a
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p2c_hs p2c_hs pressure as a function of enthalpy and entropy in subregion 2c
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.h3ab_p h3ab_p ergion 3 a b boundary for pressure/enthalpy
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3a_ph T3a_ph Region 3 a: inverse function T(p,h)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3b_ph T3b_ph Region 3 b: inverse function T(p,h)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3a_ph v3a_ph Region 3 a: inverse function v(p,h)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3b_ph v3b_ph Region 3 b: inverse function v(p,h)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3a_ps T3a_ps Region 3 a: inverse function T(p,s)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3b_ps T3b_ps Region 3 b: inverse function T(p,s)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3a_ps v3a_ps Region 3 a: inverse function v(p,s)
Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3b_ps v3b_ps Region 3 b: inverse function v(p,s)


Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g1 Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g1

Gibbs function for region 1: g(p,T)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
GibbsDerivsgdimensionless Gibbs funcion and dervatives w.r.t. pi and tau

Modelica definition

function g1 "Gibbs function for region 1: g(p,T)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature (K)";
  output Modelica.Media.Common.GibbsDerivs g 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
protected 
  Real pi1 "dimensionless pressure";
  Real tau1 "dimensionless temperature";
  Real[45] o "vector of auxiliary variables";
  Real pl "auxiliary variable";
algorithm 
  pl := min(p, data.PCRIT - 1);
  assert(p > triple.ptriple,
    "IF97 medium function g1 called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");
  assert(p <= 100.0e6,
    "IF97 medium function g1: the input pressure (= " + String(p) + " Pa) is higher than 100 Mpa");
  assert(T >= 273.15,
    "IF97 medium function g1: the temperature (= " + String(T) + " K)  is lower than 273.15 K!");
  g.p := p;
  g.T := T;
  g.R := data.RH2O;
  g.pi := p/data.PSTAR1;
  g.tau := data.TSTAR1/T;
  pi1 := 7.1000000000000 - g.pi;
  tau1 := -1.22200000000000 + g.tau;
  o[1] := tau1*tau1;
  o[2] := o[1]*o[1];
  o[3] := o[2]*o[2];
  o[4] := o[3]*tau1;
  o[5] := 1/o[4];
  o[6] := o[1]*o[2];
  o[7] := o[1]*tau1;
  o[8] := 1/o[7];
  o[9] := o[1]*o[2]*o[3];
  o[10] := 1/o[2];
  o[11] := o[2]*tau1;
  o[12] := 1/o[11];
  o[13] := o[2]*o[3];
  o[14] := 1/o[3];
  o[15] := pi1*pi1;
  o[16] := o[15]*pi1;
  o[17] := o[15]*o[15];
  o[18] := o[17]*o[17];
  o[19] := o[17]*o[18]*pi1;
  o[20] := o[15]*o[17];
  o[21] := o[3]*o[3];
  o[22] := o[21]*o[21];
  o[23] := o[22]*o[3]*tau1;
  o[24] := 1/o[23];
  o[25] := o[22]*o[3];
  o[26] := 1/o[25];
  o[27] := o[1]*o[2]*o[22]*tau1;
  o[28] := 1/o[27];
  o[29] := o[1]*o[2]*o[22];
  o[30] := 1/o[29];
  o[31] := o[1]*o[2]*o[21]*o[3]*tau1;
  o[32] := 1/o[31];
  o[33] := o[2]*o[21]*o[3]*tau1;
  o[34] := 1/o[33];
  o[35] := o[1]*o[3]*tau1;
  o[36] := 1/o[35];
  o[37] := o[1]*o[3];
  o[38] := 1/o[37];
  o[39] := 1/o[6];
  o[40] := o[1]*o[22]*o[3];
  o[41] := 1/o[40];
  o[42] := 1/o[22];
  o[43] := o[1]*o[2]*o[21]*o[3];
  o[44] := 1/o[43];
  o[45] := 1/o[13];
  g.g := pi1*(pi1*(pi1*(o[10]*(-0.000031679644845054 + o[2]*(-2.82707979853120e-6
     - 8.5205128120103e-10*o[6])) + pi1*(o[12]*(-2.24252819080000e-6 + (-6.5171222895601e-7
     - 1.43417299379240e-13*o[13])*o[7]) + pi1*(-4.0516996860117e-7*o[14]
     + o[16]*((-1.27343017416410e-9 - 1.74248712306340e-10*o[11])*o[36]
     + o[19]*(-6.8762131295531e-19*o[34] + o[15]*(1.44783078285210e-20*o[
    32] + o[20]*(2.63357816627950e-23*o[30] + pi1*(-1.19476226400710e-23*
    o[28] + pi1*(1.82280945814040e-24*o[26] - 9.3537087292458e-26*o[24]*
    pi1))))))))) + o[8]*(-0.00047184321073267 + o[7]*(-0.000300017807930260
     + (0.000047661393906987 + o[1]*(-4.4141845330846e-6 -
    7.2694996297594e-16*o[9]))*tau1))) + o[5]*(0.000283190801238040 + o[1]
    *(-0.00060706301565874 + o[6]*(-0.0189900682184190 + tau1*(-0.032529748770505
     + (-0.0218417171754140 - 0.000052838357969930*o[1])*tau1))))) + (
    0.146329712131670 + tau1*(-0.84548187169114 + tau1*(-3.7563603672040
     + tau1*(3.3855169168385 + tau1*(-0.95791963387872 + tau1*(
    0.157720385132280 + (-0.0166164171995010 + 0.00081214629983568*tau1)*
    tau1))))))/o[1];

  g.gpi := pi1*(pi1*(o[10]*(0.000095038934535162 + o[2]*(
    8.4812393955936e-6 + 2.55615384360309e-9*o[6])) + pi1*(o[12]*(
    8.9701127632000e-6 + (2.60684891582404e-6 + 5.7366919751696e-13*o[13])
    *o[7]) + pi1*(2.02584984300585e-6*o[14] + o[16]*((1.01874413933128e-8
     + 1.39398969845072e-9*o[11])*o[36] + o[19]*(1.44400475720615e-17*o[
    34] + o[15]*(-3.3300108005598e-19*o[32] + o[20]*(-7.6373766822106e-22
    *o[30] + pi1*(3.5842867920213e-22*o[28] + pi1*(-5.6507093202352e-23*o[
    26] + 2.99318679335866e-24*o[24]*pi1))))))))) + o[8]*(
    0.00094368642146534 + o[7]*(0.00060003561586052 + (-0.000095322787813974
     + o[1]*(8.8283690661692e-6 + 1.45389992595188e-15*o[9]))*tau1))) + o[
    5]*(-0.000283190801238040 + o[1]*(0.00060706301565874 + o[6]*(
    0.0189900682184190 + tau1*(0.032529748770505 + (0.0218417171754140 +
    0.000052838357969930*o[1])*tau1))));

  g.gpipi := pi1*(o[10]*(-0.000190077869070324 + o[2]*(-0.0000169624787911872
     - 5.1123076872062e-9*o[6])) + pi1*(o[12]*(-0.0000269103382896000 + (
    -7.8205467474721e-6 - 1.72100759255088e-12*o[13])*o[7]) + pi1*(-8.1033993720234e-6
    *o[14] + o[16]*((-7.1312089753190e-8 - 9.7579278891550e-9*o[11])*o[36]
     + o[19]*(-2.88800951441230e-16*o[34] + o[15]*(7.3260237612316e-18*o[
    32] + o[20]*(2.13846547101895e-20*o[30] + pi1*(-1.03944316968618e-20*
    o[28] + pi1*(1.69521279607057e-21*o[26] - 9.2788790594118e-23*o[24]*
    pi1))))))))) + o[8]*(-0.00094368642146534 + o[7]*(-0.00060003561586052
     + (0.000095322787813974 + o[1]*(-8.8283690661692e-6 -
    1.45389992595188e-15*o[9]))*tau1));

  g.gtau := pi1*(o[38]*(-0.00254871721114236 + o[1]*(0.0042494411096112
     + (0.0189900682184190 + (-0.0218417171754140 - 0.000158515073909790*
    o[1])*o[1])*o[6])) + pi1*(o[10]*(0.00141552963219801 + o[2]*(
    0.000047661393906987 + o[1]*(-0.0000132425535992538 -
    1.23581493705910e-14*o[9]))) + pi1*(o[12]*(0.000126718579380216 -
    5.1123076872062e-9*o[37]) + pi1*(o[39]*(0.0000112126409540000 + (
    1.30342445791202e-6 - 1.43417299379240e-12*o[13])*o[7]) + pi1*(
    3.2413597488094e-6*o[5] + o[16]*((1.40077319158051e-8 +
    1.04549227383804e-9*o[11])*o[45] + o[19]*(1.99410180757040e-17*o[44]
     + o[15]*(-4.4882754268415e-19*o[42] + o[20]*(-1.00075970318621e-21*o[
    28] + pi1*(4.6595728296277e-22*o[26] + pi1*(-7.2912378325616e-23*o[24]
     + 3.8350205789908e-24*o[41]*pi1))))))))))) + o[8]*(-0.292659424263340
     + tau1*(0.84548187169114 + o[1]*(3.3855169168385 + tau1*(-1.91583926775744
     + tau1*(0.47316115539684 + (-0.066465668798004 + 0.0040607314991784*
    tau1)*tau1)))));

  g.gtautau := pi1*(o[36]*(0.0254871721114236 + o[1]*(-0.033995528876889
     + (-0.037980136436838 - 0.00031703014781958*o[2])*o[6])) + pi1*(o[12]
    *(-0.0056621185287920 + o[6]*(-0.0000264851071985076 -
    1.97730389929456e-13*o[9])) + pi1*((-0.00063359289690108 -
    2.55615384360309e-8*o[37])*o[39] + pi1*(pi1*(-0.0000291722377392842*o[
    38] + o[16]*(o[19]*(-5.9823054227112e-16*o[32] + o[15]*(o[20]*(
    3.9029628424262e-20*o[26] + pi1*(-1.86382913185108e-20*o[24] + pi1*(
    2.98940751135026e-21*o[41] - (1.61070864317613e-22*pi1)/(o[1]*o[22]*o[
    3]*tau1)))) + 1.43624813658928e-17/(o[22]*tau1))) + (-1.68092782989661e-7
     - 7.3184459168663e-9*o[11])/(o[2]*o[3]*tau1))) + (-0.000067275845724000
     + (-3.9102733737361e-6 - 1.29075569441316e-11*o[13])*o[7])/(o[1]*o[2]
    *tau1))))) + o[10]*(0.87797827279002 + tau1*(-1.69096374338228 + o[7]
    *(-1.91583926775744 + tau1*(0.94632231079368 + (-0.199397006394012 +
    0.0162429259967136*tau1)*tau1))));

  g.gtaupi := o[38]*(0.00254871721114236 + o[1]*(-0.0042494411096112 + (-0.0189900682184190
     + (0.0218417171754140 + 0.000158515073909790*o[1])*o[1])*o[6])) +
    pi1*(o[10]*(-0.00283105926439602 + o[2]*(-0.000095322787813974 + o[1]
    *(0.0000264851071985076 + 2.47162987411820e-14*o[9]))) + pi1*(o[12]*(
    -0.00038015573814065 + 1.53369230616185e-8*o[37]) + pi1*(o[39]*(-0.000044850563816000
     + (-5.2136978316481e-6 + 5.7366919751696e-12*o[13])*o[7]) + pi1*(-0.0000162067987440468
    *o[5] + o[16]*((-1.12061855326441e-7 - 8.3639381907043e-9*o[11])*o[45]
     + o[19]*(-4.1876137958978e-16*o[44] + o[15]*(1.03230334817355e-17*o[
    42] + o[20]*(2.90220313924001e-20*o[28] + pi1*(-1.39787184888831e-20*
    o[26] + pi1*(2.26028372809410e-21*o[24] - 1.22720658527705e-22*o[41]*
    pi1))))))))));
end g1;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g2 Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g2

Gibbs function for region 2: g(p,T)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
GibbsDerivsgdimensionless Gibbs funcion and dervatives w.r.t. pi and tau

Modelica definition

function g2 "Gibbs function for region 2: g(p,T)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature (K)";
  output Modelica.Media.Common.GibbsDerivs g 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
protected 
  Real tau2 "dimensionless temperature";
  Real[55] o "vector of auxiliary variables";
algorithm 
  g.p := p;
  g.T := T;
  g.R := data.RH2O;
  assert(p > triple.ptriple,
    "IF97 medium function g2 called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");
  assert(p <= 100.0e6,
    "IF97 medium function g2: the input pressure (= " + String(p) + " Pa) is higher than 100 Mpa");
  assert(T >= 273.15,
    "IF97 medium function g2: the temperature (= " + String(T) + " K) is lower than 273.15 K!");
  assert(T <= 1073.15,
    "IF97 medium function g2: the input temperature (= " + String(T) + " K) is higher than the limit of 1073.15 K");
  g.pi := p/data.PSTAR2;
  g.tau := data.TSTAR2/T;
  tau2 := -0.5 + g.tau;
  o[1] := tau2*tau2;
  o[2] := o[1]*tau2;
  o[3] := -0.050325278727930*o[2];
  o[4] := -0.057581259083432 + o[3];
  o[5] := o[4]*tau2;
  o[6] := -0.045996013696365 + o[5];
  o[7] := o[6]*tau2;
  o[8] := -0.0178348622923580 + o[7];
  o[9] := o[8]*tau2;
  o[10] := o[1]*o[1];
  o[11] := o[10]*o[10];
  o[12] := o[11]*o[11];
  o[13] := o[10]*o[11]*o[12]*tau2;
  o[14] := o[1]*o[10]*tau2;
  o[15] := o[10]*o[11]*tau2;
  o[16] := o[1]*o[12]*tau2;
  o[17] := o[1]*o[11]*tau2;
  o[18] := o[1]*o[10]*o[11];
  o[19] := o[10]*o[11]*o[12];
  o[20] := o[1]*o[10];
  o[21] := g.pi*g.pi;
  o[22] := o[21]*o[21];
  o[23] := o[21]*o[22];
  o[24] := o[10]*o[12]*tau2;
  o[25] := o[12]*o[12];
  o[26] := o[11]*o[12]*o[25]*tau2;
  o[27] := o[10]*o[12];
  o[28] := o[1]*o[10]*o[11]*tau2;
  o[29] := o[10]*o[12]*o[25]*tau2;
  o[30] := o[1]*o[10]*o[25]*tau2;
  o[31] := o[1]*o[11]*o[12];
  o[32] := o[1]*o[12];
  o[33] := g.tau*g.tau;
  o[34] := o[33]*o[33];
  o[35] := -0.000053349095828174*o[13];
  o[36] := -0.087594591301146 + o[35];
  o[37] := o[2]*o[36];
  o[38] := -0.0078785554486710 + o[37];
  o[39] := o[1]*o[38];
  o[40] := -0.00037897975032630 + o[39];
  o[41] := o[40]*tau2;
  o[42] := -0.000066065283340406 + o[41];
  o[43] := o[42]*tau2;
  o[44] := 5.7870447262208e-6*tau2;
  o[45] := -0.301951672367580*o[2];
  o[46] := -0.172743777250296 + o[45];
  o[47] := o[46]*tau2;
  o[48] := -0.091992027392730 + o[47];
  o[49] := o[48]*tau2;
  o[50] := o[1]*o[11];
  o[51] := o[10]*o[11];
  o[52] := o[11]*o[12]*o[25];
  o[53] := o[10]*o[12]*o[25];
  o[54] := o[1]*o[10]*o[25];
  o[55] := o[11]*o[12]*tau2;

  g.g := g.pi*(-0.00177317424732130 + o[9] + g.pi*(tau2*(-0.000033032641670203
     + (-0.000189489875163150 + o[1]*(-0.0039392777243355 + (-0.043797295650573
     - 0.0000266745479140870*o[13])*o[2]))*tau2) + g.pi*(
    2.04817376923090e-8 + (4.3870667284435e-7 + o[1]*(-0.000032277677238570
     + (-0.00150339245421480 - 0.040668253562649*o[13])*o[2]))*tau2 + g.
    pi*(g.pi*(2.29220763376610e-6*o[14] + g.pi*((-1.67147664510610e-11 +
    o[15]*(-0.00211714723213550 - 23.8957419341040*o[16]))*o[2] + g.pi*(-5.9059564324270e-18
     + o[17]*(-1.26218088991010e-6 - 0.038946842435739*o[18]) + g.pi*(o[
    11]*(1.12562113604590e-11 - 8.2311340897998*o[19]) + g.pi*(
    1.98097128020880e-8*o[15] + g.pi*(o[10]*(1.04069652101740e-19 + (-1.02347470959290e-13
     - 1.00181793795110e-9*o[10])*o[20]) + o[23]*(o[13]*(-8.0882908646985e-11
     + 0.106930318794090*o[24]) + o[21]*(-0.33662250574171*o[26] + o[21]*
    (o[27]*(8.9185845355421e-25 + (3.06293168762320e-13 -
    4.2002467698208e-6*o[15])*o[28]) + g.pi*(-5.9056029685639e-26*o[24]
     + g.pi*(3.7826947613457e-6*o[29] + g.pi*(-1.27686089346810e-15*o[30]
     + o[31]*(7.3087610595061e-29 + o[18]*(5.5414715350778e-17 -
    9.4369707241210e-7*o[32]))*g.pi)))))))))))) + tau2*(-7.8847309559367e-10
     + (1.27907178522850e-8 + 4.8225372718507e-7*tau2)*tau2))))) + (-0.0056087911830200
     + g.tau*(0.071452738814550 + g.tau*(-0.40710498239280 + g.tau*(
    1.42408197144400 + g.tau*(-4.3839511194500 + g.tau*(-9.6927686002170
     + g.tau*(10.0866556801800 + (-0.284086326077200 + 0.0212684635330700
    *g.tau)*g.tau) + Modelica.Math.log(g.pi)))))))/(o[34]*g.tau);

  g.gpi := (1.00000000000000 + g.pi*(-0.00177317424732130 + o[9] + g.pi*(
    o[43] + g.pi*(6.1445213076927e-8 + (1.31612001853305e-6 + o[1]*(-0.000096833031715710
     + (-0.0045101773626444 - 0.122004760687947*o[13])*o[2]))*tau2 + g.pi
    *(g.pi*(0.0000114610381688305*o[14] + g.pi*((-1.00288598706366e-10 +
    o[15]*(-0.0127028833928130 - 143.374451604624*o[16]))*o[2] + g.pi*(-4.1341695026989e-17
     + o[17]*(-8.8352662293707e-6 - 0.272627897050173*o[18]) + g.pi*(o[11]
    *(9.0049690883672e-11 - 65.849072718398*o[19]) + g.pi*(
    1.78287415218792e-7*o[15] + g.pi*(o[10]*(1.04069652101740e-18 + (-1.02347470959290e-12
     - 1.00181793795110e-8*o[10])*o[20]) + o[23]*(o[13]*(-1.29412653835176e-9
     + 1.71088510070544*o[24]) + o[21]*(-6.0592051033508*o[26] + o[21]*(o[
    27]*(1.78371690710842e-23 + (6.1258633752464e-12 -
    0.000084004935396416*o[15])*o[28]) + g.pi*(-1.24017662339842e-24*o[24]
     + g.pi*(0.000083219284749605*o[29] + g.pi*(-2.93678005497663e-14*o[
    30] + o[31]*(1.75410265428146e-27 + o[18]*(1.32995316841867e-15 -
    0.0000226487297378904*o[32]))*g.pi)))))))))))) + tau2*(-3.15389238237468e-9
     + (5.1162871409140e-8 + 1.92901490874028e-6*tau2)*tau2))))))/g.pi;

  g.gpipi := (-1.00000000000000 + o[21]*(o[43] + g.pi*(
    1.22890426153854e-7 + (2.63224003706610e-6 + o[1]*(-0.000193666063431420
     + (-0.0090203547252888 - 0.244009521375894*o[13])*o[2]))*tau2 + g.pi
    *(g.pi*(0.000045844152675322*o[14] + g.pi*((-5.0144299353183e-10 + o[
    15]*(-0.063514416964065 - 716.87225802312*o[16]))*o[2] + g.pi*(-2.48050170161934e-16
     + o[17]*(-0.000053011597376224 - 1.63576738230104*o[18]) + g.pi*(o[
    11]*(6.3034783618570e-10 - 460.94350902879*o[19]) + g.pi*(
    1.42629932175034e-6*o[15] + g.pi*(o[10]*(9.3662686891566e-18 + (-9.2112723863361e-12
     - 9.0163614415599e-8*o[10])*o[20]) + o[23]*(o[13]*(-1.94118980752764e-8
     + 25.6632765105816*o[24]) + o[21]*(-103.006486756963*o[26] + o[21]*(
    o[27]*(3.3890621235060e-22 + (1.16391404129682e-10 -
    0.00159609377253190*o[15])*o[28]) + g.pi*(-2.48035324679684e-23*o[24]
     + g.pi*(0.00174760497974171*o[29] + g.pi*(-6.4609161209486e-13*o[30]
     + o[31]*(4.0344361048474e-26 + o[18]*(3.05889228736295e-14 -
    0.00052092078397148*o[32]))*g.pi)))))))))))) + tau2*(-9.4616771471240e-9
     + (1.53488614227420e-7 + o[44])*tau2)))))/o[21];

  g.gtau := (0.0280439559151000 + g.tau*(-0.285810955258200 + g.tau*(
    1.22131494717840 + g.tau*(-2.84816394288800 + g.tau*(4.3839511194500
     + o[33]*(10.0866556801800 + (-0.56817265215440 + 0.063805390599210*g.
     tau)*g.tau))))))/(o[33]*o[34]) + g.pi*(-0.0178348622923580 + o[49]
     + g.pi*(-0.000033032641670203 + (-0.00037897975032630 + o[1]*(-0.0157571108973420
     + (-0.306581069554011 - 0.00096028372490713*o[13])*o[2]))*tau2 + g.
    pi*(4.3870667284435e-7 + o[1]*(-0.000096833031715710 + (-0.0090203547252888
     - 1.42338887469272*o[13])*o[2]) + g.pi*(-7.8847309559367e-10 + g.pi*
    (0.0000160454534363627*o[20] + g.pi*(o[1]*(-5.0144299353183e-11 + o[
    15]*(-0.033874355714168 - 836.35096769364*o[16])) + g.pi*((-0.0000138839897890111
     - 0.97367106089347*o[18])*o[50] + g.pi*(o[14]*(9.0049690883672e-11
     - 296.320827232793*o[19]) + g.pi*(2.57526266427144e-7*o[51] + g.pi*(
    o[2]*(4.1627860840696e-19 + (-1.02347470959290e-12 -
    1.40254511313154e-8*o[10])*o[20]) + o[23]*(o[19]*(-2.34560435076256e-9
     + 5.3465159397045*o[24]) + o[21]*(-19.1874828272775*o[52] + o[21]*(o[
    16]*(1.78371690710842e-23 + (1.07202609066812e-11 -
    0.000201611844951398*o[15])*o[28]) + g.pi*(-1.24017662339842e-24*o[27]
     + g.pi*(0.000200482822351322*o[53] + g.pi*(-4.9797574845256e-14*o[54]
     + (1.90027787547159e-27 + o[18]*(2.21658861403112e-15 -
    0.000054734430199902*o[32]))*o[55]*g.pi)))))))))))) + (
    2.55814357045700e-8 + 1.44676118155521e-6*tau2)*tau2))));

  g.gtautau := (-0.168263735490600 + g.tau*(1.42905477629100 + g.tau*(-4.8852597887136
     + g.tau*(8.5444918286640 + g.tau*(-8.7679022389000 + o[33]*(-0.56817265215440
     + 0.127610781198420*g.tau)*g.tau)))))/(o[33]*o[34]*g.tau) + g.pi*(-0.091992027392730
     + (-0.34548755450059 - 1.50975836183790*o[2])*tau2 + g.pi*(-0.00037897975032630
     + o[1]*(-0.047271332692026 + (-1.83948641732407 - 0.033609930371750*
    o[13])*o[2]) + g.pi*((-0.000193666063431420 + (-0.045101773626444 -
    48.395221739552*o[13])*o[2])*tau2 + g.pi*(2.55814357045700e-8 +
    2.89352236311042e-6*tau2 + g.pi*(0.000096272720618176*o[10]*tau2 + g.
    pi*((-1.00288598706366e-10 + o[15]*(-0.50811533571252 -
    28435.9329015838*o[16]))*tau2 + g.pi*(o[11]*(-0.000138839897890111 -
    23.3681054614434*o[18])*tau2 + g.pi*((6.3034783618570e-10 -
    10371.2289531477*o[19])*o[20] + g.pi*(3.09031519712573e-6*o[17] + g.
    pi*(o[1]*(1.24883582522088e-18 + (-9.2112723863361e-12 -
    1.82330864707100e-7*o[10])*o[20]) + o[23]*(o[1]*o[11]*o[12]*(-6.5676921821352e-8
     + 261.979281045521*o[24])*tau2 + o[21]*(-1074.49903832754*o[1]*o[10]
    *o[12]*o[25]*tau2 + o[21]*((3.3890621235060e-22 + (
    3.6448887082716e-10 - 0.0094757567127157*o[15])*o[28])*o[32] + g.pi*(
    -2.48035324679684e-23*o[16] + g.pi*(0.0104251067622687*o[1]*o[12]*o[
    25]*tau2 + g.pi*(o[11]*o[12]*(4.7506946886790e-26 + o[18]*(
    8.6446955947214e-14 - 0.00311986252139440*o[32]))*g.pi -
    1.89230784411972e-12*o[10]*o[25]*tau2))))))))))))))));

  g.gtaupi := -0.0178348622923580 + o[49] + g.pi*(-0.000066065283340406
     + (-0.00075795950065260 + o[1]*(-0.0315142217946840 + (-0.61316213910802
     - 0.00192056744981426*o[13])*o[2]))*tau2 + g.pi*(1.31612001853305e-6
     + o[1]*(-0.000290499095147130 + (-0.0270610641758664 -
    4.2701666240781*o[13])*o[2]) + g.pi*(-3.15389238237468e-9 + g.pi*(
    0.000080227267181813*o[20] + g.pi*(o[1]*(-3.00865796119098e-10 + o[15]
    *(-0.203246134285008 - 5018.1058061618*o[16])) + g.pi*((-0.000097187928523078
     - 6.8156974262543*o[18])*o[50] + g.pi*(o[14]*(7.2039752706938e-10 -
    2370.56661786234*o[19]) + g.pi*(2.31773639784430e-6*o[51] + g.pi*(o[2]
    *(4.1627860840696e-18 + (-1.02347470959290e-11 - 1.40254511313154e-7*
    o[10])*o[20]) + o[23]*(o[19]*(-3.7529669612201e-8 + 85.544255035272*o[
    24]) + o[21]*(-345.37469089099*o[52] + o[21]*(o[16]*(
    3.5674338142168e-22 + (2.14405218133624e-10 - 0.0040322368990280*o[15])
    *o[28]) + g.pi*(-2.60437090913668e-23*o[27] + g.pi*(
    0.0044106220917291*o[53] + g.pi*(-1.14534422144089e-12*o[54] + (
    4.5606669011318e-26 + o[18]*(5.3198126736747e-14 -
    0.00131362632479764*o[32]))*o[55]*g.pi)))))))))))) + (
    1.02325742818280e-7 + o[44])*tau2)));
end g2;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g2metastable Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g2metastable

Gibbs function for metastable part of region 2: g(p,T)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
GibbsDerivsgdimensionless Gibbs funcion and dervatives w.r.t. pi and tau

Modelica definition

function g2metastable 
  "Gibbs function for metastable part of region 2: g(p,T)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature (K)";
  output Modelica.Media.Common.GibbsDerivs g 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
protected 
  Real pi "dimensionless pressure";
  Real tau "dimensionless temperature";
  Real tau2 "dimensionless temperature";
  Real[27] o "vector of auxiliary variables";
algorithm 
  assert(p > triple.ptriple,
    "IF97 medium function g2metastable called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");
  assert(p <= 100.0e6,
    "IF97 medium function g2metastable: the input pressure (= " + String(p) + " Pa) is higher than 100 Mpa");
  assert(T >= 273.15,
    "IF97 medium function g2metastable: the temperature (= " + String(T) + " K) is lower than 273.15 K!");
  assert(T <= 1073.15,
    "IF97 medium function g2metastable: the input temperature (= " + String(T) + " K) is higher than the limit of 1073.15 K");
  g.p := p;
  g.T := T;
  g.R := data.RH2O;
  g.pi := p/data.PSTAR2;
  g.tau := data.TSTAR2/T;
  tau2 := -0.5 + g.tau;
  o[1] := tau2*tau2;
  o[2] := o[1]*tau2;
  o[3] := o[1]*o[1];
  o[4] := o[1]*o[3];
  o[5] := -0.0040813178534455*o[4];
  o[6] := -0.072334555213245 + o[5];
  o[7] := o[2]*o[6];
  o[8] := -0.088223831943146 + o[7];
  o[9] := o[1]*o[8];
  o[10] := o[3]*o[3];
  o[11] := o[10]*tau2;
  o[12] := o[10]*o[3];
  o[13] := o[1]*o[3]*tau2;
  o[14] := g.tau*g.tau;
  o[15] := o[14]*o[14];
  o[16] := -0.015238081817394*o[11];
  o[17] := -0.106091843797284 + o[16];
  o[18] := o[17]*o[4];
  o[19] := 0.0040195606760414 + o[18];
  o[20] := o[19]*tau2;
  o[21] := g.pi*g.pi;
  o[22] := -0.0448944963879005*o[4];
  o[23] := -0.361672776066225 + o[22];
  o[24] := o[2]*o[23];
  o[25] := -0.176447663886292 + o[24];
  o[26] := o[25]*tau2;
  o[27] := o[3]*tau2;

  g.g := g.pi*(-0.0073362260186506 + o[9] + g.pi*(g.pi*((-0.0063498037657313
     - 0.086043093028588*o[12])*o[3] + g.pi*(o[13]*(0.007532158152277 -
    0.0079238375446139*o[2]) + o[11]*g.pi*(-0.00022888160778447 -
    0.002645650148281*tau2))) + (0.0020097803380207 + (-0.053045921898642
     - 0.007619040908697*o[11])*o[4])*tau2)) + (-0.00560879118302 + g.tau
    *(0.07145273881455 + g.tau*(-0.4071049823928 + g.tau*(1.424081971444
     + g.tau*(-4.38395111945 + g.tau*(-9.6937268393049 + g.tau*(
    10.087275970006 + (-0.2840863260772 + 0.02126846353307*g.tau)*g.tau)
     + Modelica.Math.log(g.pi)))))))/(o[15]*g.tau);

  g.gpi := (1.0 + g.pi*(-0.0073362260186506 + o[9] + g.pi*(o[20] + g.pi*(
    (-0.0190494112971939 - 0.258129279085764*o[12])*o[3] + g.pi*(o[13]*(
    0.030128632609108 - 0.0316953501784556*o[2]) + o[11]*g.pi*(-0.00114440803892235
     - 0.013228250741405*tau2))))))/g.pi;

  g.gpipi := (-1. + o[21]*(o[20] + g.pi*((-0.0380988225943878 -
    0.516258558171528*o[12])*o[3] + g.pi*(o[13]*(0.090385897827324 -
    0.0950860505353668*o[2]) + o[11]*g.pi*(-0.0045776321556894 -
    0.05291300296562*tau2)))))/o[21];

  g.gtau := (0.0280439559151 + g.tau*(-0.2858109552582 + g.tau*(
    1.2213149471784 + g.tau*(-2.848163942888 + g.tau*(4.38395111945 + o[
    14]*(10.087275970006 + (-0.5681726521544 + 0.06380539059921*g.tau)*g.
    tau))))))/(o[14]*o[15]) + g.pi*(o[26] + g.pi*(0.0020097803380207 + (-0.371321453290494
     - 0.121904654539152*o[11])*o[4] + g.pi*((-0.0253992150629252 -
    1.37668948845741*o[12])*o[2] + g.pi*((0.052725107065939 -
    0.079238375446139*o[2])*o[4] + o[10]*g.pi*(-0.00205993447006023 -
    0.02645650148281*tau2)))));

  g.gtautau := (-0.1682637354906 + g.tau*(1.429054776291 + g.tau*(-4.8852597887136
     + g.tau*(8.544491828664 + g.tau*(-8.7679022389 + o[14]*(-0.5681726521544
     + 0.12761078119842*g.tau)*g.tau)))))/(o[14]*o[15]*g.tau) + g.pi*(-0.176447663886292
     + o[2]*(-1.4466911042649 - 0.448944963879005*o[4]) + g.pi*((-2.22792871974296
     - 1.82856981808728*o[11])*o[27] + g.pi*(o[1]*(-0.0761976451887756 -
    20.6503423268611*o[12]) + g.pi*((0.316350642395634 -
    0.713145379015251*o[2])*o[27] + o[13]*g.pi*(-0.0164794757604818 -
    0.23810851334529*tau2)))));

  g.gtaupi := o[26] + g.pi*(0.0040195606760414 + (-0.742642906580988 -
    0.243809309078304*o[11])*o[4] + g.pi*((-0.0761976451887756 -
    4.13006846537222*o[12])*o[2] + g.pi*((0.210900428263756 -
    0.316953501784556*o[2])*o[4] + o[10]*g.pi*(-0.0102996723503012 -
    0.13228250741405*tau2))));
end g2metastable;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.f3 Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.f3

Helmholtz function for region 3: f(d,T)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
HelmholtzDerivsfdimensionless Helmholtz function and dervatives w.r.t. delta and tau

Modelica definition

function f3 "Helmholtz function for region 3: f(d,T)"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature (K)";
  output Modelica.Media.Common.HelmholtzDerivs f 
    "dimensionless Helmholtz function and dervatives w.r.t. delta and tau";
protected 
  Real[40] o "vector of auxiliary variables";
algorithm 
  f.T := T;
  f.d := d;
  f.R := data.RH2O;
  f.tau := data.TCRIT/T;
  f.delta := if (d == data.DCRIT and T == data.TCRIT) then 1 - Modelica.
    Constants.eps else abs(d/data.DCRIT);
  o[1] := f.tau*f.tau;
  o[2] := o[1]*o[1];
  o[3] := o[2]*f.tau;
  o[4] := o[1]*f.tau;
  o[5] := o[2]*o[2];
  o[6] := o[1]*o[5]*f.tau;
  o[7] := o[5]*f.tau;
  o[8] := -0.64207765181607*o[1];
  o[9] := 0.88521043984318 + o[8];
  o[10] := o[7]*o[9];
  o[11] := -1.15244078066810 + o[10];
  o[12] := o[11]*o[2];
  o[13] := -1.26543154777140 + o[12];
  o[14] := o[1]*o[13];
  o[15] := o[1]*o[2]*o[5]*f.tau;
  o[16] := o[2]*o[5];
  o[17] := o[1]*o[5];
  o[18] := o[5]*o[5];
  o[19] := o[1]*o[18]*o[2];
  o[20] := o[1]*o[18]*o[2]*f.tau;
  o[21] := o[18]*o[5];
  o[22] := o[1]*o[18]*o[5];
  o[23] := 0.251168168486160*o[2];
  o[24] := 0.078841073758308 + o[23];
  o[25] := o[15]*o[24];
  o[26] := -6.1005234513930 + o[25];
  o[27] := o[26]*f.tau;
  o[28] := 9.7944563083754 + o[27];
  o[29] := o[2]*o[28];
  o[30] := -1.70429417648412 + o[29];
  o[31] := o[1]*o[30];
  o[32] := f.delta*f.delta;
  o[33] := -10.9153200808732*o[1];
  o[34] := 13.2781565976477 + o[33];
  o[35] := o[34]*o[7];
  o[36] := -6.9146446840086 + o[35];
  o[37] := o[2]*o[36];
  o[38] := -2.53086309554280 + o[37];
  o[39] := o[38]*f.tau;
  o[40] := o[18]*o[5]*f.tau;

  f.f := -15.7328452902390 + f.tau*(20.9443969743070 + (-7.6867707878716
     + o[3]*(2.61859477879540 + o[4]*(-2.80807811486200 + o[1]*(
    1.20533696965170 - 0.0084566812812502*o[6]))))*f.tau) + f.delta*(o[14]
     + f.delta*(0.38493460186671 + o[1]*(-0.85214708824206 + o[2]*(
    4.8972281541877 + (-3.05026172569650 + o[15]*(0.039420536879154 +
    0.125584084243080*o[2]))*f.tau)) + f.delta*(-0.279993296987100 + o[1]
    *(1.38997995694600 + o[1]*(-2.01899150235700 + o[16]*(-0.0082147637173963
     - 0.47596035734923*o[17]))) + f.delta*(0.043984074473500 + o[1]*(-0.44476435428739
     + o[1]*(0.90572070719733 + 0.70522450087967*o[19])) + f.delta*(f.
    delta*(-0.0221754008730960 + o[1]*(0.094260751665092 +
    0.164362784479610*o[21]) + f.delta*(-0.0135033722413480*o[1] + f.
    delta*(-0.0148343453524720*o[22] + f.delta*(o[1]*(0.00057922953628084
     + 0.0032308904703711*o[21]) + f.delta*(0.000080964802996215 -
    0.000044923899061815*f.delta*o[22] - 0.000165576797950370*f.tau)))))
     + (0.107705126263320 + o[1]*(-0.32913623258954 - 0.50871062041158*o[
    20]))*f.tau))))) + 1.06580700285130*Modelica.Math.log(f.delta);

  f.fdelta := (1.06580700285130 + f.delta*(o[14] + f.delta*(
    0.76986920373342 + o[31] + f.delta*(-0.83997989096130 + o[1]*(
    4.1699398708380 + o[1]*(-6.0569745070710 + o[16]*(-0.0246442911521889
     - 1.42788107204769*o[17]))) + f.delta*(0.175936297894000 + o[1]*(-1.77905741714956
     + o[1]*(3.6228828287893 + 2.82089800351868*o[19])) + f.delta*(f.
    delta*(-0.133052405238576 + o[1]*(0.56556450999055 + 0.98617670687766
    *o[21]) + f.delta*(-0.094523605689436*o[1] + f.delta*(-0.118674762819776
    *o[22] + f.delta*(o[1]*(0.0052130658265276 + 0.0290780142333399*o[21])
     + f.delta*(0.00080964802996215 - 0.00049416288967996*f.delta*o[22]
     - 0.00165576797950370*f.tau))))) + (0.53852563131660 + o[1]*(-1.64568116294770
     - 2.54355310205790*o[20]))*f.tau))))))/f.delta;

  f.fdeltadelta := (-1.06580700285130 + o[32]*(0.76986920373342 + o[31]
     + f.delta*(-1.67995978192260 + o[1]*(8.3398797416760 + o[1]*(-12.1139490141420
     + o[16]*(-0.049288582304378 - 2.85576214409538*o[17]))) + f.delta*(
    0.52780889368200 + o[1]*(-5.3371722514487 + o[1]*(10.8686484863680 +
    8.4626940105560*o[19])) + f.delta*(f.delta*(-0.66526202619288 + o[1]*
    (2.82782254995276 + 4.9308835343883*o[21]) + f.delta*(-0.56714163413662
    *o[1] + f.delta*(-0.83072333973843*o[22] + f.delta*(o[1]*(
    0.041704526612220 + 0.232624113866719*o[21]) + f.delta*(
    0.0072868322696594 - 0.0049416288967996*f.delta*o[22] -
    0.0149019118155333*f.tau))))) + (2.15410252526640 + o[1]*(-6.5827246517908
     - 10.1742124082316*o[20]))*f.tau)))))/o[32];

  f.ftau := 20.9443969743070 + (-15.3735415757432 + o[3]*(
    18.3301634515678 + o[4]*(-28.0807811486200 + o[1]*(14.4640436358204
     - 0.194503669468755*o[6]))))*f.tau + f.delta*(o[39] + f.delta*(f.tau
    *(-1.70429417648412 + o[2]*(29.3833689251262 + (-21.3518320798755 + o[
    15]*(0.86725181134139 + 3.2651861903201*o[2]))*f.tau)) + f.delta*((
    2.77995991389200 + o[1]*(-8.0759660094280 + o[16]*(-0.131436219478341
     - 12.3749692910800*o[17])))*f.tau + f.delta*((-0.88952870857478 + o[
    1]*(3.6228828287893 + 18.3358370228714*o[19]))*f.tau + f.delta*(
    0.107705126263320 + o[1]*(-0.98740869776862 - 13.2264761307011*o[20])
     + f.delta*((0.188521503330184 + 4.2734323964699*o[21])*f.tau + f.
    delta*(-0.0270067444826960*f.tau + f.delta*(-0.38569297916427*o[40]
     + f.delta*(f.delta*(-0.000165576797950370 - 0.00116802137560719*f.
    delta*o[40]) + (0.00115845907256168 + 0.084003152229649*o[21])*f.tau)))))))));

  f.ftautau := -15.3735415757432 + o[3]*(109.980980709407 + o[4]*(-252.727030337580
     + o[1]*(159.104479994024 - 4.2790807283126*o[6]))) + f.delta*(-2.53086309554280
     + o[2]*(-34.573223420043 + (185.894192367068 - 174.645121293971*o[1])
    *o[7]) + f.delta*(-1.70429417648412 + o[2]*(146.916844625631 + (-128.110992479253
     + o[15]*(18.2122880381691 + 81.629654758002*o[2]))*f.tau) + f.delta*
    (2.77995991389200 + o[1]*(-24.2278980282840 + o[16]*(-1.97154329217511
     - 309.374232277000*o[17])) + f.delta*(-0.88952870857478 + o[1]*(
    10.8686484863680 + 458.39592557179*o[19]) + f.delta*(f.delta*(
    0.188521503330184 + 106.835809911747*o[21] + f.delta*(-0.0270067444826960
     + f.delta*(-9.6423244791068*o[21] + f.delta*(0.00115845907256168 +
    2.10007880574121*o[21] - 0.0292005343901797*o[21]*o[32])))) + (-1.97481739553724
     - 330.66190326753*o[20])*f.tau)))));

  f.fdeltatau := o[39] + f.delta*(f.tau*(-3.4085883529682 + o[2]*(
    58.766737850252 + (-42.703664159751 + o[15]*(1.73450362268278 +
    6.5303723806402*o[2]))*f.tau)) + f.delta*((8.3398797416760 + o[1]*(-24.2278980282840
     + o[16]*(-0.39430865843502 - 37.124907873240*o[17])))*f.tau + f.
    delta*((-3.5581148342991 + o[1]*(14.4915313151573 + 73.343348091486*o[
    19]))*f.tau + f.delta*(0.53852563131660 + o[1]*(-4.9370434888431 -
    66.132380653505*o[20]) + f.delta*((1.13112901998110 +
    25.6405943788192*o[21])*f.tau + f.delta*(-0.189047211378872*f.tau + f.
     delta*(-3.08554383331418*o[40] + f.delta*(f.delta*(-0.00165576797950370
     - 0.0128482351316791*f.delta*o[40]) + (0.0104261316530551 +
    0.75602837006684*o[21])*f.tau))))))));
end f3;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g5 Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g5

base function for region 5: g(p,T)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
GibbsDerivsgdimensionless Gibbs funcion and dervatives w.r.t. pi and tau

Modelica definition

function g5 "base function for region 5: g(p,T)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature (K)";
  output Modelica.Media.Common.GibbsDerivs g 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
protected 
  Real[11] o "vector of auxiliary variables";
algorithm 
  assert(p > triple.ptriple,
    "IF97 medium function g5 called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");
  assert(p <= data.PLIMIT5,
    "IF97 medium function g5: input pressure (= " + String(p) + " Pa) is higher than 10 Mpa in region 5");
  assert(T <= 2273.15,
    "IF97 medium function g5: input temperature (= " + String(T) + " K) is higher than limit of 2273.15K in region 5");
  g.p := p;
  g.T := T;
  g.R := data.RH2O;
  g.pi := p/data.PSTAR5;
  g.tau := data.TSTAR5/T;
  o[1] := g.tau*g.tau;
  o[2] := -0.0045942820899910*o[1];
  o[3] := 0.00217746787145710 + o[2];
  o[4] := o[3]*g.tau;
  o[5] := o[1]*g.tau;
  o[6] := o[1]*o[1];
  o[7] := o[6]*o[6];
  o[8] := o[7]*g.tau;
  o[9] := -7.9449656719138e-6*o[8];
  o[10] := g.pi*g.pi;
  o[11] := -0.0137828462699730*o[1];

  g.g := g.pi*(-0.000125631835895920 + o[4] + g.pi*(-3.9724828359569e-6*o[
    8] + 1.29192282897840e-7*o[5]*g.pi)) + (-0.0248051489334660 + g.tau*(
    0.36901534980333 + g.tau*(-3.11613182139250 + g.tau*(-13.1799836742010
     + (6.8540841634434 - 0.32961626538917*g.tau)*g.tau +
    Modelica.Math.log(g.pi)))))/o[5];

  g.gpi := (1.0 + g.pi*(-0.000125631835895920 + o[4] + g.pi*(o[9] +
    3.8757684869352e-7*o[5]*g.pi)))/g.pi;

  g.gpipi := (-1.00000000000000 + o[10]*(o[9] + 7.7515369738704e-7*o[5]*g.
     pi))/o[10];

  g.gtau := g.pi*(0.00217746787145710 + o[11] + g.pi*(-0.000035752345523612
    *o[7] + 3.8757684869352e-7*o[1]*g.pi)) + (0.074415446800398 + g.tau*(
    -0.73803069960666 + (3.11613182139250 + o[1]*(6.8540841634434 -
    0.65923253077834*g.tau))*g.tau))/o[6];

  g.gtautau := (-0.297661787201592 + g.tau*(2.21409209881998 + (-6.2322636427850
     - 0.65923253077834*o[5])*g.tau))/(o[6]*g.tau) + g.pi*(-0.0275656925399460
    *g.tau + g.pi*(-0.000286018764188897*o[1]*o[6]*g.tau +
    7.7515369738704e-7*g.pi*g.tau));

  g.gtaupi := 0.00217746787145710 + o[11] + g.pi*(-0.000071504691047224*o[
    7] + 1.16273054608056e-6*o[1]*g.pi);
end g5;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.gibbs Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.gibbs

Gibbs function for region 1, 2 or 5: g(p,T,region)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature (K) [K]
Integerregion IF97 region, 1, 2 or 5

Outputs

TypeNameDescription
Realgdimensionless Gibbs funcion

Modelica definition

function gibbs "Gibbs function for region 1, 2 or 5: g(p,T,region)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature (K)";
 input Integer region "IF97 region, 1, 2 or 5";
  output Real g "dimensionless Gibbs funcion";
protected 
  Modelica.Media.Common.GibbsDerivs gibbs 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
algorithm 
  assert(region == 1 or region == 2 or region == 5,
    "IF97 medium function gibbs called with wrong region (= " + String(region) + ").\n" +
    "Only regions 1, 2 or 5 are possible");
  if region
     == 1 then
    gibbs
   := g1(p,T);
  elseif 
  region == 2 then
    gibbs
   := g2(p,T);
  else
    gibbs
   := g5(p,T);
  end if;
  g :=
gibbs.g;
end gibbs;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g1pitau Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g1pitau

derivative of g w.r.t. pi and tau

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
Realpidimensionless pressure
Realtaudimensionless temperature
Realgpidimensionless dervative of Gibbs function w.r.t. pi
Realgtaudimensionless dervative of Gibbs function w.r.t. tau

Modelica definition

function g1pitau "derivative of g w.r.t. pi and tau"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature (K)";
  output Real pi "dimensionless pressure";
  output Real tau "dimensionless temperature";
  output Real gpi "dimensionless dervative of Gibbs function w.r.t. pi";
  output Real gtau "dimensionless dervative of Gibbs function w.r.t. tau";
protected 
  Real pi1 "dimensionless pressure";
  Real tau1 "dimensionless temperature";
  Real[28] o "vector of auxiliary variables";
algorithm 
  assert(p > triple.ptriple,
    "IF97 medium function g1pitau called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");
  assert(p <= 100.0e6,
    "IF97 medium function g1pitau: the input pressure (= " + String(p) + " Pa) is higher than 100 Mpa");
  assert(T >= 273.15,
    "IF97 medium function g1pitau: the temperature (= " + String(T) + " K) is lower than 273.15 K!");
  pi := p/data.PSTAR1;
  tau := data.TSTAR1/T;
  pi1 := 7.1 - pi;
  tau1 := -1.222 + tau;
  o[1] := tau1*tau1;
  o[2] := o[1]*tau1;
  o[3] := 1/o[2];
  o[4] := o[1]*o[1];
  o[5] := o[4]*o[4];
  o[6] := o[1]*o[5];
  o[7] := o[1]*o[4];
  o[8] := 1/o[4];
  o[9] := o[1]*o[4]*o[5];
  o[10] := o[4]*tau1;
  o[11] := 1/o[10];
  o[12] := o[4]*o[5];
  o[13] := o[5]*tau1;
  o[14] := 1/o[13];
  o[15] := pi1*pi1;
  o[16] := o[15]*pi1;
  o[17] := o[15]*o[15];
  o[18] := o[17]*o[17];
  o[19] := o[17]*o[18]*pi1;
  o[20] := o[15]*o[17];
  o[21] := o[5]*o[5];
  o[22] := o[21]*o[21];
  o[23] := o[22]*o[5]*tau1;
  o[24] := 1/o[23];
  o[25] := o[22]*o[5];
  o[26] := 1/o[25];
  o[27] := o[1]*o[22]*o[4]*tau1;
  o[28] := 1/o[27];
  gtau := pi1*((-0.00254871721114236 + o[1]*(0.00424944110961118 + (
    0.018990068218419 + (-0.021841717175414 - 0.00015851507390979*o[1])*o[
    1])*o[7]))/o[6] + pi1*(o[8]*(0.00141552963219801 + o[4]*(
    0.000047661393906987 + o[1]*(-0.0000132425535992538 -
    1.2358149370591e-14*o[9]))) + pi1*(o[11]*(0.000126718579380216 -
    5.11230768720618e-9*o[6]) + pi1*((0.000011212640954 + (
    1.30342445791202e-6 - 1.4341729937924e-12*o[12])*o[2])/o[7] + pi1*(
    3.24135974880936e-6*o[14] + o[16]*((1.40077319158051e-8 +
    1.04549227383804e-9*o[10])/o[12] + o[19]*(1.9941018075704e-17/(o[1]*o[
    21]*o[4]*o[5]) + o[15]*(-4.48827542684151e-19/o[22] + o[20]*(-1.00075970318621e-21
    *o[28] + pi1*(4.65957282962769e-22*o[26] + pi1*(-7.2912378325616e-23*
    o[24] + (3.83502057899078e-24*pi1)/(o[1]*o[22]*o[5])))))))))))) + o[3]
    *(-0.29265942426334 + tau1*(0.84548187169114 + o[1]*(3.3855169168385
     + tau1*(-1.91583926775744 + tau1*(0.47316115539684 + (-0.066465668798004
     + 0.0040607314991784*tau1)*tau1)))));
  gpi := pi1*(pi1*((0.000095038934535162 + o[4]*(8.4812393955936e-6 +
    2.55615384360309e-9*o[7]))*o[8] + pi1*(o[11]*(8.9701127632e-6 + (
    2.60684891582404e-6 + 5.7366919751696e-13*o[12])*o[2]) + pi1*(
    2.02584984300585e-6/o[5] + o[16]*(o[19]*(o[15]*(o[20]*(-7.63737668221055e-22
    /(o[1]*o[22]*o[4]) + pi1*(3.5842867920213e-22*o[28] + pi1*(-5.65070932023524e-23
    *o[26] + 2.99318679335866e-24*o[24]*pi1))) - 3.33001080055983e-19/(o[
    1]*o[21]*o[4]*o[5]*tau1)) + 1.44400475720615e-17/(o[21]*o[4]*o[5]*
    tau1)) + (1.01874413933128e-8 + 1.39398969845072e-9*o[10])/(o[1]*o[5]
    *tau1))))) + o[3]*(0.00094368642146534 + o[2]*(0.00060003561586052 +
    (-0.000095322787813974 + o[1]*(8.8283690661692e-6 +
    1.45389992595188e-15*o[9]))*tau1))) + o[14]*(-0.00028319080123804 + o[
    1]*(0.00060706301565874 + o[7]*(0.018990068218419 + tau1*(
    0.032529748770505 + (0.021841717175414 + 0.00005283835796993*o[1])*
    tau1))));
end g1pitau;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g2pitau Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g2pitau

derivative of g w.r.t. pi and tau

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
Realpidimensionless pressure
Realtaudimensionless temperature
Realgpidimensionless dervative of Gibbs function w.r.t. pi
Realgtaudimensionless dervative of Gibbs function w.r.t. tau

Modelica definition

function g2pitau "derivative of g w.r.t. pi and tau"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature (K)";
  output Real pi "dimensionless pressure";
  output Real tau "dimensionless temperature";
  output Real gpi "dimensionless dervative of Gibbs function w.r.t. pi";
  output Real gtau "dimensionless dervative of Gibbs function w.r.t. tau";
protected 
  Real tau2 "dimensionless temperature";
  Real[22] o "vector of auxiliary variables";
algorithm 
  assert(p > triple.ptriple,
    "IF97 medium function g2pitau called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");
  assert(p <= 100.0e6,
    "IF97 medium function g2pitau: the input pressure (= " + String(p) + " Pa) is higher than 100 Mpa");
  assert(T >= 273.15,
    "IF97 medium function g2pitau: the temperature (= " + String(T) + " K) is lower than 273.15 K!");
  assert(T <= 1073.15,
    "IF97 medium function g2pitau: the input temperature (= " + String(T) + " K) is higher than the limit of 1073.15 K");
  pi := p/data.PSTAR2;
  tau := data.TSTAR2/T;
  tau2 := -0.5 + tau;
  o[1] := tau*tau;
  o[2] := o[1]*o[1];
  o[3] := tau2*tau2;
  o[4] := o[3]*tau2;
  o[5] := o[3]*o[3];
  o[6] := o[5]*o[5];
  o[7] := o[6]*o[6];
  o[8] := o[5]*o[6]*o[7]*tau2;
  o[9] := o[3]*o[5];
  o[10] := o[5]*o[6]*tau2;
  o[11] := o[3]*o[7]*tau2;
  o[12] := o[3]*o[5]*o[6];
  o[13] := o[3]*o[5]*tau2;
  o[14] := o[5]*o[6]*o[7];
  o[15] := pi*pi;
  o[16] := o[15]*o[15];
  o[17] := o[15]*o[16];
  o[18] := o[5]*o[7]*tau2;
  o[19] := o[7]*o[7];
  o[20] := o[3]*o[5]*o[6]*tau2;
  o[21] := o[5]*o[7];
  o[22] := o[3]*o[7];
  gtau := (0.0280439559151 + tau*(-0.2858109552582 + tau*(1.2213149471784
     + tau*(-2.848163942888 + tau*(4.38395111945 + o[1]*(10.08665568018
     + (-0.5681726521544 + 0.06380539059921*tau)*tau))))))/(o[1]*o[2]) +
    pi*(-0.017834862292358 + tau2*(-0.09199202739273 + (-0.172743777250296
     - 0.30195167236758*o[4])*tau2) + pi*(-0.000033032641670203 + (-0.0003789797503263
     + o[3]*(-0.015757110897342 + o[4]*(-0.306581069554011 -
    0.000960283724907132*o[8])))*tau2 + pi*(4.3870667284435e-7 + o[3]*(-0.00009683303171571
     + o[4]*(-0.0090203547252888 - 1.42338887469272*o[8])) + pi*(-7.8847309559367e-10
     + (2.558143570457e-8 + 1.44676118155521e-6*tau2)*tau2 + pi*(
    0.0000160454534363627*o[9] + pi*((-5.0144299353183e-11 + o[10]*(-0.033874355714168
     - 836.35096769364*o[11]))*o[3] + pi*((-0.0000138839897890111 -
    0.973671060893475*o[12])*o[3]*o[6] + pi*(o[13]*(9.0049690883672e-11
     - 296.320827232793*o[14]) + pi*(2.57526266427144e-7*o[5]*o[6] + pi*(
    o[4]*(4.1627860840696e-19 + (-1.0234747095929e-12 -
    1.40254511313154e-8*o[5])*o[9]) + o[17]*(o[14]*(-2.34560435076256e-9
     + 5.3465159397045*o[18]) + o[15]*(-19.1874828272775*o[19]*o[6]*o[7]
     + o[15]*(o[11]*(1.78371690710842e-23 + (1.07202609066812e-11 -
    0.000201611844951398*o[10])*o[20]) + pi*(-1.24017662339842e-24*o[21]
     + pi*(0.000200482822351322*o[19]*o[5]*o[7] + pi*(-4.97975748452559e-14
    *o[19]*o[3]*o[5] + (1.90027787547159e-27 + o[12]*(
    2.21658861403112e-15 - 0.0000547344301999018*o[22]))*o[6]*o[7]*pi*
    tau2))))))))))))))));
  gpi := (1. + pi*(-0.0017731742473213 + tau2*(-0.017834862292358 + tau2*
    (-0.045996013696365 + (-0.057581259083432 - 0.05032527872793*o[4])*
    tau2)) + pi*(tau2*(-0.000066065283340406 + (-0.0003789797503263 + o[3]
    *(-0.007878555448671 + o[4]*(-0.087594591301146 -
    0.000053349095828174*o[8])))*tau2) + pi*(6.1445213076927e-8 + (
    1.31612001853305e-6 + o[3]*(-0.00009683303171571 + o[4]*(-0.0045101773626444
     - 0.122004760687947*o[8])))*tau2 + pi*(tau2*(-3.15389238237468e-9 +
    (5.116287140914e-8 + 1.92901490874028e-6*tau2)*tau2) + pi*(
    0.0000114610381688305*o[13] + pi*((-1.00288598706366e-10 + o[10]*(-0.012702883392813
     - 143.374451604624*o[11]))*o[4] + pi*(-4.1341695026989e-17 + (-8.8352662293707e-6
     - 0.272627897050173*o[12])*o[3]*o[6]*tau2 + pi*((9.0049690883672e-11
     - 65.8490727183984*o[14])*o[6] + pi*(1.78287415218792e-7*o[10] + pi*
    (o[5]*(1.0406965210174e-18 + (-1.0234747095929e-12 -
    1.0018179379511e-8*o[5])*o[9]) + o[17]*((-1.29412653835176e-9 +
    1.71088510070544*o[18])*o[8] + o[15]*(-6.05920510335078*o[19]*o[6]*o[
    7]*tau2 + o[15]*((1.78371690710842e-23 + (6.1258633752464e-12 -
    0.000084004935396416*o[10])*o[20])*o[21] + pi*(-1.24017662339842e-24*
    o[18] + pi*(0.0000832192847496054*o[19]*o[5]*o[7]*tau2 + pi*((
    1.75410265428146e-27 + o[12]*(1.32995316841867e-15 -
    0.0000226487297378904*o[22]))*o[3]*o[6]*o[7]*pi -
    2.93678005497663e-14*o[19]*o[3]*o[5]*tau2)))))))))))))))))/pi;
end g2pitau;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g5pitau Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.g5pitau

derivative of g w.r.t. pi and tau

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
Realpidimensionless pressure
Realtaudimensionless temperature
Realgpidimensionless dervative of Gibbs function w.r.t. pi
Realgtaudimensionless dervative of Gibbs function w.r.t. tau

Modelica definition

function g5pitau "derivative of g w.r.t. pi and tau"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature (K)";
  output Real pi "dimensionless pressure";
  output Real tau "dimensionless temperature";
  output Real gpi "dimensionless dervative of Gibbs function w.r.t. pi";
  output Real gtau "dimensionless dervative of Gibbs function w.r.t. tau";
protected 
  Real[3] o "vector of auxiliary variables";
algorithm 
  assert(p > triple.ptriple,
    "IF97 medium function g5pitau called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");
  assert(p <= data.PLIMIT5,
    "IF97 medium function g5pitau: input pressure (= " + String(p) + " Pa) is higher than 10 Mpa in region 5");
  assert(T <= 2273.15,
    "IF97 medium function g5pitau: input temperature (= " + String(T) + " K) is higher than limit of 2273.15 K in region 5");
  pi := p/data.PSTAR5;
  tau := data.TSTAR5/T;
  o[1] := tau*tau;
  o[2] := o[1]*o[1];
  o[3] := o[2]*o[2];
  gtau := pi*(0.0021774678714571 - 0.013782846269973*o[1] + pi*(-0.0000357523455236121
    *o[3] + 3.8757684869352e-7*o[1]*pi)) + (0.074415446800398 + tau*(-0.73803069960666
     + (3.1161318213925 + o[1]*(6.8540841634434 - 0.65923253077834*tau))*
    tau))/o[2];
  gpi := (1.0 + pi*(-0.00012563183589592 + (0.0021774678714571 -
    0.004594282089991*o[1])*tau + pi*(-7.9449656719138e-6*o[3]*tau +
    3.8757684869352e-7*o[1]*pi*tau)))/pi;
end g5pitau;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.f3deltatau Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.f3deltatau

1st derivatives of f w.r.t. delta and tau

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
Realdeltadimensionless density
Realtaudimensionless temperature
Realfdeltadimensionless dervative of Helmholtz function w.r.t. delta
Realftaudimensionless dervative of Helmholtz function w.r.t. tau

Modelica definition

function f3deltatau "1st derivatives of f w.r.t. delta and tau"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature (K)";
  output Real delta "dimensionless density";
  output Real tau "dimensionless temperature";
  output Real fdelta 
    "dimensionless dervative of Helmholtz function w.r.t. delta";
  output Real ftau "dimensionless dervative of Helmholtz function w.r.t. tau";
protected 
  Real[13] o "vector of auxiliary variables";
algorithm 
  tau := data.TCRIT/T;
  delta := if (d == data.DCRIT and T == data.TCRIT) then 1 + Modelica.
    Constants.eps else d/data.DCRIT;
  o[1] := tau*tau;
  o[2] := o[1]*o[1];
  o[3] := o[2]*o[2];
  o[4] := o[3]*tau;
  o[5] := o[1]*o[2]*o[3]*tau;
  o[6] := o[2]*o[3];
  o[7] := o[1]*o[3];
  o[8] := o[3]*o[3];
  o[9] := o[1]*o[2]*o[8];
  o[10] := o[1]*o[2]*o[8]*tau;
  o[11] := o[3]*o[8];
  o[12] := o[1]*o[3]*o[8];
  o[13] := o[3]*o[8]*tau;
  fdelta := (1.0658070028513 + delta*(o[1]*(-1.2654315477714 + o[2]*(-1.1524407806681
     + (0.88521043984318 - 0.64207765181607*o[1])*o[4])) + delta*(
    0.76986920373342 + o[1]*(-1.70429417648412 + o[2]*(9.7944563083754 +
    (-6.100523451393 + (0.078841073758308 + 0.25116816848616*o[2])*o[5])*
    tau)) + delta*(-0.8399798909613 + o[1]*(4.169939870838 + o[1]*(-6.056974507071
     + o[6]*(-0.0246442911521889 - 1.42788107204769*o[7]))) + delta*(
    0.175936297894 + o[1]*(-1.77905741714956 + o[1]*(3.62288282878932 +
    2.82089800351868*o[9])) + delta*(delta*(-0.133052405238576 + o[1]*(
    0.565564509990552 + 0.98617670687766*o[11]) + delta*(-0.094523605689436
    *o[1] + delta*(-0.118674762819776*o[12] + delta*(o[1]*(
    0.00521306582652756 + 0.0290780142333399*o[11]) + delta*(
    0.00080964802996215 - 0.000494162889679965*delta*o[12] -
    0.0016557679795037*tau))))) + (0.5385256313166 + o[1]*(-1.6456811629477
     - 2.5435531020579*o[10]))*tau))))))/delta;
  ftau := 20.944396974307 + tau*(-15.3735415757432 + o[2]*tau*(
    18.3301634515678 + o[1]*tau*(-28.08078114862 + o[1]*(14.4640436358204
     - 0.194503669468755*o[1]*o[3]*tau)))) + delta*((-2.5308630955428 + o[
    2]*(-6.9146446840086 + (13.2781565976477 - 10.9153200808732*o[1])*o[4]))
    *tau + delta*(tau*(-1.70429417648412 + o[2]*(29.3833689251262 + (-21.3518320798755
     + (0.867251811341388 + 3.26518619032008*o[2])*o[5])*tau)) + delta*((
    2.779959913892 + o[1]*(-8.075966009428 + o[6]*(-0.131436219478341 -
    12.37496929108*o[7])))*tau + delta*((-0.88952870857478 + o[1]*(
    3.62288282878932 + 18.3358370228714*o[9]))*tau + delta*(
    0.10770512626332 + o[1]*(-0.98740869776862 - 13.2264761307011*o[10])
     + delta*((0.188521503330184 + 4.27343239646986*o[11])*tau + delta*(-0.027006744482696
    *tau + delta*(-0.385692979164272*o[13] + delta*(delta*(-0.00016557679795037
     - 0.00116802137560719*delta*o[13]) + (0.00115845907256168 +
    0.0840031522296486*o[11])*tau)))))))));
end f3deltatau;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tph1 Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tph1

inverse function for region 1: T(p,h)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]

Outputs

TypeNameDescription
TemperatureTtemperature (K) [K]

Modelica definition

function tph1 "inverse function for region 1: T(p,h)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  output SI.Temperature T "temperature (K)";
protected 
  Real pi "dimensionless pressure";
  Real eta1 "dimensionless specific enthalpy";
  Real[3] o "vector of auxiliary variables";
algorithm 
  assert(p > triple.ptriple,
    "IF97 medium function tph1 called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");
  pi := p/data.PSTAR2;
  eta1 := h/data.HSTAR1 + 1.0;
  o[1] := eta1*eta1;
  o[2] := o[1]*o[1];
  o[3] := o[2]*o[2];
  T := -238.724899245210 - 13.3917448726020*pi + eta1*(404.21188637945 +
    43.211039183559*pi + eta1*(113.497468817180 - 54.010067170506*pi +
    eta1*(30.5358922039160*pi + eta1*(-6.5964749423638*pi + o[1]*(-5.8457616048039
     + o[2]*(pi*(0.0093965400878363 + (-0.0000258586412820730 +
    6.6456186191635e-8*pi)*pi) + o[2]*o[3]*(-0.000152854824131400 + o[1]*
    o[3]*(-1.08667076953770e-6 + pi*(1.15736475053400e-7 + pi*(-4.0644363084799e-9
     + pi*(8.0670734103027e-11 + pi*(-9.3477771213947e-13 + (
    5.8265442020601e-15 - 1.50201859535030e-17*pi)*pi))))))))))));
end tph1;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps1 Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps1

inverse function for region 1: T(p,s)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
TemperatureTtemperature (K) [K]

Modelica definition

function tps1 "inverse function for region 1: T(p,s)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  output SI.Temperature T "temperature (K)";
protected 
  constant SI.Pressure pstar=1.0e6;
  constant SI.SpecificEntropy sstar=1.0e3;
  Real pi "dimensionless pressure";
  Real sigma1 "dimensionless specific entropy";
  Real[6] o "vector of auxiliary variables";
algorithm 
  pi := p/pstar;
  assert(p > triple.ptriple,
    "IF97 medium function tps1 called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");

  sigma1 := s/sstar + 2.0;
  o[1] := sigma1*sigma1;
  o[2] := o[1]*o[1];
  o[3] := o[2]*o[2];
  o[4] := o[3]*o[3];
  o[5] := o[4]*o[4];
  o[6] := o[1]*o[2]*o[4];

  T := 174.782680583070 + sigma1*(34.806930892873 + sigma1*(
    6.5292584978455 + (0.33039981775489 + o[3]*(-1.92813829231960e-7 -
    2.49091972445730e-23*o[2]*o[4]))*sigma1)) + pi*(-0.261076364893320 +
    pi*(0.00056608900654837 + pi*(o[1]*o[3]*(2.64004413606890e-13 +
    7.8124600459723e-29*o[6]) - 3.07321999036680e-31*o[5]*pi) + sigma1*(-0.00032635483139717
     + sigma1*(0.000044778286690632 + o[1]*o[2]*(-5.1322156908507e-10 -
    4.2522657042207e-26*o[6])*sigma1))) + sigma1*(0.225929659815860 +
    sigma1*(-0.064256463395226 + sigma1*(0.0078876289270526 + o[3]*sigma1
    *(3.5672110607366e-10 + 1.73324969948950e-24*o[1]*o[4]*sigma1)))));
end tps1;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tph2 Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tph2

reverse function for region 2: T(p,h)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]

Outputs

TypeNameDescription
TemperatureTtemperature (K) [K]

Modelica definition

function tph2 "reverse function for region 2: T(p,h)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  output SI.Temperature T "temperature (K)";
protected 
  Real pi "dimensionless pressure";
  Real pi2b "dimensionless pressure";
  Real pi2c "dimensionless pressure";
  Real eta "dimensionless specific enthalpy";
  Real etabc "dimensionless specific enthalpy";
  Real eta2a "dimensionless specific enthalpy";
  Real eta2b "dimensionless specific enthalpy";
  Real eta2c "dimensionless specific enthalpy";
  Real[8] o "vector of auxiliary variables";
algorithm 
  pi := p*data.IPSTAR;
  eta := h*data.IHSTAR;
  etabc := h*1.0e-3;
  if (pi < 4.0) then
    eta2a := eta - 2.1;
    o[1] := eta2a*eta2a;
    o[2] := o[1]*o[1];
    o[3] := pi*pi;
    o[4] := o[3]*o[3];
    o[5] := o[3]*pi;
    T := 1089.89523182880 + (1.84457493557900 - 0.0061707422868339*pi)*pi
       + eta2a*(849.51654495535 - 4.1792700549624*pi + eta2a*(-107.817480918260
       + (6.2478196935812 - 0.310780466295830*pi)*pi + eta2a*(
      33.153654801263 - 17.3445631081140*pi + o[2]*(-7.4232016790248 + pi
      *(-200.581768620960 + 11.6708730771070*pi) + o[1]*(271.960654737960
      *pi + o[1]*(-455.11318285818*pi + eta2a*(1.38657242832260*o[4] + o[
      1]*o[2]*(3091.96886047550*pi + o[1]*(11.7650487243560 + o[2]*(-13551.3342407750
      *o[5] + o[2]*(-62.459855192507*o[3]*o[4]*pi + o[2]*(o[4]*(
      235988.325565140 + 7399.9835474766*pi) + o[1]*(19127.7292396600*o[3]
      *o[4] + o[1]*(o[3]*(1.28127984040460e8 - 551966.97030060*o[5]) + o[
      1]*(-9.8554909623276e8*o[3] + o[1]*(2.82245469730020e9*o[3] + o[1]*
      (o[3]*(-3.5948971410703e9 + 3.7154085996233e6*o[5]) + o[1]*pi*(
      252266.403578720 + pi*(1.72273499131970e9 + pi*(1.28487346646500e7
       + (-1.31052365450540e7 - 415351.64835634*o[3])*pi))))))))))))))))))));
  elseif (pi < (0.12809002730136e-03*etabc - 0.67955786399241)*etabc +
      0.90584278514723e3) then
    eta2b := eta - 2.6;
    pi2b := pi - 2.0;
    o[1] := pi2b*pi2b;
    o[2] := o[1]*pi2b;
    o[3] := o[1]*o[1];
    o[4] := eta2b*eta2b;
    o[5] := o[4]*o[4];
    o[6] := o[4]*o[5];
    o[7] := o[5]*o[5];
    T := 1489.50410795160 + 0.93747147377932*pi2b + eta2b*(
      743.07798314034 + o[2]*(0.000110328317899990 - 1.75652339694070e-18
      *o[1]*o[3]) + eta2b*(-97.708318797837 + pi2b*(3.3593118604916 +
      pi2b*(-0.0218107553247610 + pi2b*(0.000189552483879020 + (
      2.86402374774560e-7 - 8.1456365207833e-14*o[2])*pi2b))) + o[5]*(
      3.3809355601454*pi2b + o[4]*(-0.108297844036770*o[1] + o[5]*(
      2.47424647056740 + (0.168445396719040 + o[1]*(0.00308915411605370
       - 0.0000107798573575120*pi2b))*pi2b + o[6]*(-0.63281320016026 +
      pi2b*(0.73875745236695 + (-0.046333324635812 + o[1]*(-0.000076462712454814
       + 2.82172816350400e-7*pi2b))*pi2b) + o[6]*(1.13859521296580 + pi2b
      *(-0.47128737436186 + o[1]*(0.00135555045549490 + (
      0.0000140523928183160 + 1.27049022719450e-6*pi2b)*pi2b)) + o[5]*(-0.47811863648625
       + (0.150202731397070 + o[2]*(-0.0000310838143314340 + o[1]*(-1.10301392389090e-8
       - 2.51805456829620e-11*pi2b)))*pi2b + o[5]*o[7]*(
      0.0085208123431544 + pi2b*(-0.00217641142197500 + pi2b*(
      0.000071280351959551 + o[1]*(-1.03027382121030e-6 + (
      7.3803353468292e-8 + 8.6934156344163e-15*o[3])*pi2b))))))))))));
  else
    eta2c := eta - 1.8;
    pi2c := pi + 25.0;
    o[1] := pi2c*pi2c;
    o[2] := o[1]*o[1];
    o[3] := o[1]*o[2]*pi2c;
    o[4] := 1/o[3];
    o[5] := o[1]*o[2];
    o[6] := eta2c*eta2c;
    o[7] := o[2]*o[2];
    o[8] := o[6]*o[6];
    T := eta2c*((859777.22535580 + o[1]*(482.19755109255 +
      1.12615974072300e-12*o[5]))/o[1] + eta2c*((-5.8340131851590e11 + (
      2.08255445631710e10 + 31081.0884227140*o[2])*pi2c)/o[5] + o[6]*(o[8]
      *(o[6]*(1.23245796908320e-7*o[5] + o[6]*(-1.16069211309840e-6*o[5]
       + o[8]*(0.0000278463670885540*o[5] + (-0.00059270038474176*o[5] +
      0.00129185829918780*o[5]*o[6])*o[8]))) - 10.8429848800770*pi2c) + o[
      4]*(7.3263350902181e12 + o[7]*(3.7966001272486 + (-0.045364172676660
       - 1.78049822406860e-11*o[2])*pi2c))))) + o[4]*(-3.2368398555242e12
       + pi2c*(3.5825089945447e11 + pi2c*(-1.07830682174700e10 + o[1]*
      pi2c*(610747.83564516 + pi2c*(-25745.7236041700 + (1208.23158659360
       + 1.45591156586980e-13*o[5])*pi2c)))));
  end if;
end tph2;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2a Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2a

reverse function for region 2a: T(p,s)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
TemperatureTtemperature (K) [K]

Modelica definition

function tps2a "reverse function for region 2a: T(p,s)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  output SI.Temperature T "temperature (K)";
protected 
  Real[12] o "vector of auxiliary variables";
  constant Real IPSTAR=1.0e-6 "scaling variable";
  constant Real ISSTAR2A=1/2000.0 "scaling variable";
  Real pi "dimensionless pressure";
  Real sigma2a "dimensionless specific entropy";
algorithm 
  pi := p*IPSTAR;
  sigma2a := s*ISSTAR2A - 2.0;
  o[1] := pi^0.5;
  o[2] := sigma2a*sigma2a;
  o[3] := o[2]*o[2];
  o[4] := o[3]*o[3];
  o[5] := o[4]*o[4];
  o[6] := pi^0.25;
  o[7] := o[2]*o[4]*o[5];
  o[8] := 1/o[7];
  o[9] := o[3]*sigma2a;
  o[10] := o[2]*o[3]*sigma2a;
  o[11] := o[3]*o[4]*sigma2a;
  o[12] := o[2]*sigma2a;
  T := ((-392359.83861984 + (515265.73827270 + o[3]*(40482.443161048 + o[
    2]*o[3]*(-321.93790923902 + o[2]*(96.961424218694 - 22.8678463717730*
    sigma2a))))*sigma2a)/(o[4]*o[5]) + o[6]*((-449429.14124357 + o[3]*(-5011.8336020166
     + 0.35684463560015*o[4]*sigma2a))/(o[2]*o[5]*sigma2a) + o[6]*(o[8]*(
    44235.335848190 + o[9]*(-13673.3888117080 + o[3]*(421632.60207864 + (
    22516.9258374750 + o[10]*(474.42144865646 - 149.311307976470*sigma2a))
    *sigma2a))) + o[6]*((-197811.263204520 - 23554.3994707600*sigma2a)/(o[
    2]*o[3]*o[4]*sigma2a) + o[6]*((-19070.6163020760 + o[11]*(
    55375.669883164 + (3829.3691437363 - 603.91860580567*o[2])*o[3]))*o[8]
     + o[6]*((1936.31026203310 + o[2]*(4266.0643698610 + o[2]*o[3]*o[4]*(
    -5978.0638872718 - 704.01463926862*o[9])))/(o[2]*o[4]*o[5]*sigma2a)
     + o[1]*((338.36784107553 + o[12]*(20.8627866351870 + (
    0.033834172656196 - 0.000043124428414893*o[12])*o[3]))*sigma2a + o[6]
    *(166.537913564120 + sigma2a*(-139.862920558980 + o[3]*(-0.78849547999872
     + (0.072132411753872 + o[3]*(-0.0059754839398283 + (-0.0000121413589539040
     + 2.32270967338710e-7*o[2])*o[3]))*sigma2a)) + o[6]*(-10.5384635661940
     + o[3]*(2.07189254965020 + (-0.072193155260427 + 2.07498870811200e-7
    *o[4])*o[9]) + o[6]*(o[6]*(o[12]*(0.210375278936190 +
    0.000256812397299990*o[3]*o[4]) + (-0.0127990029337810 -
    8.2198102652018e-6*o[11])*o[6]*o[9]) + o[10]*(-0.0183406579113790 +
    2.90362723486960e-7*o[2]*o[4]*sigma2a)))))))))))/(o[1]*pi);
end tps2a;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2b Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2b

reverse function for region 2b: T(p,s)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
TemperatureTtemperature (K) [K]

Modelica definition

function tps2b "reverse function for region 2b: T(p,s)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  output SI.Temperature T "temperature (K)";
protected 
  Real[8] o "vector of auxiliary variables";
  constant Real IPSTAR=1.0e-6 "scaling variable";
  constant Real ISSTAR2B=1/785.3 "scaling variable";
  Real pi "dimensionless pressure";
  Real sigma2b "dimensionless specific entropy";
algorithm 
  pi := p*IPSTAR;
  sigma2b := 10.0 - s*ISSTAR2B;
  o[1] := pi*pi;
  o[2] := o[1]*o[1];
  o[3] := sigma2b*sigma2b;
  o[4] := o[3]*o[3];
  o[5] := o[4]*o[4];
  o[6] := o[3]*o[5]*sigma2b;
  o[7] := o[3]*o[5];
  o[8] := o[3]*sigma2b;
  T := (316876.65083497 + 20.8641758818580*o[6] + pi*(-398593.99803599 -
    21.8160585188770*o[6] + pi*(223697.851942420 + (-2784.17034458170 +
    9.9207436071480*o[7])*sigma2b + pi*(-75197.512299157 + (
    2970.86059511580 + o[7]*(-3.4406878548526 + 0.38815564249115*sigma2b))
    *sigma2b + pi*(17511.2950857500 + sigma2b*(-1423.71128544490 + (
    1.09438033641670 + 0.89971619308495*o[4])*o[4]*sigma2b) + pi*(-3375.9740098958
     + (471.62885818355 + o[4]*(-1.91882419936790 + o[8]*(
    0.41078580492196 - 0.33465378172097*sigma2b)))*sigma2b + pi*(
    1387.00347775050 + sigma2b*(-406.63326195838 + sigma2b*(
    41.727347159610 + o[3]*(2.19325494345320 + sigma2b*(-1.03200500090770
     + (0.35882943516703 + 0.0052511453726066*o[8])*sigma2b)))) + pi*(
    12.8389164507050 + sigma2b*(-2.86424372193810 + sigma2b*(
    0.56912683664855 + (-0.099962954584931 + o[4]*(-0.0032632037778459 +
    0.000233209225767230*sigma2b))*sigma2b)) + pi*(-0.153348098574500 + (
    0.0290722882399020 + 0.00037534702741167*o[4])*sigma2b + pi*(
    0.00172966917024110 + (-0.00038556050844504 - 0.000035017712292608*o[
    3])*sigma2b + pi*(-0.0000145663936314920 + 5.6420857267269e-6*sigma2b
     + pi*(4.1286150074605e-8 + (-2.06846711188240e-8 +
    1.64093936747250e-9*sigma2b)*sigma2b))))))))))))/(o[1]*o[2]);
end tps2b;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2c Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2c

reverse function for region 2c: T(p,s)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
TemperatureTtemperature (K) [K]

Modelica definition

function tps2c "reverse function for region 2c: T(p,s)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  output SI.Temperature T "temperature (K)";
protected 
  constant Real IPSTAR=1.0e-6 "scaling variable";
  constant Real ISSTAR2C=1/2925.1 "scaling variable";
  Real pi "dimensionless pressure";
  Real sigma2c "dimensionless specific entropy";
  Real[3] o "vector of auxiliary variables";
algorithm 
  pi := p*IPSTAR;
  sigma2c := 2.0 - s*ISSTAR2C;
  o[1] := pi*pi;
  o[2] := sigma2c*sigma2c;
  o[3] := o[2]*o[2];
  T := (909.68501005365 + 2404.56670884200*sigma2c + pi*(-591.62326387130
     + pi*(541.45404128074 + sigma2c*(-270.983084111920 + (
    979.76525097926 - 469.66772959435*sigma2c)*sigma2c) + pi*(
    14.3992746047230 + (-19.1042042304290 + o[2]*(5.3299167111971 -
    21.2529753759340*sigma2c))*sigma2c + pi*(-0.311473344137600 + (
    0.60334840894623 - 0.042764839702509*sigma2c)*sigma2c + pi*(
    0.0058185597255259 + (-0.0145970082847530 + 0.0056631175631027*o[3])*
    sigma2c + pi*(-0.000076155864584577 + sigma2c*(0.000224403429193320
     - 0.0000125610950134130*o[2]*sigma2c) + pi*(6.3323132660934e-7 + (-2.05419896753750e-6
     + 3.6405370390082e-8*sigma2c)*sigma2c + pi*(-2.97598977892150e-9 +
    1.01366185297630e-8*sigma2c + pi*(5.9925719692351e-12 + sigma2c*(-2.06778701051640e-11
     + o[2]*(-2.08742781818860e-11 + (1.01621668250890e-10 -
    1.64298282813470e-10*sigma2c)*sigma2c))))))))))))/o[1];
end tps2c;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2 Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tps2

reverse function for region 2: T(p,s)

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
TemperatureTtemperature (K) [K]

Modelica definition

function tps2 "reverse function for region 2: T(p,s)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  output SI.Temperature T "temperature (K)";
protected 
  Real pi "dimensionless pressure";
  constant SI.SpecificEntropy SLIMIT=5.85e3 
    "subregion boundary specific entropy between regions 2a and 2b";
algorithm 
  if p < 4.0e6 then
    T := tps2a(p, s);
  elseif s > SLIMIT then
    T := tps2b(p, s);
  else
    T := tps2c(p, s);
  end if;
end tps2;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tsat Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tsat

region 4 saturation temperature as a function of pressure

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]

Outputs

TypeNameDescription
Temperaturet_sattemperature [K]

Modelica definition

function tsat 
  "region 4 saturation temperature as a function of pressure"
  annotation(derivative=tsat_der);

  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  output SI.Temperature t_sat "temperature";
protected 
  Real pi "dimensionless pressure";
  Real[20] o "vector of auxiliary variables";
algorithm 
  assert(p > triple.ptriple,
    "IF97 medium function tsat called with too low pressure\n" +
    "p = " + String(p) + " Pa <= " + String(triple.ptriple) + " Pa (triple point pressure)");
//  assert(p <= data.PCRIT,
//    "tsat: input pressure is higher than the critical point pressure");
  pi := min(p,data.PCRIT)*data.IPSTAR;
  o[1] := pi^0.25;
  o[2] := -3.2325550322333e6*o[1];
  o[3] := pi^0.5;
  o[4] := -724213.16703206*o[3];
  o[5] := 405113.40542057 + o[2] + o[4];
  o[6] := -17.0738469400920*o[1];
  o[7] := 14.9151086135300 + o[3] + o[6];
  o[8] := -4.0*o[5]*o[7];
  o[9] := 12020.8247024700*o[1];
  o[10] := 1167.05214527670*o[3];
  o[11] := -4823.2657361591 + o[10] + o[9];
  o[12] := o[11]*o[11];
  o[13] := o[12] + o[8];
  o[14] := o[13]^0.5;
  o[15] := -o[14];
  o[16] := -12020.8247024700*o[1];
  o[17] := -1167.05214527670*o[3];
  o[18] := 4823.2657361591 + o[15] + o[16] + o[17];
  o[19] := 1/o[18];
  o[20] := 2.0*o[19]*o[5];

  t_sat := 0.5*(650.17534844798 + o[20] - (-4.0*(-0.238555575678490 +
    1300.35069689596*o[19]*o[5]) + (650.17534844798 + o[20])^2.0)^0.5);
end tsat;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.dtsatofp Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.dtsatofp

derivative of saturation temperature w.r.t. pressure

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]

Outputs

TypeNameDescription
Realdtsatderivative of T w.r.t. p [K/Pa]

Modelica definition

function dtsatofp 
  "derivative of saturation temperature w.r.t. pressure"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  output Real dtsat(unit="K/Pa") "derivative of T w.r.t. p";
protected 
  Real pi "dimensionless pressure";
  Real[49] o "vector of auxiliary variables";
algorithm 
  pi := max(Modelica.Constants.small,p*data.IPSTAR);
  o[1] := pi^0.75;
  o[2] := 1/o[1];
  o[3] := -4.268461735023*o[2];
  o[4] := sqrt(pi);
  o[5] := 1/o[4];
  o[6] := 0.5*o[5];
  o[7] := o[3] + o[6];
  o[8] := pi^0.25;
  o[9] := -3.2325550322333e6*o[8];
  o[10] := -724213.16703206*o[4];
  o[11] := 405113.40542057 + o[10] + o[9];
  o[12] := -4*o[11]*o[7];
  o[13] := -808138.758058325*o[2];
  o[14] := -362106.58351603*o[5];
  o[15] := o[13] + o[14];
  o[16] := -17.073846940092*o[8];
  o[17] := 14.91510861353 + o[16] + o[4];
  o[18] := -4*o[15]*o[17];
  o[19] := 3005.2061756175*o[2];
  o[20] := 583.52607263835*o[5];
  o[21] := o[19] + o[20];
  o[22] := 12020.82470247*o[8];
  o[23] := 1167.0521452767*o[4];
  o[24] := -4823.2657361591 + o[22] + o[23];
  o[25] := 2.0*o[21]*o[24];
  o[26] := o[12] + o[18] + o[25];
  o[27] := -4.0*o[11]*o[17];
  o[28] := o[24]*o[24];
  o[29] := o[27] + o[28];
  o[30] := sqrt(o[29]);
  o[31] := 1/o[30];
  o[32] := (-o[30]);
  o[33] := -12020.82470247*o[8];
  o[34] := -1167.0521452767*o[4];
  o[35] := 4823.2657361591 + o[32] + o[33] + o[34];
  o[36] := o[30];
  o[37] := -4823.2657361591 + o[22] + o[23] + o[36];
  o[38] := o[37]*o[37];
  o[39] := 1/o[38];
  o[40] := -1.72207339365771*o[30];
  o[41] := 21592.2055343628*o[8];
  o[42] := o[30]*o[8];
  o[43] := -8192.87114842946*o[4];
  o[44] := -0.510632954559659*o[30]*o[4];
  o[45] := -3100.02526152368*o[1];
  o[46] := pi;
  o[47] := 1295.95640782102*o[46];
  o[48] := 2862.09212505088 + o[40] + o[41] + o[42] + o[43] + o[44] + o[
    45] + o[47];
  o[49] := 1/(o[35]*o[35]);
  dtsat := data.IPSTAR*0.5*((2.0*o[15])/o[35] - 2.*o[11]*(-3005.2061756175
    *o[2] - 0.5*o[26]*o[31] - 583.52607263835*o[5])*o[49] - (
    20953.46356643991*(o[39]*(1295.95640782102 + 5398.05138359071*o[2] +
    0.25*o[2]*o[30] - 0.861036696828853*o[26]*o[31] - 0.255316477279829*o[
    26]*o[31]*o[4] - 4096.43557421473*o[5] - 0.255316477279829*o[30]*o[5]
     - 2325.01894614276/o[8] + 0.5*o[26]*o[31]*o[8]) - 2.0*(o[19] + o[20]
     + 0.5*o[26]*o[31])*o[48]*o[37]^(-3)))/sqrt(o[39]*o[48]));
end dtsatofp;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tsat_der Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.tsat_der

derivative function for tsat

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
Realder_p pressure derivatrive [Pa/s]

Outputs

TypeNameDescription
Realder_tsattemperature derivative [K/s]

Modelica definition

function tsat_der "derivative function for tsat"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input Real der_p(unit="Pa/s") "pressure derivatrive";
  output Real der_tsat(unit="K/s") "temperature derivative";
protected 
  Real dtp;
algorithm 
  dtp := dtsatofp(p);
  der_tsat := dtp*der_p;
end tsat_der;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.psat Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.psat

region 4 saturation pressure as a functionx of temperature

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
Pressurep_satpressure [Pa]

Modelica definition

function psat 
  "region 4 saturation pressure as a functionx of temperature"
  annotation(derivative=psat_der);

  extends Modelica.Icons.Function;
  input SI.Temperature T "temperature (K)";
  output SI.Pressure p_sat "pressure";
protected 
  Real[8] o "vector of auxiliary variables";
  Real Tlim=min(T, data.TCRIT);
algorithm 
  assert(T >= 273.16,
    "IF97 medium function psat: input temperature (= " + String(triple.ptriple) + " K).\n" +
    "lower than the triple point temperature 273.16 K");
  o[1] := -650.17534844798 + Tlim;
  o[2] := 1/o[1];
  o[3] := -0.238555575678490*o[2];
  o[4] := o[3] + Tlim;
  o[5] := -4823.2657361591*o[4];
  o[6] := o[4]*o[4];
  o[7] := 14.9151086135300*o[6];
  o[8] := 405113.40542057 + o[5] + o[7];
  p_sat := 16.0e6*o[8]*o[8]*o[8]*o[8]*1/(3.2325550322333e6 -
    12020.8247024700*o[4] + 17.0738469400920*o[6] + (-4.0*(-724213.16703206
     + 1167.05214527670*o[4] + o[6])*o[8] + (-3.2325550322333e6 +
    12020.8247024700*o[4] - 17.0738469400920*o[6])^2.0)^0.5)^4.0;
end psat;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.dptofT Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.dptofT

derivative of pressure w.r.t. temperature along the saturation pressure curve

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
Realdpttemperature derivative of pressure [Pa/K]

Modelica definition

function dptofT 
  "derivative of pressure w.r.t. temperature along the saturation pressure curve"

  extends Modelica.Icons.Function;
  input SI.Temperature T "temperature (K)";
  output Real dpt(unit = "Pa/K") "temperature derivative of pressure";
protected 
  Real[31] o "vector of auxiliary variables";
  Real Tlim "temeprature limited to TCRIT";
algorithm 
  Tlim := min(T, data.TCRIT);
  o[1] := -650.17534844798 + Tlim;
  o[2] := 1/o[1];
  o[3] := -0.238555575678490*o[2];
  o[4] := o[3] + Tlim;
  o[5] := -4823.2657361591*o[4];
  o[6] := o[4]*o[4];
  o[7] := 14.9151086135300*o[6];
  o[8] := 405113.40542057 + o[5] + o[7];
  o[9] := o[8]*o[8];
  o[10] := o[9]*o[9];
  o[11] := o[1]*o[1];
  o[12] := 1/o[11];
  o[13] := 0.238555575678490*o[12];
  o[14] := 1.00000000000000 + o[13];
  o[15] := 12020.8247024700*o[4];
  o[16] := -17.0738469400920*o[6];
  o[17] := -3.2325550322333e6 + o[15] + o[16];
  o[18] := -4823.2657361591*o[14];
  o[19] := 29.8302172270600*o[14]*o[4];
  o[20] := o[18] + o[19];
  o[21] := 1167.05214527670*o[4];
  o[22] := -724213.16703206 + o[21] + o[6];
  o[23] := o[17]*o[17];
  o[24] := -4.0000000000000*o[22]*o[8];
  o[25] := o[23] + o[24];
  o[26] := sqrt(o[25]);
  o[27] := -12020.8247024700*o[4];
  o[28] := 17.0738469400920*o[6];
  o[29] := 3.2325550322333e6 + o[26] + o[27] + o[28];
  o[30] := o[29]*o[29];
  o[31] := o[30]*o[30];
  dpt := 1e6*((-64.0*o[10]*(-12020.8247024700*o[14] + 34.147693880184*o[
    14]*o[4] + (0.5*(-4.0*o[20]*o[22] + 2.00000000000000*o[17]*(
    12020.8247024700*o[14] - 34.147693880184*o[14]*o[4]) - 4.0*(
    1167.05214527670*o[14] + 2.0*o[14]*o[4])*o[8]))/o[26]))/(o[29]*o[31])
     + (64.*o[20]*o[8]*o[9])/o[31]);
end dptofT;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.psat_der Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.psat_der

derivative function for psat

Information

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
TemperatureT temperature (K) [K]
Realder_T temperature derivative [K/s]

Outputs

TypeNameDescription
Realder_psatpressure [Pa/s]

Modelica definition

function psat_der "derivative function for psat"
  extends Modelica.Icons.Function;
  input SI.Temperature T "temperature (K)";
  input Real der_T(unit = "K/s") "temperature derivative";
  output Real der_psat(unit = "Pa/s") "pressure";
protected 
  Real dpt;
algorithm 
  dpt := dptofT(T);
  der_psat := dpt*der_T;
end psat_der;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p1_hs Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p1_hs

pressure as a function of ehtnalpy and entropy in region 1

Information


Equation number 1 from:
The International Association for the Properties of Water and Steam
Gaithersburg, Maryland, USA
September 2001
Supplementary Release on  Backward Equations for Pressure as a Function of Enthalpy and Entropy p(h,s) to the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of Water and Steam

  

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
SpecificEnthalpyh specific enthalpy [J/kg]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
PressurepPressure [Pa]

Modelica definition

function p1_hs 
  "pressure as a function of ehtnalpy and entropy in region 1"
  extends Modelica.Icons.Function;
  input SI.SpecificEnthalpy h "specific enthalpy";
  input SI.SpecificEntropy s "specific entropy";
  output SI.Pressure p "Pressure";
  constant Real[:] n=
    {-0.691997014660582,-0.183612548787560e2,-0.928332409297335e1,0.659639569909906e2,
     -0.162060388912024e2,0.450620017338667e3,0.854680678224170e3,0.607523214001162e4,0.326487682621856e2,
     -0.269408844582931e2,-0.319947848334300e3,-0.928354307043320e3,0.303634537455249e2,-0.650540422444146e2,
     -0.430991316516130e4,-0.747512324096068e3,0.730000345529245e3,0.114284032569021e4,-0.436407041874559e3};
  constant Real[:] I = {0,0,0,0,0,0,0,0,1,1,1,1,2,2,2,3,4,4,5};
  constant Real[:] J = {0,1,2,4,5,6,8,14,0,1,4,6,0,1,10,4,1,4,0};
  constant SI.SpecificEnthalpy hstar = 3400e3 "normalization enthalpy";
  constant SI.Pressure pstar = 100e6 "normalization pressure";
  constant SI.SpecificEntropy sstar = 7.6e3 "normalization entropy";
protected 
  Real eta = h/hstar "normalized specific enthalpy";
  Real sigma = s/sstar "normalized specific entropy";
algorithm 
  p := sum(n[i]*(eta + 0.05)^I[i]*(sigma + 0.05)^J[i] for i in 1:19)*pstar;
end p1_hs;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.h2ab_s Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.h2ab_s

boundary between regions 2a and 2b

Information


  

Equation number 2 from:
The International Association for the Properties of Water and Steam
Gaithersburg, Maryland, USA
September 2001
Supplementary Release on  Backward Equations for Pressure as a Function of Enthalpy and Entropy p(h,s) to the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of Water and Steam

  

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
SpecificEntropys Entropy [J/(kg.K)]

Outputs

TypeNameDescription
SpecificEnthalpyhEnthalpy [J/kg]

Modelica definition

function h2ab_s "boundary between regions 2a and 2b"
  extends Modelica.Icons.Function;
  output SI.SpecificEnthalpy h "Enthalpy";
  input SI.SpecificEntropy s "Entropy";
protected 
  constant Real[:] n = {-0.349898083432139e4,0.257560716905876e4,-0.421073558227969e3,0.276349063799944e2};
  constant SI.SpecificEnthalpy hstar = 1e3 "normalization enthalpy";
  constant SI.SpecificEntropy sstar = 1e3 "normalization entropy";
  Real sigma = s/sstar "normalized specific entropy";
algorithm 
  h := (n[1] + n[2]*sigma + n[3]*sigma^2 + n[4]*sigma^3)*hstar;
end h2ab_s;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p2a_hs Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p2a_hs

pressure as a function of enthalpy and entropy in subregion 2a

Information


  

Equation number 3 from:
The International Association for the Properties of Water and Steam
Gaithersburg, Maryland, USA
September 2001
Supplementary Release on  Backward Equations for Pressure as a Function of Enthalpy and Entropy p(h,s) to the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of Water and Steam

  

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
SpecificEnthalpyh specific enthalpy [J/kg]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
PressurepPressure [Pa]

Modelica definition

function p2a_hs 
  "pressure as a function of enthalpy and entropy in subregion 2a"
  extends Modelica.Icons.Function;
  input SI.SpecificEnthalpy h "specific enthalpy";
  input SI.SpecificEntropy s "specific entropy";
  output SI.Pressure p "Pressure";
  constant Real[:] n=
    {-0.182575361923032e-1,-0.125229548799536,0.592290437320145,0.604769706185122e1,
     0.238624965444474e3,-0.298639090222922e3,0.512250813040750e-1,-0.437266515606486,0.413336902999504,
     -0.516468254574773e1,-0.557014838445711e1,0.128555037824478e2,0.114144108953290e2,-0.119504225652714e3,
     -0.284777985961560e4,0.431757846408006e4,0.112894040802650e1,0.197409186206319e4,0.151612444706087e4,
     0.141324451421235e-1,0.585501282219601,-0.297258075863012e1,0.594567314847319e1,-0.623656565798905e4,
     0.965986235133332e4,0.681500934948134e1,-0.633207286824489e4,-0.558919224465760e1,0.400645798472063e-1};
  constant Real[:] I = {0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,2,2,2,3,3,3,3,3,4,5,5,6,7};
  constant Real[:] J = {1,3,6,16,20,22,0,1,2,3,5,6,10,16,20,22,3,16,20,0,2,3,6,16,16,3,16,3,1};
  constant SI.SpecificEnthalpy hstar = 4200e3 "normalization enthalpy";
  constant SI.Pressure pstar = 4e6 "normalization pressure";
  constant SI.SpecificEntropy sstar = 12e3 "normalization entropy";
protected 
  Real eta = h/hstar "normalized specific enthalpy";
  Real sigma = s/sstar "normalized specific entropy";
algorithm 
  p := sum(n[i]*(eta - 0.5)^I[i]*(sigma - 1.2)^J[i] for i in 1:29)^4*pstar;
end p2a_hs;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p2b_hs Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p2b_hs

pressure as a function of enthalpy and entropy in subregion 2a

Information


Equation number 4 from:
The International Association for the Properties of Water and Steam
Gaithersburg, Maryland, USA
September 2001
Supplementary Release on  Backward Equations for Pressure as a Function of Enthalpy and Entropy p(h,s) to the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
SpecificEnthalpyh specific enthalpy [J/kg]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
PressurepPressure [Pa]

Modelica definition

function p2b_hs 
  "pressure as a function of enthalpy and entropy in subregion 2a"
  extends Modelica.Icons.Function;
  input SI.SpecificEnthalpy h "specific enthalpy";
  input SI.SpecificEntropy s "specific entropy";
  output SI.Pressure p "Pressure";
  constant Real[:] n=
    {0.801496989929495e-1,-0.543862807146111,0.337455597421283,0.890555451157450e1,
     0.313840736431485e3,0.797367065977789,-0.121616973556240e1,0.872803386937477e1,-0.169769781757602e2,
     -0.186552827328416e3,0.951159274344237e5,-0.189168510120494e2,-0.433407037194840e4,0.543212633012715e9,
     0.144793408386013,0.128024559637516e3,-0.672309534071268e5,0.336972380095287e8,-0.586634196762720e3,
     -0.221403224769889e11,0.171606668708389e4,-0.570817595806302e9,-0.312109693178482e4,-0.207841384633010e7,
     0.305605946157786e13,0.322157004314333e4,0.326810259797295e12,-0.144104158934487e4,0.410694867802691e3,
     0.109077066873024e12,-0.247964654258893e14,0.188801906865134e10,-0.123651009018773e15};
  constant Real[:] I = {0,0,0,0,0,1,1,1,1,1,1,2,2,2,3,3,3,3,4,4,5,5,6,6,6,7,7,8,8,8,8,12,14};
  constant Real[:] J = {0,1,2,4,8,0,1,2,3,5,12,1,6,18,0,1,7,12,1,16,1,12,1,8,18,1,16,1,3,14,18,10,16};
  constant SI.SpecificEnthalpy hstar = 4100e3 "normalization enthalpy";
  constant SI.Pressure pstar = 100e6 "normalization pressure";
  constant SI.SpecificEntropy sstar = 7.9e3 "normalization entropy";
protected 
  Real eta = h/hstar "normalized specific enthalpy";
  Real sigma = s/sstar "normalized specific entropy";
algorithm 
  p := sum(n[i]*(eta - 0.6)^I[i]*(sigma - 1.01)^J[i] for i in 1:33)^4*pstar;
end p2b_hs;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p2c_hs Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.p2c_hs

pressure as a function of enthalpy and entropy in subregion 2c

Information


      

Equation number 5 from:
The International Association for the Properties of Water and Steam
Gaithersburg, Maryland, USA
September 2001
Supplementary Release on  Backward Equations for Pressure as a Function of Enthalpy and Entropy p(h,s) to the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of Water and Steam

      

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
SpecificEnthalpyh specific enthalpy [J/kg]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
PressurepPressure [Pa]

Modelica definition

function p2c_hs 
  "pressure as a function of enthalpy and entropy in subregion 2c"
    extends Modelica.Icons.Function;
    input SI.SpecificEnthalpy h "specific enthalpy";
    input SI.SpecificEntropy s "specific entropy";
    output SI.Pressure p "Pressure";
    constant Real[:] n=
      {0.112225607199012,-0.339005953606712e1,-0.320503911730094e2,-0.197597305104900e3,
       -0.407693861553446e3,0.132943775222331e5,0.170846839774007e1,0.373694198142245e2,0.358144365815434e4,
       0.423014446424664e6,-0.751071025760063e9,0.523446127607898e2,-0.228351290812417e3,-0.960652417056937e6,
       -0.807059292526074e8,0.162698017225669e13,0.772465073604171,0.463929973837746e5,-0.137317885134128e8,
       0.170470392630512e13,-0.251104628187308e14,0.317748830835520e14,0.538685623675312e2,-0.553089094625169e5,
       -0.102861522421405e7,0.204249418756234e13,0.273918446626977e9,-0.263963146312685e16,-0.107890854108088e10,
       -0.296492620980124e11,-0.111754907323424e16};
    constant Real[:] I = {0,0,0,0,0,0,1,1,1,1,1,2,2,2,2,2,3,3,3,3,3,4,5,5,5,5,6,6,10,12,16};
    constant Real[:] J = {0,1,2,3,4,8,0,2,5,8,14,2,3,7,10,18,0,5,8,16,18,18,1,4,6,14,8,18,7,7,10};
    constant SI.SpecificEnthalpy hstar = 3500e3 "normalization enthalpy";
    constant SI.Pressure pstar = 100e6 "normalization pressure";
    constant SI.SpecificEntropy sstar = 5.9e3 "normalization entropy";
protected 
    Real eta = h/hstar "normalized specific enthalpy";
    Real sigma = s/sstar "normalized specific entropy";
algorithm 
    p := sum(n[i]*(eta - 0.7)^I[i]*(sigma - 1.1)^J[i] for i in 1:31)^4*pstar;
end p2c_hs;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.h3ab_p Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.h3ab_p

ergion 3 a b boundary for pressure/enthalpy

Information


      

 Equation number 1 from:

 [1] The international Association for the Properties of Water and Steam
 Vejle, Denmark
 August 2003
 Supplementary Release on Backward Equations for the Fucnctions T(p,h), v(p,h) and T(p,s),
 v(p,s) for Region 3 of the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of
 Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep Pressure [Pa]

Outputs

TypeNameDescription
SpecificEnthalpyhEnthalpy [J/kg]

Modelica definition

function h3ab_p "ergion 3 a b boundary for pressure/enthalpy"
    extends Modelica.Icons.Function;
    output SI.SpecificEnthalpy h "Enthalpy";
    input SI.Pressure p "Pressure";
protected 
    constant Real[:] n = {0.201464004206875e4,0.374696550136983e1,-0.219921901054187e-1,0.875131686009950e-4};
    constant SI.SpecificEnthalpy hstar = 1000 "normalization enthalpy";
    constant SI.Pressure pstar = 1e6 "normalization pressure";
    Real pi = p/pstar "normalized specific pressure";

algorithm 
    h := (n[1] + n[2]*pi + n[3]*pi^2 + n[4]*pi^3)*hstar;
end h3ab_p;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3a_ph Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3a_ph

Region 3 a: inverse function T(p,h)

Information


 

 Equation number 2 from:

 [1] The international Association for the Properties of Water and Steam
 Vejle, Denmark
 August 2003
 Supplementary Release on Backward Equations for the Fucnctions T(p,h), v(p,h) and T(p,s),
 v(p,s) for Region 3 of the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of
 Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep Pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]

Outputs

TypeNameDescription
Temp_KTTemperature [K]

Modelica definition

function T3a_ph "Region 3 a: inverse function T(p,h)"
    extends Modelica.Icons.Function;
    input SI.Pressure p "Pressure";
    input SI.SpecificEnthalpy h "specific enthalpy";
    output SI.Temp_K T "Temperature";
protected 
    constant Real[:] n=
      {-0.133645667811215e-6,0.455912656802978e-5,-0.146294640700979e-4, 0.639341312970080e-2,0.372783927268847e3,
       -0.718654377460447e4,0.573494752103400e6,-0.267569329111439e7,-0.334066283302614e-4,-0.245479214069597e-1,
       0.478087847764996e2,0.764664131818904e-5,0.128350627676972e-2,0.171219081377331e-1,-0.851007304583213e1,
       -0.136513461629781e-1,-0.384460997596657e-5,0.337423807911655e-2,-0.551624873066791,0.729202277107470,
       -0.992522757376041e-2,-0.119308831407288,0.793929190615421,0.454270731799386,0.209998591259910,
       -0.642109823904738e-2,-0.235155868604540e-1,0.252233108341612e-2,-0.764885133368119e-2,0.136176427574291e-1,
       -0.133027883575669e-1};
    constant Real[:] I = {-12,-12,-12,-12,-12,-12,-12,-12,-10,-10,
                          -10,-8,-8,-8,-8,-5,-3,-2,-2,-2,-1,-1,0,0,1,3,3,4,4,10,12};
    constant Real[:] J = { 0,1,2,6,14,16,20,22,1,5,12,0,2,4,10,2,0,1,3,4,0,2,0,1,1,0,1,0,3,4,5};
    constant SI.SpecificEnthalpy hstar = 2300e3 "normalization enthalpy";
    constant SI.Pressure pstar = 100e6 "normalization pressure";
    constant SI.Temp_K Tstar = 760 "normalization temperature";
    Real pi = p/pstar "normalized specific pressure";
    Real eta = h/hstar "normalized specific enthalpy";
algorithm 
    T := sum(n[i]*(pi + 0.240)^I[i]*(eta - 0.615)^J[i] for i in 1:31)*Tstar;
end T3a_ph;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3b_ph Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3b_ph

Region 3 b: inverse function T(p,h)

Information


 

 Equation number 3 from:

 [1] The international Association for the Properties of Water and Steam
 Vejle, Denmark
 August 2003
 Supplementary Release on Backward Equations for the Fucnctions T(p,h), v(p,h) and T(p,s),
 v(p,s) for Region 3 of the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of
 Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep Pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]

Outputs

TypeNameDescription
Temp_KTTemperature [K]

Modelica definition

function T3b_ph "Region 3 b: inverse function T(p,h)"
    extends Modelica.Icons.Function;
    input SI.Pressure p "Pressure";
    input SI.SpecificEnthalpy h "specific enthalpy";
    output SI.Temp_K T "Temperature";
protected 
    constant Real[:] n=
      {0.323254573644920e-4,-0.127575556587181e-3,-0.475851877356068e-3,0.156183014181602e-2,
       0.105724860113781,-0.858514221132534e2,0.724140095480911e3,0.296475810273257e-2,-0.592721983365988e-2,
       -0.126305422818666e-1,-0.115716196364853,0.849000969739595e2,-0.108602260086615e-1,0.154304475328851e-1,
       0.750455441524466e-1,0.252520973612982e-1,-0.602507901232996e-1,-0.307622221350501e1,-0.574011959864879e-1,
       0.503471360939849e1,-0.925081888584834,0.391733882917546e1,-0.773146007130190e2,0.949308762098587e4,
       -0.141043719679409e7,0.849166230819026e7,0.861095729446704,0.323346442811720,0.873281936020439,
       -0.436653048526683,0.286596714529479,-0.131778331276228,0.676682064330275e-2};
    constant Real[:] I = {-12,-12,-10,-10,-10,-10,-10,-8,-8,-8,-8,
                          -8,-6,-6,-6,-4,-4,-3,-2,-2,-1,-1,-1,-1,-1,-1,0,0,1,3,5,6,8};
    constant Real[:] J = {0,1,0,1,5,10,12,0,1,2,4,10,0,1,2,0,1,5,0,4,2,4,6,10,14,16,0,2,1,1,1,1,1};
    constant SI.Temp_K Tstar = 860 "normalization temperature";
    constant SI.Pressure pstar = 100e6 "normalization pressure";
    constant SI.SpecificEnthalpy hstar = 2800e3 "normalization enthalpy";
    Real pi = p/pstar "normalized specific pressure";
    Real eta = h/hstar "normalized specific enthalpy";
algorithm 
    T := sum(n[i]*(pi + 0.298)^I[i]*(eta - 0.720)^J[i] for i in 1:33)*Tstar;
end T3b_ph;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3a_ph Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3a_ph

Region 3 a: inverse function v(p,h)

Information


 

 Equation number 4 from:

 [1] The international Association for the Properties of Water and Steam
 Vejle, Denmark
 August 2003
 Supplementary Release on Backward Equations for the Fucnctions T(p,h), v(p,h) and T(p,s),
 v(p,s) for Region 3 of the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of
 Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep Pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]

Outputs

TypeNameDescription
SpecificVolumevspecific volume [m3/kg]

Modelica definition

function v3a_ph "Region 3 a: inverse function v(p,h)"
    extends Modelica.Icons.Function;
    input SI.Pressure p "Pressure";
    input SI.SpecificEnthalpy h "specific enthalpy";
    output SI.SpecificVolume v "specific volume";
protected 
    constant Real[:] n=
      { 0.529944062966028e-2,-0.170099690234461,0.111323814312927e2,-0.217898123145125e4,
       -0.506061827980875e-3,0.556495239685324,-0.943672726094016e1,-0.297856807561527,0.939353943717186e2,
       0.192944939465981e-1,0.421740664704763,-0.368914126282330e7,-0.737566847600639e-2,-0.354753242424366,
       -0.199768169338727e1,0.115456297059049e1,0.568366875815960e4,0.808169540124668e-2,0.172416341519307,
       0.104270175292927e1,-0.297691372792847,0.560394465163593,0.275234661176914,-0.148347894866012,
       -0.651142513478515e-1,-0.292468715386302e1,0.664876096952665e-1,0.352335014263844e1,-0.146340792313332e-1,
       -0.224503486668184e1,0.110533464706142e1,-0.408757344495612e-1};
    constant Real[:] I = {-12,-12,-12,-12,-10,-10,-10,-8,-8,-6,
                          -6,-6,-4,-4,-3,-2,-2,-1,-1,-1,-1,0,0,1,1,1,2,2,3,4,5,8};
    constant Real[:] J = {6,8,12,18,4,7,10,5,12,3,4,22,2,3,7,3,16,0,1,2,3,0,1,0,1,2,0,2,0,2,2,2};
    constant SI.Volume vstar = 0.0028 "normalization temperature";
    constant SI.Pressure pstar = 100e6 "normalization pressure";
    constant SI.SpecificEnthalpy hstar = 2100e3 "normalization enthalpy";
    Real pi = p/pstar "normalized specific pressure";
    Real eta = h/hstar "normalized specific enthalpy";
algorithm 
    v := sum(n[i]*(pi + 0.128)^I[i]*(eta - 0.727)^J[i] for i in 1:32)*vstar;
end v3a_ph;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3b_ph Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3b_ph

Region 3 b: inverse function v(p,h)

Information


 

 Equation number 5 from:

 [1] The international Association for the Properties of Water and Steam
 Vejle, Denmark
 August 2003
 Supplementary Release on Backward Equations for the Fucnctions T(p,h), v(p,h) and T(p,s),
 v(p,s) for Region 3 of the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of
 Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep Pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]

Outputs

TypeNameDescription
SpecificVolumevspecific volume [m3/kg]

Modelica definition

function v3b_ph "Region 3 b: inverse function v(p,h)"
    extends Modelica.Icons.Function;
    input SI.Pressure p "Pressure";
    input SI.SpecificEnthalpy h "specific enthalpy";
    output SI.SpecificVolume v "specific volume";
protected 
    constant Real[:] n=
      { -0.225196934336318e-8,0.140674363313486e-7,0.233784085280560e-5,-0.331833715229001e-4,
       0.107956778514318e-2,-0.271382067378863,0.107202262490333e1,-0.853821329075382,-0.215214194340526e-4,
       0.769656088222730e-3,-0.431136580433864e-2,0.453342167309331,-0.507749535873652,-0.100475154528389e3,
       -0.219201924648793,-0.321087965668917e1,0.607567815637771e3,0.557686450685932e-3,0.187499040029550,
       0.905368030448107e-2,0.285417173048685,0.329924030996098e-1,0.239897419685483,0.482754995951394e1,
       -0.118035753702231e2,0.169490044091791,-0.179967222507787e-1,0.371810116332674e-1,-0.536288335065096e-1,
       0.160697101092520e1};
    constant Real[:] I = {-12,-12,-8,-8,-8,-8,-8,-8,-6,-6,
                          -6,-6,-6,-6,-4,-4,-4,-3,-3,-2,-2,-1,-1,-1,-1,0,1,1,2,2};
    constant Real[:] J = {0,1,0,1,3,6,7,8,0,1,2,5,6,10,3,6,10,0,2,1,2,0,1,4,5,0,0,1,2,6};
    constant SI.Volume vstar = 0.0088 "normalization temperature";
    constant SI.Pressure pstar = 100e6 "normalization pressure";
    constant SI.SpecificEnthalpy hstar = 2800e3 "normalization enthalpy";
    Real pi = p/pstar "normalized specific pressure";
    Real eta = h/hstar "normalized specific enthalpy";
algorithm 
    v := sum(n[i]*(pi + 0.0661)^I[i]*(eta - 0.720)^J[i] for i in 1:30)*vstar;
end v3b_ph;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3a_ps Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3a_ps

Region 3 a: inverse function T(p,s)

Information


 

 Equation number 6 from:

 [1] The international Association for the Properties of Water and Steam
 Vejle, Denmark
 August 2003
 Supplementary Release on Backward Equations for the Fucnctions T(p,h), v(p,h) and T(p,s),
 v(p,s) for Region 3 of the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of
 Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep Pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
Temp_KTTemperature [K]

Modelica definition

function T3a_ps "Region 3 a: inverse function T(p,s)"
    extends Modelica.Icons.Function;
    input SI.Pressure p "Pressure";
    input SI.SpecificEntropy s "specific entropy";
    output SI.Temp_K T "Temperature";
protected 
    constant Real[:] n=
      {0.150042008263875e10,-0.159397258480424e12,0.502181140217975e-3,-0.672057767855466e2,
       0.145058545404456e4,-0.823889534888890e4,-0.154852214233853,0.112305046746695e2,-0.297000213482822e2,
       0.438565132635495e11,0.137837838635464e-2,-0.297478527157462e1,0.971777947349413e13,-0.571527767052398e-4,
       0.288307949778420e5,-0.744428289262703e14,0.128017324848921e2,-0.368275545889071e3,0.664768904779177e16,
       0.449359251958880e-1,-0.422897836099655e1,-0.240614376434179,-0.474341365254924e1,0.724093999126110,
       0.923874349695897,0.399043655281015e1,0.384066651868009e-1,-0.359344365571848e-2,-0.735196448821653,
       0.188367048396131,0.141064266818704e-3,-0.257418501496337e-2,0.123220024851555e-2};
    constant Real[:] I = {-12,-12,-10,-10,-10,-10,-8,-8,
                          -8,-8,-6,-6,-6,-5,-5,-5,-4,-4,-4,-2,-2,-1,-1,0,0,0,1,2,2,3,8,8,10};
    constant Real[:] J = {28,32,4,10,12,14,5,7,8,28,2,6,32,0,14,32,6,10,36,1,4,1,6,0,1,4,0,0,3,2,0,1,2};
    constant SI.Temp_K Tstar = 760 "normalization temperature";
    constant SI.Pressure pstar = 100e6 "normalization pressure";
    constant SI.SpecificEntropy sstar = 4.4e3 "normalization entropy";
    Real pi = p/pstar "normalized specific pressure";
    Real sigma = s/sstar "normalized specific entropy";
algorithm 
    T := sum(n[i]*(pi + 0.240)^I[i]*(sigma - 0.703)^J[i] for i in 1:33)*Tstar;
end T3a_ps;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3b_ps Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.T3b_ps

Region 3 b: inverse function T(p,s)

Information


 

 Equation number 7 from:

 [1] The international Association for the Properties of Water and Steam
 Vejle, Denmark
 August 2003
 Supplementary Release on Backward Equations for the Fucnctions T(p,h), v(p,h) and T(p,s),
 v(p,s) for Region 3 of the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of
 Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep Pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
Temp_KTTemperature [K]

Modelica definition

function T3b_ps "Region 3 b: inverse function T(p,s)"
    extends Modelica.Icons.Function;
    input SI.Pressure p "Pressure";
    input SI.SpecificEntropy s "specific entropy";
    output SI.Temp_K T "Temperature";
protected 
    constant Real[:] n=
      {0.527111701601660,-0.401317830052742e2,0.153020073134484e3,-0.224799398218827e4,
       -0.193993484669048,-0.140467557893768e1,0.426799878114024e2,0.752810643416743,0.226657238616417e2,
       -0.622873556909932e3,-0.660823667935396,0.841267087271658,-0.253717501764397e2,0.485708963532948e3,
       0.880531517490555e3,0.265015592794626e7,-0.359287150025783,-0.656991567673753e3,0.241768149185367e1,
       0.856873461222588,0.655143675313458,-0.213535213206406,0.562974957606348e-2,-0.316955725450471e15,
       -0.699997000152457e-3,0.119845803210767e-1,0.193848122022095e-4,-0.215095749182309e-4};
    constant Real[:] I = {-12,-12,-12,-12,-8,-8,-8,-6,-6,-6,-5,-5,-5,-5,-5,-4,-3,-3,-2,0,2,3,4,5,6,8,12,14};
    constant Real[:] J = {1,3,4,7,0,1,3,0,2,4,0,1,2,4,6,12,1,6,2,0,1,1,0,24,0,3,1,2};
    constant SI.Temp_K Tstar = 860 "normalization temperature";
    constant SI.Pressure pstar = 100e6 "normalization pressure";
    constant SI.SpecificEntropy sstar = 5.3e3 "normalization entropy";
    Real pi = p/pstar "normalized specific pressure";
    Real sigma = s/sstar "normalized specific entropy";
algorithm 
    T := sum(n[i]*(pi + 0.760)^I[i]*(sigma - 0.818)^J[i] for i in 1:28)*Tstar;
end T3b_ps;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3a_ps Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3a_ps

Region 3 a: inverse function v(p,s)

Information


 

 Equation number 8 from:

 [1] The international Association for the Properties of Water and Steam
 Vejle, Denmark
 August 2003
 Supplementary Release on Backward Equations for the Fucnctions T(p,h), v(p,h) and T(p,s),
 v(p,s) for Region 3 of the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of
 Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep Pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
SpecificVolumevspecific volume [m3/kg]

Modelica definition

function v3a_ps "Region 3 a: inverse function v(p,s)"
    extends Modelica.Icons.Function;
    input SI.Pressure p "Pressure";
    input SI.SpecificEntropy s "specific entropy";
    output SI.SpecificVolume v "specific volume";
protected 
    constant Real[:] n=
      {0.795544074093975e2,-0.238261242984590e4,0.176813100617787e5,-0.110524727080379e-2,
       -0.153213833655326e2,0.297544599376982e3,-0.350315206871242e8,0.277513761062119,-0.523964271036888,
       -0.148011182995403e6,0.160014899374266e7,0.170802322663427e13,0.246866996006494e-3,0.165326084797980e1,
       -0.118008384666987,0.253798642355900e1,0.965127704669424,-0.282172420532826e2,0.203224612353823,
       0.110648186063513e1,0.526127948451280,0.277000018736321,0.108153340501132e1,-0.744127885357893e-1,
       0.164094443541384e-1,-0.680468275301065e-1,0.257988576101640e-1,-0.145749861944416e-3};
    constant Real[:] I = {-12,-12,-12,-10,-10,-10,-10,-8,-8,-8,-8,-6,-5,-4,-3,-3,-2,-2,-1,-1,0,0,0,1,2,4,5,6};
    constant Real[:] J = {10,12,14,4,8,10,20,5,6,14,16,28,1,5,2,4,3,8,1,2,0,1,3,0,0,2,2,0};
    constant SI.Volume vstar = 0.0028 "normalization temperature";
    constant SI.Pressure pstar = 100e6 "normalization pressure";
    constant SI.SpecificEntropy sstar = 4.4e3 "normalization entropy";
    Real pi = p/pstar "normalized specific pressure";
    Real sigma = s/sstar "normalized specific entropy";
algorithm 
    v := sum(n[i]*(pi + 0.187)^I[i]*(sigma - 0.755)^J[i] for i in 1:28)*vstar;
end v3a_ps;

Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3b_ps Modelica.Media.Water.IF97_Utilities.BaseIF97.Basic.v3b_ps

Region 3 b: inverse function v(p,s)

Information


 

 Equation number 9 from:

 [1] The international Association for the Properties of Water and Steam
 Vejle, Denmark
 August 2003
 Supplementary Release on Backward Equations for the Fucnctions T(p,h), v(p,h) and T(p,s),
 v(p,s) for Region 3 of the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of
 Water and Steam

Extends from Modelica.Icons.Function (Icon for functions).

Inputs

TypeNameDefaultDescription
Pressurep Pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]

Outputs

TypeNameDescription
SpecificVolumevspecific volume [m3/kg]

Modelica definition

function v3b_ps "Region 3 b: inverse function v(p,s)"
    extends Modelica.Icons.Function;
    input SI.Pressure p "Pressure";
    input SI.SpecificEntropy s "specific entropy";
    output SI.SpecificVolume v "specific volume";
protected 
    constant Real[:] n=
      {0.591599780322238e-4,-0.185465997137856e-2,0.104190510480013e-1,0.598647302038590e-2,
       -0.771391189901699,0.172549765557036e1,-0.467076079846526e-3,0.134533823384439e-1,-0.808094336805495e-1,
       0.508139374365767,0.128584643361683e-2,-0.163899353915435e1,0.586938199318063e1,-0.292466667918613e1,
       -0.614076301499537e-2,0.576199014049172e1,-0.121613320606788e2,0.167637540957944e1,-0.744135838773463e1,
       0.378168091437659e-1,0.401432203027688e1,0.160279837479185e2,0.317848779347728e1,-0.358362310304853e1,
       -0.115995260446827e7,0.199256573577909,-0.122270624794624,-0.191449143716586e2,-0.150448002905284e-1,
       0.146407900162154e2,-0.327477787188230e1};
    constant Real[:] I = {-12,-12,-12,-12,-12,-12,-10,-10,
                          -10,-10,-8,-5,-5,-5,-4,-4,-4,-4,-3,-2,-2,-2,-2,-2,-2,0,0,0,1,1,2};
    constant Real[:] J = {0,1,2,3,5,6,0,1,2,4,0,1,2,3,0,1,2,3,1,0,1,2,3,4,12,0,1,2,0,2,2};
    constant SI.Volume vstar = 0.0088 "normalization temperature";
    constant SI.Pressure pstar = 100e6 "normalization pressure";
    constant SI.SpecificEntropy sstar = 5.3e3 "normalization entropy";
    Real pi = p/pstar "normalized specific pressure";
    Real sigma = s/sstar "normalized specific entropy";
algorithm 
    v := sum(n[i]*(pi + 0.298)^I[i]*(sigma - 0.816)^J[i] for i in 1:31)*vstar;
end v3b_ps;

Automatically generated Fri Nov 12 16:31:39 2010.