Articles | Volume 18, issue 3
Research article
11 Mar 2024
Research article |  | 11 Mar 2024

Impact of boundary conditions on the modeled thermal regime of the Antarctic ice sheet

In-Woo Park, Emilia Kyung Jin, Mathieu Morlighem, and Kang-Kun Lee

A realistic initialization of ice flow models is critical for predicting future changes in ice sheet mass balance and their associated contribution to sea level rise. The initial thermal state of an ice sheet is particularly important, as it controls ice viscosity and basal conditions, thereby influencing the overall ice velocity. Englacial and subglacial conditions, however, remain poorly understood due to insufficient direct measurements, which complicate the initialization and validation of thermal models. Here, we investigate the impact of using different geothermal heat flux (GHF) datasets and vertical velocity profiles on the thermal state of the Antarctic ice sheet and compare our modeled temperatures to in situ measurements from 15 boreholes. We find that the temperature profile is more sensitive to vertical velocity than to GHF. The basal temperature of grounded ice and the amount of basal melting are influenced by both selection of GHF and vertical velocity. More importantly, we find that the standard approach, which consists of combining basal sliding speed and incompressibility to derive vertical velocities, provides reasonably good results in fast-flow regions (ice velocity >50 m yr−1) but performs poorly in slow-flow regions (ice velocity <50 m yr−1). Furthermore, the modeled temperature profiles in ice streams, where bed geometry is generated using a mass conservation approach, show better agreement with observed borehole temperatures compared to kriging-based bed geometry.

1 Introduction

Global warming has been responsible for rapid sea level rise from the mass loss of ice sheets and glaciers over the past few decades. The mass loss of the Antarctic ice sheet has more than tripled over the past 3 decades (IPCC AR6 Chapter 9; Fox-Kemper et al.2021). The retrograde bed slopes in deep submarine basins (Schoof2007), the intrusion of warm water in ice shelf cavities (Alley et al.2016), and the collapse of ice shelves can accelerate this mass loss (Scambos2004), especially in West Antarctica. Ice sheet models have been developed to capture these processes (e.g., Larour et al.2012b; Gillet-Chaulet et al.2012; Pollard and DeConto2012) and provide projections of future contributions of the ice sheets to sea level rise under different warming scenarios (DeConto and Pollard2016; Seroussi et al.2020). However, the uncertainty in these projections remains high partly due to poorly constrained model inputs, such as bed geometry, basal conditions, ice mechanical properties, or oversimplified parameterization of melting rates under floating ice shelves (e.g., Schlegel et al.2013; Brondex et al.2019).

A critical aspect of ice sheet models is their initial conditions. Several important properties, such as ice elevation and surface ice velocity, can be directly observed at the surface of the ice sheet, whereas observing englacial and subglacial properties, such as ice temperature and geothermal heat flux, remains particularly challenging, and direct measurements of these properties are scarce.

In order to get reasonable estimates of these englacial and subglacial fields, inversion techniques are routinely employed to estimate basal friction and ice shelf rigidity (MacAyeal1993; Khazendar et al.2007; Morlighem et al.2010; Gillet-Chaulet2020). These inverse-modeling approaches have not been applied to the ice thermal regime of the ice sheet, which remains highly uncertain despite its critical control over ice viscosity and basal friction. Critically, the geothermal heat flux (GHF) is an important parameter that affects basal temperature, water production, and ice dynamics (Pattyn et al.2008; Seroussi et al.2017; Smith-Johnsen et al.2020b), yet large uncertainties in the spatial variation and magnitude of GHFs in Antarctica still remain.

Previous studies have attempted to infer the GHFs using different methods such as a seismic model (Shapiro and Ritzwoller2004; An et al.2015), magnetic satellite data (Maule et al.2005), and a combination of seismic and magnetic satellite data (Martos et al.2017). The most accurate measurements are from in situ borehole measurements of temperature profiles that can be used to constrain the GHFs (Dahl-Jensen et al.1999; Mony et al.2020; Talalay et al.2020).

While drilling boreholes requires a lot of resources and effort, the boreholes provide critical insights into subsurface conditions and lead to a better understanding of the current subglacial and englacial environments, as well as past climate (Augustin and Antonelli2002; Motoyama2007; Slawny et al.2014; Fisher et al.2015; Priscu et al.2021; Mulvaney et al.2021; Smith et al.2021). Borehole temperature profiles can also be utilized to validate thermomechanical ice sheet models. As boreholes provide vertical temperature profiles, a one-dimensional thermal model is generally utilized to estimate GHFs (Mony et al.2020) and reconstruct past climates (Zagorodnov et al.2012; Yang et al.2018). Since one-dimensional thermal models typically neglect horizontal advection and only consider vertical advection and diffusion (Engelhardt2004a; Mony et al.2020; Talalay et al.2020), they have strong limitations and may not be applicable in regions of fast flow. The vertical velocities used in one-dimensional thermal models are generally recovered through the equation of incompressibility, assuming a stationary bed and no sliding (Hindmarsh1999). Only a handful of three-dimensional thermomechanical ice sheet models have utilized these borehole temperature profiles for validation (Joughin et al.2004; Pattyn2010; Seroussi et al.2013). Moreover, measurements of borehole temperatures in fast-flow regions remain scarce due to the technical difficulty of drilling boreholes in these regions (Engelhardt2004b; Doyle et al.2018; Anker et al.2021).

In addition to being sensitive to the GHF, the ice thermal regime is also particularly sensitive to horizontal and vertical ice velocities. While surface ice velocities can be spatially and temporally observed through satellite remote sensing (Mouginot et al.2012; Derkacheva et al.2020), englacial velocities are difficult to observe remotely. Few measurements of internal vertical ice velocities are available through direct methods, such as fiber optic instruments (Pettit et al.2011) and borehole optical televiewer (OPTV) logging (Hubbard et al.2020), or indirect methods, such as phase-sensitive radio-echo sounders (Gillet-Chaulet et al.2011; Kingslake et al.2014). Due to scarcity of internal ice velocity measurements, three-dimensional mechanical models, such as higher order (HO; Pattyn2003) and full Stokes (FS), are used to estimate internal ice velocities (Pattyn2003; Larour et al.2012b). The ice velocities from mechanical models can, in turn, be used as input variables to compute three-dimensional ice temperature.

Overall, the difficulty in estimating GHF, combined with the lack of observations of subsurface ice velocities and temperature, limits our ability to capture the thermal regime of the ice sheet and increases the uncertainty in future mass projections. Here, we perform a suite of sensitivity experiments using a three-dimensional thermomechanical model and using various GHF sources and different approaches to construct vertical ice velocities. We then compare each modeled temperature to 15 temperature profiles from in situ borehole drilling campaigns, including three boreholes located in fast-flow regions, to determine which combination of parameters best reproduces measured temperature profiles.

2 Methods

2.1 Ice flow model

We use the Ice-sheet and Sea-level System Model (ISSM) to model the stress balance and thermal state across the entire Antarctic continent (Larour et al.2012b). We rely on an anisotropic mesh with a resolution varying from 2 km in coastal regions to 40 km near ice divides, and we refine the mesh to a 2 km mesh around the locations of boreholes where temperature measurements are available. The mesh comprises a total of over 1 million prismatic elements distributed vertically over 15 layers. We use a three-dimensional higher-order (HO; Pattyn2003) model and assume that the ice viscosity follows Glen’s flow (Glen1955):

(1) μ = B 2 ε ˙ e n - 1 n ,

where B is the ice rigidity (Pa s-1/3), ε˙e is the effective strain rate (s−1), and n is Glen’s law exponent, the value of which is 3 in this study. We also utilize the Budd friction law (Budd et al.1979; Morlighem et al.2010):

(2) τ b = - α 2 N v b ,

where α is the friction coefficient (yr0.5 m−0.5), N is the effective pressure (taken here as simply ρigH+ρwgmax0,b), and vb is the basal ice velocity vector. ρi is the ice density, ρw is the water density, H is the ice thickness, and b is the bed elevation with respect to sea level. The friction coefficient under grounded ice and the ice rigidity of floating ice shelves are estimated based on an inverse method (Morlighem et al.2010). To minimize misfit between modeled and observed ice velocities, the surface ice velocity of Rignot (2017) is used. The ice rigidity under grounded ice is estimated using the temperature–rigidity relationship (Cuffey and Paterson2010, p. 72–77).

We use an enthalpy model that considers the transition between cold and temperate ice and the conservation of the total energy balance (Aschwanden et al.2012; Seroussi et al.2013; Kleiner et al.2015). Here, the enthalpy model is referred to as the thermal model and assumes that the ice is in a thermal steady state:

(3) 0 = - v E + ϕ i + k i c i ρ i E , if  E < E s k T pmp + k 0 E , if  E E s ,

where v=(vx,vy,vz) is the ice velocity vector, E is the enthalpy, ϕi is the internal deformation heat, Es is the enthalpy of pure ice, k=(1-ω)ki+ωkw is the mixture thermal conductivity (with ω representing water content), ki and kw are the thermal conductivity of pure ice and liquid water, k0 is a small positive constant (Aschwanden et al.2012), ci is the heat capacity of ice, and Tpmp is the pressure melting point of ice.

