Computational Methods in Science and Engineering

Solution of Nonlinear Stokes Equations Discretized By High-Order Finite Elements on Nonconforming and Anisotropic Meshes, with Application to Ice Sheet Dynamics


Motivated by the need for efficient and accurate simulation of the dynamics of the polar ice sheets, we design high-order finite element discretizations and scalable solvers for the solution of nonlinear incompressible Stokes equations. In particular, we focus on power-law, shear thinning rheologies commonly used in modeling ice dynamics and other geophysical flows. We use nonconforming hexahedral meshes and the conforming inf-sup stable finite element velocity-pressure pairings $\mathbb{Q}_k \times \mathbb{Q}^{disc}_{k-2}$ or $\mathbb{Q}_k \times \mathbb{P}^{disc}_{k-1}$, where $k\ge 2$ is the polynomial order of the velocity space. To solve the nonlinear equations, we propose a Newton--Krylov method with a block upper triangular preconditioner for the linearized Stokes systems. The diagonal blocks of this preconditioner are sparse approximations of the (1,1)-block and of its Schur complement. The (1,1)-block is approximated using linear finite elements based on the nodes of the high-order discretization, and the application of its inverse is approximated using algebraic multigrid with an incomplete factorization smoother. This preconditioner is designed to be efficient on anisotropic meshes, which are necessary to match the high aspect ratio domains typical for ice sheets. As part of this work, we develop and make available extensions to two libraries---a hybrid meshing scheme for the \tt p4est parallel adaptive mesh refinement library and a modified smoothed aggregation scheme for PETSc---to improve their support for solving PDEs in high aspect ratio domains. In a comprehensive numerical study, we find that our solver yields fast convergence that is independent of the element aspect ratio, the occurrence of nonconforming interfaces, and the mesh refinement and that depends only weakly on the polynomial finite element order. We simulate the ice flow in a realistic description of the Antarctic ice sheet derived from field data and study the parallel scalability of our solver for problems with up to 383 million unknowns.


  1. viscous incompressible flow
  2. nonlinear Stokes equations
  3. shear-thinning
  4. high-order finite elements
  5. preconditioning
  6. multigrid
  7. Newton--Krylov method
  8. ice sheet modeling
  9. Antarctic ice sheet

MSC codes

  1. 65N30
  2. 65N55
  3. 76A05
  4. 65F08
  5. 49M15

Get full access to this article

View all available purchase options and get full access to this article.


