Modelica.Media.Water.IF97_Utilities

Low level and utility computation for high accuracy water properties according to the IAPWS/IF97 standard

Information


      

Package description:

This package provides high accuracy physical properties for water according to the IAPWS/IF97 standard. It has been part of the ThermoFluid Modelica library and been extended, reorganized and documented to become part of the Modelica Standard library.

An important feature that distinguishes this implementation of the IF97 steam property standard is that this implementation has been explicitly designed to work well in dynamic simulations. Computational performance has been of high importance. This means that there often exist several ways to get the same result from different functions if one of the functions is called often but can be optimized for that purpose.

The original documentation of the IAPWS/IF97 steam properties can freely be distributed with computer implementations, so for curious minds the complete standard documentation is provided with the Modelica properties library. The following documents are included (in directory Modelica/Resources/Documentation/Media/Water/IF97documentation):

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 BaseIF97 Modelica Physical Property Model: the new industrial formulation IAPWS-IF97
Modelica.Media.Water.IF97_Utilities.iter iter  
Modelica.Media.Water.IF97_Utilities.waterBaseProp_ph waterBaseProp_ph intermediate property record for water
Modelica.Media.Water.IF97_Utilities.waterBaseProp_ps waterBaseProp_ps intermediate property record for water
Modelica.Media.Water.IF97_Utilities.rho_props_ps rho_props_ps density as function of pressure and specific entropy
Modelica.Media.Water.IF97_Utilities.rho_ps rho_ps density as function of pressure and specific entropy
Modelica.Media.Water.IF97_Utilities.T_props_ps T_props_ps temperature as function of pressure and specific entropy
Modelica.Media.Water.IF97_Utilities.T_ps T_ps temperature as function of pressure and specific entropy
Modelica.Media.Water.IF97_Utilities.h_props_ps h_props_ps specific enthalpy as function or pressure and temperature
Modelica.Media.Water.IF97_Utilities.h_ps h_ps specific enthalpy as function or pressure and temperature
Modelica.Media.Water.IF97_Utilities.phase_ps phase_ps phase as a function of pressure and specific entropy
Modelica.Media.Water.IF97_Utilities.phase_ph phase_ph phase as a function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.phase_dT phase_dT phase as a function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.rho_props_ph rho_props_ph density as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.rho_ph rho_ph density as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.rho_ph_der rho_ph_der derivative function of rho_ph
Modelica.Media.Water.IF97_Utilities.T_props_ph T_props_ph temperature as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.T_ph T_ph temperature as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.T_ph_der T_ph_der derivative function of T_ph
Modelica.Media.Water.IF97_Utilities.s_props_ph s_props_ph specific entropy as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.s_ph s_ph specific entropy as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.s_ph_der s_ph_der specific entropy as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.cv_props_ph cv_props_ph specific heat capacity at constant volume as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.cv_ph cv_ph specific heat capacity at constant volume as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.regionAssertReal regionAssertReal assert function for inlining
Modelica.Media.Water.IF97_Utilities.cp_props_ph cp_props_ph specific heat capacity at constant pressure as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.cp_ph cp_ph specific heat capacity at constant pressure as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.beta_props_ph beta_props_ph isobaric expansion coefficient as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.beta_ph beta_ph isobaric expansion coefficient as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.kappa_props_ph kappa_props_ph isothermal compressibility factor as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.kappa_ph kappa_ph isothermal compressibility factor as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.velocityOfSound_props_ph velocityOfSound_props_ph speed of sound as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.velocityOfSound_ph velocityOfSound_ph  
Modelica.Media.Water.IF97_Utilities.isentropicExponent_props_ph isentropicExponent_props_ph isentropic exponent as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.isentropicExponent_ph isentropicExponent_ph isentropic exponent as function of pressure and specific enthalpy
Modelica.Media.Water.IF97_Utilities.ddph_props ddph_props density derivative by pressure
Modelica.Media.Water.IF97_Utilities.ddph ddph density derivative by pressure
Modelica.Media.Water.IF97_Utilities.ddhp_props ddhp_props density derivative by specific enthalpy
Modelica.Media.Water.IF97_Utilities.ddhp ddhp density derivative by specific enthalpy
Modelica.Media.Water.IF97_Utilities.waterBaseProp_pT waterBaseProp_pT intermediate property record for water (p and T prefered states)
Modelica.Media.Water.IF97_Utilities.rho_props_pT rho_props_pT density as function or pressure and temperature
Modelica.Media.Water.IF97_Utilities.rho_pT rho_pT density as function or pressure and temperature
Modelica.Media.Water.IF97_Utilities.h_props_pT h_props_pT specific enthalpy as function or pressure and temperature
Modelica.Media.Water.IF97_Utilities.h_pT h_pT specific enthalpy as function or pressure and temperature
Modelica.Media.Water.IF97_Utilities.h_pT_der h_pT_der derivative function of h_pT
Modelica.Media.Water.IF97_Utilities.rho_pT_der rho_pT_der derivative function of rho_pT
Modelica.Media.Water.IF97_Utilities.s_props_pT s_props_pT specific entropy as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.s_pT s_pT temperature as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.cv_props_pT cv_props_pT specific heat capacity at constant volume as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.cv_pT cv_pT specific heat capacity at constant volume as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.cp_props_pT cp_props_pT specific heat capacity at constant pressure as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.cp_pT cp_pT specific heat capacity at constant pressure as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.beta_props_pT beta_props_pT isobaric expansion coefficient as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.beta_pT beta_pT isobaric expansion coefficient as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.kappa_props_pT kappa_props_pT isothermal compressibility factor as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.kappa_pT kappa_pT isothermal compressibility factor as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.velocityOfSound_props_pT velocityOfSound_props_pT speed of sound as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.velocityOfSound_pT velocityOfSound_pT speed of sound as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.isentropicExponent_props_pT isentropicExponent_props_pT isentropic exponent as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.isentropicExponent_pT isentropicExponent_pT isentropic exponent as function of pressure and temperature
Modelica.Media.Water.IF97_Utilities.waterBaseProp_dT waterBaseProp_dT intermediate property record for water (d and T prefered states)
Modelica.Media.Water.IF97_Utilities.h_props_dT h_props_dT specific enthalpy as function of density and temperature
Modelica.Media.Water.IF97_Utilities.h_dT h_dT specific enthalpy as function of density and temperature
Modelica.Media.Water.IF97_Utilities.h_dT_der h_dT_der derivative function of h_dT
Modelica.Media.Water.IF97_Utilities.p_props_dT p_props_dT pressure as function of density and temperature
Modelica.Media.Water.IF97_Utilities.p_dT p_dT pressure as function of density and temperature
Modelica.Media.Water.IF97_Utilities.p_dT_der p_dT_der derivative function of p_dT
Modelica.Media.Water.IF97_Utilities.s_props_dT s_props_dT specific entropy as function of density and temperature
Modelica.Media.Water.IF97_Utilities.s_dT s_dT temperature as function of density and temperature
Modelica.Media.Water.IF97_Utilities.cv_props_dT cv_props_dT specific heat capacity at constant volume as function of density and temperature
Modelica.Media.Water.IF97_Utilities.cv_dT cv_dT specific heat capacity at constant volume as function of density and temperature
Modelica.Media.Water.IF97_Utilities.cp_props_dT cp_props_dT specific heat capacity at constant pressure as function of density and temperature
Modelica.Media.Water.IF97_Utilities.cp_dT cp_dT specific heat capacity at constant pressure as function of density and temperature
Modelica.Media.Water.IF97_Utilities.beta_props_dT beta_props_dT isobaric expansion coefficient as function of density and temperature
Modelica.Media.Water.IF97_Utilities.beta_dT beta_dT isobaric expansion coefficient as function of density and temperature
Modelica.Media.Water.IF97_Utilities.kappa_props_dT kappa_props_dT isothermal compressibility factor as function of density and temperature
Modelica.Media.Water.IF97_Utilities.kappa_dT kappa_dT isothermal compressibility factor as function of density and temperature
Modelica.Media.Water.IF97_Utilities.velocityOfSound_props_dT velocityOfSound_props_dT speed of sound as function of density and temperature
Modelica.Media.Water.IF97_Utilities.velocityOfSound_dT velocityOfSound_dT speed of sound as function of density and temperature
Modelica.Media.Water.IF97_Utilities.isentropicExponent_props_dT isentropicExponent_props_dT isentropic exponent as function of density and temperature
Modelica.Media.Water.IF97_Utilities.isentropicExponent_dT isentropicExponent_dT isentropic exponent as function of density and temperature
Modelica.Media.Water.IF97_Utilities.ThermoFluidSpecial ThermoFluidSpecial  
Modelica.Media.Water.IF97_Utilities.hl_p hl_p compute the saturated liquid specific h(p)
Modelica.Media.Water.IF97_Utilities.hv_p hv_p compute the saturated vapour specific h(p)
Modelica.Media.Water.IF97_Utilities.sl_p sl_p compute the saturated liquid specific s(p)
Modelica.Media.Water.IF97_Utilities.sv_p sv_p compute the saturated vapour specific s(p)
Modelica.Media.Water.IF97_Utilities.rhol_T rhol_T compute the saturated liquid d(T)
Modelica.Media.Water.IF97_Utilities.rhov_T rhov_T compute the saturated vapour d(T)
Modelica.Media.Water.IF97_Utilities.rhol_p rhol_p compute the saturated liquid d(p)
Modelica.Media.Water.IF97_Utilities.rhov_p rhov_p compute the saturated vapour d(p)
Modelica.Media.Water.IF97_Utilities.dynamicViscosity dynamicViscosity compute eta(d,T) in the one-phase region
Modelica.Media.Water.IF97_Utilities.thermalConductivity thermalConductivity compute lambda(d,T,p) in the one-phase region
Modelica.Media.Water.IF97_Utilities.surfaceTension surfaceTension compute sigma(T) at saturation T
Modelica.Media.Water.IF97_Utilities.isentropicEnthalpy isentropicEnthalpy isentropic specific enthalpy from p,s (preferably use dynamicIsentropicEnthalpy in dynamic simulation!)
Modelica.Media.Water.IF97_Utilities.isentropicEnthalpy_props isentropicEnthalpy_props  
Modelica.Media.Water.IF97_Utilities.isentropicEnthalpy_der isentropicEnthalpy_der derivative of isentropic specific enthalpy from p,s
Modelica.Media.Water.IF97_Utilities.dynamicIsentropicEnthalpy dynamicIsentropicEnthalpy isentropic specific enthalpy from p,s and good guesses of d and T


Modelica.Media.Water.IF97_Utilities.iter Modelica.Media.Water.IF97_Utilities.iter

Modelica definition

replaceable record iter = BaseIF97.IterationData;

Modelica.Media.Water.IF97_Utilities.waterBaseProp_ph Modelica.Media.Water.IF97_Utilities.waterBaseProp_ph

