Add adaptive flow direction with accumulation-dependent exponent

[brendancol]

Problem

D8, Dinf, and MFD each use the same routing logic everywhere on the raster. D8 is right in channels but draws artificial parallel lines on hillslopes. MFD is right on hillslopes but smears flow in valleys where it should converge. You have to pick one and accept its tradeoffs everywhere.

Proposed algorithm

Two-pass adaptive MFD that varies the partition exponent per cell based on upstream accumulation:

1. Run standard MFD accumulation with a low exponent (e.g. p=1.0) to get an initial accumulation grid.
2. Recompute flow fractions with a per-cell exponent derived from that accumulation:

p(a) = p_min + (p_max - p_min) * min(a / a_threshold, 1.0)

• a = upstream accumulation from pass 1
• p_min (default ~1.0) = dispersion on hillslopes
• p_max (default ~8.0) = concentration in channels
• a_threshold = accumulation at which routing becomes fully D8-like

Low-accumulation cells (hilltops, hillslopes) get a low exponent and spread flow. High-accumulation cells (valleys, channels) get a high exponent and send nearly all flow to the steepest neighbor.

Implementation notes

• flow_direction_mfd.py already computes weight = (slope * contour_length) ^ p. The main change is making p a 2D array instead of a scalar.
• CPU and GPU kernels need a minor signature change to accept a p_grid alongside elevation.
• Dask: accumulation must finish before the second pass, but each pass parallelizes the same way MFD does now.

Function signature

def flow_direction_adaptive(
    agg: xr.DataArray,
    p_min: float = 1.0,
    p_max: float = 8.0,
    a_threshold: float = 100.0,
    name: str = 'flow_direction_adaptive',
    boundary: str = 'nan',
) -> xr.DataArray

Returns the same 3D (8, H, W) fractional flow array as flow_direction_mfd.

References

• Holmgren, P. (1994). Multiple flow direction algorithms for runoff modelling in grid based elevation models: an empirical evaluation. Hydrological Processes, 8(4), 327-334.
• Qin, C. et al. (2007). An adaptive approach to selecting a flow-partition exponent for a multiple-flow-direction algorithm. International Journal of Geographical Information Science, 21(4), 443-458.