6.5. Ice-thickness evolution

Several options are available for the evolution of the ice thickness, selected by the run-specs-header parameter THK_EVOL:

  • 0: No evolution of the ice thickness, kept fixed on the initial thickness.

  • 1: Free evolution of the ice thickness.

  • 2: Evolution of the ice thickness, but the ice topography (zs, zb, zl, H) is nugded towards a prescribed target with a time-dependent relaxation time read from the file TARGET_TOPO_TAU0_FILE (to be located in sico_in/general).

  • 3: Evolution of the ice thickness, but the ice topography (zs, zb, zl, H) is nugded towards a prescribed target with the constant relaxation time TARGET_TOPO_TAU0.

For the cases THK_EVOL >= 1, the ice-thickness equation

(6.10)\[\frac{\partial{}H}{\partial{}t} = -\left( \frac{\partial{}Q_x}{\partial{}x} + \frac{\partial{}Q_y}{\partial{}y} \right) + a_\mathrm{s} + a_\mathrm{b}\]

is solved, where \(H\) is the ice thickness, \((Q_x,Q_y)\) the volume flux (depth-integrated horizontal velocity), \(a_\mathrm{s}\) the surface mass balance (SMB) and \(a_\mathrm{b}\) the basal mass balance (BMB). Note that SMB and BMB are consistently counted as positive for a volume gain and negative for a volume loss.

For the cases THK_EVOL = 2, 3, nudging of the ice topography towards a prescribed target topography (e.g., a slightly smoothed present-day topography) is employed. During each time step, after solving the ice-thickness equation, the relaxation equations

(6.11)\[\frac{\partial{}h}{\partial{}t} = -\frac{h-h_\mathrm{target}}{\tau}\,, \quad \frac{\partial{}b}{\partial{}t} = -\frac{b-b_\mathrm{target}}{\tau}\,, \quad \frac{\partial{}z_\mathrm{l}}{\partial{}t} = -\frac{z_\mathrm{l}-z_\mathrm{l,target}}{\tau}\,,\]

are applied, where \(h\) is the surface topography, \(b\) the ice-base topography (\(H=h-b\)), \(z_\mathrm{l}\) the lithosphere-surface (bedrock/seabed) topography, \((\cdot)_\mathrm{target}\) the respective nudging target and \(\tau\) the relaxation time (Rueckamp et al. [52]). The nudging target is specified by the parameter TARGET_TOPO_DAT_NAME, to be assigned a text string with the file name (typically a time-slice output file from a previous simulation, see Section “Output files”). In addition, the path where this file is located must be given as an option -t /path/to/nudging/target/directory to the script sico.sh when running SICOPOLIS (see Section “How to run a simulation”).

Nudging is equivalent to applying an SMB correction, which is diagnosed by the model.

For the cases THK_EVOL >= 1, the maximum ice extent can be constrained spatially by a prescribed mask file, specified by the parameter MASK_MAXEXTENT_FILE. This file, in either NetCDF (*.nc, recommended) or ASCII (any other file extension) format, is to be located in sico_in/ant for Antarctica, sico_in/grl for Greenland, etc. The mask must be a 2D integer array that matches the horizontal grid, with values 1 for grid points allowed to glaciate, and 0 for grid points not allowed to glaciate. For NetCDF, the required variable name of the mask is mask_maxextent. If MASK_MAXEXTENT_FILE is set to 'none' or undefined, no spatial constraint is applied.

The numerical scheme for solving the ice-thickness equation can be chosen by the parameter CALCTHK:

  • 1: Explicit scheme for the diffusive SIA ice-surface equation.

  • 2: Over-implicit scheme for the diffusive SIA ice-surface equation, iterative built-in SOR solver.

  • 4: Explicit scheme for the general ice-thickness equation.

For more details on the cases CALCTHK = 1, 2, see Greve and Calov [32]. However, these cases are only applicable for SIA dynamics (see “Ice dynamics”), and they do not preserve mass well on complex bed topographies. Therefore, for most real-world problems, CALCTHK = 4 is the recommended setting.