AdvectionSolverTVD: Model horizontal advection of a given quantity using a Total Variation Diminishing (TVD) solution method#

Solve advection numerically using Total Variation Diminishing method.

class AdvectionSolverTVD(*args, **kwds)[source]#

Bases: Component

Numerical solution for advection using a Total Variation Diminishing method.

The component is restricted to regular grids (e.g., Raster or Hex). If multiple fields are advected, the advection__flux field will apply to the last one listed.

Parameters:

grid (RasterModelGrid or HexModelGrid) – A Landlab grid object.
fields_to_advect (field name or list or (n_nodes,) array (default None)) – A node field of scalar values that will be advected, or list of fields. If not given, the component creates a generic field, initialized to zeros, called advected__quantity. If list >1 element given, advection will be applied to each field in it.
advection_direction_is_steady (bool (default False)) – Indicates whether the directions of advection are expected to remain steady throughout a run. If True, some computation time is saved by calculating upwind links only once.

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> from landlab.components import AdvectionSolverTVD
>>> grid = RasterModelGrid((3, 7))
>>> s = grid.add_zeros("advected__quantity", at="node")
>>> s[9:12] = np.array([1.0, 2.0, 1.0])
>>> u = grid.add_zeros("advection__velocity", at="link")
>>> u[grid.horizontal_links] = 1.0
>>> advec = AdvectionSolverTVD(grid, fields_to_advect="advected__quantity")
>>> for _ in range(5):
...     advec.update(0.2)
>>> np.argmax(s[7:14])
4

Initialize AdvectionSolverTVD.

__init__(grid, fields_to_advect=None, advection_direction_is_steady=False)[source]#: Initialize AdvectionSolverTVD.

calc_rate_of_change_at_nodes(scalar, flux, dt, update_upwind_links=False)[source]#

Calculate time rate of change in the advected quantity at nodes.

Parameters:

scalar ((n_nodes, ) array) – Scalar at-node field of values to be advected.
dt (float) – Time-step duration. Needed to calculate the Courant number.
update_upwind_links (bool (optional; default False)) – If True, upwind links will be updated (set to True if the direction of advection is changing through time; if there are multiple advected quantities, it only needs to be True for the first one updated, and the update will be used for the others)

run_one_step(dt)[source]#

Update the solution by one time step dt.

Same as update.

Parameters:: dt (float) – Time-step duration. Needed to calculate the Courant number.

update(dt)[source]#

Update the solution by one time step dt.

Same as run_one_step.

Parameters:: dt (float) – Time-step duration. Needed to calculate the Courant number.

find_upwind_link_at_link(grid, u)[source]#

Return the upwind link at every link.

For all links, return ID of upwind link, defined based on the sign of u. If u is zero, the upwind link is found as though u were positive.

For instance (see examples below), consider a 3x4 raster grid with link numbering:

.-14-.-15-.-16-.
|    |    |    |
10  11   12   13
|    |    |    |
.--7-.--8-.--9-.
|    |    |    |
3    4    5    6
|    |    |    |
.--0-.--1-.--2-.

There are at most 7 active links (4, 5, 7, 8, 9, 11, 12). If u is positive everywhere, then the upwind links are:

.----.-14-.-15-.
|    |    |    |
3    4    5    6
|    |    |    |
.----.--7-.--8-.
|    |    |    |
|    |    |    |
|    |    |    |
.----.--0-.--1-.

If u is negative everywhere, then the upwind links are:

.-15-.-16-.----.
|    |    |    |
|    |    |    |
|    |    |    |
.--8-.--9-.----.
|    |    |    |
10  11   12   13
|    |    |    |
.--1-.--2-.----.

Parameters:

grid (RasterModelGrid or HexModelGrid) – A landlab grid.
u (float or (n_links,) ndarray) – Array of at-link values used to determine which node is upwind.

Returns:

The upwind links.

Return type:

(n_links,) ndarray of int

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> grid = RasterModelGrid((3, 4))

>>> uwl = find_upwind_link_at_link(grid, 1.0)
>>> uwl[grid.vertical_links].reshape((2, 4))
array([[-1, -1, -1, -1],
       [ 3,  4,  5,  6]])
>>> uwl[grid.horizontal_links].reshape((3, 3))
array([[-1,  0,  1],
       [-1,  7,  8],
       [-1, 14, 15]])

>>> uwl = find_upwind_link_at_link(grid, -1.0)
>>> uwl[grid.vertical_links].reshape((2, 4))
array([[10, 11, 12, 13],
       [-1, -1, -1, -1]])
>>> uwl[grid.horizontal_links].reshape((3, 3))
array([[ 1,  2, -1],
       [ 8,  9, -1],
       [15, 16, -1]])

>>> u = np.zeros(grid.number_of_links)
>>> u[4:6] = -1
>>> u[7] = -1
>>> u[8:10] = 1
>>> u[11:13] = 1
>>> u[grid.vertical_links].reshape((2, 4))
array([[ 0., -1., -1.,  0.],
       [ 0.,  1.,  1.,  0.]])
>>> u[grid.horizontal_links].reshape((3, 3))
array([[ 0.,  0.,  0.],
       [-1.,  1.,  1.],
       [ 0.,  0.,  0.]])
>>> uwl = find_upwind_link_at_link(grid, u)
>>> uwl[grid.vertical_links].reshape((2, 4))
array([[-1, 11, 12, -1],
       [ 3,  4, 5,   6]])
>>> uwl[grid.horizontal_links].reshape((3, 3))
array([[-1,  0,  1],
       [ 8,  7,  8],
       [-1, 14, 15]])

upwind_to_local_grad_ratio(grid, v, uwll, out=None)[source]#

Calculate and return ratio of upwind to local gradient in v.

Gradients are defined on links. Upwind is pre-determined via parameter uwll (upwind link at link), which can be obtained using the find_upwind_link_at_link function.

In Total Variation Diminishing (TVD) numerical schemes, this ratio is input to a flux limiter to calculate the weighting factor for higher-order vs. lower-order terms.

Parameters:

grid (RasterModelGrid or HexModelGrid) – A landlab grid.
v ((n_links,) ndarray) – Array of at-link values of which to calculate the gradient.
uwll ((n_links,) ndarray) – Array of upwind links for every link (as returned, for example, by find_upwind_link_at_link).
out ((n_links,) ndarray, optional) – If provided, place output into this array. Otherwise, create a new array.

Returns:

The ratio of the gradients. For links that have a gradient of zero, the ratio is set to one. For links that do not have an upwind link, the ratio is also set to one.

Return type:

(n_links,) ndarray of int