tross.rst - pism - [fork] customized build of PISM, the parallel ice sheet model (tillflux branch) (HTM) git clone git://src.adamsgaard.dk/pism (DIR) Log (DIR) Files (DIR) Refs (DIR) LICENSE --- tross.rst (12024B) --- 1 .. include:: ../../global.txt 2 3 .. _sec-ross: 4 5 An SSA flow model for the Ross Ice Shelf in Antarctica 6 ------------------------------------------------------ 7 8 As part of the EISMINT series of intercomparisons, MacAyeal and others :cite:`MacAyealetal` 9 successfully validated early-1990s ice shelf numerical models using velocity data for the 10 Ross ice shelf. The data were from the RIGGS survey :cite:`RIGGS2`, acquired in the period 11 1973--1978 and measured at a few hundred locations in a grid across the shelf. Substantial 12 modelling developments followed EISMINT-Ross, including inverse modeling to recover 13 depth-averaged viscosity :cite:`RommelaereMacAyeal` and parameter-sensitivity studies 14 :cite:`HumbertGreveHutter`. Previous PISM versions set up the EISMINT-Ross flow model and 15 performed the diagnostic computation, with RIGGS data for validation. 16 17 However, availability of rich new data sets for ice sheet modeling, including the ALBMAP 18 v1 :cite:`LeBrocqetal2010` ice sheet geometry, bedrock, and climate data set, and the 19 radar-derived (InSAR) MEaSUREs Antarctica Velocity Map :cite:`Rignotetal2011`, allows us to 20 use more complete, recent, and higher-resolution data for the same basic job. Furthermore 21 one can extend the diagnostic Ross ice shelf calculation both to other ice shelves around 22 Antarctica and to time-evolving ("prognostic") cases using the eigencalving 23 :cite:`Levermannetal2012` mechanisms. 24 25 The scripts in this subsection are found in directory ``examples/ross/``. In summary, the 26 script ``preprocess.py`` downloads easily-available data and builds a NetCDF input file 27 for PISM. For the diagnostic computation we document first, the script ``run_diag.sh`` (in 28 subdirectory ``examples/ross/diagnostic/``) runs PISM. The script ``plot.py`` shows a 29 comparison of observations and model results, as in :numref:`fig-rosspython`. 30 31 Preprocessing input data 32 ^^^^^^^^^^^^^^^^^^^^^^^^ 33 34 NSIDC_ requires registration to download the `MEaSUREs InSAR-Based Antarctica Ice Velocity 35 Map <measures-ross_>`_; please register, log in, and get the *first*, 900 m version 36 of it (``antarctica_ice_velocity_900m.nc``) here: 37 38 https://n5eil01u.ecs.nsidc.org/MEASURES/NSIDC-0484.001/1996.01.01/ 39 40 Put this file in ``examples/ross``. 41 42 The script ``preprocess.py`` then downloads ALBMAP using ``wget``; these two files total 43 around 350 Mb. Then it uses NCO_ to cut out the relevant portion of the grid and CDO_ to 44 conservatively-interpolate the high-resolution (900 m) velocity data onto the coarser (5 45 km) geometry grid used in ALBMAP. The script ``nc2cdo.py`` from directory ``util/``, 46 prepares the NetCDF file for the application of CDO, which requires complete projection 47 information. Run 48 49 .. code-block:: none 50 51 cd examples/ross/ 52 ./preprocess.py 53 54 to download ALBMAP and finish the pre-processing. 55 56 The NetCDF file ``Ross_combined.nc`` produced by ``preprocess.py`` contains ice thickness, 57 bed elevations, surface temperature, net accumulation, as well as latitude and longitude 58 values. All of these are typical of ice sheet modelling data, both in diagnostic runs and 59 as needed to initialize and provide boundary conditions for prognostic (evolutionary) 60 runs; see below for the prognostic case with these data. The ``_combined`` file also has 61 variables ``u_ssa_bc`` and ``v_ssa_bc`` for the boundary values used in the (diagnostic 62 and prognostic) computation of velocity. They are used at all grounded locations and at 63 ice shelf cells that are immediate neighbors of grounded ice. The variable ``bc_mask`` 64 specifies these locations. Finally the variables ``u_ssa_bc,v_ssa_bc``, which contain 65 observed values, are used after the run to compare to the computed interior velocities. 66 67 Diagnostic computation of ice shelf velocity 68 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ 69 70 The diagnostic velocity computation bootstraps from ``Ross_combined.nc`` and performs a 71 short run; in the `211\times 211` grid case we demonstrate below, the key parts of the 72 PISM command are 73 74 .. code-block:: none 75 76 pismr -i ../Ross_combined.nc -bootstrap -Mx 211 -My 211 -Mz 3 -Lz 3000 -z_spacing equal \ 77 -surface given -stress_balance ssa -energy none -no_mass -yield_stress constant -tauc 1e6 \ 78 -pik -ssa_dirichlet_bc -y 1.0 -ssa_e 0.6 -ssafd_ksp_monitor 79 80 The computational grid here is the "native" `5` km data grid used in ALBMAP. Regarding the 81 options, 82 83 - The maximum thickness of the ice is `2766` m so we choose a height for the computational 84 box large enough to contain the ice (i.e. ``-Lz 3000``). Vertical grid resolution 85 is, however, unimportant in this case because we use the SSA stress balance only, and 86 the temperature set at bootstrapping suffices to determine the ice softness; thus the 87 options ``-Mz 3 -z_spacing equal``. 88 89 - Option ``-stress_balance ssa`` selects the SSA stress balance and turns off the SIA 90 stress balance computation, since our goal is to model the ice shelf. It also side-steps 91 a technical issue: PISM uses periodic boundary conditions at domain boundaries and most 92 fields in this setup are not periodic. Turning off SIA avoids operations such as 93 differencing surface elevation across the domain edges. For a more complete solution to 94 this technical issue see section :ref:`sec-jako` about a regional model using PISM's 95 "regional mode" ``pismr -regional`` and the option :opt:`-no_model_strip`. 96 97 - Option ``-y 1.0 -no_mass -energy none`` chooses a "diagnostic" run: in absence of 98 geometry evolution and stability restrictions of the energy balance model a 99 one-year-long run will be covered by exactly one time step. 100 101 - Option ``-pik`` is equivalent to ``-cfbc -part_grid -kill_icebergs -subgl`` in this 102 non-evolving example. Note that ``-kill_icebergs`` removes effectively-detached bits of 103 ice, especially in McMurdo sound area, so that the SSA problem is well-posed for the 104 grounded-ice-sheet-connected ice shelf. 105 106 - Option :opt:`-ssa_dirichlet_bc` forces the use of fields ``u_ssa_bc, v_ssa_bc, bc_mask`` 107 described above. The field ``bc_mask`` is `1` at boundary condition locations, and `0` 108 elsewhere. For the prognostic runs below, the ice thickness is also fixed at boundary 109 condition locations, so as to prescribe ice flux as an ice shelf input. 110 111 - Options ``-yield_stress constant -tauc 1e6`` essentially just turn off the 112 grounded-ice evolving yield stress mechanism, which is inactive anyway, and force a high 113 resistance under grounded ice so it does not slide. 114 115 - Option ``-ssa_e 0.6`` is the single tuned parameter; this value gives good 116 correlation between observed and modeled velocity magnitudes. 117 118 - Option ``-ssafd_ksp_monitor`` provides feedback on the linear solver iterations 119 "underneath" the nonlinear (shear-thinning) SSA solver iteration. 120 121 There is no need to type in the above command; just run 122 123 .. code-block:: none 124 125 cd diagnostic/ 126 ./run_diag.sh 2 211 0.6 127 128 Note ``run_diag.sh`` accepts three arguments: ``run_diag.sh N Mx E`` does a run 129 with ``N`` MPI processes, an ``Mx`` by ``Mx`` grid, and option 130 ``-ssa_e E``. The choices above give a run which only takes a few seconds, and it 131 produces output file ``diag_Mx211.nc``. 132 133 There are many reasonable choices for the effective softness of an ice shelf, as ice 134 density, temperature, and the presence of fractures all influence the effective softness. 135 Using an enhancement factor :opt:`-ssa_e 0.6` acknowledges that the physical justification 136 for tuning the ice softness is uncertain. One could instead use the temperature itself or 137 the ice density\ [#]_ as tuning parameters, and these are worthwhile experiments for the 138 interested PISM user. 139 140 The script ``plot.py`` takes PISM output such as ``diag_Mx211.nc`` to produce 141 :numref:`fig-rosspython`. The run shown in the figure used an enhancement factor of 142 `0.6` as above. The thin black line outlines the floating shelf, which is the actual 143 modeling domain here. To generate this figure yourself, run 144 145 .. code-block:: none 146 147 ../plot.py diag_Mx211.nc 148 149 .. figure:: figures/ross-results.png 150 :name: fig-rosspython 151 152 *Left*: Color is speed in m/a. Arrows are observed (white) and modeled (black) 153 velocities. *Right*: Comparison between modeled and observed speeds at points plotted 154 on the left. 155 156 Extending this example to other ice shelves 157 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ 158 159 The SSA diagnostic solution described in this section can be easily applied to other ice 160 shelves in Antarctica, such as the Filchner-Ronne Ice Shelf modeled using PISM in 161 :cite:`AlbrechtLevermann2012`, for example. 162 163 Simply choose a different rectangular domain, within the area covered by the 164 whole-Antarctic data-sets used here, at the preprocessing stage. In particular you should 165 modify the lines "``ncks -O -d x1,439,649 -d y1,250,460 ...``" (for ALBMAP data) and 166 "``ncks -d x,2200,3700 -d y,3500,4700 ...``" (for MEaSUREs velocity data) in the 167 script ``examples/ross/preprocess.py``. 168 169 Prognostic modelling using eigencalving 170 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ 171 172 Next we summarize how you can create an evolving-geometry model of the Ross ice shelf with 173 constant-in-time inflow across the fixed grounding line. See ``README.md`` and 174 ``run_prog.sh`` in ``examples/ross/prognostic/``. This example also demonstrates the 175 :opt:`-calving eigen_calving` model for a moving calving front :cite:`Levermannetal2012`. 176 177 Start by running ``preprocess.py`` in ``examples/ross/`` as described above. If 178 you have already done the diagnostic example above, then this stage is complete. 179 180 Then change to the ``prognostic/`` directory and run the default example: 181 182 .. code-block:: none 183 184 cd examples/ross/prognostic/ 185 ./run_prog.sh 4 211 0.6 100 186 187 This 100 model year run on 4 processes and a 5 km grid took about forty minutes on a 2016 188 laptop. It starts with a bootstrapping stage which does a ``-y 1.0`` run, which generates 189 ``startfile_Mx211.nc``. It then re-initializes to start the prognostic run itself. See the 190 ``README.md`` for a bit more on the arguments taken by ``run_prog.sh`` and on viewing the 191 output files. 192 193 The PISM command done here is (essentially, and without showing diagnostic output choices) 194 195 .. code-block:: none 196 197 pismr -i startfile_Mx211.nc -surface given -stress_balance ssa \ 198 -yield_stress constant -tauc 1e6 -pik -ssa_dirichlet_bc -ssa_e 0.6 \ 199 -y 100 -o prog_Mx211_yr100.nc -o_size big \ 200 -calving eigen_calving,thickness_calving -eigen_calving_K 1e17 \ 201 -front_retreat_cfl -thickness_calving_threshold 50.0 202 203 Several of these options are different from those used in the diagnostic case. First, 204 while the command ``-pik`` is the same as before, now each part of its expansion, namely 205 ``-cfbc -part_grid -kill_icebergs -subgl``, is important. As the calving front evolves 206 (i.e. regardless of the calving law choices), option ``-part_grid`` moves the calving 207 front by one grid cell only when the cell is full of the ice flowing into it; see 208 :cite:`Albrechtetal2011`. The option ``-kill_icebergs`` is essential to maintain well-posedness 209 of the SSA velocity problem at each time step :cite:`Winkelmannetal2011`. See section 210 :ref:`sec-pism-pik`. 211 212 Option combination 213 214 .. code-block:: none 215 216 -calving eigen_calving,thickness_calving -eigen_calving_K 1e17 \ 217 -front_retreat_cfl -thickness_calving_threshold 50.0 218 219 specifies that ice at the calving front will be removed if either a criterion on the 220 product of principal stresses is satisfied :cite:`Levermannetal2012`, namely ``eigen_calving`` 221 with the given constant `K`, or if the ice thickness goes below the given threshold of 50 222 meters. See section :ref:`sec-calving`. 223 224 .. %FIXME Use evolving fracture density. See ``README.md``, ``preprocess_frac.py``, and 225 ``run_frac.sh`` in directory ``examples/ross/fracture_density/``. This example 226 demonstrates the fracture density transport model in :cite:`AlbrechtLevermann2012`. 227 228 .. rubric:: Footnotes 229 230 .. [#] High accumulation rates, cold firn with minimal compression, and basal freeze-on of 231 marine ice may all generate significant variation in shelf density.