The surface temperature is constrained using mean 2 m air temperature data from ERA-Interim, which assimilated the recent atmospheric conditions from 1979 to 2018 with a 0.125° × 0.125° resolution (Dee et al.2011). At the bottom, we impose a Neumann boundary condition with a heat flux from GHF and frictional heating. The basal temperature under floating ice shelves is set to the pressure melting point. An anisotropic streamline upwind Petrov–Galerkin (SUPG) method is adopted since it is more accurate than the original SUPG scheme, which is sensitive to low aspect ratios between the horizontal and vertical resolution meshes (Rückamp et al.2020). The stress balance and thermal state are closely coupled because the internal deformation and frictional heat from the stress balance affect the thermal model. In turn, the ice rigidity inferred from the thermal model influences the stress balance model. To capture this coupling and to reach thermomechanical consistency, we iterate 10 times by solving iteratively the stress balance and thermal model until it reaches convergence. The convergence is reached when the difference in mean basal temperature is lower than 0.5 °C between two consecutive iterations.

We use the surface elevation from the Reference Elevation Model of Antarctica (REMA; Howat et al.2019). The bed geometry is from BedMachine version 1 (Morlighem et al.2020), which used the mass conservation method to generate the bed geometry in fast-flow regions and streamline diffusion in slow-flow regions (Morlighem et al.2010).

2.2 Vertical velocities

We compute the thermal state of the ice sheet using two different vertical velocity profiles: (1) vertical velocity computed by solving for incompressibility while accounting for the inferred basal sliding (hereafter IVz) and (2) the equation of incompressibility of ice while not allowing basal sliding when surface ice velocities are below 10 m yr−1 (hereafter IVz-nosliding). In other words, IVz-nosliding ignores the inferred basal sliding velocities from the initial inversion and assumes that the bed is frozen when surface velocities are <10 m yr−1.

For IVz and IVz-nosliding, we recover the vertical velocity from the continuity equation as follows:

(4) v z z = v z b + b z - v x x - v y y d z .

For IVz-nosliding, we set vx,yb=0, while for IVz the basal vertical velocity is set as

(5) v z b = v x b b x + v y b b y - M ˙ b ,

where M˙b is the basal melting rate (in m yr−1 ice equivalent).

2.3 Geothermal heat flux

We compare four different geothermal flux datasets: Shapiro and Ritzwoller (2004) (SR), which used a seismic model to extrapolate heat-flow measurements; (2) Maule et al. (2005) (Maule), which used a magnetic model with satellite magnetic data; (3) An et al. (2015) (An), which used a crust–lithosphere temperature model; and (4) Martos et al. (2017) (Martos), which inferred the GHF by compiling aeromagnetic data. The mean GHF on grounded ice is 60.78 mW m−2 for SR, 65.61 mW m−2 for Maule, 54.66 mW m−2 for An, and 65.49 mW m−2 for Martos.

2.4 Borehole temperature measurements

To validate the thermal models, we compile all 15 available borehole temperature profiles listed in Table 1. The 10 boreholes in the West Antarctic Ice Sheet region are drilled at Whillans Ice Stream (WIS), Bindschadler Ice Stream (BIS), Engelhardt Ice Ridge (ER), Kamb Ice Stream (KIS), Raymond Ice Ridge (RR), Unicorn (UC), Alley Ice Stream (AIS), and Siple Dome (SD) (Engelhardt2004a) (Fig. 1b). Here, we use borehole names from Engelhardt (2004b): ER-1996-12, SD-1997-1, RR-1997-42, KIS-1996-2, KIS-2000-1,2, UC-1993-11, UC-1993-14, AIS/WIS-1991-1, AIS/WIS-1995-4,7, and BIS-1998-4,5.

Figure 1(a) Borehole locations with temperature measurements overlaid over ice velocity (Rignot2017). The dashed black box shows the location of (b). The solid black box in (a) indicates each basin from Jourdain et al. (2020), and each number indicates each basin number. We use a different symbol for each borehole based on the shape of its temperature profile (blue dots and magenta triangles indicate linear profiles and concave profiles, respectively). The gray contours indicate surface elevations, with dashed lines for every 500 m and solid lines for every 1000 m. The white-dot contours indicate regions where ice velocity is 10 m yr−1. (b) Enlargement of borehole locations at West Antarctica overlaid over the ice velocity. The borehole names are abbreviated: WIS, Whillans Ice Stream; BIS, Bindschadler Ice Stream; ER, Engelhardt Ice Ridge; KIS, Kamb Ice Stream; RR, Raymond Ice Ridge; UC, Unicorn; AIS, Alley Ice Stream; and SD, Siple Dome.

Hondoh et al. (2002)Engelhardt (2004b)Engelhardt (2004b)Yang et al. (2018)Engelhardt (2004b)Engelhardt (2004b)Slawny et al. (2014)Van Ommen et al. (1999)Zagorodnov et al. (2012)Engelhardt (2004b)Engelhardt (2004b)Engelhardt (2004b)Engelhardt (2004b)Engelhardt (2004b)Engelhardt (2004b)

Table 1Summarized information for each borehole. The observed ice velocity is from Rignot (2017), and the number in parentheses indicates error in magnitude of ice velocity. The date refers to when the boreholes were drilled.

Download Print Version | Download XLSX

Since the vertical distance between the temperature measurements along the borehole profile and the triangular mesh is not uniform, we calculate a weighted absolute misfit between the modeled and the measured temperatures (or modeled ice surface velocities) when evaluating the thermal model's performance:

(6) misfit = i = 1 n obs w i Y i mod - Y i obs ,

where nobs is the number of measured points at each borehole (or the number of observed ice velocities), i indicates the index of the specific measured elevation (or index of the ice velocity area), wi is a weight calculated from the ratio of a specific measured point’s occupying length to the total measured length (or ratio of the measured area to the total area), and Yi is the temperature (or ice velocity magnitude). The subscripts “obs” and “mod” indicate the observed and modeled variables, respectively.

Since the ice thickness and the surface temperature of the ice flow model are not always exactly consistent with the observed borehole data, we make adjustments using an exponential decaying correction following Pattyn (2010):

(7) X corr = X + X 0 - X exp - x - x w 2 + y - y w 2 σ ,

where (xw,yw) is the location of the borehole, X0 is the observed quantity, and X is the model ice thickness or surface temperature. The surface temperature at each borehole (except for SD-1997-1, RR-1997-42, UC-1993-11, AIS/WIS-1988-1, and BIS-1998-4,5) is corrected where the given climatological temperature is relatively higher than the observed surface temperature (Table 1). Xcorr is the corrected data, and σ is the radius of influence, which is set here to 50 km. The geometry from BedMachine is constrained using radar-derived ice thickness measurements, except for at Dome Fuji and Law Dome, for which few thickness measurements were available. These two locations are the only places where an ice thickness correction is applied so that the ice thickness is 3090 m at Dome Fuji and 1220 m at Law Dome.

3 Results

3.1 Model experiments

To estimate the ice temperature of the entire Antarctic continent, we perform eight different experiments by combining two different vertical velocity profiles (IVz and IVz-nosliding) and four different GHF datasets. Table 2 shows the weighted absolute misfits between the modeled and observed surface ice velocities across the entire domain. The mean ice surface velocity misfit is 12.5 m yr−1 for the IVz group and 19.5 m yr−1 for the IVz-nosliding group. The standard deviation in the ice velocity misfit is 0.09 m yr−1 for the IVz group and 0.35 m yr−1 for the IVz-nosliding group.

Table 2Experimental design for eight simulations using different vertical velocities and geothermal heat fluxes. The value between parentheses under each experiment represents the weighted absolute misfit between observed and modeled surface ice velocity across the entire domain.

Download Print Version | Download XLSX

Figure 2 displays the measured and modeled vertical profiles of the ice temperature at the 15 borehole locations (see Fig. S1). The measured vertical profiles of the borehole temperatures, marked as dashed black lines in Fig. 2, can be categorized into two groups based on temperature profile shapes. One group exhibits concave profiles, for which the vertical advection toward the bed dominates, while the other group has a more linear shape, for which vertical diffusion dominates. Dome Fuji, SD-1997-1, RR-1997-42, ER-1996-12, and Styx Glacier at slow-flow regions show diffusion-dominant temperature profiles compared to the WAIS Divide, Bruce Plateau, Law Dome, KIS-1996-2, KIS-2000-1,2, UC-1993-11, and UC-1993-14, where the advection toward the bed dominates. Note that AIS/WIS-1991-1, AIS/WIS-1995-4,7, and BIS-1998-4,5 are located in regions with comparatively high ice velocity compared to other boreholes and have concave temperature profiles. To clearly define this specific fast-flow region, we refer to AIS/WIS and BIS as the Siple Coast fast-flow region.

Figure 2Observed and modeled vertical temperature profiles from eight different experiments at 15 borehole locations. Blue and red boxes indicate slow-flow and Siple Coast fast-flow regions, respectively. The bottom elevation at each borehole is set considering the ice thickness, as listed in Table 1. An asterisk next to the borehole name indicates that the drilling reaches the bedrock.


