water mover package load

[cnicol-gwlogic]

Hi,

I'm using flopy (develop branch) to build, write and load an mf6 disu model, with water mover package. On loading mvr, I receive this error:

MFDataException: An error occurred in data element "perioddata" model "mymodel" package "mvr". The error occurred while setting array data in the "_build_recarray" method.
Additional Information:
(1) An error occurred when storing data "perioddata" in a recarray. perioddata data is a one or two dimensional list containing the variables "<pname1> <id1> <pname2> <id2> <mvrtype> <value>" (some variables may be optional, see MF6 documentation), but data "[('drn', (0,), 'maw', (0,), 'factor', 1.0), ('drn', (1,), 'maw', (0,), 'factor', 1.0), ('drn', (2,), 'maw', (0,), 'factor', 1.0)]" was supplied.
(2) Error occurred while loading data "perioddata" in block "period" from file "temp\t504_mf6disu\.\mymodel.mvr".

I have replicated the issue with a modified form of autotest t504 test_mf6disu(). See here.

If I change the mvr file perioddata id1 and id2 variables to be tuples (in the mvr modflow file...eg "(1,)" rather than "1"), flopy will load it ok. Mf6 however spits it on this as I would expect. Flopy writes the mvr correctly (no tuples), but as far as I can tell it is expecting some unusual recarray shape on load for the perioddata item.

Advice appreciated. Thanks a lot.

[christianlangevin]

Any insight from here? https://github.com/MODFLOW-USGS/modflow6/blob/4a64ea77a813036a9d86030bfd643f244f39fe0b/autotest/test_gwt_prudic2004t2fmiats.py#L387

[spaulins-usgs]

Looks like this is a bug with flopy which only occurs in the mvr package when the unstructured grid (disu) is being used. The bug has to do with the way flopy tries to determine if a particular column of data is a cellid (flopy treats cellids differently than other data). In this case flopy is incorrectly interpreting the mvr field as a cellid.

@langevin-usgs, I think the best solution is to add an optional "cellid" field to the dfns. All data fields that are a cellids would contain a line "cellid True" in the dfn. This way flopy does not have to deduce whether a particular field is a cellid.

[christianlangevin]

Don't we already call fields "cellid" that are cell ids? Do we ever have a data field with cellids, but we call the data field something else? We can add another field to the dfns, but it would be nice to avoid that unless it's absolutely necessary.

[spaulins-usgs]

We mostly call fields "cellid" with a few exceptions. The exceptions are:

OBS package "id" field
GNC package "cellidn", "cellidm", "cellidsj" fields
GWF exchange package "cellidm1", "cellidm2" fields

The bug I am trying to solve now is that both the OBS package and the MVR package have a field that starts with "id". With Flopy's current code it is interpreting starting with "id" as "might be a cellid" and then runs further checks to try and confirm this. Though now that I think about it adding a "cellid" field to the DFN will not solve resolve cellid for the OBS package, since in that case it may or may not be a cellid.

Perhaps the best rule to determine cellids is any name containing "cellid" somewhere in the name is a cellid. Then just treat the OBS package as a special case. However, the OBS package may see the same type of bug we are currently seeing in the MVR package. I still need a better rule to determine whether in the OBS package is a cellid or is something else.

[cnicol-gwlogic]

Thanks for the thinking on this. @spaulins-usgs - re: the obs issue, I think you're right that this is closely related. I (briefly) noticed that for a gwf node obs, flopy was not adding 1 to the cellid on obs file write (for this same disu model). I can replicate it in my example if you want.

Re: " I still need a better rule to determine whether in the OBS package is a cellid or is something else." - could you not parse the obstype field against a list of cellid types vs lak/sfr/whatever types? Or I guess to keep it generic / future-proof in flopy, maybe a user-defined field is needed in the obs package flagging cellid vs other index.

[spaulins-usgs]

@cnicol-gwlogic, I pushed a fix to the develop branch. This should fix the problem you are seeing with the MVR package and should also allow flopy to read/write OBS regardless of whether it is able to determine if is a cellid. However, there is still a potential problem where flopy may not always correctly determine whether an OBS is a cellid (vs a boundname). I am still working on finding the best approach to resolve this.

[cnicol-gwlogic]

@spaulins-usgs thanks so much, appreciate it. Worked perfectly FYI. Feel free to close this issue.