Skip to content

mibitrans.transport API reference

model_parent

Results

Object that holds model results and input parameters for individual runs.

Source code in mibitrans/transport/model_parent.py
class Results:
    """Object that holds model results and input parameters for individual runs."""

    def __init__(self, model):
        """Records input parameters and resulting output based given model.

        Args:
            model (Transport3D): Model object from which to initialize results. Should be child class of Transport3D.

        Properties:
            model_type (Transport3D) : Class instance of model used to generate results.
            short_description (str) : Short description of model.
            x (np.ndarray) : Numpy array with model x (longitudinal direction) discretization, corresponding to
                model_parameters, with distance in [m].
            y (np.ndarray) : Numpy array with model y (transverse horizontal direction) discretization, corresponding to
                model_parameters, with distance in [m].
            t (np.ndarray) : Numpy array with model t (time) discretization, corresponding to
                model_parameters, with time in [days].
            hydrological_parameters (HydrologicalParameters) : Dataclass holding the hydrological parameters used to
                run the model.
            attenuation_parameters (AttenuationParameters) : Dataclass holding the attenuation parameters used to run
                the model.
            source_parameters (SourceParameters) : Dataclass holding the source parameters used to run the model.
            model_parameters (ModelParameters): Dataclass holding the model parameters used to run the model.
            electron_acceptors (ElectronAcceptors): Dataclass holding the electron acceptor concentrations used to run
                the model. Only for instant reaction, None for other models.
            utilization_factor (UtilizationFactor): Dataclass holding the electron acceptor utilization factors used to
                run the model. Only for instant reaction, None for other models.
            mode (str) : Model mode of the used model. Either 'linear' or 'instant_reaction'
            rv (float) : Retarded flow velocity, as v / R [m/day].
            k_source (float) : Source depletion rate [1/days]. For infinite source mass, k_source = 0, and therefore, no
                source depletion takes place.
            c_source (np.ndarray) : Initial nett source zone concentrations. For multiple source zones, nett
                concentration in nth source zone is original concentration minus concentration in source zone n - 1. For
                instant reaction model, the biodegradation capacity is added to the outermost source zone.
            biodegradation_capacity (float) : Maximum capacity of biodegradation taking place, based on electron
                acceptor concentrations and utilization factor.
            cxyt (np.ndarray) : Three-dimensional numpy array with concentrations for all x, y and t positions. Indexed
             as cxyt[t,y,x]. In [g/m3].
            relative_cxyt (np.ndarray) : Three-dimensional numpy array with relative concentrations for all x, y and t
                positions. Compared to maximum source zone concentrations.
            cxyt_noBC (np.ndarray) : Three-dimensional numpy array with concentrations for all x, y and t of instant
                reaction models, without subtracting the biodegradation capacity, in [g/m3].
            input_parameters (dict) : Dictionary of input parameter dataclasses for the model. Does not include instant
                reaction parameters.

        Methods:
            centerline : Plot center of contaminant plume, at a specified time and y position.
            transverse : Plot concentration distribution as a line horizontal transverse to the plume extent.
            breakthrough : Plot contaminant breakthrough curve at given x and y position in model domain.
            plume_2d : Plot contaminant plume as a 2D colormesh, at a specified time.
            plume_3d : Plot contaminant plume as a 3D surface, at a specified time.
            mass_balance : Return a mass balance object with source and plume characteristics at given time(s).
        """
        self._model_type = model.__class__
        self._short_description = model.short_description
        self._x = model.x
        self._y = model.y
        self._t = model.t

        # All properties of Transport3D that are objects should be copied; if not copied, changing them in the class
        # object where they originated will also change them here, which is not the intended behaviour.
        self._hydrological_parameters = copy.copy(model.hydrological_parameters)
        self._attenuation_parameters = copy.copy(model.attenuation_parameters)
        self._source_parameters = copy.copy(model.source_parameters)
        self._model_parameters = copy.copy(model.model_parameters)
        self._electron_acceptors = copy.copy(model._electron_acceptors)
        self._utilization_factor = copy.copy(model._utilization_factor)

        self._mode = model.mode
        self._rv = model.rv
        self._k_source = model.k_source
        self._c_source = model.c_source
        self._biodegradation_capacity = model.biodegradation_capacity

        self._cxyt = model.cxyt
        self._relative_cxyt = model.relative_cxyt
        self._cxyt_noBC = model.cxyt_noBC

    @property
    def model_type(self):
        """Class object of the model that generated the results."""
        return self._model_type

    @property
    def short_description(self):
        """Short description of the model that generated the results."""
        return self._short_description

    @property
    def x(self):
        """Model x discretization array."""
        return self._x

    @property
    def y(self):
        """Model y discretization array."""
        return self._y

    @property
    def t(self):
        """Model t discretization array."""
        return self._t

    @property
    def hydrological_parameters(self):
        """Hydrological parameters of the model used for the results."""
        return self._hydrological_parameters

    @property
    def attenuation_parameters(self):
        """Attenuation parameters of the model used for the results."""
        return self._attenuation_parameters

    @property
    def source_parameters(self):
        """Source parameters of the model used for the results."""
        return self._source_parameters

    @property
    def model_parameters(self):
        """Space-time discretization parameters of the model used for the results."""
        return self._model_parameters

    @property
    def electron_acceptors(self):
        """Electron acceptor/byproduct concentrations of the model used for the results."""
        return self._electron_acceptors

    @property
    def utilization_factor(self):
        """Utilization factor of the model used for the results."""
        return self._utilization_factor

    @property
    def mode(self):
        """Model mode used for running the model."""
        return self._mode

    @property
    def rv(self):
        """Retarded flow velocity used in the model."""
        return self._rv

    @property
    def k_source(self):
        """Source depletion rate used in the model."""
        return self._k_source

    @property
    def c_source(self):
        """Nett source zone concentration used in the model."""
        return self._c_source

    @property
    def biodegradation_capacity(self):
        """Biodegradation capacity of the model used for the results. Only for instant reaction models."""
        return self._biodegradation_capacity

    @property
    def cxyt(self):
        """Modelled concentration for all x, y and t, using the input parameters present in this object."""
        return self._cxyt

    @property
    def relative_cxyt(self):
        """Modelled concentration for all x, y and t, divided by the maximum source zone concentration."""
        return self._relative_cxyt

    @property
    def cxyt_noBC(self):
        """Concentration in domain without subtracting biodegradation capacity, in the instant reaction model."""
        return self._cxyt_noBC

    @property
    def input_parameters(self):
        """Return the input arguments for the model in the form of a dictionary, based on current values."""
        return dict(
            hydrological_parameters=self.hydrological_parameters,
            attenuation_parameters=self.attenuation_parameters,
            source_parameters=self.source_parameters,
            model_parameters=self.model_parameters,
        )

    def centerline(self, y_position=0, time=None, relative_concentration=False, animate=False, **kwargs):
        """Plot center of contaminant plume of this model, at a specified time and y position.

        Args:
            y_position (float, optional): y-position across the plume (transverse horizontal direction) for the plot.
                By default, the center of the plume at y=0 is plotted.
            time (float, optional): Point of time for the plot. Will show the closest time step to given value.
                By default, last point in time is plotted.
            relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
                source zone concentrations at t=0. By default, absolute concentrations are shown.
            animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
                False.
            **kwargs : Arguments to be passed to plt.plot().

        """
        if animate:
            anim = pline.centerline(
                self,
                y_position=y_position,
                time=time,
                relative_concentration=relative_concentration,
                animate=animate,
                **kwargs,
            )
            return anim
        else:
            pline.centerline(
                self,
                y_position=y_position,
                time=time,
                relative_concentration=relative_concentration,
                animate=animate,
                **kwargs,
            )
            return None

    def transverse(self, x_position, time=None, relative_concentration=False, animate=False, **kwargs):
        """Plot concentration distribution as a line horizontal transverse to the plume extent.

        Args:
            x_position : x-position along the plume (longitudinal direction) for the plot.
            time (float): Point of time for the plot. Will show the closest time step to given value.
                By default, last point in time is plotted.
            relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
                source zone concentrations at t=0. By default, absolute concentrations are shown.
            animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
                False.
            **kwargs : Arguments to be passed to plt.plot().
        """
        if animate:
            anim = pline.transverse(
                self,
                x_position=x_position,
                time=time,
                relative_concentration=relative_concentration,
                animate=animate,
                **kwargs,
            )
            return anim
        else:
            pline.transverse(
                self,
                x_position=x_position,
                time=time,
                relative_concentration=relative_concentration,
                animate=animate,
                **kwargs,
            )
            return None

    def breakthrough(self, x_position, y_position=0, relative_concentration=False, animate=False, **kwargs):
        """Plot contaminant breakthrough curve at given x and y position in model domain.

        Args:
            x_position : x-position along the plume (longitudinal direction).
            y_position : y-position across the plume (transverse horizontal direction).
                By default, at the center of the plume (at y=0).
            relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
                source zone concentrations at t=0. By default, absolute concentrations are shown.
            animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
                False.
            **kwargs : Arguments to be passed to plt.plot().
        """
        if animate:
            anim = pline.breakthrough(
                self,
                x_position=x_position,
                y_position=y_position,
                relative_concentration=relative_concentration,
                animate=animate,
                **kwargs,
            )
            return anim
        else:
            pline.breakthrough(
                self,
                x_position=x_position,
                y_position=y_position,
                relative_concentration=relative_concentration,
                animate=animate,
                **kwargs,
            )
            return None

    def plume_2d(self, time=None, relative_concentration=False, animate=False, **kwargs):
        """Plot contaminant plume as a 2D colormesh, at a specified time.

        Args:
            time (float): Point of time for the plot. Will show the closest time step to given value.
                By default, last point in time is plotted.
            relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
                source zone concentrations at t=0. By default, absolute concentrations are shown.
            animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
                False.
            **kwargs : Arguments to be passed to plt.pcolormesh().

        Returns a matrix plot of the input plume as object.
        """
        anim = psurf.plume_2d(self, time=time, relative_concentration=relative_concentration, animate=animate, **kwargs)
        return anim

    def plume_3d(self, time=None, relative_concentration=False, animate=False, **kwargs):
        """Plot contaminant plume as a 3D surface, at a specified time.

        Args:
            time (float): Point of time for the plot. Will show the closest time step to given value.
                By default, last point in time is plotted.
            relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
                source zone concentrations at t=0. By default, absolute concentrations are shown.
            animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
                False.
            **kwargs : Arguments to be passed to plt.plot_surface().

        Returns:
            ax (matplotlib.axes._axes.Axes) : Matplotlib Axes object of plume plot.
                or if animate == True
            anim (matplotib.animation.FuncAnimation) : Matplotlib FuncAnimation object of plume plot.
        """
        ax_or_anim = psurf.plume_3d(
            self, time=time, relative_concentration=relative_concentration, animate=animate, **kwargs
        )
        return ax_or_anim

    def mass_balance(self, time="all", verbose=False):
        """Return a mass balance object with source and plume characteristics at given time(s).

        Args:
            time (float | str): Time at which to initially calculate the mass balance. Either as a value between 0 and
                model end time. Or as 'all', which will calculate mass balance attributes for each time step as arrays.
            verbose (bool, optional): Verbose mode. Defaults to False.

        Returns:
            mass_balance_object: Object of the MassBalance class. Output is accessed through object properties. Can be
                called to change the time of the mass balance.

        The mass balance object has the following properties:
            plume_mass: Mass of the contaminant plume inside the model extent, at the given time(s), in [g].
            source_mass: Mass of the contaminant source at the given time(s), in [g]. No values are given for models
                with infinite source mass.
            delta_source: Difference in mass between contaminant source at given time and source at t = 0, in [g].
            degraded_mass: Mass of plume contaminant degradation at the given time(s), compared to a model without
                degradation, in [g]. Has no value if model does not consider degradation.
            model_without_degradation: Object of model without degradation. Has no value if model does not consider
                degradation.
            instant_reaction_degraded_mass(self): Difference in plume mass instant reaction with and without
                biodegradation capacity subtracted, in [g].
            electron_acceptor_change(self): Change in electron acceptor/byproduct masses at the given time(s), in [g].
                Only for instant reaction.
        """
        return MassBalance(self, time, verbose)

