Source-transformation Differentiation of a C++-like Ice Sheet model

International audienceAlgorithmic Differentiation (AD) has become one of the most powerful tools to improve our understanding of theEarth System. If AD has been used by the ocean and atmospheric circulation modelingcommunity for almost 20 years, it is relatively new in the ice sheet modeling community. The Ice SheetSystem Model (ISSM) is a C++, object-oriented, massively parallelized, new generation ice sheet model that recentlyimplemented AD to improve its data assimilation capabilities. ISSM currently relies on Object Overloading throughADOL-C and AMPI. However, experience shows that Object Overloading AD on ISSM is significantly more memoryintensive compared to the primal code. We want to investigate other AD approaches to improve the performance ofthe AD adjoint. Yet, to our knowledge, there is no source-to-source AD tool that supports C++.To overcome this problem, we have developed a prototype of ISSM entirely in C, called Boreas, in order to testsource-to-source transformation and compare the performance of these two approaches to AD. Boreas is a clone ofISSM, the main difference with ISSM is that all the objects are converted to C-structures and some function nameshave been adapted in order to be unique, but the code architectures are identical. The programming style of Boreas isa first attempt at defining a programming style of (or a sub-language of) C++ that source-transformation AD couldhandle. We present here the first results of Source-Transformation AD of Boread with the AD tool Tapenade

Extended enthalpy formulations in the Ice-sheet and Sea-level System Model (ISSM) version 4.17: discontinuous conductivity and anisotropic streamline upwind Petrov--Galerkin (SUPG) method

The thermal state of an ice sheet is an important control on its past and future evolution. Some parts of the ice sheet may be polythermal, leading to discontinuous properties at the cold–temperate transition surface (CTS). These discontinuities require a careful treatment in ice sheet models (ISMs). Additionally, the highly anisotropic geometry of the 3D elements in ice sheet modelling poses a problem for stabilization approaches in advection-dominated problems. Here, we present extended enthalpy formulations within the finite-element Ice-Sheet and Sea-Level System model (ISSM) that show a better performance than earlier implementations. In a first polythermal-slab experiment, we found that the treatment of the discontinuous conductivities at the CTS with a geometric mean produces more accurate results compared to the arithmetic or harmonic mean. This improvement is particularly efficient when applied to coarse vertical resolutions. In a second ice dome experiment, we find that the numerical solution is sensitive to the choice of stabilization parameters in the well-established streamline upwind Petrov–Galerkin (SUPG) method. As standard literature values for the SUPG stabilization parameter do not account for the highly anisotropic geometry of the 3D elements in ice sheet modelling, we propose a novel anisotropic SUPG (ASUPG) formulation. This formulation circumvents the problem of high aspect ratio by treating the horizontal and vertical directions separately in the stabilization coefficients. The ASUPG method provides accurate results for the thermodynamic equation on geometries with very small aspect ratios like ice sheets

An approach to computing discrete adjoints for MPI-parallelized models applied to Ice Sheet System Model 4.11

Within the framework of sea-level rise projections, there is a strong need for hindcast validation of the evolution of polar ice sheets in a way that tightly matches observational records (from radar, gravity, and altimetry observations mainly). However, the computational requirements for making hindcast reconstructions possible are severe and rely mainly on the evaluation of the adjoint state of transient ice-flow models. Here, we look at the computation of adjoints in the context of the NASA/JPL/UCI Ice Sheet System Model (ISSM), written in C++ and designed for parallel execution with MPI. We present the adaptations required in the way the software is designed and written, but also generic adaptations in the tools facilitating the adjoint computations. We concentrate on the use of operator overloading coupled with the AdjoinableMPI library to achieve the adjoint computation of the ISSM. We present a comprehensive approach to (1) carry out type changing through the ISSM, hence facilitating operator overloading, (2) bind to external solvers such as MUMPS and GSL-LU, and (3) handle MPI-based parallelism to scale the capability. We demonstrate the success of the approach by computing sensitivities of hindcast metrics such as the misfit to observed records of surface altimetry on the northeastern Greenland Ice Stream, or the misfit to observed records of surface velocities on Upernavik Glacier, central West Greenland. We also provide metrics for the scalability of the approach, and the expected performance. This approach has the potential to enable a new generation of hindcast-validated projections that make full use of the wealth of datasets currently being collected, or already collected, in Greenland and Antarctica

Greenland Subglacial Discharge as a Driver of Hotspots of Increasing Coastal Chlorophyll Since the Early 2000s

