The Einstein constraints equations are of fundamental interest in mathematical physics in the study of Einstein's theory of general relativity. This coupled nonlinear elliptic system must also be solved in various forms for gravitational wave simulation. These equations have been studied intensively for half a century; they are a particular example of a ``critical exponent'' problem, further examples of which arise commonly in geometric analysis. In this lecture, we first present some new results that extend the known solution theory for the Einstein constraint equations. We then examine Galerkin and Petrov-Galerkin type methods, and develop both a priori and a posterior error estimates for particular classes of methods, including both classical finite element methods and methods in the setting of finite element exterior calculus. Our goal is to develop provably good methods which do not depend on unrealistic angle conditions to control nonlinear terms in the equations. This is joint work with a number of collaborators over several years.