IceUtility namespace : general purpose ice sheet functions common to the 'big' classes (AmrIce, JFNKSolver, etc) More...

Functions
void	computeC0 (Vector< LevelData< FArrayBox > * > &a_vectC0, const Vector< LevelData< FArrayBox > * > &a_vectC, const Vector< DisjointBoxLayout > &a_grids, const Vector< RefCountedPtr< LevelSigmaCS > > &a_coordSys, const Vector< Real > a_dx, int a_finest_level)
	Compute a C0 such that basal traction Tb = C0 u + f(C,u) More...

void	addWallDrag (FArrayBox &a_drag, const BaseFab< int > &a_mask, const FArrayBox &a_usrf, const FArrayBox &a_thk, const FArrayBox &a_topg, const FArrayBox &a_beta, const Real &a_extra, const RealVect &a_dx, const Box &a_box)

void	addThinIceDrag (FArrayBox &a_drag, const BaseFab< int > &a_mask, const FArrayBox &a_thk, const Real &a_extra, const Real &a_thin, const Box &a_box)

void	defineRHS (Vector< LevelData< FArrayBox > * > &a_rhs, const Vector< RefCountedPtr< LevelSigmaCS > > &a_CS, const Vector< DisjointBoxLayout > &a_grids, const Vector< RealVect > &a_dx)
	compute RHS for velocity field solve More...

void	applyHelmOp (LevelData< FArrayBox > &a_lapPhi, const LevelData< FArrayBox > &a_phi, const Real &a_a, const Real &a_b, const DisjointBoxLayout &a_grids, const RealVect &a_dx)
	compute the operator L phi = a * phi + b * grad^2 (phi) More...

void	applyDiv (LevelData< FArrayBox > &a_divU, const LevelData< FluxBox > &a_u, const DisjointBoxLayout &a_grids, const RealVect &a_dx)

void	applyGradSq (LevelData< FArrayBox > &a_gradPhiSq, const LevelData< FArrayBox > &a_phi, const DisjointBoxLayout &a_grids, const RealVect &a_dx)

void	computeFaceFluxCentered (LevelData< FluxBox > &a_us, const LevelData< FArrayBox > &a_u, const LevelData< FArrayBox > &a_s, const DisjointBoxLayout &a_grids)

void	computeFaceFluxUpwind (LevelData< FluxBox > &a_us, const LevelData< FluxBox > &a_u, const LevelData< FArrayBox > &a_s, const DisjointBoxLayout &a_grids)

void	computeA (LevelData< FArrayBox > &a_A, const Vector< Real > &a_sigma, const LevelSigmaCS &a_coordSys, const RateFactor *a_rateFactor, const LevelData< FArrayBox > &a_internalEnergy)
	compute cell centered rate factor A from the temperature More...

void	extrapVelocityToMargin (LevelData< FluxBox > &a_faceVel, const LevelData< FArrayBox > &a_cellVel, const LevelSigmaCS &a_coordSys)
	extrapolate face centered velocity field (usually derived by cell-to-face average) to the margins More...

