Inspect MODFLOW 6 interactively

pymf6 allows to access all MODFLOW 6 variables at run time. First we import the class MF6:

from pymf6.mf6 import MF6

We use a simple model:

nam_file = 'examples/head_controlled_well/models/pymf6/mfsim.nam'

and instantiate the class:

%xmode Verbose
Exception reporting mode: Verbose
mf6 = MF6(nam_file=nam_file)

The instance has HTML representation with some meta data:

mf6

MF6

pymf6 configuration data

pymf6 version:1.0.4.dev3+g73beafa6
xmipy version:1.2.0
ini file path:HOME/pymf6.ini
dll file path:path/to/dll/libmf6.dll
MODFLOW version:6.4.1

MODFLOW variable documentation is available

The same information is available as text:

mf6.info
========================
pymf6 configuration data
========================
pymf6 version: 1.0.4.dev3+g73beafa6
xmipy version: 1.2.0
ini file path: HOME/pymf6.ini
dll file path: path/to/dll/libmf6.dll
MODFLOW version: 6.4.1
MODFLOW variable documentation is available

MODFLOW 6 supports multiple simulations. This example only has one:

mf6.simulation
modeltype namefile modelname
gwf6 headconwell.nam HEADCONWELL

The names of the simulations are available as list of strings:

mf6.simulation.model_names
['HEADCONWELL']

The models have their own representations:

mf6.simulation.models
[Model HEADCONWELL 
 15 packages
 49 variables.]

The meta information is also available as a dictionary:

mf6.simulation.models_meta
[{'modeltype': 'gwf6',
  'namefile': 'headconwell.nam',
  'modelname': 'HEADCONWELL'}]

This example only one solution group:

mf6.simulation.solution_groups
[Solution 1 
 1 packages
 68 variables.]

This is the time discretization is:

mf6.simulation.TDIS

Package TDIS

number of variables: 19

These are the names of the variables:

mf6.simulation.TDIS.var_names
['NPER',
 'ITMUNI',
 'KPER',
 'KSTP',
 'INATS',
 'READNEWDATA',
 'ENDOFPERIOD',
 'ENDOFSIMULATION',
 'DELT',
 'PERTIM',
 'TOTIM',
 'TOTIMC',
 'DELTSAV',
 'TOTIMSAV',
 'PERTIMSAV',
 'TOTALSIMTIME',
 'PERLEN',
 'NSTP',
 'TSMULT']

There a scalar values:

mf6.simulation.TDIS.NPER

Variable NPER

value: 4
docstring: set equal to nper
docstring source: src/Timing/ats.f90 line 25
docstring: number of stress period
docstring source: src/Timing/tdis.f90 line 22

and arrays:

mf6.simulation.TDIS.PERLEN

Variable PERLEN

value: array([ 1., 10., 10., 10.])
shape: (4,)
docstring: length of each stress period
docstring source: src/Timing/tdis.f90 line 38

The docstrings are extracted from the MODFLOW 6 source code from the used MODFLOW version as displayed with mf6.info.

The instance has many variables:

len(mf6.vars)
529

Filter for all TDIS entries:

{k: v for k, v in mf6.vars.items() if 'TDIS' in k}
{'TDIS/NPER': array([4], dtype=int32),
 'TDIS/ITMUNI': array([4], dtype=int32),
 'TDIS/KPER': array([0], dtype=int32),
 'TDIS/KSTP': array([0], dtype=int32),
 'TDIS/INATS': array([0], dtype=int32),
 'TDIS/DELT': array([0.]),
 'TDIS/PERTIM': array([0.]),
 'TDIS/TOTIM': array([0.]),
 'TDIS/TOTIMC': array([0.]),
 'TDIS/DELTSAV': array([0.]),
 'TDIS/TOTIMSAV': array([0.]),
 'TDIS/PERTIMSAV': array([0.]),
 'TDIS/TOTALSIMTIME': array([31.]),
 'TDIS/PERLEN': array([ 1., 10., 10., 10.]),
 'TDIS/NSTP': array([  1, 120, 120, 120], dtype=int32),
 'TDIS/TSMULT': array([1., 1., 1., 1.])}

Show the first 50:

