/[ascend]/branches/fprops-incomp/models/johnpye/fprops/towerarray.py
ViewVC logotype

Annotation of /branches/fprops-incomp/models/johnpye/fprops/towerarray.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3470 - (hide annotations) (download) (as text)
Thu Feb 13 01:56:33 2020 UTC (7 months, 1 week ago) by jpye
File MIME type: text/x-python
File size: 12418 byte(s)
add reporting of dp

1 jpye 3463 # coding=utf-8
2     import matplotlib
3     matplotlib.use('GtkAgg')
4     import matplotlib.pyplot as pl
5     import numpy as np
6     import math
7    
8     class Label:
9     def __init__(self,ix,iy,i):
10     self.ix = ix
11     self.iy = iy
12     self.i = i
13     def __repr__(self):
14     return "{%d,%d}"%(self.ix,self,iy)
15    
16     class Node:
17     def __init__(self,x,y,lab):
18     self.x = x
19     self.y = y
20     self.lab = lab
21     self.links_to = []
22     self.links_from = []
23     def distance_from(self,cnode):
24     return math.sqrt((cnode.x-self.x)**2 + (cnode.y-self.y)**2)
25    
26     def __repr__(self):
27     return "<Node %s @%.0f,%.0f>" % (repr(self.lab),self.x, self.y)
28    
29     class Link:
30     def __init__(self,lab1,lab2,id):
31     self.lab1 = lab1
32     self.lab2 = lab2
33     self.id = id
34     def __repr__(self):
35     return "<Link %d: %s--%s>" % (self.id,repr(self.lab1),repr(self.lab2))
36 jpye 3467 def plot(self,arr,xsep,ysep):
37 jpye 3463 #print [arr[self.lab1].x,arr[self.lab2].x],[arr[self.lab1].y,arr[self.lab1].y]
38     x = arr[self.lab1].x
39     y = arr[self.lab1].y
40     dx = arr[self.lab2].x - arr[self.lab1].x
41     dy = arr[self.lab2].y - arr[self.lab1].y
42     l = math.sqrt(dx**2+dy**2)
43     hw = 0.2 * xsep
44     hl = 1 * hw
45     l1 = l + hl
46     dx *= l/l1
47     dy *= l/l1
48     pl.arrow(x,y,dx,dy,width=xsep*0.03,head_width=hw,head_length=hl,color='#ff0000')
49     pl.annotate("id=%d"%(self.id,),(x+dx/2,y+dy/2))
50     def node_from(self,arr):
51     return arr[self.lab1]
52     def node_to(self,arr):
53     return arr[self.lab2]
54     def reverse(self):
55     l = self.lab1
56     self.lab1 = self.lab2
57     self.lab2 = l
58     def length(self,arr):
59     dx = arr[self.lab2].x - arr[self.lab1].x
60     dy = arr[self.lab2].y - arr[self.lab1].y
61     l = math.sqrt(dx**2+dy**2)
62     return l
63    
64 jpye 3467 def create_field_model(nx,ny,Qdot_onecoll_MW=10.,mdot_onecoll=None
65     ,ix_pb=None,iy_pb=None,suffix=None,f=None):
66 jpye 3466 xsep = 290.
67     ysep = 240.
68 jpye 3463
69 jpye 3467 if f is None:
70     import sys
71     f = sys.stdout
72     if ix_pb is None:
73     ix_pb = int(nx)/2
74     if iy_pb is None:
75     iy_pb = int(ny)/2
76     if mdot_onecoll is None:
77     mdot_onecoll = 40 * (Qdot_onecoll_MW / 10)
78     if suffix is None:
79     suffix = "_%dx%d_%fMW"%(nx,ny,Qdot_onecoll_MW)
80    
81     ix = np.arange(0,nx)
82     x = ix*xsep;
83     iy = np.arange(0,ny)
84     y = iy*ysep;
85     IX,IY = np.meshgrid(ix,iy)
86     X,Y = np.meshgrid(x,y)
87     lab = zip(IX.ravel(),IY.ravel())
88     loc = zip(X.ravel(),Y.ravel())
89 jpye 3463
90 jpye 3467 arr = {}
91     for i,l in enumerate(lab):
92     print "node label",l
93     arr[l] = Node(x=loc[i][0],y=loc[i][1],lab=l)
94 jpye 3463
95 jpye 3467 print "arr =",arr
96 jpye 3463
97 jpye 3467 centre = arr[(ix_pb,iy_pb)]
98 jpye 3463
99 jpye 3467 # create the link network...
100     links = []
101     id = 0
102     for ix1 in ix:
103     print "col",ix1
104     for iy1 in np.arange(ny-1):
105     l = Link((ix1,iy1),(ix1,iy1+1),id); id += 1
106     print "row",iy1,"to",iy1+1,":",l
107     links.append(l)
108     if ix1 != 0:
109     print ix1
110     links.append(Link((ix1-1,iy_pb),(ix1,iy_pb),id)); id += 1
111 jpye 3463
112 jpye 3467 print links
113 jpye 3463
114 jpye 3467 f.write("""
115     MODEL layout{suffix};
116 jpye 3463 n_DI IS_A integer_constant;
117     n_DI :== {n};
118     n_PH IS_A integer_constant;
119     n_PH :== {nl};
120     n_PC ALIASES n_PH;
121     DI[0..n_DI-1] IS_A boiler_simple;
122     TH[0..n_DI-1] IS_A throttle;
123     PH[0..n_PH-1] IS_A pipe_heat_loss_insul;
124     PC[0..n_PC-1] IS_A pipe_heat_loss_insul;
125     JO[0..n_DI-1] IS_A merge_generic;
126     TE[0..n_DI-1] IS_A tee_generic;
127    
128     PC[0..n_PH-1].config, PH[0..n_PH-1].config ARE_THE_SAME;
129    
130     corr_intconv ALIASES PC[0].config.corr_intconv;
131     k_pipe ALIASES PC[0].config.k_pipe;
132     k_insul ALIASES PC[0].config.k_insul;
133     h_ext ALIASES PC[0].config.h_ext;
134     emiss ALIASES PC[0].config.emiss;
135     eps ALIASES PC[0].config.eps;
136    
137 jpye 3467 corr_intconv :== 'none';\n\n""".format(
138     n=len(arr)
139     ,nl=len(links)
140     ,suffix=suffix
141     ))
142 jpye 3463
143 jpye 3467 # reverse mappings
144 jpye 3463
145 jpye 3467 rarr = {}
146     rlinks = {}
147     for i,lab in enumerate(arr.keys()):
148     rarr[lab] = i
149     for i,l in enumerate(links):
150     rlinks[l] = i
151 jpye 3463
152 jpye 3467 # create the ASCEND code...
153 jpye 3463
154 jpye 3467 if 0:
155     f.write("\tNOTES\n");
156     # output the positions of each collector (in comments)
157     for i,lab in enumerate(arr.keys()):
158     print i
159     f.write("\t\t'location' DI[{i}].inlet {{{x},{x}}};\n".format(
160     i=i
161     ,x=arr[lab].x
162     ,y=arr[lab].y
163     ))
164     f.write("\tEND NOTES;\n\n");
165 jpye 3463
166 jpye 3467 for i,l in enumerate(links):
167     A = l.node_from(arr)
168     B = l.node_to(arr)
169     print "link",l,": A_dist = %f, B_dist = %f" % (A.distance_from(centre), B.distance_from(centre))
170     if A.distance_from(centre) > B.distance_from(centre):
171     l.reverse()
172     B.links_from.append(l)
173     A.links_to.append(l)
174     print "link reverse",l
175     else:
176     A.links_from.append(l)
177     B.links_to.append(l)
178 jpye 3463
179 jpye 3467 for i,n in enumerate(arr):
180     assert len(arr[n].links_to) <= 1
181     # print "NOTE: Too many links_to node %s: %s"%(n,arr[n].links_to)
182 jpye 3463
183 jpye 3467 f.write("\t(* links at node {i} {lab} *)\n".format(i=i,lab=n))
184 jpye 3463
185 jpye 3467 if len(arr[n].links_from) == 0:
186     f.write("\tTE[{i}] IS_REFINED_TO tee_trivial;\n\tJO[{i}] IS_REFINED_TO merge_trivial;\n".format(i=i))
187     else:
188     f.write("\tTE[{i}] IS_REFINED_TO tee_n({nout});\n\tJO[{i}] IS_REFINED_TO merge_n({nout});\n".format(
189     i=i
190     ,lab = n
191     ,nin = len(arr[n].links_to)
192     ,nout = len(arr[n].links_from)+1
193     ))
194 jpye 3463
195 jpye 3467 for j,lt in enumerate(arr[n].links_to):
196     f.write("\tTE[{i}].inlet, PC[{li}].outlet ARE_THE_SAME;\n".format(
197     i=i
198     ,li=rlinks[lt]
199     ))
200 jpye 3463
201 jpye 3467 f.write("\tTE[{i}].outlet[0], DI[{i}].inlet ARE_THE_SAME;\n".format(i=i))
202     for j,lt in enumerate(arr[n].links_from):
203     f.write("\tTE[{i}].outlet[{j}], PC[{li}].inlet ARE_THE_SAME;\n".format(
204     i=i
205     ,j=j+1
206     ,li=rlinks[lt]
207     ))
208 jpye 3463
209 jpye 3467 f.write("\tDI[{i}].outlet, TH[{i}].inlet ARE_THE_SAME;\n".format(i=i));
210     f.write("\tTH[{i}].outlet, JO[{i}].inlet[0] ARE_THE_SAME;\n".format(i=i));
211 jpye 3463
212 jpye 3467 if len(arr[n].links_from) == 0:
213     f.write("\tTH[{i}] IS_REFINED_TO throttle_trivial;\n".format(i=i));
214 jpye 3463
215 jpye 3467 for j,lt in enumerate(arr[n].links_from):
216     f.write("\tJO[{i}].inlet[{j}], PH[{li}].outlet ARE_THE_SAME;\n".format(
217     i=i
218     ,j=j+1
219     ,li=rlinks[lt]
220     ))
221 jpye 3463
222 jpye 3467 for j,lt in enumerate(arr[n].links_to):
223     f.write("\tJO[{i}].outlet, PH[{li}].inlet ARE_THE_SAME;\n".format(
224     i=i
225     ,li=rlinks[lt]
226     ))
227 jpye 3463
228    
229    
230 jpye 3467 import sys
231     print "centre.lab",centre.lab
232     print "rarr[centre.lab]",rarr[centre.lab]
233 jpye 3463
234    
235 jpye 3467 f.write("""
236 jpye 3463 inlet, outlet IS_A stream_node;
237     inlet, TE[{i}].inlet ARE_THE_SAME;
238     outlet, JO[{i}].outlet ARE_THE_SAME;
239    
240     cd ALIASES inlet.cd;
241     cd.component :== 'sodium';
242     cd.type :== 'incomp';
243    
244     METHODS
245 jpye 3467 """.format(
246     i=rarr[centre.lab]
247 jpye 3463 ))
248    
249 jpye 3467 f.write("METHOD set_pipe_lengths;\n");
250     for i,l in enumerate(links):
251     f.write("""\tFIX PH[{i}].L := {L} {{m}};
252     FIX PC[{i}].L := {L} {{m}};\n""".format(
253     i=i
254     ,L = l.length(arr)
255     ))
256     f.write("END set_pipe_lengths;\n");
257 jpye 3463
258 jpye 3467
259     f.write("""
260 jpye 3463 METHOD default_self;
261     FOR i IN [0..n_DI-1] DO
262     RUN TE[i].default_self;
263     RUN DI[i].default_self;
264     RUN TH[i].default_self;
265     RUN JO[i].default_self;
266     END FOR;
267     FOR i IN [0..n_PH-1] DO
268     RUN PH[i].default_self;
269     RUN PC[i].default_self;
270     END FOR;
271     END default_self;
272 jpye 3467 END layout{suffix};
273     """.format(
274     i=rarr[centre.lab]
275     ,suffix=suffix
276     ))
277 jpye 3463
278 jpye 3467 # note, we don't use str.format() here because lots of units with curly {}:
279     f.write("""
280 jpye 3463 UNITS
281     (* currency conversions as of 10 Feb 2020 *)
282     AUD = {USD/1.49};
283     EUR = {1.10*USD};
284     END UNITS;
285 jpye 3467 """);
286 jpye 3463
287 jpye 3467 f.write("""
288     MODEL towerarray{suffix} REFINES layout{suffix};
289     """.format(suffix=suffix))
290 jpye 3463
291 jpye 3467 f.write("""
292 jpye 3463 PC[0..n_PH-1].Vel_out, PH[0..n_PH-1].Vel_out ARE_THE_SAME;
293     Vel ALIASES PC[0].Vel_out;
294    
295     PC[0..n_PH-1].t_insul, PH[0..n_PH-1].t_insul ARE_THE_SAME;
296 jpye 3465
297     FOR i IN [0..n_PH-1] CREATE
298     (* curve fit to Sch10S pipe sizes *)
299     PC[i].t_pipe = 1{mm} * (PC[i].D/1{mm})^0.3 / 1.222;
300     PH[i].t_pipe = 1{mm} * (PH[i].D/1{mm})^0.3 / 1.222;
301     END FOR;
302     (*PC[0..n_PH-1].t_pipe, PH[0..n_PH-1].t_pipe ARE_THE_SAME;
303     t_pipe ALIASES PC[0].t_pipe;*)
304 jpye 3463 t_insul ALIASES PC[0].t_insul;
305    
306    
307     (* total materials *)
308 jpye 3469 m_pipe, m_insul, m_sodium IS_A mass;
309     V_insul IS_A volume;
310 jpye 3463 rho_pipe, rho_insul IS_A mass_density;
311     L_pipe_tot IS_A distance;
312     m_pipe = rho_pipe * SUM[PC[i].solid_pipe.V + PH[i].solid_pipe.V | i IN [0..n_PH-1]];
313 jpye 3469 V_insul = SUM[PC[i].solid_insul.V + PH[i].solid_insul.V | i IN [0..n_PH-1]];
314 jpye 3463 L_pipe_tot = SUM[PC[i].L + PH[i].L | i IN [0..n_PH-1]];
315 jpye 3469 m_sodium = SUM[PC[i].solid.V / PC[i].inlet.v + PH[i].solid.V / PH[i].inlet.v | i IN [0..n_PH-1]];
316     m_insul = rho_insul * V_insul;
317 jpye 3463
318     (* costs *)
319 jpye 3469 C_pipe, C_insul, C_inst_pipe, C_inst_insul, C_supp, C_sodium, C_tot IS_A monetary_unit;
320     c_pipe IS_A cost_per_mass;
321     c_insul IS_A cost_per_volume;
322 jpye 3463 C_pipe = c_pipe * m_pipe;
323 jpye 3469 C_insul = c_insul * V_insul;
324     C_sodium = c_sodium * m_sodium;
325 jpye 3463
326     c_inst_pipe_0, c_inst_insul_0, c_supp_0 IS_A cost_per_length;
327     c_inst_pipe_1, c_inst_insul_1 IS_A cost_per_area;
328 jpye 3469 c_sodium IS_A cost_per_mass;
329 jpye 3463
330     C_inst_pipe = SUM[PC[i].L * (c_inst_pipe_1*PC[i].D_o + c_inst_pipe_0)
331 jpye 3467 + PH[i].L * (c_inst_pipe_1*PH[i].D_o + c_inst_pipe_0) | i IN [0..n_PH-1]];
332 jpye 3463
333     C_inst_insul = SUM[PC[i].L * (c_inst_insul_1*PC[i].D_o + c_inst_insul_0)
334 jpye 3467 + PH[i].L * (c_inst_insul_1*PH[i].D_o + c_inst_insul_0) | i IN [0..n_PH-1]];
335 jpye 3463
336     C_supp = c_supp_0 * SUM[PC[i].L + PH[i].L | i IN [0..n_PH-1]];
337    
338     C_tot = C_pipe + C_insul + C_inst_pipe + C_inst_insul + C_supp;
339 jpye 3465
340 jpye 3463 T_amb IS_A temperature;
341     T_amb, PC[0..n_PH-1].T_amb, PH[0..n_PH-1].T_amb ARE_THE_SAME;
342    
343     Qdot_onecoll IS_A energy_rate;
344     Qdot_onecoll, DI[0..n_DI-1].Qdot ARE_THE_SAME;
345     DI[0..n_DI-1].eta ARE_THE_SAME; eta_DI ALIASES DI[1].eta;
346    
347     Qdot_coll_tot IS_A energy_rate;
348     Qdot_coll_tot = SUM[DI[i].Qdot | i IN [0..n_DI-1]];
349    
350     Qdot_loss_tot IS_A energy_rate;
351     Qdot_loss_tot = SUM[PC[i].Q + PH[i].Q|i IN [0..n_PH-1]];
352    
353     Qdot_net IS_A energy_rate;
354     Qdot_net = outlet.mdot*outlet.h - inlet.mdot*inlet.h;
355    
356     eta_th_array IS_A fraction;
357     eta_th_array = Qdot_net / Qdot_coll_tot;
358 jpye 3467
359 jpye 3463 (* try: same mass flow for each dish *)
360     DI[0..n_DI-1].mdot ARE_THE_SAME;
361     mdot_onecoll IS_A mass_rate;
362     mdot_onecoll, DI[0].mdot ARE_THE_SAME;
363    
364     T_in ALIASES inlet.T;
365     T_out ALIASES outlet.T;
366    
367     p_in ALIASES inlet.p;
368     p_out ALIASES outlet.p;
369 jpye 3470 dp IS_A delta_pressure;
370     dp = p_in - p_out;
371 jpye 3463
372     METHODS
373     METHOD on_load;
374     RUN default_self;
375    
376     RUN set_pipe_lengths;
377    
378 jpye 3466 FIX Qdot_onecoll := %f {MW}; (* value substituted in script! *)
379 jpye 3463 FIX eta_DI := 1;
380 jpye 3465
381 jpye 3463 FIX Vel := 15 {ft/s}; (* NAK Hbk v3 p6: 4-16" size, 15 *)
382    
383 jpye 3465 (*FIX t_pipe := 4 {mm}; (* FIXME implement equation here *)*)
384     FIX t_insul := 100 {mm};
385    
386 jpye 3463 FIX eps := 0.09 {mm};
387     FIX h_ext := 10 {W/m^2/K};
388     FIX emiss := 0.8;
389     FIX k_insul := 0.05 {W/m/K}; (* 'microporous insulation board' @ 700��C, https://is.gd/5j9Gkw *)
390     FIX k_pipe := 22.4 {W/m/K}; (* Haynes 230 @ 700��C: https://is.gd/F0RICh *)
391     FIX rho_pipe := 9 {g/cm^3}; (* Haynes 230 @ 700��C: https://is.gd/F0RICh *)
392 jpye 3469 FIX rho_insul := 380 {kg/m^3}; (* 'microporous insulation board', approx value from Felix *)
393 jpye 3463
394     FIX c_pipe := 8 {USD/kg} (* 316 stainless steel assumed *);
395 jpye 3469
396     FIX c_insul := 10000 {USD/m^3}; (* old quote from brandname product...or maybe 250 USD/kg if alibaba.com *)
397 jpye 3463 FIX c_inst_pipe_0 := 22 {AUD/m};
398     FIX c_inst_pipe_1 := 0.3734 {AUD/m/mm};
399     FIX c_inst_insul_0 := 29 {AUD/m};
400     FIX c_inst_insul_1 := 0.6257 {AUD/m/mm};
401     FIX c_supp_0 := 100 {USD/m}; (* this one's in yankee dollars *)
402 jpye 3469 FIX c_sodium := 3000 {USD/t}; (* based on a quick look at alibaba https://is.gd/pPKCYt *)
403 jpye 3465
404 jpye 3463 (*
405     installation cost for pipe, per length installed
406     0.3734 AUD/m/mm * d_o_pipe + 22 AUD/m
407    
408     installation cost for insulation, per length installed
409     0.6257 AUD/m/mm * d_o_pipe + 29 AUD/m
410    
411     anchors and supports
412     - supports with 2.5 m spacing, 175 AUD/support, plus
413     - 400 AUD/anchor, one per every second support.
414     totals to 100 USD/m.
415    
416     insulation cost
417     10000 USD/m^3 (or maybe 250....?)
418     *)
419    
420     FIX inlet.p := 50 {bar};
421 jpye 3465 FIX inlet.T := 520 {K} + 273.15 {K};
422 jpye 3466 FIX mdot_onecoll := %f {kg/s}; (* same for all dishes *)
423 jpye 3463 FIX T_amb := 300 {K};
424    
425     (* initial guesses to help the solver *)
426     PC[0..n_PH-1].inlet.h := 400 {kJ/kg};
427    
428     SOLVER QRSlv;
429     OPTION convopt 'RELNOM_SCALE';
430 jpye 3465 END on_load;
431     METHOD correct_dp_and_Tout;
432     FREE Vel;
433     FIX p_in := 5 {bar};
434     FIX p_out := 3 {bar};
435     FREE mdot_onecoll;
436     FIX T_out := 740 {K} + 273.15{K};
437     END correct_dp_and_Tout;
438 jpye 3469 METHOD show_temperatures;
439     EXTERNAL disharray_temperature_list(SELF);
440     END show_temperatures;
441 jpye 3467 """ % (Qdot_onecoll_MW,mdot_onecoll));
442 jpye 3463
443 jpye 3467 f.write("END towerarray{suffix};\n".format(suffix=suffix));
444 jpye 3463
445 jpye 3467 pl.figure()
446     pl.axes().set_aspect('equal', 'datalim')
447     pl.plot(X,Y,'bo')
448     pl.axis([x.min() - xsep/2, x.max() + xsep/2, y.min() - ysep/2, y.max() + ysep/2])
449     for lab in arr:
450     pl.annotate("%s\n%d"%(lab,rarr[lab]),(arr[lab].x,arr[lab].y),ha='right')
451 jpye 3463
452 jpye 3467 for l in links:
453     print "link",l
454     l.plot(arr,xsep,ysep)
455 jpye 3463
456 jpye 3467 if __name__=='__main__':
457 jpye 3463
458 jpye 3467 f = open("layout.a4c","w")
459     f.write("""
460     REQUIRE "johnpye/fprops/pipe.a4c";
461     """);
462 jpye 3463
463 jpye 3467 create_field_model(1,4,iy_pb=0,suffix="_chain1",f=f)
464    
465     create_field_model(5,1,ix_pb=0,Qdot_onecoll_MW=70,suffix="_chain2",f=f)
466    
467     f.close()
468    
469     pl.show()
470    
471    

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