Skip to content

Commit

Permalink
fix(binary): fix binary header information
Browse files Browse the repository at this point in the history
Replaces PR #1848
  • Loading branch information
jdhughes-usgs committed Jul 16, 2023
1 parent 9168b3a commit d6c6a8f
Show file tree
Hide file tree
Showing 4 changed files with 179 additions and 65 deletions.
187 changes: 142 additions & 45 deletions autotest/test_mf6.py
Original file line number Diff line number Diff line change
Expand Up @@ -503,40 +503,88 @@ def test_subdir(function_tmpdir):


@requires_exe("mf6")
def test_binary_write(function_tmpdir):
@pytest.mark.parametrize("layered", [True, False])
def test_binary_write(function_tmpdir, layered):
nlay, nrow, ncol = 2, 1, 10
shape2d = (nrow, ncol)
shape3d = (nlay, nrow, ncol)

# data for layers
botm = [4.0, 0.0]
strt = [5.0, 10.0]

# create binary data structured
if layered:
idomain_data = []
botm_data = []
strt_data = []
for k in range(nlay):
idomain_data.append(
{
"factor": 1.0,
"filename": f"idomain_l{k+1}.bin",
"data": 1,
"binary": True,
"iprn": 1,
}
)
botm_data.append(
{
"filename": f"botm_l{k+1}.bin",
"binary": True,
"iprn": 1,
"data": np.full(shape2d, botm[k], dtype=float),
}
)
strt_data.append(
{
"filename": f"strt_l{k+1}.bin",
"binary": True,
"iprn": 1,
"data": np.full(shape2d, strt[k], dtype=float),
}
)
else:
idomain_data = {
"filename": "idomain.bin",
"binary": True,
"iprn": 1,
"data": 1,
}
botm_data = {
"filename": "botm.bin",
"binary": True,
"iprn": 1,
"data": np.array(
[
np.full(shape2d, botm[0], dtype=float),
np.full(shape2d, botm[1], dtype=float),
]
),
}
strt_data = {
"filename": "strt.bin",
"binary": True,
"iprn": 1,
"data": np.array(
[
np.full(shape2d, strt[0], dtype=float),
np.full(shape2d, strt[1], dtype=float),
]
),
}