3.2 Borehole temperature profiles

To provide a quantitative comparison between the modeled and observed borehole temperatures, a weighted absolute misfit is calculated (Fig. 3). The average temperature misfit values for IVz and IVz-nosliding are 5.64 °C and 3.61 °C, respectively, and 1.69 °C and 2.50 °C for slow-flow and Siple Coast fast-flow regions, respectively. The temperature misfit value of the IVz-nosliding group is lower than that of the IVz group; however, the misfit temperatures in the Siple Coast fast-flow regions for IVz and IVz-nosliding are not exactly the same. The spread in misfits among the different vertical velocity schemes is larger than the one obtained when varying GHF. This shows that the difference in GHF has a limited influence on estimating the overall temperature profiles, while the choice of vertical velocities has a stronger impact. Both the IVz and IVz-nosliding groups demonstrate good performance in Siple Coast fast-flow regions, such as AIS/WIS and BIS. In the case of slow-flow regions, the thermal model's performance for the IVz-nosliding group is improved compared to the IVz group, and the model produces a reduced temperature misfit. However, none of the experiments successfully reproduce the temperature profiles at KIS boreholes, where the ice has been stagnant since around 1850 CE (Alley et al.1994; Joughin and Tulaczyk2002). This history cannot be captured by our thermal steady-state assumption. A more detailed description of misfit values for each borehole can be found in the next section.

Figure 3Weighted absolute misfit between observed and modeled borehole temperatures according to each experiment. The absolute temperature misfit is truncated over 5 °C.


First, we focus on the three borehole profiles: SD, RR, and Dome Fuji. They all have linear temperature profiles, which are rarely observed in general borehole temperature profiles. SD and RR are adjacent to each other, but measurements of borehole temperatures at RR are limited to the top few hundred meters. Dome Fuji is located in the interior of the ice sheet. For these boreholes, the IVz group does not capture the linear shape of the temperature profiles. The IVz-nosliding group at these boreholes has a misfit value within 2 °C, which is lower than that of the IVz group (Fig. 3). The basal temperatures from the IVz-nosliding group reach the pressure melting point at SD, RR (Engelhardt2004b), and Dome Fuji. In the case of An, the GHF at each borehole is 40.1 mW m−2 for Dome Fuji, 64.9 mW m−2 for SD, and 65.3 mW m−2 for RR, which is lower than the values from other GHF sources. The basal modeled temperature for An is the lowest and does not reach pressure melting point. The depth-averaged vertical velocity at Dome Fuji is −0.14 m yr−1 for IVz (where a negative value means the vector is oriented downward), which is a higher value than that of IVz-nosliding (−0.01 m yr−1) (Table 3). The depth-averaged vertical velocities of IVz at SD and RR are also higher than those of IVz-nosliding. This suggests a larger advection toward the ice sheet base in the IVz group, where downward heat advection is more dominant than the diffusion process and leads to a colder basal temperature compared to in the IVz-nosliding group.

Table 3Depth-averaged vertical velocity for each experiment at each borehole. Positive values indicate upward advection.

Download Print Version | Download XLSX

At the borehole of Styx Glacier, both IVz and IVz-nosliding groups display similar average misfit values of ∼0.64 °C and  0.40 °C, which show good agreement with the observed temperature profiles. The drilling depth of Styx Glacier is about 210.5 m (Yang et al.2018), and the ice thickness measured with a ground-penetrating radar survey is about 550 m (Hur2013). While we cannot definitively confirm the basal condition from observations, the thermal model results suggest that none of the glaciers reach the melting point.

The UC boreholes are located at an area of stagnant ice and have a relatively high basal temperature gradient compared to the other adjacent boreholes, such as the AIS/WIS boreholes (Engelhardt2004b). The mean GHF in the UC region is approximately 81.4 mW m−2 for SR, 86.5 mW m−2 for Maule, 62.8 mW m−2 for An, and 95.6 mW m−2 for Martos. The current modeled temperature profiles at UC-1993-11 and UC-1993-14 agree well with the measured temperature regardless of the choice of GHFs. The misfit value for the modeled and observed temperatures from the IVz-nosliding group is lower than that of the IVz group. In addition, the misfit of UC-1993-14 for IVz-nosliding is lower than that of UC-1993-11 (Fig. 3). UC-1993-14 is located in a slow region; however, UC-1993-11 is adjacent to the shear margin of the AIS ice stream, which induces a sharp transition in the basal velocity constraints for the IVz-nosliding group where the ice velocity crosses 10 m yr−1. While the IVz-nosliding group better captures the observed temperature profiles for UC-1993-14, this is not the case for UC-1993-11.

The modeled basal temperature at the WAIS Divide reaches the pressure melting point for only the SR and Martos IVz groups. The GHF is approximately 112.6 mW m−2 for SR and 141 mW m−2 for Martos; these values are higher than those of the other two GHF datasets, which are 60.3 mW m−2 for Maule and 68.9 mW m−2 for An. The basal melting rate of the IVz-nosliding group is 7.9 mm yr−1 for SR, 2.5 mm yr−1 for Maule, 3.4 mm yr−1 for An, and 10.9 mm yr−1 for Martos. GHF estimations in previous studies are 113.3±16.9 mW m−2 from Talalay et al. (2020) and 90.5 mW m−2 from Mony et al. (2020). The thickness at WAIS Divide is 3455 m (WAIS Divide Project Members2013). However, the drilling depth is 3405 m (Slawny et al.2014) and does not reach the bed, so we do not know the rate of basal melting. According to Talalay et al. (2020), the estimated basal temperature at WAIS Divide reaches the pressure melting point, and the basal melting rate is about 3.7±1.7 mm yr−1. All experiments show reasonably good agreement in terms of the shape of the observed borehole temperature profile at WAIS Divide regardless of the choice of GHF. The average misfit value of the borehole temperature for IVz is 2.90 °C and is better than that of IVz-nosliding (Fig. 3).

At Law Dome, the misfit between the observed and modeled temperatures is 2.9 °C and 1.5 °C for the IVz and IVz-nosliding groups, respectively (Fig. 3). A primary difference between IVz and IVz-nosliding is the depth-averaged vertical velocity, the value of which is −1.5 m yr−1 for the IVz group and −0.1 m yr−1 for the IVz-nosliding group (Table 3). In the Law Dome case, we confirm that the use of IVz-nosliding improves the model’s vertical temperature profile (Fig. 2).

The observed ice velocity at Bruce Plateau is 49 m yr−1 according to Rignot (2017), which is higher than the previously reported value of 10±4 m yr−1 (Zagorodnov et al.2012). We find that none of the modeled thermal profiles can reproduce the upper part of the observed ice temperature that captured the colder surface temperature of past climates (Zagorodnov et al.2012). The mean vertical velocity for the IVz group is −5.24 m yr−1 and is −2.71 m yr−1 for the IVz-nosliding group; these values indicate high vertical advection toward the bottom.

With the exception of ER-1997-12, neither the IVz group nor the IVz-nosliding group captures the observed temperature profiles at the KIS boreholes. All modeled temperature profiles exhibit a convex shape (Fig. 2). At ER-1997-12, the mean misfit between the modeled and observed temperature is 2.7 °C for the IVz-nosliding group and 5.7 °C for the IVz group (Fig. 3).

The AIS/WIS and BIS boreholes are located in the fast-flow regions of the Siple Coast, where the ice velocities are 365 m yr−1 for AIS/WIS-1991-1, 379 m yr−1 for AIS/WIS-1995-4,7, and 376 m yr−1 for BIS-1998-4,5 from Rignot (2017). The average misfit values of the IVz group are 1.38 °C for AIS/WIS-1988-1, 2.16 °C for AIS/WIS-1995-4,7, and 0.86 °C for BIS-1998-4,5 (Fig. 3). In these regions, both IVz and IVz-nosliding allow for basal sliding. However, there are differences in misfit values between IVz groups and IVz-nosliding groups at AIS/WIS. The reason for these differences is that the modeled ice velocities of IVz-nosliding in the AIS/WIS region are slower than the ones from IVz because it is a narrow ice stream that is influenced by the no-sliding constraint along its sides, resulting in higher misfit values for the IVz-nosliding group compared to the IVz group. The misfit between the modeled and observed temperatures at BIS is lower than that of AIS/WIS. In fast-flow regions, the advection, estimated through the stress balance of the ice and the ice incompressibility, plays a more important role in the thermal model than diffusion does. In these advection-dominated regions, the temperature is sensitive to bed geometry. A likely explanation for the difference in the misfits between the BIS and AIS/WIS regions is that the bed geometry in the BIS region is constructed using a mass conservation approach, which relies on the equation of ice incompressibility. In contrast, the bed geometry in the AIS/WIS region is constructed using the stream diffusion method, similar to kriging (Fig. S5). This suggests that enhancement in the quality of the geometry and utilizing the mass conservation method in the Siple Coast fast-flow regions would improve the estimation of the vertical velocity by the IVz equation with sliding and improve the overall performance of the thermal model. The AIS/WIS-1995-4,7 borehole is located at the center of the ice stream, whereas AIS/WIS-1988-1 is relatively near the margin of the ice stream. Although the bed geometry at AIS/WIS is constructed using the kriging method, IVz reproduces the temperature profile reasonably well at the center of fast ice flow regions.

