Dissolution of pyrite, with and without fixed oxygen fugacity

Section 14.2 of Bethke (2007) and involves studying the dissolution of pyrite, with and without fixed oxygen fugacity.

The initial composition of the water is shown in Table 1. In addition, in the initial configuration:

  • hematite is used in place of Fe and its free concentration is mol (mg)

  • O(g) is used in place of O(aq), and its fugacity is set to ;

  • the pH is 6.5.

Table 1: Initial bulk composition

SpeciesConcentration (molal)Approx conc (mg.kg)
Ca4
Mg1
Na2
HCO18
SO3
Cl (charge-balance)5
Hematite (in place of Fe) mol (free)1 mg (free)

Two cases are studied:

  • a case in which oxygen cannot or leave the system after initial equilibration

  • a case in which the oxygen fugacity is fixed at 0.2 throughout the simulation

The MOOSE input files for both these cases are very similar. A GeochemicalModelDefinition defines the basis species, the equilibrium minerals and equilibrium gases:

[UserObjects]
  [definition]
    type = GeochemicalModelDefinition
    database_file = "../../../database/moose_geochemdb.json"
    basis_species = "H2O H+ Fe++ Ca++ Mg++ Na+ HCO3- SO4-- Cl- O2(aq)"
    equilibrium_minerals = "Hematite Pyrite"
    equilibrium_gases = "O2(g)"
    piecewise_linear_interpolation = true # for comparison with GWB
  []
[]
(modules/geochemistry/test/tests/time_dependent_reactions/dissolution_pyrite_1.i)

An Executioner defines the time-stepping and end-time:

[Executioner]
  type = Transient
  dt = 1
  end_time = 10
[]
(modules/geochemistry/test/tests/time_dependent_reactions/dissolution_pyrite_1.i)

The figures below were generated using a time-step size of 0.05, while the regression-test files use a larger time-step size for efficiency.

A set of Postprocessors define the desired outputs using the AuxVariables automatically added by the TimeDependentReactionSolver:

[Postprocessors]
  [mg_Hematite]
    type = PointValue
    point = '0 0 0'
    variable = 'free_mg_Hematite'
  []
  [mg_Pyrite]
    type = PointValue
    point = '0 0 0'
    variable = 'free_mg_Pyrite'
  []
  [pH]
    type = PointValue
    point = '0 0 0'
    variable = 'pH'
  []
  [molal_CO2aq]
    type = PointValue
    point = '0 0 0'
    variable = 'molal_CO2(aq)'
  []
  [molal_HCO3-]
    type = PointValue
    point = '0 0 0'
    variable = 'molal_HCO3-'
  []
  [molal_SO4--]
    type = PointValue
    point = '0 0 0'
    variable = 'molal_SO4--'
  []
  [molal_Fe++]
    type = PointValue
    point = '0 0 0'
    variable = 'molal_Fe++'
  []
  [molal_O2aq]
    type = PointValue
    point = '0 0 0'
    variable = 'molal_O2(aq)'
  []
[]
(modules/geochemistry/test/tests/time_dependent_reactions/dissolution_pyrite_1.i)

Case 1: when oxygen cannot enter or leave the system

10mg of pyrite is slowly added to the system, without allowing any oxygen to enter or exit the system.

MOOSE input file

The TimeDependentReactionSolver involves specifying:

  • the swaps

  • the initial free mole number of Hematite, the pH via the H activity, the bulk mole composition of the aqueous species, and the O2(g) fugacity

  • the system is closed at (the default)

  • that at the activity constraints are removed from H (so the pH can vary during the simulation because an external agent is no longer fixing it) and O2(g) (meaning the external agent controlling the O2(g) fugacity is no longer adding/removing O2(g) from the aqueous solution)

  • that 1mg/second of pyrite is added to the system

