Divalent metal ions are of fundamental importance to the function and folding of nucleic acids. Divalent metal ion–nucleic acid interactions are complex in nature and include both territorial and site specific binding. Commonly employed nonbonded divalent ion models, however, are often parametrized against bulk ion properties and are subsequently utilized in biomolecular simulations without considering any data related to interactions at specific nucleic acid sites. Previously, we assessed the ability of 17 different nonbonded Mg2+ ion models to reproduce different properties of Mg2+ in aqueous solution including radial distribution functions, solvation free energies, water exchange rates, and translational diffusion coefficients. In the present work, we depart from the recently developed 12–6–4 potential models for divalent metal ions developed by Li and Merz and tune the pairwise parameters for Mg2+, Mn2+, Zn2+, and Cd2+binding dimethyl phosphate, adenosine, and guanosine in order to reproduce experimental site specific binding free energies derived from potentiometric pH titration data. We further apply these parameters to investigate a metal ion migration previously proposed to occur during the catalytic reaction of the hammerhead ribozyme. The new parameters are shown to be accurate and balanced for nucleic acid binding in comparison with available experimental data and provide an important tool for molecular dynamics and free energy simulations of nucleic acids where these ions may exhibit different binding modes.