revisions to supp tech chapter

doc/SuppTechInfo/ist-sorption.tex

@@ -1,4 +1,6 @@
-The Groundwater Transport (GWT) Model \citep{modflow6gwt} in \mf was designed to simulate a range of solute-transport processes including the sorption and desorption of solute mass within aquifer material.  For the mobile domain, sorption can be represented using a linear isotherm as well as the nonlinear Freundlich and Langmuir isotherms.  For the ``immobile'' domain, the GWT model was limited to representing linear isotherms only.  The program checked to make sure that if sorption was represented in the immobile domain that sorption in the mobile domain was represented using a linear isotherm.  As reported by \mf users, there is a need to support nonlinear sorption isotherms for mobile-immobile domain simulations.  The purpose of this chapter is to describe the mathematical approach that was used to implement the nonlinear Freundlich and Langmuir isotherms to represent sorption in the immobile domain. There is still a requirement to use consistent sorption approaches in both the mobile and immobile domains, but with these changes, any of the three sorption isotherms can now be used as long as the same isotherm is used for both the mobile and immobile domains.
+The Groundwater Transport (GWT) Model \citep{modflow6gwt} in \mf was designed to simulate a range of solute-transport processes including the sorption and desorption of solute within the mobile and immobile domains.  For the mobile domain, sorption can be represented using a linear isotherm as well as the nonlinear Freundlich and Langmuir isotherms.  The type of sorption isotherm used for the immobile domain is specified by the user in the input file for the Mobile Storage and Transfer (MST) Package.  For the immobile domain, which is activated using the Immobile Storage and Transfer (IST) Package, the GWT model was limited to representing linear isotherms only.  The program checked to make sure that if sorption was represented in the immobile domain that sorption in the mobile domain was represented using a linear isotherm.  This was the behavior of \mf up to and including version 6.5.0.
+
+As reported by \mf users, some problems may require use of nonlinear sorption isotherms for mobile-immobile domain simulations.  In response to this need, the IST Package was extended to support the nonlinear Freundlich and Langmuir isotherms, in addition to the linear isotherm.  The purpose of this chapter is to describe the mathematical approach that was used to implement these nonlinear isotherms in the IST Package. The \mf program still requires use consistent sorption approaches in both the mobile and immobile domains, but with these changes, any of the three sorption isotherms can now be used as long as the same type of isotherm is used for both the mobile and immobile domains.
 
 Chapter \ref{ch:sorption} of this document presents a revised parameterization of the mobile and immobile domains.  Included in Chapter 9 is the revised form of the partial differential equation governing solute mass within an immobile domain (equation \ref{eqn:gwtistpde}).  Equation~\ref{eqn:gwtistpde} is reproduced here as,
 
@@ -27,7 +29,7 @@
 $S_w$ is the water saturation defined as the volume of water per volume of voids ($L^3/L^3$), 
 and $C$ is the mobile-domain volumetric concentration of solute expressed as mass of dissolved solute per unit volume of mobile-domain fluid ($M/L^3$).
 
-The GWT Model documentation report \citep{modflow6gwt} describes the approach for including the effects of an immobile domain on solute transport.   The approach is described in \cite{zheng2002} and involves discretizing the solute balance equation for the immobile domain.  \cite{modflow6gwt} present a form of the discretized balance equation in their equation 7--4.  That balance equation was based on a linear sorption isotherm and included the outdated parameterization, which was revised according to the approach described in Chapter \ref{ch:sorption}.  The following is a discretized form of equation~\ref{eqn:gwtistpde2} and includes a more general representation of the sorption isotherm as well as the updated parameterization described in Chapter \ref{ch:sorption}:
+The GWT Model documentation report \citep{modflow6gwt} describes the approach for including the effects of an immobile domain on solute transport.   The approach is based on the method described in \cite{zheng2002} and involves discretizing the solute balance equation for the immobile domain.  \cite{modflow6gwt} present a form of the discretized balance equation in their equation 7--4.  That balance equation was based on a linear sorption isotherm and used the outdated parameterization, which was revised according the description in Chapter \ref{ch:sorption}.  The following discretized form of equation~\ref{eqn:gwtistpde2} includes a more general representation of the sorption isotherm, in order to support both linear and nonlinear expressions, as well as the updated parameterization described in Chapter \ref{ch:sorption},
 
 \begin{equation}
 \label{eqn:imdomain2}
