CLIMATOLOGICAL CONTROLS ON THE CHEMICAL STRATIGRAPHY OF
THE MARTIAN POLAR LAYER DEPOSITS.
E. Vos, O. Aharonson1 , Dept. of Earth & Planetary Sciences, Weizmann Institute of Science, Rehovot, Israel
(Eran.Vos@weizmann.ac.il), 1 also at: Planetary Science Institute, Tucson, Arizona, USA., N. Schörghofer, Planetary
Science Institute, Tucson, Arizona, USA., F. Forget, E. Millour, Laboratoire de Météorologie Dynamique/IPSL,
Sorbonne Université, ENS, PSL Research University, Ecole Polytechnique, CNRS, Paris, France, L. Rossi, M. Vals,
F. Montmessin, Laboratoire Atmosphères, Milieux, Observations Spatiales (LATMOS/CNRS), Paris, France.
Introduction
The growth of the North Polar Layered Deposits (NPLD)
at the expense of tropical ice is thought to have started
approximately 4.5 Myr ago [1] when the mean obliquity
dropped from ∼35◦ to ∼25◦ [2]. During this epoch, the
obliquity has varied between 15◦ and 35◦ , and the eccentricity between 0 and 0.13 [2]. Alternating layers seen
in the NPLD in images [3] and radar data [4] provide an
archive of past climate oscillations, that are widely held
to be linked to orbital variations[5, 1, 6, and many others]. Ice also exists in the subsurface; both theoretical
and observational work show that shallow ground-ice is
abundant in the mid to high latitudes [7, 8, 9] The volume
fraction and extent of the subsurface ice can alter the atmospheric energy budget by storing and releasing heat
at different times, and this in turn can influence the physical and chemical evolution of the polar caps [9, 10]. On
Earth, The isotopic signal in ice cores from Antarctica
and Greenland have revealed past climate oscillations
[11]. Mars lacks an ocean to buffer the atmospheric
isotopic ratio as on Earth, thus the isotopic composition
of the present atmosphere is controlled by the evaporative fluxes of the caps. Observations from orbiters and
ground-based telescopes measured the atmospheric D/H
ratio to be in the range of 1-8 × Vienna Standard Mean
Ocean Water (VSMOW) [12, 13, 14], with an average of
4.6±0.7 [14]. The atmospheric values provide evidence
for both net water loss in the long-term evolution of
Mars [15, 13], as well as seasonal, cloud physics [16],
and geographic source effects. Previous work [5, 17]
shows a D/H record in the polar deposits is expected,
although potentially complex to interpret than Earth’s
due to the aforementioned effects. In this work, we seek
to quantify the processes that affect the D/H ratio of
the accumulating ice on the PLDs, in order to provide
a framework for interpreting the chemical signal in a
vertical profile of polar ice.

Model
The results reported in this paper are based on simulations of the Mars LMD-GCM with the full water cycle that includes treatment of surface ice, atmospheric
vapor and ice clouds and has been described in detail

previously [18, 19, 16, 20, 17]. In addition to dust
and water tracers [16], we use the HDO tracer [21] in
both the vapor and ice phases to track the D/H ratio.
Temperature-dependent isotopic fractionation of water
occurs at condensation either in clouds or directly onto
the surface. The equilibrium solid-vapor fractionation
coefficient (α) assumed is based on lab measurements
[22]. Ice clouds are taken to be radiatively passive.
The simulations reported here have a horizontal resolution of 64 by 64, corresponding to 2.8125◦ in latitude,
and 5.625◦ in longitude, with 29 vertical layers extending to 80 km. The simulations were performed with
an assumed past ice distribution in which surface ice is
present in the tropics at locations on the eastern flanks
of the Tharsis rise where accumulation is predicted at
high obliquity [23]. The presence of such a past tropical
surface ice reservoir is supported by geologic evidence
[24]. We tested thick and thin initial tropical ice layers.
In the thick case, the tropical ice layer has a thickness
of ∼60 m Polar Equivalent Layer (PEL), defined as a
layer of the same volume, spread evenly from 81◦ to
90◦ latitude. This case is designed to reach steady state,
to investigate the fluxes and D/H ratio of ice reaching
the polar cap with a persistent tropical source. In the
thin test case, we initialize the model with a ∼0.45 m
PEL tropical ice layer, in order to investigate the different stages of the evolution of ice migration as the finite
initial source is exhausted.

