Phase equilibrium calculations with quasi-Newton methods

The phase split problem, formulated as an unconstrained minimization of the Gibbs free energy, is commonly solved by the second-order Newton method, preceded by a number of first-order successive substitutions. For difficult problems, the convergence radius of the Newton method may be small and a high number of successive substitution iterations may be required before the switch, or in repeated switch-backs. An interesting alternative is given by the quasi-Newton methods, representing a good compromise between complexity and convergence speed. The quasi-Newton BFGS method exhibits a super-linear convergence rate (in some cases without step length control) and a rank two update of the Hessian matrix guarantees a hereditary positive definiteness. In this work, a scaling methodology is proposed for finding the appropriate change of variables for phase equilibrium problems; applied to the two-phase split problem, the resulting change of variables leads to a Hessian matrix of the form H = I + D + ND, where I is the identity matrix, D is a diagonal matrix with elements vanishing at the solution, and ND is an effective low-rank matrix. The results of numerical experiments carried out on several test cases show that the BFGS method using the proposed variables is more robust and efficient than previous implementations (from the literature and open source codes). A two-parameter cubic equation of state was used in this work, but any equation of state can be used. The quasi-Newton methods are particularly suited for thermodynamic models for which the Hessian matrix is difficult or costly to obtain.