@@ -46,9 +48,9 @@
 \end{split}
 \end{equation}
 
-\noindent where $V_{cell}$ is the volume of the model cell and $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ is a sorption term that describes how the sorbate concentration changes in response to changes in the dissolved immobile domain concentration.
+\noindent where $V_{cell}$ is the volume of the model cell and $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ is a sorption term that describes the change in sorbate concentration versus change in the dissolved immobile domain concentration.
 
-The $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ term in equation~\ref{eqn:imdomain2} can be determined from the mathematical expression used for the different isotherms.   Mathematical expressions for the linear, Freundlich, and Langmuir isotherms, written in terms of immobile domain concentrations with subscript $im$, are
+The $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ term in equation~\ref{eqn:imdomain2} can be determined from the mathematical expression used for the different isotherms.   Mathematical expressions for the linear, Freundlich, and Langmuir isotherms, respectively, written in terms of immobile domain concentrations with subscript $im$, are
 
 \begin{equation}
 \label{eqn:linear}
@@ -69,7 +71,7 @@
 
 \noindent where $K_d$ is the linear distribution coefficient ($L^3/M$), $K_f$ is the Freundlich constant $(L^3 / M)^a$, $a$ is the dimensionless Freundlich exponent, $K_l$ is the Langmuir constant $(L^3 / M)$ and $\overline{S}$ is the total concentration of sorption sites available $(M/M)$.
 
-By differentiating equations \ref{eqn:linear} to \ref{eqn:langmuir}, expressions for $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ can be written as
+By differentiating equations \ref{eqn:linear} to \ref{eqn:langmuir}, corresponding expressions for $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ can be written as
 
 \begin{equation}
 \label{eqn:lineardcbardc}
@@ -90,4 +92,6 @@
 
 \noindent for the linear, Freundlich, and Langmuir isotherms, respectively.  
 
-Unlike the linear isotherm, the Freundlich and Langmuir isotherms have a nonlinear relation to $C_{im}$.  In order to calculate a value for $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ in equations~\ref{eqn:freundlichdcbardc} and \ref{eqn:langmuirdcbardc}, the $C_{im}$ value used in these equations is approximated as the arithmetic average of $C_{im}^{t}$ and the value of $C_{im}^{t + \Delta t}$ from the previous outer iteration.  Use of a previous iterate in these nonlinear calculations for $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ may cause convergence problems for some applications, compared to use of the linear expression, which does not depend on $C_{im}$.  
+Unlike the linear isotherm, the Freundlich and Langmuir isotherms have a nonlinear relation to $C_{im}$.  In order to calculate a value for $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ in equations~\ref{eqn:freundlichdcbardc} and \ref{eqn:langmuirdcbardc}, the $C_{im}$ value used in these equations is approximated as the arithmetic average of $C_{im}^{t}$ and the value of $C_{im}^{t + \Delta t}$ from the previous outer iteration.  Use of a previous iterate in these nonlinear calculations for $\frac{\partial \overline{C}_{im}}{\partial C_{im}}$ may cause convergence problems for some applications, compared to use of the linear expression, which does not depend on $C_{im}$.
+
+The remaining details of the implementation are identical to those described by \cite{modflow6gwt}.  Equation~\ref{eqn:imdomain2} is solved for $C_{im}^{t + \Delta t}$, which is then substituted into the equation for the transfer of mass between the mobile and immobile domains.  The resulting terms are then added to the system of equations and solved numerically.