General comments:

This manuscript highlights the computational efficiency of using the method of separation of variables to solve the diffusion equation rather than direct numerical methods. While this is well-known, the application to the basal temperature evolution in an ice column had previously only been done for infinitely deep ice; this manuscript considers the case of a shallow ice sheet and performs some example calculations for different initial conditions and boundary conditions.

Overall, I was somewhat disappointed that the authors did not go into more depth in analysing and describing their results, in particular exploring the wealth of curious trends shown in figure 4 - most of the paper is instead given over to a routine description of the method of separation of variables. In particular, given the stated threshold of 2km for the solution to approach the infinite depth limit, it would be nice to explore what factors set this threshold. Looking at figure 5 there seems to be a rather narrow band of depth values for which T is finite but larger than the MacAyeal solution. I think figure 4b also shows this rather sudden regime change.

Specific comments:

If Equation (6) were given as cot(L*sqrt(λ)) = βλ, there would be no need to treat β=0 as a special case.

Figure 4 - the values of the parameters held fixed are not given.

Figure 4d - interesting that T is non-montonic with L at -14°C. Why is this?

Figure 4c - this figure shows the most interesting trends, but is barely discussed in the text. Perhaps using θ_L/L as the primary variable instead would clarify the impact of the temperature gradient on the basal evolution.

Line 162 - where T saturates to above 25kyr, are we in fact in a limit where T is infinite?

Convergence towards no dependence on the detailed surface boundary conditions as L→∞ could be moved to an appendix for better flow of the manuscript.

Technical corrections:

Figure 4 colorbar caption could be oriented to match the axis label.

Recommendation:  This manuscript, in anything like its current form, does not seem to contain a publishable idea.  The most generous interpretation is that other researchers, over decades of analysis of temperature conduction in a solid rod, have failed to notice an intrinsic timescale which might relate to ice sheet binge-purge cycles.  If that is so, something this reader thinks is not true, then the way the article is written must be completely redone.  Critically, issues of incoherent definition ("potential periodicity" is here meaningless) and essentially-disregarded parameter dependence (the assumed initial basal temperature and geothermal rates are in fact dominant) must be somehow overcome.  (It would be a different paper if so.)  In any case, the many time scales potentially associated to full, physically-clear binge-purge mechanisms must be carefully considered if the claimed special time scale here is to be taken seriously.

Summary of the manuscript:  The Introduction ties binge-purge (Heinrich event) cycles to ice temperature (which is fine) and concludes by asserting that 7ka periodicity is widely used in the literature.  Section 2 sets up an initial-boundary value problem for a motionless ice column of finite length, with geothermal (Neumann) basal and Robin surface boundary conditions, and linear-in-height initial temperature.  Sections 3 and 4 sketch, with details in the Appendices, a Fourier series solution of the problem, in which (generally) the eigenvalues solve a transcendental equation requiring numerical solution.  Section 5 visualizes the temperature profiles and their time-dependence, with an emphasis on how they depend on the ice thickness L and on beta, an insulation coefficient in the surface Robin condition.  Section 6 starts by defining a certain solution time as "potential periodicity"---there is no given justification for connecting *this* solution time to periodicity!---and then illustrates and discusses dependence of this time on parameters.  Section 7 then focuses on the dependence of the time on L, as L becomes large, revealing a time 6944a in the limit.  (This value, conveniently near 7ka, entirely depends on the assumed conditions at the base, namely the initial basal temperature theta_b and the geothermal rate G/k.)  Finally the Conclusion again emphasizes the role of L.  Appendices then give details of the standard Fourier series analysis.

Major concerns:

Understanding the consequences of conservation of energy in ice sheets is a nontrivial matter, thus it is included as a 3D partial differential equation into most modern ice sheet modeling efforts, and it is important because internal energy (e.g. temperature) is tied to the long time-scales at which ice sheets change.  Because ice sheets are thin, variations in the vertical are generally larger than in the horizontal, but nonetheless the problem is advection-dominated.  In ice columns near the divide the strongest direction of ice advection is typically vertical, but over large areas of ice sheets this direction is horizontal so that column-wise temperature distributions are commonly far from what any isolated vertical-column model might generate.  Furthermore the bases of ice sheets are usually near or at the pressure-melting point.  The thermo-mechanical condition of near-basal ice can dominate overall ice sheet dynamics because the presence of pressurized liquid water facilitates ice deformation and basal sliding.  The near-basal thermal regime is dominated by geothermal flux, dissipation heat from sliding, and at times the transport of liquid water from elsewhere (e.g. ice surface or through subglacial hydrology).  Because of the strong role of liquid water, it follows that conservation of energy is a two-phase problem, thus not one which can be well-modeled by temperature alone.