Results
Condensation of water vapor leads to fractionation. Figure 1 shows the zonal mean seasonal evolution of atmospheric water vapor and the D/H ratio relative to the
initial value of the source. We plot three cases, that
are different in orbital state and ice distribution. Atmospheric humidity is significantly higher when a tropical
ice reservoir is available, and peaks near mid-summer in
each pole. In general, the atmospheric vapor is isotopically depleted because the HDO is either preferentially
in ice clouds or in surface ice that sedimented at an
earlier stage. This phenomenon is amplified when a
tropical reservoir is present [17]. As described in the
previous section, we performed simulations with initial
conditions that assume ice deposits placed in the trop-

Figure 1: Zonal mean seasonal evolution of atmospheric water vapor and its D/H ratio for three different states: Current orbit and
ice distribution (left column), current orbit and a thick tropical ice distribution (middle column), and an orbit with Lp =90◦ and a
thick tropical ice distribution (right column).

ics. Figure 2a shows The NPLD growth rate for a thick
tropical ice distribution as a function of obliquity, , and
perihelion longitude, Lp , for present-day eccentricity. In
addition to obliquity controlling the annual mean insolation as a function of latitude, and hence the sign and
magnitude of polar accumulation, there is also a strong
dependence on Lp . This phase controls the length and
intensity of the high-humidity summer season, with a
difference of a factor of ∼ 2 in accumulation rate between a prolonged mild northern summer to a short intense one. We also calculate the mean D/H ratio of ice
that condensed in the polar region and plot it in Figure 2b.
The contours show D/H values in δ D. The hydrogen isotopic anomaly depends on the orbital elements, and in
particular Lp , with a difference of more than 100h over
a precession cycle. A short, intense summer leads to
higher D/H ratio of the accumulating ice. Notably, the
polar δ D value is negative for all cases. This occurs because as vapor transports from the tropics poleward, the
water molecules first deposit in the intermediate latitudes
and preferentially condense out the heavy isotope. The
atmosphere becomes depleted in the HDO and therefore water reaching the pole is isotopically light. Under
the assumption that the majority of the tropical source
ultimately migrates to the polar cap, the depletion in
HDO of the early deposits would be compensated by
enrichment of subsequent deposits sourced from the intermediate latitudes. In order to test if the simulated ice
core is affected by the history of the migration we considered an initial ice distribution with a smaller amount
(thin tropical ice layer test case), placed in the equatorial
region. We use this approach as an approximation, as
the stages of the migration of this limited deposit over
the simulation timescale reflects those of a much larger
deposit, migrating over a much longer time (for more
details on the model assumptions and the details of the

ice migration see [17]). At first, as long as the tropical
reservoir is present, ice migrates to latitudes higher than
±45◦ , but only a small portion of the ice that sublimated
from the tropics reaches the polar caps. When ice is no
longer available in the tropics, the planet’s mean humidity drops and ice migrates from these quasi-stable regions towards higher latitudes, poleward of about ±55◦ .
As the simulation advances and more ice migrates, the
location of the source (the equator-most ice still available) migrates poleward. At the last stage, when the
boundary has migrated to ∼ ±75◦ , the ice accumulating
on the caps is sourced only from these high latitude locations (either in the same hemisphere, or the opposite
one). Thus, the upper layers of the cap are expected
to represent a deposit that migrated from high latitudes
(polewards of 70◦ ), and not directly from the tropical
reservoir. The migration stages mentioned above affect
the D/H ratio of the condensed ice as can be seen in Figure 3. Different hydrogen isotopic sections appear along
the profile due to the different migration stages and ice
source locations mentioned above. The bottom portion
of the profile (up to ∼0.2 ice fraction transported) is
highly depleted in HDO relative to the equatorial source
in both poles, due to the mid-high latitude reservoir fractionating the atmosphere as it grows. Once the tropical
source is depleted, the enriched mid-high latitude deposit becomes the source, and the growing polar ice
is therefore also enriched. As the polar ice grows at
the expense of ever-higher latitude ice, its isotope ratio
changes in accordance with the source, with a smaller
modification due to heavy isotope trapping in intermediate latitudes between the source and the pole. In both
hemispheres, and for all profiles examined, our simulations show that the upper layers of each PLD is enriched
in HDO compared to the average D/H ratio of that PLD
(dashed line). The effect is due to the HDO depletion

manuscript submitted to JGR Planets
REFERENCES

a)

b)
-120

3

-160

-160

2

-160

3

-200
-20

2

-200

0

2

-160

3
3

3

-160

-160

-120

-120

-120
2

-80

3

-80

4

-80

3

4

3

-120

-120

Figure 2: a) Simulation results of the North Polar ice flux as function of obliquity and perihelion longitude for eccentricity of

