I can't run my model in any way

[Alvin-prc]

Could you please help me see what is wrong with my model that I can't run anyway
import os
import flopy
import numpy as np
import flopy.utils.binaryfile as bf
import matplotlib.pyplot as plt

############################################################################

1. MODFLOW

############################################################################
modelname = 'mod_flow'
workspace = './model_files'
mf = flopy.modflow.Modflow(modelname, exe_name="mf2005dbl.exe", model_ws=workspace)
############################################################################

2. DIS

############################################################################
Lx = 800.0 # x 长度
Ly = 800.0 # y 长度
ztop = 30.0 # z 顶部高程
zbot = 0.0 # z 底部高程
nlay = 1 # 含水层层数
nrow = 80 # 行数
ncol = 80 # 列数
delc = Lx / ncol # x 步长
delr = Ly / nrow # y 步长
botm = np.linspace(ztop, zbot, nlay + 1) # z 每一层的底部高程
nper = 10 # 应力期数
perlen = np.array(np.ones(nper) * 180) # 应力期长度数组
nstp = 2 # 每个应力期的时间步数
steady = np.ones(nper, dtype=bool)

Unit System

itmuni = 4 # time units (1=seconds)
lenuni = 2 # lenght units (3=centimeters)
dis = flopy.modflow.ModflowDis(mf, nlay, nrow, ncol, delr=delr,
delc=delc, top=ztop, botm=botm[1:],
nper=nper, perlen=perlen, nstp=nstp, steady=steady,
itmuni=itmuni, lenuni=lenuni)
###################################################################################

3. BAS

###################################################################################
ibound = np.ones((nlay, nrow, ncol), dtype=np.int32)
strt = 100.0 * np.ones((nlay, nrow, ncol), dtype=np.float32)
bas = flopy.modflow.ModflowBas(mf, ibound=ibound, strt=strt)
###################################################################################

4. LPF

###################################################################################
laytyp = 0 # 含水层类型
vka = 1.0 # 垂直渗透系数，单位：m/d
hk = np.ones((nlay, nrow, ncol), dtype=np.float32) # 水平渗透系数，单位：m/d

hk[0, :, :] = np.loadtxt('real_K.txt', dtype=np.float32)

sy = 0.1 # 给水度，单位：1/m

ss = 1.e-4 # 单位储水量，单位：1/m

lpf = flopy.modflow.ModflowLpf(mf, laytyp=0, chani=-1, layvka=1, hk=hk, vka=1, hdry=-888.0)
#########################################################

5. CHD

#########################################################
left_bound = []
right_bound = []
left_head = 120
right_head = 100
for lay in range(nlay):
for row in range(nrow):
left_bound.append([lay, row, 0, left_head, left_head])
right_bound.append([lay, row, -1, right_head, right_head])

bound_sp = left_bound + right_bound
chd_spd = {0: bound_sp}
chd = flopy.modflow.ModflowChd(mf, stress_period_data=chd_spd)

#########################################################

6. Link to MT3DMS

#########################################################
lmt = flopy.modflow.ModflowLmt(mf, output_file_name='mt3d_link.ftl')
#########################################################

8. OC

#########################################################
stress_period_data = {}
for kper in range(nper):
for kstp in range(nstp):
stress_period_data[(kper, kstp)] = ["save head"]
oc = flopy.modflow.ModflowOc(mf, stress_period_data=stress_period_data, compact=True)
#########################################################

9. PCG

#########################################################
pcg = flopy.modflow.ModflowPcg(mf, mxiter=20, iter1=50, hclose=1e-03, rclose=1e-03, relax=1.0)

mf.write_input()

try:
os.remove(os.path.join(workspace, modelname + '.hds'))
except FileNotFoundError:
pass
success, mfoutput = mf.run_model(pause=False, report=True)
if not success:
raise Exception('MODFLOW did not terminate normally.')

[spaulins-usgs]

There is a problem with the CHD package cells you are storing in the right_bound variable.

right_bound.append([lay, row, -1, right_head, right_head])

You have column set to -1, which is not valid. If you intend to have these CHD cells on the right-hand edge of your model change that line to:

right_bound.append([lay, row, ncol - 1, right_head, right_head])

[Alvin-prc]

Thank you, you made my day!