Molecular statics calculations on gas-phase and solvated clusters and on gas-phase and solvated slabs representing aqueous species and surfaces were applied to investigate acid/base reactions on silica surfaces. Our gas-phase approach, which was previously applied to goethite, predicts a surface pK(a) of 8.5 for the reaction >SiOH-->>SiO- + H+ which is in good agreement with estimates based on potentiometric titration. However, the model gives an unrealistically large pK(a) for the reaction >SiOH2+ --> >SiOH + H+. The model dependence of this result was checked by using two different types of interaction potentials, one based on quantum mechanical calculations on H4SiO4 clusters, and another empirical model fitted to the structure and elastic properties of Lu-quartz. Because these models gave similar results, we hypothesize that the failure of the gas-phase models is due to intrinsic solvation effects not accounted for by our previously developed correlations. We tested this idea by carrying out energy minimization calculations on gas-phase clusters with one hydration shell as well as molecular dynamics simulations on fully-solvated H5SiO4+ and a fully solvated (0001) surface of beta-quartz. Though we are unable to establish a quantitative measure of the pK(a) of SiOH2 groups, the solvated systems do indicate that SiOH groups do not protonate in any of our solvated models.