Chemical model of Red Sea brine

This example closely follows Section 6.3 of Bethke (2007).

A chemical analysis of the major element composition of hot hydrothermal brines of the Red Sea is shown in Table 1. In addition, the temperature is around 60C and the pH is 5.6.

Table 1: Major element composition of Red Sea brine

SpeciesConcentration (mg.kg)
Cl156030
Na92600
Ca5150
K1870
SO840
Mg764
HCO140
Fe81
Zn5.4
F5
Ba0.9
Pb0.63
Cu0.26

MOOSE input file: no precipitation or dissolution

In order to estimate the brine's oxidation state, Bethke (2007) recommends using the mineral Sphalerite instead of primary species O2(aq), and the mineral Barite instead of primary aqueous species Ba++. Bethke (2007) uses a small initial value of g of free amounts of each of these minerals to avoid dissolution when considering precipitation/dissolution later in the calculation.

The MOOSE input file contains the usual GeochemicalModelDefinition that specifies the database file to use, and in this case the basis species and equilibrium minerals. Various minerals are included in order to make the comparison to the precipitation+dissolution case (below) clearer. The TimeIndependentReactionSolver will prevent all mineral precipitation will be forbidden, except the Sphalerite and Barite which are constrained to have a fixed number of free moles, so the additional minerals do not impact the system whatsoever, except their saturation indices will be computed after the solution is found. 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+ Na+ K+ Mg++ Ca++ Cl- SO4-- HCO3- Cu+ F- Fe++ Pb++ Zn++ O2(aq) Ba++"
    equilibrium_minerals = "Sphalerite Barite Fluorite Chalcocite Bornite Chalcopyrite Pyrite Galena Covellite"
    piecewise_linear_interpolation = true # for comparison with GWB
  []
[]
(modules/geochemistry/test/tests/equilibrium_models/red_sea_no_precip.i)

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

  • The swaps are defined so that the minerals Sphalerite and Barite can be provided with a small number of free moles.

  • The pH is fixed using the activity of H.

  • There is 1kg of solvent water.

  • 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 = "O2(aq) Ba++"
  swap_into_basis = "Sphalerite Barite"
  prevent_precipitation = "Sphalerite Barite Fluorite Chalcocite Bornite Chalcopyrite Pyrite Galena Covellite"
  charge_balance_species = "Cl-"
  constraint_species = "H2O              H+            Na+              K+               Mg++             Ca++             Cl-              SO4--            HCO3-            Cu+              F-               Fe++            Pb++              Zn++             Sphalerite   Barite"
  constraint_value = "  1.0              -5.6          5.42             0.0643           0.0423           0.173            5.89             0.0118           0.00309          5.50E-06         0.000354         0.00195          4.09E-06         0.000111         1E-11        0.5E-11"
  constraint_meaning = "kg_solvent_water log10activity bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition free_mineral free_mineral"
  constraint_unit = "   kg               dimensionless moles            moles            moles            moles            moles            moles            moles            moles            moles            moles            moles            moles            moles        moles"
  ramp_max_ionic_strength_initial = 0 # not needed in this simple example
  temperature = 60
  stoichiometric_ionic_str_using_Cl_only = true # for comparison with GWB
  mol_cutoff = 1E-7
  abs_tol = 1E-12
[]
(modules/geochemistry/test/tests/equilibrium_models/red_sea_no_precip.i)

Geochemists Workbench input file: no precipitation or dissolution

The analogous Geochemists Workbench input file for this case is

# React script that is equivalent to red_sea_no_precip.i
data = thermo.tdat verify
conductivity = conductivity-USGS.dat
temperature = 60 C
H2O          = 1 free kg
Cl-          = 5.89 mol
balance on Cl-
Na+          = 5.42 mol
Ca++         = 0.173 mol
K+           = 0.0643 mol
SO4--        = 0.0118 mol
Mg++         = 0.0423 mol
HCO3-        = 0.00309 mol
Fe++         = 0.00195 mol
Zn++         = 0.000111 mol
F-           = 0.000354 mol
swap Barite for Ba++
Barite       = 0.5e-11 free mol
Pb++         = 4.09e-06 mol
Cu+          = 5.50e-6 mol
pH           = 5.6
swap Sphalerite for O2(aq)
Sphalerite   = 1e-11 free mol
suppress ALL
unsuppress  Barite Sphalerite
printout  species = long
epsilon = 1e-12
(modules/geochemistry/test/tests/equilibrium_models/red_sea_no_precip.rea)

Results: no precipitation or dissolution

The geochemistry results match those from Geochemists Workbench exactly.

Error and charge-neutrality error

The geochemistry simulation reports an error of 3.862e-14mol, and that the charge of the solution is -1.732e-17mol.

Solution mass

The solution mass is 1.345kg.

Ionic strength and water activity

The ionic strength is greater than 3 (which is the default upper-bound for applicability of the Debye-Huckel theory) and the water activity is 0.899.

