[CHNOSZ-commits] r350 - in pkg/CHNOSZ: . R demo inst man tests/testthat
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Nov 11 05:55:33 CET 2018
Author: jedick
Date: 2018-11-11 05:55:31 +0100 (Sun, 11 Nov 2018)
New Revision: 350
Modified:
pkg/CHNOSZ/DESCRIPTION
pkg/CHNOSZ/R/diagram.R
pkg/CHNOSZ/R/nonideal.R
pkg/CHNOSZ/R/util.expression.R
pkg/CHNOSZ/demo/gold.R
pkg/CHNOSZ/inst/NEWS
pkg/CHNOSZ/man/berman.Rd
pkg/CHNOSZ/man/buffer.Rd
pkg/CHNOSZ/man/diagram.Rd
pkg/CHNOSZ/man/util.expression.Rd
pkg/CHNOSZ/tests/testthat/test-util.R
Log:
expr.species(): reorganize arguments
Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/DESCRIPTION 2018-11-11 04:55:31 UTC (rev 350)
@@ -1,6 +1,6 @@
Date: 2018-11-11
Package: CHNOSZ
-Version: 1.1.3-57
+Version: 1.1.3-58
Title: Thermodynamic Calculations and Diagrams for Geo(bio)chemistry
Authors at R: c(
person("Jeffrey", "Dick", , "j3ffdick at gmail.com", role = c("aut", "cre"),
Modified: pkg/CHNOSZ/R/diagram.R
===================================================================
--- pkg/CHNOSZ/R/diagram.R 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/R/diagram.R 2018-11-11 04:55:31 UTC (rev 350)
@@ -105,7 +105,7 @@
}
## use molality instead of activity if the affinity calculation include ionic strength 20171101
- use.molality <- "IS" %in% names(eout)
+ molality <- "IS" %in% names(eout)
## when can normalize and as.residue be used
if(normalize | as.residue) {
@@ -277,7 +277,7 @@
### 0-D diagram - bar graph of properties of species or reactions
# plot setup
- if(missing(ylab)) ylab <- axis.label(plotvar, units="", use.molality=use.molality)
+ if(missing(ylab)) ylab <- axis.label(plotvar, units="", molality=molality)
barplot(unlist(plotvals), names.arg=names, ylab=ylab, cex.names=cex.names, col=col, ...)
if(!is.null(main)) title(main=main)
@@ -288,8 +288,8 @@
if(missing(xlim)) xlim <- range(xvalues) # TODO: this is backward if the vals are not increasing
# initialize the plot
if(!add) {
- if(missing(xlab)) xlab <- axis.label(eout$vars[1], basis=eout$basis, use.molality=use.molality)
- if(missing(ylab)) ylab <- axis.label(plotvar, units="", use.molality=use.molality)
+ if(missing(xlab)) xlab <- axis.label(eout$vars[1], basis=eout$basis, molality=molality)
+ if(missing(ylab)) ylab <- axis.label(plotvar, units="", molality=molality)
# to get range for y-axis, use only those points that are in the xrange
if(is.null(ylim)) {
isx <- xvalues >= min(xlim) & xvalues <= max(xlim)
@@ -563,8 +563,8 @@
ylim <- c(ys[1], tail(ys, 1))
# initialize the plot
if(!add) {
- if(is.null(xlab)) xlab <- axis.label(eout$vars[1], basis=eout$basis, use.molality=use.molality)
- if(is.null(ylab)) ylab <- axis.label(eout$vars[2], basis=eout$basis, use.molality=use.molality)
+ if(is.null(xlab)) xlab <- axis.label(eout$vars[1], basis=eout$basis, molality=molality)
+ if(is.null(ylab)) ylab <- axis.label(eout$vars[2], basis=eout$basis, molality=molality)
if(tplot) thermo.plot.new(xlim=xlim, ylim=ylim, xlab=xlab, ylab=ylab,
cex=cex, cex.axis=cex.axis, mar=mar, yline=yline, side=side, ...)
else plot(0, 0, type="n", xlim=xlim, ylim=ylim, xlab=xlab, ylab=ylab, ...)
Modified: pkg/CHNOSZ/R/nonideal.R
===================================================================
--- pkg/CHNOSZ/R/nonideal.R 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/R/nonideal.R 2018-11-11 04:55:31 UTC (rev 350)
@@ -93,8 +93,8 @@
Z[i] <- thisZ
}
# get species formulas to assign acirc 20181105
+ formula <- get("thermo")$obigt$formula[species]
if(grepl("Bdot", method)) {
- formula <- get("thermo")$obigt$formula[species]
# "ion size paramter" taken from UT_SIZES.REF of HCh package (Shvarov and Bastrakov, 1999),
# based on Table 2.7 of Garrels and Christ, 1965
acircdat <- c("Rb+"=2.5, "Cs+"=2.5, "NH4+"=2.5, "Tl+"=2.5, "Ag+"=2.5,
@@ -109,10 +109,10 @@
"Th+4"=11, "Zr+4"=11, "Ce+4"=11, "Sn+4"=11)
acirc <- as.numeric(acircdat[formula])
acirc[is.na(acirc)] <- 4.5
- # make a message
- nZ <- sum(Z!=0)
- if(nZ > 1) message("nonideal: using ", paste(acirc[Z!=0], collapse=" "), " for ion size parameters of ", paste(formula[Z!=0], collapse=" "))
- else if(nZ==1) message("nonideal: using ", acirc[Z!=0], " for ion size parameter of ", formula[Z!=0])
+ ## make a message
+ #nZ <- sum(Z!=0)
+ #if(nZ > 1) message("nonideal: using ", paste(acirc[Z!=0], collapse=" "), " for ion size parameters of ", paste(formula[Z!=0], collapse=" "))
+ #else if(nZ==1) message("nonideal: using ", acirc[Z!=0], " for ion size parameter of ", formula[Z!=0])
# use correct units (cm) for ion size parameter
acirc <- acirc * 10^-8
} else if(grepl("bgamma", method)) {
@@ -128,7 +128,7 @@
iH <- info("H+")
ie <- info("e-")
speciesprops <- as.list(speciesprops)
- ncharged <- nneutral <- 0
+ icharged <- ineutral <- logical(length(species))
for(i in 1:length(species)) {
myprops <- speciesprops[[i]]
# to keep unit activity coefficients of the proton and electron
@@ -172,11 +172,11 @@
}
# save the calculated properties and increment progress counters
speciesprops[[i]] <- myprops
- ncharged <- ncharged + sum(didcharged)
- nneutral <- nneutral + sum(didneutral)
+ if(didcharged) icharged[i] <- TRUE
+ if(didneutral) ineutral[i] <- TRUE
}
- if(ncharged > 0) message("nonideal: calculations for ", ncharged, " charged species (", mettext(method), ")")
- if(nneutral > 0) message("nonideal: calculations for ", nneutral, " neutral species (Setchenow equation)")
+ if(sum(icharged) > 0) message("nonideal: calculations for ", paste(formula[icharged], collapse=", "), " (", mettext(method), ")")
+ if(sum(ineutral) > 0) message("nonideal: calculations for ", paste(formula[ineutral], collapse=", "), " (Setchenow equation)")
return(speciesprops)
}
Modified: pkg/CHNOSZ/R/util.expression.R
===================================================================
--- pkg/CHNOSZ/R/util.expression.R 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/R/util.expression.R 2018-11-11 04:55:31 UTC (rev 350)
@@ -5,7 +5,7 @@
## if this file is interactively sourced, the following are also needed to provide unexported functions:
#source("util.character.R")
-expr.species <- function(species, state="", log="", value=NULL, use.makeup=FALSE, use.molality=FALSE) {
+expr.species <- function(species, state = "aq", value = NULL, log = FALSE, molality = FALSE, use.state = FALSE, use.makeup = FALSE) {
# make plotting expressions for chemical formulas
# that include subscripts, superscripts (if charged)
# and optionally designations of states +/- loga or logf prefix
@@ -31,9 +31,6 @@
# 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))
expr <- substitute(a[b], list(a=expr, b=coeff))
} else {
# for charged species, don't show "Z" but do show e.g. "+2"
@@ -46,36 +43,30 @@
}
}
}
- # write a designation of physical state
- ## deprecated 20181101
- ## use the state given in log if it's a gas or neutral aqueous species
- #if(log %in% c("g", "gas")) state <- "g"
- #else if(!"Z" %in% names(elements) & !missing(log)) state <- log
- if(state != "") {
- # subscript it if we're not in a log expression
- if(log != "") expr <- substitute(a*group('(',italic(b),')'),list(a=expr, b=state))
- else expr <- substitute(a[group('(',italic(b),')')],list(a=expr, b=state))
+ # write the physical state
+ if(use.state) {
+ # subscript it if we're not giving the value
+ if(is.null(value)) expr <- substitute(a[group('(',italic(b),')')],list(a=expr, b=state))
+ else expr <- substitute(a*group('(',italic(b),')'),list(a=expr, b=state))
}
- # write logarithm of activity or fugacity
- # (or molality 20171101)
- if(log != "") {
- if(log == "aq") {
- if(use.molality) acity <- "m"
- else acity <- "a"
- } else if(log %in% c("cr", "liq", "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))
- expr <- substitute(a[b], list(a=logacity, b=expr))
- # write a value if given
+ # write a variable and value if given
+ if(!is.null(value) | log | molality) {
+ # write [logarithm of] activity or fugacity (or molality 20171101)
+ var <- "a"
+ if(molality) var <- "m"
+ if(state %in% c("g", "gas")) var <- "f"
+ expr <- substitute(italic(a)[b], list(a = var, b = expr))
+ # use the logarithm?
+ if(log) expr <- substitute(log ~ a, list(a = expr))
+ # write the value if not NULL or NA
if(!is.null(value)) {
- expr <- substitute(a==b, list(a=expr, b=value))
+ if(!is.na(value)) expr <- substitute(a == b, list(a = expr, b = value))
}
}
return(expr)
}
-expr.property <- function(property, use.molality=FALSE) {
+expr.property <- function(property, molality=FALSE) {
# a way to make expressions for various properties
# e.g. expr.property('DG0r') for standard molal Gibbs
# energy change of reaction
@@ -86,7 +77,7 @@
if(property=="logK") return(quote(log~italic(K)))
# grepl here b/c diagram() uses "loga.equil" and "loga.basis"
if(grepl("loga", property)) {
- if(use.molality) return(quote(log~italic(m)))
+ if(molality) return(quote(log~italic(m)))
else return(quote(log~italic(a)))
}
if(property=="alpha") return(quote(alpha))
@@ -160,7 +151,7 @@
return(expr)
}
-axis.label <- function(label, units=NULL, basis=get("thermo")$basis, prefix="", use.molality=FALSE) {
+axis.label <- function(label, units=NULL, basis=get("thermo")$basis, prefix="", molality=FALSE) {
# make a formatted axis label from a generic description
# it can be a chemical property, condition, or chemical activity in the system
# if the label matches one of the basis species
@@ -172,11 +163,11 @@
# 20090215: the state this basis species is in
state <- basis$state[match(label, rownames(basis))]
# get the formatted label
- desc <- expr.species(label, log=state, use.molality=use.molality)
+ desc <- expr.species(label, state = state, log = TRUE, molality = molality)
} else {
# the label is for a chemical property or condition
# make the label by putting a comma between the property and the units
- property <- expr.property(label, use.molality=use.molality)
+ property <- expr.property(label, molality=molality)
if(is.null(units)) units <- expr.units(label, prefix=prefix)
# no comma needed if there are no units
if(units=="") desc <- substitute(a, list(a=property))
@@ -187,22 +178,27 @@
}
describe.basis <- function(basis = get("thermo")$basis, ibasis = 1:nrow(basis),
- digits = 1, oneline = FALSE, use.molality = FALSE, use.pH = TRUE) {
+ digits = 1, oneline = FALSE, molality = FALSE, use.pH = TRUE) {
# make expressions for the chemical activities/fugacities of the basis species
propexpr <- valexpr <- character()
for(i in ibasis) {
- # propexpr is logarithm of activity or fugacity
- if(rownames(basis)[i]=="H+" & use.pH) thispropexpr <- "pH"
- else thispropexpr <- expr.species(rownames(basis)[i], log=basis$state[i], use.molality = use.molality)
- propexpr <- c(propexpr, thispropexpr)
if(can.be.numeric(basis$logact[i])) {
# we have an as.numeric here in case the basis$logact is character
# (by inclusion of a buffer for one of the other basis species)
- if(thispropexpr=="pH") valexpr <- c(valexpr, format(round(-as.numeric(basis$logact[i]), digits), nsmall=digits))
- else valexpr <- c(valexpr, format(round(as.numeric(basis$logact[i]), digits), nsmall=digits))
+ if(rownames(basis)[i]=="H+" & use.pH) {
+ propexpr <- c(propexpr, "pH")
+ valexpr <- c(valexpr, format(round(-as.numeric(basis$logact[i]), digits), nsmall=digits))
+ } else {
+ # propexpr is logarithm of activity or fugacity
+ propexpr <- c(propexpr, expr.species(rownames(basis)[i], state=basis$state[i], log = TRUE, molality = molality))
+ valexpr <- c(valexpr, format(round(as.numeric(basis$logact[i]), digits), nsmall=digits))
+ }
} else {
# a non-numeric value is the name of a buffer
valexpr <- c(valexpr, basis$logact[i])
+ # propexpr is pH, activity or fugacity
+ if(rownames(basis)[i]=="H+" & use.pH) propexpr <- c(propexpr, "pH")
+ else propexpr <- c(propexpr, expr.species(rownames(basis)[i], state=basis$state[i], value = NA, log = FALSE, molality = molality))
}
}
# write an equals sign between the property and value
@@ -259,8 +255,8 @@
if(i %in% iname) species <- reaction$name[i]
else {
# should the chemical formula have a state?
- if(identical(states,"all")) species <- expr.species(reaction$formula[i], state=reaction$state[i])
- else species <- expr.species(reaction$formula[i])
+ if(identical(states,"all")) species <- expr.species(reaction$formula[i], state=reaction$state[i], use.state=TRUE)
+ else species <- expr.species(reaction$formula[i], state=reaction$state[i])
}
# get the absolute value of the reaction coefficient
abscoeff <- abs(reaction$coeff[i])
@@ -287,7 +283,7 @@
}
# make formatted text for activity ratio 20170217
-ratlab <- function(ion="K+", use.molality=FALSE) {
+ratlab <- function(ion="K+", molality=FALSE) {
# the charge
Z <- makeup(ion)["Z"]
# the text for the exponent on aH+
@@ -295,8 +291,8 @@
# the expression for the ion and H+
expr.ion <- expr.species(ion)
expr.H <- expr.species("H+")
- # with use.molality, change a to m
- a <- ifelse(use.molality, "m", "a")
+ # with molality, change a to m
+ a <- ifelse(molality, "m", "a")
# the final expression
if(exp.H=="1") substitute(log~(italic(a)[expr.ion] / italic(a)[expr.H]), list(a=a, expr.ion=expr.ion, expr.H=expr.H))
else substitute(log~(italic(a)[expr.ion] / italic(a)[expr.H]^exp.H), list(a=a, expr.ion=expr.ion, expr.H=expr.H, exp.H=exp.H))
Modified: pkg/CHNOSZ/demo/gold.R
===================================================================
--- pkg/CHNOSZ/demo/gold.R 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/demo/gold.R 2018-11-11 04:55:31 UTC (rev 350)
@@ -156,11 +156,13 @@
diagram(s, add = TRUE, type = "loga.balance", lwd = 3, lty = 2)
# make legend and title
dP <- describe.property("P", 1000)
- dNaCl <- expression(NaCl == 2~mol~kg^-1)
- dK <- describe.basis(ibasis=5, use.molality=TRUE)
- legend("topleft", c(dP, dNaCl, dK), bty = "n")
- dbasis <- describe.basis(ibasis = c(9, 7, 10))
- legend(320, -3, dbasis, bty = "n")
+ dNaCl <- expression(italic(m)[NaCl] == 1.5)
+ dKCl <- expression(italic(m)[KCl] == 0.5)
+ legend("topleft", c(dP, dNaCl, dKCl), bty = "n")
+ dH2S <- describe.basis(ibasis = 7, molality=TRUE)
+ dO2 <- describe.basis(ibasis = 9)
+ dpH <- describe.basis(ibasis = 10)
+ legend(300, -3, c(dH2S, dO2, dpH), bty = "n")
title(main=("After Williams-Jones et al., 2009, Fig. 2B"), font.main = 1)
}
@@ -188,11 +190,13 @@
diagram(s, add = TRUE, type = "loga.balance", lwd = 3, lty = 2)
# make legend and title
dP <- describe.property("P", 1000)
- dNaCl <- expression(NaCl == 2~mol~kg^-1)
- dK <- describe.basis(ibasis=5, use.molality=TRUE)
- legend("topleft", c(dP, dNaCl, dK), bty = "n")
- dbasis <- describe.basis(ibasis = c(9, 7, 10))
- legend(320, -3, dbasis, bty = "n")
+ dNaCl <- expression(italic(m)[NaCl] == 1.5)
+ dKCl <- expression(italic(m)[KCl] == 0.5)
+ legend("topleft", c(dP, dNaCl, dKCl), bty = "n")
+ dH2S <- expr.species("H2S", value = 0.01, molality = TRUE)
+ dO2 <- describe.basis(ibasis = 9)
+ dpH <- describe.basis(ibasis = 10)
+ legend(300, -3, c(dH2S, dO2, dpH), bty = "n")
title(main=("After Williams-Jones et al., 2009, Fig. 2A"), font.main = 1)
}
Modified: pkg/CHNOSZ/inst/NEWS
===================================================================
--- pkg/CHNOSZ/inst/NEWS 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/inst/NEWS 2018-11-11 04:55:31 UTC (rev 350)
@@ -1,4 +1,4 @@
-CHANGES IN CHNOSZ 1.1.3-56 (2018-11-08)
+CHANGES IN CHNOSZ 1.1.3-58 (2018-11-11)
---------------------------------------
NEW FEATURES
@@ -37,6 +37,12 @@
of NaCl in water, taking account of activity coefficients and the
reaction Na+ + Cl- = NaCl(aq).
+- Add dumpdata() for returning/writing all packaged thermodynamic data
+ (including default database and optional data files).
+
+- Add demo/bison.R (average oxidation state of carbon of metagenome-
+ derived proteins in different microbial phyla at Bison Pool)
+
THERMODYNAMIC DATA
- The Berman data (Berman, 1988 and later additions) have replaced the
@@ -130,47 +136,46 @@
NA where the density of H2O is less than 0.35 g/cm3, avoiding the
output of bogus values in this region. Thanks to Evgeniy Bastrakov.
-OTHER CHANGES
+COMPUTATIONAL OPTIONS
- Add 'exceed.rhomin' argument to subcrt() and affinity() to enable
output of properties for species in the revised HKF model below 0.35
g/cm3.
-- To provide betters diagnostics for potential web apps, warning
+- In equilibrate(), accept a length > 1 'normalize' argument to
+ normalize the chemical formulas of only the selected species.
+
+- Add thermo$opt$maxcores (default 2) to specify maximum number of
+ cores for parallel calculations with palply().
+
+- Keywords in basis(): Change 'CHNOPS+' to use O2 instead of e-, and add
+ 'CHNOPSe' and 'MgCHNOPSe' for sets of basis species that have e-.
+
+- Add 'keep.duplicates' argument to thermo.refs(). Set it to TRUE to
+ output a single primary reference for each species, keeping any
+ duplicated references (but not including any secondary references in
+ thermo$obigt$ref2). Thanks to Evgeniy Bastrakov for the suggestion.
+
+USABILITY ENHANCEMENTS
+
+- To provide better diagnostics for potential web apps, warning
messages produced by subcrt() are now available in the output of
affinity(), under 'sout$warnings'.
- Change internal variable names in subcrt() for better readability
(sinfo -> ispecies, inpho -> iphases, sinph -> phasespecies).
-- Add dumpdata() for returning/writing all packaged thermodynamic data
- (including default database and optional data files).
-
- TODO: fix overly long message for info("SiO2").
-- In equilibrate(), accept a length > 1 'normalize' argument to
- normalize the chemical formulas of only the selected species.
-
- Add C implementation of counting occurrences of all letters in a
string (src/count_letters.c) to speed up operation of count.aa().
- read.fasta(): add support for file connections created using
archive::archive_read (https://github.com/jimhester/archive).
-- Add demo/bison.R (average oxidation state of carbon of metagenome-
- derived proteins in different microbial phyla at Bison Pool)
+- The arguments in expr.species() have been reorganized for more
+ flexible and concise usage.
-- Add thermo$opt$maxcores (default 2) to specify maximum number of
- cores for parallel calculations with palply().
-
-- Keywords in basis(): Change 'CHNOPS+' to use O2 instead of e-, and add
- 'CHNOPSe' and 'MgCHNOPSe' for sets of basis species that have e-.
-
-- Add 'keep.duplicates' argument to thermo.refs(). Set it to TRUE to
- output a single primary reference for each species, keeping any
- duplicated references (but not including any secondary references in
- thermo$obigt$ref2). Thanks to Evgeniy Bastrakov for the suggestion.
-
CHANGES IN CHNOSZ 1.1.3 (2017-11-13)
------------------------------------
Modified: pkg/CHNOSZ/man/berman.Rd
===================================================================
--- pkg/CHNOSZ/man/berman.Rd 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/man/berman.Rd 2018-11-11 04:55:31 UTC (rev 350)
@@ -111,7 +111,7 @@
res <- 400
a <- affinity(`K+` = c(0, 5, res), T = c(200, 650, res), P = 1000, exceed.Ttr = TRUE)
# make base plot with colors and no lines
-diagram(a, xlab = ratlab("K+", use.molality = TRUE), lty = 0, fill = "terrain")
+diagram(a, xlab = ratlab("K+", molality = TRUE), lty = 0, fill = "terrain")
# add the lines, extending into the low-density region (exceed.rhomin = TRUE)
a <- affinity(`K+` = c(0, 5, res), T = c(200, 650, res), P = 1000,
exceed.Ttr = TRUE, exceed.rhomin = TRUE)
Modified: pkg/CHNOSZ/man/buffer.Rd
===================================================================
--- pkg/CHNOSZ/man/buffer.Rd 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/man/buffer.Rd 2018-11-11 04:55:31 UTC (rev 350)
@@ -147,7 +147,7 @@
a <- affinity(O2=c(-85, -70, 4), T=c(25, 100, 4))
d <- diagram(a, type="CO2", add=TRUE, lty=2)
# add a legend
-lAC <- expr.species("CH3COOH", log="aq")
+lAC <- expr.species("CH3COOH", log=TRUE)
ltext <- c(as.expression(lAC), -3, -10)
lty <- c(NA, 1, 2)
legend("topright", legend=ltext, lty=lty, bg="white")
Modified: pkg/CHNOSZ/man/diagram.Rd
===================================================================
--- pkg/CHNOSZ/man/diagram.Rd 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/man/diagram.Rd 2018-11-11 04:55:31 UTC (rev 350)
@@ -166,7 +166,7 @@
\section{Activity Coefficients}{
The wording in this page and names of variables in functions refer exclusively to \samp{activities} of aqueous species.
However, if activity coefficients are calculated (using the \code{IS} argument in \code{\link{affinity}}), then these variables are effectively transformed to molalities (see \code{tests/testthat/} \code{test-logmolality.R}).
-So that the labels on diagrams are adjusted accordingly, \code{\link{diagram}} sets the \code{use.molality} argument of \code{\link{axis.label}} to TRUE if \code{IS} was supplied as an argument to \code{\link{affinity}}.
+So that the labels on diagrams are adjusted accordingly, \code{\link{diagram}} sets the \code{molality} argument of \code{\link{axis.label}} to TRUE if \code{IS} was supplied as an argument to \code{\link{affinity}}.
The labeling as molality takes effect even if \code{IS} is set to 0; this way, by including (or not) the \code{IS = 0} argument to \code{affinity}, the user decides whether to label aqueous species variables as molality (or activity) for calculations at zero ionic strength (where molality = activity).
}
@@ -273,9 +273,8 @@
# at 1 bar only, so we permit calculation at higher temperatures
a <- affinity(T=c(200, 900, 99), P=c(0, 9000, 101), exceed.Ttr=TRUE)
d <- diagram(a, fill=NULL)
-bexpr <- sapply(c("Al2O3", "SiO2", "H2O"), expr.species, simplify=FALSE)
-btext <- substitute(Al2O3 - SiO2 - H2O, unlist(bexpr))
-mtitle(c(as.expression(btext), "after Helgeson et al., 1978"))
+slab <- syslab(c("Al2O3", "SiO2", "H2O"))
+mtitle(c(as.expression(slab), "after Helgeson et al., 1978"))
# find the approximate position of the triple point
tp <- find.tp(d$predominant)
Ttp <- a$vals[[1]][tp[1, 2]]
Modified: pkg/CHNOSZ/man/util.expression.Rd
===================================================================
--- pkg/CHNOSZ/man/util.expression.Rd 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/man/util.expression.Rd 2018-11-11 04:55:31 UTC (rev 350)
@@ -14,28 +14,29 @@
\description{Generate expressions suitable for axis labels and plot legends describing chemical species, properties and reactions.}
\usage{
- expr.species(species, state = "", log = "", value=NULL, use.makeup=FALSE,
- use.molality = FALSE)
- expr.property(property, use.molality = FALSE)
+ expr.species(species, state = "aq", value=NULL, log=FALSE, molality=FALSE,
+ use.state=FALSE, use.makeup=FALSE)
+ expr.property(property, molality = FALSE)
expr.units(property, prefix = "", per = "mol")
axis.label(label, units = NULL, basis = get("thermo")$basis, prefix = "",
- use.molality = FALSE)
+ molality = FALSE)
describe.basis(basis = get("thermo")$basis, ibasis = 1:nrow(basis),
- digits = 1, oneline = FALSE, use.molality = FALSE, use.pH = TRUE)
+ digits = 1, oneline = FALSE, molality = FALSE, use.pH = TRUE)
describe.property(property, value, digits = 0, oneline = FALSE,
ret.val = FALSE)
describe.reaction(reaction, iname = numeric(), states = NULL)
syslab(system = c("K2O", "Al2O3", "SiO2", "H2O"), dash="\u2013")
- ratlab(ion = "K+", use.molality = FALSE)
+ ratlab(ion = "K+", molality = FALSE)
}
\arguments{
\item{species}{character, formula of a chemical species}
\item{state}{character, designation of physical state}
- \item{log}{character, designation of physical state (for logarithm of activity or fugacity)}
\item{value}{numeric, logarithm of activity or fugacity of species, or value of other property}
+ \item{log}{logical, write logarithm of activity/fugacity/molality?}
+ \item{molality}{logical, use molality (m) instead of activity (a) for aqueous species?}
+ \item{use.state}{logical, include state in expression?}
\item{use.makeup}{logical, use \code{\link{makeup}} to count the elements?}
- \item{use.molality}{logical, use molality (m) instead of activity (a) for aqueous species?}
\item{use.pH}{logical, use pH instead of log activity of H+?}
\item{property}{character, description of chemical property}
\item{prefix}{character, prefix for units}
@@ -59,10 +60,15 @@
The \code{expr.*} functions create \code{\link{expression}}s using the \code{\link{plotmath}} syntax to describe the names and states and logarithms of activity or fugacity of chemical species, conditions including temperature and pressure and chemical properties such as Gibbs energy and volume.
-\code{expr.species} takes as input the formula of a single chemical \code{species} and constructs an expression including subscripted coefficients, and a suffixed designation of physical \code{state} (italicized, in parentheses) if provided.
-If \code{log} designates a physical state (as in \code{\link{thermo}$obigt$state}), the expression includes a \samp{log} prefix, followed by \samp{f} for fugacity of gaseous species, or \samp{a} for activity of species in all other states.
+\code{expr.species} constructs a formatted expression using the formula or name of a single chemical \code{species}.
+With no other arguments, the formula is just formatted with the appropriate subscripts and superscripts.
+Providing the physical \code{state} adds a variable to the expression (\emph{a} for aqueous species and pure phases, except \emph{f} for gases).
+Set \code{molality} to TRUE to write \emph{m} instead of \emph{a} for aqueous species.
+The state itself is written in the expression if \code{use.state} is TRUE.
+If \code{log} is TRUE, the expression includes a \samp{log} prefix.
+Finally, provide a value in \code{value} to write an equation (something like logfO2 = -70), or set it to NA to only write the variable itself (e.g. logfO2).
Set \code{use.makeup} to TRUE to use \code{\link{makeup}} to parse the chemical formula.
-This can have the undesirable effect of reordering and grouping all the elements, and has been replaced with a different splitting algorithm so that coefficients and charges are sub/superscripted without affecting the intervening text.
+This was an older default action that had the undesirable effect of reordering and grouping all the elements, and has been replaced with a different splitting algorithm so that coefficients and charges are sub/superscripted without affecting the intervening text.
\code{expr.property} accepts a description in \code{property} that indicates the chemical property of interest.
Uppercase letters are italicized, and lowercase letters are italicized and subscripted.
@@ -100,7 +106,7 @@
If this matches the chemical formula of one of the basis species in the \code{basis} argument, the expression for the label is generated using \code{expr.species} with \code{log} set to the physical state of the basis species.
Otherwise, the expression is built by combining the output of \code{expr.property} with \code{expr.units} (or the value in \code{units}, if it is supplied), placing a comma between the two.
This function is used extensively in \code{\link{diagram}} and also appears in many of the examples.
-Note that \code{\link{diagram}} sets \code{use.molality} to TRUE if \code{IS} was supplied as an argument to \code{\link{affinity}}.
+Note that \code{\link{diagram}} sets \code{molality} to TRUE if \code{IS} was supplied as an argument to \code{\link{affinity}}.
\code{describe.basis} makes an expression summarizing the basis species definition (logarithms of activity or fugacity of the basis species) provided in \code{basis}; only the basis species identified by \code{ibasis} are included.
@@ -124,10 +130,10 @@
text0 <- function(...) text(..., adj=0)
# species
text0(1, 1, expr.species("CO2"))
-text0(1, 2, expr.species("CO2", state="aq"))
-text0(1, 3, expr.species("CO2", state="aq", log="aq"))
-text0(1, 4, expr.species("CO2", log="aq"))
-text0(1, 5, expr.species("CO2", log="aq", value=-3))
+text0(1, 2, expr.species("CO2", use.state=TRUE))
+text0(1, 3, expr.species("CO2", log=TRUE, use.state=TRUE))
+text0(1, 4, expr.species("CO2", log=TRUE))
+text0(1, 5, expr.species("CO2", log=TRUE, value=-3))
# properties
text0(2, 1, expr.property("A"))
text0(2, 2, expr.property("DV"))
Modified: pkg/CHNOSZ/tests/testthat/test-util.R
===================================================================
--- pkg/CHNOSZ/tests/testthat/test-util.R 2018-11-11 02:22:09 UTC (rev 349)
+++ pkg/CHNOSZ/tests/testthat/test-util.R 2018-11-11 04:55:31 UTC (rev 350)
@@ -36,10 +36,8 @@
expect_equal(as.numeric(testGHS[1, 3]), testent[2])
})
-
test_that("expr.species() produces expected errors", {
expect_error(expr.species(c("H2O", "CO2")), "more than one species")
- expect_error(expr.species("CO2", log = "aqq"), "'aqq' is not a recognized state")
})
test_that("[P|T|E].units() do not accept invalid units", {
More information about the CHNOSZ-commits
mailing list