[CHNOSZ-commits] r93 - in pkg/CHNOSZ: . R inst tests tests/testthat

Fri Oct 9 14:03:10 CEST 2015

Author: jedick
Date: 2015-10-09 14:03:09 +0200 (Fri, 09 Oct 2015)
New Revision: 93

Added:
   pkg/CHNOSZ/tests/testthat.R
Removed:
   pkg/CHNOSZ/tests/test-all.R
Modified:
   pkg/CHNOSZ/DESCRIPTION
   pkg/CHNOSZ/R/mosaic.R
   pkg/CHNOSZ/inst/NEWS
   pkg/CHNOSZ/tests/testthat/test-mosaic.R
Log:
fix calculation of total activity of basis species in mosaic(, blend=TRUE)


Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================

--- pkg/CHNOSZ/DESCRIPTION	2015-10-08 12:20:21 UTC (rev 92)
+++ pkg/CHNOSZ/DESCRIPTION	2015-10-09 12:03:09 UTC (rev 93)
@@ -1,6 +1,6 @@
-Date: 2015-10-08
+Date: 2015-10-09
 Package: CHNOSZ
-Version: 1.0.5-5
+Version: 1.0.5-6
 Title: Chemical Thermodynamics and Activity Diagrams
 Author: Jeffrey Dick
 Maintainer: Jeffrey Dick <j3ffdick at gmail.com>

Modified: pkg/CHNOSZ/R/mosaic.R
===================================================================
--- pkg/CHNOSZ/R/mosaic.R	2015-10-08 12:20:21 UTC (rev 92)
+++ pkg/CHNOSZ/R/mosaic.R	2015-10-09 12:03:09 UTC (rev 93)
@@ -59,13 +59,13 @@
     # calculate affinities using relative abundances of basis species
     e <- equilibrate(A.bases)
     # what is the total activity of the basis species?
-    loga.tot <- sum(10^unlist(e$loga.equil))
+    a.tot <- Reduce("+", lapply(e$loga.equil, function(x) 10^x))
     for(j in seq_along(affs)) {
       for(i in seq_along(A.species$values)) {
         # start with zero affinity
         if(j==1) A.species$values[[i]][] <- 0
         # add affinity scaled by __relative__ abundance of this basis species
-        A.species$values[[i]] <- A.species$values[[i]] + affs[[j]]$values[[i]] * 10^e$loga.equil[[j]]/loga.tot
+        A.species$values[[i]] <- A.species$values[[i]] + affs[[j]]$values[[i]] * 10^e$loga.equil[[j]]/a.tot
       }
     }
   } else {

Modified: pkg/CHNOSZ/inst/NEWS
===================================================================
--- pkg/CHNOSZ/inst/NEWS	2015-10-08 12:20:21 UTC (rev 92)
+++ pkg/CHNOSZ/inst/NEWS	2015-10-09 12:03:09 UTC (rev 93)
@@ -1,4 +1,4 @@
-CHANGES IN CHNOSZ 1.0.5-8 (2015-10-08)
+CHANGES IN CHNOSZ 1.0.5-6 (2015-10-09)
 --------------------------------------
 
 - Rewrite rho.IAPWS95() to be able to invert density from IAPWS-95
@@ -26,9 +26,10 @@
 - Fix bugs in which.pmax() that prevented proper assignment of
   attributes in output, and functionality for lists of length 1.
 
-- mosaic() now correctly multiplies affinities by _relative_ abundances
-  of basis species when blend=TRUE. Thanks to Grayson Boyer for the bug
-  report that led to this fix.
+- mosaic() now multiplies affinities by _relative_ abundances of basis
+  species when blend=TRUE. Also, fix bug where using blend=TRUE gave
+  wrong results due to incorrect calculation of total activity of basis
+  species. Thanks to Grayson Boyer for the bug reports.
 
 - Add grDevices, graphics, stats, utils to Imports: in DESCRIPTION, and
   add specific importsFrom lines in NAMESPACE.

Deleted: pkg/CHNOSZ/tests/test-all.R
===================================================================
--- pkg/CHNOSZ/tests/test-all.R	2015-10-08 12:20:21 UTC (rev 92)
+++ pkg/CHNOSZ/tests/test-all.R	2015-10-09 12:03:09 UTC (rev 93)
@@ -1,8 +0,0 @@
-library(testthat)
-library(CHNOSZ)
-
-# as fix for https://github.com/hadley/testthat/issues/129, https://github.com/hadley/testthat/issues/144
-# (since we use makeCluster() etc via palply)
-Sys.setenv("R_TESTS" = "") 
-
-test_check("CHNOSZ")

Modified: pkg/CHNOSZ/tests/testthat/test-mosaic.R
===================================================================
--- pkg/CHNOSZ/tests/testthat/test-mosaic.R	2015-10-08 12:20:21 UTC (rev 92)
+++ pkg/CHNOSZ/tests/testthat/test-mosaic.R	2015-10-09 12:03:09 UTC (rev 93)
@@ -14,4 +14,32 @@
   expect_equal(a$values, m2$A.species$values)
 })
 
+test_that("blend=TRUE produces reasonable values", {
+  # a more rigorous test than above. this was failing because loga.tot (actually, a.tot)
+  # was computed incorrectly, by sum'ing an unlist'ed list (the affinities of basis species)
+  # to produce a single value; corrected by using Reduce for addition of vectors/arrays in the list.
+  # example adapted from ?mosaic
+  basis(c("FeO", "SO4-2", "H2O", "H+", "e-"))
+  basis("SO4-2", -6)
+  basis("Eh", -0.15)
+  species(c("hematite", "magnetite"))
+  # the basis species we'll swap through
+  bases <- c("SO4-2", "HSO4-", "HS-", "H2S")         
+  # calculate affinities using the predominant basis species
+  pH <- c(0, 14, 29)
+  m1 <- mosaic(bases, pH=pH)
+  m2 <- mosaic(bases, pH=pH, blend=TRUE)
+  # these species have no S so the results should be the same
+  expect_equal(m1$A.species$values, m2$A.species$values)
+  species(c("pyrrhotite", "pyrite"))
+  # now with S-bearing species ...
+  m3 <- mosaic(bases, pH=pH)
+  m4 <- mosaic(bases, pH=pH, blend=TRUE)
+  # the results are different ...
+  expect_equal(sapply(m3$A.species$values, "[", 13), sapply(m4$A.species$values, "[", 13), tol=1e-1)
+  # but more similar at extreme pH values
+  expect_equal(sapply(m3$A.species$values, "[", 1), sapply(m4$A.species$values, "[", 1), tol=1e-6)
+  expect_equal(sapply(m3$A.species$values, "[", 29), sapply(m4$A.species$values, "[", 29), tol=1e-10)
+})
+
 # TODO: test that basis specifications can be exchanged between bases and bases2 without altering output

Copied: pkg/CHNOSZ/tests/testthat.R (from rev 92, pkg/CHNOSZ/tests/test-all.R)
===================================================================
--- pkg/CHNOSZ/tests/testthat.R	                        (rev 0)
+++ pkg/CHNOSZ/tests/testthat.R	2015-10-09 12:03:09 UTC (rev 93)
@@ -0,0 +1,8 @@
+library(testthat)
+library(CHNOSZ)
+
+# as fix for https://github.com/hadley/testthat/issues/129, https://github.com/hadley/testthat/issues/144
+# (since we use makeCluster() etc via palply)
+Sys.setenv("R_TESTS" = "") 
+
+test_check("CHNOSZ")