attenuation_parameters property

Attenuation parameters of the model used for the results.

biodegradation_capacity property

Biodegradation capacity of the model used for the results. Only for instant reaction models.

c_source property

Nett source zone concentration used in the model.

cxyt property

Modelled concentration for all x, y and t, using the input parameters present in this object.

cxyt_noBC property

Concentration in domain without subtracting biodegradation capacity, in the instant reaction model.

electron_acceptors property

Electron acceptor/byproduct concentrations of the model used for the results.

hydrological_parameters property

Hydrological parameters of the model used for the results.

input_parameters property

Return the input arguments for the model in the form of a dictionary, based on current values.

k_source property

Source depletion rate used in the model.

mode property

Model mode used for running the model.

model_parameters property

Space-time discretization parameters of the model used for the results.

model_type property

Class object of the model that generated the results.

relative_cxyt property

Modelled concentration for all x, y and t, divided by the maximum source zone concentration.

rv property

Retarded flow velocity used in the model.

short_description property

Short description of the model that generated the results.

source_parameters property

Source parameters of the model used for the results.

t property

Model t discretization array.

utilization_factor property

Utilization factor of the model used for the results.

x property

Model x discretization array.

y property

Model y discretization array.

__init__(model)

Records input parameters and resulting output based given model.

Parameters:

Name Type Description Default
model Transport3D

Model object from which to initialize results. Should be child class of Transport3D.

required
Properties

model_type (Transport3D) : Class instance of model used to generate results. short_description (str) : Short description of model. x (np.ndarray) : Numpy array with model x (longitudinal direction) discretization, corresponding to model_parameters, with distance in [m]. y (np.ndarray) : Numpy array with model y (transverse horizontal direction) discretization, corresponding to model_parameters, with distance in [m]. t (np.ndarray) : Numpy array with model t (time) discretization, corresponding to model_parameters, with time in [days]. hydrological_parameters (HydrologicalParameters) : Dataclass holding the hydrological parameters used to run the model. attenuation_parameters (AttenuationParameters) : Dataclass holding the attenuation parameters used to run the model. source_parameters (SourceParameters) : Dataclass holding the source parameters used to run the model. model_parameters (ModelParameters): Dataclass holding the model parameters used to run the model. electron_acceptors (ElectronAcceptors): Dataclass holding the electron acceptor concentrations used to run the model. Only for instant reaction, None for other models. utilization_factor (UtilizationFactor): Dataclass holding the electron acceptor utilization factors used to run the model. Only for instant reaction, None for other models. mode (str) : Model mode of the used model. Either ‘linear’ or ‘instant_reaction’ rv (float) : Retarded flow velocity, as v / R [m/day]. k_source (float) : Source depletion rate [1/days]. For infinite source mass, k_source = 0, and therefore, no source depletion takes place. c_source (np.ndarray) : Initial nett source zone concentrations. For multiple source zones, nett concentration in nth source zone is original concentration minus concentration in source zone n - 1. For instant reaction model, the biodegradation capacity is added to the outermost source zone. biodegradation_capacity (float) : Maximum capacity of biodegradation taking place, based on electron acceptor concentrations and utilization factor. cxyt (np.ndarray) : Three-dimensional numpy array with concentrations for all x, y and t positions. Indexed as cxyt[t,y,x]. In [g/m3]. relative_cxyt (np.ndarray) : Three-dimensional numpy array with relative concentrations for all x, y and t positions. Compared to maximum source zone concentrations. cxyt_noBC (np.ndarray) : Three-dimensional numpy array with concentrations for all x, y and t of instant reaction models, without subtracting the biodegradation capacity, in [g/m3]. input_parameters (dict) : Dictionary of input parameter dataclasses for the model. Does not include instant reaction parameters.

Functions:

Name Description
centerline

Plot center of contaminant plume, at a specified time and y position.

transverse

Plot concentration distribution as a line horizontal transverse to the plume extent.

breakthrough

Plot contaminant breakthrough curve at given x and y position in model domain.

plume_2d

Plot contaminant plume as a 2D colormesh, at a specified time.

plume_3d

Plot contaminant plume as a 3D surface, at a specified time.

mass_balance

Return a mass balance object with source and plume characteristics at given time(s).

Source code in mibitrans/transport/model_parent.py
def __init__(self, model):
    """Records input parameters and resulting output based given model.

    Args:
        model (Transport3D): Model object from which to initialize results. Should be child class of Transport3D.

    Properties:
        model_type (Transport3D) : Class instance of model used to generate results.
        short_description (str) : Short description of model.
        x (np.ndarray) : Numpy array with model x (longitudinal direction) discretization, corresponding to
            model_parameters, with distance in [m].
        y (np.ndarray) : Numpy array with model y (transverse horizontal direction) discretization, corresponding to
            model_parameters, with distance in [m].
        t (np.ndarray) : Numpy array with model t (time) discretization, corresponding to
            model_parameters, with time in [days].
        hydrological_parameters (HydrologicalParameters) : Dataclass holding the hydrological parameters used to
            run the model.
        attenuation_parameters (AttenuationParameters) : Dataclass holding the attenuation parameters used to run
            the model.
        source_parameters (SourceParameters) : Dataclass holding the source parameters used to run the model.
        model_parameters (ModelParameters): Dataclass holding the model parameters used to run the model.
        electron_acceptors (ElectronAcceptors): Dataclass holding the electron acceptor concentrations used to run
            the model. Only for instant reaction, None for other models.
        utilization_factor (UtilizationFactor): Dataclass holding the electron acceptor utilization factors used to
            run the model. Only for instant reaction, None for other models.
        mode (str) : Model mode of the used model. Either 'linear' or 'instant_reaction'
        rv (float) : Retarded flow velocity, as v / R [m/day].
        k_source (float) : Source depletion rate [1/days]. For infinite source mass, k_source = 0, and therefore, no
            source depletion takes place.
        c_source (np.ndarray) : Initial nett source zone concentrations. For multiple source zones, nett
            concentration in nth source zone is original concentration minus concentration in source zone n - 1. For
            instant reaction model, the biodegradation capacity is added to the outermost source zone.
        biodegradation_capacity (float) : Maximum capacity of biodegradation taking place, based on electron
            acceptor concentrations and utilization factor.
        cxyt (np.ndarray) : Three-dimensional numpy array with concentrations for all x, y and t positions. Indexed
         as cxyt[t,y,x]. In [g/m3].
        relative_cxyt (np.ndarray) : Three-dimensional numpy array with relative concentrations for all x, y and t
            positions. Compared to maximum source zone concentrations.
        cxyt_noBC (np.ndarray) : Three-dimensional numpy array with concentrations for all x, y and t of instant
            reaction models, without subtracting the biodegradation capacity, in [g/m3].
        input_parameters (dict) : Dictionary of input parameter dataclasses for the model. Does not include instant
            reaction parameters.

    Methods:
        centerline : Plot center of contaminant plume, at a specified time and y position.
        transverse : Plot concentration distribution as a line horizontal transverse to the plume extent.
        breakthrough : Plot contaminant breakthrough curve at given x and y position in model domain.
        plume_2d : Plot contaminant plume as a 2D colormesh, at a specified time.
        plume_3d : Plot contaminant plume as a 3D surface, at a specified time.
        mass_balance : Return a mass balance object with source and plume characteristics at given time(s).
    """
    self._model_type = model.__class__
    self._short_description = model.short_description
    self._x = model.x
    self._y = model.y
    self._t = model.t

    # All properties of Transport3D that are objects should be copied; if not copied, changing them in the class
    # object where they originated will also change them here, which is not the intended behaviour.
    self._hydrological_parameters = copy.copy(model.hydrological_parameters)
    self._attenuation_parameters = copy.copy(model.attenuation_parameters)
    self._source_parameters = copy.copy(model.source_parameters)
    self._model_parameters = copy.copy(model.model_parameters)
    self._electron_acceptors = copy.copy(model._electron_acceptors)
    self._utilization_factor = copy.copy(model._utilization_factor)

    self._mode = model.mode
    self._rv = model.rv
    self._k_source = model.k_source
    self._c_source = model.c_source
    self._biodegradation_capacity = model.biodegradation_capacity

    self._cxyt = model.cxyt
    self._relative_cxyt = model.relative_cxyt
    self._cxyt_noBC = model.cxyt_noBC

breakthrough(x_position, y_position=0, relative_concentration=False, animate=False, **kwargs)

Plot contaminant breakthrough curve at given x and y position in model domain.

Parameters:

Name Type Description Default
x_position

x-position along the plume (longitudinal direction).

required
y_position

y-position across the plume (transverse horizontal direction). By default, at the center of the plume (at y=0).

required
relative_concentration (bool, optional)

If set to True, will plot concentrations relative to maximum source zone concentrations at t=0. By default, absolute concentrations are shown.

required
animate bool

If True, animation of contaminant plume until given time is shown. Default is False.

False
**kwargs

Arguments to be passed to plt.plot().

required
Source code in mibitrans/transport/model_parent.py
def breakthrough(self, x_position, y_position=0, relative_concentration=False, animate=False, **kwargs):
    """Plot contaminant breakthrough curve at given x and y position in model domain.

    Args:
        x_position : x-position along the plume (longitudinal direction).
        y_position : y-position across the plume (transverse horizontal direction).
            By default, at the center of the plume (at y=0).
        relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
            source zone concentrations at t=0. By default, absolute concentrations are shown.
        animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
            False.
        **kwargs : Arguments to be passed to plt.plot().
    """
    if animate:
        anim = pline.breakthrough(
            self,
            x_position=x_position,
            y_position=y_position,
            relative_concentration=relative_concentration,
            animate=animate,
            **kwargs,
        )
        return anim
    else:
        pline.breakthrough(
            self,
            x_position=x_position,
            y_position=y_position,
            relative_concentration=relative_concentration,
            animate=animate,
            **kwargs,
        )
        return None

centerline(y_position=0, time=None, relative_concentration=False, animate=False, **kwargs)

Plot center of contaminant plume of this model, at a specified time and y position.

Parameters:

Name Type Description Default
y_position float

y-position across the plume (transverse horizontal direction) for the plot. By default, the center of the plume at y=0 is plotted.

0
time float

Point of time for the plot. Will show the closest time step to given value. By default, last point in time is plotted.

None
relative_concentration (bool, optional)

If set to True, will plot concentrations relative to maximum source zone concentrations at t=0. By default, absolute concentrations are shown.

required
animate bool

If True, animation of contaminant plume until given time is shown. Default is False.

False
**kwargs

Arguments to be passed to plt.plot().

