Equations of State - American Chemical Society

23 Convergence Behavior of Single-Stage Flash Calculations 1

Marinus P. W. Rijkers and Robert A. Heidemann

Downloaded by COLUMBIA UNIV on July 10, 2013 | http://pubs.acs.org Publication Date: March 24, 1986 | doi: 10.1021/bk-1986-0300.ch023

Department of Chemical and Petroleum Engineering, University of Calgary, Calgary, Alberta T2C 1R4, Canada

Flash calculation procedures have been investigated to determine when a flash routine (1) converges to the trivial solution or (2) converges very slowly or not at a l l . The two mixtures used as the basis for the study are a CO containing natural gas and a CO -H S-CH mixture that shows liquid-liquid separations. The Peng-Robinson equation was used to describe all the equilibrium phases. The computation method is the successive substitution procedure with four different approaches for updating the Κ values. 2

2

2

4

When a s i n g l e e q u a t i o n o f s t a t e i s used t o model a l l t h e e q u i l i b r i u m phases i n a f l a s h c a l c u l a t i o n , a t r i v i a l s o l u t i o n can be reached i n w h i c h t h e c a l c u l a t e d phases a r e i d e n t i c a l and have a l l t h e p r o p e r t i e s of t h e feed stream. The n e c e s s a r y c o n d i t i o n s f o r e q u i l i b r i u m , i . e . , g

±

= In f

V ±

L

- in f . = 0

;

i=l,...,N

(1)

are o b v i o u s l y s a t i s f i e d by such a s o l u t i o n . Most o f t h e p r a c t i c a l e q u a t i o n s o f s t a t e have a t most threevolume r o o t s a t any p r e s s u r e and temperature. Whenever t h r e e volume r o o t s a r e obtained f o r the feed mixture i t i s p o s s i b l e t o avoid the t r i v i a l s o l u t i o n s i n c e a " l i q u i d " f e e d and "vapor" f e e d c a n be assigned different volumes at the flash conditions and t h e f u g a c i t i e s i n t h e l i q u i d - l i k e and v a p o r - l i k e phases w i l l s u r e l y be different. The problems w i t h t r i v i a l s o l u t i o n s can occur o n l y a t Τ and Ρ c o n d i t i o n s where t h e e q u a t i o n o f s t a t e f o r t h e f e e d m i x t u r e has o n l y one p o s i t i v e r e a l volume r o o t . There have been s e v e r a l s u g g e s t i o n s made f o r o b t a i n i n g from an e q u a t i o n o f s t a t e a " l i q u i d - l i k e " o r " v a p o r - l i k e " volume, as needed, i n r e g i o n s where t h e r e i s o n l y one volume r o o t ( 1 - 4 ) . These i d e a s are reported t o be h e l p f u l i n avoiding t r i v i a l solutions i n 7

Current address: Technische Hogeschool Delft, The Netherlands. 0097-6156/ 86/ 0300-0476506.00/ 0 © 1986 American Chemical Society

In Equations of State; Chao, K., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1986.


23.

Single-Stage Flash Calculations

RIJKERS A N D H E I D E M A N N

All

s i m u l a t i o n s where many e q u i l i b r i u m c a l c u l a t i o n s a r e performed i n sequence. An a l t e r n a t i v e approach i n d e a l i n g w i t h s u c c e s s i v e e q u i l i b r i u m computations has been to use converged s o l u t i o n s a t one o r more p o i n t s t o i n i t i a t e a new c o m p u t a t i o n ( 5 - 7 ) . I t i s assumed t h a t i n i t i a t i n g near the s o l u t i o n has the advantage of r e d u c i n g the number of i t e r a t i o n s r e q u i r e d t o r e a c h convergence and of i n c r e a s i n g the p r o b a b i l i t y of convergence t o the c o r r e c t , n o n - t r i v i a l , s o l u t i o n . I t has a l s o been r e p o r t e d t h a t , i n g e n e r a l , f l a s h r o u t i n e s can converge v e r y s l o w l y i n some p a r t s o f the P-T c o - e x i s t e n c e r e g i o n . M i c h e l s e n (8,9) has p r e s e n t e d an a n a l y s i s of the c o n v e n t i o n a l s u c c e s s i v e s u b s t i t u t i o n method and has shown v e r y c l e a r l y that r a p i d convergence c o u l d not be expected a t c o n d i t i o n s c l o s e t o the c r i t i c a l p o i n t of the m i x t u r e b e i n g f l a s h e d . In t h i s p a p e r , a s t u d y has been performed t o determine the e f f e c t o f the i n i t i a t i o n p r o c e d u r e on whether o r not the t r i v i a l s o l u t i o n i s reached. Also, several alternative flash calculation a l g o r i t h m s have been examined t o see what convergence b e h a v i o r can be expected i n a l l p a r t s of the phase diagram. Initiation The c a l c u l a t i o n s performed here a r e a l l based on use of Κ v a l u e s . F o r s u b s t a n c e i , p r e s e n t i n the l i q u i d and vapor phases w i t h mole f r a c t i o n s x^, and y^, r e s p e c t i v e l y , the phase d i s t r i b u t i o n c o e f f i c i e n t s a r e d e f i n e d by (2)

K. = y./x. The K. a r e e v e n t u a l l y t o be found from the Peng-Robinson state (10).

e q u a t i o n of

1

The

trivial

s o l u t i o n i s c h a r a c t e r i z e d by K_ = l f o r a l l i .

t i a t i n g computations w i t h

Ini

v a l u e s near u n i t y would appear t o r i s k

