[CHNOSZ-commits] r566 - in pkg/CHNOSZ: . R inst man man/macros tests/testthat vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sat Jul 18 05:41:13 CEST 2020
Author: jedick
Date: 2020-07-18 05:41:12 +0200 (Sat, 18 Jul 2020)
New Revision: 566
Added:
pkg/CHNOSZ/R/combine.R
pkg/CHNOSZ/man/combine.Rd
Removed:
pkg/CHNOSZ/R/duplex.R
pkg/CHNOSZ/man/duplex.Rd
Modified:
pkg/CHNOSZ/DESCRIPTION
pkg/CHNOSZ/NAMESPACE
pkg/CHNOSZ/R/diagram.R
pkg/CHNOSZ/R/examples.R
pkg/CHNOSZ/inst/NEWS.Rd
pkg/CHNOSZ/man/macros/macros.Rd
pkg/CHNOSZ/tests/testthat/test-diagram.R
pkg/CHNOSZ/vignettes/mklinks.sh
pkg/CHNOSZ/vignettes/multi-metal.Rmd
pkg/CHNOSZ/vignettes/vig.bib
Log:
Add combine() and examples to multi-metal.Rmd
Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/DESCRIPTION 2020-07-18 03:41:12 UTC (rev 566)
@@ -1,6 +1,6 @@
-Date: 2020-07-16
+Date: 2020-07-18
Package: CHNOSZ
-Version: 1.3.6-39
+Version: 1.3.6-40
Title: Thermodynamic Calculations and Diagrams for Geochemistry
Authors at R: c(
person("Jeffrey", "Dick", , "j3ffdick at gmail.com", role = c("aut", "cre"),
Modified: pkg/CHNOSZ/NAMESPACE
===================================================================
--- pkg/CHNOSZ/NAMESPACE 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/NAMESPACE 2020-07-18 03:41:12 UTC (rev 566)
@@ -59,7 +59,7 @@
"CHNOSZ", "thermo", "reset", "OBIGT", "retrieve", "AkDi", "moles",
"lNaCl", "lS", "lT", "lP", "lTP", "lex",
# added 20200716 or later
- "duplex"
+ "duplex", "combine"
)
# Load shared objects
Copied: pkg/CHNOSZ/R/combine.R (from rev 565, pkg/CHNOSZ/R/duplex.R)
===================================================================
--- pkg/CHNOSZ/R/combine.R (rev 0)
+++ pkg/CHNOSZ/R/combine.R 2020-07-18 03:41:12 UTC (rev 566)
@@ -0,0 +1,138 @@
+# CHNOSZ/duplex.R
+# Combine diagrams for two metals
+# 20200713 first version jmd
+
+# Function to combine two diagrams (simple overlay, no interaction) 20200717
+# -- makes new "species" from all combinations of those in d1 and d2
+combine <- function(d1, d2) {
+ check_d1_d2(d1, d2)
+
+ # Index all combinations of species in d1 and d2
+ i1 <- 1:nrow(d1$species)
+ i2 <- 1:nrow(d2$species)
+ combs <- expand.grid(i1, i2)
+
+ # Get species rows for each combination
+ s1 <- d1$species[combs[, 1], ]
+ s2 <- d2$species[combs[, 2], ]
+ # Make a new species data frame
+ nbasis <- nrow(d1$basis)
+ species <- s1[, 1:nbasis] + s2[, 1:nbasis]
+ ispecies <- paste(s1$ispecies, s2$ispecies, sep = ",")
+ logact <- paste(s1$logact, s2$logact, sep = ",")
+ state <- paste(s1$state, s2$state, sep = ",")
+ # Use names from diagram()
+ if(is.expression(d1$names) & is.expression(d2$names)) {
+ name <- lapply(1:nrow(combs), function(i) bquote(.(d1$names[[combs[i, 1]]])+.(d2$names[[combs[i, 2]]])))
+ name <- unlist(lapply(name, deparse, width.cutoff = 500, control = NULL))
+ } else if(is.expression(d1$names)) {
+ name <- lapply(1:nrow(combs), function(i) bquote(.(d1$names[[combs[i, 1]]])+.(d2$names[combs[i, 2]])))
+ name <- unlist(lapply(name, deparse, width.cutoff = 500, control = NULL))
+ } else if(is.expression(d2$names)) {
+ name <- lapply(1:nrow(combs), function(i) bquote(.(d1$names[combs[i, 1]])+.(d2$names[[combs[i, 2]]])))
+ name <- unlist(lapply(name, deparse, width.cutoff = 500, control = NULL))
+ } else name <- paste(d1$names[combs[, 1]], d2$names[combs[, 2]], sep="+")
+ if(length(name) != nrow(combs)) stop("deparse()-ing expressions gives unequal length; try diagram(., format.names = FALSE)")
+ species <- cbind(species, ispecies, logact, state, name)
+
+ # Get affinities for each combination
+ v1 <- d1$values[combs[, 1]]
+ v2 <- d2$values[combs[, 2]]
+ values <- Map("+", v1, v2)
+ # Assign -Inf affinity where a species isn't predominant
+ for(i in seq_along(values)) {
+ i1 <- combs[i, 1]
+ i2 <- combs[i, 2]
+ ip1 <- d1$predominant == i1
+ ip2 <- d2$predominant == i2
+ ip12 <- ip1 & ip2
+ values[[i]][!ip12] <- -Inf
+ }
+
+ # Use d1 as a template for the new affinity object
+ anew <- d1[1:11]
+ # Insert combined results
+ anew$species <- species
+ anew$values <- values
+ # We don't have sout (results from subcrt()) for the combined "species"
+ anew$sout <- NULL
+ anew
+}
+
+# Function to make a new "affinity" object from two diagrams 20200713
+# -- uses *secondary* balancing coefficients to combine the diagrams
+duplex <- function(d1, d2, balance = NULL) {
+ check_d1_d2(d1, d2)
+
+ # Combine the species data frames
+ species <- rbind(d1$species, d2$species)
+ # Combine the sout objects (results from subcrt())
+ only2 <- !d2$sout$species$ispecies %in% d1$sout$species$ispecies
+ sout <- d1$sout
+ sout$species <- rbind(sout$species, d2$sout$species[only2, ])
+ sout$out <- c(sout$out, d2$sout$out[only2])
+ # Combine the affinity values divided by the *primary*
+ # balancing coefficients ("plotvals" from diagram())
+ values <- c(d1$plotvals, d2$plotvals)
+
+ # Use d1 as a template for the new affinity object
+ anew <- d1[1:11]
+ # Insert combined results
+ anew$species <- species
+ anew$sout <- sout
+ anew$values <- values
+
+ # Figure out the *secondary* balancing coefficients
+ n.balance <- balance(anew, balance = balance)$n.balance
+ # In the Fe-Cu-S-O-H example all the coefficients on H+ are negative
+ if(all(n.balance < 0)) n.balance <- -n.balance
+ n1 <- nrow(d1$species)
+ n.balance.1 <- n.balance[1:n1]
+ n.balance.2 <- n.balance[(n1+1):length(n.balance)]
+
+ # Make empty matrices to hold affinities and balancing coefficients
+ a1 <- d1$values[[1]]
+ a1[] <- NA
+ b2 <- a2 <- b1 <- a1
+ # Get the affinities (per mole of species, not divided by any balancing coefficients)
+ # and the secondary balancing coefficients for the predominant species in each diagram
+ p1 <- d1$predominant
+ for(ip in unique(as.vector(p1))) {
+ a1[p1 == ip] <- d1$values[[ip]][p1 == ip]
+ b1[p1 == ip] <- n.balance.1[ip]
+ }
+ p2 <- d2$predominant
+ for(ip in unique(as.vector(p2))) {
+ a2[p2 == ip] <- d2$values[[ip]][p2 == ip]
+ b2[p2 == ip] <- n.balance.2[ip]
+ }
+ # Divide the affinities by the secondary balancing coefficients
+ ab1 <- a1 / b1
+ ab2 <- a2 / b2
+ # Identify the species with the highest affinity (predominant in the *secondary* reactions)
+ i1 <- ab1 > ab2
+ # Suppress non-predominant species at each grid point
+ for(i in 1:n1) anew$values[[i]][!i1] <- -Inf
+ for(i in (n1+1):length(n.balance)) anew$values[[i]][i1] <- -Inf
+
+ anew
+
+}
+
+### unexported function ###
+
+# Check that d1 and d2 can be combined
+# Extracted from duplex() 20200717
+check_d1_d2 <- function(d1, d2) {
+ # Check that the basis species are the same
+ if(!identical(d1$basis, d2$basis)) stop("basis species in objects 'd1' and 'd2' are not identical")
+ # Check that the variables and their values are the same
+ if(!identical(d1$vars, d2$vars)) stop("variable names in objects 'd1' and 'd2' are not identical")
+ if(!identical(d1$vals, d2$vals)) stop("variable values in objects 'd1' and 'd2' are not identical")
+ # Check that T and P are the same
+ if(!identical(d1$T, d2$T)) stop("temperatures in objects 'd1' and 'd2' are not identical")
+ if(!identical(d1$P, d2$P)) stop("pressures in objects 'd1' and 'd2' are not identical")
+ # Check that we have plotvals and predominant (from diagram())
+ if(is.null(d1$plotvals) | is.null(d1$predominant)) stop("object 'd1' is missing 'plotvals' or 'predominant' components (not made by diagram()?)")
+ if(is.null(d2$plotvals) | is.null(d2$predominant)) stop("object 'd2' is missing 'plotvals' or 'predominant' components (not made by diagram()?)")
+}
Modified: pkg/CHNOSZ/R/diagram.R
===================================================================
--- pkg/CHNOSZ/R/diagram.R 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/R/diagram.R 2020-07-18 03:41:12 UTC (rev 566)
@@ -40,8 +40,10 @@
### argument handling ###
- ## check that eout is valid input
- if(!"sout" %in% names(eout)) stop("'eout' does not look like output from equil() or affinity()")
+ ## check that eout is a valid object
+ efun <- eout$fun
+ if(length(efun)==0) efun <- ""
+ if(!efun %in% c("affinity", "equilibrate", "solubility")) stop("'eout' is not the output from affinity(), equilibrate(), or solubility()")
## 'type' can be:
# 'auto' - property from affinity() (1D) or maximum affinity (affinity 2D) (aout) or loga.equil (eout)
@@ -273,13 +275,26 @@
## apply formatting to chemical formulas 20170204
if(all(grepl("_", names))) is.pname <- TRUE
if(format.names & !is.pname) {
- exprnames <- as.expression(names)
+ # check if names are a deparsed expression (used in combine()) 20200718
+ parsed <- FALSE
+ if(any(grepl("paste\\(", names))) {
+ exprnames <- parse(text = names)
+ if(length(exprnames) != length(names)) stop("parse()-ing names gives length not equal to number of names")
+ parsed <- TRUE
+ } else {
+ exprnames <- as.expression(names)
+ # get formatted chemical formulas
+ for(i in seq_along(exprnames)) exprnames[[i]] <- expr.species(exprnames[[i]])
+ }
+ # apply bold or italic
+ bold <- rep(bold, length.out = length(exprnames))
+ italic <- rep(italic, length.out = length(exprnames))
for(i in seq_along(exprnames)) {
- exprnames[[i]] <- expr.species(exprnames[[i]])
- if(bold) exprnames[[i]] <- substitute(bold(a), list(a=exprnames[[i]]))
- if(italic) exprnames[[i]] <- substitute(italic(a), list(a=exprnames[[i]]))
+ if(bold[i]) exprnames[[i]] <- substitute(bold(a), list(a=exprnames[[i]]))
+ if(italic[i]) exprnames[[i]] <- substitute(italic(a), list(a=exprnames[[i]]))
}
- names <- exprnames
+ # only use the expression if it's different from the unformatted names
+ if(parsed | !identical(as.character(exprnames), names)) names <- exprnames
}
if(nd==0) {
@@ -544,7 +559,7 @@
lapply(linesout, `length<-`, max(sapply(linesout, length)))
}
## label plot function
- plot.names <- function(out, xs, ys, xlim, ylim, names) {
+ plot.names <- function(out, xs, ys, xlim, ylim, names, srt) {
# calculate coordinates for field labels
# revisions: 20091116 for speed, 20190223 work with user-specified xlim and ylim
namesx <- namesy <- rep(NA, length(names))
@@ -564,7 +579,14 @@
namesy[i] <- mean(ysth)
}
# fields that really exist on the plot
- if(!is.null(names)) text(namesx, namesy, labels=names, cex=cex.names, col=col.names, font=font, family=family)
+ if(!is.null(names)) {
+ cex <- rep(cex.names, length.out = length(names))
+ col <- rep(col.names, length.out = length(names))
+ font <- rep(font, length.out = length(names))
+ family <- rep(family, length.out = length(names))
+ srt <- rep(srt, length.out = length(names))
+ for(i in seq_along(names)) text(namesx[i], namesy[i], labels=names[i], cex=cex[i], col=col[i], font=font[i], family=family[i], srt = srt[i])
+ }
return(list(namesx=namesx, namesy=namesy))
}
@@ -646,7 +668,7 @@
# put predominance matrix in the right order for image() etc
zs <- t(predominant[, ncol(predominant):1])
if(!is.null(fill)) fill.color(xs, ys, zs, fill, ngroups)
- pn <- plot.names(zs, xs, ys, xlim, ylim, names)
+ pn <- plot.names(zs, xs, ys, xlim, ylim, names, srt)
# only draw the lines if there is more than one field 20180923
# (to avoid warnings from contour, which seem to be associated with weird
# font metric state and subsequent errors adding e.g. subscripted text to plot)
Deleted: pkg/CHNOSZ/R/duplex.R
===================================================================
--- pkg/CHNOSZ/R/duplex.R 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/R/duplex.R 2020-07-18 03:41:12 UTC (rev 566)
@@ -1,74 +0,0 @@
-# CHNOSZ/duplex.R
-# Combine diagrams for two metals
-# 20200713 first version jmd
-
-# Function to make a new "affinity" object from two diagrams;
-# uses *secondary* balancing coefficients to combine the diagrams
-duplex <- function(d1, d2, balance = NULL) {
- # Check that the basis species are the same
- if(!identical(d1$basis, d2$basis)) stop("basis species are not identical")
- # Check that the variables and their values are the same
- if(!identical(d1$vars, d2$vars)) stop("variable names are not identical")
- if(!identical(d1$vals, d2$vals)) stop("variable values are not identical")
- # Check that T and P are the same
- if(!identical(d1$T, d2$T)) stop("temperatures are not identical")
- if(!identical(d1$P, d2$P)) stop("pressures are not identical")
- # Check that we have plotvals and predominant (from diagram())
- if(is.null(d1$plotvals) | is.null(d1$predominant)) stop("d1 is missing 'plotvals' or 'predominant' components (not made by diagram()?)")
- if(is.null(d2$plotvals) | is.null(d2$predominant)) stop("d2 is missing 'plotvals' or 'predominant' components (not made by diagram()?)")
-
- # Combine the species data frames
- species <- rbind(d1$species, d2$species)
- # Combine the sout objects (results from subcrt())
- only2 <- !d2$sout$species$ispecies %in% d1$sout$species$ispecies
- sout <- d1$sout
- sout$species <- rbind(sout$species, d2$sout$species[only2, ])
- sout$out <- c(sout$out, d2$sout$out[only2])
- # Combine the affinity values divided by the *primary*
- # balancing coefficients ("plotvals" from diagram())
- values <- c(d1$plotvals, d2$plotvals)
-
- # Use d1 as a template for the new affinity object
- anew <- d1[1:11]
- # Insert combined results
- anew$species <- species
- anew$sout <- sout
- anew$values <- values
-
- # Figure out the *secondary* balancing coefficients
- n.balance <- balance(anew, balance = balance)$n.balance
- # In the Fe-Cu-S-O-H example all the coefficients on H+ are negative
- if(all(n.balance < 0)) n.balance <- -n.balance
- n1 <- nrow(d1$species)
- n.balance.1 <- n.balance[1:n1]
- n.balance.2 <- n.balance[(n1+1):length(n.balance)]
-
- # Make empty matrices to hold affinities and balancing coefficients
- a1 <- d1$values[[1]]
- a1[] <- NA
- b2 <- a2 <- b1 <- a1
- # Get the affinities (per mole of species, not divided by any balancing coefficients)
- # and the secondary balancing coefficients for the predominant species in each diagram
- p1 <- d1$predominant
- for(ip in unique(as.vector(p1))) {
- a1[p1 == ip] <- d1$values[[ip]][p1 == ip]
- b1[p1 == ip] <- n.balance.1[ip]
- }
- p2 <- d2$predominant
- for(ip in unique(as.vector(p2))) {
- a2[p2 == ip] <- d2$values[[ip]][p2 == ip]
- b2[p2 == ip] <- n.balance.2[ip]
- }
- # Divide the affinities by the secondary balancing coefficients
- ab1 <- a1 / b1
- ab2 <- a2 / b2
- # Identify the species with the highest affinity (predominant in the *secondary* reactions)
- i1 <- ab1 > ab2
- # Suppress non-predominant species at each grid point
- for(i in 1:n1) anew$values[[i]][!i1] <- -Inf
- for(i in (n1+1):length(n.balance)) anew$values[[i]][i1] <- -Inf
-
- anew
-
-}
-
Modified: pkg/CHNOSZ/R/examples.R
===================================================================
--- pkg/CHNOSZ/R/examples.R 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/R/examples.R 2020-07-18 03:41:12 UTC (rev 566)
@@ -12,7 +12,7 @@
"hkf", "water", "IAPWS95", "subcrt", "berman",
"makeup", "basis", "swap.basis", "species", "affinity",
"solubility", "equilibrate",
- "diagram", "mosaic", "duplex",
+ "diagram", "mosaic", "combine",
"buffer", "nonideal", "NaCl",
"add.protein", "protein", "ionize.aa",
"objective", "revisit", "EOSregress", "wjd")
Modified: pkg/CHNOSZ/inst/NEWS.Rd
===================================================================
--- pkg/CHNOSZ/inst/NEWS.Rd 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/inst/NEWS.Rd 2020-07-18 03:41:12 UTC (rev 566)
@@ -6,7 +6,7 @@
\newcommand{\H2O}{\ifelse{latex}{\eqn{\mathrm{H_{2}O}}}{\ifelse{html}{\out{H<sub>2</sub>O}}{H2O}}}
\newcommand{\Hplus}{\ifelse{latex}{\eqn{\mathrm{H^{+}}}}{\ifelse{html}{\out{H<sup>+</sup>}}{H+}}}
-\section{Changes in CHNOSZ version 1.3.6-38 (2020-07-16)}{
+\section{Changes in CHNOSZ version 1.3.6-40 (2020-07-18)}{
\subsection{MAJOR CHANGES}{
\itemize{
@@ -29,16 +29,38 @@
\code{add.OBIGT()} and \code{mod.OBIGT()} replace the previous
\code{add.obigt()} and \code{mod.obigt()}.
- \item An \strong{add} argument has been added to \code{species()}. By
- default, loading new species now causes any previous species definition
- to be removed first. Therefore, \code{species(delete = TRUE)} is no
- longer needed to clear the species definition in a series of calculations
- for different species. To restore the previous behavior, where new
- species are added to an existing definiton, set \samp{add = TRUE}.
+ \item An \strong{add} argument has been added to \code{species()}. With
+ the default of \code{FALSE}, loading new species now causes any previous species
+ definition to be removed first. Therefore, \code{species(delete = TRUE)}
+ is no longer needed to clear the species definition in a series of
+ calculations for different species. To add new species to an existing
+ system, use \samp{add = TRUE}.
}
}
+ \subsection{NEW FEATURES}{
+ \itemize{
+
+ \item Add function \strong{combine()} for combining diagrams for different
+ systems (i.e., simple overlay with labels for species from both systems).
+
+ \item Add function \strong{duplex()} for making a new diagram by secondary
+ balancing between two systems.
+
+ \item Add a \strong{predominant} argument to \code{mosaic()} to use
+ previously calculated predominances of species (e.g. minerals) for the
+ changing basis species. This allows \code{mosaic()} calculations to be
+ linked in series, for instance to sequentially add metals (Fe, then Cu)
+ to a diagram.
+
+ \item Add vignette \strong{multi-metal.Rmd} for examples that use these
+ new features to make diagrams for minerals with multiple metals
+ (specifically Fe and Cu).
+
+ }
+ }
+
\subsection{CHANGES TO OBIGT DATABASE}{
\itemize{
@@ -75,9 +97,6 @@
\subsection{DEMOS AND VIGNETTES}{
\itemize{
- \item Add \samp{multi-metal.Rmd} for examples of diagrams with multiple
- metals.
-
\item Add \samp{demo/comproportionation.R}: Gibbs energy of sulfur
comproportionation, after
\href{https://doi.org/10.1111/1462-2920.14982}{Amend et al., 2020}.
@@ -121,15 +140,15 @@
\item \code{mosaic()} now allows a \strong{blend} argument of length > 1 to
apply a specific setting to each group of basis species.
- \item Add a \strong{predominant} argument to \code{mosaic()} to use
- previously calculated predominances of species (e.g. minerals) for the
- changing basis species. This allows \code{mosaic()} calculations to be
- linked in series, for instance to sequentially add metals (Fe, then Cu)
- to a diagram.
-
\item Add a \strong{bottom} argument to \code{ratlab()} to allow changing
the ion in the denominator to something other than \Hplus.
+ \item The \strong{srt} argument in \code{diagram()} now can be used to
+ rotate field labels, not only line labels. This and other arguments
+ (\strong{cex}, \strong{col}, \strong{font}, \strong{family},
+ \strong{bold}, \strong{italic}) now can have length > 1 to apply
+ different settings to each species.
+
}
}
@@ -170,6 +189,18 @@
\item TODO: OBIGT.Rmd: change "CHNOSZ" references to "OBIGT".
+ \item TODO: fix logic for negating balancing coefficients in diagram().
+
+ \item TODO: add elements from RH95.
+
+ \item TODO: check that Rd <-> vignette links work in Rstudio help viewer.
+
+ \item TODO: add check to mosaic() that 'predominant' values have the right dimensions.
+
+ \item TODO: remove stopifnot() from examples.
+
+ \item TODO: get axis labels at intervals of 5 (0, 5, 10, 15) in thermo.plot.new(c(0, 15), c(0, 15), "x", "y").
+
}
}
Copied: pkg/CHNOSZ/man/combine.Rd (from rev 565, pkg/CHNOSZ/man/duplex.Rd)
===================================================================
--- pkg/CHNOSZ/man/combine.Rd (rev 0)
+++ pkg/CHNOSZ/man/combine.Rd 2020-07-18 03:41:12 UTC (rev 566)
@@ -0,0 +1,74 @@
+\encoding{UTF-8}
+\name{combine}
+\alias{combine}
+\alias{duplex}
+\title{Combine Diagrams}
+\description{
+ Combine two diagrams for different systems by simple overlay or using a secondary balancing constraint.
+}
+
+\usage{
+ combine(d1, d2)
+ duplex(d1, d2, balance = NULL)
+}
+
+\arguments{
+ \item{d1}{list, output of \code{\link{diagram}} for first system}
+ \item{d2}{list, output of \code{diagram} for second system}
+ \item{balance}{character or numeric, specification of secondary balancing coefficients}
+}
+
+\details{
+
+These functions both make a new \code{\link{affinity}} object from two diagrams.
+In essence, \code{combine} identifies the intersection of predominance fields for all possible combinations of species (without interaction between the systems), while \code{duplex} creates a new diagram by comparing the affinities of reactions between species in both systems.
+Both functions mask the non-predominant species by assigning them -Inf values of affinity, so the result can be used to make a new diagram that shows the combined system.
+
+\code{combine} makes a simple overlay of two diagrams using new "species" generated from all combinations of those in \code{d1} and \code{d2}.
+The new names are formed from the \code{names} used in the source diagrams; for example if "Cp" and "Py" are predominant minerals at the same position in diagrams 1 and 2, the field for the combined diagram will be labeled "Cp+Py".
+The resulting affinities are simply the sum of affinities of the two species; they are assigned values of -Inf wherever one of the species is not predominant in either of the diagrams.
+
+Diagrams for different systems likely use different \emph{primary} balancing coefficients, such as balancing on different metals.
+\code{duplex} uses \emph{secondary} balancing coefficients, specified acording to \code{balance} (see \code{\link{equilibrate}} for a description of this argument), to determine the reactions between the species in the two systems.
+The affinities of these reactions are then used \emph{only} to identify the predominant species at each grid point.
+The \emph{returned} value of affinity are carried forward from those used to make the source diagrams (\samp{plotvals} in \code{d1} and \code{d2}), and therefore reflect the primary balancing coefficients.
+The returned values are assigned -Inf wherever that species is determined to not predominate according to the secondary balancing.
+
+\code{combine} yields finite values of affinity for only a single species at any grid point, the final diagram can be made with any setting of \code{balance}.
+However, for \code{duplex}, \code{balance} in the final diagram should be set to \samp{1} in order to preserve the primary balancing coefficients.
+
+}
+
+\value{
+A list object with the same structure as the output from \code{\link{affinity}}, so it can be used as input to \code{diagram}.
+}
+
+\seealso{
+A longer example is in the vignette \viglink{multi-metal}.
+}
+
+\examples{\dontshow{opar <- par(no.readonly = TRUE)}
+par(mfrow = c(2, 2))
+# Define basis species with Fe and Cu
+basis(c("Fe+2", "Cu+", "hydrogen sulfide", "oxygen", "H2O", "H+"))
+xlab <- ratlab("Fe+2", "Cu+")
+# Calculate diagram for only Fe-bearing minerals
+species(c("pyrite", "pyrrhotite", "magnetite", "hematite"))
+aFe <- affinity("Fe+2" = c(0, 12), O2 = c(-40, -16), T = 400, P = 2000)
+dFe <- diagram(aFe, xlab = xlab, main = "Fe-S-O-H")
+# Calculate diagram for only Cu-bearing minerals
+species(c("covellite", "chalcocite", "tenorite", "cuprite"))
+aCu <- affinity(aFe) # argument recall
+dCu <- diagram(aCu, xlab = xlab, main = "Cu-S-O-H")
+### combine() diagram
+ac <- combine(dFe, dCu)
+diagram(ac, xlab = xlab, main = "Cu-Fe-S-O-H with combine()")
+### duplex() diagram
+ad <- duplex(dFe, dCu)
+diagram(ad, xlab = xlab, balance = 1, main = "Cu-Fe-S-O-H with duplex()")
+db <- describe.basis(ibasis = 3)
+leg <- lex(lTP(400, 2000), db)
+legend("bottomleft", legend = leg, bty = "n")
+\dontshow{par(opar)}}
+
+\concept{Extended workflow}
Deleted: pkg/CHNOSZ/man/duplex.Rd
===================================================================
--- pkg/CHNOSZ/man/duplex.Rd 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/man/duplex.Rd 2020-07-18 03:41:12 UTC (rev 566)
@@ -1,60 +0,0 @@
-\encoding{UTF-8}
-\name{duplex}
-\alias{duplex}
-\title{Combine Diagrams for Two Metals}
-\description{
- Combine diagrams from two systems using a secondary balancing constraint.
-}
-
-\usage{
- duplex(d1, d2, balance = NULL)
-}
-
-\arguments{
- \item{d1}{list, output of \code{\link{diagram}} for first system}
- \item{d2}{list, output of \code{diagram} for second system}
- \item{balance}{character or numeric, specification of secondary balancing coefficients}
-}
-
-\details{
-
-This function makes a new \code{\link{affinity}} object from two diagrams.
-Each of the diagrams likely uses different \emph{primary} balancing coefficients (e.g., balancing on different metals).
-This function uses \emph{secondary} balancing coefficients to combine the diagrams.
-
-See \code{\link{equilibrate}} for a description of the \code{balance} argument.
-
-}
-
-\value{
-A list object with the same structure as the output from \code{\link{affinity}}, so it can be used as input to \code{diagram}.
-}
-
-\seealso{
-A longer example is in the vignette \viglink{multi-metal}.
-}
-
-\examples{
-def.par <- par(no.readonly = TRUE)
-layout(matrix(c(1, 1, 2, 2, 0, 3, 3, 0), nrow = 4))
-# Define basis species with Fe and Cu
-basis(c("Fe+2", "Cu+", "hydrogen sulfide", "oxygen", "H2O", "H+"))
-xlab <- ratlab("Fe+2", "Cu+")
-# Calculate diagram for only Fe-bearing minerals
-species(c("pyrite", "pyrrhotite", "magnetite", "hematite"))
-aFe <- affinity("Fe+2" = c(0, 12), O2 = c(-40, -16), T = 400, P = 2000)
-dFe <- diagram(aFe, xlab = xlab, main = "Fe-S-O-H")
-# Calculate diagram for only Cu-bearing minerals
-species(c("covellite", "chalcocite", "tenorite", "cuprite"))
-aCu <- affinity(aFe) # argument recall
-dCu <- diagram(aCu, xlab = xlab, main = "Cu-S-O-H")
-# Combine the diagrams
-ad <- duplex(dFe, dCu)
-diagram(ad, xlab = xlab, balance = 1, main = "Cu-Fe-S-O-H")
-db <- describe.basis(ibasis = 3)
-leg <- lex(lTP(400, 2000), db)
-legend("bottomleft", legend = leg, bty = "n")
-par(def.par)
-}
-
-\concept{Extended workflow}
Modified: pkg/CHNOSZ/man/macros/macros.Rd
===================================================================
--- pkg/CHNOSZ/man/macros/macros.Rd 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/man/macros/macros.Rd 2020-07-18 03:41:12 UTC (rev 566)
@@ -51,4 +51,4 @@
\newcommand{\alpha}{\ifelse{latex}{\eqn{\alpha}}{\ifelse{html}{\out{α}}{α}}}
% links to vignettes 20200716
-\newcommand{\viglink}{\ifelse{html}{\out{<a href="../doc/#1.html"><strong>#1.Rmd</strong></a>}}{\bold{#1.Rmd}}}
+\newcommand{\viglink}{\ifelse{html}{\out{<a href="../doc/#1.html"><strong>#1</strong></a>}}{\bold{#1.html}}}
Modified: pkg/CHNOSZ/tests/testthat/test-diagram.R
===================================================================
--- pkg/CHNOSZ/tests/testthat/test-diagram.R 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/tests/testthat/test-diagram.R 2020-07-18 03:41:12 UTC (rev 566)
@@ -1,7 +1,7 @@
context("diagram")
test_that("expected errors are produced for inconsistent arguments", {
- expect_error(diagram(list()), "'eout' does not look like output from equil\\(\\) or affinity\\(\\)")
+ expect_error(diagram(list()), "'eout' is not the output from")
basis("CHNOS")
species(c("glycine", "alanine"))
a <- affinity()
Modified: pkg/CHNOSZ/vignettes/mklinks.sh
===================================================================
--- pkg/CHNOSZ/vignettes/mklinks.sh 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/vignettes/mklinks.sh 2020-07-18 03:41:12 UTC (rev 566)
@@ -109,8 +109,10 @@
sed -i 's/add.OBIGT()/<a href="..\/html\/add.OBIGT.html">add.OBIGT()<\/a>/g' equilibrium.html
# add links to multi-metal.html 20200716
-sed -i 's/duplex()/<a href="..\/html\/duplex.html">duplex()<\/a>/g' multi-metal.html
-sed -i 's/ratlab()/<a href="..\/html\/util.expression.html">ratlab()<\/a>/g' multi-metal.html
+sed -i 's/affinity()/<a href="..\/html\/affinity.html">affinity()<\/a>/g' multi-metal.html
+sed -i 's/combine()/<a href="..\/html\/combine.html">combine()<\/a>/g' multi-metal.html
sed -i 's/diagram()/<a href="..\/html\/diagram.html">diagram()<\/a>/g' multi-metal.html
sed -i 's/mosaic()/<a href="..\/html\/mosaic.html">mosaic()<\/a>/g' multi-metal.html
sed -i 's/equilibrate()/<a href="..\/html\/equilibrate.html">equilibrate()<\/a>/g' multi-metal.html
+sed -i 's/duplex()/<a href="..\/html\/combine.html">duplex()<\/a>/g' multi-metal.html
+sed -i 's/ratlab()/<a href="..\/html\/util.expression.html">ratlab()<\/a>/g' multi-metal.html
Modified: pkg/CHNOSZ/vignettes/multi-metal.Rmd
===================================================================
--- pkg/CHNOSZ/vignettes/multi-metal.Rmd 2020-07-16 05:47:38 UTC (rev 565)
+++ pkg/CHNOSZ/vignettes/multi-metal.Rmd 2020-07-18 03:41:12 UTC (rev 566)
@@ -72,16 +72,46 @@
This vignette was compiled on `r Sys.Date()` with CHNOSZ version `r sessionInfo()$otherPkgs$CHNOSZ$Version`.
Basic diagrams in CHNOSZ are made for reactions that are *balanced on an element* (see [Equilibrium in CHNOSZ](equilibrium.html)) and therefore represent minerals or aqueous species that all have one element, often a metal, in common.
-This vignette describes some methods for constructing diagrams for elements with multiple metals.
+The package documentation has many examples of diagrams for a single metal appearing in different minerals or complexed with different ligands, but a common request is to make diagrams for multiple metals.
+This vignette describes some methods for constructing diagrams for multi-metal minerals and other multi-element systems.
The methods are **simple overlay**, **mosaic series**, and **secondary balancing**.
## Simple Overlay
Simple overlay refers to independent calculations for two different systems that are displayed on the same diagram.
-It is easy to make such a diagram, but there is no interaction between the systems.
+Although is easy to make such a diagram, there is no interaction between the systems.
-## Mosaic Series
+This example starts with a log*f*~O<sub>2</sub>~--pH base diagram for the C-O-H system then overlays a diagram for S-O-H.
+The second call to `affinity()` uses the argument recall feature, where the arguments are taken from the previous command.
+This allows calculations to be run at the same conditions for a different system.
+This feature is also used in other examples in this vignette.
+```{r overlay, echo = 1:8, eval = FALSE}
+par(mfrow = c(1, 2))
+basis("CHNOS+")
+species(c("CH4", "CO2", "HCO3-", "CO3-2"))
+aC <- affinity(pH = c(0, 14), O2 = c(-75, -60))
+dC <- diagram(aC)
+species(c("H2S", "HS-", "HSO4-", "SO4-2"))
+aS <- affinity(aC) # argument recall
+dS <- diagram(aS, add = TRUE, col = 4, col.names = 4)
+aCS <- combine(dC, dS)
+diagram(aCS)
+legend("topright", legend = lTP(25, 1))
+```
+
+The second diagram is just like the first, except the function `combine()` is used to label the fields with names of species from both systems and we add a legend to indicate the temperature and pressure.
+
+```{r overlay, echo = 9:11, results = "hide", message = FALSE, fig.width = 10, fig.height = 5, out.width = "100%"}
+```
+
+Note that these are predominance diagrams, so they show only the species with highest activity; there is in fact a distribution of activities of aqueous species that is not visible here.
+
+Tip: the names of the fields in the second diagram come from `aCS$species$name`, which are expressions made by combining `aC$names` and `aS$names`.
+If you prefer plain text names without formatting, add `format.names = FALSE` to all of the `diagram()` calls.
+
+## Mosaic Series 1
+
A mosaic diagram shows the effects of changing basis species on the stabilities of minerals.
The Fe-S-O-H system is a common example: the speciation of aqueous sulfur species affects the stabilities of iron oxides and sulfides.
Examples of mosaic diagrams with Fe or other single metals are given elsewhere.
@@ -93,7 +123,7 @@
It is important to put Cu^+^ first so that it will be used as the balance for the reactions with Cu-bearing minerals (which also have Fe).
Pyrite is chosen as the starting Fe-bearing basis species, which will be changed as indicated in `bases2`.
-```{r series1, results = "hide", message = FALSE}
+```{r series1_1, results = "hide", message = FALSE}
logaH2S <- -2
T <- 200
pH <- c(0, 12, 500)
@@ -104,7 +134,7 @@
bases2 <- c("pyrite", "pyrrhotite", "magnetite", "hematite")
```
-Now we calculate affinities for minerals in the Fe-S-O-H system using changing aqueous sulfur species in `bases1`.
+Now we calculate affinities for minerals in the Fe-S-O-H system that take account of the changing aqueous sulfur species in `bases1`.
The result is used to make different layers of the diagram:
1. Water stability region (bounded by the grey area)
@@ -111,7 +141,7 @@
2. Predominance fields for the aqueous S species (italic blue text and dashed lines)
3. Stability areas for the Fe-bearing minerals (black text and solid lines)
-```{r series2, eval = FALSE, echo = 1:6}
+```{r series1_2, eval = FALSE, echo = 1:6}
species(bases2)
mFe <- mosaic(bases1, pH = pH, O2 = O2, T = T)
diagram(mFe$A.bases, lty = 0, names = NA)
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/chnosz -r 566
More information about the CHNOSZ-commits
mailing list