[Help request] Suggestion how to keep Faults (HFB) and aquifer/aquitard parameters together during regridding

[FransRoelofsen]

From Sprint session 7/5/2026.
Scale of aquifer parameters often differs from target model scale. Faults are rasterizes on target model. This can result in edge cases where HFB is not overlapping Geohydrology. This can cause vertical leakage of groundwater near faults.

[Huite]

Do you have a screenshot or picture to illustrate the problem?

I can imagine something like: there's a clay layer, and the fault connects exactly. So the overlying aquifer is "closed off" on the bottom by the clay layer, and, say, on the left by the fault.

Then when regridding, maybe some sand conductivity gets mixed into the aggregated conductivity. The fault is drawn to some nearest cell boundary, e.g. on the left of the aggregated conductivity cell. In that case, the bottom is no longer as closed off as it was.

There is a workaround, I guess, and it looks rather tedious. I'm also unsure whether it can easily be automated.
For the regridding, the workaround is to exlude values if they lie beyond the fault. So one would have to mark cells being "left of" the fault and "right of". This should occur for both the original grid and the original lines; and the new grid and the newly snapped lines.

We have some utilities to find the affected cells. We can find which cells are connected to a particular cell edge (where the fault is snapped to). But this isn't really enough, I think.
If I sketch in my note book a bit, it seems like the important data are: which cells are "swept" by (re-)snapping the fault geometry?
This marks all cells which could potentially cause problems.

Given that we have those cells, there might be an easy solution. They should be ignored from the perspective of maintaining continuity for the fualt. The regridder methods are robust in dealing with NaN as NoData. So if we mark them and set them to NoData, then the smearing shouldn't happen.

Next question is: how to find the "swept" cells easily? This might have an easy answer. We have a pre- and post-snap geometry to work with, which both consist of edges of two vertices. The swept area for an edge is the convex hull of the four total vertices.

To illustrate it:

Let's say we have 10 m grid cells originally. The red line marks the fault on e.g. 30 meter from the center line.
In the new grid we have 100 m cells, and the center line forms a cell boundary.

The problem area is marked by sweeping the edges from the original position to the newly snapped position:

Because the cells on the right, in normal regridding operations, will consist of a mix of clay and sand, and will therefore leak quite a bit in the new 100 m discretization.
Anyway, by connecting the old and new edge vertices, we should be able to fully identify this area:

These are simple convex polygons, so we can find them with the celltree spatial index.

Next we eliminate the values in the original grid:

And then regrid to find:

Which remains hydraulically closed on the right.

I'm not 100% sure whether this is entirely robust or whether it might create NoData holes -- I don't think so, because the snapping occurs from old grid cell edge to new grid cell edge. So there should always be data available (even if it's only a single value).
At any rate, if it does occur, it should simple to treat: also regrid the with the original data, then do a combine first:

masked_data = data.where(~swept)
out1 = regridder.regrid(masked_data)
out2 = regridder.regrid(data)
out = out1.combine_first(out2)

We could test whether this is robust; if so then it's a fairly easy fix. Although it does require custom logic, since it means the NPF package has to be HFB aware (which is managable).