[CHNOSZ-commits] r633 - in pkg/CHNOSZ: . inst inst/extdata/OBIGT vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri Feb 19 08:11:24 CET 2021
Author: jedick
Date: 2021-02-19 08:11:24 +0100 (Fri, 19 Feb 2021)
New Revision: 633
Modified:
pkg/CHNOSZ/DESCRIPTION
pkg/CHNOSZ/inst/NEWS.Rd
pkg/CHNOSZ/inst/TODO
pkg/CHNOSZ/inst/extdata/OBIGT/Berman_cr.csv
pkg/CHNOSZ/vignettes/multi-metal.Rmd
Log:
OBIGT/Berman_cr.csv: Use standard abbreviations for magnetite and hematite
Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION 2020-12-21 05:18:34 UTC (rev 632)
+++ pkg/CHNOSZ/DESCRIPTION 2021-02-19 07:11:24 UTC (rev 633)
@@ -1,6 +1,6 @@
-Date: 2020-12-21
+Date: 2021-02-19
Package: CHNOSZ
-Version: 1.4.0-3
+Version: 1.4.0-4
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/inst/NEWS.Rd
===================================================================
--- pkg/CHNOSZ/inst/NEWS.Rd 2020-12-21 05:18:34 UTC (rev 632)
+++ pkg/CHNOSZ/inst/NEWS.Rd 2021-02-19 07:11:24 UTC (rev 633)
@@ -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.0-3 (2020-12-21)}{
+\section{Changes in CHNOSZ version 1.4.0-4 (2021-02-19)}{
\itemize{
@@ -26,6 +26,9 @@
can now be used to draw contours or a color image showing the activities of
the predominant species.
+ \item Use standard abbreviations for hematite (Hem) and magnetite (Mag) in
+ OBIGT/Berman_cr.csv.
+
}
}
Modified: pkg/CHNOSZ/inst/TODO
===================================================================
--- pkg/CHNOSZ/inst/TODO 2020-12-21 05:18:34 UTC (rev 632)
+++ pkg/CHNOSZ/inst/TODO 2021-02-19 07:11:24 UTC (rev 633)
@@ -35,5 +35,46 @@
- Show species name in message: mod.OBIGT(5, G = 0)
- Check that protein ionization calculations are not affected by E.units()
+library(testthat)
+reset()
+basis(c("CO2", "H2", "NH4+", "H2O", "H2S", "H+"))
+a1 <- affinity(iprotein = 1)
+E.units("J")
+a2 <- affinity(iprotein = 1)
+expect_equal(a1$values, a2$values)
- diagram(): scale predominant.values in output by balancing coefficients?
+
+- add HKF parameters from LaRowe and Amend, 2016 (https://doi.org/10.1038/ismej.2015.227)
+
+[20210115]
+
+BUG: Why doesn't this work?
+> basis("QEC")
+> subcrt("LYSC_CHICK", 1/129)
+Error in count.elements(formula) :
+ 'C2.15503876355161e-07H3.72868217191069e-07O1.72868217157562eZ-7' is not a simple chemical formula
+
+[20210116]
+
+- Add an option to disable use of ionize.aa() (it's always active if the basis species have H+)
+
+[20210122] Fix incorrect states in messages:
+> subcrt("C10H22")
+info.character: found C10H22(decane); other available phases are gas, gas, gas, liq, liq, liq
+> subcrt("C10H22", "liq")
+info.character: found C10H22(liq); other available states are liq, liq
+
+[20210201]
+
+- Move alunite from SUPCRT92.csv to inorganic_cr.csv.
+
+[20210212]
+
+- Fix formatting of messages from palply() (seen with read.fasta())
+
+[20210219]
+
+- Don't use data(thermo) in demo Shiny app (citrate)
+
+- Check all mineral abbreviations in Berman_cr.csv.
Modified: pkg/CHNOSZ/inst/extdata/OBIGT/Berman_cr.csv
===================================================================
--- pkg/CHNOSZ/inst/extdata/OBIGT/Berman_cr.csv 2020-12-21 05:18:34 UTC (rev 632)
+++ pkg/CHNOSZ/inst/extdata/OBIGT/Berman_cr.csv 2021-02-19 07:11:24 UTC (rev 633)
@@ -29,7 +29,7 @@
forsterite,Fo,Mg2SiO4,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
gehlenite,Ge,Al2Ca2SiO7,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
grossular,Grs,Ca3Al2Si3O12,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
-hematite,Hm,Fe2O3,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
+hematite,Hem,Fe2O3,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
ilmenite,Ilm,FeTiO3,cr,Ber88,Ber90.1,2017-10-03,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
jadeite,Jd,NaAlSi2O6,cr,Ber88,SHD91,2017-10-03,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
kaolinite,Kln,Al2Si2O9H4,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
@@ -37,7 +37,7 @@
lawsonite,Lw,CaAl2Si2O10H4,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
lime,Lm,CaO,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
magnesite,Mst,MgCO3,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
-magnetite,Mt,Fe3O4,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
+magnetite,Mag,Fe3O4,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
margarite,Mrg,CaAl4Si2O12H2,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
meionite,Me,Ca4Al6Si6O27C,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
merwinite,Mw,Ca3MgSi2O8,cr,Ber88,NA,2017-10-01,cal,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
Modified: pkg/CHNOSZ/vignettes/multi-metal.Rmd
===================================================================
--- pkg/CHNOSZ/vignettes/multi-metal.Rmd 2020-12-21 05:18:34 UTC (rev 632)
+++ pkg/CHNOSZ/vignettes/multi-metal.Rmd 2021-02-19 07:11:24 UTC (rev 633)
@@ -93,18 +93,21 @@
basis("CHNOS+")
species(c("CH4", "CO2", "HCO3-", "CO3-2"))
aC <- affinity(pH = c(0, 14), O2 = c(-75, -60))
-dC <- diagram(aC)
+dC <- diagram(aC, dx = c(0, 1, 0, 0), dy = c(0, 1, 0, 0))
species(c("H2S", "HS-", "HSO4-", "SO4-2"))
aS <- affinity(aC) # argument recall
-dS <- diagram(aS, add = TRUE, col = 4, col.names = 4)
+dS <- diagram(aS, add = TRUE, col = 4, col.names = 4, dx = c(0, -0.5, 0, 0))
aCS <- mash(dC, dS)
-diagram(aCS)
+srt <- c(0, 0, 90, 0, 0, 0, 90, 0, 0, 0)
+cex.names <- c(1, 1, 0.8, 1, 1, 1, 1, 1, 1, 1)
+dy <- c(0, 0, 0, -0.2, 0, 0, 0, 0, 0, 0)
+diagram(aCS, srt = srt, cex.names = cex.names, dy = dy)
legend("topright", legend = lTP(25, 1), bty = "n")
```
The second diagram is just like the first, except the function `mash()` is used to label the fields with names of species from both systems, and a legend is added to indicate the temperature and pressure.
-```{r mash, echo = 9:11, results = "hide", message = FALSE, fig.width = 10, fig.height = 5, out.width = "100%"}
+```{r mash, echo = 9:14, results = "hide", message = FALSE, fig.width = 10, fig.height = 5, out.width = "100%"}
```
Note that these are predominance diagrams, so they show only the species with highest activity; there is in fact a distribution of activities of aqueous species that is not visible here.
@@ -158,7 +161,8 @@
Fe.cr <- eVtoJ(Fe.cr)
V.cr <- eVtoJ(V.cr)
-# Gibbs energies of formation (J / mol) for aqueous species from Wagman et al.
+# Gibbs energies of formation (J / mol) for aqueous species
+# Most are from Wagman et al., 1982
Fe.aq <- 1000 * c("Fe+2" = -78.90, "Fe+3" = -4.7, "FeO2-2" = -295.3,
"FeOH+" = -277.4, "FeOH+2" = -229.41, "HFeO2-" = -377.7,
"Fe(OH)2+" = -438.0, "Fe(OH)3" = -659.3,
@@ -208,7 +212,7 @@
```
</div>
-This code produce species indices in the OBIGT database for Fe- and V-bearing aqueous species (`iFe.aq`, `iV.aq`), solids (`iFe.cr`, `iV.cr`), and bimetallic solids (`iFeV.cr`), which are used in the following diagrams.
+This code produces species indices in the OBIGT database for Fe- and V-bearing aqueous species (`iFe.aq`, `iV.aq`), solids (`iFe.cr`, `iV.cr`), and bimetallic solids (`iFeV.cr`), which are used in the following diagrams.
Now we set up the plot area and assign activities of aqueous species to 10^-5^, which is the default value for diagrams on the MP website (from the page for a material: "Generate Phase Diagram" -- "Aqueous Stability (Pourbaix)").
The following commands compute Eh-pH diagrams for the single-metal systems Fe-O-H and V-O-H.
@@ -237,16 +241,31 @@
dFeV <- diagram(aFeV, plot.it = FALSE, bold = TRUE)
# 1:1 mixture (Fe:V)
a11 <- mix(dFe, dV, dFeV, c(1, 1))
-diagram(a11, min.area = 0.01)
+# Adjust labels 20210219
+iV2O3 <- info("V2O3")
+iFeO <- info("FeO", "cr")
+iFe3V <- info("Fe3V")
+srt <- rep(0, nrow(a11$species))
+srt[a11$species$ispecies == paste(iFeO, iV2O3, sep = ",")] <- 90
+srt[a11$species$ispecies == paste(iV2O3, iFe3V, sep = ",")] <- -13
+diagram(a11, min.area = 0.01, srt = srt)
title("Fe:V = 1:1")
label.figure(lTP(25, 1), xfrac = 0.12)
# 1:3 mixture
a13 <- mix(dFe, dV, dFeV, c(1, 3))
-diagram(a13, min.area = 0.01)
+srt <- rep(0, nrow(a13$species))
+srt[a13$species$ispecies == paste(iFeO, iV2O3, sep = ",")] <- 90
+srt[a13$species$ispecies == paste(iV2O3, iFe3V, sep = ",")] <- -13
+diagram(a13, min.area = 0.01, srt = srt)
title("Fe:V = 1:3")
# 1:5 mixture
a15 <- mix(dFe, dV, dFeV, c(1, 5))
-diagram(a15, min.area = 0.01)
+iFeV3 <- info("FeV3")
+srt <- rep(0, nrow(a15$species))
+srt[a15$species$ispecies == paste(iFeO, iV2O3, sep = ",")] <- 90
+srt[a15$species$ispecies == paste(iV2O3, iFe3V, sep = ",")] <- -13
+srt[a15$species$ispecies == paste(iV2O3, iFeV3, sep = ",")] <- -13
+diagram(a15, min.area = 0.01, srt = srt)
title("Fe:V = 1:5")
```
@@ -258,7 +277,7 @@
Finally, the `diagram()`s are plotted; the `min.area` argument is used to remove labels for very small fields.
Regarding the legend, it should be noted that although the DFT calculations for solids are made for zero temperature and zero pressure [@SZS_17], the standard Gibbs energies of aqueous species [e.g. @WEP_82] are modified by a correction term so that they can be combined with DFT energies to reproduce the experimental energy for dissolution of a representative material for each metal at 25 °C and 1 bar [@PWLC12].
-```{r mixing1, echo = 16:37, message = FALSE, results = "hide", fig.width = 9, fig.height = 3, out.width = "100%", out.extra='class="full-width"', pngquant = pngquant}
+```{r mixing1, echo = 16:47, message = FALSE, results = "hide", fig.width = 9, fig.height = 3, out.width = "100%", out.extra='class="full-width"', pngquant = pngquant}
```
In these diagrams, changing the Fe:V ratio affects the fully reduced metallic species.
@@ -269,7 +288,7 @@
Let's make another diagram for the 1:1 Fe:V composition over a broader range of Eh and pH.
The diagram shows a stable assemblage of Fe~2~O~3~ with an oxidized bimetallic material, [Fe~2~V~4~O~13~](https://materialsproject.org/materials/mp-1200054/).
-```{r FeVO4, eval = FALSE, echo = 1:20}
+```{r FeVO4, eval = FALSE, echo = 1:28}
par(mfrow = c(1, 2))
# Fe-bearing species
basis(c("VO+2", "Fe+2", "H2O", "e-", "H+"))
@@ -288,7 +307,15 @@
dFeV <- diagram(aFeV, plot.it = FALSE, bold = TRUE)
# 1:1 mixture (Fe:V)
a11 <- mix(dFe, dV, dFeV, c(1, 1))
-d11 <- diagram(a11, min.area = 0.01)
+# Adjust labels 20210219
+iV2O3 <- info("V2O3")
+iFe3V <- info("Fe3V")
+iVO4m3 <- info("VO4-3")
+iFe2O3 <- info("Fe2O3")
+srt <- rep(0, nrow(a11$species))
+srt[a11$species$ispecies == paste(iV2O3, iFe3V, sep = ",")] <- -13
+srt[a11$species$ispecies == paste(iFe2O3, iVO4m3, sep = ",")] <- 90
+d11 <- diagram(a11, min.area = 0.01, srt = srt)
water.lines(d11, col = "orangered")
# Calculate affinity of FeVO4
@@ -317,7 +344,7 @@
Given the diagram for the stable Fe-, V- and bimetallic materials *mixed with the same stoichiometry* as FeVO~4~ (1:1 Fe:V), the difference between their affinities of formation and that of FeVO~4~ corresponds to the Pourbaix energy difference (-Δ*G*~pbx~).
This is plotted as a color map in the second diagram.
-```{r FeVO4, echo = 22:34, message = FALSE, results = "hide", fig.width = 10, fig.height = 5, out.width = "100%", pngquant = pngquant}
+```{r FeVO4, echo = 30:42, message = FALSE, results = "hide", fig.width = 10, fig.height = 5, out.width = "100%", pngquant = pngquant}
```
Now we locate the pH and Eh that maximize the affinity (that is, minimize Δ*G*~pbx~) of FeVO~4~ compared to the stable species.
@@ -382,11 +409,12 @@
logaH2S <- -2
T <- 200
pH <- c(0, 14, 500)
-O2 <- c(-45, -30, 500)
+O2 <- c(-46, -31, 500)
basis(c("Cu+", "pyrite", "H2S", "oxygen", "H2O", "H+"))
basis("H2S", logaH2S)
S.aq <- c("H2S", "HS-", "HSO4-", "SO4-2")
Fe.cr <- c("pyrite", "pyrrhotite", "magnetite", "hematite")
+Fe.abbrv <- c("Py", "Po", "Mag", "Hem")
```
Now we calculate affinities for minerals in the Fe-S-O-H system that take account of the changing aqueous sulfur species in `S.aq`.
@@ -399,16 +427,17 @@
```{r stack1_2, eval = FALSE, echo = 1:4}
species(Fe.cr)
mFe <- mosaic(S.aq, pH = pH, O2 = O2, T = T)
-diagram(mFe$A.bases, lty = 2, col = 4, col.names = 4, italic = TRUE)
-dFe <- diagram(mFe$A.species, add = TRUE, lwd = 2)
+diagram(mFe$A.bases, lty = 2, col = 4, col.names = 4, italic = TRUE, dx = c(0, 1, 0, 0), dy = c(-1, 0, 0, 0))
+dFe <- diagram(mFe$A.species, add = TRUE, lwd = 2, names = Fe.abbrv, dy = c(0, 0, 1, 0))
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))
mFeCu <- mosaic(list(S.aq, Fe.cr), pH = pH, O2 = O2,
T = T, stable = list(NULL, dFe$predominant))
-diagram(mFeCu$A.species, add = TRUE, col = 2, col.names = 2, bold = TRUE)
+diagram(mFeCu$A.species, add = TRUE, col = 2, col.names = 2, bold = TRUE, names = FeCu.abbrv)
col <- c("#FF8C00", rep(NA, 6))
-diagram(mFeCu$A.species, add = TRUE, col = col, lwd = 2, col.names = col, bold = TRUE)
+diagram(mFeCu$A.species, add = TRUE, col = col, lwd = 2, col.names = col, bold = TRUE, names = FeCu.abbrv)
TP <- describe.property(c("T", "P"), c(T, "Psat"))
legend("topright", TP, bty = "n")
title(paste("Cu-Fe-S-O-H; Total S =", 10^logaH2S, "m"))
@@ -425,11 +454,11 @@
After that we add the legend and title.
-```{r stack1_2, echo=5:15, results = "hide", message = FALSE, fig.width = 6, fig.height = 5, out.width = "75%", fig.align = "center", pngquant = pngquant}
+```{r stack1_2, echo=5:16, results = "hide", message = FALSE, fig.width = 6, fig.height = 5, out.width = "75%", fig.align = "center", pngquant = pngquant}
```
This diagram has a distinctive chalcopyrite "hook" surrounded by a thin bornite field.
-Only the chalcopyrite-bornite reaction in the pyrite field is shown in many published diagrams [e.g. @And75;@Gio02], but diagrams with a similar chalcopyrite wedge or hook geometry can be seen in @BBR77 and @Bri80.
+Only the chalcopyrite-bornite reaction in the pyrite field is shown in some published diagrams [e.g. @And75;@Gio02], but diagrams with a similar chalcopyrite wedge or hook geometry can be seen in @BBR77 and @Bri80.
## Mixing 2
@@ -458,7 +487,7 @@
species(Fe.cr)
mFe <- mosaic(S.aq, pH = pH, O2 = O2, T = T)
-diagram(mFe$A.bases, lty = 2, col = 4, col.names = 4, italic = TRUE)
+diagram(mFe$A.bases, lty = 2, col = 4, col.names = 4, italic = TRUE, dx = c(0, 1, 0, 0), dy = c(-1, 0, -2, 0))
names <- info(info(Fe.cr))$abbrv
dFe <- diagram(mFe$A.species, add = TRUE, names = names)
@@ -481,18 +510,21 @@
ifelse(grepl("Ccp", a$species$name), "#FF8C0088",
ifelse(grepl("Bn", a$species$name), "#DC143C88", NA)
)}
-srt <- function(a) ifelse(a$species$name %in% c("Mt+Bn", "Mt+Cct", "Mt+Ccp"), 80, 0)
+srt <- function(a) ifelse(a$species$name %in% c("Mag+Bn", "Mag+Ccp"), 80, 0)
+cex.names <- function(a) ifelse(a$species$name %in% c("Mag+Bn", "Mag+Ccp", "Mag+Cct"), 0.8, 1)
+dx <- function(a) sapply(a$species$name, switch, "Mag+Bn" = 0.15, "Mag+Cct" = 0.5, 0)
+dy <- function(a) sapply(a$species$name, switch, "Py+Bn" = -1, "Po+Bn" = -0.8, "Po+Cu" = -0.8, "Mag+Ccp" = -1, 0)
a11 <- mix(dFe, dCu, dFeCu)
-diagram(a11, fill = fill(a11), srt = srt(a11), min.area = 0.01)
+diagram(a11, fill = fill(a11), srt = srt(a11), min.area = 0.01, dx = dx(a11), dy = dy(a11), cex.names = cex.names(a11))
title("Fe:Cu = 1:1")
a21 <- mix(dFe, dCu, dFeCu, c(2, 1))
-diagram(a21, fill = fill(a21), srt = srt(a21), min.area = 0.01)
+diagram(a21, fill = fill(a21), srt = srt(a21), min.area = 0.01, dx = dx(a21), dy = dy(a21), cex.names = cex.names(a21))
title("Fe:Cu = 2:1")
a12 <- mix(dFe, dCu, dFeCu, c(1, 2))
-diagram(a12, fill = fill(a12), srt = srt(a12), min.area = 0.01)
+diagram(a12, fill = fill(a12), srt = srt(a12), min.area = 0.01, dx = dx(a12), dy = dy(a12), cex.names = cex.names(a12))
title("Fe:Cu = 1:2")
```
</div>
@@ -525,7 +557,7 @@
aFe <- affinity(S2 = c(-34, -10), O2 = c(-55, -40), T = T)
# Order species by a function of composition to make colors
oFe <- order(aFe$species$S2 - aFe$species$O2)
-fill <- terrain.colors(length(oFe), alpha = 0.3)[oFe]
+fill <- terrain.colors(length(oFe), alpha = 0.2)[oFe]
abbrv <- info(aFe$species$ispecies)$abbrv
dFe <- diagram(aFe, names = abbrv, fill = fill)
title("Fe-S-O-H")
@@ -536,9 +568,9 @@
mCu <- mosaic(bFe, S2 = c(-34, -10), O2 = c(-55, -40),
T = T, stable = list(dFe$predominant))
oCu <- order(mCu$A.species$species$S2 - mCu$A.species$species$O2)
-fill <- terrain.colors(length(oCu), alpha = 0.3)[oCu]
+fill <- terrain.colors(length(oCu), alpha = 0.2)[oCu]
abbrv <- info(mCu$A.species$species$ispecies)$abbrv
-dCu <- diagram(mCu$A.species, names = abbrv, fill = fill)
+dCu <- diagram(mCu$A.species, names = abbrv, fill = fill, dx = c(0, 0, 0, 0, 1.8))
title("Cu-Fe-S-O-H")
# Mash the diagrams and adjust labels
@@ -545,14 +577,14 @@
aFeCu <- mash(dFe, dCu)
names <- aFeCu$species$name
srt <- rep(0, length(names))
-srt[names %in% c("Mt+Cu", "Hm+Cu")] <- 90
-srt[names %in% c("Mt+Bn", "Hm+Bn")] <- 63
-srt[names %in% c("Mt+Ccp", "Hm+Ccp")] <- 66
+srt[names %in% c("Mag+Cu", "Hem+Cu")] <- 90
+srt[names %in% c("Mag+Bn", "Hem+Bn")] <- 63
+srt[names %in% c("Mag+Ccp", "Hem+Ccp")] <- 66
srt[names %in% c("Py+Bn", "Py+Cv")] <- 90
cex <- rep(1, length(names))
-cex[names %in% c("Hm+Ccp", "Hm+Cv")] <- 0.8
+cex[names %in% c("Hem+Ccp", "Hem+Cv")] <- 0.8
oFeCu <- order(aFeCu$species$S2 - aFeCu$species$O2)
-fill <- terrain.colors(length(oFeCu), alpha = 0.3)[oFeCu]
+fill <- terrain.colors(length(oFeCu), alpha = 0.2)[oFeCu]
diagram(aFeCu, cex.names = cex, srt = srt, fill = fill)
legend("topleft", legend = lTP(T, "Psat"), bg = "white")
title("Cu-Fe-S-O-H")
@@ -702,22 +734,25 @@
# Only Fe-bearing minerals
species(c("pyrite", "pyrrhotite", "magnetite", "hematite"))
aFe <- affinity("Fe+2" = c(0, 12), O2 = c(-40, -16), T = 400, P = 2000)
-dFe <- diagram(aFe, xlab = xlab)
-title(paste("Only Fe; 1° balance:", dFe$balance))
+names <- info(aFe$species$ispecies)$abbrv
+dFe <- diagram(aFe, xlab = xlab, names = names)
+title(bquote("Only Fe; 1° balance:" ~ .(expr.species(dFe$balance))))
label.plot("A")
# Only Cu-bearing minerals
species(c("covellite", "chalcocite", "tenorite", "cuprite"))
aCu <- affinity(aFe) # argument recall
-dCu <- diagram(aCu, xlab = xlab)
-title(paste("Only Cu; 1° balance:", dCu$balance))
+names <- info(aCu$species$ispecies)$abbrv
+dCu <- diagram(aCu, xlab = xlab, names = names)
+title(bquote("Only Cu; 1° balance:" ~ .(expr.species(dCu$balance))))
label.plot("B")
# Only Fe- AND Cu-bearing minerals
species(c("chalcopyrite", "bornite"))
aFeCu <- affinity(aFe)
-dFeCu <- diagram(aFeCu, xlab = xlab, balance = "H+")
-title(paste("Only Fe+Cu; 1° balance:", dFeCu$balance))
+names <- info(aFeCu$species$ispecies)$abbrv
+dFeCu <- diagram(aFeCu, xlab = xlab, balance = "H+", names = names)
+title(bquote("Only Fe+Cu; 1° balance:" ~ .(expr.species(dFeCu$balance))))
label.plot("C")
### SECONDARY balancing
@@ -724,16 +759,17 @@
# Fe- or Cu-bearing minerals
ad1 <- rebalance(dFe, dCu, balance = "H+")
-d1 <- diagram(ad1, xlab = xlab, balance = 1)
-title("Only Fe or Cu; 2° balance: H+")
+names <- info(ad1$species$ispecies)$abbrv
+d1 <- diagram(ad1, xlab = xlab, balance = 1, names = names)
+title(bquote("Only Fe or Cu; 2° balance:" ~ .(expr.species("H+"))))
label.plot("D")
# All minerals
d1$values <- c(dFe$values, dCu$values)
ad2 <- rebalance(d1, dFeCu, balance = "H+")
-abbrv <- info(ad2$species$ispecies)$abbrv
-diagram(ad2, xlab = xlab, balance = 1, names = abbrv)
-title("Fe and/or Cu; 2° balance: H+")
+names <- info(ad2$species$ispecies)$abbrv
+diagram(ad2, xlab = xlab, balance = 1, names = names)
+title(bquote("Fe and/or Cu; 2° balance:" ~ .(expr.species("H+"))))
label.plot("E")
db <- describe.basis(ibasis = 3)
@@ -750,7 +786,7 @@
We then calculate the diagrams for the primary balancing coefficients, for the groups of only Fe-, only Cu-, and only Fe+Cu-bearing minerals.
It is obvious that the first two systems are balanced on Fe and Cu, respectively, but the third has a somewhat unusual balance: H^+^.
See Reaction 4 of @MH85 for an example.
-```{r rebalance, eval = FALSE, echo = 10:30}
+```{r rebalance, eval = FALSE, echo = 10:33}
```
Now comes the secondary balancing, where all reactions, not only that between bornite and chalcopyrite, are balanced on H^+^.
@@ -760,7 +796,7 @@
Then we rebalance diagrams D and C to make the final diagram in E.
The fields in this diagram are labeled with mineral abbreviations from the OBIGT database.
-```{r rebalance, eval = FALSE, echo = 33:43}
+```{r rebalance, eval = FALSE, echo = 36:49}
```
```{r rebalance, echo = FALSE, results = "hide", message = FALSE, fig.width = 6.5, fig.height = 5, out.width = "100%", fig.align = "center", pngquant = pngquant}
@@ -789,12 +825,15 @@
species(c("pyrite", "magnetite", "hematite", "covellite", "tenorite",
"chalcopyrite", "bornite"))
a <- affinity("Cu+" = c(-8, 2, 500), "Fe+2" = c(-4, 12, 500), T = 400, P = 2000)
-d <- diagram(a, xlab = ratlab("Cu+"), ylab = ratlab("Fe+2"), balance = "O2")
-title(paste("Cu-Fe-S-O-H; 1° balance:", d$balance))
+names <- info(a$species$ispecies)$abbrv
+d <- diagram(a, xlab = ratlab("Cu+"), ylab = ratlab("Fe+2"), balance = "O2", names = names)
+title(bquote("Cu-Fe-S-O-H; 1° balance:" ~ .(expr.species(d$balance))))
# Add saturation lines
species(c("pyrrhotite", "ferrous-oxide", "chalcocite", "cuprite"))
asat <- affinity(a) # argument recall
-diagram(asat, type = "saturation", add = TRUE, lty = 2, col = 4)
+names <- asat$species$name
+names[2] <- "ferrous oxide"
+diagram(asat, type = "saturation", add = TRUE, lty = 2, col = 4, names = names)
legend("topleft", legend = lTP(400, 2000), bty = "n")
```
@@ -839,5 +878,6 @@
## Document History
* 2020-07-15 First version.
+* 2021-02-19 Improve mineral abbreviations and placement of labels.
## References
More information about the CHNOSZ-commits
mailing list