[CHNOSZ-commits] r704 - in pkg/CHNOSZ: . R inst inst/tinytest man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Tue Feb 8 12:29:21 CET 2022


Author: jedick
Date: 2022-02-08 12:29:21 +0100 (Tue, 08 Feb 2022)
New Revision: 704

Modified:
   pkg/CHNOSZ/DESCRIPTION
   pkg/CHNOSZ/R/basis.R
   pkg/CHNOSZ/R/swap.basis.R
   pkg/CHNOSZ/inst/NEWS.Rd
   pkg/CHNOSZ/inst/TODO
   pkg/CHNOSZ/inst/tinytest/test-basis.R
   pkg/CHNOSZ/man/basis.Rd
Log:
Add an 'add' argument to basis()


Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION	2022-02-08 09:42:37 UTC (rev 703)
+++ pkg/CHNOSZ/DESCRIPTION	2022-02-08 11:29:21 UTC (rev 704)
@@ -1,6 +1,6 @@
 Date: 2022-02-08
 Package: CHNOSZ
-Version: 1.4.1-29
+Version: 1.4.1-30
 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/R/basis.R
===================================================================
--- pkg/CHNOSZ/R/basis.R	2022-02-08 09:42:37 UTC (rev 703)
+++ pkg/CHNOSZ/R/basis.R	2022-02-08 11:29:21 UTC (rev 704)
@@ -1,10 +1,13 @@
 # CHNOSZ/basis.R
 # Set up the basis species of a thermodynamic system
 
