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

[noreply at r-forge.r-project.org]

Author: jedick
Date: 2019-05-05 09:45:34 +0200 (Sun, 05 May 2019)
New Revision: 460

Modified:
   pkg/CHNOSZ/DESCRIPTION
   pkg/CHNOSZ/R/equilibrate.R
   pkg/CHNOSZ/inst/NEWS
   pkg/CHNOSZ/man/equilibrate.Rd
Log:
equilibrate(): handle output from mosaic by combining activities of basis species and formed species


Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION	2019-05-05 06:49:13 UTC (rev 459)
+++ pkg/CHNOSZ/DESCRIPTION	2019-05-05 07:45:34 UTC (rev 460)
@@ -1,6 +1,6 @@
 Date: 2019-05-05
 Package: CHNOSZ
-Version: 1.3.2-7
+Version: 1.3.2-8
 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/equilibrate.R
===================================================================
--- pkg/CHNOSZ/R/equilibrate.R	2019-05-05 06:49:13 UTC (rev 459)
+++ pkg/CHNOSZ/R/equilibrate.R	2019-05-05 07:45:34 UTC (rev 460)
@@ -9,9 +9,26 @@
 equilibrate <- function(aout, balance=NULL, loga.balance=NULL, 
   ispecies=1:length(aout$values), normalize=FALSE, as.residue=FALSE,
   method=c("boltzmann", "reaction"), tol=.Machine$double.eps^0.25) {
-  ### set up calculation of equilibrium activities of species from the affinities 
-  ### of their formation reactions from basis species at known activities
+  ### calculate equilibrium activities of species from the affinities 
+  ### of their formation reactions from basis species at given activities
   ### split from diagram() 20120925 jmd
+  ## if aout is the output from mosaic(), combine the equilibrium activities of basis species
+  ## and formed species into an object that can be plotted with diagram() 20190505
+  if(aout$fun == "mosaic") {
+    # calculate equilibrium activities of species
+    if(missing(ispecies)) ispecies <- 1:length(aout$A.species$values)
+    if(missing(method)) eqc <- equilibrate(aout$A.species, balance = balance, loga.balance = loga.balance,
+      ispecies = ispecies, normalize = normalize, as.residue = as.residue, tol = tol)
+    else eqc <- equilibrate(aout$A.species, balance = balance, loga.balance = loga.balance,
+      ispecies = ispecies, normalize = normalize, as.residue = as.residue, method = method, tol = tol)
+    # make combined object for basis species and species:
+    # put together the species matrix and logarithms of equilibrium activity
+    eqc$species <- rbind(aout$E.bases[[1]]$species, eqc$species)
+    eqc$loga.equil <- c(aout$E.bases[[1]]$loga.equil, eqc$loga.equil)
+    # we also need to combine 'values' (values of affinity) because diagram() uses this to get the number of species
+    eqc$values <- c(aout$E.bases[[1]]$values, eqc$values)
+    return(eqc)
+  }
   ## number of possible species
   nspecies <- length(aout$values)
   ## get the balancing coefficients

Modified: pkg/CHNOSZ/inst/NEWS
===================================================================
--- pkg/CHNOSZ/inst/NEWS	2019-05-05 06:49:13 UTC (rev 459)
+++ pkg/CHNOSZ/inst/NEWS	2019-05-05 07:45:34 UTC (rev 460)
@@ -1,4 +1,4 @@
-CHANGES IN CHNOSZ 1.3.2-5 (2019-05-05)
+CHANGES IN CHNOSZ 1.3.2-8 (2019-05-05)
 --------------------------------------
 
 BUG FIXES
@@ -32,6 +32,13 @@
 - In mosaic(), include the result of equilibrate() for each group of
   basis species in the output ('E.bases').
 
+- equilibrate(): make it possible to process the output of mosaic().
+  This combines the equilibrium activities of the formed species with
+  the changing basis species to make an object that can be plotted with
+  diagram(). This can be useful for conserving two elements in a system
+  (one in the basis species, the other in the formed species). Thanks
+  to Kirt Robinson for the feature request and test system.
+
 - demo/aluminum.R: add calculations using Si(OH)4 from Akinfiev-Diamond
   model (SiO2 in these reactions is replaced by Si(OH)4 - 2 H2O).
 

Modified: pkg/CHNOSZ/man/equilibrate.Rd
===================================================================
--- pkg/CHNOSZ/man/equilibrate.Rd	2019-05-05 06:49:13 UTC (rev 459)
+++ pkg/CHNOSZ/man/equilibrate.Rd	2019-05-05 07:45:34 UTC (rev 460)
@@ -32,10 +32,11 @@
 \details{
 \code{equilibrate} calculates the chemical activities of species in metastable equilibrium, for constant temperature, pressure and chemical activities of basis species, using specified balancing constraints on reactions between species.
 
-It takes as input \code{aout}, the output from \code{\link{affinity}}, which may be calculated from a multidimensional grid of conditions.
+It takes as input \code{aout}, the output from \code{\link{affinity}}, giving the chemical affinities of formation reaction of each species, which may be calculated on a multidimensional grid of conditions.
+Alternatively, \code{aout} can be the output from \code{\link{mosaic}}, in which case the equilibrium activities of the formed species are calculated and combined with those of the changing basis species to make an object that can be plotted with \code{\link{diagram}}.
+
 The equilibrium chemical activities of species are calculated using either the \code{equil.reaction} or \code{equil.boltzmann} functions, the latter only if the balance is on one mole of species.
 
-\code{aout} contains the chemical affinities of formation reactions of each species of interest.
 \code{equilibrate} needs to be provided constraints on how to balance the reactions representing transformations between the species.
 \code{balance} indicates the balancing component, according to the following scheme:
 
@@ -125,6 +126,22 @@
 logK <- subcrt(c("CO2","H2O","HCO3-","H+"), c(-1, -1, 1, 1), T=25)$out$logK
 # we could decrease tolerance here by increasing res
 stopifnot(all.equal(pKa.equil, -logK, tolerance=1e-2))
+
+# equilibrate with mosaic: balancing on N and C
+loga_N <- -4
+loga_C <- -3
+basis(c("CO2", "NH3", "O2", "H2O", "H+"))
+basis("NH3", loga_N)
+species(c("acetamide", "acetic acid", "acetate"))
+# this calculates equilibrium activities of NH3 and NH4+ for given loga_N
+# and calculates the corresponding affinities of the formed species
+m <- mosaic(c("NH3", "NH4+"), pH = c(0, 14))
+# this calculated equilibrium activities of the formed species for given loga_C
+# and combines them with the activities of the changing basis species (NH3 and NH4+)
+eqc <- equilibrate(m, loga.balance = loga_C)
+diagram(eqc, ylim = c(-10, -2))
+title(main = paste("log(total N in basis species) =", loga_N,
+                   "\nlog(total C in formed species) =", loga_C), font.main = 1)
 }
 
 \references{