Specific recognition of ligands including metal ions by proteins is the key of many crucial biological functions and systems. Accurate prediction of the binding strength not only sheds light on the mechanism of the molecular recognition but also provides the most important prerequisite of drug discovery. Computational modeling of molecular binding has gained a great deal of attentions in the last few decades since the advancement of computer power and availability of high-resolution crystal structures. However there still exist two major challenges in the field of molecular modeling, i.e. sampling issue and accuracy of the models. In this work, I have dedicated to tackling these two problems with a noval polarizable force field which is believed to produce more accurate description of molecular interactions than classic non-polarizable force fields. We first developed the model for divalent cations Mg²⁺ and Ca²⁺, deriving the parameters from quantum mechanics. To understand the hydration thermodynamics of these ions we have performed molecular dynamics simulations using our AMOEBA force field. Both the water structures around ions and the solvation free energies were in great accordance with experiment data. We have also simulated and calculated the binding free energies of a series of benzamidine-like inhibitors to trypsin using explicit solvent approach by free energy perturbation. The calculated binding free energies are well within the accuracy of experimental measurement and the direction of change is predicted correctly in all cases. Finally, we computed the hydration free energies of a few organic molecules and automated the calculation procedure. / text