water_tools

The route_water routine manages the water routing. During route_water, water iteration is repeated a total of itermax times. During each of these iterations of the water routing, the following methods are called in order:

water_tools.init_water_iteration()	Init the water iteration routine.
water_tools.run_water_iteration()	Run a single iteration of travel paths for all water parcels.
water_tools.compute_free_surface()	Calculate free surface after routing all water parcels.
water_tools.finalize_water_iteration(iteration)	Finish updating flow fields.

Public API methods attached to model

The following methods are defined in the water_tools class, of the pyDeltaRCM.water_tools module.

water_tools()

class pyDeltaRCM.water_tools.water_tools

build_weight_array(array, fix_edges: bool = False, normalize: bool = False)

Weighting array of neighbors.

Create np.array((8,L,W)) of quantity a in each of the neighbors to a cell

check_for_boundary(inds)

Check whether parcels have reached the boundary.

Checks whether any parcels have reached the model boundaries. If they have, then update the information in free_surf_flag.

Parameters:

inds (ndarray) – Unraveled indicies of parcels.

check_size_of_indices_matrix(it) → None

Check if step path matrix needs to be made larger.

Initial size of self.free_surf_walk_inds is half of self.stepmax because the number of iterations doesn’t go beyond that for many timesteps.

Once it reaches it > self.stepmax/2 once, make the size self.iter for all further timesteps

compute_free_surface() → None

Calculate free surface after routing all water parcels.

This method uses the free_surf_walk_inds matrix accumulated during the routing of the water parcels (in run_water_iteration) to determine the free surface. The operations of the free surface computation are placed in a jitted function _accumulate_free_surface_walks. Following this computation, the free surface is smoothed by steps in finalize_free_surface.

See [1] and [2] for a complete description of hydrodynamic assumptions in the DeltaRCM model.

Examples

The sequence of compute_free_surface is depicted in the figures below. The first image depicts the “input” to compute_free_surface, which is the current bed elevation, and the path of each water parcel (cyan lines in right image).

(png, hires.png)

The compute_free_surface method then calls the _accumulate_free_surface_walks function to determine 1) the number of times each cell has been visited by a water parcel (sfc_visit), and 2) the total sum of expected elevations of the water surface at each cell (sfc_sum).

The output from _accumulate_free_surface_walks is then used to calculate a new stage surface (H_new) based only on the water parcel paths and expected water surface elevations, approximately as Hnew = sfc_sum / sfc_visit (“computed Hnew” in figure below)

Following this step, a correction is applied forcing the new free surface to be below sea level and above or equal to the land surface elevation over the model domain (“stage” in figure below). This surface Hnew is used in the following operation finalize_free_surface.

Finally, the updated water surface is combined with the previous timestep’s water surface and an underrelaxation coefficient omega_sfc.

(png, hires.png)

[1]

A reduced-complexity model for river delta formation – Part 1: Modeling deltas with channel dynamics, M. Liang, V. R. Voller, and C. Paola, Earth Surf. Dynam., 3, 67–86, 2015. https://doi.org/10.5194/esurf-3-67-2015

[2]

A reduced-complexity model for river delta formation – Part 2: Assessment of the flow routing scheme, M. Liang, N. Geleynse, D. A. Edmonds, and P. Passalacqua, Earth Surf. Dynam., 3, 87–104, 2015. https://doi.org/10.5194/esurf-3-87-2015

finalize_free_surface() → None

Finalize the water free surface.

This method occurs after the initial computation of the free surface, and creates the new free surface by smoothing the newly computed free surface with a jitted function (_smooth_free_surface), and then combining the old surface and new smoothed surface with underrelaxation. Finally, a flooding_correction is applied.

finalize_water_iteration(iteration: int) → None

Finish updating flow fields.

Clean up at end of water iteration

flooding_correction() → None

Flood dry cells along the shore if necessary.

Check the neighbors of all dry cells. If any dry cells have wet neighbors, check that their stage is not higher than the bed elevation of the center cell. If it is, flood the dry cell.

