In this paper, we study the scalability of the recently developed Albany/FELIX finite-element based code for the first-order Stokes momentum balance equations for ice flow. We focus our analysis on the performance of two possible preconditioners for the iterative solution of the sparse linear systems that arise from the discretization of the governing equations: (1) a preconditioner based on the incomplete LU (ILU) factorization, and (2) a newly-developed, algebraic multigrid (AMG) preconditioner based on the idea of semi-coarsening. A strong scalability study on a realistic, high resolution Greenland ice sheet problem reveals that, for a given number of processor cores, the AMG preconditioner results in faster linear solve times but the ILU preconditioner exhibits better scalability. A weak scalability study is also performed on a realistic, moderate resolution Antarctic ice sheet problem. Because a substantial fraction of the domain contains floating ice shelves, it is fundamentally different from the Greenland ice sheet problem. Here, we show that as the problem size increases, the performance of the ILU preconditioner deteriorates and the AMG preconditioner maintains scalability. This is because the linear systems are extremely ill-conditioned in the presence of floating ice shelves, which has a more negative effect on ILU preconditioning than on AMG preconditioning. Overall, we find that access to multiple preconditioning methods allows for robust and scalable solves for a range of realistic, large-scale, high-resolution ice sheet problems. Reference: Tezaur, I., R. Tuminaro, M. Perego, A. Salinger, S. Price. On the scalability of the Albany/FELIX first-order Stokes approximation ice sheet solver for large-scale simulations of the Greenland and Antarctic ice sheets, Procedia Computer Science, ICCS 2015 International Conference on Computational Science 51 (2015) 2026-2035. doi: 10.1016/j.procs.2015.05.467.