Thermodynamic phase equilibrium calculations for systems containing polymers and/or many molecular species are computationally expensive and difficult because the number of equations that need to be solved increases with the number of species present and because of the large size asymmetry that often leads to numerical instability. By taking advantage of the underlying form of recently proposed polymer equations of state (such as the statistical associating fluid theory or SAFT), modifications can be made that greatly simplify both the calculation of thermodynamic stability and phase equilibrium calculations for polydisperse polymer solutions. The proposed algorithm is applicable for polymer solutions including copolymers that may contain any number of pseudocomponents. The stability algorithm consists of the generation of trial compositions of the polymer pseudocomponents from the parameter space of a distribution function rather than from the more conventional composition space of each pseudocomponent. This, together with the proposed pressure-temperature (PP flash algorithm in which the number of equations to be solved is made independent of the number of pseudocomponents used to represent a polymer's molecular weight distribution (MWD), considerably reduces the complexity and dimensionality of the phase equilibrium problem. Calculations using this approach show promise in the ability to detect thermodynamic instability in a polydisperse polymer system, even for polymers with multimodal MWD that are in solution with multiple solvents.