Biogeochemistry with a sulfate reducer

This follows the example in Section 18.5 of Bethke (2007). Suppose that a sulfate-reducing microbe acts in the presence of acetate (CHCOO) to catalyse the reaction (1) Bethke details the initial composition of the aqueous solution for this problem. Most important for this example are the initial bulk composition of the acetate, which is moles, the initial bulk composition of sulfate, which is large, and the initial mass of the sulfate reducer, which is 0.1mg.

In this example, microbe mortality is unimportant. For each mole turnover of this reaction, Bethke estimates the microbe mass increases by g.

The rate of Eq. (1) is assumed to be of the form (2) where

  • (units: kg) is the mass of solvent water;

  • mol.mg.smol.g.day is the intrinsic rate constant, where the mass (grams) in the denominator is the mass of microbe, and the direction is forward (dissolution of acetate only, not the reverse);

  • (units: g.kg) is the microbe concentration, where the mass (kg) in the denominator is the mass of solvent water;

  • (units: mol.kg) is the molality of CHCOO (this is the free concentration, not the bulk composition);

  • mol.kg is the half-saturation constant of the acetate;

  • is the activity product of Eq. (1);

  • is the equilibrium constant of Eq. (1);

  • J.mol is an estimate of the energy captured by the microbe for each mole turnover of Eq. (1);

  • J.K.mol is the gas constant;

  • (units: K) is the temperature.

As discussed on the main biogeochemistry page, there are two ways of modelling this scenario:

  1. the microbe is treated as a new basis species, and the acetate (CHCOO) as a kinetic species that is not an equilibrium with the aqueous solution;

  2. the microbe is treated as a kinetic species, and all other species (acetate, etc) are at equilibrium

Input files and results for both of these methods are found below.

Method 1

Here, the microbe is considered as a primary species. The aqueous solution as described by Bethke contains a number of ionic species: the important thing here is the existence of Biomass1 in the basis_species list:

[UserObjects]
  [definition]
    type = GeochemicalModelDefinition
    database_file = "../../../database/moose_geochemdb.json"
    basis_species = "H2O Na+ Ca++ Fe++ Cl- SO4-- HCO3- O2(aq) H+ Biomass1"
    equilibrium_minerals = "Mackinawite" # other minerals make marginal difference
    kinetic_redox = "CH3COO-"
    kinetic_rate_descriptions = "rate_sulfate_reducer"
    piecewise_linear_interpolation = true # comparison with GWB
  []
[]
(modules/geochemistry/test/tests/kinetics/bio_sulfate_1.i)

The Biomass1 species has a molar mass of g.mol.

The acetate is considered a kinetic species. Eq. (2) is implemented using the following GeochemistryKineticRate UserObject:

[UserObjects]
  [rate_sulfate_reducer]
    type = GeochemistryKineticRate
    kinetic_species_name = "CH3COO-"
    intrinsic_rate_constant = 0.0864 # 1E-9 mol/mg/s = 0.0864 mol/g/day
    multiply_by_mass = false
    kinetic_molal_index = 1.0
    kinetic_monod_index = 1.0
    kinetic_half_saturation = 70E-6
    promoting_species_names = 'H2O Biomass1'
    promoting_indices = '1 1'
    direction = dissolution
    non_kinetic_biological_catalyst = Biomass1
    non_kinetic_biological_efficiency = 4.3
    energy_captured = 45E3
    theta = 0.2
    eta = 1
  []
[]
(modules/geochemistry/test/tests/kinetics/bio_sulfate_1.i)

Most of the values are self-explanatory. The parameter direction = dissolution because the microbe acts to convert acetate into bicarbonate, not the reverse. Notice the non_kinetic_biological_catalyst. Because one mole of Eq. (1) produces of microbe,


non_kinetic_biological_efficiency = 4.3

Consider the initial composition. Above, it is specified that the initial bulk composition of acetate is moles. It is easy to use the geochemistry module to show that the molality of acetate is 0.0008643mol.kg (for instance, remove the sulfate_reducer from bio_sulfate_2.i mentioned below and use a TimeIndependentReactionSolver). This is the molality that enters into Eq. (2). There are two possibilities:

  • initialize the kinetic acetate to 0.0008643moles. Then the rate is likely to be reasonably accurate, but there won't be as much acetate to react compared with the case when moles are used, so the final biomass will be smaller.

  • initialize the kinetic acetate to moles. Then the rate will be slightly higher than it should be, but the full amount of acetate will react, so the final biomass will be reasonably correct. This approach is used to generate the graphs below.