3.3 Subglacial conditions

Figure 4a and b show the mean and the standard deviation of the basal temperature distribution for the eight experiments. The mean basal temperature at the main ice trunk, where the ice primarily discharges into the ocean, reaches the ice pressure melting point. The standard deviation of the basal temperature is higher in the internal ice compared to in the peripheral regions. In the case of IVz-nosliding, constraining the basal velocity to zero in slow-flow regions leads to a warmer basal temperature distribution compared to in the IVz group. In slow-flow regions, the basal temperature of the IVz group shows a notable difference depending on the choice of GHF. The modeled basal temperatures in the Maule and Martos experiments, which have higher mean GHF values (Table 4), are warmer than those in the SR and An experiments, as expected. The mean GHF of An is the lowest compared to the other GHFs, and, therefore, the basal temperature at each borehole modeled with the An GHF is lower than those of the other GHFs.

Figure 4(a) Mean and (b) standard deviation of the basal temperature distribution from eight experiments. (c–j) Basal temperature distribution for each experiment. The temperature legend is truncated below −10 °C. A hatched region with a white line indicates that the basal temperature of ice reaches the pressure melting point.

All the experiments generally indicate that most of the regions experiencing basal melting are concentrated in fast-flow regions, where basal frictional heat is significant and provides enough heat for the ice to reach the pressure melting point (Fig. 4). Since IVz-nosliding displays lower vertical advection than that of IVz, the basal temperature of the IVz-nosliding group in slow-flow regions is warmer than that of the IVz group (Fig. 4c-j).

Figure 5Upper panels (a–d) are the geothermal heat flux distributions of each source. Middle panels (e–l) are the basal melting rate distributions, with the value at the bottom left indicating the total grounded ice melting volume for each experiment. The basal melting rate exceeding 50 mm yr−1 is truncated. Lower panels (m–q) are the difference in the basal melting rate between IVz-nosliding and IVz for each geothermal heat flux. The green crosses on the geothermal heat flux map indicate the borehole location. The color maps for the geothermal heat flux and the difference in basal melting rates are from Crameri et al. (2020).

The mean total grounded ice melting volume is 26.62 Gt yr−1 for the IVz group and 29.77 Gt yr−1 for the IVz-nosliding group (Table 4). The total grounded ice melting volume for the IVz-nosliding group is 3.15 Gt yr−1 higher than that of the IVz group. Compared to IVz, IVz-nosliding suggests a mean total basal melting volume increase of 1.89 Gt yr−1 (60 %) and 1.26 Gt yr−1 (40 %) in the slow-flow and fast-flow regions, respectively. The total grounded ice melting volume is proportional to the GHF magnitude. Each basin displays significant differences in terms of the grounded ice melting volume depending on the GHF source. Note that the GHF from An, which is the lowest value among all GHFs, shows the lowest total grounded ice melting volume.

Table 4Grounded ice basal melting volumes of eight experiments at each basin (Fig. 1) and the total grounded ice melting volume and the total grounded melting fraction corresponding to each experiment.

Download Print Version | Download XLSX

4 Discussion

Previous studies that have successfully reproduced borehole temperature profiles using one-dimensional thermal analytical solutions have been limited to slow-flow regions (Joughin et al.2003; Mony et al.2020; Talalay et al.2020). These studies have demonstrated good agreement between modeled and observed temperatures, which is expected given their simplicity and the tunability of the analytical solutions. One important tunable parameter is the analytical vertical velocities, which rely on ice surface mass balance (Hindmarsh1999; Joughin et al.2003; Talalay et al.2020). The choice of vertical velocity is a key factor in reproducing borehole temperature profiles. Uncertainties in the GHF have also been identified as a major factor in reproducing observed borehole temperature profiles (Talalay et al.2020; Mony et al.2020). On the other hand, some other studies have shown that uncertainties in the GHF have little influence on model performance in terms of ice dynamics (Larour et al.2012a; Smith-Johnsen et al.2020a) and simulating future projections (Schlegel et al.2018; Smith-Johnsen et al.2020b; Seroussi et al.2013). Therefore, to test other factors, such as different GHFs and vertical velocities, which may affect the calculation of borehole temperatures, we use a three-dimensional thermomechanical model in order to account for both horizontal and vertical advection. We compare our calculated temperatures to observed borehole profiles in both fast-flow and slow-flow regions.

In slow-flow regions, we find that IVz-nosliding experiments show reasonably good agreement with the observed borehole temperature profiles. However, the three-dimensional thermal model occasionally estimates convex temperature profiles, which are not consistent with the observations, such as the KIS boreholes. Compared to other boreholes, the ice velocities at KIS and ER gradually decrease from upstream to downstream and coincide with the presence of a basal ridge (Price et al.2001; Ng and Conway2004) (see also Fig. S2). In the past, the KIS and ER regions experienced faster ice flow, and the ice stream started to stagnate around 1850 CE (Alley et al.1994; Joughin and Tulaczyk2002). There are hypotheses explaining the stagnation in the KIS region: the water piracy hypothesis (Alley et al.1994) or the removal of basal water contributing to the loss of lubrication (Tulaczyk et al.2000; Bougamont et al.2003). This history results in colder temperatures in the upper part of the ice column, which contains ice that was deposited farther upstream where the surface temperature was lower than it is at the current location of the boreholes. This ice was then transported downstream to the current location (Hills et al.2023). In model experiments, Bougamont et al. (2015) revealed changes in the tributaries at KIS and ER using a plastic till deformation friction law including a simple subglacial hydrology model. In contrast, we employ the Budd-type friction law and assume the effective pressure fully connected to the ocean part, not including changes in the effective pressure. The variation in effective pressures also changed the basal ice velocity in the Budd-type friction law. In addition, a selection of other types of friction laws, including Weertman (Weertman1974), Schoof (Schoof2005), and Coulomb (Tsai et al.2015) types, also influences the initialization and future fate of ice (Brondex et al.2017, 2019). Further investigation is required, such as the application of other types of friction laws or initialization with paleo spin-up, to better understand temperature profiles.

In Siple Coast fast-flow regions, Joughin et al. (2004) utilized a thermal model with vertical velocity derived from an analytical solution, which reproduced the observed borehole temperature profile of BIS-1998-4,5 with good agreement (UpD in Joughin et al.2004). Here, we also find that the modeled temperature using a vertical velocity based on the equation of incompressibility without any constraint or tunable parameter also agrees well with the observed temperatures in this sector.

The total grounded ice melting volume for both the IVz and IVz-nosliding groups falls within the range reported by previous studies. It is lower than 65 Gt yr−1 from Pattyn (2010) and higher than 14.7 Gt yr−1 from Llubes et al. (2006). In the study by Joughin et al. (2009), a homogeneous GHF value of 70 mW m−2 was adopted, which is similar to the mean GHFs from Maule (66.95 mW m−2) and An (67.15 mW m−2) at basin 10, which includes the Pine Island and Thwaites glaciers (see basin in Fig. 1). However, the total grounded ice melting mass estimated by Joughin et al. (2009) of 5.2 Gt yr−1 is higher than that of the IVz group (average value of Maule and An) of 3.5 Gt yr−1 and the IVz-nosliding group (average value of Maule and An) of 3.8 Gt yr−1. The thermal models have been employed to explore the thermal regime of ice and to estimate basal melting rates beneath grounded ice. In the thermal model's advection term, the horizontal components of the ice velocity are estimated using the stress balance equations, whereas the vertical velocity is recovered from the ice incompressibility. Under kriging-based bed topography, the vertical velocity in fast-flow regions leads to large flux divergences (Seroussi et al.2011). In contrast, mass-conservation-based bed geometries, such as BedMachine (Morlighem et al.2017, 2020), preserve low flux divergence. We confirm that in areas where the bed geometry was inferred from mass continuity, the more accurate estimates of the vertical velocity provide a viable input for estimates of temperature profiles, for example, in the Siple Coast fast-flow regions. Additionally, we expect this study to provide a reliable understanding of temperature profiles in the other fast-flow regions generated with mass conservation. We should highlight that the good agreement between modeled and observed temperatures in fast-flow regions does not guarantee that the magnitude of basal melting volume is accurate, as it depends on both geothermal heat fluxes and frictional heat.

We find that the impact of using different GHF fields has only a modest influence on the ice temperature field and the total grounded ice basal melting volume. Under these circumstances, our results reveal that the shapes of the borehole temperature profile are less sensitive to the current estimated GHFs than previously reported. It is also worth noting that the initialization with the GHF from An results in underestimated basal temperatures and a lower total grounded ice melting volume due to an excessively low GHF value compared to other datasets.

