The incorporation of polarization in multiscale quantum-mechanics/molecular-mechanics (QM/MM) simulations is important for a variety of applications, for example, charge-transfer reactions. A recently developed formalism based on a density functional theory description of the QM region and a potential energy function for H2O molecules that includes quadrupole as well as dipole polarizability of the MM region is used to simulate liquid water and water clusters. Analysis of the energy, atomic forces, MM polarization, and structure is presented. A quantitative assessment of the QM/MM-MM/MM interaction energy differences of all possible QM/MM configurations of (H2O)n clusters shows that the interquartile range of the distributions of the QM/MM binding energies is never more than 20 meV/molecule higher or lower than the binding energies produced with either of the single-model results. Comparing these interaction energy differences with the QM/MM induction differences show that they are not systematically caused by the induced MM moments of our polarizable embedding scheme. Optimized hexamer geometries as well as the liquid water structure are shown to be improved in comparison with results obtained using point-charge based embedding models neglecting polarization.