intermediate property record for water

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase0phase: 2 for two-phase, 1 for one phase, 0 if unknown
Integerregion0if 0, do region computation, otherwise assume the region is this input

Outputs

TypeNameDescription
IF97BaseTwoPhaseauxauxiliary record

Modelica definition

function waterBaseProp_ph "intermediate property record for water"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase =  0 
    "phase: 2 for two-phase, 1 for one phase, 0 if unknown";
  input Integer region = 0 
    "if 0, do region computation, otherwise assume the region is this input";
  output Common.IF97BaseTwoPhase aux "auxiliary record";
protected 
  Common.GibbsDerivs g 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Common.HelmholtzDerivs f 
    "dimensionless Helmholtz funcion and dervatives w.r.t. delta and tau";
  Integer error "error flag for inverse iterations";
  SI.SpecificEnthalpy h_liq "liquid specific enthalpy";
  SI.Density d_liq "liquid density";
  SI.SpecificEnthalpy h_vap "vapour specific enthalpy";
  SI.Density d_vap "vapour density";
  Common.PhaseBoundaryProperties liq "phase boundary property record";
  Common.PhaseBoundaryProperties vap "phase boundary property record";
  Common.GibbsDerivs gl 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Common.GibbsDerivs gv 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Modelica.Media.Common.HelmholtzDerivs fl 
    "dimensionless Helmholtz function and dervatives w.r.t. delta and tau";
  Modelica.Media.Common.HelmholtzDerivs fv 
    "dimensionless Helmholtz function and dervatives w.r.t. delta and tau";
  SI.Temperature t1 
    "temperature at phase boundary, using inverse from region 1";
  SI.Temperature t2 
    "temperature at phase boundary, using inverse from region 2";
