landlab.components.flow_director.flow_direction_mfd¶
flow_direction_mfd.py: calculate multiple-flow-direction flow directions.
Works on both a regular or irregular grid. Also calculates flow proportions.
KRB Jan 2017
- flow_directions_mfd(elev, neighbors_at_node, links_at_node, active_link_dir_at_node, link_slope, baselevel_nodes=None, partition_method='slope')[source]¶
Find multiple-flow-direction flow directions on a grid.
Finds and returns flow directions and proportions for a given elevation grid. Each node is assigned multiple flow directions, toward all of the N neighboring nodes that are lower than it. If none of the neighboring nodes are lower, it is assigned to itself. Flow proportions can be calculated as proportional to slope (default) or proportional to the square root of slope, which is the solution to a steady kinematic wave.
- Parameters:
elev (array_like) – Elevations at nodes.
neighbors_at_node (array_like (num nodes, max neighbors at node)) – For each node, the link IDs of active links.
links_at_node (array_like (num nodes, max neighbors at node))
link_dir_at_node (array_like (num nodes, max neighbors at node)) – IDs of the head node for each link.
link_slope (array_like) – slope of each link, defined POSITIVE DOWNHILL (i.e., a negative value means the link runs uphill from the fromnode to the tonode).
baselevel_nodes (array_like, optional) – IDs of open boundary (baselevel) nodes.
partition_method (string, optional) – Method for partitioning flow. Options include ‘slope’ (default) and ‘square_root_of_slope’.
- Returns:
receivers (ndarray of size (num nodes, max neighbors at node)) – For each node, the IDs of the nodes that receive its flow. For nodes that do not direct flow to all neighbors, BAD_INDEX_VALUE is given as a placeholder. The ID of the node itself is given if no other receiver is assigned.
proportions (ndarray of size (num nodes, max neighbors at node)) – For each receiver, the proportion of flow (between 0 and 1) is given. A proportion of zero indicates that the link does not have flow along it.
slopes (ndarray of size (num nodes, max neighbors at node)) – For each node in the array
recievers
, the slope value (positive downhill) in the direction of flow. If no flow occurs (value ofrecievers
is -1), then this array is set to 0.steepest_slope (ndarray) – The slope value (positive downhill) in the direction of flow.
steepest_receiver (ndarray) – For each node, the node ID of the node connected by the steepest link. BAD_INDEX_VALUE is given if no flow emmanates from the node.
sink (ndarray) – IDs of nodes that are flow sinks (they are their own receivers)
receiver_links (ndarray of size (num nodes, max neighbors at node)) – ID of links that leads from each node to its receiver, or BAD_INDEX_VALUE if no flow occurs on this link.
steepest_link (ndarray) – For each node, the link ID of the steepest link. BAD_INDEX_VALUE is given if no flow emmanates from the node.
Examples
>>> from landlab import RasterModelGrid >>> import numpy as np >>> from landlab.components.flow_director.flow_direction_mfd import ( ... flow_directions_mfd, ... ) >>> grid = RasterModelGrid((3, 3), xy_spacing=(1, 1)) >>> elev = grid.add_field( ... "topographic__elevation", ... grid.node_x + grid.node_y, ... at="node", ... )
For the first example, we will not pass any diagonal elements to the flow direction algorithm.
>>> neighbors_at_node = grid.adjacent_nodes_at_node >>> links_at_node = grid.links_at_node >>> active_link_dir_at_node = grid.active_link_dirs_at_node >>> link_slope = np.arctan(grid.calc_grad_at_link(elev)) >>> slopes_to_neighbors_at_node = ( ... link_slope[links_at_node] * active_link_dir_at_node ... ) >>> ( ... receivers, ... proportions, ... slopes, ... steepest_slope, ... steepest_receiver, ... sink, ... receiver_links, ... steepest_link, ... ) = flow_directions_mfd( ... elev, ... neighbors_at_node, ... links_at_node, ... active_link_dir_at_node, ... link_slope, ... baselevel_nodes=None, ... partition_method="slope", ... ) >>> receivers array([[ 0, -1, -1, -1], [ 1, -1, -1, -1], [ 2, -1, -1, -1], [ 3, -1, -1, -1], [-1, -1, 3, 1], [-1, -1, 4, -1], [ 6, -1, -1, -1], [-1, -1, -1, 4], [ 8, -1, -1, -1]]) >>> proportions array([[1. , 0. , 0. , 0. ], [1. , 0. , 0. , 0. ], [1. , 0. , 0. , 0. ], [1. , 0. , 0. , 0. ], [0. , 0. , 0.5, 0.5], [0. , 0. , 1. , 0. ], [1. , 0. , 0. , 0. ], [0. , 0. , 0. , 1. ], [1. , 0. , 0. , 0. ]]) >>> proportions.sum(axis=-1) array([1., 1., 1., 1., 1., 1., 1., 1., 1.])
In the second example, we will pass diagonal elements to the flow direction algorithm.
>>> dal = grid.active_d8 >>> neighbors_at_node = np.hstack( ... (grid.adjacent_nodes_at_node, grid.diagonal_adjacent_nodes_at_node) ... ) >>> links_at_node = grid.d8s_at_node >>> active_link_dir_at_node = grid.active_d8_dirs_at_node
We need to create a list of diagonal links since it doesn’t exist.
>>> diag_links = np.sort(np.unique(grid.d8s_at_node[:, 4:])) >>> diag_links = diag_links[diag_links > 0] >>> diag_grads = np.zeros(diag_links.shape) >>> where_active_diag = dal >= diag_links.min() >>> active_diags_inds = dal[where_active_diag] - diag_links.min() >>> diag_grads = grid.calc_grad_at_diagonal(elev) >>> ortho_grads = grid.calc_grad_at_link(elev) >>> link_slope = np.hstack((np.arctan(ortho_grads), np.arctan(diag_grads))) >>> ( ... receivers, ... proportions, ... slopes, ... steepest_slope, ... steepest_receiver, ... sink, ... receiver_links, ... steepest_link, ... ) = flow_directions_mfd( ... elev, ... neighbors_at_node, ... links_at_node, ... active_link_dir_at_node, ... link_slope, ... baselevel_nodes=None, ... partition_method="slope", ... ) >>> receivers array([[ 0, -1, -1, -1, -1, -1, -1, -1], [ 1, -1, -1, -1, -1, -1, -1, -1], [ 2, -1, -1, -1, -1, -1, -1, -1], [ 3, -1, -1, -1, -1, -1, -1, -1], [-1, -1, 3, 1, -1, -1, 0, -1], [-1, -1, 4, -1, -1, -1, -1, -1], [ 6, -1, -1, -1, -1, -1, -1, -1], [-1, -1, -1, 4, -1, -1, -1, -1], [-1, -1, -1, -1, -1, -1, 4, -1]]) >>> proportions array([[1. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [1. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [1. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [1. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [0. , 0. , 0.31091174, 0.31091174, 0. , 0. , 0.37817653, 0. ], [0. , 0. , 1. , 0. , 0. , 0. , 0. , 0. ], [1. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [0. , 0. , 0. , 1. , 0. , 0. , 0. , 0. ], [0. , 0. , 0. , 0. , 0. , 0. , 1. , 0. ]]) >>> slopes array([[0. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [0. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [0. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [0. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [0. , 0. , 0.78539816, 0.78539816, 0. , 0. , 0.95531662, 0. ], [0. , 0. , 0.78539816, 0. , 0. , 0. , 0. , 0. ], [0. , 0. , 0. , 0. , 0. , 0. , 0. , 0. ], [0. , 0. , 0. , 0.78539816, 0. , 0. , 0. , 0. ], [0. , 0. , 0. , 0. , 0. , 0. , 0.95531662, 0. ]]) >>> proportions.sum(axis=-1) array([1., 1., 1., 1., 1., 1., 1., 1., 1.])