Century-scale simulations of the West Antarctic Ice Sheet

56 Pages·2015·4.86 MB·English

TheCryosphereDiscuss.,9,1887–1942,2015 D
www.the-cryosphere-discuss.net/9/1887/2015/ isc TCD
u
doi:10.5194/tcd-9-1887-2015 s
s
©Author(s)2015.CCAttribution3.0License. io 9,1887–1942,2015
n
P
a
Thisdiscussionpaperis/hasbeenunderreviewforthejournalTheCryosphere(TC). p
e Century-scale
PleaserefertothecorrespondingfinalpaperinTCifavailable. r
simulations of the
|
West Antarctic Ice
Century-scale simulations of the D
Sheet
is
c
u
response of the West Antarctic Ice Sheet s S.L.Cornfordetal.
s
io
n
to a warming climate P
a
p
e TitlePage
r
1 2 1 2 3
S. L. Cornford , D. F. Martin , A. J. Payne , E. G. Ng , A. M. Le Brocq , | Abstract Introduction
1 1 1 4,5
R. M. Gladstone , T. L. Edwards , S. R. Shannon , C. Agosta ,
M. R. van den Broeke6, H. H. Hellmer7, G. Krinner4, S. R. M. Ligtenberg6, Dis Conclusions References
c
R. Timmermann7, and D. G. Vaughan8 us Tables Figures
s
io
1CentreforPolarObservationandModelling,SchoolofGeographicalSciences,Universityof n
(cid:74) (cid:73)
P
Bristol,Bristol,BS81SS,UK a
2ComputationalResearchDivision,LawrenceBerkeleyNationalLaboratory, pe (cid:74) (cid:73)
r
Berkeley,California,USA
3Geography,CollegeofLifeandEnvironmentalSciences,UniversityofExeter, | Back Close
Exeter,EX44RJ,UK D
4UJF-Grenoble1/CNRS,LaboratoiredeGlaciologieetGéophysiquedel’Environnement isc FullScreen/Esc
u
(LGGE)UMR5183,38041Grenoble,France s
5DépartementdeGéographie,UniversitédeLiège,Liège,Belgium sio Printer-friendlyVersion
n
6InstituteforMarineandAtmosphericResearch,UtrechtUniversity,Utrecht,theNetherlands P InteractiveDiscussion
7Alfred-Wegener-InstitutfürPolarundMeeresforschung,Bussestrasse24, ap
e
27570Bremerhaven,Germany r
1887 |
8BritishAntarcticSurvey,MadingleyRoad,Cambridge,CB30ET,UK D
is
c TCD
u
Received:9February2015–Accepted:23February2015–Published:23March2015 s
s
io 9,1887–1942,2015
Correspondenceto:S.L.Cornford([email protected]) n
P
PublishedbyCopernicusPublicationsonbehalfoftheEuropeanGeosciencesUnion. a
p
e Century-scale
r
simulations of the
|
West Antarctic Ice
D
Sheet
is
c
u
s S.L.Cornfordetal.
s
io
n
P
a
p
e TitlePage
r
| Abstract Introduction
D Conclusions References
is
c
us Tables Figures
s
io
n
(cid:74) (cid:73)
P
a
p
e (cid:74) (cid:73)
r
| Back Close
D
is FullScreen/Esc
c
u
s
sio Printer-friendlyVersion
n
P InteractiveDiscussion
a
p
e
r
1888 |
Abstract D
is
c TCD
u
We use the BISICLES adaptive mesh ice sheet model to carry out one, two, and three s
s
century simulations of the fast-flowing ice streams of the West Antarctic Ice Sheet. io 9,1887–1942,2015
n
Eachofthesimulationsbeginswithageometryandvelocityclosetopresentdayobser- P
a
5 vations,andevolvesaccordingtovariationinmeteoriciceaccumulation,iceshelfmelt- pe Century-scale
r
ing,andmeshresolution.Futurechangesinaccumulationandmeltratesrangefromno simulations of the
change,throughanomaliescomputedbyatmosphereandoceanmodelsdrivenbythe | West Antarctic Ice
E1andA1Bemissionsscenarios,tospatiallyuniformmeltratesanomaliesthatremove D Sheet
is
most of the ice shelves over a few centuries. We find that variation in the resulting ice c
u
10 dreysnoalumtiiocsn,isaldthoomuignhatiecdebayccthuemcuhlaotiicoenoafffiencittisalthceonndeittiocnhsa,nigceeisnhveolflummeeltarabtoevaenfldotmaetioshn ssion S.L.Cornfordetal.
toasimilardegree.Givensufficientmeltrates,wecomputegroundinglineretreatover Pa
p
hundreds of kilometers in every major ice stream, but the ocean models do not predict e TitlePage
r
such melt rates outside of the Amundsen Sea Embayment until after 2100. Sensitiv-
| Abstract Introduction
ity to mesh resolution is spurious, and we find that sub-kilometer resolution is needed
15
along most regions of the grounding line to avoid systematic under-estimates of the Dis Conclusions References
c
retreat rate, although resolution requirements are more stringent in some regions – us Tables Figures
for example the Amundsen Sea Embayment – than others – such as the Möller and s
io
Institute ice streams. n (cid:74) (cid:73)
P
a
p
e (cid:74) (cid:73)
r
1 Introduction
20
| Back Close
The present day West Antarctic Ice Sheet (WAIS) is experiencing an imbalance be- D
is FullScreen/Esc
tweenthemassitreceivesassnowfallandthatwhichitlosesthroughdischargetothe c
u
oceans (Rignot, 2008; Pritchard et al., 2009; Shepherd et al., 2012; Mouginot et al., s
2014; Rignot et al., 2014). In several areas this has led to the persistent loss of ice sion Printer-friendlyVersion
amountingtoasignificantcontributiontosea-levelrise.Continuedaccelerationofthese P InteractiveDiscussion
25 a
p
losseswouldimplyasignificantadditionalglobalsea-levelriseincomingdecadesand e
r
1889 |
centuries.PhysicallybasedprojectionsofthecontributionoftheWAIStosealevelrise D
is
are hampered by two main factors. The first of these is the lack of a fully coupled cli- c TCD
u
s
mateandicesheetmodel,inwhichtheprincipalforcingontheicesheet(accumulation s
io 9,1887–1942,2015
at the upper surface and submarine ice shelf melt) is determined within the model. n
5 The second is the technical difficulty involved in calculating the flow of ice across the Pa
p
ice sheet’s grounding line and the consequent grounding line migration. Progress has e Century-scale
r
been made in both areas (Goldberg et al., 2012a, b; Pattyn et al., 2013; Asay-Davis simulations of the
|
etal.,2015),butthecomputationalexpenseoffully-coupledice/oceanmodelsatsuffi- West Antarctic Ice
cient resolution and for sufficient integration times remains prohibitive. D Sheet
is
c
10 We approximate full coupling between the ice sheet and the rest of the climate sys- u
s S.L.Cornfordetal.
tem by imposing combinations of published accumulation and melt-rate anomaly data s
io
on the BISICLES adaptive mesh ice sheet model (Cornford et al., 2013). Two emis- n
P
sion scenarios are included: SRES E1, a mitigation scenario in which emissions are a
p
e TitlePage
stabilized by 2050 at 500ppm CO2, and A1B, a balanced scenario close to the centre r
15 of the SRES range. These were used to drive global climate warming in the UKMO | Abstract Introduction
HadCM3 and MPI ECHAM5 global climate models, which have among the highest
D Conclusions References
skill scores in the CMIP3 model group (based on Antarctic SMB, surface air tempera- is
c
ture, mean sea level pressure, and height and temperature at 500hPa, Connolley and us Tables Figures
s
Bracegirdle, 2007). The resulting global climate projections provided boundary con- io
n
ditions to two high resolution atmosphere models: RACMO2 (Ligtenberg et al., 2013) (cid:74) (cid:73)
20 P
a
andLMDZ4(Agostaetal.,2013)andtwooceanmodels:themediumresolutionBRIOS p
e (cid:74) (cid:73)
(Bremerhaven Regional Ice–Ocean Simulations) (Hellmer et al., 2012) and the higher r
resolution FESOM (Finite-element Sea ice-Ocean Model) (Timmermann and Hellmer, | Back Close
2013),ultimatelyproviding sevensetsofaccumulationdataand eightsetsofmelt-rate
D
data. is FullScreen/Esc
25 c
u
At the same time, we examine the response of the ice sheet model to variability s
beyond the scope of the atmosphere and ocean models. The climate projections de- sio Printer-friendlyVersion
n
scribed above tend to agree on the timing and magnitude of future accumulation and P InteractiveDiscussion
a
melt-rate increases, if not the distribution. We complement them with some simplified, p
e
r
1890 |
widespread melt-rate increases, as well as projections further into the future, in order D
is
to investigate the additional response to more extreme scenarios. The century-scale c TCD
u
s
evolution of the ice sheet model is also sensitive to its present day state, especially in s
io 9,1887–1942,2015
the Amundsen Sea Embayment, and we evaluate at least a part of this sensitivity – n
P
5 which will prove to be substantial – by varying the initial accumulation rate and hence a
p
the initial thinning rate. e Century-scale
r
In summary, the aim of this paper is to consider the response of the West Antarc- simulations of the
|
tic ice streams to process-based and simplified projections of future ocean and at- West Antarctic Ice
D
mosphere warming over the 21st and 22nd centuries. We focus on West Antarctica Sheet
is
c
10 primarily because of constraints on available computational resources; however these u
s S.L.Cornfordetal.
areas are also thought to be most vulnerable to future grounding line retreat because s
io
of their deep bedrock and changes in oceanic forcing (Hellmer et al., 2012; Pritchard n
P
et al., 2012; Ross et al., 2012; Joughin et al., 2014). a
p
e TitlePage
r
| Abstract Introduction
2 Methods
D Conclusions References
is
15 2.1 Model equations cus Tables Figures
s
BISICLES employs a vertically integrated ice flow model based on Schoof and Hind- io
n
marsh (2010) which includes longitudinal and lateral stresses and a simplified treat- P (cid:74) (cid:73)
a
ment of vertical shear stress which is best suited to ice shelves and fast flowing ice p
e (cid:74) (cid:73)
r
streams. Ice is assumed to be in hydrostatic equilibrium so that given bedrock eleva-
tion b and ice thickness h the upper surface elevation s is | Back Close
20
D
(cid:20) (cid:18) ρ (cid:19) (cid:21) is FullScreen/Esc
s=max h+b, 1− i h , (1) c
u
ρ s
w sio Printer-friendlyVersion
n
in which ρ and ρ are the densities of ice and ocean water.
i w P InteractiveDiscussion
a
p
e
r
1891 |
The ice thickness h and horizontal velocity u satisfy a two-dimensional mass trans- D
is
port equation c TCD
u
s
s
∂h +∇·[uh]=M −M , (2) ion 9,1887–1942,2015
∂t s b P
a
p
and two dimensional stress-balance equation e Century-scale
r
simulations of the
∇·[φhµ¯(2(cid:15)˙ +2tr((cid:15)˙)I)]+τb=ρgh∇s, (3) |
5 i West Antarctic Ice
D
Sheet
together with lateral boundary conditions. The terms on the right hand side of Eq. (2), is
c
Ms and Mb, are accumulation and melt rates. As for Eq. (3), (cid:15)˙ is the horizontal strain- uss S.L.Cornfordetal.
rate tensor, io
n
P
1(cid:104) (cid:105) a
(cid:15)˙ = ∇u+(∇u)T (4) p
2 er TitlePage
andI istheidentitytensor.Theverticallyintegratedeffectiveviscosityφhµ¯ iscomputed | Abstract Introduction
10
from the vertically varying effective viscosity µ through D Conclusions References
is
c
(cid:90)s us Tables Figures
φhµ¯(x,y)=φ µ(x,y,z)dz, (5) sio
n
(cid:74) (cid:73)
s−h P
a
p
where µ includes a contribution from vertical shear and satisfies e (cid:74) (cid:73)
r
2µA(T)(4µ2(cid:15)˙2+|ρg(s−z)∇s|2)(n−1)/2=1, (6) | Back Close
i
D
where the flow rate exponent n=3, φ is a stiffening factor (or, equivalently, φ−n is an is FullScreen/Esc
15 c
u
enhancement factor), and A(T) depends on the ice temperature T through the Arrhe- s
nius law described by Hooke (1981), sio Printer-friendlyVersion
n
P InteractiveDiscussion
(cid:18) 3f Q (cid:19) a
A(T)=A exp − (7) pe
0 [T −T]k RT r
r
1892 |
where A =0.093Pa−3a−1, Q/R =9.48×103K, f =0.53 Kk, k =1.17 and T = D
0 r is
273.39K.Thecoefficientφisestimatedbysolvinganinverseproblem(seeSect.2.4), cu TCD
s
and it is simply a convenient way to represent several conflated factors: uncertainty in s
io 9,1887–1942,2015
both temperature T and the form of A(T), macroscopic damage, and fabric formation. n
P
Finally, the basal traction is determined by a viscous law: a
5 p
e Century-scale
(cid:40)−C|u|m−1u if ρi h>−b r simulations of the
τb= ρw , (8) | West Antarctic Ice
0 otherwise
D
Sheet
is
with m=1. Like φ, C will be determined by solving an inverse problem, as described cu
s S.L.Cornfordetal.
in Sect. 2.4.1. Our choice of a linear viscous law may well bias our results toward sio
n
excessive grounding line retreat: in previous work on Pine Island Glacier non-linear
P
laws with m<1 have led to slower rates of retreat (Joughin et al., 2010; Favier et al., a
10 p
e TitlePage
2014). r
We hold the fields C and φ constant throughout our simulations. That is not to say | Abstract Introduction
thatthesefieldsoughtnotchangeoverthecourseofoneortwocenturies;forexample
D Conclusions References
regions of damage (low φ) might well propagate with the grounding line as englacial is
c
15 stresses grow in regions previously dominated by the balance between gravitational us Tables Figures
s
and basal shear stress. Rather, we lack models of sufficient skill for the present, and io
n
anticipate incorporating progress in damage models (Borstad et al., 2012) and hydrol- P (cid:74) (cid:73)
a
ogy models (Werder et al., 2013) in future calculations. We note, however, that the p
e (cid:74) (cid:73)
maps of C and φ we use (see Sect. 2.4.4) already feature slippery beds and weak r
shear margins hundreds of kilometers upstream from the grounding line. | Back Close
20
D
2.2 Model domains and boundary conditions is FullScreen/Esc
c
u
s
Wecarriedoutcalculationsonthreerectangulardomains,showninFig.1.Thelargest sio Printer-friendlyVersion
n
ofthese(RISFRIS)coverstheRossandFilchner–Ronneiceshelvesandtheirtributary
P InteractiveDiscussion
ice streams, while two smaller domains cover the Amundsen Sea Embayment (ASE) ap
e
and Marie-Byrd Land (MBL). Each of the rectangular domains is split into an active r
25
1893 |
region Ω , where ice is permitted to flow, and a quiescent region Ω where ice is D
taken to Vbe stationary. For example, in the RISFRIS domain, Ω coveQrs the present isc TCD
V u
s
day drainage basins of the Ross and Filchner–Ronne ice shelves, and the ice shelves s
themselves, while Ω covers the Amundsen Sea Embayment, Marie-Byrd Land, and ion 9,1887–1942,2015
Q
5 partoftheAntarcticPeninsula.Likewise,insidetheASEdomainΩVspansthedrainage Pap
basin of Pine Island, Thwaites, Smith, Pope, and Kohler glaciers. This construction e Century-scale
r
assumesthattheicedivideswillnotstrayfromtheircurrentconfiguration,andsolimits simulations of the
|
us to simulations over a few centuries. West Antarctic Ice
D
Reflection boundary conditions were applied at the edge of each domain. If n is Sheet
is
10 normal to a boundary and t is parallel to it, cu
s S.L.Cornfordetal.
s
u·n=0, t·∇u·n=0, ∇h·n=0. (9) io
n
P
In practice, these boundary conditions are unimportant because of the presence of a
p
quiescentregionsandcalvingfrontsinsidethedomain.Inthequiescentregions,weset e TitlePage
r
the basal traction coefficient to a large value, C=105Pam−1a so that at the interface
| Abstract Introduction
between Ω and Ω ,
15 V Q
D Conclusions References
u≈0 and uh≈0 (10) is
c
us Tables Figures
whileatthecalvingfront(whichisfixed),weimposetheusualconditionsonthenormal s
io
and transverse stress: n
(cid:74) (cid:73)
P
1 (cid:18) ρ (cid:19) a
n·[φhµ¯(2(cid:15)˙ +2tr((cid:15)˙)I)]= 2ρig 1− ρi h2n. (11) per (cid:74) (cid:73)
w
These boundary conditions, and indeed, Eq. (11) alone for a problem whose entire | Back Close
20
boundary is a fixed calving front, are sufficient provided that h(x,y,t=0) is given and D
that the basal friction coefficient C(x,y) is non-zero in at least part of the ice sheet. isc FullScreen/Esc
u
s
2.3 Adaptive mesh refinement sio Printer-friendlyVersion
n
P InteractiveDiscussion
Fine spatial resolution, or other careful treatment, is held to be crucial when simulat- ap
e
ing grounding line migration (Vieli and Payne, 2005; Durand et al., 2009). Indeed, the r
25
1894 |
BISICLESicesheetmodelwasdesignedprimarilywiththisinmind,anddiscretizesthe D
is
stressandmassbalanceEqs.(2)and(3)onblock-structuredmeshesbuiltfromrectan- c TCD
u
gular subsets of uniform grids with resolution ∆x(cid:96), with 0≤(cid:96) ≤L and 2∆x(cid:96)+1=∆x(cid:96). ss
io 9,1887–1942,2015
Whilerestrictiveinsomesenses–allmodeldomainsmustberectangular,forexample n
P
5 –thesemesheshaveasignaladvantage:itisstraightforwardtogeneratenewmeshes a
p
as the ice sheet evolves, and to transfer the previous time-step’s ice thickness data e Century-scale
r
to the new mesh in a conservative fashion. It is also relatively easy to study conver- simulations of the
|
gence with mesh resolution by running the same experiment for successive values of West Antarctic Ice
L and check that the differences between, say, the volume above flotation calculated D Sheet
is
c
10 in each case converge at the expected rate. We include the results of such a study in u
s S.L.Cornfordetal.
Sect. 3.1.1. s
io
n
P
2.4 Model data requirements a
p
e TitlePage
r
Time-dependent simulations require initial ice thickness data h (x,y)as well as accu-
0 | Abstract Introduction
mulationratesa(x,y,t)andmeltratesM(x,y,t)forEq.(2),togetherwithabedrockele-
vationmapb(x,y),abasalfrictioncoefficientfieldC(x,y),atemperaturefieldT(x,y,z) D Conclusions References
15 is
and a stiffening factor φ(x,y) to solve Eq. (3). Bedrock elevation and initial ice thick- c
us Tables Figures
nessdatafortheRISFRISandMBLdomainsweretakenfromtheALBMAP5kmDEM s
io
(Le Brocq et al., 2010). A custom map of bedrock elevation and ice thickness set on n
(cid:74) (cid:73)
P
a 1km grid was used for the ASE domain: it is close to the more recent Bedmap2 a
p
(Fretwell et al., 2013) data, and was used before for studies of Pine Island Glacier e (cid:74) (cid:73)
20 r
(Favier et al., 2014). It was prepared in a similar manner to ALBMAP, but includes
| Back Close
extra data from high resolution airborne radar (Vaughan et al., 2006) and submarine
D
surveys (Jenkins et al., 2010). It also includes a pinning point at the tip of Thwaites is FullScreen/Esc
c
Glacier’s slower flowing eastern ice shelf, a feature that is clearly visible in the velocity u
s
25 data(Joughinetal.,2009;Rignotetal.,2011),thatcorrespondstopeakoneofthetwo sio Printer-friendlyVersion
describedinTintoandBell(2011),butisabsentinthebathymetrydata.Weraisedthe n
P InteractiveDiscussion
bathymetry by 120m to ground the ice in that region. Ice temperature data is provided a
p
e
r
1895 |
by a three-dimensional thermo-mechanical model (Pattyn, 2010) and is held fixed in D
is
time. c TCD
u
The remaining data are rather more complicated. The basal friction and stiffening ss
coefficients inside the drainage basin are estimated by solving an inverse problem, ion 9,1887–1942,2015
P
5 describedinSect.2.4.1,whiletheaccumulationandmeltratesarecomposedofinitial a
p
accumulation and melt rates a and M , described in Sect. 2.4.2, and future climate e Century-scale
0 0 r
anomalies, described in Sect. 2.5. simulations of the
|
West Antarctic Ice
2.4.1 Basal friction and stiffening coefficients D Sheet
is
c
u
We estimate basal friction and stiffening coefficient by solving an inverse problem sim- s S.L.Cornfordetal.
s
ilar to those of MacAyeal (1993), Joughin et al. (2009) and Morlighem et al. (2010), io
10 n
Broadly speaking, we choose smooth fields C(x,y) and φ(x,y) that minimise the mis- P
a
p
matchbetweenmodeledandobservedspeeds.Anonlinearconjugategradientmethod e TitlePage
r
was employed to seek a minimum of the objective function
| Abstract Introduction
J =Jm+Jp (12) D Conclusions References
is
c
15 composed from a misfit function us Tables Figures
s
(cid:90) io
J = 1 α2(x,y)(cid:0)|u|−|u |(cid:1)2dΩ (13) nP (cid:74) (cid:73)
m 2 u obs a
p
Ω e (cid:74) (cid:73)
V r
and a Tikhonov penalty function | Back Close
D
α2 (cid:90) α2 (cid:90) is FullScreen/Esc
J = C |∇C|2dΩ+ φ |∇φ|2dΩ. (14) cu
p 2 2 s
Ω Ω sio Printer-friendlyVersion
V V n
P InteractiveDiscussion
The coefficient αu2(x,y) is related to error estimates for the observed velocity and we ap
e
set it to 1 where velocity data is available and 0 elsewhere. r
20
1896 |