[CHNOSZ-commits] r639 - in pkg/CHNOSZ: . vignettes

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Wed Feb 24 10:24:38 CET 2021


Author: jedick
Date: 2021-02-24 10:24:38 +0100 (Wed, 24 Feb 2021)
New Revision: 639

Modified:
   pkg/CHNOSZ/DESCRIPTION
   pkg/CHNOSZ/vignettes/multi-metal.Rmd
Log:
multi-metal.Rmd: Add Cl to system in Mosaic Stacking 2


Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION	2021-02-23 09:17:50 UTC (rev 638)
+++ pkg/CHNOSZ/DESCRIPTION	2021-02-24 09:24:38 UTC (rev 639)
@@ -1,6 +1,6 @@
-Date: 2021-02-23
+Date: 2021-02-24
 Package: CHNOSZ
-Version: 1.4.0-8
+Version: 1.4.0-9
 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/vignettes/multi-metal.Rmd
===================================================================
--- pkg/CHNOSZ/vignettes/multi-metal.Rmd	2021-02-23 09:17:50 UTC (rev 638)
+++ pkg/CHNOSZ/vignettes/multi-metal.Rmd	2021-02-24 09:24:38 UTC (rev 639)
@@ -65,6 +65,7 @@
 if (!nzchar(Sys.which("pngquant"))) pngquant <- NULL
 ## logK with a thin space 20200627
 logK <- "log <i>K</i>"
+
 ## Resolution settings
 # Change this to TRUE to make high-resolution plots
 # (default is FALSE to save time in CRAN checks)
