/[ascend]/trunk/models/steam/dsgsat3.a4c
ViewVC logotype

Contents of /trunk/models/steam/dsgsat3.a4c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1277 - (show annotations) (download) (as text)
Sun Feb 4 07:39:31 2007 UTC (13 years, 6 months ago) by johnpye
File MIME type: text/x-ascend
File size: 6834 byte(s)
All of TestSteam passes.
Have increased zone of convergence for dsgsat3 model by changing to some central differences. Randomly.
1 REQUIRE "ivpsystem.a4l";
2 REQUIRE "atoms.a4l";
3 REQUIRE "johnpye/thermo_types.a4c";
4
5 (*
6 An attempt to model direct steam generation in pipe flow, limited to the
7 saturated regime, and with constant-valued friction factor. External heat
8 loss is also simplified.
9 *)
10 REQUIRE "steam/satsteamstream.a4c";
11
12 MODEL dsgsat3;
13 n IS_A integer_constant;
14 n :== 5;
15
16 (* temporal derivatives *)
17 drho_dt[2..n] IS_A density_rate;
18 dmdot_dt[2..n] IS_A mass_rate_rate;
19 du_dt[2..n] IS_A specific_energy_rate;
20 dTw_dt[2..n] IS_A temperature_rate;
21
22 (* wall properties *)
23 rho_w IS_A mass_density;
24 D, D_2 IS_A distance;
25 c_w IS_A specific_heat_capacity;
26 A, A_w IS_A area;
27 h_int IS_A heat_transfer_coefficient; (* internal *)
28 h_ext IS_A heat_transfer_coefficient; (* external *)
29 z_A: A = 1{PI}*D^2/4;
30 z_Aw: A_w = 1{PI}*(D_2^2 - D^2)/4;
31 dz IS_A distance;
32 L IS_A distance;
33 z_dz: dz = L / (n - 1);
34
35 (* fluid properties *)
36 node[1..n] IS_A satsteamstream;
37
38 (* flow properties *)
39 vel[1..n] IS_A speed;
40 T_w[2..n] IS_A temperature;
41 T[1..n] IS_A temperature;
42
43 (* constants, for the moment: *)
44 f IS_A positive_factor;
45 (* mu_f IS_A viscosity; *)
46 T_amb IS_A temperature;
47
48 (* system dynamics *)
49 qdot_t[2..n], qdot_l[2..n] IS_A power_per_length;
50 qdot_s IS_A power_per_length;
51
52 FOR i IN [1..n] CREATE
53 z_vel[i]: vel[i] = node[i].v*node[i].mdot/A;
54 END FOR;
55
56 (* some aliases just for easier review of the state of the model *)
57 x[1..n] IS_A fraction;
58 mdot[1..n] IS_A mass_rate;
59 p[1..n] IS_A pressure;
60 FOR i IN [1..n] CREATE
61 x[i], node[i].x ARE_THE_SAME;
62 mdot[i], node[i].mdot ARE_THE_SAME;
63 p[i], node[i].p ARE_THE_SAME;
64 T[i], node[i].T ARE_THE_SAME;
65 END FOR;
66
67 (* differential equations *)
68
69 (* use central difference for several derivatives *)
70 FOR i IN [2..n-1] CREATE
71 z_massbal[i]: A * drho_dt[i] * dz = - (node[i+1].mdot - node[i-1].mdot) / 2.0;
72
73 z_enbal[i]: dz * (qdot_t[i] - node[i].rho * A * du_dt[i]) =
74 + node[i].mdot * (node[i].u - node[i-1].u)
75 + (node[i+1].p*node[i+1].v*node[i+1].mdot - node[i-1].p*node[i-1].v*node[i-1].mdot) / 2;
76
77 z_mombal[i]: - dz/A*dmdot_dt[i] =
78 (node[i+1].p-node[i-1].p) / 2.0
79 + dz * f/D/2 * node[i].rho * vel[i]^2
80 + (node[i].rho*vel[i]^2 - node[i-1].rho*vel[i-1]^2);
81
82 z_wall[i]: rho_w*A_w*c_w*dTw_dt[i] = qdot_s - qdot_l[i] - qdot_t[i];
83 z_loss[i]: qdot_l[i] = h_ext*(1{PI}*D_2)*(T_w[i] - T_amb);
84 z_trans[i]: qdot_t[i] = h_int*(1{PI}*D) *(T_w[i] - node[i].T);
85
86 END FOR;
87
88 (* the last equations use backwards difference *)
89 FOR i IN [n] CREATE
90 z_massbal1[i]: A * drho_dt[i] * dz = - (node[i].mdot - node[i-1].mdot);
91
92 z_enbal1[i]: dz * (qdot_t[i] - node[i].rho * A * du_dt[i]) =
93 + node[i].mdot * (node[i].u - node[i-1].u)
94 + (node[i].p*node[i].v*node[i].mdot - node[i-1].p*node[i-1].v*node[i-1].mdot);
95
96 z_mombal1[i]: - dz/A*dmdot_dt[i] =
97 (node[i].p-node[i-1].p)
98 + dz * f/D/2 * node[i].rho * vel[i]^2
99 + (node[i].rho*vel[i]^2 - node[i-1].rho*vel[i-1]^2);
100
101 z_wall1[i]: rho_w*A_w*c_w*dTw_dt[i] = qdot_s - qdot_l[i] - qdot_t[i];
102 z_loss1[i]: qdot_l[i] = h_ext*(1{PI}*D_2)*(T_w[i] - T_amb);
103 z_trans1[i]: qdot_t[i] = h_int*(1{PI}*D) *(T_w[i] - node[i].T);
104
105 END FOR;
106
107 t IS_A time;
108 METHODS
109 METHOD bound_self;
110 vel[1..n].upper_bound := 100 {m/s};
111 qdot_l[2..n].lower_bound := 0 {W/m};
112 FOR i IN [1..n] DO
113 RUN node[i].bound_self;
114 END FOR;
115 END bound_self;
116 METHOD default_self;
117 D := 0.06 {m};
118 D_2 := 0.07 {m};
119 A_w := 0.25{PI}*D_2^2 -0.25{PI}*D^2;
120 A := 1 {m^2};
121 L := 50 {m};
122 T_amb := 298 {K};
123 dz := 3.08 {m};
124 f := 0.01;
125 h_ext := 10 {W/m^2/K};
126 h_int := 100 {W/m^2/K};
127 qdot_s := 700 {W/m};
128 rho_w := 1000 {kg/m^3};
129 t := 0 {s};
130 FOR i IN [1..n] DO
131 T[i] := 298 {K};
132 vel[i] := 1 {m/s};
133 node[i].v := 0.2 {L/kg};
134 node[i].rho := 6 {kg/L};
135 node[i].dp_dT := +0.5 {kPa/K};
136 node[i].p := 5 {bar};
137 RUN node[i].default_self;
138 END FOR;
139 FOR i IN [2..n] DO
140 T_w[i] := 298 {K};
141 drho_dt[i] := 0 {kg/m^3/s};
142 dmdot_dt[i] := 0 {kg/s/s};
143 du_dt[i] := 0 {kJ/kg/s};
144 dTw_dt[i] := 0 {K/s};
145 qdot_t[i] := 0 {W/m};
146 qdot_l[i] := 0 {W/m};
147 END FOR;
148 END default_self;
149 METHOD specify;
150 (* change to a proper steady-state problem, with fluid properties FREEd *)
151 FOR i IN [1..n] DO
152 RUN node[i].specify;
153 FIX dTw_dt[i]; FREE T_w[i];
154 END FOR;
155 FIX node[1].p;
156 FREE node[1].T;
157 FIX qdot_s;
158 FIX D, D_2, L;
159 FIX h_int, c_w, rho_w, h_ext;
160 FIX f;
161 (* FIX mu_f; *)
162 FIX T_amb;
163 (* fix derivatives to zero *)
164 FOR i IN [2..n] DO
165 (* FIX dmdot_dt[i]; FREE node[i].mdot; *)
166 FREE node[i].x; FIX node[i].p;
167 FIX drho_dt[i]; FREE node[i].p;
168 FIX du_dt[i]; FREE node[i].T;
169 FREE mdot[i]; FIX dmdot_dt[i];
170 END FOR;
171 END specify;
172 METHOD values;
173 L := 9 {m};
174 h_int := 100 {W/m^2/K};
175 h_ext := 20 {W/m^2/K};
176 f := 0.01;
177 node[1].mdot := 0.5 {kg/s};
178 node[1].p := 8 {bar};
179 node[1].x := 0.3;
180 qdot_s := 1000 {W/m^2} * D_2 * 10;
181 FOR i IN [2..n] DO
182 dmdot_dt[i] := 0.0 {kg/s/s};
183 du_dt[i] := 0 {kJ/kg/s};
184 END FOR;
185 END values;
186 METHOD on_load;
187 RUN configure_steady;
188 RUN ode_init;
189 END on_load;
190 (*---------------- a physically sensible steady-state configuration-----------*)
191 METHOD configure_steady;
192 RUN default_self;
193 RUN reset;
194 RUN bound_steady;
195 RUN values;
196 END configure_steady;
197 METHOD bound_steady;
198 RUN bound_self;
199 T_w[2..n].upper_bound := 1000 {K};
200 END bound_steady;
201 (*------------------------- the dynamic problem ------------------------------*)
202 METHOD configure_dynamic;
203 FOR i IN [2..n] DO
204 FREE drho_dt[i]; FIX node[i].rho;
205 FREE dmdot_dt[i]; FIX node[i].mdot;
206 FREE du_dt[i]; FIX node[i].u;
207 FREE dTw_dt[i]; FIX T_w[i];
208 FREE node[i].x;
209 FREE node[i].T;
210 END FOR;
211 t := 0 {s};
212 END configure_dynamic;
213 METHOD free_states;
214 FOR i IN [2..n] DO
215 FREE node[i].rho;
216 FREE node[i].mdot;
217 FREE node[i].u;
218 FREE T_w[i];
219 END FOR;
220 END free_states;
221 METHOD ode_init;
222 (* add the necessary meta data to allow solving with the integrator *)
223 t.ode_type := -1;
224
225 FOR i IN [2..n] DO
226 drho_dt[i].ode_id := 4*i; node[i].rho.ode_id := 4*i;
227 drho_dt[i].ode_type := 2; node[i].rho.ode_type := 1;
228
229 dmdot_dt[i].ode_id := 4*i+1; node[i].mdot.ode_id := 4*i+1;
230 dmdot_dt[i].ode_type := 2; node[i].mdot.ode_type := 1;
231
232 du_dt[i].ode_id := 4*i+2; node[i].u.ode_id := 4*i+2;
233 du_dt[i].ode_type := 2; node[i].u.ode_type := 1;
234
235 dTw_dt[i].ode_id := 4*i+3; T_w[i].ode_id := 4*i+3;
236 dTw_dt[i].ode_type := 2; T_w[i].ode_type := 1;
237 END FOR;
238
239 FOR i IN [1..n] DO
240 p[i].obs_id := 1 + 10*i;
241 x[i].obs_id := 2 + 10*i;
242 node[i].mdot.obs_id := 4 + 10*i;
243 END FOR;
244 FOR i IN [2..n] DO
245 qdot_t[i].obs_id := 3 + 10*i;
246 T_w[i].obs_id := 5 + 10*i;
247 T[i].obs_id := 6 + 10*i;
248 END FOR;
249 END ode_init;
250 METHOD fix_outlet_quality;
251 FIX x[n];
252 FREE node[1].mdot;
253 END fix_outlet_quality;
254
255 END dsgsat3;
256 ADD NOTES IN dsgsat2;
257 'QRSlv' iterationlimit {50}
258 END NOTES;

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