Phase separation has emerged as an essential concept for the spatial organization inside biological cells. However, despite the clear relevance to virtually all physiological functions, we understand surprisingly little about what phases form in a system of many interacting components, like in cells. Here we introduce a numerical method based on physical relaxation dynamics to study the coexisting phases in such systems. We use our approach to optimize interactions between components, similar to how evolution might have optimized the interactions of proteins. These evolved interactions robustly lead to a defined number of phases, despite substantial uncertainties in the initial composition, while random or designed interactions perform much worse. Moreover, the optimized interactions are robust to perturbations, and they allow fast adaption to new target phase counts. We thus show that genetically encoded interactions of proteins provide versatile control of phase behavior. The phases forming in our system are also a concrete example of a robust emergent property that does not rely on fine-tuning the parameters of individual constituents.