required
Source code in mibitrans/transport/model_parent.py
def centerline(self, y_position=0, time=None, relative_concentration=False, animate=False, **kwargs):
    """Plot center of contaminant plume of this model, at a specified time and y position.

    Args:
        y_position (float, optional): y-position across the plume (transverse horizontal direction) for the plot.
            By default, the center of the plume at y=0 is plotted.
        time (float, optional): Point of time for the plot. Will show the closest time step to given value.
            By default, last point in time is plotted.
        relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
            source zone concentrations at t=0. By default, absolute concentrations are shown.
        animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
            False.
        **kwargs : Arguments to be passed to plt.plot().

    """
    if animate:
        anim = pline.centerline(
            self,
            y_position=y_position,
            time=time,
            relative_concentration=relative_concentration,
            animate=animate,
            **kwargs,
        )
        return anim
    else:
        pline.centerline(
            self,
            y_position=y_position,
            time=time,
            relative_concentration=relative_concentration,
            animate=animate,
            **kwargs,
        )
        return None

mass_balance(time='all', verbose=False)

Return a mass balance object with source and plume characteristics at given time(s).

Parameters:

Name Type Description Default
time float | str

Time at which to initially calculate the mass balance. Either as a value between 0 and model end time. Or as ‘all’, which will calculate mass balance attributes for each time step as arrays.

'all'
verbose bool

Verbose mode. Defaults to False.

False

Returns:

Name Type Description
mass_balance_object

Object of the MassBalance class. Output is accessed through object properties. Can be called to change the time of the mass balance.

The mass balance object has the following properties

plume_mass: Mass of the contaminant plume inside the model extent, at the given time(s), in [g]. source_mass: Mass of the contaminant source at the given time(s), in [g]. No values are given for models with infinite source mass. delta_source: Difference in mass between contaminant source at given time and source at t = 0, in [g]. degraded_mass: Mass of plume contaminant degradation at the given time(s), compared to a model without degradation, in [g]. Has no value if model does not consider degradation. model_without_degradation: Object of model without degradation. Has no value if model does not consider degradation. instant_reaction_degraded_mass(self): Difference in plume mass instant reaction with and without biodegradation capacity subtracted, in [g]. electron_acceptor_change(self): Change in electron acceptor/byproduct masses at the given time(s), in [g]. Only for instant reaction.

Source code in mibitrans/transport/model_parent.py
def mass_balance(self, time="all", verbose=False):
    """Return a mass balance object with source and plume characteristics at given time(s).

    Args:
        time (float | str): Time at which to initially calculate the mass balance. Either as a value between 0 and
            model end time. Or as 'all', which will calculate mass balance attributes for each time step as arrays.
        verbose (bool, optional): Verbose mode. Defaults to False.

    Returns:
        mass_balance_object: Object of the MassBalance class. Output is accessed through object properties. Can be
            called to change the time of the mass balance.

    The mass balance object has the following properties:
        plume_mass: Mass of the contaminant plume inside the model extent, at the given time(s), in [g].
        source_mass: Mass of the contaminant source at the given time(s), in [g]. No values are given for models
            with infinite source mass.
        delta_source: Difference in mass between contaminant source at given time and source at t = 0, in [g].
        degraded_mass: Mass of plume contaminant degradation at the given time(s), compared to a model without
            degradation, in [g]. Has no value if model does not consider degradation.
        model_without_degradation: Object of model without degradation. Has no value if model does not consider
            degradation.
        instant_reaction_degraded_mass(self): Difference in plume mass instant reaction with and without
            biodegradation capacity subtracted, in [g].
        electron_acceptor_change(self): Change in electron acceptor/byproduct masses at the given time(s), in [g].
            Only for instant reaction.
    """
    return MassBalance(self, time, verbose)

plume_2d(time=None, relative_concentration=False, animate=False, **kwargs)

Plot contaminant plume as a 2D colormesh, at a specified time.

Parameters:

Name Type Description Default
time float

Point of time for the plot. Will show the closest time step to given value. By default, last point in time is plotted.

None
relative_concentration (bool, optional)

If set to True, will plot concentrations relative to maximum source zone concentrations at t=0. By default, absolute concentrations are shown.

required
animate bool

If True, animation of contaminant plume until given time is shown. Default is False.

False
**kwargs

Arguments to be passed to plt.pcolormesh().

required

Returns a matrix plot of the input plume as object.

Source code in mibitrans/transport/model_parent.py
def plume_2d(self, time=None, relative_concentration=False, animate=False, **kwargs):
    """Plot contaminant plume as a 2D colormesh, at a specified time.

    Args:
        time (float): Point of time for the plot. Will show the closest time step to given value.
            By default, last point in time is plotted.
        relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
            source zone concentrations at t=0. By default, absolute concentrations are shown.
        animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
            False.
        **kwargs : Arguments to be passed to plt.pcolormesh().

    Returns a matrix plot of the input plume as object.
    """
    anim = psurf.plume_2d(self, time=time, relative_concentration=relative_concentration, animate=animate, **kwargs)
    return anim

plume_3d(time=None, relative_concentration=False, animate=False, **kwargs)

Plot contaminant plume as a 3D surface, at a specified time.

Parameters:

Name Type Description Default
time float

Point of time for the plot. Will show the closest time step to given value. By default, last point in time is plotted.

None
relative_concentration (bool, optional)

If set to True, will plot concentrations relative to maximum source zone concentrations at t=0. By default, absolute concentrations are shown.

required
animate bool

If True, animation of contaminant plume until given time is shown. Default is False.

False
**kwargs

Arguments to be passed to plt.plot_surface().

required

Returns:

Type Description

ax (matplotlib.axes._axes.Axes) : Matplotlib Axes object of plume plot. or if animate == True

anim (matplotib.animation.FuncAnimation) : Matplotlib FuncAnimation object of plume plot.

Source code in mibitrans/transport/model_parent.py
def plume_3d(self, time=None, relative_concentration=False, animate=False, **kwargs):
    """Plot contaminant plume as a 3D surface, at a specified time.

    Args:
        time (float): Point of time for the plot. Will show the closest time step to given value.
            By default, last point in time is plotted.
        relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
            source zone concentrations at t=0. By default, absolute concentrations are shown.
        animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
            False.
        **kwargs : Arguments to be passed to plt.plot_surface().

    Returns:
        ax (matplotlib.axes._axes.Axes) : Matplotlib Axes object of plume plot.
            or if animate == True
        anim (matplotib.animation.FuncAnimation) : Matplotlib FuncAnimation object of plume plot.
    """
    ax_or_anim = psurf.plume_3d(
        self, time=time, relative_concentration=relative_concentration, animate=animate, **kwargs
    )
    return ax_or_anim

transverse(x_position, time=None, relative_concentration=False, animate=False, **kwargs)

Plot concentration distribution as a line horizontal transverse to the plume extent.

Parameters:

Name Type Description Default
x_position

x-position along the plume (longitudinal direction) for the plot.

required
time float

Point of time for the plot. Will show the closest time step to given value. By default, last point in time is plotted.

None
relative_concentration (bool, optional)

If set to True, will plot concentrations relative to maximum source zone concentrations at t=0. By default, absolute concentrations are shown.

required
animate bool

If True, animation of contaminant plume until given time is shown. Default is False.

False
**kwargs

Arguments to be passed to plt.plot().

required
Source code in mibitrans/transport/model_parent.py
def transverse(self, x_position, time=None, relative_concentration=False, animate=False, **kwargs):
    """Plot concentration distribution as a line horizontal transverse to the plume extent.

    Args:
        x_position : x-position along the plume (longitudinal direction) for the plot.
        time (float): Point of time for the plot. Will show the closest time step to given value.
            By default, last point in time is plotted.
        relative_concentration (bool, optional) : If set to True, will plot concentrations relative to maximum
            source zone concentrations at t=0. By default, absolute concentrations are shown.
        animate (bool, optional): If True, animation of contaminant plume until given time is shown. Default is
            False.
        **kwargs : Arguments to be passed to plt.plot().
    """
    if animate:
        anim = pline.transverse(
            self,
            x_position=x_position,
            time=time,
            relative_concentration=relative_concentration,
            animate=animate,
            **kwargs,
        )
        return anim
    else:
        pline.transverse(
            self,
            x_position=x_position,
            time=time,
            relative_concentration=relative_concentration,
            animate=animate,
            **kwargs,
        )
        return None

Transport3D

Bases: ABC

Parent class for all 3-dimensional analytical solutions.

