Andrew, Michael, thanks a lot for the replies.
Andrew, this is roughly what I meant by "searching around". The only
difference is that I also modify the objective function to maximize
diversity and add another constraint to make sure that all subsequent
solutions are within some neighborhood of the optimal one.
The method works OK, the only problem is performance. Basically
finding K solutions amounts to solving K MIP problems. I was hoping
that I could sort of hook up inside GLPK and keep track of solutions
(similarly to what Michael suggested). A naive thought was to use a
callback routine and store solutions when handling GLP_IBINGO events.
I guess that diversity might be really poor then (even if I relax my
requirements and ask only for a set of sub-optimal solutions, i.e.
those within a certain distance from the optimal).