PISM bed smoother, plus bed roughness parameterization, based on Schoof (2003). More...

#include <PISMBedSmoother.hh>

Public Member Functions
	PISMBedSmoother (IceGrid &g, const NCConfigVariable &conf, PetscInt MAX_GHOSTS)
virtual	~PISMBedSmoother ()
virtual PetscErrorCode	preprocess_bed (IceModelVec2S topg, PetscReal n, PetscReal lambda)
virtual PetscErrorCode	preprocess_bed (IceModelVec2S topg, PetscReal n, PetscInt Nx_in, PetscInt Ny_in)
virtual PetscErrorCode	get_smoothing_domain (PetscInt &Nx_out, PetscInt &Ny_out)
virtual PetscInt	get_max_ghosts ()
virtual PetscErrorCode	get_smoothed_thk (IceModelVec2S usurf, IceModelVec2S thk, IceModelVec2Int mask, PetscInt GHOSTS, IceModelVec2S *thksmooth)
	Computes a smoothed thickness map.
virtual PetscErrorCode	get_theta (IceModelVec2S usurf, PetscReal n, PetscInt GHOSTS, IceModelVec2S *theta)
Public Attributes
IceModelVec2S	topgsmooth
	smoothed bed elevation; publicly-available; set by calling preprocess_bed(); has ghosts with width get_max_ghosts()
Protected Member Functions
PetscErrorCode	allocate ()
PetscErrorCode	deallocate ()
PetscErrorCode	smooth_the_bed_on_proc0 ()
	Computes the smoothed bed by a simple average over a rectangle of grid points.
PetscErrorCode	compute_coefficients_on_proc0 (PetscReal n)
Protected Attributes
IceGrid &	grid
const NCConfigVariable &	config
IceModelVec2S	maxtl
IceModelVec2S	C2
IceModelVec2S	C3
IceModelVec2S	C4
PetscInt	Nx
PetscInt	Ny
PetscInt	maxGHOSTS
Vec	g2
Vec	g2natural
	global Vecs used to transfer data to/from processor 0.
VecScatter	scatter
	VecScatter used to transfer data to/from processor 0.
Vec	topgp0
	original bed elevation on processor 0
Vec	topgsmoothp0
	smoothed bed elevation on processor 0
Vec	maxtlp0
	maximum elevation at (i,j) of local topography (nearby patch)
Vec	C2p0
Vec	C3p0
Vec	C4p0

Detailed Description

PISM bed smoother, plus bed roughness parameterization, based on Schoof (2003).

This class both smooths the bed and computes coefficients for an approximation to Schoof's $\theta$ . The factor $0\le \theta \le 1$ multiplies the diffusivity in the theory of the effect of bed roughness in the SIA by Christian Schoof (2003; The effect of basal topography on ice sheet dynamics) [Schoofbasaltopg2003].

For additional information on this class see page Using Schoof's (2003) parameterized bed roughness technique in PISM.

The user of this class hands PISMBedSmoother an "original" topography, and it is preprocessed to fill the smoothed topography topgsmooth, and the coefficients in an approximation to $\theta$ . This is done by a call to preprocess_bed(). The call requires the half-width of the smoothing square (a distance in m), or the number of grid points in each direction in the smoothing rectangle, and the Glen exponent.

The call to preprocess_bed() must be repeated any time the "original" topography changes, for instance at the start of an IceModel run, or at a bed deformation step in an IceModel run.

PISMBedSmoother then provides three major functionalities, all of which must follow the call to preprocess_bed():

User accesses public IceModelVec2S topgsmooth, the smoothed bed itself.
User asks get_smoothed_thk() for gridded values of the consistent smoothed version of the ice thickness, which is the thickness corresponding to a given surface elevation and the pre-computed smoothed bed.
User asks for gridded values of $\theta(h,x,y)$ using get_theta().

Here is a basic example of the creation and usage of a PISMBedSmoother instance. Note that PISMBedSmoother will update ghosted values in topgsmooth, and other internal fields, and update them in the return fields thksmooth, theta, if asked. In IceModel::velocitySIAStaggered(), MAX_GHOSTS=2

    PISMBedSmoother smoother(grid, config, 2);
    const PetscReal n = 3.0,
                    lambda = 50.0e3;
    ierr = smoother.preprocess_bed(topg, n, lambda); CHKERRQ(ierr);
    ierr = smoother.get_smoothed_thk(usurf, thk, 1, &thksmooth); CHKERRQ(ierr);
    ierr = smoother.get_theta(usurf, n, 1, &theta); CHKERRQ(ierr);

See IceGrid and NCConfigVariable documentation for initializing grid and config. Note we assume topg, usurf, thk, thksmooth, and theta are all created IceModelVec2S instances.

Definition at line 72 of file PISMBedSmoother.hh.

Constructor & Destructor Documentation

PISMBedSmoother	(	IceGrid &	g,
		const NCConfigVariable &	conf,
		PetscInt	MAX_GHOSTS
	)

Definition at line 23 of file PISMBedSmoother.cc.

References allocate(), config, and grid.

~PISMBedSmoother ( ) [virtual]

Definition at line 41 of file PISMBedSmoother.cc.

References deallocate(), and grid.

Member Function Documentation

PetscErrorCode allocate ( ) [protected]

Definition at line 71 of file PISMBedSmoother.cc.

