affinity {CHNOSZ} | R Documentation |
Calculate the chemical affinities of formation reactions of species. Do it for single values of temperature, pressure, ionic strength and chemical activities of the basis species, or as a function of one or two of these variables. Or, return other properties including standard molal Gibbs energies of basis species and species of interest, and standard molal Gibbs energies, equilibrium constants and reaction activity products of formation reactions.
energy(what, vars, vals, lims, T=thermo$opt$Tr, P="Psat", IS=0, sout=NULL, do.phases=FALSE) energy.args(args, quiet = FALSE) affinity(..., property=NULL, sout=NULL, do.phases=FALSE, return.buffer=FALSE, balance="PBB", quiet=FALSE, iprotein=NULL, logact.protein=-3)
what |
character, name of property to calculate. |
vars |
character, names of variables over which to calculate a property. |
vals |
list of numeric, grid values for each variable. |
lims |
list of numeric, limits of the values for each variable. |
T |
numeric, temperature. Default is to take the temperature from thermo$opt$Tr , which corresponds to 25 degrees C. |
P |
numeric, pressure, or character "Psat" (default), which denotes 1 bar or the saturation vapor pressure of water above 100 degrees C (see water ). |
IS |
numeric, ionic strength default is 0 mol kg^-1. |
sout |
list, output from subcrt function. |
do.phases |
logical, allow subcrt to compute properties for stable phases? |
args |
list, defines the variables over which to calculate properties. |
quiet |
logical, print fewer messages to the screen? |
... |
numeric, zero or more named arguments, used to identify the variables of interest in the calculations. |
property |
character, denoting the property to be calculated. Default is A, for chemical affinity of formation reactions of species of interest. |
return.buffer |
logical. If TRUE , and a buffer has been associated with one or more basis species in the system, return the values of the activities of the basis species calculated using the buffer (it is not necessary in this case to have defined any species of interest). Default is FALSE . |
balance |
character. This argument is passed to buffer to identify a conserved basis species or other quantity in a chemical activity buffer. Default is PBB. |
iprotein |
numeric, indices of proteins in thermo$protein for which to calculate properties. |
logact.protein |
numeric, logarithms of activities of proteins identified in protein . |
The chemical affinity relates the driving force for a reaction to proceed in a forward or reverse direction to the standard Gibbs energy of the reaction and to the physical and chemical environment. The calculation of chemical affinity A is expressed by A=RT*ln(K/Q) (Kondepudi and Prigogine, 1998), where K denotes the equilibrium constant of the reaction and Q stands for the activity product of the species in the reaction. (The equilibrium constant is related to standard Gibbs energy of reaction, DeltaG0r, by DeltaG0r = -2.303*RT*logK, where R and T stand for, respectively, the gas constant and temperature). With the approach of a given reaction to a state of equilibrium, the chemical affinity tends toward zero, or K = Q.
In normal practice, the user calls affinity
in order to calculate the chemical affinities of formation reactions of species of interest. This function in turn calls energy.args
and energy
to perform the calculations. These functions rely on the global definitions of the basis
species and species
of interest.
energy
is the workhorse of the affinity
calculations. Given n (which can be zero, one, or more) names of basis species and/or T, P, or IS as the var
s, it calculates the property given in what
on an n-dimensional grid for each of the values (vals
) of the corresponding variable. The limits for each variable given in lims
indicate the minimum and maximum value and, if a third value is supplied, the resolution, or number of points in the given dimension. If T, P, and/or IS are not among the var
s, their constant values can be supplied in T
(in Kelvin), P
(in bar, or Psat), and IS
(in mol kg^-1). sout
, if provided, replaces the call to subcrt
which can greatly speed up the calculations if this intermediate step is stored by other functions (e.g., transfer
). do.phases
is passed to subcrt
so that the properties of stable mineral phases as a function of temperature can optionally be calculated.
energy.args
is used by affinity
to generate the argument list for energy
. affinity
passes the ...
list as args
to energy.args
, which returns an argument list appropriate for energy
. energy.args
also has the job of converting Eh to pe as a temperature-dependent function (see convert
), and converting pe and pH to logarithms of activities of the electron and protein, respectively (i.e., negating the values).
The property
argument to affinity
is analogous to the what
argument of energy
. Valid properties are A or NULL for chemical affinity, logK or logQ for logarithm of equilibrium constant and reaction activity product, or any of the properties available in subcrt
except for rho. The properties returned are those of the formation reactions of the species of interest from the basis species. It is also possible to calculate the properties of the species of interest themselves (not their formation reactions) by setting the property
to G.species, Cp.species, etc.
Zero, one, or more leading arguments to the function identify which of the chemical activities of basis species, temperature, pressure and/or ionic strength to vary. The names of each of these arguments may be the formula of any of the basis species of the system, or T, P, pe, pH, Eh, or IS (but names may not be repeated). The value of each argument is of the form c(min,max)
or c(min,max,res)
where min
and max
refer to the minimimum and maximum values of variable identified by the name of the argument, and res
denotes the resolution, or number of points along which to do the calculations (this is assigned a default value of 128 if it is missing, or 72 if T and P are both variables). For any arguments that specify the names of basis species, the numerical values refer to the logarithms of the activities of that basis species, or logarithms of fugacities if it is a gas. Unlike the energy
function, the units of T and P in affinity
are those the user has set using nuts
(on program start-up these are degrees C and bar, respectively).
If proteins are in the species definition (which are distinguished from other species by having an underscore character in the name), and the basis species are charged, and thermo$opt$ionize
is TRUE, affinity
calls ionize
to calculate the properties of charged proteins. If one or more buffers are assigned to the definition of basis
species, affinity
calls buffer
to calculate the logarithms of activities of these basis species from the buffer.
The commutative properties of the protein group additivity algorithm can be utilized in order to calculate the properties of many proteins in a fraction of the time it would take to calculate them individually. The iprotein
and logact.protein
arguments can be used to compute the chemical affinities of formation of proteins that are not in the global species
definition. The appropriate basis
species (e.g., CHNOSZ or CHNOSZ+) still must be defined prior to calling affinity
. The user can give the indices of desired proteins contained in thermo$protein
and affinity
will automatically add to the species list the amino acid residues and protein backbone group (found as “proteins” in thermo$protein
) then calculate the properties of the reactions for the residues, and add them together to get those of the indicated proteins. The iprotein
option is compatible with calculations for ionize
d proteins.
The values of chemical affinity of formation reactions of the species of interest are returned in dimensionless (base 10) units (i.e., A/2.303RT).
To use the names of charged basis species such as K+ and SO4-2 as the arguments, they should be enclosed in quotes (see the example for for aluminum speciation in diagram
).
For energy
, a list the first element of which is sout
(the intermediate results from subcrt
) and the second element of which is a
, which contains the calculated properties. The latter itself is a list, one element for each species of interest,the elements of which have dimensions that are the number of var
iables passed to the function.
For energy.arg
, a list with elements what
, vars
, vals
, lims
, T
, P
, IS
that are each of the type for the corresponding arguments in energy
.
For affinity
, a list, elements of which are sout
, property
(name of the calculated property), basis
and species
(definition of basis species and species of interest in effect at runtime), T
and P
(temperature and pressure, in the system units of Kelvin and bar, of length two (corresponding to minimum/maximum values) if either one is a variable of interest or length one otherwise), xname
(the name of the first variable of interest, "" if none is present), xlim
(if a first variable of interest is present, numeric of length 3 specifying the (minimum, maximum, resolution) of this variable), yname
(the name of the second variable of interest, "" if none present), ylim
(analogous to xlim
but for a second variable of interest), values
. The latter is itself a list, each element of which corresponds to a species of interest (names of the elements in this list are the character versions of the index of the species in thermo$obigt
) and which contains a numeric value if the number of variables of interest is zero or matrix otherwise.
If pe
or Eh
, or pH
is/are among the variables of interest, xnames
and/or ynames
become e- or H+ (respective to the property) and the minimum and maximum values in xlim
and/or ylim
adjusted accordingly (using convert
).
Amend, J. P. and Shock, E. L., 1998. Energetics of amino acid synthesis in hydrothermal ecosystems. Science, 281, 1659-1662.
Amend, J. P. and Shock, E. L., 2001. Energetics of overall metabolic reactions of thermophilic and hyperthermophilic Archaea and Bacteria. FEMS Microbiol. Rev., 25, 175-243.
Kondepudi, D. K. and Prigogine, I., 1998. Modern Thermodynamics: From Heat Engines to Dissipative Structures, John Wiley & Sons, New York, 486 p.
Schulte, M. D. and Shock, E. L., 1995. Thermodynamics of Strecker synthesis in hydrothermal systems. Orig. Life Evol. Biosph., 25, 161-173.
For prerequisites to affinity
, see basis
and species
. Some calculations may invoke ionize
and buffer
. To visualize the results, use diagram
.
## set up a system and calculate ## chemical affinities of formation reactions. basis(c('SiO2','MgO','H2O','O2'),c(-5,-5,0,999)) species(c('quartz','enstatite','forsterite')) # chemical affinities (A/2.303RT) at 25 deg C and 1 bar affinity() # at higher temperature and pressure affinity(T=500,P=2000) # ten different temperatures at one pressure affinity(T=c(500,1000,10),P=2000) # at 25 temperatures and pressures affinity(T=c(500,1000,5),P=c(1000,5000,5)) # as a function of logarithm of activity of MgO affinity(MgO=c(-10,-5,10)) ## equilibrium constants of formation reactions affinity(property='logK') # Standard molal Gibbs energies of species, # dimensionless (-G/2.303RT) affinity(property='G.species') # Standard molal Gibbs energies of species, # user dimensions (default: calories) affinity(property='G') ## Activity of glycine as a function of those of ## HCN and formaldehyde (200 C, 300 bar) ## After Schulte and Shock, 1995, Fig. 5 # we can define the basis as this: basis(c('formaldehyde','H2O','HCN','O2')) species('glycine') t <- affinity(HCHO=c(-10,-2,9),HCN=c(-18,-2,9),T=200,P=300) # that gave us *affinities* (dimensionless) for logact(glycine)=-3 # (the default). we can now find the *activities* that # correspond to affinity=0 logact.glycine <- species()$logact + t$values[[1]] # note transposition of the z-value matrix in the following command contour(x=-10:-2,y=seq(-18,-2,by=2),z=t(logact.glycine), xlab=axis.label('HCHO'),ylab=axis.label('HCN'), labcex=1,xaxs='i',yaxs='i') title(main=paste('log activity glycine, after Schulte and Shock, 1995\n', describe(x=NULL),', logaH2O = 1')) ### amino acid synthesis at low and high temperatures ### after Amend and Shock, 1998 # select the basis species and species of interest # and set their activities (first for the 18 degree C case) basis(c('H2O','CO2','NH4+','H2','H+','H2S'), rep('aq',6),log10(c(1,1e-4,5e-8,2e-9,5e-9,1e-15))) species(c('alanine','argininate','asparagine','aspartate','cysteine', 'glutamate','glutamine','glycine','histidine','isoleucine', 'leucine','lysinium','methionine','phenylalanine','proline', 'serine','threonine','tryptophan','tyrosine','valine'), log10(c(3.9,0.7,1.1,3.3,0.5,3.8,1.0,5.8,1.2,0.7, 0.8,1.0,2.8,0.5,0.5,4.6,5.8,0.6,0.9,2.8)/1e9)) Tc <- 18 T <- convert(Tc,'K') # converting A (dimensionless) to G of reaction (cal/mol) # is like converting log K to standard G of reaction AS98.18 <- convert(convert(as.numeric(affinity(T=Tc)$values),'G',T=T),'J')/1000 # the 100 degree C case Tc <- 100 T <- convert(Tc,'K') basis(c('H2O','CO2','NH4+','H2','H+','H2S'), log10(c(1,2.2e-3,2.9e-6,3.4e-4,1.9e-6,1.6e-3))) species(1:20,log10(c(2.8e-9,5.0e-10,7.9e-10,2.4e-9,3.6e-10, 2.7e-9,7.2e-10,4.2e-9,8.6e-10,5.0e-10, 5.7e-10,7.2e-10,2.0e-9,3.6e-10,3.6e-10, 3.3e-9,4.2e-9,4.3e-10,6.5e-10,2.0e-9))) AS98.100 <- convert(convert(as.numeric(affinity(T=Tc)$values),'G',T=T),'J')/1000 # the nominal carbon oxidation state Z.C <- ZC(as.character(thermo$obigt$formula[thermo$species$ispecies])) # put them together print(data.frame(T100=AS98.100,T18=AS98.18,Z.C=Z.C)) # values not exactly reproducing AS98 - different amino acid parameters ## affinities of metabolic reactions ## after Amend and Shock, 2001, Fig. 7 basis(c('CO2','H2','NH3','O2','H2S','H+')) basis(c('O2','H2'),'aq') # O2 and H2 were gas species('H2O') doplot <- function(T) { res <- 20 t <- affinity(H2=c(-10,0,res),O2=c(-10,0,res),T=T) T.K <- convert(T,'K') # temperature in Kelvin t <- convert(t$values[[1]],'G',T.K) # affinities (cal/mol) t <- convert(t,'J') # affinities (Joule) contour(x=seq(-10,0,length.out=res), y=seq(-10,0,length.out=res),z=t(t/1000), labcex=1,xlab=axis.label('H2'),ylab=axis.label('O2')) } layout(matrix(c(1,1,2,3,4,5),ncol=2,byrow=TRUE),heights=c(1,4,4)) T <- c(25,55,100,150) par(mar=c(0,0,0,0)) plot.new() text(0.5,0.1,paste(c('H2(aq) + 0.5O2(aq) = H2O(liq)\n\n', 'after Amend and Shock, 2001')),cex=2) par(mar=c(3,3,0.5,0.5),cex=1.3,mgp=c(2,1,0)) for(i in 1:length(T)) doplot(T[i]) # this is so the plots in the next examples show up OK layout(matrix(1)) ## continuation of last example, affinity calculations ## in three dimensions print(affinity(H2=c(-10,0,3),O2=c(-10,0,3),T=c(25,150,4))$values)