Ice sheet mass loss caused by dust and black carbon accumulation

Albedo is the dominant factor governing surface melt variability in the ablation area of ice sheets and glaciers. Aerosols such as mineral dust and black carbon (soot) accumulate on the ice surface and cause a darker surface and therefore a lower albedo. The darkening effect on the ice surface is currently not included in sea level projections, and the effect is unknown. We present a model framework which includes ice dynamics, aerosol transport, aerosol accumulation and the darkening effect on ice albedo and its consequences for surface melt. The model is applied to a simplified geometry resembling the conditions of the Greenland ice sheet, and it is forced by several temperature scenarios to quantify the darkening effect of aerosols on future mass loss. The effect of aerosols depends non-linearly on the temperature rise due to the feedback between aerosol accumulation and surface melt. According to our conceptual model, accounting for black carbon and dust in future projections of ice sheet changes until the year 3000 could induce an additional volume loss of 7 %. Since we have ignored some feedback processes, the impact might be even larger.


Introduction
The Greenland ice sheet contributes to sea level rise through dynamic processes and surface mass balance (SMB). After the year 2009, about 84 % mass loss of the Greenland ice sheet was due to a reduced SMB (Enderlin et al., 2014). The SMB is the net balance between snow accumulation and ablation consisting of melt and sublimation of snow and ice. Therefore calculations of SMB depend on the accuracy of snowfall and ablation. Ablation is largely controlled by nearsurface temperature and absorbed short-wave radiation. Surface albedo and the amount of incident radiation controls the portion of absorbed short-wave radiation (e.g. van Angelen et al., 2012;Bougamont et al., 2005;Tedesco et al., 2011).
In the Fifth Assessment Report (AR5; IPCC, 2013), ablation was computed by energy balance models (EBMs) and with the positive-degree-day (PDD) method, which does not explicitly include surface albedo. Regional climate models with energy balance models produced a 14-31 % higher sea level rise contribution from the Greenland ice sheet than models using PDD (Goelzer et al., 2013) because PDD can not account for the positive feedback of albedo to nearsurface temperature. Surface albedo implementations vary through the different EBMs considered in the AR5. The surface albedo is primarily determined by whether the surface consists of snow or ice and by the optical properties of the surface. All the EBMs in the AR5 used for sea level rise predictions have rather sophisticated snow albedo schemes, while the ice albedo is often constant in space and time (e.g. Mernild et al., 2010;Rae et al., 2012).
Surface albedo is determined by the optical properties of the surface, the angle of incidence of downward radiation and the ratio of direct to diffuse radiation. The optical properties of ice are mainly determined by the specific surface area of ice (that results from the combined effect of air bubbles, cracks, etc.). In addition, the optical properties are altered by the content of liquid water and impurities (Cuffey and Paterson, 2010). These impurities consist of aerosols such as mineral dust and black carbon (BC) and impurities related to biological activity. Dust makes up most of the im- purity mass on ice (Bøggild et al., 2010;Takeuchi et al., 2014), but BC has a higher effect on albedo per mass (Warren and Wiscombe, 1980). The darkening effect and dynamics of biological-related impurities remain to be quantified Yallop et al., 2012).
BC is a carbonaceous material which strongly absorbs visible light and is formed primarily in flames (Bond et al., 2013). Ice core data from Greenland revealed higher concentrations in ice dated younger than 1850 (McConnell et al., 2007) due to increased emissions after the industrial revolution. BC concentrations decreased back to almost preindustrial levels around 1950 due to successful mitigation in North America, which was found to be the source of BC in Greenland during that period. In contrast to that, BC emissions in Asia have been rising since the year 2000 (Lei et al., 2011;Lu et al., 2011).
Asia was found to be the main source of dust during the last glacial (Biscaye et al., 1997). During glacials the dust concentration in ice is 10-100 times higher than during interglacials such as the current Holocene (Steffensen, 1997). Ice from the last glacial is found in the deeper parts of the Greenland ice sheet (MacGregor et al., 2015) and therefore acts as a reservoir for large amounts of dust (Reeh et al., 1991;Bøggild et al., 1996).
Dust and BC are transported with the ice flow towards the ablation zone, where melt releases the particles. Over time the impurities from melt-out and atmospheric deposition accumulate, which enhances melt and therefore causes more particles to be released. This positive feedback is currently unquantified and not included in the AR5 sea level rise experiments.
In this study, we investigate the impact on ice sheet volume caused by the darkening effect of dust and BC on the ice surface of the Greenland ice sheet. We introduce a new model framework which includes transport, melt-out and accumulation of aerosols and its effect on ice albedo coupled to the SMB and ice dynamics. This model is then applied to a simplified geometry and forced by future climate scenarios until the year 3000.

