[CHNOSZ-commits] r151 - in pkg/CHNOSZ: . demo inst man vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Wed Feb 15 16:34:03 CET 2017
Author: jedick
Date: 2017-02-15 16:34:03 +0100 (Wed, 15 Feb 2017)
New Revision: 151
Modified:
pkg/CHNOSZ/DESCRIPTION
pkg/CHNOSZ/demo/buffer.R
pkg/CHNOSZ/inst/NEWS
pkg/CHNOSZ/man/CHNOSZ-package.Rd
pkg/CHNOSZ/vignettes/anintro.Rmd
Log:
anintro.Rmd: copyediting
Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION 2017-02-14 11:37:01 UTC (rev 150)
+++ pkg/CHNOSZ/DESCRIPTION 2017-02-15 15:34:03 UTC (rev 151)
@@ -1,6 +1,6 @@
-Date: 2017-02-14
+Date: 2017-02-15
Package: CHNOSZ
-Version: 1.0.8-40
+Version: 1.0.8-41
Title: Chemical Thermodynamics and Activity Diagrams
Author: Jeffrey Dick
Maintainer: Jeffrey Dick <j3ffdick at gmail.com>
Modified: pkg/CHNOSZ/demo/buffer.R
===================================================================
--- pkg/CHNOSZ/demo/buffer.R 2017-02-14 11:37:01 UTC (rev 150)
+++ pkg/CHNOSZ/demo/buffer.R 2017-02-15 15:34:03 UTC (rev 151)
@@ -28,6 +28,7 @@
text(a$vals[[1]][13], mean(sapply(d$plotvals, c)[13, ]), logact)
}
# add legend
-legend("topright", legend = c(describe.property("P", 300), describe.basis(ibasis=c(2,4)),
- "minerals", "formaldehyde", "HCN"), lty=c(NA, NA, NA, 1, 3, 2), col=c(NA, NA, NA, 3, 1, 1),
- bg="white", cex=0.9)
+legend("topright", legend = c("minerals", "formaldehyde", "HCN"),
+ lty=c(1, 3, 2), col=c(3, 1, 1), bg="white", cex=0.9)
+legend("bottomright", legend = c(describe.property("P", 300),
+ describe.basis(ibasis=c(2,4))), bg="white", cex=0.9)
Modified: pkg/CHNOSZ/inst/NEWS
===================================================================
--- pkg/CHNOSZ/inst/NEWS 2017-02-14 11:37:01 UTC (rev 150)
+++ pkg/CHNOSZ/inst/NEWS 2017-02-15 15:34:03 UTC (rev 151)
@@ -1,4 +1,4 @@
-CHANGES IN CHNOSZ 1.0.8-40 (2017-02-14)
+CHANGES IN CHNOSZ 1.0.8-41 (2017-02-15)
---------------------------------------
DOCUMENTATION:
@@ -6,11 +6,11 @@
- Replace anintro.Rnw (Sweave) with rewritten anintro.Rmd (knitr, Tufte
style).
-- Add eos-regress.Rmd vignette; update related functions.
+- Add eos-regress.Rmd.
-- Demos: add bugstab.R, rename buffer.R to protbuff.R, add new buffer.R,
- remove nucleobase.R, and add protein.equil.R, add.obigt.R, and
- affinity.R (examples moved from help pages).
+- Demos: Move examples from help pages to protein.equil.R, add.obigt.R,
+ and affinity.R; add bugstab.R; update buffer.R and rename old buffer.R
+ to protbuff.R; remove nucleobase.R.
NEW DATA AND FEATURES:
Modified: pkg/CHNOSZ/man/CHNOSZ-package.Rd
===================================================================
--- pkg/CHNOSZ/man/CHNOSZ-package.Rd 2017-02-14 11:37:01 UTC (rev 150)
+++ pkg/CHNOSZ/man/CHNOSZ-package.Rd 2017-02-15 15:34:03 UTC (rev 151)
@@ -55,7 +55,7 @@
\item Activity diagrams - plot the equilibrium activities at a single point (as barplots), or as a function of one (species activity diagrams) or two (predominance diagrams) variables (\code{\link{diagram}}).
- \item Buffer calculations (\bold{experimental}) - compute activities of basis species that are determined by a buffer of one or more species (e.g., pyrite-pyrrhotite-magnetite; acetic acid-CO2) (\code{\link{buffer}}).
+ \item Buffer calculations - compute activities of basis species that are determined by a buffer of one or more species (e.g., pyrite-pyrrhotite-magnetite; acetic acid-CO2) (\code{\link{buffer}}).
\item Activity statistics (\bold{experimental}) - calculate summary statistics for equilibrium activities of species (\code{\link{revisit}}).
Modified: pkg/CHNOSZ/vignettes/anintro.Rmd
===================================================================
--- pkg/CHNOSZ/vignettes/anintro.Rmd 2017-02-14 11:37:01 UTC (rev 150)
+++ pkg/CHNOSZ/vignettes/anintro.Rmd 2017-02-15 15:34:03 UTC (rev 151)
@@ -85,18 +85,18 @@
# First steps
This document introduces the usage of CHNOSZ, a package for the [R software environment](http://r-project.org).
-For more information on R, see "An Introduction to R" and other documents in [The R Manuals](http://cran.r-project.org/manuals.html).
+For more information on R, see "[An Introduction to R](http://cran.r-project.org/manuals.html)" and the [Contributed Documentation](https://cran.r-project.org/other-docs.html) for R.
CHNOSZ has been developed since 2006 to support research projects in geochemistry and compositional biology.
The package provides functions and a thermodynamic database that can be used to calculate the stoichiometric and energetic properties of reactions among minerals and organic or inorganic aqueous species.
-These same functions enable calculations of chemical affinities and metastable equilibrium distributions of proteins.
+These functions also enable calculations of chemical affinities and metastable equilibrium distributions of proteins.
## Installing and loading CHNOSZ
-After installing R, install CHNOSZ by selecting the "Install packages from CRAN" or similar menu item in the R GUI or by using the following command:
+After starting R, install CHNOSZ by selecting the "Install packages from CRAN" or similar menu item in the R GUI or by using the following command:
```{marginfigure}
-Or, install the package from a package file, which you can download from [CRAN](https://cran.r-project.org/package=CHNOSZ) or from the [CHNOSZ website](http://chnosz.net/download).
-This vignette depends on some features that are only available in the development version of CHNOSZ, which can be found at [R-Forge](https://r-forge.r-project.org/projects/chnosz/).
+Or, install the package from a package file, which you can download from [CRAN](https://cran.r-project.org/package=CHNOSZ) or (for the development version) from [R-Forge](https://r-forge.r-project.org/projects/chnosz/).
+This vignette depends on some features that are only available in the development version of CHNOSZ.
```
```{r install_CHNOSZ, eval=FALSE}
install.packages("CHNOSZ")
@@ -130,24 +130,24 @@
Many functions in CHNOSZ have no side effects.
That is, the function only returns a result; to use the result elsewhere, it can be assigned to a variable with `<-`.
-In this document (but not in the code chunks), the names of these functions are set in <span style="color:green">green text</span>.
+In this document, the names of these functions are set in <span style="color:green">green text</span> (not applicable to the code chunks).
```{marginfigure}
-When they are mentioned, names of functions in the base and recommended packages of R are said to belong to R. Example: Use R's `plot()` to plot the points.
+When they are mentioned, names of functions in the base and recommended packages of R are said to belong to R. Example: Use R's `plot()` to plot the data.
```
Major functions without side effects in CHNOSZ are:
-* <span style="color:green">`info()`</span>: search for species in the thermodynamic database
-* <span style="color:green">`subcrt()`</span>: calculate the thermodynamic properties of species and reactions
-* <span style="color:green">`affinity()`</span>: calculate the affinities of formation reactions using given chemical activities
-* <span style="color:green">`equilibrate()`</span>: calculate the equilibrium chemical activities of the species of interest
-* <span style="color:green">`diagram()`</span>: plot the results
+* <span style="color:green">`info()`</span>: search for species in the thermodynamic database;
+* <span style="color:green">`subcrt()`</span>: calculate the thermodynamic properties of species and reactions;
+* <span style="color:green">`affinity()`</span>: calculate the affinities of formation reactions using given chemical activities;
+* <span style="color:green">`equilibrate()`</span>: calculate the equilibrium chemical activities of the species of interest;
+* <span style="color:green">`diagram()`</span>: plot the results.
Some functions in CHNOSZ do have side effects: they modify the `thermo` data object in the current R session.
-In this document, the names of these functions are set in <span style="color:red">red text</span>.
-The major functions with side effects are:
+In this document, the names of these functions are set in <span style="color:red">red text</span> (not in the code chunks).
+Major functions with side effects are:
-* <span style="color:red">`basis()`</span>: set the basis species and their chemical activities
-* <span style="color:red">`species()`</span>: set the species of interest and their (non-equilibrium) chemical activities
+* <span style="color:red">`basis()`</span>: set the basis species and their chemical activities;
+* <span style="color:red">`species()`</span>: set the species of interest and their (non-equilibrium) chemical activities;
* <span style="color:red">`data(thermo)`</span>: reset the database, restoring all settings to their default values.
The following pseudocode shows a common sequence of commands.
@@ -163,7 +163,7 @@
# Thermodynamic database and chemical formulas
-While an attempt has been made to provide a primary database (`OBIGT.csv`) that is generally internally consistent, all thermodynamic data, calculations, and examples are provided *as is*.
+While an attempt has been made to provide a primary database (`OBIGT.csv`) that is generally internally consistent, all thermodynamic data, calculations, and examples are provided **as is**.
```{marginfigure}
For crucial problems, check not only the accuracy of the database, but also the *suitability of the data* for your problem.
If there is any doubt about the suitability of data, please consult the primary sources (see <span style="color:blue">?`thermo.refs`</span>).
@@ -225,7 +225,7 @@
info(88)$formula
```
-We can use <span style="color:green">`makeup()`</span> to count the elements in the formula, followed by <span style="color:green">`as.chemical.formula()`</span> to return to the formula:
+We can use <span style="color:green">`makeup()`</span> to count the elements in the formula, followed by <span style="color:green">`as.chemical.formula()`</span> to rewrite the formula on one line:
```{r makeup_88}
makeup(88)
as.chemical.formula(makeup(88))
@@ -242,7 +242,7 @@
To calculate the standard molal properties of species and reactions, use <span style="color:green">`subcrt()`</span>.
```{marginfigure}
-The inspiration for the name <span style="color:green">`subcrt()`</span>, and the source of the FORTRAN subroutine used to calculate the thermodynamic properties of H<sub>2</sub>O, is SUPCRT (Johnson et al., 1992).
+The inspiration for the name <span style="color:green">`subcrt()`</span>, and the source of the Fortran subroutine used to calculate the thermodynamic properties of H<sub>2</sub>O, is SUPCRT (Johnson et al., 1992).
```
<sup>[- at JOH92]</sup>
If no reaction coefficients are given, <span style="color:green">`subcrt()`</span> calculates the standard molal properties of invididual species:
@@ -250,13 +250,13 @@
subcrt("water")
```
-That uses the default temperature and pressure settings: equally-spaced temperature grid from 0 to 350 °C at *P*<sub>sat</sub>.
+That uses the default temperature and pressure settings, i.e. equally spaced temperature intervals from 0 to 350 °C at *P*<sub>sat</sub>.
```{marginfigure}
*P*<sub>sat</sub> is 1 bar below 100 °C, or the pressure of liquid-vapor coexistence at higher temperatures.
```
-The columns in the output are temperature, pressure, density of water, logarithm of the equilibrium constant (only meaningful for reactions; see below), standard molal Gibbs energy and enthalpy of formation from the elements, standard molal entropy, volume, and heat capacity.
+The columns in the output are temperature, pressure, density of water, logarithm of the equilibrium constant (only meaningful for reactions; [see below](#Properties-of-reactions)), standard molal Gibbs energy and enthalpy of formation from the elements, standard molal entropy, volume, and heat capacity.
```{marginfigure}
-The corresponding units are °C (`T`), bar (`P`), g cm<sup>-3</sup> (`rho`), cal mol<sup>-1</sup> (`G` and `H`), cal K<sup>-1</sup> mol<sup>-1</sup> (`S`), cm<sup>3</sup> mol<sup>-1</sup> (`V`), and cal K<sup>-1</sup> mol<sup>-1</sup> (`Cp`).
+The corresponding units are °C (`T`), bar (`P`), g cm<sup>-3</sup> (`rho`), cal mol<sup>-1</sup> (`G` and `H`), cal K<sup>-1</sup> mol<sup>-1</sup> (`S` and `Cp`), and cm<sup>3</sup> mol<sup>-1</sup> (`V`).
```
A custom temperature-pressure grid can be specified.
@@ -281,13 +281,12 @@
The default units of temperature, pressure, and energy are °C, bar, and calories.
The functions <span style="color:red">`T.units()`</span>, <span style="color:red">`P.units()`</span>, and <span style="color:red">`E.units()`</span> can be used to change the units used by various functions in CHNOSZ.
-What is the Gibbs energy (J/mol) of aqueous methane at 298.15 K and 0.1 MPa?
+What is the Gibbs energy in J/mol of aqueous methane at 298.15 K and 0.1 MPa?
```{r units_methane, message=FALSE}
T.units("K")
P.units("MPa")
E.units("J")
subcrt("methane", T = 298.15, P = 0.1)$out$methane$G
-data(thermo) ## restore default settings
```
A related function, <span style="color:green">`convert()`</span>, can be used to convert given values between units.
@@ -298,8 +297,12 @@
As expected, we get the same result from both operations.
-# Properties of reactions
+Use <span style="color:red">`data(thermo)`</span> to restore the units and all other settings for CHNOSZ to their defaults:
+```{r data_thermo}
+```
+# Properties of reactions {#Properties-of-reactions}
+
## Reaction definitions
To calculate the thermodynamic properties of reactions, give the names of species, the physical states (optional), and reaction coefficients as the arguments to <span style="color:green">`subcrt()`</span>.
@@ -351,10 +354,9 @@
## Setting the basis species
-_Basis species_ are a minimal number of chemical species that represent the compositional variation in a system.
+_Basis species_ are a minimal number of chemical species that linearly combine to give the composition of any species in the system.
The basis species are similar to thermodynamic components, but can include charged species.
-You might want to use basis species to automatically balance reactions.
-A basis setting is also required for making chemical activity diagrams.
+Basis species are used in CHNOSZ to automatically balance reactions; they are also required for making chemical activity diagrams.
Let's start with an example that doesn't work:
```{r basis_singular, error=TRUE}
@@ -383,7 +385,7 @@
subcrt(c("acetate", "methane"), c(-1, 1))$reaction
```
-We can similarly consider reactions for hydrogenotrophic methanogenesis as well as acetate oxidation (no production of methane):
+We can similarly consider reactions for hydrogenotrophic methanogenesis as well as acetate oxidation (without production of methane):
```{r subcrt_methanogenesis, message=FALSE}
acetate_oxidation <- subcrt("acetate", -1)
hydrogenotrophic <- subcrt("methane", 1)
@@ -421,6 +423,7 @@
The activity of acetate and fugacity of methane, as well as temperature and pressure, are set in the call to <span style="color:green">`subcrt()`</span>:
```{r basis_mayumi, message=FALSE, results="hide"}
E.units("J")
+basis(c("CO2", "H2", "H2O", "H+"))
basis(c("CO2", "H2"), "gas")
basis(c("H2", "pH"), c(-3.92, 7.3))
```
@@ -453,7 +456,7 @@
Now we're ready to calculate and plot the affinities.
Here, we use R's `lapply()` to list the results at two values of logarithm of fugacity of CO<sub>2</sub>.
We insert an empty reaction to get a line at zero affinity.
-`do.call(rbind, Adat)` turns the list into a data frame that can be plotted with R's `matplot()`.
+R's handy `do.call()` and `rbind()` are used to turn the list into a data frame that can be plotted with R's `matplot()`.
There, we plot the negative affinities, equal to Gibbs energy, as shown in the plot of Mayumi et al. (2013).
```{r methanogenesis_plot, fig.margin=TRUE, fig.width=4.1, fig.height=4.1, small.mar=TRUE, dpi=dpi, out.width="100%", echo=FALSE, message=FALSE, fig.cap="Gibbs energies of acetate oxidation and methanogenesis (after Mayumi et al., 2013).", cache=TRUE, pngquant=pngquant, timeit=timeit}
@@ -481,18 +484,23 @@
# Using affinity
-<span style="color:green">`affinity()`</span> offers calculations of chemical affinity of formation reactions over a configurable range of T, P, and activities of basis species.
+<span style="color:green">`affinity()`</span> offers calculations of chemical affinity of formation reactions over a configurable range of *T*, *P*, and activities of basis species.
-## Species of interest
+## Species of interest {#Species-of-interest}
By *formation reaction* is meant the stoichiometric requirements for formation of one mole of any species from the basis species.
The <span style="color:red">`species()`</span> function is used to set these *species of interest*.
Let's consider the stoichiometry of some aqueous sulfur-bearing species.
-Here we use <span style="color:red">`basis()`</span> with a keyword to identify a preset basis definition.
+Here we use <span style="color:red">`basis()`</span> with a keyword to load a preset basis definition.
```{marginfigure}
-Possible keywords are `CHNOS` (including CO<sub>2</sub>, H<sub>2</sub>O, NH<sub>3</sub>, H<sub>2</sub>S, and O<sub>2</sub>), `CHNOS+` (also including H<sup>+</sup>), `CHNOSe` (including H<sup>+</sup>, and *e*<sup>-</sup> instead of O<sub>2</sub>).
+Some available keywords are `CHNOS` (including CO<sub>2</sub>, H<sub>2</sub>O, NH<sub>3</sub>, H<sub>2</sub>S, and O<sub>2</sub>), `CHNOS+` (also including H<sup>+</sup>), and `CHNOSe` (including H<sup>+</sup>, and *e*<sup>-</sup> instead of O<sub>2</sub>).
See <span style="color:blue">`?basis`</span> for more options.
```
+```{marginfigure}
+What is `SO42-`? Is it 1 S, 4 O, and 2 -ve charges, or 1 S, 42 O, and 1 -ve charge?
+The ambiguity of a digit that could belong to the coefficient for the following charge or to that for the preceding element is why formulas in CHNOSZ are written with the number of charges after the + or - symbol.
+`SO4-2` is unambiguously parsed as 1 S, 4 O and 2 -ve charges.
+```
```{r basis_CHNOSZ, results="hide"}
basis("CHNOS+")
```
@@ -500,7 +508,9 @@
species(c("H2S", "HS-", "HSO4-", "SO4-2"))
```
-Now, we can use <span style="color:green">`affinity()`</span> to calculate the affinities of the formation reactions of each of the species:
+Aqueous species are assigned default activities of 10<sup>-3</sup> (`logact` is -3).
+Now, we can use <span style="color:green">`affinity()`</span> to calculate the affinities of the formation reactions of each of the species.
+R's `unlist()` is used here simply to turn the list of values of affinity into a numeric object that can be printed in a couple of lines (note that the names correspond to `ispecies` above):
```{marginfigure}
The values returned by <span style="color:green">`affinity()`</span> are dimensionless, i.e. *A*/2.303*RT*.
```
@@ -508,15 +518,15 @@
unlist(affinity()$values)
```
-The same result (in energetic units) could be obtained using <span style="color:green">`subcrt()`</span>, but <span style="color:green">`affinity()`</span> has the advantage of being able to perform calculations on a grid of *T*, *P*, or activities of basis species.
+The same result (but expressed in units of cal/mol or J/mol) could be obtained using <span style="color:green">`subcrt()`</span>; however, <span style="color:green">`affinity()`</span> has the advantage of being able to perform calculations on a grid of *T*, *P*, or activities of basis species.
Let's choose a set of variables commonly used in aqueous speciation diagrams: Eh and pH.
To use Eh as a variable, the electron (*e*<sup>-</sup>) should be in the basis.
-To put the electron in there, we could use a different keyword (<span style="color:red">`basis("CHNOSe")`</span>), or swap oxygen out of the existing basis:
+To put the electron in there, we can use a different keyword (<span style="color:red">`basis("CHNOSe")`</span>), or swap oxygen out of the existing basis:
```{r swap_basis}
swap.basis("O2", "e-")
```
-The <span style="color:red">`swap.basis()`</span> changed the basis species and recalculated their activities, but preserved the species of interest.
+<span style="color:red">`swap.basis()`</span> changed the basis species and recalculated their activities, but preserved the species of interest.
```{marginfigure}
That is, running <span style="color:green">`affinity()`</span>`$values` again would give the same result.
```
@@ -536,7 +546,7 @@
```{r EhpH_plot, echo=2, eval=FALSE}
```
-Note that the calculation of affinity implies a non-equilibrium reference state of equal activities of species (see above).
+Note that the calculation of affinity implies a non-equilibrium reference state of equal activities of species ([see above](#Species-of-interest)).
Generally, then, <span style="color:green">`diagram()`</span> gives a *potential diagram* because it shows regions of maximum affinity.
In systems where equilibrium is attainable, it makes sense to call this a *predominance diagram*, showing regions of maximum activity.
@@ -545,17 +555,17 @@
names = c("hydrogen sulfide", "bisulfide", "bisulfate", "sulfate"),
tplot = FALSE, main = "sulfur species, 25 °C", bty = "n")
```
-The default colors for diagrams shown on the screen uses R's `heat.colors()` palette.
-Some arguments in <span style="color:green">`diagram()`</span> can be used to control the color, labels, and lines, and title (`main`).
+The default colors for diagrams shown on the screen use R's `heat.colors()` palette.
+Some arguments in <span style="color:green">`diagram()`</span> can be used to control the color, labels, and lines, and title.
The `tplot` argument turns off plot customizations used in CHNOSZ.
Additional arguments are passed to R's plotting functions; here, we use `bty` to remove the box around the plot:
```{r EhpH_plot_color, echo=TRUE, eval=FALSE}
```
-## Mosaic diagrams
+## Mosaic diagrams {#Mosaic-diagrams}
If sulfur is in the basis species, then we should consider that its speciation is sensitive to Eh and pH, as shown in the preceding diagram.
-Mosaic diagrams (or combined diagrams), which are often shown for metal oxide, sulfide, and carbonate minerals, account for speciation of the basis species.
+Mosaic diagrams (or combined diagrams), which are often shown for oxide, sulfide, and carbonate minerals, account for speciation of the basis species.
These diagrams are made by constructing individual diagrams for the possible basis species.
The individual diagrams are then combined, each one contributing to the final diagram only in the range of stability of the corresponding basis species.
@@ -574,7 +584,7 @@
Those are temperatures in Kelvin (regardless of the <span style="color:red">`T.units()`</span>); at 200 °C we should use the second phase.
Next we define the basis, and set the activities of the H<sub>2</sub>S and Cl<sup>-</sup> basis species.
-These represent the total activity of S and Cl in the system, which are distributed among the minerals and aqueous species (i.e., not the basis species).
+These represent the total activity of S and Cl in the system, which are distributed among the minerals and aqueous species.
Three minerals and the aqueous copper chloride species are included:
```{r copper_setup, echo=TRUE, results="hide"}
basis(c("Cu", "H2S", "Cl-", "H2O", "H+", "e-"))
@@ -586,8 +596,8 @@
```
We use <span style="color:green">`mosaic()`</span> to generate and combine diagrams for each candidate basis species (H<sub>2</sub>S, HS<sup>-</sup>, HSO<sub>4</sub><sup>-</sup>, or SO<sub>4</sub><sup>-2</sup>) as a function of Eh and pH.
-The key argument is `bases`, which identifies the candidate basis species (starting with the one in the current basis).
-The other arguments, like those of <span style="color:green">`affinity()`</span>, specify the ranges of the variables; `res` indicates the grid resolution to use for each variable (higher than the default of 128).
+The key argument is `bases`, which identifies the candidate basis species, starting with the one in the current basis.
+The other arguments, like those of <span style="color:green">`affinity()`</span>, specify the ranges of the variables; `res` indicates the grid resolution to use for each variable (the default is 128).
The first call to <span style="color:green">`diagram()`</span> plots the species of interest; the second adds the predominance fields of the basis species.
Finally, <span style="color:green">`water.lines()`</span> is used to add the stability limits of water at the given temperature.
@@ -596,15 +606,16 @@
res <- 200
bases <- c("H2S", "HS-", "HSO4-", "SO4-2")
m1 <- mosaic(bases, blend = TRUE, pH = c(0, 12, res), Eh=c(-1.2, 0.75, res), T=T)
-diagram(m1$A.species, lwd = 2)
+diagram(m1$A.species, lwd = 2, fill = NA)
diagram(m1$A.bases, add = TRUE, col = "blue", col.names = "blue", lty = 2)
water.lines("pH", "Eh", T = convert(T, "K"), col = "red", lwd = 2, lty = 2)
```
-The argument `blend=TRUE` is used to combine the diagrams according to the equilibrium activities of the basis species (see below).
+The argument `blend = TRUE` is used to combine the diagrams according to the equilibrium activities of the basis species by themselves ([see below](#Equilibration)).
The smooth transitions between basis species cause the appearance of curved lines on the plot.
Without that argument, the diagrams would be combined using the dominant basis species, and all of the line segments would be straight.
+<a name="mosaicfun"></a>
We have seen the effects of speciation of S in the basis species.
However, the choice of other basis species can also affect the diagram.
For instance, we can use H<sub>2</sub> or `r o2` in place of *e*<sup>-</sup>.
@@ -626,13 +637,15 @@
mosaicfun(list(O2 = c(-70, 0, res)))
```
-## *T*, *P*, activity transects
+## *T*, *P*, activity transects {#T-P-activity-transects}
-Above, we used evenly-spaced grids of *T*, *P*, and/or chemical activities of basis species; the ranges of variables were given by two or three values (minimum, maximum, and optionally resolution).
+Above, we used evenly-spaced grids of chemical activities of basis species; the ranges of variables were given by two or three values (minimum, maximum, and optionally resolution).
<span style="color:green">`affinity()`</span> can also perform calculations along a transect, i.e. a particular path along one or more variables.
A transect is calculated when there are four or more values assigned to the variable(s).
Let's use this feature to calculate affinities (negative Gibbs energies) of methanogenesis and biosynthetic reactions in a hydrothermal system.
-Some results of mixing calculations for seawater and vent fluid from the Rainbow hydrothermal field, reported by @SC10, are included in a data file in CHNOSZ:
+
+Some results of mixing calculations for seawater and vent fluid from the Rainbow hydrothermal field, calculated using EQ3/6 by @SC10, are included in a data file in CHNOSZ.
+Reading the file with R's `read.csv()`, we set `check.names = FALSE` to preserve the `NH4+` column name (which is not a syntactically valid variable name):
```{r rainbow_data}
file <- system.file("extdata/cpetc/SC10_Rainbow.csv", package = "CHNOSZ")
rb <- read.csv(file, check.names = FALSE)
@@ -651,8 +664,8 @@
"methane", "leucine", "tryptophan", "n-nonanoic acid"), -6)
```
-Now we can calculate affinity along the transect of changing temperature and activities of five basis species.
-Each variable is given as a named argument; the name for `NH4+` must be quoted.
+Now we can calculate affinities along the transect of changing temperature and activities of five basis species.
+Each variable is given as a named argument; `NH4+` must be quoted.
```{marginfigure}
A shorter expression would use R's `do.call()` to construct the function call: `do.call(`<span style="color:green">`affinity`</span>`, as.list(rb))`.
```
@@ -661,7 +674,7 @@
That conversion requires temperature in Kelvin, which is obtained by conversion from °C.
We finish with a negation (affinity is negative Gibbs energy) and scaling from cal to kcal.
```
-Using <span style="color:green">`convert()`</span>, we also convert the result from dimensionless values (*A*/2.303*RT*) to kcal/mol.
+Using R's `lapply()` to run <span style="color:green">`convert()`</span> for each species, we convert the affinity from dimensionless values (*A*/2.303*RT*) to cal/mol.
```{r rainbow_affinity, message=FALSE}
a <- affinity(T = rb$T, CO2 = rb$CO2, H2 = rb$H2,
`NH4+` = rb$`NH4+`, H2S = rb$H2S, pH = rb$pH)
@@ -678,25 +691,25 @@
Finally, we use <span style="color:green">`diagram()`</span> to plot the results.
Although only temperature is shown on the *x* axis, pH and the activities of CO<sub>2</sub>, H<sub>2</sub>, NH<sub>4</sub><sup>+</sup>, and H<sub>2</sub>S are also varied according to the data in `rb`.
By default, <span style="color:green">`diagram()`</span> attempts to scale the affinities by dividing by the reaction coefficients of a shared basis species (in this case, CO<sub>2</sub>).
-To override that behavior, we set `balance = 1` to plot the affinities of the formation reactions as written (per mole of the species being formed).
+To override that behavior, we set `balance = 1` to plot the affinities of the formation reactions as written (per mole of the product species).
```{r rainbow_diagram, eval=FALSE}
```
-## Buffers (experimental feature)
+## Buffers
-There is one other feature of <span style="color:green">`affinity()`</span> to mention here.
+There is one other feature of <span style="color:green">`affinity()`</span> that can be mentioned here.
Can we go the other direction: calculate the activities of basis species from the activities of the species of interest?
This question relates to the concept of chemical activity buffers.
In CHNOSZ there are two ways to perform buffer calculations:
-1. Assign the name of a buffer (listed in `thermo$buffer`) to the basis species
-* more versatile (multiple activities can be buffered, e.g. both S<sub>2</sub> and O<sub>2</sub> by pyrite-pyrrhotite-magnetite)
-* the buffers are active in calculations of affinity of other species
-* use <span style="color:red">`mod.buffer()`</span> to change or add buffers in `thermo$buffer`
-* <span style="color:blue">`demo(buffer)`</span> uses it for mineral buffers (solid lines)
-2. Use the `what` argument of <span style="color:green">`diagram()`</span> to solve for the activity of the indicated basis species
-* more convenient (the buffers come from the currently defined species of interest), but only a single basis species can be buffered, and it's not used in the calculation of affinity
-* <span style="color:blue">`demo(buffer)`</span> uses it for aqueous organic species as buffers (dotted and dashed lines)
+1. Assign the name of a buffer (listed in `thermo$buffer`) to the basis species:
+* more versatile (multiple activities can be buffered, e.g. both S<sub>2</sub> and O<sub>2</sub> by pyrite-pyrrhotite-magnetite);
+* the buffers are active in calculations of affinity of other species;
+* use <span style="color:red">`mod.buffer()`</span> to change or add buffers in `thermo$buffer`;
+* <span style="color:blue">`demo(buffer)`</span> uses it for mineral buffers (solid lines).
+2. Use the `what` argument of <span style="color:green">`diagram()`</span> to solve for the activity of the indicated basis species:
+* more convenient (the buffers come from the currently defined species of interest), but only a single basis species can be buffered, and it's not used in the calculation of affinity;
+* <span style="color:blue">`demo(buffer)`</span> uses it for aqueous organic species as buffers (dotted and dashed lines).
As an example of method 1, let's look at the pyrite-pyrrhotite-magnetite (PPM) buffer at 300 °C.
```{marginfigure}
@@ -712,11 +725,13 @@
```{marginfigure}
The affinity of formation of pyrite happens to be zero because it is identical to one of the selected basis species.
```
-```{r PPM_affinity, message=FALSE}
+```{r PPM_affinity, message=FALSE, echo=1}
unlist(affinity(T = 300, P = 100)$values)
+## 2031 1999 2036
+## 0 0 0
```
-We use <span style="color:red">`mod.buffer()`</span> to choose the `cr2` phase of pyrrhotite, which is stable at this temperature.
+We use <span style="color:red">`mod.buffer()`</span> to choose the `cr2` phase of pyrrhotite, which is stable at this temperature ([see above](#Mosaic-diagrams) for how to get this information for minerals with phase transitions).
Then, we set up H<sub>2</sub>S and `r o2` to be buffered by PPM, and inspect their buffered activities:
```{r PPM_setup, results="hide"}
mod.buffer("PPM", "pyrrhotite", "cr2")
@@ -726,38 +741,39 @@
unlist(affinity(T = 300, P = 100, return.buffer = TRUE)[1:3])
```
-```{r demo_buffer, fig.margin=TRUE, fig.width=4, fig.height=4, small.mar=TRUE, dpi=dpi, out.width="100%", message=FALSE, echo=FALSE, cache=TRUE, fig.cap="Values of log<i>f</i><sub>H<sub>2</sub></sub> corresponding to mineral buffers or to given activities of aqueous species.", pngquant=pngquant, timeit=timeit}
+<!-- put demo(buffer) here for appealing placement on page -->
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/chnosz -r 151
More information about the CHNOSZ-commits
mailing list