# binary data that does not vary by layers
top_data = {
"filename": "top.bin",
"binary": True,
"iprn": 0,
"iprn": 1,
"data": 10.0,
}
botm_data = {
"filename": "botm.bin",
"binary": True,
"iprn": 0,
"data": np.array(
[
np.full(shape2d, 4.0, dtype=float),
np.full(shape2d, 0.0, dtype=float),
]
),
}
strt_data = {
"filename": "strt.bin",
"binary": True,
"iprn": 0,
"data": np.full(shape3d, 10.0, dtype=float),
}
rch_data = {
0: {
"filename": "recharge.bin",
"binary": True,
"iprn": 0,
"iprn": 1,
"data": 0.000001,
},
}
Expand All @@ -548,7 +596,7 @@ def test_binary_write(function_tmpdir):
0: {
"filename": "chd.bin",
"binary": True,
"iprn": 0,
"iprn": 1,
"data": chd_data,
},
}
Expand All @@ -566,14 +614,15 @@ def test_binary_write(function_tmpdir):
delc=1.0,
top=top_data,
botm=botm_data,
idomain=idomain_data,
)
ModflowGwfnpf(
gwf,
icelltype=1,
)
ModflowGwfic(
gwf,
strt=10.0,
strt=strt_data,
)
ModflowGwfchd(
gwf,
Expand All @@ -588,8 +637,8 @@ def test_binary_write(function_tmpdir):


@requires_exe("mf6")
@pytest.mark.skip(reason="todo:: after flopy binary fix.")
def test_vor_binary_write(function_tmpdir):
@pytest.mark.parametrize("layered", [True, False])
def test_vor_binary_write(function_tmpdir, layered):
# build voronoi grid
boundary = [(0.0, 0.0), (0.0, 1.0), (10.0, 1.0), (10.0, 0.0)]
triangle_ws = function_tmpdir / "triangle"
Expand All @@ -606,37 +655,84 @@ def test_vor_binary_write(function_tmpdir):

# problem dimensions
nlay = 2
shape3d = (nlay, vor.ncpl)

# data for layers
botm = [4.0, 0.0]
strt = [5.0, 10.0]

# build binary data
if layered:
idomain_data = []
botm_data = []
strt_data = []
for k in range(nlay):
idomain_data.append(
{
"factor": 1.0,
"filename": f"idomain_l{k + 1}.bin",
"data": 1,
"binary": True,
"iprn": 1,
}
)
botm_data.append(
{
"filename": f"botm_l{k + 1}.bin",
"binary": True,
"iprn": 1,
"data": np.full(vor.ncpl, botm[k], dtype=float),
}
)
strt_data.append(
{
"filename": f"strt_l{k + 1}.bin",
"binary": True,
"iprn": 1,
"data": np.full(vor.ncpl, strt[k], dtype=float),
}
)
else:
idomain_data = {
"filename": "idomain.bin",
"binary": True,
"iprn": 1,
"data": 1,
}
botm_data = {
"filename": "botm.bin",
"binary": True,
"iprn": 1,
"data": np.array(
[
np.full(vor.ncpl, botm[0], dtype=float),
np.full(vor.ncpl, botm[1], dtype=float),
]
),
}
strt_data = {
"filename": "strt.bin",
"binary": True,
"iprn": 1,
"data": np.array(
[
np.full(vor.ncpl, strt[0], dtype=float),
np.full(vor.ncpl, strt[1], dtype=float),
]
),
}

# binary data that does not vary by layers
top_data = {
"filename": "top.bin",
"binary": True,
"iprn": 0,
"iprn": 1,
"data": 10.0,
}
botm_data = {
"filename": "botm.bin",
"binary": True,
"iprn": 0,
"data": np.array(
[
np.full(vor.ncpl, 4.0, dtype=float),
np.full(vor.ncpl, 0.0, dtype=float),
]
),
}
strt_data = {
"filename": "strt.bin",
"binary": True,
"iprn": 0,
"data": np.full(shape3d, 10.0, dtype=float),
}
rch_data = {
0: {
"filename": "recharge.bin",
"binary": True,
"iprn": 0,
"iprn": 1,
"data": np.full(vor.ncpl, 0.000001, dtype=float), # 0.000001,
},
}
Expand Down Expand Up @@ -668,6 +764,7 @@ def test_vor_binary_write(function_tmpdir):
cell2d=vor.get_disv_gridprops()["cell2d"],
top=top_data,
botm=botm_data,
idomain=idomain_data,
xorigin=0.0,
yorigin=0.0,
)
Expand Down
44 changes: 30 additions & 14 deletions flopy/mf6/data/mffileaccess.py
Original file line number Diff line number Diff line change
Expand Up @@ -212,6 +212,8 @@ def write_binary_file(
):
data = self._resolve_cellid_numbers_to_file(data)
fd = self._open_ext_file(fname, binary=True, write=True)
if data.size == modelgrid.nnodes:
write_multi_layer = False
if write_multi_layer:
for layer, value in enumerate(data):
self._write_layer(
Expand Down Expand Up @@ -252,7 +254,7 @@ def _write_layer(
ilay=None,
):
header_data = self._get_header(
modelgrid, modeltime, stress_period, precision, text, fname, ilay
modelgrid, modeltime, stress_period, precision, text, fname, ilay=ilay, data=data,
)
header_data.tofile(fd)
data.tofile(fd)
Expand All @@ -266,6 +268,7 @@ def _get_header(
text,
fname,
ilay=None,
data=None,
):
# handle dis (row, col, lay), disv (ncpl, lay), and disu (nodes) cases
if modelgrid is not None and modeltime is not None:
Expand All @@ -274,13 +277,18 @@ def _get_header(
if ilay is None:
ilay = modelgrid.nlay
if modelgrid.grid_type == "structured":
m1, m2, m3 = modelgrid.ncol, modelgrid.nrow, ilay
if data is not None:
shape3d = modelgrid.nlay * modelgrid.nrow * modelgrid.ncol
if data.size == shape3d:
m1, m2, m3 = shape3d, 1, 1
return BinaryHeader.create(
bintype="vardis",
precision=precision,
text=text,
nrow=modelgrid.nrow,
ncol=modelgrid.ncol,
ilay=ilay,
m1=m1,
m2=m2,
m3-m3,
pertim=pertim,
totim=totim,
kstp=1,
Expand All @@ -289,24 +297,30 @@ def _get_header(
elif modelgrid.grid_type == "vertex":
if ilay is None:
ilay = modelgrid.nlay
m1, m2, m3 = modelgrid.ncpl, 1, ilay
if data is not None:
shape3d = modelgrid.nlay * modelgrid.ncpl
if data.size == shape3d:
m1, m2, m3 = shape3d, 1, 1
return BinaryHeader.create(
bintype="vardisv",
precision=precision,
text=text,
ncpl=modelgrid.ncpl,
ilay=ilay,
m3=1,
m1=m1,
m2=m2,
m3=m3,
pertim=pertim,
totim=totim,
kstp=1,
kper=stress_period,
)
elif modelgrid.grid_type == "unstructured":
m1, m2, m3 = modelgrid.nnodes, 1, 1
return BinaryHeader.create(
bintype="vardisu",
precision=precision,
text=text,
nodes=modelgrid.nnodes,
m1=m1,
m2=1,
m3=1,
pertim=pertim,
Expand All @@ -317,13 +331,14 @@ def _get_header(
else:
if ilay is None:
ilay = 1
m1, m2, m3 = 1, 1, ilay
header = BinaryHeader.create(
bintype="vardis",
precision=precision,
text=text,
nrow=1,
ncol=1,
ilay=ilay,
m1=m1,
m2=m2,
m3=m3,
pertim=pertim,
totim=totim,
kstp=1,
Expand All @@ -339,14 +354,15 @@ def _get_header(
"binary file {}.".format(fname)
)
else:
m1, m2, m3 = 1, 1, 1
pertim = np.float64(1.0)
header = BinaryHeader.create(
bintype="vardis",
precision=precision,
text=text,
nrow=1,
ncol=1,
ilay=1,
m1=m1,
m2=m2,
m3=m3,
pertim=pertim,
totim=pertim,
kstp=1,
Expand Down
1 change: 1 addition & 0 deletions flopy/utils/binaryfile.py
Original file line number Diff line number Diff line change
Expand Up @@ -194,6 +194,7 @@ def set_values(self, **kwargs):
"ilay",
"ncpl",
"nodes",
"m1",
"m2",
"m3",
]
Expand Down
Loading

0 comments on commit d6c6a8f

Please sign in to comment.