[CHNOSZ-commits] r158 - in pkg/CHNOSZ: . R demo inst man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri Feb 17 14:49:46 CET 2017
Author: jedick
Date: 2017-02-17 14:49:45 +0100 (Fri, 17 Feb 2017)
New Revision: 158
Added:
pkg/CHNOSZ/demo/activity_ratios.R
Modified:
pkg/CHNOSZ/DESCRIPTION
pkg/CHNOSZ/R/diagram.R
pkg/CHNOSZ/R/examples.R
pkg/CHNOSZ/R/util.expression.R
pkg/CHNOSZ/demo/00Index
pkg/CHNOSZ/inst/NEWS
pkg/CHNOSZ/man/diagram.Rd
pkg/CHNOSZ/man/examples.Rd
pkg/CHNOSZ/man/protein.Rd
pkg/CHNOSZ/man/util.expression.Rd
Log:
Add demo activity_ratios.R
Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/DESCRIPTION 2017-02-17 13:49:45 UTC (rev 158)
@@ -1,6 +1,6 @@
Date: 2017-02-17
Package: CHNOSZ
-Version: 1.0.8-47
+Version: 1.0.8-48
Title: Chemical Thermodynamics and Activity Diagrams
Author: Jeffrey Dick
Maintainer: Jeffrey Dick <j3ffdick at gmail.com>
Modified: pkg/CHNOSZ/R/diagram.R
===================================================================
--- pkg/CHNOSZ/R/diagram.R 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/R/diagram.R 2017-02-17 13:49:45 UTC (rev 158)
@@ -12,7 +12,7 @@
what="loga.equil", alpha=FALSE, normalize=FALSE, as.residue=FALSE, balance=NULL,
groups=as.list(1:length(eout$values)), xrange=NULL,
# plot dimensions
- mar=NULL, yline=par("mgp")[1]+0.7, side=1:4,
+ mar=NULL, yline=par("mgp")[1]+0.5, side=1:4,
# axes
ylog=TRUE, xlim=NULL, ylim=NULL, xlab=NULL, ylab=NULL,
# sizes
Modified: pkg/CHNOSZ/R/examples.R
===================================================================
--- pkg/CHNOSZ/R/examples.R 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/R/examples.R 2017-02-17 13:49:45 UTC (rev 158)
@@ -5,7 +5,7 @@
examples <- function(do.png=FALSE) {
# run all the examples in CHNOSZ documentation
.ptime <- proc.time()
- topics <- c("CHNOSZ-package", "thermo", "sideeffects", "examples",
+ topics <- c("thermo", "sideeffects", "examples",
"util.args", "util.array", "util.blast", "util.character",
"util.data", "util.expression", "util.fasta", "util.formula", "util.matrix",
"util.misc", "util.program",
@@ -31,7 +31,7 @@
demos <- function(which=c("sources", "protein.equil", "add.obigt", "affinity", "NaCl", "density",
"ORP", "revisit", "findit", "ionize", "buffer", "protbuff", "yeastgfp", "mosaic",
- "copper", "solubility", "wjd", "dehydration", "bugstab"), to.file=FALSE) {
+ "copper", "solubility", "wjd", "dehydration", "bugstab", "activity_ratios"), to.file=FALSE) {
# run one or more demos from CHNOSZ with ask=FALSE, and return the value of the last one
for(i in 1:length(which)) {
# say something so the user sees where we are
Modified: pkg/CHNOSZ/R/util.expression.R
===================================================================
--- pkg/CHNOSZ/R/util.expression.R 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/R/util.expression.R 2017-02-17 13:49:45 UTC (rev 158)
@@ -21,9 +21,10 @@
# recover the coefficient
if(elements[i]==1) coeff <- "" else coeff <- elements[i]
# append the coefficient
- # subscripts within subscripts (log) are too small
- if(log != "") expr <- substitute(a*b, list(a=expr, b=coeff))
- else expr <- substitute(a[b], list(a=expr, b=coeff))
+ ## subscripts within subscripts (log) are too small
+ #if(log != "") expr <- substitute(a*b, list(a=expr, b=coeff))
+ #else expr <- substitute(a[b], list(a=expr, b=coeff))
+ expr <- substitute(a[b], list(a=expr, b=coeff))
} else {
# for charged species, don't show "Z" but do show e.g. "+2"
coeff <- elements[i]
@@ -49,7 +50,7 @@
if(log %in% c("aq", "cr", "liq", "cr1", "cr2", "cr3", "cr4")) acity <- "a"
else if(log %in% c("g", "gas")) acity <- "f"
else stop(paste("'", log, "' is not a recognized state", sep=""))
- logacity <- substitute(log*italic(a), list(a=acity))
+ logacity <- substitute(log~italic(a), list(a=acity))
expr <- substitute(a[b], list(a=logacity, b=expr))
# write a value if given
if(!is.null(value)) {
@@ -248,3 +249,26 @@
return(desc)
}
+# make formatted text for activity ratio 20170217
+ratlab <- function(ion="K+") {
+ # the charge
+ Z <- makeup(ion)["Z"]
+ # the text for the exponent on aH+
+ exp.H <- as.character(Z)
+ # the expression for the ion and H+
+ expr.ion <- expr.species(ion)
+ expr.H <- expr.species("H+")
+ # the final expression
+ if(exp.H=="1") substitute(log~(italic(a)[expr.ion] / italic(a)[expr.H]), list(expr.ion=expr.ion, expr.H=expr.H))
+ else substitute(log~(italic(a)[expr.ion] / italic(a)[expr.H]^exp.H), list(expr.ion=expr.ion, expr.H=expr.H, exp.H=exp.H))
+}
+
+# make formatted text for thermodynamic system 20170217
+syslab <- function(system = c("K2O", "Al2O3", "SiO2", "H2O")) {
+ for(i in seq_along(system)) {
+ expr <- expr.species(system[i])
+ # use en dash here
+ if(i==1) lab <- expr else lab <- substitute(a*"\u2013"*b, list(a=lab, b=expr))
+ }
+ lab
+}
Modified: pkg/CHNOSZ/demo/00Index
===================================================================
--- pkg/CHNOSZ/demo/00Index 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/demo/00Index 2017-02-17 13:49:45 UTC (rev 158)
@@ -17,3 +17,4 @@
wjd run.wjd() with proteins: cell periphery of yeast
dehydration log K of dehydration reactions; SVG file contains tooltips and links
bugstab formation potential of microbial proteins in colorectal cancer
+activity_ratios mineral stability plots with activity ratios on the axes
Added: pkg/CHNOSZ/demo/activity_ratios.R
===================================================================
--- pkg/CHNOSZ/demo/activity_ratios.R (rev 0)
+++ pkg/CHNOSZ/demo/activity_ratios.R 2017-02-17 13:49:45 UTC (rev 158)
@@ -0,0 +1,61 @@
+## Equilibrium activity diagrams for minerals using activity ratios as variables
+## These are made with pH = 0 (activity of H+ = 1), so (activity of the ion) is equal to
+## (activity of the ion) / [(activity of H+) ^ (charge of the ion)]
+
+par(mfrow = c(2, 2))
+res <- 200
+fill <- "heat"
+
+## K2O-Al2O3-SiO2-H2O, 25 degree C, 1 bar
+## Steinmann et al., 1994 (http://ccm.geoscienceworld.org/content/42/2/197)
+## Garrels and Christ, p. 361 (http://www.worldcat.org/oclc/517586)
+## https://wwwbrr.cr.usgs.gov/projects/GWC_coupled/phreeqc/html/final-75.html
+# use add.obigt to get H4SiO4
+add.obigt()
+basis(c("Al+3", "H4SiO4", "K+", "H2O", "H+", "O2"))
+species(c("gibbsite", "kaolinite", "pyrophyllite", "k-feldspar", "muscovite"))
+a <- affinity(H4SiO4 = c(-6, 0, res), `K+` = c(0, 8, res))
+diagram(a, ylab = ratlab("K+"), fill = fill)
+title(main = syslab(c("K2O", "Al2O3", "SiO2", "H2O")))
+legend("topright", describe.property(c("T", "P"), c(25, 1)), bty = "n")
+# reset database
+data(thermo)
+
+## H2O-CaO-MgO-SiO2 at 300 degree C and 1000 bar
+## Helgeson et al., 1969, p. 136 (http://www.worldcat.org/oclc/902423149)
+## Bowers et al., 1984, p. 246 (http://www.worldcat.org/oclc/224591948)
+basis(c("H2O", "Ca+2", "Mg+2", "SiO2", "O2", "H+"))
+species(c("quartz", "talc", "chrysotile", "forsterite", "monticellite",
+ "merwinite", "wollastonite", "diopside", "tremolite"))
+# calculate the chemical affinities of formation reactions
+a <- affinity("Mg+2" = c(4, 9, res), "Ca+2" = c(5, 14, res), T = 300, P = 1000)
+diagram(a, xlab = ratlab("Mg+2"), ylab = ratlab("Ca+2"), fill = fill)
+title(main = syslab(c("H2O", "CaO", "MgO", "SiO2")))
+legend("bottomright", describe.property(c("T", "P"), c(300, 1000)), bty = "n")
+# note: Bowers et al. use more complicated variables
+# (involving the hydration numbers of H2O and the ion)
+# with accordingly different axis ranges
+
+## MgO-CaO-SiO2-H2O at 300 degree C and Psat
+## Russell et al., 2010 (http://dx.doi.org/10.1111/j.1472-4669.2010.00249.x)
+basis(c("Mg+2", "Ca+2", "SiO2", "H2O", "O2", "H+"))
+species(c("brucite", "chrysotile", "talc", "tremolite", "diopside", "akermanite"))
+a <- affinity(SiO2 = c(-10, 0, res), `Ca+2` = c(0, 20, res), T = 300)
+diagram(a, ylab = ratlab("Ca+2"), fill = fill)
+title(main = syslab(c("MgO", "CaO", "SiO2", "H2O")))
+legend("topright", describe.property(c("T", "P"), c(300, 85.84)), bty = "n")
+
+## CaO-MgO-SiO2-H2O and
+## CaO-Al2O3-MgO-SiO2-H2O at 300 degree C and 500 bar
+## Bach and Klein, 2009 (http://dx.doi.org/10.1016/j.lithos.2008.10.022)
+basis(c("Ca+2", "Al+3", "Mg+2", "SiO2", "H2O", "O2", "H+"))
+species(c("clinochlore,14a", "clinozoisite", "prehnite", "grossular"))
+a <- affinity(SiO2 = c(-5, 0, res), `Ca+2` = c(6, 11, res), T = 300, P = 500)
+diagram(a, ylab = ratlab("Ca+2"), balance = "Al+3", fill = fill)
+# (Hmmm... where is clinochlore? it doesn't appear on our diagram)
+species(delete = TRUE)
+species(c("brucite", "chrysotile", "talc", "tremolite", "diopside"))
+a <- affinity(SiO2 = c(-5, 0, res), `Ca+2` = c(6, 11, res), T = 300, P = 500)
+diagram(a, add = TRUE, col = "blue", col.names = "blue")
+title(main = syslab(c("CaO", "Al2O3", "MgO", "SiO2", "H2O")))
+legend("topright", describe.property(c("T", "P"), c(300, 500)), bty = "n")
Modified: pkg/CHNOSZ/inst/NEWS
===================================================================
--- pkg/CHNOSZ/inst/NEWS 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/inst/NEWS 2017-02-17 13:49:45 UTC (rev 158)
@@ -1,31 +1,47 @@
-CHANGES IN CHNOSZ 1.0.8-47 (2017-02-17)
+CHANGES IN CHNOSZ 1.0.8-48 (2017-02-17)
---------------------------------------
DOCUMENTATION:
-- Replace anintro.Rnw (Sweave) with rewritten anintro.Rmd (knitr, Tufte
- style).
+- Replace anintro.Rnw (using Sweave) with rewritten anintro.Rmd (using
+ knitr with Tufte style).
-- Add eos-regress.Rmd.
+- Add eos-regress.Rmd (Regressing thermodynamic data).
-- 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 demos: bugstab.R (potential diagrams for microbial proteins in
+ colorectal cancer), activity_ratios.R (mineral stability diagrams
+ with activity ratios on the axes).
-NEW DATA AND FEATURES:
+- Converted demos: move examples from help pages to protein.equil.R,
+ add.obigt.R, and affinity.R.
-- OBIGT.csv: Add data for citric acid cycle metabolites from Canovas and
- Shock, 2016.
+- Updated and removed demos: update buffer.R (logfH2 set by mineral
+ buffers and aqueous species) and rename old buffer.R to protbuff.R
+ (proteins as buffers); remove nucleobase.R.
-- extdata: Add fasta/rubisco.fasta and cpetc/*.csv files SS98_Fig5a,
- SS98_Fig5b, PM90 (all used by anintro.Rmd).
+CHANGES TO diagram():
+- For the `groups` argument, activities are multiplied by the balance
+ coefficients before the summation.
+
+- Draw smoother boundary lines (with contour(); `lty` is usable), at the
+ expense of not being precisely aligned with the grid used for plotting
+ colors. For the old behavior, use `dotted=0`.
+
+- Default of `legend.x` changed to NA (instead of making a legend, put
+ labels near the lines).
+
+- Field labels: apply formatting (subscripts and superscripts) to names
+ that parse as chemical formulas; use `format.names=FALSE` to disable.
+
+- Chemical formulas that are part of the axis labels (log activity or
+ fugacity) now have formatting applied.
+
+NEW FEATURES:
+
- Add ZC.col() for generating a red-grey-blue color scale from
numeric values.
-- For the `groups` argument in diagram(), activities are multiplied
- by the balance coefficients before the summation.
-
- Add "QEC" as a keyword for preset species in basis() (glutamine,
glutamic acid, cysteine, H2O, oxygen).
@@ -36,13 +52,6 @@
(accession, name, organism) into columns of the output data frame
(protein, abbrv, organism).
-- diagram() draws smoother boundary lines (with contour(); `lty` is
- usable), at the expense of not being precisely aligned with the grid
- used for plotting colors. For the old behavior, use `dotted=0`.
-
-- diagram() applies formatting (subscripts and superscripts) to names
- that parse as chemical formulas; use `format.names=FALSE` to disable.
-
- Add `...` argument to thermo.plot.new() (additional arguments for
par()).
@@ -50,7 +59,7 @@
marks and lines; this is used in diagram() to redraw the axes on
filled diagrams.
-- Add `...` argument to read.aaa() (additional arguments for
+- Add `...` argument to read.aa() (additional arguments for
read.csv()).
- seq2aa() removes newlines and whitespace before counting the
@@ -59,9 +68,17 @@
- read.fasta(): change argument `i` to `iseq`; this is used to select
particular sequences to read from the file.
-- diagram(): default of `legend.x` changed to NA (instead of making a
- legend, put labels near the lines).
+- Add syslab() and ratlab() to make formatted expressions for
+ sets of thermodynamic components and activity ratios.
+NEW DATA:
+
+- OBIGT.csv: Add data for citric acid cycle metabolites from Canovas and
+ Shock, 2016.
+
+- extdata: Add fasta/rubisco.fasta and cpetc/*.csv files SS98_Fig5a,
+ SS98_Fig5b, PM90 (all used by anintro.Rmd).
+
CLEANUP AND BUG FIXES:
- subcrt() returns `loggam` using the common logarithm; add
Modified: pkg/CHNOSZ/man/diagram.Rd
===================================================================
--- pkg/CHNOSZ/man/diagram.Rd 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/man/diagram.Rd 2017-02-17 13:49:45 UTC (rev 158)
@@ -10,7 +10,7 @@
\usage{
diagram(eout, what = "loga.equil", alpha = FALSE, normalize = FALSE,
as.residue = FALSE, balance=NULL, groups=as.list(1:length(eout$values)),
- xrange=NULL, mar=NULL, yline=par("mgp")[1]+0.7, side=1:4,
+ xrange=NULL, mar=NULL, yline=par("mgp")[1]+0.5, side=1:4,
ylog=TRUE, xlim=NULL, ylim=NULL, xlab=NULL, ylab=NULL,
cex=par("cex"), cex.names=1, cex.axis=par("cex"),
lty=NULL, lwd=par("lwd"), dotted=NULL,
@@ -199,28 +199,6 @@
title(main=paste("Fe-S-O, 200 degrees C, 1 bar",
"After Helgeson, 1970", sep="\n"))
-## Mineral equilibrium activity diagram
-## After Bowers et al., 1984, p. 246
-# Consider the system denoted by BJH84 as
-# HCl-H2O-CaO-CO2-MgO-(SiO2) at 300 deg C and 1000 bar
-# HCl, CO2, O2 have zero stoichiometric coeffs in the species
-# Ca+2, Mg+2 are going to be on the diagram
-# SiO2 will be conserved (the immobile component)
-# (try changing any of the 999's, the diagram will be unaffected)
-basis(c("HCl","H2O","Ca+2","CO2","Mg+2","SiO2","O2","H+"),
- c(999,0,999,999,999,999,999,-7))
-# we have to tell CHNOSZ which species to include; there are
-# others that could be in this system
-species(c("quartz","talc","forsterite","tremolite","diopside",
- "wollastonite","monticellite","merwinite"))
-# calculate the chemical affinities of formation reactions
-a <- affinity("Mg+2"=c(-12,-4),"Ca+2"=c(-8,0),T=300,P=1000)
-diagram(a)
-title(main=paste("HCl-H2O-CaO-CO2-MgO-(SiO2) at 300 degrees C,",
- "1000 bar and pH=7. After Bowers et al., 1984", sep="\n"))
-# note: Bowers et al. use different variables (e.g. (a_Ca+2)/(a_H+)^2),
-# so the plot here is not an exact reproduction
-
## Temperature-Pressure: kayanite-sillimanite-andalusite
# cf. Fig. 49 of Helgeson et al., 1978
# this is a system of one component (Al2SiO5), but we need the same
@@ -250,8 +228,6 @@
Aksu, S. and Doyle, F. M. (2001) Electrochemistry of copper in aqueous glycine solutions. \emph{J. Electrochem. Soc.} \bold{148}, B51--B57. \url{http://dx.doi.org/10.1149/1.1344532}
- Bowers, T. S., Jackson, K. J. and Helgeson, H. C., 1984. \emph{Equilibrium Activity Diagrams for Coexisting Minerals and Aqueous Solutions at Pressures and Temperatures to 5 kb and 600} \eqn{^{\circ}}{degrees }\emph{C}. Springer-Verlag, Heidelberg, 397 p. \url{http://www.worldcat.org/oclc/224591948}
-
Helgeson, H. C. (1970) A chemical and thermodynamic model of ore deposition in hydrothermal systems. \emph{Mineral. Soc. Amer. Spec. Pap.} \bold{3}, 155--186. \url{http://www.worldcat.org/oclc/583263}
Helgeson, H. C., Delany, J. M., Nesbitt, H. W. and Bird, D. K. (1978) Summary and critique of the thermodynamic properties of rock-forming minerals. \emph{Am. J. Sci.} \bold{278-A}, 1--229. \url{http://www.worldcat.org/oclc/13594862}
Modified: pkg/CHNOSZ/man/examples.Rd
===================================================================
--- pkg/CHNOSZ/man/examples.Rd 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/man/examples.Rd 2017-02-17 13:49:45 UTC (rev 158)
@@ -16,7 +16,8 @@
examples(do.png = FALSE)
demos(which = c("sources", "protein.equil", "add.obigt", "affinity", "NaCl", "density",
"ORP", "revisit", "findit", "ionize", "buffer", "protbuff", "yeastgfp", "mosaic",
- "copper", "solubility", "wjd", "dehydration", "bugstab"), to.file=FALSE)
+ "copper", "solubility", "wjd", "dehydration", "bugstab", "activity_ratios"),
+ to.file=FALSE)
}
\details{
@@ -50,6 +51,7 @@
\code{wjd} \tab run.wjd() with proteins: cell periphery of yeast \cr
\code{dehydration} \tab log K of dehydration reactions; SVG file contains tooltips and links \cr
\code{bugstab} \tab formation potential of microbial proteins in colorectal cancer (Dick, 2016) \cr
+ \code{activity_ratios} \tab mineral stability plots with activity ratios on the axes \cr
}
}
Modified: pkg/CHNOSZ/man/protein.Rd
===================================================================
--- pkg/CHNOSZ/man/protein.Rd 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/man/protein.Rd 2017-02-17 13:49:45 UTC (rev 158)
@@ -58,7 +58,7 @@
a <- affinity(O2=c(-100, -65))
# try normalize=FALSE to make Fig. 5a in the paper
e <- equilibrate(a, normalize=TRUE)
-d <- diagram(e, ylim=c(-5, -1), names=organisms)
+d <- diagram(e, ylim=c(-5, -1), names=organisms, format.names=FALSE)
# add water stability line
abline(v=-83.1, lty=2)
title(main="Surface-layer proteins, after Dick, 2008")
Modified: pkg/CHNOSZ/man/util.expression.Rd
===================================================================
--- pkg/CHNOSZ/man/util.expression.Rd 2017-02-17 02:51:34 UTC (rev 157)
+++ pkg/CHNOSZ/man/util.expression.Rd 2017-02-17 13:49:45 UTC (rev 158)
@@ -7,6 +7,8 @@
\alias{describe.basis}
\alias{describe.property}
\alias{describe.reaction}
+\alias{syslab}
+\alias{ratlab}
\title{Functions to Express Chemical Formulas and Properties}
\description{Generate expressions suitable for axis labels and plot legends describing chemical species, properties and reactions.}
@@ -20,6 +22,8 @@
describe.property(property, value, digits = 1, oneline = FALSE,
ret.val = FALSE)
describe.reaction(reaction, iname = numeric(), states = NULL)
+ syslab(system = c("K2O", "Al2O3", "SiO2", "H2O"))
+ ratlab(ion = "K+")
}
\arguments{
@@ -40,6 +44,8 @@
\item{reaction}{data frame, definition of reaction}
\item{iname}{numeric, show names instead of formulas for these species}
\item{states}{character, if \samp{all}, show states for all species}
+ \item{system}{character, thermodynamic components}
+ \item{ion}{character, an ion}
}
\details{
@@ -85,8 +91,16 @@
\code{describe.reaction} makes an expression summarizing a chemical reaction. The \code{reaction} data frame can be generated using \code{\link{subcrt}}. Based on the sign of their reaction coefficients, species are placed on the reactant (left) or product (right) side of the reaction, where the species with their coefficients are separated by plus signs; the two sides of the reaction are separated by an equals sign. Coefficients equal to 1 are not shown. Chemical formulas of species include a designation of physical state if \code{states} is \samp{all}. Names of species (as provided in \code{reaction}) are shown instead of chemical formulas for the species identified by \code{iname}.
+\code{syslab} formats the given thermodynamic components (using \code{expr.species}) and adds intervening en dashes.
+
+\code{ratlab} produces a expression for the activity ratio, viz. (activity of the ion) / [(activity of H+) ^ (charge of the ion)].
+
}
+\seealso{
+\code{demo(activity_ratio)} for examples of \code{syslab} and \code{ratlab}.
+}
+
\examples{\dontshow{data(thermo)}
## show descriptions of species and properties on a plot
plot(0, 0, xlim=c(1,5), ylim=c(1,5), xlab="function", ylab="example")
More information about the CHNOSZ-commits
mailing list