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

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

Parent Directory Parent Directory | Revision Log Revision Log


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

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