Class soilSPL#

class eroder.soilSPL.soilSPL(*args, **kwargs)[source]#

The class computes river incision expressed using a stream power formulation function of river discharge and slope also accounting for soil production.

A non-linear diffusion of soil based on soil thickness is also implemented in this class.

If the user has turned-on the sedimentation capability, this class will solve implicitly the stream power formulation accounting for a sediment transport/deposition term (Yuan et al, 2019).

Methods

diffuseSoil()

For river-transported sediments reaching the marine realm, this function computes the related marine deposition diffusion.

erodepSPLsoil()

Modified stream power law model used to represent erosion by rivers also taking into account the role played by sediments in modulating erosion and deposition rate, considering non-linear slope dependency and accounting for soil production.

updateSoilThickness([deposition])

Updates soil thickness through time from the increment in self.tmp.

Initialise

__init__(*args, **kwargs)

Initialisation of soilSPL class.

Public Methods

updateSoilThickness([deposition])

Updates soil thickness through time from the increment in self.tmp.

erodepSPLsoil()

Modified stream power law model used to represent erosion by rivers also taking into account the role played by sediments in modulating erosion and deposition rate, considering non-linear slope dependency and accounting for soil production.

diffuseSoil()

For river-transported sediments reaching the marine realm, this function computes the related marine deposition diffusion.

Private Methods

_form_residual_soil(snes, h, F)

The nonlinear system (SNES) at each time step is solved iteratively by assessing the residual of the SPL equation accounting for erosion, deposition (transport-limited) and soil production.

_monitorsoil(snes, its, norm)

Non-linear SPL with soil production solver convergence evaluation.

_build_soil_snes([primary])

Construct (and return) a soil-SPL SNES solver and its residual vector.

_soilErodibility(PA, surfK)

(Re)build the dt-scaled stream-power erodibility coefficients Kbr (bedrock) and K_soil (soil layer) for the CURRENT self.dt and self.hOldArray.

_soilSNESsolve(x, guess)

Solve the soil-SPL SNES into x from initial guess guess: the primary solver, then the complementary quasi-Newton / ngmres fallback on failure (same primary→fallback net as the flow KSP).

_adaptiveSubstepSoil(x, PA, surfK)

Adaptive sub-stepping for the soil-SPL solve when the full-Δt SNES diverges.

_solveSoil()

Solves the non-linear stream power law for the transport limited and soil case.

_getEroDepRateSoil()

This function computes erosion deposition rates in metres per year and associated soil evolution.

_evalFunctionSoil(ts, t, x, xdot, f)

The non-linear system for soil diffusion is solved iteratively using PETSc time stepping and SNES solution and is based on Rosenbrock W-scheme (rosw).

_evalJacobianSoil(ts, t, x, xdot, a, A, B)

The non-linear system for soil diffusion is solved iteratively using PETSc time stepping and SNES solution and is based on Rosenbrock W-scheme (rosw).

_evalSolutionSoil(t, x)

Evaluate the initial solution of the SNES system.

Public functions#

soilSPL.updateSoilThickness(deposition=True)[source]#

Updates soil thickness through time from the increment in self.tmp.

Lsoil is a subaerial regolith cover, so the increment is not retained at subaqueous nodes (marine or ponded lake) — see _subaqueousMask — and is frozen under ice (_iceFrozenMask). These gates run in every mode, so callers invoke this unconditionally.

deposition distinguishes the two kinds of increment:

  • deposition=True (lake/pit, marine): in regolith mode the deposit is stratigraphy, not soil, so the increment is skipped (the gates still run — a cell newly ponded by this step’s deposition is re-zeroed consistently with lumped mode). In lumped mode it is added.

  • deposition=False (soil creep / hillslope transport): always added — creep moves the weathering regolith itself, in both modes.

soilSPL.erodepSPLsoil()[source]#

Modified stream power law model used to represent erosion by rivers also taking into account the role played by sediments in modulating erosion and deposition rate, considering non-linear slope dependency and accounting for soil production.

It calls the private function _getEroDepRateSoil described above. Once erosion/deposition rates have been calculated, the function computes local thicknesses and soil evolution for the considered time step and update local elevation and cumulative erosion, deposition values.

soilSPL.diffuseSoil()[source]#

For river-transported sediments reaching the marine realm, this function computes the related marine deposition diffusion. It is based on a non-linear diffusion approach.

