/[ascend]/trunk/models/test/bug567/brayton.a4c
ViewVC logotype

Contents of /trunk/models/test/bug567/brayton.a4c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2687 - (show annotations) (download) (as text)
Fri Mar 1 02:15:09 2013 UTC (11 years, 7 months ago) by jpye
File MIME type: text/x-ascend
File size: 14808 byte(s)
Moving reproducible version of bug 567 into models/test directory.
Valgrind test shows invalid memory access in RemoveRelation of data freed in DestroyInstanceParts.
1 (* ASCEND modelling environment
2 Copyright (C) 2007, 2008, 2009, 2010 Carnegie Mellon University
3
4 This program is free software; you can redistribute it and/or modify
5 it under the terms of the GNU General Public License as published by
6 the Free Software Foundation; either version 2, or (at your option)
7 any later version.
8
9 This program is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 GNU General Public License for more details.
13
14 You should have received a copy of the GNU General Public License
15 along with this program. If not, see <http://www.gnu.org/licenses/>.
16 *)(*
17 This file contains some models of Brayton engines and associated cycles,
18 following the development of Çengel & Boles 'Thermodynamcs: An Engineering
19 Approach, 6th Ed, McGraw-Hill, 2008.
20
21 Author: John Pye
22 *)
23
24 REQUIRE "atoms.a4l";
25 REQUIRE "johnpye/thermo_types.a4c";
26 REQUIRE "johnpye/airprops.a4c";
27 IMPORT "sensitivity/solve";
28
29 (* first some models of air as an ideal gas *)
30
31 MODEL ideal_gas_base;
32 M IS_A molar_weight_constant;
33 c_p IS_A specific_heat_capacity;
34 s IS_A specific_entropy;
35 h IS_A specific_enthalpy;
36 v IS_A specific_volume;
37 T IS_A temperature;
38 p IS_A pressure;
39 R IS_A specific_gas_constant;
40
41 R :== 1{GAS_C} / M;
42 p * v = R * T;
43
44 METHODS
45 METHOD bound_self;
46 s.lower_bound := -5 {kJ/kg/K};
47 END bound_self;
48 END ideal_gas_base;
49
50 (*
51 Ideal air assuming ideal gas and constant cp.
52 *)
53 MODEL simple_ideal_air
54 REFINES ideal_gas_base;
55
56 M :== 28.958600656 {kg/kmol};
57
58 c_p = 1.005 {kJ/kg/K};
59
60 T_ref IS_A temperature_constant;
61 p_ref IS_A pressure_constant;
62
63 s = c_p * ln(T / T_ref) - R * ln(p / p_ref);
64 h = c_p * (T - T_ref);
65
66 T_ref :== 273.15 {K};
67 p_ref :== 1 {bar};
68
69 METHODS
70 METHOD on_load;
71 RUN ClearAll;
72 RUN bound_self;
73 FIX T, p;
74 T := 300 {K};
75 p := 1 {bar};
76 END on_load;
77 END simple_ideal_air;
78
79 (*
80 Ideal air, using a quartic polynomial for c_p as given in Moran & Shapiro
81 4th Ed.
82 *)
83 MODEL ideal_air_poly REFINES ideal_gas_base;
84 M :== 28.958600656 {kg/kmol};
85
86 a[0..4] IS_A real_constant;
87 a[0] :== 3.653;
88 a[1] :== -1.337e-3{1/K};
89 a[2] :== 3.294e-6{1/K^2};
90 a[3] :== -1.913e-9{1/K^3};
91 a[4] :== 0.2763e-12{1/K^4};
92
93 T_ref IS_A temperature_constant;
94 p_ref IS_A pressure_constant;
95 h_ref IS_A real_constant;
96 s_ref IS_A real_constant;
97
98 T_ref :== 300 {K};
99 p_ref :== 1 {bar};
100 h_ref :== -4.40119 {kJ/kg};
101 s_ref :== 0. {kJ/kg/K};
102
103 c_p * M = 1{GAS_C} * SUM[a[i]*T^i | i IN [0..4]];
104
105 (h - h_ref) * M = 1{GAS_C} * SUM[a[i]/(i+1) * T^(i+1) | i IN[0..4]];
106
107 s0 IS_A specific_entropy;
108 s0 = R * (a[0]*ln(T/T_ref) + SUM[a[i]/i * (T - T_ref)^i | i IN[1..4]]) + 1.294559 {kJ/kg/K} + 0.38191663487 {kJ/kg/K};
109
110 s = s0 - R * ln(p/p_ref);
111
112 METHODS
113 METHOD bound_self;
114 RUN ideal_gas_base::bound_self;
115 s0.lower_bound := -1e20 {kJ/kg/K};
116 END bound_self;
117 METHOD on_load;
118 RUN ClearAll;
119 RUN bound_self;
120 FIX T, p;
121 T := 200 {K};
122 p := 1 {bar};
123 END on_load;
124 END ideal_air_poly;
125
126 IMPORT "johnpye/datareader/datareader";
127 MODEL drconf;
128 filename IS_A symbol_constant;
129 format IS_A symbol_constant;
130 parameters IS_A symbol_constant;
131 END drconf;
132
133 MODEL ideal_air_table REFINES ideal_gas_base;
134 M :== 28.958600656 {kg/kmol};
135
136 dataparams IS_A drconf;
137 dataparams.filename :== 'johnpye/idealair.csv';
138 dataparams.format :== 'CSV';
139 dataparams.parameters :== '1,4,6';
140
141 s0 IS_A specific_entropy;
142 u IS_A specific_energy;
143 p_ref IS_A pressure;
144
145 data: datareader(
146 T : INPUT;
147 u, s0 :OUTPUT;
148 dataparams : DATA
149 );
150
151 h = u + R * T;
152
153 s = s0 - R * ln(p/p_ref);
154 END ideal_air_table;
155
156 (*
157 Thermo properties
158 *)
159 MODEL air_state;
160 air IS_A ideal_air_poly;
161 p ALIASES air.p;
162 T ALIASES air.T;
163 h ALIASES air.h;
164 s ALIASES air.s;
165 v ALIASES air.v;
166 METHODS
167 METHOD default;
168 p := 10{bar};
169 p.nominal := 42 {bar};
170 h := 2000 {kJ/kg};
171
172 T := 400 {K};
173 v.nominal := 10 {L/kg};
174 s := 4 {kJ/kg/K};
175 END default;
176 METHOD solve;
177 EXTERNAL do_solve(SELF);
178 END solve;
179 METHOD on_load;
180 RUN default_all;
181 FIX p, h;
182 END on_load;
183 END air_state;
184
185
186 (* a simple connector that includes calculation of steam properties *)
187 MODEL air_node;
188 state IS_A air_state;
189 p ALIASES state.p;
190 h ALIASES state.h;
191 v ALIASES state.v;
192 T ALIASES state.T;
193 s ALIASES state.s;
194 mdot IS_A mass_rate;
195 METHODS
196 METHOD default;
197 mdot.nominal := 2 {kg/s};
198 END default;
199 METHOD solve;
200 EXTERNAL do_solve(SELF);
201 END solve;
202 METHOD on_load;
203 RUN default; RUN reset; RUN values;
204 FIX p,h;
205 END on_load;
206 END air_node;
207
208 MODEL air_equipment;
209 inlet "in: inlet air stream" IS_A air_node;
210 outlet "out: outlet air stream" IS_A air_node;
211
212 inlet.mdot, outlet.mdot ARE_THE_SAME;
213 mdot ALIASES inlet.mdot;
214 END air_equipment;
215
216
217 MODEL compressor REFINES air_equipment;
218 NOTES
219 'block' SELF {Simple model of a compressor using isentropic efficiency}
220 END NOTES;
221
222 dp IS_A delta_pressure;
223 inlet.p + dp = outlet.p;
224
225 outlet_is IS_A air_state;
226 outlet_is.p, outlet.p ARE_THE_SAME;
227
228 outlet_is.s, inlet.s ARE_THE_SAME;
229 eta IS_A fraction;
230
231 r IS_A factor;
232 r * inlet.p = outlet.p;
233
234 eta_eq:eta * (inlet.h - outlet.h) = (inlet.h - outlet_is.h);
235
236 (* work done on the environment, will be negative *)
237 Wdot IS_A energy_rate;
238 Wdot_eq:Wdot * eta = mdot * (inlet.h - outlet_is.h);
239
240 w IS_A specific_energy;
241 w_eq:w = eta * (outlet.h - inlet.h);
242
243 END compressor;
244
245 MODEL compressor_test REFINES compressor;
246 (* no equations here *)
247 METHODS
248 METHOD on_load;
249 FIX inlet.T;
250 FIX inlet.p;
251
252 inlet.T := 300 {K};
253 inlet.p := 1 {bar};
254
255 FIX r;
256 FIX eta;
257 FIX mdot;
258
259 r := 8;
260 eta := 0.8;
261 mdot := 1 {kg/s};
262 END on_load;
263 END compressor_test;
264
265
266
267 MODEL gas_turbine REFINES air_equipment;
268 NOTES
269 'block' SELF {Simple turbine model}
270 END NOTES;
271
272 dp IS_A delta_pressure;
273 inlet.p + dp = outlet.p;
274
275 outlet_is IS_A air_state;
276 outlet_is.p, outlet.p ARE_THE_SAME;
277 outlet_is.s, inlet.s ARE_THE_SAME;
278
279 eta IS_A fraction;
280 eta_eq:eta * (inlet.h - outlet_is.h) = (inlet.h - outlet.h);
281
282 (* work done on the environment, will be positive *)
283 Wdot IS_A energy_rate;
284 Wedot_eq:Wdot = mdot * (inlet.h - outlet.h);
285
286 w IS_A specific_energy;
287 w_eq:w = inlet.h - outlet.h;
288
289 r IS_A factor;
290 r * outlet.p = inlet.p;
291
292 END gas_turbine;
293
294 MODEL gas_turbine_test REFINES gas_turbine;
295 (* no equations here *)
296 METHODS
297 METHOD on_load;
298 FIX inlet.p;
299 FIX inlet.T;
300 FIX outlet.p;
301 FIX eta;
302 FIX mdot;
303
304 inlet.p := 15 {bar};
305 inlet.T := 1200 {K};
306 outlet.p := 1 {bar};
307 eta := 0.85;
308 mdot := 1 {kg/s};
309 END on_load;
310 END gas_turbine_test;
311
312
313
314
315 (*
316 simple model assumes no pressure drop, but heating losses due to
317 flue gas temperature
318 *)
319 MODEL combustor REFINES air_equipment;
320 NOTES
321 'block' SELF {Simple combustion chamber model}
322 END NOTES;
323
324 inlet.p, outlet.p ARE_THE_SAME;
325 Qdot_fuel IS_A energy_rate;
326 Qdot IS_A energy_rate;
327
328 Qdot = mdot * (outlet.h - inlet.h);
329
330 eta IS_A fraction;
331 Qdot = eta * Qdot_fuel;
332
333 qdot_fuel IS_A specific_energy_rate;
334 qdot_fuel * mdot = Qdot_fuel;
335 END combustor;
336
337 MODEL combustor_test REFINES combustor;
338 (* nothing here *)
339 METHODS
340 METHOD on_load;
341 FIX inlet.p;
342 FIX inlet.T;
343 FIX eta;
344 FIX outlet.T;
345 FIX mdot;
346
347 inlet.p := 15 {bar};
348 inlet.T := 500 {K};
349
350 eta := 0.8;
351 outlet.T := 700 {K};
352 mdot := 1 {kg/s};
353 END on_load;
354 END combustor_test;
355
356
357
358 (*
359 this is really simple (fluid props permitting): just work out the heat
360 that must be expelled to get the gas down to a specified temperature
361 *)
362 MODEL dissipator REFINES air_equipment;
363 NOTES
364 'block' SELF {Simple condenser model}
365 END NOTES;
366
367 inlet.p, outlet.p ARE_THE_SAME;
368 Qdot IS_A energy_rate;
369
370 Qdot = mdot * (outlet.h - inlet.h);
371
372 END dissipator;
373
374 MODEL dissipator_test REFINES dissipator;
375 (* nothing here *)
376 METHODS
377 METHOD on_load;
378 FIX inlet.p, inlet.T;
379 FIX outlet.T;
380 FIX mdot;
381
382 inlet.p := 1 {bar};
383 inlet.T := 500 {K};
384 outlet.T := 300 {K};
385 mdot := 1 {kg/s};
386 END on_load;
387 END dissipator_test;
388
389
390 MODEL brayton;
391 NOTES
392 'description' SELF {
393 This is a model of a simple Brayton cycle with
394 irreversible compressor (eta=0.8) and turbine (eta=0.85) operating
395 between 300 K and 1300 K, with a compression ratio of 8 and an
396 assumed inlet pressure of 1 bar. Based on examples 9-5 and 9-6 from
397 Çengel & Boles, 'Thermodynamics: An Engineering Approach', 6th Ed,
398 McGraw-Hill, 2008}
399 END NOTES;
400
401 CO IS_A compressor;
402 TU IS_A gas_turbine;
403 BU IS_A combustor;
404 DI IS_A dissipator;
405
406 CO.outlet, BU.inlet ARE_THE_SAME;
407 BU.outlet, TU.inlet ARE_THE_SAME;
408 TU.outlet, DI.inlet ARE_THE_SAME;
409 DI.outlet, CO.inlet ARE_THE_SAME;
410 braytonpressureratio IS_A positive_factor;
411 braytonpressureratio * CO.inlet.p = TU.outlet.p;
412
413 Wdot_CO ALIASES CO.Wdot;
414 Wdot_TU ALIASES TU.Wdot;
415 Wdot IS_A energy_rate;
416 Wdot = Wdot_CO + Wdot_TU;
417
418 Qdot_BU ALIASES BU.Qdot;
419 Qdot_DI ALIASES DI.Qdot;
420
421 Qdot IS_A energy_rate;
422 Qdot = Qdot_BU + Qdot_DI;
423
424 Edot IS_A energy_rate;
425 Edot = Wdot - Qdot;
426
427 eta IS_A fraction;
428 eta = Wdot / Qdot_BU;
429
430 r_bw IS_A factor;
431 r_bw = -Wdot_CO / Wdot_TU;
432
433 state[1..4] IS_A air_node;
434 state[1], CO.inlet ARE_THE_SAME;
435 state[2], BU.inlet ARE_THE_SAME;
436 state[3], TU.inlet ARE_THE_SAME;
437 state[4], DI.inlet ARE_THE_SAME;
438
439 eta_TU ALIASES TU.eta;
440 eta_CO ALIASES CO.eta;
441
442 METHODS
443 METHOD on_load;
444 FIX CO.eta, TU.eta;
445 CO.eta := 0.8;
446 TU.eta := 0.85;
447 FIX CO.inlet.T, TU.inlet.T;
448 CO.inlet.T := 300 {K};
449 TU.inlet.T := 1300 {K};
450 FIX CO.r;
451 CO.r := 8;
452 FIX CO.inlet.p;
453 CO.inlet.p := 1 {bar};
454 FIX CO.inlet.mdot;
455 CO.inlet.mdot := 1 {kg/s};
456 FIX BU.eta;
457 BU.eta := 1;
458 END on_load;
459 END brayton;
460
461
462 (*
463 Regenerator: air-to-air heat exchanger
464
465 Assumption: fluid on both sides have the same c_p.
466 *)
467 MODEL regenerator REFINES air_equipment;
468 inlet_hot, outlet_hot IS_A air_node;
469
470 inlet.p, outlet.p ARE_THE_SAME;
471 inlet_hot.p, outlet_hot.p ARE_THE_SAME;
472
473 inlet_hot.mdot, outlet_hot.mdot ARE_THE_SAME;
474 mdot_hot ALIASES inlet_hot.mdot;
475
476 (* for perfect eps=1 case: inlet_hot.T, outlet.T ARE_THE_SAME;*)
477
478 epsilon IS_A fraction;
479
480 Qdot IS_A energy_rate;
481 mdot_min IS_A mass_rate;
482 mdot_min = inlet.mdot + 0.5*(inlet.mdot - inlet_hot.mdot + abs(inlet.mdot - inlet_hot.mdot));
483
484 Qdot = epsilon * mdot_min * (inlet_hot.h - inlet.h);
485 outlet.h = inlet.h + Qdot/inlet.mdot;
486 outlet_hot.h = inlet_hot.h - Qdot/inlet_hot.mdot;
487 END regenerator;
488
489 MODEL regenerator_test REFINES regenerator;
490 METHODS
491 METHOD on_load;
492 FIX inlet.mdot, inlet.p, inlet.T;
493 FIX inlet_hot.mdot, inlet_hot.p, inlet_hot.T;
494 inlet.mdot := 1 {kg/s};
495 inlet.p := 1 {bar};
496 inlet.T := 300 {K};
497 inlet_hot.mdot := 1.05 {kg/s};
498 inlet_hot.p := 15 {bar};
499 inlet_hot.T := 500 {K};
500 FIX epsilon;
501 epsilon := 0.8;
502 END on_load;
503 END regenerator_test;
504
505
506
507 MODEL brayton_regenerative;
508 NOTES
509 'description' SELF {
510 This is a model of a regenerative Brayton cycle with
511 irreversible compressor (eta=0.8) and turbine (eta=0.85) operating
512 between 300 K and 1300 K, with a compression ratio of 8 and an
513 assumed inlet pressure of 1 bar. The regenerator effectiveness is
514 0.8.
515
516 Based on example 9-7 from Çengel & Boles, 'Thermodynamics: An
517 Engineering Approach', 6th Ed, McGraw-Hill, 2008}
518 END NOTES;
519
520 CO IS_A compressor;
521 TU IS_A gas_turbine;
522 BU IS_A combustor;
523 DI IS_A dissipator;
524 RE IS_A regenerator;
525
526 CO.outlet, RE.inlet ARE_THE_SAME;
527 RE.outlet, BU.inlet ARE_THE_SAME;
528 BU.outlet, TU.inlet ARE_THE_SAME;
529 TU.outlet, RE.inlet_hot ARE_THE_SAME;
530 RE.outlet_hot, DI.inlet ARE_THE_SAME;
531 DI.outlet, CO.inlet ARE_THE_SAME;
532
533 Wdot_CO ALIASES CO.Wdot;
534 Wdot_TU ALIASES TU.Wdot;
535 Wdot IS_A energy_rate;
536 Wdot = Wdot_CO + Wdot_TU;
537
538 Qdot_BU ALIASES BU.Qdot;
539 Qdot_DI ALIASES DI.Qdot;
540
541 Qdot IS_A energy_rate;
542 Qdot = Qdot_BU + Qdot_DI;
543
544 Edot IS_A energy_rate;
545 Edot = Wdot - Qdot;
546
547 eta IS_A factor;
548 eta = Wdot / Qdot_BU;
549
550 r_bw IS_A factor;
551 r_bw = -Wdot_CO / Wdot_TU;
552
553 Qdot_RE ALIASES RE.Qdot;
554
555 eta_TU ALIASES TU.eta;
556 eta_CO ALIASES CO.eta;
557 epsilon_RE ALIASES RE.epsilon;
558
559 braytonpressureratio ALIASES CO.r;
560 METHODS
561 METHOD on_load;
562 FIX CO.eta, TU.eta;
563 CO.eta := 0.88;
564 TU.eta := 0.85;
565 FIX CO.inlet.T, TU.inlet.T;
566 CO.inlet.T := 300 {K};
567 TU.inlet.T := 1300 {K};
568 FIX CO.r;
569 CO.r := 4.5;
570 FIX CO.inlet.p;
571 CO.inlet.p := 1 {bar};
572 FIX CO.inlet.mdot;
573 CO.inlet.mdot := 1 {kg/s};
574 FIX BU.eta;
575 BU.eta := 1;
576 FIX RE.epsilon;
577 RE.epsilon := 0.8;
578 END on_load;
579 END brayton_regenerative;
580
581
582
583
584 MODEL brayton_intercool_reheat_regen;
585 NOTES
586 'description' SELF {
587 This is a model of a Brayton cycle with intercooling, reheating, and
588 regeneration.
589
590 It has an irreversible compressor (eta=0.8) and turbine (eta=0.85)
591 and operates between 300 K and 1300 K, with a compression ratio of 8
592 and an assumed inlet pressure of 1 bar. The regenerator
593 effectiveness is 0.8.
594
595 In adding the intercooling and reheating stages, we assume an
596 intermediate pressure that results in two compression stages of
597 equal pressure ratio, and two turbine stages of equal pressure
598 ratio.
599
600 Based on example 9-8 from Çengel & Boles, 'Thermodynamics: An
601 Engineering Approach', 6th Ed, McGraw-Hill, 2008}
602 END NOTES;
603
604 CO1, CO2 IS_A compressor;
605 TU1, TU2 IS_A gas_turbine;
606 BU IS_A combustor;
607 DI IS_A dissipator;
608 RE IS_A regenerator;
609 IC IS_A dissipator;
610 RH IS_A combustor;
611
612 CO1.outlet, IC.inlet ARE_THE_SAME;
613 IC.outlet, CO2.inlet ARE_THE_SAME;
614 CO2.outlet, RE.inlet ARE_THE_SAME;
615 RE.outlet, BU.inlet ARE_THE_SAME;
616 BU.outlet, TU1.inlet ARE_THE_SAME;
617 TU1.outlet, RH.inlet ARE_THE_SAME;
618 RH.outlet, TU2.inlet ARE_THE_SAME;
619 TU2.outlet, RE.inlet_hot ARE_THE_SAME;
620 RE.outlet_hot, DI.inlet ARE_THE_SAME;
621 DI.outlet, CO1.inlet ARE_THE_SAME;
622
623 Wdot_CO1 ALIASES CO1.Wdot;
624 Wdot_CO2 ALIASES CO2.Wdot;
625 Wdot_TU1 ALIASES TU1.Wdot;
626 Wdot_TU2 ALIASES TU2.Wdot;
627
628 Wdot_CO, Wdot_TU, Wdot IS_A energy_rate;
629 Wdot_CO = Wdot_CO1 + Wdot_CO2;
630 Wdot_TU = Wdot_TU1 + Wdot_TU2;
631 Wdot = Wdot_CO + Wdot_TU;
632
633 Qdot_BU ALIASES BU.Qdot;
634 Qdot_DI ALIASES DI.Qdot;
635 Qdot_IC ALIASES IC.Qdot;
636 Qdot_RH ALIASES RH.Qdot;
637
638 Qdot IS_A energy_rate;
639 Qdot = Qdot_BU + Qdot_DI + Qdot_IC + Qdot_RH;
640
641 Edot IS_A energy_rate;
642 Edot = Wdot - Qdot;
643
644 eta IS_A factor;
645 eta = Wdot / Qdot_BU;
646
647 CO1.r = CO2.r;
648 TU1.r = TU2.r;
649
650 RH.outlet.T = BU.outlet.T;
651 IC.outlet.T = DI.outlet.T;
652
653 r IS_A factor;
654 r = CO2.outlet.p / CO1.inlet.p;
655
656 r_bw IS_A factor;
657 r_bw = -Wdot_CO / Wdot_TU;
658
659 Qdot_RE ALIASES RE.Qdot;
660
661 eta_TU1 ALIASES TU1.eta;
662 eta_TU2 ALIASES TU2.eta;
663 eta_CO1 ALIASES CO1.eta;
664 eta_CO2 ALIASES CO2.eta;
665 epsilon_RE ALIASES RE.epsilon;
666 METHODS
667 METHOD on_load;
668 FIX CO1.eta, CO2.eta, TU1.eta, TU2.eta;
669 CO1.eta := 0.8;
670 CO2.eta := 0.8;
671 TU1.eta := 0.85;
672 TU2.eta := 0.85;
673 FIX CO1.inlet.T, TU1.inlet.T;
674 CO1.inlet.T := 300 {K};
675 TU1.inlet.T := 1300 {K};
676 FIX r;
677 r := 8;
678 FIX CO1.inlet.p;
679 CO1.inlet.p := 1 {bar};
680 FIX CO1.inlet.mdot;
681 CO1.inlet.mdot := 1 {kg/s};
682 FIX BU.eta, RH.eta;
683 BU.eta := 1;
684 RH.eta := 1;
685 FIX RE.epsilon;
686 RE.epsilon := 0.8;
687 END on_load;
688 END brayton_intercool_reheat_regen;

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