The current manuscript considers none of these realities, nor does it provide this reader any insights about ice sheet thermodynamics.  Instead it examines a conduction-only isolated column model.  Within this narrow, unpromising model it proceeds to ignore the dominant parameter dependencies and instead extract a special 7ka time scale, a time scale for temperature change at the base, by surreptitiously fixing some dominant, but unexamined, values.  Then it confusingly discusses dependence on less-dominant parameters, especially ice thickness L and surface conduction beta, simultaneously arguing that L is important and irrelevant.

Thus the manuscript first fails to consider the actual thermodynamics of ice sheets, and then it makes unreasonable claims for the relevance of its very-simplified model.  An extremely well-trod mathematical analysis, namely Fourier series applied to conduction in an interval, a problem already addressed by Fourier and Kelvin, is offered as new and insightful, which it is not.  The modeled time evolution of a column's basal temperature profile simply does depend strongly on the column thickness L, despite the "strongly dependent" claim in the abstract (line 5).  The particular 7ka time scale revealed herein, and unconvincingly tied to binge-purge oscillations and Heinrich events, actually does have strong dependence on particular basal parameters in the model, namely the assumed geothermal flux rate and initial basal temperature.  However, this special time scale would in any case be destroyed by any (here missing) advection mechanism including sliding, critical to any serious discussion of binge-purge.

A key sentence (lines 138-140) is that "We further calculate the time required for the column base to reach the melting point ..., hereinafter referred as potential periodicity".  There is no offered justification for why this solution time is a "periodicity" for anything!  Indeed binge-purge is a periodic mechanism, one of great interest and importance, but there is not even an attempt to explain why this time is related to the desired periodicity.

This "potential periodicity" time is completely dependent on a parameter which is completely arbitrary, namely theta_b = -10C as the starting point at time 0.  It also depends strongly on the geothermal flux rate, which is known to vary substantially over a continent.  (Geothermal flux rates are available for modern North America and thus could be used to explore this parameter dependence.)  As shown in Figure 4(d), stably across a broad range of ice thicknesses L, variation of theta_b from -15C to -5C implies "potential periodicity" which ranges from about 4ka to about 20ka.  Lines 161-162 actually mention this but the rest of the manuscript drops it: "the potential periodicity appears to be rather sensitive to the initial basal temperature, rapidly saturating to values above 25 kyr for theta_b < -11C".  Attempting to interpret time scales as depending on L seems to deliberately ignore that they depend much more strongly on an uninspected parameters theta_b and G/k.  Possibly theta_b should be regarded here as a proxy for the coldness of the cold part of the atmospheric-driver temperature cycles, but (as far as I can tell) even this is not argued.

Finally I want to describe two important figures, so as to illustrate the inappropriateness of the manuscript's analysis.

Figure 4:  What the parts of this Figure actually show, though this is ignored, is that the strongest dependence of the "potential periodicity" time is on the geothermal flux rate and the initial basal temperature.  The discussion of dependence on air temperature and ice thickness is mostly a distraction.

Figure 5:  Here is my attempt to say what is shown in this Figure; note that Figure 2b in particular supports my interpretation.  A geothermal rate and ice conductivity are fixed, giving a fixed value G/k.  An initial basal temperature (theta_b) is fixed, most likely as -10C consistently with Figures 2a and 3, though its value is unstated.  Then the time for the base to warm to 0C (the mis-named "potential periodicity") is shown as a function of ice thickness L.  Different surface boundary condition treatments give several curves, but for L > 2.5km they all coincide at a time about 7ka.  I observe that the explanation for this value of 7ka is actually quite clear!  Namely, as long as the top of the ice is far away, the chosen values of the initial basal temperature and the geothermal flux rate will determine the time taken for the base of the ice to warm up to 0C; this is a balance of upward conduction with the delivered heat in the time interval.  Thus the special value 7ka is actually (and strongly, and entirely as L goes to infinity) a function of theta_b, G, and k, which were all fixed at certain values for no stated reason.  This dependence should be examined, but instead the paper looks elsewhere, at L and beta, and then it spins the results as related to Heinrich events.