tBedSmoother.hh - pism - [fork] customized build of PISM, the parallel ice sheet model (tillflux branch)
HTML git clone git://src.adamsgaard.dk/pism
DIR Log
DIR Files
DIR Refs
DIR LICENSE
---
tBedSmoother.hh (4958B)
---
1 // Copyright (C) 2010, 2011, 2013, 2014, 2015, 2016, 2017 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/iceModelVec.hh"
25 #include "pism/util/ConfigInterface.hh"
26
27 namespace pism {
28
29 class IceGrid;
30 class Config;
31 class IceModelVec2CellType;
32
33 namespace stressbalance {
34
35 //! PISM bed smoother, plus bed roughness parameterization, based on Schoof (2003).
36 /*!
37 This class both smooths the bed and computes coefficients for an approximation
38 to Schoof's \f$\theta\f$. The factor \f$0\le \theta \le 1\f$ multiplies the diffusivity
39 in the theory of the effect of bed roughness in the SIA by Christian Schoof
40 (2003; *The effect of basal topography on ice sheet dynamics*) [\ref
41 Schoofbasaltopg2003].
42
43 For additional information on this class see page \ref bedrough.
44
45 The user of this class hands BedSmoother an "original" topography, and it
46 is preprocessed to fill the smoothed topography `topgsmooth`, and the
47 coefficients in an approximation to \f$\theta\f$. This is done by a call to
48 `preprocess_bed()`. The call requires the half-width of the smoothing square
49 (a distance in m), or the number of grid points in each direction in the
50 smoothing rectangle, and the Glen exponent.
51
52 The call to `preprocess_bed()` <b>must be repeated</b> any time the "original"
53 topography changes, for instance at the start of an IceModel run, or at a bed
54 deformation step in an IceModel run.
55
56 BedSmoother then provides three major functionalities, all of which \e must
57 \e follow the call to `preprocess_bed()`:
58 -# User accesses public IceModelVec2S `topgsmooth`, the smoothed bed itself.
59 -# User asks `get_smoothed_thk()` for gridded values of the consistent smoothed
60 version of the ice thickness, which is the thickness corresponding to a given
61 surface elevation and the pre-computed smoothed bed.
62 -# User asks for gridded values of \f$\theta(h,x,y)\f$ using `get_theta()`.
63
64 Here is a basic example of the creation and usage of a BedSmoother instance.
65 Note that BedSmoother will update ghosted values in `topgsmooth`, and other
66 internal fields, and update them in the return fields `thksmooth`, `theta`,
67 if asked. In IceModel::velocitySIAStaggered()
68 \code
69 BedSmoother smoother(grid, 2);
70 const double n = 3.0,
71 lambda = 50.0e3;
72 ierr = smoother.preprocess_bed(topg, n, lambda); CHKERRQ(ierr);
73 ierr = smoother.get_smoothed_thk(usurf, thk, 1, &thksmooth); CHKERRQ(ierr);
74 ierr = smoother.get_theta(usurf, n, 1, &theta); CHKERRQ(ierr);
75 \endcode
76 See IceGrid documentation for initializing `grid`. Note we assume
77 `topg`, `usurf`, `thk`, `thksmooth`, and `theta` are all created
78 IceModelVec2S instances.
79 */
80 class BedSmoother {
81 public:
82 BedSmoother(IceGrid::ConstPtr g, int MAX_GHOSTS);
83 virtual ~BedSmoother();
84
85 virtual void preprocess_bed(const IceModelVec2S &topg);
86
87 virtual void smoothed_thk(const IceModelVec2S &usurf,
88 const IceModelVec2S &thk,
89 const IceModelVec2CellType &mask,
90 IceModelVec2S &thksmooth) const;
91
92 virtual void theta(const IceModelVec2S &usurf, IceModelVec2S &result) const;
93
94 const IceModelVec2S& smoothed_bed() const;
95 protected:
96 //! smoothed bed elevation; set by calling preprocess_bed()
97 IceModelVec2S m_topgsmooth;
98
99 IceGrid::ConstPtr m_grid;
100 const Config::ConstPtr m_config;
101 IceModelVec2S m_maxtl, m_C2, m_C3, m_C4;
102
103 int m_Nx, m_Ny; //!< number of grid points to smooth over; e.g.
104 //!i=-Nx,-Nx+1,...,-1,0,1,...,Nx-1,Nx; note Nx>=1 and Ny>=1
105 //!always, unless lambda<=0
106
107 double m_Glen_exponent, m_smoothing_range;
108
109 petsc::Vec::Ptr m_topgp0, //!< original bed elevation on processor 0
110 m_topgsmoothp0, //!< smoothed bed elevation on processor 0
111 m_maxtlp0, //!< maximum elevation at (i,j) of local topography (nearby patch)
112 m_C2p0, m_C3p0, m_C4p0;
113
114 virtual void preprocess_bed(const IceModelVec2S &topg,
115 unsigned int Nx_in, unsigned int Ny_in);
116
117 void smooth_the_bed_on_proc0();
118 void compute_coefficients_on_proc0();
119 };
120
121 } // end of namespace stressbalance
122 } // end of namespace pism
123
124 #endif // __BedSmoother_hh
125