References C2, C2p0, C3, C3p0, C4, C4p0, g2, g2natural, grid, maxGHOSTS, maxtl, maxtlp0, scatter, topgp0, topgsmooth, and topgsmoothp0.

Referenced by PISMBedSmoother().

PetscErrorCode compute_coefficients_on_proc0 ( PetscReal n ) [protected]

Definition at line 252 of file PISMBedSmoother.cc.

References b0, C2p0, C3p0, C4p0, count, grid, maxtlp0, n, Nx, Ny, topgp0, and topgsmoothp0.

Referenced by preprocess_bed().

PetscErrorCode deallocate ( ) [protected]

Definition at line 120 of file PISMBedSmoother.cc.

References C2p0, C3p0, C4p0, g2, g2natural, maxtlp0, scatter, topgp0, and topgsmoothp0.

Referenced by ~PISMBedSmoother().

virtual PetscInt get_max_ghosts ( ) [inline, virtual]

Definition at line 83 of file PISMBedSmoother.hh.

References maxGHOSTS.

PetscErrorCode get_smoothed_thk	(	IceModelVec2S	usurf,
		IceModelVec2S	thk,
		IceModelVec2Int	mask,
		PetscInt	GHOSTS,
		IceModelVec2S *	thksmooth
	)		`[virtual]`

Computes a smoothed thickness map.

The result thksmooth is the difference between the given upper surface elevation (usurf) and the stored smoothed bed topography (topgsmooth), except where the given original thickness (thk) is zero. In places where the original thickness is zero, the result thksmooth is also set to zero.

Ghosted values are updated directly and no communication occurs. In fact, we assume usurf, thk, and thksmooth all have stencil width at least equal to GHOSTS. We check whether topgsmooth, which has stencil width maxGHOSTS, has at least GHOSTS stencil width, and throw an error if not.

Call preprocess_bed() first.

Definition at line 333 of file PISMBedSmoother.cc.

References grid, maxGHOSTS, maxtl, and topgsmooth.

Referenced by SIAFD::compute_diffusive_flux(), SIAFD::compute_diffusivity(), SIAFD::compute_I(), and SIAFD::compute_sigma().

PetscErrorCode get_smoothing_domain	(	PetscInt &	Nx_out,
		PetscInt &	Ny_out
	)		`[virtual]`

Call preprocess_bed() first.

Definition at line 210 of file PISMBedSmoother.cc.

References Nx, and Ny.

Referenced by main().

PetscErrorCode get_theta	(	IceModelVec2S	usurf,
		PetscReal	n,
		PetscInt	GHOSTS,
		IceModelVec2S *	theta
	)		`[virtual]`

Implements the strategy for computing $\theta(h,x,y)$ from previously- stored coefficients, described on page Using Schoof's (2003) parameterized bed roughness technique in PISM and in [Schoofbasaltopg2003].

Specifically, $\theta = \omega^{-n}$ where $\omega$ is a local average of a rational function of surface elevation, approximated here by a Taylor polynomial:

$\omega = \fint \left(1 - \frac{\tilde b(x_1,x_2,\xi_1,\xi_2)}{H} \right)^{-(n+2)/n}\,d\xi_1\,d\xi_2 \approx 1 + C_2 H^{-2} + C_3 H^{-3} + C_4 H^{-4}$

where $h =$ usurf, $H = h -$ topgsmooth and $\tilde b$ is the local bed topography, a function with mean zero. The coefficients $C_2,C_3,C_4$ , which depend on $x,y$ , are precomputed by preprocess_bed().

Ghosted values are updated directly and no communication occurs. In fact, we assume usurf and theta have stencil width at least equal to GHOSTS. We check whether topgsmooth, which has stencil width maxGHOSTS, has at least GHOSTS stencil width, and throw an error if not.

Call preprocess_bed() first.

Definition at line 412 of file PISMBedSmoother.cc.

References C2, C3, C4, grid, exactV::H(), maxGHOSTS, maxtl, Nx, Ny, and topgsmooth.

Referenced by SIAFD::compute_diffusive_flux(), and main().

PetscErrorCode preprocess_bed	(	IceModelVec2S	topg,
		PetscReal	n,
		PetscInt	Nx_in,
		PetscInt	Ny_in
	)		`[virtual]`

Inputs Nx,Ny gives half-width in number of grid points, over which to do the average.

Definition at line 174 of file PISMBedSmoother.cc.

References C2, C2p0, C3, C3p0, C4, C4p0, compute_coefficients_on_proc0(), g2, g2natural, grid, maxtl, maxtlp0, Nx, Ny, scatter, smooth_the_bed_on_proc0(), topgp0, topgsmooth, and topgsmoothp0.

PetscErrorCode preprocess_bed	(	IceModelVec2S	topg,
		PetscReal	n,
		PetscReal	lambda
	)		`[virtual]`

Input lambda gives physical half-width (in m) of square over which to do the average. Only square smoothing domains are allowed with this call, which is the default case.

Definition at line 143 of file PISMBedSmoother.cc.

References grid, Nx, Ny, and topgsmooth.

Referenced by main(), and SIAFD::update().

PetscErrorCode smooth_the_bed_on_proc0 ( ) [protected]

Computes the smoothed bed by a simple average over a rectangle of grid points.

Definition at line 218 of file PISMBedSmoother.cc.

References b0, count, grid, Nx, Ny, topgp0, and topgsmoothp0.

Referenced by preprocess_bed().