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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 923 - (show annotations) (download) (as text)
Thu Nov 9 10:54:14 2006 UTC (15 years, 7 months ago) by johnpye
File MIME type: text/x-ascend
File size: 4744 byte(s)
Moving this file from my (separate) dsg-transient project into ASCEND.
1 REQUIRE "ivpsystem.a4l";
2 REQUIRE "atoms.a4l";
3 REQUIRE "johnpye/thermo_types.a4c";
4
5 IMPORT "johnpye/extpy/extpy";
6 IMPORT "johnpye/solve";
7 IMPORT "johnpye/solvernotes";
8
9 (*
10 This model uses some ASCEND models from the freesteam library. See
11 http://freesteam.sf.net/ for more information. This model doesn't actually
12 require compiled binaries of freesteam, so you can just download the .a4c
13 files if you wish.
14 *)
15 REQUIRE "satsteamstream.a4c";
16
17
18 MODEL dsgsat2;
19 n IS_A integer_constant;
20 n :== 3;
21
22 (* temporal derivatives *)
23 drho_dt[2..n] IS_A density_rate;
24 dmdot_dt[2..n] IS_A mass_rate_rate;
25 drhou_dt[2..n] IS_A power_per_volume;
26 dTw_dt[2..n] IS_A temperature_rate;
27
28 (* wall properties *)
29 rho_w IS_A mass_density;
30 D, D_2 IS_A distance;
31 c_w IS_A specific_heat_capacity;
32 A, A_w IS_A area;
33 h_int IS_A heat_transfer_coefficient; (* internal *)
34 h_ext IS_A heat_transfer_coefficient; (* external *)
35 A = 1{PI}*D^2/4;
36 A_w = 1{PI}*(D_2^2 - D^2)/4;
37 dz IS_A distance;
38 L IS_A distance;
39 dz = L / n;
40
41 (* fluid properties *)
42 node[1..n] IS_A satsteamstream;
43
44 (* flow properties *)
45 vel[1..n] IS_A speed;
46 T_w[1..n] IS_A temperature;
47
48 (* constants, for the moment: *)
49 f IS_A positive_factor;
50 mu_f IS_A viscosity;
51 T_amb IS_A temperature;
52
53 (* system dynamics *)
54 qdot_t[2..n], qdot_l[2..n] IS_A power_per_length;
55 qdot_s IS_A power_per_length;
56 rhou[1..n] IS_A energy_per_volume;
57
58 FOR i IN [1..n] CREATE
59 vel[i] = node[i].v*node[i].mdot/A;
60 rhou[i] = node[i].rho * node[i].u;
61 END FOR;
62
63 (* some aliases just for easier review of the state of the model *)
64 x[1..n] IS_A fraction;
65 mdot[1..n] IS_A mass_rate;
66 p[1..n] IS_A pressure;
67 FOR i IN [1..n] CREATE
68 x[i], node[i].x ARE_THE_SAME;
69 mdot[i], node[i].mdot ARE_THE_SAME;
70 p[i], node[i].p ARE_THE_SAME;
71 END FOR;
72
73 (* differential equations *)
74 FOR i IN [2..n] CREATE
75 A * drho_dt[i] = - (node[i].mdot - node[i-1].mdot)/dz;
76 1/A*dmdot_dt[i] = -(node[i].p-node[i-1].p)/dz - f/D/2*node[i].rho*node[i].v^2* (node[i].rho*vel[i]^2 - node[i-1].rho*vel[i-1]^2)/dz;
77 A * drhou_dt[i] = qdot_t[i] - (node[i].Hdot - node[i-1].Hdot)/dz;
78 rho_w*A_w*c_w*dTw_dt[i] = qdot_s - qdot_l[i] - qdot_t[i];
79 qdot_l[i] = h_ext*(1{PI}*D_2)*(T_w[i] - T_amb);
80 qdot_t[i] = h_int*(1{PI}*D) *(T_w[i] - node[i].T);
81 END FOR;
82
83 t IS_A time;
84 METHODS
85 METHOD specify;
86 RUN node[1].specify;
87 FIX qdot_s;
88 FIX D, D_2, L;
89 FIX h_int, c_w, rho_w, h_ext;
90 FIX f, mu_f;
91 FIX T_amb;
92 (* fix derivatives to zero *)
93 FOR i IN [2..n] DO
94 FIX drho_dt[i]; FREE node[i].rho;
95 FIX dmdot_dt[i]; FREE node[i].mdot;
96 FIX drhou_dt[i]; FREE rhou[i];
97 FIX dTw_dt[i]; FREE T_w[i];
98 END FOR;
99 (* FIX node[3].rho; *)
100
101 END specify;
102 METHOD values;
103 node[1].T := 400 {K};
104 node[1].x := 0.1;
105 qdot_s := 0 {W/m};
106 D := 60 {mm}; D_2 := 70 {mm};
107 L := 100 {m};
108 A_w := 1{PI}*D_2^2;
109 h_int := 10 {W/m^2/K}; c_w := 0.47 {J/g/K}; rho_w := 7.8 {g/cm^3}; h_ext := 10 {W/m^2/K};
110 f := 0.005; mu_f := 4.5e-5 {Pa*s};
111 T_amb := 300 {K};
112 FOR i IN [2..n] DO
113 drho_dt[i] := 0 {kg/m^3/s};
114 dmdot_dt[i] := 0 {kg/s/s};
115 drhou_dt[i] := 0 {kJ/m^3/s};
116 dTw_dt[i] := 0 {K/s};
117 END FOR;
118 END values;
119 METHOD on_load;
120 RUN specify;
121 FOR i IN [1..n] DO
122 RUN node[i].bound_self;
123 END FOR;
124 RUN values;
125 RUN solvernotes;
126 RUN solve; (* wish I knew why it needs three goes *)
127 RUN ode_init;
128 qdot_s := 10 {W/m};
129 RUN solve; (* after fixing the states and freeing the derivatives *)
130
131 END on_load;
132 METHOD ode_init;
133 (* get the model into the required state for solving as ODE *)
134 t.ode_type := -1;
135 t := 0 {s};
136 FOR i IN [2..n] DO
137 FREE drho_dt[i]; FIX node[i].rho;
138 FREE dmdot_dt[i]; FIX node[i].mdot;
139 FREE drhou_dt[i]; FIX rhou[i];
140 FREE dTw_dt[i]; FIX T_w[i];
141 END FOR;
142
143 FOR i IN [2..n] DO
144 drho_dt[i].ode_id := 4*i; node[i].rho.ode_id := 4*i;
145 drho_dt[i].ode_type := 2; node[i].rho.ode_type := 1;
146
147 dmdot_dt[i].ode_id := 4*i+1; node[i].mdot.ode_id := 4*i+1;
148 dmdot_dt[i].ode_type := 2; node[i].mdot.ode_type := 1;
149
150 drhou_dt[i].ode_id := 4*i+2; rhou[i].ode_id := 4*i+2;
151 drhou_dt[i].ode_type := 2; rhou[i].ode_type := 1;
152
153 dTw_dt[i].ode_id := 4*i+3; T_w[i].ode_id := 4*i+3;
154 dTw_dt[i].ode_type := 2; T_w[i].ode_type := 1;
155
156 (*
157 p[i].obs_id := 4*i;
158 x[i].obs_id := 4*i+1;
159 qdot_t[i].obs_id := 4*i+2;
160 T_w[i].obs_id := 4*i+3;
161 *)
162 END FOR;
163
164 FOR i IN [1,n] DO
165 p[i].obs_id := 4*i;
166 x[i].obs_id := 4*i+1;
167 END FOR;
168 FOR i IN [2,n] DO
169 (* qdot_t[i].obs_id := 4*i+2; *)
170 T_w[i].obs_id := 4*i+3;
171 END FOR;
172 END ode_init;
173
174 METHOD fix_outlet_quality;
175 FIX x[n];
176 FREE node[1].mdot;
177 END fix_outlet_quality;
178
179 METHOD solvernotes;
180 EXTERNAL solvernotes(SELF);
181 END setup_solver;
182
183 METHOD solve;
184 EXTERNAL solve(SELF);
185 END solve;
186
187 END dsgsat2;
188 ADD NOTES IN dsgsat2;
189 'QRSlv' iterationlimit {50}
190 END NOTES;

Properties

Name Value
svn:executable *

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