We present a specially designed evolutionary algorithm for the prediction of surface reconstructions. This new technique allows one to automatically explore all the low-energy configurations with variable surface atoms and variable surface unit cells through the whole chemical potential range. The power of evolutionary search is demonstrated by the efficient identification of diamond 2×1 (100) and 2×1 (111) surfaces with a fixed number of surface atom and a fixed cell size. With further variation of surface unit cells, we study the reconstructions of the polar surface MgO (111). Experiment has detected an oxygen trimer (ozone) motif (Plass et al, 1998). We predict a new version of this motif which can be thermodynamically stable at extreme oxygen rich condition. Finally, we perform a variable stoichiometry search for a complex ternary system: semi-polar GaN (1011) with and without adsorbed oxygen. The search yields a non-intuitive reconstruction based on N3-trimers. These examples demonstrate that an automated scheme to explore the energy landscape of surfaces will improve our understanding of surface reconstructions. The method presented in this report can be generally applied to binary and multi-component systems.