The IVz-nosliding experiment has the advantage of better simulating the vertical temperature profiles in slow-flow regions compared to IVz. However, it tends to produce large discrepancies between modeled and observed surface ice velocities (Fig. S3). For instance, the An IVz-nosliding thermal model experiments exhibit the largest misfits in ice velocity among all the experiments, as the lowest value of average GHF leads to relatively high ice rigidity that perturbs the ice flow in the slow-flow regions. In contrast, the IVz experiment shows relatively smaller misfit values in surface ice velocity because sliding compensates for the underestimated internal deformations in the slow-flow region. In general, we find that IVz leads to a higher depth-averaged ice rigidity compared to IVz-nosliding in slow regions due to presence of colder ice temperatures (Fig. S4). Higher ice rigidity causes ice to deform less vertically, through vertical shear, and the surface ice velocity with no sliding cannot reproduce the observed surface velocities. In other words, the surface ice velocity of IVz-nosliding shows a larger ice velocity misfit compared to that of the IVz group, because the basal velocities are constrained to zero and cannot compensate for the high velocity misfit. Furthermore, the adoption of no sliding in specific regions results in a sharp transition zone in ice rigidity, B. This occurs because the basal velocity near the transition zone does not smoothly change from no sliding to sliding (Fig. S4). Therefore, additional work is required to address and resolve the smooth transition between no sliding and sliding.

In slow-flow regions, a competition between vertical diffusion and advection determines the shape of the temperature profiles and the bottom temperatures. In the IVz experiments, the boundary condition for basal vertical velocity is recovered with the gradient of the bed geometry and the basal melting volume. This approach provides relatively high vertical velocities in slow-flow regions. The vertical velocities are not always in agreement with the analytical expression of vertical velocities assuming a stationary bed and no sliding. As the depth-averaged vertical velocity of IVz is higher than that of IVz-nosliding, cold surface temperatures can be more effectively transferred deeper into the ice column.

The surface temperature of ice would be one of the factors in considering the boundary conditions of a thermal model. While ERA5 (Hersbach et al.2023), RACMO2.3p2 forced with ERA5 (van Wessem et al.2023), and MERRA-2 (Global Modeling and Assimilation Office (GMAO)2015) are the recent reanalysis datasets, they display some discrepancies between the climatological mean 2 m air temperature (1980–2018) and the observed surface temperature at each borehole (Fig. S6). For the comparison with a different version of ECMWF (European Centre for Medium-Range Weather Forecasts) reanalysis data, we perform experiments in the same manner, utilizing 2 m air temperature from ERA5. These results display no significant differences compared to experiments using ERA-Interim (Fig. S7). However, in the case of SD, RR, and AIS/WIS (only for the IVz-nosliding case), they display slight discrepancies in surface temperature, leading to shifts in the modeled temperature profiles. In fact, the improvement in surface temperature and the accurate correction would bring the modeled temperatures into closer agreement with observations.

Finally, borehole temperatures have a long-term memory of past climate air temperatures and are a good proxy for reconstruction over a few hundred years or longer using inverse modeling (Nagornov et al.2001; Zagorodnov et al.2012). This history is not accounted for in this study, as we assumed thermal steady state using current climatological information. Despite this strong limitation, we find that this approach provides temperature profiles that are in good agreement with observations.

5 Conclusions

In this study, we used a three-dimensional thermomechanical model of Antarctica with different sources of GHF and vertical velocity fields to reproduce different thermal states of the Antarctic ice sheet, and we compared the results to 15 in situ measured borehole temperature profiles in slow-flow and fast-flow regions. Comparing the modeled to the measured borehole temperature profiles, we confirm that the vertical ice velocity based on the equation of incompressibility (IVz) is suitable for fast-flow regions, such as BIS, where the bed geometry is constructed using the mass conservation method, while an IVz that ignores basal sliding (IVz-nosliding) performs better in slow-flow regions. Our results show that the vertical temperature profile is more sensitive to the vertical velocity. In addition, the basal conditions, such as temperature and melting rate, are sensitive to both GHF and the vertical velocity field. The total grounded ice melting volume and basal temperature are proportional to the magnitude of the average GHF values for the same vertical velocity method. Finally, constraining the basal velocity to zero in slow-flow regions is a reasonable assumption and leads to a more realistic temperature profile.

Code and data availability

ISSM is open source and can be downloaded at (last access: 20 February 2024; Larour et al.2012b). The Law Dome temperature profile by Van Ommen et al. (1999) is available at (van Ommen2023). The Dome Fuji temperature profile is available in Hondoh et al. (2002). The Styx Glacier borehole temperature profile by Yang et al. (2018) can be obtained by personal communication with the author. The Bruce Plateau temperature profile is available in Zagorodnov et al. (2012). The WAIS Divide borehole temperature profile by Cuffey and Clow (2014) is available at The SD, RR, UC, ER, KIS, AIS/WIS, and BIS borehole temperature profiles by Engelhardt (2004b) are available at (Engelhardt2013). The GHF map by Shapiro and Ritzwoller (2004) and Maule et al. (2005) is available from ALBMAP v1.0 (; Le Brocq et al.2010). The GHF map by An et al. (2015) is available at The GHF map by Martos et al. (2017) is available at (Martos2017). The 2 m air temperature by Dee et al. (2011) is available at The ice velocity map by Rignot (2017) is available at The bed geometry of Antarctica by Morlighem et al. (2020) is available at (Morlighem2019). The surface elevation map by Howat et al. (2019) is available at The results of the gridded basal temperature field are available from the KDPC (Korea Polar Data Center) (; Park et al.2023).


The supplement related to this article is available online at:

Author contributions

IWP designed and performed experiments with inputs from EKJ, MM, and KKL. All authors participated in the writing of the article.

Competing interests

The contact author has declared that none of the authors has any competing interests.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors.


The authors are grateful to the editor, Benjamin Smith, for handling our article and to the two reviewers, Tyler Pelle and an anonymous reviewer, for their helpful comments on the article.

Financial support

This research has been supported by the Korea Institute of Marine Science and Technology Promotion (KIMST) funded by the Ministry of Oceans and Fisheries (grant nos. RS‐2023‐00256677 and PM24020). Mathieu Morlighem has been supported by the PROPHET project, a component of the International Thwaites Glacier Collaboration (ITGC) funded by the National Science Foundation (NSF; grant no. 1739031).

Review statement

This paper was edited by Benjamin Smith and reviewed by Tyler Pelle and one anonymous referee.


Alley, K. E., Scambos, T. A., Siegfried, M. R., and Fricker, H. A.: Impacts of Warm Water on Antarctic Ice Shelf Stability through Basal Channel Formation, Nat. Geosci., 9, 290–293,, 2016. a

Alley, R. B., Anandakrishnan, S., Bentley, C. R., and Lord, N.: A Water-Piracy Hypothesis for the Stagnation of Ice Stream C, Antarctica, Ann. Glaciol., 20, 187–194,, 1994. a, b, c

An, M., Wiens, D. A., Zhao, Y., Feng, M., Nyblade, A., Kanao, M., Li, Y., Maggi, A., and Lévêque, J.: Temperature, Lithosphere‐asthenosphere Boundary, and Heat Flux beneath the Antarctic Plate Inferred from Seismic Velocities, J. Geophys. Res.-Sol. Ea., 120, 8720–8742, 2015 (data available at:, last access: 20 February 2024). a, b, c

Anker, P. G. D., Makinson, K., Nicholls, K. W., and Smith, A. M.: The BEAMISH Hot Water Drill System and Its Use on the Rutford Ice Stream, Antarctica, Ann. Glaciol., 62, 233–249,, 2021. a

Aschwanden, A., Bueler, E., Khroulev, C., and Blatter, H.: An Enthalpy Formulation for Glaciers and Ice Sheets, J. Glaciol., 58, 441–457,, 2012. a, b

Augustin, L. and Antonelli, A.: The EPICA deep drilling program, Mem. Natl. Inst. Polar Res., Special Issue, 226–244, 2002. a

Bougamont, M., Tulaczyk, S., and Joughin, I.: Response of Subglacial Sediments to Basal Freeze‐on 2. Application in Numerical Modeling of the Recent Stoppage of Ice Stream C, West Antarctica, J. Geophys. Res.-Sol. Ea., 108, 2223,, 2003. a

Bougamont, M., Christoffersen, P., Price, S. F., Fricker, H. A., Tulaczyk, S., and Carter, S. P.: Reactivation of Kamb Ice Stream Tributaries Triggers Century‐scale Reorganization of Siple Coast Ice Flow in West Antarctica, Geophys. Res. Lett., 42, 8471–8480,, 2015. a

Brondex, J., Gagliardini, O., Gillet-Chaulet, F., and Durand, G.: Sensitivity of Grounding Line Dynamics to the Choice of the Friction Law, J. Glaciol., 63, 854–866,, 2017. a

Brondex, J., Gillet-Chaulet, F., and Gagliardini, O.: Sensitivity of centennial mass loss projections of the Amundsen basin to the friction law, The Cryosphere, 13, 177–195,, 2019. a, b

Budd, W. F., Keage, P. L., and Blundy, N. A.: Empirical Studies of Ice Sliding, J. Glaciol., 23, 157–170, 1979. a

