1

我有一个带输出的工作模型。我试图利用 flopy 来询问结果。具体来说,我希望能够(即时)为特定单元生成区域预算信息。我当然可以不厌其烦地构建我感兴趣的区域文件,但我希望能够根据需要提取结果,而无需经过中间步骤。

我正在寻找一个可以在 Python 中使用的工作流,以仅在 Python 控制台中为特定单元格从 CBB(已生成)中提取净流量信息(宁愿避免生成区域文件、导入然后提取的传统样式) .

编辑

试图让该实用程序完全正常工作时遇到了问题。这是最新的尝试。无法判断 flopy 是否在 CBB 本身上窒息,或者我提供的区域是否给它带来了问题。

该模型为 7 层,196 行,241 列。我正在尝试提取 Layer=3、Row=58、Col=30。我为区域创建的列表对象在这里:

zon_lst = []
for lay in range(7):
    for row in range(196):
        for col in range(241):
            if lay == 2 and row == 57 and col == 29:
                zon_lst.append(1)
            else:
                zon_lst.append(0)

然后我使用以下内容将其写入文件:

def chunk_list(alist, n):
    for i in range(0, len(alist), n):
        yield alist[i:i + n]

def zon_gen(mylst, rows, cols, file_cols):
    # Assumes integers for zonation file
    frmt = '{:>3}'
    list_by_lay = chunk_list(mylst, rows * cols)
    astr = '\n'.join(['\n'.join([''.join([frmt.format(i) for i in seq]) for seq in chunk_list(layer, file_cols)]) for layer in list_by_lay])
    return astr

zon_str = zon_gen(zon_lst, 196, 241, 10)
with open('t_26.zon', 'w+') as file:
    file.write('# For Scoping Work\n')
    file.write('1\n')
    file.write('t_26\n')
    file.write('INTERNAL 1 (10I3)  1\n')
    file.write(zon_str)

然后我为区域预算类/方法构建我的 modflow 模型:

import flopy
mf = flopy.modflow.Modflow(modelname='t_26_scope', version='mf2k')
zon = flopy.modflow.ModflowZon.load(r"t_26.zon",mf,nrow=196, ncol=241)
zb = flopy.utils.zonbud.ZoneBudget(r"P2Rv8.2_1000yr.cbb", zon)

所有这些都运行良好,直到最后一个命令,我得到以下错误:

Traceback (most recent call last):
  File "<input>", line 1, in <module>
  File "C:\ProgramData\Miniconda3\envs\ca_cie\lib\site-packages\flopy\utils\zonbud.py", line 53, in __init__
    self.cbc = CellBudgetFile(cbc_file)
  File "C:\ProgramData\Miniconda3\envs\ca_cie\lib\site-packages\flopy\utils\binaryfile.py", line 618, in __init__
    self._build_index()
  File "C:\ProgramData\Miniconda3\envs\ca_cie\lib\site-packages\flopy\utils\binaryfile.py", line 708, in _build_index
    self._skip_record(header)
  File "C:\ProgramData\Miniconda3\envs\ca_cie\lib\site-packages\flopy\utils\binaryfile.py", line 768, in _skip_record
    self.file.seek(nbytes, 1)
OSError: [Errno 22] Invalid argument

请记住,我仍然对跳过文件写入部分感兴趣。我把它包括在内是为了展示我到目前为止的工作

4

2 回答 2

0

您可以使用ZoneBudget. flopy它接受一个 numpy 数组作为输入;read_zbarray和伴随的write_zbarray实用程序被包括在内是为了向后兼容。如果您遇到问题,请在此处发布最低工作示例,我会看看是否可以帮助您。

于 2019-10-29T13:35:53.397 回答
0

ZoneBudget 类所需的矩阵必须具有以下维度(通用术语):

ndarray_dim = lays * rows * cols

在我放在这里的示例问题中,对应的尺寸是:7 * 196 * 241,结果是:

len(my_array)
>>> 7
len(my_array[0])
>>> 196
len(my_array[0][0])
>>> 241

正如所指出的,我需要找出逐个单元预算二进制文件的二进制读取器无法读取文件的原因,因此我将准备作为后续讨论的另一个讨论。感谢您幽默@Jason Bellino 的问题!

附带说明一下,区域文件生成器不会像 flopy 所期望的那样生成块样式文件。我在这里添加了我的更新代码:

def chunk_list(alist, n):
for i in range(0, len(alist), n):
    yield alist[i:i + n]


def block_gen(mylst, rows, cols, file_cols, field_len):
    # mylst is a 1D array whose length matches that of the number of cells in the model
    # rows, cols: These are the total rows, columns; respectively
    # Assumes integers for zonation file
    frmt = '{:>%d}' % field_len
    zon_str = ''
    for lay in chunk_list(mylst, (rows * cols)):
        zon_str += 'INTERNAL          ({:}I{})\n'.format(file_cols, field_len)
        for block in chunk_list(lay, cols):
            for line in chunk_list(block, file_cols):
                zon_str += ''.join([frmt.format(cell) for cell in line]) + '\n'
    return zon_str


def write_zb_zonefile(filepath, lays, rows, cols, zon_str):
    with open(filepath, 'w+') as file:
        file.write('{:<6}{:<6}{:<6}\n'.format(lays, rows, cols))
        file.write(zon_str)
    return

其中,如果您使用生成的文件,flopy.utils.zonbud.read_zbarray(<filepath>)您将获得一个可用的 zone_dict 以应用于flopy.utils.zonbud.ZoneBudget课程。

因为我想避免每次需要新的区域预算时都必须写出一个文件,所以我修改了block_gen()函数并创建了一个将生成可在ZoneBudget类中直接使用的内容的函数:

from numpy import array
def arr_block_gen(mylst, rows, cols):
    # mylst: 1D array containing zone definition for every cell in the model
    # rows: total number of rows in the model
    # cols: total number of columns in the model
    zon_arr = []
    lay_idx = 0
    for lay in chunk_list(mylst, (rows * cols)):
        row_idx = 0
        zon_arr.append([])
        for block in chunk_list(lay, cols):
            zon_arr[lay_idx].append(block)
            row_idx += 1
        lay_idx += 1
    zon_arr = array(zon_arr)
    return zon_arr

UPDATE 原来该方法被炸毁的原因是由于该方法的默认行为是将 CBB 文件路径传递给 binaryfile 读取器类,其默认值是以单精度读取 CBB 文件。因此,我通过使用二进制文件读取器创建一个 CBB 对象(传递“双精度”作为精度)然后显式传递 CBB 对象而不是像我在问题中显示的文件来解决这个问题。

于 2019-10-29T20:14:42.930 回答