[TimeDependentReactionSolver]
  model_definition = definition
  geochemistry_reactor_name = reactor
  swap_out_of_basis = "O2(aq) Fe++"
  swap_into_basis = "O2(g) Hematite"
  charge_balance_species = "Cl-"
  constraint_species = "H2O              Hematite     H+            Ca++             Mg++             Na+              HCO3-            SO4--            Cl-              O2(g)"
  constraint_value = "  1.0              1            -6.5          4                1                2                18               3                5                0.2"
  constraint_meaning = "kg_solvent_water free_mineral log10activity bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition fugacity"
  constraint_unit = "   kg               mg           dimensionless mg               mg               mg               mg               mg               mg               dimensionless"
  remove_fixed_activity_name = "H+ O2(g)"
  remove_fixed_activity_time = '0  0'
  source_species_names = "Pyrite"
  source_species_rates = 8.336E-6 # = 1mg(pyrite)/second, 1mg(pyrite) = 8.34E-6
  ramp_max_ionic_strength_initial = 0 # not needed in this simple problem
  stoichiometric_ionic_str_using_Cl_only = true # for comparison with GWB
  abs_tol = 1E-13
  execute_console_output_on = '' # only CSV output is required
[]
(modules/geochemistry/test/tests/time_dependent_reactions/dissolution_pyrite_1.i)

GWB input file

The equivalent Geochemists Workbench input file is

# React script that is equivalent to dissolution_pyrite_1.i
Ca++ = 9.98e-5 mol
Mg++ = 4.114e-5 mol
Na+ = 8.7e-5 mol
HCO3- = 2.95e-4 mol
SO4-- = 3.123e-5 mol
Cl- = 1.4103e-5 mol
balance on Cl-
H+ = 6.5 pH
swap O2(g) for O2(aq)
f O2(g) = 0.2
swap Hematite for Fe++
Hematite = 6.26e-6 free mol
epsilon = 1e-13
react 8.336e-5 mol Pyrite
go
(modules/geochemistry/test/tests/time_dependent_reactions/dissolution_pyrite_1.rea)

Results

The following graphs demonstrate agreement between the geochemistry module and the GWB software, which may also be compared with Bethke (2007) Figure 14.2 and 14.3.

Figure 1: Precipitated volumes as pyrite is reacted in a system that is closed to oxygen

Figure 2: Solution pH as pyrite is reacted in a system that is closed to oxygen

Figure 3: Species concentrations as pyrite is reacted in a system that is closed to oxygen

Case 2: fixed oxygen fugacity

1000mg of pyrite is slowly added to the system, with fixed fugacity of oxygen (so that oxygen can enter or leave the system)

MOOSE input file

The TimeDependentReactionSolver is almost identical to Case 1. The differences are:

  • that at the activity constraints are removed from H, but not O2(g)

  • that 1000mg/second of pyrite is added to the system

[TimeDependentReactionSolver]
  model_definition = definition
  geochemistry_reactor_name = reactor
  swap_out_of_basis = "O2(aq) Fe++"
  swap_into_basis = "O2(g) Hematite"
  charge_balance_species = "Cl-"
  constraint_species = "H2O              Hematite     H+            Ca++             Mg++             Na+              HCO3-            SO4--            Cl-              O2(g)"
  constraint_value = "  1.0              1            -6.5          4                1                2                18               3                5                0.2"
  constraint_meaning = "kg_solvent_water free_mineral log10activity bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition fugacity"
  constraint_unit = "   kg               mg           dimensionless mg               mg               mg               mg               mg               mg               dimensionless"
  remove_fixed_activity_name = "H+"
  remove_fixed_activity_time = '0 '
  source_species_names = "Pyrite"
  source_species_rates = 8.336E-4 # = 0.1g(pyrite)/second over 10 seconds, 1g(pyrite) = 8.34E-3
  ramp_max_ionic_strength_initial = 0 # not needed in this simple problem
  stoichiometric_ionic_str_using_Cl_only = true # for comparison with GWB
  abs_tol = 1E-13
  execute_console_output_on = '' # only CSV output is required
[]
(modules/geochemistry/test/tests/time_dependent_reactions/dissolution_pyrite_2.i)

Results

The following graphs demonstrate agreement between the geochemistry module and the GWB software, which may also be compared with Bethke (2007) Figure 14.4.

Figure 4: Precipitated volumes as pyrite is reacted in a system that is open to oxygen

Figure 5: Solution pH as pyrite is reacted in a system that is open to oxygen

The two cases are clearly dramatically different. In Case 1, about 8mg of pyrite dissolves, precipitating hematite, before hematite suddenly dissolves and no further pyrite dissolution occurs. In Case 2, pyrite dissolution and hematite precipitation continues indefinitely in a rather acidic solution.

References

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