Add spatial hash describe function

[wyatt-fluidnumerics]

Description

Currently, it is difficult to get a sense for many of the dimensions of the spatial hash table that gets built over unstructured, and many structured grids. This makes understanding the performance speed up due to spatial hashing, or potential bugs in the spatial hashing such as that addressed by PR #2780, difficult. This PR adds functionality for .describe() to the SpatialHash class. The function is used exclusively for debugging or logging and provides useful statistics regarding the size of the hash grid and the ratio of hash grid cells to mesh grid faces.

Checklist

• This PR targets the correct branch (main for normal development, v3-support for v3 support)

AI Disclosure

• This PR contains AI-generated content.
  • I have tested any AI-generated content in my PR.
  • I take responsibility for any AI-generated content in my PR.
  • Describe how you used it (e.g., by pasting your prompt): I used Claude Code to generate an initial implementation of the spatialhash_describe function. I then made significant modifications to the resulting code to improve the utility of the reported statistics.

[erikvansebille]

Thanks for this PR, @wyatt-fluidnumerics. I like it a lot!

See below a screenshot for one of my Delft3D curvilinear meshes

Very useful information; as I now wonder why the max faces per occupied cell is more than 175k(!)

See below some suggestions/comments on this PR itself

[erikvansebille]

Very useful information; as I now wonder why the max faces per occupied cell is more than 175k(!)

I dug a bit deeper, and it turns out that counts.max() is for those hash cells where keys == 0.

I think/assume(?) that these are the hash cells that fall outside the model grid (i.e. no mesh faces occupy that hash cell), which makes sense because my grid is highly curved (see screenshot below)

So the question now is; should we only consider the counts where keys > 0 in this describe method?

If I do counts = hash_table["counts"][hash_table["keys"] > 0] on the grid below I get a max count of 44, which is much more reasonable than the 175k hash cells with keys==0!

[wyatt-fluidnumerics]

Interesting, 44 certainly seems a lot more reasonable, let me dig into this some myself as well. My understanding was that keys only stores hash cells that contain faces from the model grid's mesh, but this could be wrong.

[erikvansebille]

If hash_table["keys"] == 0 indeed represents all hash cells that don't contain any mesh face, then we could add that number as a separate line in the describe(). For example as "Hash cells without mesh faces"

And note that the computation of n_occupied_cells = int(hash_table["keys"].size) should then probably also ignore keys == 0

[wyatt-fluidnumerics]

I am pretty sure that every entry in keys corresponds to an occupied hash cell.
https://github.com/Parcels-code/Parcels/blob/main/src/parcels/_core/spatialhash.py#L339-L366
The keys for the hash table are constructed by indexing into an array of sorted morton codes at the location where that hash cell's entries in faces begin.

[erikvansebille]

OK, but if I do np.isnan(fieldset_wind.models[0].grid._spatialhash._xlow).sum() I get 175698, which is exactly the same number as the counts[0]...

So, could these be all the grid cells with NaN lon/lat locations from Delft3D?

[wyatt-fluidnumerics]

This would make sense, the actual meaning of key == 0 is the cell occupies the (0, 0, 0) cell on the unit cube the hash grid creates. So if NaN lon/lat get masked to 0 then I believe that would be expected. I don't totally see how this leads to an enormous number of faces.

[wyatt-fluidnumerics]

I think I found where the problem occurs. The code snippet below prints 0. So if you have a NaN value then it gets clipped to (0, 0, 0) coord equivalent int in quantize_coordinates

  import numpy as np
  with np.errstate(invalid="ignore"):
      xq = np.clip(np.array([np.nan]) * 1023, 0, 1023).astype(np.uint32)

  print(xq)

Normally I think this would produce a RuntimeWarning but those are silenced for the block of code.

[erikvansebille]

PR looks great; thanks!