This manuscript, "Nunataks as barriers to ice flow: implications for palaeo ice-sheet reconstructions", presents a suite of ice sheet model simulations to investigate the flow of ice around nunataks. This problem is of interest for two reasons which are foci of this paper: (1) exposure and thinning history of ice sheets are generally measured from samples collected at Nunataks, and (2) Nunataks are features that may not be resolved in coarse-resolution simulations of ice flow evolution. Several ensembles of simulations are presented and discussed which specifically interrogate these questions using idealized configuration of a 2D SSA ice sheet model.

Overall, the core scientific idea of this paper is very interesting and novel. The results showing variations in exposure history on different sides of Nunataks are important and likely to have first-order implications for how an entire subfield (cosmogenic exposure age dating) interprets their results. I have some questions about the details of the modeling and some suggestions with respect to the presentation, but given a suitable revised manuscript, I see no reason why this work shouldn't eventually be published in The Cryosphere.

Major issues:

1. This study uses an SSA model with a minimum ice thickness to simulate shallow ice flow around high-relief obstacles. In principle, I can see how this model is appropriate, but it is a fairly non-standard use of SSA. For one, it is unclear to what extent SSA applies for ice that changes from kilometer-scale thickness to meter scale thickness. Further, even if the exposure threshold you use is above the minimum imposed modeled thickness, the use of a minimum in of itself will potentially affect how ice flow occurs near the ice margin on the Nunatak (i.e. constant thickness at the minimum over space vs. gradually decreasing thickness over space). There are also missing details regarding the stress and flux boundary conditions at these sorts of margins. At a minimum, the model descriptions needs more details and justification for why we should believe SSA in such a situation (as opposed to, say, Full Stokes). It may also help to investigate how the results (particular in terms of simulated exposure timing) depend on the choices of minimum model ice thickness and the boundary conditions therein.

2. I come away from reading this study mostly convinced that there are potential issues with the way cosmogenic exposure ages are interpreted. However, I would also have liked to see a bit more concrete solutions/suggestions to solve this problem. Is the answer that every time someone measures a cosmogenic exposure ages, they will have to do detailed high-resolution modeling of ice flow in the region from that time period to place the specific measurements in context? This seems like quite a lot of work. Two less onerous possibilities (only one of which is currently mentioned in the manuscript) intrigue me: (1) can already-measured exposure ages be corrected using the results of this study? What would be necessary to do so? Does one need multiple ages from different parts of the same Nunatak to do such a correction? (2) It seems like you suggest in the discussion that sampling Nunatak at locations that are perpendicular to flow (i.e. on the "sides") would mitigate the age bias somewhat. However, most of your discussion of your experiments was focused on the difference between upstream and downstream. Overall, I think the focus of the discussion needs to be shifted towards providing potential ways forward on a solution, particularly for those who may not have the capability to do detailed simulations in support for exposure ages.

3. I found the text of section 3 very dense with descriptions of results. Rather than re-describe a lot details that are in the figures, making more effort to synthesize results and describe their physical meaning would be welcome. Also, generally shortening section 3 would help as well.

Minor issues:

Throughout: new results should be described using the present tense, whereas prior studies should be described using past tense. There is some inconsistency throughout on this.

Line 15: while rapid ice flow through outlet glaciers...

Line 15: what do you mean by "alleviated the differential response" - this is confusing

Line 23-43: These two paragraphs are quite general and can probably be shortened to one paragraphs

Line 27: ice loss from marine-terminating outlet glaciers

Line 40: drag is exerted on fjords?

Line 56: setting where rock samples are acquired

Line 69: that there is no systematic approach to selecting the sampling

Figure 2: is the coloring of the bars necessary? The information seems redundant with the radial axis.

Figure 2: later in the discussion you seemingly suggest that you might also have information on ages for each of these samples. If so, it would seem to make sense (or at least be interesting) to incorporate this extra data in this figure. Even if there is sampling bias, it would be interesting to show if there was a systematic difference in exposure ages in upstream and downstream samples. However, if this data is not available, thats not an issue.