Source code in mibitrans/transport/model_parent.py
class Transport3D(ABC):
    """Parent class for all 3-dimensional analytical solutions."""

    def __init__(
        self, hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose=False
    ):
        """Initialize parent class object.

        Args:
            hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters) : Dataclass object containing
                hydrological parameters from HydrologicalParameters.
            attenuation_parameters (mibitrans.data.read.AttenuationParameters) : Dataclass object containing adsorption,
                degradation and diffusion parameters from AttenuationParameters.
            source_parameters (mibitrans.data.read.SourceParameters) : Dataclass object containing source parameters
                from SourceParameters.
            model_parameters (mibitrans.data.read.ModelParameters) : Dataclass object containing model parameters from
                ModelParameters.
            verbose (bool, optional): Verbose mode. Defaults to False.
        """
        # Check if input arguments are of the correct dataclass
        for key, value in locals().items():
            if key not in ["self", "verbose"]:
                self._check_input_dataclasses(key, value)

        self._hyd_pars = copy.copy(hydrological_parameters)
        self._att_pars = copy.copy(attenuation_parameters)
        self._src_pars = copy.copy(source_parameters)
        self._mod_pars = copy.copy(model_parameters)
        self._decay_rate = self._att_pars.decay_rate

        self.verbose = verbose
        self._mode = "linear"
        self._electron_acceptors = None
        self._utilization_factor = None
        self.biodegradation_capacity = None
        self.cxyt_noBC = None
        self._pre_run_initialization_parameters()

    @property
    def input_parameters(self):
        """Return the input arguments for the model in the form of a dictionary, based on current values."""
        return dict(
            hydrological_parameters=self.hydrological_parameters,
            attenuation_parameters=self.attenuation_parameters,
            source_parameters=self.source_parameters,
            model_parameters=self.model_parameters,
            verbose=self.verbose,
        )

    @property
    def hydrological_parameters(self):
        """Rename to shorthand form of hydrological_parameters inside class for ease of use."""
        return self._hyd_pars

    @hydrological_parameters.setter
    def hydrological_parameters(self, value):
        self._check_input_dataclasses("hydrological_parameters", value)
        self._hyd_pars = copy.copy(value)

    @property
    def attenuation_parameters(self):
        """Rename to shorthand form of attenuation_parameters inside class for ease of use."""
        return self._att_pars

    @attenuation_parameters.setter
    def attenuation_parameters(self, value):
        self._check_input_dataclasses("attenuation_parameters", value)
        self._att_pars = copy.copy(value)

    @property
    def source_parameters(self):
        """Rename to shorthand form of source_parameters inside class for ease of use."""
        return self._src_pars

    @source_parameters.setter
    def source_parameters(self, value):
        self._check_input_dataclasses("source_parameters", value)
        self._src_pars = copy.copy(value)

    @property
    def model_parameters(self):
        """Rename to shorthand form of model_parameters inside class for ease of use."""
        return self._mod_pars

    @model_parameters.setter
    def model_parameters(self, value):
        self._check_input_dataclasses("model_parameters", value)
        self._mod_pars = copy.copy(value)

    @property
    def mode(self):
        """Model mode property. Either 'linear' or 'instant_reaction'."""
        return self._mode

    @mode.setter
    def mode(self, value):
        match value:
            case "linear" | "linear decay" | "linear_decay" | 0:
                self._mode = "linear"
            case "instant" | "instant_reaction" | "instant reaction" | 1:
                if self._electron_acceptors is None or self._utilization_factor is None:
                    raise ValueError(
                        "Model mode was set to 'instant reaction', but electron acceptor parameters are "
                        "missing. Use the instant_reaction method to supply the electron acceptor "
                        "concentrations."
                    )
                self._mode = "instant_reaction"
            case _:
                warnings.warn(f"Mode '{value}' not recognized. Defaulting to 'linear' instead.", UserWarning)
                self._mode = "linear"

    @property
    def electron_acceptors(self):
        """Return dictionary of electron acceptor parameters."""
        return self._electron_acceptors.dictionary

    @property
    def utilization_factor(self):
        """Return dictionary of utilization factor property."""
        return self._utilization_factor.dictionary

    @property
    def relative_cxyt(self):
        """Compute relative concentration c(x,y,t)/c0, where c0 is the maximum source zone concentration at t=0."""
        maximum_concentration = np.max(self.source_parameters.source_zone_concentration)
        relative_cxyt = self.cxyt / maximum_concentration
        return relative_cxyt

    @property
    @abstractmethod
    def short_description(self):
        """Short string describing model type."""
        pass

    @abstractmethod
    def run(self):
        """Method that runs the model and ensures that initialisation is performed."""
        pass

    @abstractmethod
    def sample(self, x_position, y_position, t_position):
        """Method that calculates concentration at single, specified location in model domain."""
        pass

    @abstractmethod
    def _calculate_concentration_for_all_xyt(self) -> np.ndarray:
        """Method that calculates and return concentration array for all model x, y and t."""
        pass

    def _pre_run_initialization_parameters(self):
        """Parameter initialization for model."""
        # One-dimensional model domain arrays
        self.x = np.arange(0, self._mod_pars.model_length + self._mod_pars.dx, self._mod_pars.dx)
        self.y = self._calculate_y_discretization()
        self.t = np.arange(self._mod_pars.dt, self._mod_pars.model_time + self._mod_pars.dt, self._mod_pars.dt)

        # Three-dimensional model domain arrays
        self.xxx = self.x[None, None, :]
        self.yyy = self.y[None, :, None]
        self.ttt = self.t[:, None, None]

        if (
            self._att_pars.bulk_density is not None
            and self._att_pars.partition_coefficient is not None
            and self._att_pars.fraction_organic_carbon is not None
        ):
            self._att_pars.calculate_retardation(self._hyd_pars.porosity)

        self.rv = self._hyd_pars.velocity / self._att_pars.retardation

        # cxyt is concentration output array
        self.cxyt = np.zeros((len(self.t), len(self.y), len(self.x)))

        # Calculate retardation if not already specified in adsorption_parameters
        self.k_source = self._calculate_source_decay()
        self.y_source = self._src_pars.source_zone_boundary
        # Subtract outer source zones from inner source zones
        self.c_source = self._src_pars.source_zone_concentration.copy()
        self.c_source[:-1] = self.c_source[:-1] - self.c_source[1:]
        if self._mode == "instant_reaction":
            self.c_source[-1] += self.biodegradation_capacity
            self._decay_rate = 0
        else:
            self._decay_rate = self._att_pars.decay_rate

    def _calculate_source_decay(self):
        """Calculate source decay/depletion."""
        if self._src_pars.total_mass != np.inf:
            Q, c0_avg = calculate_discharge_and_average_source_zone_concentration(self)
            k_source = Q * c0_avg / self._src_pars.total_mass
        # If source mass is not a float, it is an infinite source, therefore, no source decay takes place.
        else:
            k_source = 0

        return k_source

    def _check_input_dataclasses(self, key, value):
        """Check if input parameters are the correct dataclasses. Raise an error if not."""
        dataclass_dict = {
            "hydrological_parameters": mibitrans.data.parameters.HydrologicalParameters,
            "attenuation_parameters": mibitrans.data.parameters.AttenuationParameters,
            "source_parameters": mibitrans.data.parameters.SourceParameters,
            "model_parameters": mibitrans.data.parameters.ModelParameters,
        }

        if not isinstance(value, dataclass_dict[key]):
            raise TypeError(f"Input argument {key} should be {dataclass_dict[key]}, but is {type(value)} instead.")

    def _calculate_y_discretization(self):
        """Calculate y-direction discretization."""
        if self._mod_pars.model_width >= 2 * self._src_pars.source_zone_boundary[-1]:
            y = np.arange(
                -self._mod_pars.model_width / 2, self._mod_pars.model_width / 2 + self._mod_pars.dy, self._mod_pars.dy
            )
        else:
            y = np.arange(
                -self._src_pars.source_zone_boundary[-1],
                self._src_pars.source_zone_boundary[-1] + self._mod_pars.dy,
                self._mod_pars.dy,
            )
            warnings.warn(
                "Source zone boundary is larger than model width. Model width adjusted to fit entire source zone."
            )
        return y

    def _calculate_biodegradation_capacity(self):
        """Determine biodegradation capacity based on electron acceptor concentrations and utilization factor."""
        biodegradation_capacity = 0
        for key, item in self._utilization_factor.dictionary.items():
            biodegradation_capacity += getattr(self._electron_acceptors, util_to_conc_name[key]) / item

        return biodegradation_capacity

    def instant_reaction(
        self,
        electron_acceptors: list | np.ndarray | dict | ElectronAcceptors,
        utilization_factor: list | np.ndarray | dict | UtilizationFactor = UtilizationFactor(
            util_oxygen=3.14, util_nitrate=4.9, util_ferrous_iron=21.8, util_sulfate=4.7, util_methane=0.78
        ),
    ):
        """Enable and set up parameters for instant reaction model.

        Instant reaction model assumes that biodegradation is an instantaneous process compared to the groundwater flow
        velocity. The biodegradation is assumed to be governed by the availability of electron acceptors, and quantified
        using  stoichiometric relations from the degradation reactions. Considered are concentrations of acceptors
        Oxygen, Nitrate and Sulfate, and reduced species Ferrous Iron and Methane.

        Args:
            electron_acceptors (ElectronAcceptors): ElectronAcceptor dataclass containing electron acceptor
                concentrations. Alternatively provided as list, numpy array or dictionary corresponding with
                delta_oxygen, delta_nitrate, ferrous_iron, delta_sulfate and methane. For more information, see
                documentation for ElectronAcceptors.
            utilization_factor (UtilizationFactor, optional): UtilizationFactor dataclass containing electron acceptor
                utilization factors. Alternatively provided as list, numpy array or dictionary corresponding with
                information, see documentation of UtilizationFactor. By default, electron acceptor utilization factors
                for a BTEX mixture are used, based on values by Wiedemeier et al. (1995).
        """
        self._electron_acceptors, self._utilization_factor = _check_instant_reaction_acceptor_input(
            electron_acceptors, utilization_factor
        )
        self._mode = "instant_reaction"
        self.biodegradation_capacity = self._calculate_biodegradation_capacity()
        self.cxyt_noBC = 0
        self._pre_run_initialization_parameters()

    def _check_model_mode_before_run(self):
        self._pre_run_initialization_parameters()
        if self._mode == "linear":
            if self.biodegradation_capacity is not None:
                warnings.warn(
                    "Instant reaction parameters are present while model mode is linear. "
                    "Make sure that this is indeed the desired model."
                )
        if self._mode == "instant_reaction":
            if self.biodegradation_capacity is None:
                raise ValueError(
                    "Instant reaction parameters are not present. "
                    "Please provide them with the 'instant_reaction' class method."
                )

attenuation_parameters property writable

Rename to shorthand form of attenuation_parameters inside class for ease of use.

electron_acceptors property

Return dictionary of electron acceptor parameters.

hydrological_parameters property writable

Rename to shorthand form of hydrological_parameters inside class for ease of use.

input_parameters property

Return the input arguments for the model in the form of a dictionary, based on current values.

mode property writable

Model mode property. Either ‘linear’ or ‘instant_reaction’.

model_parameters property writable

Rename to shorthand form of model_parameters inside class for ease of use.

relative_cxyt property

Compute relative concentration c(x,y,t)/c0, where c0 is the maximum source zone concentration at t=0.

short_description abstractmethod property

Short string describing model type.

source_parameters property writable

Rename to shorthand form of source_parameters inside class for ease of use.

utilization_factor property

Return dictionary of utilization factor property.

__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose=False)

Initialize parent class object.

Parameters:

Name Type Description Default
hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters)

Dataclass object containing hydrological parameters from HydrologicalParameters.

required
attenuation_parameters (mibitrans.data.read.AttenuationParameters)

Dataclass object containing adsorption, degradation and diffusion parameters from AttenuationParameters.

required
source_parameters (mibitrans.data.read.SourceParameters)

Dataclass object containing source parameters from SourceParameters.

required
model_parameters (mibitrans.data.read.ModelParameters)

Dataclass object containing model parameters from ModelParameters.

required
verbose bool

Verbose mode. Defaults to False.

False
Source code in mibitrans/transport/model_parent.py
def __init__(
    self, hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose=False
):
    """Initialize parent class object.

    Args:
        hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters) : Dataclass object containing
            hydrological parameters from HydrologicalParameters.
        attenuation_parameters (mibitrans.data.read.AttenuationParameters) : Dataclass object containing adsorption,
            degradation and diffusion parameters from AttenuationParameters.
        source_parameters (mibitrans.data.read.SourceParameters) : Dataclass object containing source parameters
            from SourceParameters.
        model_parameters (mibitrans.data.read.ModelParameters) : Dataclass object containing model parameters from
            ModelParameters.
        verbose (bool, optional): Verbose mode. Defaults to False.
    """
    # Check if input arguments are of the correct dataclass
    for key, value in locals().items():
        if key not in ["self", "verbose"]:
            self._check_input_dataclasses(key, value)

    self._hyd_pars = copy.copy(hydrological_parameters)
    self._att_pars = copy.copy(attenuation_parameters)
    self._src_pars = copy.copy(source_parameters)
    self._mod_pars = copy.copy(model_parameters)
    self._decay_rate = self._att_pars.decay_rate

    self.verbose = verbose
    self._mode = "linear"
    self._electron_acceptors = None
    self._utilization_factor = None
    self.biodegradation_capacity = None
    self.cxyt_noBC = None
    self._pre_run_initialization_parameters()

instant_reaction(electron_acceptors, utilization_factor=UtilizationFactor(util_oxygen=3.14, util_nitrate=4.9, util_ferrous_iron=21.8, util_sulfate=4.7, util_methane=0.78))

Enable and set up parameters for instant reaction model.

Instant reaction model assumes that biodegradation is an instantaneous process compared to the groundwater flow velocity. The biodegradation is assumed to be governed by the availability of electron acceptors, and quantified using stoichiometric relations from the degradation reactions. Considered are concentrations of acceptors Oxygen, Nitrate and Sulfate, and reduced species Ferrous Iron and Methane.

Parameters:

Name Type Description Default
electron_acceptors ElectronAcceptors

ElectronAcceptor dataclass containing electron acceptor concentrations. Alternatively provided as list, numpy array or dictionary corresponding with delta_oxygen, delta_nitrate, ferrous_iron, delta_sulfate and methane. For more information, see documentation for ElectronAcceptors.

required
utilization_factor UtilizationFactor

UtilizationFactor dataclass containing electron acceptor utilization factors. Alternatively provided as list, numpy array or dictionary corresponding with information, see documentation of UtilizationFactor. By default, electron acceptor utilization factors for a BTEX mixture are used, based on values by Wiedemeier et al. (1995).

