Process simulation using continuation method in complex domain

Computers

them.

Engng

Vol. 22, Suppl.. pp. S94kS946.1998

0 1998 Elsevier Science Ltd. All rights reserved

Pergamon PII: SOO98-1354(98)00186-O

F’rinted in Great Britain 0098-1354/98 $19.00 + 0.00

Process Simulation using Continuation Method in Complex Domain bY

Farhang. Jalali Chem. Eng. Dept.,, Shiraz University,Shiraz, IRAN Abstract

In recent years there has been active interest in process simulationin the complex domain. In this work a new procedure regartig fixed-point continuation method in complex domain using b&cation theory has been proposed. Based on this method, homotopy branches in real space are connected to each other through biication branches in complex space. Therefore, all of the solution branches are found by just one initial guess. There are some examplespresented to show how it can be done in solving nonlinear sets of solutions, specially in phase equilibrium problems which contain the most noulineariries in process simulations.The results are quit interesting and it can be extended to the flowsheetcalculations. 0 1998 Elsevier Science Ltd. All rights reserved. Introduction

Over the past three decades, many computer-aidedmethods for process design and control in &emical industry have evolved. Today, chemical processes are considered with an increasing number of complex models. Among the fundamental tasks in cumputer-aideddesign of chemical processes is that of solving the largescale nonlinear algebraic equation systems involved in steady-state flowsheet modeling, Zimey and Stadtherr (1988). and Seader et al. (1990). There are two reasons for this interest. The first one comes &II the broadly recognized need to automate the task of numerical modeling. The second reason is that a reliable method for solving the required nonlinear equation systems still needs to be developed. Several d.ifIlcultiesiu canverging iterative calculationshave been experienced using available software. The programming problems form a major subset of the field of mathematical programming. In particular, problems related to chemical process design and control can often be formulated as nonlinear optimization problems. These problems usually have nonconvexities that imply difficulties in determining a global solution. Chemical engineering examples of such problems are mostly in reactor ne&vork synthesis, phase equilibrium, phase and chemical reaction equilibrium, heat exchanger network design. batch processes,and optimal design of distillation. Much of the design literature is related to solving the material and energy balances, phase equilibrium, transport equations, and chemical kinetic expressions. For steady-state analysis, the resulting nonlinear equation system, i.e. (1) f(x,g)=O is solved, for the most part using Newton-basedmethods. Here, f is a vector function of vector x, and p is a vector of parameters. Nonlinear equations, however, have many roots in general, and it is highly desirable in many chemical engineering problems to find all real solutions. Methods for finding additional roots are discussed by Allgower and Georg (1980). By far, the most succes&l of these approaches, which however is limited to systems of polynomial equations and does not apply to general systemsof nonlinear equations, is the global homotopy method as discussedby Garcia and Zangwill(19791,Watson (1986), and Li (1987). On the other hand the behavior of iterative map8 of complex variables have been studied in mathematical physics for many years. From these studies many important contributions to the theory of nonlinear systems such as bifurcation theory have come to see. Furthermore, many of these concepts have found applications in fluid mechanics, chemically reacting systems, process control, and other areas of science and engineering. But, very little has had a direct impact on the development of methods of solving sets of nonlinear algebraic equations. There are very few studies that deal with multivariable maps in the complex domain. Recently, Lucia and Xu (1992) have studied the behavior of several commonly used nonlinear equation-solvingmethods in process simulation. The objective of this work is to study the behavior of continuation method using Bifurcation theory in the complex domain in the context of multivariablephase equilibrium which has the most nonlinear nature in algebraic systems. It is shown that the non-convergentbehavior in R space can be transferred to complex domain and then come back to the real space to get the real solutions s943

s944

European

Symposmm

on Computer

Aided Process

Engineering--X

Continuation Method Conthmation method is a major advance in the solution of nonlinear systems of equations. It is 3 global method that has the capability of finding all the roots to a nonlinear system. Continuation generates all of the possible solutions from which physically meaningful solutions can be selected. The importance of the continuation technique comes from the disadvantage of Newton’s method, which can generate a singularity, or near-singularity of the Jacobean matrix and, thus, prevent convergence. In highly nonlinear systems of equations, a large or even global convergence domain is often desirable. This is especially true when the algorithm is part of a large simulation system where a good initial estimate may be difficult to provide and failure to converge would be resulted. A method that allows for a greatly expanded convergence domain relative to Newton’s method is differential-arclength homotopy continuation. Waybum and Seader (1987) discuss the conditions under which homotopy continuation is guaranteed to converge, but can not guarantee to find all of the solutions. Because of the features mentioned above, the application of homotopy continuation to engineering is receiving much attention. Waybum and Seader (1987) give an excellent review of continuation methods and their application to separation problems. A system of nonlinear algebraic equations can be written as: f@=O (2) where E Rn-+ Rn, the continuously differentiable function; and x is the variable vector. The homotopy functions provide a smooth transition between an approximation to the solution (often linear or nearly linear) and the true solution. Therefore, they gradually mtroduce the complex nonlinearities. In many formulations, the homotopy function is defined as: H(x , t) = t f(x) + (1 - t) g(x) (3) where f(x) is the system of equations to be solved and g(x) is the simple system of equations for which a solution is known. Here, t is a scalar homotopy parameter, which is gradually varied Tom 0 to 1 as the path, is tracked horn the staring point to a solution. Bifurcation

Theory

Bifurcation theory deals with the analysis of branch points of nonlinear equations in Banach space. Banach space is a real or complex vector space in which a vector has a non-negative length (or norm), and in which every Cauchy sequence converges to a point of the space. The key point in bifurcation theory is the Implicit function theorem (Golubitsky and Schaeffer, 1985). By exploiting this theorem, one can determine if a set of equations has a unique solutiou. Complex computation To show the complex computation, the following example is solved using a homotopy function. Consider the following cubic equation, presented by Choi and Book (1989). This problem can not reach the solution in real space with u=2 as an initial guess using the homotopy continuation method. However, by extending the computation into complex space, all of the three solutions to the cubic equation are found. The original equation is: -u ( u - l)( u + 1) = 0 (4) Using a complex variable, the equation becomes -z ( z - l)( z + 1) = 0 (5) where z=x+iy (6) To solve the eqation, it is split into two equations; one is the real part of the equaticm and the other one is the imaginary part of the equation as follows: F (1) = x (-x2 + 39 + 1) = 0 F (2) = Y (y2 -3x2+1)=0 Using the fixed-point homotopy fuuction, H (1) = t F(1) + (1-t)( x - xo) H(2)=tF(2)+(1-t)(y-y”) The computation is started from x=2, and y=O. Figure1 shows the result.

(7) (8) (9) (10)

S945

European Symposium on Computer Aided Process Engineering-8

4326

. .. .

. . . ‘._,, ‘.._....._.. . . . ...‘. ,....’ _,._........

‘5 l8

. .

i,

6 os

?.

%.I L

.i _:’ .: ;:

-2 -3 -4 .’ -5. -5

* -4

J

-3

-2

-1

0

1

2

3

4

5

Applkation to phase equilibrium

The basis of most chemical process design problems, such as distillatiou cohmm design cx altemative separation systems, lies in the fuudamentalproblem of phase equilibrium. Given a set of feed compositions for a multi-component system at a given temperamre and pressure, the state of phases existing at expdlibri~ as well as the compositiou and quantity of each is determined. Using the aiteriou of equilibrium for a system under constraiuts of fixed pressum and temperature, the Gibbs fmxtion has to reach its minimum value. The formulationcan be writteu iu texmsof Lagrange multipliers to transfann the problem to an unconstrained form. This problem is difkeutiated with respect to the variable iu the fun&m, thm a system ofnonliuear algebraic functions is resulted. Two examples of the above algorithm as applied to phase equilibrium iu noukleal systems are presented to show the abilities of the propoxd method. The lirst example comes &II the paper by E&auk and Barrufet (1988). The mathematical model for Gibbs free energy for a binary mixture usiug Marguk activity coe!Xcieutmodel is used in which the variables cau be cousidered as complex oues. III this work, we first calculate in real space with the initial guess x1 =0.8. TheresultisshowuiuFigure2 Asambe seen &om the figure only one solution GUIbe found at t=l. In tie next step we tried the problem in both real and complex domain. The computatiouis done with the same initial guess (x=0.8, y=O.O)and all of the three solutions are traced by homotopy brauchtx iu both spaces. The result is shown in Figure 3. III the second exampk another phase equilibrium problem (derived &om S. M. Walas, 1985) is tried with the same activity coefficieut. Iu this example the problem is solved in the real danain with initial value of 0.5 and only one solution is reached by fixed-point homotopy brauch, Figure 4. The problem is tried again in the real and complex space with the same initial value (x=0.5, ~0.0) and all the three solutio&aregaineQ Figure5 _ Rgure 3.Homology Path (or Canplex space 1 ___.“.

...

.._

..

. .

-

0.2

1.5 hcmotq

2 Parameter

2.5

3

I. __I 3.!5

L..-.__I_*

-4

-3

, -2

-1 hcmci~

LL

1

0 Parameler

r’......,.,.I.

1

2

4

3

European

S946

rigure 4-l iomdopy --------I-.

0.75r’-1--

Symposium

on Computer

Path 101 Real Space

Alded Process

1 _ ~_.

----

07.

.-l

Figure 5.Homolo~~y Path lor Complex Space __~~__ ~_~__~__ I _ _,,_,._ _.

OY-

0.65 -

..,

,.... ., I..

OS-

0.6

; 0.7 -

/--

5 LxO.S-

:i

0.5

I z,,,_...

I.......

2 p 0.45

.P G f 0.4

x.0.55 s t z

Engmeering--P

m g 03

0.4 0.35

. .. .. .. .. .. . .

. . . . . . . . . . . . . . . . . . . . .._

: I

0.2

0.3 0.25’ 0

0.5

1

-I 1.5 homolopy Parameter

2

2.5

I 3

2

fi 3 homolopy Parameler

4

5

In this research, chemical process simulation were presented using homotopy continuation method in complex plane using Bifurcation branches as connecting tool between two real spaces. For the above problems, fixed-point homotopy was applied because they are so sensitive to the initial guess. The flxedpoint bomotopy can have a scaling problem, and in tbe above examples, scaling is an important factor in the computationsthat should be considered. Complexitlcation can occur if the equations related to the design can not be transferred into complex variables, of course in that case, there are some modification available to overcome the difficulties. In order to prevent a singtdarity during the computation,it is possibleto increase the precision of the iteration, so the calculation can pass through the singular points. Acknowledgment

F. Jalali would like to appreciate the University Council of Research, Shiraz university, for supporting this work under the grant awarded to the author. AIlgowcr, E., and K. Georg, “Simplicial and Continuation Methods for Approximating Fixed Points and Solutionsto Systemsof Equations,”SIAMRev., 22.28 (1980). Choi. S.H., and N.L. Book, “An EvaIuation of the Starting Point Criteria for the Global Fixed-Point Homotopy Continuation Method,” Department of Chemical Engineering, University of Missouri-Rolla (1989). Eubank P.T., and M.A. Barrufet,” A Simple Algorithm For Calculation of Phase Separation”, Chem. Eng. Edu., 3641, (winter 1988). Garcia, C.B., and W.I. Zangwiil, ‘Determining All Solutions to Certain Systems of Nonlinear Equations,” Math. Gpns. Res., 4, 1 (1979). Li, T.Y., “SolvingPolynomialSystems,”Math. Intelligence,9(3), 33 (1987). Lucia, A., and J. Xu, “ChemicalProcess Gptimization Using Newton-likeMethods,”Comput. Chem. Eng., 14,2.119 (1992). Seader, J.D., M. Kuno, W.J. Lin, S.A. Johnson, K. Unsworth, and J.W. Wiskin, “Mapped Continuation Methods for Computing all Solutions to General Systems of Nonlinear Equations,”ComputersChem. Eng., 14(l), 71 (1990). Walas. S. M.,”Phase Equilibrii in ChemicalEngineering”, Butterworth,Stoneham. MA (1985) Watson, L.T., “Engineering Applications of the Chow-York Algorithm,” Appl. Math. Compum, 9, 111 (1986). Waybum, T.L., and J.D. Seader, “HomotopyContinuation Methods for Computer-AidedProcess Design,” Comput. Chem. Eng., 11(1).7(1987). Zimey, S.E., and M.A. Stadtherr. “Computational Experiments in Equation-based Chemical Process Flowsheeting,”Comput. Chem. Eng., 12(12), 1171 (1988).

-