get_inlet_weights_water(**kwargs)

Get weight for inlet cells for water parcels.

This method determines the weights describing which inlet cells water parcels are fed into the domain (where the cells are self.inlet). Internally, pyDeltaRCM.shared_tools._get_inlet_weights is called, returning a balanced weighting across all inlet cells.

Note

Reimplement this method in custom subclasses as needed.

get_water_weight_array() → None

Get step direction weights for each cell.

This method is called once, before parcels are stepped, because the weights do not change during the stepping of parcels.

See Hydrodynamics for a description of the model design and equations/algorithms for how water weights are determined.

Note

No computation actually occurs in this method. Internally, the method calls the jitted _get_water_weight_array, which in turn loops through all of the cells in the model domain and determines the water weight for that cell with get_weight_sfc_int and _get_weight_at_cell_water.

Examples

The following figure shows several examples of locations within the model domain, and the corresponding water routing weights determined for that location.

(png, hires.png)

get_wet_mask_nh()

Get wet mask.

Returns np.array((8,L,W)), for each neighbor around a cell with 1 if the neighbor is wet and 0 if dry

init_water_iteration() → None

Init the water iteration routine.

route_water() → None

Water routing main method.

This is the main method for water routing in the model. It is called once per update() call.

Internally, this method calls:

• init_water_iteration

• run_water_iteration

• compute_free_surface

• finalize_water_iteration

run_water_iteration() → None

Run a single iteration of travel paths for all water parcels.

Runs all water parcels (Np_water parcels) for stepmax steps, or until the parcels reach a boundary.

All parcels are processed in parallel, taking one step for each loop of the while loop.

Example

On the initial delta surface, see how ten selected parcels are routed through the domain:

(png, hires.png)

update_Q(dist, current_inds, next_inds, astep, jstep, istep, update_current: bool = False, update_next: bool = False) → None

Update discharge field values.

Method is called after one set of water parcel steps.

update_flow_field(iteration: int) → None

Update water discharge.

Update the water discharge field after one set of water parcels iteration.

update_velocity_field() → None

Update flow velocity fields.

Update the flow velocity fields after one set of water parcels iteration.

water_tools helper functions

The following routines are jitted for speed. They generally take a large number of arrays and constants and return a new array(s) to continue with the model progression within the methods defined above.

pyDeltaRCM.water_tools._get_weight_at_cell_water(ind, weight_sfc, weight_int, depth_nbrs, mod_water_weight, ct_nbrs, dry_depth: float, gamma: float, theta: float)

Compute water weights for a given cell.

This is a jitted function called by _get_water_weight_array().

pyDeltaRCM.water_tools._choose_next_directions(inds: ndarray, water_weights: ndarray) → ndarray

Get new cell locations, based on water weights.

Algorithm is to:

1. loop through each parcel, which is described by a pair in the inds array.

2. determine the water weights for that location (from water_weights)

3. choose a new cell based on the probabilities of the weights (using the random_pick function)

Parameters:

• inds (ndarray) – Current unraveled indices of the parcels. (N,) ndarray containing the unraveled indices.

• water_weights (ndarray) – Weights of every water cell. (LxW, 9) ndarray, uses unraveled indicies along 0th dimension; 9 cells represent self and 8 neighboring cells.

Returns:

next_direction – The direction to move towards the new cell for water parcels, relative to the current location. I.e., this is the D8 direction the parcel is going to travel in the next stage, pyDeltaRCM.shared_tools._calculate_new_ind.

Return type:

ndarray

pyDeltaRCM.water_tools._calculate_new_inds(current_inds: ndarray, new_direction, ravel_walk)

Calculate the new location (current_inds) of parcels.

Use the information of the current parcel (current_inds) in conjunction with the D8 direction the parcel needs to travel (new_direction) to determine the new current_inds of each parcel.

In implementation, we use the flattened ravel_walk array, but the result is identical to unraveling the index, adding iwalk and jwalk to the location and then raveling the index back.

