Skip to content

Commit

Permalink
fix(mf6model.write_input): handle cases of single model or multi-mode…
Browse files Browse the repository at this point in the history 
…l LGR simulation by working with the simulation-level model dictionary to implement model write functionality that is not included in the simulation-level write_simulation(). Previously, BCs were presumably not getting removed from inactive cells, auxillary header information such as modflow-setup and model version were not getting written to child model package files, and child model SFR packagedata was not getting written to an external file.
  • Loading branch information
@aleaf
aleaf committed Aug 27, 2024
1 parent 85b7816 commit 105f42d
Show file tree
Hide file tree
Showing 2 changed files with 67 additions and 53 deletions.
108 changes: 57 additions & 51 deletions mfsetup/mf6model.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -970,66 +970,72 @@ def write_input(self):
"""
# prior to writing output
# remove any BCs in inactive cells
pckgs = ['chd', 'drn', 'ghb', 'riv', 'wel']
for pckg in pckgs:
package_instance = getattr(self, pckg.lower(), None)
if package_instance is not None:
external_files = self.cfg[pckg.lower()]['stress_period_data']
remove_inactive_bcs(package_instance,
external_files=external_files)
if hasattr(self, 'obs'):
for obs_package_instance in self.obs:
remove_inactive_obs(obs_package_instance)

# write the model with flopy
# but skip the sfr package
# by monkey-patching the write method
def skip_write(**kwargs):
pass
if hasattr(self, 'sfr'):
self.sfr.write = skip_write
# handle cases of single model or multi-model LGR simulation
# by working with the simulation-level model dictionary
for model_name, model in self.simulation.model_dict.items():
pckgs = ['chd', 'drn', 'ghb', 'riv', 'wel']
for pckg in pckgs:
package_instance = getattr(model, pckg.lower(), None)
if package_instance is not None:
external_files = model.cfg[pckg.lower()]['stress_period_data']
remove_inactive_bcs(package_instance,
external_files=external_files)
if hasattr(model, 'obs'):
for obs_package_instance in model.obs:
remove_inactive_obs(obs_package_instance)

# write the model with flopy
# but skip the sfr package
# by monkey-patching the write method
def skip_write(**kwargs):
pass
if hasattr(model, 'sfr'):
model.sfr.write = skip_write
self.simulation.write_simulation()

# write the sfr package with SFRmaker
if 'SFR' in ' '.join(self.get_package_list()):
options = []
for k, b in self.cfg['sfr']['options'].items():
if k == 'mover':
if 'mvr' not in self.simulation.package_key_dict:
continue
options.append(k)
if 'save_flows' in options:
budget_fileout = '{}.{}'.format(self.name,
self.cfg['sfr']['budget_fileout'])
stage_fileout = '{}.{}'.format(self.name,
self.cfg['sfr']['stage_fileout'])
options.append('budget fileout {}'.format(budget_fileout))
options.append('stage fileout {}'.format(stage_fileout))
if len(self.sfrdata.observations) > 0:
options.append('obs6 filein {}.{}'.format(self.name,
self.cfg['sfr']['obs6_filein_fmt'])
)
self.sfrdata.write_package(idomain=self.idomain,
version='mf6',
options=options,
external_files_path=self.external_path
)
# post-flopy write actions
for model_name, model in self.simulation.model_dict.items():
# write the sfr package with SFRmaker
if 'SFR' in ' '.join(model.get_package_list()):
options = []
for k, b in model.cfg['sfr']['options'].items():
if k == 'mover':
if 'mvr' not in model.simulation.package_key_dict:
continue
options.append(k)
if 'save_flows' in options:
budget_fileout = '{}.{}'.format(model_name,
model.cfg['sfr']['budget_fileout'])
stage_fileout = '{}.{}'.format(model_name,
model.cfg['sfr']['stage_fileout'])
options.append('budget fileout {}'.format(budget_fileout))
options.append('stage fileout {}'.format(stage_fileout))
if len(model.sfrdata.observations) > 0:
options.append('obs6 filein {}.{}'.format(model_name,
model.cfg['sfr']['obs6_filein_fmt'])
)
model.sfrdata.write_package(idomain=model.idomain,
version='mf6',
options=options,
external_files_path=model.external_path
)
# add version info to package file headers
files = [model.namefile]
files += [p.filename for p in model.packagelist]
files += [p[0].filename for k, p in model.simulation.package_key_dict.items()]
for f in files:
add_version_to_fileheader(f, model_info=model.header)

if not model.cfg['mfsetup_options']['keep_original_arrays']:
shutil.rmtree(model.tmpdir)

# label stress periods in tdis file with comments
self.perioddata.sort_values(by='per', inplace=True)
add_date_comments_to_tdis(self.simulation.tdis.filename,
self.perioddata.start_datetime,
self.perioddata.end_datetime
)
# add version info to file headers
files = [self.namefile]
files += [p.filename for p in self.packagelist]
files += [p[0].filename for k, p in self.simulation.package_key_dict.items()]
for f in files:
add_version_to_fileheader(f, model_info=self.header)

if not self.cfg['mfsetup_options']['keep_original_arrays']:
shutil.rmtree(self.tmpdir)



@staticmethod
Expand Down
12 changes: 10 additions & 2 deletions mfsetup/tests/test_lgr.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -191,6 +191,8 @@ def test_lgr_grid_setup(lgr_grid_spacing, layer_refinement_input,
m.modelgrid.write_shapefile(outfolder / 'pleasant_lgr_parent_grid.shp')
inset_model.modelgrid.write_shapefile(outfolder / 'pleasant_lgr_inset_grid.shp')

# convert layer_refinement to a list
# (mimicking what is done internally in modflow-setup)
if np.isscalar(layer_refinement):
layer_refinement = np.array([1] * m.modelgrid.nlay)
elif isinstance(layer_refinement, dict):
Expand Down Expand Up @@ -289,13 +291,19 @@ def test_lgr_parent_bcs_in_lgr_area(pleasant_vertical_lgr_setup_from_yaml):

def test_mover_get_sfr_package_connections(pleasant_lgr_setup_from_yaml):
m = pleasant_lgr_setup_from_yaml
inset_model = m.inset['plsnt_lgr_inset']

# check that packagedata for both models was written to external files
assert Path(m.external_path, f"{m.name}_packagedata.dat").exists()
assert Path(inset_model.external_path, f"{inset_model.name}_packagedata.dat").exists()

parent_reach_data = m.sfrdata.reach_data
inset_reach_data = m.inset['plsnt_lgr_inset'].sfrdata.reach_data
inset_reach_data = inset_model.sfrdata.reach_data
to_inset, to_parent = get_sfr_package_connections(parent_reach_data, inset_reach_data,
distance_threshold=200)
assert len(to_inset) == 0
# verify that the last reaches in the two segments are keys
last_reaches = m.inset['plsnt_lgr_inset'].sfrdata.reach_data.groupby('iseg').last().rno
last_reaches = inset_model.sfrdata.reach_data.groupby('iseg').last().rno
assert not any(set(to_parent.keys()).difference(last_reaches))
# verify that the first reaches are headwaters
outreaches = set(m.sfrdata.reach_data.outreach)
Expand Down

0 comments on commit 105f42d

Please sign in to comment.