Figure
3. polar
a) Simulation
results
of the
North
Polar ice flux
0.093. b) North
ice accumulation
D/H ratio
for the
same configuration
as in as
a. function of obliquity noted as
.

✏ and Lp for eccentricity of 0.093 with a thick tropical reservoir. A short and more intense sum-

213

humidity. We tested how the position and distribution
the subsurface
affects
simulated
polar profile.
mer leads to higher accumulation rate up to obliquityof⇠27,
beyond ice
which
it the
decreases
rapidly.
We ran two sets of simulations, one with ground-ice at
a depth of 15 cm extending to latitude 60◦ similar to at
b) North polar ice accumulation D/H ratio for the same
orbital configuration and ice distribupresent, and one where the thermal inertia in the subsurface was set to mimic the presence of an ice sheet
tion as in a. A short, intense summer leads to higheratD/H
ratio.
isotopic
is mostly
a depth
of 5 The
cm that
extendslayering
to latitude
30◦ in both
hemispheres, similar to past conditions. The difference
determined by Lp .
between the two test cases as function of Lp is shown
in Figure 4, where (a) shows the change in flux to the
. NPLD, and (b) the change in the D/H ratio of the condensed ice. Both the flux and the D/H ratio are higher
in all Lp tested for the first case (margin at 60◦ latitude).
This phenomenon occurs because when the subsurface
extend further towards the equator and closer to the surface it acts to reduce the difference between seasons;
on Lp . This phase controls the length and intensity
of the
summer
sea-leads to
as shown
in high-humidity
Figure 2 a more intense
summer
higher accumulation and a larger D/H ratio. This result
son, with a di↵erence of a factor of ⇠2 in accumulation
rate
between
a prolonged
mild
supports the
claim
that the subsurface
extent,
depth, and
Figure 3: North and South Polar D/H profiles as a function
ice
volume
(not
shown
here)
will
affect
the
physical
and
◦
of the fraction
of ice transported
for obliquity
of 25
, eccennorthern
summer
to a◦ short
intense
one.
Higher chemical
eccentricity
simulations
show ice
this
ef- The restratigraphy
of a simulated
core.
tricity of 0.093 and Lp of 90 and an initial thin tropical ice
sults shown in this abstract and in previous published
distribution.
The data and
pointsthe
are separated
one Marscan
year,reachwork
[17] highlight
the fact
in order to interpret
fect
is amplified
NPLDbygrowth
6 mm/Mars
year.
Atthat
obliquities
higher the
and the dashed line is the mean D/H ratio of each final ice cap.
chemical stratigraphy of the PLDs a model that takes
into account
both
the atmosphereand
evolution
over seathan 25 the pattern changes, the dependence becomes
more
complicated,
the mean
sonal timescales, and the ice reservoirs evolution over
timescales,
is needed.
In theThis
conference
of polar ice that accumulates
early along
witha lower
lat-peak millennial
accumulation
rate no longer
shows
single
in Lp over
a precession
cycle.
oc- we
will discuss the development of such a model.
itude enriched deposits. The simulations also show that
the southern
ice with
deposits
enrichedhumidity
in deuterium
com- at moderate obliquities there is continucurs
because
a are
tropical
source
pared to the north. Subsurface ice can affect the NPLD
References
physical stratigraphy via several mechanisms. It acts as
ous
polar accumulation throughout the year. At obliquities ⇠30 and above, polar ice
a water source, as well as a heat reservoir altering the
[1] B. Levrard et al. “Recent formation and evolution of
energy budget. Previous modeling work [9] that is supbegins
to observations
sublimate [7]
during
season; at even
higher
obliquities,
NPLD
fluxfrom
northern
Martian
polar layeredthe
deposits
as inferred
ported by
showedthe
thatsummer
the subsurface
a Global Climate Model”. In: J. Geophys. Res-Planet.
ice distribution can extend to ∼30◦ or retreat to ∼ 60◦ ,
isdepending
at a state
of orbital
net loss.
Overall,andthe
simulations show112.E6
the growth
(2007). rate of the north poon the
configuration
atmospheric

214

lar cap is 1-6 mm/Mars year, while loss rates can reach 100 mm/Mars year. Levrard et

215

al. (2007) found similar loss rates. The results in Figure 3a imply that the stratigraphy

216

in the cap is controlled not just by obliquity cycles, but also strongly, by precession cy-

217

cles.

204

205

206

207

208

209

210

211

212

218

To study the hydrogen isotopic signal expected in the Martian polar cap, we an-

REFERENCES

b)

North Pole Flux (PEL\Mars year)

a)

