PISM, A Parallel Ice Sheet Model
stable v2.1-1-g6902d5502 committed by Ed Bueler on 2023-12-20 08:38:27 -0800
|
◆ compute_diffusivity()
Compute the SIA diffusivity. If full_update, also store delta on the staggered grid. Recall that \( Q = -D \nabla h \) is the diffusive flux in the mass-continuity equation \[ \frac{\partial H}{\partial t} = M - S - \nabla \cdot (Q + \mathbf{U}_b H),\] where \(h\) is the ice surface elevation, \(M\) is the top surface accumulation/ablation rate, \(S\) is the basal mass balance and \(\mathbf{U}_b\) is the thickness-advective (in PISM: usually SSA) ice velocity. Recall also that at any particular point in the map-plane (i.e. if \(x\) and \(y\) are fixed) \[ D = 2\int_b^h F(z)P(z)(h-z)dz, \] where \(F(z)\) is a constitutive function and \(P(z)\) is the pressure at a level \(z\). By defining \[ \delta(z) = 2F(z)P(z) \] one can write \[D = \int_b^h\delta(z)(h-z)dz. \] The advantage is that it is then possible to avoid re-evaluating \(F(z)\) (which is computationally expensive) in the horizontal ice velocity (see compute_3d_horizontal_velocity()) computation. This method computes \(D\) and stores \(\delta\) in delta[0,1] if full_update is true. The trapezoidal rule is used to approximate the integral.
Definition at line 529 of file SIAFD.cc. References pism::Geometry::cell_type, pism::ParallelSection::check(), pism::Time::current(), pism::ParallelSection::failed(), pism::RuntimeError::formatted(), pism::array::Array3D::get_column(), pism::GlobalMax(), pism::GlobalSum(), pism::Geometry::ice_surface_elevation, pism::Geometry::ice_thickness, interglacial(), pism::k, m_bed_smoother, pism::Component::m_config, pism::stressbalance::SSB_Modifier::m_D_max, m_delta_0, m_delta_1, m_e_factor, m_e_factor_interglacial, pism::stressbalance::SSB_Modifier::m_EC, pism::stressbalance::SSB_Modifier::m_flow_law, pism::Component::m_grid, pism::Component::m_log, m_seconds_per_year, m_work_2d_0, m_work_2d_1, pism::array::max(), PISM_ERROR_LOCATION, pism::array::Array::set(), pism::array::Array3D::set_column(), pism::stressbalance::BedSmoother::smoothed_thk(), pism::array::Array::stencil_width(), pism::stressbalance::BedSmoother::theta(), pism::Component::time(), pism::grid::X_PERIODIC, and pism::grid::Y_PERIODIC. Referenced by update(). |