The development of scalable robust solvers for unstructured finite element applications related to viscous flow problems in earth sciences is an active research area. Solving high-resolution convection problems with order of magnitude 108 degrees of freedom requires solvers that scale well, with respect to both the number of degrees of freedom as well as having optimal parallel scaling characteristics on computer clusters. We investigate the use of a smoothed aggregation (SA) algebraic multigrid (AMG)-type solution strategy to construct efficient preconditioners for the Stokes equation. We integrate AMG in our solver scheme as a preconditioner to the conjugate gradient method (CG) used during the construction of a block triangular preconditioner (BTR) to the Stokes equation, accelerating the convergence rate of the generalized conjugate residual method (GCR). We abbreviate this procedure as BTA-GCR. For our experiments, we use unstructured grids with quadratic finite elements, making the model flexible with respect to geometry and topology and O(h3) accurate. We find that AMG-type methods scale linearly (O(n)), with respect to the number of degrees of freedom, n. Although not all parts of AMG have preferred parallel scaling characteristics, we show that it is possible to tune AMG, resulting in parallel scaling characteristics that we consider optimal, for our experiments with up to 100 million degrees of freedom. Furthermore, AMG-type methods are shown to be robust methods, allowing us to solve very ill-conditioned systems resulting from strongly varying material properties over short distances in the model interior.