We report the first large-scale statistical study of very high-lying eigenmodes ͑quantum states͒ of the mushroom billiard proposed by L. A. Bunimovich ͓Chaos 11, 802 ͑2001͔͒. The phase space of this mixed system is unusual in that it has a single regular region and a single chaotic region, and no KAM hierarchy. We verify Percival's conjecture to high accuracy ͑1.7%͒. We propose a model for dynamical tunneling and show that it predicts well the chaotic components of predominantly regular modes. Our model explains our observed density of such superpositions dying as E −1/3 ͑E is the eigenvalue͒. We compare eigenvalue spacing distributions against Random Matrix Theory expectations, using 16 000 odd modes ͑an order of magnitude more than any existing study͒. We outline new variants of mesh-free boundary collocation methods which enable us to achieve high accuracy and high mode numbers ͑ϳ10 5 ͒ orders of magnitude faster than with competing methods. © 2007 American Institute of Physics. ͓DOI: 10.1063/1.2816946͔ Quantum chaos is the study of the quantum (wave) properties of Hamiltonian systems whose classical (ray) dynamics is chaotic. Billiards are some of the simplest and most studied examples; physically their wave analogs are vibrating membranes, quantum, electromagnetic, or acoustic cavities. They continue to provide a wealth of theoretical challenges. In particular "mixed" systems, where ray phase space has both regular and chaotic regions (the generic case), are difficult to analyze. Six years ago Bunimovich described 1 a mushroom billiard with simple mixed dynamics free of the usual island hierarchies of Kolmogorov-Arnold-Moser (KAM). He concluded by anticipating the growth of "quantum mushrooms;" it is this gardening task that we achieve here, by developing advanced numerical methods to collect an unprecedented large number n of eigenmodes (much higher than competing numerics 2 or microwave studies 3 ). Since uncertainties scale as n −1Õ2 , a large n is vital for accurate spectral statistics and for studying the semiclassical (high eigenvalue) limit. We address three main issues: (i) The conjecture of Percival 4 that semiclassically modes live exclusively in invariant (regular or chaotic) regions, and occur in proportion to the phase space volumes. (ii) The mechanism for dynamical tunneling, or quantum coupling between classically isolated phase space regions. (iii) The distribution of spacings of nearest-neighbor eigenvalues, about which recent questions have been raised. 3 We show many pictures of modes, including the boundary phase space (the so-called Husimi function).