In this thesis we consider a fully coupled model which aims at reproducing some qualitative features of the electro-mechanical activity of the heart. The models used to describe both the electrical and mechanical activities are relatively simple. However, coupling them together can give rise to numerical instabilities or incorrect predictions. After having introduced each of the sub-models of the fully coupled system we perform some numerical experiments to draw some insights on the numerical approximation of this problem. Firstly we focus on the numerical approximation of the Aliev-Panfilov model, which controls the electrical activation of the muscle. We verify that different approaches can be followed to solve such a problem by the finite element method reducing the computational effort. However each approach can lead to inaccurate predictions of the front velocity. Then we suggest also two numerical schemes for time integration particularly suited for PDEs such as the Aliev-Panfilov model: the operator splitting method and the Runge-Kutta-Chebyshev method. When considering the fully coupled problem, we examine two ways of reducing the computational cost: treating some of the coupling terms explicitly or solving the linearised system iteratively. We verify that with the first choice we can experience numerical instabilities depending on the numerical scheme used for time integration. On the other hand, when solving the linearised system iteratively, key points to solve the problem efficiently are the choice of an adaptive stopping criterion and a good preconditioner. From the numerical experiments performed we conclude that the coupling between the active stress and the mechanics is very influential on the stability of the system and on the convergence of the residuals.