Accepted Manuscript The Formation and Stability of Buried Polar CO2 Deposits on Mars Curtis V. Manning, Carver Bierson, Nathaniel E. Putzig, Christopher P. McKay PII: DOI: Reference: S0019-1035(18)30236-7 https://doi.org/10.1016/j.icarus.2018.07.021 YICAR 12970 To appear in: Icarus Received date: Revised date: Accepted date: 6 April 2018 22 June 2018 26 July 2018 Please cite this article as: Curtis V. Manning, Carver Bierson, Nathaniel E. Putzig, Christopher P. McKay, The Formation and Stability of Buried Polar CO2 Deposits on Mars, Icarus (2018), doi: https://doi.org/10.1016/j.icarus.2018.07.021 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. ACCEPTED MANUSCRIPT 1 HIGHLIGHTS ? The Clancy effect allows burial of CO2 deposits to preserve them during higher obliquity phases. ? Further stabilization of buried CO2 occurs at depth by close-off of pore space. AC CE PT ED M AN US CR IP T ? The deepest CO2 deposits approach the triple point temp. making basal melting possible. ACCEPTED MANUSCRIPT The Formation and Stability of Buried Polar CO2 Deposits on Mars AC CE PT ED M AN US CR IP T In the shelter of some of those scarps and troughs, deposits have been found that are thought to consist of solid CO2 covered by a layer of insulating water ice (Phillips et al., 2011; Bierson et al., 2016 ? hereafter, B2016). Curtis V. Manninga , Carver Biersonb , Nathaniel Recent radar observations of the southern polar cap with SHARAD, the Shallow Radar sounder, on the E. Putzigc and Christopher P. McKaya a NASA-Ames Research Center, MS 245-3, Moffett Field, CA Mars Reconnaissance Observer (MRO) have revealed 94035-1000, USA radar ?reflection-free zones? (RFZ) that require mab U.C. Santa Cruz, 1156 High St. Santa Cruz, CA, 95064 c terial with a dielectric constant close to that of solid Planetary Science Institute, 1546 Cole Blvd., Lakewood, CO CO2 in order for the depth analysis to render under80401 email@example.com lying water-ice layers flat in radargrams (Phillips et al., 2011), in continuity with abutting layers. Further Draft version August 18, 2018 investigations show the presence of sublimation pits in the unit known as Aa3 of the residual cap (Kolb and Tanaka, 2006; Tanaka et al., 2007) that overlies ABSTRACT reflection-free zones (RFZ), that pure Shallow Radar soundings of the south polar layered depositsthe (SPLD) have revealed buried CO2implying ice CO fills those RFZs (Phillips et al., 2011; B2016). with as many as three distinct layers separated by thinner layers of 2water ice (Bierson et al., 2016). This These deposits more than double the atmospheric layering is suggestive of formation by cyclical processes such as obliquity, eccentricity and advance of of COto on Mars (Putzig et 2 have perihelion. If obliquity is the main driver, then there were manyexchangeable opportunities inventory for CO2 deposits al., 2018), and are divided into layers, bounded by formed within the last 3.5-4 Myr. To persist, however, CO2 deposits must be followed by the deposition thinner layers of water ice. Up to three CO2 layof an insulating layer of porous water ice that can seal in the CO2 to survive the obliquity maxima ers are found in a given location. that follow. We suggest that the existing deposits were formed within the last 350 kyr. Each obliq- The maximum thickness of theindeposits is a little over one kilometer uity minimum was quickly followed by a period in which perihelion occurred the northern summer (B2016). The layering is suggestive (southern winter) solstice. Under these conditions, an enhanced deposition of low porosity, fine-grained of depositional cycles, likely by the oscillawater ice could occur on the winter pole. A similar series of depositions and dominated burials by water icequasi-periodic could tions of the obliquity (B2016). If so, then the Laskar have occurred between 2.75 and 2.2 million years ago, but these deposits were unlikely to have survived et al. (2004) obliquity calculations suggest that there the many high obliquity swings that occurred between 2.2 Ma and 400 ka. The formation and burial were upwards of twenty opportunities for the formaof CO2 ice could also have occurred in similar fashion earlier during the last ? 100 Myr in which the tion of deposits during flux the may last be ? 3.5 southern polar layered deposits were formed. Although at current times the geothermal too to 4 Myr of low average obliquity, as shown in Figure low to cause basal melting of CO2 deposits, in earlier times, deposits formed in deeper scarps and valleys 1. theBycyclic pattern in these deposits suggests could have led to basal CO2 melting and sequestration into the While regolith. extension, basal melting the deposits may be tied to the ? 125 and sequestration in the regolith may provide an explanation for the fate of the thick CO2 atmosphere kyr obliquity cycle, they may also be dependent on other factors, (pCO2 > ? 0.5 bar) implied by climate models of early Mars. such as eccentricity or the seasonal timing of perihelion. In particular, the survival of a CO2 deposit Keywords depends on their being sealed by a blanket of insuMars, polar caps, basal melting, obliquity cycles lating porous water ice before the subsequent rise in the obliquity causes the ambient surface temperature 1. INTRODUCTION to rise above the frost point of CO2 . Another reThe Martian south polar ice cap consists of 3 comquirement for deposit stability is that the water-ice ponents, the seasonal cap, the residual cap, and the deposit should be deep enough that the thermal grasouth polar layered deposits (SPLD). The seasonal dient would accelerate the sintering process near the cap is composed CO2 ice less than about 1 meter bottom of the deposit to the extent that the ice grains thick. The residual cap is composed of CO2 which compress to the point that the pores are cut off from survives the summer, and is about 10 meters thick. each other so that CO2 vapor cannot diffuse out. The SPLD is a more extensive deposit of water ice The temporal variations of obliquity, eccentricity with a small amount of dust, and some layers with a and perihelion advance, each with different cadences, high fraction of dust. It is approximately 3.5 km thick may result in CO2 and water-ice deposition cycles at its thickest (Tanaka, Kolb and Fortezzo, 2007). that, for a given low-obliquity excursion, may or may The layered deposits have scarps, gullies and troughs, not be synchronized in a way to form and seal a formed during periods of strong ablation (Milkovich CO2 ice deposit. The availability of the Laskar et and Plaut, 2008), although it remains possible, by al. (2004) calculations of obliquity and orbital paanalogy to the northern cap, that the chasmata and rameters of Mars for the last 20 Myr allows a detroughs may be constructional features formed during tailed assessment of opportunities for the formation accretion (Holt et al., 2010; Smith and Holt, 2010)). ACCEPTED MANUSCRIPT 2 2. SHARAD DATA AND CO2 DEPOSITS 3. BURIAL AND PRESERVATION OF CO2 DEPOSITS We have noted that the existence of today?s buried CO2 deposits in the SPLD depend on the subsequent burial by an insulating boundary layer while the surface CO2 cap is still stable to ambient temperatures. The conditions for transport of water southward are sensitive to the orbital parameters. According to Clancy et al. (1996), water ice could accumulate in the south when perihelion occurs in the northern summer at rates approaching a solid ice-equivalent of up to 10 mm per year. If we assume an average rate half of that, then during a ? 10 kyr period during which perihelion occurs in the northern summer (southern winter), an equivalent of 50 m of zeroporosity ice could accumulate. However, some of this layer would be lost during the subsequent southern summer, during which ablation would dominate in the south. Given a 30 m deposit (the average of BL1 or BL2 ) of pure ice after an obliquity cycle, then before compaction the high porosity snow could have initially been ? 100 m thick if it had a porosity of 70%. However, many processes could have driven this porosity down, including sintering and compaction due to the overburden of the seasonal cap. The current top boundary layer, BL3 , is less than 20 m thick, some of which is mostly buried by unit Aa4 (Tanaka et al., 2007). However, it is very doubtful that BL3 is as compacted as the boundary layers below it, and so is likely to have a lower water content. CE PT ED M AN US In B2016, the Shallow Radar instrument on the Mars Reconnaissance Observer (MRO) was used to collect a total of 429 radar scans of the SPLD. These were used to analyze the layered deposits at latitudes north of 87? S (see B2016 Fig. 1a). The RFZs are detected within the residual ice cap, and taken to be CO2 . In Table 1, we reproduce data from Table 1 of B2016. In the first column are shown the layer designations (naming conventions of units set by Tanaka et al., 2007), in column two are the compositions. These units are divided vertically into three sections separated by horizontal lines. The top section is the residual cap. The middle section is what could be called the upper SPLD composed of the alternating layers of RFZs and boundary layers, which are numbered from the bottom up. Below the bottom line is the lower SPLD with two units of water ice that extend down to the regolith. In column 3 are posted the averaged thicknesses in meters of each layer. We note that the round numbers presented in Table 1 mask irregularities in the thickness of the RFZs and boundary layers, and not all layers are present at each location. As a result, adding the layer thicknesses from the top to the bottom of layer Aa3a , one finds a thickness of 1.2 km, although the deepest part of the CO2 deposit is only a little over 1 km in depth (B2016). current cap (Nye et al., 2000). On current Mars, CO2 deposits need to be contained, or else they will flow on geologically short time scales. Radargrams presented in B2016 show that the deposits are exposed in scarps and troughs. Radargrams also show the details of the layering of the water-ice cap that contain the deposits. Below some of the CO2 deposits, horizontal layering in the water ice can be seen. Elsewhere, there are apparent signs of differential settling and locations of unconformity in the water-ice layering, often aligning with the scarp slopes that contain the deposits, extending into lower water-ice units. CR IP T and burial of CO2 deposits on Mars. Below, we present the B2016 data on RFZs (�, followed, in � by a discussion of the factors affecting the deposition and preservation of buried CO2 deposits in recent, and the less-recent past. In � we touch on the possibility of basal melting of CO2 ice. Basal melting is strongly limited by the geothermal flux, the depth of the deposit and, to a lesser extent, by the conditions within the top boundary layer (Mellon, 1996). We discuss our results and present conclusions in � AC Table 1: Layering data (B2016) Unit name Composition Unit thickness (m) Aa4 CO2 /H2 O < ? 10 BL3 H2 O < 20 Aa3c CO2 300 BL2 H2 O 40 Aa3b CO2 < 700 BL1 H2 O 20 Aa3a CO2 200 Aa2 H2 O 300 Aa1 H2 O 3500 Because of the unusually low viscosity of CO2 , a cap of pure CO2 could not support the profile of the 3.1. The Clancy effect(1996) Currently, Mars is near aphelion in the northern late spring to summer season, Lp = 251? . Under these conditions, with perihelion near the southern summer solstice, the transport of water vapor is generally northward (Clancy et al., 1996). However, about 25 kyr ago, the longitude of perihelion was in the northern summer (Lp = 71? ), and when this occurs, conditions are optimal for southward transport during the southern winter. Under such conditions, especially when the eccentricity is large, temperatures rise, water sublimes from the northern cap, and moves toward the equator. When perihelion occurs in ACCEPTED MANUSCRIPT 3 CR IP T high obliquity, a switch to low obliquity would likely provide a substantial source for both poles, regardless of eccentricity. Breaks in a continuity of episodes of delivery can also occur during periods of high polar insolation (e.g., Putzig et al., 2009). The Clancy effect implies that the accumulation of water-ice layers on the north and south poplar layered deposits should not have a direct correspondence with each other since the timing of periods in which net accumulation could be expected are separated by ? 25 kyr, by which time eccentricity and obliquity will have changed significantly. However, the major depositional breaks between sets of deposits (e.g., Milkovich and Plaut, 2008) very likely correspond between north and south, being sensitive both to high obliquity and or low eccentricity. If a correspondence between these ablative breaks that define the sets of layers is made, then the variation of layers within the respective sets (both north and south) can be modeled if these sets are properly correlated with the Laskar et al., (2004) obliquity and orbital data. Then the details of the corresponding sets could provide insight into the details of accumulation and ablation on the respective polar layered deposits. According to Levrard et al., (2007), ice migrates to the poles either from equatorial regions or the higher latitudes of either hemisphere. It is likely that at 5 Ma, before which the obliquity was high, the equatorial regions would have had a significant amount of water while the higher latitudes would have had little. After millions of years of low obliquity, however, the low-latitude source of water ice has probably been depleted. Generally, therefore, the rate of growth of the polar caps during a long-term period of low average obliquity should decline in time, with seasonal exchanges driven primarily by the Clancy effect. Judging from the estimated cratering age of the SPLD, there must have been a significant period of low obliquity, ending between about 30 and 100 or so Myr (Koutnik, Byrne and Murray, 2002) in which much of the equatorial ice was transported to the poles. Following that there must have been a period of high obliquity to dessicate the poles, leaving a thick lag at Promethei Lingula, followed, during the last 5 Myr, by a new low-obliquity stage. It would appear that the 300 m of water ice accumulated in the lower SPLD was added during the last 5 Myr since it shows no sign of having had an extended period of ablation. While we can know little of the obliquity and orbital parameters beyond 10 or 20 Myr, we have detailed information for the last 10 Myr. This can be used to diagnose the forces which formed and buried the CO2 deposits in the upper section of the SPLD. AC CE PT ED M AN US the northern summer, atmospheric water vapor concentrations can be an order of magnitude greater than at aphelion (Clancy, 1996). Since the temperature in the northern hemisphere is from 20 to 25 K greater under these conditions than at aphelion, the altitude at which the saturation of water vapor would occur is 25 or 30 km. This vapor can then be caught up in the ascending branch of the Hadley circulation and taken southward. When perihelion occurs in the southern summer, conditions are reversed; temperatures in the northern summer are lower, and the saturation of water vapor occurs at altitudes of about 10 km, below the reach of the Hadley circulation. The water vapor nucleates on dust, freezes and gravitationally settles. Thus, perihelion advance can very broadly control the growth and ablation of the respective ice caps. Simulations of the south polar ice cap by Montmessin et al. (2007) confirm this view for both phases of the perihelion cycle. An exception to this pattern may occur when global dust storms cause the atmosphere to warm during the northern aphelion, such as in 1976 and 1978 (during the Viking mission), to the point that higher concentrations of water occurred, allowing the vapor to rise to high altitudes and become entrained in the Hadley circulation (Clancy et al., 1996). The precession-induced variability in water vapor transport is called the ?Clancy effect? (e.g., Montmessin et al., 2007). Because the perihelion advance is a 51,000 year cycle, it reverses direction every 25,500 yr. High eccentricities enhance the net yearly transport of water since it enhances the north/south temperature asymmetry, helping to magnify the Clancy effect. According to the modeling work of Levrard et al. (2007; Fig. 2) higher obliquities may stifle net yearly transport since summer ablation at the aphelion would increase more rapidly with obliquity than the accumulation rate at perihelion. Under higher obliquity water tends to move toward the lower latitudes. However, for moderately high obliquity, this effect may be may be effectively counter balanced if the eccentricity is high since, if aphelion is in the southern summer, a large eccentricity can reduce polar insolation to a level normally delivered only at a lower obliquity. While the short-term variability of the eccentricity has a 25,500 kyr period, there is a longer 1.8 Myr period over which the eccentricity moves from a high of ' 1.2 to zero and back again (Ward, 1979). When the eccentricity is near zero, there is a cessation of accumulation at the poles (Clancy et al., 1996), producing a break in a series of layers, as seen in the NPLD. However, this conclusion must depend on there being only a small rate of delivery of water from the low latitudes, as is the case today, but after a long period of 3.2. The recent past ACCEPTED MANUSCRIPT 4 3.2.1. The last 520 kyr AC CE PT ED M AN US To give a picture of the conditions during which the CO2 deposits on the SPLD were formed, we present, in Figure 2, the recent historical trends in orbital parameters and obliquity (Laskar et al., 2004) in correspondence with our simulated evolution of CO2 inventories over the last 520 kyr. To chart the evolution of CO2 inventories we use the Mars Evolution Code (MEC; Manning et al., 2006). MEC is a yearly averaged energy balance model of CO2 evolution to which the many sources and sinks are fully coupled, and driven by the modeled obliquities of Laskar et al. (2004). The upper three panels of Fig. 2 show the trends in eccentricity, longitude of perihelion, and obliquity, respectively. The horizontal line in the obliquity plot shows the approximate obliquity (? ' 26? ), at which a surface CO2 ice cap would sublime away. The bottom panel of Fig. 2 shows simulation results, initialized so that the current atmospheric pressure (orange ?+? sign) can be reached. The running atmospheric partial pressure of CO2 , PCO2 , is in black, the regolith content is in red, and the surface CO2 ice deposit is shown in blue. The calculation of polar temperatures in MEC are based on a model of Marinova et al. (2005), developed to account for the meridional advection of heat from equator to pole. This is in turn based on an analytical approximation of an advection parameter, ?, by Gierasch and Toon (1973). Meridional adevection of heat plays an important role in determining the obliquity at which the atmosphere may collapse. However, recent general climate model studies by Soto et al. (2015) claim that models of meridional transport, such as Gierasch and Toon, (2015), Marinova et al (2005), and Manning et al. (2006; see Eq. A32), may overestimate the magnitude of heat transport, causing an underestimation of the obliquity below which atmospheric collapse may be triggered. Polar temperatures are calculated using yearly averaged solar insolation, greenhouse warming and meridional advection. The latter is based on the equator to pole thermal gradient (Haberle et al., 1994). The calibration procedure for MEC is to adopt an average polar albedo, then adjust the magnitude of the advection parameter to achieve an atmospheric pressure of ? 6.5 mbars under the current obliquity, as well as to arrange for the disappearance of the cap when the obliquity rises to approximately 27.5? (Francois et al., 1990). We have found that reducing the advection parameter to one tenth of that used in our study can be compensated by a relatively small reduction of the albedo (i.e., Ap ' 0.700 reduced to Ap = 0.686. Results show that for atmospheric pressures visited by the simulation in Fig. 2, the obliquities at which the atmosphere both collapses and re-sublimes would be increased by less than one degree of arc for an advection parameter reduced by a factor of ten. This difference has an insignificant affect our results, primarily because the response of calibration to a reduction in the advection parameter is to lower the albedo, which tends to stabilize pole temperatures. Also, for smaller advection parameter values, greenhouse warming dominates the advective heat transport, so that further reduction of the advection parameter produces smaller changes. However, for significantly higher atmospheric pressures (e.g., > ? 100 mbars), the reduction of advective transport can have a significant effect on both collapse and resublimation obliquities. Note that the Manning et al. (2006) model includes an active regolith which releases CO2 when the atmosphere collapses (Zent & Quinn, 1995; Toon et al., 1980). This serves to make CO2 caps more massive than would be the case when ignoring the regolith content since the decline in atmospheric pressure causes the regolith to outgas. Our parameterization of the regolith is presented in Manning et al. (2006), Fig. 2 of Appendix 6, where the current regolith content was assumed to be 50 mbars. However, according to Zent and Quinn (1995), without a very deep powdered regolith (e.g., 100 m), the most likely exchangeable inventory is between about 1 and 4 kPa. We scaled the regolith content to 60% of our previous value so that the adsorbed inventory would be ? 30 mbars at the current averaged temperature and pressure of Mars (see Eq. A21 of Manning et al., 2006). With an active regolith, if a 12 mbar atmosphere were to collapse to PCO2 = 1.0 mbar, we find that the CO2 ice cap thereby formed is equivalent to about 26 mbar atmospheric. This only includes 11 mbars from the atmosphere, so the other 15 mbars must come from the regolith. This result implies that the regolith strongly buffers the atmospheric pressure. Thus, with an active regolith, if the current inventory of frozen CO2 within the RFZs (? 7 mbars; Putzig et al., 2018) were to sublime, the atmosphere (? 6.5 mbar) would pick up only about 2.9 mbar of CO2 making it 9.4 mbars while the regolith would pick up about 4.1 mbars. The vertical lines in Fig. 2 are aligned with particular longitudes of perihelion; red (dot-dashed line) for LP = 90? , and blue (short-dashed lines) for the flanking values, LP = 0? and 180? . These are intended to help the eye in judging the synchronization of the southward transport of water with the CR IP T The obliquity and orbital parameter simulations of Laskar et al. (2002; 2004) provide an opportunity to study the effects of obliquity, perihelion advance and eccentricity on the residual cap. ACCEPTED MANUSCRIPT 5 CR IP T B2016 as the top and bottom of the water-ice layer. However, the analysis above suggests an alternative interpretation. We propose these two reflections may be two distinct water-ice layers with thicknesses below the ? 15 meter resolution of SHARAD separated by a thin CO2 layer of about 10 meters. If true, it would help remove any ambiguity associated with the suggested B2016 obliquity associations of the respective CO2 deposits. Looking a little further back in time, we note that in Fig. 2, the obliquity minimum at ? 450 ka has no corresponding deposit in the SPLD under our suggested interpretation, although perihelion occurs in the northern summer within the time that the CO2 ice cap is still massive. We find two possible reasons why this unit might be missing. First, the eccentricity is relatively low and declining while perihelion is in the southern winter, suggesting that the Clancy effect could have been weak, producing a relatively thin insulating layer of snow. Another factor could be that, following the burial of the deposit, the obliquity reaches about 30? . A thin water-ice cap and a relatively strong period of ablation could therefore explain its absence. Simulations with the LMD Martian Global Climate Model (Forget et al., 1999; Montmessin et al., 2004) by Levrard et al., (2007) of the north polar layered deposits (NPLD; Figure 1), suggest that at this obliquity, ? = 30? , the ablation of a couple of millimeters per Mars year is likely. In principle this result is applicable to the south over a perihelion cycle, since the two caps appear to be in a near depositional equilibrium. We note that the ablation of two millimeters per Mars year, summed over 10 kyr of ablation is still the equivalent of about ten meters of zero porosity ice lost over a perihelion cycle. Considering a possibly thin water-ice layer, this may be significant. It is similarly notable that the deposit at 320 ka, here assumed to be Aa3a , also has the perihelion rather optimally placed relative to the obliquity minimum. Nevertheless, this deposit has a rather spotty presence on the SPLD. As can be seen in B2016 Fig. 1b, RFZs are generally localized in three separate locations, low on the cap, high on the cap, and in between. The deposit of greatest area, layer Aa3a is low on the cap, but it is relatively thin. The CO2 deposits half way up the cap have no representative of unit Aa3a (just Aa3b and Aa3c ). There is only a small unit Aa3a deposit high on the cap. We note, however, that the high-obliquity period that followed (240 ? t ? 300 ka) barely exceeds the 26? limit of surface CO2 stability, suggesting that ablation should not have been a serious threat to its stability. On the other hand, the eccentricity during southern winter solstice is relatively low and declining during the pe- AC CE PT ED M AN US formation of the buried CO2 caps. Horizontal bars are placed between these vertical lines at an ice cap pressure-equivalent at the southern winter solstice for each potential CO2 ice burial. In the calculation of the evolution of CO2 shown in Fig. 2, it was assumed that 10% of frozen CO2 was buried each time the CO2 cap forms. In defense of this, we note that existing buried CO2 is only found in sheltered regions of the SPLD that amount to about 10% of its surface area, suggesting that the more exposed CO2 sublimes when obliquities rise. This is a convenient way to simulate the progressive burial of some of the CO2 cap. In actuality, burial may depend on the seasonal dust storms, trends in eccentricity with time, and the availability of deep sheltered spaces for CO2 deposits, as well as the variable longitude of perihelion, eccentricity and obliquity (Clancy et al., 1996). The CO2 ice deposits represented in Fig. 2 are surface caps since MEC does not currently model their burial. In B2016 it was suggested that the CO2 deposits in the SPLD may have formed during the last 340 kyr; unit Aa3a at 340 ka, unit Aa3b at 200 ka and unit Aa3c at 22 ka. Figure 2 shows that in each case, the minima in obliquity are followed by a period in which perihelion occurs in the northern summer, as we are led to believe should be the case for CO2 burial (Clancy et al., 1996). To these three deposits, we now suggest another; in Figure 3 we show an enlargement of the radargram from the upper-right corner of B2016 Fig. 2a, showing that the boundary layer, BL2 , shows two separate reflections, with a RFZ in between. If this deposit (let us provisionally call it Aa3b0 ) is indeed a CO2 deposit, then there must be an obliquity minimum following the 200 ka unit Aa3b , and perihelion in northern summer must quickly follow that. Reference to Fig. 2 shows that all conditions are met, as there is a small obliquity minimum at ? 125 ka, followed shortly afterwards by perihelion in the northern summer. The possible bifurcation of BL2 can be seen at other locations of B2016 Fig. 2a and 2b that show BL2 , but this double layer is not apparent in other boundary layers. In their Fig. 2b (upper right), this same locale (see the paths of the two radar scans in Fig. 1b of B2016) has an irregular BL2 with the two reflections erratically distributed, as though the depth of the boundary layer varied within the footprint of the radargram. We must note, however, that the eccentricity is low during the ? 10 kyr period at ? 125 ka in which BL2b would be expected to have formed, suggesting a possibly thin water-ice layer. However, the obliquity that follows this period barely exceeds 26? , thus minimizing ablation from the following summers. The two reflections of BL2 were interpreted by ACCEPTED MANUSCRIPT 6 CR IP T of 6.5 mbar (dotted) and 12 mbar (solid line), and the figure shows the effect both on the atmosphere and ? regolith at obliquities ? > ? 26 . The situation is in fact quite comparable to the last 350 kyr, including the presence of what could plausibly produce a double water-ice layer with CO2 in between. In the comparison, it would also be necessary to explain the absence of the first possible deposit at ? 2.7 ka. However, the high obliquities near t = 1 and 2 Ma augur for high rates of ablation of the upper boundary layer and the sublimation of any near-surface CO2 deposits. We suggest these deposits could not have survived the much more intensive insolation that was to follow. 3.2.4. Opportunities for basal melting in the last ? 100 Myr Although the above analysis only pertains to the last ? 3.5 Myr, the southern polar layered deposits are thought to be many millions of years old, judging from crater retention ages (Koutnik, Byrne & Murray, 2002). It is reasonable to consider whether, during the extended period of their construction, there could have been opportunities for similar atmospheric collapses, with CO2 deposits in crevasses and scarps, buried under snow, similar to the deposits in the SPLD. Indeed, shallow radar studies have found numerous RFZs within the lower section of the SPLD (Phillips et al., 2011, supplementary material). These RFZs are not confirmed as CO2 since we do not yet have a good means to establish their dielectric constant since no favorable layering geometries have been found. The stratigraphic study of Milkovich and Plaut (2008) finds that the southern polar layered deposits are composed of three major sequences or units. The second unit, P LL, which is exposed in Promethei Lingula, consists of a sequence of sets of layers that are interrupted by unconformities suggestive of periods of ablative climate change. This sequence ends with a strong period of ablation that produced deep scarps and the curvilinear valleys (Milkovich and Plaut 2008). PLL has a domed shape of height ? 3 km above the southern plains (Tanaka et al. 2007; Milkovich and Plaut, 2008). According to Koutnik et al. (2002), the crater-retention age of this unit is between 30 and 100 Myr of age. We suggest that during the construction of P LL and of the following ? 300 meters of layered deposits upon which the lower SPLD is placed, there may have been many opportunities for the deep burial of CO2 deposits some of which could have resulted in the basal melting of CO2 . ED M AN US riod in which the Clancy effect would otherwise have been active, suggesting weak water transport. We suggest that other, perhaps mesoscale factors, may have affected the distribution of CO2 deposits, or subsequent water-ice deposits. We note also that the top boundary layer, BL3 , is thin (< 20 m), especially considering that it is unlikely to have been compacted as much as the lower boundary layers. The buried boundary layers are likely to have been sintered at depth (a temperaturesensitive process), and crushed by the overburden of hundreds of meters of CO2 and water ice. The sublimation pits in BL3 (B2016) may be a sign of the permeability of the layer, suggesting that not all of unit Aa3c will long endure. While the deposits Aa3b (200 ka), Aa3c (22 ka), and the proposed deposit at Aa3b0 (125 ka), seem to obey the general criteria for CO2 deposits outlined earlier, the obliquity fluctuation at 450 ka, and to a lesser extent, 340 ka, show a lack of proportionality in their resulting deposits with respect to obliquity minimum and perihelion placement. The variability of the Aa3a deposit size with position on the SPLD suggests that the value and trend of eccentricity may be a significant factor in determining the thickness of the boundary layer. To help disentangle the relative strengths of obliquity and eccentricity and perihelion in determining water transport will take a more detailed model of atmospheric evolution such as a zonal/seasonal model or a GCM. 3.2.2. Obliquity variations between 520 ka and 1 Ma AC CE PT For obliquity variations between 520 ka and 1.0 Ma, there are few situations in which the proper relative positioning of perihelion and obliquity minima occur. One case at ? 800 ka finds the relative perihelion and obliquity minimum ideally placed for a strong Clancy effect, but evidence for that deposit is lacking in the SHARAD data. We note that the lowobliquity excursion at ? 800 ka was followed by a high-obliquity excursion in excess of ? 35? (see Fig. 1). According to Figure 1 of Levrard et al. (2007), in the northern summer season, ablation can far exceed the winter accumulation rate at such high obliquity, suggesting that a hypothetical deposit of water ice would not survive a subsequent obliquity maximum of ? ? 35? , allowing the just-buried CO2 deposit to sublime. Again, this model of the NPLD ought to, under averaged conditions, accurately approximate the response on the south pole. 3.2.3. Opportunities at 2.75 ? t ? 2.2 Ma In Figure 4, we explore another period of time in which the variations are minimal ? those between 2.2 Ma and 2.75 Ma. We ran the MEC at initial pressures 4. THE THERMAL STRUCTURE OF CO2 SLABS ACCEPTED MANUSCRIPT 7 ?T = H ?z, k (1) 4.1. The effective thermal conductivity of the upper boundary layer Although we are neglecting the water ice of the inner boundary layers, we cannot neglect the thermal effect of the insulating upper boundary layer. Mars conditions are ideal for forming a very low thermal conductivity upper boundary layer; moisture delivered in the Hadley circulation could condense around sub-micron sized dust as nucleation sites, freeze, and fall to the residual cap, forming dry, particulate snow. On the other hand, it is theoretically possible that saturated air advected into the area close to the residual cap could be directly deposited on the residual cap, forming something akin to hoarfrost, which could have a significantly greater thermal conductivity. Under the conditions of a relatively rapid deposit delivered by Hadley circulation, we assume that a dry snow is the more likely initial result. We assume the boundary layers are composed of frozen water that was condensed around dust grains as nucleation sites. We assume the snow particles have a diameter of approximately 50 祄. Terrestrial observations of the thermal conductivity of dry snow suggest a lower bound of about k = 0.078 W m?1 K?1 (e.g., Mellon, Ferguson and Putzig, 2008). However, when the mean free path of molecules is greater than the pore sizes, the thermal conductivity is controlled by the atmospheric pressure. For dust-nucleated snow with a size of ? 50 祄, the thermal conductivity is reduced by over an order of magnitude at current Martian pressures (Mellon, McKay, and Grant, 2015; Edgett & Christensen, 1991), suggesting ?1 ?1 that kH2 O < ? 0.0078 W m K . However, the sintering process, principally produced by grain boundary diffusion and sublimation-condensation, can change the effective grain size at temperature-sensitive time scales, causing an increased surface area of contact between grains, and a higher thermal conductivity. The process of sintering can fuse individual grains together, increasing the thermal conductivity of the material. In an experimental study by Kaempfer and Schneebeli (2007), they recorded the advance of the sintering process at a range of temperatures over periods of up to a year. They used a parameter called the specific surface area (SSA), which is the surface H = 20 (� mW m?2 , M AN US where ?z is the slab thickness in meters, and H is the geothermal flux (Mellon, 1996). The canonical value of H is 30 mW m?2 (Clifford, 1987; Heldmann et al., 2005), however recent studies with SHARAD suggest that there has been very little isostatic compensation (less than 100 m) in the Gemina Lingula lobe of the NPLD (Phillips et al,., 2008), or the main lobe of the NPLD (Putzig et al., 2009), implying an elastic thickness greater than about 300 km (Phillips et al., 2008; supplementary material). This appears to limit the value of the geothermal flux at the north pole to H ? 13 mW m?2 (Dehant et al., 2012). However, stagnant lid models suggest that there is significant variability in the thermal flux over the planet, some of which suggest a south polar value near H ? 20 mW m?2 (Dehant et al., 2012). Thus, we consider that the south polar geothermal flux is in the range, perspective, but this is not particularly important for understanding the issues of basal melting. Because boundary layers, BL1 and BL2 , are relatively thin and their thermal conductivity is likely to be high, the temperature difference due to their presence would add to only a little over 1 K at the bottom of the CO2 slab. In the following analysis, we simplify the problem by considering a pure CO2 slab of arbitrary thickness. CR IP T Here we examine the possibility of basal melting, and under some conditions, sequestration into the regolith. In equilibrium, the temperature shift across a layer of thermal conductivity, k, is (2) 93.4 W m?1 K?1 , T CE kCO2 (T ) = PT ED although the global average appears close to 30 mW m?2 . Please note that, for the convenience of presentation, we have given H in units mW m?2 , but it must be converted into W m?2 our equations. Below, we study the thermal structure of CO2 slabs. The thermal conductivity of CO2 is a function of temperature, well-modeled by the following function, (3) AC where T is the temperature in Kelvins (Kravchenko and Krupskii 1986). For instance, a moderately thick upper boundary layer could result in a temperature T = 165 K at the top of the CO2 slab, for which the thermal conductivity would be k = 0.566 W m?1 K?1 , while at the triple point of CO2 , T ' 217 K, k = 0.43 W m?1 K?1 . Substituting Eq. 3 into the differential form of Eq. 1, one finds the following integrable expression, dT H = dz, (4) T 93.4 For the melting of a thick CO2 slab to occur, it must reach a temperature of 216.67 K at its base. We round this up to 217 K for our calculations. The fact that the CO2 ice is interspersed with layers of water ice is interesting from the stratigraphic ACCEPTED MANUSCRIPT 8 CR IP T and the thermal gradient within the top boundary layer would accelerate the diffusion of CO2 from the deposit. Such conditions could lead to a permanent feature being formed, such as the sublimation pits seen overlying the RFZs today (Phillips et al., 2011). Finally we suggest that seasonal frost could affect the surface of the top boundary layer by penetrating and binding the snow particles together, increasing conductivity at the surface if the boundary layer is not covered with seasonal/permanent CO2 ice. However, we believe this is unlikely to penetrate to a significant portion of the top 5 meters, so we neglect this effect. If we assume that both BL1 and BL2 have a porosity near zero (note that at only 90 meters under a CO2 cap on Mars, the pressure would exceed 5 bars), then we may get an estimate of their thickness when they were first formed. The lowest boundary layer, BL1 , which is observed to be ? 20 m (Table 1), would be 67 meters if it originally had a porosity of 70%. However, as we have seen, this maximal porosity may not last long. If we accept that BL2 is actually two layers of water ice, each 15 meters (and CO2 ice of thickness 10 m), then this would convert to 45 meters each. However, these differences in thickness would not significantly affect the temperture jump since, as we suggest above, below about 10 m of the top waterice layer, the conductivity would increase to the point that it would not significantly impede thermal conduction. The southern cap is said to currently have a permanent cap of CO2 , less than about 10 m thick overlying the SPLD, which can survive the summer due to its high albedo. In this case, the surface temperature would remain near the current frost point of CO2 , about 150 K. Under a higher obliquity, this CO2 surface cap would most likely disappear during the summer season, and seasonal changes of surface temperature would accelerate the sintering process. The seasonal skin depth would vary from about 0.54 m (for k = 0.01 W m?1 K?1 and ? = 300 kgm?3 ) and 1.1 m (for k = 0.1 W m?1 K?1 and ? = 700 kgm?3 ). But during the last 400 kyr, obliquities were low, and a permanent CO2 cap could have persisted. We summarize these results as follows. For the current Mars, we judge that the top boundary layer can be decomposed into three layers. Starting from the top, a 5 meter layer of snow, then a 5-meter layer of moderately sintered snow, and below that, a high-conductivity layer, which contributes little to the temperature jump between the top and bottom of the top boundary layer. We assume a geothermal flux, H = 20 mWm?2 (Eq. 2). For the top layer,with k = 0.008 W m?1 K?1 , the thermal gradient is dT /dz ? H/k = 2.5 Km?1 , or 12.5 K AC CE PT ED M AN US area of snow accessible to gases per unit volume. Experiments measured the deformation and densification of grains of snow as small as 50祄 for a range of temperatures. We found that measured rates were proportional to the vapor pressure over ice. At 219 K they report that the decline of the SSA is about 19% in a year. By calculating the partial pressure over water ice at temperatures down to 150 K, we were able to approximate the time scale of sintering as a function of temperature. We found that at 150 K, it would take about 4 � 105 years for the SSA to decline by 20%, roughly three obliquity cycles, and hence a negligible rate for a south polar surface boundary layer. On the other hand, at Martian atmospheric pressures, newly formed snow at this temperature would have a very ?1 ?1 low thermal conductivity (kH2 O < ? 0.008 W m K ). In this case, there is a substantial thermal gradient dT /dz ' 2.5 K m?1 at equilibrium for a geothermal flux of H = 20 mW m?2 . Under these conditions, if the surface temperature were 150 K, at a depth of five meters, the temperature would rise to ? 163 K at current Martian atmospheric pressures. At 163 K, the rate of the sintering process is about 1.5 orders of magnitude greater than that at 150 K, so that the SSA would decline by 20% in only about 15 kyr, roughly the length of the period in which perihelion aligns with northern summer. Thus, with the Clancy effect, the great majority of the upper boundary layer would be deposited within this sintering time scale. Below 5 meters, however, the thermal conductivity of the top boundary layer could be significantly affected. We suggest that between 5 and 10 meters depth, the thermal conductivity could be increased to k ' 0.02 W m?1 K ?1 , only moderately higher than at 5 meters depth since our initial assumption is that the boundary layer for unit Aa3c formed about 22 ka, not much longer than the sintering time-scale of 15 kyr at that depth. Below 10 meters depth, we expect that the effects of the overburden pressure on the density and the thermal conductivity, would become relevant, accelerating grain boundary diffusion, sublimationcondensation, and densification. We suggest that ?1 ?1 the thermal conductivity k > ? 1.0 W m K , and that at some point, the sintering process will cause the pores to close, preventing the diffusive venting of CO2 . This pore blockage may be a crucial aspect of the stability of buried CO2 deposits which, on Earth, occurs between snow densities of 820 and 840 kgm?3 (Herron and Langway, 1980). Sintering would have to proceed to this point in a time-scale less than about 10 kyr since obliquity is generally rising after the CO2 deposit is formed, and rising temperatures ACCEPTED MANUSCRIPT 9 CR IP T would just reach the triple point temperature is, z ? = 1265 meters for H = 20 mW m?2 . The extremal values for H = 20 (� mW m?2 , and from Eq. 5, yield are z ? = 936 m and z ? = 1816 m for the upper and lower bounds of H, respectively. We note here that, on the basis of this, basal melting is a reasonable possibility during the previous ? 5 My of growth of the residual cap. Note, however, that the thermal gradient in the CO2 slab is, dz/dT = k/H, CO2 ?1 ?1 where we let kef f = 0.5 Wm K , (the value of k at 187 K). We find dz/dT = 25 mK?1 for the central value, 33 mK?1 for H = 15 mWm?2 , and 20 mK?1 for H = 25 mWm?2 . Thus, if our estimate of the temperature at the bottom of the upper boundary layer was high by, say, 5 K, then z ? would be greater by approximately five times the above values ? 125 meters for H = 20mWm?2 . 4.2.2. Deposits with z 0 > z ? AN US over 5 meters. For the zone (k = 0.02 W m?1 K?1 ), 5 ? z ? 10 m, dT /dz ' 1.0 Km?1 , or 5 K over the next 5 meters. For the remainder of the top boundary layer (k = 1.0 W m?1 K?1 ), we estimate dT /dz = 0.02 Km?1 , or dz/dT = 50 mK ?1 . At this thermal conductivity, over 20 meters there would be a thermal jump of 0.4 K. Summing these temperatures we find that the temperature rise from the top of BL3 to its bottom would be about 18 K. For a surface temperature of 150 K, this would imply a temperature of 168 K at the top of the unit, Aa3c . Applying the extrema of our estimate of the geothermal flux at the south pole (see Eq. 2), we find that the temperature at the top of the CO2 slab would be approximately CO2 Ttop = 168 K. It is, however, likely that in time, such as the 20 kyr since we believe unit Aa3c formed, sintering will advance and raise the thermal conductivity of the top boundary layer. Thus, we suggest that 168 K is likely an upper limit. We believe that Ttop = 163 K could reasonably represent the minimum. An average of the extremes gives, CO2 Ttop = 165.5 �4 K, (5) M where we use the extremal values of H (Eq. 1) to generate the errors. 4.2. The thermal structure of a CO2 slab CE PT ED Here we treat the CO2 deposit as a single slab. We have a top temperature, Ttop = 165.5 � 3.4 K, and a bottom temperature Tbott ? 217 K, where the upper limit is the triple point of CO2 . As a simplification, we neglect the boundary layers that divide the CO2 slab, since they have negligible impact on the thermal gradient of the CO2 . We consider two cases, first of which is a CO2 slab whose thickness is less than or equal to that at which can melting can be approached. Following that, we consider the melting rate with a still deeper deposit. AC 4.2.1. Depth of melting Let us consider a single slab of thickness z ? such that the temperature at the bottom just reaches the triple point of CO2 without causing any melting. In this case, the slab would transmit the full value of the geothermal flux, H. Integrating Eq. 3, and solving for thickness, z ? , we have, Tbott 93.4 ? . (6) ln z = H Ttop Substituting the central value of the temperature, Ttop = 165.5 K, from Eq. 5, and Tbott = 217 K into Eq. 6, the thickness of a CO2 slab that For a CO2 slab of thickness z 0 > z ? , we expect basal melting to occur. In this case, the geothermal flux through the slab at equilibrium, H 0 , would be less than the current average value, assumed to be H = 20 (� mW m2 , since the thermal gradient at equilibrium would be more gradual, ?T kef f , (7) z0 where ?T is the difference between the triple point and the temperature, Ttop , at the top of the CO2 slab. We let kef f = 0.5 Wm?1 K?1 , as above. Since the heat flux that is devoted to melting CO2 ice is H ? H 0 , the melting rate can be expressed as follows, ?T kef f 1 d2 MCO2 H ? H0 1 = = ? , (8) dA dt ?f us H 0 ?f us H 0 z ? z 0 H0 = where the heat of fusion for CO2 is ?f us H 0 = 8347 J mol?1 or 1.897 � 105 J kg?1 . By dividing Eq. 8 by the density of CO2 (? 1600 kg m?3 ), we have the vertical melting rate at depth z 0 > z ? . Results for a range of z 0 under the three conditions discussed above, are shown in Figure 5, where respective line types are listed in Column 3 of Table 2. Melting rates are given in millimeters per year. In reality, the uneven topography containing the RFZs guarantees that each deep deposit would have a distribution of melting rates, so that the weighted mean melting rates of a deposit would likely be significantly smaller than at the deepest parts. Table 2 Minimal slab thickness for basal melting CO2 H (mW m?2 ) ?T Ttop z ? (m) L-type 15 12.1 162.1 1782 long-dash 20 15.5 165.5 1195 short-dash 25 18.9 168.9 868 solid ACCEPTED MANUSCRIPT 10 AC CE PT ED M AN US Our analysis of CO2 ice deposit survival is based on observations of B2016, the analytical work of Clancy et al. (1996), and simulations from the Manning et al. (2006) model of atmospheric evolution. This analysis shows the consistency of the identifications of the obliquity minima in the last 350 Ma with the layers suggested by B2016. It should be noted, however, that it is still possible that the series at 2.6 ? t ? 2.3 Ma could theoretically be the source of the deposits, but many of the the series of high-obliquity excursions between 2.2Ma and 500 ka would likely have ablated and sublimed the deposits, leaving the scarps exposed for the later series of obliquity minima at t < ? 400 ka. In addition, the idea of the bifurcated boundary layer, BL2 , suggests that the essentially zero porosity 40 m layer suggested by B2016 analysis is actually two layers of water ice, perhaps 15 m thick each, plus one layer of CO2 about 10 m thick. In this case, we would have three boundary layers of thickness ? < 20 m composed of very low porosity ice (not counting BL3 ). Correcting for 70% porosity when they were the dry snow was first formed, they would have been from 45 to 67 meters thick. While the Clancy effect correctly predicts the preservation of the CO2 deposits we do find, it may be less efficient in predicting the absence of the one at 450 ka, and the relatively less-robust deposit of 350 ka. It may be that these failures are due to the relatively low and declining eccentricity during the perihelion period at northern summer (southern winter) solstice, which would generally reduce the net transport of water. However, the moderately high obliquity maximum (? ' 30? ) following deposition of the 400 ka layer may also have helped to erode what may have been an already thin layer of snow. The cause of the variability of the deposit size of unit Aa3c with altitude on the cap cannot be fully explained by eccentricity or other orbital factors since the obliquity maximum that followed its emplacement barely exceeded 26? . For this, we may have to consider mesoscale effects. The analysis of the thermal structure of the top boundary layer suggests that most top boundary layers of water ice, when first formed, could have been significantly thicker than the current one, perhaps 40 to 75 meters thick, given an initial 70% porosity. However, the time scale for densification may be rapid (e.g., ? ? 14 kyr) due to temperature-accelerated sintering and compression in the lower half of the layer. For CO2 deposits, the thermal gradient is very sensitive to the geothermal flux. However, for the central value, H = 20 mWm?2 , the pont where basal melting begins is only about 200 meters below the current deepest CO2 deposit (B2016). For depths in excess of that at which temperatures just reach the triple point of CO2 , the temperature gradient of the slab (and thermal flux) is less, with the excess devoted to basal melting. Although the limits on H suggest that it is not likely that basal melting is currently occurring, it is not excluded. Several million years ago, during the construction of the SPLD, there could have been troughs and valleys deep enough to accommodate significant basal melting. The opportunity for sequestration into the regolith would be limited by whether the deposits touched or approached the underlying regolith. In that case, the exchangeable inventory of CO2 on Mars could have been significantly diminished by earlier episodes of the formation of deep frozen CO2 deposits, especially at a time when the geothermal flux was greater, as noted by Kurahashi-Nakamura and Tajika(2006) CR IP T 5. SUMMARY AND CONCLUSIONS Conclusions We have analyzed the conditions under which buried CO2 ice caps could occur, considering Shallow Radar probes of existing deposits, their thermal structure, and the possibility of basal melting of CO2 within the context of a declining atmospheric pressure and geothermal flux. Based on this we reach the following conclusions: 1. The Clancy effect, in which the alignment of perihelion with northern summer when the obliquity has reached a minimum and a CO2 deposit has been formed, can explain the burial and stabilization of the CO2 layers in the SPLD. 2. A second important stabilization factor for CO2 deposits is the sintering and densification of the top boundary layer with depth to a density greater than the ? 800 kg m?3 pore cut-off limit, which seals in the CO2 . This allows the deposit to survive the following high-obliquity phase of quasi-periodic oscillations, when the CO2 ice would otherwise be unstable. 3. The depth at which basal melting may occur under the current geothermal flux is possibly only about 200 meters deeper than the deepest of the current deposits, opening up the possibility of basal melting and sequestration of CO2 , both in the recent and distant past. This has implications for the evolution of the Martian atmospheric pressure with time. REFERENCES Bierson, C. J., Phillips, R. J., Smith, I. B., Wood, S. E., Putzig, N. E., Nunes, D. and Byrne, S., 2016. (B2016). Stratigraphy and evolution of the buried CO2 deposit in the Martian south polar cap. 2016 Geophys. Res. Lett., 43,4172-4179. ACCEPTED MANUSCRIPT 11 CR IP T Mellon, M. T. and Fergason, R. L. and Putzig, N. E., 2008. The thermal inertia of the surface of Mars, in The Martian Surface: Composition, Mineralogy and Physical Properties, ed. by J.F. Bell III (Cambridge University Press, Cambridge, 2008). Mellon, M. T. and McKay, C. P. and Grant, J. A., 2015. The Thermal Conductivity of Soils with Bimodal Grain-Sizes at Mars Pressures. LPSC 46, 2837. Milkovich, S. M. and Plaut, J. J., 2008. Martian South Polar Layered Deposit stratigraphy and implications for accumulation history. J. Geophys. Res., 113, E06007. Montmessin, F. and Forget, F. and Rannou, P. and Cabane, M. and Haberle, R. M., 2004. Origin and role of water ice clouds in the Martian water cycle as inferred from a general circulation model. J. Geophys. Res.(planets), 109, E18, doi:10.1029/2004JE002284 Montmessin, F. and Haberle, R. M. and Forget, F. and Langevin, Y. and Clancy, R. T. and Bibring, J.-P., (2007). On the origin of perennial water ice at the south pole of Mars: A precession-controlled mechanism?. J. Geophys. Res., 112, E08S17,10.1029/2007JE002902, http://adsabs.harvard.edu/abs/2007JGRE..112.8S17M Nye, J. F., Durham, W. B., Schenk, P. M. and Moore, J. M., 2000. The Instability of a South Polar Cap on Mars Composed of Carbon Dioxide. Icarus, 144, 449-455. Phillips, R. J. and Zuber, M. T. and Smrekar, S. E. and 24 coauthors, 2008. Mars North Polar Deposits: Stratigraphy, Age, and Geodynamical Response. Science, 320, 1182, doi:10.1126/science.1157546. Phillips, R. J. and Davis, B. J. and Tanaka, K. L., and 15 coauthors 2011. Massive CO2 ice deposits sequestered in the south polar layered deposits of Mars. Science, 332, 838841. Putzig, N.E. and Phillips, R. J. and Campbell, B. A. and Holt, J. W. and Plaut, J. J. and Carter, L. M. and Egan, A. F. and Bernardini, F. and Safaeinili, A. and Seu, R., 2009. Subsurface structure of Planum Boreum from Mars Reconnaissance Orbiter Shallow Radar soundings. Icarus, 204, 443-457. doi:10.1016/j.icarus.2009.07.034 Putzig, N. E. and Smith, I. B. and Perry, M. R. and Foss, F. J. and Campbell, B. A. and Phillips, R. J. and Seu, R., 2018. Three-dimensional radar imaging of structures and craters in the Martian polar caps. Icarus, 308, 138-147. doi:10.1016/j.icarus.2017.09.023. and heat transfer through the snow on the ice of Beaufort Sea. J. Geophys. Res., bf 107, 8043. Smith, I. B., and Holt, J. W., (2010). Onset and migration of spiral troughs on Mars revealed by orbital radar, Nature, 465(7297), 450-453. Tanaka, K. L. and Kolb, E. J. and Fortezzo, C., 2007. Recent Advances in the Stratigraphy of the Polar Regions of Mars. Seventh Int. Conf. on Mars, Pasadena, Ca., July 7-13. Toon, O. B. and Pollack, J. B. and Ward, W. and Burns, J. A. and Bilski, K., 1980. The astronomical theory of climatic change on Mars. Icarus, 44, 552-607. Ward, W. R., 1979. Present obliquity oscillations of Mars: Fourth-order accuracy in orbital e and i. J. Geophys. Res., 84, 237-241. Zent, A. P. and Quinn, R. C., 1995. Simultaneous adsorption of CO2 and H2 O under Mars-like conditions and application to the evolution of the Martian climate. J. Geophys. Res., 100, 5341-5349. AC CE PT ED M AN US Clancy, R. T. and Grossman, A. W. and Wolff, M. J. and James, P. B. and Rudy, D. J. and Billawala, Y. N. and Sandor, B. J. and Lee, S. W. and Muhleman, D. O., 1996. Water Vapor Saturation at Low Altitudes around Mars Aphelion: A Key to Mars Climate? Icarus, 122,36-62. Clifford, S. M., 1987. Polar basal melting on Mars. J. Geophys. Res., 92, 9135-9152. Edgett, K. S. and Christensen, P. R., 1991. The particle size of Martian aeolian dunes. J. Geophys. Res., 96, 22. Forget, F. and Hourdin, F. and Fournier, R. and Hourdin, C. and Talagrand, O. and Collins, M. and Lewis, S. R. and Read, P. L. and Huot, J.-P., 1999. Improved general circulation models of the Martian atmosphere from the surface to above 80 km. J. Geophys. Res., 104, 24155-24176, doi:10.1029/1999JE001025 Heldmann, J. L., Toon, O. B., Pollard, W. H., Mellon, M. T., Pitlick, J., McKay, C. P. and Andersen, D. T., 2005. Formation of Martian gullies by the action of liquid water flowing under current Martian environmental conditions. J. Geophys. Res., 110, E05004. Herron, M. M. and Langway, Jr., C. C., 1980. Firn densification: an empirical model. JGlac, 25, 373-385.doi: 10.1017/S0022143000015239 Holt, J. W., Fishbaugh, K. E., Byrne, S., Christian, S., Tanaka, K., Russell, P. S., Herkenhoff, K. E., Safaeinili, A., Putzig, N. E., and Phillips, R. J. (2010). The construction of Chasma Boreale on Mars, Nature, 465(7297), 446-449, doi:10.1038/nature09050. Kaempfer, T. U. and Schneebeli, M., 2007. Observation of isothermal metamorphism of new snow and interpretation as a sintering process. J. Geophys. Res., 112, D24101. Kolb, E. J., and Tanaka, K. L. (2006), Accumulation and erosion of south polar layered deposits in the Promethei Lingula region, Planum Australe, Mars, Mars, 2, doi:10.1555/mars.2006.0001. Koutnik, M. and Byrne, S. and Murray, B., 2002. South Polar Layered Deposits of Mars: The cratering record. J. Geophys. Res., 107, 10-1. Kravchenko, Yu. G. and Krypskii, I. N., 1986. Thermal conductivity of solid N2 O and CO2 . Sov. J. Low Temp. Phys. 12(1), 79-83. Kurahashi-Nakamura, T. and Tajika, E., 2006. Atmospheric collapse and transport of carbon dioxide into the subsurface on early Mars. Geophys. Res. Lett., 33.18, 205. Laskar, J., Correia, A. C. M., Gastineau, M., Joutel, F., Levrard, B. and Robutel, P., 2004. Long term evolution and chaotic diffusion of the insolation quantities of Mars. Icarus, 170, 343-364. Laskar, J., Levrard, B., and Mustard, J. F., 2002. Orbital forcing of the martian polar layered deposits. 2002 Nature, 419, 375?377. Manning, C. V., McKay, C. P. and Zahnle, K. J., 2006. Thick and thin models of the evolution of carbon dioxide on Mars. Icarus, 180, 38-59. Mellon, M. T., (1996). Limits on the CO2 content of the Martian polar deposits. Icarus, 124, 268?279. ACCEPTED MANUSCRIPT CR IP T 12 AN US Fig. 1.? The Laskar et al. (2002) values for the obliquity during the last 4 Myr. The ?plus? sign marks the present obliquity. Red (short-dashed) line indicates the approximate obliquity at which a 12 mbar atmosphere would collapse. Green (long dashed)is the collapse obliquity for a 140 mbar atmosphere. The blue solid line shows the obliquity at which a perennial CO2 ice cap becomes unstable. The positions of the colored lines are based on the model of Manning et al., 2006 (see below). M 0.15 0.1 0.05 ED 0 300 200 100 0 PT 30 AC CE 20 Atmosphere 100 Regolith 10 1 -0.5 -0.4 -0.3 -0.2 -0.1 0 Time (Myr) Fig. 2.? The last 0.52 Myr of Mars evolution. From top down, panels show the eccentricity, the angle of perihelion, Lp , and the obliquity, all from from Laskar et al. (2004), and on the bottom, a simulation of the atmospheric evolution model MEC (Manning et al., 2006) showing the trends of atmospheric pressure, in black, the regolith content, in red, and the CO2 ice in blue. The current pCO2 is denoted by the orange ?+? sign. Vertical lines are at Lp values 0? (blue; short-dash), 90? (red; dot-dash), and 180? (blue; short-dash) intended to guide the eye for alignment with obliquity and CO2 ice cap. The RFZ unit names corresponding to buried deposits are shown in the bottom panel, and in Table 1. See text for discussion. ACCEPTED MANUSCRIPT CR IP T 13 M AN US Fig. 3.? Radargram from B2016 Fig. 2a (SHARAD observation 5968-01). The maximum depth of RFZs is a little over 1 km. Note the bifurcation of BL2 , suggesting that the little obliquity dip at ? 125 ka has produced some buried CO2 . ED 0.15 0.1 0.05 0 300 PT 200 100 0 AC CE 30 20 Atmosphere 100 Regolith 10 1 -2.7 -2.6 -2.5 -2.4 -2.3 -2.2 Time (Myr) Fig. 4.? Same as Figure 2 except that we explore the time range, 2.75 Ma < t < 2.20 Ma, when obliquity variations are small (see Fig. 1). We show results for two different atmospheric pressures, PCO2 = 6.5 mbar and 12 mbar; the traces corresponding to a maximum 6.5 mbar atmospheric inventory have dotted lines. ACCEPTED MANUSCRIPT PT ED M AN US CR IP T 14 AC CE Fig. 5.? The vertical melting rates derived from Eq. 8 for H = 15, 20, and 25 mW m?2 (long-dash, short-dash and solid lines, respectively). The dotted lines around the short-dashed lines are the result of applying the maxima and minima of Ttop of Eq. 5 to Eq. 6 and then to Eq. 8. The ordinate gives vertical melting rate in mm/yr. To transform vertical melting rates to areal melting rates, multiply by ?CO2 ' 1600 kg m?3 .