pH, pe and Eh

The pH is 5.6 (as specified) the pe is -1.67, and Eh = -0.111V.

Aqueous species distribution

The geochemistry output matches Bethke (2007) (and the GWB software) who writes that the molalities of the most abundant species results are as shown in Table 2.

Table 2: Calculated molalities, activity coefficients and activities of the most abundant species in Red Sea water

SpeciesMolality (mol.kg)Activity coeffloga
Cl5.1820.61250.5017
Na4.860.70360.5341
NaCl0.5511.0-0.2587
CaCl0.12770.7036-1.047
K0.6125-1.438
Ca0.1941-2.064
MgCl0.7036-1.756
Mg0.2895-2.310
NaSO0.7036-2.325
SO0.0985-3.579
CO(aq)1.0-2.743

Mineral saturation indices

The saturation indices of the equilibrium solution in Table 2 are greater than 0 for a number of minerals: bornite, chalcopyrite, chalcocite, pyrite, fluorite, galena and covellite. Both GWB and geochemistry produce identical results

Mineral precipitation and dissolution

The results are slightly different when minerals are allowed to precipitate and the Sphalerite and Barite are allowed to dissolve.

MOOSE input file: allowing precipitation or dissolution

The MOOSE input file is almost identical to the one above

  • the prevent_precipitation option is removed

  • a bulk number of moles are provided to Sphalerite and Barite so that they can dissolve. If a free number of moles had been specified then geochemistry assumes that moles of the mineral are added or removed by an external agent in order to keep the free number as specified. Hence, specifying free mole numbers prevents any dissolution. The bulk number of moles is set to the amount predicted by red_sea_no_precip.i

[TimeIndependentReactionSolver]
  model_definition = definition
  swap_out_of_basis = "O2(aq) Ba++"
  swap_into_basis = "Sphalerite Barite"
  charge_balance_species = "Cl-"
  constraint_species = "H2O              H+            Na+              K+               Mg++             Ca++             Cl-              SO4--            HCO3-            Cu+              F-               Fe++            Pb++              Zn++             Sphalerite       Barite"
  constraint_value = "  1.0              -5.6          5.42             0.0643           0.0423           0.173            5.89             0.0118           0.00309          5.50E-06         0.000354         0.00195          4.09E-06         0.000111         5.87E-8          9.772E-6"
  constraint_meaning = "kg_solvent_water log10activity bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition bulk_composition"
  constraint_unit = "   kg               dimensionless moles            moles            moles            moles            moles            moles            moles            moles            moles            moles            moles            moles            moles            moles"
  ramp_max_ionic_strength_initial = 0
  temperature = 60
  stoichiometric_ionic_str_using_Cl_only = true # for comparison with GWB
  mol_cutoff = 1E-7
  abs_tol = 1E-12
[]
(modules/geochemistry/test/tests/equilibrium_models/red_sea_precip.i)

Geochemists Workbench input file: allowing precipitation or dissolution

The analogous Geochemists Workbench input file for this case is

# React script that is equivalent to red_sea_precip.i
data = thermo.tdat verify
conductivity = conductivity-USGS.dat
temperature = 60 C
H2O          = 1 free kg
Cl-          = 5.89 mol
balance on Cl-
Na+          = 5.42 mol
Ca++         = 0.173 mol
K+           = 0.0643 mol
SO4--        = 0.0118 mol
Mg++         = 0.0423 mol
HCO3-        = 0.00309 mol
Fe++         = 0.00195 mol
Zn++         = 0.000111 mol
F-           = 0.000354 mol
swap Barite for Ba++
Barite       = 0.5e-11 free mol
Pb++         = 4.09e-06 mol
Cu+          = 5.50e-6 mol
pH           = 5.6
swap Sphalerite for O2(aq)
Sphalerite   = 1e-11 free mol
printout  species = long
epsilon = 1e-14
(modules/geochemistry/test/tests/equilibrium_models/red_sea_precip.rea)
commentnote

In this input file, Geochemists Workbench differs subtly from the geochemistry module. A free number of mineral moles is specified in Geochemists Workbench, which then proceeds to solve the system without precipitation or dissolution, then fixes the bulk mole number of Sphalerite and Barite, then continues the solve allowing precipitation and dissolution. In contrast in geochemistry, if a free number of mineral moles is specified, it is assumed that the condition is meant to hold throughout the entire solve process. Hence, in geochemistry, a bulk number of moles is specified (alternately, a TimeDependentReactionSolver could be used to replicate GWB's procedure).

Results: allowing precipitation and dissolution

Allowing the minerals to precipitate, both codes and Bethke (2007) predicts that the Sphalerite dissolves and that 3 minerals precipitate in small quantities, as shown in Table 3

Table 3: Calculated final mass of each precipitate in the stable phase assemblage for the Red Sea brine.

MineralMass (g)
Fluorite
Chalcocite
Barite

References

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