Skip to content
This repository
Browse code

derivatives along saturation boundary

  • Loading branch information...
commit d9215d417c19b3dd185e0880a090747fbe86255f 1 parent 294e829
Matthis Thorade authored
129 Examples/Validate_Derivatives_SaturationBoundary.mo
... ... @@ -0,0 +1,129 @@
  1 +within HelmholtzMedia.Examples;
  2 +model Validate_Derivatives_SaturationBoundary
  3 + "compare analytical derivatives to numerical derivatives"
  4 +
  5 + package medium = HelmholtzFluids.R134a;
  6 + // choose subcritical T
  7 + parameter medium.Temperature T=298.15;
  8 +
  9 + medium.SaturationProperties sat;
  10 + medium.SaturationProperties sat_Tplus;
  11 + medium.SaturationProperties sat_Tminus;
  12 + medium.SaturationProperties sat_pplus;
  13 + medium.SaturationProperties sat_pminus;
  14 +
  15 + // Enthalpy derivatives
  16 + // medium.Types.DerEnthalpyByDensity dhdT_analytical;
  17 + // medium.Types.DerEnthalpyByDensity dhdT_numerical;
  18 + medium.Types.DerEnthalpyByTemperature dhT_liq_analytical;
  19 + medium.Types.DerEnthalpyByTemperature dhT_liq_numerical;
  20 + // Energy derivatives
  21 + // medium.Types.DerEnergyByDensity dudT_analytical;
  22 + // medium.Types.DerEnergyByDensity dudT_numerical;
  23 + // medium.Types.DerEnergyByTemperature duTd_analytical;
  24 + // medium.Types.DerEnergyByTemperature duTd_numerical;
  25 + // Entropy derivatives
  26 + // medium.Types.DerEntropyByDensity dsdT_analytical;
  27 + // medium.Types.DerEntropyByDensity dsdT_numerical;
  28 + // medium.Types.DerEntropyByTemperature dsTd_analytical;
  29 + // medium.Types.DerEntropyByTemperature dsTd_numerical;
  30 + // Gibbs derivatives
  31 + // medium.Types.DerEnergyByDensity dgdT_analytical;
  32 + // medium.Types.DerEnergyByDensity dgdT_numerical;
  33 + // medium.Types.DerEnergyByTemperature dgTd_analytical;
  34 + // medium.Types.DerEnergyByTemperature dgTd_numerical;
  35 + // Density derivatives
  36 + medium.DerDensityByTemperature ddT_liq_analytical;
  37 + medium.DerDensityByTemperature ddT_liq_numerical;
  38 + medium.DerDensityByTemperature ddT_vap_analytical;
  39 + medium.DerDensityByTemperature ddT_vap_numerical;
  40 + medium.DerDensityByPressure ddp_liq_analytical;
  41 + medium.DerDensityByPressure ddp_liq_numerical;
  42 + medium.DerDensityByPressure ddp_vap_analytical;
  43 + medium.DerDensityByPressure ddp_vap_numerical;
  44 +
  45 +equation
  46 + sat=medium.setSat_T(T=T);
  47 + sat_Tplus = medium.setSat_T(T=1.0001*T);
  48 + sat_Tminus = medium.setSat_T(T=0.9999*T);
  49 + sat_pplus = medium.setSat_p(p=1.0001*sat.psat);
  50 + sat_pminus = medium.setSat_p(p=0.9999*sat.psat);
  51 +
  52 + Modelica.Utilities.Streams.print("====|====|====|====|====|====|====|====|====|====|====|====|====|====|====|====|");
  53 +
  54 + Modelica.Utilities.Streams.print("Density");
  55 + // check (dd/dT)@liq
  56 + ddT_liq_analytical = medium.density_derT_p(state=sat.liq) +medium.density_derp_T(state=sat.liq)*medium.saturationPressure_derT(T=sat.Tsat, sat=sat);
  57 + ddT_liq_numerical = (sat_Tplus.liq.d - sat_Tminus.liq.d)/(sat_Tplus.liq.T - sat_Tminus.liq.T);
  58 + Modelica.Utilities.Streams.print("(dd/dT)@liq analytical= " + String(ddT_liq_analytical));
  59 + Modelica.Utilities.Streams.print("(dd/dT)@liq numerical= " + String(ddT_liq_numerical));
  60 + // check (dd/dT)@vap
  61 + ddT_vap_analytical = medium.density_derT_p(state=sat.vap) +medium.density_derp_T(state=sat.vap)*medium.saturationPressure_derT(T=sat.Tsat, sat=sat);
  62 + ddT_vap_numerical = (sat_Tplus.vap.d - sat_Tminus.vap.d)/(sat_Tplus.vap.T - sat_Tminus.vap.T);
  63 + Modelica.Utilities.Streams.print("(dd/dT)@vap analytical= " + String(ddT_vap_analytical));
  64 + Modelica.Utilities.Streams.print("(dd/dT)@vap numerical= " + String(ddT_vap_numerical));
  65 + // check (dd/dp)@liq
  66 + ddp_liq_analytical = medium.density_derp_T(state=sat.liq) +medium.density_derT_p(state=sat.liq)*medium.saturationTemperature_derp(p=sat.psat, sat=sat);
  67 + ddp_liq_numerical = (sat_pplus.liq.d - sat_pminus.liq.d)/(sat_pplus.liq.p - sat_pminus.liq.p);
  68 + Modelica.Utilities.Streams.print("(dd/dp)@liq analytical= " + String(ddp_liq_analytical));
  69 + Modelica.Utilities.Streams.print("(dd/dp)@liq numerical= " + String(ddp_liq_numerical));
  70 + // check (dd/dp)@vap
  71 + ddp_vap_analytical = medium.density_derp_T(state=sat.vap) +medium.density_derT_p(state=sat.vap)*medium.saturationTemperature_derp(p=sat.psat, sat=sat);
  72 + ddp_vap_numerical = (sat_pplus.vap.d - sat_pminus.vap.d)/(sat_pplus.vap.p - sat_pminus.vap.p);
  73 + Modelica.Utilities.Streams.print("(dd/dp)@vap analytical= " + String(ddp_vap_analytical));
  74 + Modelica.Utilities.Streams.print("(dd/dp)@vap numerical= " + String(ddp_vap_numerical));
  75 +
  76 + Modelica.Utilities.Streams.print(" ");
  77 + // Modelica.Utilities.Streams.print("Enthalpy");
  78 + // check (dh/dd)@T=const
  79 + // dhdT_analytical = medium.specificEnthalpy_derd_T(state=state, f=f);
  80 + // dhdT_numerical = (d_plus.h-d_minus.h)/(d_plus.d-d_minus.d);
  81 + // Modelica.Utilities.Streams.print("(dh/dd)@T=const analytical= " + String(dhdT_analytical));
  82 + // Modelica.Utilities.Streams.print("(dh/dd)@T=const numerical= " + String(dhdT_numerical));
  83 + // check (dh/dT)@d=const
  84 + dhT_liq_analytical = medium.specificEnthalpy_derT_d(state=sat.liq);
  85 + dhT_liq_numerical = (sat_Tplus.liq.h-sat_Tminus.liq.h)/(sat_Tplus.liq.T-sat_Tminus.liq.T);
  86 + Modelica.Utilities.Streams.print("(dh/dT)@liq analytical= " + String(dhT_liq_analytical));
  87 + Modelica.Utilities.Streams.print("(dh/dT)@liq numerical= " + String(dhT_liq_numerical));
  88 +
  89 + // Modelica.Utilities.Streams.print(" ");
  90 + // Modelica.Utilities.Streams.print("internal Energy");
  91 + // check (du/dd)@T=const
  92 + // dudT_analytical = f.R*T/d*f.tau*f.delta*f.rtd;
  93 + // dudT_numerical = (d_plus.u-d_minus.u)/(d_plus.d-d_minus.d);
  94 + // Modelica.Utilities.Streams.print("(du/dd)@T=const analytical= " + String(dudT_analytical));
  95 + // Modelica.Utilities.Streams.print("(du/dd)@T=const numerical= " + String(dudT_numerical));
  96 + // check (du/dT)@d=const
  97 + // duTd_analytical = medium.specificHeatCapacityCv(state=state, f=f);
  98 + // duTd_numerical = (T_plus.u-T_minus.u)/(T_plus.T-T_minus.T);
  99 + // Modelica.Utilities.Streams.print("(du/dT)@d=const analytical= " + String(duTd_analytical));
  100 + // Modelica.Utilities.Streams.print("(du/dT)@d=const numerical= " + String(duTd_numerical));
  101 +
  102 + // Modelica.Utilities.Streams.print(" ");
  103 + // Modelica.Utilities.Streams.print("Entropy");
  104 + // check (ds/dd)@T=const
  105 + // dsdT_analytical = f.R/d*(-(1+f.delta*f.rd)+(0+f.tau*f.delta*f.rtd));
  106 + // dsdT_numerical = (d_plus.s-d_minus.s)/(d_plus.d-d_minus.d);
  107 + // Modelica.Utilities.Streams.print("(ds/dd)@T=const analytical= " + String(dsdT_analytical));
  108 + // Modelica.Utilities.Streams.print("(ds/dd)@T=const numerical= " + String(dsdT_numerical));
  109 + // check (ds/dT)@d=const
  110 + // dsTd_analytical = f.R/T*(-f.tau^2*(f.itt+f.rtt));
  111 + // dsTd_numerical = (T_plus.s-T_minus.s)/(T_plus.T-T_minus.T);
  112 + // Modelica.Utilities.Streams.print("(ds/dT)@d=const analytical= " + String(dsTd_analytical));
  113 + // Modelica.Utilities.Streams.print("(ds/dT)@d=const numerical= " + String(dsTd_numerical));
  114 +
  115 + // Modelica.Utilities.Streams.print(" ");
  116 + // Modelica.Utilities.Streams.print("Gibbs energy");
  117 + // check (dg/dd)@T=const
  118 + // dgdT_analytical = f.R*T/d*(1+2*f.delta*f.rd + f.delta^2*f.rdd);
  119 + // dgdT_numerical = ((d_plus.h-d_plus.T*d_plus.s)-(d_minus.h-d_minus.T*d_minus.s))/(d_plus.d-d_minus.d);
  120 + // Modelica.Utilities.Streams.print("(dg/dd)@T=const analytical= " + String(dgdT_analytical));
  121 + // Modelica.Utilities.Streams.print("(dg/dd)@T=const numerical= " + String(dgdT_numerical));
  122 + // check (dg/dT)@d=const
  123 + // dgTd_analytical = f.R*(f.i+f.r + 1+f.delta*f.rd -f.tau*(f.it+f.rt) - f.tau*f.delta*f.rtd);
  124 + // dgTd_numerical = ((T_plus.h-T_plus.T*T_plus.s)-(T_minus.h-T_minus.T*T_minus.s))/(T_plus.T-T_minus.T);
  125 + // Modelica.Utilities.Streams.print("(dg/dT)@d=const analytical= " + String(dgTd_analytical));
  126 + // Modelica.Utilities.Streams.print("(dg/dT)@d=const numerical= " + String(dgTd_numerical));
  127 +
  128 + annotation (experiment(NumberOfIntervals=1), __Dymola_experimentSetupOutput);
  129 +end Validate_Derivatives_SaturationBoundary;
