This study is aimed at developing an explicit solvent model for highly charged species such as polyoxometalates. The model includes solvent molecules in the first solvation shell explicitly, and long-range bulk effects and counter-ions as a set of single point charges. The model strongly stabilizes the electronic structure of the Keggin anion. The energies of the Kohn-Sham orbitals obtained using our model lie very close to those computed using the COSMO continuum solvent model; moreover, the total solvation energy evaluated with our model compares well to the value calculated by COSMO.