reset()
## list the buffers
thermo()$buffer
# another way to do it, for a specific buffer
print(mod.buffer("PPM"))
## buffer made of one species
# calculate the activity of CO2 in equilibrium with
# (a buffer made of) acetic acid at a given activity
basis("CHNOS")
basis("CO2", "AC")
# what activity of acetic acid are we using?
print(mod.buffer("AC"))
# return the activity of CO2
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
affinity(O2 = c(-85, -70, 4), T = c(25, 100, 4), return.buffer = TRUE)
# 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)
## 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))
# 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
# setup H2S and O2 to be buffered by 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)
# 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
basis("CHNOSe")
# here we will buffer the activity of the electron by O2
mod.buffer("O2", "O2", "gas", 999)
basis("e-", "O2")
# start our plot, then loop over values of logfO2
thermo.plot.new(xlim = c(0, 14), ylim = c(-0.8, 1.2),
xlab = "pH",ylab = axis.label("Eh"))
# the upper and lower lines correspond to the upper
# and lower stability limits of water
logfO2 <- c(0, -20, -40, -60, -83.1)
for(i in 1:5) {
# update the logarithm of fugacity (logact) of O2 in the buffer
mod.buffer("O2", "O2", "gas", logfO2[i])
# get the values of the logarithm of activity of the electron
a <- affinity(pH = c(0, 14, 15), return.buffer = TRUE)
# convert values of pe (-logact of the electron) to Eh
Eh <- convert(-as.numeric(a$`e-`), "Eh")
lines(seq(0, 14, length.out = 15), Eh)
# add some labels
text(seq(0, 14, length.out = 15)[i*2+2], Eh[i*2+2],
paste("logfO2 =", logfO2[i]))
}
title(main = paste("Relation between logfO2(g), Eh and pH at\n",
"25 degC and 1 bar. After Garrels, 1960"))
## buffer made of two species
# conditions for metastable equilibrium among
# CO2 and acetic acid. note their starting activities:
print(mod.buffer("CO2-AC"))
basis("CHNOS")
basis("O2", "CO2-AC")
affinity(return.buffer = TRUE) # logfO2 = -75.94248
basis("CO2", 123) # what the buffer reactions are balanced on
affinity(return.buffer = TRUE) # unchanged
# consider more oxidizing conditions
mod.buffer("CO2-AC", logact = c(0, -10))
affinity(return.buffer = TRUE)
Run the code above in your browser using DataLab