View source: R/plainMethods_stoiBuild_stoiCheck.r
| stoiCheck | R Documentation |
Validates the stoichiometry matrix by checking for conservation of mass (more typically conservation of moles).
stoiCheck(stoi, comp, env = globalenv(), zero = .Machine$double.eps * 2)
stoi |
Stoichiometry matrix either in evaluated
( |
comp |
Matrix defining the elemental composition of compounds.
Column names of |
env |
An environment or list supplying constants, functions, and
operators needed to evaluate expressions in |
zero |
A number close to zero. If the absolute result value of a mass balance computation is less than this, the result is set to 0 (exactly). |
A numeric matrix with the following properties:
There is one row for each process, thus the number and names of rows
are the same as in stoi.
There is one column per checked element, hence column names are
equal to the row names of comp.
A matrix element at position [i,k] represent the mass balance
for the process in row i with respect to the element in column
k. A value of zero indicates a close balance. Positive (negative)
values indicate that mass is gained (lost) in the respective process.
David Kneis david.kneis@tu-dresden.de
Use stoiCreate to create a stoichiometry matrix
from a set of reactions in common notation.
# Eq. 1 and 2 are from Soetaert et al. (1996), Geochimica et Cosmochimica
# Acta, 60 (6), 1019-1040. 'OM' is organic matter. Constants 'nc' and 'pc'
# represent the nitrogen/carbon and phosphorus/carbon ratio, respectively.
reactions <- c(
oxicDegrad= "OM + O2 -> CO2 + nc * NH3 + pc * H3PO4 + H2O",
denitrific= "OM + 0.8*HNO3 -> CO2 + nc*NH3 + 0.4*N2 + pc*H3PO4 + 1.4*H2O",
dissPhosp1= "H3PO4 <-> H + H2PO4",
dissPhosp2= "H2PO4 <-> H + HPO4"
)
# Non-evaluated stoichiometry matrix
stoi <- stoiCreate(reactions, toRight="_f", toLeft="_b")
# Parameters ('nc' and 'pc' according to Redfield ratio)
pars <- list(nc=16/106, pc=1/106)
# Elemental composition
comp <- rbind(
OM= c(C=1, N="nc", P="pc", H="2 + 3*nc + 3*pc"),
O2= c(C=0, N=0, P=0, H=0),
CO2= c(C=1, N=0, P=0, H=0),
NH3= c(C=0, N=1, P=0, H=3),
H3PO4= c(C=0, N=0, P=1, H=3),
H2PO4= c(C=0, N=0, P=1, H=2),
HPO4= c(C=0, N=0, P=1, H=1),
H2O= c(C=0, N=0, P=0, H=2),
HNO3= c(C=0, N=1, P=0, H=1),
N2= c(C=0, N=2, P=0, H=0),
H= c(C=0, N=0, P=0, H=1)
)
# We need the transposed form
comp <- t(comp)
# Mass balance check
bal <- stoiCheck(stoi, comp=comp, env=pars)
print(bal)
print(colSums(bal) == 0)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.