stoiCheck: Validation of a Stoichiometry Matrix

In dkneis/rodeo: A Code Generator for ODE-Based Models

Validation of a Stoichiometry Matrix

Description

Validates the stoichiometry matrix by checking for conservation of mass (more typically conservation of moles).

Usage

stoiCheck(stoi, comp, env = globalenv(), zero = .Machine$double.eps * 2)

Arguments

stoi	Stoichiometry matrix either in evaluated (numeric) or non-evaluated (character) form. A suitable matrix can be created with stoiCreate, for example.
comp	Matrix defining the elemental composition of compounds. Column names of comp need to match column names of stoi (but additional columns are allowed and columns can be in different order). There must be one row per element whose balance is to be checked and the elements' names must appear as row names. The elements of the matrix specify how much of an element is contained in a certain amount of a compound. Typically, these are molar ratios. If one works with mass ratios (not being a good idea), the information in stoi must be based on mass concentrations as well. The elements of comp are treated as mathematical expressions. Any variables, functions, or operators needed to evaluate those expressions must be provided by the specified environment env.
env	An environment or list supplying constants, functions, and operators needed to evaluate expressions in comp or stoi.
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).

Value

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.

Author(s)

David Kneis [e-mail]

See Also

Use stoiCreate to create a stoichiometry matrix from a set of reactions in common notation.

Examples

# 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)

dkneis/rodeo documentation built on Jan. 4, 2024, 2:18 p.m.