Crameri, F., Shephard, G. E., and Heron, P. J.: The Misuse of Colour in Science Communication, Nat. Commun., 11, 5444,, 2020. a

Cuffey, K. M. and Clow, G. D.: Temperature Profile of the West Antarctic Ice Sheet Divide Deep Borehole, U.S. Antarctic Program (USAP) Data Center [data set],, 2014. a

Cuffey, K. M. and Paterson, W. S. B.: The Physics of Glaciers, Academic Press, ISBN 0-08-091912-X, 2010. a

Dahl-Jensen, D., Morgan, V. I., and Elcheikh, A.: Monte Carlo Inverse Modelling of the Law Dome (Antarctica) Temperature Profile, Ann. Glaciol., 29, 145–150, 1999. a

DeConto, R. M. and Pollard, D.: Contribution of Antarctica to Past and Future Sea-Level Rise, Nature, 531, 591–597,, 2016. a

Dee, D. P., Uppala, S. M., Simmons, A. J., Berrisford, P., Poli, P., Kobayashi, S., Andrae, U., Balmaseda, M. A., Balsamo, G., Bauer, P., Bechtold, P., Beljaars, A. C. M., van de Berg, L., Bidlot, J., Bormann, N., Delsol, C., Dragani, R., Fuentes, M., Geer, A. J., Haimberger, L., Healy, S. B., Hersbach, H., Hólm, E. V., Isaksen, L., Kållberg, P., Köhler, M., Matricardi, M., McNally, A. P., Monge-Sanz, B. M., Morcrette, J.-J., Park, B.-K., Peubey, C., de Rosnay, P., Tavolato, C., Thépaut, J.-N., and Vitart, F.: The ERA-Interim Reanalysis: Configuration and Performance of the Data Assimilation System, Q. J. Roy. Meteor. Soc., 137, 553–597,, 2011 (data available at:, last access: 20 February 2024). a, b

Derkacheva, A., Mouginot, J., Millan, R., Maier, N., and Gillet-Chaulet, F.: Data Reduction Using Statistical and Regression Approaches for Ice Velocity Derived by Landsat-8, Sentinel-1 and Sentinel-2, Remote Sens., 12, 1935,, 2020. a

Doyle, S. H., Hubbard, B., Christoffersen, P., Young, T. J., Hofstede, C., Bougamont, M., Box, J. E., and Hubbard, A.: Physical Conditions of Fast Glacier Flow: 1. Measurements From Boreholes Drilled to the Bed of Store Glacier, West Greenland, J. Geophys. Res.-Earth Surf., 123, 324–348,, 2018. a

Engelhardt, H.: Ice Temperature and High Geothermal Flux at Siple Dome, West Antarctica, from Borehole Measurements, J. Glaciol., 50, 251–256,, 2004a. a, b

Engelhardt, H.: Thermal Regime and Dynamics of the West Antarctic Ice Sheet, Ann. Glaciol., 39, 85–92,, 2004b. a, b, c, d, e, f, g, h, i, j, k, l, m, n, o

Engelhardt, H.: Temperature of the West Antarctic Ice Sheet” U.S. Antarctic Program (USAP) Data Center [data set],, 2013. a

Fisher, A. T., Mankoff, K. D., Tulaczyk, S. M., Tyler, S. W., Foley, N., and and the WISSARD Science Team: High Geothermal Heat Flux Measured below the West Antarctic Ice Sheet, Sci. Adv., 1, e1500093,, 2015. a

Fox-Kemper, B., Hewitt, H. T., Xiao, C., Aðalgeirsdóttir, G., Drijfhout, S. S., Edwards, T. L., Golledge, N. R., Hemer, M., Kopp, R. E., Krinner, G., Mix, A., Notz, D., Nowicki, S., Nurhati, I. S., Ruiz, L., Sallée, J.-B., Slangen, A. B. A., and Yu, Y.: Ocean, Cryosphere and Sea Level Change, in: Climate Change 2021: The Physical Science Basis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, edited by: Masson-Delmotte, V., Zhai, P., Pirani, A., Connors, S. L., Péan, C., Berger, S., Caud, N., Chen, Y., Goldfarb, L., Gomis, M. I., Huang, M., Leitzell, K., Lonnoy, E., Matthews, J. B. R., Maycock, T. K., Waterfield, T., Yelekçi, O., Yu, R., and Zhou, B., Cambridge University Press, Cambridge, United Kingdom and New York, NY, USA, 1211–1362,, 2021. a

Gillet-Chaulet, F.: Assimilation of surface observations in a transient marine ice sheet model using an ensemble Kalman filter, The Cryosphere, 14, 811–832,, 2020. a

Gillet-Chaulet, F., Hindmarsh, R. C. A., Corr, H. F. J., King, E. C., and Jenkins, A.: In-Situ Quantification of Ice Rheology and Direct Measurement of the Raymond Effect at Summit, Greenland Using a Phase-Sensitive Radar, Geophys. Res. Lett., 38, L24503,, 2011. a

Gillet-Chaulet, F., Gagliardini, O., Seddik, H., Nodet, M., Durand, G., Ritz, C., Zwinger, T., Greve, R., and Vaughan, D. G.: Greenland ice sheet contribution to sea-level rise from a new-generation ice-sheet model, The Cryosphere, 6, 1561–1576,, 2012. a

Glen, J. W.: The Creep of Polycrystalline Ice, P. Roy. Soc. Lond. A, 228, 519–538, 1955. a

Global Modeling and Assimilation Office (GMAO): MERRA-2 instM_2d_asm_Nx: 2d,Monthly Mean, Single-Level, Assimilation, Single-Level Diagnostics V5.12.4, Greenbelt, MD, USA: Goddard Space Flight Center Distributed Active Archive Center (GSFC DAAC) [data set],, 2015. a

Hersbach, H., Bell, B., Berrisford, P., Biavati, G., Horányi, A., Muñoz Sabater, J., Nicolas, J., Peubey, C., Rozum, I., Schepers, D., Simmons, A., Soci, C., Dee, D., and Thépaut, J.-N.: ERA5 Monthly Averaged Data on Single Levels from 1940 to Present, Copernicus Climate Change Service (C3S) Climate Data Store (CDS) [data set],, 2023. a

Hills, B. H., Christianson, K., Jacobel, R. W., Conway, H., and Pettersson, R.: Radar Attenuation Demonstrates Advective Cooling in the Siple Coast Ice Streams, J. Glaciol., 69, 566–576,, 2023. a

Hindmarsh, R. C.: On the Numerical Computation of Temperature in an Ice Sheet, J. Glaciol., 45, 568–574, 1999. a, b

Hondoh, T., Shoji, H., Watanabe, O., Salamatin, A. N., and Lipenkov, V. Y.: Depth–Age and Temperature Prediction at Dome Fuji Station, East Antarctica, Ann. Glaciol., 35, 384–390, 2002. a, b

Howat, I. M., Porter, C., Smith, B. E., Noh, M.-J., and Morin, P.: The Reference Elevation Model of Antarctica, The Cryosphere, 13, 665–674,, 2019 (data available at:, last access: 20 February 2024). a, b

Hubbard, B., Philippe, M., Pattyn, F., Drews, R., Young, T. J., Bruyninx, C., Bergeot, N., Fjøsne, K., and Tison, J.-L.: High-Resolution Distributed Vertical Strain and Velocity from Repeat Borehole Logging by Optical Televiewer: Derwael Ice Rise, Antarctica, J. Glaciol., 66, 523–529,, 2020. a

Hur, S. D.: Development of Core Technology for Ice Core Drilling and Ice Core Bank, Korea Polar Research Institute, Incheon, 398, 2013. a, b

Joughin, I. and Tulaczyk, S.: Positive Mass Balance of the Ross Ice Streams, West Antarctica, Science, 295, 476–480,, 2002. a, b

Joughin, I., Tulaczyk, S., MacAyeal, D. R., and Engelhardt, H.: Melting and Freezing beneath the Ross Ice Streams, Antarctica, J. Glaciol., 50, 96–108,, 2004. a, b, c

Joughin, I., Tulaczyk, S., Bamber, J. L., Blankenship, D., Holt, J. W., Scambos, T., and Vaughan, D. G.: Basal Conditions for Pine Island and Thwaites Glaciers, West Antarctica, Determined Using Satellite and Airborne Data, J. Glaciol., 55, 245–257,, 2009. a, b

Joughin, I. R., Tulaczyk, S., and Engelhardt, H. F.: Basal Melt beneath Whillans Ice Stream and Ice Streams A and C, West Antarctica, Ann. Glaciol., 36, 257–262,, 2003. a, b

Jourdain, N. C., Asay-Davis, X., Hattermann, T., Straneo, F., Seroussi, H., Little, C. M., and Nowicki, S.: A protocol for calculating basal melt rates in the ISMIP6 Antarctic ice sheet projections, The Cryosphere, 14, 3111–3134,, 2020. a

