The electrical activity of the heart may be modeled with a system of partial differential equations (PDEs) known as the bidomain model. Computer simulations based on these equations may become a helpful tool to understand the relationship between changes in the electrical field and various heart diseases. Because of the rapid variations in the electrical field, sufficiently accurate simulations require a fine-scale discretization of the equations. For realistic geometries this leads to a large number of grid points and consequently large linear systems to be solved for each time step. In this paper, we present a fully coupled discretization of the bidomain model, leading to a block structured linear system. We take advantage of the block structure to construct an efficient preconditioner for the linear system, by combining multigrid with an operator splitting technique.