General comments

This study addresses a longstanding question in ice sheet modeling:  Are there regions where full-Stokes (FS) models are much more accurate than the Blatter-Pattyn (BP) and other higher-order (HO) approximations, and therefore FS models are needed for accurate sea-level projections?  By implementing FS and BP-like solvers using the COMSOL finite element package, the authors have developed a useful tool for exploring differences between FS and BP, without complications due to differing numerics.  They apply these solvers to the central part of the Northeast Greenland Ice Stream (NEGIS) at a wide range of resolutions (0.1 km to 12.8 km) and find modest differences between FS and BP.  As expected, the differences are greatest at fine resolution and in regions with a high slip ratio, high aspect ratio, and/or rough topography.  Differences are small at resolutions of 1 km or coarser, as typically used in ice sheet models.  Since these differences are small compared to other uncertainties in ice sheet modeling, the authors conclude that FS models are not urgently needed for sea-level projections.

The question is an important one and is not settled.  Morlighem et al. (2010) argued that when inverting for basal drag coefficients in ISSM, bridging effects near the grounding line of Pine Island Glacier (PIG) make it essential to use FS models.  Other authors, including Nowicki and Wingham (2008) and Durand et al. (2009), have made similar arguments.  Today, many ice sheet models with HO approximations are being used for Greenland and Antarctic projections on century time scales, while FS models are used sparingly because of their complexity and expense.  As far as I know, the projections based on FS models, including Elmer/Ice, are not dramatically different from projections based on HO models.  I am inclined to agree with Ruckamp et al. that FS–BP differences are secondary compared to other uncertainties, but it is still important to explore these differences systematically.

The authors take a step in this direction with their detailed NEGIS analysis, but this step does not go far in settling the question.  They focus on part of a single ice stream without addressing other regions, such as PIG, where stress-balance terms neglected by BP could be important.  As the authors acknowledge, their analysis does not include prognostic simulations or thermomechanical coupling, and they have tested only one basal sliding law.  So the analysis does not strongly support the conclusions.

I would not expect these complex issues to be resolved in one paper, but I would like to see a broader analysis.  For example, the authors could simulate an entire ice sheet in a transient run of a few years at the highest affordable resolution.  Or they might look at regions where others have argued that FS models are needed.  Perhaps this is impractical for the COMSOL solver and would require a different model such as ISSM.  But I would like a sense that having developed the COMSOL tool, the authors have pushed it as far as they can.  I encourage them to consider other applications that would strengthen their main argument.

Also, the paper needs significant editing for correct English.

Specific comments

Abstract, l. 12:  “severe impacts on internal layers of ice sheets”.   The discussion of englacial advection is limited and does not justify this strong and rather vague statement in the Abstract.

l. 15: “no simplification”. This is a bit strong; all model equation sets have some simplifications (e.g., isotropic, temperature-dependent flow factors in FS models).

l. 26: Please cite the papers by Blatter and Pattyn here rather than below at l. 62. Also, I think “severe” should be “several”.  (Maybe this was also true in l. 12?)

l. 30: “only one contribution”. I think this is the Elmer/Ice model?  Please include the model name and a citation.

l. 39: “with the analytical solution”. I don’t think the ISMIP-HOM experiments have analytical solutions.  The FS results serve as a benchmark to which HO models are compared.

l. 40: “is prohibited”. This is not the right word.  Maybe “is  not possible”?

l. 43: “huge additional computational amount”. Is it possible to give an order-of-magnitude estimate of the additional cost?

ll. 45–55: This paragraph mentions some studies that used FS, but it would be helpful to the reader to give more details, for example the argument of Morlighem et al. (2010) that FS is needed near Antarctic grounding lines. Of the three difficulties listed here for FS–BP comparisons, only (i) is addressed by the present study.  I wish the authors had pressed their analysis further; see the general comments.

l. 52: “which resolves, e.g. bed differently”. The meaning here is not clear.

l. 60: “certainly because FS and higher-order models are too expensive for these long-time integrations”. Some HO models, especially depth-integrated models (e.g., Goldberg 2011), are, in fact, practical for long time integrations.

l. 67: “Therefore…” This is not a strong argument for focusing on a subset of NEGIS.  Does this region have characteristics that could make it especially challenging for BP models as opposed to FS?  The text mentions bed topography, and later slip and aspect ratios, but it is unclear that central NEGIS is an appropriate analog for, say, the PIG grounding line or for glaciers with rougher topography.