convergence t o the t r i v i a l s o l u t i o n . A w i d e l y used p r o c e d u r e i s t o i n i t i a t e from R a o u l t ' s l a w , K^=P^ /P where P_^ i s the vapor p r e s s u r e of s u b s t a n c e i a t the e q u i l i b r i u m t e m p e r a t u r e . A p r o p o s a l a s c r i b e d t o G.M. W i l s o n p e r m i t s e s t i m a t i o n of t h e s e R a o u l t ' s Law Κ v a l u e s from the pure component c r i t i c a l p r o p e r t i e s and a c e n t r i c f a c t o r s , i . e . S

S

Κ:. = ^

in

in(? C

/P) + 5.37 i

(1 + ω.)(1-Τ 1

/Τ)

(3)

i

I t would be i m p o s s i b l e t o examine a l l p o s s i b l e ways f o r i n i t i a t i n g the K^. We have chosen t o v a r y the i n i t i a l Κ v a l u e s from u n i t y through the W i l s o n v a l u e s a c c o r d i n g t o

in

Κ

= α

where α i s a s c a l a r .

£n The o v e r l i n e s a r e meant t o show t h a t in Κ and


478

in

EQUATIONS OF STATE: THEORIES AND APPLICATIONS

are vector q u a n t i t i e s .

The Κ v a l u e s o b t a i n e d from (4) a r e f

u n i t y when α = 0 and a r e the W i l s o n R a o u l t s law v a l u e s when α= 1. We have found i t u s e f u l t o v a r y α i n t h e i n t e r v a l 0 < α £ 1.5 and t o observe t h e converged f l a s h c a l c u l a t i o n s t o see f o r what v a l u e s o f α the t r i v i a l s o l u t i o n i s o b t a i n e d . Flash

Calculations

G i v e n t h e Κ v a l u e s i t i s p o s s i b l e t o compute how a g i v e n f e e d m i x t u r e s e p a r a t e s i n t o two phases. The mass b a l a n c e e q u a t i o n s a r e o r g a n i z e d as d e s c r i b e d , f o r i n s t a n c e , by N u l l ( 1 1 ) , i . e . , f o r one mole o f feed w i t h mole f r a c t i o n s z^, (1-ν)χ + V K x Downloaded by COLUMBIA UNIV on July 10, 2013 | http://pubs.acs.org Publication Date: March 24, 1986 | doi: 10.1021/bk-1986-0300.ch023

±

±

= z

±

(5)

±9

hence,

x

« z /[l + (K -1)V]

(6)

and

y

- Κ ζ / [ l + (K -1)V]

(7)

±

±

The l i q u i d and vapor mole f r a c t i o n s must sum t o u n i t y , o r h(v) - I (y.-χ,.) - I ( K i ) r

Z j

/[i+(K i)v] r

=o

(8)

F i n a l l y , e q u a t i o n (8) i s an e q u a t i o n i n one unknown, V, which can be s o l v e d by t h e Newton-Raphson p r o c e d u r e . The r e s u l t i n g v a l u e o f V, when s u b s t i t u t e d i n t o e q u a t i o n s (6) and ( 7 ) , g i v e s t h e e q u i l i b r i u m phase c o m p o s i t i o n s . The b e h a v i o r o f e q u a t i o n (8) has been a n a l y z e d by Nghiem e t a l . (12). They p o i n t o u t t h a t h(V) = 0 has one and o n l y one m e a n i n g f u l r o o t , 0 < V < 1, i f and o n l y i f h(0) £ 0 and h ( l ) < 0. I f e i t h e r o f these i n e q u a l i t i e s i s v i o l a t e d , t h e m i x t u r e w i l l be i n o n l y one phase. In p a r t i c u l a r , i f h ( 0 ) < 0, then V=0 and the m i x t u r e i s a l l i n the l i q u i d phase w i t h mole f r a c t i o n s x^ = ζ . The c a l c u l a t e d mole f r a c t i o n s i n the vapor phase, y^ = K^x^, w i l l n o t sum t o u n i t y . I f f u r t h e r c a l c u l a t i o n s a r e r e q u i r e d , i t i s necessary t o normalize; i.e. to set y, = K.x./ I K.x. i i ι J jJ S i m i l a r l y , i f h ( l ) > 0, then V = l , y mole f r a c t i o n s a r e found from x

i

=

(9) = z^ and m e a n i n g f u l

7

> 0 c o r r e s p o n d s t o the r e g i o n where T


the e q u a t i o n of s t a t e has t h r e e volume r o o t s . The boundary curve of the r e g i o n has been c o n s t r u c t e d by f i x i n g the volume, c a l c u l a t i n g the temperature t h a t makes (3Ρ/3ν)^ = 0 and then c a l c u l a t i n g the p r e s s u r e from the e q u a t i o n s of s t a t e e v a l u a t e d a t the volume and temperature. The cusp i s the " m e c h a n i c a l c r i t i c a l p o i n t " f o r the m i x t u r e and i s w e l l i n s i d e the two phase r e g i o n . The d i f f u s i o n a l s t a b i l i t y l i m i t , d e t ( B ) = 0, was e v a l u a t e d i n the same manner as the l i m i t of m e c h a n i c a l s t a b i l i t y and as d e s c r i b e d by Heidemann and K h a l i l ( 1 9 ) . The c r i t i c a l p o i n t of the m i x t u r e l i e s on t h i s c u r v e . For the n a t u r a l gas m i x t u r e , the two b o u n d a r i e s f o r the unstable region enclose o n l y a s m a l l p a r t of the two-phase r e g i o n . The u n s t a b l e r e g i o n f o r the CH^ - l ^ S - CO^ m i x t u r e i n c l u d e s most of the phase diagram. I n s i d e these u n s t a b l e r e g i o n s s o l u t i o n i s impossible with successive s u b s t i t u t i o n .

