source: branches/gui/sample/miscellaneous/tenprobs/prob09.mso @ 697

#*---------------------------------------------------------------------
* 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.
*
*----------------------------------------------------------------------
* 9. Gas phase catalytic reactor
*----------------------------------------------------------------------
*
*   Description:
*      This problem is part of a collection of 10 representative
*       problems in Chemical Engineering for solution by numerical methods
*       developed for Cutlip (1998).
*
*   Subject:
*       * Reaction Engineering
*
*       Concepts utilized:
*               Design of a gas phase catalytic reactor with pressure drop for
*       a first order reversible gas phase reaction.
*
*       Numerical method:
*               * Simultaneous ODEs with known boundary conditions
*
*   Reference:
*               * CUTLIP et al. A collection of 10 numerical problems in 
*       chemical engineering solved by various mathematical software 
*       packages. Comp. Appl. in Eng. Education. v. 6, 169-180, 1998. 
*       * More informations and a detailed description of all problems
*       is available online in http://www.polymath-software.com/ASEE
* 
*----------------------------------------------------------------------
* Author: Rodolfo Rodrigues
* GIMSCOP/UFRGS - Group of Integration, Modeling, Simulation, 
*                                       Control, and Optimization of Processes
* $Id$
*--------------------------------------------------------------------*#
using "types";



Model stream
        PARAMETERS
outer NComp     as Integer              (Brief="Number of components", Lower=1);
        cp(NComp)       as cp_mol;
        
        SET
        cp = [40, 80]*'J/mol/K';
        
        VARIABLES
        C(NComp)        as conc_mol             (Brief="Molar concentration", Lower=0, DisplayUnit='mol/l');
        T                       as temperature;
        P                       as pressure;
end



FlowSheet reactor
        PARAMETERS
        NComp           as Integer              (Brief="Number of components", Lower=1);
        NReac           as Integer              (Brief="Number of reactions");
        
        Tref            as temperature  (Brief="Reference temperature", Default=450);
        Ta                      as temperature;
        
        ko(NReac)       as Real                 (Brief="Frequency factor at Treff", Unit='l^2/kg/min/mol');
        Ko(NReac)       as Real                 (Unit='l/mol');
        E(NReac)        as energy_mol   (Brief="Activation energy");
        R                       as Real                 (Brief="Universal gas constant", Unit='J/mol/K', Default=8.314);
        
        DH(NReac)       as enth_mol;
        U                       as Real                 (Unit='J/kg/K/min');
        alpha           as Real                 (Unit='1/kg');
        
        
        VARIABLES
        Inlet           as stream;
        Outlet          as stream;
        
        Fo(NComp)       as flow_mol             (Brief="Inlet molar flow");
        X                       as fraction             (Brief="Molar conversion of A");
        y                       as fraction             (Brief="Normalized pressure");
        Tn                      as Real                 (Brief="Temperature by 1000", Lower=0.273);
        
        
        r(NReac)        as Real                 (Brief="Mass rate of reaction", Unit='mol/kg/min');
        k(NReac)        as Real                 (Brief="Specific rate of reaction", Unit='l^2/kg/min/mol');
        K(NReac)        as Real                 (Brief="Equilibrium constant", Unit='l/mol');
        
        W                       as mass                 (Brief="Catalytic weight");
        
        
        EQUATIONS
        "Change time in W"
        W = time*'kg/s';
        
        "Mole balance"
        diff(X) = -r(1)/Fo(1)*'kg/s';
        
        "Rate of reaction"
        -r = k*(Outlet.C(1)^2 - Outlet.C(2)/K);
        
        
        "Specific rate of reaction"
        k = ko*exp(E/R*(1/Tref - 1/Outlet.T));
        
        "Equilibrium constant"
        K = Ko*exp(DH/R*(1/Tref - 1/Outlet.T));
        
        
        "Molar concentration of A"
        Outlet.C(1) = Inlet.C(1)*(1 - X)/(1 - 0.5*X)*y*(Inlet.T/Outlet.T);
        
        "Molar concentration of C"
        Outlet.C(2) = 0.5*Inlet.C(1)*X/(1 - 0.5*X)*y*(Inlet.T/Outlet.T);
        
        
        "Pressure drop"
        diff(y) = -alpha*(1 - 0.5*X)/2/y*(Outlet.T/Inlet.T)*'kg/s';
        
        "Energy balance"
        diff(Outlet.T) = (U*(Ta - Outlet.T) + r*DH)/(Fo(1)*Inlet.cp(1))*'kg/s';
        
        
        "Normalized pressure"
        y = Outlet.P/Inlet.P;
        
        "Temperature/1000"
        Tn = Outlet.T/1e3/'K';
        
        
        SET
        NComp = 2; # A, and C
        
        DH(1) = -4e4*'J/mol';
        ko(1) = 0.5*'l^2/kg/min/mol'; # at Tref
        Ko(1) = 2.5e4*'l/mol'; # at Tref
        E(1)  = 4.18e4*'J/mol';
        
        Ta = 500*'K';
        U = 0.8*'J/kg/K/min';
        alpha = 0.015/'kg';

        
        SPECIFY
        Fo = [5, 0]*'mol/min';
        Inlet.C = [0.271, 0]*'mol/l'; 
        Inlet.T = 450*'K';
        Inlet.P = 10*'atm';

        
        INITIAL
        "Molar conversion"
        X = 0;
        
        "Drop pressure"
        y = 1;
        
        "Temperature"
        Outlet.T = Inlet.T;


        OPTIONS
        TimeStart = 0;
        TimeStep = 0.25;
        TimeEnd = 20;
        
        TimeUnit = 's';
end