In dynamic-R/RTM: Reaction Transport Modelling: Course Tutorials and Exercises
Aim
In parts I and II of the reader series "Local equilibrium chemistry in R", we have shown how to approach pH modeling in scenarios where proton dynamics is affected by 1 or 3 fast reversible reactions, respectively. The reactions corresponded to dissociation of ammonium or dissociation of carbonic acid and self-ionization of water. Here, we will combine what we have learned in those two parts and show how to model pH changes in a system where all of these fast reversible reactions occur simultaneously. An important step will be the generalization of the concept of alkalinity to include contributions from all charged species (carbonate, bicarbonate, ammonium, hydroxide, and protons). We will illustrate the model by studying pH dynamics in water induced by organic matter degradation. Additionally, we will revisit the ammonia degassing model originally developed in part I of this reader series to illustrate the pH buffering provided by the dissolved carbonate species in the system.
Approach
The aproach followed in this reader is based on a set of tasks that you should be able to solve on your own if you have understood the material presented in parts I and II of this reader series. Most of the derivations and R-coding ideas presented in parts I and II will be directly reused here. We encourage you to first try to solve the tasks on your own. If you feel you are stuck, you can look up the solutions provided at the end of this document.
Model formulation
We consider here four fast reversible reactions:
$$
\begin{array}{rl}
Reaction~1!!: & H_2CO_3 \rightleftarrows HCO_3^- + H^+ \[4mm]
Reaction~2!!: & HCO_3^- \rightleftarrows CO_3^{2-} + H^+ \[4mm]
Reaction~3!!: & NH_4^+ \rightleftarrows NH_3 + H^+ \[4mm]
Reaction~4!!: & H_2O \rightleftarrows OH^- + H^+
\end{array}
$$
Additionally, we consider two slow irreversible processes, one that adds $H_2CO_3$ to the system at a rate $R_C$, and one that adds $NH_3$ to the system at a rate $R_N$. We assume that the forward and backward reactions 1--4 are elementary reactions. We denote the corresponding rate constants as $k_{if}$ (units of $s^{-1}$) and $k_{ib}$ (units of $(mol~L^{-1})^{-1}~s^{-1}$) for the reaction $i$, where $i=1, \dots ,4$.
Task 1: Mass balance equations
- Formulate the mass balance equations considering the reactions 1--4 and the processes with rate $R_C$ and $R_N$.
You should arrive at 8 differential equations for the following 8 state variables: $[H_2CO_3]$, $[HCO_3^-]$, $[CO_3^{2-}]$, $[NH_4^+]$, $[NH_3]$, $[H_2O]$, $[OH^-]$ and $[H^+]$.
Task 2: Equilibrium conditions
Assume that the rates of the slow processes are zero ($R_C=R_N=0$). Using the mass balance equations formulated in Task 1,
- derive relationships between the equilibrium concentrations of the state variables. Hint: consider that the time-derivatives of the state variables are zero when the system is in equilibrium.
You should arrive at four equations that relate the equilibrium concentrations of [$H_2CO_3$], [$HCO_3^-$], [$CO_3^{2-}$], $[NH_4^+]$, $[NH_3]$, $[OH^-]$, and [$H^+$]. To simplify the relationships, define the ratios between the forward and backward rate constants as the so-called dissociation constants: $K_1 = k_{1f}/k_{1b}$, $K_2 = k_{2f}/k_{2b}$, $K_n = k_{3f}/k_{3b}$, and $K_w' = k_{4f}/k_{4b}$.
Task 3: Lump-sum species
- Show that if you define the lump-sum species according to
$$
\begin{array}{rcl}
DIC & = & H_2CO_3 + HCO_3^- + CO_3^{2-}, \[4mm]
NH_x & = & NH_4^+ + NH_3, \[4mm]
H_xO & = & H_2O + OH^-, \[4mm]
ALK & = & HCO_3^- + 2\cdot CO_3^{2-} + OH^- - NH_4^+ - H^+,
\end{array}
$$
their time derivatives are only affected by the slow processes (rates $R_C$ and $R_N$) but not by the reversible reactions 1--4.
Note that these new state variables correspond to the total concentration of dissolved inorganic carbon ($DIC$), total concentration of dissolved ammonia ($NH_x$), total concentration of hydroxide ($H_xO$), and the excess concentration of negatively charged species ($ALK$). Also note that $ALK$ introduced here is a further generalization of the alkalinity introduced in part II of this reader series. This generalization can proceed further to include other rapidly dissociating chemical species in the system (e.g., borate, silicate, phosphate, etc.) and ultimately results in the so-called total alkalinity. However, we will limit our analysis here to the contributions from carbonate, ammonium and water.
Task 4: From lump-sum species to original species
Assume that the system is in an equilibrium and that the concentrations $[DIC]$, $[NH_x]$, $[H_xO$], $[ALK]$ and $[H^+]$ are known. Note that $pH = -\log_{10}{[H^+]}$, where $[H^+]$ is expressed in $mol~L^{-1}$ of solution. Using the equilibrium relationships found in Task 2,
-
derive expressions for $[H_2CO_3]$, $[HCO_3^-]$, $[CO_3^{2-}]$, $[NH_4^+]$, $[NH_3]$, and $[OH^-]$ as a function of $[DIC]$, $[NH_x]$, $[H_xO]$ and $[H^+]$.
-
find a relationship between $[ALK]$, $[DIC]$, $[NH_x]$, $[H_xO]$ and $[H^+]$ from which one of the five quantities can be calculated provided that the other four are known.
A similar task has been solved in parts I and II of this reader series. Nevertheless, it is useful to revise the steps and check the formulas again.
Task 5: Application to organic matter mineralization
A possible process that adds $H_2CO_3$ and $NH_3$ to the system is the oxic mineralization of organic matter (OM), where OM is represented by the chemical formula $CH_2O(NH_3){NC}$. Assume that the stoichiometry of this process is
$$
\quad CH_2O(NH_3){NC} + O_2 \rightarrow H_2CO_3 + {NC} \cdot NH_3
$$
- Update your model developed in part II of this reader series so that it will allow you to predict the dynamics of the carbonate species, ammonia species, water species, alkalinity and $pH$ as a result of OM mineralization.
Use the following assumptions:
-
OM mineralization proceeds according to the first-order kinetics with respect to OM. Thus, consider the following rate expressions: $R_C = r_C\cdot [OM]$ and $R_N = NC\cdot R_C$, where $r_C$ is the rate constant ($s^{-1}$) and $[OM]$ is the OM concentration ($mol~C~L^{-1}$). Assume values of $r_C=0.1~d^{-1}$ and $NC=16/106$, the latter corresponding to the Redfield ratio of N:C in organic matter.
-
OM mineralization is a slow process relative to the reversible reactions 1--4. Thus, you can apply the "local equilibrium assumption".
Solve the model for 10 days. For a temperature of 25 $^\circ$C, salinity 0 and pressure of 1 bar, assume the following values of the equilibrium constants, as calculated using the R-package AquaEnv (Hofmann et al., 2010):
require(AquaEnv)
# water properties
a <- aquaenv(S=0, t=25, k1k2="millero")
rho <- a$density # density, g/L
MW <- 18.0153 # molar weight, g/mol
cHxO <- rho/MW # total H2O concentration, mol/L
# dissociation constants (converted from mol/kg-soln -> mol/L)
K1 <- a$K_CO2 * 1e-3 * rho # mol/L
K2 <- a$K_HCO3 * 1e-3 * rho # mol/L
Kn <- a$K_NH4 * 1e-3 * rho # mol/L
KW <- a$K_W * 1e-6 * rho^2 # (mol/L)^2
KW.prime <- KW / cHxO # mol/L
# process properties
rC <- 0.1 # d-1
NC <- 16/106 # molN/molC
knitr::kable(cbind(c("$K_1$", "$K_2$", "$K_n$", "$K_w'$", "$K_w$"),
format(c(K1, K2, Kn, KW.prime, KW),digits=3),
c(rep("$mol~L^{-1}$",4), "$mol^2~L^{-2}$")))
For the initial conditions, assume the following values:
knitr::kable(cbind(c("$pH$", "$[OM]$", "$[DIC]$", "$[NH_x]$", "$[H_xO]$"),
format(c(7.9, 100, 2000, 10, cHxO), digits=2),
c("--", rep("$\\mu mol~L^{-1}$",3), "$mol~L^{-1}$")))
Task 6: Does the formulation of slow reversible processes matter?
In the literature, processes such as OM mineralization are often formulated in different ways, depending on the preference or perspective. For example, the stoichiometry of OM mineralization can be represented by the following reaction
$$
\quad CH_2O(NH_3)_{NC} + O_2 \rightarrow HCO_3^- + H^- + {NC} \cdot NH_3.
$$
instead of the one shown above.
- Show that this alternative formulation does not affect the mass balance equations for the lump-sum variables $DIC$, $ALK$, $NH_x$ and $H_xO$.
Task 7: Application to ammonia degassing
The model developed in Task 5 should be equivalent to the ammonia degassing model developed in part I of this reader series in the limit of $DIC\rightarrow 0$ and $H_xO \rightarrow 0$.
- Check that this is true by comparing the results of the two models.
In this comparison, assume the rate constant of $\lambda = 1~d^{-1}$ for ammonia degassing, initial $pH=8$, and initial $[NH_x] = 1~mmol~L^{-1}$.
After this comparison,
- apply the model to study the magnitude of changes in water $pH$ associated with ammonia degassing at different $DIC$ concentrations ($0$, $1~mmol~L^{-1}$ and $2~mmol~L^{-1}$).
This application will illustrate the $pH$ buffering role of $DIC$ in natural waters.
dynamic-R/RTM documentation built on Feb. 28, 2025, 1:23 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Aim
In parts I and II of the reader series "Local equilibrium chemistry in R", we have shown how to approach pH modeling in scenarios where proton dynamics is affected by 1 or 3 fast reversible reactions, respectively. The reactions corresponded to dissociation of ammonium or dissociation of carbonic acid and self-ionization of water. Here, we will combine what we have learned in those two parts and show how to model pH changes in a system where all of these fast reversible reactions occur simultaneously. An important step will be the generalization of the concept of alkalinity to include contributions from all charged species (carbonate, bicarbonate, ammonium, hydroxide, and protons). We will illustrate the model by studying pH dynamics in water induced by organic matter degradation. Additionally, we will revisit the ammonia degassing model originally developed in part I of this reader series to illustrate the pH buffering provided by the dissolved carbonate species in the system.
Approach
The aproach followed in this reader is based on a set of tasks that you should be able to solve on your own if you have understood the material presented in parts I and II of this reader series. Most of the derivations and R-coding ideas presented in parts I and II will be directly reused here. We encourage you to first try to solve the tasks on your own. If you feel you are stuck, you can look up the solutions provided at the end of this document.
Model formulation
We consider here four fast reversible reactions:
$$ \begin{array}{rl} Reaction~1!!: & H_2CO_3 \rightleftarrows HCO_3^- + H^+ \[4mm] Reaction~2!!: & HCO_3^- \rightleftarrows CO_3^{2-} + H^+ \[4mm] Reaction~3!!: & NH_4^+ \rightleftarrows NH_3 + H^+ \[4mm] Reaction~4!!: & H_2O \rightleftarrows OH^- + H^+ \end{array} $$
Additionally, we consider two slow irreversible processes, one that adds $H_2CO_3$ to the system at a rate $R_C$, and one that adds $NH_3$ to the system at a rate $R_N$. We assume that the forward and backward reactions 1--4 are elementary reactions. We denote the corresponding rate constants as $k_{if}$ (units of $s^{-1}$) and $k_{ib}$ (units of $(mol~L^{-1})^{-1}~s^{-1}$) for the reaction $i$, where $i=1, \dots ,4$.
Task 1: Mass balance equations
- Formulate the mass balance equations considering the reactions 1--4 and the processes with rate $R_C$ and $R_N$.
You should arrive at 8 differential equations for the following 8 state variables: $[H_2CO_3]$, $[HCO_3^-]$, $[CO_3^{2-}]$, $[NH_4^+]$, $[NH_3]$, $[H_2O]$, $[OH^-]$ and $[H^+]$.
Task 2: Equilibrium conditions
Assume that the rates of the slow processes are zero ($R_C=R_N=0$). Using the mass balance equations formulated in Task 1,
- derive relationships between the equilibrium concentrations of the state variables. Hint: consider that the time-derivatives of the state variables are zero when the system is in equilibrium.
You should arrive at four equations that relate the equilibrium concentrations of [$H_2CO_3$], [$HCO_3^-$], [$CO_3^{2-}$], $[NH_4^+]$, $[NH_3]$, $[OH^-]$, and [$H^+$]. To simplify the relationships, define the ratios between the forward and backward rate constants as the so-called dissociation constants: $K_1 = k_{1f}/k_{1b}$, $K_2 = k_{2f}/k_{2b}$, $K_n = k_{3f}/k_{3b}$, and $K_w' = k_{4f}/k_{4b}$.
Task 3: Lump-sum species
- Show that if you define the lump-sum species according to $$ \begin{array}{rcl} DIC & = & H_2CO_3 + HCO_3^- + CO_3^{2-}, \[4mm] NH_x & = & NH_4^+ + NH_3, \[4mm] H_xO & = & H_2O + OH^-, \[4mm] ALK & = & HCO_3^- + 2\cdot CO_3^{2-} + OH^- - NH_4^+ - H^+, \end{array} $$ their time derivatives are only affected by the slow processes (rates $R_C$ and $R_N$) but not by the reversible reactions 1--4.
Note that these new state variables correspond to the total concentration of dissolved inorganic carbon ($DIC$), total concentration of dissolved ammonia ($NH_x$), total concentration of hydroxide ($H_xO$), and the excess concentration of negatively charged species ($ALK$). Also note that $ALK$ introduced here is a further generalization of the alkalinity introduced in part II of this reader series. This generalization can proceed further to include other rapidly dissociating chemical species in the system (e.g., borate, silicate, phosphate, etc.) and ultimately results in the so-called total alkalinity. However, we will limit our analysis here to the contributions from carbonate, ammonium and water.
Task 4: From lump-sum species to original species
Assume that the system is in an equilibrium and that the concentrations $[DIC]$, $[NH_x]$, $[H_xO$], $[ALK]$ and $[H^+]$ are known. Note that $pH = -\log_{10}{[H^+]}$, where $[H^+]$ is expressed in $mol~L^{-1}$ of solution. Using the equilibrium relationships found in Task 2,
-
derive expressions for $[H_2CO_3]$, $[HCO_3^-]$, $[CO_3^{2-}]$, $[NH_4^+]$, $[NH_3]$, and $[OH^-]$ as a function of $[DIC]$, $[NH_x]$, $[H_xO]$ and $[H^+]$.
-
find a relationship between $[ALK]$, $[DIC]$, $[NH_x]$, $[H_xO]$ and $[H^+]$ from which one of the five quantities can be calculated provided that the other four are known.
A similar task has been solved in parts I and II of this reader series. Nevertheless, it is useful to revise the steps and check the formulas again.
Task 5: Application to organic matter mineralization
A possible process that adds $H_2CO_3$ and $NH_3$ to the system is the oxic mineralization of organic matter (OM), where OM is represented by the chemical formula $CH_2O(NH_3){NC}$. Assume that the stoichiometry of this process is $$ \quad CH_2O(NH_3){NC} + O_2 \rightarrow H_2CO_3 + {NC} \cdot NH_3 $$
- Update your model developed in part II of this reader series so that it will allow you to predict the dynamics of the carbonate species, ammonia species, water species, alkalinity and $pH$ as a result of OM mineralization.
Use the following assumptions:
-
OM mineralization proceeds according to the first-order kinetics with respect to OM. Thus, consider the following rate expressions: $R_C = r_C\cdot [OM]$ and $R_N = NC\cdot R_C$, where $r_C$ is the rate constant ($s^{-1}$) and $[OM]$ is the OM concentration ($mol~C~L^{-1}$). Assume values of $r_C=0.1~d^{-1}$ and $NC=16/106$, the latter corresponding to the Redfield ratio of N:C in organic matter.
-
OM mineralization is a slow process relative to the reversible reactions 1--4. Thus, you can apply the "local equilibrium assumption".
Solve the model for 10 days. For a temperature of 25 $^\circ$C, salinity 0 and pressure of 1 bar, assume the following values of the equilibrium constants, as calculated using the R-package AquaEnv (Hofmann et al., 2010):
require(AquaEnv) # water properties a <- aquaenv(S=0, t=25, k1k2="millero") rho <- a$density # density, g/L MW <- 18.0153 # molar weight, g/mol cHxO <- rho/MW # total H2O concentration, mol/L # dissociation constants (converted from mol/kg-soln -> mol/L) K1 <- a$K_CO2 * 1e-3 * rho # mol/L K2 <- a$K_HCO3 * 1e-3 * rho # mol/L Kn <- a$K_NH4 * 1e-3 * rho # mol/L KW <- a$K_W * 1e-6 * rho^2 # (mol/L)^2 KW.prime <- KW / cHxO # mol/L # process properties rC <- 0.1 # d-1 NC <- 16/106 # molN/molC
knitr::kable(cbind(c("$K_1$", "$K_2$", "$K_n$", "$K_w'$", "$K_w$"), format(c(K1, K2, Kn, KW.prime, KW),digits=3), c(rep("$mol~L^{-1}$",4), "$mol^2~L^{-2}$")))
For the initial conditions, assume the following values:
knitr::kable(cbind(c("$pH$", "$[OM]$", "$[DIC]$", "$[NH_x]$", "$[H_xO]$"), format(c(7.9, 100, 2000, 10, cHxO), digits=2), c("--", rep("$\\mu mol~L^{-1}$",3), "$mol~L^{-1}$")))
Task 6: Does the formulation of slow reversible processes matter?
In the literature, processes such as OM mineralization are often formulated in different ways, depending on the preference or perspective. For example, the stoichiometry of OM mineralization can be represented by the following reaction $$ \quad CH_2O(NH_3)_{NC} + O_2 \rightarrow HCO_3^- + H^- + {NC} \cdot NH_3. $$ instead of the one shown above.
- Show that this alternative formulation does not affect the mass balance equations for the lump-sum variables $DIC$, $ALK$, $NH_x$ and $H_xO$.
Task 7: Application to ammonia degassing
The model developed in Task 5 should be equivalent to the ammonia degassing model developed in part I of this reader series in the limit of $DIC\rightarrow 0$ and $H_xO \rightarrow 0$.
- Check that this is true by comparing the results of the two models.
In this comparison, assume the rate constant of $\lambda = 1~d^{-1}$ for ammonia degassing, initial $pH=8$, and initial $[NH_x] = 1~mmol~L^{-1}$.
After this comparison,
- apply the model to study the magnitude of changes in water $pH$ associated with ammonia degassing at different $DIC$ concentrations ($0$, $1~mmol~L^{-1}$ and $2~mmol~L^{-1}$).
This application will illustrate the $pH$ buffering role of $DIC$ in natural waters.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.