A comprehensive thermodynamic model based on the electrolyte non-random two-liquid (eNRTL) theory is developed for aqueous Na+-K+-Mg2+-SO42- quaternary system. The model accounts for the liquid phase nonideality by two binary interaction parameters per water-electrolyte and electrolyte-electrolyte pair. The temperature dependence of the binary interaction parameters is further correlated with three temperature coefficients associated with a Gibbs Helmholtz type expression. We obtain the binary interaction parameters for the H2O:(Mg2+ SO42-) pair, the (Na+ SO42-):(Mg2+ SO42-) pair and the (K+SO42-):(Mg2+ SO42-) pair by regressing appropriate thermodynamic data. Integrated with other eNRTL binary interaction parameters in the literature, the eNRTL model accurately represents various thermodynamic properties of aqueous Na+ -K+ -Mg2+ -SO42- quaternary system and its subsystems with temperatures up to 473.15 K and electrolyte concentrations up to saturation. (C) 2015 Elsevier B.V. All rights reserved.