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

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Mon Jul 27 13:44:09 CEST 2020


Author: jedick
Date: 2020-07-27 13:44:09 +0200 (Mon, 27 Jul 2020)
New Revision: 589

Added:
   pkg/CHNOSZ/demo/berman.R
Modified:
   pkg/CHNOSZ/DESCRIPTION
   pkg/CHNOSZ/R/examples.R
   pkg/CHNOSZ/demo/00Index
   pkg/CHNOSZ/demo/AkDi.R
   pkg/CHNOSZ/demo/DEW.R
   pkg/CHNOSZ/demo/NaCl.R
   pkg/CHNOSZ/inst/NEWS.Rd
   pkg/CHNOSZ/man/affinity.Rd
   pkg/CHNOSZ/man/berman.Rd
   pkg/CHNOSZ/man/examples.Rd
   pkg/CHNOSZ/man/mosaic.Rd
   pkg/CHNOSZ/man/subcrt.Rd
   pkg/CHNOSZ/man/swap.basis.Rd
Log:
Suppress warnings about low H2O density in examples and demos


Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/DESCRIPTION	2020-07-27 11:44:09 UTC (rev 589)
@@ -1,6 +1,6 @@
 Date: 2020-07-27
 Package: CHNOSZ
-Version: 1.3.6-62
+Version: 1.3.6-63
 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/examples.R