Line 78: To perform these tests

Line 104: how is the ice front simulated?

Line 118: I felt at this point like it would be helpful to have a figure just showing the bed geometries used in the simulations because it was hard to visualize from the text alone

Line 121: idealized nunatak dimensions

Line 122: are constrained by 33 real nunatak topographic profiles from Antarctica

Line 128: is there any expectation that your results are dependent on the sliding law, or coefficient chosen? Would more rapid sliding give a drastically different result?

Equation 2: is Lx = L*x or L_x?

Line 145-146: confusing sentence

Section 2.5: here and elsewhere you say that typical paleo-ice sheet mode configurations are >5 km. For ice sheet wide simulations, sure, this is true, but I think for individual catchments, there are plenty of people doing ~1 km and higher resolution simulations, even for paleo simulations. Its pretty well within the realm of current SSA models to simulate such resolution on a spatially limited basis (even without adaptive meshing or other sophisticated techniques) for 10's of kyr.

Line 190-193: I'm a bit confused because the ice intersects "above" the lowest ice surface elevation both upstream and downstream. Either this is a mistake, or more explanation is needed because it is counterintuitive that ice thickness would be higher downstream.

Line 199-201: confusing sentence

Line 213: this connection between slope and exposure is not very clear. Could be described more

Figure 4: it is not clear what all the different lines meaning

Figure 4: why does it appear like the "idealized" nunatak (brown in panels as and b) is bumpy? I'm guessing a plotting artifact, but its a bit distracting

Figure 6: difference with respect to what?

Line 287: with-width

Line 294: response to what?

Line 299: In their fjord experiment, Frank

Line 309: between ice thickness on the up and downstream sides o the nunatak increases...

Line 327-328: would this add or subtract to the effect from flow around the nunatak? Could the whole exposure age discrepancy be compensated by this, or made even more severe? Seems like its worthy of a test in your model.

Line 338: where does the 400 m figure come from?

Line 361-364: I referenced this above, but this point seems very important and could do with more references to results, since you should be able to show this in your simulations

In this manuscript, the authors apply the ice sheet model Ua to an idealized, rectangular geometry including Gaussian nunataks, and investigate how the ice flow and geometry react to these barriers. They report that, for conditions mimicking those typical for Antarctica, nunataks can influence the geometry significantly, such that the ice thickens upstream and thins downstream of the obstacle. Details vary depending on scenarios and mesh resolution. The results are relevant because they point out that the interpretation of cosmogenic exposure dating in the vicinity of nunataks in terms of past ice-thickness changes is far from being straightforward and can be grossly wrong if the flow disturbance caused by the nunataks is ignored.

This is a novel contribution with important implications for reconstructions of past ice-sheet geometry. The presentation is generally good. Referee #1 has already pointed out a number of detailed issues, which I am not going to repeat. However, I would like to raise two major issues that deserve some attention:

(1) I have some doubts whether the shelfy-stream approximation (SStA) is appropriate for the tackled problem. SStA is SSA (shallow-shelf approximation) with basal drag, thus it assumes plug flow and a shallow geometry without sharp gradients. Especially when zooming in high resolution to the immediate vicinity of nunataks, these conditions are certainly violated, as the flow regime is strongly 3D and steep surface gradients occur. I think that a proper simulation of this behaviour requires nothing less than full-Stokes dynamics. There is probably not much to be done about the shortcoming at this stage, but at the minimum, the limitation should be discussed very openly and clearly.

BTW, which boundary condition is applied at the nunataks? I suppose they are treated like the ice base?

(2) What I am also missing is a discussion about the potential impact of crevassing in the vicinity of nunataks. I would expect that crevasses often occur around such flow obstacles, with the consequence that the ice is effectively softened and the large-scale flow less disturbed than under the assumption of undamaged ice that has to find its way around the obstacle. This effect has the potential to alter/weaken the described influence on the ice surface significantly. At least, a qualitative discussion about it deems appropriate. Even better would be to investigate the effect by running some tests with an assumed softening around the nunataks.

Please find review comments in the attached document.