Figure 4: a) North polar flux as function of Lp for two subsurface ice cases, the first assumes ice from the pole to latitude 60◦ at a
depth of 15 cm, and the second, to latitude 30◦ at depth of 5 cm. b) same as a, only the D/H ratio of the condensed ice.
.
[2]

[3]

[4]

[5]
[6]

[7]

[8]

[9]

[10]

[11]

[12]

[13]

J. Laskar et al. “A long-term numerical solution for the
insolation quantities of the Earth”. In: Astron. Astrophys. 428.1 (2004), pp. 261–285.
K. E. Fishbaugh et al. “Martian north polar layered deposits stratigraphy: Implications for accumulation rates
and flow”. In: J. Geophys. Res-Planet. 111.E6 (2006).
Isaac B. Smith et al. “An ice age recorded in the polar
deposits of Mars”. In: Science 352.6289 (2016), pp. 1075–
1078.
E. Vos et al. “Dynamic and isotopic evolution of ice
reservoirs on Mars”. In: Icarus 324 (2019), pp. 1–7.
Isaac B. Smith et al. “The Holy Grail: A road map for
unlocking the climate record stored within Mars’ polar
layered deposits”. In: Planet. Space. Sci 184 (2020),
p. 104841.
G.A. Morgan et al. “Availability of subsurface waterice resources in the northern mid-latitudes of Mars”.
In: Nature Astronomy 5 (2021), pp. 230–236.
Norbert Schorghofer et al. “Stability and exchange of
subsurface ice on Mars”. In: J. Geophys. Res-Planet.
110.E5 (2005).
E. Vos et al. “The Role of Mid-Latitude Ground Ice
in North Polar Layered Deposits Formation”. In: Mars
Polar Conference abstract (2020).
Robert M. Haberle et al. “The effect of ground ice on
the Martian seasonal CO2 cycle”. In: Planet. Space. Sci
56.2 (2008), pp. 251–255.
R. Alley, ed. The Two-Mile Time Machine: Ice Cores,
Abrupt Climate Change, and Our Future. Princeton University Press, 2000.
B.M. Jakosky et al. “Loss of the Martian atmosphere to
space: Present-day loss rates determined from MAVEN
observations and integrated loss through time”. In: Icarus
315 (2018), pp. 146–157.
G. L. Villanueva et al. “Strong water isotopic anomalies
in the martian atmosphere: Probing current and ancient
reservoirs”. In: Science 348.6231 (2015), pp. 218–221.

[14]

[15]

[16]

[17]

[18]

[19]

[20]

[21]

[22]

[23]

[24]

V. Krasnopolsky. “Variations of the HDO/H2O ratio in
the martian atmosphere and loss of water from Mars”.
In: Icarus 257 (2015), pp. 377–386.
N. R. Alsaeed et al. “Mars Water and D/H Evolution
From 3.3 Ga to Present”. In: J. Geophys. Res-Planet.
124.12 (2019), pp. 3344–3353.
F. Montmessin et al. “Modeling the annual cycle of
HDO in the Martian atmosphere”. In: J. Geophys. ResPlanet. 110.E3 (2005).
E. Vos et al. “Stratigraphic and Isotopic Evolution of
the Martian Polar Caps From Paleo-Climate Models”.
In: Journal of Geophysical Research: Planets 127.3 (2022),
e2021JE007115.
François Forget et al. “Improved general circulation models of the Martian atmosphere from the surface to above
80 km”. In: J. Geophys. Res-Planet. 104.E10 (1999).
F. Montmessin et al. “Origin and role of water ice clouds
in the Martian water cycle as inferred from a general
circulation model”. In: J. Geophys. Res-Planet. 109.E10
(2004).
J.-B. Madeleine et al. “Revisiting the radiative impact
of dust on Mars using the LMD Global Climate Model”.
In: J. Geophys. Res-Planet. 116.E11, E11010 (2011).
Loı̈c Rossi et al. “The effect of the Martian 2018 global
dust storm on HDO as predicted by a Mars Global Climate Model”. In: Geophys. Res. Lett. (2021).
L. Merlivat et al. “Fractionnement isotopique lors des
changements d‘état solide-vapeur et liquide-vapeur de
l’eau à des températures inférieures à 0◦ C”. In: Tellus
19.1 (1967).
F. Forget et al. “Formation of Glaciers on Mars by Atmospheric Precipitation at High Obliquity”. In: Science
311 (2006), pp. 368–371.
J. W. Head et al. “Tropical to mid-latitude snow and ice
accumulation, flow and glaciation on Mars”. In: Nature
434 (2005), pp. 346–351.