mf6rtm.simulation package

Submodules

mf6rtm.simulation.discretization module

This module contains the Discretization functions for various ModFlow grid types.

mf6rtm.simulation.discretization.grid_dimensions(modflow_api)

Return the grid dimensions for the MODFLOW 6 model.

Parameters:

modflow_api (Mf6API) – MODFLOW 6 API instance to query the grid from.

Returns:

(nlay, nrow, ncol) for a DIS grid or (nlay, ncpl) for a DISV grid.

Return type:

tuple[int, ...]

Raises:

ValueError – If the grid type is not supported (only DIS and DISV are).

mf6rtm.simulation.discretization.total_cells_in_grid(modflow_api)

Return the total number of cells in the MODFLOW 6 grid.

Parameters:

modflow_api (Mf6API) – MODFLOW 6 API instance to query the grid from.

Returns:

Product of the grid dimensions (nlay * nrow * ncol for DIS, nlay * ncpl for DISV).

Return type:

int

mf6rtm.simulation.mf6api module

MODFLOW 6 API wrapper for mf6rtm.

Defines Mf6API, a subclass of modflowapi.ModflowApi that loads the coupled flow/transport simulation and drives its solve loop inside the mf6rtm coupling.

class mf6rtm.simulation.mf6api.Mf6API(wd, dll)

Bases: ModflowApi

MODFLOW 6 API wrapper driving the flow and transport solve.

Extends modflowapi.ModflowApi to load the coupled simulation and run the per-timestep GWF/GWT solve loop that Mf6RTM alternates with the PhreeqcRM reaction step.

Parameters:

  • wd (os.PathLike) – Working directory containing the MODFLOW 6 simulation.

  • dll (os.PathLike) – Path to the MODFLOW 6 shared library (libmf6).

property grid_type: str

Return the grid type of the MODFLOW 6 model.

mf6rtm.simulation.phreeqcbmi module

PHREEQC BMI wrapper for mf6rtm.

Defines PhreeqcBMI, a thin subclass of phreeqcrm.BMIPhreeqcRM that adapts the PhreeqcRM Basic Model Interface for use inside the mf6rtm coupling loop.

class mf6rtm.simulation.phreeqcbmi.PhreeqcBMI(yaml='mf6rtm.yaml')

Bases: BMIPhreeqcRM

Basic Model Interface wrapper around PhreeqcRM.

Extends phreeqcrm.BMIPhreeqcRM with helpers used by Mf6RTM to run the geochemical reaction step: scalar setters, selected-output handling and per-timestep solution updates.

Parameters:

yaml (str, optional) – Path to the PhreeqcRM YAML configuration file used to initialize the chemistry. Default is "mf6rtm.yaml".

get_grid_to_map()

Get the PhreeqcRM grid-to-map mapping.

Returns:

The grid-to-map mapping returned by GetGridToMap.

Return type:

array-like

set_scalar(var_name, value)

Set a scalar PhreeqcRM BMI variable.

Parameters:

  • var_name (str) – Name of the BMI variable to set.

  • value (Any) – Scalar value to assign; cast to the variable’s declared dtype.

Raises:

ValueError – If var_name does not refer to a scalar (dimension != 1).

mf6rtm.simulation.solver module

The solver module provides the Mf6RTM class that couples modflowapi and phreeqcrm, along with functions to run the coupled simulations.

class mf6rtm.simulation.solver.Mf6RTM(wd, mf6api, phreeqcbmi, output_format=None)

Bases: object

Reactive transport loop coupling MODFLOW 6 and PhreeqcRM.

Drives the sequential non-iterative (SNIA) operator-splitting loop: each time step MODFLOW 6 solves flow and conservative transport, the resulting concentrations are handed to PhreeqcRM to solve chemistry, and the reacted concentrations are written back to the transport models.

Instances are normally created via prep_to_run() and run with solve() rather than constructed directly. See __init__ for the full list of parameters and attributes.

get_saturation_from_mf6()

Get the cell saturation array from MODFLOW 6.

Reads the saturation for the transport models and stores it on the PhreeqcRM interface (sat_now).

Returns:

Saturation values for the grid; identical across all components.

Return type:

dict[Any, ndarray]

get_time_units_from_mf6()

Get the simulation time units from MODFLOW 6.

Returns:

