An algorithm for computing a few of the smallest (or largest) eigenvalues and associated eigenvectors of the large sparse generalized eigenvalue problem $Ax = \lambda Bx$ is presented. The matrices A and B are assumed to be symmetric, and haphazardly sparse, with B being positive definite. The problem is treated as one of constrained optimization and an inverse iteration is developed which requires the solution of linear algebraic systems only to the accuracy demanded by a given subspace. The rate of convergence of the method is established, and a technique for improving it is discussed. Numerical experiments and comparisons with other methods are presented.