Model framework and set-up
Aerosols reach the ablation zone by four different mechanisms ( Fig. 1): via direct atmospheric deposition from local sources (k I ), by large-scale transport (k II ) or by melt-out of englacial particles which have been transported via ice flow (k III ). Local production (k IV ) is only relevant for microbiological activity and is 0 for dust and BC considered in this study. The atmospheric contribution (k II ) and input from the surrounding tundra (k I ) are prescribed and constant in the simulations. Melt-out of aerosols (k III ) requires, besides the SMB, the near-surface concentration of englacial aerosols ([ι n,englacial )]). This is calculated via a tracer transport module which depends on the aerosol time series, velocities and ice sheet dimensions (Fig. 2). Aerosol accumulation, ice albedo and surface mass balance are presented in detail in Goelles and Bøggild (2015), and therefore their description is kept to a minimum in the following sections. The common quantities and physical constants are listed in Table 1. The model is realised in Mathematica (version 10; Wolfram Research, Inc., 2014) using self-coded solvers for the differential equations.

Ice dynamics
We employed the shallow ice approximation in plane strain, that is, only the vertical x-z plane is considered. The corresponding finite difference discretisation is described in detail by Greve and Blatter (2009). The ice thickness equa- Specific surface area of ice 2 cm 2 g −1 tion is solved by an implicit scheme (e.g Greve and Calov, 2002). A terrain-following coordinate transformation is used in order to explicitly compute surface values. These are required for the tracer transport in order to calculate the englacial aerosol concentration close to the surface and ultimately ice albedo values. Basal sliding is implemented with a Weertman-type sliding law including sub-melt sliding (Greve, 2005).

Tracer transport
Dust and BC concentrations are indirectly derived from the time of deposition t d and the aerosol time series (see next section). Besides the time t d , also the x coordinate of deposition x d is calculated via a semi-Lagrangian transport scheme based on Clarke and Marshall (2002). The scheme provides both t d and x d as functions of (t, x, z). It has advanced stability compared to Eulerian advection while still using a regular grid. The scheme is common in atmospheric models and has been used before in ice sheet models for δ 18 O transport (Clarke and Marshall, 2002;Clarke et al., 2005;Lhomme et al., 2005;Lhomme, 2004;Goelles et al., 2014). Especially the indirect method is suited for discontinuous data such as aerosol concentration while still using a coarse vertical resolution. At each grid point a back trajectory is calculated based on the velocity field from the ice dynamics module. Here we use the first-order backtracking (Clarke and Marshall, 2002;Clarke et al., 2005;Lhomme et al., 2005;Lhomme, 2004) for simplicity and efficiency, knowing that this method has its limitations close to the margin (Goelles et al., 2014) (see discussion in Sect. 5.2).

Aerosol time series
The aerosol concentration [ι n,englacial ] of every grid point can be derived from the aerosol time series f (t d ) and the deposi- Data from ice cores with an established age-depth model are required for the right-hand side of Eq. (1). A comparison between the NGRIP and NEEM ice cores (Bigler, 2012) of Greenland showed very similar results in the two cores although they are more then 350 km apart. Therefore, we use only the information of one ice core for the whole ice sheet, and f in Eq.
(1) has no explicit dependence on x. The aerosol with two additional trajectories of the future. In the period before data are available the BC concentration is set to 0. (b) The dust concentration from NGRIP (Ruth, 2007), where the missing Holocene period (after −8 ka) is assigned a constant value of 20 ng g −1 , and a value of 0 is assumed after the year 2000.
time series in Fig. 3 are extended beyond the data range in order to allow simulations until the year 3000.

Aerosol accumulation
The accumulation module handles both bulk BC and dust content of the snowpack and the mass per area on the ice surface. It is described in detail in Goelles and Bøggild (2015). Aerosols accumulate in the snowpack from sources k I and k II and are released onto the ice surface as the snow disappears at the start of the melt season. As ice is exposed, the sources k I and k II , in addition to aerosols from melt-out (k III ), directly contribute to the aerosol amount on the ice surface. The accumulation is counteracted by a daily reduction fraction r ice of the order of 1 per mille per day when ice is exposed. In the snowpack all aerosols are conserved (r snow = 0).