UtilizationFactor(util_oxygen=3.14, util_nitrate=4.9, util_ferrous_iron=21.8, util_sulfate=4.7, util_methane=0.78)
Source code in mibitrans/transport/model_parent.py
def instant_reaction(
    self,
    electron_acceptors: list | np.ndarray | dict | ElectronAcceptors,
    utilization_factor: list | np.ndarray | dict | UtilizationFactor = UtilizationFactor(
        util_oxygen=3.14, util_nitrate=4.9, util_ferrous_iron=21.8, util_sulfate=4.7, util_methane=0.78
    ),
):
    """Enable and set up parameters for instant reaction model.

    Instant reaction model assumes that biodegradation is an instantaneous process compared to the groundwater flow
    velocity. The biodegradation is assumed to be governed by the availability of electron acceptors, and quantified
    using  stoichiometric relations from the degradation reactions. Considered are concentrations of acceptors
    Oxygen, Nitrate and Sulfate, and reduced species Ferrous Iron and Methane.

    Args:
        electron_acceptors (ElectronAcceptors): ElectronAcceptor dataclass containing electron acceptor
            concentrations. Alternatively provided as list, numpy array or dictionary corresponding with
            delta_oxygen, delta_nitrate, ferrous_iron, delta_sulfate and methane. For more information, see
            documentation for ElectronAcceptors.
        utilization_factor (UtilizationFactor, optional): UtilizationFactor dataclass containing electron acceptor
            utilization factors. Alternatively provided as list, numpy array or dictionary corresponding with
            information, see documentation of UtilizationFactor. By default, electron acceptor utilization factors
            for a BTEX mixture are used, based on values by Wiedemeier et al. (1995).
    """
    self._electron_acceptors, self._utilization_factor = _check_instant_reaction_acceptor_input(
        electron_acceptors, utilization_factor
    )
    self._mode = "instant_reaction"
    self.biodegradation_capacity = self._calculate_biodegradation_capacity()
    self.cxyt_noBC = 0
    self._pre_run_initialization_parameters()

run() abstractmethod

Method that runs the model and ensures that initialisation is performed.

Source code in mibitrans/transport/model_parent.py
@abstractmethod
def run(self):
    """Method that runs the model and ensures that initialisation is performed."""
    pass

sample(x_position, y_position, t_position) abstractmethod

Method that calculates concentration at single, specified location in model domain.

Source code in mibitrans/transport/model_parent.py
@abstractmethod
def sample(self, x_position, y_position, t_position):
    """Method that calculates concentration at single, specified location in model domain."""
    pass

models

Anatrans

Bases: Transport3D

Model class using an analytical solution based on Bear (1979), Domenico (1987) & Newell et al. (1997).

Under the assumption that C(x,y,z,t) = C(x,t) * C(y,t) * C(z,t), the 3D ADE can be broken up in three separate differential equations which can be solved individually. For C(x,t) the solution is given in Bear (1979), C(y,t) and C(z,t) can be derived from Crank (1975). The equation used for Anatrans is the combination of these solutions, with addition of source depletion, source superposition and instant reaction model, described in Newell et al. (1997) and implemented in the BIOSCREEN screening model. The solution of Newell et al. (1997) is based of the Domenico (1987) solution, a truncated version of the equation described above, which introduces an error with a size dependent on the ratio of flow velocity and longitudinal dispersivity. Anatrans instead uses the fully untruncated version.

Bear, J. 1979. Hydraulics of Ground Water. New York: McGraw-Hill.

Crank, J. 1975. The mathematics of Diffusion. New York: Oxford University Press.

Domenico, P. A. (1987). An analytical model for multidimensional transport of a decaying contaminant species. Journal of Hydrology, 91(1-2), 49-58.

Newell, C. J., McLeod, R. K., & Gonzales, J. R. (1997). BIOSCREEN natural attenuation decision support system version 1.4 revisions, Tech. rep., U.S. EPA.

Source code in mibitrans/transport/models.py
class Anatrans(Transport3D):
    """Model class using an analytical solution based on Bear (1979), Domenico (1987) & Newell et al. (1997).

    Under the assumption that C(x,y,z,t) = C(x,t) * C(y,t) * C(z,t), the 3D ADE can be broken up in three separate
    differential equations which can be solved individually. For C(x,t) the solution is given in Bear (1979), C(y,t) and
    C(z,t) can be derived from Crank (1975). The equation used for Anatrans is the combination of these solutions, with
    addition of source depletion, source superposition and instant reaction model, described in Newell et al. (1997) and
    implemented in the BIOSCREEN screening model. The solution of Newell et al. (1997) is based of the Domenico (1987)
    solution, a truncated version of the equation described above, which introduces an error with a size dependent on
    the ratio of flow velocity and longitudinal dispersivity. Anatrans instead uses the fully untruncated version.

    Bear, J. 1979. Hydraulics of Ground Water. New York: McGraw-Hill.

    Crank, J. 1975. The mathematics of Diffusion. New York: Oxford University Press.

    Domenico, P. A. (1987). An analytical model for multidimensional transport of a decaying contaminant species.
    Journal of Hydrology, 91(1-2), 49-58.

    Newell, C. J., McLeod, R. K., & Gonzales, J. R. (1997). BIOSCREEN natural attenuation decision support
    system version 1.4 revisions, Tech. rep., U.S. EPA.
    """

    def __init__(
        self,
        hydrological_parameters,
        attenuation_parameters,
        source_parameters,
        model_parameters,
        verbose=False,
    ):
        """Initialize model object.

        Args:
            hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters) : Dataclass object containing
                hydrological parameters from HydrologicalParameters.
            attenuation_parameters (mibitrans.data.read.AttenuationParameters) : Dataclass object containing adsorption,
                degradation and diffusion parameters from AttenuationParameters.
            source_parameters (mibitrans.data.read.SourceParameters) : Dataclass object containing source parameters
                from SourceParameters.
            model_parameters (mibitrans.data.read.ModelParameters) : Dataclass object containing model parameters from
                ModelParameters.
            verbose (bool, optional): Verbose mode. Defaults to False.

        Attributes:
            mode (str) : Current model mode. Is 'linear' by default. Once instant reaction parameters are provided. Use
                this attribute to switch between 'linear' and 'instant_reaction' models.
            cxyt (np.ndarray) : Output array containing concentrations in model domain, in [g/m^3]. Indexed as [t,y,x]
            relative_cxyt (np.ndarray) : Output array with concentrations in model domain, divided by the maximum source
                zone concentration at t=0. Indexed as [t,y,x].
            x (np.ndarray) : Discretized model x-dimension, in [m].
            y (np.ndarray) : Discretized model y-dimension, in [y].
            t (np.ndarray) : Discretized model t-dimension, in [days].
            c_source (np.ndarray) : Nett source zone concentrations, accounting for source superposition, in [g/m^3].
            vr (float) : Retarded groundwater flow velocity, in [m/d].
            k_source (float) : Source zone decay rate, in [1/days].
            biodegradation_capacity (float) : Maximum capacity of biodegradation given provided electron acceptor
            concentrations, in [g/m^3].

        Methods:
            run : Run model with current parameters, returns Results object.
            sample : Calculate concentration at any given position and point in time.
            instant_reaction : Activate the instant reaction model by providing electron acceptor concentrations. And
                optionally electron acceptor utilization factors. Switch between model modes by using the mode
                attribute.
            centerline : Plot center of contaminant plume of this model, at a specified time and y position.
            transverse : Plot concentration distribution as a line horizontal transverse to the plume extent.
            breakthrough : Plot contaminant breakthrough curve at given x and y position in model domain.
            plume_2d : Plot contaminant plume as a 2D colormesh, at a specified time.
            plume_3d : Plot contaminant plume as a 3D surface, at a specified time.

        Raises:
            TypeError : If input is not of the correct Dataclass.

        """
        super().__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose)
        if self._att_pars.diffusion != 0:
            warnings.warn("Domenico model does not consider molecular diffusion.", UserWarning)

    @property
    def short_description(self):
        """Short description of model type."""
        return "Anatrans model"

    def run(self):
        """Calculate the concentration for all discretized x, y and t using the analytical transport model."""
        self._check_model_mode_before_run()
        with np.errstate(divide="ignore", invalid="ignore"):
            self.cxyt = self._calculate_concentration_for_all_xyt(self.xxx, self.yyy, self.ttt)
        return Results(self)

    def sample(self, x_position, y_position, time):
        """Give concentration at any given position and point in time.

        Args:
            x_position (float): x position in domain extent [m].
            y_position (float): y position in domain extent [m].
            time (float): time for which concentration is sampled [days].

        Returns:
            concentration (float): concentration at given position and point in time [g/m^3].

        """
        for par, value in locals().items():
            if par != "self":
                validate_input_values(par, value)

        self._pre_run_initialization_parameters()

        if self.mode == "instant_reaction":
            save_c_noBC = self.cxyt_noBC.copy()
        x = np.array([x_position])
        y = np.array([y_position])
        t = np.array([time])
        concentration = self._calculate_concentration_for_all_xyt(x, y, t)[0]
        if self.mode == "instant_reaction":
            self.cxyt_noBC = save_c_noBC
        return concentration

    def _equation_term_x(self, xxx, ttt, decay_sqrt):
        return np.exp(xxx * (1 - decay_sqrt) / (self._hyd_pars.alpha_x * 2)) * erfc(
            (xxx - self.rv * ttt * decay_sqrt) / (2 * np.sqrt(self._hyd_pars.alpha_x * self.rv * ttt))
        )

    def _equation_term_additional_x(self, xxx, ttt, decay_sqrt):
        erfc_inner = (xxx + decay_sqrt * self.rv * ttt) / (2 * np.sqrt(self._hyd_pars.alpha_x * self.rv * ttt))
        # Additional term is prone to overflow of exp and underflow of erfc under certain parameter combinations.
        # To decrease cases, used erfcx. Where erfcx(a) = exp(a**2)*erfc(a) -> exp(b)*erfc(a) = exp(b - a**2) * erfcx(a)
        term = np.exp(xxx * (1 + decay_sqrt) * (1 / 2) / self._hyd_pars.alpha_x - erfc_inner**2) * erfcx(erfc_inner)
        return term

    def _equation_term_z(self, xxx):
        inner_term = self._src_pars.depth / (2 * np.sqrt(self._hyd_pars.alpha_z * xxx))
        return erf(inner_term) - erf(-inner_term)

    def _equation_term_source_decay(self, xxx, ttt):
        term = np.exp(-self.k_source * (ttt - xxx / self.rv))
        # Term can be max 1; can not have 'generation' of solute ahead of advection.
        return np.where(term > 1, 1, term)

    def _equation_term_y(self, i, xxx, yyy):
        div_term = 2 * np.sqrt(self._hyd_pars.alpha_y * xxx)
        term = erf((yyy + self.y_source[i]) / div_term) - erf((yyy - self.y_source[i]) / div_term)
        term[np.isnan(term)] = 0
        return term

    def _calculate_concentration_for_all_xyt(self, xxx, yyy, ttt):
        cxyt = 0
        decay_sqrt = np.sqrt(1 + 4 * self._decay_rate * self._hyd_pars.alpha_x / self.rv)
        x_term = self._equation_term_x(xxx, ttt, decay_sqrt)
        additional_x = self._equation_term_additional_x(xxx, ttt, decay_sqrt)
        z_term = self._equation_term_z(xxx)
        source_decay = self._equation_term_source_decay(xxx, ttt)
        for i in range(len(self.c_source)):
            y_term = self._equation_term_y(i, xxx, yyy)
            cxyt_step = 1 / 8 * self.c_source[i] * source_decay * (x_term + additional_x) * y_term * z_term
            cxyt += cxyt_step
        if self._mode == "instant_reaction":
            self.cxyt_noBC = cxyt.copy()
            cxyt -= self.biodegradation_capacity
            cxyt = np.where(cxyt < 0, 0, cxyt)
        self.has_run = True
        return cxyt

