PISM, A Parallel Ice Sheet Model  stable v2.0.6 committed by Constantine Khrulev on 2023-01-23 15:14:38 -0900

## ◆ compute_vertical_velocity()

 void pism::stressbalance::StressBalance::compute_vertical_velocity ( const IceModelVec2CellType & mask, const IceModelVec3 & u, const IceModelVec3 & v, const IceModelVec2S * basal_melt_rate, IceModelVec3 & result )
protectedvirtual

Compute vertical velocity using incompressibility of the ice.

The vertical velocity $$w(x,y,z,t)$$ is the velocity relative to the location of the base of the ice column. That is, the vertical velocity computed here is identified as $$\tilde w(x,y,s,t)$$ in the page []vertchange.

Thus $$w<0$$ here means that that that part of the ice is getting closer to the base of the ice, and so on. The slope of the bed (i.e. relative to the geoid) and/or the motion of the bed (i.e. from bed deformation) do not affect the vertical velocity.

In fact the following statement is exactly true if the basal melt rate is zero: the vertical velocity at a point in the ice is positive (negative) if and only if the average horizontal divergence of the horizontal velocity, in the portion of the ice column below that point, is negative (positive). In particular, because $$z=0$$ is the location of the base of the ice always, the only way to have $$w(x,y,0,t) \ne 0$$ is to have a basal melt rate.

Incompressibility itself says

$\nabla\cdot\mathbf{U} + \frac{\partial w}{\partial z} = 0.$

This is immediately equivalent to the integral

$w(x,y,z,t) = - \int_{b(x,y,t)}^{z} \nabla\cdot\mathbf{U}\,d\zeta + w_b(x,y,t).$

Here the value $$w_b(x,y,t)$$ is either zero or the negative of the basal melt rate according to the value of the flag geometry.update.use_basal_melt_rate.

The vertical integral is computed by the trapezoid rule.

Definition at line 280 of file StressBalance.cc.

Referenced by update().