Ice albedo
Ice albedo is determined via the specific surface area of icê S and the reductions caused by dust and BC, as described in Goelles and Bøggild (2015). The dust concentration is converted into an equivalent BC concentration, and the albedo reduction is calculated with the parameterisation dα ι,1 by Gardner and Sharp (2010). A conversion from aerosol mass per area to ppm is required as the parameterisation is formulated in terms of BC concentration. This is done via the effective depth d eff representing the absorption length in ice.
Only an active fraction F of the aerosols influences the ice albedo, while the remainder is concealed in cryoconite holes.

Surface mass balance
The surface mass balance is calculated with a simplified energy-balance model (e.g. Oerlemans, 2001) which is optimised for Greenland (Robinson et al., 2010(Robinson et al., , 2011Goelles and Bøggild, 2015). It also includes a simple snowpack model and refreezing. If the snow depth exceeds 5 m, the ice equivalent part is added to the ice thickness H (x, t). The snow depth is also required for the surface albedo to distinguish between ice and snow. The snow albedo is divided into the two values α snow, dry and α snow, wet . These values are kept constant throughout the experiments in order to separate the ice albedo effect.

EISMINT and RCPs
The experiments are based on the set-up of the European Ice Sheet Modelling Initiative (EISMINT Phase 2; Payne et al., 2000) combined with AR5 temperature projections for Representative Concentration Pathways (RCP; IPCC, 2013), with standard parameters listed in Table 1. These parameters are common for all simulations unless stated otherwise. The EISMINT boundary conditions are symmetrical around x su and roughly mimic conditions of a west-east cross section of Greenland. The simple, symmetric geometry is ideally suited as a test case for new methods with the additional benefit of ease of interpretation. In Goelles and Bøggild (2015) we introduced an annual temperature parameterisation with a trapezoidal shape: where T + is the mean summer temperature between t ,start and t ,end . Based on the EISMINT set-up the temperature distribution in x direction is where T +,su is the minimum temperature located in at the summit and S T is the gradient with horizontal distance. Temperature evolution scenarios are based on long-term projections until 2300 for the global annual mean surface temperature (Fig. 4). The RCPs are named after their expected radiative forcing in the year 2100 compared to preindustrial values. The four scenarios are RCP2.6, which is a mitigation scenario leading to low forcing, two stabilising scenarios (RCP4.5 and RCP6.0) and one scenario with very high concentrations (RCP8.5). The temperature anomalies are added each year to every grid cell.
The precipitation (see Fig. S1a in the Supplement) is kept on the same level for all RCP scenarios and corresponds to 0.5 m a −1 ice equivalent in the accumulation zone, except for one experiment with an increase of 20 % in precipitation (RCP4.5 + 20 % precip †).
The insolation at the top of the atmosphere (see Fig. S1b) is based on the calculations by Liou (2002) for 67 • N, in the centre-south of the Greenland ice sheet for which the EIS-MINT boundary conditions are roughly representative.

Spin-up
A spin-up model run acts as the common starting point for all experiments. The ice sheet after the spin-up is in steady state with the dust and BC concentrations computed over the whole domain. The spin-up starts at −30 ka and ends in the year 2000 with a resolution in x of 25 km and 20 layers in z direction. During the spin-up, the ice albedo is based on a specific surface area of 2 cm 2 g −1 of ice of a density of 880 kg m −3 and above (Dadic et al., 2013), equivalent to a clean ice albedo of 0.53.

Spin-up
It takes about 10 000 years for the spin-up to reach equilibrium (see Fig. S2).
The simulated englacial aerosol concentration, depositional x and time are shown in Fig. 5. Of those quantities, the depositional time t d (panel d) is the most important one as the englacial BC and dust concentrations are inferred from it (Eq. 1). Ice older than −20 ka can only be found at the very bottom and not close to the surface, which indicates that the spin-up is sufficiently long. Figure 6 shows surface values of the same quantities as in Fig. 5. At each time step, the values are joined together, which results in graphs showing the surface evolution of the respective quantities over the whole period. The ice sheet expands in the first few thousand years until a dynamic equilibrium is reached, as can be seen in the widening of the surface values in the upper sections of the graphs.
The high BC concentrations during the period 1850-1950 are indicated in light green in Figs. 5d and 6d. This period is short compared to the whole spin-up and is visible as a thin light-green stripe close to the surface in the cross section (Fig. 5b) and in the bottom part of the time series plot (Fig. 6b). Most of the ice was deposited before 1850, and therefore most parts of the ice have a 0 BC concentration.
The sequence of peaks in the dust time series can be seen in Fig. 6a as a sequence of red stripes, indicating high concentrations. The stripes are first visible in the accumulation area in the centre and later in the ablation zone on the sides. The low Holocene values dominate after −10 ka, and dust concentrations are only high in the outermost grid cells on each site. This is the region where ice from the last glacial reaches the surface (Fig. 6d).

Ice sheet volume evolution
Results of the experiments driven by the four RCP scenarios with (dashed lines) and without (solid lines) aerosols from the year 2000 to 3000 are shown in Fig. 7. In general, the ice sheet is thinner everywhere in the year 3000 when dust and BC are considered (panel b). The main height differences (panel c) are located at the margin. This is a result of the lower annual SMB (panel d), which in general is lower in the ablation zone when aerosols are considered. Panel e shows the volume change relative to the constant climate/noaerosols run. The inset shows the period with transient temperatures until 2300 in more detail. Panel f displays the volume change of each RCP scenario when aerosols are considered relative to the no-aerosols run. The higher the temperatures, the more influential are the aerosols.
After the year 2300, the temperatures remain unchanged at the 2300 level, and the ice dynamics still responds to the increased melt at the margin, and in addition aerosols continue to accumulate and fade away.
In the year 3000, the ice sheet volume is 1.46 % smaller for RCP4.5 and 7.61 % smaller for RCP8.5 when aerosols are considered.

RCP4.5 in more detail
The evolution of the RCP4.5 simulation with aerosols (RCP4.5 aerosols) is now presented in more detail. Figure   at the very margin are caused by dust (see Fig. 9a and c).
The low values below 0.40 close the accumulation zone are caused by BC, as the younger ice there contains BC from the period with higher concentrations. Until 2300 black carbon lowers the ice albedo of up to three grid points, and later on only one is affected. The main part of the ablation zone in Fig. 9 has low values of aerosol concentrations. This is because the englacial BC is 0 when the age of the ice is younger than at the start of the time series but older than glacial ice with high dust concentrations. Therefore, the BC content at the surface is only due to the atmospheric deposition which is low (Goelles and Bøggild, 2015). The surface amount of BC accumulates to just 200 ng m −2 . The surface content of dust contains, besides the atmospheric signal, melt-out of the low englacial concentration during the Holocene of 20 ng g −1 and accumulates to 0.4 g m −2 . The combined effect of dust and BC causes a constant albedo reduction of about 1.0 %. This causes a lower SMB in the central part of the ablation zone when aerosols are considered (Fig. 7).
Irregular peaks in the aerosol concentration can be seen in the englacial concentrations in panels a and b. Those peaks cause also higher amounts of aerosols at the surface due to the slow loss of aerosols at the surface. Nevertheless, over millennial time scales, the prolonged residence time of some decades is comparably short. This can be seen for example in Fig. 9b and d after the year 2550. Figure 10 shows the volume change relative to the RCP4.5 simulation without aerosols, and how the results depend on different parameter values. The blue thick line is the same as in Fig. 7f. The simulations marked with † are compared to another reference. Run RCP4.5 mound aerosol † is compared to a run with mound topography and RCP4.5 conditions but no aerosols. RCP4.5 + 20 % precip † is compared to a simulation with 20 % more precipitation and RCP4.5 conditions but no aerosols. Simulation RCP4.5 aerosols † is compared to a reference for which the ice albedo already considers the 1.0 % reduction mentioned above. Experiment RCP4.5 mound aerosol † uses 500 m mounds from the EISMINT (experiment K) set-up, and the volume change is calculated from a spin-up with these mounds. The effect of aerosols is slightly lower during the simulation but very close to the original simulation in the year 3000.
The highest impact of aerosols was reached with an active fraction F of 0.8, closely followed by an effective depth d eff of 3 m. Both parameters determine how much aerosols influence the ice albedo.
Considering only BC (RCP4.5 BC only) leads to 0.42 % and dust alone to 1.16 % additional ice sheet volume loss in the year 3000 compared to the simulation RCP4.5 without aerosols. Therefore, dust is responsible for the major part of ice sheet loss.

The effect of aerosols on ice volume
The additional ice loss caused by aerosols increases nonlinearly with temperature ( Fig. 7e and f). The non-linearity is caused by the relationship between BC concentration and albedo reduction (Fig. S4)  between ice melt and aerosol accumulation. We did not include other feedbacks such as the elevation feedback and the albedo feedback to temperature. There is a weak negative elevation feedback included due to the parameterisation of the transmissivity (Robinson et al., 2010). Otherwise, the elevation and albedo feedbacks are positive, which further enhances melt. Due to this, and since we used global temperature anomalies for the simulations, the strength of the aerosol effect found in this study is expected to be a lower boundary. The two BC scenarios (RCP4.5 BC future0 and RCP4.5 BC future3) result in a similar evolution as the standard setup. This is a consequence of the time lag between aerosol deposition and melt-out at the ablation zone. Only the highest grid point in the ablation zone is affected by aerosols deposited after 2000 (see Fig. 9f). Therefore the effect of different BC and also dust emissions from the year 2000 onwards is limited until the year 3000. Nevertheless, BC and dust deposition has a direct effect on snow albedo and a long-term effect on ice albedo.
Ice albedo is sensitive to the small addition of impurities when the ice is clean (Fig. S4a). If the same amount of impurities is added to an already dark ice surface, the additional effect is weaker as compared to a clean surface. This could lead to an overestimation of the effect, whereas the missing feedbacks can lead to an underestimation; which effects dominate requires additional research. Nevertheless, the nonlinear nature of the effect and the amplification due to rising temperatures presented in this study remain.

Assumptions, simplifications and uncertainties
The magnitude and non-linear dependence on future temperature of the aerosol effect were shown with a new model framework and applied to a simplified geometry experiment. Here we discuss if and how the assumptions and simplifications could have lead to an overestimation of the effect.
The experiments were performed using a flat bedrock topography; however, the effect of 500 m high mounds was also investigated (Fig. 10, run RCP4.5 mound aerosol †). The effect was sometimes weaker but very similar at the end of the simulation. The topography of the bedrock below the Greenland ice sheet is more complicated than sinusoidal mounds (Bamber et al., 2013); nevertheless, the effect remained even though the timing when aerosols emerge is influenced by the topography.
The ice dynamics module was compared to an analytical solution (Cuffey and Paterson, 2010). Close to the ice margin the shallow-ice approximation employed can be violated, and thus velocities can be computed incorrectly when the surface slope is steep (e.g. Ahlkrona et al., 2013;Greve and Blatter, 2009). This typically occurs in areas where the ice sheet reaches the coast and calving is present. However, in our case, we have pronounced ablation zones on both sides of the two-dimensional domain (Fig. S3). Therefore, the simulated ice sheet never reaches the boundaries of the domain and maximal surface slopes at the ice margin are just around 1 • , which is well within acceptable limits for the validity of the shallow-ice approximation.
The tracer transport, which is used to derive the englacial aerosol concentration, influences the rate of melt-out and therefore the amount of accumulated aerosols on the ice surface. The accuracy of the module depends on the velocities and the numerics of the transport scheme. Here we use a first-order scheme which was found to deliver different results in the ablation zone (Goelles et al., 2014). Nevertheless the calculated ages are comparable to radio stratigraphy from Greenland (MacGregor et al., 2015) as the horizon of the boundary between Holocene and glacial ice is at similar depths.
The aerosol time series of mineral dust and BC have a more direct effect on the results. They determine the englacial concentration, which directly influences the ice albedo and determines the amount of melt-out. We assumed that the aerosol concentration solely depends on the age of the ice because distant ice cores correlate well (Bigler, 2012).
However, if this were not the case, and aerosols were deposited in smaller patches rather than uniformly, then they would also be released in smaller areas in the ablation zone. Owing to the reduced response when aerosols are added to an already dark surface (Fig. S4), the effect of aerosols would then be smaller.
Nevertheless, close to the ELA the surrounding tundra could contribute locally to the englacial dust concentration (Wientjes et al., 2012). As mentioned in Goelles and Bøggild (2015) there might be a "threshold elevation" up to which local dust contributes. BC concentration peaks due to forest fires (Fig. 3) are only available between 1788 and 2000. Also, dust ice core data of the Holocene are not available. Since the central part of the ablation zone has its origin in the Holocene, the englacial aerosols consist of only dust which has a low and constant value. This is the cause for the slight decrease of albedo in the centre part of the ablation zone.
The aerosol accumulation determines the amount of aerosols which darken the ice surface and therefore plays an important role in the overall estimation of the effect. The direct input from the atmosphere (k II ) was constant in this study, and the input from the tundra was 0. This is most likely true for BC as local sources on Greenland are negligible. Below the "threshold elevation", dust from the surrounding tundra might contribute significantly. The dust amount at the margin is already large due to melt-out of glacial dust. Therefore, an additional amount of dust from the tundra has a smaller effect as the surface is already dark (Fig. S4). Beside the input, the daily reduction on ice r ice determines the surface amount. This parameter was found to be of the order of 1 per mille per day, which is comparable with measurements (Goelles and Bøggild, 2015). Under warming scenarios, the daily reduction might increase with increased surface run-off due to meltwater and rainfall, which would have a stabilising effect that is currently not captured.
Under recent conditions, the amount of meltwater run-off is typically 1 to 2 magnitudes higher than summer rainfall. Therefore, even with more rainfall during summer in the future, the potential increase of the daily reduction will most likely be determined by the meltwater run-off.
The relationship of aerosols and ice albedo depends on the specific surface area of ice, englacial concentration, the amount at the surface, the conversion to BC concentration and the darkening parameterisation as well as the active fraction F. The conversion to concentration is necessary because of the formulation of the darkening parameterisation and depends on the effective depth d eff . The value is based on the absorption length in ice which varies greatly with wavelength (Warren and Brandt, 2008). The conversion of mass per area to concentration (ppm) would not be necessary if an ice albedo reduction parameterisation based on mass of BC per area existed. Also, the active fraction F is similarly powerful as the effective depth and not well constrained. The parameter lumps together all surface processes which keep aerosols from darkening the ice surface. This part of the model is based on a lower level of understanding and yet causes a wide spread in the outcome; therefore it requires more attention (see Goelles and Bøggild, 2015;Hodson, 2014). Nevertheless, even at a low active fraction or a high effective depth the effect remains (Fig. 10).
The surface mass balance determines the snow depth, which governs how long ice is exposed and the amount of outcropping aerosols (k III ). The time span of ice exposure is important for the actual effect of ice albedo and, secondly, for the period of aerosol reduction (r ice,n ). Especially the longwave radiation coefficient λ is very influential as it scales the temperature dependence of melt (Fig. S4). Nevertheless, tests with data from western Greenland showed results comparable to observations (Goelles and Bøggild, 2015).
Overall the magnitude of the effect of aerosols depends on all the parts discussed, but even though limitations exist the effect is intuitively understandable. More ice melt causes an increase of outcropping aerosols, which have a long residence time at the surface. These aerosols darken the ice surface, which further enhances ice melt. The presented model components capture the main effect, while the exact timing when and where aerosols emerge is harder to achieve.

Conclusions
We tested the ice volume response to darker ice in the ablation zone caused by accumulation of dust and black carbon. We introduced a new model framework which includes advection, melt-out and accumulation of aerosols and its dark-ening effect on the surface mass balance. The response of the ice volume to the aerosols depends non-linearly on the future temperature because of a positive feedback between ice melt and aerosol melt-out which is disproportionally larger in warmer climate scenarios. The exact timing when and where darkening occurs is difficult to achieve; however, the overall effect is captured with the presented framework.
In order to isolate the effect of ice albedo, we kept the snow albedo unaffected by impurities. A 0.01 lower fresh snow albedo causes a decrease of 27 Gt a −1 of the total surface mass balance for the Greenland ice sheet (Dumont et al., 2014). Therefore the combined effect of impurities on snow and ice albedo is significantly larger then the ice albedo effect alone.
The effect of black carbon from the industrial revolution was significant and will have an effect for a long time to come. Currently, Asian black carbon emissions are rising which could lead, besides to the darkening of snow, to darker ice via direct deposition and over long time scales via ice flow.
The presented principles are not just limited to Greenland in the future but could also be applied to palaeo-climatic studies, detailed studies for alpine glaciers or the termination of the little ice age (Painter et al., 2013).
We investigated the response under the RCP4.5 scenario in more detail and tested the sensitivity to several parameters. Considering that the temperatures in the Arctic rise higher than the global mean and the elevation and ice albedo feedbacks were not considered, the estimated effect can be seen as a lower boundary estimate. The presented simulations were based on a simplified geometry in two dimensions; nevertheless the age structure and overall system resembles the Greenland ice sheet. The presented simulations should not be seen as forecasts but emphasise the importance of the effect.
The Supplement related to this article is available online at doi:10.5194/tc-9-1845-2015-supplement.