Khazendar, A., Rignot, E., and Larour, E.: Larsen B Ice Shelf Rheology Preceding Its Disintegration Inferred by a Control Method, Geophys. Res. Lett., 34, L19503,, 2007. a

Kingslake, J., Hindmarsh, R. C., Ahalgeirsdóttir, G., Conway, H., Corr, H. F., Gillet-Chaulet, F., Martín, C., King, E. C., Mulvaney, R., and Pritchard, H. D.: Full-Depth Englacial Vertical Ice Sheet Velocities Measured Using Phase-Sensitive Radar, J. Geophys. Res.-Earth Surf., 119, 2604–2618, 2014. a

Kleiner, T., Rückamp, M., Bondzio, J. H., and Humbert, A.: Enthalpy benchmark experiments for numerical ice sheet models, The Cryosphere, 9, 217–228,, 2015. a

Larour, E., Morlighem, M., Seroussi, H., Schiermeier, J., and Rignot, E.: Ice Flow Sensitivity to Geothermal Heat Flux of Pine Island Glacier, Antarctica, J. Geophys. Res.-Earth Surf., 117, F04023,, 2012a. a

Larour, E., Seroussi, H., Morlighem, M., and Rignot, E.: Continental Scale, High Order, High Spatial Resolution, Ice Sheet Modeling Using the Ice Sheet System Model (ISSM), J. Geophys. Res.-Earth Surf., 117, F01022,, 2012b (data available at:, last access: 20 February 2024). a, b, c, d

Le Brocq, A. M., Payne, A. J., and Vieli, A.: Antarctic dataset in NetCDF format, PANGAEA [data set],, 2010. a

Llubes, M., Lanseau, C., and Rémy, F.: Relations between Basal Condition, Subglacial Hydrological Networks and Geothermal Flux in Antarctica, Earth Planet. Sc. Lett., 241, 655–662,, 2006. a

MacAyeal, D. R.: A Tutorial on the Use of Control Methods in Ice-Sheet Modeling, J. Glaciol., 39, 91–98, 1993. a

Martos, Y. M.: Antarctic geothermal heat flux distribution and estimated Curie Depths, links to gridded files, PANGAEA [data set],, 2017. a

Martos, Y. M., Catalán, M., Jordan, T. A., Golynsky, A., Golynsky, D., Eagles, G., and Vaughan, D. G.: Heat Flux Distribution of Antarctica Unveiled, Geophys. Res. Lett., 44, 11417–11426,, 2017. a, b, c

Maule, C. F., Purucker, M. E., Olsen, N., and Mosegaard, K.: Heat Flux Anomalies in Antarctica Revealed by Satellite Magnetic Data, Science, 309, 464–467,, 2005. a, b, c

Mony, L., Roberts, J. L., and Halpin, J. A.: Inferring Geothermal Heat Flux from an Ice-Borehole Temperature Profile at Law Dome, East Antarctica, J. Glaciol., 66, 509–519,, 2020. a, b, c, d, e, f

Morgan, V. I., Wookey, C. W., Li, J., Van Ommen, T. D., Skinner, W., and Fitzpatrick, M. F.: Site Information and Initial Results from Deep Ice Drilling on Law Dome, Antarctica, J. Glaciol., 43, 3–10, 1997. a

Morlighem, M.: MEaSUREs BedMachine Antarctica, Version 1, Boulder, Colorado USA, NASA National Snow and Ice Data Center Distributed Active Archive Center [data set],, 2019. a

Morlighem, M., Rignot, E., Seroussi, H., Larour, E., Ben Dhia, H., and Aubry, D.: 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,, 2010. a, b, c, d

Morlighem, M., Williams, C. N., Rignot, E., An, L., Arndt, J. E., Bamber, J. L., Catania, G., Chauché, N., Dowdeswell, J. A., Dorschel, B., Fenty, I., Hogan, K., Howat, I., Hubbard, A., Jakobsson, M., Jordan, T. M., Kjeldsen, K. K., Millan, R., Mayer, L., Mouginot, J., Noël, B. P. Y., O'Cofaigh, C., Palmer, S., Rysgaard, S., Seroussi, H., Siegert, M. J., Slabon, P., Straneo, F., van den Broeke, M. R., Weinrebe, W., Wood, M., and Zinglersen, K. B.: BedMachine v3: Complete Bed Topography and Ocean Bathymetry Mapping of Greenland From Multibeam Echo Sounding Combined With Mass Conservation, Geophys. Res. Lett., 44, 11051–11061,, 2017. a

Morlighem, M., Rignot, E., Binder, T., Blankenship, D., Drews, R., Eagles, G., Eisen, O., Ferraccioli, F., Forsberg, R., Fretwell, P., Goel, V., Greenbaum, J. S., Gudmundsson, H., Guo, J., Helm, V., Hofstede, C., Howat, I., Humbert, A., Jokat, W., Karlsson, N. B., Lee, W. S., Matsuoka, K., Millan, R., Mouginot, J., Paden, J., Pattyn, F., Roberts, J., Rosier, S., Ruppel, A., Seroussi, H., Smith, E. C., Steinhage, D., Sun, B., van den Broeke, M. R., van Ommen, T. D., van Wessem, M., and Young, D. A.: Deep Glacial Troughs and Stabilizing Ridges Unveiled beneath the Margins of the Antarctic Ice Sheet, Nat. Geosci., 13, 132–137,, 2020. a, b, c

Motoyama, H.: The Second Deep Ice Coring Project at Dome Fuji, Antarctica, Sci. Dril., 5, 41–43,, 2007. a

Mouginot, J., Scheuchl, B., and Rignot, E.: Mapping of Ice Motion in Antarctica Using Synthetic-Aperture Radar Data, Remote Sens., 4, 2753–2767, 2012. a

Mulvaney, R., Rix, J., Polfrey, S., Grieman, M., Martìn, C., Nehrbass-Ahles, C., Rowell, I., Tuckwell, R., and Wolff, E.: Ice Drilling on Skytrain Ice Rise and Sherman Island, Antarctica, Ann. Glaciol., 62, 311–323,, 2021. a

Nagornov, O. V., Konovalov, Y. V., Zagorodnov, V. S., and Thompson, L. G.: Reconstruction of the Surface Temperature of Arctic Glaciers from the Data of Temperature Measurements in Wells, Journal of Engineering Physics and Thermophysics, 74, 253–265, 2001. a

Ng, F. and Conway, H.: Fast-Flow Signature in the Stagnated Kamb Ice Stream, West Antarctica, Geology, 32, 481,, 2004. a

Park, I.-W., Jin, E., Morlighem, M., and Lee, K.-K.: The basal temperature data for the entire Antarctic region derived from the thermal model of the ISSM, Korea Polar Research Institute [data set],, 2023. a

Parrenin, F., Dreyfus, G., Durand, G., Fujita, S., Gagliardini, O., Gillet, F., Jouzel, J., Kawamura, K., Lhomme, N., Masson-Delmotte, V., Ritz, C., Schwander, J., Shoji, H., Uemura, R., Watanabe, O., and Yoshida, N.: 1-D-ice flow modelling at EPICA Dome C and Dome Fuji, East Antarctica, Clim. Past, 3, 243–259,, 2007. a

Pattyn, F.: A New Three‐dimensional Higher‐order Thermomechanical Ice Sheet Model: Basic Sensitivity, Ice Stream Development, and Ice Flow across Subglacial Lakes, J. Geophys. Res.-Sol. Ea., 108, 2382,, 2003. a, b, c

Pattyn, F.: Antarctic Subglacial Conditions Inferred from a Hybrid Ice Sheet/Ice Stream Model, Earth Planet. Sc. Lett., 295, 451–461,, 2010. a, b, c

Pattyn, F., Perichon, L., Aschwanden, A., Breuer, B., de Smedt, B., Gagliardini, O., Gudmundsson, G. H., Hindmarsh, R. C. A., Hubbard, A., Johnson, J. V., Kleiner, T., Konovalov, Y., Martin, C., Payne, A. J., Pollard, D., Price, S., Rückamp, M., Saito, F., Souček, O., Sugiyama, S., and Zwinger, T.: Benchmark experiments for higher-order and full-Stokes ice sheet models (ISMIP–HOM), The Cryosphere, 2, 95–108,, 2008. a

Pettit, E. C., Waddington, E. D., Harrison, W. D., Thorsteinsson, T., Elsberg, D., Morack, J., and Zumberge, M. A.: The Crossover Stress, Anisotropy and the Ice Flow Law at Siple Dome, West Antarctica, J. Glaciol., 57, 39–52,, 2011. a

Pollard, D. and DeConto, R. M.: Description of a hybrid ice sheet-shelf model, and application to Antarctica, Geosci. Model Dev., 5, 1273–1295,, 2012. a

Price, S. F., Bindschadler, R. A., Hulbe, C. L., and Joughin, I. R.: Post-Stagnation Behavior in the Upstream Regions of Ice Stream C, West Antarctica, J. Glaciol., 47, 283–294,, 2001. a