===================================================================
--- pkg/CHNOSZ/R/examples.R	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/R/examples.R	2020-07-27 11:44:09 UTC (rev 589)
@@ -33,7 +33,7 @@
 demos <- function(which=c("sources", "protein.equil", "affinity", "NaCl", "density", 
   "ORP", "findit", "ionize", "buffer", "protbuff", "glycinate",
   "mosaic", "copper", "arsenic", "solubility", "gold", "contour", "sphalerite",
-  "Shh", "saturation", "adenine", "DEW", "lambda", "TCA", "aluminum",
+  "Shh", "saturation", "adenine", "DEW", "lambda", "berman", "TCA", "aluminum",
   "AkDi", "comproportionation"), save.png=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)) {

Modified: pkg/CHNOSZ/demo/00Index
===================================================================
--- pkg/CHNOSZ/demo/00Index	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/demo/00Index	2020-07-27 11:44:09 UTC (rev 589)
@@ -22,6 +22,7 @@
 adenine         HKF parameters regressed from heat capacity and volume of aqueous adenine
 DEW             Deep Earth Water (DEW) model for high pressures
 lambda          Thermodynamic properties of lambda transition in quartz
+berman          Mineral stabilities in Berman and Helgeson et al. datasets
 TCA             Standard Gibbs energies of steps of the citric acid cycle
 aluminum        Reactions involving Al-bearing minerals
 carboxylase     Rank abundance distribution for RuBisCO and acetyl-CoA carboxylase

Modified: pkg/CHNOSZ/demo/AkDi.R
===================================================================
--- pkg/CHNOSZ/demo/AkDi.R	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/demo/AkDi.R	2020-07-27 11:44:09 UTC (rev 589)
@@ -10,7 +10,7 @@
   # use alternative parameters for H2S? (AD03 Table 1)
   if(altH2S) mod.OBIGT("H2S", state="aq", a=-11.2303, b=12.6104, c=-0.2102)
   # get properties of aq - gas reaction
-  sres <- subcrt(c(name, name), c("aq", "gas"), c(-1, 1), T = T, P = P)
+  sres <- suppressWarnings(subcrt(c(name, name), c("aq", "gas"), c(-1, 1), T = T, P = P))
   # calculate natural logarithm of Henry's constant in mole-fraction units
   ln_KH <- log(1000/18.0153) + log(10) * sres$out$logK
   # plot with units of reciprocal temperature (1000/K)

Modified: pkg/CHNOSZ/demo/DEW.R
===================================================================
--- pkg/CHNOSZ/demo/DEW.R	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/demo/DEW.R	2020-07-27 11:44:09 UTC (rev 589)
@@ -27,9 +27,9 @@
 PT20.0 <- data.frame(P=20000, T=seq(200, 800, 10))
 PT <- rbind(PT0.5, PT1.0, PT2.0, PT5.0, PT10.0, PT20.0)
 # reaction 1: quartz = SiO2(aq) [equivalent to quartz + 3 H2O = Si(OH)4]
-SiO2_logK <- subcrt(c("quartz", "SiO2"), c("cr", "aq"), c(-1, 1), P=PT$P, T=PT$T)$out$logK
+SiO2_logK <- suppressWarnings(subcrt(c("quartz", "SiO2"), c("cr", "aq"), c(-1, 1), P=PT$P, T=PT$T))$out$logK
 # reaction 2: 2 quartz = Si2O4(aq) [equivalent to 2 quartz + 3 H2O = Si2O(OH)6]
-Si2O4_logK <- subcrt(c("quartz", "Si2O4"), c("cr", "aq"), c(-2, 1), P=PT$P, T=PT$T)$out$logK
+Si2O4_logK <- suppressWarnings(subcrt(c("quartz", "Si2O4"), c("cr", "aq"), c(-2, 1), P=PT$P, T=PT$T))$out$logK
 # plot the sum of molalities (== activities) for each pressure
 plot(c(200, 1000), c(-2.5, 0.5), type="n", xlab=axis.label("T"), ylab="log molality")
 for(P in unique(PT$P)) {
@@ -208,3 +208,6 @@
 stopifnot(maxdiff(sres$out[[1]]$loggam, loggamma) < 0.023)
 # if m_star in nonideal() was zero, we could decrease the tolerance here
 #stopifnot(maxdiff(sres$out[[1]]$loggam, loggamma) < 0.004)
+
+# reset OBIGT database
+reset()

Modified: pkg/CHNOSZ/demo/NaCl.R
===================================================================
--- pkg/CHNOSZ/demo/NaCl.R	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/demo/NaCl.R	2020-07-27 11:44:09 UTC (rev 589)
@@ -35,7 +35,7 @@
 coeffs <- c(-1, 1, 1)
 logK <- numeric()
 for(i in 1:length(T)) {
-  s <- subcrt(species, coeffs, T=T[[i]], P=P[[i]])
+  s <- suppressWarnings(subcrt(species, coeffs, T=T[[i]], P=P[[i]]))
   if(i==2) lty <- 3 else lty <- 1
   lines(s$out$T, s$out$logK, lty=lty)
   # keep the calculated values for each experimental condition (excluding Psat)

Added: pkg/CHNOSZ/demo/berman.R
===================================================================
--- pkg/CHNOSZ/demo/berman.R	                        (rev 0)
+++ pkg/CHNOSZ/demo/berman.R	2020-07-27 11:44:09 UTC (rev 589)
@@ -0,0 +1,56 @@
+# CHNOSZ/demo/berman.R
+# moved from berman.Rd 20200727
+### compare mineral stabilities in the Berman and Helgeson datasets
+### on a T - log(K+/H+) diagram, after Sverjensky et al., 1991
+### (doi:10.1016/0016-7037(91)90157-Z)
+library(CHNOSZ)
+
+## set up the system: basis species
+basis(c("K+", "Al+3", "quartz", "H2O", "O2", "H+"))
+# use pH = 0 so that aK+ = aK+/aH+
+basis("pH", 0)
+# load the species
+species(c("K-feldspar", "muscovite", "kaolinite",
+          "pyrophyllite", "andalusite"), "cr")
+
+## start with the data from Helgeson et al., 1978
+add.OBIGT("SUPCRT92")
+# calculate affinities in aK+ - temperature space
+# exceed.Tr: enable calculations above stated temperature limit of pyrophyllite
+res <- 400
+a <- suppressWarnings(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+", 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)
+diagram(a, add = TRUE, names = FALSE, col = 2, lwd = 1.5, lty = 2)
+# the list of references:
+ref1 <- thermo.refs(species()$ispecies)$key
+
+## now use the (default) data from Berman, 1988
+# this resets the thermodynamic database
+# without affecting the basis and species settings
+OBIGT()
+# we can check that we have Berman's quartz
+# and not coesite or some other phase of SiO2
+iSiO2 <- rownames(basis()) == "SiO2"
+stopifnot(info(basis()$ispecies[iSiO2])$name == "quartz")
+# Berman's dataset doesn't have the upper temperature limits,
+# so we don't need exceed.Ttr here
+a <- affinity(`K+` = c(0, 5, res), T = c(200, 650, res), P = 1000, exceed.rhomin = TRUE)
+diagram(a, add = TRUE, names = FALSE, col = 4, lwd = 1.5)
+# the list of references:
+ref2 <- thermo.refs(species()$ispecies)$key
+ref2 <- paste(ref2, collapse = ", ")
+# add legend and title
+legend("top", "low-density region", text.font = 3, bty = "n")
+legend("topleft", describe.property(c("P", "IS"), c(1000, 1)), bty = "n", inset = c(0, 0.1))
+legend("topleft", c("Helgeson et al., 1978", "Berman, 1988 and\nSverjensky et al., 1991"),
+       lty = c(2, 1), lwd = 1.5, col = c(2, 4), bty = "n", inset = c(0.25, 0.08))
+title(main = syslab(c("K2O", "Al2O3", "SiO2", "H2O", "HCl")), line = 1.8)
+title(main = "After Sverjensky et al., 1991",
+      line = 0.3, font.main = 1)
+
+# cleanup for next example
+reset()

Modified: pkg/CHNOSZ/inst/NEWS.Rd
===================================================================
--- pkg/CHNOSZ/inst/NEWS.Rd	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/inst/NEWS.Rd	2020-07-27 11:44:09 UTC (rev 589)
@@ -9,7 +9,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.3.6-61 (2020-07-27)}{
+\section{Changes in CHNOSZ version 1.3.6-62 (2020-07-27)}{
 
   \subsection{MAJOR CHANGES}{
     \itemize{
@@ -132,6 +132,8 @@
       and \samp{bugstab.R} demo in the \pkg{JMDplots} package
       (\url{https://github.com/jedick/JMDplots}).
 
+      \item New demo \samp{berman.R}, extracted from \samp{berman.Rd}.
+
     }
   }
 

Modified: pkg/CHNOSZ/man/affinity.Rd
===================================================================
--- pkg/CHNOSZ/man/affinity.Rd	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/man/affinity.Rd	2020-07-27 11:44:09 UTC (rev 589)
@@ -86,8 +86,9 @@
 affinity()
 # at higher temperature and pressure
 affinity(T=500, P=2000)
-# at 25 temperatures and pressures
-affinity(T=c(500, 1000, 5), P=c(1000, 5000, 5))
+# at 25 temperatures and pressures,
+# some are in the low-density region so we suppress warnings
+suppressWarnings(affinity(T=c(500, 1000, 5), P=c(1000, 5000, 5)))
 # equilibrium constants of formation reactions
 affinity(property="logK")
 # standard molal Gibbs energies of species,

Modified: pkg/CHNOSZ/man/berman.Rd
===================================================================
--- pkg/CHNOSZ/man/berman.Rd	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/man/berman.Rd	2020-07-27 11:44:09 UTC (rev 589)
@@ -94,56 +94,6 @@
 DGrxn <- G_Cs - G_aQz
 stopifnot(all(abs(DGrxn) < 100))
 
-### compare mineral stabilities in the Berman and Helgeson datasets
-### on a T - log(K+/H+) diagram, after Sverjensky et al., 1991
-### (doi:10.1016/0016-7037(91)90157-Z)
-## set up the system: basis species
-basis(c("K+", "Al+3", "quartz", "H2O", "O2", "H+"))
-# use pH = 0 so that aK+ = aK+/aH+
-basis("pH", 0)
-# load the species
-species(c("K-feldspar", "muscovite", "kaolinite",
-          "pyrophyllite", "andalusite"), "cr")
-## start with the data from Helgeson et al., 1978
-add.OBIGT("SUPCRT92")
-# calculate affinities in aK+ - temperature space
-# exceed.Tr: enable calculations above stated temperature limit of pyrophyllite
-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+", 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)
-diagram(a, add = TRUE, names = FALSE, col = "red", lwd = 1.5)
-# the list of references:
-ref1 <- thermo.refs(species()$ispecies)$key
-## now use the (default) data from Berman, 1988
-# this resets the thermodynamic database
-# without affecting the basis and species settings
-OBIGT()
-# we can check that we have Berman's quartz
-# and not coesite or some other phase of SiO2
-iSiO2 <- rownames(basis()) == "SiO2"
-stopifnot(info(basis()$ispecies[iSiO2])$name == "quartz")
-# Berman's dataset doesn't have the upper temperature limits,
-# so we don't need exceed.Ttr here
-a <- affinity(`K+` = c(0, 5, res), T = c(200, 650, res), P = 1000, exceed.rhomin = TRUE)
-diagram(a, add = TRUE, names = FALSE, col = "blue", lwd = 1.5)
-# the list of references:
-ref2 <- thermo.refs(species()$ispecies)$key
-ref2 <- paste(ref2, collapse = ", ")
-# add legend and title
-legend("top", "low-density region", text.font = 3, bty = "n")
-legend("topleft", describe.property(c("P", "IS"), c(1000, 1)), bty = "n")
-legend("left", c(ref1, ref2),
-       lty = c(1, 1), lwd = 1.5, col = c(2, 4), bty = "n")
-title(main = syslab(c("K2O", "Al2O3", "SiO2", "H2O", "HCl")), line = 1.8)
-title(main = "Helgeson and Berman minerals, after Sverjensky et al., 1991",
-      line = 0.3, font.main = 1)
-# cleanup for next example
-reset()
-
 # make a P-T diagram for SiO2 minerals (Ber88 Fig. 4)
 basis(c("SiO2", "O2"), c("cr", "gas"))
 species(c("quartz", "quartz,beta", "coesite"), "cr")

Modified: pkg/CHNOSZ/man/examples.Rd
===================================================================
--- pkg/CHNOSZ/man/examples.Rd	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/man/examples.Rd	2020-07-27 11:44:09 UTC (rev 589)
@@ -17,7 +17,7 @@
     "density", "ORP", "findit", "ionize", "buffer", "protbuff",
     "glycinate", "mosaic", "copper", "arsenic", "solubility", "gold",
     "contour", "sphalerite", "Shh", "saturation",
-    "adenine", "DEW", "lambda", "TCA", "aluminum", "AkDi",
+    "adenine", "DEW", "lambda", "berman", "TCA", "aluminum", "AkDi",
     "comproportionation"),
     save.png=FALSE)
 }
@@ -55,6 +55,7 @@
     \code{adenine} \tab HKF regression of heat capacity and volume of aqueous adenine (Lowe et al., 2017) \cr
     \code{DEW} \tab Deep Earth Water (DEW) model for high pressures (Sverjensky et al., 2014a and 2014b) \cr
     \code{lambda} \tab Effects of lambda transition on thermodynamic properties of quartz (Berman, 1988) \cr
+    \code{berman} \tab ineral stabilities in Berman and Helgeson et al. datasets (Sverjensky et al., 1991) \cr
     \code{TCA} \tab Standard Gibbs energies of the tricarboxylic (citric) acid cycle (Canovas and Shock, 2016) \cr
     \code{aluminum} \tab Reactions involving Al-bearing minerals (Zimmer et al., 2016; Tutolo et al., 2014) \cr
     \code{carboxylase} \tab Rank abundance distribution for RuBisCO and acetyl-CoA carboxylase \cr

Modified: pkg/CHNOSZ/man/mosaic.Rd
===================================================================
--- pkg/CHNOSZ/man/mosaic.Rd	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/man/mosaic.Rd	2020-07-27 11:44:09 UTC (rev 589)
@@ -73,7 +73,7 @@
 basis(c("FeO", "SO4-2", "H2O", "H+", "e-"))
 basis("SO4-2", -6)
 species(c("Fe+2", "Fe+3"), -6)
-species(c("pyrrhotite", "pyrite", "hematite", "magnetite"))
+species(c("pyrrhotite", "pyrite", "hematite", "magnetite"), add = TRUE)
 # the basis species we'll swap through
 bases <- c("SO4-2", "HSO4-", "HS-", "H2S")
 # calculate affinities using the relative abundances of the basis species

Modified: pkg/CHNOSZ/man/subcrt.Rd
===================================================================
--- pkg/CHNOSZ/man/subcrt.Rd	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/man/subcrt.Rd	2020-07-27 11:44:09 UTC (rev 589)
@@ -171,10 +171,8 @@
 # above the T, P limits for the H2O equations of state
 subcrt("alanine", T=c(2250, 2251), P=c(30000, 30001), grid="T")
 # Psat is not defined above the critical point
-\dontrun{
-## (also gives a warning)
-subcrt("alanine", T=seq(0, 5000, by=1000))
-}
+## this also gives a warning
+suppressWarnings(subcrt("alanine", T=seq(0, 5000, by=1000)))
 
 ## minerals with phase transitions
 # compare calculated values of heat capacity of iron with
@@ -200,14 +198,13 @@
 P <- seq(500, 5000, 500)
 # this is like the temperature specification used 
 # in the example by Johnson et al., 1992
-# T <- seq(0, 1000, 100)
-# we use this one to avoid warnings at 0 deg C, 5000 bar
-T <- c(2, seq(100, 1000, 100))
-s <- subcrt(c("andradite", "carbon dioxide", "H2S", "Cu+", "quartz", 
+T <- seq(0, 1000, 100)
+s <- suppressWarnings(subcrt(
+  c("andradite", "carbon dioxide", "H2S", "Cu+", "quartz", 
   "calcite", "chalcopyrite", "H+", "H2O"), 
   coeff=c(-1, -3, -4, -2, 3, 3, 2, 2, 3),
   state=c("cr", "g", "aq", "aq", "cr", "cr", "cr", "aq", "liq"),
-  P=P, T=T, grid="P")
+  P=P, T=T, grid="P"))
 # The results are not identical to SUPCRT92, as CHNOSZ has updated
 # parameters for species e.g. Cu+ from Shock et al., 1997.
 # Check the calculated phase transitions for chalcopyrite

Modified: pkg/CHNOSZ/man/swap.basis.Rd
===================================================================
--- pkg/CHNOSZ/man/swap.basis.Rd	2020-07-27 04:02:39 UTC (rev 588)
+++ pkg/CHNOSZ/man/swap.basis.Rd	2020-07-27 11:44:09 UTC (rev 589)
@@ -57,15 +57,15 @@
 # we have returned to starting point
 stopifnot(all.equal(b1$logact, b3$logact))
 
-## demonstrating the interconvertibility between 
-## chemical potentials of elements and logarithms 
-## of activities of basis species at high temperature
-#basis("CHNOS+")
-#bl1 <- basis()$logact
-#emu <- element.mu(T=100)
-#bl2 <- basis.logact(emu, T=100)
-## note that basis.logact produces a named array
-#stopifnot(all.equal(bl1, as.numeric(bl2)))
+# demonstrating the interconvertibility between 
+# chemical potentials of elements and logarithms 
+# of activities of basis species at high temperature
+basis("CHNOS+")
+bl1 <- basis()$logact
+emu <- element.mu(T=100)
+bl2 <- basis.logact(emu, T=100)
+# note that basis.logact produces a named array
+stopifnot(all.equal(bl1, as.numeric(bl2)))
 
 ## swapping basis species while species are defined
 ## and using numeric species indices



More information about the CHNOSZ-commits mailing list