2

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

Community
  • 1
  • 1
Jacob
  • 58
  • 9
  • This appears to be an issue with the binary file reader; try to loading your cbb file using the CellBudgetFile class directly and see if you get the same error. As a side note, the ZoneBudget class was not intended for use with the ModflowZon package. It appears the object returned from that module is an instance of a ModflowZon package as opposed to a numpy array. – Jason Bellino Oct 29 '19 at 17:12
  • I did find that the error is related to the binaryfile reader class. It cannot read the CBB that my model has generated. However, using the zonebudget utility by the USGS, I have been able to read the same CBB (just tested it). What is the format of the numpy array expected by the zonebudget reader? That's all I'm really after at this point. I'll have to find the bug in the reader and submit that as another question/error report on GitHub. – Jacob Oct 29 '19 at 18:17
  • 1
    The flopy ZoneBudget class will accept a 2-D array of dimensions nrow * ncol, or a 3-D array of dimensions nlay * nrow * ncol; where nlay, nrow, ncol are the total number of layers, rows, and columns in the model, respectively. If a 2-D array is used, ZoneBudget will use the same zone distribution for each layer in the binary output. – Jason Bellino Oct 29 '19 at 21:09
  • On another note, the error that I was running into is that my CBB file is in double precision. If you use the flopy.utils.zonbud.ZoneBudget() class as it is, by default the class will specific "single precision" for the CBB file read. After fixing this and applying my answer I was able to generate zonebudget results on the fly. – Jacob Oct 30 '19 at 21:17
  • Do you think it would it be helpful for users to have the option to specify precision when passing a filename string or is the current ability to pass a CellBudgetFile object directly enough? One would have to know whether to pass a precision keyword arg either way, so I was leaning toward leaving it as is. – Jason Bellino Oct 31 '19 at 13:55
  • Either way for me at this point, though from a desire to see flopy continue to mature, an error message would be as helpful as having an explicit option in my opinion. Either I'd see it as an option to read my CBB file with double/single precision (default as single), or I'd have an error message ask if I'm using a double-precision file. My personal preference would be to have what you say: an argument in the method that allows me to explicitly define the precision of the CBB file. – Jacob Nov 11 '19 at 23:52

2 Answers2

0

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.

Jason Bellino
  • 494
  • 1
  • 6
  • 18
0

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.

Jacob
  • 58
  • 9