[CHNOSZ-commits] r593 - in pkg/CHNOSZ: . man tests/testthat vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Aug 11 02:57:48 CEST 2020
Author: jedick
Date: 2020-08-11 02:57:48 +0200 (Tue, 11 Aug 2020)
New Revision: 593
Modified:
pkg/CHNOSZ/DESCRIPTION
pkg/CHNOSZ/man/berman.Rd
pkg/CHNOSZ/man/buffer.Rd
pkg/CHNOSZ/man/diagram.Rd
pkg/CHNOSZ/tests/testthat/test-diagram.R
pkg/CHNOSZ/vignettes/multi-metal.Rmd
Log:
Remove stopifnot() from examples, round 2
Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION 2020-07-28 12:06:41 UTC (rev 592)
+++ pkg/CHNOSZ/DESCRIPTION 2020-08-11 00:57:48 UTC (rev 593)
@@ -1,6 +1,6 @@
-Date: 2020-07-28
+Date: 2020-08-11
Package: CHNOSZ
-Version: 1.3.6-66
+Version: 1.3.6-67
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/man/berman.Rd
===================================================================
--- pkg/CHNOSZ/man/berman.Rd 2020-07-28 12:06:41 UTC (rev 592)
+++ pkg/CHNOSZ/man/berman.Rd 2020-08-11 00:57:48 UTC (rev 593)
@@ -92,7 +92,7 @@
# that is close to the univariant curve (Ber88 Fig. 4),
# so the difference in G is close to 0
DGrxn <- G_Cs - G_aQz
-stopifnot(all(abs(DGrxn) < 100))
+all(abs(DGrxn) < 100) # TRUE
# make a P-T diagram for SiO2 minerals (Ber88 Fig. 4)
basis(c("SiO2", "O2"), c("cr", "gas"))
@@ -113,10 +113,10 @@
thermo("opt$Berman" = file.path(datadir, "BA96_berman.csv"))
Gex_BA96 <- subcrt(species, coeffs, T=seq(600, 1500, 50), P=1)$out$G
# Ber88 is lower than BA96 at low T
-stopifnot((Gex_BA96 - Gex_Ber88)[1] > 0)
+(Gex_BA96 - Gex_Ber88)[1] > 0 # TRUE
# the curves cross at about 725 deg C (BA96 Fig. 8)
# (actually, in our calculation they cross closer to 800 deg C)
-stopifnot(T[which.min(abs(Gex_BA96 - Gex_Ber88))] == 800)
+T[which.min(abs(Gex_BA96 - Gex_Ber88))] # 800
# reset the database (thermo()$opt$E.units, thermo()$OBIGT, and thermo()$opt$Berman)
reset()
}
Modified: pkg/CHNOSZ/man/buffer.Rd
===================================================================
--- pkg/CHNOSZ/man/buffer.Rd 2020-07-28 12:06:41 UTC (rev 592)
+++ pkg/CHNOSZ/man/buffer.Rd 2020-08-11 00:57:48 UTC (rev 593)
@@ -42,6 +42,10 @@
This function is compatible with systems of proteins, but note that for buffers \emph{made} of proteins the buffer calculations presently use whole protein formulas (instead of residue equivalents) and consider nonionized proteins only.
}
+\seealso{
+ \code{link{diagram}} with \code{type} set to the name of a basis species solves for the activity of the basis species.
+}
+
\examples{
\dontshow{reset()}
## list the buffers
@@ -57,8 +61,7 @@
# what activity of acetic acid are we using?
print(mod.buffer("AC"))
# return the activity of CO2
-(logaCO2 <- affinity(return.buffer=TRUE)$CO2)
-stopifnot(all.equal(logaCO2, -7.05752136))
+affinity(return.buffer=TRUE)$CO2 # -7.057521
# as a function of oxygen fugacity
affinity(O2=c(-85,-70,4),return.buffer=TRUE)
# as a function of logfO2 and temperature
@@ -66,27 +69,25 @@
# change the activity of species in the buffer
mod.buffer("AC",logact=-10)
affinity(O2=c(-85,-70,4),T=c(25,100,4),return.buffer=TRUE)
-# see below for a different strategy using the
-# 'type' argument of diagram
## buffer made of three species
## Pyrite-Pyrrhotite-Magnetite (PPM)
# specify basis species and initial activities
-basis(c("FeS2","H2S","O2","H2O"),c(0,-10,-50,0))
+basis(c("FeS2", "H2S", "O2", "H2O"), c(0, -10, -50, 0))
# note that the affinity of formation of pyrite,
# which corresponds to FeS2 in the basis, is zero
-species(c("pyrite","pyrrhotite","magnetite"))
-affinity(T=c(200,400,11),P=2000)$values
+species(c("pyrite", "pyrrhotite", "magnetite"))
+affinity(T = c(200, 400, 11), P = 2000)$values
# setup H2S and O2 to be buffered by PPM
-basis(c("H2S","O2"),c("PPM","PPM"))
+basis(c("H2S", "O2"), c("PPM", "PPM"))
# inspect values of H2S activity and O2 fugacity
-affinity(T=c(200, 400, 11), P=2000, return.buffer=TRUE, exceed.Ttr=TRUE)
-# now, the affinities of formation reactions of
-# species in the buffer are all equal to zero
-print(a <- affinity(T=c(200, 400, 11), P=2000,
- exceed.Ttr=TRUE)$values)
-for(i in 1:length(a)) stopifnot(isTRUE(
- all.equal(as.numeric(a[[i]]),rep(0,length(a[[i]])))))
+affinity(T = c(200, 400, 11), P = 2000, return.buffer = TRUE, exceed.Ttr = TRUE)
+# calculate affinities of formation reactions of species in the buffer
+a <- affinity(T = c(200, 400, 11), P = 2000, exceed.Ttr = TRUE)$values
+# the affinities for species in the buffer are all equal to zero
+all.equal(as.numeric(a[[1]]), rep(0, 11)) # TRUE
+all.equal(as.numeric(a[[2]]), rep(0, 11)) # TRUE
+all.equal(as.numeric(a[[3]]), rep(0, 11)) # TRUE
## buffer made of one species: show values of logfO2 on an
## Eh-pH diagram; after Garrels, 1960, Figure 6
@@ -127,34 +128,6 @@
# consider more oxidizing conditions
mod.buffer("CO2-AC",logact=c(0,-10))
affinity(return.buffer=TRUE)
-
-# one can solve for the logarithm of activity of a
-# basis species using the 'type' argument of diagram
-basis("CHNOS")
-basis("CO2", 999)
-species("acetic acid", -3)
-a <- affinity(O2=c(-85, -70, 4), T=c(25, 100, 4))
-# write a title with formulas and subscripts
-lCO2 <- axis.label("CO2")
-main <- substitute(a~~b~~c,list(a=lCO2, b="buffered by",
- c="acetic acid"))
-d <- diagram(a, type="CO2", main=main)
-species(1, -10)
-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=TRUE)
-ltext <- c(as.expression(lAC), -3, -10)
-lty <- c(NA, 1, 2)
-legend("topright", legend=ltext, lty=lty, bg="white")
-# do return.buffer and diagram(type=...) give the same results?
-and <- as.numeric(d$plotvals[[1]])
-basis("CO2", "AC")
-mod.buffer("AC", logact=-10)
-a.buffer <- affinity(O2=c(-85, -70, 4), T=c(25, 100, 4),
- return.buffer=TRUE)
-ana <- as.numeric(unlist(a.buffer[[1]]))
-stopifnot(all.equal(ana, and))
}
\references{
Modified: pkg/CHNOSZ/man/diagram.Rd
===================================================================
--- pkg/CHNOSZ/man/diagram.Rd 2020-07-28 12:06:41 UTC (rev 592)
+++ pkg/CHNOSZ/man/diagram.Rd 2020-08-11 00:57:48 UTC (rev 593)
@@ -241,8 +241,6 @@
"goethite", "melanterite", "pyrite"))
a <- affinity(pH=c(-1, 4, 256), pe=c(-5, 23, 256))
d <- diagram(a, main="Fe-S-O-H, after Majzlan et al., 2006")
-# the first four species show up in order near pe=15
-stopifnot(all.equal(unique(d$predominant[, 183]), 1:4))
water.lines(d, lwd=2)
text(3, 22, describe.basis(thermo()$basis[2:3,], digits=2, oneline=TRUE))
text(3, 21, describe.property(c("T", "P"), c(25, 1), oneline=TRUE))
@@ -281,10 +279,6 @@
Ttp <- a$vals[[1]][tp[1, 2]]
Ptp <- rev(a$vals[[2]])[tp[1, 1]]
points(Ttp, Ptp, pch=10, cex=5)
-# some testing of the overall geometry
-stopifnot(species()$name[d$predominant[1, 1]]=="andalusite")
-stopifnot(species()$name[d$predominant[1, 101]]=="kyanite")
-stopifnot(species()$name[d$predominant[99, 101]]=="sillimanite")
}
\references{
Modified: pkg/CHNOSZ/tests/testthat/test-diagram.R
===================================================================
--- pkg/CHNOSZ/tests/testthat/test-diagram.R 2020-07-28 12:06:41 UTC (rev 592)
+++ pkg/CHNOSZ/tests/testthat/test-diagram.R 2020-08-11 00:57:48 UTC (rev 593)
@@ -77,14 +77,6 @@
expect_equal(d3$predominant, d4$predominant)
})
-test_that("diagram() handles 2D plots with different x and y resolution and warns for >1 species in contour plot", {
- basis("CHNOS")
- species(c("alanine", "glycine", "serine", "methionine"))
- a <- affinity(T=c(0, 200, 6), O2=c(-90, -60, 5))
- # now the warning is invokes next to the actual plotting, so no warning is produced with plot.it=FALSE 20180315
- #expect_warning(diagram(a, type="CO2", plot.it=FALSE), "showing only first species in 2-D property diagram")
-})
-
test_that("NaN values from equilibrate() are preserved (as NA in predominance calculation)", {
# example provided by Grayson Boyer 20170411
basis(c("H2", "O2", "CO2"), c(-7.19, -60, -2.65))
@@ -114,3 +106,34 @@
# d <- diagram(a, fill="heat", xlim=c(-50, -20), ylim=c(-90, -50), plot.it=FALSE)
# expect_equal(sum(is.na(d$namesx)), 2)
#})
+
+test_that("P-T diagram has expected geometry", {
+ # modifified from kayanite-sillimanite-andalusite example in ?diagram 20200811
+ basis(c("corundum", "quartz", "oxygen"))
+ species(c("kyanite", "sillimanite", "andalusite"))
+ a <- affinity(T = c(200, 900, 50), P = c(0, 9000, 51), exceed.Ttr = TRUE)
+ d <- diagram(a, plot.it = FALSE)
+ expect_equal(species()$name[d$predominant[1, 1]], "andalusite")
+ expect_equal(species()$name[d$predominant[1, 51]], "kyanite")
+ expect_equal(species()$name[d$predominant[50, 51]], "sillimanite")
+})
+
+test_that("diagram(type = .) and affinity(return.buffer = TRUE) give the same results", {
+ # extracted from ?buffer 20200811
+ O2 <- c(-85, -70, 4)
+ T <- c(25, 100, 4)
+ basis("CHNOS")
+ basis("CO2", 999)
+ species("acetic acid", -3)
+ species(1, -10)
+ a <- affinity(O2 = O2, T = T)
+ d <- diagram(a, type = "CO2", plot.it = FALSE)
+ and <- as.numeric(d$plotvals[[1]])
+ # now do the calculation with affinity(return.buffer = TRUE)
+ basis("CO2", "AC")
+ mod.buffer("AC", logact = -10)
+ a.buffer <- affinity(O2 = O2, T = T, return.buffer = TRUE)
+ ana <- as.numeric(unlist(a.buffer[[1]]))
+ expect_equal(ana, and)
+})
+
Modified: pkg/CHNOSZ/vignettes/multi-metal.Rmd
===================================================================
--- pkg/CHNOSZ/vignettes/multi-metal.Rmd 2020-07-28 12:06:41 UTC (rev 592)
+++ pkg/CHNOSZ/vignettes/multi-metal.Rmd 2020-08-11 00:57:48 UTC (rev 593)
@@ -835,7 +835,7 @@
```
The diagram shows the expected ionization of acetic acid and NH~3~ at different pHs.
-The appearance of acetamide (CH~3~CONH~2~) is a consequence of the interaction between the N-bearing and C-bearing species, and is analogous to the formation of a bimetallic complex.
+The predicted appearance of acetamide (CH~3~CONH~2~) is a consequence of the interaction between the N-bearing and C-bearing species, and is analogous to the formation of a bimetallic complex.
*Thanks to Kirt Robinson for the feature request and test case used in this example.*
More information about the CHNOSZ-commits
mailing list