the

trivial

I n i t i a t i o n and the T r i v i a l S o l u t i o n For the two m i x t u r e s of F i g u r e s 1 and 2 a l a r g e number of f l a s h c a l c u l a t i o n s have been performed. For the m i x t u r e of F i g u r e 1, the two-phase r e g i o n was i n t e r s e c t e d a t p r e s s u r e s of 20, 65 and 70 bar and f l a s h c a l c u l a t i o n s were performed w i t h v a r y i n g i n i t i a l Κ v a l u e s over the whole temperature range w i t h i n the two-phase r e g i o n . For the m i x t u r e of F i g u r e 2, e q u i v a l e n t c a l c u l a t i o n s were performed at 40 and 100 b a r . F i g u r e s 3 through 6 and F i g u r e s 7 and 8 c o n t a i n the r e s u l t s f o r the two m i x t u r e s a t the v a r i o u s p r e s s u r e s . The i n i t i a l Κ v a l u e s were varied, as described earlier, by letting the "Phase D i s t r i b u t i o n C o e f f i c i e n t I n i t i a t i o n Parameter", a, take on v a l u e s between z e r o and 1.5. The v a l u e of α i s shown on the o r d i n a t e of the F i g u r e s and the temperature i s shown on the a b s c i s s a . The a r e a s of the Τ - α p l a n e are l a b e l l e d t o show the s o l u t i o n reached. For the f i r s t m i x t u r e a t 20 b a r , the two-phase r e g i o n extends from 170 to 254 K. Between 170 and 194 Κ the m i x t u r e i s e i t h e r mechanically or d i f f u s i o n a l l y unstable or both. As e x p e c t e d , the c o n v e n t i o n a l s u c c e s s i v e s u b s t i t u t i o n a l g o r i t h m d i d not converge t o the t r i v i a l s o l u t i o n w i t h i n t h i s temperature i n t e r v a l f o r any α v a l u e . T h i s i s i n d i c a t e d by the c o r r e s p o n d i n g c r o s s - h a t c h e d r e g i o n of F i g u r e 3. A l s o shown i n F i g u r e 3 i s the r e g i o n of low α v a l u e s , o u t s i d e the u n s t a b l e temperature i n t e r v a l , where COSS converged t o the trivial solution. To a v o i d the t r i v i a l s o l u t i o n i t i s only n e c e s s a r y t o have α l a r g e enough. Near the phase boundary, 254 K, l a r g e r α values are r e q u i r e d . However, no t r i v i a l s o l u t i o n s were found a t any temperature so l o n g as α was g r e a t e r than about 0.45. I n p a r t i c u l a r , i f the W i l s o n Κ v a l u e s are used (a = 1 ) , the t r i v i a l s o l u t i o n would n e v e r be reached a t 20 bar and any temperature.



484

EQUATIONS OF

STATE: T H E O R I E S A N D

APPLICATIONS

