Contents lists available atScienceDirect
Earth
and
Planetary
Science
Letters
www.elsevier.com/locate/epsl
Effect
of
GIA
models
with
3D
composite
mantle
viscosity
on
GRACE
mass
balance
estimates
for
Antarctica
Wouter van
der
Wal
a,
∗
,
Pippa
L. Whitehouse
b,
Ernst
J.O. Schrama
aaAerospace Engineering, Delft University of Technology, Kluyverweg 1, 2629 HS, Delft, Netherlands
bDepartment of Geography, Durham University, Lower Mountjoy, South Road, Durham, DH1 3LE, United Kingdom
a
r
t
i
c
l
e
i
n
f
o
a
b
s
t
r
a
c
t
Article history:
Received5September2014
Receivedinrevisedform31December2014 Accepted1January2015
Availableonline28January2015 Editor:Y.Ricard Keywords: glacialrebound mantlerheology viscosity time-variablegravity GRACE Antarctica
Seismic data indicate that there are large viscosity variations in the mantle beneath Antarctica. ConsiderationofsuchvariationswouldaffectpredictionsofmodelsofGlacialIsostaticAdjustment(GIA), whichareusedtocorrectsatellitemeasurementsoficemasschange.However,mostGIAmodelsused forthatpurpose haveassumedthemantletobeuniformlystratifiedintermsofviscosity.Thegoalof thisstudyistoestimatetheeffectoflateralvariationsinviscosityonAntarcticmassbalanceestimates derivedfromtheGravityRecoveryandClimateExperiment(GRACE)data.Tothisend,recently-developed globalGIAmodelsbasedonlateralvariationsinmantletemperaturearetunedtofitconstraintsinthe northernhemisphereandthencomparedtoGPS-derivedupliftratesinAntarctica.
We findthatthesemodelscan provideabetterfitto GPSupliftrates inAntarcticathanexistingGIA models witharadially-varying(1D)rheology.When3Dviscosity modelsincombinationwithspecific ice loadinghistoriesare usedtocorrectGRACEmeasurements,masslossinAntarcticaissmallerthan previouslyfoundforthesameiceloadinghistoriesandtheirpreferred1Dviscosityprofiles.Thevariation in mass balanceestimatesarising fromusing differentplausible realizationsof3D viscosity amounts to 20 Gt/yr forthe ICE-5G ice model and 16 Gt/yr forthe W12a ice model; thesevalues are larger thanthe GRACEmeasurementerror,butsmaller thanthe variationarisingfromunknown icehistory. Whilethereexist1DEarthmodelsthatcanreproducethetotalmassbalanceestimatesderivedusing3D Earthmodels,thespatialpatternofgravityratescanbesignificantlyaffectedby3Dviscosityinaway thatcannotbereproducedbyGIAmodelswith1Dviscosity.Asanexample,modelswith1Dviscosity alwayspredict maximum gravityrates intheRoss Seaforthe ICE-5G icemodel,however,for oneof thethreepreferred 3Dmodelsthemaximum(forthesameicemodel)isfoundneartheWeddell Sea. Thisdemonstratesthat3Dvariationsinviscosityaffectthesensitivityofpresent-dayupliftandgravity ratestochangesinthetimingoftheicehistory.Inparticular,lowviscosities(<1019Pa s)foundinWest Antarcticamakethemantleverysensitivetorecentchangesinicethickness.
©2015TheAuthors.PublishedbyElsevierB.V.ThisisanopenaccessarticleundertheCCBY-NC-ND license(http://creativecommons.org/licenses/by-nc-nd/4.0/).
1. Introduction
Measurementsof time-variable gravity from the GRACE satel-lite mission show continuous decrease of mass over Antarctica since the GRACE launch in 2002 (Velicogna and Wahr, 2006; King et al., 2012). A large part of the gravity change reflects mass redistribution in the solid Earthas the viscous mantle re-sponds to past changes in ice load, a process known as Glacial IsostaticAdjustment(GIA).Inordertodeterminepresenticemass change in Antarctica GRACE measurements have to be corrected forGIA,either:(i)byemployingageophysicalmodelforGIA(e.g.,
*
Correspondingauthor.Tel.:+31152782086.E-mail address:w.vanderwal@tudelft.nl(W. vanderWal).
Velicogna and Wahr, 2006; Chen et al., 2009; King et al., 2012), or (ii) by employing other datasets withdifferentsensitivities to GIA and icemelt such as GPS or satellite altimetry (Wahr et al., 1995, 2000; Riva et al., 2009; Wu et al., 2010; Wang et al., 2013; Sasgen et al., 2013; Gunter et al., 2014).Method (ii)hasthe advan-tagethatitisnotnecessarytorelyongeophysicalmodelsofGIA. However, itrequiresaccurate knowledgeoffirncompactionto be able to relate volume changes measured by satellite altimetry to masschangesmeasuredbyGRACE.Also,satellitealtimetry, specifi-callyICESat,suffersfrominhomogeneoustemporalandspatial cov-erage, cloud cover, detector saturation andinter-campaign biases (Shuman et al., 2006; Shepherd et al., 2012). Inversion of space-geodetic data(Wu et al., 2010; Sasgen et al., 2013) is sensitiveto data distribution, andspurious signals can be generatedin areas withfewer data. Finally,method (ii)doesnot makeuse ofmany http://dx.doi.org/10.1016/j.epsl.2015.01.001
key constraintson the GIAprocess, such as historic sea-level in-dicators,geomorphological,geologicalandglaciologicalconstraints ontheshapeandthicknessoftheicesheet,andknowledgeofthe interioroftheEarthbelowAntarctica.
Consideringtheadvantagesanddisadvantages ofeachmethod, thereismerit inpursuingbothmethodsinparallelforestimating present-dayice massbalance. Inthis studywe selectmethod (i) andfocus on the unknown structure of the Earth and how this affectspredictedgravitychangesduetoGIA.
MostGRACEmassbalanceestimatesforAntarcticathatrelyon method(i) haveassumed that deformation inthe Earth’smantle canbeparameterizedusingaviscositydistributionthatonlyvaries with depth (e.g., Velicogna and Wahr, 2006; Chen et al., 2009; King et al., 2012). In the following,this parameterizationwill be referred to as 1D viscosity. However, from surface wave data it isclear that the mantle is very different beneathEast andWest Antarctica(Ritzwoller et al., 2001; Danesi and Morelli, 2001). As-sumingthatdifferencesinseismicvelocitiesstemfromdifferences in mantle temperature, Kaufmann et al. (2005) suggest that a largedifferenceinmantleviscosity existsbetweenEastandWest Antarctica.Suchdifferencescouldhavealargeeffecton(regional) massbalanceestimatesforAntarcticaandtheinterpretationofGPS andaltimetrydatabecausepresent-dayupliftratesaresensitiveto thelocalviscousrelaxationtime.Belowwediscussprevious stud-iesandopenquestions.
Kaufmann et al. (2005) showed that the inclusion of 3D vis-cosity within a GIA model results in an uplift rate pattern that is similar to a 1D viscosity model. A et al. (2013) have com-putedgravityratesforacompressibleGIAmodelwith3D viscos-ity and found that the effect on mass change estimates is only mildlydifferent comparedwitha 1D model.Differencesbetween 1D and3D modelstherefore seemto be smallerthan the uncer-tainty resulting from poor knowledge of the ice loading history. However, Kaufmann et al. (2005) and A et al. (2013) only con-sidered one 3D viscosity distribution for Antarctica, while there areactuallymanyunknownsassociatedwithderiving3Dviscosity variationsfromseismicinformation(e.g.Ivins and Sammis, 1995; Trampert and van der Hilst, 2005).Inordertofullyinvestigatethe effectof3DviscosityinGIAmodels,theuncertaintyinproducing 3D viscosity mapsshould be considered, includingthe effects of differentflowlaws.Inourapproachwetakeintoaccountthetwo maintypesofdeformationinthemantle(diffusionanddislocation creep)in a so-calledcomposite rheology (Gasperini et al., 1992; van der Wal et al., 2010). Due to the difficulty of modeling the gradient in Earth structure that exists between East and West Antarctica,regional GIAmodels, which adopta 1D viscosity pro-file,havebeenusedtoinvestigatetheGIAsignalinspecificregions ofAntarctica(e.g.,Ivins et al., 2011;Nield et al., 2012, 2014).Such modelsarelikelytocontinuetobeusedbecausetheycanachieve thenecessaryspatial resolution forstudyingthe Earth’sresponse tochangesiniceloadingonaregionalscale.Therefore,oneofour aimsisalsotoproducearangeofviscositymapsthatcanbeused inregionalGIAstudiesthatadopta1Dviscosityprofile.
ClearlyGIAmodelswith3Dviscosityaremorecomputationally expensivethan GIA modelswith 1D viscosity, andthereforeit is importanttodeterminewhetheritisnecessarytouse3D models tocorrectGRACE mass balanceestimates, orwhether1D models aresufficient (see Section 3.5). Maybethe rangeofmass balance estimatesproducedusingasuitablywiderangeof1DGIAmodels containsthemassbalanceestimatethatwouldbeproducedusing a 3D Earthmodel, ormaybe there are importantregional differ-encesbetween themass change predictedby a 1D modeland a 3Dmodel.
Insummarytheresearchquestionstobeansweredinthisstudy are:
1. WhatistheeffectofusingGIAmodels withdifferent3D vis-cosity distributions on Antarctic mass balance estimates de-rivedfromGRACE?
2. Whatrangeoflateral(effective)variationsinviscositycanbe expectedbeneathAntarctica?
3. Can thegravityratepatternfromGIAmodelswith3D rheol-ogybereproducedbyaGIAmodelwith1Drheology? Inthis studythefree-air gravityanomaly rateis computedat the Earth’s surface, and this will be referred to in the following assimplythegravityrate.Section2describesthemostimportant featuresofthenumericalGIAmodelandiceloadinghistoriesused. Thissection isfollowedby thepresentationofviscositymapsfor the preferred3D modelsand a comparisonof modelpredictions with GPS data in Antarctica. Finally,mass balance estimates and acomparisonofpredictions from1Dand3DGIAmodelsare pre-sented.
2. Methods
2.1. Finite-elementmodel
TheGIAmodelisbasedonthecommercialfinite-element soft-ware ABAQUSTM, following Wu (2004). Elements have a 2◦
×
2◦ resolution at the surface and, as described in that paper, self-consistent sea levels and self-gravitation are included, but not compressibility,geocentermotionandshorelinemigration.Viscous parameters,asdescribedbelow,aredefinedforlayerswith bound-aries at 35, 70, 120, 230, 400, 670, 1170 and 3480 km depth. Elasticparametersareasinvan der Wal et al. (2013)with bound-ariestakenatthemajorseismicdiscontinuities at400,670,1170 and3480kmdepthaswellasat120kmdepth.Densityand rigid-ity for each layer are obtained by volume-averaging layers from thePREMmodel(Dziewonski and Anderson, 1981) withsmall ad-justments in order to better match densityjumps. The model is extended to include the two main types of deformation in the mantle:diffusion creepanddislocation creep(van der Wal et al., 2010).Here we useacomposite rheology(Gasperini et al., 1992; van der Wal et al., 2010) basedontheflowlawsfordiffusionand dislocation creep in olivine. We assume that olivine is the main mantle material and consider variations in grain size and water content.Byvaryingtheseparametersweintroducelargevariations intheviscosities thatare derivedfromthermalanomalies.In ad-dition,becausestrainratefordislocationcreepdependsonstress, effectiveviscosityvarieswithstressandhencewithtime. Includ-ing dislocationcreep resultsinsmallpresent-dayuplift-rates,but this ispartly countered by using a combinationof diffusion and dislocationcreep(van der Wal et al., 2010).Individualstrain componentsarecalculatedas(van der Wal et al., 2013):
ε
=
Bdiffqt
+
Bdislqnt
,
(1)whereBdiff andBdislarecreepparameterscomputedfromtheflow lawfordiffusionanddislocationcreep,respectively,t istime,n is
the stress exponent,and q is the von Mises stress q
=
3 2
σ
i jσ
i j withσ
i j anelementofthedeviatoricstresstensor.Above400km, where olivine isthe main mantlematerial, the olivine flow laws from Hirth and Kohlstedt (2003) are used to compute Bdiff and Bdisl:B
=
Ad−pf H2Oreαϕe−E+P V
R T
,
(2)inwhich A and
α
areconstants,d isthegrainsize, f H2O iswatercontent,
ϕ
ismeltfraction, E isactivationenergy,P ispressure,Visactivationvolume, R isthegasconstant,T isabsolute tempera-ture,andp andr arethegrainsizeandwaterfugacityexponents,
respectively.Ofthese, E, V , p,r and A aretakenfromHirth and Kohlstedt (2003)foreitherdiffusionordislocationcreep.The pres-sure as a function of depth is calculated by assuming that the pressure gradient is equal to 0.033 GPa/km (Keary et al., 2009). Grainsize, watercontentandtemperatureare unknownandwill be varied asdescribed later, while melt contentis setto zero.It hasbeenshownthatinahigh-temperatureregionsuchasIceland, meltcontentasmodeled byEq.(2)hasarelativelysmallinfluence on effective viscosity compared to grain size and water content (Barnhoorn et al., 2011a). Effective viscositycan becalculated by (van der Wal et al., 2013):
η
eff=
1 3Bdiff
+
3Bdislqn−1.
(3)Thetop35kmoftheEarthareassumedtobenon-viscous.Below that,theeffectiveviscositydetermines whetheran elementis re-spondingviscouslyornot,andhencewhetheritcanbeconsidered partofthelithosphere.FortheEarthlayersbelow400km,values for Bdiff and Bdisl are assumed to vary onlyradially becausethe olivineflowlawsofHirth and Kohlstedt (2003)donotholdinthis region.Bdiff andBdisl valuesforthesedepthsaretakenfroma3D GIAmodelthathasbeentunedto fitarangeofrelative sea-level data(van der Wal et al., 2010).
Temperatureisderivedintwodifferentways,fromsurfaceheat flow data(labeled HF) andfroma global seismicmodel (labeled SEIS).Themostimportantstepsinobtainingthetemperaturemaps are described in the following, butmore information is given in
van der Wal et al. (2013). In the HFapproach surface heat flow mapsareusedfromShapiro and Ritzwoller (2004).They extrapo-latedheatflowdatafromPollack et al. (1993)toareaswherenone was available. The extrapolationis basedon a shearwave veloc-itymodelandassumesathermallyhomogeneouscrust.Geotherms arecomputedbyintegratingtheequationfor1Dsteady-stateheat transfer,assumingconstantheatgeneration.
Because very few heat flow measurements exist for Antarc-tica, standard deviations in inferred heatflow there, as derived by Shapiro and Ritzwoller (2004), are large.Instead ofusing this standarddeviationtodetermineuncertainties associatedwiththe temperaturedistribution,weuseasecond,independent,approach toobtain an estimate oftemperature. Inthe SEISapproach later-allyvaryingvelocityanomaliesfromGrand (2002)areconvertedto temperatureusingthedepth-dependenttemperaturederivative of seismicwavevelocitiesgiveninKarato (2008).Here,itisassumed thatall seismicanomalies areduetothermalanomalies,whilein realitychemicalheterogeneityhasaninfluence.Intheupper man-tle the effect ofchemical heterogeneity is probablysmaller than theeffectofthermalanomalies(Cammarano et al., 2011),but nev-erthelessit caninfluence GIApredictions (Wu et al., 2013). Thus SEIS temperature estimates are an upper bound for lateral vari-ations in temperatures. Large differences exist between different tomographymodels (Schaeffer and Lebedev, in press), and hence will lead to differencesin thermal maps. Our use of two differ-ent methods to obtain thermal anomalies captures some of the variationarising fromuncertainty inthe approaches,but the un-certainty arisingfromusing differentseismictomographymodels is an important target for future work. Another interesting ap-proach is that of Priestley and McKenzie (2013), who estimate viscositydirectly fromshearwave velocity models and geophysi-calandpetrologicaldata.
We found an error in the calculation of the SEIS tempera-ture model in van der Wal et al. (2013), which resulted in the temperatures being too high at shallow depths and too low for deeperlayers. Theeffecton sea-levelcurvesissmallforthebest fitting models, butuplift rates were affected more therefore the recalculatedratesare showninSection 3.1.The twomethods for
computingtemperatures resultinmarkedlydifferenttemperature distributionswithSEIShavinglowertemperaturesthanHF.
ThefinaltwoparametersthatmustbedefinedwithintheEarth model are grain size and watercontent. Grain size is varied be-tween 1,4 and10 mm,which istherange foundforkimberlites andperidotites(Dijkstra et al., 2002).Watercontentisvaried be-tween afullywet(1000 ppmH2O)anda fullydrystate. Varying
themantletemperature(SEISandHF),grainsize(1/4/10 mm)and watercontent (wet/dry)results inatotal of2
×
3×
2=
12 com-binations of mantleparameters that are investigatedforeach ice loadinghistory(seethenext section).Effective viscosities are calculated for each ice–Earth model combination, andwill vary over both space andtime. Variations inspacearecausedbyspatialvariationsintemperatureandstress. Variationsintimeareduetothenon-linearpartofthecomposite rheologyflowlaw(secondtermofEq.(1)),whichhasbeenshown to affectviscosity by two orders ofmagnitude during theglacial cycle, neglectingtheinfluenceofbackgroundstress(Barnhoorn et al., 2011b).
2.2. Icemodels
We use two different ice loading histories for Antarctica: ICE-5Gv1.2(Peltier, 2004;referredtoasICE-5G)andW12a( White-house et al., 2012a, 2012b), both of which have previously been usedtocorrectGRACEmeasurementsforAntarcticGIAeffects.The Antarctic componentof themore recentICE-6G model was pub-lished during preparation ofthis manuscript (Argus et al., 2014) but the ice-loading history associated with this model was not available forcomparison withearliermodelsatthe time. We as-sumethat ICE-5GandW12aspanreasonable possibleiceloading variations (see alsoIvins et al., 2013,Fig. 2) and donot use the
Ivins et al. (2013) ice-loading history in order tolimit our com-putational effort. An importantlimitation ofall threeice-loading historiesisthattheyaretunedtofittosea-levelorupliftdata as-sumingalaterallyhomogeneousEarth.Infuturework3Dviscosity shouldbeconsideredwhendevelopingice-loadinghistories.
TheICE-5GandW12a modelshavebeeninterpolatedontothe 2◦
×
2◦equiangulargridofthefinite-elementmodel,andice thick-ness changes are defined at1000 yr intervalsbetween 20 kaBP andthepresent.Priorto20kaBPicethicknessisassumedto in-crease linearly over 90 ka. The main differences between W12a andICE-5Gare:– W12a incorporates a larger number of palaeo ice thickness constraints,derivedfromexposureagedating,whichwerenot availablewhenICE-5Gwasdeveloped
– W12awasdevelopedusinganumericalice-sheetmodelwhile ICE-5Gwasdirectlytunedtofitfieldobservations
– W12a makes use ofnear-field relative sea-level datato fine-tunethemodelwhereasICE-5Gistunedusingaglobalrelative sea-leveldataset
Asaresultofthesedifferences,thetotalmeltwatercontribution fromAntarctica sincetheLastGlacial Maximum(LGM)issmaller in the W12a model than the ICE-5G model.Both models infact define globalicethicknesschangesthroughoutthelastglacial cy-cle,butoutsideAntarcticaW12aisidenticaltotheICE-5Gloading history.
These two ice models are used to solve the sea-level equa-tion (Farrell and Clark, 1976) andhencedeterminegravitationally self-consistentglobalvariationsinrelativesealevelandEarth de-formationthroughoutthelast glacial cycleat1000yr timesteps. ForICE-5G,present-dayupliftratesareobtainedbynumerical dif-ferencingofthepredictedsolidearthdisplacement1000yrbefore andafterthepresent.Hencetheratesarecentered onpresent.For
Fig. 1. MaximumupliftratefordifferentmodelsinNorthAmerica(a)andScandinavia(b).Thegreybarindicatesthemaximumobservedupliftratewithonestandard deviationinNorthAmericaandScandinavia,accordingto
Sella
etal. (2007)andLidberg
etal. (2007),respectively.For(a)theICE-5Gicemodelisusedandfor(b)anice modelisusedthatwasdevelopedindependentlyfromGIAobservationsandmantleviscosity,seevan
derWaletal. (2013).the W12a model the derivative is calculated over a different in-tervalbecausethereareicethicknesschangesupto500yrbefore presentwhichwouldresultinlargeelasticeffectscontributingto theupliftrateifitwas calculatedinthesamewayasforICE-5G. Theuplift rateforthe W12amodel isthereforecalculatedasthe differenceindisplacementbetweenpresentand100yrinthe fu-ture.Thisonlyrequirestheadditionofoneextratimestepinthe computation and is found to be sufficiently accurate. The differ-encebetweenratescenteredat500andat50yrinthefuturein termsofupliftrateisatmost1.1 mm/yrinareasofmaximum up-liftrateandlessthana fewtenthsofmm/yroutsidethoseareas. Itfollowsthatratescenteredat50yrinthefuturewilldifferfrom ratescenteredatpresentbymuchlessthanthisamount.
3. Resultsanddiscussion
Thepreferred3D GIAmodels,basedon comparisonwith con-straintsonnorthernhemisphereGIA,arepresentedinSection3.1. Section 3.2presents mapsof effective viscosityfor thepreferred 3D GIA models. Section 3.3 compares uplift rates from the 3D modelswithGPS-measuredupliftratesinAntarctica.Theeffecton GRACEmassbalanceestimatesisdiscussed inSection 3.4.Finally, comparisons betweengravity rates from 1D and 3D GIA models aremadeinSection3.5.
3.1.Preferred3DGIAmodels
The123Dmodelspredictverydifferentupliftratesandgravity rates.Inordertoselectwhichofthemodelscanresultinrealistic upliftrates,we comparemodeloutput toobservations inregions whereGIA upliftrates area clearly observed,i.e.Scandinavia and NorthAmerica(Fig. 1).Intheabsenceofotherinformationonflow lawparameters suitable for Antarctica, we assume that the flow lawparameters thatresultinrealistic upliftratesinthenorthern hemispherealsoresultinreasonable upliftratesinAntarctica.For thiscomparison the ICE-5G model is used; the W12a icemodel wouldgive nearly identical results because its loading history is identicalto ICE-5Gin thenorthern hemisphere andthe effectof usingadifferentAntarcticmelthistoryonupliftratesinthe north-ern hemisphere is negligible. All the models predict uplift rates thataretoolowinScandinaviaandNorthAmericawhenICE-5Gis used,a knownresultformodelswithnon-linearrheology and,to alesserextent,compositerheology(van der Wal et al., 2010).
Thismisfit decreases in Scandinavia when the loadhistory is derived from a paleo ice height model that is not based on an earthmodelwithMaxwellrheologyand1Dviscosity(van der Wal et al., 2013),however,themaximumobservedupliftratesarestill not reproduced in this case (Fig. 1b). It appears from the figure
that under certain conditionsan increase in grain size could in-creasethe predictedupliftratebutthis was notfound to be the case in previous work (van der Wal et al., 2013). Therefore, for thisstudywe simplyselectthemodelsthat yieldthehighest up-liftrateseventhoughtheyarestillsomewhatbelowthemeasured maximumupliftrate. ForNorthAmericaandScandinaviathebest model is a dryrheology with 10 mm grain size combined with temperature model HF (labeled HF10D). If the SEIS temperature model is accepted, a model with dry rheology and 4 mm grain size(labeled S4D)givesthelargestupliftrates.
The model predictions were also compared with relative sea-leveldatainFennoscandia(van der Wal et al., 2013).Inthatcase thebest modelisbased onthe HFtemperaturemodelin combi-nation with a wet rheology and 10 mm grain size. Because the HF10W modelgivesa very poorfitto upliftrates inScandinavia (Fig. 1b)weinsteadadopttheS10Wmodelasourthirdpreferred model.This isa reasonable trade-off sincethe S10Wmodelonly gives a slightly worse fit to the relative sea-level data than the HF10W model.The three modelsHF10D, S4D,andS10W willbe used toinvestigatetherange inpredictionsone can getfrom3D models.
3.2. Effectiveviscositymaps
MapsofeffectiveviscosityinAntarcticaareplottedforthefirst preferredmodel(HF10D) inFig. 2usingW12a attime 14ka be-forepresent.RecallthataccordingtoEq.(3)effectiveviscosityisa function ofthe(vonMises) stress whichalsodependsonthe ice model.Aviscositymap foricemodelICE-5G,aswell asfora dif-ferentepoch (present) isprovidedin supplementarymaterial A.1. These additional viscositymaps are similar to Fig. 2, but viscos-itycanbeuptotwoordersofmagnitudesmallerduetolargerice thicknesschanges prescribed by the ICE-5G icemodeland up to oneorderofmagnitudelargerduetothelowerstressesatpresent comparedto14kabeforepresent.Morediscussionisprovidedin thesupplementary materialA.1andan extensiveanalysisof tem-poralchangesinviscosityispresentedinBarnhoorn et al. (2011b). In Fig. 2 it can be seen that, at a depth of 52 km, viscosity islow(<1020Pa s)inWestAntarcticaandthenorthern Antarctic Peninsula,whileEastAntarcticahashighviscosity(>1024Pa s).At 95km depthfor modelHF10D, viscosityin mostofEast Antarc-ticaistoohighforanyviscousdeformation,withtheexceptionof the coastalparts ofDronningMaudLand. At a depthof 145 km, viscosities in West Antarctica increase to around 1021 Pa s while
viscosities around much of coastal East Antarctica approach the samevalue.Finally,at200kmdepthviscositiesarelowenoughfor viscousdeformationtooccurbeneaththewholeofEastAntarctica. Lateralvariations atthisdeptharesmallbecausethetemperature
Fig. 2. Effectiveviscosityat4depthsforpreferredmodelHF10D,asderivedusing theW12aiceloadinghistory.DMdenotesDronningMaudLand.
Fig. 3. Effectiveviscosityat4depthsforpreferredmodelS10W,asderivedusing theW12aiceloadinghistory.
isdetermined moreby themantleadiabatthanthe value of sur-faceheatflow.
Fig. 3 shows the viscosity for model S10W, which provides
the second best-fit to historic sea levels in Scandinavia and sec-ond best fit to uplift rates in Antarctica (next section). For this model temperatures are derived from seismic velocity anoma-lies and the flow laws are those for wet olivine (see details in
Fig. 4. ViscosityrangesbelowAntarcticaasafunctionofdepth,averagedoverall timesteps,formodelsbasedontheW12aiceloadinghistory.
van der Wal et al., 2013).Because oflowtemperaturesatshallow depths (see Section 2.1), effective viscosities down to 95 km are large enough that no viscous deformation will occur across the wholeofAntarctica,exceptinthewesternRossSea,inthismodel. At 145 km depth the viscosity beneath the Antarctica Peninsula andbeneaththeRossIceShelfislowerthan1018Pa s,which cor-respondstorelaxationtimesontheorderofdecades.Thisis some-what belowestimatesina recentstudyinthe northernAntarctic Peninsulawhichfoundthatviscositiesof1018Pa s arerequiredin
order to matchobserved upliftrates following the2002 breakup of the Larsen B Ice Shelf (Nield et al., 2014). However, we note thattransientcreepmaybeinoperationoverthesetimescales,as suggested by experimental data (Faul and Jackson, 2005). Such a process isnotconsidered hereorbyNield et al. (2014),although thestress-dependenceinEq.(1)makestheviscosityweakly time-dependentinthisstudy.At200kmdepth,viscosityinsomecoastal regions of East Antarctica, e.g. Dronning Maud Land, drops to 1021Pa s orbelowforthefirsttime.
Therelationshipsbetweenviscosityvariationsanddepthbelow Antarctica are shown in Fig. 4. Below 200 km the model based on heatflow showsa viscosityslightlyincreasing withdepth,its value closeto that of a two-layer approximation of VM2, where VM2 is the viscosityprofile that is used to constructthe ICE-5G ice loadinghistory(Peltier, 2004). Bycomparing theSEIS models it can be seen that the effect of larger grain size, which acts to increaseviscosity,canbemorethancompensatedbyhavingawet insteadofadryrheology.
Both the wet rheology of model S10W and the small grain size inmodel S4Dresultin verylow viscosities in theflow laws of Hirth and Kohlstedt (2003), and consequently too small up-lift rates in Fennoscandia and North America (Fig. 1). However, whencombinedwiththeICE-5Gice-loadinghistory,thesemodels canreproduceGPS-measuredupliftratesinAntarcticabetterthan models witha1D rheology,aswillbe showninthenext section. Whetherthisimprovedfitisbecausethe3D modelsbetter repre-sent the Earthstructurebeneath Antarctica,orwhether errorsin the iceand Earthmodels cancelout, can onlybe determined by independently testing the accuracy of the ice models such as in
Whitehouse et al. (2012a) andArgus et al. (2014).While3D rhe-ology can onlybe constrainedvia GIA-relatedobservations ifthe icemodelisconstrainedindependentlyusingiceextentdata, fur-therevidenceoflowviscositiescomesfromxenolithswitholivine samples withsmallgrain size (0.1-4 mm)andhydrousminerals; seethecompilationinsupplementarymaterialA.2.
3.3. ComparisonwithGPSupliftrates
The 3D modelpredictions are comparedwithGPS uplift rates fromArgus et al. (2014).Theelasticupliftcorrectioninthatstudy
Table 1
Misfitbetweenmodeled upliftratesandselectedGPSupliftratesfrom
Argus
etal. (2014)withelasticratecorrectionsfromThomas
etal. (2011).Ice model W12a ICE-5G
Earth model Whitehouse et al. (2012b) HF10D S4D S10W VM2Peltier (2004) HF10D S4D S10W
Misfit (Eq.(4)) 0.91 0.74 1.2 2.6 1.3 0.83 0.61 0.63
relieson a GIA model and only corrects forlong-wavelength ef-fects. Therefore here the modeled elastic uplift correction from
Thomas et al. (2011)isusedforallstations.Ratesonthenorthern Antarctic Peninsula are presumed to reflect mostly elastic uplift (Thomas et al., 2011) sotheyarenotconsideredinthisstudy, nei-ther are stations for which the time series is shorter than 5 yr. Atotalof23stations passthesecriteria. TheGPS uplift ratesare giveninITRF2008,whilethe originofthemodelreferenceframe istheinstantaneous centerofmassoftheEarth(CM).Adriftcan exist between ITRF2008 and the CM frame. Such a drift mani-festsitselfmostlyasabiasbetweenmodeledandmeasureduplift rates. For thisreason we only consider uplift rates relative to a specificsite(asinvan der Wal et al., 2011);inthiscase,thesite withthesmallestmovement(maximummodeled upliftrateof1.2 mm/yr,acrossallmodels)whichalsohappenstobethesitewith thelongesttimeseries(Mawson). Themodeled upliftrateatthis siteissubtractedfromallmodeledupliftrates.Thisprocedurealso largelyremoves theeffect ofrotational feedbackto GIA whichis presentinthemeasuredupliftratesbutabsentfromthemodeled upliftrates.
Misfits between modeled and observed uplift rates are com-putedaccordingtothefollowingdefinition:
χ
2=
1 N N i=1 oi−
piσ
i 2,
(4)where N is the numberof observations (22), oi are the elastic-corrected relative uplift rate observations, pi are the predicted relativeupliftratesfromthemodelsinterpolatedattheGPSsites and
σi
are the standard deviations from Argus et al. (2014), not including the error in the elastic correction. Misfits are listed in Table 1. In there, each ice model is combined with the three preferred 3D Earth models. For reference the table also shows results derived from the published uplift rates for each ice model, which were produced using a specific 1D earth model. For ICE-5G this is the VM2 viscosity profile; uplift rates for ICE-5Gv1.3 in combination with VM2 L90 are taken fromhttp://www.atmosp.physics.utoronto.ca/~peltier/data.php (last ac-cessedon August2014).FortheW12a modeltheoptimum Earth model is constrained by relative sea-level data, and uplift rates aretakenfromWhitehouse et al. (2012b).Thecombinationofthe W12a ice model and the earth model parameters as derived in
Whitehouse et al. (2012b)willbereferredtoastheWhitehouse et al. (2012b)model,withpredictionstakenfromthatpaper.Wenote that both thesepublished models were derived using a spectral GIA model, and therefore differenceswith the 3D model results maybe duein partto the useof a finiteelement model inthis study.
Toinvestigatethiswereproduced upliftratesfromthe White-house et al. (2012b) modelwith thefinite element model,using thesameelasticandviscous profileasWhitehouse et al. (2012b)
(supplementarymaterialA.4).Upliftratesfromtheoriginal White-house et al. (2012b) model result in a misfit of 0.91, while the finite-element reproduction of this model gives a slightly larger misfitof0.95,mainlyduetothesmoothingoftheiceloadinthe lower-resolutionfinite-element model.Thissuggeststhatthe two computational methods give comparable results. In addition, the factthat afinite-elementmodelwith3D Earthstructure (HF10D) is able to produce smaller misfits than the finite-element ver-sion of the 1D Whitehouse et al. (2012b) model indicates that
Fig. 5. Uplift ratemapsfor the W12aice loadinghistory.Left:3Dearth model HF10D,right:upliftratesfrom
Whitehouse
etal. (2012b).(Forinterpretationofthe referencestocolorinthisfigure,thereaderisreferredtothewebversionofthis article.)the improvementis mostlikely dueto the imposed 3D viscosity variations.
FortheICE-5Gmodel,all3Dmodelsresultinanimproved mis-fit.A histogramof thedifferencesin supplementary material A.3 showsthat thismainlyarisesduetothe reductionof previously-large uplift rates at a few sites. The 3D model that gave the second-highestupliftratesinthenorthernhemisphere(S4D)leads to the best fit for the ICE-5G model in Antarctica. It is possible that the rheology of model S4D better reflects the rheology in Antarctica compared to models HF10D and S10W, but the small differenceinmisfitisunlikelytobe significantinthepresenceof othermodelerrors.
Uplift rates for the case when W12a is combined with the HF10Dearthmodelareplottedtogetherwiththeupliftratesfrom
Whitehouse et al. (2012b) in Fig. 5. The colored dots depict the elastic-corrected GPS uplift ratesfrom Argus et al. (2014).Uplift rates forthe3D compositerheology model are smaller,asfound previouslyforScandinaviaandNorthAmerica(van der Wal et al., 2013).Smallerverticalmotioncanresultbothfromlowerthan av-erageviscosityleadingtofastrelaxation,asseenintheAmundsen SeaSector, aswell ashigherthan averageviscosity, asfound be-lowcentralEastAntarctica(seeFig. 2)wheresubsidenceratesare reduced.
Whilethe smallerupliftratesinthenorthern hemisphere un-derpredict observed uplift rates,the introduction of3D structure and composite rheology improves the fit to observed rates in Antarctica.OneexplanationisthefactthatinAntarcticacomposite rheologydoesnotreduce upliftratesasmuchasitdoesinNorth America. ThemaximumupliftrateformodelHF10D inAntarctica is larger than theVM2 uplift ratethere, while inNorth America andScandinaviatheHF10D upliftratesare belowtheVM2uplift rates(seevan der Wal et al., 2013,Fig.12).However,another ex-planation could be that larger variations in Earthstructure exist beneathAntarctica which increasesthe influenceof 3D rheology. Note that the uncertainty in ice models in Antarctica is larger, thereforeimprovementsinfitarelesssignificant.
3.4. Massbalanceestimates
ToobtainmassbalanceestimatesfromGRACE,we useRelease 5monthlygravityfieldsfromtheCenterforSpaceResearch
span-Table 2
MassbalanceestimatesforGIAmodelswithvaryingEarthmodelparametersanddifferenticemodels.Theerrorinthetrendderivedfromcalibratedstandarddeviationsof theGRACEmonthlygravityfieldsis9.2Gt/yrforAntarctica,3.2Gt/yrforWestAntarctica,1.9Gt/yrfortheAntarcticPeninsulaand6.7Gt/yrforEastAntarctica.Themodels inrows5to7are1DmodelsthatbestapproximatetheHF10D,S4DandS10Wmodels,respectively;theyarediscussedfurtherinSection3.5.
GIAEarthmodel Icemodel Icemasschange (Gt/yr)
All Antarctica West Antarctica Antarctic Peninsula East Antarctica
HF10D ICE-5G −166 −147 −36 17 S4D ICE-5G −154 −141 −35 23 S10W ICE-5G −146 −135 −34 23 VM2 ICE-5G −178 −153 −41 15 1D/HF10D ICE-5G −163 −145 −35 17 1D/S4D ICE-5G −160 −144 −34 19 1D/S10W ICE-5G −146 −136 −33 23 HF10D W12a −55 −122 −30 97 S4D W12a −48 −120 −25 97 S10W W12a −39 −117 −19 97
Whitehouse et al. (2012b) W12a −98 −146 −39 87
ningFebruary 2003 to June 2013. The procedure is described in
Schrama et al. (2014);abriefsummaryisprovidedinthe follow-ing. The C20 coefficient is replaced by the values from Satellite
LaserRanging (SLR) ranging (Cheng et al., 2013) and continental waterstorage changes are accountedfor usingthe GLDASmodel (Rodell et al., 2004).TheGRACEdataareinvertedforwater equiv-alentheightin10 242globallydistributedmascons.Errorsare de-terminedbypropagatingcalibratedstandarddeviationsofthesets ofmonthlycoefficientsintothetrendestimate.Loadingatdegree 1resultingingeocentermotionisnotdirectlyobservedbyGRACE butcanbe estimatedindirectlyby assuming that massloss from icesheetsandglaciersaddstotheoceans.Themethodisdescribed inSchrama et al. (2014)whereitisfoundthatAntarcticareceives acorrectionof33Gt/yr.
Themassbalance estimatesfromGRACE arecorrected forGIA toyieldestimatesoficemassbalanceinAntarctica.OurGIA mod-elsdonotincludetheeffectofachangeintherotationalpotential asaresultofachangeinthe momentofinertiaduetomass re-distribution(Milne and Mitrovica, 1998).Toestimatetheeffectwe replaced theStokes coefficientsfor degree 2order 1with values computedusingaGIAmodelthatdoesincluderotationalfeedback (Peltier et al., 2012) andfound a small increase inmass balance estimatesof3Gt/yr.
The massbalance derived using differentGIA models is sum-marized in Table 2. We use the three preferred models as an indicationofthepossiblespreadinmassbalanceestimates result-ing fromunknown 3D rheology.This isnot an uncertainty range intheformal sense,aswe didnotinvestigatevariationinall pa-rametersinthe3Drheology,andnostatisticaltestwithrespectto thedataisperformed.Moreover,theicemodelsneglectedcoupling withthe3DEarthrheology,sothetruespreadfrom3D rheology couldbelarger.Stillwethinkthatitisinsightfultoseetheimpact ofthreerealizationsof3Drheology thatcanprovideareasonable fitto GIAobservations onthe northernhemisphere aswell asin Antarctica.
ThemassbalanceestimatesinTable 2aredifferentfrom previ-ousestimatespartlyduetoanaccelerationinicemelt(Velicogna, 2009) andanomaloussnowfallin2009acrossDronningMaudLand (Boening et al., 2012).FortheICE-5Gmodel,theuseofthethree different 3D viscosity models results in a variation in mass bal-ance of 20 Gt/yr, with mostof the variation in West Antarctica, anda mean value ofmass lossthat is 23 Gt/yrsmaller thanfor ICE-5G/VM2. For W12a the mass balance estimates vary by 16 Gt/yr,andthemeanmasslossis51Gt/yrsmallerthanwhenthe (1D) Whitehouse et al. (2012b) model is used to correctGRACE. Some of thisdiscrepancymay be accountedfor by model differ-ences: whenthe finite-element modelis run using the 1D Earth model parameters of Whitehouse et al. (2012b) it underpredicts
the signalduetoGIAby 5Gt/yrforthewholeofAntarctica,and by 14 Gt/yrforWestAntarctica(seesupplementary materialA.4), however,evenafteraccountingforthis,thestepfrom1Dto3D re-sults is still significant.About half of thedifference between the 1D and3DmodelsiscomingfromWestAntarcticawherethe3D modelspredictsmallerupliftratesthanthe1DWhitehousemodel (see Fig. 5).ForEastAntarctica,greater ratesoficemassgain are estimatedwhen3Dmodelsareused.
From Table 2, the maximum difference between two models
that use the same earth model but different ice models is 111 Gt/yr. Thus, the range in mass balance estimates resulting from variations in 3D viscosity (across our three preferred models) is less than theuncertainty caused by variations in theice loading history(asalsofoundbyA et al., 2013),butlargerthanthe uncer-tainty intheGRACE-derivedtrends.Therangedueto3Dviscosity variations iscomparabletotheuncertaintyin1DGIAmodels de-rivedinKing et al. (2012)andIvins et al. (2013),18and13 Gt/yr respectively, butit issmallerthanthe110Gt/yrinBarletta et al. (2008).
3.5. Approximatinga3DGIAmodelwitha1DGIAmodel
Wewishtodeterminewhethergravityratepredictionsderived using a GIA model with 3D viscosity variations may be well-approximated by a GIA model with 1D viscosity variations, for the purpose of computing mass balance estimates for Antarctica asa whole andregionally. Toinvestigatethisrequires knowledge ofwhich1Dviscosityprofilecorrespondsbesttoacertain3D vis-cosity distribution.However, sensitivityof the gravity rateto 3D variations inviscositydependsonthesizeoftheloadandonthe viscosityitself,whicharenotwell known.Therefore,itisdifficult toaveragethe3D viscositystructureinawaythatrepresentsthe truesensitivityoftheloadingprocess.Toadd tothat,viscosityin ourmodelisalsoafunctionoftime.Therefore,weopttocompute gravity ratesforarangeof1D models thatadoptdifferentupper andlowermantleviscosityvalues.Wethencomputethemisfit be-tween these1D modelsanda3D model,andusethemodelwith smallestmisfitasthebest1Dapproximationtothe3Dmodel.For thistesttheICE-5Gicemodelisselectedastheloadinghistory.We expect similar conclusions forothericemodels, butthatremains to be investigated. The gravity ratesof 1D models are computed usingthespectral modelofvan der Wal et al. (2011).Differences inspatialresolutionbetweentheFEmodelandthespectralmodel lead to small differencesin predictions, asdemonstrated in sup-plementarymaterial A.4.Misfitiscomputedbetweenthe1D and 3Dmodelsforall2◦
×
2◦ gridcellswithinthelandareaof Antarc-tica,accountingforthereductioninareatowardsthepole.The3D modelsandthe1DmodelsthatbestapproximatethemareshownFig. 6. Free-airgravity anomalyrates forthe 3Dmodels fromSection 3.1(right column)andthe1Dmodelsthatbestapproximatethem(leftcolumn).Maximum sphericalharmonicdegreeusedinbothmodelsis 90.Upperandlowermantle vis-cosities(UM/LM)areasindicatedinthefiguretitles(times1020 Pa s),aswellas
lithospherethickness(Li).‘S’denotesthelocationofSipleDome.
inFig. 6andmassbalanceestimatesforthebest-fitting1D mod-els are added to Table 2 as rows 5 to 7. The estimates for the 1D models differ from the 3D model by at most6 Gt/yr, which iswithin therangeofpreviously-computeddifferencesbetweena finite-elementmodelandaspectral1Dmodel.
Thequestion remains asto whetherthere are regional differ-encesbetweena 3D modelandthe1D model that best approxi-matesit.Ina1Dmodelthegravityratepatterntendstoresemble the distribution oftotal ice thicknesschange, even though small differencesinthetimingofmelt mightexistfromonelocationto another.Fora3D modelthepatternoftotalicethicknesschange andthepatternofpresent-daygravityratescan beverydifferent fromeachother.
Studyingthe spatial patternin Fig. 6,it can be seen that the gravityratepattern of3D models HF10D andS4Dare quite well
approximated bya 1D model; themaximumgravity rateinboth the1Dand3DmodelsisfoundatSipleDomewheremostofthe ice thickness change since LGM took place according to ICE-5G. However, modelS10W predicts the maximum gravity rateto be intheWeddellSea,whileitscorresponding1Dmodelpredictsthe maximum gravity rate to be at Siple Dome. Indeed, we verified that for the ICE-5G ice model (Peltier, 2004) thelocation of the predictedmaximumupliftratefor1D modelsthatsamplealarge rangeofupper/lower mantleviscositycombinationsalways corre-spondstotheSipleDome.ThisalsoholdstruefortheIJ05model (compare Fig. 4 ofIvins and James, 2005totheir Fig. 2) andthe IJ05_R2 model (compare Figs. 4 and 5 of Ivins et al., 2013 with their Fig. 3d).Notealso that theHF10D modelchanges the loca-tionofthemaximumupliftratefortheWhitehouse et al. (2012b)
modelfromtheWeddellSeatoneartheRossSea(Fig. 5).
4. Conclusions
From a set of GIA models with 3D viscosity, three preferred models were selected which provided thebest fit to either GPS-observeduplift ratesinthe northernhemisphere,orrelative sea-leveldatainFennoscandia.Allthreemodelsincludeviscosity pro-files that vary by several orders of magnitude within the upper mantle.Themodelthatpredictsthebest-fittingupliftratesinthe northern hemisphere is based on flow laws for dry olivine with large grain size (10 mm). This results in viscosity values below 1019 Pa s for parts of West Antarctica at 95 km depth, increas-ingtoalmost1022Pa s at300 kmdepth.Viscositiesthatare even lower are obtained for an alternative model with a wet olivine rheology.Althoughusingmineralflowlawstocomputeviscosities is uncertain,therheological parameters forlow viscosities are in agreementwithxenolithfindingsinAntarcticawhichrevealgrain sizessmallerthan1mmandwhichshowthepresenceofhydrous mineralsinmantlerocks(SupplementarymaterialA.2).
Usingthe3D viscositymodelstocorrectGRACEdata(February 2003–June2013)forGIAeffectsresultsinAntarcticmassbalance estimates of
−
146 to−
166 Gt/yr for the ICE-5G icemodel and−
39 to−
55 Gt/yr fortheW12amodel.Thesevaluesareless neg-ativethanearlierestimatesbasedon1D modelsforthesameice loading histories. It is possible to finda 1D modelthat approxi-matesthegravityratesfromeachofthe3Dmodelsforthepurpose ofGRACEmassbalanceestimates.However,estimatesbasedon3D models are outside the confidenceintervals for earlierpublished mass balance estimates based on 1D GIA models. The variation aroundthemeanresultingfromtheintroductionofarangeof3D viscositymodels is10 and8 Gt/yrforICE-5GandW12a, respec-tively.Thereducedmeanicemelt estimatesaswellasthe varia-tionaroundthemeanindicatesthat 3Dviscositycansignificantly affectmassbalanceestimates.Inpractice,uncertainties associated with3Dviscosityarelikelytobe evengreaterifthetrade-off be-tweeniceloadingand3DEarthrheology weretakenintoaccount duringdevelopmentoftheice-loadinghistory.Forone3D modelwithwetrheology (andlowereffective vis-cosity)thepredictedspatialpatternofgravityrateswasmarkedly different to the patterns produced by the closest 1D model ap-proximation.For example,the location ofthe largestgravity rate forthis3Dmodel(basedonICE-5G)nolongercorrespondstothe location ofgreatesticethicknesschangesincethe LGM,asisthe case for1D models. This demonstrates that future mass balance studiesthatuseGRACEtodeterminethespatialdistributionofice mass changewill benefitfromthe useofmore realistic viscosity distributionswithinAntarcticGIAmodels.Italsoindicatesthatice loadinghistoriesthat havebeentuned tofitGIAobservations us-ing1Dviscosityprofilesmaybeinerror.
Acknowledgements
Wethankthetwo reviewersfortheircommentswhichhelped us to improve the manuscript. We gratefully acknowledge Mer-ijn Logtestijn for compiling the xenolith studies in Supplemen-tarymaterialA.2andMattKingforvaluablediscussions regarding theGPS data.WWisfunded bythe NetherlandsOrganisation for Scientific Research(NWOGrant Number863.11.015).PLWis sup-portedbyaNERCIndependentResearchFellowship(grantnumber NE/K009958/1).
Appendix A. Supplementarymaterial
Supplementarymaterialrelatedtothisarticlecanbefound on-lineathttp://dx.doi.org/10.1016/j.epsl.2015.01.001.
References
A,G.,Wahr,J.,Zhong,S.,2013.Computationsoftheviscoelasticresponseofa3-D compressibleEarthtosurfaceloading:anapplicationtoGlacialIsostatic Adjust-mentinAntarcticaandCanada.Geophys.J.Int. 192(2),557–572.
Argus,D.F.,Peltier,W.R.,Drummond,R.,Moore,A.W.,2014.TheAntarctica compo-nentofpostglacialreboundmodelICE-6G_C(VM5a)basedonGPSpositioning, exposureagedatingoficethicknesses,andrelativesealevelhistories.Geophys. J.Int. 198(1),537–563.
Barletta,V.R.,Sabadini,R.,Bordoni,A.,2008.IsolatingthePGRsignalintheGRACE data:impactonmassbalanceestimatesinAntarcticaandGreenland.Geophys. J.Int. 172(1),18–30.
Barnhoorn,A.,vanderWal, W.,Drury,M.R.,2011a.Upper mantleviscosityand lithosphericthicknessunderIceland.J.Geodyn. 52(3–4),260–270.
Barnhoorn,A.,vanderWal,W.,Vermeersen,L.L.A.,Drury,M.R.,2011b.Lateral, ra-dial,and temporalvariations in uppermantle viscosityand rheology under Scandinavia.Geochem.Geophys.Geosyst. 12,Q01007.
Boening, C., Lebsock, M., Landerer, F.,Stephens, G., 2012. Snowfall-drivenmass changeontheEastAntarcticicesheet.Geophys.Res.Lett. 39(21),L21501. Cammarano,F.,Tackley,P.,Boschi,L.,2011.Seismic,petrologicalandgeodynamical
constraintsonthermalandcompositionalstructureoftheuppermantle:global thermochemicalmodels.Geophys.J.Int. 187(3),1301–1318.
Chen,J.L.,Wilson,C.R.,Blankenship,D.,Tapley,B.D.,2009.AcceleratedAntarcticice lossfromsatellitegravitymeasurements.Nat.Geosci. 2(12),859–862.
Cheng,M.,Tapley,B.D.,Ries,J.C.,2013.DecelerationintheEarth’soblateness.J. Geo-phys.Res.,SolidEarth 118,740–747.
http://dx.doi.org/10.1002/jgrb.50058
.Danesi,S.,Morelli,A.,2001.StructureoftheuppermantleundertheAntarcticPlate fromsurfacewavetomography.Geophys.Res.Lett. 28(23),4395–4398.
Dijkstra,A.H.,Drury, M.R.,Frijhoff,R.M.,2002.Microstructuresandlatticefabrics inthe Hiltimantlesection (OmanOphiolite):evidencefor shearlocalization andmeltweakeninginthecrust–mantletransitionzone?J.Geophys.Res. 107 (2270).
http://dx.doi.org/10.1029/2001JB000458
.Dziewonski,A.M.,Anderson, D.L., 1981.PreliminaryreferenceEarthmodel.Phys. EarthPlanet.Inter. 25(4),297–356.
Farrell,W.E.,Clark,J.A.,1976.Onpostglacialsealevel.Geophys.J.R.Astron.Soc. 46, 647–667.
Faul,U.H.,Jackson,I.,2005.Theseismologicalsignatureoftemperatureandgrain sizevariationsintheuppermantle.EarthPlanet.Sci.Lett. 234(1–2),119–134. Gasperini, P., Yuen, D.A., Sabadini, R., 1992. Postglacial rebound with a
non-newtonianuppermantleandanewtonianlowermantlerheology.Geophys.Res. Lett. 19(16),1711–1714.
Grand,S.P.,2002.Mantleshear-wavetomographyandthefateofsubductedslabs. Philos.Trans.R.Soc.Lond.A 360,2475–2491.
Gunter,B.C.,Didova,O.,Riva,R.E.M.,Ligtenberg,S.R.M.,Lenaerts,J.T.M.,King,M.A., vanden Broeke, M.R., Urban, T., 2014. Empirical estimation of present-day Antarcticglacialisostaticadjustmentand icemasschange.Cryosphere 8 (2), 743–760.
Hirth,G.,Kohlstedt,D.,2003.RheologyoftheUppermantleandthemantlewedge: aviewfromtheexperimentalists.In:Eiler,J.(Ed.),InsidetheSubduction Fac-tory.In:GeophysicalMonographSeries,vol. 138.AmericanGeophysicalUnion, Washington,DC,pp. 83–105.
Ivins,E.R.,James,T.S., 2005.Antarcticglacialisostaticadjustment:anew assess-ment.Antarct.Sci. 17(4),537–549.
Ivins,E.R.,James,T.S., Wahr,J.,Schrama,E.J.O.,Landerer,F.W.,Simon,K.M.,2013. AntarcticcontributiontosealevelriseobservedbyGRACEwithimprovedGIA correction.J.Geophys.Res.,SolidEarth 118(6),3126–3141.
Ivins,E.R.,Sammis,C.G.,1995.Onlateralviscositycontrastinthemantleandthe rheologyoflow-frequencygeodynamics.Geophys.J.Int. 123(2),305–322. Ivins,E.R.,Watkins,M.M.,Yuan,D.-N.,Dietrich,R.,Casassa,G.,Rülke,A.,2011.
On-landicelossandglacialisostaticadjustmentattheDrakePassage:2003–2009. J.Geophys.Res.,SolidEarth 116(B2),B02403.
Karato,S.-i.,2008.DeformationofEarthMaterials.CambridgeUniversityPress, Cam-bridge,UK.463 pp.
Kaufmann,G.,Wu,P.,Ivins,E.R.,2005.Lateralviscosityvariationsbeneath Antarc-ticaandtheirimplicationsonregionalreboundmotionsandseismotectonics. J. Geodyn. 39(2),165–181.
Keary, P., Klepeis, K.A., Vine, F.J., 2009. Global Tectonics, 3rd edition. Wiley– Blackwell,WestSussex,UK.463 pp.
King, M.A., Bingham, R.J., Moore,P., Whitehouse, P.L., Bentley,M.J., Milne,G.A., 2012.Lowersatellite-gravimetryestimatesofAntarcticsea-levelcontribution. Nature 491(7425),586–589.
Lidberg, M.,Johansson,J.M.,Scherneck,H.-G., Davis,J.L., 2007.Animproved and extendedGPS-derived3Dvelocityfieldoftheglacialisostaticadjustment(GIA) inFennoscandia.J.Geod. 81,213–230.
Milne,G.A.,Mitrovica,J.X.,1998.Postglacialsea-levelchangeonarotatingEarth. Geophys.J.Int. 133,1–19.
Nield,G.A.,Barletta,V.R.,Bordoni,A.,King,M.A.,Whitehouse,P.L.,Clarke,P.J., Do-mack,E.,Scambos,T.A.,Berthier,E.,2014.Rapidbedrockupliftinthe Antarc-ticPeninsulaexplainedbyviscoelasticresponsetorecenticeunloading.Earth Planet.Sci.Lett. 397,32–41.
Nield,G.A.,Whitehouse,P.L.,King,M.A.,Clarke,P.J.,Bentley,M.J.,2012.Increased iceloadingintheAntarcticPeninsulasincethe1850sanditseffectonglacial isostaticadjustment.Geophys.Res.Lett. 39(17),L17504.
Peltier,W.R.,2004.Globalglacialisostasyandthesurfaceoftheice-ageEarth:the ICE-5G(VM2)modelandGRACE.Annu.Rev.EarthPlanet.Sci. 32,111–149. Peltier,W.R.,Drummond,R.,Roy,K.,2012.CommentonOceanmassfromGRACE
andglacialisostaticadjustmentbyD.P. Chambersetal.J.Geophys.Res.,Solid Earth 117(B11),B11403.
Pollack, H.N.,Hurter,S.J.,Johnson,J.R.,1993.Heatflow fromtheEarth’sinterior: analysisoftheglobaldataset.Rev.Geophys. 31,267–280.
Priestley,K.,McKenzie,D.,2013.Therelationshipbetweenshearwavevelocity, tem-perature,attenuation andviscosity inthe shallowpartofthe mantle. Earth Planet.Sci.Lett. 381,78–91.
Ritzwoller,M.H.,Shapiro,N.M.,Levshin,A.L.,Leahy,G.M.,2001.Crustalandupper mantlestructurebeneathAntarcticaandsurroundingoceans.J.Geophys.Res., SolidEarth 106(B12),30645–30670.
Riva,R.E.M.,Gunter,B.C.,Urban,T.J.,Vermeersen,L.L.A.,Lindenbergh,R.C.,Helsen, M.M.,Bamber,J.L.,vandeWal,R.S.W.,vandenBroeke,M.R.,Schutz,B.E.,2009. GlacialIsostaticAdjustmentoverAntarcticafromcombinedICESatandGRACE satellitedata.EarthPlanet.Sci.Lett. 288(3–4),516–523.
Rodell,M.,Houser,P.R.,Jambor,U.,Gottschalck,J.,Mitchell,K.,Meng,C.-J.,Arsenault, K.,Cosgrove,B.,Radakovich,J.,Bosilovich,M.,Entin,J.K.,Walker,J.P.,Lohmann, D.,Toll,D.,2004.Thegloballanddataassimilationsystem.Bull.Am.Meteorol. Soc. 85(3),381–394.
Sasgen,I.,Konrad,H.,Ivins,E.R.,Vanden Broeke,M.R.,Bamber,J.L., Martinec,Z., Klemann,V.,2013.Antarcticice-massbalance2003to2012:regional reanal-ysis ofGRACEsatellitegravimetrymeasurementswith improvedestimateof glacial-isostatic adjustmentbased onGPS upliftrates.TheCryosphere 7 (5), 1499–1512.
Schaeffer,A.J.,Lebedev,S.,inpress.Globalheterogeneityofthelithosphereand un-derlying mantle:aseismologicalappraisalbased onmultimodesurface-wave dispersionanalysis,shear-velocitytomography,andtectonicregionalization.In: TheEarth’sHeterogeneousMantle.Springer,
http://homepages.dias.ie/~aschaeff/
papers/Schaeffer+Lebedev_Springer_2013_small.pdf.
Schrama, E.J.O., Wouters, B., Rietbroek,R., 2014. A masconapproach to assess icesheetandglaciermassbalancesandtheiruncertaintiesfromGRACEdata. J. Geophys. Res., Solid Earth 119 (7), 6048–6066. http://dx.doi.org/10.1002/ 2013JB010923.
Sella,G.F.,Stein,S.,Dixon,T.H.,Craymer,M.,James,T.S.,Mazzotti,S.,Dokka,R.K., 2007.Observationofglacialisostaticadjustmentin“stable”NorthAmericawith GPS.Geophys.Res.Lett. 34,L02306.
http://dx.doi.org/10.1029/2006GL027081
.Shapiro,N.M.,Ritzwoller,M.H.,2004.Inferringsurfaceheatfluxdistributionsguided byaglobalseismicmodel:particularapplicationtoAntarctica.EarthPlanet.Sci. Lett. 223(1–2),213–224.
Shepherd, A.,et al.,2012. A reconciledestimateofice-sheetmassbalance. Sci-ence 338(6111),1183–1189.
Shuman,C.A.,Zwally,H.J.,Schutz,B.E.,Brenner,A.C.,DiMarzio,J.P.,Suchdeo,V.P., Fricker,H.A.,2006.ICESatAntarcticelevation data:preliminaryprecisionand accuracyassessment.Geophys.Res.Lett. 33(7),L07501.
Thomas,I.D.,King,M.A.,Bentley,M.J.,Whitehouse,P.L.,Penna,N.T.,Williams,S.D.P., Riva,R.E.M.,Lavallee,D.A.,Clarke,P.J.,King,E.C.,Hindmarsh,R.C.A.,Koivula,H., 2011.WidespreadlowratesofAntarcticglacialisostaticadjustmentrevealedby GPSobservations.Geophys.Res.Lett. 38(22),L22302.
Trampert,J.,vanderHilst,R.D.,2005.Towardsaquantitativeinterpretationofglobal seismictomography. In:vanderHilst, R.D.,et al. (Eds.),Earth’sDeep Man-tle:Structure,CompositionandEvolution.In:GeophysicalMonograph,vol. 160. AmericanGeophysicalUnion,pp. 47–62.
vanderWal,W.,Barnhoorn,A.,Stocchi,P.,Gradmann,S.,Wu,P.,Drury,M., Ver-meersen,B.,2013.Glacialisostaticadjustmentmodelwithcomposite3-DEarth rheologyforFennoscandia.Geophys.J.Int. 194(1),61–77.
vanderWal, W.,Kurtenbach,E.,Kusche, J.,Vermeersen, L.L.A.,2011. Radialand tangentialgravityratesfromGRACEinareasofglacialisostaticadjustment. Geo-phys.J.Int. 187,797–812.
vanderWal, W.,Wu,P.,Wang,H.,Sideris,M.G.,2010.Sealevelsandupliftrate fromcompositerheologyinglacialisostaticadjustmentmodeling.J.Geodyn. 50 (1),38–48.
Velicogna,I.,2009.IncreasingratesoficemasslossfromtheGreenlandand Antarc-ticicesheetsrevealedbyGRACE.Geophys.Res.Lett. 36(19),L19503. Velicogna,I.,Wahr,J.,2006.Measurementsoftime-variablegravityshowmassloss
inAntarctica.Science 311(5768),1754–1756.
Wahr,J., DaZhong, H., Trupin, A.,1995. Predictions of verticaluplift caused by changingpolaricevolumesonaviscoelasticEarth.Geophys.Res.Lett. 22(8), 977–980.
http://dx.doi.org/10.1029/94GL02840
.Wahr, J., Wingham, D., Bentley, C., 2000. A method of combining ICESat and GRACEsatellitedatatoconstrainAntarcticmassbalance.J.Geophys.Res.,Solid Earth 105(B7),16279–16294.
Wang,H.,Jia,L.,Steffen,H.,Wu,P.,Jiang,L.,Hsu,H. ,Xiang,L.,Wang,Z.,Hu,B.,
2013.IncreasedwaterstorageinNorthAmericaandScandinaviafromGRACE gravitydata.Nat.Geosci. 6(1),38–42.
Whitehouse,P.L.,Bentley,M.J.,LeBrocq,A.M.,2012a.Adeglacialmodelfor Antarc-tica: geologicalconstraintsand glaciological modelling asa basisfor anew modelofAntarcticglacialisostaticadjustment.Quat.Sci.Rev. 32,1–24. Whitehouse,P.L.,Bentley,M.J.,Milne,G.A.,King,M.A.,Thomas,I.D.,2012b.Anew
glacialisostaticadjustment modelfor Antarctica:calibratedandtestedusing observationsofrelativesea-levelchangeandpresent-dayupliftrates.Geophys. J.Int. 190(3),1464–1482.
Wu,P.,2004.Usingcommercialfiniteelementpackagesforthestudyofearth de-formations,sealevelsandthestateofstress.Geophys.J.Int. 158(2),401–408. Wu,X.,Heflin,M.B.,Schotman,H.,Vermeersen,L.L.A.,Dong,D.,Gross,R.S.,Ivins,
E.R.,Moore,A.W.,Owen,S.E.,2010.Simultaneousestimationofglobal present-daywatertransportandglacialisostaticadjustment.Nat.Geosci. 3(9),642–646. Wu,P.,Wang,H.,Steffen,H.,2013.Theroleofthermaleffectonmantleseismic anomaliesunderLaurentiaandFennoscandiafromobservationsofGlacial Iso-staticAdjustment.Geophys.J.Int. 192(1),7–17.
