The spectral element method is used to discretize self-adjoint elliptic equations in three-dimensional domains. The domain is decomposed into hexahedral elements, and in each of the elements the discretization space is the set of polynomials of degree N in each variable. A conforming Galerkin formulation is used, the corresponding integrals are computed approximately with Gauss-Lobatto-Legendre (GLL) quadrature rules of order N, and a Lagrange interpolation basis associated with the GLL nodes is used. Fast methods are developed for solving the resulting linear system by the preconditioned conjugate gradient method. The conforming finite element space on the GLL mesh, consisting of piecewise Q(1) or P-1 functions, produces a stiffness matrix K-h that is known to be spectrally equivalent to the spectral element stiffness matrix K-N. K-h is replaced by a preconditioner Kh which is well adapted to parallel computer architectures. The preconditioned operator is then <(K)over tilde (-1)(h)> K-N. Techniques for nonregular meshes are developed, which make it possible to estimate the condition number of <(K)over tilde (-1)(h)> K-N, where <(K)over tilde (h)> is a standard finite element preconditioner of K-h, based on the GLL mesh. Two finite element-based preconditioners: the wirebasket method of Smith and the overlapping Schwarz algorithm for the spectral element method are given ay examples of the use of these tools. Numerical experiments performed by Pahl are briefly discussed to illustrate the efficiency of these methods in two dimensions.34624822502