F i g u r e 4 shows convergence b e h a v i o r w i t h COSS a p p l i e d t o the f i r s t m i x t u r e a t 65 b a r . The two-phase r e g i o n i s c o n t a i n e d between 215 and 248 K. (The p r e s s u r e i s above the " m e c h a n i c a l c r i t i c a l p r e s s u r e " , hence the m i x t u r e i s not m e c h a n i c a l l y u n s t a b l e a t any temperature.) Figure 4 i n d i c a t e s that t r i v i a l s o l u t i o n s were reached a t temperatures between the phase boundary and the s t a b i l i t y l i m i t when the Κ v a l u e s were i n i t i a t e d too near u n i t y (a below about 0.6). At temperatures g r e a t e r than 218 K, the t r i v i a l s o l u t i o n was avoided f o r α > 0.35. C e r t a i n l y , no trivial s o l u t i o n s were encountered when COSS was i n i t i a t e d a t the W i l s o n Κ v a l u e s (a = 1 ) . A l t h o u g h convergence t o the t r i v i a l s o l u t i o n i s i m p o s s i b l e i n the u n s t a b l e r e g i o n , the number of i t e r a t i o n s r e q u i r e d c o u l d become large. In F i g u r e 4, open boxes are used to i n d i c a t e p o i n t s where convergence was not a c h i e v e d . Convergence f a i l u r e i s d i s c u s s e d i n more d e t a i l l a t e r . The c a l c u l a t i o n s f o r m i x t u r e I a t 65 bar were repeated u s i n g the f i r s t a c c e l e r a t e d s u c c e s s i v e s u b s t i t u t i o n p r o c e d u r e (ACSS1). The convergence b e h a v i o r i s shown i n F i g u r e 5. I t w i l l be seen t h a t ACSSl and COSS behave s i m i l a r l y as r e g a r d s the tendency t o converge to the t r i v i a l s o l u t i o n . S i n c e b o t h p r o c e d u r e s can be regarded as steep descent a l g o r i t h m s , n e i t h e r can reach the t r i v i a l s o l u t i o n i n the u n s t a b l e temperature i n t e r v a l . A s i d e from some i s o l a t e d p o i n t s , ACSSl and COSS converge t o the c o r r e c t s o l u t i o n a t the same v a l u e s of a. F i g u r e 6 shows the r e s u l t s of ACSSl f l a s h c a l c u l a t i o n s f o r m i x t u r e I at 70 b a r . The i s o b a r at t h i s p r e s s u r e does not i n t e r s e c t the u n s t a b l e r e g i o n . P o t e n t i a l l y , t h e r e f o r e , the t r i v i a l s o l u t i o n c o u l d be reached a t any temperature w i t h i n the two-phase r e g i o n between about 221 and 244 K. The f i g u r e shows, however, t h a t i f α i s l a r g e enough, ( i . e . , α > 0.6) the t r i v i a l s o l u t i o n i s never r e a c h e d . The r e s u l t s are of the same type as seen i n F i g u r e s 3-5 and demonstrate t h a t i n i t i a t i n g w i t h the W i l s o n Κ values i s sufficient to avoid t r i v i a l s o l u t i o n s even a l o n g a difficult retrograde isobar. Convergence r e s u l t s f o r the m i x t u r e of F i g u r e 2 a t 40 and 100 b a r , u s i n g COSS, are shown i n F i g u r e s 7 and 8. The i s o b a r at 40 bar extends from a low-temperature l i q u i d - l i q u i d r e g i o n , through a v a p o r - l i q u i d r e g i o n , t o the phase boundary a t about 287 K. The m i x t u r e i s t h e r m o d y n a m i c a l l y u n s t a b l e a t temperatures below 259 K, t h e r e f o r e i t i s i m p o s s i b l e f o r COSS to converge t o the trivial s o l u t i o n at any temperature below 259 K, r e g a r d l e s s of the i n i t i a l Κ v a l u e s used. The p l o t i n F i g u r e 7 c o n t a i n s o n l y a v e r y s m a l l r e g i o n (which i s not c r o s s - h a t c h e d ) bounded by 259 < Τ < 287 Κ and α < 0.3 i n w h i c h t r i v i a l s o l u t i o n s were found. As shown i n F i g u r e 7, two d i f f e r e n t two-phase s o l u t i o n s were p o s s i b l e i n the temperature i n t e r v a l 188 < Τ < 202 Κ. Depending on the i n i t i a l Κ v a l u e s , i t was p o s s i b l e t o r e a c h e i t h e r a v a p o r - l i q u i d solution or a l i q u i d - l i q u i d solution. There i s o b v i o u s l y a t h r e e - p h a s e r e g i o n i n t h i s v i c i n i t y and i t s e x i s t e n c e c o m p l i c a t e s f l a s h c a l c u l a t i o n s when the number of phases s p e c i f i e d i s l e s s than the number of phases w h i c h s h o u l d a c t u a l l y e x i s t a t e q u i l i b r i u m ( a c c o r d i n g t o the model). No attempt was made t o l o c a t e the three-phase r e g i o n . I t s h o u l d be n o t e d , however, t h a t convergence




!

!

160.0

180.0

F i g u r e 3.


I

200.0

I

220.0

TEMPERATURE (K)

485

1

»

240.0

260.0

Region o f S o l u t i o n f o r M i x t u r e I a t 20 Bar

TEMPERATURE (K) F i g u r e 4.

Region o f S o l u t i o n f o r M i x t u r e I a t 65 Bar Using COSS


EQUATIONS O F STATE: THEORIES A N D APPLICATIONS


486

210.0

215.0

220.0

225.0

230.0

235.0

240.0

245.0

250.0

255.0

TEMPERATURE (K) F i g u r e 5.

Region o f S o l u t i o n f o r M i x t u r e I a t 65 Bar U s i n g


ACSS1


23.


Τ 180.0

F i g u r e 7.

ι 180.0

F i g u r e 8.

I 200.0

487


I 220.0

I I 240.0 260.0 TEMPERATURE (K)

I 280.0

Region o f S o l u t i o n f o r M i x t u r e II

ι

ι

200.0

220.0

ι

1

240.0 260.0 TEMPERATURE (K)

I 300.0

r 320.0

a t 40 Bar

1

I

»

280.0

300.0

320.0

Region o f S o l u t i o n f o r M i x t u r e II

a t 100 Bar


488



