By combining ab initio results for the electronic structure and phonon spectrum with the group theory, we establish the origin of the Verwey transition in Fe3O4. Two primary order parameters with X3 and ∆5 symmetries are identified. They induce the phase transformation from the hightemperature cubic to the low-temperature monoclinic structure. The on-site Coulomb interaction U between 3d electrons at Fe ions plays a crucial role in this transition -it amplifies the coupling of phonons to conduction electrons and thus opens a gap at the Fermi energy. Published in Phys. Rev. Lett. 97, 156402 (2006). PACS numbers: 71.30.+h, 64.70.Kb, 75.50.Gg The discovery of the remarkable discontinuous drop of the electrical conductivity by two orders of magnitude in magnetite (Fe 3 O 4 ) below the Verwey transition (VT) at T V = 122 K [1] triggered intensive studies of its microscopic origin which have been continued for more than six decades. Verwey suggested explaining it by a metalinsulator transition due to the reduction of electron mobility caused by ordering of Fe 3+ and Fe 2+ ions below T V . In spite of great experimental effort, however, the existence of ionic-like charge ordering at T < T V could not be confirmed, and the origin of the VT in magnetite remains still puzzling [2].Magnetite is a ferrimagnetic spinel with anomalously high critical temperature T c ≃ 860 K. Hence, it is viewed as an ideal candidate for room temperature spintronic applications. It crystalizes in the inverse spinel cubic structure, Fd3m, with two types of Fe sites: the tetrahedral A sites and the octahedral B ones, see Ref. There are also other arguments against simple ionic mechanism of the VT. Firstly, the replacement of oxygen O 16 by O 18 isotope increases T V by a few degrees [9], indicating that the transition cannot be purely electronic.Secondly, signals in diffuse neutron scattering observed above T V at k ∆ = (0, 0, 1 2 ), halfway between Γ and X points with k X = (0, 0, 1), and equivalent points of the reciprocal lattice [10] (in units of 2π a , where a is the cubic lattice constant), change into Bragg peaks below T V . The intensity of critical scattering was succesfully calculated on the basis of the transverse acoustic (TA) phonon mode with the ∆ 5 symmetry. Further neutron studies discovered diffuse scattering with large intensities close to Γ and X points, and the dominant component of X 3 symmetry [11]. Other observations from Raman [12], X-ray [13] and nuclear inelastic scattering measurements [14] suggest that phonons play an essential role in the VT. But phonons alone cannot explain the metal-insulator nature of the transition either. Actually, if they alone were responsible, a phonon soft mode would have been found for the cubic phase, while the present ab initio calculations do not imply such a soft mode behavior.The purpose of this Letter is to resolve the above controversies and to demonstrate a cooperative scenario of the VT. Following the pioneering work of Ihle and Lorenz [15], who considered weak intersite Coulom...