Greenland Ice Sheet
Present
This configuration has been developped by F. GilletChaulet (IGEGrenoble). Please report any bug, missing information etc… and contribute to improve this configuration by adding new features!!
Summary:
The aim of this configuration is to provide an accurate representation of the current fast flow dynamics and run century scale forecast simulations of the Greenland Ice Sheet.
The mesh resolution is adapted using present day observations of the surface velocity and ice thickness to equidistribute the error on the ice fluxes. This allow to individually represent all the ice streams that are at least few km in width. In addition the mesh is refined near the margins to capture the retreat of the terrestrial margins in future simulations.
The model can be initialized with basal conditions that have been optimized using the observed surface velocities.
This configuration shows a good scalabity for parallel computing and allows to run 100 years forecast simulations in only few hours.
Below we give details on the configuration and provide links to usefull datasets to model the Greenland Ice Sheet with this configuration.
Configuration details:
ELMER
MESH

Adapting the mesh resolution: Elmer/Ice includes routines for mesh adaptation using
Mmg. See the
documentation.
DATA
Below are the links to standard datasets for the Greenland Ice Sheet.
See [ELMER_TRUNK]/elmerice/IceSheet/Greenland/DATA/GET_DATA.sh to automatically download and process (when required) the following datasets.


SMB forcing: MAR v3.5.2. Get the yearly outputs projected on the 5km grid.

Processing required: None
details: This dataset contains the slip coefficient and the (ssa) model velocity optimized to fit the observed velocity (MEaSUREs Multiyear Greenland Ice Sheet Velocity Mosaic, Version 1) projected on a regular grid with 250m spacing.

Processing required: None
details: This dataset contains the vertically averaged viscosity obtained after a paleospinup with the SICOPOLIS model (Seddik et al., J. Glaciol, 2012)
SSA CONFIGURATION
See [ELMER_TRUNK]/elmerice/IceSheet/Greenland/SSA for the configuration files.
Summary:
Physics: Shallow Stream approximation (SSA)
Thermal state: constant (initial viscosity field – see DATA)
Basal Condition: linear friction law (can easily be changed to nonlinear Weertman); Constant friction coefficient.
Calving: fixed front position
Terrestrial margins: retreat/advance (up to the domain limits)
Surface Mass Balance: SMB interpolated from RCM output.
Basal Mass balance: None
Description:
Parameters:
SSA.IN contains level0 parameters that can be adjusted by beginners; e.g. name of the simulation, time steps, Output intervals, path to initialisation data sets. SSA.IN is included at the top of the .sif configurations files and are in MATC format (comment lines start with !, parameter names start with $).
[ELMER_TRUNK]/elmerice/IceSheet/Parameters/Physical_Params.IN: contains the values of usual physical parameters (sea level elevation, gravity, glen index etc…) and can be included in SSA.IN.
Advanced users can adjust the configuration .sif files for advanced tuning of the simulations
There is no hard coded unit system in Elmer. Units must be consistant for all the parameters. Default values provided use : a, m and MPa
Initialisation:
INIT.sif is used to get the initial conditions.
Algorithm:
Initialisation of
H and
bedrock:
GridDataReader Solver read the variables
thickness and
bed in the
$TOPOGRAPHY_DATA netcdf file.
Initialisation of
Zs,
Zb and
GroundedMask:
The Flotation solver applies the floatation criteria to initialise the surface elevation, the bottom surface elevation and the mask that define floating and grounded nodes.
Initialisation of
smb:
GridDataReader Solver read the variable
smb_ice (the surface mass balance in ice equivalent) in the
$SMB_DATA netcdf file. This is a constant surface mass balance in m/a (ice equivalent) used to run, e.g., a control experiment with constant forcing.
Initialisation of the ice viscosity
Mu:
GridDataReader Solver read the variable
MuMean (the vertically averaged viscosity for the SSA; must be consistent with the value of the Glen Index) in the
$VISCOSITY_DATA netcdf file.
Initialisation of the friction parameter
slc:
GridDataReader Solver read the variable
C1 in the
$SLIP_DATA netcdf file. Must be consistent with the basal friction law (linear by default).
Initialisation of the observed velocity field
Vobs:
GridDataReader Solver read the variable
vnorm in the
$VELOCITY_DATA netcdf file. This is only for postprocessing purpose and has no influence on the model results.
Remarks:
All netcdf files must contains the dimensions and associated variables x and y. This can be adjusted in the .sif files.
The GridDataReader Solver does not check the variable attribute “_FillValue”; so filled values are used by the interpolation scheme and must have a physical meaning.
Simulation:
SSA.sif controls the simulation.
Algorithm:
Read the initialisation result file and initialise time=0
Do t=1, $Intervals
time=time+$dt
if [(
$FORCING_Update) == “always”] .AND. [ MOD( t1, (1/
$dt) ) = 0 ]; Update
SMB :
GridDataReader Solver read the variable
smb_ice (the surface mass balance in ice equivalent) in the
$FORCING_DATA netcdf file. Read the time index (floor(time)+1); i.e. the smb is updated every year, and the simulation time is converted to an integer to read the corresponding variable value.
The Flotation solver applies the floatation criteria to update the surface elevation, the bottom surface elevation and the mask.
The SSA solver compute the ice flow velocities according to current geometry and boundary conditions.