l. 74: “ratio between basal sliding and gravitationally driven flow”. I think this should be “basal sliding and internal deformation”.

l. 100: “usually includes the effective pressure”. I suggest “often includes”.  Weertman-type power laws, for example, do not include effective pressure.

ll. 130–140. I like this approach to coding FS and BP in a way that isolates the discarded stress terms and minimizes differences in model numerics.

l. 164: “limiting case”. I think what is meant here is that each solver has a range of resolutions where it is applicable, and 0.4 km lies in the overlap region.  Please clarify; it would help to state the range of applicability for each solver.

l. 188: “On purpose…”. Why was a region chosen far upstream of the grounding line, given that grounding lines are important for sea-level projections and might be regions where FS–BP differences are large?

l. 193: Why a Budd-like friction law? Most ice sheet models are now using a power law, a Coulomb law, or some hybrid combination.  See, e.g., Sect. 2.1 of Asay-Davis et al. (2016).  In general, the analysis would be more compelling if it included more than one basal friction law.

l. 209: What is meant by a symmetry BC? Which field or fields are symmetric across the boundary?

l. 212: The wide range of resolutions is a strong point of this study.

l. 220: “very extreme”. I’m not sure why the authors chose some unrealistic values.  I would suggest three values: E = 1, plus values that are on the low side and the high side but still physically plausible (say, 0.5 and 3).  Below, when there are FS–BS differences with E = 0.1 or E = 6, it is hard to know whether to take those differences seriously.

l. 256: “up to 43 m/a compared to FS”. To give readers a sense of the percentage error, please state the FS value.

l. 260: “for very soft ice…”. Since E = 6 may be unrealistic, the significance of these differences is unclear.  Similarly for the stiff-ice case, l. 273.  See the comment above.

l. 273: “Maximum differences…”. Is this in relative rather than absolute terms?

ll. 315ff: I don’t think Sect. 5.5 adds much to the paper. I would guess that the FS–BP differences are small compared to the vertical diffusion that would be associated with remapping variables onto a modest number of layers, but this could only be shown in a prognostic run.  I would drop this section if it isn’t possible to say more.

l. 330: “whether FS or BP-like is located below or above each other”. Please clarify what this means.

l. 351: “particularly in extreme cases…” See comments above about unrealistic E values.

l. 359: “It might be favorable…” I agree with the statement, but it is not strongly supported by the single example.  See general comments.

l. 361: “our simulations are not prognostic”. As stated above, this is an important limitation.  Analyzing a prognostic problem, if possible, would strengthen the paper.

l. 367: The authors cite Morlighem et al. (2010), but that paper draws a different conclusion (that FS models are essential). Does the NEGIS analysis cast any doubt on the Morlighem conclusions?

l. 375: “there are indications that small initial differences become much larger over long time integrations.” It is good to acknowledge a study’s limitations, but this is another example of the analysis being too limited to draw broad conclusions.

l. 382: “The model disagreements still tend to diverge below…”. This wording is unclear.  Maybe “The models still do not agree at…”

l. 387: “a view on particle pathways…” Again, I think a diagnostic run is not sufficient to draw conclusions on englacial layering.

l. 391: “FS will start to matter…”. This has been shown only for regions that are like NEGIS in relevant ways, where there are no differences introduced by long-term advection, thermodynamic evolution, or grounding lines. 

l. 395: “the use of FS seems not an urgent issue.” See the general comments.

Tab. 1:  This is a short list of constants.  Are there any others?

Fig. 4:  The inset panels in the upper left of each panel are not described in the caption and are hard to read.  Possibly expand to full panels with a separate caption.

Fig. 5:  The axis labels are hard to read.

The figures labeled as A2 to A10 are not related to Appendix A.  These should be included instead as supplementary material.

Technical corrections

The paper contains many grammatical errors and uses of non-standard English.  This is a partial list.

Title:  Usually I have seen “Northeast” rather than “North East” for NEGIS

l. 2: “its owing” is not idiomatic

l. 3: “consequences caused by” is redundant

l. 7: “increases” -> “increase”

l. 15: “is by using” -> “is given by”

l. 26: “computational” -> “computationally”

l. 32: “fidelity to accurately simulate” is awkward. Maybe “ability to accurately simulate”

l. 43: “computational amount consumed by” -> “computational cost”

