modflowpy / flopy

A Python package to create, run, and post-process MODFLOW-based models.
https://flopy.readthedocs.io
Other
517 stars 313 forks source link

bug: `flopy.modpath.ParticleData` does not accept regular numpy array as input #1628

Closed mbakker7 closed 1 year ago

mbakker7 commented 1 year ago

The flopy.modpath.ParticleData class does not take a regular numpy array with three columns as input.

To reproduce:

import numpy as np
import flopy as fp

partlocs = np.array([(1, 1, 1),
                     (1, 1, 2)])

fp.modpath.ParticleData(partlocs=partlocs, structured=True)

The error message is ValueError: could not broadcast input array from shape (2,3) into shape (2,). The bug seems to be in mp7particledata.py on line 168, where the dtype of the input array is changed:

        elif isinstance(partlocs, np.ndarray):
            dtypein = partlocs.dtype
            if dtypein != dtype:
                partlocs = np.array(partlocs, dtype=dtype)

This conversion to a new dtype does not give the desired result as can be seen:

newdtype = np.dtype([('k', '<i4'), ('i', '<i4'), ('j', '<i4')])
np.array(partlocs, dtype=newdtype)

gives

array([[(1, 1, 1), (1, 1, 1), (1, 1, 1)],
       [(1, 1, 1), (1, 1, 1), (2, 2, 2)]],
      dtype=[('k', '<i4'), ('i', '<i4'), ('j', '<i4')])

The correct conversion requires something (arguably ugly) as

from numpy.lib.recfunctions import unstructured_to_structured
partlocs = unstructured_to_structured(partlocs, dtype=newdtype)

which gives the desired

array([(1, 1, 1), (1, 1, 2)],
      dtype=[('k', '<i4'), ('i', '<i4'), ('j', '<i4')])

Incidentally, it appears that flopy.modpath.ParticleData also doesn't accept a list of lists such as

partlocs = [[1, 1, 1], [1, 1, 2]]

which gives the same error. A list of tuples works fine, however.

wpbonelli commented 1 year ago

Thank you @mbakker7 the unstructured_to_structured solution was added in #1752, this should now work correctly