Mammalian metallothioneins (MTs) comprise a Zn(3)Cys(9) cluster in the beta domain and a Zn(4)Cys(11) cluster in the alpha domain. They play a crucial role in storing and donating Zn2+ ions to target metalloproteins and have been implicated in several diseases, thus understanding how MTs release Zn2+ is of widespread interest. In this work, we present a strategy to compute the free energy for releasing Zn2+ from MTs using a combination of classical molecular dynamics (MD) simulations, quantum-mechanics/molecular-mechanics (QM/MM) minimizations, and continuum dielectric calculations. The methodology is shown to reproduce the experimental observations that (1) the Zn-binding sites do not have equal Zn2+ affinity and (2) the isolated beta domain is thermodynamically less stable and releases Zn2+ faster with oxidizing agents than the isolated a domain. It was used to compute the free energies for Zn2+ release from the metal cluster in the absence and presence of the protein matrix (protein architecture and coupled proteinwater interactions) to yield the respective disulfide-bonded product. The results show the importance of the protein matrix as well as protein dynamics and coupled conformational changes in accounting for the differential Zn2+-releasing propensity of the two domains with oxidizing agents.