short_description property

Short description of model type.

__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose=False)

Initialize model object.

Parameters:

Name Type Description Default
hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters)

Dataclass object containing hydrological parameters from HydrologicalParameters.

required
attenuation_parameters (mibitrans.data.read.AttenuationParameters)

Dataclass object containing adsorption, degradation and diffusion parameters from AttenuationParameters.

required
source_parameters (mibitrans.data.read.SourceParameters)

Dataclass object containing source parameters from SourceParameters.

required
model_parameters (mibitrans.data.read.ModelParameters)

Dataclass object containing model parameters from ModelParameters.

required
verbose bool

Verbose mode. Defaults to False.

False

Attributes:

Name Type Description
mode str)

Current model mode. Is ‘linear’ by default. Once instant reaction parameters are provided. Use this attribute to switch between ‘linear’ and ‘instant_reaction’ models.

cxyt np.ndarray)

Output array containing concentrations in model domain, in [g/m^3]. Indexed as [t,y,x]

relative_cxyt np.ndarray)

Output array with concentrations in model domain, divided by the maximum source zone concentration at t=0. Indexed as [t,y,x].

x np.ndarray)

Discretized model x-dimension, in [m].

y np.ndarray)

Discretized model y-dimension, in [y].

t np.ndarray)

Discretized model t-dimension, in [days].

c_source np.ndarray)

Nett source zone concentrations, accounting for source superposition, in [g/m^3].

vr float)

Retarded groundwater flow velocity, in [m/d].

k_source float)

Source zone decay rate, in [1/days].

biodegradation_capacity float)

Maximum capacity of biodegradation given provided electron acceptor

Functions:

Name Description
run

Run model with current parameters, returns Results object.

sample

Calculate concentration at any given position and point in time.

instant_reaction

Activate the instant reaction model by providing electron acceptor concentrations. And optionally electron acceptor utilization factors. Switch between model modes by using the mode attribute.

centerline

Plot center of contaminant plume of this model, at a specified time and y position.

transverse

Plot concentration distribution as a line horizontal transverse to the plume extent.

breakthrough

Plot contaminant breakthrough curve at given x and y position in model domain.

plume_2d

Plot contaminant plume as a 2D colormesh, at a specified time.

plume_3d

Plot contaminant plume as a 3D surface, at a specified time.

Raises:

Type Description
TypeError

If input is not of the correct Dataclass.

Source code in mibitrans/transport/models.py
def __init__(
    self,
    hydrological_parameters,
    attenuation_parameters,
    source_parameters,
    model_parameters,
    verbose=False,
):
    """Initialize model object.

    Args:
        hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters) : Dataclass object containing
            hydrological parameters from HydrologicalParameters.
        attenuation_parameters (mibitrans.data.read.AttenuationParameters) : Dataclass object containing adsorption,
            degradation and diffusion parameters from AttenuationParameters.
        source_parameters (mibitrans.data.read.SourceParameters) : Dataclass object containing source parameters
            from SourceParameters.
        model_parameters (mibitrans.data.read.ModelParameters) : Dataclass object containing model parameters from
            ModelParameters.
        verbose (bool, optional): Verbose mode. Defaults to False.

    Attributes:
        mode (str) : Current model mode. Is 'linear' by default. Once instant reaction parameters are provided. Use
            this attribute to switch between 'linear' and 'instant_reaction' models.
        cxyt (np.ndarray) : Output array containing concentrations in model domain, in [g/m^3]. Indexed as [t,y,x]
        relative_cxyt (np.ndarray) : Output array with concentrations in model domain, divided by the maximum source
            zone concentration at t=0. Indexed as [t,y,x].
        x (np.ndarray) : Discretized model x-dimension, in [m].
        y (np.ndarray) : Discretized model y-dimension, in [y].
        t (np.ndarray) : Discretized model t-dimension, in [days].
        c_source (np.ndarray) : Nett source zone concentrations, accounting for source superposition, in [g/m^3].
        vr (float) : Retarded groundwater flow velocity, in [m/d].
        k_source (float) : Source zone decay rate, in [1/days].
        biodegradation_capacity (float) : Maximum capacity of biodegradation given provided electron acceptor
        concentrations, in [g/m^3].

    Methods:
        run : Run model with current parameters, returns Results object.
        sample : Calculate concentration at any given position and point in time.
        instant_reaction : Activate the instant reaction model by providing electron acceptor concentrations. And
            optionally electron acceptor utilization factors. Switch between model modes by using the mode
            attribute.
        centerline : Plot center of contaminant plume of this model, at a specified time and y position.
        transverse : Plot concentration distribution as a line horizontal transverse to the plume extent.
        breakthrough : Plot contaminant breakthrough curve at given x and y position in model domain.
        plume_2d : Plot contaminant plume as a 2D colormesh, at a specified time.
        plume_3d : Plot contaminant plume as a 3D surface, at a specified time.

    Raises:
        TypeError : If input is not of the correct Dataclass.

    """
    super().__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose)
    if self._att_pars.diffusion != 0:
        warnings.warn("Domenico model does not consider molecular diffusion.", UserWarning)

run()

Calculate the concentration for all discretized x, y and t using the analytical transport model.

Source code in mibitrans/transport/models.py
def run(self):
    """Calculate the concentration for all discretized x, y and t using the analytical transport model."""
    self._check_model_mode_before_run()
    with np.errstate(divide="ignore", invalid="ignore"):
        self.cxyt = self._calculate_concentration_for_all_xyt(self.xxx, self.yyy, self.ttt)
    return Results(self)

sample(x_position, y_position, time)

Give concentration at any given position and point in time.

Parameters:

Name Type Description Default
x_position float

x position in domain extent [m].

required
y_position float

y position in domain extent [m].

required
time float

time for which concentration is sampled [days].

required

Returns:

Name Type Description
concentration float

concentration at given position and point in time [g/m^3].

Source code in mibitrans/transport/models.py
def sample(self, x_position, y_position, time):
    """Give concentration at any given position and point in time.

    Args:
        x_position (float): x position in domain extent [m].
        y_position (float): y position in domain extent [m].
        time (float): time for which concentration is sampled [days].

    Returns:
        concentration (float): concentration at given position and point in time [g/m^3].

    """
    for par, value in locals().items():
        if par != "self":
            validate_input_values(par, value)

    self._pre_run_initialization_parameters()

    if self.mode == "instant_reaction":
        save_c_noBC = self.cxyt_noBC.copy()
    x = np.array([x_position])
    y = np.array([y_position])
    t = np.array([time])
    concentration = self._calculate_concentration_for_all_xyt(x, y, t)[0]
    if self.mode == "instant_reaction":
        self.cxyt_noBC = save_c_noBC
    return concentration

Bioscreen

Bases: Anatrans

Model class using the analytical solution implemented in the BIOSCREEN screening model, Newell et al. (1997).

This model is an exact implementation of the transport equations implemented in the BIOSCREEN screening model of Newell et al. (1997), which is based on the Domenico (1987) analytical model. Using a truncated version of the equation used in the Anatrans model. This model is implemented as a method of comparison with the original BIOSCREEN software. And is included for legacy reasons, since it is the first model implemented in the mibitrans package, serving as a basis for the other models. However, caution should be taken when using this model, since a varying error is introduced by using the truncated analytical solution. The error is most prominent for shorter times and distances from the source, and depends on the ratio of flow velocity and longitudinal dispersivity. For modelling, the Anatrans (untruncated approximate solution) and Mibitrans (exact analytical solution) models are recommended instead.

Domenico, P. A. (1987). An analytical model for multidimensional transport of a decaying contaminant species. Journal of Hydrology, 91(1-2), 49-58.

Newell, C. J., McLeod, R. K., & Gonzales, J. R. (1997). BIOSCREEN natural attenuation decision support system version 1.4 revisions, Tech. rep., U.S. EPA.

Source code in mibitrans/transport/models.py
class Bioscreen(Anatrans):
    """Model class using the analytical solution implemented in the BIOSCREEN screening model, Newell et al. (1997).

    This model is an exact implementation of the transport equations implemented in the BIOSCREEN screening model of
    Newell et al. (1997), which is based on the Domenico (1987) analytical model. Using a truncated version of the
    equation used in the Anatrans model. This model is implemented as a method of comparison with the original BIOSCREEN
    software. And is included for legacy reasons, since it is the first model implemented in the mibitrans package,
    serving as a basis for the other models. However, caution should be taken when using this model, since a varying
    error is introduced by using the truncated analytical solution. The error is most prominent for shorter times and
    distances from the source, and depends on the ratio of flow velocity and longitudinal dispersivity. For modelling,
    the Anatrans (untruncated approximate solution) and Mibitrans (exact analytical solution) models are recommended
    instead.

    Domenico, P. A. (1987). An analytical model for multidimensional transport of a decaying contaminant species.
    Journal of Hydrology, 91(1-2), 49-58.

    Newell, C. J., McLeod, R. K., & Gonzales, J. R. (1997). BIOSCREEN natural attenuation decision support
    system version 1.4 revisions, Tech. rep., U.S. EPA.
    """

    def __init__(
        self,
        hydrological_parameters,
        attenuation_parameters,
        source_parameters,
        model_parameters,
        verbose=False,
    ):
        """Initialize model object.

        Args:
            hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters) : Dataclass object containing
                hydrological parameters from HydrologicalParameters.
            attenuation_parameters (mibitrans.data.read.AttenuationParameters) : Dataclass object containing adsorption,
                degradation and diffusion parameters from AttenuationParameters.
            source_parameters (mibitrans.data.read.SourceParameters) : Dataclass object containing source parameters
                from SourceParameters.
            model_parameters (mibitrans.data.read.ModelParameters) : Dataclass object containing model parameters from
                ModelParameters.
            verbose (bool, optional): Verbose mode. Defaults to False.

        Attributes:
            mode (str) : Current model mode. Is 'linear' by default. Once instant reaction parameters are provided. Use
                this attribute to switch between 'linear' and 'instant_reaction' models.
            cxyt (np.ndarray) : Output array containing concentrations in model domain, in [g/m^3]. Indexed as [t,y,x]
            relative_cxyt (np.ndarray) : Output array with concentrations in model domain, divided by the maximum source
                zone concentration at t=0. Indexed as [t,y,x].
            x (np.ndarray) : Discretized model x-dimension, in [m].
            y (np.ndarray) : Discretized model y-dimension, in [y].
            t (np.ndarray) : Discretized model t-dimension, in [days].
            c_source (np.ndarray) : Nett source zone concentrations, accounting for source superposition, in [g/m^3].
            vr (float) : Retarded groundwater flow velocity, in [m/d].
            k_source (float) : Source zone decay rate, in [1/days].
            biodegradation_capacity (float) : Maximum capacity of biodegradation given provided electron acceptor
            concentrations, in [g/m^3].

        Methods:
            run : Run model with current parameters, returns Results object.
            sample : Calculate concentration at any given position and point in time.
            instant_reaction : Activate the instant reaction model by providing electron acceptor concentrations. And
                optionally electron acceptor utilization factors. Switch between model modes by using the mode
                attribute.

        Raises:
            TypeError : If input is not of the correct Dataclass.

        """
        super().__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose)

    @property
    def short_description(self):
        """Short description of model type."""
        return "Bioscreen model"

    def _calculate_concentration_for_all_xyt(self, xxx, yyy, ttt):
        # Difference with the Anatrans solution is the lack of additional term.
        cxyt = 0
        with np.errstate(divide="ignore", invalid="ignore"):
            decay_sqrt = np.sqrt(1 + 4 * self._decay_rate * self._hyd_pars.alpha_x / self.rv)
            x_term = self._equation_term_x(xxx, ttt, decay_sqrt)
            z_term = self._equation_term_z(xxx)
            source_decay = self._equation_term_source_decay(xxx, ttt)
            for i in range(len(self.c_source)):
                y_term = self._equation_term_y(i, xxx, yyy)
                cxyt_step = 1 / 8 * self.c_source[i] * source_decay * x_term * y_term * z_term
                cxyt += cxyt_step
        if self._mode == "instant_reaction":
            self.cxyt_noBC = cxyt.copy()
            cxyt -= self.biodegradation_capacity
            cxyt = np.where(cxyt < 0, 0, cxyt)
        self.has_run = True
        return cxyt