2  Examples/Validate_Derivatives_SinglePhase.mo
@@ -91,7 +91,7 @@ equation
91 91 Modelica.Utilities.Streams.print("(dg/dd)@T=const analytical= " + String(dgdT_analytical));
92 92 Modelica.Utilities.Streams.print("(dg/dd)@T=const numerical= " + String(dgdT_numerical));
93 93 // check (dg/dT)@d=const
94   - dgTd_analytical = f.R*(f.i+f.r+1+f.delta*f.rd -f.tau*(f.it+f.rt) - f.tau*f.delta*f.rtd);
  94 + dgTd_analytical = f.R*(f.i+f.r + 1+f.delta*f.rd -f.tau*(f.it+f.rt) - f.tau*f.delta*f.rtd);
95 95 dgTd_numerical = ((T_plus.h-T_plus.T*T_plus.s)-(T_minus.h-T_minus.T*T_minus.s))/(T_plus.T-T_minus.T);
96 96 Modelica.Utilities.Streams.print("(dg/dT)@d=const analytical= " + String(dgTd_analytical));
97 97 Modelica.Utilities.Streams.print("(dg/dT)@d=const numerical= " + String(dgTd_numerical));
1  Examples/package.order
@@ -14,5 +14,6 @@ Test_inverse
14 14 Test_AncillaryFunctions
15 15 Validate_HelmholtzDerivatives
16 16 Validate_Derivatives_SinglePhase
  17 +Validate_Derivatives_SaturationBoundary
