source: mso/sample/stage_separators/sample_column.mso @ 112

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

Add perturbation in the Column_ctrl example.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 12.4 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 112 2007-01-15 01:01:46Z 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 CalcObject(Brief="Physical Properties",File="vrpp");
35        NComp   as Integer;
36
37        SET
38        PP.Components = [ "isobutane", "benzene"];
39        PP.LiquidModel = "PR";
40        PP.VapourModel = "PR";
41        NComp = PP.NumberOfComponents;
42       
43        DEVICES
44        sec as Section_Column;
45        feed as stream_therm;
46        reb as stream_therm;
47        cond as stream_therm;
48        zero as stream;
49       
50        CONNECTIONS
51        feed to sec.trays(2).Inlet;
52        zero to sec.trays(1).Inlet;
53        reb to sec.trays(2).InletV;
54        cond to sec.trays(1).InletL;
55       
56        SPECIFY
57        feed.F = 113.4 * "kmol/h";
58        feed.T = 291 * "K";
59        feed.P = 168.3 * "kPa";
60        feed.z = [0.5, 0.5];
61        feed.v = 0;
62
63        zero.F = 0 * "kmol/h";
64        zero.T = 300 * "K";
65        zero.P = 1 * "atm";
66        zero.z = [0.5, 0.5];
67        zero.v = 0;
68        zero.h = 0 * "J/mol";
69       
70        cond.F = 68 * "kmol/h";
71        cond.P = 150 * "kPa";
72        cond.T = 281.75 * "K";
73        cond.v = 0.0;
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        reb.v = 1.0;
81       
82        sec.trays.Emv = 1;
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 *"K";
98        sec.trays.Level = 0.9 * sec.trays.hw;
99        sec.trays.OutletL.z(1) = 0.5;
100       
101        OPTIONS
102        #relativeAccuracy = 1e-3;
103        NLASolver = "sundials";
104        time = [0:10:1000];
105       
106# After running few seconds of transient the steady-state
107# can be obtained by using the results from that transient:
108
109        #guessFile="SectionColumn_Test_with2tray.rlt";
110        #mode = "steady";
111end
112
113# column section with 8 trays
114FlowSheet SectionColumn_Test_with8tray
115        PARAMETERS
116        PP      as CalcObject(Brief="Physical Properties",File="vrpp");
117        NComp   as Integer;
118
119        SET
120        PP.Components = [ "isobutane", "benzene"];
121        PP.LiquidModel = "PR";
122        PP.VapourModel = "PR";
123        PP.Derivatives = 0;
124        NComp = PP.NumberOfComponents;
125       
126        DEVICES
127        sec as Section_Column;
128        feed as stream_therm;
129        reb as stream_therm;
130        cond as stream_therm;
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        feed.v = 0;
148
149        zero.F = 0 * "kmol/h";
150        zero.T = 300 * "K";
151        zero.P = 1 * "atm";
152        zero.z = [0.5, 0.5];
153        zero.v = 0;
154        zero.h = 0 * "J/mol";
155       
156        cond.F = 68 * "kmol/h";
157        cond.P = 150 * "kPa";
158        cond.T = 281.75 * "K";
159        cond.v = 0.0;
160        cond.z = [0.6664, 0.3336];
161
162        reb.F = 153 * "kmol/h";
163        reb.P = 185 * "kPa";
164        reb.T = 328.12 * "K";
165        reb.z = [0.001848, 0.9982];
166        reb.v = 1.0;
167       
168        sec.trays.Emv = 1;
169
170        SET
171        sec.NTrays = 8;
172        #COLUMN
173        sec.trays.V = 4 * "ft^3";
174        sec.trays.Ah = 0.394 * "ft^2";
175        sec.trays.lw = 20.94 * "in";
176        sec.trays.hw = 0.125 * "ft";
177        sec.trays.Q = 0 * "kW";
178        sec.trays.beta = 0.6;
179        sec.trays.alfa = 4;
180        sec.trays.Ap = 3.94 * "ft^2";
181
182        INITIAL
183        sec.trays.OutletL.T = [290:(330-290)/(sec.NTrays-1):330] *"K";
184        sec.trays.Level = 0.3 * sec.trays.hw;
185        sec.trays.OutletL.z(1) = 0.5;
186       
187        OPTIONS
188        relativeAccuracy = 1e-5;
189        time = [0:1:100];
190        #guessFile="SectionColumn_Test_with8tray.rlt";
191        #mode = "steady";
192end
193
194
195FlowSheet Distillation_kettle_cond_Test
196        PARAMETERS
197        PP      as CalcObject(Brief="Physical Properties",File="vrpp");
198        NComp   as Integer;
199       
200        VARIABLES
201        Qc as heat_rate (Brief="Heat rate removed from condenser");
202        Qr as heat_rate (Brief="Heat rate supplied to reboiler");
203       
204        SET
205        PP.Components = [ "isobutane", "n-pentane", "propylene",
206                "benzene", "isobutene" ];
207        PP.LiquidModel = "PR";
208        PP.VapourModel = "PR";
209        PP.Derivatives = 1;
210        NComp = PP.NumberOfComponents;
211
212        DEVICES
213        col as Distillation_kettle_cond;
214        feed as streamTP;
215        zero as stream;
216       
217        CONNECTIONS
218        feed to col.trays(5).Inlet;
219        zero to col.reb.Inlet;
220        zero to col.trays([1:4]).Inlet;
221        zero to col.trays([6:col.NTrays]).Inlet;
222        Qc to col.cond.Q;
223        Qr to col.reb.Q;
224       
225        SPECIFY
226        feed.F = 113.4 * "kmol/h";
227        feed.T = 291 * "K";
228        feed.P = 168.3 * "kPa";
229        feed.z = 1/NComp;
230       
231        zero.F = 0 * "kmol/h";
232        zero.T = 300 * "K";
233        zero.P = 1 * "atm";
234        zero.z = 1/NComp;
235        zero.v = 0;
236        zero.h = 0 * "J/mol";
237       
238        col.sptop.Outlet2.F = 85 * "kmol/h";
239        col.reb.OutletL.F = 28.4 * "kmol/h";
240        col.sptop.frac = 0.444445;
241        col.cond.OutletV.F = 0 * "kmol/h";
242        Qr = 3.7743e6 * "kJ/h";
243        Qc = -3.71e6 * "kJ/h";
244        col.pump1.dP = 16 * "kPa";
245        col.trays.Emv = 1;
246       
247        SET
248        col.NTrays = 8;
249        col.cond.V = 2 * "m^3";
250        col.cond.Across = 1 * "m^2";
251        col.trays.V = 4 * "ft^3";
252        col.trays.Ah = 0.394 * "ft^2";
253        col.trays.lw = 20.94 * "in";
254        col.trays.hw = 0.125 * "ft";
255        col.trays.Q = 0 * "kW";
256        col.trays.beta = 0.6;
257        col.trays.alfa = 4;
258        col.trays.Ap = 3.94 * "ft^2";
259        col.reb.V = 2 * "m^3";
260        col.reb.Across = 1 * "m^2";
261       
262        INITIAL
263        # condenser
264        col.cond.OutletL.T = 260 *"K";
265        col.cond.Level = 1 * "m";
266        col.cond.OutletL.z([1:4]) = [0.65, 0.05, 0.01, 0.01];
267
268        # reboiler
269        col.reb.OutletL.T = 330 *"K";
270        col.reb.Level = 1 * "m";
271        col.reb.OutletL.z([1:4]) = [0.1, 0.7, 0.01, 0.01];
272
273        # column trays
274        col.trays.OutletL.T = [290:(330-290)/(col.NTrays-1):330] * "K";
275        col.trays.Level = 1.2 * col.trays.hw;
276        col.trays.OutletL.z([1:4]) = [0.5, 0.05, 0.01, 0.01];
277
278        OPTIONS
279        relativeAccuracy = 1e-3;
280        time = [0:0.01:1, 2:50];
281        #guessFile="Distillation_kettle_cond_Test.rlt";
282        #mode = "steady";       
283end
284
285FlowSheet Column_ctrl
286        PARAMETERS
287        PP      as CalcObject(Brief="Physical Properties",File="vrpp");
288        NComp   as Integer;
289       
290        Qcmin as heat_rate (Brief="Minimum Condenser Heat supplied");
291        Qcmax as heat_rate (Brief="Maximum Condenser Heat supplied");
292        Qrmin as heat_rate (Brief="Minimum Reboiler Heat supplied");
293        Qrmax as heat_rate (Brief="Maximum Reboiler Heat supplied");
294        Frmin as flow_mol (Brief="Minimum bottom flow rate");
295        Frmax as flow_mol (Brief="Maximum bottom flow rate");
296        Fcmin as flow_mol (Brief="Minimum reflux flow rate");
297        Fcmax as flow_mol (Brief="Maximum reflux flow rate");
298        Hmint as length (Brief="Minimum liquid level in top tank");
299        Hmaxt as length (Brief="Maximum liquid level in top tank");
300    Hminb as length (Brief="Minimum liquid level in reboiler");
301        Hmaxb as length (Brief="Maximum liquid level in reboiler");
302        Pmax as pressure (Brief="Maximum column pressure");
303        Pmin as pressure (Brief="Minimum column pressure");
304        Tmax as temperature (Brief="Maximum column temperature");
305        Tmin as temperature (Brief="Minimum column temperature");
306
307        VARIABLES
308        Qc as heat_rate (Brief="Heat rate removed from condenser");
309        Qr as heat_rate (Brief="Heat rate supplied to reboiler");
310        Had_top as Real (Brief="Dimensionless condenser level");
311        Had_bot as Real (Brief="Dimensionless reboiler level");
312        Pad as Real (Brief="Dimensionless pressure");
313        Tad as Real (Brief="Dimensionless temperature");
314        RR      as positive (Brief="Reflux ratio");
315
316        SET
317        PP.Components = [ "isobutane", "n-pentane", "propylene",
318                "benzene", "isobutene" ];
319        PP.LiquidModel = "PR";
320        PP.VapourModel = "PR";
321        PP.Derivatives = 1;
322        NComp = PP.NumberOfComponents;
323
324        DEVICES
325        col as Distillation_kettle_cond;
326        feed as streamTP;
327        zero as stream;
328        TCcond as PIDIncr_Ideal_AW;
329        LCtop as PIDIncr_Ideal_AW;
330        LCbot as PIDIncr_Ideal_AW;
331        PC as PIDIncr_Ideal_AW;
332
333        CONNECTIONS
334        feed to col.trays(5).Inlet;
335        zero to col.reb.Inlet;
336        zero to col.trays([1:4]).Inlet;
337        zero to col.trays([6:col.NTrays]).Inlet;
338        Qc to col.cond.Q;
339        Qr to col.reb.Q;
340
341        EQUATIONS
342   "Temperature Controller"
343        TCcond.Parameters.tau = 0*"s"; 
344        TCcond.Parameters.tauSet = 0*"s";       
345        TCcond.Parameters.alpha = 0.3;
346        TCcond.Parameters.bias = 0.5;   
347        TCcond.Parameters.gamma = 1;
348        TCcond.Parameters.beta = 1;
349        TCcond.Options.action = 1;
350        TCcond.Options.clip = 1;
351        TCcond.Options.autoMan = 0;
352        TCcond.Parameters.intTime = 60*"s";
353        TCcond.Parameters.gain = 0.6;
354        TCcond.Parameters.derivTime = 1*"s";
355        TCcond.Ports.setPoint = ((15+273.15) * "K" - Tmin)/(Tmax-Tmin);
356        TCcond.Ports.input = Tad;
357        Tad = (col.cond.OutletL.T-Tmin)/(Tmax-Tmin);
358        Qc = Qcmin+(Qcmax-Qcmin)*TCcond.Ports.output;   
359
360        "Pressure Controller"
361        PC.Parameters.tau = 0*"s";     
362        PC.Parameters.tauSet = 0*"s";   
363        PC.Parameters.alpha = 0.3;
364        PC.Parameters.bias = 0;
365        PC.Parameters.gamma = 1;
366        PC.Parameters.beta = 1;
367        PC.Options.action = -1;
368        PC.Options.clip = 1;
369        PC.Options.autoMan = 0;
370        PC.Parameters.intTime = 50*"s";
371        PC.Parameters.gain = 0.5;
372        PC.Parameters.derivTime = 1*"s";
373        PC.Ports.setPoint = (2.0*"bar"-Pmin)/(Pmax-Pmin);
374        PC.Ports.input = Pad;
375        Pad = (col.cond.OutletV.P-Pmin)/(Pmax-Pmin);
376        col.cond.OutletV.F = (Fcmin+(Fcmax-Fcmin)*PC.Ports.output);     
377       
378        "Ttop Level Controller"
379        LCtop.Parameters.tau = 0*"s";   
380        LCtop.Parameters.tauSet = 0*"s";       
381        LCtop.Parameters.alpha = 0.3;
382        LCtop.Parameters.bias = 0.5;   
383        LCtop.Parameters.gamma = 1;
384        LCtop.Parameters.beta = 1;
385        LCtop.Options.action = -1;
386        LCtop.Options.clip = 1;
387        LCtop.Options.autoMan = 0;
388        LCtop.Parameters.intTime = 10*"s";
389        LCtop.Parameters.gain = 1;
390        LCtop.Parameters.derivTime = 1*"s";
391        LCtop.Ports.setPoint = (1.0 * "m" - Hmint)/(Hmaxt-Hmint);
392        LCtop.Ports.input = Had_top;
393        Had_top = (col.cond.Level-Hmint)/(Hmaxt-Hmint);
394        col.sptop.Outlet1.F = Fcmin + (Fcmax-Fcmin) * LCtop.Ports.output;
395
396        "Tbottom Level Controller"
397        LCbot.Parameters.tau = 0*"s";   
398        LCbot.Parameters.tauSet = 0*"s";       
399        LCbot.Parameters.alpha = 0.3;
400        LCbot.Parameters.bias = 0.5;   
401        LCbot.Parameters.gamma = 1;
402        LCbot.Parameters.beta = 1;
403        LCbot.Options.action = -1;
404        LCbot.Options.clip = 1;
405        LCbot.Options.autoMan = 0;
406        LCbot.Parameters.intTime = 100*"s";
407        LCbot.Parameters.gain = 1;
408        LCbot.Parameters.derivTime = 1*"s";
409        LCbot.Ports.setPoint = (1.0 * "m" - Hminb)/(Hmaxb-Hminb);
410        LCbot.Ports.input = Had_bot;
411        Had_bot = (col.reb.Level-Hminb)/(Hmaxb-Hminb);
412        col.reb.OutletL.F = Frmin + (Frmax-Frmin) * LCbot.Ports.output;
413
414        RR * (col.cond.OutletV.F + col.sptop.Outlet1.F) = col.sptop.Outlet2.F;
415       
416        if time < 1 * "h" then
417                col.sptop.Outlet2.F = 70 * "kmol/h"; # reflux
418        else
419                col.sptop.Outlet2.F = 85 * "kmol/h"; # reflux
420        end
421
422        SPECIFY
423        feed.F = 113.4 * "kmol/h";
424        feed.T = 291 * "K";
425        feed.P = 168.3 * "kPa";
426        feed.z = 1/NComp;
427       
428        zero.F = 0 * "kmol/h";
429        zero.T = 300 * "K";
430        zero.P = 1 * "atm";
431        zero.z = 1/NComp;
432        zero.v = 0;
433        zero.h = 0 * "J/mol";
434       
435        Qr = 3e6 * "kJ/h";
436        col.pump1.dP = 16 * "kPa";
437        col.trays.Emv = 1;
438
439        SET
440        col.NTrays = 8;
441        col.cond.V = 2 * "m^3";
442        col.cond.Across = 1 * "m^2";
443        col.trays.V = 4 * "ft^3";
444        col.trays.Ah = 0.394 * "ft^2";
445        col.trays.lw = 20.94 * "in";
446        col.trays.hw = 0.125 * "ft";
447        col.trays.Q = 0 * "kW";
448        col.trays.beta = 0.6;
449        col.trays.alfa = 4;
450        col.trays.Ap = 3.94 * "ft^2";
451        col.reb.V = 2 * "m^3";
452        col.reb.Across = 1 * "m^2";
453
454        Qrmax = 5e6 * "kJ/h";
455        Qrmin = 1e6 * "kJ/h";
456        Frmin = 0 * "kmol/h";
457        Frmax = 60 * "kmol/h";
458        Fcmin = 0 * "kmol/h";
459        Fcmax = 120 * "kmol/h";
460        Hmint = 0 * "m";
461        Hmaxt = 2 * "m";
462        Hminb = 0 * "m";
463        Hmaxb = 2 * "m";
464        Pmin = 0.5 * "bar";
465        Pmax = 4 * "bar";
466        Qcmax = -5e5 * "kJ/h";
467        Qcmin = -5e6 * "kJ/h";
468        Tmax = (30+273.15) * "K";
469        Tmin = (-20+273.15) * "K";
470       
471        INITIAL
472        # condenser
473        col.cond.OutletL.T = 260 *"K";
474        col.cond.Level = 1 * "m";
475        col.cond.OutletL.z([1:4]) = [0.2, 0.2, 0.4, 0.05];
476
477        # reboiler
478        col.reb.OutletL.T = 350 *"K";
479        col.reb.Level = 1 * "m";
480        col.reb.OutletL.z([1:4]) = [0.1, 0.4, 0.1, 0.3];
481       
482        # column trays
483        col.trays.OutletL.T = [290:(330-290)/(col.NTrays-1):330] * "K";
484        col.trays.Level = 1.2 * col.trays.hw;
485        col.trays.OutletL.z([1:4]) = [0.15, 0.3, 0.25, 0.2];
486
487        OPTIONS
488        relativeAccuracy = 1e-3;
489        time = [0:0.01:0.1, 0.11:0.01:2]*"h";
490        #initialFile = "Column_ctrl.rlt";
491        #guessFile = "Column_ctrl.rlt";
492        #mode = "steady";
493end
494
Note: See TracBrowser for help on using the repository browser.