thermodynamic model has been developed for calculating the phase behavior and caloric properties of rare earth elements in multicomponent sulfate solutions over wide ranges of temperature (up to 300 degrees C) at concentrations extending to solid-liquid saturation. The model has been constructed using the previously developed Mixed-Solvent Electrolyte (MSE) framework, which combines the standard-state thermochemical properties of aqueous and solid species with an ion-interaction formulation for the excess Gibbs energy. The model accurately reproduces the properties of binary aqueous rare earth sulfate solutions and multicomponent systems that include sulfuric acid and/or sodium sulfate. Solid-liquid phase diagrams have been constructed to elucidate the solubility of both stable and metastable solid phases. Regularities in the solubility behavior have been identified as a function of the cation radius for rare earth sulfates and double sodium-rare earth sulfate salts. The solubility of rare earth sulfates at ambient conditions shows an s-shape pattern with a local maximum at Pr, followed by a plateau between Sm and Tb and then an increase starting from Tb. On the other hand, the double salts REE2Na2(SO4)(4) show a solubility minimum at Pr, followed by a strong increase in solubility for heavier rare earths. For the double sulfate salts for which no experimental data are available, solubilities have been estimated using an analysis of trends in thermochemical properties across the lanthanide series. The model is part of a systematic effort to develop a comprehensive treatment of the properties of rare earth salt solutions to facilitate the design and optimization of processes for the recovery of rare earth elements. (C) 2018 Elsevier Ltd.