A bound-constrained formulation for complex solution phase minimization

Abstract. Prediction of mineral phase assemblages is essential to better understand the dynamics of the solid Earth, such as metamorphic processes, magmatism and the formation of mineral ore deposits. While recently developed thermodynamic databases allow the prediction of stable phase mineral assemblages for an increasing range of pressure, temperature and compositional spaces, the increasing complexity of these databases results in a significant increase of computational cost, hindering our ability to perform realistic models of reactive fluid/magma transport. Presently, prediction of stable phase equilibrium in complex systems is therefore largely limited by how efficiently single phase minimization can be performed, as more than 75 % of the total computational time is generally dedicated to individual solution phase minimization. This limitation becomes critical for non-ideal solution phase models that involve both a large number of chemical components, and mixing on a large number of sites, resulting in many inequality constraints of the form 0≤xlM≤1, where xlM is the fraction of element l mixing on site M. Here, we present a general reformulation of complex non-ideal solution phases from the thermodynamic database of Holland et al. (2018), which comprises equations of state for multiple mineral solid solutions appearing in magmatic systems, as well as multicomponent silicate melt and aqueous fluid phases. Using a nullspace approach, non-linear inequality constraints governing the site fractions are transformed into equality constraints, and the resulting problem is turned into an bound-constrained optimization problem, subsequently optimized using efficient gradient-based methods. To test our formulation, we apply it to several equations of state for solution phases known for their complexity and compare the results of our approach against classical optimization algorithms supporting inequality constraints. We find that the the Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithm yields by far the best performance and stability with respect to the other investigated methods, improving the minimization time of individual solution phase by a factor ≥10. We estimate that our new approach can improve the computational time of stable phase equilibrium by a factor ≥5, thus potentially allowing to model realistic reactive fluid/magmatic systems by directly integrating phase equilibrium calculations in multiphase thermomechanical codes.