Review of Obase et al. (TC, 2022 submitted)

This manuscript presents simulations of a 1D age and temperature model, mainly for the Dome Fuji ice core and region, but also for the EDC ice core.

The main aim of the manuscript, as the title reads, is to infer potential old ice drilling sites in the Dome Fuji region.

The manuscript is generally of excellent quality. It is precise and reads well.

However, I have a few suggestions for the authors which could further improve the relevance of the manuscript.

I let the authors decide if they want to include these suggestions in their simulations, or simply discuss them in the discussion and outlook sections.

Main comments:

- The model is interesting since it is a transient model, while other models used for the same purpose were steady (or pseudo-steady).

However, the authors do not use the full power of this transient aspect of the model, since they fixed the ice thickness.

As the authors wrote, the ice thickness is a primary parameter controlling the basal melting/temperature and therefore basal ice age.

Therefore, a glacial-interglacial ice thickness change of 200 m can have an important impact on the simulations.

- The authors find a shift between observed and simulated temperature profile near the bed.

They reckon that this is due to polythermal ice, but there is another explanation.

Indeed, the pressure melting point is not so well known.

Appart from pressure, it also depend on the impurities and air content of the ice.

Catherine Ritz discussed that in a thesis 30 years ago, and this discussion is still relevant I think.

- The 1D simulations for EDC are not discussed as much as the simulations for DF.

I understand the authors have deeper interests for Dome Fuji, but I think it could make the manuscript more valuable if the EDC case is discussed more.

For example, I would have been interested by a graph showing the basal melting variations at EDC with time.

- On the contrary, I did not find the simulation along the DF-NDF transect so interesting.

To make it really interesting, it would have been necessary to invert the parameters (in particular accu and GHF) to fit the observed isochrones.

There is no reason to assume accu varies linearly and GHF is constant.

- The thermal parameters of the ice (conductivity, heat capacity) are not so well known.

There are several parametrizations.

Conductivity also depend on the fabric, which makes it even more challenging to estimate them.

I think a discussion on these different parametrizations would have been valuable.

- Catherine Ritz showed a long time ago that it is best to simulate the temperature variations in the bedrock.

Indeed, temperature waves propagates in the upper continental crust and the geothermal flux at the ice-bedrock interface cannot be assumed constant with time.

- It would have been interesting to make a Monte-Carlo simulation for DF and EDC to see which sets of parameters are acceptable.

Here, the parameters are changed one after the other but there are probably covariances.

Minor comments:

- l. 77-80: What is important for applying a 1D ice flow model is not the value of the horizontal velocity, but how the ice flow parameters (e.g., ice thickness) varies upstream.

For example, an ice flow line can be 100 km long with a surface velocity of 1 m/yr, and a 1D model could still be appropriate if everything is constant upstream.

- l. 97-98: Parrenin et al. (2017) did not exactly assume that basal melting was constant.

They used the pseudo-steady assumtion, which states that temporal variations in basal melting are the same than temporal variations of surface accumulation rate.

- eq. (2) and l. 130-131: I think there is an inconsistency here.

As eq. (2) is written, a positive value of Mb means ice refreezing, not melting.

- l. 133-134 and eq. (3): This equation was first formulated in Parrenin and Hindmarsh (2007).

- l. 145-146: Regarding the neglecting of heat production, I think it could justified by the small ice deformation near a dome (very low horizontal shear which is the dominant factor elsewhere).

- l. 146-148: There are different parametrizations of ice conductivity and thermal capacity (see comment above).

These are not discussed here, but I reckon they can have an important effect.

- l. 148: Is it not 917 kg/m^3 the standard value for ice density? (note the wrong unit in the manuscript).

- eq. (5) and (6) assumes a constant geothermal heat flux, which is not the case since heat waves propagate in the upper continental crust (see comment above).

- l. 164-166: I don't understand this sentence here. The formulation of the model does not allow for polythermal ice, so there is no reason to decrease the vertical resolution.

- l. 203-205: Could you please write the equation relating SAT and accu?

Is it the saturation vapor pressure relationship?

- l. 216: I find it a shame that the ice thickness is fixed despite the model being transient (see comment above).

- l. 228-230: It would have been possible to initialize the age and temperature profile with steady profile, instead of constant values, for a faster convergence.

- l. 247-248: The obs-model temperature shift near the bed is probably due to the formulation of the pressure melting point (see comment above).

- l. 261-263: For sure! Without basal melting, the age is infinite at the base.

- l. 296: Parrenin et al. (2017) also estimated the GHF at EDC (Figure 5c), but the value (~60 mW/m^2) is far higher than what you obtained here.

