Loading [MathJax]/extensions/tex2jax.js
PISM, A Parallel Ice Sheet Model 2.2.1-cd005eec8 committed by Constantine Khrulev on 2025-03-07
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
BedSmoother.hh
Go to the documentation of this file.
1// Copyright (C) 2010, 2011, 2013, 2014, 2015, 2016, 2017, 2020, 2021, 2022 Ed Bueler and Constantine Khroulev
2//
3// This file is part of PISM.
4//
5// PISM is free software; you can redistribute it and/or modify it under the
6// terms of the GNU General Public License as published by the Free Software
7// Foundation; either version 3 of the License, or (at your option) any later
8// version.
9//
10// PISM is distributed in the hope that it will be useful, but WITHOUT ANY
11// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
12// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
13// details.
14//
15// You should have received a copy of the GNU General Public License
16// along with PISM; if not, write to the Free Software
17// Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
18
19#ifndef __BedSmoother_hh
20#define __BedSmoother_hh
21
22#include <petsc.h>
23
24#include "pism/util/array/Scalar.hh"
25#include "pism/util/ConfigInterface.hh"
26
27namespace pism {
28
29class Grid;
30class Config;
31
32namespace array {
33class CellType1;
34class CellType2;
35class CellType;
36} // end of namespace array
37
38namespace stressbalance {
39
40//! PISM bed smoother, plus bed roughness parameterization, based on Schoof (2003).
41/*!
42 This class both smooths the bed and computes coefficients for an approximation
43 to Schoof's \f$\theta\f$. The factor \f$0\le \theta \le 1\f$ multiplies the diffusivity
44 in the theory of the effect of bed roughness in the SIA by Christian Schoof
45 (2003; *The effect of basal topography on ice sheet dynamics*) [\ref
46 Schoofbasaltopg2003].
47
48 For additional information on this class see page \ref bedrough.
49
50 The user of this class hands BedSmoother an "original" topography, and it
51 is preprocessed to fill the smoothed topography `topgsmooth`, and the
52 coefficients in an approximation to \f$\theta\f$. This is done by a call to
53 `preprocess_bed()`. The call requires the half-width of the smoothing square
54 (a distance in m), or the number of grid points in each direction in the
55 smoothing rectangle, and the Glen exponent.
56
57 The call to `preprocess_bed()` *must be repeated* any time the "original"
58 topography changes, for instance at the start of an IceModel run, or at a bed
59 deformation step in an IceModel run.
60
61 BedSmoother then provides three major functionalities, all of which \e must
62 \e follow the call to `preprocess_bed()`:
63 -# User accesses public array::Scalar `topgsmooth`, the smoothed bed itself.
64 -# User asks `get_smoothed_thk()` for gridded values of the consistent smoothed
65 version of the ice thickness, which is the thickness corresponding to a given
66 surface elevation and the pre-computed smoothed bed.
67 -# User asks for gridded values of \f$\theta(h,x,y)\f$ using `get_theta()`.
68
69 Here is a basic example of the creation and usage of a BedSmoother instance.
70 Note that BedSmoother will update ghosted values in `topgsmooth`, and other
71 internal fields, and update them in the return fields `thksmooth`, `theta`,
72 if asked. In IceModel::velocitySIAStaggered()
73 \code
74 BedSmoother smoother(grid, 2);
75 const double n = 3.0,
76 lambda = 50.0e3;
77 ierr = smoother.preprocess_bed(topg, n, lambda); CHKERRQ(ierr);
78 ierr = smoother.get_smoothed_thk(usurf, thk, 1, &thksmooth); CHKERRQ(ierr);
79 ierr = smoother.get_theta(usurf, n, 1, &theta); CHKERRQ(ierr);
80 \endcode
81 See Grid documentation for initializing `grid`. Note we assume
82 `topg`, `usurf`, `thk`, `thksmooth`, and `theta` are all created
83 array::Scalar instances.
84*/
86public:
87 BedSmoother(std::shared_ptr<const Grid> g);
88 virtual ~BedSmoother() = default;
89
90 void preprocess_bed(const array::Scalar &topg);
91
92 void smoothed_thk(const array::Scalar &usurf,
93 const array::Scalar &thk,
94 const array::CellType2 &mask,
95 array::Scalar &thksmooth) const;
96
97 void theta(const array::Scalar &usurf, array::Scalar &result) const;
98
99 const array::Scalar& smoothed_bed() const;
100protected:
101 std::shared_ptr<const Grid> m_grid;
103
104 //! smoothed bed elevation; set by calling preprocess_bed()
106
108
109 /* number of grid points to smooth over; e.g. i=-Nx,-Nx+1,...,-1,0,1,...,Nx-1,Nx; note
110 Nx>=1 and Ny>=1 always, unless lambda<=0
111 */
112 int m_Nx, m_Ny;
113
115
116 //! original bed elevation on processor 0
117 std::shared_ptr<petsc::Vec> m_topgp0;
118 //! smoothed bed elevation on processor 0
119 std::shared_ptr<petsc::Vec> m_topgsmoothp0;
120 //! maximum elevation at (i,j) of local topography (nearby patch)
121 std::shared_ptr<petsc::Vec> m_maxtlp0, m_C2p0, m_C3p0, m_C4p0;
122
123 void preprocess_bed(const array::Scalar &topg,
124 unsigned int Nx_in, unsigned int Ny_in);
125
128};
129
130} // end of namespace stressbalance
131} // end of namespace pism
132
133#endif // __BedSmoother_hh
134
std::shared_ptr< const Config > ConstPtr
std::shared_ptr< petsc::Vec > m_C3p0
void smoothed_thk(const array::Scalar &usurf, const array::Scalar &thk, const array::CellType2 &mask, array::Scalar &thksmooth) const
Computes a smoothed thickness map.
std::shared_ptr< petsc::Vec > m_topgp0
original bed elevation on processor 0
array::Scalar2 m_topgsmooth
smoothed bed elevation; set by calling preprocess_bed()
std::shared_ptr< petsc::Vec > m_maxtlp0
maximum elevation at (i,j) of local topography (nearby patch)
void smooth_the_bed_on_proc0()
Computes the smoothed bed by a simple average over a rectangle of grid points.
const array::Scalar & smoothed_bed() const
void theta(const array::Scalar &usurf, array::Scalar &result) const
const Config::ConstPtr m_config
std::shared_ptr< petsc::Vec > m_topgsmoothp0
smoothed bed elevation on processor 0
std::shared_ptr< petsc::Vec > m_C4p0
std::shared_ptr< const Grid > m_grid
void preprocess_bed(const array::Scalar &topg)
std::shared_ptr< petsc::Vec > m_C2p0
PISM bed smoother, plus bed roughness parameterization, based on Schoof (2003).
static const double g
Definition exactTestP.cc:36