Introduction to Matlab - PowerPoint by tyndale

VIEWS: 36 PAGES: 16

									 Differential-Algebraic Equation Systems
1. DAE systems
2. Binary flash unit
3. Exothermic CSTR
4. Matlab example
                        DAE Systems
   General model structure
       dy
           f ( x, y , z ) y ( x0 )  y 0 y  R n
       dx
       0  g( x, y, z ) 0  g( x0 , y 0 , z 0 ) x  R m

    » Differential-algebraic equation (DAE) system
    » Comprised of differential (x) and algebraic variables (z)
    » Initial condition for algebraic variables must be computed
   Numerical solution
    » Considerably more difficult than ODE systems due to
      presence of algebraic equations (AEs)
    » Need specialized techniques
    » Difficulty determined by DAE system index
                     DAE System Index
   Definition
    » The minimum number of times the algebraic equations
      must differentiated to recover the algebraic variables
    » Most codes developed for index 0 & index 1 systems
   Index 0 DAE systems
           0  g( x, y, z )  z  Γ( x, y )

     dy
         f ( x, y, z )  f [ x, y, Γ( x, y )]  F( x, y )
     dx
                  Index 2 DAE Systems

                               g g dy g dz
      0  g ( x, y , z )  0         
                               x y dx z dx
           1                             1
dz     g   g g dy      g              g g            
        x y dx    z 
                                                f ( x, y , z )   F ( x, y , z )
                                                x y            
dx     z                                                  

                    dy
                        f ( x, y , z )
                    dx
            Implicit Function Theorem
   Problem
                                                z h ( x , y )
                            m  n 1
    0  g( x, y, z) g : R              R   m
                                                                0  g[ x, y, h( x, y )]

   Assumptions
                            g       
      g (a, b, c)  0 rank  (a, b, c)  m
                            z       

   Result – there is a neighborhood of a & b where:
                     h ( a, b )  c
       h( x, y )  
                     g[ x, y, h( x, y )]  0
   Implication – the function g(x,y,z) can be solved
    numerically for z
                                  Binary Flash Unit
                  Schematic diagram                             Assumptions
                                                                  »   Saturated liquid feed
                                        V, y
                                                   Vapor          »   Perfect mixing
                                                                  »   Negligible heat losses
                                                                  »   Negligible vapor holdup
           F, xF                                                  »   Negligible energy
Feed                                                                  accumulation
                                               Q                  »   Constant heat of
                             H                         Heat
                                                                      vaporization
                                                                  »   Constant relative
                                 L, x
                                          Liquid                      volatility


                  Vapor-liquid equilibrium
                              x
                      y
                         1  (  1) x
                    Model Formulation
      Mass balance
         dH
             F V  L
         dt
      Component balance

                    d ( Hx)
                             FxF  Vy  Lx
                       dt

        d ( Hx)    dH    dx                     dx
                x    H     x( F  V  L)  H
           dt      dt    dt                     dt

  dx
