We have carried out both nb initio molecular orbital theory and density functional theory studies of mechanisms of proton transfer reactions in multiple hydrogen bond systems using formamidine and its mono-, di-, and trihydrated complexes as model systems. The highest level of nb initio theory, namely CCSD(T)/6-31G(d,p) at the optimized MP2 geometries, predicts the tautomerization in formamidine to have a high classical barrier of 48.8 kcal/mol. Adding one, two, or three waters to form cyclic hydrogen bond clusters stabilizes the transition state assisting proton transfer via a concerted mechanism and reduces this barrier to 21.9, 20.0, or 23.7 kcal/mol, respectively. Compared to our best ab initio CCSD(T)/6-3 IG(d,p)//MP2 results, we found that, among the local DFT JMW and VWN and nonlocal DFT B-LYP, B-P86, B3-P86, BH&H-LYP, and B3-LYP methods, only the hybrid BH&H-LYP method is capable of predicting the structure and energetic information of both the minimum energy and transition structures at a comparable accuracy with the MP2 level. We also found that using numerical atomic orbital or DFT-based Gaussian-type-orbital (GTO) basis sets yields slightly more accurate DFT results than using an HF-based GTO basis set at the 6-31G(d,p) level.