-
Notifications
You must be signed in to change notification settings - Fork 125
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
fix(ObsModule): implemented non-advancing output (#55)
Modified the pbservation process to use non-advancing output instead of fixed length strings when writing ascii output. The previous use of fixed length strings resulted in truncation of ascii observation output when the product of user-specified digits + 7 and the number of observations exceeded 5000. Closes #54
- Loading branch information
1 parent
df05c2e
commit cfe6e02
Showing
10 changed files
with
287 additions
and
35 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,244 @@ | ||
import os | ||
import sys | ||
import numpy as np | ||
|
||
try: | ||
import flopy | ||
except: | ||
msg = 'Error. FloPy package is not available.\n' | ||
msg += 'Try installing using the following command:\n' | ||
msg += ' pip install flopy' | ||
raise Exception(msg) | ||
|
||
from framework import testing_framework | ||
from simulation import Simulation | ||
|
||
ex = ['utl03_obs'] | ||
exdirs = [] | ||
for s in ex: | ||
exdirs.append(os.path.join('temp', s)) | ||
ddir = 'data' | ||
|
||
# temporal discretization | ||
nper = 2 | ||
tdis_rc = [(1., 1, 1.), | ||
(1., 1, 1.)] | ||
|
||
# spatial discretization data | ||
nlay, nrow, ncol = 1, 21, 21 | ||
top = 200. | ||
botm = 0. | ||
|
||
# hydraulic properties | ||
laytyp = 1 | ||
hk = 1. | ||
delr = delc = 100. | ||
strt = botm + 20. | ||
|
||
# constant head data | ||
chdlocl = [(0, i, 0) for i in range(nrow)] | ||
chdlocr = [(0, i, ncol - 1) for i in range(nrow)] | ||
chdl = 100. | ||
chdr = 75. | ||
|
||
c60 = [] | ||
c61 = [] | ||
for loc in chdlocl: | ||
c60.append([loc, chdl]) | ||
c61.append([loc, chdl + 10]) | ||
for loc in chdlocr: | ||
c60.append([loc, chdr]) | ||
c61.append([loc, chdr - 10]) | ||
cd6 = {0: c60, 1: c61} | ||
|
||
# gwf obs | ||
obs_data0 = [('h{:04d}'.format(i + 1), 'HEAD', (0, 10, 10)) for i in | ||
range(1000)] | ||
obs_data1 = [('h{:04d}'.format(i + 1001), 'HEAD', (0, 1, 1)) for i in | ||
range(737)] | ||
|
||
# solver data | ||
nouter, ninner = 100, 300 | ||
hclose, rclose, relax = 1e-6, 0.01, 1. | ||
|
||
|
||
def build_mf6(idx, ws, binaryobs=True): | ||
name = ex[idx] | ||
|
||
# build MODFLOW 6 files | ||
sim = flopy.mf6.MFSimulation(sim_name=name, version='mf6', | ||
exe_name='mf6', | ||
sim_ws=ws) | ||
# create tdis package | ||
flopy.mf6.ModflowTdis(sim, time_units='DAYS', | ||
nper=nper, perioddata=tdis_rc) | ||
|
||
# create gwf model | ||
gwf = flopy.mf6.ModflowGwf(sim, modelname=name, | ||
model_nam_file='{}.nam'.format(name), | ||
save_flows=True) | ||
|
||
# create iterative model solution and register the gwf model with it | ||
ims = flopy.mf6.ModflowIms(sim, print_option='SUMMARY', | ||
outer_hclose=hclose, | ||
outer_maximum=nouter, | ||
under_relaxation='none', | ||
backtracking_number=0, | ||
inner_maximum=ninner, | ||
inner_hclose=hclose, | ||
rcloserecord=rclose, | ||
linear_acceleration='CG', | ||
scaling_method='NONE', | ||
reordering_method='NONE', | ||
relaxation_factor=relax, | ||
number_orthogonalizations=7) | ||
sim.register_ims_package(ims, [gwf.name]) | ||
|
||
flopy.mf6.ModflowGwfdis(gwf, nlay=nlay, nrow=nrow, ncol=ncol, | ||
delr=delr, delc=delc, | ||
top=top, botm=botm, | ||
idomain=1, | ||
fname='{}.dis'.format(name)) | ||
|
||
# initial conditions | ||
flopy.mf6.ModflowGwfic(gwf, strt=strt, | ||
fname='{}.ic'.format(name)) | ||
|
||
# node property flow | ||
flopy.mf6.ModflowGwfnpf(gwf, icelltype=1, k=hk) | ||
|
||
# gwf head observations | ||
if binaryobs: | ||
obs_recarray = {'head0.obs.bsv': obs_data0, | ||
'head1.obs.bsv': obs_data1} | ||
else: | ||
obs_recarray = {'head0.obs.csv': obs_data0, | ||
'head1.obs.csv': obs_data1} | ||
|
||
o = flopy.mf6.ModflowUtlobs(gwf, pname='head_obs', | ||
fname='{}.obs'.format(name), | ||
digits=10, print_input=True, | ||
continuous=obs_recarray) | ||
|
||
# chd files | ||
flopy.mf6.modflow.ModflowGwfchd(gwf, stress_period_data=cd6) | ||
|
||
# output control | ||
flopy.mf6.ModflowGwfoc(gwf, | ||
budget_filerecord='{}.cbc'.format(name), | ||
head_filerecord='{}.hds'.format(name), | ||
headprintrecord=[('COLUMNS', 10, 'WIDTH', 15, | ||
'DIGITS', 6, 'GENERAL')], | ||
saverecord=[('HEAD', 'LAST')], | ||
printrecord=[('HEAD', 'LAST'), | ||
('BUDGET', 'LAST')]) | ||
|
||
return sim | ||
|
||
|
||
def get_model(idx, dir): | ||
ws = dir | ||
# build mf6 with ascii observation output | ||
sim = build_mf6(idx, ws, binaryobs=False) | ||
|
||
# build mf6 with binary observation output | ||
wsc = os.path.join(ws, 'mf6') | ||
mc = build_mf6(idx, wsc, binaryobs=True) | ||
|
||
return sim, mc | ||
|
||
|
||
def build_models(): | ||
for idx, dir in enumerate(exdirs): | ||
sim, mc = get_model(idx, dir) | ||
sim.write_simulation() | ||
mc.write_simulation() | ||
hack_binary_obs(idx, dir) | ||
return | ||
|
||
|
||
def hack_binary_obs(idx, dir): | ||
name = ex[idx] | ||
ws = dir | ||
wsc = os.path.join(ws, 'mf6') | ||
fname = name + '.obs' | ||
fpth = os.path.join(wsc, fname) | ||
with open(fpth, 'r') as f: | ||
lines = f.readlines() | ||
with open(fpth, 'w') as f: | ||
for line in lines: | ||
line = line.rstrip() | ||
if 'BEGIN continuous FILEOUT' in line: | ||
line += ' BINARY' | ||
f.write('{}\n'.format(line)) | ||
f.close() | ||
return | ||
|
||
def eval_obs(sim): | ||
print('evaluating observations...') | ||
|
||
# get results from the observation files | ||
pth = sim.simpath | ||
files = [fn for fn in os.listdir(pth) if '.csv' in fn] | ||
for file in files: | ||
pth0 = os.path.join(pth, file) | ||
pth1 = os.path.join(pth, 'mf6', file.replace('.csv', '.bsv')) | ||
d0 = flopy.utils.Mf6Obs(pth0, isBinary=False).get_data() | ||
d1 = flopy.utils.Mf6Obs(pth1, isBinary=True).get_data() | ||
names0 = d0.dtype.names | ||
names1 = d1.dtype.names | ||
msg = 'The number of columns ({}) '.format(len(names0)) + \ | ||
'in {} '.format(pth0) + 'is not equal to ' + \ | ||
'the number of columns ({}) '.format(len(names1)) + \ | ||
'in {}.'.format(pth1) | ||
assert len(names0) == len(names1), msg | ||
msg = 'The number of rows ({}) '.format(d0.shape[0]) + \ | ||
'in {} '.format(pth0) + 'is not equal to ' + \ | ||
'the number of rows ({}) '.format(d1.shape[0]) + \ | ||
'in {}.'.format(pth1) | ||
assert d0.shape[0] == d1.shape[0], msg | ||
for name in names0: | ||
msg = "The values for column '{}' ".format(name) + \ | ||
'are not within 1e-5 of each other' | ||
assert np.allclose(d0[name], d1[name], rtol=1e-5), msg | ||
|
||
return | ||
|
||
|
||
|
||
# - No need to change any code below | ||
def test_mf6model(): | ||
# initialize testing framework | ||
test = testing_framework() | ||
|
||
# build the models | ||
build_models() | ||
|
||
# run the test models | ||
for dir in exdirs: | ||
yield test.run_mf6, Simulation(dir, exfunc=eval_obs) | ||
|
||
return | ||
|
||
|
||
def main(): | ||
# initialize testing framework | ||
test = testing_framework() | ||
|
||
# build the models | ||
build_models() | ||
|
||
# run the test models | ||
for dir in exdirs: | ||
sim = Simulation(dir, exfunc=eval_obs) | ||
test.run_mf6(sim) | ||
|
||
return | ||
|
||
|
||
if __name__ == "__main__": | ||
# print message | ||
print('standalone run of {}'.format(os.path.basename(__file__))) | ||
|
||
# run main routine | ||
main() |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -1,3 +1,3 @@ | ||
\newcommand{\modflowversion}{mf6.0.3.30} | ||
\newcommand{\modflowdate}{October 29, 2018} | ||
\newcommand{\modflowversion}{mf6.0.3.31} | ||
\newcommand{\modflowdate}{October 30, 2018} | ||
\newcommand{\currentmodflowversion}{Version \modflowversion---\modflowdate} |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Oops, something went wrong.