l. 55: “processes and interactions that make it difficult”

l. 68: “higher variable” -> “highly variable”

l. 74: “the enhancement factor” -> “an enhancement factor”

l. 83: “isometric” -> “isotropic”

l. 99: “outwards of” -> “out of”

l. 109: “explained below” -> “as explained below”

l. 110: “are forming” -> “form”

l. 122: “simplifications … reduces” -> “simplifications … reduce”

l. 126: “are following” -> “follow”

l. 147: “computational amount” -> “computational cost”

l. 188: “upstream from the grounding line”, “downstream from the ice divide”

l. 193: “friction type” -> “friction law”?

l. 232: “towards” -> “and”, “consumed computational resources” -> “computational costs”

l. 236: “compared exemplary” -> “compared”

l. 244: “excepted” -> “expected”

l. 264: “unveils” -> “shows” or “reveals”

l. 267: “show increasing trends” -> “increase”

l. 293: “leveled out” -> “compensated”

l. 309: “less” -> “fewer”

l. 378: “alleviate” -> “allow” or “enable”

l. 379: “issues” -> “differences”

l. 393: “uncertainties” -> “uncertainties in”

Tab. 2 caption:  “Number for” -> “Number of”,  “exemplary listed” -> “listed”

References

Davis, X. S., et al., Experimental design for three interrelated marine ice sheet and ocean model intercomparison projects: MISMIP v. 3 (MISMIP+), ISOMIP v. 2 (ISOMIP+) and MISOMIP v. 1 (MISOMIP1) (2016), Geosci. Model Dev., 9, 2471–2497, doi:10.5194/gmd-9-2471-2016.

Durand, G., O. Gagliardini, B. de Fleurian, T. Zwinger, and E. Le Meur (2009), Marine ice sheet dynamics: Hysteresis and neutral equilibrium, J. Geophys. Res., 114, F03009, doi:10.1029/2008JF001170.

Goldberg, D., A variationally derived, depth-integrated approximation to a higher-order glaciological flow model  (2011), J. Glaciol, 57, 157–170.

Morlighem, M., E. Rignot, H. Seroussi, E. Larour, H. Ben Dhia, and D. Aubry (2010), Spatial patterns of basal drag inferred using control methods from a full‐Stokes and simpler models for Pine Island Glacier, West Antarctica, Geophys. Res. Lett., 37, L14502, doi:10.1029/2010GL043853.

Nowicki, S. M. J., and D. J. Wingham (2008), Conditions for a steady ice sheet–ice shelf junction, Earth Plan. Res. Lett., 265, 246-255.

General comment

The paper addresses the question whether the computationally expensive full Stokes equations are really necessary, or if the cheaper Blatter-Pattyn model is sufficient. They do so by simulating a part of the NEGIS ice stream with both models using the same code (COMSOL) and as similar numerical discretization as possible. The results are then compared using a series of different measures. The experiments show that unless the ice is very stiff, the difference in velocity field amounts to a few percent, but that there is potentially a larger difference in internal ice deformation, which may have implications for paleo reconstructions.

The topic is a very important one and this type of paper is needed. I also appreciate that the authors measure the error in several different ways. However, the study is limited, and I am afraid that the reader might draw the overhasty conclusion that full Stokes is not needed based on these results, while there might be other situations where full Stokes is important. The main limitation as I see it is:

• The lack of a grounding line experiment: I understand that it might be difficult to implement grounding line migration in a commercial software, but I think at the very least the authors can look at a time-independent grounding line problem, comparing velocity fields and buoyancy balance.

Also, I think the authors should consider the following issues:

• The lack of time dependence: The study does not include time dependent simulations. The impact of the model differences on surface evolution can be measured without actually running a surface evolution, I think it would suffice to look at how the velocity field would change the surface after one time step. However, I wonder if it could be a problem that the initial surface has not been relaxed. Are we looking at an artificial initial shock transient? Is that relevant? I would like to see at least a discussion justifying the lack of relaxation.

• The boundary conditions are retrieved by inverting with BP (and a different code). What does it mean for FS that the boundary conditions are consistent with another model? It would be interesting to see another set of experiments, where the slip coefficients are retrieved using FS (e.g. with Elmer) and are then used for both FS and BP. Also, would inverting at a higher resolution make a difference? Would there be high frequency effects that FS would pick up?

