Mask topsystem packages where IBOUND == -1

[JoerivanEngelen]

In issue #1766 it turned out that iMOD5 filters away RIV & DRN cells where IBOUND == -1 (and thus CHD cells are created), as these create shortcircuits in model cells, potentially leading to unnecessary stability issues.

We can approach this in two ways:

• Only implement this in the from_imod5 logic. We can then relatively easily use IBOUND == -1 as mask.
• Generalize this: Create a method on the groundwatermodel flow object to filter all GHB, RIV, & DRN cells where a CHD cell is located. MODFLOW6 only accepts one CHD package, so we don't require any logic to merge CHD packages together for this masking. This can then also be called when calling Modflow6Simulation.clip_box with states_for_boundary.

I prefer the latter option.