The MODFLOW 6 TDIS time units (e.g. "days").

Return type:

str

is_reactive_tstep()

Check whether the current timestep should be reactive.

Returns:

True if the current timestep should be reactive, based on the run configuration; False otherwise.

Return type:

bool

print_warning_user_active()

Print a warning if reaction timing is set to ‘user’.

set_emulator_output_add_variables()

Add emulator target and feature variables to the output.

Updates selected_output with variables defined in the configuration. Defaults to empty lists if not provided.

selected_output.target_var

Target variables for emulator training.

Type:

list of str

selected_output.feat_var

Feature variables for emulator training.

Type:

list of str

Return type:

None

set_emulator_training()

Configure emulator training output.

Reads emulator_training_data from the configuration. If enabled, sets up emulator output variables; otherwise disables training data.

ml_output

Whether emulator training data output is enabled.

Type:

bool

Return type:

None

set_kiter()

Increment and return the coupling iteration counter.

Initializes kiter to 0 on the first call, then increments it on each subsequent call.

Returns:

The current iteration counter.

Return type:

int

set_min_concentration(min_concentration)

Set the floor value for non-charge concentrations passed to PhreeqcRM.

Parameters:

min_concentration (float | None) – Any non-charge concentration below this value will be replaced with this value before calling SetConcentrations. Use a small positive number (e.g. 1e-30) to avoid PHREEQC divide-by-zero errors. Pass None to disable filtering entirely.

Return type:

None

set_time_conversion()

Set the time conversion factor from the MODFLOW 6 time units.

Return type:

None

solve()

Run the coupled MODFLOW 6 / PhreeqcRM time-stepping loop.

Returns:

True if the run finished and both interfaces were finalized successfully.

Return type:

bool

mf6rtm.simulation.solver.check_config_file(wd)

Load and validate the mf6rtm TOML configuration for a run.

Parameters:

wd (PathLike) – Model working directory; must contain mf6rtm.toml.

Returns:

The loaded, validated configuration.

Return type:

tuple[PathLike, PathLike]

Raises:

AssertionError – If mf6rtm.toml is not found in wd.

mf6rtm.simulation.solver.check_nam_files(wd)

Check that the required MODFLOW 6 name file is present.

Parameters:

wd (PathLike) – Model working directory to search for .nam files.

Returns:

Path to mfsim.nam in the working directory.

Return type:

tuple[PathLike, PathLike]

Raises:

AssertionError – If mfsim.nam is not found in wd.

mf6rtm.simulation.solver.get_conc_change_mask(ck, ci, threshold=1e-15)

Get the active/inactive cell mask for transport-driven concentration change.

Informs PhreeqcRM which cells to update: a cell is active (1) when any component’s relative concentration change exceeds threshold.

Parameters:

  • ck (ndarray[tuple[int, int], float64]) – Current-iteration concentrations, 2D array (ncomps, nxyz).

  • ci (ndarray[tuple[int, int], float64]) – Previous-iteration concentrations, 2D array (ncomps, nxyz).

  • threshold (float) – Relative-change threshold below which a cell is considered unchanged. Default is 1e-15 (float64 relative precision).

Returns:

1D mask of size nxyz: 0 for inactive cells, 1 for cells to update.

Return type:

ndarray[float64]

mf6rtm.simulation.solver.get_inactive_idx(arr, val=1e+30)

Get the indices of values greater than or equal to a threshold.

Parameters:

  • arr (ndarray) – Input array to test.

  • val (float) – Threshold marking inactive/dry cells. Default is 1e30.

Returns:

Indices of elements greater than or equal to val.

Return type:

list of int

mf6rtm.simulation.solver.get_less_than_zero_idx(arr)

Get the indices of all negative values in an array.

Parameters:

arr (array-like) – Input array to test.

Returns:

The result of numpy.where(arr < 0): index arrays for elements less than zero.

Return type:

tuple of numpy.ndarray

mf6rtm.simulation.solver.initialize_interfaces(wd, nthread=1, libname=None, output_format=None)

Initialize the MODFLOW 6 and PhreeqcRM interfaces.

Parameters:

  • wd (PathLike) – Model working directory.

  • nthread (int) – Number of threads for PhreeqcRM; values greater than 1 are written into the YAML file. Default is 1.

  • libname (Path | None) – Path to the MODFLOW 6 shared library. Default is None.

  • output_format (str) – Selected-output format passed to the Mf6RTM instance. Default is None.