dict(list(mf6.vars.items())[:50])
{'TDIS/NPER': array([4], dtype=int32),
 'TDIS/ITMUNI': array([4], dtype=int32),
 'TDIS/KPER': array([0], dtype=int32),
 'TDIS/KSTP': array([0], dtype=int32),
 'TDIS/INATS': array([0], dtype=int32),
 'TDIS/DELT': array([0.]),
 'TDIS/PERTIM': array([0.]),
 'TDIS/TOTIM': array([0.]),
 'TDIS/TOTIMC': array([0.]),
 'TDIS/DELTSAV': array([0.]),
 'TDIS/TOTIMSAV': array([0.]),
 'TDIS/PERTIMSAV': array([0.]),
 'TDIS/TOTALSIMTIME': array([31.]),
 'TDIS/PERLEN': array([ 1., 10., 10., 10.]),
 'TDIS/NSTP': array([  1, 120, 120, 120], dtype=int32),
 'TDIS/TSMULT': array([1., 1., 1., 1.]),
 'HEADCONWELL/ID': array([1], dtype=int32),
 'HEADCONWELL/IOUT': array([1005], dtype=int32),
 'HEADCONWELL/INEWTON': array([0], dtype=int32),
 'HEADCONWELL/IPRPAK': array([0], dtype=int32),
 'HEADCONWELL/IPRFLOW': array([0], dtype=int32),
 'HEADCONWELL/IPAKCB': array([-1], dtype=int32),
 'HEADCONWELL/IDSOLN': array([1], dtype=int32),
 'HEADCONWELL/NEQ': array([100], dtype=int32),
 'HEADCONWELL/NJA': array([460], dtype=int32),
 'HEADCONWELL/ICNVG': array([0], dtype=int32),
 'HEADCONWELL/MOFFSET': array([0], dtype=int32),
 'HEADCONWELL/INIC': array([1007], dtype=int32),
 'HEADCONWELL/INOC': array([1012], dtype=int32),
 'HEADCONWELL/INNPF': array([1008], dtype=int32),
 'HEADCONWELL/INBUY': array([0], dtype=int32),
 'HEADCONWELL/INVSC': array([0], dtype=int32),
 'HEADCONWELL/INSTO': array([1009], dtype=int32),
 'HEADCONWELL/INCSUB': array([0], dtype=int32),
 'HEADCONWELL/INMVR': array([0], dtype=int32),
 'HEADCONWELL/INHFB': array([0], dtype=int32),
 'HEADCONWELL/INGNC': array([0], dtype=int32),
 'HEADCONWELL/INOBS': array([0], dtype=int32),
 'HEADCONWELL/ISS': array([0], dtype=int32),
 'HEADCONWELL/INEWTONUR': array([0], dtype=int32),
 'HEADCONWELL/DIS/INUNIT': array([1006], dtype=int32),
 'HEADCONWELL/DIS/IOUT': array([1005], dtype=int32),
 'HEADCONWELL/DIS/NODES': array([100], dtype=int32),
 'HEADCONWELL/DIS/NODESUSER': array([100], dtype=int32),
 'HEADCONWELL/DIS/NDIM': array([3], dtype=int32),
 'HEADCONWELL/DIS/ICONDIR': array([1], dtype=int32),
 'HEADCONWELL/DIS/NOGRB': array([0], dtype=int32),
 'HEADCONWELL/DIS/XORIGIN': array([0.]),
 'HEADCONWELL/DIS/YORIGIN': array([0.]),
 'HEADCONWELL/DIS/ANGROT': array([0.])}

Many values are still zero. We can get the time steps:

steps = mf6.steps()

We can do single time steps:

next(steps)
0.0

Now, the values have changed:

dict(list(mf6.vars.items())[:5])
{'TDIS/NPER': array([4], dtype=int32),
 'TDIS/ITMUNI': array([4], dtype=int32),
 'TDIS/KPER': array([1], dtype=int32),
 'TDIS/KSTP': array([1], dtype=int32),
 'TDIS/INATS': array([0], dtype=int32)}

We can get access to well rates:

wel = mf6.vars['HEADCONWELL/WEL_0/BOUND']

This model has only one well:

wel
array([[0.]])

The first stress period is steady state hand has not well extraction rate. Doing the next steps will go the next stress period with a rate:

next(steps)
1.0
wel
array([[-0.05]])

Fast forward to next stress period, i.e. to day 11:

mf6.simulation.TDIS.PERLEN

Variable PERLEN

value: array([ 1., 10., 10., 10.])
shape: (4,)
docstring: length of each stress period
docstring source: src/Timing/tdis.f90 line 38
for step in steps:
    if step > 11:
        break

Now the value has changed:

wel
array([[-0.5]])