We have critically examined and compared various ways to obtain standard harmonic molecular mechanics (MM) force-field parameters for metal sites in proteins, using the 12 most common Zn2+ sites as test cases. We show that the parametrization of metal sites is hard to treat with automatic methods. The choice of method is a compromise between speed and accuracy and therefore depends on the intended use of the parameters. If the metal site is not of central interest in the investigation, for example, a structural metal far from the active site, a simple and fast parametrization is normally enough, using either a nonbonded model with restraints or a bonded parametrization based on the method of Seminario. On the other hand, if the metal site is of central interest in the investigation, a more accurate method is needed to give quantitative results, for example, the method by Norrby and Liljefors. The former methods are semiautomatic and can be performed in seconds, once a quantum mechanical (QM) geometry optimization and frequency calculation has been performed, whereas the latter method typically takes several days and requires significant human intervention. All approaches require a careful selection of the atom types used. For a nonbonded model, standard atom types can be used, whereas for a bonded model, it is normally wise to use special atom types for each metal ligand. For accurate results, new atom types for all atoms in the metal site can be used. Atomic charges should also be considered. Typically, QM restrained electrostatic potential charges are accurate and easy to obtain once the QM calculation is performed, and they allow for charge transfer within the complex. For negatively charged complexes, it should be checked that hydrogen atoms of the ligands get proper charges. Finally, water ligands pose severe problems for bonded models in force fields that ignore nonbonded interactions for atoms separated by two bonds. Complexes with a single water ligand can normally be accurately treated with a bonded potential, once it is ensured that the water H atoms have nonzero Lennard-Jones parameters. However, for metal sites with several water molecules, a nonbonded model with restraints (taken from the QM calculations) is more stable.