-basis <- function(species=NULL, state=NULL, logact=NULL, delete=FALSE) {
+basis <- function(species = NULL, state = NULL, logact = NULL, delete = FALSE, add = FALSE) {
+
+  # Get the current thermo and basis settings
   thermo <- get("thermo", CHNOSZ)
   oldbasis <- thermo$basis
-  ## Delete the basis species if requested
+
+  # Delete the basis species if requested
   if(delete | identical(species, "")) {
     thermo$basis <- NULL
     thermo$species <- NULL
@@ -11,21 +14,26 @@
     assign("thermo", thermo, CHNOSZ)
     return(invisible(oldbasis))
   }
-  ## Return the basis definition if requested
+
+  # Return the basis definition if requested
   if(is.null(species)) return(oldbasis)
-  ## From now on we need something to work with
-  if(length(species)==0) stop("species argument is empty")
+
+  # From now on we need something to work with
+  if(length(species) == 0) stop("species argument is empty")
+
   # Is the species one of the preset keywords?
   if(species[1] %in% preset.basis()) return(preset.basis(species[1]))
   # The species names/formulas have to be unique
-  if(!length(unique(species))==length(species)) stop("species names are not unique")
-  ## Processing 'state' and 'logact' arguments
+  if(!length(unique(species)) == length(species)) stop("species names are not unique")
+
+  # Processing 'state' and 'logact' arguments
   # They should be same length as species
-  if(!is.null(state)) state <- rep(state, length.out=length(species))
-  if(!is.null(logact)) logact <- rep(logact, length.out=length(species))
+  if(!is.null(state)) state <- rep(state, length.out = length(species))
+  if(!is.null(logact)) logact <- rep(logact, length.out = length(species))
+
   # Results should be identical for
-  # basis(c('H2O','CO2','H2'), rep('aq',3), c(0,-3,-3))
-  # basis(c('H2O','CO2','H2'), c(0,-3,-3), rep('aq',3))
+  # basis(c("H2O", "CO2", "H2"), rep("aq", 3), c(0, -3, -3))
+  # basis(c("H2O", "CO2", "H2"), c(0, -3, -3), rep("aq", 3))
   # First of all, do we have a third argument?
   if(!is.null(logact)) {
     # Does the 3rd argument look like states?
@@ -42,32 +50,38 @@
       state <- NULL
     }
   }
-  ## Processing 'species' argument
+
+  # Processing 'species' argument
   # pH transformation
   if("pH" %in% species) {
-    logact[species=="pH"] <- -logact[species=="pH"]
-    if(!is.null(logact)) species[species=="pH"] <- "H+"
+    logact[species == "pH"] <- -logact[species == "pH"]
+    if(!is.null(logact)) species[species == "pH"] <- "H+"
   }
   # Eh and pe transformations
   if("pe" %in% species) {
-    logact[species=="pe"] <- -logact[species=="pe"]
-    if(!is.null(logact)) species[species=="pe"] <- "e-"
+    logact[species == "pe"] <- -logact[species == "pe"]
+    if(!is.null(logact)) species[species == "pe"] <- "e-"
   }
   if("Eh" %in% species) {
     # 20090209 Should be careful with this conversion as it's only for 25 degC
     # To be sure, just don't call species("Eh")
-    if(!is.null(logact)) logact[species=="Eh"] <- -convert(logact[species=="Eh"],"pe")
-    species[species=="Eh"] <- "e-"
+    if(!is.null(logact)) logact[species == "Eh"] <- -convert(logact[species == "Eh"],"pe")
+    species[species == "Eh"] <- "e-"
   }
-  ## If all species are in the existing basis definition, 
-  ## *and* at least one of state or logact is not NULL
-  ## modify the states and/or logacts of the existing basis species
+
+  # If all species are in the existing basis definition, 
+  # *and* at least one of state or logact is not NULL,
+  # *and* add is not TRUE,
+  # then modify the states and/or logacts of the existing basis species
   if(all(species %in% rownames(oldbasis)) | all(species %in% oldbasis$ispecies)) 
     if(!is.null(state) | !is.null(logact))
-      return(mod.basis(species, state, logact))
-  ## We're on to making a new basis definition
-  # use default logacts if they aren't present
+      if(!add)
+        return(mod.basis(species, state, logact))
+
+  # If we get to this point, we're making a new basis definition or adding to an existing one
+  # Use default logacts if they aren't present
   if(is.null(logact)) logact <- rep(0, length(species))
+
   # If species argument is numeric, it's species indices
   if(is.numeric(species[1])) {
     ispecies <- species
@@ -82,9 +96,45 @@
     # We don't accept any of those
     if(is.list(ispecies)) ina <- ina | sapply(ispecies,length) > 1
   }
-  if(any(ina)) stop(paste("species not available:", paste(species[ina], "(", state[ina], ")", sep="", collapse=" ")))
+  if(any(ina)) stop(paste("species not available:", paste(species[ina], "(", state[ina], ")", sep = "", collapse = " ")))
+
+  # Are we adding a species to an existing basis definition? 20220208
+  if(add) {
+    # Check that the species aren't already in the basis definition 20220208
+    ihave <- ispecies %in% oldbasis$ispecies
+    if(any(ihave)) {
+      if(length(ihave) == 1) stop(paste("this species is already in the basis definition:", species[ihave]))
+      stop(paste("these species are already in the basis definition:", paste(species[ihave], collapse = " ")))
+    }
+    # Append the new basis species
+    ispecies <- c(oldbasis$ispecies, ispecies)
+    logact <- c(oldbasis$logact, logact)
+  }
   # Load new basis species
-  return(put.basis(ispecies, logact))
+  newthermo <- put.basis(ispecies, logact)
+  newbasis <- newthermo$basis
+
+  if(add) {
+    oldspecies <- thermo$species
+    # If there are formed species, we need to update their formation reactions
+    if(!is.null(oldspecies)) {
+      nold <- nrow(oldbasis)
+      nnew <- nrow(newbasis)
+      nspecies <- nrow(oldspecies)
+      # Stoichiometric coefficients have zero for the added basis species
+      stoich <- cbind(oldspecies[, 1:nold], matrix(0, nspecies, nnew - nold))
+      colnames(stoich) <- rownames(newbasis)
+      # Complete species data frame
+      newspecies <- cbind(stoich, oldspecies[, -(1:nold)])
+      newthermo$species <- newspecies
+      assign("thermo", newthermo, CHNOSZ)
+    }
+  } else {
+    # Remove the formed species since there's no guarantee the new basis includes all their elements
+    species(delete = TRUE)
+  }
+
+  newbasis
 }
 
 ### unexported functions ###
