- Email: [email protected]
Sliding Motion of Discontinuous Dynamical Systems Described by Differential Algebraic Equations
IFAC
Copyright «:> IFAC Fault Detection, Supervision and Safety of Technical Processes, Washington, D.C., USA, 2003
c
0
C>
Publications www.elsevier.comllocatelifac
SLIDING MOTION OF DISCONTINUOUS DYNAMICAL SYSTEMS DESCRIBED BY DIFFERENTIAL ALGEBRAIC EQUATIONS Jyoti Agrawal· Kannan M. Moudgalya .,1 Amiya K. Pani··
• Department of Chemical Engineering •• Department of Mathematics Indian Institute of Technology Bombay, India
Abstract: Starting from the equations describing the equivalent dynamics of ordinary differential equations (ODEs), an equivalent dynamical differential algebraic equation (DAE) system is derived. It is shown that a procedure similar to the one available for ODEs can be used to construct an equivalent DAE system as well. The proposed method is shown to work with a class of discontinuous DAE systems typified by an ideal gas-liquid system and a soft drink manufacturing process. Copyright © 2003 !FAC Keywords: Sliding mode, Equivalent dynamics, Hybrid system, Singular system, Multiple forcing functions, Discontinuous system, Soft drink process.
1. INTRODUCTION
reasonably well understood phenomenon for systems made up of ODEs. Although discontinuous systems made of DA& also are known to possess this property (Moudgalya and Jaguste 2001) , there exists no mathematical treatment of these systems.
Simulation is a popular method to assess the safety of processes. This approach becomes even more important in difficult processes, such as those, for example, described by DA& and discontinuities.
Most of the time, DAEs are converted into ODEs for analysis and numerical solution. Brenan et al. (1996) argue why it is better to treat the DAEs directly. Vora et al. (2002) demonstrate with an example that it is possible for a simple controller designed using a DAE control theory to work better than an "advanced" controller designed through the ODE route. Moreover, at times it may be difficult to convert the DAEs into 00&.
Discontinuous dynamical systems occur reasonably frequently in real life, see for example, the works of Mosterman et al. (1997), Moudgalya and Jaguste (2001) and Agrawal et al. (2002). These authors articulate the need to understand such systems and to solve them efficiently. Dynamical systems are generally modelled by DA&. In view of this discussion, it is not surprising that discontinuous DAE systems occur reasonably often.
This study proposes a method for treating discontinuous DAEs in the DAE framework itself. This is illustrated with a simple gas-liquid system (Moudgalya and Ryali 2(01) and a more complicated soft drink process.
Discontinuous systems often exhibit the property of sliding, a state, in which the resultant of the forcing functions is used to drive the system. This is also known as equivalent dynamics. This is a 1 Corresponding author. Fr rnail : [email protected]
735
r~ [FL-kLZ(;;$t -Po.,)]'
2. AN IDEAL GAS-LIQUID SYSTEM IN SLIDING MOTION Consider the ideal gas-liquid system studied by Moudgalya and Ryali (2001), with two phases in the tank and one outlet (see Fig. 1).
(10)
with the switching junction tp(Mc , Md > 0 and when the gas comes out, it is
Ga~
Liquid
(11)
Gas uid
Li'1
Gas Liquid -------
with tp(Mc, ML ) < is given by
--+-----i
o. The switching function
(12)
Fig. 1. A simple discontinuous system
It is shown by Moudgalya and Ryali (2001) that this system reaches the sliding mode, characterised by the following equivalent dynamics (Filippov 1988),
In the left figure, as the outlet tube dips into it, the liquid comes out. As the gas is sent in but not removed, the pressure builds up and hence the liquid outflow increases with an eventual dropping of liquid level in the tank. When the level falls below the dip tube, the gas sta.rts coming out. If on the other hand, the liquid level is below the dip tube, as in the right figure, gas only flows out; the liquid level rises and eventually reaches the dip tube. There are two DAE models that describe these two situations.
:t [~~]
=
10,
(13)
10 =
or + (1 - o:)r, 0 ~
0=
('\/I{J,j-) (Vtp, 1-) - (vC{', 1+) ,
jo
= [c 2 - ~IMc] ,
Cl
= V _ Vd
C2
= (:;
kcxR:T
(1) (2)
~ 1, (14)
(15)
+
> 0,
(16)
~Lx + Pout) kcx.
They have shown that the behaviour of the equivalent system (16) is in exact agreement with the numerical solution of discontinuous system, (9) .
(3)
(4) (gas model)
dMc
--=F.c- G dt ' dML -F
--;:u- -
0
with the following specific values,
(liquid model)
Else If &PL < Vd
tp
3. TOWARDS EQUIVALENT DYNAMICS FOR DAES IN SLIDING MOTION
(5)
(6)
L,
ML) ,
IlJ'+ O=V- ( M c P PL 0= G - kcx(P - Pout).
Although the DAEs can be converted into ODEs in the case of a simple system, such as the one discussed above, this may not be possible in a general case. An example in case is the High Density Polyethylene reactor (Moudgalya and Jaguste 2001) described by two DAE models consisting of 15 and 12 equations. This problem could become worse if more chemical components are present. In situations like these, we may not be able to convert the DAEs into ODEs by simple substitution. For the reasons discussed in Section 1 also, we may want to handle the DAEs directly. These reasons have motivated us to consider the construction of equivalent dynamics described by DAEs.
(7) (8)
The system keeps switching between these two models. The symbols are the same as the ones used by Agrawal et al. (2002) . Through substitution, a system of ODEs with two forcing functions is arrived at .
d [Mc] dt ML =j,
(9)
where, j takes one of the two expressions j+ or
1- · When the liquid flows out of the system, the
We begin with the equivalent dynamic description available from the ODE approach . Substituting
forcing function is
736
for j+ and f-, respectively from (10) and (11) and making use of (14), (13) becomes dM
c -----;It
= o.Fc
+ (1
(17)
The mathematical formulation of this procedure is straight forward : (1) Form an equivalent dynamical system by combining ODEs only. (2) Augment the above system with algebraie constraints. In the next section, we present a Filippov like procedure to combine the entire DAE sets.
(18)
4. EQUIVALENT DAES IN SLIDING
- o.)Fc
_ (1 _ o.)kc x
(Mc~
_ Pout) ,
V-~
PL
L
dM -dt
= o.FL + (1
KT o.kLX (Mc M - Pout ) V -
7ft
- o.)FL .
We would like to find out whether the Filippov's method of multiplying the entire set of ODEs by 0. can be extended to the entire set of DAEs as well. That is, we would like to multiply each equation in (1) to (4) by 0. and add to (1 0.) times the corresponding equation in the gas model. One can easily check that this procedure will give rise to (24) to (26) . But there is a problem with the equations for outflow rates, (4) and (8) : for liquid outflow rate denoted by (4), there is no corresponding equation in the gas model. Similarly, there is no equation for gas outflow in the liquid model.
From (3) or (7) we get,
P= McRT
(19)
V-ML.' PL
and from (4), (8) we get, L = kLX(P - Pout), G = kcx(P - Pout),
(20) (21)
using which, (17) and (18) become dM
c = o.Fc + (1 -----;It dM
-----;ItL
= o.FL
-
- o.)Fc - (1 - o.)G,
o.L + (1 - o.)FL ·
(22)
Thus we can extend Filippov's procedure to the entire DAE set only if we have equations that describe the chosen phenomenon in both the models. This motivates us to define the missing equations. As no liquid flows out of the system during the gas model, one may be tempted to use L = 0 as the liquid flow equation in the gas model. Unfortunately, such a combination would produce an equation different from (4) and hence will violate the derivation of equations (24) - (28) : recall that to derive this equivalent dynamics, we have used (4).
(23)
While the substitutions indicated by (20) and (21) occur only once, that indicated by (19) has been used twice. This observation will be used in the subsequent discussion . We rewrite the above equations and order them in a way suitable for further discussion: dMc (24) ( i t = o.Fc + (1 - o.)(Fc - G), dML
(it
= o.(FL -
L)
+ (1 -
o.)FL'
KT+ML) O=V- ( M c P PL 0= L - kLx(P - Poud, 0= G - kcx(P - Poud.
,
(25)
Then we ask the question, "With what equation can we combine (4) so that it comes out unchanged?" The answer is obvious: combine it with an equation of the form 0 = O. Similarly, (8) can be combined with the equation 0 = O.
(26) (27)
(28)
This set of DAEs describes the equivalent dynamics of the gas-liquid system in sliding motion. It is identical to the one obtained by Moudgalya and Ryali (2001). This can be verified by straight forward substitution of algebraic variables in differential equations.
We present the mathematical formulation of this second approach now. Suppose the following discontinuous DAEs are in a sliding state.
The above set of equations can also be derived directly from (1) - (8) , as described next . If we add 0. times the equation set (1) - (2) with (I - 0.) times the set (5) - (6), we get the differential equation set (24) - (25) . We augment this with all the algebraic constraints present in (1) - (8). Note that the volume constraint given by (3) and (7) indicates the same equation and hence needs to be included only once. This corresponds to (19) being used twice in the substitution process that has helped arrive at (22) and (23).
where, j, yE !Rn; g, z E !Rm and
dy
dt =
fey , z),
(29)
0= g(y, z),
(30)
cp(y, z) > 0, cp(y, z) < 0, f{J(Y, z)
> 0,
cp(y, z) < O.
(31) (32)
With the usual vector notation (33)
737
we get the following equivalent dynamics:
while the model is to be integrated do if the acid model is applicable then repeat if gas model was used previously then re-initialize the acid model end if integrate the acid model until the acid model is invalid else repeat if acid model was used before then re-initialize the gas model end if integrate the gas model until the gas model is invalid end if end while
(34)
[;] = Fo where,
Fo
= o:F+ + (1 -
0=
o:)P-, 0:::; 0: :::; 1 (35) ("VC{', (1- , g-)) (\lC{', (j-,g-)) _ ("VC{', (j+,g+)) . (36)
One can see this procedure to be identical to the one proposed by Filippov for ODEs. Note that this method of combining the equations is not generic because it requires domain knowledge. For example, we have to combine the liquid flow equation with the corresponding liquid flow equation only. Unfortunately, it is not clear how to overcome this requirement . It is noteworthy, however, that Filippov's approach also has this restriction implicitly. For example, we combine the two equations corresponding to dMc/dt only - we do not combine the equation corresponding to dMc/dt in one model with that corresponding to dML/dt in the other.
Fig. 2. Pseudo-code to integrate soft drink process Table 2. Soft drink process parameters Parameter PL pa V Vd
T Po",
5. APPLICATION
x kL kG
In order to test this equivalent DAE dynamic approach, we have modelled a realistic DAE system, possibly applicable in soft drink manufacturing. It is well known that soft drinks have CO 2 and sugar syrup. Let CO 2 and syrup be fed into a mixing tank and withdrawn through a single outlet. As in the ideal gas liquid system, the model that is to be used to describe the gas outflow is referred to as the gas model, while the corresponding model for acid outflow is referred to as the acid model, both consisting of DAEs.
Ft F2 Ite (T
Value 50 mol/I 16mol/I 10 I 2.25 I 293 K 1 atm 1.0 2.51/atm/s 3.0I/atm/s 5.5 mol/s 7.5 mol/s 0.433/400 I/mole/s 1640 atm/mole frac
Table 3. Initial conditions and DASSL parameters Variable Model
Value
MJ M2
0.72 mol 95.0 mol 0.0 mol
M3
The gas and the acid models have different number of variables. This is because the acid phase consists of H2COa , H 20 and dissolved CO 2 while the gas phase has only CO2 . Acid model is described by 3 ODEs and 6 algebraic equations and gas model is described by 3 ODEs and 4 algebraic equations. Each model of the system is shown in Table 1. The solvability of the acid and the gas models is tested using the incidence matrix analysis proposed by Moudgalya (2001). Both the models are found to be index-l DAEs. The discontinuous system is integrated using the state of the art DAE integrator, DASSL (Brenan et al. 1996), using the algorithm in Fig. 2. The parameters in Table 2 are used in this simulation. In Table 2, initial conditions (left) and DASSL parameters (right) are given. The results are plotted in Fig. 3 and 4 with dotted lines. Sliding motion is observed in this simulation.
Gas
Parameter rtol atol ~t
Value 10-' 10- 4 10- 4
We next solve this system using the proposed equivalent dynamic approach. One can check that these DAEs cannot be reduced to ODEs through substitution. The missing equations in columns 2 and 3 of Table 1 are represented by an equation of the form 0 = o. The equivalent DAE dynamic approach is then applied to this system: find 0: times the equations in column 2 and 1 - 0: times the corresponding equations in column 3 to arrive at the following equations: dM,
-;It
= o:(G - LdF, - G - r,
dM2 -dt = -o:TLJ7.~ + F2 dM
-dt-a = -o:La + r '
738
r,
(37) (38)
(39)
Table 1. Discontinuous DAE model for CO2 mlxmg in soft-drink process with rp= & PL
+ &P.
- Vd·
Model
Acid Model (F
C02 balance
aMI = FI _ LI - r d:U-2 - - = F2 - L2 - r d'l:l3 -;it =-L3+ r 0= Ll + L2 + L3 - kLX(P - P ou,) 0M2 _ ( L2 ) LI + L2 + L3 - Ml + ~ + M3 03 _ ( L3 ) LI + L2 + L3 - Ml + M2 + M3
H20 balance H2C03 balance Acid outflow
H20 in outllow H2C03 in outflow Gas outflow Gas-liquid equilibrium
C02 in the reactor Volume constraint
0=
M2
MI
+ M2 + M3
_ (
0= M +~: +
dt 0=0
0=0 0=0 0= G - kex(P - P ou,)
0= P _
UMI Ml +M2+M3 0= Ml - (Mg + Md (MIRT M2 M~) O=V- -p+;:-+-;:-
+ L2 + L3
+~:
M3 - (Ll + 0= Ll + L2 + L3 - kLx(P - Pout)' 0= G - kax(P - Pout)' 0= Ml - (Mg + M/), 0- p_ (FM/ M/+M2+M3'
UMI
M/+M2+M3 0= MJ - (Mg + Md O=V- (MIRT + M2 P
+ M3)
PT.
p,
), (40)
LJ '
2.5
(41) M, (mol) U
(42) (43) (44)
0.:)
(45)
... .. . . . .. . .. .. .
0 0
lOO
and
0= V _
d:U- 1 --=F2 -r d'lr3 --=r
0=0 0= P _
L2
Ll
Gas Model (F- , '" < 0) aMI = FI -G _ r
, '" > 0)
..
ID
30
30
10
30
30
50
eo
50
eo
go
(M1KI' + M2 + M3) , P PL Pa
80
(46)
70 M, (mol)
where, Q can be found on substitution of F+ and F- from Table 1 in (36) as
+
F2 Pa r(PL - Pa) Q= ~~-~~-~
eo 50 40 30
(47)
20
L 2 Pa + L 3 PL The equivalent DAE system is also of index 1. It is integrated using DASSL with initial conditions obtained at the beginning of the sliding motion, which is found during simulation of discontinuous system. These results are plotted with solid lines in Fig. 3 and 4: in the former, the total holdup of the three components inside the reactor (C0 2 , H 2 0, H2C03) are plotted; in the latter, from the top, fraction of the gas dissolved in liquid, acid concentration and C02 concentration, respectively, in the liquid phase are plotted.
10 0
40
90
80 10
eo M,(mol)
50
.a 30 30
10 0
0
Fig. 3. Dotted line: discontinuous system. Solid line: equivalent dynamics.
6. CONCLUSIONS
Simulation of a discontinuous DAE system is expensive because of the frequent initialisations and the associated startup problems. Moreover, one j·s also constrained to take small steps because of difficulties, such as, discontinuity sticking (Moudgalya and Ryali 2001) . The equivalent dynamic formulation overcomes these difficulties: the proposed method is several orders of magnitude more efficient than the discontinuous method .
It is shown that an equivalent dynamic description consisting of DAEs can be obtained for a class of discontinuous DAEs. Two methods have been proposed towards this end . The first method combines only the ODEs and augments them with algebraic constraints. The second method is more Filippov like in the sense that the algebraic equations also participate in the combining step.
739
0 .11
r-----,-----.---,.-----,,-----,..---,
0.1
M/ P
0.09
Pout
0.08
R
0.07 M,(mol) 0.06
r
0.0.'>
T
0.04
~t
0.03
V Vd x
0.02 L _ - - '_ _-'-_ _....i-_ _L-_-:'::--_-:: o 10 20 30 40 50 60
Moles of CO2 in liquid phase mole Pressure atm Outlet pressure atm Gas constant I atm/mole/K H 2 C03 formation rate = Kc M'vM , mole/I/s Temperature in the reactor K integration interval s Volume of the reactor I Volume below the dip tube I Valve opening
. - - - - - . - - - . - - - - . -... ------r... . .... ... \ --r - - - - _'
:
. .... ... .
50
"------
Greek symbols .'
40 %H,CO, 30 20
PL Pa
10
Kc
OL----'lO---20L--~30---40L--~50--~60
Henry's constant for CO2 atm/mole frac Molar density of water moles/I Molar density of acid moles/I Rate constant I/moljs Switch function
200~----.,.---_.----.-----,---_,--___, 18
REFERENCES
16
.. .>
Agrawal, J ., K. M. Moudgalya and A. K. Pani (2002). An efficient integration algorithm for a class of discontinuous dynamical systems in sliding motion. In: Proceedings of American Control Conference. Anchorage, Alaska. pp. 698-703. Brenan, K. E., S. C. Campbell and L. R. Petzold (1996) . Numerical Solution of Initial- Value Problems in Differential Algebraic Equations. SIAM. Philadelphia. Filippov, A. F. (1988). Differential Equations with Discontinuous Righthand Sides. Mathematics and Its Applications (Soviet Series). Kluwer Academic Publishers. Dordrecht. Mosterman, P. J., F. Zhao and G . Biswas (1997) . Sliding mode model semantics and simulation for hybrid systems. In: Hybrid Systems. pp. 218-237. Moudgalya, K. M. and S. Jaguste (2001) . A class of discontinuous dynamical systems H. An industrial slurry high density polyethylene reactor. Chemical Engineering Science 56, 3611362l. Moudgalya, K. M. and V. Ryali (2001) . A class of discontinuous dynamical systems I. An ideal gas - liquid system. Chemical Engineering Science 56, 3595-3609. Moudgalya, Kannan M. (2001) . A class of discontinuous dynamical systems Ill. Degrees of freedom analysis. Chemical Engineering Science 56(18) , 5443-5447. Vora, P., K. M. Moudgalya and A. K. Pani (2002). Control of a high index DAE system through a linear control law. In: Proceedings of American Control Conference. Anchorage, Alaska. pp. 465-470.
12
%CO,
10
6
Fig. 4. Dotted line: discontinuous system. Solid line: equivalent dynamics. NOMENCLATURE Roman Symbols
F± Fl
F2 Fe FL f± G
g± kL ke L
LI
L2 L3
Me MI
M2 M3 Mg
Augmented forcing functions for DAEs Flow-rate of CO2 in feed mole/s Flow-rate of water in feed mole/s Gas feed rate, mol/s Liquid feed rate, mol/s Forcing function Gas outflow rate, moljs and subscript to indicate gas Vector of algebraic constraints Valve coefficient for gas flow Valve coefficient for liquid flow Liquid outflow rate, mol/s and subscript to indicate liquid CO2 outflow in acid model kmole/hr Water outflow in acid model mole/s H2 C03 outflow in acid model mole/s Gas holdup in the system, mol Moles of CO2 in reactor mole Moles of water in reactor mole Moles of H 2 C03 in reactor mole Moles of CO2 in gas phase mole
740
Copyright «:> IFAC Fault Detection, Supervision and Safety of Technical Processes, Washington, D.C., USA, 2003
c
0
C>
Publications www.elsevier.comllocatelifac
SLIDING MOTION OF DISCONTINUOUS DYNAMICAL SYSTEMS DESCRIBED BY DIFFERENTIAL ALGEBRAIC EQUATIONS Jyoti Agrawal· Kannan M. Moudgalya .,1 Amiya K. Pani··
• Department of Chemical Engineering •• Department of Mathematics Indian Institute of Technology Bombay, India
Abstract: Starting from the equations describing the equivalent dynamics of ordinary differential equations (ODEs), an equivalent dynamical differential algebraic equation (DAE) system is derived. It is shown that a procedure similar to the one available for ODEs can be used to construct an equivalent DAE system as well. The proposed method is shown to work with a class of discontinuous DAE systems typified by an ideal gas-liquid system and a soft drink manufacturing process. Copyright © 2003 !FAC Keywords: Sliding mode, Equivalent dynamics, Hybrid system, Singular system, Multiple forcing functions, Discontinuous system, Soft drink process.
1. INTRODUCTION
reasonably well understood phenomenon for systems made up of ODEs. Although discontinuous systems made of DA& also are known to possess this property (Moudgalya and Jaguste 2001) , there exists no mathematical treatment of these systems.
Simulation is a popular method to assess the safety of processes. This approach becomes even more important in difficult processes, such as those, for example, described by DA& and discontinuities.
Most of the time, DAEs are converted into ODEs for analysis and numerical solution. Brenan et al. (1996) argue why it is better to treat the DAEs directly. Vora et al. (2002) demonstrate with an example that it is possible for a simple controller designed using a DAE control theory to work better than an "advanced" controller designed through the ODE route. Moreover, at times it may be difficult to convert the DAEs into 00&.
Discontinuous dynamical systems occur reasonably frequently in real life, see for example, the works of Mosterman et al. (1997), Moudgalya and Jaguste (2001) and Agrawal et al. (2002). These authors articulate the need to understand such systems and to solve them efficiently. Dynamical systems are generally modelled by DA&. In view of this discussion, it is not surprising that discontinuous DAE systems occur reasonably often.
This study proposes a method for treating discontinuous DAEs in the DAE framework itself. This is illustrated with a simple gas-liquid system (Moudgalya and Ryali 2(01) and a more complicated soft drink process.
Discontinuous systems often exhibit the property of sliding, a state, in which the resultant of the forcing functions is used to drive the system. This is also known as equivalent dynamics. This is a 1 Corresponding author. Fr rnail : [email protected]
735
r~ [FL-kLZ(;;$t -Po.,)]'
2. AN IDEAL GAS-LIQUID SYSTEM IN SLIDING MOTION Consider the ideal gas-liquid system studied by Moudgalya and Ryali (2001), with two phases in the tank and one outlet (see Fig. 1).
(10)
with the switching junction tp(Mc , Md > 0 and when the gas comes out, it is
Ga~
Liquid
(11)
Gas uid
Li'1
Gas Liquid -------
with tp(Mc, ML ) < is given by
--+-----i
o. The switching function
(12)
Fig. 1. A simple discontinuous system
It is shown by Moudgalya and Ryali (2001) that this system reaches the sliding mode, characterised by the following equivalent dynamics (Filippov 1988),
In the left figure, as the outlet tube dips into it, the liquid comes out. As the gas is sent in but not removed, the pressure builds up and hence the liquid outflow increases with an eventual dropping of liquid level in the tank. When the level falls below the dip tube, the gas sta.rts coming out. If on the other hand, the liquid level is below the dip tube, as in the right figure, gas only flows out; the liquid level rises and eventually reaches the dip tube. There are two DAE models that describe these two situations.
:t [~~]
=
10,
(13)
10 =
or + (1 - o:)r, 0 ~
0=
('\/I{J,j-) (Vtp, 1-) - (vC{', 1+) ,
jo
= [c 2 - ~IMc] ,
Cl
= V _ Vd
C2
= (:;
kcxR:T
(1) (2)
~ 1, (14)
(15)
+
> 0,
(16)
~Lx + Pout) kcx.
They have shown that the behaviour of the equivalent system (16) is in exact agreement with the numerical solution of discontinuous system, (9) .
(3)
(4) (gas model)
dMc
--=F.c- G dt ' dML -F
--;:u- -
0
with the following specific values,
(liquid model)
Else If &PL < Vd
tp
3. TOWARDS EQUIVALENT DYNAMICS FOR DAES IN SLIDING MOTION
(5)
(6)
L,
ML) ,
IlJ'+ O=V- ( M c P PL 0= G - kcx(P - Pout).
Although the DAEs can be converted into ODEs in the case of a simple system, such as the one discussed above, this may not be possible in a general case. An example in case is the High Density Polyethylene reactor (Moudgalya and Jaguste 2001) described by two DAE models consisting of 15 and 12 equations. This problem could become worse if more chemical components are present. In situations like these, we may not be able to convert the DAEs into ODEs by simple substitution. For the reasons discussed in Section 1 also, we may want to handle the DAEs directly. These reasons have motivated us to consider the construction of equivalent dynamics described by DAEs.
(7) (8)
The system keeps switching between these two models. The symbols are the same as the ones used by Agrawal et al. (2002) . Through substitution, a system of ODEs with two forcing functions is arrived at .
d [Mc] dt ML =j,
(9)
where, j takes one of the two expressions j+ or
1- · When the liquid flows out of the system, the
We begin with the equivalent dynamic description available from the ODE approach . Substituting
forcing function is
736
for j+ and f-, respectively from (10) and (11) and making use of (14), (13) becomes dM
c -----;It
= o.Fc
+ (1
(17)
The mathematical formulation of this procedure is straight forward : (1) Form an equivalent dynamical system by combining ODEs only. (2) Augment the above system with algebraie constraints. In the next section, we present a Filippov like procedure to combine the entire DAE sets.
(18)
4. EQUIVALENT DAES IN SLIDING
- o.)Fc
_ (1 _ o.)kc x
(Mc~
_ Pout) ,
V-~
PL
L
dM -dt
= o.FL + (1
KT o.kLX (Mc M - Pout ) V -
7ft
- o.)FL .
We would like to find out whether the Filippov's method of multiplying the entire set of ODEs by 0. can be extended to the entire set of DAEs as well. That is, we would like to multiply each equation in (1) to (4) by 0. and add to (1 0.) times the corresponding equation in the gas model. One can easily check that this procedure will give rise to (24) to (26) . But there is a problem with the equations for outflow rates, (4) and (8) : for liquid outflow rate denoted by (4), there is no corresponding equation in the gas model. Similarly, there is no equation for gas outflow in the liquid model.
From (3) or (7) we get,
P= McRT
(19)
V-ML.' PL
and from (4), (8) we get, L = kLX(P - Pout), G = kcx(P - Pout),
(20) (21)
using which, (17) and (18) become dM
c = o.Fc + (1 -----;It dM
-----;ItL
= o.FL
-
- o.)Fc - (1 - o.)G,
o.L + (1 - o.)FL ·
(22)
Thus we can extend Filippov's procedure to the entire DAE set only if we have equations that describe the chosen phenomenon in both the models. This motivates us to define the missing equations. As no liquid flows out of the system during the gas model, one may be tempted to use L = 0 as the liquid flow equation in the gas model. Unfortunately, such a combination would produce an equation different from (4) and hence will violate the derivation of equations (24) - (28) : recall that to derive this equivalent dynamics, we have used (4).
(23)
While the substitutions indicated by (20) and (21) occur only once, that indicated by (19) has been used twice. This observation will be used in the subsequent discussion . We rewrite the above equations and order them in a way suitable for further discussion: dMc (24) ( i t = o.Fc + (1 - o.)(Fc - G), dML
(it
= o.(FL -
L)
+ (1 -
o.)FL'
KT+ML) O=V- ( M c P PL 0= L - kLx(P - Poud, 0= G - kcx(P - Poud.
,
(25)
Then we ask the question, "With what equation can we combine (4) so that it comes out unchanged?" The answer is obvious: combine it with an equation of the form 0 = O. Similarly, (8) can be combined with the equation 0 = O.
(26) (27)
(28)
This set of DAEs describes the equivalent dynamics of the gas-liquid system in sliding motion. It is identical to the one obtained by Moudgalya and Ryali (2001). This can be verified by straight forward substitution of algebraic variables in differential equations.
We present the mathematical formulation of this second approach now. Suppose the following discontinuous DAEs are in a sliding state.
The above set of equations can also be derived directly from (1) - (8) , as described next . If we add 0. times the equation set (1) - (2) with (I - 0.) times the set (5) - (6), we get the differential equation set (24) - (25) . We augment this with all the algebraic constraints present in (1) - (8). Note that the volume constraint given by (3) and (7) indicates the same equation and hence needs to be included only once. This corresponds to (19) being used twice in the substitution process that has helped arrive at (22) and (23).
where, j, yE !Rn; g, z E !Rm and
dy
dt =
fey , z),
(29)
0= g(y, z),
(30)
cp(y, z) > 0, cp(y, z) < 0, f{J(Y, z)
> 0,
cp(y, z) < O.
(31) (32)
With the usual vector notation (33)
737
we get the following equivalent dynamics:
while the model is to be integrated do if the acid model is applicable then repeat if gas model was used previously then re-initialize the acid model end if integrate the acid model until the acid model is invalid else repeat if acid model was used before then re-initialize the gas model end if integrate the gas model until the gas model is invalid end if end while
(34)
[;] = Fo where,
Fo
= o:F+ + (1 -
0=
o:)P-, 0:::; 0: :::; 1 (35) ("VC{', (1- , g-)) (\lC{', (j-,g-)) _ ("VC{', (j+,g+)) . (36)
One can see this procedure to be identical to the one proposed by Filippov for ODEs. Note that this method of combining the equations is not generic because it requires domain knowledge. For example, we have to combine the liquid flow equation with the corresponding liquid flow equation only. Unfortunately, it is not clear how to overcome this requirement . It is noteworthy, however, that Filippov's approach also has this restriction implicitly. For example, we combine the two equations corresponding to dMc/dt only - we do not combine the equation corresponding to dMc/dt in one model with that corresponding to dML/dt in the other.
Fig. 2. Pseudo-code to integrate soft drink process Table 2. Soft drink process parameters Parameter PL pa V Vd
T Po",
5. APPLICATION
x kL kG
In order to test this equivalent DAE dynamic approach, we have modelled a realistic DAE system, possibly applicable in soft drink manufacturing. It is well known that soft drinks have CO 2 and sugar syrup. Let CO 2 and syrup be fed into a mixing tank and withdrawn through a single outlet. As in the ideal gas liquid system, the model that is to be used to describe the gas outflow is referred to as the gas model, while the corresponding model for acid outflow is referred to as the acid model, both consisting of DAEs.
Ft F2 Ite (T
Value 50 mol/I 16mol/I 10 I 2.25 I 293 K 1 atm 1.0 2.51/atm/s 3.0I/atm/s 5.5 mol/s 7.5 mol/s 0.433/400 I/mole/s 1640 atm/mole frac
Table 3. Initial conditions and DASSL parameters Variable Model
Value
MJ M2
0.72 mol 95.0 mol 0.0 mol
M3
The gas and the acid models have different number of variables. This is because the acid phase consists of H2COa , H 20 and dissolved CO 2 while the gas phase has only CO2 . Acid model is described by 3 ODEs and 6 algebraic equations and gas model is described by 3 ODEs and 4 algebraic equations. Each model of the system is shown in Table 1. The solvability of the acid and the gas models is tested using the incidence matrix analysis proposed by Moudgalya (2001). Both the models are found to be index-l DAEs. The discontinuous system is integrated using the state of the art DAE integrator, DASSL (Brenan et al. 1996), using the algorithm in Fig. 2. The parameters in Table 2 are used in this simulation. In Table 2, initial conditions (left) and DASSL parameters (right) are given. The results are plotted in Fig. 3 and 4 with dotted lines. Sliding motion is observed in this simulation.
Gas
Parameter rtol atol ~t
Value 10-' 10- 4 10- 4
We next solve this system using the proposed equivalent dynamic approach. One can check that these DAEs cannot be reduced to ODEs through substitution. The missing equations in columns 2 and 3 of Table 1 are represented by an equation of the form 0 = o. The equivalent DAE dynamic approach is then applied to this system: find 0: times the equations in column 2 and 1 - 0: times the corresponding equations in column 3 to arrive at the following equations: dM,
-;It
= o:(G - LdF, - G - r,
dM2 -dt = -o:TLJ7.~ + F2 dM
-dt-a = -o:La + r '
738
r,
(37) (38)
(39)
Table 1. Discontinuous DAE model for CO2 mlxmg in soft-drink process with rp= & PL
+ &P.
- Vd·
Model
Acid Model (F
C02 balance
aMI = FI _ LI - r d:U-2 - - = F2 - L2 - r d'l:l3 -;it =-L3+ r 0= Ll + L2 + L3 - kLX(P - P ou,) 0M2 _ ( L2 ) LI + L2 + L3 - Ml + ~ + M3 03 _ ( L3 ) LI + L2 + L3 - Ml + M2 + M3
H20 balance H2C03 balance Acid outflow
H20 in outllow H2C03 in outflow Gas outflow Gas-liquid equilibrium
C02 in the reactor Volume constraint
0=
M2
MI
+ M2 + M3
_ (
0= M +~: +
dt 0=0
0=0 0=0 0= G - kex(P - P ou,)
0= P _
UMI Ml +M2+M3 0= Ml - (Mg + Md (MIRT M2 M~) O=V- -p+;:-+-;:-
+ L2 + L3
+~:
M3 - (Ll + 0= Ll + L2 + L3 - kLx(P - Pout)' 0= G - kax(P - Pout)' 0= Ml - (Mg + M/), 0- p_ (FM/ M/+M2+M3'
UMI
M/+M2+M3 0= MJ - (Mg + Md O=V- (MIRT + M2 P
+ M3)
PT.
p,
), (40)
LJ '
2.5
(41) M, (mol) U
(42) (43) (44)
0.:)
(45)
... .. . . . .. . .. .. .
0 0
lOO
and
0= V _
d:U- 1 --=F2 -r d'lr3 --=r
0=0 0= P _
L2
Ll
Gas Model (F- , '" < 0) aMI = FI -G _ r
, '" > 0)
..
ID
30
30
10
30
30
50
eo
50
eo
go
(M1KI' + M2 + M3) , P PL Pa
80
(46)
70 M, (mol)
where, Q can be found on substitution of F+ and F- from Table 1 in (36) as
+
F2 Pa r(PL - Pa) Q= ~~-~~-~
eo 50 40 30
(47)
20
L 2 Pa + L 3 PL The equivalent DAE system is also of index 1. It is integrated using DASSL with initial conditions obtained at the beginning of the sliding motion, which is found during simulation of discontinuous system. These results are plotted with solid lines in Fig. 3 and 4: in the former, the total holdup of the three components inside the reactor (C0 2 , H 2 0, H2C03) are plotted; in the latter, from the top, fraction of the gas dissolved in liquid, acid concentration and C02 concentration, respectively, in the liquid phase are plotted.
10 0
40
90
80 10
eo M,(mol)
50
.a 30 30
10 0
0
Fig. 3. Dotted line: discontinuous system. Solid line: equivalent dynamics.
6. CONCLUSIONS
Simulation of a discontinuous DAE system is expensive because of the frequent initialisations and the associated startup problems. Moreover, one j·s also constrained to take small steps because of difficulties, such as, discontinuity sticking (Moudgalya and Ryali 2001) . The equivalent dynamic formulation overcomes these difficulties: the proposed method is several orders of magnitude more efficient than the discontinuous method .
It is shown that an equivalent dynamic description consisting of DAEs can be obtained for a class of discontinuous DAEs. Two methods have been proposed towards this end . The first method combines only the ODEs and augments them with algebraic constraints. The second method is more Filippov like in the sense that the algebraic equations also participate in the combining step.
739
0 .11
r-----,-----.---,.-----,,-----,..---,
0.1
M/ P
0.09
Pout
0.08
R
0.07 M,(mol) 0.06
r
0.0.'>
T
0.04
~t
0.03
V Vd x
0.02 L _ - - '_ _-'-_ _....i-_ _L-_-:'::--_-:: o 10 20 30 40 50 60
Moles of CO2 in liquid phase mole Pressure atm Outlet pressure atm Gas constant I atm/mole/K H 2 C03 formation rate = Kc M'vM , mole/I/s Temperature in the reactor K integration interval s Volume of the reactor I Volume below the dip tube I Valve opening
. - - - - - . - - - . - - - - . -... ------r... . .... ... \ --r - - - - _'
:
. .... ... .
50
"------
Greek symbols .'
40 %H,CO, 30 20
PL Pa
10
Kc
OL----'lO---20L--~30---40L--~50--~60
Henry's constant for CO2 atm/mole frac Molar density of water moles/I Molar density of acid moles/I Rate constant I/moljs Switch function
200~----.,.---_.----.-----,---_,--___, 18
REFERENCES
16
.. .>
Agrawal, J ., K. M. Moudgalya and A. K. Pani (2002). An efficient integration algorithm for a class of discontinuous dynamical systems in sliding motion. In: Proceedings of American Control Conference. Anchorage, Alaska. pp. 698-703. Brenan, K. E., S. C. Campbell and L. R. Petzold (1996) . Numerical Solution of Initial- Value Problems in Differential Algebraic Equations. SIAM. Philadelphia. Filippov, A. F. (1988). Differential Equations with Discontinuous Righthand Sides. Mathematics and Its Applications (Soviet Series). Kluwer Academic Publishers. Dordrecht. Mosterman, P. J., F. Zhao and G . Biswas (1997) . Sliding mode model semantics and simulation for hybrid systems. In: Hybrid Systems. pp. 218-237. Moudgalya, K. M. and S. Jaguste (2001) . A class of discontinuous dynamical systems H. An industrial slurry high density polyethylene reactor. Chemical Engineering Science 56, 3611362l. Moudgalya, K. M. and V. Ryali (2001) . A class of discontinuous dynamical systems I. An ideal gas - liquid system. Chemical Engineering Science 56, 3595-3609. Moudgalya, Kannan M. (2001) . A class of discontinuous dynamical systems Ill. Degrees of freedom analysis. Chemical Engineering Science 56(18) , 5443-5447. Vora, P., K. M. Moudgalya and A. K. Pani (2002). Control of a high index DAE system through a linear control law. In: Proceedings of American Control Conference. Anchorage, Alaska. pp. 465-470.
12
%CO,
10
6
Fig. 4. Dotted line: discontinuous system. Solid line: equivalent dynamics. NOMENCLATURE Roman Symbols
F± Fl
F2 Fe FL f± G
g± kL ke L
LI
L2 L3
Me MI
M2 M3 Mg
Augmented forcing functions for DAEs Flow-rate of CO2 in feed mole/s Flow-rate of water in feed mole/s Gas feed rate, mol/s Liquid feed rate, mol/s Forcing function Gas outflow rate, moljs and subscript to indicate gas Vector of algebraic constraints Valve coefficient for gas flow Valve coefficient for liquid flow Liquid outflow rate, mol/s and subscript to indicate liquid CO2 outflow in acid model kmole/hr Water outflow in acid model mole/s H2 C03 outflow in acid model mole/s Gas holdup in the system, mol Moles of CO2 in reactor mole Moles of water in reactor mole Moles of H 2 C03 in reactor mole Moles of CO2 in gas phase mole
740