H     ( F  V  L) x  FxF  Vy  Lx  F ( xF  x)  V ( y  x)
  dt
              Model Formulation cont.
   Steady-state energy balance
                       Q
      VH v  Q  V 
                      H v
   Index 0 DAE model

       dx                  Q
     H     F ( xF  x)        ( y  x)   dx
       dt                 H                    f ( x, y )
                        x                 dt
          0 y                             0  g ( x, y )
                   1  (  1) x
                          Exothermic CSTR

                                            Reaction: A  B
                                            Eliminate assumption of
                                             constant cooling jacket
                                             temperature
                                            Approximate cooling jacket
                                             as well-mixed system



    Mass and energy balances on tank

                  dC A
                V       q(C Ai  C A )  Vk0 exp( E / RT )C A
                   dt
         dT
    VCp     qC p (Ti  T )  (H )Vk0 exp( E / RT )C A  UA(Tc  T )
         dt
                 Exothermic CSTR cont.
   Energy balance on cooling jacket
               dTc
      cVcC pc       c qcC pc (Tci  Tc )  UA(Tc  T )
                dt

   Assume small thermal capacitance
                 dTc
      cVcC pc        0  c qcC pc (Tci  Tc )  UA(Tc  T )
                  dt
   Index 0 DAE model
                   dC A
                 V       q(C Ai  C A )  Vk0 exp( E / RT )c A
                    dt
         dT
    VCp     qC p (Ti  T )  (H )Vk0 exp( E / RT )c A  UA(Tc  T )
         dt
                  0   c qcC pc (Tci  Tc )  UA(Tc  T )
                   Matlab Mass Matrix
    Mass matrix representation
                  dy
       M ( x, y )     f ( x, y )
                  dx

    CSTR example
         1 0 0
     M  0 1 0 
               
         0 0 0 
               
                              dy1 
                              dx 
                              dy   1
                                       f ( x, y1 , y2 , y3 ) 
     ( x, y )  f ( x, y )   2    f 2 ( x, y1 , y2 , y3 )
  dy
M
  dx                          dx                           
                              0   f 3 ( x, y1 , y2 , y3 ) 
                                                             
                             
                                  
                                   
                           M-file cstr_dae.m

function f = cstr_dae(y,qc)

q = 1;           %   m^3/h
Cai = 10;        %   kmol/m^3
Ti = 298;        %   K
ko = 3.493e7;    %   h^-1
E = 11843;       %   kcal/kmol
H = -5960;       %   kcal/kmol
rhoCp = 500;     %   kcal/m^3/K
UA = 150;        %   kcal/H/K
R = 1.987;       %   kcal/kmol/K
V = 1;           %   m^3
Tci = Ti;

Ca = y(1);   T   = y(2);    Tc = y(3);

f = [
    q/V*(Cai - Ca) - ko*exp(-E/R/T)*Ca;
    q/V*(Ti-T) - H/rhoCp*ko*exp(-E/R/T)*Ca + UA/V/rhoCp*(Tc-T);
    rhoCp*qc*(Tci-Tc)-UA*(Tc-T)
];
             CSTR Model Solution
 Find a steady-state solution
>> qc = 1;
>> f = @(y,qc) cstr_dae(y,qc);
>> yss = fsolve(f,[9 310 300],[],qc)
yss = 8.2931 314.5318 301.8150
 Integrate model equations from steady state

>> M = [1 0 0; 0 1 0; 0 0 0];
>> df = @(t,y,qc) cstr_dae(y,qc);
>> [t1,y1]=ode15s(df,[0,50],yss,[],qc)
>> ax1=plotyy(t1,y1(:,1),t1,y1(:,2));
>> xlabel('Time [h]');
>> ylabel(ax1(1),'Concentration [kmol/m^{3}]');
>> ylabel(ax1(2),'Temperature [K]');
                                            Matlab Results for qc = 1

                           8.2931                                                           314.5318
Concentration [kmol/m 3]




                                                                                                       Temperature [K]
                           8.2931                                                           314.5318




                           8.2931                                                           314.5318
                                    0   5    10   15   20      25      30   35   40   45   50
                                                            Time [h]
           CSTR Model Solution cont.
   Integrate model equations for different qc
>> qc=0.5;
>> [t2,y2]=ode15s(df,[0,50],yss,[],qc)
>> ax2=plotyy(t2,y2(:,1),t2,y2(:,2));
>> xlabel('Time [h]');
>> ylabel(ax2(1),'Concentration [kmol/m^{3}]');
>> ylabel(ax2(2),'Temperature [K]');
                                        Matlab Results for qc = 0.5
                           10                                                           400




                           8                                                            380
Concentration [kmol/m 3]




                                                                                              Temperature [K]
                           6                                                            360




                           4                                                            340




                           2                                                            320




                           0                                                            300
                                0   5    10   15   20      25      30   35   40   45   50
                                                        Time [h]

								
To top