@@ -126,9 +176,6 @@
   # Ready to assign to the global thermo object
   thermo$basis <- comp
   assign("thermo", thermo, CHNOSZ)
-  # Remove the species since there's no guarantee the new basis includes all their elements
-  species(delete=TRUE)
-  return(thermo$basis)
 }
 
 # Modify the states or logact values in the existing basis definition

Modified: pkg/CHNOSZ/R/swap.basis.R
===================================================================
--- pkg/CHNOSZ/R/swap.basis.R	2022-02-08 09:42:37 UTC (rev 703)
+++ pkg/CHNOSZ/R/swap.basis.R	2022-02-08 11:29:21 UTC (rev 704)
@@ -93,7 +93,7 @@
   # try to load the new basis species
   ispecies <- oldbasis$ispecies
   ispecies[ib] <- ispecies2
-  newbasis <- put.basis(ispecies)
+  newbasis <- put.basis(ispecies)$basis
   # now deal with the activities
   if(!all(can.be.numeric(oldbasis$logact))) {
     # if there are any buffers, just set the old activities

Modified: pkg/CHNOSZ/inst/NEWS.Rd
===================================================================
--- pkg/CHNOSZ/inst/NEWS.Rd	2022-02-08 09:42:37 UTC (rev 703)
+++ pkg/CHNOSZ/inst/NEWS.Rd	2022-02-08 11:29:21 UTC (rev 704)
@@ -10,7 +10,7 @@
 \newcommand{\s}{\ifelse{latex}{\eqn{_{#1}}}{\ifelse{html}{\out{<sub>#1</sub>}}{#1}}}
 \newcommand{\S}{\ifelse{latex}{\eqn{^{#1}}}{\ifelse{html}{\out{<sup>#1</sup>}}{^#1}}}
 
-\section{Changes in CHNOSZ version 1.4.1-29 (2022-02-08)}{
+\section{Changes in CHNOSZ version 1.4.1-30 (2022-02-08)}{
 
   \subsection{PLANNED API CHANGE}{
     \itemize{
@@ -75,6 +75,11 @@
   \subsection{OTHER CHANGES}{
     \itemize{
 
+      \item Tests are now run using the \CRANpkg{tinytest} package.
+
+      \item Add an \strong{add} argument to \code{basis()} to allow adding a
+      species to an existing set of basis species.
+
       \item The \code{AkDi()} function has been renamed to \code{AD()}, and all
       variables and data files likewise use the acronym \acronym{AD}. In
       particular, the Akinfiev-Diamond model is activated for an aqueous
@@ -83,8 +88,6 @@
       \item Names of functions, variables, and data files now use capitalized
       \samp{Berman}, not \samp{berman}.
 
-      \item Tests are now run using the \CRANpkg{tinytest} package.
-
       \item In \file{vignettes/multi-metal.Rmd} (\dQuote{Diagrams with multiple
         metals}), add a link to the associated paper
       (\href{https://doi.org/10.1016/j.acags.2021.100059}{Dick, 2021}).

Modified: pkg/CHNOSZ/inst/TODO
===================================================================
--- pkg/CHNOSZ/inst/TODO	2022-02-08 09:42:37 UTC (rev 703)
+++ pkg/CHNOSZ/inst/TODO	2022-02-08 11:29:21 UTC (rev 704)
@@ -39,10 +39,6 @@
 
 - Merge SK95.csv into SLOP98.csv
 
-[20210722]
-
-- Add pyrobitumen properties from Helgeson et al. (2009)
-
 [20220203]
 
 - Allow numeric values for bases in mosaic() (test-mix.R)

Modified: pkg/CHNOSZ/inst/tinytest/test-basis.R
===================================================================
--- pkg/CHNOSZ/inst/tinytest/test-basis.R	2022-02-08 09:42:37 UTC (rev 703)
+++ pkg/CHNOSZ/inst/tinytest/test-basis.R	2022-02-08 11:29:21 UTC (rev 704)
@@ -37,3 +37,19 @@
 basis("CHNOS+")  # first basis species is CO2(aq)
 expect_equal(basis("CO2", "gas")$state[1], "gas", info = info)
 expect_equal(basis("CO2", "aq")$state[1], "aq", info = info)
+
+# 20220208
+info <- "Adding basis species to an existing definition works"
+basis(c("iron", "oxygen", "H2O", "H+"), c(0, -80, 1, -7))
+species(c("Fe", "Fe+2", "Fe+2"), c(0, -6, -6))
+a1 <- affinity(pH = c(0, 14))
+expect_error(basis(c("iron"), add = TRUE), "this species is already in the basis definition")
+expect_error(basis(c("iron", "oxygen", "H2O", "H+"), add = TRUE), "these species are already in the basis definition")
+expect_error(basis("Fe+2", add = TRUE), "the number of basis species is greater than the number of elements and charge")
+expect_silent(newbasis <- basis(c("H2S", "Cl-"), c(-5, -3), add = TRUE))
+expect_equal(newbasis$logact[5:6], c(-5, -3))
+newspecies <- species()
+expect_equal(colnames(newspecies)[5:6], c("H2S", "Cl-"))
+expect_equal(newspecies$H2S, numeric(3))
+a2 <- affinity(pH = c(0, 14))
+expect_identical(a1$values, a2$values)

Modified: pkg/CHNOSZ/man/basis.Rd
===================================================================
--- pkg/CHNOSZ/man/basis.Rd	2022-02-08 09:42:37 UTC (rev 703)
+++ pkg/CHNOSZ/man/basis.Rd	2022-02-08 11:29:21 UTC (rev 704)
@@ -7,7 +7,7 @@
 }
 
 \usage{
-  basis(species = NULL, state = NULL, logact = NULL, delete = FALSE)
+  basis(species = NULL, state = NULL, logact = NULL, delete = FALSE, add = FALSE)
 }
 
 \arguments{
@@ -14,7 +14,8 @@
   \item{species}{character, names or formulas of species, or numeric, indices of species}
   \item{state}{character, physical states or names of buffers}
   \item{logact}{numeric, logarithms of activities or fugacities}
-  \item{delete}{logical, delete the current basis species definition?}
+  \item{delete}{logical, delete the current basis definition?}
+  \item{add}{logical, add species to the current basis definition?}
 }
 
 \details{
@@ -32,8 +33,12 @@
 In case \samp{pH}, \samp{pe} or \samp{Eh} is named, the logarithm of activity of the basis species is converted from these values.
 For example, a value of 7 for pH is stored as a logarithm of activity of -7.
 
-Whenever \code{basis} is called with NULL values of both \code{state} and \code{logact}, the new set of species, if they are a valid basis set, completely replaces any existing basis definition.
-If this occurs, any existing species definition (created by the \code{species} function) is deleted.
+If \code{add} is TRUE, then the function attempts to \emph{add} the indicated \code{species} to the basis definition.
+This only works if the enlarged set of species is a valid basis set as described above.
+If the formed \code{\link{species}} are currently defined, their formation reactions are modified accordingly (with zeroes for the newly added basis species).
+
+If \code{add} is FALSE, and if \code{basis} is called with NULL values of both \code{state} and \code{logact}, the new set of species, if they are a valid basis set, completely replaces any existing basis definition.
+When this occurs, any existing species definition (created by the \code{species} function) is deleted.
 Call \code{basis} with \code{delete} set to TRUE or \code{species} set to \samp{""} to clear the basis definition and that of the \code{\link{species}}, if present.
 
 If the value of \code{basis} is one of the keywords in the following table, the corresponding set of basis species is loaded, and their activities are given preset values.



More information about the CHNOSZ-commits mailing list