| EOSregress {CHNOSZ} | R Documentation |
Fit experimental volumes and heat capacities using regression equations. Possible models include the Helgeson-Kirkham-Flowers (HKF) equations of state, or other equations defined using any combination of terms derived from the temperature, pressure and thermodynamic and electrostatic properties of water.
EOSregress(exptdata, var = "", T.max = 9999, ...)
EOSvar(var, T, P, ...)
EOScalc(coefficients, T, P, ...)
EOSplot(exptdata, var = NULL, T.max = 9999, T.plot = NULL,
fun.legend = "topleft", coefficients = NULL, add = FALSE,
lty = par("lty"), col=par("col"), ...)
EOSlab(var, coeff = "")
EOScoeffs(species, property, P=1)
Cp_s_var(T = 298.15, P = 1, omega.PrTr = 0, Z = 0)
V_s_var(T = 298.15, P = 1, omega.PrTr = 0, Z = 0)
exptdata |
dataframe, experimental data |
var |
character, name(s) of variables in the regression equations |
T.max |
numeric, maximum temperature for regression, in degrees Kelvin |
T |
numeric, temperature in Kelvin |
P |
numeric, pressure in bars |
... |
arguments specifying additional dependencies of the regression variables |
T.plot |
numeric, upper limit of temperature range to plot |
fun.legend |
character, where to place legend on plot |
coefficients |
dataframe, coefficients to use to make line on plot |
add |
logical, add lines to an existing plot? |
lty |
line style |
col |
color of lines |
coeff |
numeric, value of equation of state parameter for plot legend |
species |
character, name of aqueous species |
property |
character, ‘Cp’ or ‘V’ |
omega.PrTr |
numeric, value of omega at reference T and P |
Z |
numeric, charge |
EOSregress uses a linear model (lm) to regress the experimental heat capacity or volume data in exptdata, which is a data frame with columns ‘T’ (temperature in degrees Kelvin), ‘P’ (pressure in bars), and ‘Cp’ or ‘V’ (heat capacity in cal/mol.K or volume in cm3/mol).
The ‘Cp’ or ‘V’ data must be in the third column.
Only data below the temperature of T.max are included in the regression.
The regression formula is specified by a vector of names in var.
The names of the variables can be any combination of the following (listed in the order of search): variables listed in the following table, any available property of water (e.g. ‘V’, ‘alpha’, ‘QBorn’), or the name of a function that can be found using get in the default environment.
Examples of the latter are Cp_s_var, V_s_var, or functions defined by the user in the global environment; the arguments of these functions must include, but are not limited to, T and P.
T | T (temperature) |
P | P (pressure) |
TTheta | (T-\Theta) (\Theta = 228 K) |
invTTheta | 1/(T-\Theta) |
TTheta2 | (T-\Theta)^2 |
invTTheta2 | 1/(T-\Theta)^2 |
invPPsi | 1/(P+\Psi) (\Psi = 2600 bar) |
invPPsiTTheta | 1/((P+\Psi)(T-\Theta)) |
TXBorn | TX (temperature times X Born function) |
drho.dT | d\rho/dT (temperature derivative of density of water) |
V.kT | V\kappa_T (volume times isothermal compressibility of water)
|
EOSvar calculates the value of the variable named var (defined as described above) at the specified T (temperature in degrees Kelvin) and P (pressure in bars).
This function is used by EOSregress to get the values of the variables used in the regression.
EOScalc calculates the predicted heat capacities or volumes using coefficients provided by the result of EOSregress, at the temperatures and pressures specified by T and P.
EOSplot takes a table of data in exptdata, runs EOSregress and EOScalc and plots the results.
The experimental data are plotted as points, and the calculated values as a smooth line.
The point symbols are filled circles where the calculated value is within 10% of the experimental value; open circles otherwise.
EOSlab produces labels for the variables listed above that can be used as.expressions in plots.
The value of coeff is prefixed to the name of the variable (using substitute, with a multiplication symbol).
For the properties listed in the table above, and selected properties listed in water, the label is formatted using plotmath expressions (e.g., with italicized symbols and Greek letters).
If var is a user-defined function, the function can be given a ‘label’ attribute to provide plotmath-style formatting; in this case the appropriate multiplication or division symbol should be specified (see example below).
EOScoeffs retrieves coefficients in the Helgeson-Kirkham-Flowers equations from the thermodynamic database (thermo$OBIGT) for the given aqueous species.
If the property is ‘Cp’, the resulting data frame has column names of ‘(Intercept)’, ‘invTTheta2’ and ‘TX’, respectively holding the coefficients c_1, c_2 and \omega in the equation Cp^\circ = c_1 + c_2/(T-\Theta)^2 + {\omega}TX.
If the property is ‘V’, the data frame has column names of ‘(Intercept)’, ‘invTTheta’ and ‘Q’, respectively holding the coefficients \sigma, \xi and \omega in V^\circ = \sigma + \xi/(T-\Theta) - {\omega}Q.
Here, \sigma and \xi are calculated from a1, a2, a3 and a4 in thermo()$OBIGT at the pressure indicated by P (default 1 bar).
The original motivation for writing these functions was to explore alternatives or possible modifications to the revised Helgeson-Kirkham-Flowers equations applied to aqueous nonelectrolytes.
As pointed out by Schulte et al., 2001, the functional forms of the equations do not permit retrieving values of the solvation parameter (\omega) that closely represent the observed trends in both heat capacity and volume at high temperatures (above ca. 200 °C).
The examples below assume that the \omega parameter in the HKF functions is a constant (does not depend on T and P), as is appropriate for nonelectrolytes.
For charged species, the variables Cp_s_var and V_s_var can be used in the regressions.
They correspond to the solvation contribution to heat capacity or volume, respectively, in the HKF EOS, divided by the value of \omega at the reference temperature and pressure.
Because these variables are themselves a function of omega.PrTr, an iterative procedure is needed to perform the regression.
Note that variables QBorn and V_s_var are both negated, so that the value of \omega has its proper sign in the corresponding equations.
For EOSregress, an object of class “lm”.
EOSvar and EOScalc both return numeric values.
EOScoeffs returns a data frame.
Hnědkovský, L. and Wood, R. H. (1997) Apparent molar heat capacities of aqueous solutions of CH4, CO2, H2S, and NH3 at temperatures from 304 K to 704 K at a pressure of 28 MPa. J. Chem. Thermodyn. 29, 731–747. doi:10.1006/jcht.1997.0192
Schulte, M. D., Shock, E. L. and Wood, R. H. (1995) The temperature dependence of the standard-state thermodynamic properties of aqueous nonelectrolytes. Geochim. Cosmochim. Acta 65, 3919–3930. doi:10.1016/S0016-7037(01)00717-7
The vignette eos-regress has more references and examples, including an iterative method to retrieve omega.PrTr.
## Fit experimental heat capacities of CH4
## using revised Helgeson-Kirkham-Flowers equations
# Read the data from Hnedkovsky and Wood, 1997
f <- system.file("extdata/cpetc/HW97_Cp.csv", package = "CHNOSZ")
d <- read.csv(f)
# Use data for CH4
d <- d[d$species == "CH4", ]
d <- d[, -1]
# Convert J to cal and MPa to bar
d$Cp <- convert(d$Cp, "cal")
d$P <- convert(d$P, "bar")
# Specify the terms in the HKF equations
var <- c("invTTheta2", "TXBorn")
# Perform regression, with a temperature limit
EOSlm <- EOSregress(d, var, T.max = 600)
# Calculate the Cp at some temperature and pressure
EOScalc(EOSlm$coefficients, 298.15, 1)
# Get the database values of c1, c2 and omega for CH4(aq)
CH4coeffs <- EOScoeffs("CH4", "Cp")
## Make plots comparing the regressions
## with the accepted EOS parameters of CH4
opar <- par(mfrow = c(2,2))
EOSplot(d, T.max = 600)
title("Cp of CH4(aq), fit to 600 K")
legend("bottomleft", pch = 1, legend = "Hnedkovsky and Wood, 1997")
EOSplot(d, coefficients = CH4coeffs)
title("Cp from EOS parameters in database")
EOSplot(d, T.max = 600, T.plot = 600)
title("Cp fit to 600 K, plot to 600 K")
EOSplot(d, coefficients = CH4coeffs, T.plot = 600)
title("Cp from EOS parameters in database")
par(opar)
# Continuing from above, with user-defined variables
Theta <- 228 # K
invTTTheta3 <- function(T, P) (2*T) / (T-T*Theta) ^ 3
invTX <- function(T, P) 1 / T * water("XBorn", T = T, P = P)[,1]
# Print the calculated values of invTTTheta3
EOSvar("invTTTheta3", d$T, d$P)
# Use invTTTheta and invTX in a regression
var <- c("invTTTheta3", "invTX")
EOSregress(d, var)
# Give them a "label" attribute for use in the legend
attr(invTTTheta3, "label") <-
quote(phantom()%*%2 * italic(T) / (italic(T) - italic(T) * Theta) ^ 3)
attr(invTX, "label") <- quote(phantom() / italic(T * X))
# Uncomment the following to make the plot
#EOSplot(d, var)
## Model experimental volumes of CH4
## using HKF equation and an exploratory one
f <- system.file("extdata/cpetc/HWM96_V.csv", package = "CHNOSZ")
d <- read.csv(f)
# Use data for CH4 near 280 bar
d <- d[d$species == "CH4", ]
d <- d[, -1]
d <- d[abs(d$P - 28) < 0.1, ]
d$P <- convert(d$P, "bar")
# The HKF equation
varHKF <- c("invTTheta", "QBorn")
# alpha is the expansivity coefficient of water
varal <- c("invTTheta", "alpha")
opar <- par(mfrow = c(2, 2))
# For both HKF and the expansivity equation,
# we'll fit up to a temperature limit
EOSplot(d, varHKF, T.max = 663, T.plot = 625)
legend("bottomright", pch = 1, legend = "Hnedkovsky et al., 1996")
title("V of CH4(aq), HKF equation")
EOSplot(d, varal, T.max = 663, T.plot = 625)
title("V of CH4(aq), expansivity equation")
EOSplot(d, varHKF, T.max = 663)
title("V of CH4(aq), HKF equation")
EOSplot(d, varal, T.max = 663)
title("V of CH4(aq), expansivity equation")
par(opar)
# Note that the volume regression using the HKF gives
# a result for omega (coefficient on Q) that is
# not consistent with the high-T heat capacities