short_description property

Short description of model type.

__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose=False)

Initialize model object.

Parameters:

Name Type Description Default
hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters)

Dataclass object containing hydrological parameters from HydrologicalParameters.

required
attenuation_parameters (mibitrans.data.read.AttenuationParameters)

Dataclass object containing adsorption, degradation and diffusion parameters from AttenuationParameters.

required
source_parameters (mibitrans.data.read.SourceParameters)

Dataclass object containing source parameters from SourceParameters.

required
model_parameters (mibitrans.data.read.ModelParameters)

Dataclass object containing model parameters from ModelParameters.

required
verbose bool

Verbose mode. Defaults to False.

False

Attributes:

Name Type Description
mode str)

Current model mode. Is ‘linear’ by default. Once instant reaction parameters are provided. Use this attribute to switch between ‘linear’ and ‘instant_reaction’ models.

cxyt np.ndarray)

Output array containing concentrations in model domain, in [g/m^3]. Indexed as [t,y,x]

relative_cxyt np.ndarray)

Output array with concentrations in model domain, divided by the maximum source zone concentration at t=0. Indexed as [t,y,x].

x np.ndarray)

Discretized model x-dimension, in [m].

y np.ndarray)

Discretized model y-dimension, in [y].

t np.ndarray)

Discretized model t-dimension, in [days].

c_source np.ndarray)

Nett source zone concentrations, accounting for source superposition, in [g/m^3].

vr float)

Retarded groundwater flow velocity, in [m/d].

k_source float)

Source zone decay rate, in [1/days].

biodegradation_capacity float)

Maximum capacity of biodegradation given provided electron acceptor

Functions:

Name Description
run

Run model with current parameters, returns Results object.

sample

Calculate concentration at any given position and point in time.

instant_reaction

Activate the instant reaction model by providing electron acceptor concentrations. And optionally electron acceptor utilization factors. Switch between model modes by using the mode attribute.

Raises:

Type Description
TypeError

If input is not of the correct Dataclass.

Source code in mibitrans/transport/models.py
def __init__(
    self,
    hydrological_parameters,
    attenuation_parameters,
    source_parameters,
    model_parameters,
    verbose=False,
):
    """Initialize model object.

    Args:
        hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters) : Dataclass object containing
            hydrological parameters from HydrologicalParameters.
        attenuation_parameters (mibitrans.data.read.AttenuationParameters) : Dataclass object containing adsorption,
            degradation and diffusion parameters from AttenuationParameters.
        source_parameters (mibitrans.data.read.SourceParameters) : Dataclass object containing source parameters
            from SourceParameters.
        model_parameters (mibitrans.data.read.ModelParameters) : Dataclass object containing model parameters from
            ModelParameters.
        verbose (bool, optional): Verbose mode. Defaults to False.

    Attributes:
        mode (str) : Current model mode. Is 'linear' by default. Once instant reaction parameters are provided. Use
            this attribute to switch between 'linear' and 'instant_reaction' models.
        cxyt (np.ndarray) : Output array containing concentrations in model domain, in [g/m^3]. Indexed as [t,y,x]
        relative_cxyt (np.ndarray) : Output array with concentrations in model domain, divided by the maximum source
            zone concentration at t=0. Indexed as [t,y,x].
        x (np.ndarray) : Discretized model x-dimension, in [m].
        y (np.ndarray) : Discretized model y-dimension, in [y].
        t (np.ndarray) : Discretized model t-dimension, in [days].
        c_source (np.ndarray) : Nett source zone concentrations, accounting for source superposition, in [g/m^3].
        vr (float) : Retarded groundwater flow velocity, in [m/d].
        k_source (float) : Source zone decay rate, in [1/days].
        biodegradation_capacity (float) : Maximum capacity of biodegradation given provided electron acceptor
        concentrations, in [g/m^3].

    Methods:
        run : Run model with current parameters, returns Results object.
        sample : Calculate concentration at any given position and point in time.
        instant_reaction : Activate the instant reaction model by providing electron acceptor concentrations. And
            optionally electron acceptor utilization factors. Switch between model modes by using the mode
            attribute.

    Raises:
        TypeError : If input is not of the correct Dataclass.

    """
    super().__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose)

Mibitrans

Bases: Transport3D

Model class using an exact analytical solution as described in Karanovic (2007), based on Wexler (1992).

Karanovic (2007) implemented the Wexler (1992) exact analytical solution in the Excel based BIOSCREEN-AT, and added source depletion, akin to that implemented in its predecessor BIOSCREEN by Newell et al. (1997). The Mibitrans model allows for the same method as used in BIOSCREEN-AT, but expands it by allowing multiple source zones (by means of superposition) and including the instant reaction model. These were present in the original BIOSCREEN, but not reimplemented in BIOSCREEN-AT. Using a single source zone in this model, and not using the instant reaction option will make the Mibitrans solution resolve to the equation described in Karanovic (2007). Which in turn resolves to the Wexler (1992) solution if source depletion is disabled.

Karanovic, M., Neville, C. J., & Andrews, C. B. (2007). BIOSCREEN‐AT: BIOSCREEN with an exact analytical solution. Groundwater, 45(2), 242-245.

Newell, C. J., McLeod, R. K., & Gonzales, J. R. (1997). BIOSCREEN natural attenuation decision support system version 1.4 revisions, Tech. rep., U.S. EPA.

Wexler, E. J. (1992). Analytical solutions for one-, two-, and three-dimensional solute transport in ground-water systems with uniform flow. US Government Printing Office.

Source code in mibitrans/transport/models.py
class Mibitrans(Transport3D):
    """Model class using an exact analytical solution as described in Karanovic (2007), based on Wexler (1992).

    Karanovic (2007) implemented the Wexler (1992) exact analytical solution in the Excel based BIOSCREEN-AT, and added
    source depletion, akin to that implemented in its predecessor BIOSCREEN by Newell et al. (1997). The Mibitrans model
    allows for the same method as used in BIOSCREEN-AT, but expands it by allowing multiple source zones (by means of
    superposition) and including the instant reaction model. These were present in the original BIOSCREEN, but not
    reimplemented in BIOSCREEN-AT. Using a single source zone in this model, and not using the instant reaction option
    will make the Mibitrans solution resolve to the equation described in Karanovic (2007). Which in turn resolves to
    the Wexler (1992) solution if source depletion is disabled.

    Karanovic, M., Neville, C. J., & Andrews, C. B. (2007). BIOSCREEN‐AT: BIOSCREEN with an exact analytical solution.
    Groundwater, 45(2), 242-245.

    Newell, C. J., McLeod, R. K., & Gonzales, J. R. (1997). BIOSCREEN natural attenuation decision support
    system version 1.4 revisions, Tech. rep., U.S. EPA.

    Wexler, E. J. (1992). Analytical solutions for one-, two-, and three-dimensional solute transport in ground-water
    systems with uniform flow. US Government Printing Office.
    """

    def __init__(
        self,
        hydrological_parameters,
        attenuation_parameters,
        source_parameters,
        model_parameters,
        verbose=False,
    ):
        """Initialize model object.

        Args:
            hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters) : Dataclass object containing
                hydrological parameters from HydrologicalParameters.
            attenuation_parameters (mibitrans.data.read.AttenuationParameters) : Dataclass object containing adsorption,
                degradation and diffusion parameters from AttenuationParameters.
            source_parameters (mibitrans.data.read.SourceParameters) : Dataclass object containing source parameters
                from SourceParameters.
            model_parameters (mibitrans.data.read.ModelParameters) : Dataclass object containing model parameters from
                ModelParameters.
            verbose (bool, optional): Verbose mode. Defaults to False.

        Attributes:
            mode (str) : Current model mode. Is 'linear' by default. Once instant reaction parameters are provided. Use
                this attribute to switch between 'linear' and 'instant_reaction' models.
            cxyt (np.ndarray) : Output array containing concentrations in model domain, in [g/m^3]. Indexed as [t,y,x]
            relative_cxyt (np.ndarray) : Output array with concentrations in model domain, divided by the maximum source
                zone concentration at t=0. Indexed as [t,y,x].
            x (np.ndarray) : Discretized model x-dimension, in [m].
            y (np.ndarray) : Discretized model y-dimension, in [y].
            t (np.ndarray) : Discretized model t-dimension, in [days].
            c_source (np.ndarray) : Nett source zone concentrations, accounting for source superposition, in [g/m^3].
            vr (float) : Retarded groundwater flow velocity, in [m/d].
            k_source (float) : Source zone decay rate, in [1/days].
            biodegradation_capacity (float) : Maximum capacity of biodegradation given provided electron acceptor
            concentrations, in [g/m^3].

        Methods:
            run : Run model with current parameters, returns Results object.
            sample : Calculate concentration at any given position and point in time.
            instant_reaction : Activate the instant reaction model by providing electron acceptor concentrations. And
                optionally electron acceptor utilization factors. Switch between model modes by using the mode
                attribute.

        Raises:
            TypeError : If input is not of the correct Dataclass.

        """
        super().__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose)

    @property
    def short_description(self):
        """Return short description of model type."""
        if self.biodegradation_capacity:
            return "Mibitrans Instant Reaction"
        else:
            return "Mibitrans Linear"

    def run(self):
        """Calculate the concentration for all discretized x, y and t using the analytical transport model."""
        self._check_model_mode_before_run()
        with np.errstate(divide="ignore", invalid="ignore"):
            self.cxyt = self._calculate_concentration_for_all_xyt()
        return Results(self)

    def sample(self, x_position, y_position, time):
        """Give concentration at any given position and point in time.

        Args:
            x_position (float): x position in domain extent [m].
            y_position (float): y position in domain extent [m].
            time (float): time for which concentration is sampled [days].

        Returns:
            concentration (float): concentration at given position and point in time [g/m^3].

        """
        # Different sample method than parent class, as field-wide calculations use array indices
        for par, value in locals().items():
            if par != "self":
                validate_input_values(par, value)

        self._pre_run_initialization_parameters()

        def integrand(t, sz):
            div_term = 2 * np.sqrt(self.disp_y * t**4)
            inner_term = self._src_pars.depth / (2 * np.sqrt(self.disp_z * t**4))
            integrand_results = (
                1
                / (t**3)
                * (
                    np.exp(
                        (-self.k_source - self._decay_rate) * t**4
                        - (x_position - self.rv * t**4) ** 2 / (4 * self.disp_x * t**4)
                    )
                    * (
                        erfc((y_position - self.y_source[sz]) / div_term)
                        - erfc((y_position + self.y_source[sz]) / div_term)
                    )
                    * (erfc(-inner_term) - erfc(inner_term))
                )
            )
            return integrand_results

        conc_array = np.zeros(len(self.c_source))
        error_array = np.zeros(len(self.c_source))
        time = time ** (1 / 4)
        with np.errstate(divide="ignore", invalid="ignore"):
            for sz in range(len(self.c_source)):
                integral_term, error = quad(integrand, 0, time, limit=10000, args=(sz,))
                source_term = (
                    self.c_source[sz] * x_position / (8 * np.sqrt(np.pi * self.disp_x)) * np.exp(-self.k_source * time)
                )
                conc_array[sz] = 4 * integral_term * source_term
                error_array[sz] = error
            concentration = np.sum(conc_array)
            if self._mode == "instant_reaction":
                concentration -= self.biodegradation_capacity
                if concentration < 0:
                    concentration = 0
        return concentration

    def _pre_run_initialization_parameters(self):
        super()._pre_run_initialization_parameters()
        self.disp_x = self._hyd_pars.alpha_x * self.rv + self._att_pars.diffusion
        self.disp_y = self._hyd_pars.alpha_y * self.rv + self._att_pars.diffusion
        self.disp_z = self._hyd_pars.alpha_z * self.rv + self._att_pars.diffusion
        # self.integral_term = np.zeros(self.ttt.shape)
        # Stores integral error for each time step and source zone
        self.error_size = np.zeros((len(self._src_pars.source_zone_boundary), len(self.t)))

    def _calculate_concentration_for_all_xyt(self):
        cxyt = self.cxyt.copy()
        for sz in range(len(self.c_source)):
            integral_sum = self._equation_term_integral(sz)
            source_term = self._equation_term_source(sz)
            cxyt[:, :, 1:] += integral_sum[:, :, 1:] * source_term
            # If x=0, equation resolves to c=0, therefore, x=0 needs to be evaluated separately
            cxyt[:, :, 0] += self._equation_term_source_x_is_zero(sz)[:, :, 0]
        if self._mode == "instant_reaction":
            self.cxyt_noBC = cxyt.copy()
            cxyt -= self.biodegradation_capacity
            cxyt = np.where(cxyt < 0, 0, cxyt)
        return cxyt

    def _equation_term_integral(self, sz):
        integral_term = np.zeros(self.cxyt.shape)
        for j in range(len(self.t)):
            if self.verbose:
                print("integrating for source zone ", sz, " and t =", self.t[j], "days")
            if j == 0:
                lower_bound = 0
            else:
                lower_bound = self.t[j - 1]
            upper_bound = self.t[j]
            integral_term[j, :, 1:], self.error_size[sz, j] = quad_vec(
                self._equation_integrand, lower_bound, upper_bound, limit=10000 // len(self.t), args=(sz,)
            )
        integral_sum = np.cumsum(integral_term, axis=0)
        return integral_sum

    def _equation_integrand(self, t, sz):
        term = 1 / (t ** (3 / 2)) * self._equation_term_x(t) * self._equation_term_y(t, sz) * self._equation_term_z(t)
        term[np.isnan(term)] = 0
        return term

    def _equation_term_x(self, t):
        term = np.exp(
            (-self.k_source - self._decay_rate) * t - (self.xxx[:, :, 1:] - self.rv * t) ** 2 / (4 * self.disp_x * t)
        )
        term[np.isnan(term)] = 0
        return term

    def _equation_term_y(self, t, sz):
        div_term = 2 * np.sqrt(self.disp_y * t)
        term = erfc((self.yyy - self.y_source[sz]) / div_term) - erfc((self.yyy + self.y_source[sz]) / div_term)
        term[np.isnan(term)] = 0
        return term

    def _equation_term_z(self, t):
        if t == 0 or self.disp_z == 0:
            inner_term = 2
        else:
            inner_term = self._src_pars.depth / (2 * np.sqrt(self.disp_z * t))
        return erfc(-inner_term) - erfc(inner_term)

    def _equation_term_source(self, sz):
        return (
            self.c_source[sz]
            * self.xxx[:, :, 1:]
            / (8 * np.sqrt(np.pi * self.disp_x))
            * np.exp(-self.k_source * self.ttt)
        )

    def _equation_term_source_x_is_zero(self, sz):
        # Select y-positions of current source zone
        zone_location = np.where(abs(self.yyy) <= self.y_source[sz], 1, 0)
        return self.c_source[sz] * zone_location * np.exp(-self.k_source * self.ttt)