17 18 Validate_Derivatives_TwoPhase
18 19 Tests
2  Interfaces/PartialHelmholtzMedium/Types/DerEnergyByPressure.mo
... ... @@ -0,0 +1,2 @@
  1 +within HelmholtzMedia.Interfaces.PartialHelmholtzMedium.Types;
  2 +type DerEnergyByPressure = Real (final unit="(J/kg)/Pa");
2  Interfaces/PartialHelmholtzMedium/Types/DerEnthalpyByPressure.mo
... ... @@ -0,0 +1,2 @@
  1 +within HelmholtzMedia.Interfaces.PartialHelmholtzMedium.Types;
  2 +type DerEnthalpyByPressure = Real (final unit="(J/kg)/Pa");
2  Interfaces/PartialHelmholtzMedium/Types/DerEntropyByPressure.mo
... ... @@ -0,0 +1,2 @@
  1 +within HelmholtzMedia.Interfaces.PartialHelmholtzMedium.Types;
  2 +type DerEntropyByPressure = Real (final unit="(J/kg.K)/Pa");
10 Interfaces/PartialHelmholtzMedium/Types/package.mo
... ... @@ -1,5 +1,11 @@
1 1 within HelmholtzMedia.Interfaces.PartialHelmholtzMedium;
2 2 package Types
3   -
4   -
  3 +//
  4 + //
  5 + //
  6 + //
  7 + //
  8 + //
  9 + //
  10 + //
5 11 end Types;
3  Interfaces/PartialHelmholtzMedium/Types/package.order
@@ -7,8 +7,11 @@ DerPressureByDensity
7 7 DerPressureByTemperature
8 8 DerEnthalpyByDensity
9 9 DerEnthalpyByTemperature
  10 +DerEnthalpyByPressure
10 11 DerEnergyByDensity
11 12 DerEnergyByTemperature
  13 +DerEnergyByPressure
12 14 DerEntropyByDensity
13 15 DerEntropyByTemperature
  16 +DerEntropyByPressure
14 17 JouleThomsonCoefficient
3  Interfaces/PartialHelmholtzMedium/package.mo
@@ -112,7 +112,8 @@ protected
112 112 T=T,
113 113 phase=1);
114 114 sat.psat := sat.liq.p;
115   - // else: do nothing
  115 + // elseif (phase==1) then
  116 + // print("sat_of_T_EoS has been called from single-phase state","printlog.txt");
116 117 end if;
117 118
118 119 annotation (Documentation(info="<html>

0 comments on commit d9215d4

Please sign in to comment.
Something went wrong with that request. Please try again.