- l. 299: I find this paragraph a bit short (see comment above).

- l. 317-326: It would have been interesting to make a Monte-Carlo simulation to see which sets of parameters are acceptable (see comment above).

- l. 364-379: This is a very simplified transect simulation (see comment above).

- l. 424-425: Instead of using polythermal ice, use a different parametrization of pressure melting point.

- l. 460-462: Of course, Lilien et al. find different results since they simulated BELDC and not EDC, with a very different ice thickness.

"A one-dimensional temperature and age modeling study for selecting the drill site of the oldest ice core around Dom

e Fuji, Antarctica" by Obase et al. details experiments utilising one dimensional age and temperature modeling of the Antarctic ice sheet.  The validity of the model is demonstrated by comparisons with ice-core based age reconstructions and temperature measurements at both Dome Fuji and the EPICA Dome C core sites.  Parameter sensitivity and selection studies for the Dome Fuji region are then conducted, and finally the optimised model applied to a ground-based radar survey in the region, and the simulated age horizons compared to isochrones from the radar survey.

Overall, the paper is well written and easy to follow and is worthy of publication in the Cryosphere, after minor re visions as detailed below.

Minor issues

============

L83-85 The logic here isn't quite correct.  While a lower accumulation rate is necessary to increase the number of years in a given thickness of ice, a lower accumulation rate will also reduce the vertical advection of cold from the surface down into the interior of the ice sheet, therefore increasing the temperature of the ice.  So accumulation rate plays a dual and potentially competing role, but in terms of basal melt rates, lower accumulation is not necessarily a good thing.

L104-106 Parrenin et al 2017 (doi:10.5194/tc-11-2427-2017) applied a time varying rate factor to both the accumulation and melt rates in there 1-D modelling around EDC.  This rate factor was based on variations from the EDC ice core for the last 800ka and was constant before 800ka.

L139-140 Need to make it clear that "m" is Fischer et al's equivalent to "p".  Suggest re-wording from "in the case of m=0.5 in their study" to "where their parameter m fulfils a similar role to p in this study, the case of m=0.5"

L141-142 m=0.5 is only smaller than p=3 for zeta<0.3.  Suggest re-wording from "with a smaller vertical velocity, particularly near the base of the ice" -> "with a smaller vertical velocity in the lower approximately third of the ice" or "with a smaller vertical velocity near the base of the ice"

L161 Are you really calculating the temperature gradient at ice-bedrock interface using a central difference?  If so you would need to be modelling the temperature down into the bedrock.  If you are doing this, you should mention that the thermal domain extends down into the bedrock and give the boundary conditions at the bottom of the rock domain.  If you are only modelling the thermal domain in the ice, then you must be using a one-sided difference discretization at the ice-bedrock interface

L224-226 I think that you have swapped around your "above" and "below" in this sentence.  Surely the age modelling based on orbital tuning of the gas record is for the oldest, and therefore the deepest, part of the ice core, and the matching with AICC2012 is for the younger and shallower part of the core.

L247-248 If the simulated temperatures are colder, especially in the middle of the ice column, this suggests that the downward advection of surface cold is probably too large, indicating that the p value might not be optimal.  It might be worth adding a sentence here outlining this.

L272-273 Your estimate of an annual layer thickness of 0.1mm (Figure 6b, dark blue line) is for a GHF of 52 mW/m^2. You state on lines 250-251 that there has been no melt for a GHF of 52 mW/m2, therefore the age will be greater than 1.5Ma.  At a minimum, you need to delete "of 1.5 MA BP ice" on line 272 because you don't know the age in this case.

L302-314 It is somewhat ambiguous as to what you mean by "different amplitude of temperature changes", especially given your comment on lines 308-209 "because mean temperature over the glacial cycles increased if we reduce a small temperature amplitude of glacial-interglacial cycles."  Presumably, this means that you have kept the interglacial temperatures unchanged and increased the glacial temperatures to change the "amplitude of the changes".  If thiis is the case you should state this somewhere in Section 3.4

L317-326 You might also want to mention that the GHF may vary over the spatial scale of the radar survey, (e.g. Carson et al 2013, doi:10.1144/jgs2013-030), especially given the sensitivity to GHF that you mention on line 276

L348-349 is the impact of the spatial distribution of SMB minor because 1) the sensitivity to SMB is low and/or 2) the spatial variability of SMB is low?

L390 For the radar transect between DF and NDF, while the old ice occurs "where the ice is thin", this is at the expense of the age resolution.  It would be good to add some words to point that out.