ind_tuple = shared_tools.custom_unravel(ind, domain_shape)
new_ind = (ind_tuple[0] + jwalk[newd],
           ind_tuple[1] + iwalk[newd])
new_inds[p] = shared_tools.custom_ravel(new_ind, domain_shape)

pyDeltaRCM.water_tools._check_for_loops(free_surf_walk_inds: ndarray, new_inds, _step: int, L0: int, CTR: int, stage_above_SL: ndarray)

Check for loops in water parcel pathways.

Look for looping random walks. I.e., this function checks for where a parcel will return on its new_inds to somewhere it has already been in free_surf_walk_inds. If the loop is found, the parcel is relocated along the mean transport vector of the parcel, which is computed as the vector from the cell (0, CTR) to the new location in new_inds.

This implementation of loop checking will relocate any parcel that has looped, but only disqualifies a parcel p from contributing to the free surface in _accumulate_free_surf_walks (i.e., looped[p] == 1) if the stage at the looped location is above the sea level in the domain.

Parameters:

• free_surf_walk_inds – Array recording the walk of parcels. Shape is (:obj:`Np_water, …)`, where the second dimension will depend on the step number, but records each step of the parcel. Each element in the array records the flat index into the domain.

• new_inds – Array recording the new index for each water parcel, if the step is taken. Shape is (Np_water, 1), with each element recording the flat index into the domain shape.

• _step – Step number of water parcels.

• L0 – Domain shape parameter, number of cells inlet length.

• CTR – Domain shape parameter, index along inlet wall making the center of the domain. I.e., (0, CTR) is the midpoint across the inlet, along the inlet domain edge.

• stage_above_SL – Water surface elevation minus the domain sea level.

Returns:

• new_inds – An updated array of parcel indicies, where the index of a parcel has been changed, if and only if, that parcel was looped.

• looped – A binary integer array indicating whether a parcel was determined to have been looped, and should be disqualified from the free surface computation.

Examples

The following shows an example of how water parcels that looped along their paths would be relocated. Note than in this example, the parcels are artificially forced to loop, just for the sake of demonstration.

(png, hires.png)

pyDeltaRCM.water_tools._update_dirQfield(qfield, dist, inds, astep, dirstep)

Update unit vector of water flux in x or y.

pyDeltaRCM.water_tools._update_absQfield(qfield, dist, inds, astep, Qp_water, dx)

Update norm of water flux vector.

pyDeltaRCM.water_tools._accumulate_free_surface_walks(free_surf_walk_inds: ndarray, free_surf_flag: ndarray, cell_type: ndarray, uw: ndarray, ux: ndarray, uy: ndarray, depth: ndarray, dx: float, u0: float, h0: float, H_SL: float, S0: float)

Accumulate the free surface by walking parcel paths.

This routine comprises the hydrodynamic physics-based computations.

Algorithm is to:

1. loop through every parcel’s directed random walk in series.

2. for a parcel’s walk, unravel the indices and determine whether the parcel should contribute to the free surface. Parcels are considered contributors if they have reached the ocean and if they are not looped pathways.

3. then, we begin at the downstream end of the parcel’s walk and iterate up-walk until, determining the Hnew for each location. Downstream of the shoreline-ocean boundary, the water surface elevation is set to the sea level. Upstream of the shoreline-ocean boundary, the water surface is determined according to the land-slope (S0) and the parcel pathway.

4. repeat from 2, for each parcel.

Examples

The following shows an example of the walk of a few water parcels, along with the resultant computed water surface.

(png, hires.png)

pyDeltaRCM.water_tools._smooth_free_surface(Hin, cell_type, Nsmooth: int, Csmooth: int)

Smooth the free surface.

Parameters:

• Hin – Stage input to the smoothing (i.e., the old stage).

• cell_type – Type of each cell.

• Nsmooth – Number of times to smooth the free surface (i.e., stage)

• Csmooth – Underrelaxation coefficient for smoothing iterations.

Examples

The following shows an example of the water surface smoothing process.

(png, hires.png)