Returns:

The coupled reactive transport object wrapping both interfaces.

Return type:

Mf6API

mf6rtm.simulation.solver.longest_common_substring(strings)

Find the longest common substring shared by a list of strings.

Used to find the common “stem” of the GWT model names for matching with PhreeqcRM components.

Parameters:

strings (list of str) – Strings to compare.

Returns:

The longest substring common to all strings, or an empty string if strings is empty or no common substring exists.

Return type:

str

mf6rtm.simulation.solver.mrbeaker()

Render the Mr. Beaker mascot image as ASCII art.

Returns:

Multi-line ASCII-art representation of the Mr. Beaker image.

Return type:

str

mf6rtm.simulation.solver.prep_to_run(wd, libname=None)

Prepares the model to run by checking if the model directory (wd) contains the necessary files and returns the path to the yaml file (phreeqcrm) and the dll file (mf6 api)

Parameters:

wd (PathLike) – The path to the working directory of model directory

Returns:

The path to the phreeqcrm model file (yaml) and the path to the MODFLOW 6 dll (associated with mf6api).

Return type:

tuple[PathLike, PathLike]

mf6rtm.simulation.solver.run_cmd(cwd=None)

Console entrypoint compatibility wrapper.

When used as a console script the entrypoint calls mf6rtm:run_cmd with no arguments. Allow cwd to be optional and default to the current working directory.

Return type:

None

mf6rtm.simulation.solver.set_nthread_yaml(yamlfile, nthread=1)

Set the number of threads in the PhreeqcRM YAML file.

Parameters:

  • yamlfile (PathLike) – Path to the PhreeqcRM YAML file to modify in place.

  • nthread (int) – Number of threads to write to the nthreads entry. Default is 1.

Return type:

None

mf6rtm.simulation.solver.solve(wd, reactive=None, nthread=1, libname=None, output_format=None, **mf6rtm_kwargs)

Prepare the coupled interfaces and run the reactive transport solve.

Parameters:

  • wd (PathLike) – Model working directory.

  • reactive (bool | None) – If given, override the reactive mode taken from the configuration. When it differs from the configured value the mode is switched (and selected output is left to MODFLOW 6 for conservative runs). Default is None.

  • nthread (int) – Number of threads for PhreeqcRM. Default is 1.

  • libname (Path | None) – Path to the MODFLOW 6 shared library. Default is None.

  • output_format (str) – Selected-output format (e.g. "csv" or "hdf5"). Default is None.

  • **mf6rtm_kwargs – Additional attributes to set on the Mf6RTM instance before solving; each key must be an existing attribute.

Returns:

True if the simulation finished successfully.

Return type:

bool

Raises:

AttributeError – If a keyword in mf6rtm_kwargs is not an attribute of Mf6RTM.

Module contents

Main module to manage the APIs PHREEQCRM and MODFLOW6 API, and solve the reactive transport loop

class mf6rtm.simulation.Mf6API(wd, dll)

Bases: ModflowApi

MODFLOW 6 API wrapper driving the flow and transport solve.

Extends modflowapi.ModflowApi to load the coupled simulation and run the per-timestep GWF/GWT solve loop that Mf6RTM alternates with the PhreeqcRM reaction step.

Parameters:

  • wd (os.PathLike) – Working directory containing the MODFLOW 6 simulation.

  • dll (os.PathLike) – Path to the MODFLOW 6 shared library (libmf6).

property grid_type: str

Return the grid type of the MODFLOW 6 model.

class mf6rtm.simulation.Mf6RTM(wd, mf6api, phreeqcbmi, output_format=None)

Bases: object

Reactive transport loop coupling MODFLOW 6 and PhreeqcRM.

Drives the sequential non-iterative (SNIA) operator-splitting loop: each time step MODFLOW 6 solves flow and conservative transport, the resulting concentrations are handed to PhreeqcRM to solve chemistry, and the reacted concentrations are written back to the transport models.

Instances are normally created via prep_to_run() and run with solve() rather than constructed directly. See __init__ for the full list of parameters and attributes.

get_saturation_from_mf6()

Get the cell saturation array from MODFLOW 6.