Department of Defense World Geodetic System 1984: Its definition and relationships with local geodetic systems, Tech. rep. NGA TR8350.2, U.S. National Geospatial-Intelligence Agency, Springfield, VA, 1984. Revised in 2004.
M. F. Adams, A parallel maximal independent set algorithm, in Proceedings of the 5th Copper Mountain Conference on Iterative Methods, Boulder, CO, 1998.
M. Ainsworth and P. Coggins, The stability of mixed hp-finite element methods for Stokes flow on high aspect ratio elements, SIAM J. Numer. Anal., 38 (2000), pp. 1721--1761.
S. Balay, J. Brown, K. Buschelman, V. Eijkhout, W. D. Gropp, D. Kaushik, M. G. Knepley, L. C. McInnes, B. F. Smith, and H. Zhang, PETSc users manual, Tech. rep. ANL-95/11, Revision 3.3, Argonne National Laboratory, Lemont, IL, 2012.
H. Blatter, Velocity and stress fields in grounded glaciers: A simple algorithm for including deviatoric stress gradients, J. Glaciology, 41 (1995), pp. 333--344.
J. Brown, Efficient nonlinear solvers for nodal high-order finite elements in 3D, J. Sci. Comput., 45 (2010), pp. 48--63.
J. Brown, B. Smith, and A. Ahmadia, Achieving textbook multigrid efficiency for hydrostatic ice sheet flow, SIAM J. Sci. Comput., 35 (2013), pp. B359--B375.
C. Burstedde, O. Ghattas, G. Stadler, T. Tu, and L. C. Wilcox, Parallel scalable adjoint-based adaptive solution for variable-viscosity Stokes flows, Comput. Methods Appl. Mech. Engrg., 198 (2009), pp. 1691--1700.
C. Burstedde, L. C. Wilcox, and O. Ghattas, p4est: Scalable algorithms for parallel adaptive mesh refinement on forests of octrees, SIAM J. Sci. Comput., 33 (2011), pp. 1103--1133.
M. Cai, A. Nonaka, J. B. Bell, B. E. Griffith, and A. Donev, Efficient variable-coefficient finite-volume Stokes solvers, Commun. Comput. Phys., 16 (2014), pp. 1263--1297.
C. Canuto, P. Gervasio, and A. Quarteroni, Finite-element preconditioning of G-NI spectral methods, SIAM J. Sci. Comput., 31 (2010), pp. 4422--4451.
E. Chow, R. D. Falgout, J. J. Hu, R. S. Tuminaro, and U. M. Yang, A survey of parallelization techniques for multigrid solvers, in Parallel Processing for Scientific Computing, M. A. Heroux, P. Raghavan, and H. D. Simon, eds., SIAM, Philadelphia, 2006, pp. 179--201.
S. L. Cornford, D. F. Martin, D. T. Graves, D. F. Ranken, A. M. Le Brocq, R. M. Gladstone, A. J. Payne, E. G. Ng, and W. H. Lipscomb, Adaptive mesh, finite volume modeling of marine ice sheets, J. Comput. Phys., 232 (2012), pp. 529--549.
H. Damanik, J. Hron, A. Ouazzi, and S. Turek, A monolithic FEM-multigrid solver for non-isothermal incompressible flow on general meshes, J. Comput. Phys., 228 (2009), pp. 3869--3881.
M. O. Deville and E. H. Mund, Finite-element preconditioning for pseudospectral solutions of elliptic problems, SIAM J. Sci. Statist. Comput., 11 (1990), pp. 311--342.
M. O. Deville, P. F. Fischer, and E. H. Mund, High-Order Methods for Incompressible Fluid Flow, Cambridge Monogr. Appl. Comput. Math. 9, Cambridge University Press, Cambridge, UK, 2002.
M. Durufle, P. Grob, and P. Joly, Influence of Gauss and Gauss-Lobatto quadrature rules on the accuracy of a quadrilateral finite element method in the time domain, Numer. Methods Partial Differential Equations, 25 (2009), pp. 526--551.
S. C. Eisenstat and H. F. Walker, Choosing the forcing terms in an inexact Newton method, SIAM J. Sci. Comput., 17 (1996), pp. 16--32.
H. Elman, V. E. Howle, J. Shadid, R. Shuttleworth, and R. Tuminaro, Block preconditioners based on approximate commutators, SIAM J. Sci. Comput., 27 (2006), pp. 1651--1668.
H. C. Elman, D. J. Silvester, and A. J. Wathen, Finite Elements and Fast Iterative Solvers with applications in incompressible fluid dynamics, Oxford University Press, Oxford, UK, 2005.
C. Förste, R. Schmidt, R. Stubenvoll, F. Flechtner, U. Meyer, R. König, H. Neumayer, R. Biancale, J.-M. Lemoine, S. Bruinsma, et al., The Geoforschungszentrum Potsdam/groupe de Recherche de Geodesie Spatiale satellite-only and combined gravity field models: EIGEN-GL04S1 and EIGEN-GL04C, J. Geodesy, 82 (2008), pp. 331--346.
M. Furuichi, D. A. May, and P. J. Tackley, Development of a Stokes flow solver robust to large viscosity jumps using a Schur complement approach with mixed precision arithmetic, J. Comput. Phys., 230 (2011), pp. 8835--8851.
O. Gagliardini and T. Zwinger, The ISMIP-HOM benchmark experiments performed using the finite-element code Elmer, The Cryosphere, 2 (2008), pp. 67--76.
M. W. Gee, J. J. Hu, and R. S. Tuminaro, A new smoothed aggregation multigrid method for anisotropic problems, Numer. Linear Algebra Appl., 16 (2009), pp. 19--37.
T. Geenen, M. ur Rehman, S. P. MacLachlan, G. Segal, C. Vuik, A. P. van den Berg, and W. Spakman, Scalable robust solvers for unstructured FE geodynamic modeling applications: Solving the Stokes equation for models with large localized viscosity contrasts, Geochem. Geophys. Geosyst., 10 (2009), Q09002.
R. Glowinski and J. Xu, Numerical Methods for Non-Newtonian Fluids: Special Volume, Handbook Numer. Anal. 16, North--Holland, Amsterdam, 2011.
P. P. Grinevich and M. A. Olshanskii, An iterative method for the Stokes-type problem with variable viscosity, SIAM J. Sci. Comput., 31 (2009), pp. 3959--3978.
V. Heuveline and F. Schieweck, On the inf-sup condition for higher order mixed fem on meshes with hanging nodes, ESAIM Math. Model. Numer. Anal., 41 (2007), pp. 1--20.
J. Heys, T. Manteuffel, S. McCormick, and L. Olson, Algebraic multigrid for higher-order finite elements, J. Comput. Phys., 204 (2005), pp. 520--532.
R. C. A. Hindmarsh, A numerical comparison of approximations to the Stokes equations used in ice sheet and glacier modeling, J. Geophys. Res., 109 (2004), F01012.
K. Hutter, Theoretical Glaciology, Mathematical Approaches to Geophysics, D. Reidel, Dordrecht, The Netherlands, 1983.
T. Isaac, C. Burstedde, and O. Ghattas, Low-cost parallel algorithms for 2:1 octree balance, in Proceedings of the 26th IEEE International Parallel and Distributed Processing Symposium, IEEE, Washington, DC, 2012.
T. Isaac, N. Petra, G. Stadler, and O. Ghattas, Scalable and efficient algorithms for the propagation of uncertainty from data through inference to prediction for large-scale problems, with application to flow of the Antarctic ice sheet, J. Comput. Phys., 296 (2015), pp. 348--368.
A. Janka, Smoothed aggregation multigrid for a Stokes problem, Comput. Vis. Sci., 11 (2008), pp. 169--180.
G. Jouvet and J. Rappaz, Analysis and finite element approximation of a nonlinear stationary Stokes problem arising in glaciology, Adv. Numer. Anal., 2011 (2012), 164581.
Y. Kanarska, A. Shchepetkin, and J. McWilliams, Algorithm for non-hydrostatic dynamics in the regional oceanic modeling system, Ocean Model., 18 (2007), pp. 143--174.
E. Karer and J. Kraus, Algebraic multigrid for finite element elasticity equations: Determination of nodal dependence via edge-matrices and two-level convergence, Internat. J. Numer. Methods Engrg., 83 (2010), pp. 642--670.
C. Karney, Algorithms for geodesics, J. Geodesy, 87 (2013), pp. 43--55.
S. Kim, Piecewise bilinear preconditioning of high-order finite element methods, Electron. Trans. Numer. Anal., 26 (2007), pp. 228--242.
A. M. Le Brocq, A. J. Payne, and A. Vieli, An improved Antarctic dataset for high resolution numerical ice sheet models (ALBMAP v1), Earth Syst. Sci. Data, 2 (2010), pp. 247--260.
W. Leng, L. Ju, M. Gunzburger, S. Price, and T. Ringler, A parallel high-order accurate finite element nonlinear Stokes ice-sheet model and benchmark experiments, J. Geophys. Res., 117 (2012), F01001.
T. A. Manteuffel, An incomplete factorization technique for positive definite linear systems, Math. Comp., 34 (1980), pp. 473--497.
D. A. May and L. Moresi, Preconditioned iterative methods for Stokes flow problems arising in computational geodynamics, Phys. Earth Planetary Interiors, 171 (2008), pp. 33--47.
D. McKenzie, The generation and compaction of partially molten rock, J. Petrology, 25 (1984), pp. 713--765.
G. A. Meehl, T. F. Stocker, W. D. Collins, A. T. Friedlingstein, A. T. Gaye, J. M. Gregory, A. Kitoh, R. Knutti, J. M. Murphy, A. Noda, S. C. B. Raper, I. G. Watterson, A. J. Weaver, and Z.-C. Zhao, Global climate projections, in Climate Change 2007: The Physical Science Basis. Contribution of Working Group I to the Fourth Assessment Report of the Intergovernmental Panel on Climate Change, Cambridge University Press, Cambridge, UK, 2007, pp. 747--845.
M. Morlighem, H. Seroussi, E. Larour, and E. Rignot, Inversion of basal friction in Antarctica using exact and incomplete adjoints of a higher-order model, J. Geophys. Res. Earth Surface, 118 (2013), pp. 1746--1753.
J. Nocedal and S. J. Wright, Numerical Optimization, Springer-Verlag, Berlin, Heidelberg, New York, 1999.
L. Olson, Algebraic multigrid preconditioning of high-order spectral elements for elliptic problems on a simplicial mesh, SIAM J. Sci. Comput., 29 (2007), pp. 2189--2209.
F. Pattyn, L. Perichon, A. Aschwanden, B. Breuer, B. de Smedt, O. Gagliardini, G. H. Gudmundsson, R. C. A. Hindmarsh, A. Hubbard, J. V. Johnson, T. Kleiner, Y. Konovalov, C. Martin, A. J. Payne, D. Pollard, S. Price, M. Ruckamp, F. Saito, O. Soucek, S. Sugiyama, and T. Zwinger, Benchmark experiments for higher-order and full-Stokes ice sheet models (ISMIP-HOM), The Cryosphere, 2 (2008), pp. 95--108.
N. Petra, H. Zhu, G. Stadler, T. J. R. Hughes, and O. Ghattas, An inexact Gauss-Newton method for inversion of basal sliding and rheology parameters in a nonlinear Stokes ice sheet model, J. Glaciology, 58 (2012), pp. 889--903.
Y. Saad, A flexible inner-outer preconditioned GMRES algorithm, SIAM J. Sci. Comput., 14 (1993), pp. 461--469.
Y. Saad, Iterative Methods for Sparse Linear Systems, 2nd ed., SIAM, Philadelphia, 2003.
G. Schubert, D. Turcotte, and P. Olsen, Mantle Convection in the Earth and Planets, Cambridge University Press, Cambridge, UK, 2001.
H. Seroussi, H. B. E. N. Dhia, M. Morlighem, E. Larour, E. Rignot, and D. Aubry, Coupling ice flow models of varying orders of complexity with the Tiling method, J. Glaciology, 58 (2012), pp. 776--786.
J. R. Shewchuk, Triangle: Engineering a 2D quality mesh generator and Delaunay triangulator, in Applied Computational Geometry: Towards Geometric Engineering, M. C. Lin and D. Manocha, eds., Lecture Notes in Comput. Sci. 1148, Springer-Verlag, Berlin, 1996, pp. 203--222.
I. K. Tezaur, M. Perego, A. G. Salinger, R. S. Tuminaro, and S. F. Price, Albany/felix: A parallel, scalable and robust, finite element, first-order Stokes approximation ice sheet solver built for advanced analysis, Geosci. Model Development, 8 (2015), pp. 1197--1220.
A. Toselli and C. Schwab, Mixed hp-finite element approximations on geometric edge and boundary layer meshes in three dimensions, Numer. Math., 94 (2003), pp. 771--801.
U. Trottenberg, C. Oosterlee, and A. Schüller, Multigrid, Academic Press, London, 2001.
P. Vaněk, M. Brezina, and J. Mandel, Convergence of algebraic multigrid based on smoothed aggregation, Numer. Math., 88 (2001), pp. 559--579.
P. Vaněk, J. Mandel, and M. Brezina, Algebraic multigrid by smoothed aggregation for second and fourth order elliptic problems, Computing, 56 (1996), pp. 179--196.
A. Zmitrowicz, Glaciers and laws of friction and sliding, Acta Mechanica, 166 (2003), pp. 185--206.

Information & Authors


Published In

cover image SIAM Journal on Scientific Computing
SIAM Journal on Scientific Computing
Pages: B804 - B833
ISSN (online): 1095-7197


Submitted: 25 June 2014
Accepted: 18 August 2015
Published online: 3 November 2015


  1. viscous incompressible flow
  2. nonlinear Stokes equations
  3. shear-thinning
  4. high-order finite elements
  5. preconditioning
  6. multigrid
  7. Newton--Krylov method
  8. ice sheet modeling
  9. Antarctic ice sheet

MSC codes

  1. 65N30
  2. 65N55
  3. 76A05
  4. 65F08
  5. 49M15



Metrics & Citations



If you have the appropriate software installed, you can download article citation data to the citation manager of your choice. Simply select your manager software from the list below and click Download.

Cited By







Copy the content Link

Share with email

Email a colleague

Share on social media

The SIAM Publications Library now uses SIAM Single Sign-On for individuals. If you do not have existing SIAM credentials, create your SIAM account