PISM, A Parallel Ice Sheet Model  stable v2.1-1-g6902d5502 committed by Ed Bueler on 2023-12-20 08:38:27 -0800

◆ compute_surface_gradient()

void pism::stressbalance::SIAFD_Regional::compute_surface_gradient ( const Inputs inputs,
array::Staggered1 h_x,
array::Staggered1 h_y 
)
privatevirtual

Compute the ice surface gradient for the SIA.

There are three methods for computing the surface gradient. Which method is controlled by configuration parameter sia.surface_gradient_method which can have values haseloff, mahaffy, or eta.

The most traditional method is to directly differentiate the surface elevation \(h\) by the Mahaffy method [Mahaffy]. The haseloff method, suggested by Marianne Haseloff, modifies the Mahaffy method only where ice-free adjacent bedrock points are above the ice surface, and in those cases the returned gradient component is zero.

The alternative method, when sia.surface_gradient_method = eta, transforms the thickness to something more regular and differentiates that. We get back to the gradient of the surface by applying the chain rule. In particular, as shown in [Vazquezetal2003] for the flat bed and \(n=3\) case, if we define

\[\eta = H^{(2n+2)/n}\]

then \(\eta\) is more regular near the margin than \(H\). So we compute the surface gradient by

\[\nabla h = \frac{n}{(2n+2)} \eta^{(-n-2)/(2n+2)} \nabla \eta + \nabla b,\]

recalling that \(h = H + b\). This method is only applied when \(\eta > 0\) at a given point; otherwise \(\nabla h = \nabla b\).

In all cases we are computing the gradient by finite differences onto a staggered grid. In the method with \(\eta\) we apply centered differences using (roughly) the same method for \(\eta\) and \(b\) that applies directly to the surface elevation \(h\) in the mahaffy and haseloff methods.

Parameters
[out]h_xthe X-component of the surface gradient, on the staggered grid
[out]h_ythe Y-component of the surface gradient, on the staggered grid

Reimplemented from pism::stressbalance::SIAFD.

Definition at line 42 of file SIAFD_Regional.cc.

References pism::array::Array2D< T >::box(), pism::Geometry::cell_type, pism::stressbalance::SIAFD::compute_surface_gradient(), pism::stressbalance::Inputs::geometry, pism::Component::m_grid, m_h_x_no_model, m_h_y_no_model, pism::stressbalance::Inputs::no_model_mask, pism::stressbalance::Inputs::no_model_surface_elevation, and pism::stressbalance::SIAFD::surface_gradient_haseloff().