Reads the saturation for the transport models and stores it on the PhreeqcRM interface (sat_now).

Returns:

Saturation values for the grid; identical across all components.

Return type:

dict[Any, ndarray]

get_time_units_from_mf6()

Get the simulation time units from MODFLOW 6.

Returns:

The MODFLOW 6 TDIS time units (e.g. "days").

Return type:

str

is_reactive_tstep()

Check whether the current timestep should be reactive.

Returns:

True if the current timestep should be reactive, based on the run configuration; False otherwise.

Return type:

bool

print_warning_user_active()

Print a warning if reaction timing is set to ‘user’.

set_emulator_output_add_variables()

Add emulator target and feature variables to the output.

Updates selected_output with variables defined in the configuration. Defaults to empty lists if not provided.

selected_output.target_var

Target variables for emulator training.

Type:

list of str

selected_output.feat_var

Feature variables for emulator training.

Type:

list of str

Return type:

None

set_emulator_training()

Configure emulator training output.

Reads emulator_training_data from the configuration. If enabled, sets up emulator output variables; otherwise disables training data.

ml_output

Whether emulator training data output is enabled.

Type:

bool

Return type:

None

set_kiter()

Increment and return the coupling iteration counter.

Initializes kiter to 0 on the first call, then increments it on each subsequent call.

Returns:

The current iteration counter.

Return type:

int

set_min_concentration(min_concentration)

Set the floor value for non-charge concentrations passed to PhreeqcRM.

Parameters:

min_concentration (float | None) – Any non-charge concentration below this value will be replaced with this value before calling SetConcentrations. Use a small positive number (e.g. 1e-30) to avoid PHREEQC divide-by-zero errors. Pass None to disable filtering entirely.

Return type:

None

set_time_conversion()

Set the time conversion factor from the MODFLOW 6 time units.

Return type:

None

solve()

Run the coupled MODFLOW 6 / PhreeqcRM time-stepping loop.

Returns:

True if the run finished and both interfaces were finalized successfully.

Return type:

bool

class mf6rtm.simulation.PhreeqcBMI(yaml='mf6rtm.yaml')

Bases: BMIPhreeqcRM

Basic Model Interface wrapper around PhreeqcRM.

Extends phreeqcrm.BMIPhreeqcRM with helpers used by Mf6RTM to run the geochemical reaction step: scalar setters, selected-output handling and per-timestep solution updates.

Parameters:

yaml (str, optional) – Path to the PhreeqcRM YAML configuration file used to initialize the chemistry. Default is "mf6rtm.yaml".

get_grid_to_map()

Get the PhreeqcRM grid-to-map mapping.

Returns:

The grid-to-map mapping returned by GetGridToMap.

Return type:

array-like

set_scalar(var_name, value)

Set a scalar PhreeqcRM BMI variable.

Parameters:

  • var_name (str) – Name of the BMI variable to set.

  • value (Any) – Scalar value to assign; cast to the variable’s declared dtype.

Raises:

ValueError – If var_name does not refer to a scalar (dimension != 1).

mf6rtm.simulation.run_cmd(cwd=None)

Console entrypoint compatibility wrapper.

When used as a console script the entrypoint calls mf6rtm:run_cmd with no arguments. Allow cwd to be optional and default to the current working directory.

Return type:

None

mf6rtm.simulation.solve(wd, reactive=None, nthread=1, libname=None, output_format=None, **mf6rtm_kwargs)

Prepare the coupled interfaces and run the reactive transport solve.

Parameters:

  • wd (PathLike) – Model working directory.

  • reactive (bool | None) – If given, override the reactive mode taken from the configuration. When it differs from the configured value the mode is switched (and selected output is left to MODFLOW 6 for conservative runs). Default is None.

  • nthread (int) – Number of threads for PhreeqcRM. Default is 1.

  • libname (Path | None) – Path to the MODFLOW 6 shared library. Default is None.

  • output_format (str) – Selected-output format (e.g. "csv" or "hdf5"). Default is None.

  • **mf6rtm_kwargs – Additional attributes to set on the Mf6RTM instance before solving; each key must be an existing attribute.

Returns:

True if the simulation finished successfully.

Return type:

bool

Raises:

AttributeError – If a keyword in mf6rtm_kwargs is not an attribute of Mf6RTM.