Subglacial discharge emerging from the base of Greenland\u27s marine-terminating glaciers drives upwelling of nutrient-rich bottom waters to the euphotic zone, which can fuel nitrate-limited phytoplankton growth. Here, we use buoyant plume theory to quantify this subglacial discharge-driven nutrient supply on a pan-Greenland scale. The modeled nitrate fluxes were concentrated in a few critical systems, with half of the total modeled nitrate flux anomaly occurring at just 14% of marine-terminating glaciers. Increasing subglacial discharge fluxes results in elevated nitrate fluxes, with the largest flux occurring at Jakobshavn Isbræ in Disko Bay, where subglacial discharge is largest. Subglacial discharge and nitrate flux anomaly also account for significant temporal variability in summer satellite chlorophyll a (Chl) within 50 km of Greenland\u27s coast, particularly in some regions in central west and northwest Greenland

A JavaScript API for the Ice Sheet System Model: towards on online interactive model for the Cryosphere Community

Abstract. Earth System Models (ESMs) are becoming increasingly complex, requiring extensive knowledge and experience to deploy and use in an efficient manner. They run on high-performance architectures that are significantly different from the everyday environments that scientists use to pre and post-process results (i.e. MATLAB, Python). This results in models that are hard to use for non specialists, and that are increasingly specific in their application. It also makes them relatively inaccessible to the wider science community, not to mention to the general public. Here, we present a new software/model paradigm that attempts to bridge the gap between the science community and the complexity of ESMs, by developing a new JavaScript Application Program Interface (API) for the Ice Sheet System Model (ISSM). The aforementioned API allows Cryosphere Scientists to run ISSM on the client-side of a webpage, within the JavaScript environment. When combined with a Web server running ISSM (using a Python API), it enables the serving of ISSM computations in an easy and straightforward way. The deep integration and similarities between all the APIs in ISSM (MATLAB, Python, and now JavaScript) significantly shortens and simplifies the turnaround of state-of-the-art science runs and their use by the larger community. We demonstrate our approach via a new Virtual Earth System Laboratory (VESL) Web site

Drivers of Change of Thwaites Glacier, West Antarctica, Between 1995 and 2015

We run several transient numerical simulations applying these three perturbations individually. Our results show that ocean-induced ice-shelf thinning generates most of the observed grounding line retreat, inland speed-up, and mass loss, in agreement with previous work. We improve the agreement with observed inland speed-up and thinning by prescribing changes in ice-shelf geometry and a reduction in basal traction over areas that became ungrounded since 1995, suggesting that shelf breakups and thinning-induced reduction in basal traction play a critical role on Thwaites's dynamics, as pointed out by previous studies. These findings suggest that modeling Thwaites's future requires reliable ocean-induced melt estimates in models that respond accurately to downstream perturbations

Implementation and performance of adaptive mesh refinement in the Ice Sheet System Model (ISSM v4.14)

Accurate projections of the evolution of ice sheets in a changing climate require a fine mesh/grid resolution in ice sheet models to correctly capture fundamental physical processes, such as the evolution of the grounding line, the region where grounded ice starts to float. The evolution of the grounding line indeed plays a major role in ice sheet dynamics, as it is a fundamental control on marine ice sheet stability. Numerical modeling of a grounding line requires significant computational resources since the accuracy of its position depends on grid or mesh resolution. A technique that improves accuracy with reduced computational cost is the adaptive mesh refinement (AMR) approach. We present here the implementation of the AMR technique in the finite element Ice Sheet System Model (ISSM) to simulate grounding line dynamics under two different benchmarks: MISMIP3d and MISMIP+. We test different refinement criteria: (a) distance around the grounding line, (b) a posteriori error estimator, the Zienkiewicz-Zhu (ZZ) error estimator, and (c) different combinations of (a) and (b). In both benchmarks, the ZZ error estimator presents high values around the grounding line. In the MISMIP+ setup, this estimator also presents high values in the grounded part of the ice sheet, following the complex shape of the bedrock geometry. The ZZ estimator helps guide the refinement procedure such that AMR performance is improved. Our results show that computational time with AMR depends on the required accuracy, but in all cases, it is significantly shorter than for uniformly refined meshes. We conclude that AMR without an associated error estimator should be avoided, especially for real glaciers that have a complex bed geometry.121215232CONSELHO NACIONAL DE DESENVOLVIMENTO CIENTÍFICO E TECNOLÓGICO - CNPQ140186/2015-