The key issue for disulfide bond engineering is to select the most appropriate location in the protein. By surveying the structure of experimentally engineered disulfide bonds, we found about half of them that have geometry incompatible with any native disulfide bond geometry. To improve the current prediction methods that tend to apply either ideal geometrical or energetical criteria to single three-dimensional structures, we have combined a novel computational protocol with the usage of multiple protein structures to take into account protein backbone flexibility. The multiple structures can be selected from either independently determined crystal structures for identical proteins, models of nuclear magnetic resonance experiments, or crystal structures of homology-related proteins. We have validated our approach by comparing the predictions with known disulfide bonds. The accuracy of prediction for native disulfide bonds reaches 99.6%. In a more stringent test on the reported engineered disulfide bonds, we have obtained a success rate of 93%. Our protocol also determines the oxido-reduction state of a predicted disulfide bond and the corresponding mutational cost. From the energy ranking, the user can easily choose top predicted sites for mutagenesis experiments. Our method provides information about local stability of the engineered disulfide bond surroundings.