This example illustrate the inaccuracies brought about by using the simple approach of treating the biomass as an equilibrium species that is generated by a kinetic reaction.

Method 2

Here, the microbe is considered as a kinetic species that is generated via Eq. (1). The other species (acetate, sulfate, etc) are considered primary species in equilibrium in the aqueous solution. The standard MOOSE geochemical database contains the required information:


    "sulfate_reducer": {
      "species": {
        "CH3COO-": -1.0,
        "SO4--": -1.0,
        "HCO3-": 2.0,
        "HS-": 1.0
      },
      "molar volume": 1,
      "molecular weight": 1E3,
	"logk": [8.502, 8.404, 8.451, 8.657, 9.023, 9.457, 9.917, 10.31]
    },

If this entry did not appear in the database, it could be manually entered, but the equilibrium constants would have to be derived from those already in the database (or from experimental measurements). This is easily done by using the GeochemicalModelInterrogator:

[GeochemicalModelInterrogator]
  model_definition = definition
  swap_out_of_basis = 'O2(aq) HS-'
  swap_into_basis = 'HS- CH4(aq)'
  temperature = 300.0
[]
(modules/geochemistry/test/tests/kinetics/bio_sulfate_0.i)

and the appropriate model definition:

[UserObjects]
  [definition]
    type = GeochemicalModelDefinition
    database_file = "../../../database/moose_geochemdb.json"
    basis_species = "H2O HCO3- SO4-- H+ O2(aq)"
    piecewise_linear_interpolation = true
  []
[]
(modules/geochemistry/test/tests/kinetics/bio_sulfate_0.i)

After the entry exists in the database, the model description starts with defining the species present:

[UserObjects]
  [definition]
    type = GeochemicalModelDefinition
    database_file = "../../../database/moose_geochemdb.json"
    basis_species = "H2O Na+ Ca++ Fe++ Cl- SO4-- HCO3- O2(aq) H+ CH3COO-"
    equilibrium_minerals = "Mackinawite" # other minerals make marginal difference
    kinetic_minerals = "sulfate_reducer"
    kinetic_rate_descriptions = "rate_sulfate_reducer"
    piecewise_linear_interpolation = true # comparison with GWB
  []
[]
(modules/geochemistry/test/tests/kinetics/bio_sulfate_2.i)

The sulfate_reducer species has a molar mass of g.mol in the database file. It is unimportant that this is definitely incorrect: it is only important to remember this value when entering numerical values for the kinetic rate, below.

Eq. (2) is implemented using the following GeochemistryKineticRate UserObject:

[UserObjects]
  [rate_sulfate_reducer]
    type = GeochemistryKineticRate
    kinetic_species_name = "sulfate_reducer"
    intrinsic_rate_constant = 0.0864 # 1E-9 mol/mg/s = 0.0864 mol/g/day
    multiply_by_mass = true
    promoting_species_names = 'CH3COO-'
    promoting_indices = 1
    promoting_monod_indices = 1
    promoting_half_saturation = 70E-6
    direction = both
    kinetic_biological_efficiency = 4.3E-3
    energy_captured = 45E3
    theta = 0.2
    eta = 1
  []
[]
(modules/geochemistry/test/tests/kinetics/bio_sulfate_2.i)

Most of the values are self-explanatory. Note:

  • multiply_by_mass = true means the rate is multiplied by the mass (in grams) of the sulfate_reducer, as desired

  • there is no need to use H2O in the promoting_species because is just the mass of the sulfate_reducer in grams

  • kinetic_biological_efficiency = 4.3E-3 because this is the number of moles of sulfate_reducer that is created for one mole Eq. (1) turnover: recall that its molar mass is g.mol

Results

As expected, the two methods differ marginally: see Figure 1. Method 1 is slightly inaccurate because the full amount of acetate is used in the rate, instead of just the free molality as in method 2. Method 2 is initialized with the same bulk composition ( moles) but some of this is bound into secondary species, which slowly release the acetate into solution in order to react via Eq. (1).

Figure 1: Results of the biologically-catalysed sulfate reduction

References

  1. Craig M. Bethke. Geochemical and Biogeochemical Reaction Modeling. Cambridge University Press, 2 edition, 2007. doi:10.1017/CBO9780511619670.[BibTeX]