algorithm 
  aux.region := if region == 0 then 
    (if phase == 2 then 4 else BaseIF97.Regions.region_ph(p=p,h= h,phase= phase)) else region;
  aux.phase := if phase <> 0 then phase else if aux.region == 4 then 2 else 1;
  aux.p := max(p,611.657);
  aux.h := max(h,1e3);
  aux.R := BaseIF97.data.RH2O;
  if (aux.region == 1) then
    aux.T := BaseIF97.Basic.tph1(aux.p, aux.h);
    g := BaseIF97.Basic.g1(p, aux.T);
    aux.s := aux.R*(g.tau*g.gtau - g.g);
    aux.rho := p/(aux.R*aux.T*g.pi*g.gpi);
    aux.vt := aux.R/p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*aux.T/(p*p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
    aux.x := 0.0;
    aux.dpT := -aux.vt/aux.vp;
  elseif (aux.region == 2) then
    aux.T := BaseIF97.Basic.tph2(aux.p, aux.h);
    g := BaseIF97.Basic.g2(p, aux.T);
    aux.s := aux.R*(g.tau*g.gtau - g.g);
    aux.rho := p/(aux.R*aux.T*g.pi*g.gpi);
    aux.vt := aux.R/p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*aux.T/(p*p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
    aux.x := 1.0;
    aux.dpT := -aux.vt/aux.vp;
  elseif (aux.region == 3) then
    (aux.rho,aux.T,error) := BaseIF97.Inverses.dtofph3(p=aux.p,h=aux.h,delp= 1.0e-7,delh=
            1.0e-6);
    f := BaseIF97.Basic.f3(aux.rho, aux.T);
    aux.h := aux.R*aux.T*(f.tau*f.ftau + f.delta*f.fdelta);
    aux.s := aux.R*(f.tau*f.ftau - f.f);
    aux.pd := aux.R*aux.T*f.delta*(2.0*f.fdelta + f.delta*f.fdeltadelta);
    aux.pt := aux.R*aux.rho*f.delta*(f.fdelta - f.tau*f.fdeltatau);
    aux.cv := abs(aux.R*(-f.tau*f.tau*f.ftautau)) 
      "can be close to neg. infinity near critical point";
    aux.cp := (aux.rho*aux.rho*aux.pd*aux.cv + aux.T*aux.pt*aux.pt)/(aux.rho*aux.rho*aux.pd);
    aux.x := 0.0;
    aux.dpT := aux.pt; /*safety against div-by-0 in initialization*/
  elseif (aux.region == 4) then
    h_liq := hl_p(p);
    h_vap := hv_p(p);
    aux.x := if (h_vap <> h_liq) then (h - h_liq)/(h_vap - h_liq) else 1.0;
    if p < BaseIF97.data.PLIMIT4A then
      t1:= BaseIF97.Basic.tph1(aux.p, h_liq);
      t2:= BaseIF97.Basic.tph2(aux.p, h_vap);
      gl := BaseIF97.Basic.g1(aux.p, t1);
      gv := BaseIF97.Basic.g2(aux.p, t2);
      liq := Common.gibbsToBoundaryProps(gl);
      vap := Common.gibbsToBoundaryProps(gv);
      aux.T := t1 + aux.x*(t2-t1);
    else
      aux.T := BaseIF97.Basic.tsat(aux.p); // how to avoid ?
      d_liq:= rhol_T(aux.T);
      d_vap:= rhov_T(aux.T);
      fl := BaseIF97.Basic.f3(d_liq, aux.T);
      fv := BaseIF97.Basic.f3(d_vap, aux.T);
      liq := Common.helmholtzToBoundaryProps(fl);
      vap := Common.helmholtzToBoundaryProps(fv);
      //  aux.dpT := BaseIF97.Basic.dptofT(aux.T);
    end if;
    aux.dpT := if (liq.d <> vap.d) then (vap.s - liq.s)*liq.d*vap.d/(liq.d - vap.d) else BaseIF97.Basic.dptofT(aux.T);
    aux.s := liq.s + aux.x*(vap.s - liq.s);
    aux.rho := liq.d*vap.d/(vap.d + aux.x*(liq.d - vap.d));
    aux.cv := Common.cv2Phase(liq, vap, aux.x, aux.T, p);
    aux.cp := liq.cp + aux.x*(vap.cp - liq.cp);
    aux.pt := liq.pt + aux.x*(vap.pt - liq.pt);
    aux.pd := liq.pd + aux.x*(vap.pd - liq.pd);
  elseif (aux.region == 5) then
    (aux.T,error) := BaseIF97.Inverses.tofph5(p=aux.p,h= aux.h,reldh= 1.0e-7);
    assert(error == 0, "error in inverse iteration of steam tables");
    g := BaseIF97.Basic.g5(aux.p, aux.T);
    aux.s := aux.R*(g.tau*g.gtau - g.g);
    aux.rho := p/(aux.R*aux.T*g.pi*g.gpi);
    aux.vt := aux.R/p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*aux.T/(p*p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
    aux.dpT := -aux.vt/aux.vp;
  else
    assert(false, "error in region computation of IF97 steam tables"
    + "(p = " + String(p) + ", h = " + String(h) + ")");
  end if;
end waterBaseProp_ph;

Modelica.Media.Water.IF97_Utilities.waterBaseProp_ps Modelica.Media.Water.IF97_Utilities.waterBaseProp_ps

intermediate property record for water

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]
Integerphase0phase: 2 for two-phase, 1 for one phase, 0 if unknown
Integerregion0if 0, do region computation, otherwise assume the region is this input

Outputs

TypeNameDescription
IF97BaseTwoPhaseauxauxiliary record

Modelica definition

function waterBaseProp_ps "intermediate property record for water"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Integer phase = 0 
    "phase: 2 for two-phase, 1 for one phase, 0 if unknown";
  input Integer region = 0 
    "if 0, do region computation, otherwise assume the region is this input";
  output Common.IF97BaseTwoPhase aux "auxiliary record";
protected 
  Common.GibbsDerivs g 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Common.HelmholtzDerivs f 
    "dimensionless Helmholtz funcion and dervatives w.r.t. delta and tau";
  Integer error "error flag for inverse iterations";
  SI.SpecificEntropy s_liq "liquid specific entropy";
  SI.Density d_liq "liquid density";
  SI.SpecificEntropy s_vap "vapour specific entropy";
  SI.Density d_vap "vapour density";
  Common.PhaseBoundaryProperties liq "phase boundary property record";
  Common.PhaseBoundaryProperties vap "phase boundary property record";
  Common.GibbsDerivs gl 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Common.GibbsDerivs gv 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Modelica.Media.Common.HelmholtzDerivs fl 
    "dimensionless Helmholtz function and dervatives w.r.t. delta and tau";
  Modelica.Media.Common.HelmholtzDerivs fv 
    "dimensionless Helmholtz function and dervatives w.r.t. delta and tau";
  SI.Temperature t1 
    "temperature at phase boundary, using inverse from region 1";
  SI.Temperature t2 
    "temperature at phase boundary, using inverse from region 2";
algorithm 
  aux.region := if region == 0 then 
    (if phase == 2 then 4 else BaseIF97.Regions.region_ps(p=p,s=s,phase=phase)) else region;
  aux.phase := if phase <> 0 then phase else if aux.region == 4 then 2 else 1;
  aux.p := p;
  aux.s := s;
  aux.R := BaseIF97.data.RH2O;
  if (aux.region == 1) then
    aux.T := BaseIF97.Basic.tps1(p, s);
    g := BaseIF97.Basic.g1(p, aux.T);
    aux.h := aux.R*aux.T*g.tau*g.gtau;
    aux.rho := p/(aux.R*aux.T*g.pi*g.gpi);
    aux.vt := aux.R/p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*aux.T/(p*p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
    aux.x := 0.0;
  elseif (aux.region == 2) then
    aux.T := BaseIF97.Basic.tps2(p, s);
    g := BaseIF97.Basic.g2(p, aux.T);
    aux.h := aux.R*aux.T*g.tau*g.gtau;
    aux.rho := p/(aux.R*aux.T*g.pi*g.gpi);
    aux.vt := aux.R/p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*aux.T/(p*p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
    aux.x := 1.0;
  elseif (aux.region == 3) then
    (aux.rho,aux.T,error) := BaseIF97.Inverses.dtofps3(p=p,s=s,delp=1.0e-7,dels=
      1.0e-6);
    f := BaseIF97.Basic.f3(aux.rho, aux.T);
    aux.h := aux.R*aux.T*(f.tau*f.ftau + f.delta*f.fdelta);
    aux.s := aux.R*(f.tau*f.ftau - f.f);
    aux.pd := aux.R*aux.T*f.delta*(2.0*f.fdelta + f.delta*f.fdeltadelta);
    aux.pt := aux.R*aux.rho*f.delta*(f.fdelta - f.tau*f.fdeltatau);
    aux.cv := aux.R*(-f.tau*f.tau*f.ftautau);
    aux.cp := (aux.rho*aux.rho*aux.pd*aux.cv + aux.T*aux.pt*aux.pt)/(aux.rho*aux.rho*aux.pd);
    aux.x := 0.0;
  elseif (aux.region == 4) then
    s_liq := BaseIF97.Regions.sl_p(p);
    s_vap := BaseIF97.Regions.sv_p(p);
    aux.x := if (s_vap <> s_liq) then (s - s_liq)/(s_vap - s_liq) else 1.0;
    if p < BaseIF97.data.PLIMIT4A then
      t1 := BaseIF97.Basic.tps1(p, s_liq);
      t2 := BaseIF97.Basic.tps2(p, s_vap);
      gl := BaseIF97.Basic.g1(p, t1);
      gv := BaseIF97.Basic.g2(p, t2);
      liq := Common.gibbsToBoundaryProps(gl);
      vap := Common.gibbsToBoundaryProps(gv);
      aux.T := t1 + aux.x*(t2 - t1);
    else
      aux.T := BaseIF97.Basic.tsat(p);
      d_liq := rhol_T(aux.T);
      d_vap := rhov_T(aux.T);
      fl := BaseIF97.Basic.f3(d_liq, aux.T);
      fv := BaseIF97.Basic.f3(d_vap, aux.T);
      liq := Common.helmholtzToBoundaryProps(fl);
      vap := Common.helmholtzToBoundaryProps(fv);
    end if;
    aux.dpT := if (liq.d <> vap.d) then (vap.s - liq.s)*liq.d*vap.d/(liq.d - vap.d) else 
         BaseIF97.Basic.dptofT(aux.T);
    aux.h := liq.h + aux.x*(vap.h - liq.h);
    aux.rho := liq.d*vap.d/(vap.d + aux.x*(liq.d - vap.d));
    aux.cv := Common.cv2Phase(liq, vap, aux.x, aux.T, p);
    aux.cp := liq.cp + aux.x*(vap.cp - liq.cp);
    aux.pt := liq.pt + aux.x*(vap.pt - liq.pt);
    aux.pd := liq.pd + aux.x*(vap.pd - liq.pd);
  elseif (aux.region == 5) then
    (aux.T,error) := BaseIF97.Inverses.tofps5(p=p,s=s,relds= 1.0e-7);
    assert(error == 0, "error in inverse iteration of steam tables");
    g := BaseIF97.Basic.g5(p, aux.T);
    aux.h := aux.R*aux.T*g.tau*g.gtau;
    aux.rho := p/(aux.R*aux.T*g.pi*g.gpi);
    aux.vt := aux.R/p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*aux.T/(p*p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
  else
    assert(false, "error in region computation of IF97 steam tables"
    + "(p = " + String(p) + ", s = " + String(s) + ")");
  end if;
end waterBaseProp_ps;

Modelica.Media.Water.IF97_Utilities.rho_props_ps Modelica.Media.Water.IF97_Utilities.rho_props_ps

density as function of pressure and specific entropy

Information

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

Inputs

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

Outputs

TypeNameDescription
Densityrhodensity [kg/m3]

Modelica definition

function rho_props_ps 
  "density as function of pressure and specific entropy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Common.IF97BaseTwoPhase properties "auxiliary record";
  output SI.Density rho "density";
algorithm 
  rho := properties.rho;
end rho_props_ps;

Modelica.Media.Water.IF97_Utilities.rho_ps Modelica.Media.Water.IF97_Utilities.rho_ps

density as function of pressure and specific entropy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Densityrhodensity [kg/m3]

Modelica definition

function rho_ps 
  "density as function of pressure and specific entropy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.Density rho "density";
algorithm 
  rho := rho_props_ps(p, s, waterBaseProp_ps(p, s, phase, region));
end rho_ps;

Modelica.Media.Water.IF97_Utilities.T_props_ps Modelica.Media.Water.IF97_Utilities.T_props_ps

temperature as function of pressure and specific entropy

Information

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

Inputs

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

Outputs

TypeNameDescription
TemperatureTtemperature [K]

Modelica definition

function T_props_ps 
  "temperature as function of pressure and specific entropy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Common.IF97BaseTwoPhase properties "auxiliary record";
  output SI.Temperature T "temperature";
algorithm 
  T := properties.T;
end T_props_ps;

Modelica.Media.Water.IF97_Utilities.T_ps Modelica.Media.Water.IF97_Utilities.T_ps

temperature as function of pressure and specific entropy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
TemperatureTTemperature [K]

Modelica definition

function T_ps 
  "temperature as function of pressure and specific entropy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.Temperature T "Temperature";
algorithm 
  T := T_props_ps(p, s, waterBaseProp_ps(p, s, phase, region));
end T_ps;

Modelica.Media.Water.IF97_Utilities.h_props_ps Modelica.Media.Water.IF97_Utilities.h_props_ps

specific enthalpy as function or pressure and temperature

Information

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

Inputs

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

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function h_props_ps 
  "specific enthalpy as function or pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificEnthalpy h "specific enthalpy";
algorithm 
  h := aux.h;
end h_props_ps;

Modelica.Media.Water.IF97_Utilities.h_ps Modelica.Media.Water.IF97_Utilities.h_ps

specific enthalpy as function or pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function h_ps 
  "specific enthalpy as function or pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificEnthalpy h "specific enthalpy";
algorithm 
  h := h_props_ps(p, s, waterBaseProp_ps(p, s, phase, region));
end h_ps;

Modelica.Media.Water.IF97_Utilities.phase_ps Modelica.Media.Water.IF97_Utilities.phase_ps

phase as a function of pressure and specific entropy

Information

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

Inputs

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

Outputs

TypeNameDescription
Integerphasetrue if in liquid or gas or supercritical region

Modelica definition

function phase_ps 
  "phase as a function of  pressure and specific entropy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  output Integer phase "true if in liquid or gas or supercritical region";
algorithm 
  phase := if ((s < sl_p(p) or s > sv_p(p)) or p > BaseIF97.data.PCRIT) then 1 else 2;
end phase_ps;

Modelica.Media.Water.IF97_Utilities.phase_ph Modelica.Media.Water.IF97_Utilities.phase_ph

phase as a function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
Integerphasetrue if in liquid or gas or supercritical region

Modelica definition

function phase_ph 
  "phase as a function of  pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  output Integer phase "true if in liquid or gas or supercritical region";
algorithm 
  phase := if ((h < hl_p(p) or h > hv_p(p)) or p > BaseIF97.data.PCRIT) then 1 else 2;
end phase_ph;

Modelica.Media.Water.IF97_Utilities.phase_dT Modelica.Media.Water.IF97_Utilities.phase_dT

phase as a function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityrho density [kg/m3]
TemperatureT temperature [K]

Outputs

TypeNameDescription
Integerphasetrue if in liquid or gas or supercritical region

Modelica definition

function phase_dT "phase as a function of  pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Density rho "density";
  input SI.Temperature T "temperature";
  output Integer phase "true if in liquid or gas or supercritical region";
algorithm 
  phase := if not ((rho < rhol_T(T) and rho > rhov_T(T)) and T < BaseIF97.
    data.TCRIT) then 1 else 2;
end phase_dT;

Modelica.Media.Water.IF97_Utilities.rho_props_ph Modelica.Media.Water.IF97_Utilities.rho_props_ph

density as function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
Densityrhodensity [kg/m3]

Modelica definition

function rho_props_ph 
  "density as function of pressure and specific enthalpy"
  annotation(derivative=rho_ph_der);
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase properties "auxiliary record";
  output SI.Density rho "density";
algorithm 
  rho := properties.rho;
end rho_props_ph;

Modelica.Media.Water.IF97_Utilities.rho_ph Modelica.Media.Water.IF97_Utilities.rho_ph

density as function of pressure and specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Densityrhodensity [kg/m3]

Modelica definition

function rho_ph 
  "density as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.Density rho "density";
algorithm 
  rho := rho_props_ph(p, h, waterBaseProp_ph(p, h, phase, region));
end rho_ph;

Modelica.Media.Water.IF97_Utilities.rho_ph_der Modelica.Media.Water.IF97_Utilities.rho_ph_der

derivative function of rho_ph

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
IF97BaseTwoPhaseaux auxiliary record
Realp_der derivative of pressure
Realh_der derivative of specific enthalpy

Outputs

TypeNameDescription
Realrho_derderivative of density

Modelica definition

function rho_ph_der "derivative function of rho_ph"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  input Real p_der "derivative of pressure";
  input Real h_der "derivative of specific enthalpy";
  output Real rho_der "derivative of density";
algorithm 
  if (aux.region == 4) then
    rho_der := (aux.rho*(aux.rho*aux.cv/aux.dpT + 1.0)/(aux.dpT*aux.T))*p_der
       + (-aux.rho*aux.rho/(aux.dpT*aux.T))*h_der;
  elseif (aux.region == 3) then
    rho_der := ((aux.rho*(aux.cv*aux.rho + aux.pt))/(aux.rho*aux.rho*aux.pd*
      aux.cv + aux.T*aux.pt*aux.pt))*p_der + (-aux.rho*aux.rho*aux.pt/(aux.
      rho*aux.rho*aux.pd*aux.cv + aux.T*aux.pt*aux.pt))*h_der;
  else
    //regions 1,2,5
    rho_der := (-aux.rho*aux.rho*(aux.vp*aux.cp - aux.vt/aux.rho + aux.T*aux.
      vt*aux.vt)/aux.cp)*p_der + (-aux.rho*aux.rho*aux.vt/(aux.cp))*h_der;
  end if;
end rho_ph_der;

Modelica.Media.Water.IF97_Utilities.T_props_ph Modelica.Media.Water.IF97_Utilities.T_props_ph

temperature as function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
TemperatureTtemperature [K]

Modelica definition

function T_props_ph 
  "temperature as function of pressure and specific enthalpy"
  annotation(derivative=T_ph_der);
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase properties "auxiliary record";
  output SI.Temperature T "temperature";
algorithm 
  T := properties.T;
end T_props_ph;

Modelica.Media.Water.IF97_Utilities.T_ph Modelica.Media.Water.IF97_Utilities.T_ph

temperature as function of pressure and specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
TemperatureTTemperature [K]

Modelica definition

function T_ph 
  "temperature as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.Temperature T "Temperature";
algorithm 
  T := T_props_ph(p, h, waterBaseProp_ph(p, h, phase, region));
end T_ph;

Modelica.Media.Water.IF97_Utilities.T_ph_der Modelica.Media.Water.IF97_Utilities.T_ph_der

derivative function of T_ph

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
IF97BaseTwoPhaseaux auxiliary record
Realp_der derivative of pressure
Realh_der derivative of specific enthalpy

Outputs

TypeNameDescription
RealT_derderivative of temperature

Modelica definition

function T_ph_der "derivative function of T_ph"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  input Real p_der "derivative of pressure";
  input Real h_der "derivative of specific enthalpy";
  output Real T_der "derivative of temperature";
algorithm 
  if (aux.region == 4) then
    T_der := 1/aux.dpT*p_der;
  elseif (aux.region == 3) then
    T_der := ((-aux.rho*aux.pd + aux.T*aux.pt)/(aux.rho*aux.rho*aux.pd*aux.cv
       + aux.T*aux.pt*aux.pt))*p_der + ((aux.rho*aux.rho*aux.pd)/(aux.rho*aux.
       rho*aux.pd*aux.cv + aux.T*aux.pt*aux.pt))*h_der;
  else
    //regions 1,2 or 5
    T_der := ((-1/aux.rho + aux.T*aux.vt)/aux.cp)*p_der + (1/aux.cp)*h_der;
  end if;
end T_ph_der;

Modelica.Media.Water.IF97_Utilities.s_props_ph Modelica.Media.Water.IF97_Utilities.s_props_ph

specific entropy as function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
SpecificEntropysspecific entropy [J/(kg.K)]

Modelica definition

function s_props_ph 
  "specific entropy as function of pressure and specific enthalpy"
  annotation(derivative=s_ph_der);
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase properties "auxiliary record";
  output SI.SpecificEntropy s "specific entropy";
algorithm 
  s := properties.s;
end s_props_ph;

Modelica.Media.Water.IF97_Utilities.s_ph Modelica.Media.Water.IF97_Utilities.s_ph

specific entropy as function of pressure and specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificEntropysspecific entropy [J/(kg.K)]

Modelica definition

function s_ph 
  "specific entropy as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase =   0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificEntropy s "specific entropy";
algorithm 
  s := s_props_ph(p, h, waterBaseProp_ph(p, h, phase, region));
end s_ph;

Modelica.Media.Water.IF97_Utilities.s_ph_der Modelica.Media.Water.IF97_Utilities.s_ph_der

specific entropy as function of pressure and specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
IF97BaseTwoPhaseaux auxiliary record
Realp_der derivative of pressure
Realh_der derivative of specific enthalpy

Outputs

TypeNameDescription
Reals_derderivative of entropy

Modelica definition

function s_ph_der 
  "specific entropy as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  input Real p_der "derivative of pressure";
  input Real h_der "derivative of specific enthalpy";
  output Real s_der "derivative of entropy";
algorithm 
  s_der := -1/(aux.rho*aux.T)*p_der + 1/aux.T*h_der;
end s_ph_der;

Modelica.Media.Water.IF97_Utilities.cv_props_ph Modelica.Media.Water.IF97_Utilities.cv_props_ph

specific heat capacity at constant volume as function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
SpecificHeatCapacitycvspecific heat capacity [J/(kg.K)]

Modelica definition

function cv_props_ph 
  "specific heat capacity at constant volume as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificHeatCapacity cv "specific heat capacity";
algorithm 
  cv := aux.cv;
end cv_props_ph;

Modelica.Media.Water.IF97_Utilities.cv_ph Modelica.Media.Water.IF97_Utilities.cv_ph

specific heat capacity at constant volume as function of pressure and specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificHeatCapacitycvspecific heat capacity [J/(kg.K)]

Modelica definition

function cv_ph 
  "specific heat capacity at constant volume as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificHeatCapacity cv "specific heat capacity";
algorithm 
  cv := cv_props_ph(p, h, waterBaseProp_ph(p, h, phase, region));
end cv_ph;

Modelica.Media.Water.IF97_Utilities.regionAssertReal Modelica.Media.Water.IF97_Utilities.regionAssertReal

assert function for inlining

Information

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

Inputs

TypeNameDefaultDescription
Booleancheck condition to check

Outputs

TypeNameDescription
Realdummydummy output

Modelica definition

function regionAssertReal "assert function for inlining"
  extends Modelica.Icons.Function;
  input Boolean check "condition to check";
  output Real dummy "dummy output";
algorithm 
  assert(check, "this function can not be called with two-phase inputs!");
end regionAssertReal;

Modelica.Media.Water.IF97_Utilities.cp_props_ph Modelica.Media.Water.IF97_Utilities.cp_props_ph

specific heat capacity at constant pressure as function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
SpecificHeatCapacitycpspecific heat capacity [J/(kg.K)]

Modelica definition

function cp_props_ph 
  "specific heat capacity at constant pressure as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificHeatCapacity cp "specific heat capacity";
algorithm 
  cp := aux.cp;
end cp_props_ph;

Modelica.Media.Water.IF97_Utilities.cp_ph Modelica.Media.Water.IF97_Utilities.cp_ph

specific heat capacity at constant pressure as function of pressure and specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificHeatCapacitycpspecific heat capacity [J/(kg.K)]

Modelica definition

function cp_ph 
  "specific heat capacity at constant pressure as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificHeatCapacity cp "specific heat capacity";
algorithm 
  cp := cp_props_ph(p, h, waterBaseProp_ph(p, h, phase, region));
end cp_ph;

Modelica.Media.Water.IF97_Utilities.beta_props_ph Modelica.Media.Water.IF97_Utilities.beta_props_ph

isobaric expansion coefficient as function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
RelativePressureCoefficientbetaisobaric expansion coefficient [1/K]

Modelica definition

function beta_props_ph 
  "isobaric expansion coefficient as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.RelativePressureCoefficient beta "isobaric expansion coefficient";
algorithm 
  beta := if aux.region == 3 or aux.region == 4 then 
    aux.pt/(aux.rho*aux.pd) else 
    aux.vt*aux.rho;
end beta_props_ph;

Modelica.Media.Water.IF97_Utilities.beta_ph Modelica.Media.Water.IF97_Utilities.beta_ph

isobaric expansion coefficient as function of pressure and specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
RelativePressureCoefficientbetaisobaric expansion coefficient [1/K]

Modelica definition

function beta_ph 
  "isobaric expansion coefficient as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.RelativePressureCoefficient beta "isobaric expansion coefficient";
algorithm 
  beta := beta_props_ph(p, h, waterBaseProp_ph(p, h, phase, region));
end beta_ph;

Modelica.Media.Water.IF97_Utilities.kappa_props_ph Modelica.Media.Water.IF97_Utilities.kappa_props_ph

isothermal compressibility factor as function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
IsothermalCompressibilitykappaisothermal compressibility factor [1/Pa]

Modelica definition

function kappa_props_ph 
  "isothermal compressibility factor as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.IsothermalCompressibility kappa "isothermal compressibility factor";
algorithm 
  kappa := if aux.region == 3 or aux.region == 4 then 
    1/(aux.rho*aux.pd) else -aux.vp*aux.rho;
end kappa_props_ph;

Modelica.Media.Water.IF97_Utilities.kappa_ph Modelica.Media.Water.IF97_Utilities.kappa_ph

isothermal compressibility factor as function of pressure and specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
IsothermalCompressibilitykappaisothermal compressibility factor [1/Pa]

Modelica definition

function kappa_ph 
  "isothermal compressibility factor as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.IsothermalCompressibility kappa "isothermal compressibility factor";
algorithm 
  kappa := kappa_props_ph(p, h, waterBaseProp_ph(p, h, phase, region));
end kappa_ph;

Modelica.Media.Water.IF97_Utilities.velocityOfSound_props_ph Modelica.Media.Water.IF97_Utilities.velocityOfSound_props_ph

speed of sound as function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
Velocityv_soundspeed of sound [m/s]

Modelica definition

function velocityOfSound_props_ph 
  "speed of sound as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.Velocity v_sound "speed of sound";
algorithm 
  // dp/drho at constant s
  v_sound := if aux.region == 3 then sqrt(max(0,(aux.pd*aux.rho*aux.rho*aux.cv + aux.pt*aux.pt*aux.T)/(aux.rho*aux.rho*aux.cv))) else 
    if aux.region == 4 then 
    sqrt(max(0,1/((aux.rho*(aux.rho*aux.cv/aux.dpT + 1.0)/(aux.dpT*aux.T)) - 1/aux.rho*aux.rho*aux.rho/(aux.dpT*aux.T)))) else 
         sqrt(max(0,-aux.cp/(aux.rho*aux.rho*(aux.vp*aux.cp+aux.vt*aux.vt*aux.T))));
end velocityOfSound_props_ph;

Modelica.Media.Water.IF97_Utilities.velocityOfSound_ph Modelica.Media.Water.IF97_Utilities.velocityOfSound_ph

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Velocityv_soundspeed of sound [m/s]

Modelica definition

function velocityOfSound_ph
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.Velocity v_sound "speed of sound";
algorithm 
  v_sound := velocityOfSound_props_ph(p, h, waterBaseProp_ph(p, h, phase, region));
end velocityOfSound_ph;

Modelica.Media.Water.IF97_Utilities.isentropicExponent_props_ph Modelica.Media.Water.IF97_Utilities.isentropicExponent_props_ph

isentropic exponent as function of pressure and specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
Realgammaisentropic exponent

Modelica definition

function isentropicExponent_props_ph 
  "isentropic exponent as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output Real gamma "isentropic exponent";
algorithm 
  gamma := if aux.region == 3 then 1/(aux.rho*p)*((aux.pd*aux.cv*aux.rho*aux.rho + aux.pt*aux.pt*aux.T)/(aux.cv)) else 
         if aux.region == 4 then 1/(aux.rho*p)*aux.dpT*aux.dpT*aux.T/aux.cv else 
    -1/(aux.rho*aux.p)*aux.cp/(aux.vp*aux.cp + aux.vt*aux.vt*aux.T);
end isentropicExponent_props_ph;

Modelica.Media.Water.IF97_Utilities.isentropicExponent_ph Modelica.Media.Water.IF97_Utilities.isentropicExponent_ph

isentropic exponent as function of pressure and specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Realgammaisentropic exponent

Modelica definition

function isentropicExponent_ph 
  "isentropic exponent as function of pressure and specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase =   0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output Real gamma "isentropic exponent";
algorithm 
  gamma := isentropicExponent_props_ph(p, h, waterBaseProp_ph(p, h, phase, region));
end isentropicExponent_ph;

Modelica.Media.Water.IF97_Utilities.ddph_props Modelica.Media.Water.IF97_Utilities.ddph_props

density derivative by pressure

Information

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

Inputs

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

Outputs

TypeNameDescription
DerDensityByPressureddphdensity derivative by pressure [s2/m2]

Modelica definition

function ddph_props "density derivative by pressure"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.DerDensityByPressure ddph "density derivative by pressure";
algorithm 
  ddph := if aux.region == 3 then 
    ((aux.rho*(aux.cv*aux.rho + aux.pt))/(aux.rho*aux.rho*aux.pd*aux.cv + aux.T*aux.pt*aux.pt)) else 
    if aux.region == 4 then  (aux.rho*(aux.rho*aux.cv/aux.dpT + 1.0)/(aux.dpT*aux.T)) else 
         (-aux.rho*aux.rho*(aux.vp*aux.cp - aux.vt/aux.rho + aux.T*aux.vt*aux.vt)/aux.cp);
end ddph_props;

Modelica.Media.Water.IF97_Utilities.ddph Modelica.Media.Water.IF97_Utilities.ddph

density derivative by pressure

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
DerDensityByPressureddphdensity derivative by pressure [s2/m2]

Modelica definition

function ddph "density derivative by pressure"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.DerDensityByPressure ddph "density derivative by pressure";
algorithm 
  ddph := ddph_props(p, h, waterBaseProp_ph(p, h, phase, region));
end ddph;

Modelica.Media.Water.IF97_Utilities.ddhp_props Modelica.Media.Water.IF97_Utilities.ddhp_props

density derivative by specific enthalpy

Information

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

Inputs

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

Outputs

TypeNameDescription
DerDensityByEnthalpyddhpdensity derivative by specific enthalpy [kg.s2/m5]

Modelica definition

function ddhp_props "density derivative by specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.DerDensityByEnthalpy ddhp "density derivative by specific enthalpy";
algorithm 
  ddhp := if aux.region == 3 then 
    -aux.rho*aux.rho*aux.pt/(aux.rho*aux.rho*aux.pd*aux.cv + aux.T*aux.pt*aux.pt) else 
    if aux.region == 4 then -aux.rho*aux.rho/(aux.dpT*aux.T) else 
         -aux.rho*aux.rho*aux.vt/(aux.cp);
end ddhp_props;

Modelica.Media.Water.IF97_Utilities.ddhp Modelica.Media.Water.IF97_Utilities.ddhp

density derivative by specific enthalpy

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEnthalpyh specific enthalpy [J/kg]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
DerDensityByEnthalpyddhpdensity derivative by specific enthalpy [kg.s2/m5]

Modelica definition

function ddhp "density derivative by specific enthalpy"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEnthalpy h "specific enthalpy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.DerDensityByEnthalpy ddhp "density derivative by specific enthalpy";
algorithm 
  ddhp := ddhp_props(p, h, waterBaseProp_ph(p, h, phase, region));
end ddhp;

Modelica.Media.Water.IF97_Utilities.waterBaseProp_pT Modelica.Media.Water.IF97_Utilities.waterBaseProp_pT

intermediate property record for water (p and T prefered states)

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
Integerregion0if 0, do region computation, otherwise assume the region is this input

Outputs

TypeNameDescription
IF97BaseTwoPhaseauxauxiliary record

Modelica definition

function waterBaseProp_pT 
  "intermediate property record for water (p and T prefered states)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Integer region = 0 
    "if 0, do region computation, otherwise assume the region is this input";
  output Common.IF97BaseTwoPhase aux "auxiliary record";
protected 
  Common.GibbsDerivs g 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Common.HelmholtzDerivs f 
    "dimensionless Helmholtz funcion and dervatives w.r.t. delta and tau";
  Integer error "error flag for inverse iterations";
algorithm 
  aux.phase := 1;
  aux.region := if region == 0 then BaseIF97.Regions.region_pT(p=p,T= T) else region;
  aux.R := BaseIF97.data.RH2O;
  aux.p := p;
  aux.T := T;
  if (aux.region == 1) then
    g := BaseIF97.Basic.g1(p, T);
    aux.h := aux.R*aux.T*g.tau*g.gtau;
    aux.s := aux.R*(g.tau*g.gtau - g.g);
    aux.rho := p/(aux.R*T*g.pi*g.gpi);
    aux.vt := aux.R/p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*T/(p*p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
    aux.x := 0.0;
  elseif (aux.region == 2) then
    g := BaseIF97.Basic.g2(p, T);
    aux.h := aux.R*aux.T*g.tau*g.gtau;
    aux.s := aux.R*(g.tau*g.gtau - g.g);
    aux.rho := p/(aux.R*T*g.pi*g.gpi);
    aux.vt := aux.R/p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*T/(p*p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
    aux.x := 1.0;
  elseif (aux.region == 3) then
    (aux.rho,error) := BaseIF97.Inverses.dofpt3(p=p,T= T,delp= 1.0e-7);
    f := BaseIF97.Basic.f3(aux.rho, T);
    aux.h := aux.R*T*(f.tau*f.ftau + f.delta*f.fdelta);
    aux.s := aux.R*(f.tau*f.ftau - f.f);
    aux.pd := aux.R*T*f.delta*(2.0*f.fdelta + f.delta*f.fdeltadelta);
    aux.pt := aux.R*aux.rho*f.delta*(f.fdelta - f.tau*f.fdeltatau);
    aux.cv := aux.R*(-f.tau*f.tau*f.ftautau);
    aux.x := 0.0;
  elseif (aux.region == 5) then
    g := BaseIF97.Basic.g5(p, T);
    aux.h := aux.R*aux.T*g.tau*g.gtau;
    aux.s := aux.R*(g.tau*g.gtau - g.g);
    aux.rho := p/(aux.R*T*g.pi*g.gpi);
    aux.vt := aux.R/p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*T/(p*p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
  else
    assert(false, "error in region computation of IF97 steam tables"
     + "(p = " + String(p) + ", T = " + String(T) + ")");
  end if;
end waterBaseProp_pT;

Modelica.Media.Water.IF97_Utilities.rho_props_pT Modelica.Media.Water.IF97_Utilities.rho_props_pT

density as function or pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
Densityrhodensity [kg/m3]

Modelica definition

function rho_props_pT 
  "density as function or pressure and temperature"
  annotation(derivative=rho_pT_der);
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.Density rho "density";
algorithm 
  rho := aux.rho;
end rho_props_pT;

Modelica.Media.Water.IF97_Utilities.rho_pT Modelica.Media.Water.IF97_Utilities.rho_pT

density as function or pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Densityrhodensity [kg/m3]

Modelica definition

function rho_pT "density as function or pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.Density rho "density";
algorithm 
  rho := rho_props_pT(p, T, waterBaseProp_pT(p, T, region));
end rho_pT;

Modelica.Media.Water.IF97_Utilities.h_props_pT Modelica.Media.Water.IF97_Utilities.h_props_pT

specific enthalpy as function or pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function h_props_pT 
  "specific enthalpy as function or pressure and temperature"
  annotation(derivative=h_pT_der);
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificEnthalpy h "specific enthalpy";
algorithm 
  h := aux.h;
end h_props_pT;

Modelica.Media.Water.IF97_Utilities.h_pT Modelica.Media.Water.IF97_Utilities.h_pT

specific enthalpy as function or pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT Temperature [K]
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function h_pT 
  "specific enthalpy as function or pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "Temperature";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificEnthalpy h "specific enthalpy";
algorithm 
  h := h_props_pT(p, T, waterBaseProp_pT(p, T, region));
end h_pT;

Modelica.Media.Water.IF97_Utilities.h_pT_der Modelica.Media.Water.IF97_Utilities.h_pT_der

derivative function of h_pT

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record
Realp_der derivative of pressure
RealT_der derivative of temperature

Outputs

TypeNameDescription
Realh_derderivative of specific enthalpy

Modelica definition

function h_pT_der "derivative function of h_pT"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  input Real p_der "derivative of pressure";
  input Real T_der "derivative of temperature";
  output Real h_der "derivative of specific enthalpy";
algorithm 
  if (aux.region == 3) then
    h_der := ((-aux.rho*aux.pd + T*aux.pt)/(aux.rho*aux.rho*aux.pd))*p_der +
      ((aux.rho*aux.rho*aux.pd*aux.cv + aux.T*aux.pt*aux.pt)/(aux.rho*aux.rho
      *aux.pd))*T_der;
  else
    //regions 1,2 or 5
    h_der := (1/aux.rho - aux.T*aux.vt)*p_der + aux.cp*T_der;
  end if;
end h_pT_der;

Modelica.Media.Water.IF97_Utilities.rho_pT_der Modelica.Media.Water.IF97_Utilities.rho_pT_der

derivative function of rho_pT

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record
Realp_der derivative of pressure
RealT_der derivative of temperature

Outputs

TypeNameDescription
Realrho_derderivative of density

Modelica definition

function rho_pT_der "derivative function of rho_pT"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  input Real p_der "derivative of pressure";
  input Real T_der "derivative of temperature";
  output Real rho_der "derivative of density";
algorithm 
  if (aux.region == 3) then
    rho_der := (1/aux.pd)*p_der - (aux.pt/aux.pd)*T_der;
  else
    //regions 1,2 or 5
    rho_der := (-aux.rho*aux.rho*aux.vp)*p_der + (-aux.rho*aux.rho*aux.vt)*
      T_der;
  end if;
end rho_pT_der;

Modelica.Media.Water.IF97_Utilities.s_props_pT Modelica.Media.Water.IF97_Utilities.s_props_pT

specific entropy as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
SpecificEntropysspecific entropy [J/(kg.K)]

Modelica definition

function s_props_pT 
  "specific entropy as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificEntropy s "specific entropy";
algorithm 
  s := aux.s;
end s_props_pT;

Modelica.Media.Water.IF97_Utilities.s_pT Modelica.Media.Water.IF97_Utilities.s_pT

temperature as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificEntropysspecific entropy [J/(kg.K)]

Modelica definition

function s_pT "temperature as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificEntropy s "specific entropy";
algorithm 
  s := s_props_pT(p, T, waterBaseProp_pT(p, T, region));
end s_pT;

Modelica.Media.Water.IF97_Utilities.cv_props_pT Modelica.Media.Water.IF97_Utilities.cv_props_pT

specific heat capacity at constant volume as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
SpecificHeatCapacitycvspecific heat capacity [J/(kg.K)]

Modelica definition

function cv_props_pT 
  "specific heat capacity at constant volume as function of pressure and temperature"

  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificHeatCapacity cv "specific heat capacity";
algorithm 
  cv := aux.cv;
end cv_props_pT;

Modelica.Media.Water.IF97_Utilities.cv_pT Modelica.Media.Water.IF97_Utilities.cv_pT

specific heat capacity at constant volume as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificHeatCapacitycvspecific heat capacity [J/(kg.K)]

Modelica definition

function cv_pT 
  "specific heat capacity at constant volume as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificHeatCapacity cv "specific heat capacity";
algorithm 
  cv := cv_props_pT(p, T, waterBaseProp_pT(p, T, region));
end cv_pT;

Modelica.Media.Water.IF97_Utilities.cp_props_pT Modelica.Media.Water.IF97_Utilities.cp_props_pT

specific heat capacity at constant pressure as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
SpecificHeatCapacitycpspecific heat capacity [J/(kg.K)]

Modelica definition

function cp_props_pT 
  "specific heat capacity at constant pressure as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificHeatCapacity cp "specific heat capacity";
algorithm 
  cp := if aux.region == 3 then 
    (aux.rho*aux.rho*aux.pd*aux.cv + aux.T*aux.pt*aux.pt)/(aux.rho*aux.rho*aux.pd) else 
    aux.cp;
end cp_props_pT;

Modelica.Media.Water.IF97_Utilities.cp_pT Modelica.Media.Water.IF97_Utilities.cp_pT

specific heat capacity at constant pressure as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificHeatCapacitycpspecific heat capacity [J/(kg.K)]

Modelica definition

function cp_pT 
  "specific heat capacity at constant pressure as function of pressure and temperature"

  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificHeatCapacity cp "specific heat capacity";
algorithm 
  cp := cp_props_pT(p, T, waterBaseProp_pT(p, T, region));
end cp_pT;

Modelica.Media.Water.IF97_Utilities.beta_props_pT Modelica.Media.Water.IF97_Utilities.beta_props_pT

isobaric expansion coefficient as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
RelativePressureCoefficientbetaisobaric expansion coefficient [1/K]

Modelica definition

function beta_props_pT 
  "isobaric expansion coefficient as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.RelativePressureCoefficient beta "isobaric expansion coefficient";
algorithm 
  beta := if aux.region == 3 then 
    aux.pt/(aux.rho*aux.pd) else 
    aux.vt*aux.rho;
end beta_props_pT;

Modelica.Media.Water.IF97_Utilities.beta_pT Modelica.Media.Water.IF97_Utilities.beta_pT

isobaric expansion coefficient as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
RelativePressureCoefficientbetaisobaric expansion coefficient [1/K]

Modelica definition

function beta_pT 
  "isobaric expansion coefficient as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.RelativePressureCoefficient beta "isobaric expansion coefficient";
algorithm 
  beta := beta_props_pT(p, T, waterBaseProp_pT(p, T, region));
end beta_pT;

Modelica.Media.Water.IF97_Utilities.kappa_props_pT Modelica.Media.Water.IF97_Utilities.kappa_props_pT

isothermal compressibility factor as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
IsothermalCompressibilitykappaisothermal compressibility factor [1/Pa]

Modelica definition

function kappa_props_pT 
  "isothermal compressibility factor as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.IsothermalCompressibility kappa "isothermal compressibility factor";
algorithm 
  kappa := if aux.region == 3 then 
    1/(aux.rho*aux.pd) else -aux.vp*aux.rho;
end kappa_props_pT;

Modelica.Media.Water.IF97_Utilities.kappa_pT Modelica.Media.Water.IF97_Utilities.kappa_pT

isothermal compressibility factor as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
IsothermalCompressibilitykappaisothermal compressibility factor [1/Pa]

Modelica definition

function kappa_pT 
  "isothermal compressibility factor as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.IsothermalCompressibility kappa "isothermal compressibility factor";
algorithm 
  kappa := kappa_props_pT(p, T, waterBaseProp_pT(p, T, region));
end kappa_pT;

Modelica.Media.Water.IF97_Utilities.velocityOfSound_props_pT Modelica.Media.Water.IF97_Utilities.velocityOfSound_props_pT

speed of sound as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
Velocityv_soundspeed of sound [m/s]

Modelica definition

function velocityOfSound_props_pT 
  "speed of sound as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.Velocity v_sound "speed of sound";
algorithm 
  // dp/drho at constant s
  v_sound := if aux.region == 3 then sqrt(max(0,(aux.pd*aux.rho*aux.rho*aux.cv + aux.pt*aux.pt*aux.T)/(aux.rho*aux.rho*aux.cv))) else 
    sqrt(max(0,-aux.cp/(aux.rho*aux.rho*(aux.vp*aux.cp+aux.vt*aux.vt*aux.T))));
end velocityOfSound_props_pT;

Modelica.Media.Water.IF97_Utilities.velocityOfSound_pT Modelica.Media.Water.IF97_Utilities.velocityOfSound_pT

speed of sound as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Velocityv_soundspeed of sound [m/s]

Modelica definition

function velocityOfSound_pT 
  "speed of sound as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.Velocity v_sound "speed of sound";
algorithm 
  v_sound := velocityOfSound_props_pT(p, T, waterBaseProp_pT(p, T, region));
end velocityOfSound_pT;

Modelica.Media.Water.IF97_Utilities.isentropicExponent_props_pT Modelica.Media.Water.IF97_Utilities.isentropicExponent_props_pT

isentropic exponent as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
Realgammaisentropic exponent

Modelica definition

function isentropicExponent_props_pT 
  "isentropic exponent as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output Real gamma "isentropic exponent";
algorithm 
  gamma := if aux.region == 3 then 1/(aux.rho*p)*((aux.pd*aux.cv*aux.rho*aux.rho + aux.pt*aux.pt*aux.T)/(aux.cv)) else 
    -1/(aux.rho*aux.p)*aux.cp/(aux.vp*aux.cp + aux.vt*aux.vt*aux.T);
end isentropicExponent_props_pT;

Modelica.Media.Water.IF97_Utilities.isentropicExponent_pT Modelica.Media.Water.IF97_Utilities.isentropicExponent_pT

isentropic exponent as function of pressure and temperature

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
TemperatureT temperature [K]
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Realgammaisentropic exponent

Modelica definition

function isentropicExponent_pT 
  "isentropic exponent as function of pressure and temperature"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.Temperature T "temperature";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output Real gamma "isentropic exponent";
algorithm 
  gamma := isentropicExponent_props_pT(p, T, waterBaseProp_pT(p, T, region));
end isentropicExponent_pT;

Modelica.Media.Water.IF97_Utilities.waterBaseProp_dT Modelica.Media.Water.IF97_Utilities.waterBaseProp_dT

intermediate property record for water (d and T prefered states)

Information

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

Inputs

TypeNameDefaultDescription
Densityrho density [kg/m3]
TemperatureT temperature [K]
Integerphase0phase: 2 for two-phase, 1 for one phase, 0 if unknown
Integerregion0if 0, do region computation, otherwise assume the region is this input

Outputs

TypeNameDescription
IF97BaseTwoPhaseauxauxiliary record

Modelica definition

function waterBaseProp_dT 
  "intermediate property record for water (d and T prefered states)"
  extends Modelica.Icons.Function;
  input SI.Density rho "density";
  input SI.Temperature T "temperature";
  input Integer phase =  0 
    "phase: 2 for two-phase, 1 for one phase, 0 if unknown";
  input Integer region = 0 
    "if 0, do region computation, otherwise assume the region is this input";
  output Common.IF97BaseTwoPhase aux "auxiliary record";
protected 
  SI.SpecificEnthalpy h_liq "liquid specific enthalpy";
  SI.Density d_liq "liquid density";
  SI.SpecificEnthalpy h_vap "vapour specific enthalpy";
  SI.Density d_vap "vapour density";
  Common.GibbsDerivs g 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Common.HelmholtzDerivs f 
    "dimensionless Helmholtz funcion and dervatives w.r.t. delta and tau";
  Modelica.Media.Common.PhaseBoundaryProperties liq 
    "phase boundary property record";
  Modelica.Media.Common.PhaseBoundaryProperties vap 
    "phase boundary property record";
  Modelica.Media.Common.GibbsDerivs gl 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Modelica.Media.Common.GibbsDerivs gv 
    "dimensionless Gibbs funcion and dervatives w.r.t. pi and tau";
  Modelica.Media.Common.HelmholtzDerivs fl 
    "dimensionless Helmholtz function and dervatives w.r.t. delta and tau";
  Modelica.Media.Common.HelmholtzDerivs fv 
    "dimensionless Helmholtz function and dervatives w.r.t. delta and tau";
  Integer error "error flag for inverse iterations";
algorithm 
  aux.region := if region == 0 then 
    (if phase == 2 then 4 else BaseIF97.Regions.region_dT(d=rho,T= T,phase= phase)) else region;
  aux.phase := if aux.region == 4 then 2 else 1;
  aux.R := BaseIF97.data.RH2O;
  aux.rho := rho;
  aux.T := T;
  if (aux.region == 1) then
    (aux.p,error) := BaseIF97.Inverses.pofdt125(d=rho,T= T,reldd= 1.0e-8,region=
             1);
    g := BaseIF97.Basic.g1(aux.p, T);
    aux.h := aux.R*aux.T*g.tau*g.gtau;
    aux.s := aux.R*(g.tau*g.gtau - g.g);
    aux.rho := aux.p/(aux.R*T*g.pi*g.gpi);
    aux.vt := aux.R/aux.p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*T/(aux.p*aux.p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
    aux.x := 0.0;
  elseif (aux.region == 2) then
    (aux.p,error) := BaseIF97.Inverses.pofdt125(d=rho,T= T,reldd= 1.0e-8,region=
             2);
    g := BaseIF97.Basic.g2(aux.p, T);
    aux.h := aux.R*aux.T*g.tau*g.gtau;
    aux.s := aux.R*(g.tau*g.gtau - g.g);
    aux.rho := aux.p/(aux.R*T*g.pi*g.gpi);
    aux.vt := aux.R/aux.p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*T/(aux.p*aux.p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
    aux.x := 1.0;
  elseif (aux.region == 3) then
    f := BaseIF97.Basic.f3(rho, T);
    aux.p := aux.R*rho*T*f.delta*f.fdelta;
    aux.h := aux.R*T*(f.tau*f.ftau + f.delta*f.fdelta);
    aux.s := aux.R*(f.tau*f.ftau - f.f);
    aux.pd := aux.R*T*f.delta*(2.0*f.fdelta + f.delta*f.fdeltadelta);
    aux.pt := aux.R*rho*f.delta*(f.fdelta - f.tau*f.fdeltatau);
    aux.cp := (aux.rho*aux.rho*aux.pd*aux.cv + aux.T*aux.pt*aux.pt)/(aux.rho*aux.rho*aux.pd);
    aux.cv := aux.R*(-f.tau*f.tau*f.ftautau);
    aux.x := 0.0;
  elseif (aux.region == 4) then
    aux.p := BaseIF97.Basic.psat(T);
    d_liq := rhol_T(T);
    d_vap := rhov_T(T);
    h_liq := hl_p(aux.p);
    h_vap := hv_p(aux.p);
    aux.x := if (d_vap <> d_liq) then (1/rho - 1/d_liq)/(1/d_vap - 1/d_liq) else 
    1.0;
    aux.h := h_liq + aux.x*(h_vap - h_liq);
    if T < BaseIF97.data.TLIMIT1 then
      gl := BaseIF97.Basic.g1(aux.p, T);
      gv := BaseIF97.Basic.g2(aux.p, T);
      liq := Common.gibbsToBoundaryProps(gl);
      vap := Common.gibbsToBoundaryProps(gv);
    else
      fl := BaseIF97.Basic.f3(d_liq, T);
      fv := BaseIF97.Basic.f3(d_vap, T);
      liq := Common.helmholtzToBoundaryProps(fl);
      vap := Common.helmholtzToBoundaryProps(fv);
    end if;
    aux.dpT := if (liq.d <> vap.d) then (vap.s - liq.s)*liq.d*vap.d/(liq.d - vap.d) else BaseIF97.Basic.dptofT(aux.T);
    aux.s := liq.s + aux.x*(vap.s - liq.s);
    aux.cv := Common.cv2Phase(liq, vap, aux.x, aux.T, aux.p);
    aux.cp := liq.cp + aux.x*(vap.cp - liq.cp);
    aux.pt := liq.pt + aux.x*(vap.pt - liq.pt);
    aux.pd := liq.pd + aux.x*(vap.pd - liq.pd);
  elseif (aux.region == 5) then
    (aux.p,error) := BaseIF97.Inverses.pofdt125(d=rho,T= T,reldd= 1.0e-8,region=5);
    g := BaseIF97.Basic.g2(aux.p, T);
    aux.h := aux.R*aux.T*g.tau*g.gtau;
    aux.s := aux.R*(g.tau*g.gtau - g.g);
    aux.rho := aux.p/(aux.R*T*g.pi*g.gpi);
    aux.vt := aux.R/aux.p*(g.pi*g.gpi - g.tau*g.pi*g.gtaupi);
    aux.vp := aux.R*T/(aux.p*aux.p)*g.pi*g.pi*g.gpipi;
    aux.cp := -aux.R*g.tau*g.tau*g.gtautau;
    aux.cv := aux.R*(-g.tau*g.tau*g.gtautau + ((g.gpi - g.tau*g.gtaupi)*(g.gpi - g.tau*g.gtaupi)/g.gpipi));
  else
    assert(false, "error in region computation of IF97 steam tables"
     + "(rho = " + String(rho) + ", T = " + String(T) + ")");
  end if;
end waterBaseProp_dT;

Modelica.Media.Water.IF97_Utilities.h_props_dT Modelica.Media.Water.IF97_Utilities.h_props_dT

specific enthalpy as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT Temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function h_props_dT 
  "specific enthalpy as function of density and temperature"
  annotation(derivative=h_dT_der);
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "Temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificEnthalpy h "specific enthalpy";
algorithm 
  h := aux.h;
end h_props_dT;

Modelica.Media.Water.IF97_Utilities.h_dT Modelica.Media.Water.IF97_Utilities.h_dT

specific enthalpy as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT Temperature [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function h_dT 
  "specific enthalpy as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "Temperature";
  input Integer phase =  0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificEnthalpy h "specific enthalpy";
algorithm 
  h := h_props_dT(d, T, waterBaseProp_dT(d, T, phase, region));
end h_dT;

Modelica.Media.Water.IF97_Utilities.h_dT_der Modelica.Media.Water.IF97_Utilities.h_dT_der

derivative function of h_dT

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record
Reald_der derivative of density
RealT_der derivative of temperature

Outputs

TypeNameDescription
Realh_derderivative of specific enthalpy

Modelica definition

function h_dT_der "derivative function of h_dT"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  input Real d_der "derivative of density";
  input Real T_der "derivative of temperature";
  output Real h_der "derivative of specific enthalpy";
algorithm 
  if (aux.region == 3) then
    h_der := ((-d*aux.pd + T*aux.pt)/(d*d))*d_der + ((aux.cv*d + aux.pt)/d)*
      T_der;
  elseif (aux.region == 4) then
    h_der := T*aux.dpT/(d*d)*d_der + ((aux.cv*d + aux.dpT)/d)*T_der;
  else
    //regions 1,2 or 5
    h_der := (-(-1/d + T*aux.vt)/(d*d*aux.vp))*d_der + ((aux.vp*aux.cp - aux.
      vt/d + T*aux.vt*aux.vt)/aux.vp)*T_der;
  end if;
end h_dT_der;

Modelica.Media.Water.IF97_Utilities.p_props_dT Modelica.Media.Water.IF97_Utilities.p_props_dT

pressure as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT Temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
Pressureppressure [Pa]

Modelica definition

function p_props_dT "pressure as function of density and temperature"
  annotation(derivative=p_dT_der);
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "Temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.Pressure p "pressure";
algorithm 
  p := aux.p;
end p_props_dT;

Modelica.Media.Water.IF97_Utilities.p_dT Modelica.Media.Water.IF97_Utilities.p_dT

pressure as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT Temperature [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Pressureppressure [Pa]

Modelica definition

function p_dT "pressure as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "Temperature";
  input Integer phase =  0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.Pressure p "pressure";
algorithm 
  p := p_props_dT(d, T, waterBaseProp_dT(d, T, phase, region));
end p_dT;

Modelica.Media.Water.IF97_Utilities.p_dT_der Modelica.Media.Water.IF97_Utilities.p_dT_der

derivative function of p_dT

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record
Reald_der derivative of density
RealT_der derivative of temperature

Outputs

TypeNameDescription
Realp_derderivative of pressure

Modelica definition

function p_dT_der "derivative function of p_dT"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  input Real d_der "derivative of density";
  input Real T_der "derivative of temperature";
  output Real p_der "derivative of pressure";
algorithm 
  if (aux.region == 3) then
    p_der := aux.pd*d_der + aux.pt*T_der;
  elseif (aux.region == 4) then
    p_der := aux.dpT*T_der;
    /*density derivative is 0.0*/
  else
    //regions 1,2 or 5
    p_der := (-1/(d*d*aux.vp))*d_der + (-aux.vt/aux.vp)*T_der;
  end if;
end p_dT_der;

Modelica.Media.Water.IF97_Utilities.s_props_dT Modelica.Media.Water.IF97_Utilities.s_props_dT

specific entropy as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT Temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
SpecificEntropysspecific entropy [J/(kg.K)]

Modelica definition

function s_props_dT 
  "specific entropy as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "Temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificEntropy s "specific entropy";
algorithm 
  s := aux.s;
end s_props_dT;

Modelica.Media.Water.IF97_Utilities.s_dT Modelica.Media.Water.IF97_Utilities.s_dT

temperature as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT Temperature [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificEntropysspecific entropy [J/(kg.K)]

Modelica definition

function s_dT "temperature as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "Temperature";
  input Integer phase =  0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificEntropy s "specific entropy";
algorithm 
  s := s_props_dT(d, T, waterBaseProp_dT(d, T, phase, region));
end s_dT;

Modelica.Media.Water.IF97_Utilities.cv_props_dT Modelica.Media.Water.IF97_Utilities.cv_props_dT

specific heat capacity at constant volume as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
SpecificHeatCapacitycvspecific heat capacity [J/(kg.K)]

Modelica definition

function cv_props_dT 
  "specific heat capacity at constant volume as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificHeatCapacity cv "specific heat capacity";
algorithm 
  cv := aux.cv;
end cv_props_dT;

Modelica.Media.Water.IF97_Utilities.cv_dT Modelica.Media.Water.IF97_Utilities.cv_dT

specific heat capacity at constant volume as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificHeatCapacitycvspecific heat capacity [J/(kg.K)]

Modelica definition

function cv_dT 
  "specific heat capacity at constant volume as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Integer phase =  0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificHeatCapacity cv "specific heat capacity";
algorithm 
  cv := cv_props_dT(d, T, waterBaseProp_dT(d, T, phase, region));
end cv_dT;

Modelica.Media.Water.IF97_Utilities.cp_props_dT Modelica.Media.Water.IF97_Utilities.cp_props_dT

specific heat capacity at constant pressure as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
SpecificHeatCapacitycpspecific heat capacity [J/(kg.K)]

Modelica definition

function cp_props_dT 
  "specific heat capacity at constant pressure as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificHeatCapacity cp "specific heat capacity";
algorithm 
  cp := aux.cp;
end cp_props_dT;

Modelica.Media.Water.IF97_Utilities.cp_dT Modelica.Media.Water.IF97_Utilities.cp_dT

specific heat capacity at constant pressure as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificHeatCapacitycpspecific heat capacity [J/(kg.K)]

Modelica definition

function cp_dT 
  "specific heat capacity at constant pressure as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Integer phase =  0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificHeatCapacity cp "specific heat capacity";
algorithm 
  cp := cp_props_dT(d, T, waterBaseProp_dT(d, T, phase, region));
end cp_dT;

Modelica.Media.Water.IF97_Utilities.beta_props_dT Modelica.Media.Water.IF97_Utilities.beta_props_dT

isobaric expansion coefficient as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
RelativePressureCoefficientbetaisobaric expansion coefficient [1/K]

Modelica definition

function beta_props_dT 
  "isobaric expansion coefficient as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.RelativePressureCoefficient beta "isobaric expansion coefficient";
algorithm 
  beta := if aux.region == 3 or aux.region == 4 then 
    aux.pt/(aux.rho*aux.pd) else 
    aux.vt*aux.rho;
end beta_props_dT;

Modelica.Media.Water.IF97_Utilities.beta_dT Modelica.Media.Water.IF97_Utilities.beta_dT

isobaric expansion coefficient as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
RelativePressureCoefficientbetaisobaric expansion coefficient [1/K]

Modelica definition

function beta_dT 
  "isobaric expansion coefficient as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Integer phase =  0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.RelativePressureCoefficient beta "isobaric expansion coefficient";
algorithm 
  beta := beta_props_dT(d, T, waterBaseProp_dT(d, T, phase, region));
end beta_dT;

Modelica.Media.Water.IF97_Utilities.kappa_props_dT Modelica.Media.Water.IF97_Utilities.kappa_props_dT

isothermal compressibility factor as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
IsothermalCompressibilitykappaisothermal compressibility factor [1/Pa]

Modelica definition

function kappa_props_dT 
  "isothermal compressibility factor as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.IsothermalCompressibility kappa "isothermal compressibility factor";
algorithm 
  kappa := if aux.region == 3 or aux.region == 4 then 
    1/(aux.rho*aux.pd) else -aux.vp*aux.rho;
end kappa_props_dT;

Modelica.Media.Water.IF97_Utilities.kappa_dT Modelica.Media.Water.IF97_Utilities.kappa_dT

isothermal compressibility factor as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
IsothermalCompressibilitykappaisothermal compressibility factor [1/Pa]

Modelica definition

function kappa_dT 
  "isothermal compressibility factor as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Integer phase =  0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.IsothermalCompressibility kappa "isothermal compressibility factor";
algorithm 
  kappa := kappa_props_dT(d, T, waterBaseProp_dT(d, T, phase, region));
end kappa_dT;

Modelica.Media.Water.IF97_Utilities.velocityOfSound_props_dT Modelica.Media.Water.IF97_Utilities.velocityOfSound_props_dT

speed of sound as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
Velocityv_soundspeed of sound [m/s]

Modelica definition

function velocityOfSound_props_dT 
  "speed of sound as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.Velocity v_sound "speed of sound";
algorithm 
  // dp/drho at constant s
  v_sound := if aux.region == 3 then sqrt(max(0,((aux.pd*aux.rho*aux.rho*aux.cv + aux.pt*aux.pt*aux.T)/(aux.rho*aux.rho*aux.cv)))) else 
    if aux.region == 4 then 
    sqrt(max(0,1/((aux.rho*(aux.rho*aux.cv/aux.dpT + 1.0)/(aux.dpT*aux.T)) - 1/aux.rho*aux.rho*aux.rho/(aux.dpT*aux.T)))) else 
         sqrt(max(0,-aux.cp/(aux.rho*aux.rho*(aux.vp*aux.cp+aux.vt*aux.vt*aux.T))));
end velocityOfSound_props_dT;

Modelica.Media.Water.IF97_Utilities.velocityOfSound_dT Modelica.Media.Water.IF97_Utilities.velocityOfSound_dT

speed of sound as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Velocityv_soundspeed of sound [m/s]

Modelica definition

function velocityOfSound_dT 
  "speed of sound as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Integer phase =  0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.Velocity v_sound "speed of sound";
algorithm 
  v_sound := velocityOfSound_props_dT(d, T, waterBaseProp_dT(d, T, phase, region));
end velocityOfSound_dT;

Modelica.Media.Water.IF97_Utilities.isentropicExponent_props_dT Modelica.Media.Water.IF97_Utilities.isentropicExponent_props_dT

isentropic exponent as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
IF97BaseTwoPhaseaux auxiliary record

Outputs

TypeNameDescription
Realgammaisentropic exponent

Modelica definition

function isentropicExponent_props_dT 
  "isentropic exponent as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output Real gamma "isentropic exponent";
algorithm 
  gamma := if aux.region == 3 then 1/(aux.rho*aux.p)*((aux.pd*aux.cv*aux.rho*aux.rho + aux.pt*aux.pt*aux.T)/(aux.cv)) else 
         if aux.region == 4 then 1/(aux.rho*aux.p)*aux.dpT*aux.dpT*aux.T/aux.cv else 
    -1/(aux.rho*aux.p)*aux.cp/(aux.vp*aux.cp + aux.vt*aux.vt*aux.T);
end isentropicExponent_props_dT;

Modelica.Media.Water.IF97_Utilities.isentropicExponent_dT Modelica.Media.Water.IF97_Utilities.isentropicExponent_dT

isentropic exponent as function of density and temperature

Information

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

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
Realgammaisentropic exponent

Modelica definition

function isentropicExponent_dT 
  "isentropic exponent as function of density and temperature"
  extends Modelica.Icons.Function;
  input SI.Density d "density";
  input SI.Temperature T "temperature";
  input Integer phase =  0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region =  0 
    "if 0, region is unknown, otherwise known and this input";
  output Real gamma "isentropic exponent";
algorithm 
  gamma := isentropicExponent_props_dT(d, T, waterBaseProp_dT(d, T, phase, region));
end isentropicExponent_dT;

Modelica.Media.Water.IF97_Utilities.hl_p Modelica.Media.Water.IF97_Utilities.hl_p

compute the saturated liquid specific h(p)

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function hl_p = BaseIF97.Regions.hl_p 
  "compute the saturated liquid specific h(p)";

Modelica.Media.Water.IF97_Utilities.hv_p Modelica.Media.Water.IF97_Utilities.hv_p

compute the saturated vapour specific h(p)

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function hv_p = BaseIF97.Regions.hv_p 
  "compute the saturated vapour specific h(p)";

Modelica.Media.Water.IF97_Utilities.sl_p Modelica.Media.Water.IF97_Utilities.sl_p

compute the saturated liquid specific s(p)

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]

Outputs

TypeNameDescription
SpecificEntropysspecific entropy [J/(kg.K)]

Modelica definition

function sl_p = BaseIF97.Regions.sl_p 
  "compute the saturated liquid specific s(p)";

Modelica.Media.Water.IF97_Utilities.sv_p Modelica.Media.Water.IF97_Utilities.sv_p

compute the saturated vapour specific s(p)

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]

Outputs

TypeNameDescription
SpecificEntropysspecific entropy [J/(kg.K)]

Modelica definition

function sv_p = BaseIF97.Regions.sv_p 
  "compute the saturated vapour specific s(p)";

Modelica.Media.Water.IF97_Utilities.rhol_T Modelica.Media.Water.IF97_Utilities.rhol_T

compute the saturated liquid d(T)

Inputs

TypeNameDefaultDescription
TemperatureT temperature [K]

Outputs

TypeNameDescription
Densityddensity of water at the boiling point [kg/m3]

Modelica definition

function rhol_T = BaseIF97.Regions.rhol_T "compute the saturated liquid d(T)";

Modelica.Media.Water.IF97_Utilities.rhov_T Modelica.Media.Water.IF97_Utilities.rhov_T

compute the saturated vapour d(T)

Inputs

TypeNameDefaultDescription
TemperatureT temperature [K]

Outputs

TypeNameDescription
Densityddensity of steam at the condensation point [kg/m3]

Modelica definition

function rhov_T = BaseIF97.Regions.rhov_T "compute the saturated vapour d(T)";

Modelica.Media.Water.IF97_Utilities.rhol_p Modelica.Media.Water.IF97_Utilities.rhol_p

compute the saturated liquid d(p)

Inputs

TypeNameDefaultDescription
Pressurep saturation pressure [Pa]

Outputs

TypeNameDescription
Densityrhodensity of steam at the condensation point [kg/m3]

Modelica definition

function rhol_p = BaseIF97.Regions.rhol_p "compute the saturated liquid d(p)";

Modelica.Media.Water.IF97_Utilities.rhov_p Modelica.Media.Water.IF97_Utilities.rhov_p

compute the saturated vapour d(p)

Inputs

TypeNameDefaultDescription
Pressurep saturation pressure [Pa]

Outputs

TypeNameDescription
Densityrhodensity of steam at the condensation point [kg/m3]

Modelica definition

function rhov_p = BaseIF97.Regions.rhov_p "compute the saturated vapour d(p)";

Modelica.Media.Water.IF97_Utilities.dynamicViscosity Modelica.Media.Water.IF97_Utilities.dynamicViscosity

compute eta(d,T) in the one-phase region

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature (K) [K]
Pressurep pressure (only needed for region of validity) [Pa]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known

Outputs

TypeNameDescription
DynamicViscosityetadynamic viscosity [Pa.s]

Modelica definition

function dynamicViscosity = BaseIF97.Transport.visc_dTp 
  "compute eta(d,T) in the one-phase region";

Modelica.Media.Water.IF97_Utilities.thermalConductivity Modelica.Media.Water.IF97_Utilities.thermalConductivity

compute lambda(d,T,p) in the one-phase region

Inputs

TypeNameDefaultDescription
Densityd density [kg/m3]
TemperatureT temperature (K) [K]
Pressurep pressure [Pa]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
BooleanindustrialMethodtrueif true, the industrial method is used, otherwise the scientific one

Outputs

TypeNameDescription
ThermalConductivitylambdathermal conductivity [W/(m.K)]

Modelica definition

function thermalConductivity = BaseIF97.Transport.cond_dTp 
  "compute lambda(d,T,p) in the one-phase region";

Modelica.Media.Water.IF97_Utilities.surfaceTension Modelica.Media.Water.IF97_Utilities.surfaceTension

compute sigma(T) at saturation T

Inputs

TypeNameDefaultDescription
TemperatureT temperature (K) [K]

Outputs

TypeNameDescription
SurfaceTensionsigmasurface tension in SI units [N/m]

Modelica definition

function surfaceTension = BaseIF97.Transport.surfaceTension 
  "compute sigma(T) at saturation T";

Modelica.Media.Water.IF97_Utilities.isentropicEnthalpy Modelica.Media.Water.IF97_Utilities.isentropicEnthalpy

isentropic specific enthalpy from p,s (preferably use dynamicIsentropicEnthalpy in dynamic simulation!)

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known
Integerregion0if 0, region is unknown, otherwise known and this input

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function isentropicEnthalpy 
  "isentropic specific enthalpy from p,s (preferably use dynamicIsentropicEnthalpy in dynamic simulation!)"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  input Integer region = 0 
    "if 0, region is unknown, otherwise known and this input";
  output SI.SpecificEnthalpy h "specific enthalpy";
algorithm 
  h := isentropicEnthalpy_props(p, s, waterBaseProp_ps(p, s, phase, region));
end isentropicEnthalpy;

Modelica.Media.Water.IF97_Utilities.isentropicEnthalpy_props Modelica.Media.Water.IF97_Utilities.isentropicEnthalpy_props

Information

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

Inputs

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

Outputs

TypeNameDescription
SpecificEnthalpyhisentropic enthalpay [J/kg]

Modelica definition

function isentropicEnthalpy_props
  annotation(derivative=isentropicEnthalpy_der);
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  output SI.SpecificEnthalpy h "isentropic enthalpay";
algorithm 
  h := aux.h;
end isentropicEnthalpy_props;

Modelica.Media.Water.IF97_Utilities.isentropicEnthalpy_der Modelica.Media.Water.IF97_Utilities.isentropicEnthalpy_der

derivative of isentropic specific enthalpy from p,s

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]
IF97BaseTwoPhaseaux auxiliary record
Realp_der pressure derivative
Reals_der entropy derivative

Outputs

TypeNameDescription
Realh_derspecific enthalpy derivative

Modelica definition

function isentropicEnthalpy_der 
  "derivative of isentropic specific enthalpy from p,s"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input Common.IF97BaseTwoPhase aux "auxiliary record";
  input Real p_der "pressure derivative";
  input Real s_der "entropy derivative";
  output Real h_der "specific enthalpy derivative";
algorithm 
  h_der := 1/aux.rho*p_der + aux.T*s_der;
end isentropicEnthalpy_der;

Modelica.Media.Water.IF97_Utilities.dynamicIsentropicEnthalpy Modelica.Media.Water.IF97_Utilities.dynamicIsentropicEnthalpy

isentropic specific enthalpy from p,s and good guesses of d and T

Information

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

Inputs

TypeNameDefaultDescription
Pressurep pressure [Pa]
SpecificEntropys specific entropy [J/(kg.K)]
Densitydguess good guess density, e.g., from adjacent volume [kg/m3]
TemperatureTguess good guess temperature, e.g., from adjacent volume [K]
Integerphase02 for two-phase, 1 for one-phase, 0 if not known

Outputs

TypeNameDescription
SpecificEnthalpyhspecific enthalpy [J/kg]

Modelica definition

function dynamicIsentropicEnthalpy 
  "isentropic specific enthalpy from p,s and good guesses of d and T"
  extends Modelica.Icons.Function;
  input SI.Pressure p "pressure";
  input SI.SpecificEntropy s "specific entropy";
  input SI.Density dguess "good guess density, e.g., from adjacent volume";
  input SI.Temperature Tguess 
    "good guess temperature, e.g., from adjacent volume";
  input Integer phase = 0 "2 for two-phase, 1 for one-phase, 0 if not known";
  output SI.SpecificEnthalpy h "specific enthalpy";
algorithm 
 h := BaseIF97.Isentropic.water_hisentropic_dyn(p,s,dguess,Tguess,0);
end dynamicIsentropicEnthalpy;

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