For postprocessing: Compute global metrics (Volume, Volume above floatation, fluxes, etc…). Additionnal reduction operations (e.g. min/max of dh/dt) can be performed with the Elmer SaveScalars solver.
If [ MOD( t1, (1/$OutPutIntervals) ) = 0 ]; Save simulation results in Elmer .result format and postprocessing .vtu format
end
Remarks:
The time stepping strategy is set to the Backward Differences Formulae (BDF) of order 2.
The bottleneck of the simulation is the resolution of the linear system of equations for the force balance. Here by default we use an iterative method (BiCGStab) with a preconditioning by ILU of order 2. Here the simulation is not aborted if the iterative method does not converge to the desired accuracy; However, changing the conditions of the test (Mesh, internal parameters) may deteriorate the convergence of the iterative method and may require to change the preconditioning or the linear solver method; see ElmerSolver manual for available options. A direct solver (mumps or cpardiso in parallel) will always converge but usually at a higher numerical cost.
Futur developments:
Here is a (nonordered) lists of possible developments to improve the current configuration; Contributions welcome:
Automatic postprocessing tools to generate standard plots
Implement subice shelf melt parameterisations
Thermomechanical coupling
Calving and moving calving fronts
Coupling basal friction and subglacial water pressure
Parameterise basal friction to extend the configuration to areas that are not currently icecovered
…
RESULTS
Example of results obtained with this configuration:
Surface velocities:
Numerical performances: Speedup and walltime achieved to run a 100 year simulation with constant forcing (dt=5e3 a):
Ice thickness variations: after 100 years for a control experiment and for a simulation forced with MAR (v3.5) for rcp 8.5 (with CanESM2); 2D views and associated histograms.
Volume change and fluxes for a control experiment (Mean (19791999) of MAR forced by ERA Interim) and for simulations with MAR forced by CanESM2 and MIROC5 under rcp8.5.
BASAL FRICTION OPTIMISATION
See [ELMER_TRUNK]/elmerice/IceSheet/Greenland/SSA_FRICTION_INIT for the configuration files.
Summary:
Optimise the basla friction field by minimising where:
is the basal friction coefficient
: is the mismatch between the model and observed velocity
: is the mismatch between the model and observed flux divergence
: is a regularisation term
The output contain all the fields required to run projections in the SSA directory.
This can be used to optimise the basal friction field for your own mesh or data sets.
Description:
Parameters:
OPTIM_BETA.IN contains level0 parameters that can be adjusted by beginners; e.g. name of the simulation, maximul number of iterations, Output intervals, path to initialisation data sets. OPTIM_BETA.IN is included at the top of the .sif configurations files and uses MATC format (comment lines start with !, parameter names start with $).
[ELMER_TRUNK]/elmerice/IceSheet/Parameters/Physical_Params.IN: contains the values of usual physical parameters (sea level elevation, gravity, glen index etc…) and can be included in OPTIM_BETA.IN.
Advanced users can adjust the configuration .sif files for advanced tuning of the optimisation
There is no hard coded unit system in Elmer. Units must be consistant for all the parameters. Default values provided use : a, m and MPa
Initialisation:
INIT_OPTIM_BETA.sif is used to get the initial conditions.
OPTIMISATION:
OPTIM_BETA.sif controls the optimisation.
The bottleneck of the simulation is the resolution of the linear system of equations for the force balance. Here by default we use an iterative method (BiCGStab) with a preconditioning by ILU of order 2. Here the simulation is not aborted if the iterative method does not converge to the desired accuracy; However, changing the conditions of the test (Mesh, internal parameters) may deteriorate the convergence of the iterative method and may require to change the preconditioning or the linear solver method; see ElmerSolver manual for available options. A direct solver (mumps or cpardiso in parallel) will always converge but usually at a higher numerical cost.