#*-------------------------------------------------------------------
* EMSO Model Library (EML) Copyright (C) 2004 - 2007 ALSOC.
*
* This LIBRARY is free software; you can distribute it and/or modify
* it under the therms of the ALSOC FREE LICENSE as available at
* http://www.enq.ufrgs.br/alsoc.
*
* EMSO Copyright (C) 2004 - 2007 ALSOC, original code
* from http://www.rps.eng.br Copyright (C) 2002-2004.
* All rights reserved.
*
* EMSO is distributed under the therms of the ALSOC LICENSE as
* available at http://www.enq.ufrgs.br/alsoc.
*
*--------------------------------------------------------------------
* Sample file for for a high-index optimal control problem.
*--------------------------------------------------------------------
* Author: Rafael de Pelegrini Soares
* $Id$
*--------------------------------------------------------------------*#
using "types";
Model FlashRaoult
ATTRIBUTES
Info = "
This is a very simple (wrong) model with dynamics only on the energy.
It should be used for ilustration purposes only.";
PARAMETERS
NComp as Integer;
# Antoine constants
A(NComp) as Real;
B(NComp) as Real;
C(NComp) as Real;
Cv as Real(Unit = 'J/mol/K');
DHvap(NComp) as energy_mol;
VARIABLES
F as flow_mol;
L as flow_mol;
V as flow_mol;
z(NComp) as fraction;
x(NComp) as fraction;
y(NComp) as fraction;
n(NComp) as mol;
nt as mol;
T as temperature;
P as pressure;
Psat(NComp) as pressure;
Q as power;
E as energy;
EQUATIONS
"Component Molar Balance"
diff(n) = F*z - (L*x + V*y);
nt = sum(n);
x = n/nt;
"Energy Balance"
diff(E) = Q - V*sum(DHvap*y);
"Internal energy"
E = nt*Cv*T;
"Raoult's Law"
P*y = Psat*x;
"Antoine for Vapour Pressure"
ln(Psat/'kPa') = A - B/(T/'K'- 273.15 + C);
"Molar Fraction sum"
sum(y) = 1;
end
FlowSheet FlashRaoultTest
DEVICES
fl as FlashRaoult;
SET
fl.NComp = 3;
# Antoine constants (Acetone, Acetonitrile, Nitromethane)
fl.A = [14.31, 14.89, 14.75];
fl.B = [2756, 3413, 3331];
fl.C = [228, 250, 227];
fl.Cv = 30 * 'J/mol/K';
fl.DHvap = [10, 20, 30] * 'kJ/mol';
EQUATIONS
# Disturb
#if time < 0.5 * 'h' then
# fl.z = [0.45, 0.35, 0.2];
#else
# fl.z = [0.55, 0.25, 0.2];
#end
SPECIFY
# Feed condition
fl.F = 1 * 'kmol/h';
# Steady-state feed
#fl.z = [0.45, 0.35, 0.2];
# Disturb on feed composition
fl.z = [0.55, 0.25, 0.2];
# Desired production of y1 (index 2 - will determine the "perfect" heat profile)
fl.V = 0.1 * 'kmol/h';
fl.y(1) = 0.7;
fl.nt = 1 * 'kmol';
# Default specification (index 1 - heat is given)
#fl.V = 0.1 * 'kmol/h';
#fl.Q = 0.37 * 'kW';
#fl.nt = 1 * 'kmol';
# Fixed Temperature (index 2 - will determine the "perfect" heat profile)
#fl.T = (90+273.15) * 'K';
#fl.V = 0.1 * 'kmol/h';
#fl.nt = 1 * 'kmol';
# Fixed Pressure (index 2 - will determine the "perfect" heat profile)
#fl.P = 1.2 * 'atm';
#fl.V = 0.1 * 'kmol/h';
#fl.nt = 1 * 'kmol';
INITIAL
#fl.T = (80+273.15) * 'K';
#fl.nt = 1 * 'kmol';
# Steady-state composition with feed = [0.45, 0.35, 0.2];
fl.x(1) = 0.4222;
fl.x(2) = 0.3622;
# steady state compositions
# diff(fl.n(1:2)) = 0 * 'mol/s';
OPTIONS
TimeEnd = 4;
TimeStep = 0.05;
TimeUnit = 'h';
#DAESolver(File="dasslc"); # slow integration
DAESolver(File="mebdf"); # much faster
Dynamic = true;
end