The phase field methods have been employed to study the surface roughening of heteroepitaxial films with cubic elastic anisotropy in three dimensions during the annealing process. And the stress field for an arbitrarily-shaped surface has been simultaneously solved using the phase field microelasticity model. Our simulations reveal that the self-organized surface morphologies are strongly dependent on the elastic anisotropy. It is interesting to note that one can obtain favorable surface morphologies via selection of appropriate epitaxial orientations combined with tuning the strength of elastic anisotropy or modifying the asymmetric strains induced by lattice mismatches between the epitaxial materials and substrates. (c) 2005 Elsevier B.V. All rights reserved.