to a n o n - t r i v i a l s o l u t i o n was o b t a i n e d a t a l l temperatures and a t a l l α v a l u e s i n t h e v i c i n i t y o f t h e t h r e e phase r e g i o n . F i g u r e 8 shows convergence p r o p e r t i e s o f t h e m i x t u r e o f F i g u r e 2 a t 100 b a r . The p r e s s u r e i s j u s t below t h e d i p i n t h e phase boundary curve a t about 235 Κ so t h a t t h e m i x t u r e remains i n two phases from v e r y low temperatures t o about 309 K. There i s a temperature i n t e r v a l from about 228 t o 249 Κ where t h e m i x t u r e i s outside the l i m i t of d i f f u s i o n a l s t a b i l i t y . L i k e w i s e , from about 297 t o 309 Κ t h e m i x t u r e i s o u t s i d e t h e l i m i t o f s t a b i l i t y . Only i n these two s e p a r a t e i n t e r v a l s a r e t r i v i a l s o l u t i o n s p o s s i b l e when flash c a l c u l a t i o n s a r e performed using COSS. These trivial s o l u t i o n s a r e avoided i f α i s g r e a t e r than about 0.45, as shown i n F i g u r e 8. A l s o i n d i c a t e d i n F i g u r e 8 a r e a few p o i n t s where convergence was n o t a c h i e v e d . These p o i n t s a r e a l l i n t h e i n t e r m e d i a t e s t a b l e temperature i n t e r v a l w i t h v e r y l o w α v a l u e s . Slow Convergence Problems In some r e g i o n s o f t h e phase diagrams t h e s u c c e s s i v e s u b s t i t u t i o n a l g o r i t h m s converged v e r y s l o w l y . I n o r d e r t o reduce computation time a d u a l a l g o r i t h m was used. The computations were begun w i t h the s u c c e s s i v e s u b s t i t u t i o n p r o c e d u r e , b u t i f convergence was n o t a c h i e v e d w i t h i n a r e a s o n a b l e number o f i t e r a t i o n s t h e Newton-Raphson a l g o r i t h m f o r u p d a t i n g t h e Κ v a l u e s was begun. T h i s procedure d i d n o t change t h e c o n c l u s i o n s as t o whether o r not t h e t r i v i a l s o l u t i o n was u l t i m a t e l y t o be reached so l o n g as t h e f l a s h c a l c u l a t i o n was w e l l on i t s way t o convergence and so l o n g as t h e r e was t o be a s p l i t i n t o two phases o f f i n i t e amount. The second c o n d i t i o n i s n e c e s s a r y because t h e J a c o b i a n m a t r i x of e q u a t i o n (26) i s u n d e f i n e d when V = 0 o r V = 1. The problem a r i ses w i t h m a t r i x IJ * g i v e n by e q u a t i o n (22) . The d u a l a l g o r i t h m improved convergence b e h a v i o r c o n s i d e r a b l y . In p r e p a r i n g F i g u r e 4, 3283 independent f l a s h c a l c u l a t i o n s were required. W i t h COSS a l o n e , 238 o f these c a l c u l a t i o n s d i d n o t converge w i t h i n 100 i t e r a t i o n s . W i t h t h e d u a l a l g o r i t h m o n l y 27 poorly i n i t i a t e d c a l c u l a t i o n s were l e f t unconverged. Figure 7 r e q u i r e d 3179 f l a s h c a l c u l a t i o n s , 382 o f w h i c h d i d n o t converge w i t h COSS a l o n e . A l l c a l c u l a t i o n s converged w i t h t h e d u a l a l g o r i t h m . The reason f o r convergence f a i l u r e even o f t h e Newton-Raphson p r o c e d u r e r e q u i r e s some d i s c u s s i o n . What was a c h i e v e d i n these cases was an o s c i l l a t o r y b e h a v i o r between a s o l u t i o n l i k e t h e t r i v i a l s o l u t i o n and one l i k e t h e c o r r e c t s o l u t i o n . These f a i l u r e s were g e n e r a l l y a s s o c i a t e d w i t h i n i t i a l guesses f o r t h e Κ v a l u e s near u n i t y ( a near z e r o ) . Speed o f Convergence The speeds o f convergence o f f o u r d i f f e r e n t a l g o r i t h m s f o r u p d a t i n g the Κ v a l u e s a r e compared i n F i g u r e s 10 through 13. The m i x t u r e o f F i g u r e 1 i s t h e b a s i s f o r t h e comparison. F l a s h c a l c u l a t i o n s were performed a t 930 p o i n t s throughout t h e two-phase r e g i o n f o r t h i s m i x t u r e and a t each p o i n t t h e c a l c u l a t i o n was i n i t i a t e d w i t h t h e W i l s o n Κ v a l u e s (a = 1 ) .



23.

RUKERS AND HEIDEMANN

Single-Stage

Flash

Calculations

489

