tssafd-cfbc.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
---
tssafd-cfbc.rst (9030B)
---
1 .. include:: ../global.txt
2
3 Calving front stress boundary condition
4 =======================================
5
6 .. contents::
7
8 .. only:: html
9
10 .. math::
11
12 \newcommand{\diff}[2]{\frac{\partial #1}{\partial #2}}
13 \newcommand{\n}{\mathbf{n}}
14 \newcommand{\nx}{\n_{x}}
15 \newcommand{\ny}{\n_{y}}
16 \newcommand{\nz}{\n_{z}}
17 \newcommand{\psw}{p_{\text{ocean}}}
18 \newcommand{\pmelange}{p_{\text{melange}}}
19 \newcommand{\pice}{p_{\text{ice}}}
20 \newcommand{\rhoi}{\rho_{\text{ice}}}
21 \newcommand{\rhosw}{\rho_{\text{ocean}}}
22 \newcommand{\zs}{z_{\text{s}}}
23 \newcommand{\td}[1]{t^{D}_{#1}}
24 \newcommand{\D}{\displaystyle}
25 \newcommand{\dx}{\Delta x}
26 \newcommand{\dy}{\Delta y}
27 \newcommand{\viscosity}{\nu}
28 \newcommand{\partI}{(2\tilde N_{xx} + \tilde N_{yy})}
29 \newcommand{\partII}{(\tilde N_{xy})}
30
31 .. raw:: latex
32
33 \providecommand{\diff}[2]{\frac{\partial #1}{\partial #2}}
34 \providecommand{\n}{\mathbf{n}}
35 \providecommand{\nx}{\n_{x}}
36 \providecommand{\ny}{\n_{y}}
37 \providecommand{\nz}{\n_{z}}
38 \providecommand{\psw}{p_{\text{ocean}}}
39 \providecommand{\pmelange}{p_{\text{melange}}}
40 \providecommand{\pice}{p_{\text{ice}}}
41 \providecommand{\rhoi}{\rho_{\text{ice}}}
42 \providecommand{\rhosw}{\rho_{\text{ocean}}}
43 \providecommand{\zs}{z_{\text{s}}}
44 \providecommand{\td}[1]{t^{D}_{#1}}
45 \providecommand{\D}{\displaystyle}
46 \providecommand{\dx}{\Delta x}
47 \providecommand{\dy}{\Delta y}
48 \providecommand{\viscosity}{\nu}
49 \providecommand{\partI}{(2\tilde N_{xx} + \tilde N_{yy})}
50 \providecommand{\partII}{(\tilde N_{xy})}
51
52 Notation
53 --------
54
55 .. csv-table::
56 :name: tab-cfbc-notation
57 :header: Variable, Meaning
58
59 `h`, ice top surface elevation
60 `b`, ice bottom surface elevation
61 `H = h - b`, ice thickness
62 `\zs`, sea level elevation
63 `g`, acceleration due to gravity
64 `\rhoi`, ice density
65 `\rhosw`, sea water density
66 `\viscosity`, vertically-averaged viscosity of ice
67 `\n`, normal vector
68 `B(T)`, ice hardness
69 `\pice`, pressure due to the weight of a column of ice
70 `\psw`, pressure due to the weight of a column of seawater
71 `D`, strain rate tensor
72 `d_{e}`, effective strain rate
73 `t`, Cauchy stress tensor
74 `t^{D}`, deviatoric stress tensor; note `\td{ij} = t_{ij} + p \delta_{ij}`
75
76 Calving front stress boundary condition
77 ---------------------------------------
78
79 In the 3D case the calving front stress boundary condition (:cite:`GreveBlatter2009`, equation
80 (6.19)) reads
81
82 .. math::
83
84 \left.t\right|_{\text{cf}} \cdot \n = -\psw \n.
85
86 Expanded in component form, and evaluating the left-hand side at the calving front and
87 assuming that the calving front face is vertical (`\nz = 0`), this gives
88
89 .. math::
90
91 \begin{array}{rcrcl}
92 (\td{xx} - p) \nx &+& \td{xy} \ny &=& -\psw \nx,\\
93 \td{xy} \nx &+& (\td{yy} - p) \ny &=& -\psw \ny,\\
94 \td{xz} \nx &+& \td{yz} \ny &=& 0.
95 \end{array}
96
97
98 Because we seek boundary conditions for the SSA stress balance, in which the
99 vertically-integrated forms of the stresses `\td{xx},\td{xy},\td{yy}` are balanced
100 separately from the shear stresses `\td{xz},\td{yz}`, the third of the above equations can
101 be omitted from the remaining analysis.
102
103 Let `\pice=\rhoi g (h-z)`. In the hydrostatic approximation, `t_{zz}=-\pice`
104 (:cite:`GreveBlatter2009`, equation (5.59)). Next, since `\td{}` has trace zero,
105
106 .. math::
107
108 p &= p - \td{xx} - \td{yy} - \td{zz}
109
110 &= - t_{zz} - \td{xx} - \td{yy}
111
112 &= \pice - \td{xx} - \td{yy}.
113
114 Thus
115
116 .. math::
117 :label: ssafd-cfbc-1
118
119 \begin{array}{rcrcl}
120 (2\td{xx} + \td{yy}) \nx &+& \td{xy} \ny &=& (\pice - \psw) \nx,\\
121 \td{xy} \nx &+& (2\td{yy} + \td{xx}) \ny &=& (\pice - \psw) \ny.\\
122 \end{array}
123
124 Now, using the "viscosity form" of the flow law
125
126 .. math::
127
128 \td{} &= 2\viscosity\, D,
129
130 \viscosity &= \frac 12 B(T) d_{e}^{1/n-1}
131
132 and the fact that in the shallow shelf approximation `u` and `v` are depth-independent, we
133 can define the membrane stress `N` as the vertically-integrated deviatoric stress
134
135 .. math::
136
137 N_{ij} = \int_{b}^{h} t^{D}_{ij} dz = 2\, \bar \viscosity\, H\, D_{ij}.
138
139 Here `\bar \viscosity` is the vertically-averaged effective viscosity.
140
141 Integrating :eq:`ssafd-cfbc-1` with respect to `z`, we get
142
143 .. math::
144 :label: ssafd-cfbc-3
145
146 \begin{array}{rcrcl}
147 (2N_{xx} + N_{yy}) \nx &+& N_{xy} \ny &=& \displaystyle \int_{b}^{h}(\pice - \psw) dz\, \nx,\\
148 N_{xy} \nx &+& (2N_{yy} + N_{xx}) \ny &=& \displaystyle \int_{b}^{h}(\pice - \psw) dz\, \ny.
149 \end{array}
150
151
152 Shallow shelf approximation
153 ---------------------------
154
155 The shallow shelf approximation written in terms of `N_{ij}` is
156
157 .. math::
158 :label: ssafd-cfbc-2
159
160 \begin{array}{rcrcl}
161 \D \diff{}{x}(2N_{xx} + N_{yy}) &+& \D \diff{N_{xy}}{y} &=& \D \rho g H \diff{h}{x},\\
162 \D \diff{N_{xy}}{x} &+& \D \diff{}{y}(2N_{yy} + N_{xx}) &=& \D \rho g H \diff{h}{y}.
163 \end{array}
164
165 Implementing the calving front boundary condition
166 -------------------------------------------------
167
168 We use centered finite difference approximations in the discretization of the SSA
169 :eq:`ssafd-cfbc-2`. Consider the first equation:
170
171 .. math::
172 :label: ssafd-cfbc-4
173
174 \diff{}{x}(2N_{xx} + N_{yy}) + \D \diff{N_{xy}}{y} = \D \rho g H \diff{h}{x}.
175
176 Let `\tilde F` be an approximation of `F` using a finite-difference scheme. Then the first
177 SSA equation is approximated by
178
179 .. math::
180
181 \frac1{\dx}\left(\partI_{i+\frac12,j} - \partI_{i-\frac12,j}\right) +
182 \frac1{\dy}\left(\partII_{i,j+\frac12} - \partII_{i,j-\frac12}\right) =
183 \rho g H \frac{h_{i+\frac12,j} - h_{i-\frac12,j}}{\dx}.
184
185
186 Now, assume that the cell boundary (face) at `i+\frac12,j` is at the
187 calving front. Then `\n = (1,0)` and from :eq:`ssafd-cfbc-3` we have
188
189 .. math::
190 :label: ssafd-cfbc-vertintbdry
191
192 2N_{xx} + N_{yy} = \int_{b}^{h}(\pice - \psw) dz.
193
194 We call the right-hand side of :eq:`ssafd-cfbc-vertintbdry` the "pressure difference term."
195
196 In forming the matrix approximation of the SSA :cite:`BBssasliding`,
197 :cite:`Winkelmannetal2011`, instead of assembling a matrix row corresponding to the
198 generic equation :eq:`ssafd-cfbc-4` we use
199
200 .. math::
201
202 \frac1{\dx}\left(\left[\int_{b}^{h}(\pice - \psw) dz\right]_{i+\frac12,j} - \partI_{i-\frac12,j}\right) +
203 \frac1{\dy}\left(\partII_{i,j+\frac12} - \partII_{i,j-\frac12}\right) =
204 \rho g H \frac{h_{i+\frac12,j} - h_{i-\frac12,j}}{\dx}.
205
206 Moving terms that do not depend on the velocity field to the
207 right-hand side, we get
208
209 .. math::
210
211 \frac1{\dx}\left(- \partI_{i-\frac12,j}\right) +
212 \frac1{\dy}\left(\partII_{i,j+\frac12} - \partII_{i,j-\frac12}\right) =
213 \rho g H \frac{h_{i+\frac12,j} - h_{i-\frac12,j}}{\dx} + \left[\frac{\int_{b}^{h}(\psw - \pice) dz}{\dx}\right]_{i+\frac12,j}.
214
215
216 The second equation and other cases (`\n = (-1,0)`, etc) are treated similarly.
217
218 Evaluating the "pressure difference term"
219 -----------------------------------------
220
221 For `z \in [b, h]` the modeled pressures are
222
223 .. math::
224
225 \begin{array}{ll}
226 \psw &=
227 \begin{cases}
228 0, & z > \zs,\\
229 \rhosw\, g (\zs - z), & z \le \zs,
230 \end{cases}\\
231 \pice &= \rhoi\, g (h - z).
232 \end{array}
233
234 Depending on the local geometry `b` is either prescribed (grounded case) or is a function
235 of the ice thickness and sea level elevation.
236
237 Floating case
238 ^^^^^^^^^^^^^
239
240 Using the flotation thickness relation `\rhoi H = \rhosw (\zs - b)`, which applies to
241 floating ice given that `b` here denotes the base elevation of the floating ice, we have
242
243 .. math::
244
245 \int_{b}^{h}(\pice - \psw) dz =
246 \frac{1}{2}\, g\, \rhoi \left(H^2 - \frac{\rhoi}{\rhosw} H^2\right)
247
248 Grounded below sea level
249 ^^^^^^^^^^^^^^^^^^^^^^^^
250
251 Because `b` denotes both the base elevation of the grounded ice, and the bedrock
252 elevation, here `\rhoi H \ge \rhosw (\zs - b)`. The integral simplifies to
253
254 .. math::
255
256 \int_{b}^{h}(\pice - \psw) dz =
257 \frac{1}{2}\, g\, \rhoi \left(H^2-\frac{\rhosw}{\rhoi}\, \left(\zs-b\right)^2\right)
258
259
260 This, in fact, applies in both floating and grounded-below-sea-level cases.
261
262 Grounded above sea level
263 ^^^^^^^^^^^^^^^^^^^^^^^^
264
265 In this case `\psw = 0`, so
266
267 .. math::
268
269 \int_{b}^{h}(\pice - \psw) dz = \frac{1}{2}\, g\, \rhoi\, H^2
270
271
272 Modeling melange back-pressure
273 ------------------------------
274
275 Let `\pmelange` be the additional melange back-pressure. Then `\psw \le \psw + \pmelange
276 \le \pice`. Put another way,
277
278 .. math::
279 :label: ssafd-cfbc-13
280
281 0 \le \pmelange \le \pice - \psw.
282
283 Let `\lambda` be the "melange back-pressure fraction" (or "relative melange pressure")
284 ranging from `0` to `1`, so that
285
286 .. math::
287 :label: ssafd-cfbc-14
288
289 \pmelange = \lambda \cdot (\pice - \psw).
290
291 Then the modified pressure difference term is
292
293 .. math::
294 :label: ssafd-cfbc-15
295
296 \int_{b}^{h}(\pice - (\psw + \pmelange)) dz &= \int_{b}^{h}(\pice - (\psw + \lambda(\pice - \psw)))\, dz
297
298 &= (1 - \lambda) \int_{b}^{h} (\pice - \psw)\, dz.