short_description property

Return short description of model type.

__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose=False)

Initialize model object.

Parameters:

Name Type Description Default
hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters)

Dataclass object containing hydrological parameters from HydrologicalParameters.

required
attenuation_parameters (mibitrans.data.read.AttenuationParameters)

Dataclass object containing adsorption, degradation and diffusion parameters from AttenuationParameters.

required
source_parameters (mibitrans.data.read.SourceParameters)

Dataclass object containing source parameters from SourceParameters.

required
model_parameters (mibitrans.data.read.ModelParameters)

Dataclass object containing model parameters from ModelParameters.

required
verbose bool

Verbose mode. Defaults to False.

False

Attributes:

Name Type Description
mode str)

Current model mode. Is ‘linear’ by default. Once instant reaction parameters are provided. Use this attribute to switch between ‘linear’ and ‘instant_reaction’ models.

cxyt np.ndarray)

Output array containing concentrations in model domain, in [g/m^3]. Indexed as [t,y,x]

relative_cxyt np.ndarray)

Output array with concentrations in model domain, divided by the maximum source zone concentration at t=0. Indexed as [t,y,x].

x np.ndarray)

Discretized model x-dimension, in [m].

y np.ndarray)

Discretized model y-dimension, in [y].

t np.ndarray)

Discretized model t-dimension, in [days].

c_source np.ndarray)

Nett source zone concentrations, accounting for source superposition, in [g/m^3].

vr float)

Retarded groundwater flow velocity, in [m/d].

k_source float)

Source zone decay rate, in [1/days].

biodegradation_capacity float)

Maximum capacity of biodegradation given provided electron acceptor

Functions:

Name Description
run

Run model with current parameters, returns Results object.

sample

Calculate concentration at any given position and point in time.

instant_reaction

Activate the instant reaction model by providing electron acceptor concentrations. And optionally electron acceptor utilization factors. Switch between model modes by using the mode attribute.

Raises:

Type Description
TypeError

If input is not of the correct Dataclass.

Source code in mibitrans/transport/models.py
def __init__(
    self,
    hydrological_parameters,
    attenuation_parameters,
    source_parameters,
    model_parameters,
    verbose=False,
):
    """Initialize model object.

    Args:
        hydrological_parameters (mibitrans.data.parameters.HydrologicalParameters) : Dataclass object containing
            hydrological parameters from HydrologicalParameters.
        attenuation_parameters (mibitrans.data.read.AttenuationParameters) : Dataclass object containing adsorption,
            degradation and diffusion parameters from AttenuationParameters.
        source_parameters (mibitrans.data.read.SourceParameters) : Dataclass object containing source parameters
            from SourceParameters.
        model_parameters (mibitrans.data.read.ModelParameters) : Dataclass object containing model parameters from
            ModelParameters.
        verbose (bool, optional): Verbose mode. Defaults to False.

    Attributes:
        mode (str) : Current model mode. Is 'linear' by default. Once instant reaction parameters are provided. Use
            this attribute to switch between 'linear' and 'instant_reaction' models.
        cxyt (np.ndarray) : Output array containing concentrations in model domain, in [g/m^3]. Indexed as [t,y,x]
        relative_cxyt (np.ndarray) : Output array with concentrations in model domain, divided by the maximum source
            zone concentration at t=0. Indexed as [t,y,x].
        x (np.ndarray) : Discretized model x-dimension, in [m].
        y (np.ndarray) : Discretized model y-dimension, in [y].
        t (np.ndarray) : Discretized model t-dimension, in [days].
        c_source (np.ndarray) : Nett source zone concentrations, accounting for source superposition, in [g/m^3].
        vr (float) : Retarded groundwater flow velocity, in [m/d].
        k_source (float) : Source zone decay rate, in [1/days].
        biodegradation_capacity (float) : Maximum capacity of biodegradation given provided electron acceptor
        concentrations, in [g/m^3].

    Methods:
        run : Run model with current parameters, returns Results object.
        sample : Calculate concentration at any given position and point in time.
        instant_reaction : Activate the instant reaction model by providing electron acceptor concentrations. And
            optionally electron acceptor utilization factors. Switch between model modes by using the mode
            attribute.

    Raises:
        TypeError : If input is not of the correct Dataclass.

    """
    super().__init__(hydrological_parameters, attenuation_parameters, source_parameters, model_parameters, verbose)

run()

Calculate the concentration for all discretized x, y and t using the analytical transport model.

Source code in mibitrans/transport/models.py
def run(self):
    """Calculate the concentration for all discretized x, y and t using the analytical transport model."""
    self._check_model_mode_before_run()
    with np.errstate(divide="ignore", invalid="ignore"):
        self.cxyt = self._calculate_concentration_for_all_xyt()
    return Results(self)

sample(x_position, y_position, time)

Give concentration at any given position and point in time.

Parameters:

Name Type Description Default
x_position float

x position in domain extent [m].

required
y_position float

y position in domain extent [m].

required
time float

time for which concentration is sampled [days].

required

Returns:

Name Type Description
concentration float

concentration at given position and point in time [g/m^3].

Source code in mibitrans/transport/models.py
def sample(self, x_position, y_position, time):
    """Give concentration at any given position and point in time.

    Args:
        x_position (float): x position in domain extent [m].
        y_position (float): y position in domain extent [m].
        time (float): time for which concentration is sampled [days].

    Returns:
        concentration (float): concentration at given position and point in time [g/m^3].

    """
    # Different sample method than parent class, as field-wide calculations use array indices
    for par, value in locals().items():
        if par != "self":
            validate_input_values(par, value)

    self._pre_run_initialization_parameters()

    def integrand(t, sz):
        div_term = 2 * np.sqrt(self.disp_y * t**4)
        inner_term = self._src_pars.depth / (2 * np.sqrt(self.disp_z * t**4))
        integrand_results = (
            1
            / (t**3)
            * (
                np.exp(
                    (-self.k_source - self._decay_rate) * t**4
                    - (x_position - self.rv * t**4) ** 2 / (4 * self.disp_x * t**4)
                )
                * (
                    erfc((y_position - self.y_source[sz]) / div_term)
                    - erfc((y_position + self.y_source[sz]) / div_term)
                )
                * (erfc(-inner_term) - erfc(inner_term))
            )
        )
        return integrand_results

    conc_array = np.zeros(len(self.c_source))
    error_array = np.zeros(len(self.c_source))
    time = time ** (1 / 4)
    with np.errstate(divide="ignore", invalid="ignore"):
        for sz in range(len(self.c_source)):
            integral_term, error = quad(integrand, 0, time, limit=10000, args=(sz,))
            source_term = (
                self.c_source[sz] * x_position / (8 * np.sqrt(np.pi * self.disp_x)) * np.exp(-self.k_source * time)
            )
            conc_array[sz] = 4 * integral_term * source_term
            error_array[sz] = error
        concentration = np.sum(conc_array)
        if self._mode == "instant_reaction":
            concentration -= self.biodegradation_capacity
            if concentration < 0:
                concentration = 0
    return concentration