F i g u r e 9 s h o w s t h e COSS r e s u l t s . A t e a c h o f t h e 930 p o i n t s t h e number o f iterations to convergence has been inserted in the diagram. Through most o f the two-phase region, convergence was o b t a i n e d i n 25 i t e r a t i o n s o r f e w e r a n d t h e a v e r a g e i t e r a t i o n c o u n t t o c o n v e r g e n c e was j u s t o v e r 1 5 . A t o n l y t h r e e p o i n t s o f t h e 930 was c o n v e r g e n c e u n a t t a i n a b l e i n 100 i t e r a t i o n s a n d t h e s e p o i n t s a r e very c l o s e to the c r i t i c a l p o i n t of the m i x t u r e . F i g u r e s 10 a n d 11 show i t e r a t i o n c o u n t s f o r t h e A C S S 1 a n d A C S S 2 algorithms, respectively. These two algorithms perform quite c o m p a r a b l y , w i t h a v e r a g e i t e r a t i o n s t o convergence o f 9.84 and 9 . 8 5 . ACSS1 is apparently superior to ACSS2 on the grounds that convergence was o b t a i n e d w i t h i n 100 i t e r a t i o n s a t a l l 9 3 0 p o i n t s w i t h ACSS1 and ACSS2 f a i l e d t o c o n v e r g e a t one p o i n t . T h e t o t a l e x e c u t i o n t i m e a t 9 3 0 p o i n t s was 1 8 1 . 8 , 1 4 1 . 6 , a n d 1 6 6 . 3 s e c o n d s f o r COSS, ACSS1 and A C S S 2 , r e s p e c t i v e l y . On t h i s a c c o u n t A C S S 1 a p p e a r s a l s o t o be s u p e r i o r t o A C S S 2 a n d t o s u c c e s s i v e s u b s t i t u t i o n w i t h o u t a c c e l e r a t i o n (COSS). F i g u r e 12 shows t h e r e s u l t s o b t a i n e d w h e n t h e Newton-Raphson a l g o r i t h m was employed to update Κ values. As b e f o r e , these calculations were a l l initiated with the Wilson Κ values. A s t r i k i n g r e s u l t i s t h a t t h e a v e r a g e number o f i t e r a t i o n s p e r c o r r e c t s o l u t i o n i s o n l y 4 . 7 1 , l e s s than h a l f the average of the other three procedures. However, the iteration count g i v e s only a partial p i c t u r e o f t h e e f f i c i e n c y o f t h e a l g o r i t h m a n d two n e g a t i v e factors must be c o n s i d e r e d . The f i r s t n e g a t i v e f e a t u r e i s t h a t t h e N-R u p d a t i n g p r o c e d u r e f o r t h e Κ v a l u e s f a i l e d t o r e a c h t h e c o r r e c t s o l u t i o n a t 40 o f 9 3 0 points. T h i s p o o r c o n v e r g e n c e b e h a v i o r was a r e s u l t o f a n u n d e f i n e d J a c o b i a n at the phase b o u n d a r i e s . I f e q u a t i o n s ( 6 ) a n d (7) d i d n o t y i e l d a two p h a s e s o l u t i o n , s m a l l q u a n t i t i e s o f L o r V w e r e i n t r o duced. C o n s e q u e n t l y the x^ and y^ were e v a l u a t e d from e q u a t i o n s (6) o r (7) a n d w e r e n o r m a l i z e d b y a p p l y i n g e q u a t i o n s ( 9 ) o r ( 1 0 ) . This approach c l e a r l y d i d not solve the problem: At 25 p o i n t s the t r i v i a l s o l u t i o n w a s r e a c h e d a n d a t t h e r e m a i n i n g 15 p o i n t s the procedure d i d not converge. The second point to consider is that each N-R iteration r e q u i r e s c o n s i d e r a b l e m a t r i x m a n i p u l a t i o n and i s much more c o n s u m i n g of computer time than a s u c c e s s i v e substitution iteration. The total time required for the 930 p o i n t s o f Figure 12 i s 424.3 s e c o n d s ; m o r e t h a n t w i c e t h e t i m e r e q u i r e d f o r t h e COSS p r o c e d u r e . The f a i l u r e t o c o n v e r g e i s no d o u b t a more s e r i o u s c o n c e r n t h a n t h e computer t i m e consumed i n t h e N-R p r o c e d u r e . Computing time c e r t a i n l y depends on t h e s k i l l o f t h e programmer and no d o u b t a N-R r o u t i n e c a n be w r i t t e n w h i c h w i l l be more e f f i c i e n t t h a n o u r s . Discussion All examples presented here support the conclusion that the successive s u b s t i t u t i o n a l g o r i t h m cannot converge to the trivial s o l u t i o n , h o w e v e r i n i t i a t e d , a t p o i n t s where t h e homogeneous phase is unstable. This property i s due t o t h e n a t u r e o f successive s u b s t i t u t i o n as a steep d e s c e n t p r o c e d u r e f o r m i n i m i z i n g the G i b b s free energy.


EQUATIONS O F STATE: THEORIES A N D APPLICATIONS

490

80.0

70.0 4

Cl C2 C3 nC4 r»C5 nC6

-

0.8608 0.0247 0.0067 0.0045 0.0024 0.0009

ALGORITHM USED : COSS MO. Of CALC. EXEC: 930 NO. Of CALC CONV.: 927 N O . or com. soot: 9 2 7 TOT. NO. ITERS. : 14.705 m. ΙΤ/CORftSOLM. : 15.5394 EXEC TIC s : «1787 COEXISTENCE

50.0

152 4336 32 28252321 20«17 I «42363229262321 2 0 « 17 « 16 V «35322926242220» « 17 « IS V

CURVE

DET(B) = 0 (dP/âv) = o

4 Ο

CALCULATION

NOT

J 21 21 2 0 2 0 » « 1 7 » 19 β Η 13 13 12 Ο η H 1 JO» «««««17«««Μ131312121111111 1««««««17«««1413131212111111»11 «17T7171717«1S15H1313121211 H » » » 10» c

CONVERGED

si-!*

40.0

17171

*»»»*»»»»***»»

7 1 5 1 3 » » « « « 1 5 1 4 13 13121211 11 « 10 « r « M M « « « « « M 1 3 l 3 1 2 1 2 1 1 11 « « » B«13H«««HH13131211 H H « « »~ ---- - JH13««eHH13l312l11111»» ' ' f- à Ζ Ί ï * » « « « » » »

30.0

i

s

I

H

1 4

1 4

β

1 5

H

H

fl

, ï 5 î î > 2 ^ Η Η β β β ΐ 4 « β β ι ι η β β » 9 5 5 1 t ï « e i l 11 » » « 9 9 1 » À5 S 2 5 £A * « « « H 1 3 1 3 1 2 Î 1 Î 1 » » 999881 S 17 « « « 14 M « / » 17 W W « 1 3 l 3 l 2 t t » » » 999888 f~21 2 0 » « 517 14 « *2021£Ϊ* « 8 8 98 97 87 87 87 8 7 7 6 7 6 ^ 202\ « V T «71«3H11 2 2ΤΙ1 «2 1101 » 9 8 8 7 7 7 7 7 6 6 6 6 10222221 2 0 » « 17 « « 2 2 2 2 « * 14 1312 11 Ο 9 9 7 7 7 7 6 6 6 6 6 6 6 24232220» « « Τ7 » 2 3 2 3 2 2 · β 43Ί2 ΤΙ « 9 9 8 8 7 6 6 6 6 6 6 6 5 5 5 Ι 2 6 2 3 2 1 2 0 « « « « « 2 4 2 9 2 4 » « 1 3 1 % if » 9 9 8 8 7 7 6 6 6 6 5 5 5 5 5 2321 » « « « « » 21 32322217 14 12 H » Α / 9 8 8 7 7 7 7 6 g ^7 7 6 6 6 6 6 6 5 5 5 5 5 5 5 5 17 17 » » 2031 442317 13 12 « τ 240.0 180.0 200.0 220.0 160.0


