[CHNOSZ-commits] r684 - in pkg/CHNOSZ: . R inst inst/tinytest man tests vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon Jan 31 14:27:27 CET 2022
Author: jedick
Date: 2022-01-31 14:27:27 +0100 (Mon, 31 Jan 2022)
New Revision: 684
Added:
pkg/CHNOSZ/inst/tinytest/
pkg/CHNOSZ/tests/tinytest.R
Removed:
pkg/CHNOSZ/R/util.test.R
pkg/CHNOSZ/man/util.test.Rd
pkg/CHNOSZ/tests/testthat.R
pkg/CHNOSZ/tests/testthat/
Modified:
pkg/CHNOSZ/DESCRIPTION
pkg/CHNOSZ/NAMESPACE
pkg/CHNOSZ/R/diagram.R
pkg/CHNOSZ/R/util.character.R
pkg/CHNOSZ/R/util.formula.R
pkg/CHNOSZ/inst/CHECKLIST
pkg/CHNOSZ/inst/NEWS.Rd
pkg/CHNOSZ/inst/tinytest/test-DEW.R
pkg/CHNOSZ/inst/tinytest/test-EOSregress.R
pkg/CHNOSZ/inst/tinytest/test-IAPWS95.R
pkg/CHNOSZ/inst/tinytest/test-add.protein.R
pkg/CHNOSZ/inst/tinytest/test-affinity.R
pkg/CHNOSZ/inst/tinytest/test-basis.R
pkg/CHNOSZ/inst/tinytest/test-berman.R
pkg/CHNOSZ/inst/tinytest/test-buffer.R
pkg/CHNOSZ/inst/tinytest/test-diagram.R
pkg/CHNOSZ/inst/tinytest/test-eos.R
pkg/CHNOSZ/inst/tinytest/test-equilibrate.R
pkg/CHNOSZ/inst/tinytest/test-findit.R
pkg/CHNOSZ/inst/tinytest/test-info.R
pkg/CHNOSZ/inst/tinytest/test-ionize.aa.R
pkg/CHNOSZ/inst/tinytest/test-ionize.aa_pK.R
pkg/CHNOSZ/inst/tinytest/test-logmolality.R
pkg/CHNOSZ/inst/tinytest/test-makeup.R
pkg/CHNOSZ/inst/tinytest/test-mix.R
pkg/CHNOSZ/inst/tinytest/test-mod.OBIGT.R
pkg/CHNOSZ/inst/tinytest/test-mosaic.R
pkg/CHNOSZ/inst/tinytest/test-nonideal.R
pkg/CHNOSZ/inst/tinytest/test-objective.R
pkg/CHNOSZ/inst/tinytest/test-protein.info.R
pkg/CHNOSZ/inst/tinytest/test-recalculate.R
pkg/CHNOSZ/inst/tinytest/test-retrieve.R
pkg/CHNOSZ/inst/tinytest/test-revisit.R
pkg/CHNOSZ/inst/tinytest/test-solubility.R
pkg/CHNOSZ/inst/tinytest/test-species.R
pkg/CHNOSZ/inst/tinytest/test-subcrt.R
pkg/CHNOSZ/inst/tinytest/test-swap.basis.R
pkg/CHNOSZ/inst/tinytest/test-thermo.R
pkg/CHNOSZ/inst/tinytest/test-util.R
pkg/CHNOSZ/inst/tinytest/test-util.affinity.R
pkg/CHNOSZ/inst/tinytest/test-util.data.R
pkg/CHNOSZ/inst/tinytest/test-util.list.R
pkg/CHNOSZ/inst/tinytest/test-util.program.R
pkg/CHNOSZ/inst/tinytest/test-util.seq.R
pkg/CHNOSZ/inst/tinytest/test-water.R
pkg/CHNOSZ/inst/tinytest/test-water.lines.R
pkg/CHNOSZ/man/CHNOSZ-package.Rd
pkg/CHNOSZ/man/NaCl.Rd
pkg/CHNOSZ/man/diagram.Rd
pkg/CHNOSZ/man/nonideal.Rd
pkg/CHNOSZ/man/thermo.Rd
pkg/CHNOSZ/vignettes/anintro.Rmd
Log:
Use tinytest to run tests
Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/DESCRIPTION 2022-01-31 13:27:27 UTC (rev 684)
@@ -1,6 +1,6 @@
-Date: 2021-09-27
+Date: 2022-01-29
Package: CHNOSZ
-Version: 1.4.1-10
+Version: 1.4.1-11
Title: Thermodynamic Calculations and Diagrams for Geochemistry
Authors at R: c(
person("Jeffrey", "Dick", , "j3ffdick at gmail.com", role = c("aut", "cre"),
@@ -9,7 +9,7 @@
Author: Jeffrey Dick [aut, cre] (0000-0002-0687-5890)
Maintainer: Jeffrey Dick <j3ffdick at gmail.com>
Depends: R (>= 3.1.0)
-Suggests: testthat, knitr, rmarkdown, tufte
+Suggests: tinytest, knitr, rmarkdown, tufte
Imports: grDevices, graphics, stats, utils
Description: An integrated set of tools for thermodynamic calculations in
aqueous geochemistry and geobiochemistry. Functions are provided for writing
Modified: pkg/CHNOSZ/NAMESPACE
===================================================================
--- pkg/CHNOSZ/NAMESPACE 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/NAMESPACE 2022-01-31 13:27:27 UTC (rev 684)
@@ -51,7 +51,7 @@
# added 20170301 or later
"GHS_Tr", "calculateDensity", "calculateGibbsOfWater",
"calculateEpsilon", "calculateQ", "water.DEW", "berman",
- "maxdiff", "expect_maxdiff", "bgamma",
+ "bgamma",
# added 20171121 or later
"dumpdata", "thermo.axis", "solubility", "NaCl",
# added 20190213 or later
Modified: pkg/CHNOSZ/R/diagram.R
===================================================================
--- pkg/CHNOSZ/R/diagram.R 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/R/diagram.R 2022-01-31 13:27:27 UTC (rev 684)
@@ -73,7 +73,7 @@
} else if(type %in% rownames(eout$basis)) {
# to calculate the loga of basis species at equilibrium
if(!missing(groups)) stop("can't plot equilibrium activities of basis species for grouped species")
- if(isTRUE(alpha) | is.character(alpha)) stop("equilibrium activities of basis species not available with alpha=TRUE")
+ if(isTRUE(alpha) | is.character(alpha)) stop("equilibrium activities of basis species not available with alpha = TRUE")
plot.loga.basis <- TRUE
} else if(type=="loga.equil" & !"loga.equil" %in% names(eout)) stop("'eout' is not the output from equil()")
else if(!type %in% c("loga.equil", "loga.balance")) stop(type, " is not a valid diagram type")
Modified: pkg/CHNOSZ/R/util.character.R
===================================================================
--- pkg/CHNOSZ/R/util.character.R 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/R/util.character.R 2022-01-31 13:27:27 UTC (rev 684)
@@ -58,23 +58,7 @@
if(length(x) == 0) FALSE else
if(length(x) > 1) as.logical(sapply(x, can.be.numeric)) else {
if(is.numeric(x)) TRUE else
- if(!is.na(as.numeric.nowarn(x))) TRUE else
+ if(!is.na(suppressWarnings(as.numeric(x)))) TRUE else
if(x %in% c('.','+','-')) TRUE else FALSE
}
}
-
-# something like R's as.numeric(), but without the "NAs introduced by coercion" warnings
-# (needed because testthat somehow detects the warnings suppressed by suppressWarnings) 20170427
-as.numeric.nowarn <- function(x) {
- if(length(x) == 0) numeric() else
- if(length(x) > 1) sapply(x, as.numeric.nowarn) else
- # http://stackoverflow.com/questions/12643009/regular-expression-for-floating-point-numbers
- if(grepl("^[+-]?([0-9]*[.])?[0-9]+$", x)) as.numeric(x) else NA_real_
-}
-
-# convert to integer without NA coercion warnings
-as.integer.nowarn <- function(x) {
- if(length(x) == 0) integer() else
- if(length(x) > 1) sapply(x, as.integer.nowarn) else
- if(grepl("[^0-9]", x)) NA_integer_ else as.integer(x)
-}
Modified: pkg/CHNOSZ/R/util.formula.R
===================================================================
--- pkg/CHNOSZ/R/util.formula.R 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/R/util.formula.R 2022-01-31 13:27:27 UTC (rev 684)
@@ -169,7 +169,7 @@
# return the values in the argument, or chemical formula(s)
# for values that are species indices
# for numeric values, get the formulas from those rownumbers of thermo()$OBIGT
- i <- as.integer.nowarn(formula)
+ i <- suppressWarnings(as.integer(formula))
# we can't have more than the number of rows in thermo()$OBIGT
thermo <- get("thermo", CHNOSZ)
iover <- i > nrow(thermo$OBIGT)
Deleted: pkg/CHNOSZ/R/util.test.R
===================================================================
--- pkg/CHNOSZ/R/util.test.R 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/R/util.test.R 2022-01-31 13:27:27 UTC (rev 684)
@@ -1,25 +0,0 @@
-# CHNOSZ/util.test.R
-# functions for writing tests 20171005
-
-# the maximum (absolute) pairwise difference between x and y
-maxdiff <- function(x, y) max(abs(y - x))
-
-# modelled after the "expect_" functions in testthat ...
-# is the maximum of the pairwise differences between two objects less than some value?
-expect_maxdiff <- function(object, expected, maxdiff = 0) {
- if(!"testthat" %in% row.names(installed.packages())) {
- stop("please install the 'testthat' package to use this function")
- } else {
- # get the names for the object and expected values
- # we use a double substitute here because
- # as.character(substitute(x$a$b)) == c("$", "x$a", "b")
- # as.character(substitute(substitute(x$a$b)))[2] == "x$a$b" --> better for printing
- lab_act <- as.character(substitute(substitute(object)))[2]
- lab_exp <- as.character(substitute(substitute(expected)))[2]
- truemd <- maxdiff(object, expected)
- failmsg <- sprintf("maxdiff(%s, %s) not less than %s", lab_act, lab_exp, maxdiff)
- truemsg <- sprintf("actual value: %s", truemd)
- testthat::expect(maxdiff(object, expected) <= maxdiff, sprintf("%s\n%s", failmsg, truemsg))
- invisible(object)
- }
-}
Modified: pkg/CHNOSZ/inst/CHECKLIST
===================================================================
--- pkg/CHNOSZ/inst/CHECKLIST 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/inst/CHECKLIST 2022-01-31 13:27:27 UTC (rev 684)
@@ -5,10 +5,8 @@
- Run examples() and demos() and inspect their output (especially plots)
-- Install the package with --install-tests and run
- library(CHNOSZ); testthat::test_package("CHNOSZ", reporter = "progress")
- - runs *all* tests (long ones are skipped under R CMD check --as-cran)
- - finds warnings from tests that aren't shown in R CMD check
+- Run the package test with:
+ suppressMessages(tinytest::test_package("CHNOSZ", at_home = TRUE))
- Check that uniprot.aa() works with current UniProt web pages
Modified: pkg/CHNOSZ/inst/NEWS.Rd
===================================================================
--- pkg/CHNOSZ/inst/NEWS.Rd 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/inst/NEWS.Rd 2022-01-31 13:27:27 UTC (rev 684)
@@ -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-9 (2021-09-27)}{
+\section{Changes in CHNOSZ version 1.4.1-11 (2022-01-30)}{
\subsection{OBIGT DATABASE}{
\itemize{
@@ -17,10 +17,11 @@
\item Move Ar, Xe, CH\s{4}, and CO\s{2} from organic_aq.csv to inorganic_aq.csv.
- \item Remove OldAA.csv (superseded thermodynamic parameters for amino
- acids). This file is now available in the JMDplots package.
+ \item Remove \samp{OldAA.csv} (superseded thermodynamic parameters for
+ amino acids). This file is now available in the
+ \href{https://github.com/jedick/JMDplots}{JMDplots} package.
- \item Fix formula of CaCl2 in DEW.csv. Thanks to Grayson Boyer.
+ \item Fix formula of CaCl\s{2} in \samp{DEW.csv}. Thanks to Grayson Boyer.
\item Add H\s{2}WO\s{4}(aq) and Cp coefficients of scheelite (CaWO\s{4})
from \href{https://doi.org/10.1016/j.chemgeo.2021.120488}{Liu et al.
@@ -29,6 +30,16 @@
}
}
+ \subsection{TESTING}{
+ \itemize{
+
+ \item Tests are now run using the \CRANpkg{tinytest} package.
+
+ \item \code{maxdiff()} and \code{expect_maxdiff()} have been removed.
+
+ }
+ }
+
\subsection{OTHER CHANGES}{
\itemize{
Modified: pkg/CHNOSZ/inst/tinytest/test-DEW.R
===================================================================
--- pkg/CHNOSZ/tests/testthat/test-DEW.R 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/inst/tinytest/test-DEW.R 2022-01-31 13:27:27 UTC (rev 684)
@@ -1,161 +1,148 @@
-context("DEW")
+# Load default settings for CHNOSZ
+reset()
-# since this is alphabetically the first test,
-# we need to load the 'thermo' object (for running tests in R CMD check)
-suppressMessages(reset())
-# get properties of water from DEW implementation in CHNOSZ
+# Get properties of water from DEW implementation in CHNOSZ
water("DEW")
-# use DEW species parameters in OBIGT database
+# Use DEW species parameters in OBIGT database
add.OBIGT("DEW")
-test_that("density of water is calculated correctly", {
- pressure <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
- temperature <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
- # density from R translation of DEW macro functions
- RDensity <- calculateDensity(pressure, temperature)
- # density from DEW spreadsheet
- DEWDensity <- c(1.108200, 0.597623, 1.196591, 0.798331, 1.321050, 1.000735, 1.578116, 1.287663)
- expect_equal(RDensity, DEWDensity, tolerance=1e-6)
-})
+info <- "Density of water is calculated correctly"
+pressure <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
+temperature <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
+# Density from R translation of DEW macro functions
+RDensity <- calculateDensity(pressure, temperature)
+# Density from DEW spreadsheet
+DEWDensity <- c(1.108200, 0.597623, 1.196591, 0.798331, 1.321050, 1.000735, 1.578116, 1.287663)
+expect_equal(RDensity, DEWDensity, tolerance = 1e-6, info = info)
-test_that("Gibbs energy of water is calculated correctly", {
- pressure <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
- temperature <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
- # Gibbs energies from R translation of DEW macro functions
- RGibbs <- calculateGibbsOfWater(pressure, temperature)
- # Gibbs energies from DEW spreadsheet
- DEWGibbs <- c(-56019.85419280258, -84262.028821198, -54155.004480575895, -81210.38766217149,
- -50735.122222685815, -76433.07602205424, -41823.26077175943, -65187.48113532527)
- expect_equal(RGibbs, DEWGibbs)
-})
+info <- "Gibbs energy of water is calculated correctly"
+pressure <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
+temperature <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
+# Gibbs energies from R translation of DEW macro functions
+RGibbs <- calculateGibbsOfWater(pressure, temperature)
+# Gibbs energies from DEW spreadsheet
+DEWGibbs <- c(-56019.85419280258, -84262.028821198, -54155.004480575895, -81210.38766217149,
+ -50735.122222685815, -76433.07602205424, -41823.26077175943, -65187.48113532527)
+expect_equal(RGibbs, DEWGibbs, info = info)
-test_that("dielectric constant of water is calculated correctly", {
- pressure <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
- temperature <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
- # epsilon from R translation of DEW macro functions
- Repsilon <- calculateEpsilon(calculateDensity(pressure, temperature), temperature)
- # epsilon from DEW spreadsheet
- DEWepsilon <- c(65.63571, 6.10465, 72.40050, 8.97800, 82.16244, 12.13131, 103.12897, 16.97266)
- expect_equal(Repsilon, DEWepsilon, tolerance=1e-7)
-})
+info <- "Dielectric constant of water is calculated correctly"
+pressure <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
+temperature <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
+# epsilon from R translation of DEW macro functions
+Repsilon <- calculateEpsilon(calculateDensity(pressure, temperature), temperature)
+# epsilon from DEW spreadsheet
+DEWepsilon <- c(65.63571, 6.10465, 72.40050, 8.97800, 82.16244, 12.13131, 103.12897, 16.97266)
+expect_equal(Repsilon, DEWepsilon, tolerance = 1e-7, info = info)
-test_that("Born coefficient Q is calculated correctly", {
- pressure <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
- temperature <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
- # Q from R translation of DEW macro functions
- RQ <- calculateQ(calculateDensity(pressure, temperature), temperature)
- # Q from DEW spreadsheet
- DEWQ <- c(0.32319817, 14.50286092, 0.19453478, 3.12650897,
- 0.10918151, 0.87729257, 0.05068788, 0.20640645) / 1e6
- expect_equal(RQ, DEWQ)
-})
+info <- "Born coefficient Q is calculated correctly"
+pressure <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
+temperature <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
+# Q from R translation of DEW macro functions
+RQ <- calculateQ(calculateDensity(pressure, temperature), temperature)
+# Q from DEW spreadsheet
+DEWQ <- c(0.32319817, 14.50286092, 0.19453478, 3.12650897,
+ 0.10918151, 0.87729257, 0.05068788, 0.20640645) / 1e6
+expect_equal(RQ, DEWQ, info = info)
-test_that("g function is calculated correctly", {
- pressure <- c(1000, 1000, 5000, 5000, 10000)
- temperature <- c(100, 1000, 100, 1000, 100)
- # note that values returned for alpha, daldT, beta are NA
- w <- water(c("rho", "alpha", "daldT", "beta", "Psat"), T=convert(temperature, "K"), P=pressure)
- # g from CHNOSZ functions
- Rg <- gfun(w$rho/1000, temperature, pressure, w$alpha, w$daldT, w$beta)$g
- # g from R translation of DEW macro functions (not used in CHNOSZ)
- DEWg <- calculateG(pressure, temperature, w$rho/1000)
- expect_equal(Rg, DEWg)
-})
+info <- "g function is calculated correctly"
+pressure <- c(1000, 1000, 5000, 5000, 10000)
+temperature <- c(100, 1000, 100, 1000, 100)
+# note that values returned for alpha, daldT, beta are NA
+w <- water(c("rho", "alpha", "daldT", "beta", "Psat"), T=convert(temperature, "K"), P=pressure)
+# g from CHNOSZ functions
+Rg <- CHNOSZ:::gfun(w$rho/1000, temperature, pressure, w$alpha, w$daldT, w$beta)$g
+# g from R translation of DEW macro functions (not used in CHNOSZ)
+DEWg <- CHNOSZ:::calculateG(pressure, temperature, w$rho/1000)
+expect_equal(Rg, DEWg, info = info)
-test_that("Gibbs energies of species are calculated correctly", {
- P <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
- T <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
- RG_HCl <- subcrt("HCl", P=P, T=T)$out[[1]]$G
- DEWG_HCl <- c(-28784.99, -58496.85, -26520.94, -55276.92, -21928.89, -50337.19, -8014.34, -36746.87)
- expect_equal(RG_HCl, DEWG_HCl, tolerance = 1e-5)
- RG_Cl <- subcrt("Cl-", P=P, T=T)$out[[1]]$G
- DEWG_Cl <- c(-30054.59, -22839.35, -27910.68, -28094.07, -23568.45, -27959.67, -10443.07, -18744.93)
- expect_equal(RG_Cl, DEWG_Cl, tolerance = 1e-7)
-})
+info <- "Gibbs energies of species are calculated correctly"
+P <- c(5000, 5000, 10000, 10000, 20000, 20000, 50000, 50000)
+T <- c(100, 1000, 100, 1000, 100, 1000, 100, 1000)
+RG_HCl <- subcrt("HCl", P=P, T=T)$out[[1]]$G
+DEWG_HCl <- c(-28784.99, -58496.85, -26520.94, -55276.92, -21928.89, -50337.19, -8014.34, -36746.87)
+expect_equal(RG_HCl, DEWG_HCl, tolerance = 1e-5, info = info)
+RG_Cl <- subcrt("Cl-", P=P, T=T)$out[[1]]$G
+DEWG_Cl <- c(-30054.59, -22839.35, -27910.68, -28094.07, -23568.45, -27959.67, -10443.07, -18744.93)
+expect_equal(RG_Cl, DEWG_Cl, tolerance = 1e-7, info = info)
-test_that("Delta G, logK, and Delta V of reactions are calculated correctly", {
- # These are reactions corresponding to Fig. 1b of Sverjensky et al., 2014 (Nat. Geosci.).
- # The properties are calculated using parameters from the DEW spreadsheet,
- # which are not necessarily identical those that were used for the paper.
- T <- 600
- P <- c(5000, 50000)
- R1 <- subcrt(c("H2O", "CO2", "HCO3-", "H+"), c(-1, -1, 1, 1), T=T, P=P)$out
- R2 <- subcrt(c("HCO3-", "CO3-2", "H+"), c(-1, 1, 1), T=T, P=P)$out
- R3 <- subcrt(c("acetic acid", "acetate", "H+"), c(-1, 1, 1), T=T, P=P)$out
- R4 <- subcrt(c("H2O", "CO2", "acetic acid", "oxygen"), c(-2, -2, 1, 2), T=T, P=P)$out
- R5 <- subcrt(c("oxygen", "CH4", "acetic acid", "H2O"), c(-2, -2, 1, 2), T=T, P=P)$out
+info <- "Delta G, logK, and Delta V of reactions are calculated correctly"
+# These are reactions corresponding to Fig. 1b of Sverjensky et al., 2014 (Nat. Geosci.).
+# The properties are calculated using parameters from the DEW spreadsheet,
+# which are not necessarily identical those that were used for the paper.
+T <- 600
+P <- c(5000, 50000)
+R1 <- subcrt(c("H2O", "CO2", "HCO3-", "H+"), c(-1, -1, 1, 1), T=T, P=P)$out
+R2 <- subcrt(c("HCO3-", "CO3-2", "H+"), c(-1, 1, 1), T=T, P=P)$out
+R3 <- subcrt(c("acetic acid", "acetate", "H+"), c(-1, 1, 1), T=T, P=P)$out
+R4 <- subcrt(c("H2O", "CO2", "acetic acid", "oxygen"), c(-2, -2, 1, 2), T=T, P=P)$out
+R5 <- subcrt(c("oxygen", "CH4", "acetic acid", "H2O"), c(-2, -2, 1, 2), T=T, P=P)$out
- # Delta G values calculated using the DEW spreadsheet (May 2017 version)
- DEW_DG <- c(38267.404507442814, 14893.170655988564, # R1
- 41407.05995898576, 21347.599525026497, # R2
- 28109.598104640143, 16112.928184675075, # R3
- 186960.22705581048, 133640.9631638353, # R4
- -141552.60059404257, -134279.54670605875) # R5
- # the aqueous-only reactions
- expect_equal(c(R1$G, R2$G, R3$G), DEW_DG[1:6], tolerance=1e-7)
- # note that there is a small error for oxygen in the DEW spreadsheet (wrong c parameter),
- # but even so, these tests have a lower tolerance because of the larger magnitude of the difference
- expect_equal(c(R4$G, R5$G), DEW_DG[7:10], tolerance=1e-9)
+# Delta G values calculated using the DEW spreadsheet (May 2017 version)
+DEW_DG <- c(38267.404507442814, 14893.170655988564, # R1
+ 41407.05995898576, 21347.599525026497, # R2
+ 28109.598104640143, 16112.928184675075, # R3
+ 186960.22705581048, 133640.9631638353, # R4
+ -141552.60059404257, -134279.54670605875) # R5
+# the aqueous-only reactions
+expect_equal(c(R1$G, R2$G, R3$G), DEW_DG[1:6], tolerance = 1e-7, info = info)
+# note that there is a small error for oxygen in the DEW spreadsheet (wrong c parameter),
+# but even so, these tests have a lower tolerance because of the larger magnitude of the difference
+expect_equal(c(R4$G, R5$G), DEW_DG[7:10], tolerance = 1e-9, info = info)
- # logK values calculated using the DEW spreadsheet
- DEW_logK <- c(-9.58455651110442, -3.7301833667366027,
- -10.370923015131565, -5.346776889042665,
- -7.040405143911882, -4.035687100632826,
- -46.826558649850625, -33.47207316851283,
- 35.45364304557209, 33.632014510923746)
- expect_equal(c(R1$logK, R2$logK, R3$logK, R4$logK, R5$logK), DEW_logK, tolerance=1e-3)
+# logK values calculated using the DEW spreadsheet
+DEW_logK <- c(-9.58455651110442, -3.7301833667366027,
+ -10.370923015131565, -5.346776889042665,
+ -7.040405143911882, -4.035687100632826,
+ -46.826558649850625, -33.47207316851283,
+ 35.45364304557209, 33.632014510923746)
+expect_equal(c(R1$logK, R2$logK, R3$logK, R4$logK, R5$logK), DEW_logK, tolerance = 1e-3, info = info)
- # Delta V values calculated using the DEW spreadsheet
- DEW_DV <- c(-45.26925983499276, -14.640599169742725,
- -47.95180846799733, -9.469432706749927,
- -27.042087540311922, -6.836267057722694,
- -18.1550937649195, 5.513800665649967,
- -37.37077435045512, -45.08570662275889)
- # for the aqueous species we're getting very close results
- # (at P=5000 bar this depends on calculating drhodP -> beta -> dgdP -> dwdP -> V correctly, which is not tested above)
- expect_equal(c(R1$V, R2$V, R3$V), DEW_DV[1:6], tolerance=1e-15)
- # TODO: why does DEW spreadsheet use V (O2,g) == 24.465?
- #expect_equal(c(R4$V, R5$V), DEW_DV[7:10])
-})
+# Delta V values calculated using the DEW spreadsheet
+DEW_DV <- c(-45.26925983499276, -14.640599169742725,
+ -47.95180846799733, -9.469432706749927,
+ -27.042087540311922, -6.836267057722694,
+ -18.1550937649195, 5.513800665649967,
+ -37.37077435045512, -45.08570662275889)
+# for the aqueous species we're getting very close results
+# (at P=5000 bar this depends on calculating drhodP -> beta -> dgdP -> dwdP -> V correctly, which is not tested above)
+expect_equal(c(R1$V, R2$V, R3$V), DEW_DV[1:6], tolerance = 1e-15, info = info)
+# TODO: why does DEW spreadsheet use V (O2,g) == 24.465?
+#expect_equal(c(R4$V, R5$V), DEW_DV[7:10])
-test_that("calculated logK values are consistent with Extended Deep Earth Water paper", {
- # Reference logK values are from Appendix D of Huang and Sverjensky, 2019 (doi:10.1016/j.gca.2019.03.027)
- # Select T and P for comparisons
- T <- c(300, 1100)
- P <- c(10000, 60000)
- # Adjust calculated logKs for different values of the
- # gas constant used in the DEW spreadsheet and CHNOSZ
- RoverR <- 1.9858775 / 1.9872
- # Calculate logK for each reaction
- logK1 <- subcrt(c("H2O", "CO2", "H2CO3"), c(-1, -1, 1), T = T, P = P)$out$logK / RoverR
- logK2 <- subcrt(c("AlO2(SiO2)-", "H+", "Al+3", "H2O", "SiO2"), c(-1, -4, 1, 2, 1), T = T, P = P)$out$logK / RoverR
- logK3 <- subcrt(c("Ca(HCO3)+", "Ca+2", "HCO3-"), c(-1, 1, 1), T = T, P = P)$out$logK / RoverR
- logK4 <- subcrt(c("Ca(HCOO)+", "Ca+2", "HCOO-"), c(-1, 1, 1), T = T, P = P)$out$logK / RoverR
- logK5 <- subcrt(c("Ca(HSiO3)+", "H+", "Ca+2", "H2O", "SiO2"), c(-1, -1, 1, 1, 1), T = T, P = P)$out$logK / RoverR
- logK6 <- subcrt(c("Fe(HCOO)+", "Fe+2", "HCOO-"), c(-1, 1, 1), T = T, P = P)$out$logK / RoverR
- logK7 <- subcrt(c("Fe(HSiO3)+", "H+", "Fe+2", "H2O", "SiO2"), c(-1, -1, 1, 1, 1), T = T, P = P)$out$logK / RoverR
- logK8 <- subcrt(c("H2O", "SiO2", "HSiO3-", "H+"), c(-1, -1, 1, 1), T = T, P = P)$out$logK / RoverR
- logK9 <- subcrt(c("MgO", "H+", "Mg+2", "H2O"), c(-1, -2, 1, 1), T = T, P = P)$out$logK / RoverR
- logK10 <- subcrt(c("Mg(SiO2)(HCO3)+", "Mg+2", "SiO2", "HCO3-"), c(-1, 1, 1, 1), T = T, P = P)$out$logK / RoverR
- logK11 <- subcrt(c("NaHCO3", "Na+", "HCO3-"), c(-1, 1, 1), T = T, P = P)$out$logK / RoverR
- logK12 <- subcrt(c("Si3O6", "SiO2"), c(-1, 3), T = T, P = P)$out$logK / RoverR
- # Make the comparisons:
- # - Reference logK values are from Appendix D of Huang and Sverjensky, 2019 (doi:10.1016/j.gca.2019.03.027)
- # - Tolerance (tol) is added if there are differences between the reference and calculated values after rounding
- # (scale = 1 is used to compare absolute differences)
- expect_equal(round(logK1, 2), c(-0.76, 1.34))
- expect_equal(round(logK2, 4), c( 6.7755, 0.9413), tol = 0.003, scale = 1)
- expect_equal(round(logK3, 4), c(-2.3759, -5.6840))
- expect_equal(round(logK4, 4), c(-2.2837, -5.5822), tol = 0.0001, scale = 1)
- expect_equal(round(logK5, 4), c( 4.0349, 0.1797), tol = 0.002, scale = 1)
- expect_equal(round(logK6, 4), c(-7.5354, -8.0238))
- expect_equal(round(logK7, 4), c(-0.6883, -1.6363), tol = 0.002, scale = 1)
- expect_equal(round(logK8, 4), c(-7.0651, -5.5067), tol = 0.002, scale = 1)
- expect_equal(round(logK9, 4), c( 8.2759, 4.3493), tol = 0.002, scale = 1)
- expect_equal(round(logK10, 4), c(-6.8106, -8.5888), tol = 0.0001, scale = 1)
- expect_equal(round(logK11, 4), c(-0.2447, -2.9235))
- expect_equal(round(logK12, 4), c( 3.3283, 0.4527))
-})
-
-# Restore the default water model and OBIGT database
-suppressMessages(reset())
+info <- "Calculated logK values are consistent with Extended Deep Earth Water paper"
+# Reference logK values are from Appendix D of Huang and Sverjensky, 2019 (doi:10.1016/j.gca.2019.03.027)
+# Select T and P for comparisons
+T <- c(300, 1100)
+P <- c(10000, 60000)
+# Adjust calculated logKs for different values of the
+# gas constant used in the DEW spreadsheet and CHNOSZ
+RoverR <- 1.9858775 / 1.9872
+# Calculate logK for each reaction
+logK1 <- subcrt(c("H2O", "CO2", "H2CO3"), c(-1, -1, 1), T = T, P = P)$out$logK / RoverR
+logK2 <- subcrt(c("AlO2(SiO2)-", "H+", "Al+3", "H2O", "SiO2"), c(-1, -4, 1, 2, 1), T = T, P = P)$out$logK / RoverR
+logK3 <- subcrt(c("Ca(HCO3)+", "Ca+2", "HCO3-"), c(-1, 1, 1), T = T, P = P)$out$logK / RoverR
+logK4 <- subcrt(c("Ca(HCOO)+", "Ca+2", "HCOO-"), c(-1, 1, 1), T = T, P = P)$out$logK / RoverR
+logK5 <- subcrt(c("Ca(HSiO3)+", "H+", "Ca+2", "H2O", "SiO2"), c(-1, -1, 1, 1, 1), T = T, P = P)$out$logK / RoverR
+logK6 <- subcrt(c("Fe(HCOO)+", "Fe+2", "HCOO-"), c(-1, 1, 1), T = T, P = P)$out$logK / RoverR
+logK7 <- subcrt(c("Fe(HSiO3)+", "H+", "Fe+2", "H2O", "SiO2"), c(-1, -1, 1, 1, 1), T = T, P = P)$out$logK / RoverR
+logK8 <- subcrt(c("H2O", "SiO2", "HSiO3-", "H+"), c(-1, -1, 1, 1), T = T, P = P)$out$logK / RoverR
+logK9 <- subcrt(c("MgO", "H+", "Mg+2", "H2O"), c(-1, -2, 1, 1), T = T, P = P)$out$logK / RoverR
+logK10 <- subcrt(c("Mg(SiO2)(HCO3)+", "Mg+2", "SiO2", "HCO3-"), c(-1, 1, 1, 1), T = T, P = P)$out$logK / RoverR
+logK11 <- subcrt(c("NaHCO3", "Na+", "HCO3-"), c(-1, 1, 1), T = T, P = P)$out$logK / RoverR
+logK12 <- subcrt(c("Si3O6", "SiO2"), c(-1, 3), T = T, P = P)$out$logK / RoverR
+# Make the comparisons:
+# - Reference logK values are from Appendix D of Huang and Sverjensky, 2019 (doi:10.1016/j.gca.2019.03.027)
+# - Tolerance (tol) is added if there are differences between the reference and calculated values after rounding
+# (scale = 1 is used to compare absolute differences)
+expect_equal(round(logK1, 2), c(-0.76, 1.34), info = info)
+expect_equal(round(logK2, 4), c( 6.7755, 0.9413), tol = 0.003, scale = 1, info = info)
+expect_equal(round(logK3, 4), c(-2.3759, -5.6840), info = info)
+expect_equal(round(logK4, 4), c(-2.2837, -5.5822), tol = 0.0001, scale = 1, info = info)
+expect_equal(round(logK5, 4), c( 4.0349, 0.1797), tol = 0.002, scale = 1, info = info)
+expect_equal(round(logK6, 4), c(-7.5354, -8.0238), info = info)
+expect_equal(round(logK7, 4), c(-0.6883, -1.6363), tol = 0.002, scale = 1, info = info)
+expect_equal(round(logK8, 4), c(-7.0651, -5.5067), tol = 0.002, scale = 1, info = info)
+expect_equal(round(logK9, 4), c( 8.2759, 4.3493), tol = 0.002, scale = 1, info = info)
+expect_equal(round(logK10, 4), c(-6.8106, -8.5888), tol = 0.0001, scale = 1, info = info)
+expect_equal(round(logK11, 4), c(-0.2447, -2.9235), info = info)
+expect_equal(round(logK12, 4), c( 3.3283, 0.4527), info = info)
Modified: pkg/CHNOSZ/inst/tinytest/test-EOSregress.R
===================================================================
--- pkg/CHNOSZ/tests/testthat/test-EOSregress.R 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/inst/tinytest/test-EOSregress.R 2022-01-31 13:27:27 UTC (rev 684)
@@ -1,44 +1,43 @@
-context("EOSregress")
+# Load default settings for CHNOSZ
+reset()
-test_that("EOSvar stops with unknown variables", {
- expect_error(EOSvar("TX", T=25, P=1), "can't find a variable named TX")
- # why can't the test find these?
- #TX <- 2
- #expect_error(EOSvar("TX", T=25, P=1), "an object named TX is not a function")
- #TX <- function(T) 2
- #expect_error(EOSvar("TX", T=25, P=1), "the arguments of TX\\(\\) are not T, P")
-})
+info <- "EOSvar stops with unknown variables"
+expect_error(EOSvar("TX", T=25, P=1), "can't find a variable named TX", info = info)
+# why can't the test find these?
+#TX <- 2
+#expect_error(EOSvar("TX", T=25, P=1), "an object named TX is not a function")
+#TX <- function(T) 2
+#expect_error(EOSvar("TX", T=25, P=1), "the arguments of TX\\(\\) are not T, P")
-test_that("regressions return known HKF parameters (neutral species)", {
- # regress computed values of heat capacity and volume of CH4(aq)
- # calculated from HKF parameters on a T-P grid
- T <- convert(seq(0, 350, 25), "K")
- P <- seq(200, 1000, 100)
- # convert=FALSE means that temperature has units of K
- CH4.prop <- subcrt("CH4", T=T, P=P, grid="T", convert=FALSE)$out[[1]]
- # terms in the HKF equations for Cp
- Cp.var <- c("invTTheta2", "TXBorn")
- # get coefficients in Cp regression
- Cp.lm <- EOSregress(CH4.prop[, c("T", "P", "Cp")], Cp.var)
- Cp.coeff <- Cp.lm$coefficients
- # terms in the HKF equations for V
- V.var <- c("invPPsi", "invTTheta", "invPPsiTTheta", "QBorn")
- # get coefficients in V regression
- V.lm <- EOSregress(CH4.prop[, c("T", "P", "V")], V.var)
- # use same units as HKF: convert from cm3.bar to calories (divide by 41.84)
- V.coeff <- convert(V.lm$coefficients, "calories")
- ## the tests: did we get the HKF parameters that are in the database?
- CH4.par <- info(info("CH4"))
- # c1 and c2
- expect_equivalent(Cp.coeff[1], CH4.par$c1)
- expect_equivalent(Cp.coeff[2], CH4.par$c2)
- # omega (from Cp)
- expect_equivalent(Cp.coeff[3], CH4.par$omega)
- # a1, a2, a3 and a4
- expect_equivalent(V.coeff[1], CH4.par$a1)
- expect_equivalent(V.coeff[2], CH4.par$a2)
- expect_equivalent(V.coeff[3], CH4.par$a3)
- expect_equivalent(V.coeff[4], CH4.par$a4)
- # omega (from V)
- expect_equivalent(V.coeff[5], CH4.par$omega)
-})
+info <- "Regressions return known HKF parameters (neutral species)"
+# regress computed values of heat capacity and volume of CH4(aq)
+# calculated from HKF parameters on a T-P grid
+T <- convert(seq(0, 350, 25), "K")
+P <- seq(200, 1000, 100)
+# convert=FALSE means that temperature has units of K
+CH4.prop <- subcrt("CH4", T=T, P=P, grid="T", convert=FALSE)$out[[1]]
+# terms in the HKF equations for Cp
+Cp.var <- c("invTTheta2", "TXBorn")
+# get coefficients in Cp regression
+Cp.lm <- EOSregress(CH4.prop[, c("T", "P", "Cp")], Cp.var)
+Cp.coeff <- Cp.lm$coefficients
+# terms in the HKF equations for V
+V.var <- c("invPPsi", "invTTheta", "invPPsiTTheta", "QBorn")
+# get coefficients in V regression
+V.lm <- EOSregress(CH4.prop[, c("T", "P", "V")], V.var)
+# use same units as HKF: convert from cm3.bar to calories (divide by 41.84)
+V.coeff <- convert(V.lm$coefficients, "calories")
+## the tests: did we get the HKF parameters that are in the database?
+CH4.par <- info(info("CH4"))
+# c1 and c2
+expect_equivalent(Cp.coeff[1], CH4.par$c1, info = info)
+expect_equivalent(Cp.coeff[2], CH4.par$c2, info = info)
+# omega (from Cp)
+expect_equivalent(Cp.coeff[3], CH4.par$omega, info = info)
+# a1, a2, a3 and a4
+expect_equivalent(V.coeff[1], CH4.par$a1, info = info)
+expect_equivalent(V.coeff[2], CH4.par$a2, info = info)
+expect_equivalent(V.coeff[3], CH4.par$a3, info = info)
+expect_equivalent(V.coeff[4], CH4.par$a4, info = info)
+# omega (from V)
+expect_equivalent(V.coeff[5], CH4.par$omega, info = info)
Modified: pkg/CHNOSZ/inst/tinytest/test-IAPWS95.R
===================================================================
--- pkg/CHNOSZ/tests/testthat/test-IAPWS95.R 2022-01-29 11:25:59 UTC (rev 683)
+++ pkg/CHNOSZ/inst/tinytest/test-IAPWS95.R 2022-01-31 13:27:27 UTC (rev 684)
@@ -1,88 +1,86 @@
-context("IAPWS95")
+# Load default settings for CHNOSZ
+reset()
-test_that("calculations of Helmholtz free energy and its derivatives are consistent with reference cases", {
- ## reference values of these terms are listed Table 6.6 of Wagner and Pruss, 2002
- p <- c("phi", "phi.delta", "phi.delta.delta", "phi.tau", "phi.tau.tau", "phi.delta.tau")
- # reference values for case 1: T=500 K, rho=838.025 kg m-3
- idealgas.ref.1 <- c(0.204797734e1, 0.384236747, -0.147637878, 0.904611106e1, -0.193249185e1, 0)
- residual.ref.1 <- c(-0.342693206e1, -0.364366650, 0.856063701, -0.581403435e1, -0.223440737e1, -0.112176915e1)
- # reference values for case 1: T=647 K, rho=358 kg m-3
- idealgas.ref.2 <- c(-0.156319605e1, 0.899441341, -0.808994726, 0.980343918e1, -0.343316334e1, 0)
- residual.ref.2 <- c(-0.121202657e1, -0.714012024, 0.475730696, -0.321722501e1, -0.996029507e1, -0.133214720e1)
- ## set up the problem
- # critical point constants
- T.critical <- 647.096 # K
- rho.critical <- 322 # kg m-3
- # T and rho for cases 1 and 2
- T <- c(500, 647)
- rho <- c(838.025, 358)
- # delta and tau for cases 1 and 2
- delta <- rho / rho.critical
- tau <- T.critical / T
- # calculated ideal gas and residual parts for case 1
- idealgas.calc.1 <- sapply(p, IAPWS95.idealgas, delta[1], tau[1])
- residual.calc.1 <- sapply(p, IAPWS95.residual, delta[1], tau[1])
- # calculated ideal gas and residual parts for case 2
- idealgas.calc.2 <- sapply(p, IAPWS95.idealgas, delta[2], tau[2])
- residual.calc.2 <- sapply(p, IAPWS95.residual, delta[2], tau[2])
- ## perform the tests
- # we almost get away without increasing the tolerance in any test ...
- expect_equal(idealgas.calc.1, idealgas.ref.1, check.attributes=FALSE)
- expect_equal(residual.calc.1, residual.ref.1, check.attributes=FALSE)
- expect_equal(idealgas.calc.2, idealgas.ref.2, check.attributes=FALSE)
- # ... however an offset is apparent in the value of the residual phi.delta.delta for case 2
- expect_equal(residual.calc.2, residual.ref.2, check.attributes=FALSE, tolerance=1e-5)
-})
+info <- "Calculations of Helmholtz free energy and its derivatives are consistent with reference cases"
+## reference values of these terms are listed Table 6.6 of Wagner and Pruss, 2002
+p <- c("phi", "phi.delta", "phi.delta.delta", "phi.tau", "phi.tau.tau", "phi.delta.tau")
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/chnosz -r 684
More information about the CHNOSZ-commits
mailing list