Chemical model of seawater

MOOSE input file: no precipitation

The MOOSE input file contains the usual GeochemicalModelDefinition that specifies the database file to use, and in this case the basis species, equilibrium minerals and equilibrium gases. The flag piecewise_linear_interpolation = true in order to compare with the Geochemists Workbench result

[UserObjects]
  [definition]
    type = GeochemicalModelDefinition
    database_file = "../../../database/moose_geochemdb.json"
    basis_species = "H2O H+ Cl- Na+ SO4-- Mg++ Ca++ K+ HCO3- SiO2(aq) O2(aq)"
    equilibrium_minerals = "Antigorite Tremolite Talc Chrysotile Sepiolite Anthophyllite Dolomite Dolomite-ord Huntite Dolomite-dis Magnesite Calcite Aragonite Quartz"
    equilibrium_gases = "O2(g) CO2(g)"
    piecewise_linear_interpolation = true # for comparison with GWB
  []
[]

To instruct MOOSE to find the equilibrium configuration, a TimeIndependentReactionSolver is used:

• The swaps are defined.

• The fugacity of the gases is fixed as defined above.

• The bulk mole number of the aqueous species is also fixed appropriately. The numbers are different than the concentration in mg.kg given in the above table, and may be worked out using the TDS.

• The prevent_precipitation input prevents any minerals from precipitating when finding the equilibrium configuration, even if their saturation indices are positive.

• The other flags enable an accurate comparison with the Geochemists Workbench software.

[TimeIndependentReactionSolver]
  model_definition = definition
  swap_out_of_basis = "H+     O2(aq)"
  swap_into_basis = "  CO2(g) O2(g)"
  charge_balance_species = "Cl-" # this means the bulk moles of Cl- will not be exactly as set below
  constraint_species = "H2O              CO2(g)    O2(g)         Cl-              Na+              SO4--            Mg++             Ca++             K+               HCO3-            SiO2(aq)"
  constraint_value = "  1.0              0.0003162 0.2           0.566            0.485            0.0292           0.055            0.0106           0.0106           0.00241          0.000103"
  constraint_meaning = "kg_solvent_water fugacity  fugacity      bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition"
  constraint_unit = "   kg          dimensionless  dimensionless moles            moles            moles            moles            moles            moles            moles            moles"
  prevent_precipitation = "Antigorite Tremolite Talc Chrysotile Sepiolite Anthophyllite Dolomite Dolomite-ord Huntite Dolomite-dis Magnesite Calcite Aragonite Quartz"
  ramp_max_ionic_strength_initial = 0 # not needed in this simple problem
  stoichiometric_ionic_str_using_Cl_only = true # for comparison with GWB
  mol_cutoff = 1E-5
  abs_tol = 1E-15
  precision = 7
[]

MOOSE input file: with precipitation

The MOOSE input file is very similar to the no-precipitation case. There are two differences:

• The prevent_precipitation input is changed.

• The swaps are different, as is the initial condition for MgCO. This is discussed in Bethke (2007).

[TimeIndependentReactionSolver]
  model_definition = definition
  swap_out_of_basis = "H+   "
  swap_into_basis = "  MgCO3"
  charge_balance_species = "Cl-" # this means the bulk moles of Cl- will not be exactly as set below
  constraint_species = "H2O              MgCO3            O2(aq)             Cl-              Na+              SO4--            Mg++             Ca++             K+               HCO3-            SiO2(aq)"
# to obtain the constraint on MgCO3: (1) run seawater_no_precip.i to obtain the free molality of MgCO3; (2) then running seawater_no_precip.i with MgCO3 in the basis (in place of H+) with that free molality, to obtain the bulk mole number
  constraint_value = "  1.0              0.0001959        0.2151E-3          0.566            0.485            0.0292           0.055            0.0106           0.0106           0.00241          0.000103"
  constraint_meaning = "kg_solvent_water bulk_composition free_concentration bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition"
  constraint_unit = "   kg               moles            molal              moles            moles            moles            moles            moles            moles            moles            moles"
  prevent_precipitation = "Dolomite-dis Dolomite-ord"
  ramp_max_ionic_strength_initial = 0 # not needed in this simple problem
  stoichiometric_ionic_str_using_Cl_only = true # for comparison with GWB
  mol_cutoff = 1E-5
  abs_tol = 1E-15
[]

[UserObjects]
  [definition]
    type = GeochemicalModelDefinition
    database_file = "../../../database/moose_geochemdb.json"
    basis_species = "H2O H+ Cl- Na+ SO4-- Mg++ Ca++ K+ HCO3- SiO2(aq) O2(aq)"
    equilibrium_minerals = "Antigorite Tremolite Talc Chrysotile Sepiolite Anthophyllite Dolomite Dolomite-ord Huntite Dolomite-dis Magnesite Calcite Aragonite Quartz"
    piecewise_linear_interpolation = true # for comparison with GWB
  []
[]

Geochemists Workbench input files

The Geochemists Workbench input file for the precipitation case is:

# React script that is analogous to the seawater_precip.i MOOSE input file
data = thermo.tdat verify
conductivity = conductivity-USGS.dat
temperature = 25 C
H2O          = 1 free kg
Cl-          = 0.566 mol
balance on Cl-
Na+          = 0.485 mol
SO4--        = 0.0292 mol
Mg++         = 0.055 mol
Ca++         = 0.0106 mol
K+           = 0.0106 mol
HCO3-        = 0.00241 mol
SiO2(aq)     = 0.000103 mol
O2(aq)       = 0.0002151 free molal
swap MgCO3 for H+
MgCO3        = 0.0001959 mol
printout  species = long
epsilon = 1e-15
go

The non-precipitation case is similar (see the seawater_no_precip.rea file).

Results

The geochemistry results mirror those from Geochemists Workbench exactly.

References