source: branches/gui/eml/stage_separators/tray.mso @ 919

Last change on this file since 919 was 901, checked in by mamuller, 14 years ago

added composition indicator on the PackedColumn? model
File size: 28.3 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* Author: Paula B. Staudt
17* $Id: tray.mso 522 2008-05-21 23:21:12Z arge $
18*--------------------------------------------------------------------*#
19
20using "streams";
21
22Model tray
23
24ATTRIBUTES
25        Pallete         = false;
26        Icon            = "icon/Tray";
27        Brief           = "Complete model of a column tray.";
28        Info            =
29"== Assumptions ==
30* both phases (liquid and vapour) exists all the time;
31* thermodymanic equilibrium with Murphree plate efficiency;
32* no entrainment of liquid or vapour phase;
33* no weeping;
34* the dymanics in the downcomer are neglected.
35
36== Options ==
37You can choose the equation for the liquid outlet flow and the vapour
38inlet flow calculation through the VapourFlowModel and LiquidFlowModel
39switchers.
40
41== References ==
42* ELGUE, S.; PRAT, L.; CABASSUD, M.; LANN, J. L.; CéZERAC, J. Dynamic models for start-up operations of batch distillation columns with experimental validation. Computers and Chemical Engineering, v. 28, p. 2735-2747, 2004.
43* FEEHERY, W. F. Dynamic Optimization with Path Constraints. Tese (Doutorado) - Massachusetts Institute of Technology, June 1998.
44* KLINGBERG, A. Modeling and Optimization of Batch Distillation. Dissertação (Mestrado) - Department of Automatic Control, Lund Institute of Technology, Lund, Sweden, fev. 2000.
45* OLSEN, I.; ENDRESTOL, G. O.; SIRA, T. A rigorous and efficient distillation column model for engineering and training simulators. Computers and Chemical Engineering,v. 21, n. Suppl, p. S193-S198, 1997.
46* REEPMEYER, F.; REPKE, J.-U.; WOZNY, G. Analysis of the start-up process for reactive distillation. Chemical Engineering Technology, v. 26, p. 81-86, 2003.
47* ROFFEL, B.; BETLEM, B.; RUIJTER, J. de. First principles dynamic modeling and multivariable control of a cryogenic distillation column process. Computers and Chemical Engineering, v. 24, p. 111-123, 2000.
48* WANG, L.; LI, P.; WOZNY, G.; WANG, S. A start-up model for simulation of batch distillation starting from a cold state. Computers and Chemical Engineering, v. 27, p.1485-1497, 2003.
49";     
50       
51PARAMETERS
52outer PP                as Plugin               (Brief = "External Physical Properties", Type="PP");
53outer NComp     as Integer              (Brief="Number of components");
54
55        Mw(NComp)                                                       as molweight            (Brief="Component Mol Weight",Hidden=true);
56        Gconst                                                          as acceleration         (Brief="Gravity Acceleration",Default=9.81,Hidden=true);
57        zero_flow                                                       as flow_mol             (Brief = "Stream Flow closed",Default = 0, Hidden=true);
58        low_flow                                                        as flow_mol             (Brief = "Low stream Flow",Default = 1E-6, Hidden=true);
59        Pi                                                                      as constant             (Brief="Pi Number",Default=3.14159265, Symbol = "\pi",Hidden=true);
60       
61        TrayDiameter_                   as length               (Brief="Tray Diameter",Default=1.600);
62        TraySpacing_                    as length               (Brief="Tray Spacing",Default=0.600);
63        Fraction_HoleArea_              as fraction     (Brief="Fraction of the active area that is occupied by the holes with respect to the total tray area",Default=0.10);
64        Fraction_DowncomerArea_ as fraction     (Brief="Fraction of the downcomer area with respect to the total tray area",Default=0.20);
65        WeirLength_                             as length               (Brief="Weir length", Default = 1);
66        WeirHeight_                     as length               (Brief="Weir height", Default= 0.05);
67        TrayLiquidPasses_               as positive     (Brief="Number of liquid passes in the tray", Lower = 1,Default=1);
68        HeatSupply_                     as heat_rate    (Brief="Rate of heat supply",Default = 0);
69        AerationFraction_               as Real                 (Brief="Aeration fraction", Default = 1);
70        DryPdropCoeff_                  as Real                 (Brief="Dry pressure drop coefficient", Default= 0.60);
71        MurphreeEff_                    as Real                 (Brief="Murphree efficiency for All Trays",Lower=0.01,Upper=1);
72
73        PlateArea_                              as area                 (Brief="Plate area = Atray - Adowncomer",Protected=true);
74        TrayVolume_                             as volume               (Brief="Total Volume of the tray",Protected=true);
75        HolesArea_                              as area                 (Brief="Total holes area",Protected=true);
76       
77        FeeheryCoeff    as Real         (Brief="Feeherys correlation coefficient", Unit='1/m^4', Default=1,Hidden=true);
78        ElgueCoeff              as Real         (Brief="Elgues correlation coefficient", Unit='kg/m/mol^2', Default=1,Hidden=true);
79        OlsenCoeff              as Real         (Brief="Olsens correlation coefficient", Default=1,Hidden=true);
80       
81        VapourFlow      as Switcher     (Brief="Flag for Vapour Flow condition",Valid = ["on", "off"], Default = "off",Hidden=true);
82        LiquidFlow      as Switcher     (Brief="Flag for Liquid Flow condition",Valid = ["on", "off"], Default = "off",Hidden=true);
83       
84VARIABLES
85
86        Inlet                           as stream                       (Brief="Feed stream", Hidden=true, PosX=0, PosY=0.4932, Symbol="_{in}");
87        LiquidSideStream        as liquid_stream        (Brief="liquid Sidestream", Hidden=true, Symbol="_{outL}");
88        VapourSideStream        as vapour_stream        (Brief="vapour Sidestream", Hidden=true, Symbol="_{outV}");
89
90in      InletLiquid     as stream                       (Brief="Inlet liquid stream", PosX=0.5195, PosY=0, Symbol="_{inL}");
91in      InletVapour     as stream                       (Brief="Inlet vapour stream", PosX=0.4994, PosY=1, Symbol="_{inV}");
92out     OutletLiquid    as liquid_stream        (Brief="Outlet liquid stream", PosX=0.8277, PosY=1, Symbol="_{outL}");
93out     OutletVapour    as vapour_stream        (Brief="Outlet vapour stream", PosX=0.8043, PosY=0, Symbol="_{outV}");
94
95        LFlowModel      as positive     (Brief="Flag for Liquid Flow Model",Lower = 1, Default = 1 , Hidden=true);
96        VFlowModel      as positive     (Brief="Flag for Vapour Flow Model",Lower = 1, Default = 1 , Hidden=true);
97
98        M(NComp)                        as mol                  (Brief="Molar Holdup in the tray");
99        ML                                      as mol                  (Brief="Molar liquid holdup");
100        MV                                      as mol                  (Brief="Molar vapour holdup");
101        E                                       as energy               (Brief="Total Energy Holdup on tray");
102        vL                                      as volume_mol   (Brief="Liquid Molar Volume");
103        vV                                      as volume_mol   (Brief="Vapour Molar volume");
104        Level                           as length               (Brief="Height of clear liquid on plate");
105        yideal(NComp)           as fraction;
106        rhoL                            as dens_mass    (Brief="Mass Density of liquid phase");
107        rhoV                            as dens_mass    (Brief="Mass Density of vapour phase");
108
109SET
110
111        Mw = PP.MolecularWeight();
112        zero_flow = 0 * 'kmol/h';
113        low_flow = 1E-6 * 'kmol/h';
114
115        PlateArea_ = 0.25*Pi*(TrayDiameter_^2)*(1-Fraction_DowncomerArea_);
116        TrayVolume_ = 0.25*Pi*(TrayDiameter_^2)*TraySpacing_;
117        HolesArea_ = 0.25*Pi*(TrayDiameter_^2)*Fraction_HoleArea_;
118
119EQUATIONS
120
121# LiquidFlow and VapourFlow equations need to be linerized to avoid indetermination !
122switch LiquidFlow
123
124        case "on":
125                       
126                        if LFlowModel equal 1 then
127                                "Francis Equation"
128                                OutletLiquid.F*vL = 1.84*'1/s'*WeirLength_*((Level-(AerationFraction_*WeirHeight_))/(AerationFraction_))^2;
129                       
130                        else if LFlowModel equal 2 then
131                                "Wang_Fl"
132                                OutletLiquid.F*vL = 1.84*'m^0.5/s'*WeirLength_*((Level-(AerationFraction_*WeirHeight_))/(AerationFraction_))^1.5;
133                       
134                        else if LFlowModel equal 3 then
135                                "Olsen"
136                                OutletLiquid.F / 'mol/s'= WeirLength_*TrayLiquidPasses_*rhoL/sum(Mw*OutletVapour.z)/(0.665*OlsenCoeff)^1.5 * ((ML*sum(Mw*OutletLiquid.z)/rhoL/PlateArea_)-WeirHeight_)^1.5 * 'm^0.5/mol';
137                       
138                        else if LFlowModel equal 4 then 
139                                "Feehery_Fl"
140                                OutletLiquid.F = WeirLength_*rhoL/sum(Mw*OutletLiquid.z) * ((Level-WeirHeight_)/750/'mm')^1.5 * 'm^2/s';
141                       
142                        else   
143                                "Roffel_Fl"
144                                OutletLiquid.F = 2/3*rhoL/sum(Mw*OutletLiquid.z)*WeirLength_*(ML*sum(Mw*OutletLiquid.z)/(PlateArea_*1.3)/rhoL)^1.5*sqrt(2*Gconst/
145                                                        (2*(1 - 0.3593/'Pa^0.0888545'*abs(OutletVapour.F*sum(Mw*OutletVapour.z)/(PlateArea_*1.3)/sqrt(rhoV))^0.177709)-1)); #/'(kg/m)^0.0888545/s^0.177709';
146end
147end
148end
149end
150               
151                when Level < (AerationFraction_ *WeirHeight_) switchto "off";
152               
153                case "off":
154               
155                "Low level"
156                OutletLiquid.F = zero_flow;
157               
158                when Level > (AerationFraction_ * WeirHeight_) switchto "on";
159               
160end
161       
162switch VapourFlow
163       
164        case "on":
165                       
166                        if VFlowModel equal 1 then
167                                "Reepmeyer"
168                                InletVapour.F*vV = sqrt((InletVapour.P - OutletVapour.P)/(rhoV*DryPdropCoeff_))*HolesArea_;
169                       
170                        else if VFlowModel equal 2 then
171                                "Feehery_Fv"
172                                InletVapour.F = rhoV/PlateArea_/FeeheryCoeff/sum(Mw*OutletVapour.z) * sqrt(((InletVapour.P - OutletVapour.P)-(rhoV*Gconst*ML*vL/PlateArea_))/rhoV);
173                       
174                        else if VFlowModel equal 3 then
175                                "Roffel_Fv"
176                                InletVapour.F^1.08 * 0.0013 * 'kg/m/mol^1.08/s^0.92*1e5' = (InletVapour.P - OutletVapour.P)*1e5 - (AerationFraction_*sum(M*Mw)/(PlateArea_*1.3)*Gconst*1e5) * (rhoV*HolesArea_/sum(Mw*OutletVapour.z))^1.08 * 'm^1.08/mol^1.08';
177                       
178                        else if VFlowModel equal 4  then
179                                "Klingberg"
180                                InletVapour.F * vV = PlateArea_ * sqrt(((InletVapour.P - OutletVapour.P)-rhoL*Gconst*Level)/rhoV);
181                       
182                        else if VFlowModel equal 5 then
183                                "Wang_Fv"
184                                InletVapour.F * vV = PlateArea_ * sqrt(((InletVapour.P - OutletVapour.P)-rhoL*Gconst*Level)/rhoV*DryPdropCoeff_);
185                               
186                        else
187                                "Elgue"
188                                InletVapour.F  = sqrt((InletVapour.P - OutletVapour.P)/ElgueCoeff);
189end
190end
191end
192end
193end
194       
195        when InletVapour.F < low_flow switchto "off";
196
197        case "off":
198                InletVapour.F = zero_flow;
199               
200        when InletVapour.P > OutletVapour.P switchto "on";
201
202end
203
204"Murphree Efficiency"
205        OutletVapour.z =  MurphreeEff_ * (yideal - InletVapour.z) + InletVapour.z;
206       
207"Energy Balance"
208        diff(E) = ( Inlet.F*Inlet.h + InletLiquid.F*InletLiquid.h + InletVapour.F*InletVapour.h- OutletLiquid.F*OutletLiquid.h - OutletVapour.F*OutletVapour.h
209        - VapourSideStream.F*VapourSideStream.h - LiquidSideStream.F*LiquidSideStream.h + HeatSupply_ );
210
211"Energy Holdup"
212        E = ML*OutletLiquid.h + MV*OutletVapour.h - OutletLiquid.P*TrayVolume_;
213
214"Geometry Constraint"
215        TrayVolume_ = ML* vL + MV*vV;
216
217"Level of clear liquid over the weir"
218        Level = ML*vL/PlateArea_;
219
220"Component Molar Balance"
221        diff(M)=Inlet.F*Inlet.z + InletLiquid.F*InletLiquid.z + InletVapour.F*InletVapour.z- OutletLiquid.F*OutletLiquid.z - OutletVapour.F*OutletVapour.z-
222        LiquidSideStream.F*LiquidSideStream.z-VapourSideStream.F*VapourSideStream.z;
223
224"Molar Holdup"
225        M = ML*OutletLiquid.z + MV*OutletVapour.z;
226
227"Mol fraction normalisation"
228        sum(OutletLiquid.z)= 1.0;
229
230"Mol fraction constraint"
231        sum(OutletLiquid.z)= sum(OutletVapour.z);
232
233"Liquid Volume"
234        vL = PP.LiquidVolume(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z);
235
236"Vapour Volume"
237        vV = PP.VapourVolume(OutletVapour.T, OutletVapour.P, OutletVapour.z);
238
239"Liquid Density"
240        rhoL = PP.LiquidDensity(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z);
241       
242"Vapour Density"
243        rhoV = PP.VapourDensity(InletVapour.T, InletVapour.P, InletVapour.z);
244
245"Chemical Equilibrium"
246        PP.LiquidFugacityCoefficient(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z)*OutletLiquid.z = PP.VapourFugacityCoefficient(OutletVapour.T, OutletVapour.P, yideal)*yideal;
247
248"Thermal Equilibrium"
249        OutletVapour.T = OutletLiquid.T;
250
251"Mechanical Equilibrium"
252        OutletVapour.P = OutletLiquid.P;
253       
254"Thermal Equilibrium Vapour Side Stream"
255        OutletVapour.T = VapourSideStream.T;
256
257"Thermal Equilibrium Liquid Side Stream"
258        OutletLiquid.T = LiquidSideStream.T;
259
260"Mechanical Equilibrium Vapour Side Stream"
261        OutletVapour.P= VapourSideStream.P;
262
263"Mechanical Equilibrium Liquid Side Stream"
264        OutletLiquid.P = LiquidSideStream.P;
265
266"Composition Liquid Side Stream"
267        OutletLiquid.z= LiquidSideStream.z;
268       
269"Composition Vapour Side Stream"
270        OutletVapour.z= VapourSideStream.z;
271
272end
273
274#*-------------------------------------------------------------------
275* Model of a tray with reaction
276*-------------------------------------------------------------------*#
277Model trayReac
278        ATTRIBUTES
279        Pallete         = false;
280        Icon            = "icon/Tray";
281        Brief           = "Model of a tray with reaction.";
282        Info            =
283"== Assumptions ==
284* both phases (liquid and vapour) exists all the time;
285* thermodymanic equilibrium with Murphree plate efficiency;
286* no entrainment of liquid or vapour phase;
287* no weeping;
288* the dymanics in the downcomer are neglected.
289       
290== Specify ==
291* the Feed stream;
292* the Liquid inlet stream;
293* the Vapour inlet stream;
294* the Vapour outlet flow (OutletVapour.F);
295* the reaction related variables.
296       
297== Initial ==
298* the plate temperature (OutletLiquid.T)
299* the liquid height (Level) OR the liquid flow OutletLiquid.F
300* (NoComps - 1) OutletLiquid compositions
301";
302
303PARAMETERS
304
305        outer PP                        as Plugin(Type="PP");
306        outer NComp     as Integer;
307       
308VARIABLES
309
310        Inlet                                           as stream                               (Brief="Feed stream", PosX=0, PosY=0.4932, Symbol="_{in}");
311        LiquidSideStream        as liquid_stream        (Brief="liquid Sidestream", Hidden=true, Symbol="_{outL}");
312        VapourSideStream as vapour_stream       (Brief="vapour Sidestream", Hidden=true, Symbol="_{outV}");
313       
314
315in              InletLiquid             as      stream                          (Brief="Inlet liquid stream", PosX=0.5195, PosY=0, Symbol="_{inL}");
316in              InletVapour             as      stream                          (Brief="Inlet vapour stream", PosX=0.4994, PosY=1, Symbol="_{inV}");
317out     OutletLiquid    as      liquid_stream           (Brief="Outlet liquid stream", PosX=0.8277, PosY=1, Symbol="_{outL}");
318out     OutletVapour    as      vapour_stream   (Brief="Outlet vapour stream", PosX=0.8043, PosY=0, Symbol="_{outV}");
319
320        yideal(NComp)   as fraction;
321
322        M(NComp)                as mol                          (Brief="Molar Holdup in the tray");
323        ML                                              as mol                          (Brief="Molar liquid holdup");
324        MV                                              as mol                          (Brief="Molar vapour holdup");
325        E                                               as energy                       (Brief="Total Energy Holdup on tray");
326        vL                                              as volume_mol   (Brief="Liquid Molar Volume");
327        vV                                              as volume_mol   (Brief="Vapour Molar volume");
328        Level                                   as length                       (Brief="Height of clear liquid on plate");
329        Vol                                             as volume;
330       
331        rhoL                    as dens_mass;
332        rhoV                    as dens_mass;
333        r3                              as reaction_mol         (Brief = "Reaction resulting ethyl acetate", DisplayUnit = 'mol/l/s');
334        C(NComp)        as conc_mol             (Brief = "Molar concentration", Lower = -1);
335       
336EQUATIONS
337
338"Molar Concentration"
339        OutletLiquid.z = vL * C;
340       
341"Reaction"
342        r3 = exp(-7150*'K'/OutletLiquid.T)*(4.85e4*C(1)*C(2) - 1.23e4*C(3)*C(4))*'l/mol/s';
343       
344"Molar Holdup"
345        M = ML*OutletLiquid.z + MV*OutletVapour.z;
346
347"Thermal Equilibrium Vapour Side Stream"
348        OutletVapour.T = VapourSideStream.T;
349
350"Thermal Equilibrium Liquid Side Stream"
351        OutletLiquid.T = LiquidSideStream.T;
352
353"Mechanical Equilibrium Vapour Side Stream"
354        OutletVapour.P= VapourSideStream.P;
355
356"Mechanical Equilibrium Liquid Side Stream"
357        OutletLiquid.P = LiquidSideStream.P;
358
359"Composition Liquid Side Stream"
360        OutletLiquid.z= LiquidSideStream.z;
361       
362"Composition Vapour Side Stream"
363        OutletVapour.z= VapourSideStream.z;
364
365"Mol fraction normalisation"
366        sum(OutletLiquid.z)= 1.0;
367
368"Liquid Volume"
369        vL = PP.LiquidVolume(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z);
370
371"Vapour Volume"
372        vV = PP.VapourVolume(OutletVapour.T, OutletVapour.P, OutletVapour.z);
373
374"Thermal Equilibrium"
375        OutletVapour.T = OutletLiquid.T;
376
377"Mechanical Equilibrium"
378        OutletVapour.P = OutletLiquid.P;
379
380        Vol = ML*vL;
381       
382"Liquid Density"
383        rhoL = PP.LiquidDensity(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z);
384
385"Vapour Density"
386        rhoV = PP.VapourDensity(InletVapour.T, InletVapour.P, InletVapour.z);
387
388"Chemical Equilibrium"
389        PP.LiquidFugacityCoefficient(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z)*OutletLiquid.z =   PP.VapourFugacityCoefficient(OutletVapour.T, OutletVapour.P, yideal)*yideal;
390
391        sum(OutletLiquid.z)= sum(OutletVapour.z);
392
393end
394
395#*-------------------------------------
396* Model of a packed column stage
397-------------------------------------*#
398Model packedStage
399
400ATTRIBUTES
401        Pallete         = false;
402        Brief           = "Complete model of a packed column stage.";
403        Info            =
404"== Specify ==
405* the Feed stream
406* the Liquid inlet stream
407* the Vapour inlet stream
408* the stage pressure drop (deltaP)
409       
410== Initial ==
411* the plate temperature (OutletLiquid.T)
412* the liquid molar holdup ML
413* (NoComps - 1) OutletLiquid compositions
414";
415
416PARAMETERS
417
418outer PP        as Plugin       (Brief = "External Physical Properties", Type="PP");
419outer NComp as Integer  (Brief = "Number Of Components");
420
421        LiquidResistanceCoeff   as positive     (Brief="Resistance coefficient on the liquid load", Default=1,Hidden=true);
422        AreaPerPackingVolume    as Real                 (Brief="surface area per packing volume", Unit='m^2/m^3',Hidden=true);
423        ColumnInternalDiameter  as length               (Brief="Column diameter",Hidden=true);
424        PackingVoidFraction             as Real                 (Brief="Void fraction of packing, (m^3 void space/m^3 packed bed)",Hidden=true);
425        HeightOfPacking                 as length               (Brief="Height of packing",Hidden=true);
426        Number_Stages                   as Integer              (Brief="Number of Stages", Default=3,Hidden=true);
427        HeatOnStage                     as heat_rate    (Brief="Rate of heat supply",Hidden=true);
428
429        HETP                    as length               (Brief="The Height Equivalent to a Theoretical Plate",Hidden=true);
430        ColumnArea              as area                 (Brief="Column Sectional Cross Area",Hidden=true);
431        V                               as volume               (Brief="Total Volume of the tray",Hidden=true);
432        Pi                              as constant     (Brief="Pi Number",Default=3.14159265, Symbol = "\pi",Hidden=true);
433        Gconst                  as acceleration (Brief="Gravity Acceleration",Default=9.81,Hidden=true);
434       
435        low_flow                as flow_mol     (Brief ="Low Flow",Default = 1E-6, Hidden=true);
436        low_pressure    as pressure     (Brief ="Low Pressure",Default = 1E-6, Hidden=true);
437        zero_flow               as flow_mol     (Brief ="No Flow",Default = 0, Hidden=true);
438
439        Mw(NComp)       as molweight    (Brief = "Component Mol Weight",Hidden=true);
440        VapourFlow  as Switcher         (Brief = "Vapour Flow", Valid = ["on", "off"], Default = "on",Hidden=true);
441
442SET
443        Mw = PP.MolecularWeight();
444       
445        ColumnArea      = 0.25*Pi*ColumnInternalDiameter^2;
446        HETP            = HeightOfPacking/Number_Stages;       
447        V                       = HETP * ColumnArea;
448       
449        low_pressure = 1E-4 * 'atm';
450        low_flow = 1E-6 * 'kmol/h';
451        zero_flow = 0 * 'kmol/h';
452
453VARIABLES
454
455        Inlet                   as stream                       (Brief="Feed stream", Symbol="_{in}",Protected=true);
456in      InletLiquid     as stream                       (Brief="Inlet liquid stream",  Symbol="_{inL}",Protected=true);
457in      InletVapour     as stream                       (Brief="Inlet vapour stream",  Symbol="_{inV}",Protected=true);
458out     OutletLiquid    as liquid_stream        (Brief="Outlet liquid stream", Symbol="_{outL}",Protected=true);
459out     OutletVapour    as vapour_stream        (Brief="Outlet vapour stream", Symbol="_{outV}",Protected=true);
460
461        M(NComp)        as mol                  (Brief="Molar Holdup in the tray", Default=0.01, Lower=-0.000001, Upper=100,Protected=true);
462        ML                      as mol                  (Brief="Molar liquid holdup", Default=0.01, Lower=0, Upper=100,Protected=true);
463        MV                      as mol                  (Brief="Molar vapour holdup", Default=0.01, Lower=0, Upper=100,Protected=true);
464        E                       as energy               (Brief="Total Energy Holdup on tray", Default=-500,Protected=true);
465        vL                      as volume_mol   (Brief="Liquid Molar Volume",Protected=true);
466        vV                      as volume_mol   (Brief="Vapour Molar volume",Protected=true);
467        miL             as viscosity    (Brief="Liquid dynamic viscosity", DisplayUnit='kg/m/s',Protected=true);
468        rhoL            as dens_mass    (Brief="Liquid mass density",Protected=true);
469        rhoV            as dens_mass    (Brief="Vapour mass density",Protected=true);
470        uL                      as velocity     (Brief="volume flow rate of liquid, m^3/m^2/s", Lower=-10, Upper=1000,Protected=true);
471        uV                      as velocity     (Brief="volume flow rate of vapor, m^3/m^2/s", Lower=-10, Upper=1000,Protected=true);
472        Al                      as area                 (Brief="Area occupied by the liquid", Default=0.001, Upper=10,Protected=true);
473        hl                      as positive     (Brief="Column holdup", Unit='m^3/m^3', Default=0.01,Upper=10,Protected=true);
474        deltaP          as pressure     (Brief="Stage Pressure drop",Protected=true);
475
476EQUATIONS
477
478switch VapourFlow
479       
480        case "on":
481"Pressure drop and Vapor flow, Billet (4-58)"
482        deltaP/HETP  = LiquidResistanceCoeff *( 0.5*AreaPerPackingVolume + 2/ColumnInternalDiameter) * 1/((PackingVoidFraction-hl)^3) * (uV^2) *rhoV;
483       
484        when InletVapour.F < low_flow switchto "off";
485       
486        case "off":
487"Vapour Flow"
488        InletVapour.F = zero_flow;
489       
490        when deltaP > low_pressure switchto "on";
491
492end
493
494"Energy Balance"
495        diff(E) = (Inlet.F*Inlet.h + InletLiquid.F*InletLiquid.h + InletVapour.F*InletVapour.h- OutletLiquid.F*OutletLiquid.h
496        - OutletVapour.F*OutletVapour.h + HeatOnStage );
497
498"Energy Holdup"
499        E = ML*OutletLiquid.h + MV*OutletVapour.h - OutletLiquid.P*V;
500
501"Geometry Constraint"
502        V*PackingVoidFraction= ML*vL + MV*vV;
503
504"Volume flow rate of vapor, m^3/m^2/s"
505        uV * (V*PackingVoidFraction/HETP - Al) = InletVapour.F * vV;
506
507"Liquid holdup"
508        hl*V*PackingVoidFraction = ML*vL;
509
510"Liquid velocity as a function of liquid holdup, Billet (4-27)"
511        hl^3 = (12/Gconst) * AreaPerPackingVolume^2 * (miL/rhoL) * uL;
512
513"Area occupied by the liquid"
514        Al = ML*vL/HETP;
515
516"Component Molar Balance"
517        diff(M)=Inlet.F*Inlet.z + InletLiquid.F*InletLiquid.z + InletVapour.F*InletVapour.z- OutletLiquid.F*OutletLiquid.z - OutletVapour.F*OutletVapour.z;
518
519"Molar Holdup"
520        M = ML*OutletLiquid.z + MV*OutletVapour.z;
521
522"Mol Fraction Normalisation"
523        sum(OutletLiquid.z)= 1.0;
524
525"Mol Fraction Constraint"
526        sum(OutletLiquid.z)= sum(OutletVapour.z);
527       
528"Liquid Volume"
529        vL = PP.LiquidVolume(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z);
530
531"Vapour Volume"
532        vV = PP.VapourVolume(OutletVapour.T, OutletVapour.P, OutletVapour.z);
533       
534"Chemical Equilibrium"
535        PP.LiquidFugacityCoefficient(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z)*OutletLiquid.z =   PP.VapourFugacityCoefficient(OutletVapour.T, OutletVapour.P, OutletVapour.z)*OutletVapour.z;
536       
537"Thermal Equilibrium"
538        OutletVapour.T = OutletLiquid.T;
539       
540"Mechanical Equilibrium"
541        OutletLiquid.P = OutletVapour.P;
542       
543"Liquid Density"
544        rhoL = PP.LiquidDensity(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z);
545
546"Vapour Density"
547        rhoV = PP.VapourDensity(InletVapour.T, InletVapour.P, InletVapour.z);
548
549"Liquid viscosity"
550        miL = PP.LiquidViscosity(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z);
551
552"Volume flow rate of liquid, m^3/m^2/s"
553        uL * Al = OutletLiquid.F * vL;
554       
555"Pressure Drop"
556        deltaP = InletVapour.P - OutletVapour.P;       
557
558end
559
560#*-------------------------------------
561* Nonequilibrium Model
562-------------------------------------*
563Model interfaceTeste
564       
565        ATTRIBUTES
566        Pallete         = false;
567        Icon            = "icon/Tray";
568        Brief           = "Descrition of variables of the equilibrium interface.";
569        Info            =
570"This model contains only the variables of the equilibrium interface.";
571
572        PARAMETERS
573outer PP as Plugin(Brief = "External Physical Properties", Type="PP");
574outer NComp as Integer;
575outer NC1 as Integer;
576       
577        VARIABLES
578        NL(NComp) as flow_mol_delta     (Brief = "Stream Molar Rate on Liquid Phase");
579        NV(NComp) as flow_mol_delta     (Brief = "Stream Molar Rate on Vapour Phase");
580        T as temperature                (Brief = "Stream Temperature");
581        P as pressure                   (Brief = "Stream Pressure");
582        x(NComp) as fraction    (Brief = "Stream Molar Fraction on Liquid Phase");
583        y(NComp) as fraction    (Brief = "Stream Molar Fraction on Vapour Phase");
584        a as area                           (Brief = "Interface Area");
585        htL as heat_trans_coeff (Brief = "Heat Transference Coefficient on Liquid Phase");
586        htV as heat_trans_coeff (Brief = "Heat Transference Coefficient on Vapour Phase");     
587        E_liq as heat_rate      (Brief = "Liquid Energy Rate at interface");
588    E_vap as heat_rate      (Brief = "Vapour Energy Rate at interface");       
589        hL as enth_mol          (Brief = "Liquid Molar Enthalpy");
590        hV as enth_mol          (Brief = "Vapour Molar Enthalpy");
591        kL(NC1,NC1) as velocity (Brief = "Mass Transfer Coefficients");
592        kV(NC1,NC1) as velocity (Brief = "Mass Transfer Coefficients");
593       
594        EQUATIONS
595        "Liquid Enthalpy"
596        hL = PP.LiquidEnthalpy(T, P, x);
597       
598        "Vapour Enthalpy"
599        hV = PP.VapourEnthalpy(T, P, y);
600
601end
602
603Model trayRateBasicTeste
604        ATTRIBUTES
605        Pallete         = false;
606        Icon            = "icon/Tray";
607        Brief           = "Basic equations of a tray rate column model.";
608        Info            =
609"This model contains only the main equations of a column tray nonequilibrium model without
610the hidraulic equations.
611       
612== Assumptions ==
613* both phases (liquid and vapour) exists all the time;
614* no entrainment of liquid or vapour phase;
615* no weeping;
616* the dymanics in the downcomer are neglected.
617";
618       
619        PARAMETERS
620outer PP as Plugin(Brief = "External Physical Properties", Type="PP");
621outer NComp as Integer;
622    NC1 as Integer;
623        V as volume(Brief="Total Volume of the tray");
624        Q as heat_rate (Brief="Rate of heat supply");
625        Ap as area (Brief="Plate area = Atray - Adowncomer");
626       
627        VARIABLES
628in      Inlet as stream (Brief="Feed stream", PosX=0, PosY=0.4932, Symbol="_{in}");
629in      InletFV as stream (Brief="Feed stream", PosX=0, PosY=0.4932, Symbol="_{in}");
630in      InletLiquid as stream (Brief="Inlet liquid stream", PosX=0.5195, PosY=0, Symbol="_{inL}");
631in      InletVapour as stream (Brief="Inlet vapour stream", PosX=0.4994, PosY=1, Symbol="_{inV}");
632out     OutletLiquid as liquid_stream (Brief="Outlet liquid stream", PosX=0.8277, PosY=1, Symbol="_{outL}");
633out     OutletVapour as vapour_stream (Brief="Outlet vapour stream", PosX=0.8043, PosY=0, Symbol="_{outV}");
634
635        M_liq(NComp) as mol (Brief="Liquid Molar Holdup in the tray");
636        M_vap(NComp) as mol (Brief="Vapour Molar Holdup in the tray");
637        ML as mol (Brief="Molar liquid holdup");
638        MV as mol (Brief="Molar vapour holdup");
639        E_liq as energy (Brief="Total Liquid Energy Holdup on tray");
640        E_vap as energy (Brief="Total Vapour Energy Holdup on tray");
641        vL as volume_mol (Brief="Liquid Molar Volume");
642        vV as volume_mol (Brief="Vapour Molar volume");
643        Level as length (Brief="Height of clear liquid on plate");
644        interf as interfaceTeste;       
645
646        SET   
647        NC1=NComp-1;
648
649        EQUATIONS
650        "Component Molar Balance"
651        diff(M_liq)=Inlet.F*Inlet.z + InletLiquid.F*InletLiquid.z
652        - OutletLiquid.F*OutletLiquid.z + interf.NL;
653       
654        diff(M_vap)=InletFV.F*InletFV.z + InletVapour.F*InletVapour.z
655        - OutletVapour.F*OutletVapour.z - interf.NV;
656       
657        "Energy Balance"
658        diff(E_liq) = Inlet.F*Inlet.h + InletLiquid.F*InletLiquid.h
659                - OutletLiquid.F*OutletLiquid.h  + Q + interf.E_liq;
660       
661        diff(E_vap) = InletFV.F*InletFV.h + InletVapour.F*InletVapour.h
662                - OutletVapour.F*OutletVapour.h  - interf.E_vap;
663       
664        "Molar Holdup"
665        M_liq = ML*OutletLiquid.z;
666       
667        M_vap = MV*OutletVapour.z;
668       
669        "Energy Holdup"
670        E_liq = ML*(OutletLiquid.h - OutletLiquid.P*vL);
671       
672        E_vap = MV*(OutletVapour.h - OutletVapour.P*vV);
673       
674        "Energy Rate through the interface"
675        interf.E_liq = interf.htL*interf.a*(interf.T-OutletLiquid.T)+sum(interf.NL)*interf.hL; 
676       
677        interf.E_vap = interf.htV*interf.a*(OutletVapour.T-interf.T)+sum(interf.NV)*interf.hV;
678       
679        "Mass Conservation"
680        interf.NL = interf.NV;
681       
682        "Energy Conservation"
683        interf.E_liq = interf.E_vap;
684       
685        "Mol fraction normalisation"
686        sum(OutletLiquid.z)= 1.0;
687        sum(OutletLiquid.z)= sum(OutletVapour.z);
688        sum(interf.x)=1.0;
689        sum(interf.x)=sum(interf.y);
690       
691        "Liquid Volume"
692        vL = PP.LiquidVolume(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z);
693        "Vapour Volume"
694        vV = PP.VapourVolume(OutletVapour.T, OutletVapour.P, OutletVapour.z);
695       
696        "Chemical Equilibrium"
697        PP.LiquidFugacityCoefficient(interf.T, interf.P, interf.x)*interf.x =
698                PP.VapourFugacityCoefficient(interf.T, interf.P, interf.y)*interf.y;
699
700        "Geometry Constraint"
701        V = ML*vL + MV*vV;
702       
703        "Level of clear liquid over the weir"
704        Level = ML*vL/Ap;
705
706        "Total Mass Transfer Rates"
707        interf.NL(1:NC1)=interf.a*sumt(interf.kL*(interf.x(1:NC1)-OutletLiquid.z(1:NC1)))/vL+
708                OutletLiquid.z(1:NC1)*sum(interf.NL);
709
710#       interf.NL(1:NC1)=0.01*'kmol/s';
711       
712        interf.NV(1:NC1)=interf.a*sumt(interf.kV*(OutletVapour.z(1:NC1)-interf.y(1:NC1)))/vV+
713                OutletVapour.z(1:NC1)*sum(interf.NV);
714
715        "Mechanical Equilibrium"
716        OutletVapour.P = OutletLiquid.P;
717        interf.P=OutletLiquid.P;
718end
719
720Model trayRateTeste as trayRateBasicTeste
721        ATTRIBUTES
722        Pallete         = false;
723        Icon            = "icon/Tray";
724        Brief           = "Complete rate model of a column tray.";
725        Info            =
726"== Specify ==
727* the Feed stream
728* the Liquid inlet stream
729* the Vapour inlet stream
730* the Vapour outlet flow (OutletVapour.F)
731       
732== Initial ==
733* the plate temperature of both phases (OutletLiquid.T and OutletVapour.T)
734* the liquid height (Level) OR the liquid flow holdup (ML)
735* the vapor holdup (MV)
736* (NoComps - 1) OutletLiquid compositions
737";
738
739        PARAMETERS
740        Ah as area (Brief="Total holes area");
741        lw as length (Brief="Weir length");
742        g as acceleration (Default=9.81);
743        hw as length (Brief="Weir height");
744        beta as fraction (Brief="Aeration fraction");
745        alfa as fraction (Brief="Dry pressure drop coefficient");
746       
747        VapourFlow as Switcher(Valid = ["on", "off"], Default = "on");
748        LiquidFlow as Switcher(Valid = ["on", "off"], Default = "on");
749       
750        VARIABLES
751        rhoL as dens_mass;
752        rhoV as dens_mass;
753
754        EQUATIONS
755        "Liquid Density"
756        rhoL = PP.LiquidDensity(OutletLiquid.T, OutletLiquid.P, OutletLiquid.z);
757        "Vapour Density"
758        rhoV = PP.VapourDensity(InletVapour.T, InletVapour.P, InletVapour.z);
759
760        switch LiquidFlow
761                case "on":
762                "Francis Equation"
763#               OutletLiquid.F*vL = 1.84*'m^0.5/s'*lw*((Level-(beta*hw))/(beta))^1.5;
764                OutletLiquid.F*vL = 1.84*'1/s'*lw*((Level-(beta*hw))/(beta))^2;
765                when Level < (beta * hw) switchto "off";
766               
767                case "off":
768                "Low level"
769                OutletLiquid.F = 0 * 'mol/h';
770                when Level > (beta * hw) + 1e-6*'m' switchto "on";
771        end
772
773        switch VapourFlow
774                case "on":
775                InletVapour.F*vV = sqrt((InletVapour.P - OutletVapour.P)/(rhoV*alfa))*Ah;
776                when InletVapour.F < 1e-6 * 'kmol/h' switchto "off";
777               
778                case "off":
779                InletVapour.F = 0 * 'mol/s';
780                when InletVapour.P > OutletVapour.P + Level*g*rhoL + 1e-1 * 'atm' switchto "on";
781        end     
782end
783
784*#
Note: See TracBrowser for help on using the repository browser.