source: trunk/sample/stage_separators/sample_column.mso @ 690

Last change on this file since 690 was 650, checked in by Argimiro Resende Secchi, 15 years ago

Fixing convergence problems!

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 22.6 KB
Line 
1#*-------------------------------------------------------------------
2* EMSO Model Library (EML) Copyright (C) 2004 - 2007 ALSOC.
3*
4* This LIBRARY is free software; you can distribute it and/or modify
5* it under the therms of the ALSOC FREE LICENSE as available at
6* http://www.enq.ufrgs.br/alsoc.
7*
8* EMSO Copyright (C) 2004 - 2007 ALSOC, original code
9* from http://www.rps.eng.br Copyright (C) 2002-2004.
10* All rights reserved.
11*
12* EMSO is distributed under the therms of the ALSOC LICENSE as
13* available at http://www.enq.ufrgs.br/alsoc.
14*
15*--------------------------------------------------------------------
16* Sample file for column model
17*--------------------------------------------------------------------
18*
19* This sample file needs VRTherm DEMO(www.vrtech.com.br) to run
20* SectionColumn_Test and needs VRTherm full to run the distillation
21* column flowsheet.
22*
23*----------------------------------------------------------------------
24* Author: Paula B. Staudt
25* $Id: sample_column.mso 650 2008-10-01 01:00:48Z arge $
26*--------------------------------------------------------------------*#
27
28using "stage_separators/column";
29using "controllers/PIDIncr";
30
31# column section with 2 trays
32FlowSheet SectionColumn_Test_with2tray
33        PARAMETERS
34        PP      as Plugin(Brief="Physical Properties",
35                Type="PP",
36                Components = [ "isobutane", "benzene"],
37                LiquidModel = "PR",
38                VapourModel = "PR"
39        );
40        NComp   as Integer;
41
42        SET
43        NComp = PP.NumberOfComponents;
44       
45        DEVICES
46        sec as Section_Column;
47        feed as liquid_stream;
48        reb as vapour_stream;
49        cond as liquid_stream;
50        zero as stream;
51       
52        CONNECTIONS
53        feed to sec.trays(2).Inlet;
54        zero to sec.trays(1).Inlet;
55        reb to sec.trays(2).InletV;
56        cond to sec.trays(1).InletL;
57       
58        SPECIFY
59        feed.F = 113.4 * 'kmol/h';
60        feed.T = 291 * 'K';
61        feed.P = 168.3 * 'kPa';
62        feed.z = [0.5, 0.5];
63
64        zero.F = 0 * 'kmol/h';
65        zero.T = 300 * 'K';
66        zero.P = 1 * 'atm';
67        zero.z = [0.5, 0.5];
68        zero.v = 0;
69        zero.h = 0 * 'J/mol';
70       
71        cond.F = 68 * 'kmol/h';
72        cond.P = 150 * 'kPa';
73        cond.T = 281.75 * 'K';
74        cond.z = [0.6664, 0.3336];
75
76#       reb.F = 153 * 'kmol/h';
77        reb.P = 185 * 'kPa';
78        reb.T = 328.12 * 'K';
79        reb.z = [0.001848, 0.9982];
80       
81        sec.trays.Emv = 1;
82        sec.trays(1).OutletV.F = 150 * 'kmol/h';
83
84        SET
85        sec.NTrays = 2;
86        #COLUMN
87        sec.trays.V = 4 * 'ft^3';
88        sec.trays.Ah = 0.394 * 'ft^2';
89        sec.trays.lw = 20.94 * 'in';
90        sec.trays.hw = 0.125 * 'ft';
91        sec.trays.Q = 0 * 'kW';
92        sec.trays.beta = 0.6;
93        sec.trays.alfa = 4;
94        sec.trays.Ap = 3.94 * 'ft^2';
95
96        INITIAL
97        sec.trays.OutletL.T = [290 300] *'K';
98        sec.trays.Level = 0.3 * sec.trays.hw;
99        sec.trays.OutletL.z(1) = 0.5;
100       
101        OPTIONS
102        TimeStep = 10;
103        TimeEnd = 1000;
104       
105# After running few seconds of transient the steady-state
106# can be obtained by using the results from that transient:
107
108        #GuessFile="SectionColumn_Test_with2tray.rlt";
109        #Dynamic = false;
110end
111
112# column section with 8 trays
113FlowSheet SectionColumn_Test_with8tray
114        PARAMETERS
115        PP      as Plugin(Brief="Physical Properties",
116                Type="PP",
117                Components = [ "isobutane", "benzene"],
118                LiquidModel = "PR",
119                VapourModel = "PR"
120        );
121        NComp   as Integer;
122
123        SET
124        NComp = PP.NumberOfComponents;
125       
126        DEVICES
127        sec as Section_Column;
128        feed as liquid_stream;
129        reb as vapour_stream;
130        cond as liquid_stream;
131        zero as stream;
132       
133        CONNECTIONS
134        feed to sec.trays(5).Inlet;
135       
136        zero to sec.trays([1:4]).Inlet;
137        zero to sec.trays([6:8]).Inlet;
138       
139        reb to sec.trays(8).InletV;
140        cond to sec.trays(1).InletL;
141       
142        SPECIFY
143        feed.F = 113.4 * 'kmol/h';
144        feed.T = 291 * 'K';
145        feed.P = 168.3 * 'kPa';
146        feed.z = [0.5, 0.5];
147
148        zero.F = 0 * 'kmol/h';
149        zero.T = 300 * 'K';
150        zero.P = 1 * 'atm';
151        zero.z = [0.5, 0.5];
152        zero.v = 0;
153        zero.h = 0 * 'J/mol';
154       
155        cond.F = 68 * 'kmol/h';
156        cond.P = 150 * 'kPa';
157        cond.T = 281.75 * 'K';
158        cond.z = [0.6664, 0.3336];
159
160#       reb.F = 153 * 'kmol/h';
161        reb.P = 185 * 'kPa';
162        reb.T = 328.12 * 'K';
163        reb.z = [0.001848, 0.9982];
164       
165        sec.trays.Emv = 1;
166        sec.trays(1).OutletV.F = 150 * 'kmol/h';
167
168        SET
169        sec.NTrays = 8;
170        #COLUMN
171        sec.trays.V = 4 * 'ft^3';
172        sec.trays.Ah = 0.394 * 'ft^2';
173        sec.trays.lw = 20.94 * 'in';
174        sec.trays.hw = 0.125 * 'ft';
175        sec.trays.Q = 0 * 'kW';
176        sec.trays.beta = 0.6;
177        sec.trays.alfa = 4;
178        sec.trays.Ap = 3.94 * 'ft^2';
179
180        INITIAL
181        sec.trays.OutletL.T = [290:(300-290)/(sec.NTrays-1):300] *'K';
182        sec.trays.Level = 0.3 * sec.trays.hw;
183        sec.trays.OutletL.z(1) = 0.5;
184       
185        OPTIONS
186        TimeStep = 1;
187        TimeEnd = 100;
188        #GuessFile="SectionColumn_Test_with8tray.rlt";
189        #Dynamic = false;
190end
191
192
193FlowSheet Distillation_kettle_cond_Test
194        PARAMETERS
195        PP      as Plugin(Brief="Physical Properties",
196                Type="PP",
197                Components = [ "isobutane", "n-pentane", "propylene",
198                "benzene", "isobutene" ],
199                LiquidModel = "PR",
200                VapourModel = "PR"
201        );
202        NComp   as Integer;
203       
204        VARIABLES
205        Qc as energy_source (Brief="Heat rate removed from condenser");
206        Qr as energy_source (Brief="Heat rate supplied to reboiler");
207       
208        SET
209        NComp = PP.NumberOfComponents;
210
211        DEVICES
212        col as Distillation_kettle_cond;
213        feed as source;
214        zero as stream;
215       
216        CONNECTIONS
217        feed.Outlet to col.trays(5).Inlet;
218        zero to col.reb.Inlet;
219        zero to col.trays([1:4]).Inlet;
220        zero to col.trays([6:col.NTrays]).Inlet;
221        Qc.OutletQ to col.cond.InletQ;
222        Qr.OutletQ to col.reb.InletQ;
223       
224        SPECIFY
225        feed.F = 113.4 * 'kmol/h';
226        feed.T = 291 * 'K';
227        feed.P = 168.3 * 'kPa';
228        feed.Composition = 1/NComp;
229       
230        zero.F = 0 * 'kmol/h';
231        zero.T = 300 * 'K';
232        zero.P = 1 * 'atm';
233        zero.z = 1/NComp;
234        zero.v = 0;
235        zero.h = 0 * 'J/mol';
236       
237        col.sptop.Outlet2.F = 85 * 'kmol/h';
238        col.reb.OutletL.F = 28.4 * 'kmol/h';
239        col.sptop.frac = 0.444445;
240        col.cond.OutletV.F = 0 * 'kmol/h';
241        Qr.OutletQ.Q = 3.7743e6 * 'kJ/h';
242        Qc.OutletQ.Q = -3.71e6 * 'kJ/h';
243        col.pump1.dP = 16 * 'kPa';
244        col.trays.Emv = 1;
245       
246        col.alfaTopo = 2;
247       
248        SET
249        col.NTrays = 8;
250        col.cond.V = 2 * 'm^3';
251        col.cond.Across = 1 * 'm^2';
252        col.trays.V = 4 * 'ft^3';
253        col.trays.Ah = 0.394 * 'ft^2';
254        col.trays.lw = 20.94 * 'in';
255        col.trays.hw = 0.125 * 'ft';
256        col.trays.Q = 0 * 'kW';
257        col.trays.beta = 0.6;
258        col.trays.alfa = 4;
259        col.trays.Ap = 3.94 * 'ft^2';
260        col.reb.V = 2 * 'm^3';
261        col.reb.Across = 1 * 'm^2';
262       
263        INITIAL
264        # condenser
265        col.cond.OutletL.T = 260 *'K';
266        col.cond.Level = 1 * 'm';
267        col.cond.OutletL.z([1:4]) = [0.65, 0.05, 0.01, 0.01];
268
269        # reboiler
270        col.reb.OutletL.T = 350 *'K';
271        col.reb.Level = 1 * 'm';
272        col.reb.OutletL.z([1:4]) = [0.1, 0.7, 0.01, 0.01];
273
274        # column trays
275        col.trays.OutletL.T = [290:(330-290)/(col.NTrays-1):330] * 'K';
276        col.trays.Level = 1.2 * col.trays.hw;
277        col.trays.OutletL.z([1:4]) = [0.5, 0.05, 0.01, 0.01];
278
279        OPTIONS
280        TimeStep = 0.1;
281        TimeEnd = 50;
282        #GuessFile="Distillation_kettle_cond_Test.rlt";
283        #Dynamic = false;       
284end
285
286
287FlowSheet Column_ctrl
288        PARAMETERS
289        PP      as Plugin(Brief="Physical Properties",
290                Type="PP",
291                Components = [ "isobutane", "n-pentane", "propylene",
292                        "benzene", "isobutene" ],
293                LiquidModel = "PR",
294                VapourModel = "PR"
295        );
296        NComp   as Integer;
297       
298        Qcmin as heat_rate (Brief="Minimum Condenser Heat supplied");
299        Qcmax as heat_rate (Brief="Maximum Condenser Heat supplied");
300        Qrmin as heat_rate (Brief="Minimum Reboiler Heat supplied");
301        Qrmax as heat_rate (Brief="Maximum Reboiler Heat supplied");
302        Frmin as flow_mol (Brief="Minimum bottom flow rate");
303        Frmax as flow_mol (Brief="Maximum bottom flow rate");
304        Fcmin as flow_mol (Brief="Minimum reflux flow rate");
305        Fcmax as flow_mol (Brief="Maximum reflux flow rate");
306        Hmint as length (Brief="Minimum liquid level in top tank");
307        Hmaxt as length (Brief="Maximum liquid level in top tank");
308    Hminb as length (Brief="Minimum liquid level in reboiler");
309        Hmaxb as length (Brief="Maximum liquid level in reboiler");
310        Pmax as pressure (Brief="Maximum column pressure");
311        Pmin as pressure (Brief="Minimum column pressure");
312        Tmax as temperature (Brief="Maximum column temperature");
313        Tmin as temperature (Brief="Minimum column temperature");
314
315        VARIABLES
316        Qc as energy_source (Brief="Heat rate removed from condenser");
317        Qr as energy_source (Brief="Heat rate supplied to reboiler");
318        Had_top as Real (Brief="Dimensionless condenser level");
319        Had_bot as Real (Brief="Dimensionless reboiler level");
320        Pad as Real (Brief="Dimensionless pressure");
321        Tad as Real (Brief="Dimensionless temperature");
322        RR      as positive (Brief="Reflux ratio");
323
324        SET
325        NComp = PP.NumberOfComponents;
326
327        DEVICES
328        col as Distillation_kettle_cond;
329        feed as source;
330        zero as stream;
331        TCcond as PIDIncr;
332        LCtop as PIDIncr;
333        LCbot as PIDIncr;
334        PC as PIDIncr;
335
336        CONNECTIONS
337        feed.Outlet to col.trays(5).Inlet;
338        zero to col.reb.Inlet;
339        zero to col.trays([1:4]).Inlet;
340        zero to col.trays([6:col.NTrays]).Inlet;
341        Qc.OutletQ to col.cond.InletQ;
342        Qr.OutletQ to col.reb.InletQ;
343
344        EQUATIONS
345   "Temperature Controller"
346        TCcond.Parameters.tau = 0*'s'; 
347        TCcond.Parameters.tauSet = 0*'s';       
348        TCcond.Parameters.alpha = 0.3;
349        TCcond.Parameters.bias = 0.5;   
350        TCcond.Parameters.gamma = 1;
351        TCcond.Parameters.beta = 1;
352        TCcond.Options.action = 1;
353        TCcond.Options.clip = 1;
354        TCcond.Options.autoMan = 0;
355        TCcond.Parameters.intTime = 60*'s';
356        TCcond.Parameters.gain = 0.6;
357        TCcond.Parameters.derivTime = 1*'s';
358        TCcond.Ports.setPoint = ((15+273.15) * 'K' - Tmin)/(Tmax-Tmin);
359        TCcond.Ports.input = Tad;
360        Tad = (col.cond.OutletL.T-Tmin)/(Tmax-Tmin);
361        Qc.OutletQ.Q = Qcmin+(Qcmax-Qcmin)*TCcond.Ports.output;
362
363        "Pressure Controller"
364        PC.Parameters.tau = 0*'s';     
365        PC.Parameters.tauSet = 0*'s';   
366        PC.Parameters.alpha = 0.3;
367        PC.Parameters.bias = 0;
368        PC.Parameters.gamma = 1;
369        PC.Parameters.beta = 1;
370        PC.Options.action = -1;
371        PC.Options.clip = 1;
372        PC.Options.autoMan = 0;
373        PC.Parameters.intTime = 50*'s';
374        PC.Parameters.gain = 0.5;
375        PC.Parameters.derivTime = 1*'s';
376        PC.Ports.setPoint = (4.0*'bar'-Pmin)/(Pmax-Pmin);
377        PC.Ports.input = Pad;
378        Pad = (col.cond.OutletV.P-Pmin)/(Pmax-Pmin);
379        col.cond.OutletV.F = (Fcmin+(Fcmax-Fcmin)*PC.Ports.output);     
380       
381        "Ttop Level Controller"
382        LCtop.Parameters.tau = 0*'s';   
383        LCtop.Parameters.tauSet = 0*'s';       
384        LCtop.Parameters.alpha = 0.3;
385        LCtop.Parameters.bias = 0.5;   
386        LCtop.Parameters.gamma = 1;
387        LCtop.Parameters.beta = 1;
388        LCtop.Options.action = -1;
389        LCtop.Options.clip = 1;
390        LCtop.Options.autoMan = 0;
391        LCtop.Parameters.intTime = 10*'s';
392        LCtop.Parameters.gain = 1;
393        LCtop.Parameters.derivTime = 0*'s';
394        LCtop.Ports.setPoint = (1.0 * 'm' - Hmint)/(Hmaxt-Hmint);
395        LCtop.Ports.input = Had_top;
396        Had_top = (col.cond.Level-Hmint)/(Hmaxt-Hmint);
397        col.sptop.Outlet1.F = Fcmin + (Fcmax-Fcmin) * LCtop.Ports.output;
398
399        "Tbottom Level Controller"
400        LCbot.Parameters.tau = 0*'s';   
401        LCbot.Parameters.tauSet = 0*'s';       
402        LCbot.Parameters.alpha = 0.3;
403        LCbot.Parameters.bias = 0.5;   
404        LCbot.Parameters.gamma = 1;
405        LCbot.Parameters.beta = 1;
406        LCbot.Options.action = -1;
407        LCbot.Options.clip = 1;
408        LCbot.Options.autoMan = 0;
409        LCbot.Parameters.intTime = 100*'s';
410        LCbot.Parameters.gain = 1;
411        LCbot.Parameters.derivTime = 0*'s';
412        LCbot.Ports.setPoint = (1.0 * 'm' - Hminb)/(Hmaxb-Hminb);
413        LCbot.Ports.input = Had_bot;
414        Had_bot = (col.reb.Level-Hminb)/(Hmaxb-Hminb);
415        col.reb.OutletL.F = Frmin + (Frmax-Frmin) * LCbot.Ports.output;
416
417        RR * (col.cond.OutletV.F + col.sptop.Outlet1.F) = col.sptop.Outlet2.F;
418       
419        if time < 1 * 'h' then
420                col.sptop.Outlet2.F = 75 * 'kmol/h'; # reflux
421        else
422                col.sptop.Outlet2.F = 85 * 'kmol/h'; # reflux
423        end
424
425        SPECIFY
426        feed.F = 113.4 * 'kmol/h';
427        feed.T = 291 * 'K';
428        feed.P = 5 * 'bar';
429        feed.Composition = 1/NComp;
430       
431        zero.F = 0 * 'kmol/h';
432        zero.T = 300 * 'K';
433        zero.P = 1 * 'atm';
434        zero.z = 1/NComp;
435        zero.v = 0;
436        zero.h = 0 * 'J/mol';
437       
438        Qr.OutletQ.Q = 4e6 * 'kJ/h';
439        col.pump1.dP = 16 * 'kPa';
440        col.trays.Emv = 1;
441        col.alfaTopo = 2;
442
443        SET
444        col.NTrays = 8;
445        col.cond.V = 2 * 'm^3';
446        col.cond.Across = 1 * 'm^2';
447        col.trays.V = 4 * 'ft^3';
448        col.trays.Ah = 0.394 * 'ft^2';
449        col.trays.lw = 20.94 * 'in';
450        col.trays.hw = 0.125 * 'ft';
451        col.trays.Q = 0 * 'kW';
452        col.trays.beta = 0.6;
453        col.trays.alfa = 4;
454        col.trays.Ap = 3.94 * 'ft^2';
455        col.reb.V = 2 * 'm^3';
456        col.reb.Across = 1 * 'm^2';
457
458        # Controllers type
459        TCcond.PID_Select = "Ideal_AW";
460        PC.PID_Select = "Ideal_AW";
461        LCtop.PID_Select = "Ideal_AW";
462        LCbot.PID_Select = "Ideal_AW";
463
464        Qrmax = 5e6 * 'kJ/h';
465        Qrmin = 1e6 * 'kJ/h';
466        Frmin = 0 * 'kmol/h';
467        Frmax = 60 * 'kmol/h';
468        Fcmin = 0 * 'kmol/h';
469        Fcmax = 120 * 'kmol/h';
470        Hmint = 0 * 'm';
471        Hmaxt = 2 * 'm';
472        Hminb = 0 * 'm';
473        Hmaxb = 2 * 'm';
474        Pmin = 0.5 * 'bar';
475        Pmax = 6 * 'bar';
476        Qcmax = -5e5 * 'kJ/h';
477        Qcmin = -5e6 * 'kJ/h';
478        Tmax = (30+273.15) * 'K';
479        Tmin = (-20+273.15) * 'K';
480       
481        INITIAL
482        # condenser
483        col.cond.OutletL.T = 260 *'K';
484        col.cond.Level = 1 * 'm';
485        col.cond.OutletL.z([1:4]) = [0.2, 0.2, 0.4, 0.05];
486
487        # reboiler
488        col.reb.OutletL.T = 350 *'K';
489        col.reb.Level = 1 * 'm';
490        col.reb.OutletL.z([1:4]) = [0.1, 0.4, 0.1, 0.3];
491
492        # column trays
493        col.trays.OutletL.T = [290:(330-290)/(col.NTrays-1):330] * 'K';
494        col.trays.Level = 1.2 * col.trays.hw;
495        col.trays.OutletL.z([1:4]) = [0.15, 0.3, 0.25, 0.2];
496
497        OPTIONS
498        TimeStep = 0.1;
499        TimeEnd = 5;
500        TimeUnit = 'h';
501        InitialFile = "Column_ctrl.rlt";
502        DAESolver(File="dassl");
503        #GuessFile = "Column_ctrl.rlt";
504        #Dynamic = false;
505end
506
507# packed column section with 8 trays
508FlowSheet PackedSectionColumn
509        PARAMETERS
510        PP      as Plugin(Brief="Physical Properties",
511                Type="PP",
512                Components = [ "isobutane", "n-pentane", "propylene",
513                "benzene", "isobutene" ],
514                LiquidModel = "PR",
515                VapourModel = "PR"
516        );
517        NComp   as Integer;
518
519        SET
520        NComp = PP.NumberOfComponents;
521       
522        DEVICES
523        sec as Packed_Section_Column;
524        feed as liquid_stream;
525        reb as vapour_stream;
526        cond as liquid_stream;
527        zero as stream;
528       
529        CONNECTIONS
530        feed to sec.stage(5).Inlet;
531        zero to sec.stage([1:4]).Inlet;
532        zero to sec.stage([6:sec.NStages]).Inlet;
533        reb to sec.stage(8).InletV;
534        cond to sec.stage(1).InletL;
535       
536        SPECIFY
537        feed.F = 113.4 * 'kmol/h';
538        feed.T = 291 * 'K';
539        feed.P = 1.66 * 'atm';
540        feed.z = [0.2, 0.2, 0.2, 0.2, 0.2];
541
542        zero.F = 0 * 'kmol/h';
543        zero.T = 300 * 'K';
544        zero.P = 1 * 'atm';
545        zero.z = [0.2, 0.2, 0.2, 0.2, 0.2];
546        zero.v = 0;
547        zero.h = 0 * 'J/mol';
548       
549        cond.F = 85 * 'kmol/h';
550        cond.P = 2.33 * 'atm';
551        cond.T = 283.5 * 'K';
552        cond.z = [0.599, 0.044, 0.035, 0.007, 0.315];
553
554        reb.F = 137.57 * 'kmol/h';
555        reb.P = 2.46 * 'atm';
556        reb.T = 325 * 'K';
557        reb.z = [0.16, 0.542, 0.013, 0.008, 0.277];
558
559        sec.dP = 0.008 * 'atm'; #queda de pressao na secao!
560
561        SET
562        sec.H = 4 * 'm';#altura do recheio
563        sec.NStages = 8; #numero de estagios teoricos do recheio
564        sec.stage.Q = 0 * 'kW';
565        sec.stage.d = 1.009 * 'm';
566        sec.stage.Cpo = 0.763;
567        sec.stage.e = 0.951;
568        sec.stage.a = 112.6 * 'm^2/m^3';
569        sec.stage.Qsil = 0.1;#coeficiente de resistencia do recheio
570
571        INITIAL
572        sec.stage.OutletL.T =[283:(325-283)/(sec.NStages-1):325] *'K';
573        sec.stage.ML = 0.01 * 'kmol';
574        sec.stage.OutletL.z([1:4]) = [0.3, 0.2, 0.002, 0.1];
575
576        OPTIONS
577        DAESolver(File="dassl");
578        NLASolver(File="nlasolver");
579        TimeStep = 10;
580        TimeEnd = 100;
581end
582
583FlowSheet Packed_kettle_cond_Test
584        PARAMETERS
585        PP      as Plugin(Brief="Physical Properties",
586                Type="PP",
587                Components = [ "isobutane", "n-pentane", "propylene",
588                "benzene", "isobutene" ],
589                LiquidModel = "PR",
590                VapourModel = "PR"
591        );
592        NComp   as Integer;
593       
594        VARIABLES
595        Qc as energy_source (Brief="Heat rate removed from condenser");
596        Qr as energy_source (Brief="Heat rate supplied to reboiler");
597       
598        SET
599        NComp = PP.NumberOfComponents;
600
601        DEVICES
602        col as PackedDistillation_kettle_cond;
603        feed as source;
604        zero as stream;
605       
606        CONNECTIONS
607        feed.Outlet to col.stage(5).Inlet;
608        zero to col.reb.Inlet;
609        zero to col.stage([1:4]).Inlet;
610        zero to col.stage([6:col.NStages]).Inlet;
611        Qc.OutletQ to col.cond.InletQ;
612        Qr.OutletQ to col.reb.InletQ;
613       
614        SPECIFY
615        feed.F = 113.4 * 'kmol/h';
616        feed.T = 291 * 'K';
617        feed.P = 168.3 * 'kPa';
618        feed.Composition = 1/NComp;
619       
620        zero.F = 0 * 'kmol/h';
621        zero.T = 300 * 'K';
622        zero.P = 1 * 'atm';
623        zero.z = 1/NComp;
624        zero.v = 0;
625        zero.h = 0 * 'J/mol';
626       
627        col.sptop.Outlet2.F = 85 * 'kmol/h';
628        col.reb.OutletL.F = 28.4 * 'kmol/h';
629        col.sptop.frac = 0.444445;
630        col.cond.OutletV.F = 0 * 'kmol/h';
631        Qr.OutletQ.Q = 3.7743e6 * 'kJ/h';
632        Qc.OutletQ.Q = -3.71e6 * 'kJ/h';
633        col.pump1.dP = 16 * 'kPa';
634
635        col.dP = 0.024 * 'atm'; #queda de pressao entre o topo e o fundo da coluna.
636
637        SET
638        col.H = 3.5 * 'm'; # altura de recheio
639        col.NStages = 8; # numero de estagios teoricos do recheio
640        col.stage.Q = 0 * 'kW';
641        col.stage.d = 2.24 * 'ft';
642        col.stage.Cpo = 0.763;
643        col.stage.e = 0.951;
644        col.stage.a = 112.6 * 'm^2/m^3';
645
646        col.cond.V = 2 * 'm^3';
647        col.cond.Across = 1 * 'm^2';
648        col.reb.V = 2 * 'm^3';
649        col.reb.Across = 1 * 'm^2';
650        col.stage.Qsil = 2; # coeficiente de resistencia do recheio
651       
652        INITIAL
653        # condenser
654        col.cond.OutletL.T = 260 *'K';
655        col.cond.Level = 1 * 'm';
656        col.cond.OutletL.z([1:4]) = [0.65, 0.05, 0.01, 0.01];
657
658        # reboiler
659        col.reb.OutletL.T = 330 *'K';
660        col.reb.Level = 1 * 'm';
661        col.reb.OutletL.z([1:4]) = [0.1, 0.7, 0.01, 0.01];
662
663        # column trays
664        col.stage.OutletL.T = [290:(330-290)/(col.NStages-1):330] * 'K';
665        col.stage.ML = 0.1 * 'kmol';
666        col.stage.OutletL.z([1:4]) = [0.15, 0.5, 0.001, 0.1];
667
668        OPTIONS
669        TimeStep = 10;
670        TimeEnd = 500;
671end
672
673FlowSheet PackedColumn_ctrl
674        PARAMETERS
675        PP      as Plugin(Brief="Physical Properties",
676                Type="PP",
677                Components = [ "isobutane", "n-pentane", "propylene",
678                        "benzene", "isobutene" ],
679                LiquidModel = "PR",
680                VapourModel = "PR"
681        );
682        NComp   as Integer;
683       
684        Qcmin as heat_rate (Brief="Minimum Condenser Heat supplied");
685        Qcmax as heat_rate (Brief="Maximum Condenser Heat supplied");
686        Qrmin as heat_rate (Brief="Minimum Reboiler Heat supplied");
687        Qrmax as heat_rate (Brief="Maximum Reboiler Heat supplied");
688        Frmin as flow_mol (Brief="Minimum bottom flow rate");
689        Frmax as flow_mol (Brief="Maximum bottom flow rate");
690        Fcmin as flow_mol (Brief="Minimum reflux flow rate");
691        Fcmax as flow_mol (Brief="Maximum reflux flow rate");
692        Hmint as length (Brief="Minimum liquid level in top tank");
693        Hmaxt as length (Brief="Maximum liquid level in top tank");
694    Hminb as length (Brief="Minimum liquid level in reboiler");
695        Hmaxb as length (Brief="Maximum liquid level in reboiler");
696        Pmax as pressure (Brief="Maximum column pressure");
697        Pmin as pressure (Brief="Minimum column pressure");
698        Tmax as temperature (Brief="Maximum column temperature");
699        Tmin as temperature (Brief="Minimum column temperature");
700
701        VARIABLES
702        Qc as energy_source (Brief="Heat rate removed from condenser");
703        Qr as energy_source (Brief="Heat rate supplied to reboiler");
704        Had_top as Real (Brief="Dimensionless condenser level");
705        Had_bot as Real (Brief="Dimensionless reboiler level");
706        Pad as Real (Brief="Dimensionless pressure");
707        Tad as Real (Brief="Dimensionless temperature");
708        RR      as positive (Brief="Reflux ratio");
709
710        SET
711        NComp = PP.NumberOfComponents;
712
713        DEVICES
714        col as PackedDistillation_kettle_cond;
715        feed as source;
716        zero as stream;
717        TCcond as PIDIncr;
718        LCtop as PIDIncr;
719        LCbot as PIDIncr;
720        PC as PIDIncr;
721
722        CONNECTIONS
723        feed.Outlet to col.stage(5).Inlet;
724        zero to col.reb.Inlet;
725        zero to col.stage([1:4]).Inlet;
726        zero to col.stage([6:col.NStages]).Inlet;
727        Qc.OutletQ to col.cond.InletQ;
728        Qr.OutletQ to col.reb.InletQ;
729
730        EQUATIONS
731   "Temperature Controller"
732        TCcond.Parameters.tau = 0*'s'; 
733        TCcond.Parameters.tauSet = 0*'s';       
734        TCcond.Parameters.alpha = 0.3;
735        TCcond.Parameters.bias = 0.5;   
736        TCcond.Parameters.gamma = 1;
737        TCcond.Parameters.beta = 1;
738        TCcond.Options.action = 1;
739        TCcond.Options.clip = 1;
740        TCcond.Options.autoMan = 0;
741        TCcond.Parameters.intTime = 60*'s';
742        TCcond.Parameters.gain = 0.6;
743        TCcond.Parameters.derivTime = 1*'s';
744        TCcond.Ports.setPoint = ((15+273.15) * 'K' - Tmin)/(Tmax-Tmin);
745        TCcond.Ports.input = Tad;
746        Tad = (col.cond.OutletL.T-Tmin)/(Tmax-Tmin);
747        Qc.OutletQ.Q = Qcmin+(Qcmax-Qcmin)*TCcond.Ports.output;
748
749        "Pressure Controller"
750        PC.Parameters.tau = 0*'s';     
751        PC.Parameters.tauSet = 0*'s';   
752        PC.Parameters.alpha = 0.3;
753        PC.Parameters.bias = 0;
754        PC.Parameters.gamma = 1;
755        PC.Parameters.beta = 1;
756        PC.Options.action = -1;
757        PC.Options.clip = 1;
758        PC.Options.autoMan = 0;
759        PC.Parameters.intTime = 50*'s';
760        PC.Parameters.gain = 0.5;
761        PC.Parameters.derivTime = 1*'s';
762        PC.Ports.setPoint = (4.0*'bar'-Pmin)/(Pmax-Pmin);
763        PC.Ports.input = Pad;
764        Pad = (col.cond.OutletV.P-Pmin)/(Pmax-Pmin);
765        col.cond.OutletV.F = (Fcmin+(Fcmax-Fcmin)*PC.Ports.output);     
766       
767        "Ttop Level Controller"
768        LCtop.Parameters.tau = 0*'s';   
769        LCtop.Parameters.tauSet = 0*'s';       
770        LCtop.Parameters.alpha = 0.3;
771        LCtop.Parameters.bias = 0.5;   
772        LCtop.Parameters.gamma = 1;
773        LCtop.Parameters.beta = 1;
774        LCtop.Options.action = -1;
775        LCtop.Options.clip = 1;
776        LCtop.Options.autoMan = 0;
777        LCtop.Parameters.intTime = 10*'s';
778        LCtop.Parameters.gain = 1;
779        LCtop.Parameters.derivTime = 0*'s';
780        LCtop.Ports.setPoint = (1.0 * 'm' - Hmint)/(Hmaxt-Hmint);
781        LCtop.Ports.input = Had_top;
782        Had_top = (col.cond.Level-Hmint)/(Hmaxt-Hmint);
783        col.sptop.Outlet1.F = Fcmin + (Fcmax-Fcmin) * LCtop.Ports.output;
784
785        "Tbottom Level Controller"
786        LCbot.Parameters.tau = 0*'s';   
787        LCbot.Parameters.tauSet = 0*'s';       
788        LCbot.Parameters.alpha = 0.3;
789        LCbot.Parameters.bias = 0.5;   
790        LCbot.Parameters.gamma = 1;
791        LCbot.Parameters.beta = 1;
792        LCbot.Options.action = -1;
793        LCbot.Options.clip = 1;
794        LCbot.Options.autoMan = 0;
795        LCbot.Parameters.intTime = 100*'s';
796        LCbot.Parameters.gain = 1;
797        LCbot.Parameters.derivTime = 0*'s';
798        LCbot.Ports.setPoint = (1.0 * 'm' - Hminb)/(Hmaxb-Hminb);
799        LCbot.Ports.input = Had_bot;
800        Had_bot = (col.reb.Level-Hminb)/(Hmaxb-Hminb);
801        col.reb.OutletL.F = Frmin + (Frmax-Frmin) * LCbot.Ports.output;
802
803        RR * (col.cond.OutletV.F + col.sptop.Outlet1.F) = col.sptop.Outlet2.F;
804       
805        if time < 1 * 'h' then
806                col.sptop.Outlet2.F = 75 * 'kmol/h'; # reflux
807        else
808                col.sptop.Outlet2.F = 85 * 'kmol/h'; # reflux
809        end
810
811        SPECIFY
812        feed.F = 113.4 * 'kmol/h';
813        feed.T = 291 * 'K';
814        feed.P = 5 * 'bar';
815        feed.Composition = 1/NComp;
816       
817        zero.F = 0 * 'kmol/h';
818        zero.T = 300 * 'K';
819        zero.P = 1 * 'atm';
820        zero.z = 1/NComp;
821        zero.v = 0;
822        zero.h = 0 * 'J/mol';
823       
824        Qr.OutletQ.Q = 4e6 * 'kJ/h';
825        col.pump1.dP = 16 * 'kPa';
826       
827        col.dP = 0.024  * 'atm'; #queda de pressao entre o topo e o fundo da coluna
828
829        SET
830        col.H = 3.5 * 'm'; #altura do recheio
831        col.NStages = 8; #numero de estagios teoricos do recheio
832        col.stage.Q = 0 * 'kW';
833        col.stage.d = 2.24 * 'ft';
834        col.stage.Cpo = 0.763;
835        col.stage.e = 0.951;
836        col.stage.a = 112.6 * 'm^2/m^3';
837
838        col.cond.V = 2 * 'm^3';
839        col.cond.Across = 1 * 'm^2';
840        col.reb.V = 2 * 'm^3';
841        col.reb.Across = 1 * 'm^2';
842
843        col.stage.Qsil = 2; #coeficiente de resistencia do recheio
844
845        # Controllers type
846        TCcond.PID_Select = "Ideal_AW";
847        PC.PID_Select = "Ideal_AW";
848        LCtop.PID_Select = "Ideal_AW";
849        LCbot.PID_Select = "Ideal_AW";
850
851        Qrmax = 5e6 * 'kJ/h';
852        Qrmin = 1e6 * 'kJ/h';
853        Frmin = 0 * 'kmol/h';
854        Frmax = 60 * 'kmol/h';
855        Fcmin = 0 * 'kmol/h';
856        Fcmax = 120 * 'kmol/h';
857        Hmint = 0 * 'm';
858        Hmaxt = 2 * 'm';
859        Hminb = 0 * 'm';
860        Hmaxb = 2 * 'm';
861        Pmin = 0.5 * 'bar';
862        Pmax = 6 * 'bar';
863        Qcmax = -5e5 * 'kJ/h';
864        Qcmin = -5e6 * 'kJ/h';
865        Tmax = (30+273.15) * 'K';
866        Tmin = (-20+273.15) * 'K';
867       
868        INITIAL
869        # condenser
870        col.cond.OutletL.T = 260 *'K';
871        col.cond.Level = 1 * 'm';
872        col.cond.OutletL.z([1:4]) = [0.2, 0.2, 0.4, 0.05];
873
874        # reboiler
875        col.reb.OutletL.T = 350 *'K';
876        col.reb.Level = 1 * 'm';
877        col.reb.OutletL.z([1:4]) = [0.1, 0.4, 0.1, 0.3];
878
879        # column trays
880        col.stage.OutletL.T = [290:(330-290)/(col.NStages-1):330] * 'K';
881        col.stage.ML = 0.1 * 'kmol';
882        col.stage.OutletL.z([1:4]) = [0.15, 0.3, 0.25, 0.2];
883
884        OPTIONS
885        TimeStep = 0.1;
886        TimeEnd = 5;
887        TimeUnit = 'h';
888        InitialFile = "Packed_Column_ctrl.rlt";
889        #GuessFile = "Packed_Column_ctrl.rlt";
890        #Dynamic = false;
891        DAESolver(File="dassl");
892end
Note: See TracBrowser for help on using the repository browser.