poet/bench/barite/README.org

Description of barite benchmark

Quick start

mpirun -np 4 ./poet barite.R barite_results
mpirun -np 4 ./poet --interp  barite_interp_eval.R barite_results

List of Files

• barite.R: POET input script for a 20x20 simulation grid
• barite_interp_eval.R: POET input script for a 400x200 simulation grid
• db_barite.dat: PHREEQC database containing the kinetic expressions for barite and celestite, stripped down from phreeqc.dat
• barite.pqi: PHREEQC input script defining the chemical system

Chemical system

The benchmark depicts an isotherm porous system at 25 °C where pure water is initially at equilibrium with celestite (strontium sulfate; brute formula: SrSO₄). A solution containing only dissolved Ba²⁺ and Cl^- diffuses into the system causing celestite dissolution. The increased concentration of dissolved sulfate SO₄^2- induces precipitation of barite (barium sulfate; brute formula: BaSO₄^2-). The overall reaction can be written:

Ba²⁺ + celestite → barite + Sr²⁺

Both celestite dissolution and barite precipitation are calculated using a kinetics rate law based on transition state theory:

rate = -S_m k_barite (1-SR_m)

where the reaction rate has units mol/s, S_m (m²) is the reactive surface area, k (mol/m²/s) is the kinetic coefficient, and SR is the saturation ratio, i.e., the ratio of the ion activity product of the reacting species and the solubility constant, calculated internally by PHREEQC from the speciated solution.

For barite, the reaction rate is computed as sum of two mechanisms, r_/acid/ and r_/neutral/:

rate_barite = S_barite (k_/acid/ + k_/neutral/) * (1 - SR_barite)

where:

k_/acid/ = 10^-6.9 e^{-30800 / R} ⋅ act(H⁺)^0.22

k_/neutral/ = 10^-7.9 e^{-30800 / R}

R (8.314462 J K^-1 mol^-1) is the gas constant.

For celestite the kinetic law considers only the acidic mechanism and reads:

rate_celestite = S_celestite 10^-5.66 e^{-23800 / R} ⋅ act(H⁺)^0.109 ⋅ (1 - SR_celestite)

The kinetic rates as implemented in the db_barite.dat file accepts one parameter which represents reactive surface area in m². For the benchmarks the surface areas are set to

• S_barite: 50 m²
• S_celestite: 10 m²

A starting seed for barite is given at 0.001 mol in each domain element.

TODO Nucleation

Geochemical benchmark inspired by Tranter et al. (2021) without nucleation.

POET simulations

Currently these benchmarks are pure diffusion simulations. There are 7 transported species: H, O, Charge, Ba, Cl, S(6), Sr.

barite.R

• Grid discretization: square domain of 1 ⋅ 1 m² discretized in 20x20 cells
• Boundary conditions: E, S and W sides of the domain are closed; the N boundary has a fixed concentration (Dirichlet) of 0.1 molal BaCl₂
• Diffusion coefficients: isotropic homogeneous α = 1E-06
• Time steps & iterations: 20 iteration with Δ t = 250 s
• DHT parameters:

H	10
O	10
Charge	3
Ba	5
Cl	5
S(6)	5
Sr	5

barite_interp_eval.R

• Grid discretization: rectangular domain of 40 ⋅ 20 m² discretized in 400 ⋅ 200 cells
• Boundary conditions: all boundaries are closed. The center of the domain at indeces (200, 100) has fixed concentration of 0.1 molal of BaCl₂
• Diffusion coefficients: isotropic homogeneous α = 1E-06
• Time steps & iterations: 200 iterations with Δ t = 250 s
• PHT parameters:

H	10
O	10
Charge	3
Ba	5
Cl	5
S(6)	5
Sr	5

References

• Tranter, Wetzel, De Lucia and Kühn (2021): Reactive transport model of kinetically controlled celestite to barite replacement, Advances in Geosciences, 1, 1–9, https://doi.org/10.5194/adgeo-1-1-2021