Icosahedral quasicrystals are the most important and thermodynamically stable in all about 200 kinds of quasicrystals currently observed. Beyond the scope of classical elasticity, apart from a phonon displacement field, there is a phason displacement field in the elasticity of the quasicrystal, which induces an important effect on the mechanical properties of the material and makes an analytical solution difficult to obtain. In this paper, a finite element algorithm for the static elasticity of icosahedral quasicrystals is developed by transforming the elastic boundary value problem of the icosahedral quasicrystals into an equivalent variational problem. Analytical and numerical solutions for an icosahedral Al-Pd-Mn quasicrystal cuboid subjected to a uniaxial tension with different phonon-phason coupling parameters are given to verify the validity of the numerical approach. A comparison between the analytical and numerical solutions of the specimen demonstrates the accuracy and efficiency of the present algorithm. Finally, in order to reveal the fracture behavior of the icosahedral Al-Pd-Mn quasicrystal, a cracked specimen with a finite size of matter is investigated, both with and without phonon-phason coupling. Meanwhile, the geometry factors are calculated, including the stress intensity factor and the crack opening displacement for the finite-size specimen. Computational results reveal the importance of phonon-phason coupling effect on the icosahedral Al-Pd-Mn quasicrystal. Furthermore, the finite element procedure can be used to solve more complicated boundary value problems.