Phase stability analysis is a crucial step in the determination of multiphase equilibrium. This analysis by the tangent plane distance (TPD) minimization is a well-known technique, as well as the difficulties in providing guarantees that the global minimum has been found. On this regard, interval methods are powerful tools since they provide such guarantees. In this work, an interval Newton method plus generalized bisection, based on a modified affine arithmetic, is used to reliably find all possible stationary points of the TPD function. Additionally, an improved convergence test is suggested as well as a special treatment for mole fraction weighted averages. Several mixtures with up to 5 components, including LLE island type ternary systems, were studied. Both activity coefficient models and cubic equations of state were considered. For all the cases tested, the proposed modified affine arithmetic method was superior to other interval-based methods. (C) 2013 Elsevier Ltd. All rights reserved.