tREADME.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
       ---
       tREADME.rst (7979B)
       ---
            1 .. default-role:: math
            2 .. |prs| replace:: *Phillips et al*
            3 .. |chs| replace:: cryo-hydrologic system
            4 
            5 One-column cryo-hydrologic warming setup
            6 ========================================
            7 
            8 This directory contains a set up of the one-column cryo-hydrologic warming experiment
            9 inspired by **Phillips, T. and Rajaram, H. and Steffen, K.**, *Cryo-hydrologic warming: A
           10 potential mechanism for rapid thermal response of ice sheets*,
           11 (https://doi.org/10.1029/2010GL044397).
           12 
           13 Summary
           14 -------
           15 
           16 In summary: in ablation areas the water generated during the melt season drains to the
           17 base through the channels in the |chs|. The presence of water
           18 in the |chs| leads to warming of the ice adjacent to the channels.
           19 
           20 The authors assume that the additional heat flux into the ice due to this mechanism is
           21 proportional to the difference in temperature between the |chs| and the ice next to it:
           22 
           23 .. math::
           24 
           25    Q = \frac{k_{i}}{R^2} (T_{CH} - T_i)
           26 
           27 Here `k_i` is the thermal conductivity of ice and `R` is the spacing between the
           28 channels in the |chs| (a poorly-constrained model parameter).
           29 
           30 |prs| model `T_i` using a familiar temperature-based energy balance equation
           31 with a parameterization of cooling due to horizontal advection. The flux `Q` defined above
           32 appears as an additional source term:
           33 
           34 .. math::
           35 
           36    \frac{\partial \rho_i C_i T_i}{\partial t} + \frac{\partial (\rho_i C_i w T_i)}{\partial z}
           37    - \frac{\partial}{\partial z}\left(k_i \frac{\partial T_i}{\partial z}\right)
           38    = - \frac{\partial \rho_i C_i u T_i}{\partial x} + Q_{s} + Q
           39 
           40 Here `\rho_i` is the density of ice, `u` and `w` are horizontal and vertical velocity
           41 components, `C_i` is the specific heat capacity of ice, and `Q_s` is the strain heating.
           42 
           43 The temperature in the |chs| is modeled using an enthalpy-based energy equation omitting
           44 advection terms:
           45 
           46 .. math::
           47 
           48    \frac{\partial \bar{\rho H}}{\partial t}
           49    - \frac{\partial}{\partial z}\left(\bar{k} \frac{\partial T_{CH}}{\partial z}\right)
           50    = - Q
           51 
           52 Here the enthalpy is defined by
           53 
           54 .. math::
           55    :name: enthalpy
           56 
           57    \rho H = (1 - \phi_w) \rho_i C_i T_{CH} + \phi_w (\rho_w C_w T_{CH} + L),
           58 
           59 where `\phi_w` is the water fraction, `L` is the latent heat of fusion, `\rho_w`, is the
           60 water density, and `C_w` is the specific heat capacity of fresh water.
           61 
           62 The thermal conductivity of the mixture, `\bar k`, is defined by
           63 
           64 .. math::
           65 
           66    \bar k = (1 - \phi_w) k_i + \phi_w k_w.
           67 
           68 The authors use Dirichlet (temperature) boundary conditions at the top boundary (assuming
           69 that the temperature at the top of the snow or ice is equal to the near-surface air
           70 temperature, but never exceeds the pressure melting point) and Dirichlet or Neumann
           71 (mentioned but not specified in the paper) at the bottom boundary.
           72 
           73 It is also mentioned that they account for *snow cover, which insulates the ice surface
           74 from cold winter temperatures* and *temperature dependence of the thermal properties was
           75 incorporated in the model*, but no details are provided.
           76 
           77 During the melt season (when the water is present) the temperature in the |chs| is set to
           78 the pressure-melting and the water fraction is set to a predetermined constant (|prs| use
           79 `0.005`, i.e. one half of a percent) while during the winter the |chs| is allowed to cool.
           80 
           81 Comparing to PISM's energy conservation model
           82 ---------------------------------------------
           83 
           84 - PISM's definition of enthalpy is equivalent to the one above, with some caveats. (The
           85   equation in |prs| has a couple of issues that I don't need to go into here.)
           86 
           87   Note that in PISM's definition the temperature never exceeds `T_m(p)`, the
           88   pressure-melting point, while in |prs| `T_{CH} > T_m(p)` is allowed, but `T_{CH}` never
           89   exceeds `T_m(p)` due to the setup they use.
           90 
           91 - PISM's energy conservation model includes all the mechanisms modeled by equations above
           92   *except* for the parameterization of advection-driven cooling not necessary in a 3D model.
           93 
           94 - PISM's enthalpy-based energy conservation model can be used to model *both* the ice and
           95   the |chs| columns.
           96 
           97 - PISM assumes that the specific heat capacity of ice is constant. The thermal
           98   conductivity can be a function of temperature, but does not depend on the water fraction
           99   (unlike `\bar k` above). Note, though, that for water fractions at or below `0.005` the
          100   value of `\bar k` is very close to `k_i`.
          101 
          102 With default parameter values PISM "ignores" thermal conductivity of temperate ice
          103 *without* completely removing it: this can be thought of as a *regularization*.
          104 Specifically, the thermal conductivity of temperate ice is set to a given fraction of
          105 `k_i` (see ``energy.enthalpy.temperate_ice_thermal_conductivity_ratio``).
          106 
          107 This one-column setup uses ``energy.ch_warming.temperate_ice_thermal_conductivity_ratio``
          108 (equal to `1`) to get us close to the equation modeling |chs| enthalpy evolution in
          109 *Phillips et al*.
          110 
          111 Technical details
          112 -----------------
          113 
          114 This setup models a 200-meter deep ice column using an equally-spaced vertical grid with a
          115 1-meter resolution (201 grid points) and takes equal time steps 1 day long each.
          116 
          117 The flux `Q` from the |chs| into the ice is added to the deformational strain heating term
          118 in PISM's energy balance equation.
          119 
          120 Ice velocity components are set to zero, eliminating advection terms and resulting in zero
          121 strain heating due to ice flow.
          122 
          123 We use a Dirichlet boundary condition at the top surface
          124 
          125 .. math::
          126 
          127    \begin{aligned}
          128      H_{\text{surface}} &=
          129      \begin{cases}
          130      H(T_{\text{surface}}), & T_{\text{surface}} < T_m, \\
          131      H(T_m), & T_{\text{surface}} \ge T_m,
          132      \end{cases}\\
          133      T_{\text{surface}} &= T_{\text{mean}} + A\cdot \cos(2 \pi t - \phi_0)
          134    \end{aligned}
          135 
          136 with time `t` in years. The values of `A` (`6 K`) and `T_{\text{mean}}` (`-5^{\circ} C`)
          137 are selected to get a 9-week-long melt season.
          138 
          139 We use a Neumann boundary condition at the bottom surface:
          140 
          141 .. math::
          142 
          143    \left.\frac{\partial H}{\partial z}\right|_{z=0} = - \frac{G}{k_i},
          144 
          145 where `G` is the upward geothermal heat flux, `W / m^2`. This setup uses `G = 0`.
          146 
          147 Given the absence of ice flow, we use the steady-state enthalpy distribution corresponding
          148 to the mean-annual temperature at the surface and the constant heat flux at the base of
          149 the ice:
          150 
          151 .. math::
          152 
          153    \begin{aligned}
          154      H_{\text{steady state}} &= H(T_{\text{steady state}}, p)\\
          155      T_{\text{steady state}} &= (h - z) \frac{G}{k_i}\\
          156    \end{aligned}
          157 
          158 where `h` is ice thickness and `p` is the hydrostatic pressure. We assume that the mean
          159 annual air temperature and ice thickness are low enough so that this does not result in a
          160 layer of temperate ice near the base.
          161 
          162 At the beginning of each time step the enthalpy in the |chs| is set to `H(T_m(p),
          163 \omega_0)` if `T_{\text{surface}} \ge T_m` (here `\omega_0 = 0.005`), simulating the
          164 presence of liquid water in the |chs| during the melt season.
          165 
          166 .. note::
          167 
          168    One could use the surface mass balance to detect the melt season instead:
          169    `T_{\text{surface}}` may not capture the daily temperature variability.
          170 
          171 Results
          172 -------
          173 
          174 This annual air temperature cycle results in a 9-week-long melt season.
          175 
          176 .. figure:: air-temperature.png
          177 
          178    Air temperature
          179 
          180 .. figure:: ice-temperature.png
          181 
          182    Ice temperature in the top 15 meters of ice. Note the gradual warming with almost no
          183    inter-annual variability at the depth of 15 meters.
          184 
          185 .. figure:: ice-temperature-curves.png
          186 
          187    Ice temperature evolution at different depths.
          188 
          189 .. figure:: ch-temperature.png
          190 
          191    Temperature in the |chs|. The end of the melt season can be seen
          192    clearly.
          193 
          194 .. figure:: ch-water-fraction.png
          195 
          196    Water fraction in the |chs|. Note that at higher depths the water fraction stays above
          197    zero all the way through the winter, so the temperature in the |chs| never drops below
          198    pressure-melting.
          199 
          200 Next steps
          201 ----------
          202 
          203 - Adjust the value of `R` and observe the changes in the evolution of ice temperature.
          204 - Compare to **Hills et al**, *Processes influencing near-surface heat transfer in
          205   Greenland’s ablation zone*, (https://doi.org/10.5194/tc-2018-51).