In this talk we present progress towards an efficient and accurate computational model for non-Newtonain glacier flow. As a part this, we generalize Glen's Law which is a widely used power law model for the material. A common criticism of the standard law is it has difficulty matching data in both low and high stress regimes simultaneously. Our generalization of the law largely overcomes this difficulty and furthermore does not suffer from numerical difficulties that have been observed in the literature for standard Glen's law flow. In regards to the flow model itself, we focus on a stress-velocity-pressure formulation that resembles a Stokes system with a nonlinear viscosity. The discretization we use employs a least squares finite element principle to minimize the error in an appropriate norm. This produces symmetric positive definite linear systems that are amenable to efficient multigrid solvers. Aside from an overview of the model, we also present some preliminary numerical results for 2-D test problems and discuss current work.