• Some tests regarding the numerics are missing. Since quite some effort is taken to treat the models with similar numerics, I would like to see some tests or discussion convincing the reader that the discretization does indeed not impact the result, since the interest is in quite small velocity differences. In particular, the inf-sup stabilization parameter may not be the same for FS and BP (it is not clear if the same stabilization is used for the BP system, but I assume so), and the problem could, if you are unlucky, be sensitive to this. Either change to Taylor-Hood elements or check that varying the stabilization parameter for the inf-sup stabilization does not alter results. Also, the element aspect ratio varies in the experiments, as the number of vertical layers are constant. Will the numerical errors of FS and BP behave the same when element aspect ratio changes? Perhaps this is not relevant but if so, a comment on why should be included.

• The idea of studying internal layers (section 5.5) is nice, but this part of the study is too limited.

Specific comments

Line 30 - worth to mention that that FS simulation was with coarse resolution

Line 76 - Consider changing the title “field equations” to e.g. “the full stokes equations”

Section 2.3 Explain to the reader in what situation all of the neglected stress components are important. Which are important for shearing margins, which are important at grounding line, etc.

Equation 9: Why breaking out 1/2 but not 1/4?

Section 2.3: Comment on that normally one would manipulate the system in the style of this page: http://websrv.cs.umt.edu/isis/index.php/Blatter-Pattyn_model   , maybe write down the “normal” BP system so that the reader more easily understands what you mean by “BP-like” in the next section

Line 112: Add a reference for the difficulties of saddle point problems

Line 137: I worry there could still be numerical issues with how the saddle point FS system

is solved 

Section 3.2: Is the BP-like system symmetric? Does that matter for the linear solver you use?

Line 149-161: shorten this paragraph

Line 180: The paragraph starting here can be clarified, especially for readers who does not have ISMIP-HOM details fresh in mind. Also I think it is the first time the abbreviation HO appears.

Line 193: Comment on why you choose this sliding law, and mention already here that k is to be found with inversion

Figure 1: To me it seems the solutions in Experiment C does not seem to agree well with previous exercises. Comment on this.

Line 207-211: The boundary conditions are not completely clear to me, please write them down as equations.

Line 215: Did you experiment with sensitivity also with respect to vertical resolution? If not please do. Note that the element aspect ratio will change if you only change horizontal resolution, probably it is not an issue here but in worst case it can impact numerics.

Line 222 - 225: Mention that SICOPOLIS use SIA/SSA

Line 226-228: This is an important point to at least discuss in your study. You find a friction coefficient that is consistent with the BP-like model, but use it also for the FS model. 

Line 230-233: This paragraph seems a little bit out of place, and the table could be moved to the appendix

Section 5.2: Here I would appreciate a discussion relating the differences to the missing stress components, or perhaps you can just mention that it will come in section 5.4

Line 273: Mention that a discussion on why stiff ice is more sensitive will come in section 5.4

Section 5.4: I appreciate this section!

Line 302: Write vb/vs first, to be consistent with the order in line 300

Line 300 - 310: This is a good experiment. However the figure (Figure 7) is hard to read. Also, comment on how you think the fact the elements are flatter for high aspect ratios impact the result, or why they don’t impact the result. 

Line 311-314: This is an important check. I think not all readers will understand why you look at the vertical velocity, add a sentence to explain. Perhaps even better, would be too look at how much the surface would move in one time step given this velocity field (this should be easy to compute, you don’t have to actually move the surface)

Section 5.5: I like the idea of looking at internal layers, but this part of the study is quite incomplete

Section 5.6: This section fits better after section 5.3

Minor language/esthetics comments:

Line 26 - “Although BP neglects severe..” - is severe the right word to use? 

Line 48 - “different results as simpler models”, “as” -> “compared to?”

Line 56 - The sentence starting with “Beside the..” is a bit akward

Line 61 - “consistent analysis” - change for “consistent numerical experiments”?

Line 95 - The sentence “Boundary condition ...is traction free” does not seem to be grammatically correct

Line 99 - “v_b the velocity” -> “v_b is the velocity”?

Line 171: The sentence about the mmr contribution is a little bit confusing

Line 236: This sentence is hard to read

Figure 4: The scatter plots are a very small

Line 330: This sentence is unclear, adding a “particle” or “layer” would help

Figure 6: Increase the font

Line 395: “Seems not an urgent issue” -> “does not seem to be an urgent issue”?