8 chirality interactions in cold atom optical lattices has been proposed $^{38}$. Our model (8) is achieved at second order of the per- turbation series. Higher order terms become trunca- tion errors but may be controlled by small parameters λ$_{x,y,z}$/J$_{cluster}$ ∼ √ | J$_{x,y,z}$ | /J$_{cluster}$ . V. CONCLUSIONS. We constructed the exactly solvable Kitaev honeycomb model 1 as the exact low energy effective Hamiltonian of a spin-1/2 model [equations (8) or (9)] with spin-rotation and time reversal symmetry. The spin in Kitaev model is represented as the pseudo-spin in the two-fold degenerate spin singlet subspace of a cluster of four antiferromag- netically coupled spin-1/2 moments. The physical spin model is a honeycomb lattice of such four-spin clusters, with certain inter-cluster interactions. The machinery for the exact mapping to pseudo-spin Hamiltonian was developed (see e.g. TABLE I), which is quite general and can be used to construct other interesting (exactly solvable) spin-1/2 models from spin rotation invariant systems. In this construction the pseudo-spin correlations in the Kitaev model will be mapped to dimer or spin-chirality correlations in the physical spin system. The correspond- ing picture of the fractionalized Majorana fermion exci- tations and Ising vortices still remain to be clarified. This exact construction contains high order physical spin interactions, which is undesirable for practical im- plementation. We described two possible approaches to reduce this problem: generating the high order spin in- teractions by perturbative expansion of the coupling to optical phonon, or the magnetic coupling between clus- ters. This perturbative construction will introduce trun- cation error of perturbation series, which may be con- trolled by small expansion parameters. Whether these constructions can be experimentally engineered is how- ever beyond the scope of this study. It is conceivable that other perturbative expansion can also generate these high order spin interactions, but this possibility will be left for future works. Acknowledgments The author thanks Ashvin Vishwanath, Yong-Baek Kim and Arun Paramekanti for inspiring discussions, and Todadri Senthil for critical comments. The author is sup- ported by the MIT Pappalardo Fellowship in Physics. Appendix A: Coupling between Distortions of a Tetrahedron and the Pseudo-spins In this Appendix we reproduce from Ref. 35 the cou- plings of all tetrahedron distortion modes to the spin Now we can use TABLE I to convert the above cou- plings into pseudo-spin. It is easy to see that f A and f T$_{2}$ 1 , 2 , 3 are all zero when converted to pseudo-spins, namely projected to the physical spin singlet sector. But f E 1 = ( P$_{14}$ + P$_{23}$ + P$_{24}$ + P$_{13}$ − 2 P$_{12}$ − 2 P$_{34}$ ) / $_{(4}$√$_{3) =}$ − $_{(}$√$_{3}$$_{/}$$_{2)}$$_{τ}$ x and f E 2 = ( P$_{24}$ + P$_{13}$ P$_{14}$ P$_{23}$ ) / $_{4 = (}$√$_{3}$$_{/}$$_{2)}$$_{τ}$ $^{y}$. This has already been noted by Tchernyshyov et al. $^{28}$, only the E modes can lift the degeneracy of the physical spin singlet ground states of the tetrahedron. Therefore the general spin lattice coupling is the form of (12) given in the main text. f A = √ 2 / 3 ( S$_{1}$ · S$_{2}$ + S$_{3}$ · S$_{4}$ + S$_{1}$ · S$_{3}$ + S$_{2}$ · S$_{4}$ + S$_{1}$ · S$_{4}$ + S$_{2}$ · S$_{3}$ ) . J $^{′}$( Q $^{A}$f $^{A}$+ Q E 1 f E 1 + Q E 2 f E 2 + Q T$_{2}$ 1 f T$_{2}$ 1 + Q T$_{2}$ 2 f T$_{2}$ 2 + Q T$_{2}$ 3 f T$_{2}$ 3 ) H$_{cluster}$$_{,}$ SL = ( J$_{cluster}$/ 2)( ∑ ℓ S$_{ℓ}$ ) 2 + J ′ ∑ ℓ<m δr$_{ℓm}$ ( S$_{ℓ}$ · S$_{m}$ ) (A1) system. And convert them to pseudo-spin notation in the physical spin singlet sector. Consider a general small distortion of the tetrahedron, the spin Hamiltonian becomes where δr$_{ℓm}$ is the change of bond length between spins ℓ and m , J ′ is the derivative of J$_{cluster}$ with respect to bond length. There are six orthogonal distortion modes of the tetra- hedron [TABLE 1.1 in Ref. $^{35}$]. One of the modes A is the trivial representation of the tetrahedral group T$_{d}$ ; two E modes form the two dimensional irreducible representa- tion of T$_{d}$ ; and three T$_{2}$ modes form the three dimen- sional irreducible representation. The E modes are also illustrated in FIG. 3. The generic couplings in (A1) [second term] can be converted to couplings to these orthogonal modes, where Q are generalized coordinates of the corresponding modes, functions f can be read off from TABLE 1.2 of Ref. $^{35}$. For the A mode, δr$_{ℓm}$ = √ 2 / 3 Q $^{A}$, so f A is The functions f E 1 , 2 for the E modes have been given before but are reproduced here, The functions f T$_{2}$ 1 , 2 , 3 for the T$_{2}$ modes are f T$_{2}$ 1 = ( S$_{2}$ · S$_{3}$ − S$_{1}$ · S$_{4}$ ) , f T$_{2}$ 2 = ( S$_{1}$ S$_{3}$ S$_{2}$ S$_{4}$ ) , f T$_{2}$ 3 = ( S$_{1}$ · S$_{2}$ − S$_{3}$ · S$_{4}$ ) f E 2 = (1 / 2)( S$_{2}$ · S$_{4}$ + S$_{1}$ · S$_{3}$ − S$_{1}$ · S$_{4}$ − S$_{2}$ · S$_{3}$ ) , f E 1 = √ 1 / 12( S$_{1}$ · S$_{4}$ + S$_{2}$ · S$_{3}$ + S$_{2}$ · S$_{4}$ + S$_{1}$ · S$_{3}$ − 2 S$_{1}$ · S$_{2}$ − 2 S$_{3}$ · S$_{4}$ ) .