What is the workflow for creating a zone for zone_budget work on the fly?

Question

I have a working model with output. I'm trying to take advantage of flopy in interrogating the results. In specific, I would like to be able to generate (on the fly) zone budget information for specific cells. I can of course go to the trouble of building the zone files which I am interested in, but I was hoping to be able to extract the results as I need them without having to go through that intermediate step.

I'm looking for a workflow that I can use in Python to extract net flow information from the CBB (already generated) for specific cells only in the Python console (would rather avoid the legacy style of generating zone files, importing, then extracting).

EDIT

Ran into problems trying to get the utility to work at all. Here's the latest attempt. Can't tell if flopy is choking on the CBB itself or if the zone I've provided is giving it problems.

The model is 7 layers, 196 Rows, 241 Columns. I am trying to extract Layer=3, Row=58, Col=30. The list object I created for the zone is here:

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)

I then wrote this to a file using the following:

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)

Then I build my modflow model for the zone budget class/methods:

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)

All of this runs fine up until the very last command, where I get the following error:

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

Please bare in mind that I am still interested in skipping the file writing portion. I included it to show my work up to this point

Answer 1

You can get extract budget information using the ZoneBudget class of flopy. It accepts a numpy array as input; the read_zbarray and write_zbarray companion utilities were included for backwards compatibility. If you run into a problem post a minimum working example here and I'll see if I can help you out.

Answer 2

The matrix required by the ZoneBudget class must have the following dimensions (in generic terms):

ndarray_dim = lays * rows * cols

In the example problem I put here, the corresponding dimensions are: 7 * 196 * 241, resulting in:

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

As was pointed out, I will need to track down the reason why the binary reader for the cell-by-cell budget binary can't read the file, so I'll prepare that as another discussion down the road. Thanks for humoring the question @Jason Bellino!

On a side note, the zone file generator doesn't generate a block-style file like what flopy is expecting. I've added my updated code here:

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

Which, if you use the resulting file with the flopy.utils.zonbud.read_zbarray(<filepath>) you will get a usable zone_dict to apply in the flopy.utils.zonbud.ZoneBudget class.

Because I wanted to avoid having to write a file out every time I wanted a new zone budget, I modified the block_gen() function and created one that will generate something that is directly usable in the ZoneBudget class:

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 Turns out the reason the method bombed out was due to the fact that the default behavior of this method is to pass the CBB file path to the binaryfile reader class, whose default is to read CBB files as single-precision. So, I worked around this by creating a CBB object using the binaryfile reader (passing "double" as the precision) and then explicitly passed the CBB object rather than the file like I show in my question.