Skip to content
New issue

Have a question about this project? # for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “#”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? # to your account

bug: modelgrid coord info dropped when using a DISV grid #2283

Closed
cneyens opened this issue Aug 7, 2024 · 2 comments · Fixed by #2284
Closed

bug: modelgrid coord info dropped when using a DISV grid #2283

cneyens opened this issue Aug 7, 2024 · 2 comments · Fixed by #2284
Labels
Milestone

Comments

@cneyens
Copy link

cneyens commented Aug 7, 2024

Describe the bug
When obtaining the modelgrid object from a model with a DISV package, the coord info present in the DISV object is not passed properly to the modelgrid object. For models with DIS, this works as expected; I haven't checked for DISU grids.

To Reproduce
Example model with a triangular grid:

Example code
import os
import numpy as np
import flopy
from flopy.utils.triangle import Triangle

# create triangular grid
triangle_ws = 'triangle'
if not os.path.exists(triangle_ws):
  os.mkdir(triangle_ws)
  
active_area = [(0, 0), (0, 1000), (1000, 1000), (1000, 0)]
tri = Triangle(model_ws=triangle_ws, angle=30)
tri.add_polygon(active_area)
tri.add_region((1, 1), maximum_area=50**2)

tri.build()

# build vertex grid object
vgrid = flopy.discretization.VertexGrid(
  vertices=tri.get_vertices(),
  cell2d=tri.get_cell2d(),
  xoff=199000,
  yoff=215500,
  crs=31370,
  angrot=30,
)

# coord info is set (also correct when using vgrid.set_coord_info()
print(vgrid)

# create MODFLOW 6 model
sim = flopy.mf6.MFSimulation(sim_name='prj-test', sim_ws='model')
tdis = flopy.mf6.ModflowTdis(sim)
ims = flopy.mf6.ModflowIms(sim)

gwf = flopy.mf6.ModflowGwf(sim, modelname='gwf')
disv = flopy.mf6.ModflowGwfdisv(
  gwf,
  xorigin=vgrid.xoffset,
  yorigin=vgrid.yoffset,
  angrot=vgrid.angrot, # no CRS info can be set in DISV
  nlay=1,
  top=0.0,
  botm=-10.0,
  ncpl=vgrid.ncpl,
  nvert=vgrid.nvert,
  cell2d=vgrid.cell2d,
  vertices=tri.get_vertices() # this is not stored in the Vertex grid object?
)

# coord info is dropped
# also the case when loading model from disk
print(gwf.modelgrid)

# disv object has coord info set properly though
print([disv.xorigin.get_data(), disv.yorigin.get_data(), disv.angrot.get_data()])

Expected behavior
I expect the modelgrid object to get the coord info from the DISV object when calling it through gwf.modelgrid, as is done for a DIS model. Pretty easy to get around this using set_coord_info() based on the DISV properties, but perhaps there's more to it.

Desktop (please complete the following information):

  • OS: Windows 10
  • Flopy v3.8.0
@cneyens cneyens added the bug label Aug 7, 2024
@martclanor
Copy link
Contributor

This seems to be caused by this block of code in flopy/mf6/mfmodel.py:575 (https://github.com/modflowpy/flopy/blob/develop/flopy/mf6/mfmodel.py#L575)

if self.get_grid_type() != DiscretizationType.DISV:
    # get coordinate data from dis file
    xorig = dis.xorigin.get_data()
    yorig = dis.yorigin.get_data()
    angrot = dis.angrot.get_data()
else:
    xorig = self._modelgrid.xoffset
    yorig = self._modelgrid.yoffset
    angrot = self._modelgrid.angrot

I am tempted to just remove this exception to the DISV case, but I'm not sure yet what it is intended for.

@langevin-usgs
Copy link
Contributor

I can't think of a reason why this was handled this way, except that xorigin, yorigin, and angrot weren't originally included with DISV and DISU.

image

# for free to join this conversation on GitHub. Already have an account? # to comment
Labels
Projects
None yet
Development

Successfully merging a pull request may close this issue.

4 participants