void	computeFaceVelocity (LevelData< FluxBox > &a_faceVelAdvection, LevelData< FluxBox > &a_faceVelTotal, LevelData< FluxBox > &a_faceDiffusivity, LevelData< FArrayBox > &a_cellDiffusivity, #if BISICLES_Z==BISICLES_LAYERED LevelData< FluxBox > &a_layerXYFaceXYVel, LevelData< FArrayBox > &a_layerSFaceXYVel, #endif const LevelData< FArrayBox > &a_velocity, const LevelSigmaCS &a_coordSys, const IceThicknessIBC a_bc, const LevelData< FArrayBox > &a_A, #if BISICLES_Z==BISICLES_LAYERED const LevelData< FArrayBox > &a_sA, const LevelData< FArrayBox > &a_bA, #endif const LevelData< FArrayBox > a_crseVelocity, const LevelData< FArrayBox > a_crseDiffusivity, int a_nRefCrse, const ConstitutiveRelation a_constitutiveRelation, bool a_additionalVelocity, bool a_implicitDiffusion)
	compute face-centered velocity and thickness diffusion coefficients from cell-centered velocity More...

void	computeSigmaVelocity (LevelData< FArrayBox > &a_uSigma, const LevelData< FluxBox > &a_layerThicknessFlux, const LevelData< FArrayBox > &a_layerSFaceXYVel, const LevelData< FArrayBox > &a_dHdt, const DisjointBoxLayout a_grid, const LevelData< FArrayBox > &a_surfaceThicknessSource, const LevelData< FArrayBox > &a_basalThicknessSource, const Vector< Real > &a_dSigma, const RealVect &a_dx, const Real &a_dt)
	compute the cross layer velocity u^sigma (the contravariant component) More...

int	eliminateFastIce (Vector< RefCountedPtr< LevelSigmaCS > > &a_coordSys, Vector< LevelData< FArrayBox > * > &a_vel, Vector< LevelData< FArrayBox > * > &a_calvedIce, Vector< LevelData< FArrayBox > * > &a_addedIce, Vector< LevelData< FArrayBox > * > &a_removedIce, const Vector< DisjointBoxLayout > &a_grids, const Vector< ProblemDomain > &a_domain, const Vector< int > &a_refRatio, Real a_crseDx, int a_finestLevel, int a_maxIter, Real a_thinIceTol, Real a_fastIceTol, bool a_edgeOnly, int a_verbosity=0)
	Identify regions of fast ice and eliminate them. return the total number of cells eliminated. More...

void	eliminateRemoteIce (Vector< RefCountedPtr< LevelSigmaCS > > &a_coordSys, Vector< LevelData< FArrayBox > * > &a_vel, Vector< LevelData< FArrayBox > * > &a_calvedIce, Vector< LevelData< FArrayBox > * > &a_addedIce, Vector< LevelData< FArrayBox > * > &a_removedIce, const Vector< DisjointBoxLayout > &a_grids, const Vector< ProblemDomain > &a_domain, const Vector< int > &a_refRatio, Real a_crseDx, int a_finestLevel, int a_maxIter, Real a_tol, int a_verbosity=0)

void	multiplyByGroundedFraction (LevelData< FArrayBox > &a_u, const LevelSigmaCS &a_coords, const DisjointBoxLayout &a_grids, int a_subdivision)
	subgrid grounding line interpolation : multiply a_u by the grounded portion of each cell More...

void	setFloatingBasalFriction (LevelData< FArrayBox > &a_C, const LevelSigmaCS &a_coords, const DisjointBoxLayout &a_grids)
	set C = 0 in floating region, potentially using subgrid interpolation More...

Detailed Description

IceUtility namespace : general purpose ice sheet functions common to the 'big' classes (AmrIce, JFNKSolver, etc)

Function Documentation

◆ addThinIceDrag()

void IceUtility::addThinIceDrag	(	FArrayBox &	a_drag,
		const BaseFab< int > &	a_mask,
		const FArrayBox &	a_thk,
		const Real &	a_extra,
		const Real &	a_thin,
		const Box &	a_box
	)

Referenced by computeC0().

◆ addWallDrag()

void IceUtility::addWallDrag	(	FArrayBox &	a_drag,
		const BaseFab< int > &	a_mask,
		const FArrayBox &	a_usrf,
		const FArrayBox &	a_thk,
		const FArrayBox &	a_topg,
		const FArrayBox &	a_beta,
		const Real &	a_extra,
		const RealVect &	a_dx,
		const Box &	a_box
	)

Referenced by computeC0(), and InverseVerticallyIntegratedVelocitySolver::nDoF().

◆ applyDiv()

void IceUtility::applyDiv	(	LevelData< FArrayBox > &	a_divU,
		const LevelData< FluxBox > &	a_u,
		const DisjointBoxLayout &	a_grids,
		const RealVect &	a_dx
	)

apply cell-centred operator div(phi) to face-centred u. Assumes that ghost cells have been set

References computeFaceFluxCentered().

Referenced by applyHelmOp().

◆ applyGradSq()

void IceUtility::applyGradSq	(	LevelData< FArrayBox > &	a_gradPhiSq,
		const LevelData< FArrayBox > &	a_phi,
		const DisjointBoxLayout &	a_grids,
		const RealVect &	a_dx
	)

apply cell-centred operator [grad(phi)]^2 to face-centred u. Assumes that ghost cells have been set

References computeA().

Referenced by computeFaceFluxUpwind(), and InverseVerticallyIntegratedVelocitySolver::nDoF().

◆ applyHelmOp()

void IceUtility::applyHelmOp	(	LevelData< FArrayBox > &	a_lapPhi,
		const LevelData< FArrayBox > &	a_phi,
		const Real &	a_a,
		const Real &	a_b,
		const DisjointBoxLayout &	a_grids,
		const RealVect &	a_dx
	)

compute the operator L phi = a * phi + b * grad^2 (phi)

References applyDiv().

Referenced by defineRHS(), and InverseVerticallyIntegratedVelocitySolver::nDoF().

◆ computeA()

void IceUtility::computeA	(	LevelData< FArrayBox > &	a_A,
		const Vector< Real > &	a_sigma,
		const LevelSigmaCS &	a_coordSys,
		const RateFactor *	a_rateFactor,
		const LevelData< FArrayBox > &	a_internalEnergy
	)

compute cell centered rate factor A from the temperature

References RateFactor::computeA(), IceThermodynamics::decomposeInternalEnergy(), LevelSigmaCS::getH(), LevelSigmaCS::gravity(), LevelSigmaCS::grids(), LevelSigmaCS::iceDensity(), and IceThermodynamics::icepmeltfactor().

Referenced by applyGradSq(), AmrIce::computeA(), PatersonRateFactor::setA0(), AmrIce::setBasalFriction(), and RateFactor::~RateFactor().

◆ computeC0()

void IceUtility::computeC0	(	Vector< LevelData< FArrayBox > * > &	a_vectC0,
		const Vector< LevelData< FArrayBox > * > &	a_vectC,
		const Vector< DisjointBoxLayout > &	a_grids,
		const Vector< RefCountedPtr< LevelSigmaCS > > &	a_coordSys,
		const Vector< Real >	a_dx,
		int	a_finest_level
	)

Compute a C0 such that basal traction Tb = C0 u + f(C,u)

Basal traction Tb = C0 u + f(C,u)

f(C,u) is the conventional part of the basal traction.

the C0 u term is used to add a linear drag, and for now is used only for 'wall drag', which allows for a crude representation of ice friction along rocky walls. Two Parmparse parameters, wall_drag.basic and wall_drag.extra, determine the calculation.

if (wall_drag.basic) is true, then C0 = (C + wall_drag.extra) u

For backward compatibility, wall_drag.basic == amr.WallDrag and wall_drag.extra = amr.wallDragExtra

References addThinIceDrag(), addWallDrag(), LevelSigmaCS::getFloatingMask(), LevelSigmaCS::getH(), LevelSigmaCS::getSurfaceHeight(), and LevelSigmaCS::getTopography().

Referenced by AmrIce::setBasalFriction().

◆ computeFaceFluxCentered()

void IceUtility::computeFaceFluxCentered	(	LevelData< FluxBox > &	a_us,
		const LevelData< FArrayBox > &	a_u,
		const LevelData< FArrayBox > &	a_s,
		const DisjointBoxLayout &	a_grids
	)

compute face-centered us given cell-centered vector u and scalar s Assumes that boundary values have been set. USES CENTERED SCHEME

References computeFaceFluxUpwind().

Referenced by applyDiv().

◆ computeFaceFluxUpwind()

void IceUtility::computeFaceFluxUpwind	(	LevelData< FluxBox > &	a_us,
		const LevelData< FluxBox > &	a_u,
		const LevelData< FArrayBox > &	a_s,
		const DisjointBoxLayout &	a_grids
	)

compute face-centered us given cell-centered vector u and scalar s Assumes that boundary values have been set. USES FIRST ORDER UPWIND SCHEME

References applyGradSq().

Referenced by computeFaceFluxCentered().

◆ computeFaceVelocity()

void IceUtility::computeFaceVelocity	(	LevelData< FluxBox > &	a_faceVelAdvection,
		LevelData< FluxBox > &	a_faceVelTotal,
		LevelData< FluxBox > &	a_faceDiffusivity,
		LevelData< FArrayBox > &	a_cellDiffusivity,
		#if	BISICLES_Z = `= BISICLES_LAYERED LevelData<FluxBox>& a_layerXYFaceXYVel`,
		LevelData< FArrayBox > &	a_layerSFaceXYVel,
		#endif const LevelData< FArrayBox > &	a_velocity,
		const LevelSigmaCS &	a_coordSys,
		const IceThicknessIBC *	a_bc,
		const LevelData< FArrayBox > &	a_A,
		#if	BISICLES_Z = `= BISICLES_LAYERED const LevelData<FArrayBox>& a_sA`,
		const LevelData< FArrayBox > &	a_bA,
		#endif const LevelData< FArrayBox > *	a_crseVelocity,
		const LevelData< FArrayBox > *	a_crseDiffusivity,
		int	a_nRefCrse,
		const ConstitutiveRelation *	a_constitutiveRelation,
		bool	a_additionalVelocity,
		bool	a_implicitDiffusion
	)

compute face-centered velocity and thickness diffusion coefficients from cell-centered velocity

References computeSigmaVelocity(), extrapVelocityToMargin(), IceThicknessIBC::modifyFaceVelocity(), and L1L2ConstitutiveRelation::modifyTransportCoefficients().

Referenced by AmrIce::computeFaceVelocity(), InverseVerticallyIntegratedVelocitySolver::computeObjectiveAndGradient(), and extrapVelocityToMargin().

◆ computeSigmaVelocity()

void IceUtility::computeSigmaVelocity	(	LevelData< FArrayBox > &	a_uSigma,
		const LevelData< FluxBox > &	a_layerThicknessFlux,
		const LevelData< FArrayBox > &	a_layerSFaceXYVel,
		const LevelData< FArrayBox > &	a_dHdt,
		const DisjointBoxLayout	a_grid,
		const LevelData< FArrayBox > &	a_surfaceThicknessSource,
		const LevelData< FArrayBox > &	a_basalThicknessSource,
		const Vector< Real > &	a_dSigma,
		const RealVect &	a_dx,
		const Real &	a_dt
	)

compute the cross layer velocity u^sigma (the contravariant component)

References eliminateFastIce().

Referenced by computeFaceVelocity(), and AmrIce::updateInternalEnergy().

◆ defineRHS()

void IceUtility::defineRHS	(	Vector< LevelData< FArrayBox > * > &	a_rhs,
		const Vector< RefCountedPtr< LevelSigmaCS > > &	a_CS,
		const Vector< DisjointBoxLayout > &	a_grids,
		const Vector< RealVect > &	a_dx
	)

compute RHS for velocity field solve

References LevelSigmaCS::anyFloating(), applyHelmOp(), LevelSigmaCS::getFloatingMask(), LevelSigmaCS::getGradSurface(), LevelSigmaCS::getH(), LevelSigmaCS::getSurfaceHeight(), LevelSigmaCS::gravity(), and LevelSigmaCS::iceDensity().

Referenced by AmrIce::defineVelRHS(), and JFNKSolver::solve().

◆ eliminateFastIce()

int IceUtility::eliminateFastIce	(	Vector< RefCountedPtr< LevelSigmaCS > > &	a_coordSys,
		Vector< LevelData< FArrayBox > * > &	a_vel,
		Vector< LevelData< FArrayBox > * > &	a_calvedIce,
		Vector< LevelData< FArrayBox > * > &	a_addedIce,
		Vector< LevelData< FArrayBox > * > &	a_removedIce,
		const Vector< DisjointBoxLayout > &	a_grids,
		const Vector< ProblemDomain > &	a_domain,
		const Vector< int > &	a_refRatio,
		Real	a_crseDx,
		int	a_finestLevel,
		int	a_maxIter,
		Real	a_thinIceTol,
		Real	a_fastIceTol,
		bool	a_edgeOnly,
		int	a_verbosity = `0`
	)

Identify regions of fast ice and eliminate them. return the total number of cells eliminated.

Identify regions of fast ice and eliminate them.

implies a call to IceUtility::eliminateRemoteIce iff any fast ice is removed

References eliminateRemoteIce(), LevelSigmaCS::getFloatingMask(), LevelSigmaCS::getH(), and CalvingModel::updateCalvedIce().

Referenced by computeSigmaVelocity(), JFNKSolver::get_evil_configuration(), and JFNKSolver::solve().

◆ eliminateRemoteIce()

void IceUtility::eliminateRemoteIce	(	Vector< RefCountedPtr< LevelSigmaCS > > &	a_coordSys,
		Vector< LevelData< FArrayBox > * > &	a_vel,
		Vector< LevelData< FArrayBox > * > &	a_calvedIce,
		Vector< LevelData< FArrayBox > * > &	a_addedIce,
		Vector< LevelData< FArrayBox > * > &	a_removedIce,
		const Vector< DisjointBoxLayout > &	a_grids,
		const Vector< ProblemDomain > &	a_domain,
		const Vector< int > &	a_refRatio,
		Real	a_crseDx,
		int	a_finestLevel,
		int	a_maxIter,
		Real	a_tol,
		int	a_verbosity = `0`
	)

Identify regions of floating ice that are remote from grounded ice and eliminate them.

Identify regions of floating ice that are remote from grounded ice and eliminate them. Regions of floating ice unconnected to land lead to an ill-posed problem, with a zero basal traction coefficient and Neumann boundary conditions. Here, we attempt to identify them by carrying out a procedure that in the worst case can be O(N^2).

References LevelSigmaCS::getFloatingMask(), LevelSigmaCS::getH(), multiplyByGroundedFraction(), LevelSigmaCS::recomputeGeometry(), and CalvingModel::updateCalvedIce().

Referenced by eliminateFastIce(), and AmrIce::eliminateRemoteIce().

◆ extrapVelocityToMargin()

void IceUtility::extrapVelocityToMargin	(	LevelData< FluxBox > &	a_faceVel,
		const LevelData< FArrayBox > &	a_cellVel,
		const LevelSigmaCS &	a_coordSys
	)

extrapolate face centered velocity field (usually derived by cell-to-face average) to the margins

References computeFaceVelocity(), LevelSigmaCS::getH(), LevelSigmaCS::getSurfaceHeight(), and LevelSigmaCS::getTopography().

Referenced by computeFaceVelocity().

◆ multiplyByGroundedFraction()

void IceUtility::multiplyByGroundedFraction	(	LevelData< FArrayBox > &	a_u,
		const LevelSigmaCS &	a_coords,
		const DisjointBoxLayout &	a_grids,
		int	a_subdivision
	)

subgrid grounding line interpolation : multiply a_u by the grounded portion of each cell

multiply a_u by the grounded portion of each cell

For a_subdivision > 0, the grounded portion of a cell is computed by (1) approximating the thickness above/below flotation in each quarter of each cell with a bilinear formula h(x,y) (2) integrating the Heaviside function (h > 0)?1:0 over each quarter using the midpoint rule in [2^a_n]^2 subdivisions (making 4*[2^a_n]^2 subcells)

a_subdivision = 2 (subcell size = dx/8) is probably plenty, the cost is obviously O(a_subdivision^2)

References LevelSigmaCS::getH(), LevelSigmaCS::getTopography(), LevelSigmaCS::iceDensity(), column_thermodynamics::rhoi, column_thermodynamics::rhoo, LevelSigmaCS::seaLevel(), setFloatingBasalFriction(), and LevelSigmaCS::waterDensity().

Referenced by eliminateRemoteIce(), and setFloatingBasalFriction().

◆ setFloatingBasalFriction()

void IceUtility::setFloatingBasalFriction	(	LevelData< FArrayBox > &	a_C,
		const LevelSigmaCS &	a_coords,
		const DisjointBoxLayout &	a_grids
	)

set C = 0 in floating region, potentially using subgrid interpolation

set C = 0 in floating region

If nRef <= 0, C is set to zero according to the mask in a_levelCS

Otherwise Hab, the thickness above/below flotation (> 0 on grounded ice, < 0 in the shelf) is computed, and then interpolated (using the bilinear formula : this is not the usual conservative Chombo formula. C is then multplied by the integral of (Hab > 0)?1:0 over the each cell.

References LevelSigmaCS::anyFloating(), LevelSigmaCS::getFloatingMask(), and multiplyByGroundedFraction().

Referenced by multiplyByGroundedFraction(), InverseVerticallyIntegratedVelocitySolver::nDoF(), and AmrIce::setBasalFriction().

Functions