A comprehensive thermodynamic model based on electrolyte non-random two-liquid (eNRTL) equation is developed for aqueous Na+-K+-Cl--SO42- quaternary system. To account for the composition dependence of the solution nonideality, the eNRTL model requires two binary interaction parameters for each of the four electrolyte-water and four electrolyte-electrolyte pairs. We obtain the eNRTL binary parameters for (K+SO42-):H2O pair, (K+SO42-):(Ne+SO42-) pair and (K+ Cl-):(K+ SO42-) pair by regressing literature data. The eNRTL binary parameters for the other electrolyte-water and electrolyte-electrolyte pairs are retrieved from the literature. To account for the temperature dependence, each of the eNRTL binary parameter is correlated with a Gibbs-Helmholtz type expression with up to three temperature coefficients. The eNRTL model accurately represents various thermodynamic properties of the aqueous 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.