From aa4b6672730618c1deb02bd272bc5daf683875df Mon Sep 17 00:00:00 2001 From: ConstableCatnip <38230539+ConstableCatnip@users.noreply.github.com> Date: Thu, 12 Dec 2019 09:38:37 -0700 Subject: [PATCH] fix(upw): keyword for vertical K (#746) --- autotest/t003_test.py | 15 ++++++++++ examples/data/parameters/twrip_upw.nam | 14 +++++++++ examples/data/parameters/twrip_upw.nwt | 2 ++ examples/data/parameters/twrip_upw.pval | 6 ++++ examples/data/parameters/twrip_upw.upw | 35 ++++++++++++++++++++++ examples/data/parameters/twrip_upw.zon | 20 +++++++++++++ flopy/modflow/mfhob.py | 4 +-- flopy/modflow/mfupw.py | 40 ++++++++++++++++--------- 8 files changed, 120 insertions(+), 16 deletions(-) create mode 100644 examples/data/parameters/twrip_upw.nam create mode 100644 examples/data/parameters/twrip_upw.nwt create mode 100644 examples/data/parameters/twrip_upw.pval create mode 100644 examples/data/parameters/twrip_upw.upw create mode 100644 examples/data/parameters/twrip_upw.zon diff --git a/autotest/t003_test.py b/autotest/t003_test.py index e547ad2147..714fecb116 100644 --- a/autotest/t003_test.py +++ b/autotest/t003_test.py @@ -44,6 +44,20 @@ def test_loadtwrip(): assert ml.load_fail is False return + + +def test_loadtwrip_upw(): + cwd = os.getcwd() + pth = os.path.join('..', 'examples', 'data', 'parameters') + assert (os.path.isdir(pth)) + os.chdir(pth) + namefile = 'twrip_upw.nam' + ml = flopy.modflow.Modflow.load(namefile, verbose=True) + os.chdir(cwd) + assert isinstance(ml, flopy.modflow.Modflow) + assert ml.load_fail is False + + return def test_loadoc(): @@ -116,3 +130,4 @@ def test_load_nam_mt_nonexistant_file(): test_loadfreyberg() test_loadoahu() test_loadtwrip() + test_loadtwrip_upw() diff --git a/examples/data/parameters/twrip_upw.nam b/examples/data/parameters/twrip_upw.nam new file mode 100644 index 0000000000..235b0b0b14 --- /dev/null +++ b/examples/data/parameters/twrip_upw.nam @@ -0,0 +1,14 @@ +LIST 2 twrip_upw.lst +BAS6 3 twrip.ba6 +PVAL 4 twrip_upw.pval +UPW 11 twrip_upw.upw +WEL 12 twrip.wel +DRN 13 twrip.drn +RCH 18 twrip.rch +NWT 19 twrip_upw.nwt +OC 22 twrip.oc +MULT 8 twrip.mlt +ZONE 9 twrip_upw.zon +DIS 10 twrip.dis +DATA 101 twrip.ib1.ref +DATA 102 rchzones.ref diff --git a/examples/data/parameters/twrip_upw.nwt b/examples/data/parameters/twrip_upw.nwt new file mode 100644 index 0000000000..9a037ce764 --- /dev/null +++ b/examples/data/parameters/twrip_upw.nwt @@ -0,0 +1,2 @@ +# NWT package for MODFLOW-NWT, generated by Flopy. + 1.000e-02 5.000e+02 100 1.000e-05 1 0 0 SIMPLE diff --git a/examples/data/parameters/twrip_upw.pval b/examples/data/parameters/twrip_upw.pval new file mode 100644 index 0000000000..9fd4ea76c0 --- /dev/null +++ b/examples/data/parameters/twrip_upw.pval @@ -0,0 +1,6 @@ +# Example pval +4 +HK1_z1 1.0E-3 +HK1_z2 1.0E-2 +VK1_z1 5.0E-4 +VK1_z2 5.0E-3 diff --git a/examples/data/parameters/twrip_upw.upw b/examples/data/parameters/twrip_upw.upw new file mode 100644 index 0000000000..bddc04f100 --- /dev/null +++ b/examples/data/parameters/twrip_upw.upw @@ -0,0 +1,35 @@ +# UPW Package input data for example problem +0 1.00E+30 10 0 IUPWCB HDRY NPUPW IPHDRY + 1 0 0 + 0 0 0 + 1. 1. 1. + 0 0 0 + 0 0 0 +HK1_z1 HK 1.0E-3 1 +1 NONE upwzones 1 +HK1_z2 HK 1.0E-3 1 +1 NONE upwzones 2 +HK2 HK 1.0E-4 1 +2 NONE ALL +HK3 HK 2.0E-4 1 +3 NONE ALL +VK1_z1 VK 1.0E-4 1 +1 NONE upwzones 1 +VK1_z2 VK 1.0E-4 1 +1 NONE upwzones 2 +VK2 VK 1.0E-5 1 +2 NONE ALL +VK3 VK 1.0E-5 1 +3 NONE ALL +VKCB1 VKCB 1.0 1 +1 MULT1 ALL +VKCB2 VKCB 0.5 1 +2 MULT1 ALL +0 HK layer 1 +0 VK layer 1 +0 VKCB layer 1 +0 HK layer 2 +0 VK layer 2 +0 VKCB layer 2 +0 HK ayer 3 +0 VK layer 3 diff --git a/examples/data/parameters/twrip_upw.zon b/examples/data/parameters/twrip_upw.zon new file mode 100644 index 0000000000..83b3328c7c --- /dev/null +++ b/examples/data/parameters/twrip_upw.zon @@ -0,0 +1,20 @@ +2 +RCHZONES +EXTERNAL 102 1 (free) 1 +upwzones +INTERNAL 1 (15I1) -1 +111111111111111 +111111111111111 +111111111111111 +111111111111111 +111111111111111 +111111111111111 +111111111111111 +111111111111111 +111111111111111 +111111111111111 +111111112211111 +111111122221111 +111111111222111 +111111111122222 +111111111111222 diff --git a/flopy/modflow/mfhob.py b/flopy/modflow/mfhob.py index 7c1ee6058a..8d3034efae 100755 --- a/flopy/modflow/mfhob.py +++ b/flopy/modflow/mfhob.py @@ -577,9 +577,9 @@ def __init__(self, model, tomulth=1., obsname='HOBS', for key, value in self.mlay.items(): tot += value if not (np.isclose(tot, 1.0, rtol=0)): - msg = 'sum of dataset 4 proportions must equal 1.0 - ' + \ + msg = ('sum of dataset 4 proportions must equal 1.0 - ' + \ 'sum of dataset 4 proportions = {tot} for ' + \ - 'observation name {obsname}.'.format( + 'observation name {obsname}.').format( tot=tot, obsname=self.obsname) raise ValueError(msg) diff --git a/flopy/modflow/mfupw.py b/flopy/modflow/mfupw.py index c0893e9334..23f2446012 100644 --- a/flopy/modflow/mfupw.py +++ b/flopy/modflow/mfupw.py @@ -251,11 +251,12 @@ def write_file(self, check=True, f=None): # Item 0: text f_upw.write('{}\n'.format(self.heading)) # Item 1: IBCFCB, HDRY, NPLPF - f_upw.write('{0:10d}{1:10.3G}{2:10d}{3:10d}{4:s}\n'.format(self.ipakcb, - self.hdry, - self.npupw, - self.iphdry, - self.options)) + f_upw.write('{0:10d}{1:10.3G}{2:10d}{3:10d}{4:s}\n' + .format(self.ipakcb, + self.hdry, + self.npupw, + self.iphdry, + self.options)) # LAYTYP array f_upw.write(self.laytyp.string) # LAYAVG array @@ -276,7 +277,7 @@ def write_file(self, check=True, f=None): if self.chani[k] < 1: f_upw.write(self.hani[k].get_file_entry()) f_upw.write(self.vka[k].get_file_entry()) - if transient == True: + if transient: f_upw.write(self.ss[k].get_file_entry()) if self.laytyp[k] != 0: f_upw.write(self.sy[k].get_file_entry()) @@ -389,8 +390,7 @@ def load(f, model, ext_unit_dict=None, check=True): laywet = np.empty((nlay,), dtype=np.int32) laywet = read1d(f, laywet) - # Item 7: WETFCT, IWETIT, IHDWET - wetfct, iwetit, ihdwet = None, None, None + # check that LAYWET is 0 for all layers iwetdry = laywet.sum() if iwetdry > 0: raise Exception('LAYWET should be 0 for UPW') @@ -408,7 +408,10 @@ def load(f, model, ext_unit_dict=None, check=True): ss = [0] * nlay sy = [0] * nlay vkcb = [0] * nlay + # load by layer for k in range(nlay): + + # hk if model.verbose: print(' loading hk layer {0:3d}...'.format(k + 1)) if 'hk' not in par_types: @@ -419,6 +422,8 @@ def load(f, model, ext_unit_dict=None, check=True): t = mfpar.parameter_fill(model, (nrow, ncol), 'hk', parm_dict, findlayer=k) hk[k] = t + + # hani if chani[k] < 1: if model.verbose: print(' loading hani layer {0:3d}...'.format(k + 1)) @@ -430,23 +435,26 @@ def load(f, model, ext_unit_dict=None, check=True): t = mfpar.parameter_fill(model, (nrow, ncol), 'hani', parm_dict, findlayer=k) hani[k] = t + + # vka if model.verbose: print(' loading vka layer {0:3d}...'.format(k + 1)) + key = 'vk' + if layvka[k] != 0: + key = 'vani' if 'vk' not in par_types and 'vani' not in par_types: - key = 'vka' - if layvka[k] != 0: - key = 'vani' t = Util2d.load(f, model, (nrow, ncol), np.float32, key, ext_unit_dict) else: line = f.readline() - key = 'vka' - if 'vani' in par_types: - key = 'vani' t = mfpar.parameter_fill(model, (nrow, ncol), key, parm_dict, findlayer=k) vka[k] = t + + # storage properties if transient: + + # ss if model.verbose: print(' loading ss layer {0:3d}...'.format(k + 1)) if 'ss' not in par_types: @@ -457,6 +465,8 @@ def load(f, model, ext_unit_dict=None, check=True): t = mfpar.parameter_fill(model, (nrow, ncol), 'ss', parm_dict, findlayer=k) ss[k] = t + + # sy if laytyp[k] != 0: if model.verbose: print(' loading sy layer {0:3d}...'.format(k + 1)) @@ -469,6 +479,8 @@ def load(f, model, ext_unit_dict=None, check=True): t = mfpar.parameter_fill(model, (nrow, ncol), 'sy', parm_dict, findlayer=k) sy[k] = t + + # vkcb if model.get_package('DIS').laycbd[k] > 0: if model.verbose: print(' loading vkcb layer {0:3d}...'.format(k + 1))