First of all, i would like to thank the authors for responding comprehensively and thoughtfully to all comments by the referees and the editor. I apologize as well that it has taken so long for me to make my decision; it took quite a while to carefully read the new manuscript and 40-page response letter.



I think the work the authors have done to shorten the manuscript and to improve it based on suggestions and comments of the referees (where they see appropriate) is commendable. The manuscript is much more streamlined now and in reading it one is not taken on too many digressions. I think the manuscript is nearly ready for acceptance, but I would like to see just a few minor changes, listed below.



Regarding a physical system that receives quite strong attention by the scientific community, and about which very little is known, there are of course very strong opinions on how to represent the system in a physical model. In the discourse between author and referee i see a number of issues, the most prominent being the representation of channelized flow, but there are other issues as well. My feeling is that the authors received the referee comments as demands the model be modified to represent additional physical processes, leading to a lengthy response to certain comments but (on some issues) few changes to the manuscript (or to the model, but this of course would not be expected). As this paper will not go back to the referees, I feel I should respond on their behalf.



I believe that the reviewers simply wanted more explicit mention of what is and is not represented by the model — and this is understandable because, while the assumptions and limitations of the model are quite clear to a reader who “speaks the language” of the manuscript, there is likely to be a large contingent of people who don’t speak the language using the PISM hydrology model. Already PISM has a wide user community, and I expect it to grow substantially with the release of its hydrological component — and based on what I perceive to be a very large scientific community of researchers interested in subglacial hydrology, a very small portion of whom are familiar with numerical solution of PDEs, there is an ever-increasing need for this type of explicit clarity. Just as the authors are concerned about misconceptions when non-expert users write and publish high-impact papers on the result of this hydrological model - and they will - the reviewers are concerned too, and I share this concern.

Thus, I would like two things made more explicit:

1) In the paragraph beginning at line 95, the rationale for not including conduits is reasonable from a mathematical standpoint; however, it should be acknowledged that channelized flow does occur in Greenland, and is a potentially (i say “potentially” vet an confident it is more of a “definitely”) important physical process to subglacial hydrologic transport that is not captured by this model. That applied mathematics has not yet produced a coarse-grained continuum representation of conduits does not mean it is not important.

2) As I see it, eq 32(b) and eq 18 (and the "cap" for W_til) could be included with the ice dynamics component of PISM, and eqs 32(a,c) could be solved horribly incorrectly, and this would not impact ice dynamics in any way. This is what I believe reviewer 1 meant by 1-way coupling. If I am wrong I must apologize, but I must be convinced. The authors emphasize that the velocity solve is nonlocal and so yield stress does not determine basal stress locally; this is accepted and besides to point I am making. I am not saying this is *wrong*, or that hydrology can cavity size *should* affect basal stress (though in reality i believe it could, see below), but I insist this is briefly pointed out in section 4 or 5.

Minor issues:

l157: $m$ is described as being either due to basal melt or to drainage from the surface. But in 2.3 I might be wrong but I believe it is referred to as coming exclusively from basal melt -- otherwise it is implied that surface drainage goes directly into the till -- which I think is in contrast with the Tsai and Rice (2010) view of the process and others as well. Shouldn't drainage from the surface be a direct source to the hydrological system in addition to \rho_w(\partial W_{til}/\partial t - melt/rho_w)?

l532 harder, not hardest (which sounds awkward)

l533 confusing sentence; suggest to lead with the numerical advantage and follow up with the drawback, e.g. "By contrast, if \phi_0 is larger then the numerical solution of equation (31) is easier; but local changes in subglacial pressure P are damped and spread at the speed of influence to other parts of the connected subglacial hydrologic system."

Finally, I make a few comments regarding other issues I saw in the reviewer-author discourse. Action is not required, so to speak, but i urge the authors to at least consider them.

1) Reviewer 1 commented on the odd combination of a deforming-till basal stress parameterization with a cavity opening/closure hydrologic system evolution. The author interpreted this to mean that cavitation does not occur on deformable beds, and maybe that is what was meant. But in their response the authors cited Schoof (2007). In this paper, the ice-cavity interface is treated as shear stress-free, and I do not know the literature well enough but presumably other models of deformable bed cavitation use the same treatment. I would think this would mean, as with hard-bed cavitation, that the effective coarse-grained basal stress would decrease with some measure of cavity side. Such a dependence is not present in the yield stress formulation used in this model.

2) Reviewer 3 is right about the "corner case" mentioned, though he might not have worded it clearly enough. In the extreme case, if an infinite slab of ice of uniform thickness rested on a uniform layer of water of pressure that exactly balanced the overburden, i don't see what would be driving the (infinitely wide) cavity to close. Replace "infinite" with "very wide" and this is the case the reviewer referred to. Furthermore, the authors say that "none of the models in the literature toe the creep closure rate to the roughness scale". Well, perhaps they do not call it "roughness scale", but this is really just semantics. Plenty of models account for the geometric arrangement of clasts and obstacles in the closure rate, for instance Creyts and Schoof (2009) and Kyrke-Smith, Katz and Fowler (2015).

3) As the authors point out in the manuscript and response to referees, the use of the fixed drainage term C_d is a straightforward extension of Tulaczyk et al's UPB model. It is not the only possible one, however; Van der Wel et al 2013 is cited solely as an example of one-horizontal-dimension drainage, but it is not acknowledged that their model is *another* possible extension, which takes account of vertical transport, something which this model does not do. Furthermore, I don't think anything is straightforward about this drainage term. In Tulaczyk et al 2000b (p485) it is stated that the mechanism of drainage is far from clear and while a constant parameter for drainage is assumed for simplicity, a mechanism that depends on void ratio (and by extension till interstitial water pressure via your eq 19) cannot be discounted.