When optimizing the design of multi-energy systems, the operation strategy and the part-load behavior of the units must be considered in the optimization model, which therefore must be formulated as a two-stage problem. In order to guarantee computational tractability, the operation problem is solved for a limited set of typical and extreme periods. The selection of these periods is an important aspect of the design methodology, as the selection and sizing of the units is carried out on the basis of their optimal operation in the selected periods. This work proposes a novel Mixed Integer Linear Program clustering model, named k-MILP, devised to find at the same time the most representative days of the year and the extreme days. k-MILP allows controlling the features of the selected typical and extreme days and setting a maximum deviation tolerance on the integral of the load duration curves. The novel approach is tested on the design of two different multi-energy systems (a multiple-site university Campus and a single building) and compared with the two well-known clustering techniques k-means and k-medoids. Results show that k-MILP leads to a better representation of both typical and extreme operating conditions guiding towards more efficient and reliable designs. (C) 2019 Elsevier Ltd. All rights reserved.