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

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

updates

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