Abstract
We have implemented the Jacobian-free Newton-Krylov (JFNK) method
for solving the rst-order ice sheet momentum equation in order to improve
the numerical performance of the Community Ice Sheet Model (CISM), the
land ice component of the Community Earth System Model (CESM). Our
JFNK implementation is based on signi cant re-use of existing code. For example, our physics-based preconditioner uses the original Picard linear solver
in CISM. For several test cases spanning a range of geometries and boundary conditions, our JFNK implementation is 1.84-3.62 times more efficient than the standard Picard solver in CISM. Importantly, this computational
gain of JFNK over the Picard solver increases when re ning the grid. Global convergence of the JFNK solver has been signi cantly improved by rescaling
the equation for the basal boundary condition and through the use of an inexact Newton method. While a diverse set of test cases show that our JFNK
implementation is usually robust, for some problems it may fail to converge
with increasing resolution (as does the Picard solver). Globalization through
parameter continuation did not remedy this problem and future work to improve robustness will explore a combination of Picard and JFNK and the use of homotopy methods.