L466-470 The model-data discrepancy at 14-18 km from DF corresponds with a relatively cold ice-bedrock interface (Figure 15).  This suggests that perhaps the estimated GHF of 55 mW/m^2 is too low locally, leading to cold ice with little/no basl melt and therefore vertical velocities that are too low.  This is consistent with the model estimating ages that are too shallow.  Such fine spatial scale GHF variations have been noted elsewhere in Antarctic, (see comment above for lines 317-326).

L485-487 See comment above for L272-273

Specific edits

==============

L2 "around" -> "near"

L29-30 This sentence could do with a reference, perhaps something like Shakun et al 2015, doi 10.1016/j.epsl.2015.05.042

L41 "critically scientific challenges" -> "critical scientific challenge"

L59 "in the south" -> "to the south"

L63 it is unusual to talk about an "areal extent", i.e. an area and then give its size in units of length ("50km") rather than area.

L63 "NDF" has not be defined

L78 "Horizontal velocity" -> "Horizontal surface velocity"

L81 "experiments" -> "simulations"

L95-96 "convey the information of surface temperature" -> "advect and diffuse the surface temperature"

L124-124 "zeta=s/H" -> "zeta=z/H"

L131 "ablation" -> "basal melt"

L138 delete "induce"

L145-146 define "T" from equation 4

L159 "335,000 J kg^-1" -> "335 kJ kg^-1"

L242 Even though the section heading mentions "DF" it would be worth making it clear in the opening sentence.  Suggest changing "temperature profiles" -> "DF temperature profiles"

L261 for clarity, suggest changing "reconstructed profiles" to "ice core based reconstructed profiles"

L268 suggest either deleting "as an indicator of old ice" or changing "as an indicator of old ice" -> "as an indicator of sufficient resolution for dating ice based on chemical and isotopic methods"

L289 "Table 2" -> "Table 1"

L330 the results in section 3 included varying GHF, so therefore you need to delete "other"

L382-383 change "using seven colored lines" -> "for seven selected ages"

Figure 2 caption : "Equation [1]" should be "Equation [3]"

Figure 15 caption : need to include what "p" and GHF values are used for this experiment.  Presumably p=3 and GHF=55 mW/m^2

Obase et al. present results for a 1D ice and heat flow model. The goal is to inform site selection for a new core site near Dome Fuji, targeting ice older than the ~700 ka limit of the previous core. The goals of the paper are to: 1) identify parameter combinations that approximately match the Dome Fuji depth-age and borehole temperature relationships and thus can be used for predicting depth-age relationships in the vicinity; 2) identify the primary constraints on the basal ages, which they determine is ice thickness; and 3) apply the model to the radar line that stretches from the previous ice core site to a potential new site, North Dome Fuji.

I am providing only a brief review because I am concerned about the treatment of the basal thermal state in the model. In Figure 5, a change in the geothermal flux of 5 mW m-2 (from 55 to 60 mW m-2) yields a change in the average melt rate of ~2.5 mm/yr (from my eyeballing of the averages). This is too large. It should be about 0.5 mm/yr since 1 mW m-2 can melt approximately 0.1 mm/yr of ice. The caluclation is below:

the melt rate (M) equals the geothermal flux (G) divided by the latent heat (L) and the density of ice (ρ)

M = G / L / ρ = 0.001 (W/m2) / 334000 (J/kg) / 917 (kg/m3)

 So I’m confused why the values in Figure 5 change so much for the modest increase in geothermal flux. I checked this with a model run of my own transient 1D ice and heat flow model with forcings for EDC based on AICC2012. The attached figure shows that modeled melt rate agrees with the calculation above – each 1 mW m-2 of excess geothermal flux causes approximately 0.1 mm/yr of basal melting.

I wonder if the Obase model has a problem with the basal boundary. It sounds like the temperature gradient is being set directly as the ice-rock boundary, instead of in the bedrock well below.  

Unfortunately, the basal melt rate is the controlling factor on the depth-age, such that an error would affect the entire manuscript. I am not sure, but it looks like this problem is also affecting the depth-age relationship in Figure 6.

I initially wondering if there was some nonlinearity model that would amplify the basal melt rate in response to a change in geothermal flux. The basal melt rate affects the vertical velocity. But this has the impact of steepening the basal temperature gradient, allowing more of the heat to be conducted away rather than used to melt basal ice. So that works in the opposite direction. And the model run I performed suggests that there is not a significant non-linearity.

The manuscript addresses an interesting problem of calculating the temporal variations in the basal melt rate and the impact on the depth-age relationship. However, I think the authors need to provide further support that they are calculating the basal melt rate accurately before the remainder of the manuscript is evaluated.