The development of highly accurate and computationally efficient modeling software based on Gibbs energy minimization (GEM) makes it possible to thermodynamically simulate geochemically realistic subsurface fluid-rock interaction processes. This involves consideration of non-ideal multicomponent-multiphase systems that include dilute to concentrated aqueous electrolyte solutions, mineral solid solutions, supercritical fluids, silicate and metal melts, and sorption and ion exchange phases. Predicting the stability and thermodynamic properties of non-ideal solution phases over wide ranges of pressures and temperatures requires that theoretically sound and sufficiently accurate equation of state and activity models are used within the GEM framework. The variety of such models calls for a novel, flexible, and computationally efficient code architecture that supports a wide range of models of non-ideal mixing with different mathematical structures and input data. Here, we introduce the TSolMod C++ class library for equation of state and activity models, implemented within the GEMS3K solver of geochemical equilibria as part of the GEM-Selektor code package (http://gems.web.psi.ch). Essential features of the TSolMod library include a generic and flexible model parameter setup, computationally efficient data exchange with the GEM algorithm, and a straightforward extensibility with any new models of mixing. The current version of TSolMod features a comprehensive selection of fluid, gas, liquid, and solid solution models of interest for geochemical, petrological, material science, and chemical engineering applications.