Gael-de-Sailly
/
flopy
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

fix(upw): keyword for vertical K (#746)

develop
ConstableCatnip 2019-12-12 09:38:37 -07:00 committed by langevin-usgs
parent f83f382deb
commit aa4b667273
8 changed files with 120 additions and 16 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

15
autotest/t003_test.py
View File

@ -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 @@ if __name__ == '__main__':
test_loadfreyberg()
test_loadoahu()
test_loadtwrip()
test_loadtwrip_upw()

14
examples/data/parameters/twrip_upw.nam Normal file
View File

@ -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

2
examples/data/parameters/twrip_upw.nwt Normal file
View File

@ -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

6
examples/data/parameters/twrip_upw.pval Normal file
View File

@ -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

35
examples/data/parameters/twrip_upw.upw Normal file
View File

@ -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

20
examples/data/parameters/twrip_upw.zon Normal file
View File

@ -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

4
flopy/modflow/mfhob.py
View File

@ -577,9 +577,9 @@ class HeadObservation(object):
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)

40
flopy/modflow/mfupw.py
View File

@ -251,11 +251,12 @@ class ModflowUpw(Package):
# 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 @@ class ModflowUpw(Package):
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 @@ class ModflowUpw(Package):
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 @@ class ModflowUpw(Package):
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 @@ class ModflowUpw(Package):
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 @@ class ModflowUpw(Package):
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 @@ class ModflowUpw(Package):
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 @@ class ModflowUpw(Package):
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))