Priscu, J. C., Kalin, J., Winans, J., Campbell, T., Siegfried, M. R., Skidmore, M., Dore, J. E., Leventer, A., Harwood, D. M., Duling, D., Zook, R., Burnett, J., Gibson, D., Krula, E., Mironov, A., McManis, J., Roberts, G., Rosenheim, B. E., Christner, B. C., Kasic, K., Fricker, H. A., Lyons, W. B., Barker, J., Bowling, M., Collins, B., Davis, C., Gagnon, A., Gardner, C., Gustafson, C., Kim, O.-S., Li, W., Michaud, A., Patterson, M. O., Tranter, M., Venturelli, R., Vick-Majors, T., Elsworth, C., and The SALSA Science Team: Scientific Access into Mercer Subglacial Lake: Scientific Objectives, Drilling Operations and Initial Observations, Ann. Glaciol., 62, 340–352,, 2021. a

Rignot, E., Mouginot, J., and Scheuchl, B.: MEaSUREs InSAR-Based Antarctica Ice Velocity Map, Version 2, Boulder, Colorado USA. NASA National Snow and Ice Data Center Distributed Active Archive Center [data set],, 2017. a, b, c, d, e, f

Rückamp, M., Humbert, A., Kleiner, T., Morlighem, M., and Seroussi, H.: Extended enthalpy formulations in the Ice-sheet and Sea-level System Model (ISSM) version 4.17: discontinuous conductivity and anisotropic streamline upwind Petrov–Galerkin (SUPG) method, Geosci. Model Dev., 13, 4491–4501,, 2020. a

Scambos, T. A.: Glacier Acceleration and Thinning after Ice Shelf Collapse in the Larsen B Embayment, Antarctica, Geophys. Res. Lett., 31, L18402,, 2004. a

Schlegel, N. J., Larour, E., Seroussi, H., Morlighem, M., and Box, J. E.: Decadal‐scale Sensitivity of Northeast Greenland Ice Flow to Errors in Surface Mass Balance Using ISSM, J. Geophys. Res.-Earth Surf., 118, 667–680, 2013. a

Schlegel, N.-J., Seroussi, H., Schodlok, M. P., Larour, E. Y., Boening, C., Limonadi, D., Watkins, M. M., Morlighem, M., and van den Broeke, M. R.: Exploration of Antarctic Ice Sheet 100-year contribution to sea level rise and associated model uncertainties using the ISSM framework, The Cryosphere, 12, 3511–3534,, 2018. a

Schoof, C.: The Effect of Cavitation on Glacier Sliding, P. Roy. Soc. A, 461, 609–627,, 2005. a

Schoof, C.: Marine Ice-Sheet Dynamics. Part 1. The Case of Rapid Sliding, J. Fluid Mech., 573, 27–55, 2007. a

Seroussi, H., Morlighem, M., Rignot, E., Larour, E., Aubry, D., Ben Dhia, H., and Kristensen, S. S.: Ice Flux Divergence Anomalies on 79north Glacier, Greenland, Geophys. Res. Lett., 38, 2011GL047338,, 2011. a

Seroussi, H., Morlighem, M., Rignot, E., Khazendar, A., Larour, E., and Mouginot, J.: Dependence of Century-Scale Projections of the Greenland Ice Sheet on Its Thermal Regime, J. Glaciol., 59, 1024–1034,, 2013. a, b, c

Seroussi, H., Ivins, E. R., Wiens, D. A., and Bondzio, J.: Influence of a West Antarctic Mantle Plume on Ice Sheet Basal Conditions, J. Geophys. Res.-Sol. Ea., 122, 7127–7155,, 2017. a

Seroussi, H., Nowicki, S., Payne, A. J., Goelzer, H., Lipscomb, W. H., Abe-Ouchi, A., Agosta, C., Albrecht, T., Asay-Davis, X., Barthel, A., Calov, R., Cullather, R., Dumas, C., Galton-Fenzi, B. K., Gladstone, R., Golledge, N. R., Gregory, J. M., Greve, R., Hattermann, T., Hoffman, M. J., Humbert, A., Huybrechts, P., Jourdain, N. C., Kleiner, T., Larour, E., Leguy, G. R., Lowry, D. P., Little, C. M., Morlighem, M., Pattyn, F., Pelle, T., Price, S. F., Quiquet, A., Reese, R., Schlegel, N.-J., Shepherd, A., Simon, E., Smith, R. S., Straneo, F., Sun, S., Trusel, L. D., Van Breedam, J., van de Wal, R. S. W., Winkelmann, R., Zhao, C., Zhang, T., and Zwinger, T.: ISMIP6 Antarctica: a multi-model ensemble of the Antarctic ice sheet evolution over the 21st century, The Cryosphere, 14, 3033–3070,, 2020. a

Shapiro, N. M. and Ritzwoller, M. H.: Inferring Surface Heat Flux Distributions Guided by a Global Seismic Model: Particular Application to Antarctica, Earth Planet. Sc. Lett., 223, 213–224, 2004. a, b, c

Slawny, K. R., Johnson, J. A., Mortensen, N. B., Gibson, C. J., Goetz, J. J., Shturmakov, A. J., Lebar, D. A., and Wendricks, A. W.: Production Drilling at WAIS Divide, Ann. Glaciol., 55, 147–155,, 2014. a, b, c, d

Smith, A. M., Anker, P. G. D., Nicholls, K. W., Makinson, K., Murray, T., Rios-Costas, S., Brisbourne, A. M., Hodgson, D. A., Schlegel, R., and Anandakrishnan, S.: Ice Stream Subglacial Access for Ice-Sheet History and Fast Ice Flow: The BEAMISH Project on Rutford Ice Stream, West Antarctica and Initial Results on Basal Conditions, Ann. Glaciol., 62, 203–211,, 2021. a

Smith-Johnsen, S., de Fleurian, B., Schlegel, N., Seroussi, H., and Nisancioglu, K.: Exceptionally high heat flux needed to sustain the Northeast Greenland Ice Stream, The Cryosphere, 14, 841–854,, 2020a. a

Smith-Johnsen, S., Schlegel, N., de Fleurian, B., and Nisancioglu, K. H.: Sensitivity of the Northeast Greenland Ice Stream to Geothermal Heat, J. Geophys. Res.-Earth Surf., 125, e2019JF005252,, 2020b. a, b

Talalay, P., Li, Y., Augustin, L., Clow, G. D., Hong, J., Lefebvre, E., Markov, A., Motoyama, H., and Ritz, C.: Geothermal heat flux from measured temperature profiles in deep ice boreholes in Antarctica, The Cryosphere, 14, 4021–4037,, 2020. a, b, c, d, e, f, g

Tsai, V. C., Stewart, A. L., and Thompson, A. F.: Marine Ice-Sheet Profiles and Stability under Coulomb Basal Conditions, J. Glaciol., 61, 205–215,, 2015. a

Tulaczyk, S., Kamb, W. B., and Engelhardt, H. F.: Basal Mechanics of Ice Stream B, West Antarctica: 2. Undrained Plastic Bed Model, J. Geophys. Res.-Sol. Ea., 105, 483–494,, 2000. a

van Ommen, T.: Ice Core Borehole Temperatures, Law Dome 1987, Ver. 2, Australian Antarctic Data Centre [data set],, 2023. a

Van Ommen, T. D., Morgan, V. I., Jacka, T. H., Woon, S., and Elcheikh, A.: Near-Surface Temperatures in the Dome Summit South (Law Dome, East Antarctica) Borehole, Ann. Glaciol., 29, 141–144, 1999. a, b

van Wessem, J. M., van de Berg, W. J., and van den Broeke, M. R.: Data Set: Monthly Averaged RACMO2.3p2 Variables (1979–2022); Antarctica, Zenodo [data set],, 2023. a

WAIS Divide Project Members: Onset of Deglacial Warming in West Antarctica Driven by Local Orbital Forcing, Nature, 500, 440–444,, 2013. a, b

Weertman, J.: Stability of the Junction of an Ice Sheet and an Ice Shelf, J. Glaciol., 13, 3–11,, 1974. a

Yang, J.-W., Han, Y., Orsi, A. J., Kim, S.-J., Han, H., Ryu, Y., Jang, Y., Moon, J., Choi, T., Hur, S. D., and Ahn, J.: Surface Temperature in Twentieth Century at the Styx Glacier, Northern Victoria Land, Antarctica, From Borehole Thermometry, Geophys. Res. Lett., 45, 9834–9842,, 2018. a, b, c, d, e

Zagorodnov, V., Nagornov, O., Scambos, T. A., Muto, A., Mosley-Thompson, E., Pettit, E. C., and Tyuflin, S.: Borehole temperatures reveal details of 20th century warming at Bruce Plateau, Antarctic Peninsula, The Cryosphere, 6, 675–686,, 2012. a, b, c, d, e, f, g

Short summary
This study conducted 3D thermodynamic ice sheet model experiments, and modeled temperatures were compared with 15 observed borehole temperature profiles. We found that using incompressibility of ice without sliding agrees well with observed temperature profiles in slow-flow regions, while incorporating sliding in fast-flow regions captures observed temperature profiles. Also, the choice of vertical velocity scheme has a greater impact on the shape of the modeled temperature profile.