\[\frac{\partial h}{\partial t}= \nabla \cdot \left( C_d \times (1.0 - e^{-h_s/H_0} \nabla h \right)\]

It calls the following private functions:

  • _evalFunctionSoil

  • _evalJacobianSoil

  • _evalSolutionSoil

Note

PETSc SNES and time stepping TS approaches are used to solve the non-linear equation above over the considered time step.

Private functions#

soilSPL._form_residual_soil(snes, h, F)[source]#

The nonlinear system (SNES) at each time step is solved iteratively by assessing the residual of the SPL equation accounting for erosion, deposition (transport-limited) and soil production.

soilSPL._monitorsoil(snes, its, norm)[source]#

Non-linear SPL with soil production solver convergence evaluation.

soilSPL._build_soil_snes(primary=True)[source]#

Construct (and return) a soil-SPL SNES solver and its residual vector.

Two configurations are available:

  • primary (primary=True) – a Nonlinear GMRES accelerator (ngmres) right-preconditioned by Nonlinear Richardson (nrichardson). The Richardson sweep is what actually applies the Krylov solve (cg) and the multigrid preconditioner (HYPRE BoomerAMG by default, or self.soil_pc); a bare ngmres ignores the KSP/PC entirely, which is why the previous setup stalled to SNES_DIVERGED_MAX_IT on the stiff soil-production residual.

  • fallback (primary=False) – a limited-memory quasi-Newton solver (qn, L-BFGS) with a matrix-free critical-point line search. It builds an approximate Jacobian from secant updates (no analytic Jacobian required) and is markedly more robust for stiff problems; it runs only when the primary solver fails to converge.

Each SNES gets its own PETSc options prefix so per-solver options (e.g. the line-search type) do not leak into the model’s other SNES/KSP objects.

Parameters:

primary – select the primary (True) or fallback (False) solver.

Returns:

the configured (SNES, residual Vec) pair.

soilSPL._soilErodibility(PA, surfK)[source]#

(Re)build the dt-scaled stream-power erodibility coefficients Kbr (bedrock) and K_soil (soil layer) for the CURRENT self.dt and self.hOldArray. The erosion limiter follows the current elevation (drop to the receiver), so this is recomputed per sub-step by the adaptive sub-stepper as well as once for the full step.

soilSPL._soilSNESsolve(x, guess)[source]#

Solve the soil-SPL SNES into x from initial guess guess: the primary solver, then the complementary quasi-Newton / ngmres fallback on failure (same primary→fallback net as the flow KSP). Returns the final converged reason (>= 0 converged, < 0 both failed). Prints only the informational “fallback recovered” line; the caller decides what to do on an outright failure (sub-step, then discard).

soilSPL._adaptiveSubstepSoil(x, PA, surfK)[source]#

Adaptive sub-stepping for the soil-SPL solve when the full-Δt SNES diverges. The divergence is a stiffness problem: a de-armoured node that has just captured a large discharge demands an enormous single-step incision (Kbr·Sⁿ with Kbr Δt·A^m), so the slope term overshoots and the residual grows. Because Kbr Δt, splitting the step into N increments of Δt/N cuts the demanded per-solve incision by N until the residual is smooth enough to converge; summing the increments recovers the FULL step’s erosion (unlike simply discarding).

Escalates N over _soilSubstepN (default 4, 8, 16); accepts the first N whose every sub-step converges. If none does, the step is discarded (revert to the prior elevation) as a last resort. On success the step-start bookkeeping (hOld/hOldArray/soilH/dt) is restored and one final residual evaluation rebuilds nsoilH / _soilDepoGrowth consistently with the full step (used by the soil update). Returns the final reason (0 sub-stepped OK, -1 discarded).

soilSPL._solveSoil()[source]#

Solves the non-linear stream power law for the transport limited and soil case. This calls the following private function:

  • _form_residual_soil

Note

PETSc SNES approach is used to solve the nonlinear equation without forming an analytic Jacobian. The primary solver is a Nonlinear GMRES accelerator (ngmres) right-preconditioned by Nonlinear Richardson (nrichardson), whose Krylov solve uses a Preconditioned Conjugate Gradient (cg) method with a multi-grid preconditioner (HYPRE BoomerAMG by default). If the primary solver stalls on the stiff soil-production residual, a limited-memory quasi-Newton fallback (qn, L-BFGS) with a critical-point line search is used. The iteration budget, tolerances and preconditioner are configurable through the YAML soil block (maxIter, rtol, atol, pcType).

soilSPL._getEroDepRateSoil()[source]#

This function computes erosion deposition rates in metres per year and associated soil evolution. This is done on the filled elevation.

The approach is based on BasicHySa governing equations from Terrainbento (as described in Appendix B20 from Barnhart et al. (2019)).

Note

The approach uses a continuous layer of soil-alluvium, which influences both hillslope and river-induced erosion. It relies on the SPACE algorithm of Shobe et al. (2017).

soilSPL._evalFunctionSoil(ts, t, x, xdot, f)[source]#

The non-linear system for soil diffusion is solved iteratively using PETSc time stepping and SNES solution and is based on Rosenbrock W-scheme (rosw).

Here again, we evaluate the residual function on a DMPlex for an implicit time-stepping method.

soilSPL._evalJacobianSoil(ts, t, x, xdot, a, A, B)[source]#

The non-linear system for soil diffusion is solved iteratively using PETSc time stepping and SNES solution and is based on Rosenbrock W-scheme (rosw).

Here, we define the Jacobian matrix A and the preconditioner matrix B on a DMPlex.

soilSPL._evalSolutionSoil(t, x)[source]#

Evaluate the initial solution of the SNES system.