ΰ

20.0

1

M

e

e

e

e

M

U

1 3

1 7 1 8 1 3

10.0

4

9

8

8

7

1 7

#

260.0

TEMPERATURE (Κ) Figure Wilson

9. Iteration Κ Values

Counts

f o r Mixture

I.

COSS

Initiated

with

80.0

70.0 Η

cokrosmoN C02 - O.ttOO CI - 0.8606 C2 - 0.0247 C3 - 0.0067 nC4 - 0.0045 nC5 - 0.0024 nC8 - 0.0009

ALGORITHM USED : ACSS1 NO. Or CALC EXEC: 930 NO. Or CALC CONV.: 930 NO. ( Τ COfM. SOUL: 930 TOT. NO. ITERS. : 9.«2 Mf. IT/C0RR.SOLR : 9.8409 EXEC TIC · : 141596 . COEXISTENCE

50.0

4

• (dP/dv) 40.0

30.0

20.0

0.0

-I

CURVE

• DET(B)=0 =0

. n a n e n ctn n β » 9 β β β 9tt»121t131311»» 9 8 8 88 .. _ 8 1 1 1 3 1 2 1 3 1 3 1 1 1 1 » » 9 8 8 8 8 « 1 2 » 913H121211»»» 988888 *2»111112»1111»tt 99888887 » « 911 H » tt tt » 9 8 8 8 8 8 8 7 7 _ . 911 tf JO » 9 1 1 1 1 1 1 1 1 1 1 » » 9 8 8 8 8 8 7 7 7 99 9 H J D / 9 » » H H 1 2 1 I » » » 9 8 8 8 8 7 7 7 7 9 9»»/M0»»H1211tttttt 9 9 8 8 8 7 7 7 7 7 , ~ . . . » » » .V· 9 » 11 12 11 » « » 9 9 8 8 8 7 7 7 7 6 6 "9 9H tt » 9 9/fTi 1 t H 1 2 » 1 1 » t t 9 9 8 8 8 7 7 7 7 6 6 6 9 9 9«/l» 11 11 12 H 11 tt tt 9 9 8 8 7 7 7 7 6 6 6 6 6 I 9 911 W « H « « « β » 9 9 8 8 7 7 7 6 6 6 6 6 6 6 8 » 9 9 911 13 flAl 12 11 « « " 9 9 8 7 7 7 6 6 6 6 6 6 6 9*5" ^ ^811 » » 9 » 9 911 13 14/13/12 » 9 » 9 8 7 7 7 7 7 6 6 6 6 6 6 9 9 9 9 7 7 7 6 6 6 6 6 6 6 5 5 S S 5 S 6 l l » « « 1 1 1 2 » 912 14 14 I f f 11 9 9 9 9 7 6 6 6 6 6 6 6 9 9 9 9 9 9 " «1312121211 H » 1 2 « « M 1 2 21 H1 W I 1 7 7 7 7 6 6 6 6 6 6 9 9 9 9 9 9 9 9 6 6 6 6 9 9 9 9 9 9 9 9 ' " 1 2 M 1 3 1 3 1 3 » 1 4 12 WT7M1211 ~ 6 9 9 9 9 9 9 4 4 4 H 1 3 1 2 » H It II 1 2 H » « M 1 2 » 8 8 8 1/7 7 6 6 6 6 6 n n n n n « e i 3 i i 9 9 8 8 7 6 6 6/6/« 6 6 5 9 5 9 9 4 4 4 4 4 4 I 240.0 160.0 180.0 200.0 220.0

260.0

TEMPERATURE (K) Figure Wilson

10. I t e r a t i o n Κ Values

Counts

f o r Mixture

I.

ACSS1

Initiated


with

23.


COMPOS

70.0

J

-

Flash

491

Calculations

ALCORfTHM USED : ACSS2 MO. OF CALC. EXEC: 930 NO OF CALC. CONV.. 929

mow

C 0 2 - o.ttoo C I - αββΟβ

C2 C3 nC4 nC3 nC6

Single-Stage

0 0247 C 0067 30045 0 0024 00009

NO. UT CORR. SOLK.: 929

TOT. NO. ITERS. : 9.255 M. IT/CORK.SOLK : 9.8547 EXEC. T I C > : «6.336

COEXISTENCE

CURVE

• DET(B)-0

{dP/dv)--o n 1012 9 1 0 » 12 10tttttt 12 9 9 9 9

40.OJ tu H « « jf « βΗ η « 9yn m ' 9 1 4 9«Î1 tt/« 9 ' 8 8tt η τιtt» » tt ^ 9 9 9 « T1tt« Λ » 8 99 9 T i n « 9 « 8 9 913 910 « 9 « « '910 9 9 913 12 910 9 « « /9A3 « « -^8 913 13 12 8W 9 9 9T1/iZtt 9 9 „ 8 913 1 3 Π 8 9 9 9 912 1 2 V 9 9 9 ^9 810 12 12 tt 8 9 910 13 13 12tt/ 9 9 8 8 913 U 10 10 9 9 910 913\ Hn 12 1210 tt«fett y « 99888 8 ilO 9tt 10 8 β 9 913 13 10 ' ' 9 8 8 7 I ΠH tttt tt'tt 9 9TC tt 8 8 8tt14Ul3T1T1tttt • 10 10 A / 8 7 7 6 6 5 9 9T0 1 3 « 1 4 T 1 1 2 « 1 0 9 9 88 7 7/β/β 6 6 5 5

20.0


tttt

r

0.0

160 Figure Wilson

11. I t e r a t i o n Κ Values

COMPOSITION

C02 - 0.1000 C1 - 0.8606 C2 - 0.0247 C3 - 0.0067 nC4 - 0.0045 nC5 - 0.0024 nC8 — 0.0009

180.0

200.0 TEMPERATURE (K)

Counts

f o r Mixture

220.0 I.

ACSS2

240.0 Initiated

260.0 with

ALGORITHM USED : NR NO. Or CALC. EXEC: 930 NO. Or CALC CONV.: 9B NO. Or CORK. SOLK: 890 TOT. NO. ITERS. : 4,661 m. IT/CORR.SOLN. : 4.7135

EXEC TME ·

COEXISTENCE

: 424.301

CURVE

DET(B)=0

(dp/dv) = o 1—Ι — l

1

CALCULATION CONVERGED

NOT

C>3 C A L C . G I V E S T R I V I A L

SOL'N

Λ 5 Χ 5>f

\x\w\fi\ 44 4 4 -4 4 4 4 4 4 4 55554455 44445557

\\\m\\ S 5 5 S 5fil5 5 5 5 5 5 5 5^5 5 3 5 6 6 5 5 f i 555 5 5 5 5 Afi5 5 4 5 5 5 4/4/• • • • 180.0 200.0 TEMPERATURE

Figure

12. I t e r a t i o n

with Wilson

Counts f o r M i x t u r e

220.0

260.0

240.0

(K) I.Newton-Raphson

Initiated

Κ Values



492


Resorting to Newton-Raphson updating of Κ values, however, amounts to giving up the minimization character of the algorithm in favor of using an equation solving technique. The N-R scheme can (it has been demonstrated) converge to the trivial solution even at unstable points. In the examples it appears that convergence to trivial solutions is unlikely if the Κ values are initiated appropriately. For the mixture of Figure 1, initiating with the Wilson Κ values was adequate to avoid convergence to the trivial solution anywhere inside the two-phase region. This is true even though the region of instability is a very small part of the two-phase region. Initiating ^ith Κ values near unity can result in convergence to trivial solutions. This implies that there may be some risks involved in using Κ values at a converged solution to initiate a flash calculation at a nearby point, particularly near the critical. The results regarding speed of convergence indicate a definite advantage in accelerating the successive substitution algorithm. The first procedure, ACSS1, appears sufficient. The Newton-Raphson procedure used here certainly reduces the iteration count dramatically but has associated with it an increased risk that the trivial solution or non-convergence will be reached. Also, it is not certain that a reduced number of iterations implies a reduction in computer time. Other types of flash calculations than the constant Τ and Ρ flash examined in this manuscript are important. The occurrence of trivial solutions in isenthalpic flash routines, for instance, requires further study. Acknowledgments This work was supported by the National Science and Engineering Research Council of Canada. Literature Cited 1.

Poling, B.E.; Grens, E.A.; Prausnitz, J.M. Ind. Engg. Chem. Proc. Des. Dev. 1981, 20, 127. 2. Gundersen, T. Comp. Chem. Engg 1982, 6, 245. 3. Coward, I.; Gayle, S.E.; Webb, D.R. Trans. I. Chem. E. 1978, 56, 19. 4. Veeranna, D.; Rihani, D.N. Fluid Phase Equil 1984, 16, 41. 5. Asselineau, L . ; Bogdanic, G.; Vidal, J . Fluid Phase Equil 1979, 3, 273. 6. Michelsen, M.L. Fluid Phase Equil 1980, 4, 1. 7. Mehra, R.K.; Heidemann, R.A.; Aziz, K. Soc. Pet. Engrs J 1982, 22, 61. 8. Michelsen, M.L. Fluid Phase Equil 1982, 9, 1. 9. Michelsen, M.L. Fluid Phase Equil 1982, 9, 21. 10. Peng, D-Y.; Robinson, D.B. Ind. Engg. Chem. Fundamentals 1976, 15, 59. 11. Null, H.R. "Phase Equilibrium in Process Design"; Wiley Interscience Inc.: New York, 1970.


23.

12. 13. 14. 15. 16. 17. 18.


19.

RIJKERS AND HEIDEMANN

Single-St age Flash Calculations

493

Nghiem, L.X.; Aziz, K.; Li, Y.K. Soc. Pet. Engrs J1 1983, 23, 521. Anderson, T . F . ; Prausnitz, J.M. Ind. Engg. Chem. Proc. Des. Dev. 1980, 19, 9. Mehra, R.K; Heidemann, R.A.; Aziz, K. Can. J . Chem. Eng. 1983, 61, 590. Nghiem, L.X.; Heidemann, R.A. In 2nd European Symposium on Enhanced Oil Recovery 1982, p.303 Ed. Techniqs., Paris, 1982. Zeleznik, F . J . J . Assoc. Comp. Mach. 1968, 15, 265. Rijkers, M.P.W. M.Sc. Thesis, The University of Calgary, Calgary, Alberta, Canada, 1985. Amundson, N.R. "Mathematical Methods in Chemical Engineering"; Prentice-Hall, Inc.: Englewood Cliffs, N . J . , 1966. Heidemann, R.A.; Khalil, A.M. AIChE J . 1980, 26, 769.

RECEIVED November 8, 1985


Equations of State - American Chemical Society

Recommend Documents