optimize_splitting_mask(active_only=True) produces sub-models too large to run

[jdhughes-dev]

optimize_splitting_mask(active_only=True) partitions the active cells, then
assigns the inactive cells to whichever part is nearest
(model_splitter.py:741-744, NearestNeighbors(n_neighbors=1)). Proximity ignores
the shape of the result, so a coastal part is given a tongue of ocean and the
bounding box of the resulting sub-model grows far beyond the cells it solves.
The sub-models it produces can be too large to run.

For a 210,771,902 cell CONUS model split into 256 parts, active_only=True
partitions the active cells more compactly than active_only=False, and still
produces the larger model:

	sum of sub-model grids	fill
bounding boxes of the active cells alone, active_only=True	205,700,249	60.7%
written sub-models, active_only=False	311,005,557	40.2%
written sub-models, active_only=True	381,616,033	32.7%

The first row is what the partition would cost if the inactive cells were
assigned without inflating the boxes; the third is what is written. The
difference is entirely the nearest neighbor fill.

MODFLOW 6 could not run the active_only=True sub-models. It was OOM killed
during GWF-GWF exchange setup at about 28 GB per rank, against 7 GB per rank for
the active_only=False split of the same model at the same rank count, and
spreading the ranks over twice as many nodes only moved the failure later.
Active cell counts are identical in both, so the extra memory is padding.

Assigning inactive cells so as to minimize the growth of each part's bounding
box, rather than by proximity, would keep the sub-models near the 205,700,249
cells the partition implies. Cells that fall outside every part's active
bounding box are the ones that need a rule.

Worth noting the trade in the meantime: active_only=False gives compact
bounding boxes but active footprints that are not contiguous, because METIS is
free to route a part through zero weight inactive cells to stay connected. In
this model 14 of the 256 sub-models were in two or more geographic pieces, one
of them 1.86% full with a 934 MB DIS file. active_only=True fixes the
contiguity and loses the compactness.

Versions: flopy 3.10.0 (develop at 7732f39 is unchanged in these paths),
python 3.12.3, numpy 1.26.4, MODFLOW 6.8.0.

Related to #2813, which is what makes working around this expensive.

[jlarsen-usgs]

@jdhughes-dev

I think we could have a default "no-model/no-data" mapping value that we set for using active_only=True. This would reduce the complexity of the code, preserve compactness, and preserve contiguity.

[jdhughes-dev]

Correcting the closing suggestion above: minimizing bounding box growth cannot
reach the 205,700,249 cells that the active partition implies, because most of
the inactive cells are nowhere near it.

For the same model and partition:

	cells
sum of per-part bounding boxes of the active cells	205,700,249
inactive cells falling inside one of those boxes	9,557,062
inactive cells falling outside every one of them	76,330,257

The 76,330,257 cells are ocean, Canada and Mexico, 36.2% of the grid. Every one
of them has to be given to some part, and each one that lands outside a part's
box extends it. A greedy assignment would land somewhere between 205,700,249 and
the 381,616,033 written today, and choosing well is closer to bin packing than
to a proximity rule.

The more general point is that contiguity and low padding are in tension for a
structured sub-model: a contiguous active footprint is elongated, and the
bounding box of an elongated region carries more inactive cells than that of a
compact one. active_only=False buys compact boxes with non-contiguous
footprints; active_only=True buys contiguity with loose boxes. Neither gives
both.

[jdhughes-dev]

See #2816 for an approach that avoids the tension entirely: DISV sub-models hold only the cells they solve, so there is no bounding box to inflate.