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

Last change on this file since 676 was 676, checked in by gerson bicca, 14 years ago

updates in column model

File size: 22.7 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 trayBasic
23        ATTRIBUTES
24        Pallete         = false;
25        Icon            = "icon/Tray";
26        Brief           = "Basic equations of a tray column model.";
27        Info            =
28"This model contains only the main equations of a column tray equilibrium model without
29the hidraulic equations.
30       
31== Assumptions ==
32* both phases (liquid and vapour) exists all the time;
33* thermodymanic equilibrium with Murphree plate efficiency;
34* no entrainment of liquid or vapour phase;
35* no weeping;
36* the dymanics in the downcomer are neglected.
37";
38       
39PARAMETERS
40outer PP                        as Plugin               (Brief = "External Physical Properties", Type="PP");
41outer NComp     as Integer;
42
43VARIABLES
44
45        Inlet                   as stream                               (Brief="Feed stream", Hidden=true, PosX=0, PosY=0.4932, Symbol="_{in}");
46        LiquidSideStream        as liquid_stream        (Brief="liquid Sidestream", Hidden=true, Symbol="_{outL}");
47        VapourSideStream        as vapour_stream        (Brief="vapour Sidestream", Hidden=true, Symbol="_{outV}");
48       
49in      InletL                  as stream                               (Brief="Inlet liquid stream", PosX=0.5195, PosY=0, Symbol="_{inL}");
50in      InletV                  as stream                               (Brief="Inlet vapour stream", PosX=0.4994, PosY=1, Symbol="_{inV}");
51out     OutletL         as liquid_stream        (Brief="Outlet liquid stream", PosX=0.8277, PosY=1, Symbol="_{outL}");
52out     OutletV         as vapour_stream        (Brief="Outlet vapour stream", PosX=0.8043, PosY=0, Symbol="_{outV}");
53
54
55        M(NComp)                        as mol                          (Brief="Molar Holdup in the tray");
56        ML                                              as mol                          (Brief="Molar liquid holdup");
57        MV                                              as mol                          (Brief="Molar vapour holdup");
58        E                                                       as energy                       (Brief="Total Energy Holdup on tray");
59        vL                                              as volume_mol   (Brief="Liquid Molar Volume");
60        vV                                              as volume_mol   (Brief="Vapour Molar volume");
61        Level                                   as length                       (Brief="Height of clear liquid on plate");
62        yideal(NComp)           as fraction;
63       
64       
65EQUATIONS
66"Component Molar Balance"
67        diff(M)=Inlet.F*Inlet.z + InletL.F*InletL.z + InletV.F*InletV.z- OutletL.F*OutletL.z - OutletV.F*OutletV.z-
68        LiquidSideStream.F*LiquidSideStream.z-VapourSideStream.F*VapourSideStream.z;
69
70"Molar Holdup"
71        M = ML*OutletL.z + MV*OutletV.z;
72
73"Mol fraction normalisation"
74        sum(OutletL.z)= 1.0;
75        sum(OutletL.z)= sum(OutletV.z);
76
77"Liquid Volume"
78        vL = PP.LiquidVolume(OutletL.T, OutletL.P, OutletL.z);
79
80"Vapour Volume"
81        vV = PP.VapourVolume(OutletV.T, OutletV.P, OutletV.z);
82
83"Chemical Equilibrium"
84        PP.LiquidFugacityCoefficient(OutletL.T, OutletL.P, OutletL.z)*OutletL.z = PP.VapourFugacityCoefficient(OutletV.T, OutletV.P, yideal)*yideal;
85
86"Thermal Equilibrium"
87        OutletV.T = OutletL.T;
88
89"Mechanical Equilibrium"
90        OutletV.P = OutletL.P;
91       
92"Thermal Equilibrium Vapour Side Stream"
93        OutletV.T = VapourSideStream.T;
94
95"Thermal Equilibrium Liquid Side Stream"
96        OutletL.T = LiquidSideStream.T;
97
98"Mechanical Equilibrium Vapour Side Stream"
99        OutletV.P= VapourSideStream.P;
100
101"Mechanical Equilibrium Liquid Side Stream"
102        OutletL.P = LiquidSideStream.P;
103
104"Composition Liquid Side Stream"
105        OutletL.z= LiquidSideStream.z;
106       
107"Composition Vapour Side Stream"
108        OutletV.z= VapourSideStream.z;
109
110end
111
112Model tray as trayBasic
113        ATTRIBUTES
114        Pallete         = false;
115        Icon            = "icon/Tray";
116        Brief           = "Complete model of a column tray.";
117        Info            =
118"== Specify ==
119* the Feed stream
120* the Liquid inlet stream
121* the Vapour inlet stream
122* the Vapour outlet flow (OutletV.F)
123       
124== Initial ==
125* the plate temperature (OutletL.T)
126* the liquid height (Level) OR the liquid flow OutletL.F
127* (NoComps - 1) OutletL compositions
128
129== Options ==
130You can choose the equation for the liquid outlet flow and the vapour
131inlet flow calculation through the VapourFlowModel and LiquidFlowModel
132switchers.
133
134== References ==
135* 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.
136* FEEHERY, W. F. Dynamic Optimization with Path Constraints. Tese (Doutorado) - Massachusetts Institute of Technology, June 1998.
137* KLINGBERG, A. Modeling and Optimization of Batch Distillation. Dissertação (Mestrado) - Department of Automatic Control, Lund Institute of Technology, Lund, Sweden, fev. 2000.
138* 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.
139* 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.
140* 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.
141* 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.
142";     
143
144VARIABLES
145        rhoL as dens_mass;
146        rhoV as dens_mass;
147
148        btemp as Real (Brief="Temporary variable of Roffels liquid flow equation");
149       
150EQUATIONS
151
152"Liquid Density"
153        rhoL = PP.LiquidDensity(OutletL.T, OutletL.P, OutletL.z);
154       
155"Vapour Density"
156        rhoV = PP.VapourDensity(InletV.T, InletV.P, InletV.z);
157
158end
159
160#*-------------------------------------------------------------------
161* Model of a tray with reaction
162*-------------------------------------------------------------------*#
163Model trayReactTeste
164        ATTRIBUTES
165        Pallete         = false;
166        Icon            = "icon/Tray";
167        Brief           = "Model of a tray with reaction.";
168        Info            =
169"== Assumptions ==
170* both phases (liquid and vapour) exists all the time;
171* thermodymanic equilibrium with Murphree plate efficiency;
172* no entrainment of liquid or vapour phase;
173* no weeping;
174* the dymanics in the downcomer are neglected.
175       
176== Specify ==
177* the Feed stream;
178* the Liquid inlet stream;
179* the Vapour inlet stream;
180* the Vapour outlet flow (OutletV.F);
181* the reaction related variables.
182       
183== Initial ==
184* the plate temperature (OutletL.T)
185* the liquid height (Level) OR the liquid flow OutletL.F
186* (NoComps - 1) OutletL compositions
187";
188
189        PARAMETERS
190        outer PP as Plugin(Type="PP");
191        outer NComp as Integer;
192        V as volume(Brief="Total Volume of the tray");
193        Q as power (Brief="Rate of heat supply");
194        Ap as area (Brief="Plate area = Atray - Adowncomer");
195       
196        Ah as area (Brief="Total holes area");
197        lw as length (Brief="Weir length");
198        g as acceleration (Default=9.81);
199        hw as length (Brief="Weir height");
200        beta as fraction (Brief="Aeration fraction");
201        alfa as fraction (Brief="Dry pressure drop coefficient");
202
203        stoic(NComp) as Real(Brief="Stoichiometric matrix");
204        Hr as energy_mol;
205        Pstartup as pressure;
206       
207        VapourFlow as Switcher(Valid = ["on", "off"], Default = "off");
208        LiquidFlow as Switcher(Valid = ["on", "off"], Default = "off");
209       
210        VARIABLES
211in      Inlet as stream (Brief="Feed stream", PosX=0, PosY=0.4932, Symbol="_{in}");
212in      InletL as stream (Brief="Inlet liquid stream", PosX=0.5195, PosY=0, Symbol="_{inL}");
213in      InletV as stream (Brief="Inlet vapour stream", PosX=0.4994, PosY=1, Symbol="_{inV}");
214out     OutletL as liquid_stream (Brief="Outlet liquid stream", PosX=0.8277, PosY=1, Symbol="_{outL}");
215out     OutletV as vapour_stream (Brief="Outlet vapour stream", PosX=0.8043, PosY=0, Symbol="_{outV}");
216
217        yideal(NComp) as fraction;
218        Emv as Real (Brief = "Murphree efficiency");
219
220        M(NComp) as mol (Brief="Molar Holdup in the tray");
221        ML as mol (Brief="Molar liquid holdup");
222        MV as mol (Brief="Molar vapour holdup");
223        E as energy (Brief="Total Energy Holdup on tray");
224        vL as volume_mol (Brief="Liquid Molar Volume");
225        vV as volume_mol (Brief="Vapour Molar volume");
226        Level as length (Brief="Height of clear liquid on plate");
227        Vol as volume;
228       
229        rhoL as dens_mass;
230        rhoV as dens_mass;
231        r3 as reaction_mol (Brief = "Reaction resulting ethyl acetate", DisplayUnit = 'mol/l/s');
232        C(NComp) as conc_mol (Brief = "Molar concentration", Lower = -1); #, Unit = "mol/l");
233       
234        EQUATIONS
235        "Molar Concentration"
236        OutletL.z = vL * C;
237       
238        "Reaction"
239        r3 = exp(-7150*'K'/OutletL.T)*(4.85e4*C(1)*C(2) - 1.23e4*C(3)*C(4))*'l/mol/s';
240       
241        "Component Molar Balance"
242        diff(M)=Inlet.F*Inlet.z + InletL.F*InletL.z + InletV.F*InletV.z
243                - OutletL.F*OutletL.z - OutletV.F*OutletV.z + stoic*r3*ML*vL;
244       
245        "Energy Balance"
246        diff(E) = ( Inlet.F*Inlet.h + InletL.F*InletL.h + InletV.F*InletV.h
247                - OutletL.F*OutletL.h - OutletV.F*OutletV.h + Q ) + Hr * r3 * vL*ML;
248       
249        "Molar Holdup"
250        M = ML*OutletL.z + MV*OutletV.z;
251       
252        "Energy Holdup"
253        E = ML*OutletL.h + MV*OutletV.h - OutletL.P*V;
254       
255        "Mol fraction normalisation"
256        sum(OutletL.z)= 1.0;
257       
258        "Liquid Volume"
259        vL = PP.LiquidVolume(OutletL.T, OutletL.P, OutletL.z);
260        "Vapour Volume"
261        vV = PP.VapourVolume(OutletV.T, OutletV.P, OutletV.z);
262
263        "Thermal Equilibrium"
264        OutletV.T = OutletL.T;
265       
266        "Mechanical Equilibrium"
267        OutletV.P = OutletL.P;
268       
269        "Level of clear liquid over the weir"
270        Level = ML*vL/Ap;
271
272        Vol = ML*vL;
273       
274        "Liquid Density"
275        rhoL = PP.LiquidDensity(OutletL.T, OutletL.P, OutletL.z);
276        "Vapour Density"
277        rhoV = PP.VapourDensity(InletV.T, InletV.P, InletV.z);
278
279        switch LiquidFlow
280                case "on":
281                "Francis Equation"
282                OutletL.F*vL = 1.84*'1/s'*lw*((Level-(beta*hw)+1e-6*'m')/(beta))^2;
283                when Level < (beta * hw) switchto "off";
284               
285                case "off":
286                "Low level"
287                OutletL.F = 0 * 'mol/h';
288                when Level > (beta * hw) + 1e-6*'m' switchto "on";
289        end
290
291        switch VapourFlow
292                case "on":
293                #InletV.P = OutletV.P + Level*g*rhoL + rhoV*alfa*(InletV.F*vV/Ah)^2;
294                InletV.F*vV = sqrt((InletV.P - OutletV.P - Level*g*rhoL + 1e-8 * 'atm')/(rhoV*alfa))*Ah;
295                when InletV.P < OutletV.P + Level*g*rhoL switchto "off";
296               
297                case "off":
298                InletV.F = 0 * 'mol/s';
299                when InletV.P > OutletV.P + Level*g*rhoL + 3e-2 * 'atm' switchto "on";
300                #when InletV.P > OutletV.P + Level*beta*g*rhoL + 1e-2 * 'atm' switchto "on";
301        end
302
303        "Chemical Equilibrium"
304        PP.LiquidFugacityCoefficient(OutletL.T, OutletL.P, OutletL.z)*OutletL.z =
305                PP.VapourFugacityCoefficient(OutletV.T, OutletV.P, yideal)*yideal;
306       
307        OutletV.z = Emv * (yideal - InletV.z) + InletV.z;
308       
309        sum(OutletL.z)= sum(OutletV.z);
310       
311        "Geometry Constraint"
312        V = ML* vL + MV*vV;
313end
314
315#*-------------------------------------
316* Model of a packed column stage
317-------------------------------------*#
318Model packedStageTeste
319        ATTRIBUTES
320        Pallete         = false;
321        Icon            = "icon/PackedStage";
322        Brief           = "Complete model of a packed column stage.";
323        Info            =
324"== Specify ==
325* the Feed stream
326* the Liquid inlet stream
327* the Vapour inlet stream
328* the stage pressure drop (deltaP)
329       
330== Initial ==
331* the plate temperature (OutletL.T)
332* the liquid molar holdup ML
333* (NoComps - 1) OutletL compositions
334";     
335       
336        PARAMETERS
337outer PP as Plugin(Brief = "External Physical Properties", Type="PP");
338outer NComp as Integer;
339        PPwater as Plugin(Brief="Physical Properties",
340                Type="PP",
341                Components = [ "water" ],
342                LiquidModel = "PR",
343                VapourModel = "PR"
344        );
345
346        V as volume(Brief="Total Volume of the tray");
347        Q as heat_rate (Brief="Rate of heat supply");
348        d as length (Brief="Column diameter"); 
349
350        a as Real (Brief="surface area per packing volume", Unit='m^2/m^3');
351        g as acceleration;
352        e as Real (Brief="Void fraction of packing, m^3/m^3");
353        Cpo as Real (Brief="Constant for resitance equation"); # Billet and Schultes, 1999.
354        Mw(NComp)       as molweight    (Brief = "Component Mol Weight");
355        hs as length (Brief="Height of the packing stage");
356        Qsil as positive (Brief="Resistance coefficient on the liquid load", Default=1);
357
358        VARIABLES
359in      Inlet as stream (Brief="Feed stream", PosX=0, PosY=0.4932, Symbol="_{in}");
360in      InletL as stream (Brief="Inlet liquid stream", PosX=0.5195, PosY=0, Symbol="_{inL}");
361in      InletV as stream (Brief="Inlet vapour stream", PosX=0.4994, PosY=1, Symbol="_{inV}");
362out     OutletL as liquid_stream (Brief="Outlet liquid stream", PosX=0.8277, PosY=1, Symbol="_{outL}");
363out     OutletV as vapour_stream (Brief="Outlet vapour stream", PosX=0.8043, PosY=0, Symbol="_{outV}");
364
365        M(NComp) as mol (Brief="Molar Holdup in the tray", Default=0.01, Lower=0, Upper=100);
366        ML as mol (Brief="Molar liquid holdup", Default=0.01, Lower=0, Upper=100);
367        MV as mol (Brief="Molar vapour holdup", Default=0.01, Lower=0, Upper=100);
368        E as energy (Brief="Total Energy Holdup on tray", Default=-500);
369        vL as volume_mol (Brief="Liquid Molar Volume");
370        vV as volume_mol (Brief="Vapour Molar volume");
371       
372        miL as viscosity (Brief="Liquid dynamic viscosity", DisplayUnit='kg/m/s');
373        miV as viscosity (Brief="Vapor dynamic viscosity", DisplayUnit='kg/m/s');
374        rhoL as dens_mass;
375        rhoV as dens_mass;
376       
377        deltaP as pressure;
378       
379        uL as velocity (Brief="volume flow rate of liquid, m^3/m^2/s", Lower=-10, Upper=100);
380        uV as velocity (Brief="volume flow rate of vapor, m^3/m^2/s", Lower=-10, Upper=100);
381        dp as length (Brief="Particle diameter", Default=1e-3, Lower=0, Upper=10);
382        invK as positive (Brief="Wall factor", Default=1, Upper=10);
383        Rev as Real (Brief="Reynolds number of the vapor stream", Default=4000);
384        Al as area (Brief="Area occupied by the liquid", Default=0.001, Upper=1);
385        hl as positive (Brief="Column holdup", Unit='m^3/m^3', Default=0.01,Upper=10);
386
387        SET
388        Mw = PP.MolecularWeight();
389
390        EQUATIONS
391        "Component Molar Balance"
392        diff(M)=Inlet.F*Inlet.z + InletL.F*InletL.z + InletV.F*InletV.z
393                - OutletL.F*OutletL.z - OutletV.F*OutletV.z;
394
395        "Energy Balance"
396        diff(E) = ( Inlet.F*Inlet.h + InletL.F*InletL.h + InletV.F*InletV.h
397                - OutletL.F*OutletL.h - OutletV.F*OutletV.h + Q );
398       
399        "Molar Holdup"
400        M = ML*OutletL.z + MV*OutletV.z;
401       
402        "Energy Holdup"
403        E = ML*OutletL.h + MV*OutletV.h - OutletL.P*V;
404       
405        "Mol fraction normalisation"
406        sum(OutletL.z)= 1.0;
407       
408        "Liquid Volume"
409        vL = PP.LiquidVolume(OutletL.T, OutletL.P, OutletL.z);
410        "Vapour Volume"
411        vV = PP.VapourVolume(OutletV.T, OutletV.P, OutletV.z);
412       
413        "Chemical Equilibrium"
414        PP.LiquidFugacityCoefficient(OutletL.T, OutletL.P, OutletL.z)*OutletL.z =
415                PP.VapourFugacityCoefficient(OutletV.T, OutletV.P, OutletV.z)*OutletV.z;
416       
417        "Thermal Equilibrium"
418        OutletV.T = OutletL.T;
419       
420        "Mechanical Equilibrium"
421        OutletL.P = OutletV.P;
422       
423        "Geometry Constraint"
424        V*e = ML*vL + MV*vV;
425       
426        "Liquid Density"
427        rhoL = PP.LiquidDensity(OutletL.T, OutletL.P, OutletL.z);
428        "Vapour Density"
429        rhoV = PP.VapourDensity(InletV.T, InletV.P, InletV.z);
430        "Liquid viscosity"
431        miL = PP.LiquidViscosity(OutletL.T, OutletL.P, OutletL.z);
432        "Vapour viscosity"
433        miV = PP.VapourViscosity(InletV.T, InletV.P, InletV.z);
434
435        "Area occupied by the liquid"
436        Al = ML*vL/hs;
437
438        "Volume flow rate of liquid, m^3/m^2/s"
439        uL * Al = OutletL.F * vL;
440        "Volume flow rate of vapor, m^3/m^2/s"
441        uV * ((d^2*3.14159/4)*e - Al) = OutletV.F * vV;
442       
443        "Liquid holdup"
444        hl = ML*vL/V/e;
445       
446        "Particle diameter"
447        dp = 6 * (1-e)/a;
448       
449        "Wall Factor"
450        invK = (1 + (2*dp/(3*d*(1-e))));
451       
452        "Reynolds number of the vapor stream"
453        Rev*invK = dp*uV*rhoV / (miV*(1-e));
454       
455        deltaP = InletV.P - OutletV.P;
456       
457        "Pressure drop and Vapor flow"
458        deltaP/hs  = Qsil*a*uV^2*rhoV*invK / (2*(e-hl)^3);
459
460        "Liquid holdup"
461        hl = (12*miL*a^2*uL/rhoL/g)^1/3;
462end
463
464#*-------------------------------------
465* Nonequilibrium Model
466-------------------------------------*#
467Model interfaceTeste
468       
469        ATTRIBUTES
470        Pallete         = false;
471        Icon            = "icon/Tray";
472        Brief           = "Descrition of variables of the equilibrium interface.";
473        Info            =
474"This model contains only the variables of the equilibrium interface.";
475
476        PARAMETERS
477outer PP as Plugin(Brief = "External Physical Properties", Type="PP");
478outer NComp as Integer;
479outer NC1 as Integer;
480       
481        VARIABLES
482        NL(NComp) as flow_mol_delta     (Brief = "Stream Molar Rate on Liquid Phase");
483        NV(NComp) as flow_mol_delta     (Brief = "Stream Molar Rate on Vapour Phase");
484        T as temperature                (Brief = "Stream Temperature");
485        P as pressure                   (Brief = "Stream Pressure");
486        x(NComp) as fraction    (Brief = "Stream Molar Fraction on Liquid Phase");
487        y(NComp) as fraction    (Brief = "Stream Molar Fraction on Vapour Phase");
488        a as area                           (Brief = "Interface Area");
489        htL as heat_trans_coeff (Brief = "Heat Transference Coefficient on Liquid Phase");
490        htV as heat_trans_coeff (Brief = "Heat Transference Coefficient on Vapour Phase");     
491        E_liq as heat_rate      (Brief = "Liquid Energy Rate at interface");
492    E_vap as heat_rate      (Brief = "Vapour Energy Rate at interface");       
493        hL as enth_mol          (Brief = "Liquid Molar Enthalpy");
494        hV as enth_mol          (Brief = "Vapour Molar Enthalpy");
495        kL(NC1,NC1) as velocity (Brief = "Mass Transfer Coefficients");
496        kV(NC1,NC1) as velocity (Brief = "Mass Transfer Coefficients");
497       
498        EQUATIONS
499        "Liquid Enthalpy"
500        hL = PP.LiquidEnthalpy(T, P, x);
501       
502        "Vapour Enthalpy"
503        hV = PP.VapourEnthalpy(T, P, y);
504
505end
506
507Model trayRateBasicTeste
508        ATTRIBUTES
509        Pallete         = false;
510        Icon            = "icon/Tray";
511        Brief           = "Basic equations of a tray rate column model.";
512        Info            =
513"This model contains only the main equations of a column tray nonequilibrium model without
514the hidraulic equations.
515       
516== Assumptions ==
517* both phases (liquid and vapour) exists all the time;
518* no entrainment of liquid or vapour phase;
519* no weeping;
520* the dymanics in the downcomer are neglected.
521";
522       
523        PARAMETERS
524outer PP as Plugin(Brief = "External Physical Properties", Type="PP");
525outer NComp as Integer;
526    NC1 as Integer;
527        V as volume(Brief="Total Volume of the tray");
528        Q as heat_rate (Brief="Rate of heat supply");
529        Ap as area (Brief="Plate area = Atray - Adowncomer");
530       
531        VARIABLES
532in      Inlet as stream (Brief="Feed stream", PosX=0, PosY=0.4932, Symbol="_{in}");
533in      InletFV as stream (Brief="Feed stream", PosX=0, PosY=0.4932, Symbol="_{in}");
534in      InletL as stream (Brief="Inlet liquid stream", PosX=0.5195, PosY=0, Symbol="_{inL}");
535in      InletV as stream (Brief="Inlet vapour stream", PosX=0.4994, PosY=1, Symbol="_{inV}");
536out     OutletL as liquid_stream (Brief="Outlet liquid stream", PosX=0.8277, PosY=1, Symbol="_{outL}");
537out     OutletV as vapour_stream (Brief="Outlet vapour stream", PosX=0.8043, PosY=0, Symbol="_{outV}");
538
539        M_liq(NComp) as mol (Brief="Liquid Molar Holdup in the tray");
540        M_vap(NComp) as mol (Brief="Vapour Molar Holdup in the tray");
541        ML as mol (Brief="Molar liquid holdup");
542        MV as mol (Brief="Molar vapour holdup");
543        E_liq as energy (Brief="Total Liquid Energy Holdup on tray");
544        E_vap as energy (Brief="Total Vapour Energy Holdup on tray");
545        vL as volume_mol (Brief="Liquid Molar Volume");
546        vV as volume_mol (Brief="Vapour Molar volume");
547        Level as length (Brief="Height of clear liquid on plate");
548        interf as interfaceTeste;       
549
550        SET   
551        NC1=NComp-1;
552
553        EQUATIONS
554        "Component Molar Balance"
555        diff(M_liq)=Inlet.F*Inlet.z + InletL.F*InletL.z
556        - OutletL.F*OutletL.z + interf.NL;
557       
558        diff(M_vap)=InletFV.F*InletFV.z + InletV.F*InletV.z
559        - OutletV.F*OutletV.z - interf.NV;
560       
561        "Energy Balance"
562        diff(E_liq) = Inlet.F*Inlet.h + InletL.F*InletL.h
563                - OutletL.F*OutletL.h  + Q + interf.E_liq;
564       
565        diff(E_vap) = InletFV.F*InletFV.h + InletV.F*InletV.h
566                - OutletV.F*OutletV.h  - interf.E_vap;
567       
568        "Molar Holdup"
569        M_liq = ML*OutletL.z;
570       
571        M_vap = MV*OutletV.z;
572       
573        "Energy Holdup"
574        E_liq = ML*(OutletL.h - OutletL.P*vL);
575       
576        E_vap = MV*(OutletV.h - OutletV.P*vV);
577       
578        "Energy Rate through the interface"
579        interf.E_liq = interf.htL*interf.a*(interf.T-OutletL.T)+sum(interf.NL)*interf.hL;       
580       
581        interf.E_vap = interf.htV*interf.a*(OutletV.T-interf.T)+sum(interf.NV)*interf.hV;
582       
583        "Mass Conservation"
584        interf.NL = interf.NV;
585       
586        "Energy Conservation"
587        interf.E_liq = interf.E_vap;
588       
589        "Mol fraction normalisation"
590        sum(OutletL.z)= 1.0;
591        sum(OutletL.z)= sum(OutletV.z);
592        sum(interf.x)=1.0;
593        sum(interf.x)=sum(interf.y);
594       
595        "Liquid Volume"
596        vL = PP.LiquidVolume(OutletL.T, OutletL.P, OutletL.z);
597        "Vapour Volume"
598        vV = PP.VapourVolume(OutletV.T, OutletV.P, OutletV.z);
599       
600        "Chemical Equilibrium"
601        PP.LiquidFugacityCoefficient(interf.T, interf.P, interf.x)*interf.x =
602                PP.VapourFugacityCoefficient(interf.T, interf.P, interf.y)*interf.y;
603
604        "Geometry Constraint"
605        V = ML*vL + MV*vV;
606       
607        "Level of clear liquid over the weir"
608        Level = ML*vL/Ap;
609
610        "Total Mass Transfer Rates"
611        interf.NL(1:NC1)=interf.a*sumt(interf.kL*(interf.x(1:NC1)-OutletL.z(1:NC1)))/vL+
612                OutletL.z(1:NC1)*sum(interf.NL);
613
614#       interf.NL(1:NC1)=0.01*'kmol/s';
615       
616        interf.NV(1:NC1)=interf.a*sumt(interf.kV*(OutletV.z(1:NC1)-interf.y(1:NC1)))/vV+
617                OutletV.z(1:NC1)*sum(interf.NV);
618
619        "Mechanical Equilibrium"
620        OutletV.P = OutletL.P;
621        interf.P=OutletL.P;
622end
623
624Model trayRateTeste as trayRateBasicTeste
625        ATTRIBUTES
626        Pallete         = false;
627        Icon            = "icon/Tray";
628        Brief           = "Complete rate model of a column tray.";
629        Info            =
630"== Specify ==
631* the Feed stream
632* the Liquid inlet stream
633* the Vapour inlet stream
634* the Vapour outlet flow (OutletV.F)
635       
636== Initial ==
637* the plate temperature of both phases (OutletL.T and OutletV.T)
638* the liquid height (Level) OR the liquid flow holdup (ML)
639* the vapor holdup (MV)
640* (NoComps - 1) OutletL compositions
641";
642
643        PARAMETERS
644        Ah as area (Brief="Total holes area");
645        lw as length (Brief="Weir length");
646        g as acceleration (Default=9.81);
647        hw as length (Brief="Weir height");
648        beta as fraction (Brief="Aeration fraction");
649        alfa as fraction (Brief="Dry pressure drop coefficient");
650       
651        VapourFlow as Switcher(Valid = ["on", "off"], Default = "on");
652        LiquidFlow as Switcher(Valid = ["on", "off"], Default = "on");
653       
654        VARIABLES
655        rhoL as dens_mass;
656        rhoV as dens_mass;
657
658        EQUATIONS
659        "Liquid Density"
660        rhoL = PP.LiquidDensity(OutletL.T, OutletL.P, OutletL.z);
661        "Vapour Density"
662        rhoV = PP.VapourDensity(InletV.T, InletV.P, InletV.z);
663
664        switch LiquidFlow
665                case "on":
666                "Francis Equation"
667#               OutletL.F*vL = 1.84*'m^0.5/s'*lw*((Level-(beta*hw))/(beta))^1.5;
668                OutletL.F*vL = 1.84*'1/s'*lw*((Level-(beta*hw))/(beta))^2;
669                when Level < (beta * hw) switchto "off";
670               
671                case "off":
672                "Low level"
673                OutletL.F = 0 * 'mol/h';
674                when Level > (beta * hw) + 1e-6*'m' switchto "on";
675        end
676
677        switch VapourFlow
678                case "on":
679                InletV.F*vV = sqrt((InletV.P - OutletV.P)/(rhoV*alfa))*Ah;
680                when InletV.F < 1e-6 * 'kmol/h' switchto "off";
681               
682                case "off":
683                InletV.F = 0 * 'mol/s';
684                when InletV.P > OutletV.P + Level*g*rhoL + 1e-1 * 'atm' switchto "on";
685        end     
686end
Note: See TracBrowser for help on using the repository browser.