/[ascend]/trunk/models/johnpye/iapws95.a4c
ViewVC logotype

Contents of /trunk/models/johnpye/iapws95.a4c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 112 - (show annotations) (download) (as text)
Wed Dec 14 06:25:28 2005 UTC (19 years, 4 months ago) by johnpye
File MIME type: text/x-ascend
File size: 19760 byte(s)
Moving saturation curves into a separate file.
Adding a simple loop simulation using ideal-steam and fixed-density water.
1 REQUIRE "system.a4l";
2 REQUIRE "atoms.a4l";
3 REQUIRE "johnpye/thermo_types.a4c";
4 REQUIRE "johnpye/iapws_sat_curves.a4c";
5
6 (*************************************************************************
7
8 The thermodynamic properties of water calculated with the
9 IAPWS95 equations. Variables and (example/possible) units:
10
11 T Temperature, K
12 rho Density, kg/m^3
13 p Pressure, MPa
14 u Specific internal energy, kJ/kg
15 h Specific enthalpy, kJ/kg
16 s Specific entropy, kJ/(kg*K)
17 cv Isochoric specific heat, kJ/(kg*K)
18 cp Isobaric specific heat, kJ/(kg*K)
19 w Speed of sound, m/s
20
21 References:
22
23 [1] The International Association for the Properties of
24 Water and Steam, "Release on the IAPWS Formulation 1995
25 for the Thermodynamic Properties of Ordinary Water
26 Substance for General and Scientific Use", dated
27 September 1996, Fredericia, Denmark. See the "Releases"
28 section of the website http://www.iapws.org/.
29
30 [2] NIST Chemistry Webbook:
31 http://webbook.nist.gov/chemistry/fluid/
32
33 ----------------------------------------------------------------------
34
35 freesteam-ascend - IAPWS-95 steam library for ASCEND
36 Copyright (C) John Pye 2005
37 derived from work by Don Peterson for freesteam, (C) 2004.
38
39 This program is free software; you can redistribute it
40 and/or modify it under the terms of the GNU General Public
41 License as published by the Free Software Foundation; either
42 version 2 of the License, or (at your option) any later
43 version.
44
45 This program is distributed in the hope that it will be
46 useful, but WITHOUT ANY WARRANTY; without even the implied
47 warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
48 PURPOSE. See the GNU General Public License for more
49 details.
50
51 You should have received a copy of the GNU General Public
52 License along with this program; if not, write to the Free
53 Software Foundation, Inc., 59 Temple Place, Suite 330,
54 Boston, MA 02111-1307 USA
55
56 *)
57
58 MODEL iapws95_1phase REFINES thermo_state;
59
60 delta IS_A positive_variable;
61 tau IS_A positive_variable;
62
63 (*-------------- CONSTANTS ---------------*)
64 rhoc IS_A mass_density_constant;
65 Tc IS_A temperature_constant;
66
67 rhoc "density of water at the critical point"
68 :== 322 {kg/m^3};
69
70 Tc "temperature of water at the critical point"
71 :== 647.096 {K};
72
73 R IS_A specific_gas_constant;
74 R "specific gas constant for water"
75 :== 0.46151805 {kJ/kg/K};
76
77 tau = Tc / T;
78 delta = rho / rhoc;
79
80
81 range_0 IS_A set OF integer_constant;
82 range_0 :== [1..8];
83
84 range_01 IS_A set OF integer_constant;
85 range_01 :== [4..8];
86
87 range_r1 IS_A set OF integer_constant;
88 range_r1 :== [1..7];
89
90 range_r2 IS_A set OF integer_constant;
91 range_r2 :== [8..51];
92
93 range_r3 IS_A set OF integer_constant;
94 range_r3 :== [52..54];
95
96 range_r4 IS_A set OF integer_constant;
97 range_r4 :== [55..56];
98
99 n0[range_0] IS_A real_constant;
100 n0[1] :== -8.32044648201;
101
102 n0[2] :== 6.6832105268;
103 n0[3] :== 3.00632;
104 n0[4] :== 0.012436;
105
106 n0[5] :== 0.97315;
107 n0[6] :== 1.27950;
108 n0[7] :== 0.96956;
109 n0[8] :== 0.24873;
110
111 gamma0[range_01] IS_A real_constant;
112 gamma0[4] :== 1.28728967;
113 gamma0[5] :== 3.53734222;
114 gamma0[6] :== 7.74073708;
115 gamma0[7] :== 9.24437796;
116 gamma0[8] :== 27.5075105;
117
118 n[1..56] IS_A real_constant;
119 n[1] :== 0.12533547935523e-1; n[2] :== 0.78957634722828e+1; n[3] :== -0.87803203303561e+1; n[4] :== 0.31802509345418e+0;
120 n[5] :== -0.26145533859358e+0; n[6] :== -0.78199751687981e-2; n[7] :== 0.88089493102134e-2;
121 n[8] :== -0.66856572307965e+0; n[9] :== 0.20433810950965e+0; n[10] :== -0.66212605039687e-4;
122 n[11] :== -0.19232721156002e+0; n[12] :== -0.25709043003438e+0; n[13] :== 0.16074868486251e+0; n[14] :== -0.40092828925807e-1;
123 n[15] :== 0.39343422603254e-6; n[16] :== -0.75941377088144e-5; n[17] :== 0.56250979351888e-3;
124 n[18] :== -0.15608652257135e-4; n[19] :== 0.11537996422951e-8; n[20] :== 0.36582165144204e-6;
125 n[21] :== -0.13251180074668e-11; n[22] :== -0.62639586912454e-9; n[23] :== -0.10793600908932e+0; n[24] :== 0.17611491008752e-1;
126 n[25] :== 0.22132295167546e+0; n[26] :== -0.40247669763528e+0; n[27] :== 0.58083399985759e+0;
127 n[28] :== 0.49969146990806e-2; n[29] :== -0.31358700712549e-1; n[30] :== -0.74315929710341e+0;
128 n[31] :== 0.47807329915480e+0; n[32] :== 0.20527940895948e-1; n[33] :== -0.13636435110343e+0; n[34] :== 0.14180634400617e-1;
129 n[35] :== 0.83326504880713e-2; n[36] :== -0.29052336009585e-1; n[37] :== 0.38615085574206e-1;
130 n[38] :== -0.20393486513704e-1; n[39] :== -0.16554050063734e-2; n[40] :== 0.19955571979541e-2;
131 n[41] :== 0.15870308324157e-3; n[42] :== -0.16388568342530e-4; n[43] :== 0.43613615723811e-1; n[44] :== 0.34994005463765e-1;
132 n[45] :== -0.76788197844621e-1; n[46] :== 0.22446277332006e-1; n[47] :== -0.62689710414685e-4;
133 n[48] :== -0.55711118565645e-9; n[49] :== -0.19905718354408e+0; n[50] :== 0.31777497330738e+0;
134 n[51] :== -0.11841182425981e+0; n[52] :== -0.31306260323435e+2; n[53] :== 0.31546140237781e+2; n[54] :== -0.25213154341695e+4;
135 n[55] :== -0.14874640856724e+0; n[56] :== 0.31806110878444e+0;
136
137 c[1..51] IS_A integer_constant;
138 c[1] :== 0; c[2] :== 0; c[3] :== 0; c[4] :== 0; c[5] :== 0;
139 c[6] :== 0; c[7] :== 0; c[8] :== 1; c[9] :== 1; c[10] :== 1;
140 c[11] :== 1; c[12] :== 1; c[13] :== 1; c[14] :== 1; c[15] :== 1;
141 c[16] :== 1; c[17] :== 1; c[18] :== 1; c[19] :== 1; c[20] :== 1;
142 c[21] :== 1; c[22] :== 1; c[23] :== 2; c[24] :== 2; c[25] :== 2;
143 c[26] :== 2; c[27] :== 2; c[28] :== 2; c[29] :== 2; c[30] :== 2;
144 c[31] :== 2; c[32] :== 2; c[33] :== 2; c[34] :== 2; c[35] :== 2;
145 c[36] :== 2; c[37] :== 2; c[38] :== 2; c[39] :== 2; c[40] :== 2;
146 c[41] :== 2; c[42] :== 2; c[43] :== 3; c[44] :== 3; c[45] :== 3;
147 c[46] :== 3; c[47] :== 4; c[48] :== 6; c[49] :== 6; c[50] :== 6;
148 c[51] :== 6;
149
150 d[1..54] IS_A integer_constant;
151 d[1] :== 1; d[2] :== 1; d[3] :== 1; d[4] :== 2; d[5] :== 2;
152 d[6] :== 3; d[7] :== 4; d[8] :== 1; d[9] :== 1; d[10] :== 1;
153 d[11] :== 2; d[12] :== 2; d[13] :== 3; d[14] :== 4; d[15] :== 4;
154 d[16] :== 5; d[17] :== 7; d[18] :== 9; d[19] :== 10; d[20] :== 11;
155 d[21] :== 13; d[22] :== 15; d[23] :== 1; d[24] :== 2; d[25] :== 2;
156 d[26] :== 2; d[27] :== 3; d[28] :== 4; d[29] :== 4; d[30] :== 4;
157 d[31] :== 5; d[32] :== 6; d[33] :== 6; d[34] :== 7; d[35] :== 9;
158 d[36] :== 9; d[37] :== 9; d[38] :== 9; d[39] :== 9; d[40] :== 10;
159 d[41] :== 10; d[42] :== 12; d[43] :== 3; d[44] :== 4; d[45] :== 4;
160 d[46] :== 5; d[47] :== 14; d[48] :== 3; d[49] :== 6; d[50] :== 6;
161 d[51] :== 6; d[52] :== 3; d[53] :== 3; d[54] :== 3;
162
163 t[1..54] IS_A real_constant;
164 t[1] :== -0.5; t[2] :== 0.875; t[3] :== 1; t[4] :== 0.5; t[5] :== 0.75;
165 t[6] :== 0.375; t[7] :== 1; t[8] :== 4; t[9] :== 6; t[10] :== 12;
166 t[11] :== 1; t[12] :== 5; t[13] :== 4; t[14] :== 2; t[15] :== 13;
167 t[16] :== 9; t[17] :== 3; t[18] :== 4; t[19] :== 11; t[20] :== 4;
168 t[21] :== 13; t[22] :== 1; t[23] :== 7; t[24] :== 1; t[25] :== 9;
169 t[26] :== 10; t[27] :== 10; t[28] :== 3; t[29] :== 7; t[30] :== 10;
170 t[31] :== 10; t[32] :== 6; t[33] :== 10; t[34] :== 10; t[35] :== 1;
171 t[36] :== 2; t[37] :== 3; t[38] :== 4; t[39] :== 8; t[40] :== 6;
172 t[41] :== 9; t[42] :== 8; t[43] :== 16; t[44] :== 22; t[45] :== 23;
173 t[46] :== 23; t[47] :== 10; t[48] :== 50; t[49] :== 44; t[50] :== 46;
174 t[51] :== 50; t[52] :== 0; t[53] :== 1; t[54] :== 4;
175
176 (* Correlation parameters *)
177
178 (* TODO convert from C to ASCEND arrays? Note trickiness with 0- and 1-based array indices. *)
179
180 a[range_r4] IS_A real_constant;
181 a[55]:==3.5;
182 a[56]:==3.5;
183
184 b[range_r4] IS_A real_constant;
185 b[55]:== 0.85;
186 b[56]:== 0.95;
187
188 A[range_r4] IS_A real_constant;
189 A[55]:==0.32;
190 A[56]:==0.32;
191
192 B[range_r4] IS_A real_constant;
193 B[55]:==0.2;
194 B[56]:==0.2;
195
196 C[range_r4] IS_A real_constant;
197 C[55]:==28;
198 C[56]:==32;
199
200 D[range_r4] IS_A real_constant;
201 D[55]:==700;
202 D[56]:==800;
203
204 beta_r4[range_r4] IS_A real_constant;
205 beta_r4[55]:==0.3;
206 beta_r4[56]:==0.3;
207
208 alpha[range_r3] IS_A integer_constant;
209 alpha[52]:==20;
210 alpha[53]:==20;
211 alpha[54]:==20;
212
213 beta[range_r3] IS_A real_constant;
214 beta[52]:==150;
215 beta[53]:==150;
216 beta[54]:==250;
217
218 gamma[range_r3] IS_A real_constant;
219 gamma[52]:==1.21;
220 gamma[53]:==1.21;
221 gamma[54]:==1.25;
222
223 (*--------------- DIMENSIONLESS PARTIAL DERIVATIVES ---------------- *)
224
225 (*------------ IDEAL PARTS ------------*)
226
227 phi0 IS_A factor;
228 phi0_expr: phi0 =
229 SUM[ n0[i]*ln(1-exp(-tau*gamma0[i])) | i IN [range_01] ]
230 + ln(delta) + n0[1] + n0[2]*tau + n0[3]*ln(tau);
231
232 phi0delta IS_A factor;
233 phi0delta_expr: phi0delta = 1.0/delta;
234
235 phi0deltadelta IS_A factor;
236 phi0deltadelta_expr: phi0deltadelta =
237 -1.0/(delta*delta);
238
239
240 phi0tau IS_A factor;
241 phi0tau_expr: phi0tau =
242 n0[2] + n0[3]/tau
243 + SUM[ n0[i]*gamma0[i]*(1/(1-exp(-tau*gamma0[i])) - 1) | i IN [range_01] ];
244
245 phi0deltatau IS_A real_constant;
246 phi0deltatau :== 0.0;
247
248 phi0tautau IS_A factor;
249 phi0tautau_expr: phi0tautau
250 = -n0[3] / tau^2
251 - SUM [ n0[i] * gamma0[i]^2 * exp(-gamma0[i] * tau) / ( 1 - exp(-gamma0[i] * tau) )^2 | i IN range_01 ];
252
253 (*----------- 'REAL' PARTS -- CLOSE YOUR EYES -----------*)
254
255 d1 IS_A factor;
256 d1_expr: d1 = delta - 1;
257
258 t1 IS_A factor;
259 t1_expr: t1 = tau -1;
260
261 r3_b1[range_r3] IS_A factor;
262 FOR i IN range_r3 CREATE
263 r3_b1[i] = -alpha[i]*d1^2 - beta[i]* (tau - gamma[i])^2;
264 END FOR;
265
266 PSI[range_r4] IS_A factor;
267 theta[range_r4] IS_A factor;
268 DELTA[range_r4] IS_A factor;
269 dDELTA_ddelta[range_r4] IS_A factor;
270 dPSI_ddelta[range_r4] IS_A factor;
271 dDELTAbi_ddelta[range_r4] IS_A factor;
272 dDELTAbi_dtau[range_r4] IS_A factor;
273 dPSI_dtau[range_r4] IS_A factor;
274 d2DELTA_ddelta2[range_r4] IS_A factor;
275 d2DELTAbi_ddelta2[range_r4] IS_A factor;
276 d2PSI_ddelta2[range_r4] IS_A factor;
277 d2DELTAbi_dtau2[range_r4] IS_A factor;
278 d2PSI_dtau2[range_r4] IS_A factor;
279 d2PSI_ddeltadtau[range_r4] IS_A factor;
280 d2DELTAbi_ddeltadtau[range_r4] IS_A factor;
281
282 FOR i IN range_r4 CREATE
283 PSI[i] = exp(-C[i]*d1^2 - D[i]*t1^2);
284 (*theta_expr:*) theta[i] = 1 - tau + A[i] * (d1^2)^( 1/(2 * beta_r4[i]) );
285 (*DELTA_expr:*) DELTA[i] = theta[i]^2 + B[i] * (d1^2)^a[i];
286
287 (*dDELTA_ddelta_expr:*) dDELTA_ddelta[i] = d1*(A[i]*theta[i]*2/beta_r4[i]*
288 (d1^2)^(1/(2*beta_r4[i]) - 1) + 2*B[i]*a[i]*
289 (d1^2)^(a[i] - 1) );
290
291 (*dPSI_ddelta_expr:*) dPSI_ddelta[i] = -2*C[i]*d1*PSI[i];
292
293 (*dDELTAbi_ddelta_expr:*) dDELTAbi_ddelta[i] = b[i] * DELTA[i]^(b[i]-1) * dDELTA_ddelta[i];
294
295 (*dDELTAbi_dtau_expr:*) dDELTAbi_dtau[i] = -2 * theta[i] * b[i] * DELTA[i]^(b[i]-1);
296 (*dPSI_dtau_expr:*) dPSI_dtau[i] = -2 * D[i] * t1 * PSI[i];
297
298 (*d2PSI_ddelta2_expr:*) d2PSI_ddelta2[i] = (2 * C[i] * d1^2 - 1) * 2 * C[i] * PSI[i];
299
300 (*d2DELTA_ddelta2_expr:*) d2DELTA_ddelta2[i] = 1/d1*dDELTA_ddelta[i] + d1^2*(4*B[i]*a[i]*
301 (a[i]-1)*(d1^2)^(a[i]-2) + 2*A[i]^2*
302 (1/(beta_r4[i]^2))*((d1^2) ^ (1/(2*beta_r4[i])-1))^2) +
303 A[i]*theta[i]*4/beta_r4[i]*(1/(2*beta_r4[i]) - 1)*
304 (d1^2)^(1/(2*beta_r4[i]) - 2);
305
306 (*d2DELTAbi_ddelta2_expr:*) d2DELTAbi_ddelta2[i] =b[i] * (
307 DELTA[i]^(b[i]-1) *d2DELTA_ddelta2[i]
308 + (b[i]-1) * DELTA[i]^(b[i]-2) * dDELTA_ddelta[i]^2 );
309
310 (*d2DELTAbi_dtau2_expr:*) d2DELTAbi_dtau2[i] = 2 * b[i] * DELTA[i]^(b[i]-1) +
311 4 * theta[i]^2 * b[i]*(b[i]-1) * DELTA[i]^(b[i]-2);
312
313 (*d2PSI_dtau2_expr:*) d2PSI_dtau2[i] = (2 * D[i] * t1^2 - 1) * 2 * D[i] * PSI[i];
314
315 (*d2PSI_ddeltadtau_expr:*) d2PSI_ddeltadtau[i] = 4 * C[i] * D[i] * d1 * t1 * PSI[i];
316
317 (*d2DELTAbi_ddeltadtau_expr:*) d2DELTAbi_ddeltadtau[i] =
318 - A[i] * b[i]* (2/beta_r4[i]) * DELTA[i]^(b[i]-1) * d1 * (d1^2)^(1/(2*beta_r4[i])-1)
319 - 2 * theta[i] * b[i] *(b[i]-1) * DELTA[i]^(b[i]-2) * dDELTA_ddelta[i];
320
321 END FOR;
322
323 phir_r2[range_r2] IS_A factor;
324 FOR i IN range_r2 CREATE
325 phir_r2[i] = n[i] * delta^d[i] * tau^t[i] *
326 exp(-delta^c[i]);
327 END FOR;
328
329 phir IS_A factor;
330 phir_expr: phir
331 =
332 SUM[ n[i] * delta^d[i] * tau^t[i] | i IN [range_r1] ]
333 + SUM[ phir_r2[i] | i IN [range_r2] ]
334 + SUM[ n[i] * delta^d[i] * tau^t[i] * exp(
335 -alpha[i]*d1^2 - beta[i]*
336 (tau - gamma[i])*(tau - gamma[i])
337 ) | i IN [range_r3] ]
338 + SUM[ n[i] * DELTA[i]^b[i] * delta * PSI[i] | i IN [range_r4] ];
339
340
341 phirdelta_r2[range_r2] IS_A factor;
342 FOR i IN range_r2 CREATE
343 phirdelta_r2[i] = n[i] * exp(-delta^c[i]) * (delta^(d[i]-1) *
344 tau^t[i] * (d[i] - c[i] * delta^c[i]) );
345 END FOR;
346
347 phirdelta IS_A factor;
348 phirdelta_expr: phirdelta =
349 SUM[ n[i] * d[i] * delta^(d[i] - 1) * tau^t[i] | i IN range_r1 ]
350 + SUM[ phirdelta_r2[i] | i IN range_r2 ]
351 + SUM[ n[i]*delta^d[i] * tau^t[i] * exp( r3_b1[i] ) *
352 (d[i]/delta - 2*alpha[i]*d1) | i IN range_r3 ]
353 + SUM[ n[i] * (
354 DELTA[i]^b[i] * (PSI[i] + delta * dPSI_ddelta[i] )
355 + dDELTAbi_ddelta[i] * delta * PSI[i] ) | i IN range_r4 ];
356
357 phirtau_r2[range_r2] IS_A factor;
358 FOR i IN range_r2 CREATE
359 phirtau_r2[i] = n[i]*t[i]*delta^d[i]*tau^(t[i]-1) * exp(-delta^c[i]);
360 END FOR;
361
362 phirtau IS_A factor;
363 phirtau_expr: phirtau =
364 SUM[ n[i] * t[i] * delta^d[i] * tau^(t[i]-1) | i IN range_r1 ]
365 + SUM[ phirtau_r2[i] | i IN range_r2 ]
366 + SUM[ n[i] * delta^d[i] * tau^t[i] * exp( r3_b1[i] )*
367 (t[i]/tau - 2*beta[i]*(tau - gamma[i])) | i IN range_r3 ]
368 + SUM[ n[i]*delta*(dDELTAbi_dtau[i] * PSI[i] +
369 DELTA[i]^b[i]*dPSI_dtau[i]) | i IN range_r4 ];
370
371 phirdeltadelta_r2[range_r2] IS_A factor;
372 FOR i IN range_r2 CREATE
373 phirdeltadelta_r2[i] =
374 n[i] * exp(-delta^c[i]) * (
375 delta^(d[i]-2) *
376 tau^t[i] * (
377 ( d[i] - c[i]*delta^c[i] )*(d[i]- 1 - c[i] * delta^c[i] )
378 - c[i]^2 * delta^c[i]
379 )
380 );
381 END FOR;
382
383 phirdeltadelta IS_A factor;
384 phirdeltadelta_expr: phirdeltadelta =
385 SUM[
386 n[i] * d[i] * (d[i] - 1) * delta^(d[i]-2) * tau^t[i] | i IN range_r1 ]
387 + SUM [phirdeltadelta_r2[i] | i IN range_r2 ]
388 + SUM [
389 n[i]* tau^t[i] * exp( r3_b1[i] ) * (
390 - 2 * alpha[i] * delta^d[i]
391 + 4 * alpha[i]*alpha[i] * delta^d[i] * d1^2 (* d1 = d-1 = d-eps *)
392 - 4 * d[i] * alpha[i] * delta^(d[i]-1) * d1
393 + d[i] * (d[i]-1) * delta^(d[i]-1)
394 ) | i IN range_r3 ]
395 + SUM [n[i]*(
396 DELTA[i]^b[i] * (2*dPSI_ddelta[i] + delta*d2PSI_ddelta2[i] )
397 + 2*dDELTAbi_ddelta[i] * (PSI[i] + delta*dPSI_ddelta[i])
398 + d2DELTAbi_ddelta2[i] * delta * PSI[i] ) | i IN range_r4 ];
399
400
401
402 phirtautau_r2[range_r2] IS_A factor;
403 FOR i IN range_r2 CREATE
404 phirtautau_r2[i] = n[i]*t[i]*(t[i]-1)* delta^d[i] *
405 tau^(t[i]-2) * exp(- delta^c[i] );
406 END FOR;
407
408 phirtautau IS_A factor;
409 phirtautau_expr: phirtautau =
410 SUM[ n[i] * t[i] * (t[i]-1) * delta^d[i] *
411 tau^(t[i]-2) | i IN range_r1 ]
412 + SUM[ phirtautau_r2[i] | i IN range_r2 ]
413 + SUM[ n[i] * delta^d[i] * tau^t[i] * exp( r3_b1[i] ) *
414 ( (t[i]/tau - 2*beta[i]* (tau - gamma[i]) )^2 -
415 t[i]/tau^2 - 2*beta[i] ) | i IN range_r3 ]
416 + SUM[ n[i] * delta * (d2DELTAbi_dtau2[i] * PSI[i] +
417 2 * dDELTAbi_dtau[i] * dPSI_dtau[i]
418 + DELTA[i]^b[i] * d2PSI_dtau2[i] ) | i IN range_r4 ];
419
420
421 phirdeltatau_r2[range_r2] IS_A factor;
422 FOR i IN range_r2 CREATE
423 phirdeltatau_r2[i] =
424 n[i]*t[i] * delta^(d[i]-1) * tau^(t[i]-1)
425 * (d[i] - c[i] * delta^c[i]) * exp(-delta^c[i]);
426 END FOR;
427
428 phirdeltatau IS_A factor;
429 phirdeltatau_expr: phirdeltatau =
430 SUM[ n[i]*d[i]*t[i] * delta^(d[i]-1) * tau^(t[i]-1) | i IN range_r1 ]
431
432 + SUM[ phirdeltatau_r2[i] | i IN range_r2 ]
433 + SUM[ n[i] * delta^d[i] * tau^t[i] * exp( r3_b1[i] ) *
434 (d[i]/delta - 2*alpha[i]*d1)
435 * (t[i]/tau - 2*beta[i]* (tau - gamma[i]) ) | i IN range_r3 ]
436 + SUM[ n[i]*( DELTA[i]^b[i] * (dPSI_dtau[i] + delta *
437 d2PSI_ddeltadtau[i] ) + delta * dDELTAbi_ddelta[i] * dPSI_dtau[i] +
438 dDELTAbi_dtau[i] * (PSI[i] + delta * dPSI_ddelta[i] ) +
439 d2DELTAbi_ddeltadtau[i] * delta * PSI[i] ) | i IN range_r4 ];
440
441 (*--------- THERMO PROPERTY RELATIONS ----------- *)
442
443 pressure: p
444 = rho * R * T * (1 + delta*phirdelta);
445
446 internal_energy: u
447 = R * T * tau * (phi0tau + phirtau);
448
449 enthalpy: h
450 = R * T * (1 + tau*(phi0tau + phirtau) + delta*phirdelta);
451
452 entropy: s
453 = R * (tau*(phi0tau + phirtau) - phi0 - phir);
454
455 c_isochoric: cv
456 = - R * tau^2 * (phi0tautau + phirtautau);
457
458 c_isobaric: cp
459 = - R * (
460 tau^2 * (phi0tautau + phirtautau)
461 + ( ( 1 + delta*phirdelta - delta*tau*phirdeltatau )^2 )
462 / ( 1 + 2*delta*phirdelta + delta^2 * phirdeltadelta )
463 );
464
465 spd_sound: w^2
466 = R * T * ( 1 + 2*delta*phirdelta + delta^2 * phirdeltadelta -
467 ( 1 + delta*phirdelta - delta*tau*phirdeltatau )^2
468 / ( tau^2 * (phi0tautau + phirtautau) )
469 );
470
471 METHODS
472 METHOD default_self;
473 RUN ClearAll;
474 RUN specify;
475 RUN values;
476 END default_self;
477
478 METHOD specify;
479 T.fixed := TRUE;
480 rho.fixed := TRUE;
481 END specify;
482
483 METHOD values;
484 (* these are the test values from page 14 of the IAPWS-95 release *)
485 T := 500 {K};
486 rho := 838.025 {kg/m^3};
487 END values;
488
489 END iapws95_1phase;
490
491 MODEL iapws95_2phase REFINES thermo_state;
492 Sf IS_A iapws95_1phase;
493 Sg IS_A iapws95_1phase;
494 sat IS_A iapws_sat_curves;
495 sat.rhof, Sf.rho ARE_THE_SAME;
496 sat.rhog, Sg.rho ARE_THE_SAME;
497 Sf.T, Sg.T, sat.T ARE_THE_SAME;
498 rhof ALIASES Sf.rho;
499 rhog ALIASES Sg.rho;
500 T, Sf.T ARE_THE_SAME;
501 p, Sf.p ARE_THE_SAME;
502
503 x IS_A fraction;
504
505 rhog*rhof = rho*rhog + rho*(rhof - rhog) *x;
506
507 u = Sf.u + (Sg.u -Sf.u) *x;
508
509 h = Sf.h + (Sg.h -Sf.h) *x;
510 s = Sf.s + (Sg.s -Sf.s) *x;
511 cp = Sf.cp + (Sg.cp -Sf.cp) *x;
512 cv = Sf.cv + (Sg.cv -Sf.cv) *x;
513 w = Sf.w + (Sg.w - Sf.w) *x; (* check this *)
514 METHODS
515
516 METHOD default_self;
517 RUN ClearAll;
518 RUN specify;
519 RUN values;
520 END default_self;
521
522 METHOD specify;
523 T.fixed := TRUE;
524 x.fixed := TRUE;
525 END specify;
526
527 METHOD values;
528 T := 390 {K};
529 x := 0.5;
530 END values;
531
532 END iapws95_2phase;
533
534 MODEL iapws95 REFINES thermo_state;
535 S1 IS_A iapws95_1phase;
536 S2 IS_A iapws95_2phase;
537 sat IS_A iapws_sat_curves;
538
539 T, sat.T, S1.T, S2.T ARE_THE_SAME;
540 rho, S1.rho, S2.rho ARE_THE_SAME;
541 h, S1.h, S2.h ARE_THE_SAME;
542 u, S1.u, S2.u ARE_THE_SAME;
543 p, S1.p, S2.p ARE_THE_SAME;
544 s, S1.s, S2.s ARE_THE_SAME;
545 cp, S1.cp, S2.cp ARE_THE_SAME;
546 cv, S1.cv, S2.cv ARE_THE_SAME;
547 w, S1.w, S2.w ARE_THE_SAME;
548
549 CONDITIONAL
550 sub_test: rho >= sat.rhof;
551 super_test: rho <= sat.rhog;
552 sat1_test: rho < sat.rhof;
553 sat2_test: rho > sat.rhog;
554 END CONDITIONAL;
555
556 issat1 IS_A boolean_var;
557 sat1_expr: issat1 == SATISFIED(sat1_test);
558 issat2 IS_A boolean_var;
559 sat2_expr: issat1 == SATISFIED(sat2_test);
560 issat IS_A boolean_var;
561 saturation_expr: issat == issat1 OR issat2;
562 issub IS_A boolean_var;
563 subcooled_expr: issub == SATISFIED(sub_test);
564 issuper IS_A boolean_var;
565 superheated_expr: issuper == SATISFIED(super_test);
566
567 saturation_when: WHEN (issat)
568 CASE TRUE:
569 USE S2;
570 CASE FALSE:
571 USE S1;
572 END WHEN;
573
574 METHODS
575 METHOD default_self;
576 RUN ClearAll;
577 RUN specify;
578 RUN values;
579 END default_self;
580
581 METHOD specify;
582 T.fixed := TRUE;
583 rho.fixed := TRUE;
584 END specify;
585
586 METHOD values;
587 (* fixed vars *)
588 T := 300 {K};
589 rho := 100 {kg/m^3};
590 (* initial values *)
591 issat := TRUE;
592 issuper := FALSE;
593 issub := FALSE;
594 p := 1 {atm};
595 h := 400 {kJ/kg};
596 END values;
597
598 END iapws95;

john.pye@anu.edu.au
ViewVC Help
Powered by ViewVC 1.1.22