@@ -498,47 +499,67 @@
 <div id="D-stack2" style="display: none">
 ```{r stack2, eval = FALSE}
 # Define system
-logaH2S <- -2
-T <- 200
 pH <- c(0, 14, res2)
 O2 <- c(-48, -33, res2)
-basis(c("Cu+", "pyrite", "H2S", "oxygen", "H2O", "H+"))
-basis("H2S", logaH2S)
+T <- 200
+logmS <- -2
+m_NaCl <- 0.1
+logm_aq <- -6 # for both Fe- and Cu-bearing aq species
+# Basis species
 S.aq <- c("H2S", "HS-", "HSO4-", "SO4-2")
+# Minerals
 Fe.cr <- c("pyrite", "pyrrhotite", "magnetite", "hematite")
 Fe.abbrv <- c("Py", "Po", "Mag", "Hem")
+FeCu.cr <- c("chalcopyrite", "bornite")
+Cu.cr <- c("copper", "cuprite", "tenorite", "chalcocite", "covellite")
+FeCu.abbrv <- c("Ccp", "Bn", "Cu", "Cpr", "Tnr", "Cct", "Cv")
+# Aqueous species
+iFe.aq <- retrieve("Fe", c("S", "O", "H", "Cl"), "aq")
+Fe.aq <- info(iFe.aq)$name
+iCu.aq <- retrieve("Cu", c("S", "O", "H", "Cl"), "aq")
+Cu.aq <- info(iCu.aq)$name
+# Expressions for making the legend
+TPexpr <- describe.property(c("T", "P"), c(T, "Psat"))
+Sexpr <- as.expression(bquote(sum(S) == .(10^logmS)*m))
+NaClexpr <- as.expression(bquote(NaCl == .(m_NaCl)*m))
+aqexpr <- as.expression(bquote("("*aq*")"[italic(i)] == 10^.(logm_aq)*m))
 
+# Setup basis species
+basis(c("Cu+", "pyrite", "H2S", "oxygen", "H2O", "H+", "Cl-"))
+basis("H2S", logmS)
+nacl <- NaCl(T = T, P = "Psat", m_tot = m_NaCl)
+basis("Cl-", log10(nacl$m_Cl))
+# Use values for FeSO4(aq) from Table S2 of St Clair et al. (2019)
+mod.OBIGT("FeSO4", G = -203129, H = -238307)
 # Fe-bearing minerals
 species(Fe.cr)
 # Add aqueous species 20210220
-iFe.aq <- retrieve("Fe", c("S", "O", "H"), "aq")
-Fe.aq <- info(iFe.aq)$name
-logaFe <- -6
-species(iFe.aq, logaFe, add = TRUE)
-mFe <- mosaic(S.aq, pH = pH, O2 = O2, T = T)
+species(iFe.aq, logm_aq, add = TRUE)
+mFe <- mosaic(S.aq, pH = pH, O2 = O2, T = T, IS = nacl$IS)
 # Start plot with just the fields for transparency effect
 dFe <- diagram(mFe$A.species, lwd = 0, names = FALSE)
 
 # Cu-bearing minerals
-FeCu.cr <- c("chalcopyrite", "bornite")
-Cu.cr <- c("copper", "cuprite", "tenorite", "chalcocite", "covellite")
-FeCu.abbrv <- c("Ccp", "Bn", "Cu", "Cpr", "Tnr", "Cct", "Cv")
 species(c(FeCu.cr, Cu.cr))
 # Add aqueous species 20210220
-iCu.aq <- retrieve("Cu", c("S", "O", "H"), "aq")
-Cu.aq <- info(iCu.aq)$name
-logaCu <- -6
-species(iCu.aq, logaCu, add = TRUE)
+species(iCu.aq, logm_aq, add = TRUE)
 # TODO: limitation in mosaic() when using both solid and aqueous basis species (i.e. c(Fe.cr, Fe.aq)):
 #   Only the activity of the first-defined basis species is used for all basis species.
-#   The first-defined basis species is pyrite (logact = 0), but we want logact = -4 for the aq Fe species.
+#   The first-defined basis species is pyrite (logact = 0), but logact < 0 for the aq Fe species.
 # Workaround: Adjust standard Gibbs energies of aq Fe species to virtually change their activities
-DG <- convert(-logaFe, "G", T = convert(T, "K"))
+DG <- convert(-logm_aq, "G", T = convert(T, "K"))
 G.orig <- info(iFe.aq)$G
 G.new <- G.orig + DG
 mod.OBIGT(iFe.aq, G = G.new)
-mFeCu <- mosaic(list(S.aq, c(Fe.cr, Fe.aq)), pH = pH, O2 = O2, T = T, stable = list(NULL, dFe$predominant))
 
+#mFeCu <- mosaic(list(S.aq, c(Fe.cr, Fe.aq)), pH = pH, O2 = O2, T = T, stable = list(NULL, dFe$predominant))
+# Use only predominant species as basis species (to speed up calculation) 20210224
+predom <- dFe$predominant
+ipredom <- sort(unique(as.numeric(predom)))
+for(i in seq_along(ipredom)) predom[dFe$predominant == ipredom[i]] <- i
+Fe.predom <- c(Fe.cr, Fe.aq)[ipredom]
+mFeCu <- mosaic(list(S.aq, Fe.predom), pH = pH, O2 = O2, T = T, stable = list(NULL, predom), IS = nacl$IS)
+
 # Adjust labels
 bold <- c(rep(TRUE, length(FeCu.abbrv)), rep(FALSE, length(Cu.aq)))
 names <- c(FeCu.abbrv, Cu.aq)
@@ -545,22 +566,24 @@
 srt <- dx <- dy <- rep(0, length(names))
 cex <- rep(1, length(names))
 dx[names == "Cu"] <- -1.5
-dx[names == "Bn"] <- 1.3
-dx[names == "CuHS"] <- 2
+dx[names == "Bn"] <- 1.4
+dx[names == "CuHS"] <- 1
 dx[names == "Cu+"] <- -0.5
 dy[names == "Cu"] <- 3
-dy[names == "Cct"] <- 1
+dx[names == "Cct"] <- -2
+dy[names == "Cct"] <- 4
 dy[names == "CuHS"] <- 1
-dy[names == "Bn"] <- -0.5
+dy[names == "Bn"] <- -0.9
+dx[names == "CuCl2-"] <- -1
+dy[names == "CuCl2-"] <- 2
 cex[names == "Bn"] <- 0.8
-srt[names == "Bn"] <- 83
+srt[names == "Bn"] <- 85
 diagram(mFeCu$A.species, add = TRUE, col = 2, col.names = 2, names = names, bold = bold, dx = dx, dy = dy, cex.names = cex, srt = srt)
 # Add second Cu label
-text(12.5, -47, "Cu", col = 2, font = 2)
-
+text(12.3, -47, "Cu", col = 2, font = 2)
 # Highlight Bn-Ccp reaction
 col.names <- col <- rep(NA, nrow(mFeCu$A.species$species))
-col[2] <- "#FF8C00"
+col[3] <- "#FF8C00"
 col.names[1] <- "#FF8C00"
 diagram(mFeCu$A.species, add = TRUE, col = col, lwd = 2, col.names = col.names, bold = TRUE, names = names, fill = NA)
 
@@ -570,14 +593,20 @@
 names <- c(Fe.abbrv, Fe.aq)
 srt <- dx <- dy <- rep(0, length(names))
 cex <- rep(1, length(names))
-dy[names == "Hem"] <- -1.5
-dy[names == "Mag"] <- -0.5
-dx[names == "Mag"] <- 0.2
+dy[names == "Hem"] <- 0.5
+dy[names == "Mag"] <- -0.8
+dx[names == "Hem"] <- -0.5
+dx[names == "Mag"] <- 0.25
+dx[names == "Fe+2"] <- 0.5
+dx[names == "FeO2-"] <- 1
+dy[names == "FeO2-"] <- -3
+dx[names == "HFeO2-"] <- -0.5
+dy[names == "HFeO2-"] <- 1
 srt[names == "Mag"] <- 83
+srt[names == "FeSO4"] <- 90
 diagram(mFe$A.species, add = TRUE, lwd = 2, names = names, bold = bold, dx = dx, dy = dy, cex.names = cex, srt = srt, fill = NA)
-TPS <- c(describe.property(c("T", "P"), c(T, "Psat")), expression(sum(S) == 0.01*m))
-legend("topright", TPS, bty = "n")
-title("Cu-Fe-S-O-H (minerals and aqueous species)", font.main = 1)
+legend("topright", c(TPexpr, Sexpr, NaClexpr, aqexpr), bty = "n")
+title("Cu-Fe-S-O-H-Cl (minerals and aqueous species)", font.main = 1)
 
 # Restore default OBIGT database
 OBIGT()



More information about the CHNOSZ-commits mailing list