FloPy is a Python package to create, run, and post-process MODFLOW-based models
66
This module provides comprehensive support for MT3DMS and MT3D-USGS transport models with advection, dispersion, reactions, and specialized packages for complex transport scenarios. MT3DMS is the multi-species transport model for simulating contaminant fate and transport in groundwater systems.
Main MT3DMS model class that serves as a container for transport packages.
class Mt3dms:
"""MT3DMS transport model container and execution manager"""
def __init__(
self,
modelname: str = 'mt3dtest',
namefile_ext: str = 'nam',
modflowmodel: object = None,
ftlfilename: str = 'mt3d_link.ftl',
ftlfree: bool = False,
version: str = 'mt3dms',
exe_name: str = 'mt3dms5',
structured: bool = True,
listunit: int = 16,
ftlunit: int = 10,
model_ws: str = '.',
external_path: str = None,
verbose: bool = False,
load: bool = True,
silent: int = 0,
**kwargs
): ...
def write_input(self, SelPackList: list = None, check: bool = True) -> None:
"""Write MT3DMS input files"""
...
def run_model(
self,
silent: bool = False,
pause: bool = False,
report: bool = False,
normal_msg: str = 'program completed'
) -> tuple[bool, list[str]]:
"""Run the MT3DMS model and return success status and output"""
...
@classmethod
def load(
cls,
f: str,
version: str = 'mt3dms',
exe_name: str = 'mt3dms5',
verbose: bool = False,
model_ws: str = '.',
load_only: list = None,
forgive: bool = False,
check: bool = True
) -> 'Mt3dms':
"""Load existing MT3DMS model from files"""
...
def add_package(self, p: object) -> None:
"""Add package to model"""
...
def remove_package(self, pname: str) -> None:
"""Remove package from model"""
...
def get_package(self, name: str) -> object:
"""Get package by name"""
...
@property
def modelgrid(self) -> object:
"""Model grid object"""
...
@property
def packages(self) -> list[object]:
"""List of model packages"""
...
@property
def ncomp(self) -> int:
"""Number of chemical components"""
...
@property
def nlay(self) -> int:
"""Number of layers"""
...
@property
def nrow(self) -> int:
"""Number of rows"""
...
@property
def ncol(self) -> int:
"""Number of columns"""
...Basic transport package that defines fundamental transport parameters and component setup.
class Mt3dBtn:
"""MT3DMS basic transport package"""
def __init__(
self,
model: Mt3dms,
nlay: int = 1,
nrow: int = 2,
ncol: int = 2,
nper: int = 1,
ncomp: int = 1,
mcomp: int = 1,
tunit: str = 'D',
lunit: str = 'M',
munit: str = 'KG',
laycon: ArrayData = 0,
delr: ArrayData = 1.0,
delc: ArrayData = 1.0,
htop: ArrayData = 1.0,
dz: ArrayData = 1.0,
prsity: ArrayData = 0.30,
icbund: ArrayData = 1,
sconc: ArrayData = 0.0,
cinact: float = -1e30,
thkmin: float = 0.01,
ifmtcn: int = 0,
ifmtnp: int = 0,
ifmtrf: int = 0,
ifmtdp: int = 0,
savucn: bool = True,
nprs: int = 0,
timprs: ArrayData = None,
obs: ArrayData = None,
nprobs: int = 1,
chkmas: bool = True,
nprmas: int = 1,
dt0: float = 0,
mxstrn: int = 50000,
ttsmult: float = 1.0,
ttsmax: float = 0,
species_names: list[str] = None,
extension: str = 'btn',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...Parameters:
ncomp (int): Number of mobile speciesmcomp (int): Total number of species (mobile + immobile)tunit (str): Time unit ('S', 'M', 'H', 'D', 'Y')lunit (str): Length unit ('CM', 'M', 'KM', 'FT', 'MI')munit (str): Mass unit ('G', 'KG', 'LB', 'TON')prsity (ArrayData): Effective porosityicbund (ArrayData): Boundary condition arraysconc (ArrayData): Starting concentrationAdvection package for simulating solute transport by groundwater flow.
class Mt3dAdv:
"""MT3DMS advection package"""
def __init__(
self,
model: Mt3dms,
mixelm: int = 3,
percel: float = 0.75,
mxpart: int = 800000,
nadvfd: int = 1,
itrack: int = 3,
wd: float = 0.5,
dceps: float = 1e-5,
nplane: int = 2,
npl: int = 10,
nph: int = 40,
npmin: int = 5,
npmax: int = 80,
nlsink: int = None,
npsink: int = None,
dchmoc: float = 1e-3,
extension: str = 'adv',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...Parameters:
mixelm (int): Advection solution method
percel (float): Courant number for advectionmxpart (int): Maximum number of particlesDispersion package for modeling mechanical dispersion and molecular diffusion.
class Mt3dDsp:
"""MT3DMS dispersion package"""
def __init__(
self,
model: Mt3dms,
al: ArrayData = 0.01,
trpt: ArrayData = 0.1,
trpv: ArrayData = 0.01,
dmcoef: ArrayData = 0.0,
multiDiff: bool = False,
extension: str = 'dsp',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...Parameters:
al (ArrayData): Longitudinal dispersivitytrpt (ArrayData): Ratio of horizontal transverse to longitudinal dispersivitytrpv (ArrayData): Ratio of vertical transverse to longitudinal dispersivitydmcoef (ArrayData): Molecular diffusion coefficientSource/sink mixing package for handling mass loading from boundary conditions.
class Mt3dSsm:
"""MT3DMS source/sink mixing package"""
def __init__(
self,
model: Mt3dms,
crch: ArrayData = None,
cevt: ArrayData = None,
mxss: int = None,
stress_period_data: dict = None,
dtype: np.dtype = None,
extension: str = 'ssm',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...Parameters:
crch (ArrayData): Concentration in rechargecevt (ArrayData): Concentration in evapotranspirationstress_period_data (dict): Source/sink data by stress period
css: Source/sink concentrationitype: Source/sink typeChemical reactions package for simulating decay, sorption, and other reactions.
class Mt3dRct:
"""MT3DMS chemical reactions package"""
def __init__(
self,
model: Mt3dms,
isothm: int = 0,
ireact: int = 0,
igetsc: int = 0,
rhob: ArrayData = 1.8,
prsity2: ArrayData = None,
srconc: ArrayData = None,
sp1: ArrayData = None,
sp2: ArrayData = None,
rc1: ArrayData = None,
rc2: ArrayData = None,
extension: str = 'rct',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...Parameters:
isothm (int): Sorption isotherm type
ireact (int): Kinetic rate reaction typerhob (ArrayData): Bulk densitysp1 (ArrayData): First sorption parametersp2 (ArrayData): Second sorption parameterrc1 (ArrayData): First-order reaction rateGeneralized conjugate gradient solver for the matrix equations.
class Mt3dGcg:
"""MT3DMS generalized conjugate gradient solver"""
def __init__(
self,
model: Mt3dms,
mxiter: int = 1,
iter1: int = 50,
isolve: int = 3,
ncrs: int = 0,
accl: float = 1.0,
cclose: float = 1e-5,
iprgcg: int = 0,
extension: str = 'gcg',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...Parameters:
mxiter (int): Maximum outer iterationsiter1 (int): Maximum inner iterationsisolve (int): Solver type
cclose (float): Closure criterionaccl (float): Acceleration parameterTransport observation package for model calibration and analysis.
class Mt3dTob:
"""MT3DMS transport observations package"""
def __init__(
self,
model: Mt3dms,
outnam: str = 'mt3d_obs.out',
CScale: float = 1.0,
FluxGroups: int = 0,
FScale: float = 1.0,
iOutFlux: int = 0,
obs_data: list = None,
extension: str = 'tob',
no_print: bool = False,
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...Streamflow transport package for stream-aquifer solute exchange.
class Mt3dSft:
"""MT3DMS streamflow transport package"""
def __init__(
self,
model: Mt3dms,
nsfinit: int = 0,
mxsfbc: int = 0,
icbcsf: int = 0,
ioutobs: int = 0,
ietsfr: int = 0,
isfsolv: int = 1,
cclosesf: float = 1e-6,
mxitersf: int = 10,
crntsf: float = 1.0,
iprtxmd: int = 0,
coldsf: ArrayData = 0.0,
dispsf: ArrayData = 0.0,
nobssf: int = 0,
obs_sf: list = None,
sf_stress_period_data: dict = None,
extension: str = 'sft',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...Lake transport package for lake-aquifer solute exchange.
class Mt3dLkt:
"""MT3DMS lake transport package"""
def __init__(
self,
model: Mt3dms,
nlkinit: int = 0,
mxlkbc: int = 0,
icbclk: int = 0,
ietlak: int = 0,
coldlk: ArrayData = 0.0,
lk_stress_period_data: dict = None,
extension: str = 'lkt',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...Unsaturated zone transport package.
class Mt3dUzt:
"""MT3DMS unsaturated zone transport package"""
def __init__(
self,
model: Mt3dms,
iuzfbnd: ArrayData = None,
cuzinf: ArrayData = 0.0,
cuzet: ArrayData = 0.0,
cgwet: ArrayData = 0.0,
extension: str = 'uzt',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...PHC3D reactions package for advanced geochemical modeling.
class Mt3dPhc:
"""MT3DMS PHC3D reactions package"""
def __init__(
self,
model: Mt3dms,
os: int = 2,
temp: ArrayData = 25.0,
asbin: ArrayData = None,
Prime: ArrayData = None,
extension: str = 'phc',
unitnumber: int = None,
filenames: str = None,
**kwargs
): ...import flopy
import numpy as np
# Create MODFLOW model first (required for transport)
mf = flopy.modflow.Modflow(modelname='transport_base')
# Basic MODFLOW setup
nlay, nrow, ncol = 1, 50, 100
dis = flopy.modflow.ModflowDis(
mf, nlay=nlay, nrow=nrow, ncol=ncol,
delr=10.0, delc=10.0, top=10.0, botm=0.0
)
bas = flopy.modflow.ModflowBas(mf, ibound=1, strt=5.0)
lpf = flopy.modflow.ModflowLpf(mf, hk=1.0, sy=0.1, ss=1e-5)
pcg = flopy.modflow.ModflowPcg(mf)
oc = flopy.modflow.ModflowOc(mf)
# Write and run MODFLOW
mf.write_input()
success, buff = mf.run_model()
# Create MT3DMS model
mt = flopy.mt3d.Mt3dms(
modelname='transport',
modflowmodel=mf,
ftlfilename='mt3d_link.ftl'
)
# Basic transport package
btn = flopy.mt3d.Mt3dBtn(
mt,
ncomp=1, # Single species
prsity=0.25,
icbund=1, # Active transport cells
sconc=0.0, # Starting concentration
timprs=np.array([100.0, 200.0, 500.0, 1000.0]), # Output times
nprs=4,
chkmas=True,
nprmas=1
)
# Advection package
adv = flopy.mt3d.Mt3dAdv(
mt,
mixelm=3, # TVD scheme
percel=0.75,
nadvfd=1
)
# Dispersion package
dsp = flopy.mt3d.Mt3dDsp(
mt,
al=0.5, # Longitudinal dispersivity
trpt=0.1, # Transverse/longitudinal ratio
dmcoef=1e-9 # Molecular diffusion
)
# Source/sink mixing - contaminant source
ssm_data = {}
# Injection well at center with contamination
ssm_data[0] = [
[0, 25, 50, 1000.0, 2] # Layer, row, col, concentration, type (2=well)
]
ssm = flopy.mt3d.Mt3dSsm(mt, stress_period_data=ssm_data)
# Solver
gcg = flopy.mt3d.Mt3dGcg(mt, mxiter=1, iter1=50, cclose=1e-5)
# Write and run MT3DMS
mt.write_input()
success, buff = mt.run_model()import flopy
import numpy as np
# Create base MODFLOW model (abbreviated)
mf = flopy.modflow.Modflow(modelname='reactive_transport')
# ... MODFLOW setup code ...
# Multi-species MT3DMS model
mt = flopy.mt3d.Mt3dms(
modelname='reactive_mt',
modflowmodel=mf
)
# Multiple species setup
ncomp = 3 # Parent compound + 2 decay products
species_names = ['PCE', 'TCE', 'DCE']
btn = flopy.mt3d.Mt3dBtn(
mt,
ncomp=ncomp,
mcomp=ncomp,
species_names=species_names,
prsity=0.30,
icbund=1,
sconc=[[0.0] * ncol * nrow] * ncomp, # Zero initial concentration
timprs=np.logspace(0, 3, 20), # Logarithmic time spacing
savucn=True
)
# Advection
adv = flopy.mt3d.Mt3dAdv(
mt,
mixelm=0, # MOC for accuracy with reactions
percel=0.5,
mxpart=100000
)
# Dispersion - different for each species
al = [1.0, 0.8, 0.6] # Decreasing dispersivity with molecular size
dsp = flopy.mt3d.Mt3dDsp(
mt,
al=al,
trpt=0.1,
dmcoef=[1e-9, 1.2e-9, 1.5e-9] # Increasing diffusion
)
# Reactions - sequential decay chain
# PCE -> TCE -> DCE
rct = flopy.mt3d.Mt3dRct(
mt,
isothm=1, # Linear sorption
ireact=1, # First-order decay
rhob=1.8, # Bulk density
sp1=[[0.5, 0.3, 0.1]], # Kd values (decreasing sorption)
rc1=[
[0.01, 0.0, 0.0], # PCE decay rate
[0.0, 0.02, 0.0], # TCE decay rate
[0.0, 0.0, 0.05] # DCE decay rate
],
rc2=[
[0.0, 0.01, 0.0], # PCE -> TCE yield
[0.0, 0.0, 0.02], # TCE -> DCE yield
[0.0, 0.0, 0.0] # DCE -> end products
]
)
# Source with declining concentration
ssm_data = {}
for kper in range(10): # 10 stress periods
# Exponentially declining source
conc = 1000.0 * np.exp(-0.1 * kper)
ssm_data[kper] = [
[0, 10, 20, conc, 2, 0.0, 0.0] # Only PCE injected
]
ssm = flopy.mt3d.Mt3dSsm(mt, stress_period_data=ssm_data)
# Solver
gcg = flopy.mt3d.Mt3dGcg(
mt,
mxiter=1,
iter1=100,
isolve=3,
cclose=1e-6
)
# Observations at multiple locations
obs_data = []
# Monitoring wells at different distances
for i, dist in enumerate([50, 100, 200]):
row = 25
col = 20 + dist // 10
obs_data.append(['MW_{}'.format(i+1), 0, row, col])
tob = flopy.mt3d.Mt3dTob(
mt,
obs_data=obs_data,
outnam='reactive_obs.out'
)
# Write and run
mt.write_input()
success, buff = mt.run_model()import flopy
import numpy as np
# MODFLOW model with SFR package (abbreviated setup)
mf = flopy.modflow.Modflow(modelname='stream_transport')
# ... discretization, basic, lpf packages ...
# Add SFR package for stream routing
# Stream reaches along model diagonal
reach_data = []
segment_data = {}
for i in range(50): # 50 reaches
row = i
col = i
reach_data.append([
1, i+1, row, col, 5.0, 1.0, 0.001, 2.0, 1.0, 0.03, 0.5, 0.5
])
segment_data[0] = {
1: [1, 1, 0, 100.0, 0.0, 0.0, 0.0, 'stream_seg']
}
sfr = flopy.modflow.ModflowSfr2(
mf,
nstrm=50,
nss=1,
reach_data=np.array(reach_data, dtype=flopy.modflow.ModflowSfr2.get_default_reach_dtype()),
segment_data=segment_data
)
# Run MODFLOW
mf.write_input()
success, buff = mf.run_model()
# MT3DMS with stream transport
mt = flopy.mt3d.Mt3dms(
modelname='stream_mt',
modflowmodel=mf
)
# Basic transport
btn = flopy.mt3d.Mt3dBtn(
mt,
ncomp=1,
prsity=0.25,
sconc=0.0,
timprs=np.array([10.0, 30.0, 100.0, 365.0])
)
# Advection and dispersion
adv = flopy.mt3d.Mt3dAdv(mt, mixelm=3)
dsp = flopy.mt3d.Mt3dDsp(mt, al=2.0, trpt=0.1)
# Stream transport package
sf_data = {}
# Initial stream concentration and dispersion
nsfinit = 50 # Number of stream reaches
coldsf = np.ones(nsfinit) * 500.0 # Initial stream concentration
dispsf = np.ones(nsfinit) * 10.0 # Stream dispersion coefficient
# Stream boundary conditions by stress period
sf_stress_data = {}
for kper in range(4):
# Upstream inflow concentration varies seasonally
if kper == 0: # Spring - high concentration
inflow_conc = 800.0
elif kper == 1: # Summer - medium concentration
inflow_conc = 400.0
elif kper == 2: # Fall - low concentration
inflow_conc = 200.0
else: # Winter - medium concentration
inflow_conc = 500.0
sf_stress_data[kper] = [
[1, 1, inflow_conc, 2] # Segment, reach, concentration, type
]
sft = flopy.mt3d.Mt3dSft(
mt,
nsfinit=nsfinit,
mxsfbc=20,
coldsf=coldsf,
dispsf=dispsf,
sf_stress_period_data=sf_stress_data
)
# Source/sink mixing for groundwater-stream interaction
ssm_data = {}
# Point source in aquifer
ssm_data[0] = [
[0, 10, 10, 2000.0, 2] # High concentration injection
]
ssm = flopy.mt3d.Mt3dSsm(mt, stress_period_data=ssm_data)
# Solver
gcg = flopy.mt3d.Mt3dGcg(mt)
# Run transport model
mt.write_input()
success, buff = mt.run_model()import flopy
import numpy as np
# Complex reactive transport with PHC3D
mf = flopy.modflow.Modflow(modelname='geochemical_model')
# ... MODFLOW setup ...
# Multi-component reactive system
mt = flopy.mt3d.Mt3dms(modelname='geochemical_mt', modflowmodel=mf)
# Multiple chemical species
ncomp = 6 # Multiple interacting species
species_names = ['Ca', 'CO3', 'H', 'CaCO3', 'CaHCO3', 'HCO3']
btn = flopy.mt3d.Mt3dBtn(
mt,
ncomp=ncomp,
mcomp=ncomp,
species_names=species_names,
prsity=0.25,
# Initial concentrations from equilibrium calculations
sconc=[
[1.0e-3] * (nrow * ncol), # Ca concentration
[1.0e-4] * (nrow * ncol), # CO3 concentration
[1.0e-7] * (nrow * ncol), # H concentration (pH 7)
[0.0] * (nrow * ncol), # CaCO3 (mineral)
[1.0e-4] * (nrow * ncol), # CaHCO3 complex
[1.0e-3] * (nrow * ncol) # HCO3 concentration
]
)
# Transport processes
adv = flopy.mt3d.Mt3dAdv(mt, mixelm=3)
dsp = flopy.mt3d.Mt3dDsp(mt, al=1.0, trpt=0.1)
# PHC3D for geochemical reactions
# Temperature field
temp = np.ones((nlay, nrow, ncol)) * 25.0 # 25°C
# Specify which species are in equilibrium vs kinetic
os_flag = 2 # Mixed equilibrium-kinetic system
phc = flopy.mt3d.Mt3dPhc(
mt,
os=os_flag,
temp=temp
)
# Source with varying chemistry
ssm_data = {}
for kper in range(5):
# Injection with different pH values
ca_conc = 2.0e-3
co3_conc = 5.0e-4
h_conc = 1.0e-6 * (10 ** (kper - 2)) # pH 4-8 range
ssm_data[kper] = [
[0, 25, 50, ca_conc, 2, co3_conc, h_conc, 0.0, 0.0, 0.0]
]
ssm = flopy.mt3d.Mt3dSsm(mt, stress_period_data=ssm_data)
# Specialized solver for reactive transport
gcg = flopy.mt3d.Mt3dGcg(
mt,
mxiter=1,
iter1=200, # More iterations for convergence
isolve=3,
cclose=1e-7 # Tighter convergence
)
# Write and run
mt.write_input()
success, buff = mt.run_model()# Transport-specific types
ConcentrationArray = Union[float, list[float], np.ndarray]
ReactionParameter = Union[float, list[float], np.ndarray]
SpeciesData = list[ConcentrationArray]
# Boundary condition types
SourceSinkData = list[list[Union[int, float]]]
StreamTransportData = dict[int, list[list[Union[int, float]]]]
ObservationData = list[list[Union[str, int, float]]]
# Reaction types
IsothermType = Literal[0, 1, 2, 3, 4] # None, Linear, Freundlich, Langmuir, Kinetic
ReactionType = int
SolverType = Literal[1, 2, 3] # Jacobi, SSOR, MIC
# Advection methods
AdvectionMethod = Literal[-1, 0, 1, 2, 3] # FD, MOC, MMOC, HMOC, TVD
# Time and mass units
TimeUnit = Literal['S', 'M', 'H', 'D', 'Y']
LengthUnit = Literal['CM', 'M', 'KM', 'FT', 'MI']
MassUnit = Literal['G', 'KG', 'LB', 'TON']
# File format types
FileFormat = Literal['formatted', 'unformatted']
OutputFormat = intThis comprehensive documentation covers the complete MT3DMS transport modeling API including all packages, reaction types, and usage patterns. The examples demonstrate basic to advanced transport scenarios including multi-species systems, reactive transport, and coupled surface water-groundwater transport.
Install with Tessl CLI
npx tessl i tessl/pypi-flopydocs
evals
scenario-1
scenario-2
scenario-3
scenario-4
scenario-5
scenario-6
scenario-7
scenario-8
scenario-9
scenario-10