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.
Species | Concentration (mg.kg) |
---|---|
Cl | 156030 |
Na | 92600 |
Ca | 5150 |
K | 1870 |
SO | 840 |
Mg | 764 |
HCO | 140 |
Fe | 81 |
Zn | 5.4 |
F | 5 |
Ba | 0.9 |
Pb | 0.63 |
Cu | 0.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.
Species | Molality (mol.kg) | Activity coeff | loga |
---|---|---|---|
Cl | 5.182 | 0.6125 | 0.5017 |
Na | 4.86 | 0.7036 | 0.5341 |
NaCl | 0.551 | 1.0 | -0.2587 |
CaCl | 0.1277 | 0.7036 | -1.047 |
K | 0.6125 | -1.438 | |
Ca | 0.1941 | -2.064 | |
MgCl | 0.7036 | -1.756 | |
Mg | 0.2895 | -2.310 | |
NaSO | 0.7036 | -2.325 | |
SO | 0.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 removeda 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 byred_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)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
Mineral | Mass (g) |
---|---|
Fluorite | |
Chalcocite | |
Barite |
References
- Craig M. Bethke.
Geochemical and Biogeochemical Reaction Modeling.
Cambridge University Press, 2 edition, 2007.
doi:10.1017/CBO9780511619670.[BibTeX]