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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1383 - (hide annotations) (download) (as text)
Fri Apr 6 10:50:41 2007 UTC (15 years, 4 months ago) by jpye
File MIME type: text/x-ascend
File size: 8818 byte(s)
Substituted stationary momentum in dsgsat3.a4c
1 johnpye 1179 REQUIRE "ivpsystem.a4l";
2     REQUIRE "atoms.a4l";
3     REQUIRE "johnpye/thermo_types.a4c";
4    
5     (*
6 johnpye 1275 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 johnpye 1179 *)
10     REQUIRE "steam/satsteamstream.a4c";
11    
12 jpye 1331 IMPORT "johnpye/extpy/extpy";
13     IMPORT "johnpye/roots";
14    
15 johnpye 1179 MODEL dsgsat3;
16     n IS_A integer_constant;
17 jpye 1383 n :== 30;(* with L = 10m: 5,6,7,8,9,10,11 *)
18 johnpye 1305 (* with L = 5m: 2,3,4,5,7,9,11,12,13,1415,16 *)
19 johnpye 1179
20 johnpye 1305 dz IS_A real_constant;
21     L IS_A real_constant;
22 jpye 1383 L :== 11 {m};
23 johnpye 1305 dz :== L / (n-1);
24    
25     nodes,butfirst1,upwind4,central IS_A set OF integer_constant;
26     nodes :== [1..n];
27     butfirst1 :== nodes - [1];
28     upwind4 :== nodes - [1,2,n];
29     central :== nodes - [1,n];
30    
31 johnpye 1179 (* temporal derivatives *)
32 johnpye 1305 drho_dt[butfirst1] IS_A density_rate;
33 jpye 1383 (* dmdot_dt[butfirst1] IS_A mass_rate_rate; *)
34 johnpye 1305 du_dt[butfirst1] IS_A specific_energy_rate;
35     dTw_dt[butfirst1] IS_A temperature_rate;
36 johnpye 1179
37     (* wall properties *)
38     rho_w IS_A mass_density;
39     D, D_2 IS_A distance;
40     c_w IS_A specific_heat_capacity;
41     A, A_w IS_A area;
42     h_int IS_A heat_transfer_coefficient; (* internal *)
43     h_ext IS_A heat_transfer_coefficient; (* external *)
44     z_A: A = 1{PI}*D^2/4;
45     z_Aw: A_w = 1{PI}*(D_2^2 - D^2)/4;
46    
47     (* fluid properties *)
48 johnpye 1305 node[nodes] IS_A satsteamstream;
49 johnpye 1179
50     (* flow properties *)
51 johnpye 1305 vel[nodes] IS_A speed;
52     T_w[butfirst1] IS_A temperature;
53     T[nodes] IS_A temperature;
54 johnpye 1179
55 johnpye 1305 (* constant, for the moment: *)
56 johnpye 1179 f IS_A positive_factor;
57     (* mu_f IS_A viscosity; *)
58     T_amb IS_A temperature;
59    
60     (* system dynamics *)
61 johnpye 1305 qdot_t[butfirst1], qdot_l[butfirst1] IS_A power_per_length;
62 johnpye 1179 qdot_s IS_A power_per_length;
63    
64     (* some aliases just for easier review of the state of the model *)
65 johnpye 1305 x[nodes] IS_A fraction;
66     mdot[nodes] IS_A mass_rate;
67     p[nodes] IS_A pressure;
68     rho[nodes] IS_A mass_density;
69     u[nodes] IS_A specific_energy;
70     v[nodes] IS_A specific_volume;
71 jpye 1383
72     (* Qdot_t IS_A energy_rate;
73     Qdot_t = SUM[qdot_t[i] | i IN butfirst1]; *)
74    
75 johnpye 1305 FOR i IN nodes CREATE
76 jpye 1383 z_vel[i]: vel[i] = v[i]*mdot[i]/A;
77     END FOR;
78    
79     FOR i IN nodes CREATE
80 johnpye 1305 x[i], node[i].x ARE_THE_SAME;
81 johnpye 1179 mdot[i], node[i].mdot ARE_THE_SAME;
82 johnpye 1305 p[i], node[i].p ARE_THE_SAME;
83     T[i], node[i].T ARE_THE_SAME;
84     rho[i], node[i].rho ARE_THE_SAME;
85     u[i], node[i].u ARE_THE_SAME;
86     v[i], node[i].v ARE_THE_SAME;
87 johnpye 1179 END FOR;
88    
89 johnpye 1312 en_upwind4,en_central,mom_upwind4,mom_central,mass_upwind4,mass_central IS_A set OF integer_constant;
90    
91 johnpye 1305 (* mass conservation *)
92 johnpye 1312 mass_upwind4 :== upwind4;
93 jpye 1328 mass_central :== central - mass_upwind4;
94 johnpye 1312 FOR i IN mass_upwind4 CREATE (* 4-pt upwind biased *)
95 johnpye 1305 z_massbal1[i]: A * drho_dt[i] * dz =
96     - (mdot[i+1] + 6.*mdot[i] - 3.*mdot[i-1] - 2.*mdot[i-2]) / 6.;
97     END FOR;
98 johnpye 1312 FOR i IN mass_central CREATE
99 johnpye 1305 z_massbal2[i]: A * drho_dt[i] * dz =
100     - (mdot[i+1] - mdot[i-1]) / 2.;
101     END FOR;
102 johnpye 1312 FOR i IN butfirst1 - mass_upwind4 - mass_central CREATE
103 johnpye 1305 z_massbal[i]: A * drho_dt[i] * dz = - (mdot[i] - mdot[i-1]);
104     END FOR;
105 johnpye 1277
106 johnpye 1305 (* energy conservation *)
107 johnpye 1312 en_upwind4 :== [];
108 jpye 1328 en_central :== central - en_upwind4;
109 johnpye 1312 FOR i IN en_upwind4 CREATE
110 johnpye 1305 z_enbal2[i]: dz * (qdot_t[i] - rho[i] * A * du_dt[i]) =
111     + mdot[i] * (node[i+1].u + 6.*u[i] - 3.*u[i-1] - 2.*u[i-1]) / 6.
112     + (p[i+1]*node[i+1].v*mdot[i+1] - p[i-1]*v[i-1]*mdot[i-1]) / 2.;
113 johnpye 1179 END FOR;
114 johnpye 1312 FOR i IN en_central CREATE
115 johnpye 1305 z_enbal1[i]: dz * (qdot_t[i] - rho[i] * A * du_dt[i]) =
116     + mdot[i] * (u[i] - u[i-1]) (* NOTE: not central *)
117     + (p[i+1]*v[i+1]*mdot[i+1] - p[i-1]*v[i-1]*mdot[i-1]) / 2.;
118     END FOR;
119 johnpye 1312 FOR i IN butfirst1 - en_upwind4 - en_central CREATE
120 johnpye 1305 z_enbal[i]: dz * (qdot_t[i] - rho[i] * A * du_dt[i]) =
121     + mdot[i] * (u[i] - u[i-1])
122     + (p[i]*v[i]*mdot[i] - p[i-1]*v[i-1]*mdot[i-1]);
123     END FOR;
124 johnpye 1277
125 jpye 1383 (* stationary momentum *)
126     FOR i IN butfirst1 CREATE
127     z_mombal[i]: p[i] = p[i-1] - dz * f/D/2 * rho[i] * vel[i]^2;
128     END FOR;
129 johnpye 1312 mom_upwind4 :== [];
130 jpye 1383 mom_central :== [];
131     (*
132     * momentum conservation *
133 johnpye 1312 FOR i IN mom_upwind4 CREATE
134 johnpye 1305 z_mombal2[i]: - dz/A * dmdot_dt[i]
135 jpye 1383 = (p[i]-p[i-1]) * backdiff for pressure *
136 johnpye 1305 + dz * f/D/2 * rho[i] * vel[i]^2
137     + (rho[i+1]*vel[i+1]^2 + 6.*rho[i]*vel[i]^2 - 3.*rho[i-1]*vel[i-1]^2 - 2.*rho[i-2]*vel[i-2]^2) / 6.;
138 johnpye 1278 END FOR;
139 johnpye 1312 FOR i IN mom_central CREATE
140 johnpye 1305 z_mombal1[i]: - dz/A * dmdot_dt[i]
141     = (p[i+1]-p[i-1]) / 2.
142     + dz * f/D/2 * rho[i] * vel[i]^2
143     + (rho[i+1]*vel[i+1]^2 - rho[i-1]*vel[i-1]^2) / 2.;
144     END FOR;
145 johnpye 1312 FOR i IN butfirst1 - mom_upwind4 - mom_central CREATE
146 johnpye 1305 z_mombal[i]: - dz/A * dmdot_dt[i]
147     = (p[i]-p[i-1])
148     + dz * f/D/2 * rho[i] * vel[i]^2
149     + (rho[i]*vel[i]^2 - rho[i-1]*vel[i-1]^2);
150     END FOR;
151 jpye 1383 *)
152 johnpye 1277
153 johnpye 1305 (* internal/external convection, and thermal mass of wall -- no spatial derivs here *)
154     FOR i IN butfirst1 CREATE
155 johnpye 1278 z_wall[i]: rho_w*A_w*c_w*dTw_dt[i] = qdot_s - qdot_l[i] - qdot_t[i];
156     z_loss[i]: qdot_l[i] = h_ext*(1{PI}*D_2)*(T_w[i] - T_amb);
157 johnpye 1305 z_trans[i]: qdot_t[i] = h_int*(1{PI}*D) *(T_w[i] - T[i]);
158 johnpye 1277 END FOR;
159    
160 johnpye 1179 t IS_A time;
161 jpye 1383
162 johnpye 1179 METHODS
163 jpye 1383
164 johnpye 1278 METHOD bound_self;
165 johnpye 1305 vel[nodes].upper_bound := 100 {m/s};
166     qdot_l[butfirst1].lower_bound := 0 {W/m};
167     FOR i IN nodes DO
168 johnpye 1278 RUN node[i].bound_self;
169     END FOR;
170     END bound_self;
171 johnpye 1305 METHOD default;
172     (* these are initial guesses only; fixed parameters are overwritten by 'values' below *)
173 johnpye 1278 t := 0 {s};
174 johnpye 1305 FOR i IN nodes DO
175 johnpye 1312 T[i] := 298 {K};
176     vel[i] := 1 {m/s};
177 johnpye 1278 RUN node[i].default_self;
178     END FOR;
179 johnpye 1305 FOR i IN butfirst1 DO
180 johnpye 1278 drho_dt[i] := 0 {kg/m^3/s};
181 jpye 1383 (* dmdot_dt[i] := 0 {kg/s/s}; *)
182 johnpye 1278 du_dt[i] := 0 {kJ/kg/s};
183     dTw_dt[i] := 0 {K/s};
184     qdot_t[i] := 0 {W/m};
185     qdot_l[i] := 0 {W/m};
186 johnpye 1312 x[i] := x[1];
187 johnpye 1278 END FOR;
188 johnpye 1305 END default;
189 johnpye 1278 METHOD specify;
190     (* change to a proper steady-state problem, with fluid properties FREEd *)
191 johnpye 1305 FOR i IN nodes DO
192 johnpye 1278 RUN node[i].specify;
193     FIX dTw_dt[i]; FREE T_w[i];
194     END FOR;
195 johnpye 1305 FIX p[1];
196     FREE T[1];
197 johnpye 1278 FIX qdot_s;
198 johnpye 1305 FIX D, D_2;
199 johnpye 1278 FIX h_int, c_w, rho_w, h_ext;
200     FIX f;
201     (* FIX mu_f; *)
202     FIX T_amb;
203     (* fix derivatives to zero *)
204 johnpye 1311 FOR i IN butfirst1 DO
205 jpye 1383 FREE x[i]; FIX p[i]; FREE node[i].mdot;
206 johnpye 1305 FIX drho_dt[i]; FREE p[i];
207     FIX du_dt[i]; FREE T[i];
208 jpye 1383 (* FREE mdot[i]; FIX dmdot_dt[i]; *)
209 johnpye 1278 END FOR;
210     END specify;
211     METHOD values;
212 johnpye 1313 D := 0.06 {m};
213     D_2 := 0.07 {m};
214     T_amb := 298 {K};
215 jpye 1383 h_int := 1000 {W/m^2/K};
216     h_ext := 5 {W/m^2/K};
217 johnpye 1307 f := 0.03;
218 johnpye 1305 mdot[1] := 0.26 {kg/s};
219     p[1] := 10 {bar};
220 johnpye 1307 x[1] := 0.23;
221 johnpye 1313 rho_w := 1000 {kg/m^3};
222 jpye 1383 qdot_s := 100 {W/m^2} * D_2 * 10 * 11{m} / L;
223 johnpye 1311 FOR i IN butfirst1 DO
224 johnpye 1313 T_w[i] := 298 {K};
225 jpye 1383 (* dmdot_dt[i] := 0.0 {kg/s/s};*)
226 johnpye 1278 du_dt[i] := 0 {kJ/kg/s};
227 johnpye 1305 v[i] := 0.2 {L/kg};
228     rho[i] := 6 {kg/L};
229 johnpye 1278 node[i].dp_dT := +0.5 {kPa/K};
230 johnpye 1305 p[i] := 5 {bar};
231 johnpye 1278 END FOR;
232     END values;
233     METHOD on_load;
234     RUN configure_steady;
235     RUN ode_init;
236     END on_load;
237     (*---------------- a physically sensible steady-state configuration-----------*)
238     METHOD configure_steady;
239 jpye 1383 EXTERNAL defaultself_visit_childatoms(SELF);
240     EXTERNAL defaultself_visit_submodels(SELF);
241 johnpye 1278 RUN ClearAll;
242     RUN specify;
243     RUN bound_steady;
244     RUN values;
245     END configure_steady;
246     METHOD bound_steady;
247     RUN bound_self;
248 johnpye 1311 T_w[butfirst1].upper_bound := 1000 {K};
249 johnpye 1278 END bound_steady;
250     (*------------------------- the dynamic problem ------------------------------*)
251     METHOD configure_dynamic;
252 johnpye 1311 FOR i IN butfirst1 DO
253 johnpye 1305 FREE drho_dt[i]; FIX rho[i];
254 jpye 1383 (* FREE dmdot_dt[i]; FIX mdot[i]; *)
255 johnpye 1305 FREE du_dt[i]; FIX u[i];
256 johnpye 1278 FREE dTw_dt[i]; FIX T_w[i];
257 johnpye 1305 FREE x[i];
258     FREE T[i];
259 johnpye 1278 END FOR;
260     t := 0 {s};
261     END configure_dynamic;
262     METHOD free_states;
263 johnpye 1311 FOR i IN butfirst1 DO
264 johnpye 1305 FREE rho[i];
265 jpye 1383 (* FREE mdot[i]; *)
266 johnpye 1305 FREE u[i];
267 johnpye 1278 FREE T_w[i];
268     END FOR;
269     END free_states;
270     METHOD ode_init;
271     (* add the necessary meta data to allow solving with the integrator *)
272     t.ode_type := -1;
273 johnpye 1179
274 johnpye 1311 FOR i IN butfirst1 DO
275 johnpye 1305 drho_dt[i].ode_id := 4*i; rho[i].ode_id := 4*i;
276     drho_dt[i].ode_type := 2; rho[i].ode_type := 1;
277 johnpye 1179
278 jpye 1383 (* dmdot_dt[i].ode_id := 4*i+1; mdot[i].ode_id := 4*i+1;
279     dmdot_dt[i].ode_type := 2; mdot[i].ode_type := 1;*)
280 johnpye 1278
281 johnpye 1305 du_dt[i].ode_id := 4*i+2; u[i].ode_id := 4*i+2;
282     du_dt[i].ode_type := 2; u[i].ode_type := 1;
283 johnpye 1179
284 johnpye 1278 dTw_dt[i].ode_id := 4*i+3; T_w[i].ode_id := 4*i+3;
285     dTw_dt[i].ode_type := 2; T_w[i].ode_type := 1;
286     END FOR;
287 johnpye 1179
288 johnpye 1305 FOR i IN nodes DO
289 johnpye 1278 (* p[i].obs_id := 1 + 10*i; *)
290     (* x[i].obs_id := 2 + 10*i; *)
291 johnpye 1305 (* mdot[i].obs_id := 4 + 10*i; *)
292 johnpye 1279 (* node[i].h.obs_id := 3 + 10*i; *)
293 johnpye 1278 END FOR;
294 johnpye 1311 FOR i IN butfirst1 DO
295 johnpye 1278 (* qdot_t[i].obs_id := 3 + 10*i; *)
296     (* T_w[i].obs_id := 5 + 10*i; *)
297     (* T[i].obs_id := 6 + 10*i;*)
298     END FOR;
299 johnpye 1279
300 johnpye 1305 mdot[n].obs_id := 1;
301     x[n].obs_id := 1;
302 johnpye 1279 p[n].obs_id := 1;
303     vel[n].obs_id := 1;
304 jpye 1383 T_w[n].obs_id := 2;
305 johnpye 1279
306 johnpye 1278 END ode_init;
307     METHOD fix_outlet_quality;
308     FIX x[n];
309 johnpye 1305 FREE mdot[1];
310 johnpye 1278 END fix_outlet_quality;
311 johnpye 1179
312 jpye 1331 METHOD roots;
313     EXTERNAL roots(SELF);
314     END roots;
315    
316 johnpye 1179 END dsgsat3;
317     ADD NOTES IN dsgsat2;
318     'QRSlv' iterationlimit {50}
319     END NOTES;

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