the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
The Newton solver with step size control is faster than the Picard iteration in simulating ice flow (FEniCSfullStokes v1.1.0)
Angelika Humbert
Thomas Slawig
Abstract. Solving the momentum balance is the computationally expensive part of simulating the evolution of ice sheets. The momentum balance is described by the nonlinear fullStokes equations. As a nonlinear problem, they are solved iteratively. We solve these equations with Newton's method. We obtain global superlinear convergence by using a step size control. For the step size control, we need a minimization problem. Solving the fullStokes equations is equivalent to minimizing a specific convex function. We use the Armijo and the exact step sizes for Newton's method. Additionally, we use the exact step sizes for the Picard iteration. Finally, we compare the Picard iteration and the variants of Newton's method in two benchmark experiments, called ISMIPHOM experiments A and B. These experiments consist of a more realistic domain and are designed to test the quality of ice models. We obtain that Newton's method and the Picard iteration with exact step sizes greatly reduce the necessary number of iterations.
 Preprint
(1229 KB)  Metadata XML
 BibTeX
 EndNote
Niko Schmidt et al.
Status: open (until 02 Nov 2023)

RC1: 'Comment on egusphere20231569', Anonymous Referee #1, 27 Sep 2023
reply
Summary: The paper implements a Newton method with stepsize control for solving a full Stokes model describing the dynamics of glacial ice. In idealised diagnostic numerical experiments two types of step size control are compared with each other and to Picard iterations with and without stepsize control.
General comments:
The chosen topic is important and the paper is well structured and seemingly without any obvious technical mistakes. My main concern is about the novelty of this paper, for the following reasons:
 The Newton method is faster than Picard is expected and Newtons method is already used in ice sheet models extensively, even with step size control (e.g. in Elmer/Ice). It is true that Newton iterations do not converge for all cases especially without running a few initial Picard iterations, but for the numerical experiments chosen in this paper my guess is that Newton works in most existing codes already. Perhaps exact step size is new but I think in that case a brief summary of the most common step size control methods in the most used ice sheet models should be provided. Furthermore, a significant superiority of exact stepsizes compared to Armijo is not clear from this paper.
 In comparison to the previous manuscript which is referred to in this paper: Schmidt, 2013 https://arxiv.org/pdf/2307.02930.pdf In Schmidt, 2013 more mathematical details and proves are given, but the method and experiments seems identical. Is Schmidt, 2013 indented to be published?
Nevertheless, I think the paper can be published after major revisions. My suggestions for addressing my above concerns are:
 If Schmidt, 2013 https://arxiv.org/pdf/2307.02930.pdf is not already submitted somewhere, I think the mathematical proofs here are novel enough for publication. I am not sure if the Cryosphere is the best journal for this, but possibly the mathematical proofs could be explained and summarized in the main text, and the proofs themselves be kept in the appendix.
 Otherwise, a solution which probably fits the Cryosphere better, is to extend the numerical experiments and literature study so that the importance to glaciology is clearer. If the authors can show that it is possible to use Newtons method in realistic simulations in an effective way without initial Picard iterations, I think it is a major contribution to the community. Specifically the following points should be included:
 In the introduction, go through the most important ice sheet models (especially those using the full Stokes model) and review what method they are using to resolve the nonlinearity of the problem
 Extend the numerical experiments. If the aim of this paper is to show that some method is better than another and could be useful for an ice sheet modeller, it must be for a relevant example. Most importantly, a timedependancy (free surface movement) should be introduced. This is important since the first timestep (the only one studied here) usually behaves quite differently in terms of number of iterations, as compared to later timesteps when free surface has relaxed. In a real application, there is also a contact problem between the ice and bedrock, which introduces an extra nonlinearity. The more realistic example (3D, variable bedrock) the better. Also, I think CPUtimes should be measured rather than only number of iterations.
If the original paper Schmidt, 2013 https://arxiv.org/pdf/2307.02930.pdf is kept as a separate paper, I think still some of the mathematical details must be added to this paper. I am lacking some mathematical justification of some of the key points of the paper, for instance the justification of the extra diffusion term in the equations, and explanation behind Algorithm 3.
Specific comments:
The abstract and introduction could flow better. The sentences could be longer or better connected.
L 4 and 5: I think these sentences can be removed from the abstract “For the step size control, we need a minimization problem. Minimizing a specific convex function is equivalent to solving the fullStokes equations.”
L 18: nonlinear > nonlinearly?
L 18: make it clearer that the stationary variation of these equations is the standard one
L 56: These spaces are for linear Stokes, see e.g. the work by Belenki et al: https://www.jstor.org/stable/41582741 for the appropriate spaces
L 60: Why do you add a diffusion term? This is nonstandard and not included in other ice sheet models, which might make the reader wonder if these results are applicable to their model. I suggest trying to get rid of this term.
L 71: I would not say that Picard is standard. Include some references for ice sheet models that use Picard rather than Newton. Most probably use a combination of the two.
L 9293: It’s true that in Hirn \mu_0 is 0 (and even delta=0 !) . This should be discussed  why is it not possible to set them to zero in this work?
Algorithm 3: what is the stopping criteria, i.e. how many forloop iterations do you do?
L 104: For the reader of the cryosphere, I think you should introduce the two stepsize control methods with a bit more words, discuss pros and cons etc, rather than just stating the algorithm.
Equation 1011: Does the fact that the terms are smaller really mean that one can surely say that the solution is not impacted significantly, especially in a timedependent simulation it is not clear to me.
Figure 2: I think it is misrepresentative to compare to models that are not full stokes, please only include full stokes models for the dashed lines.
About the experiments: I strongly recommend that you add a more complex example, something with more complex geometry and timedependency. A timedependent simulation of the arolla glacier maybe.
About the experiments: Add plots with relative error vs CPUtime, to show if the stepsize control increases the work per iteration significantly or not
The outlook: This is a little bit informally written.
Citation: https://doi.org/10.5194/egusphere20231569RC1
Niko Schmidt et al.
Niko Schmidt et al.
Viewed
HTML  XML  Total  BibTeX  EndNote  

124  32  7  163  3  3 
 HTML: 124
 PDF: 32
 XML: 7
 Total: 163
 BibTeX: 3
 EndNote: 3
Viewed (geographical distribution)
Country  #  Views  % 

Total:  0 
HTML:  0 
PDF:  0 
XML:  0 
 1