Effect of calcite dissolution on exchangeable cations#

Written by Svetlana Kyas (ETH Zurich) on Mar 10th, 2022

In this tutorial, we study the effect of calcite dissolution on determination of CaX2 using 1-molal NH4-acetate. We will vary the acidity of the solution from 9.3 to 7.0 by adding NH4+ to the porewater.

As all the tutorials above, the numerical experiment start from importing necessary python-packages, initializing database, and defining aqueous and exchange phases to form chemical system.

from reaktoro import *
import numpy as np
import pandas as pd

# Load Phreeqc database
db = PhreeqcDatabase("phreeqc.dat")

# Define an aqueous phase
solution = AqueousPhase(speciate("H O C Na N Cl Ca Mg"))
solution.setActivityModel(chain(
    ActivityModelHKF(),
    ActivityModelDrummond("CO2")
))

Since we are mainly concerned with the extraction of Ca+2, we only include two ion exchange species into ion exchange phase. Below, we also define a single mineral phase represented by calcite.

# Define an ion exchange phase
exchange = IonExchangePhase("NaX CaX2")
exchange.setActivityModel(ActivityModelIonExchangeGainesThomas())

# Define an ion exchange phase
mineral = MineralPhase("Calcite")

After all the phases are prepared, the chemical system as well all other building blocks of equilibrium simulation, i.e., equilibrium solver and aqueous properties, can be initialized.

# Create the chemical system
system = ChemicalSystem(db, solution, exchange, mineral)

# Define equilibrium solver and equilibrate given initial state with input conditions
solver = EquilibriumSolver(system)

# Define ion exchange properties and aqueous properties
aqprops = AqueousProps(system)

Bellow, we define the range of NH4+ amount as well as the array with NaX values in the initial chemical state. We also collect the parameters and arrays of the data we are trying to study.

# Sampling arrays of NH4 ions' amounts and NaX exchange species
num_ph = 31
num_exchangers = 4
mols_NH4 = np.linspace(0, 60.0, num=num_ph)
mols_NaX = 1e-3 * np.array([0.0, 2.5, 12.5, 25])

# Output dataframe
data = pd.DataFrame(columns=["amount_NaX", "amount_NH4", "pH", "I",
                             "amount_Ca", "amount_CaX2", "delta_Calcite"])

Next, we perform a sequence of equilibrium calculations by varying amounts of NH4-acetate and NaX. At the end of each calculation, we’ll extract the following properties from the computed chemical state:

  • pH,

  • ionic strength I,

  • amount of Ca+2,

  • amount of CaX2,

  • amount of Calcite.

During each equilibrium calculation, the following two major equilbrium calculations will be defining the solubility of calcite and to formation of CaX2: $$

(5)#\[\begin{alignat}{4} {\rm CaCO{_3}} + {\rm H_2O} & \rightleftharpoons {\rm Ca}^{+2} + {\rm HCO}_3^- + {\rm OH}^- \\ \tfrac{1}{2} {\rm Ca^{+2}} + {\rm NaX} & \rightleftharpoons {\rm Na}^+ + \tfrac{1}{2} {\rm CaX}_2\\ \end{alignat}\]

$$

for mol_NaX in mols_NaX:

    for mol_NH4 in mols_NH4:

        # Initial amount of Calcite
        m0Calcite = 10.0

        # Define initial chemical state
        state = ChemicalState(system)
        state.setTemperature(25.0, "celsius")
        state.setPressure(1.0, "atm")
        # Seawater
        state.set("H2O" , 1.0    , "kg")
        state.set("Na+" , 1.10   , "mmol")
        state.set("Mg+2", 0.48   , "mmol")
        state.set("Ca+2", 1.90   , "mmol")
        state.set("NH4+", mol_NH4, "mmol")
        # Ammonia
        state.set("Calcite", m0Calcite, "mol")

        # Equilibrate chemical state corresponding to the seawater
        res = solver.solve(state)
        # Stop if the equilibration did not converge or failed
        if not res.optima.succeeded: continue

        # Update aqueous properties and evaluate pH
        aqprops.update(state)
        pH = float(aqprops.pH())

        # Exchanger
        state.set("NaX", mol_NaX, "mol")

        # Equilibrate the seawater with carbonates
        res = solver.solve(state)
        # Stop if the equilibration did not converge or failed
        if not res.optima.succeeded: continue

        # Update aqueous properties to evaluate ionic strength
        aqprops.update(state)
        chemprops = state.props()

        # Collect the value to be added to the dataframe in the following order
        # "amount_NaX", "amount_NH4", "amount_pH", "I", "amount_Ca", "amount_CaX2", "delta_Calcite"
        data.loc[len(data)] = [mol_NaX, mol_NH4, pH, float(aqprops.ionicStrength()),
                              float(chemprops.elementAmountInPhase("Ca", "AqueousPhase")), float(state.speciesAmount("CaX2")),
                              m0Calcite - float(state.speciesAmount("Calcite"))]

To visualize the obtained results, we use bokeh plotting library.

from bokeh.plotting import figure, show
from bokeh.models import HoverTool
from bokeh.io import output_notebook
from bokeh.models import ColumnDataSource
output_notebook()
Loading BokehJS ...

First, we plot the dependence of the dissolved calcite on the levels of pH and for different initial values of NaX providing the exchanger X-. We see that the solubility of calcite grows with an increase in acidity.

hovertool = HoverTool()
hovertool.tooltips = [("pH", "@pH"),
                      ("delta(Calcite)", "@delta_Calcite mol"),
                      ("amount(NaX)", "@amount_NaX mol")]

p = figure(
    title="DEPENDENCE OF CALCITE SOLUBILITY ON PH",
    x_axis_label=r'PH [-]',
    y_axis_label='AMOUNT OF DISSOLVED CALCITE [MOL]',
    sizing_mode="scale_width",
    plot_height=300)

p.add_tools(hovertool)

colors = ['teal', 'darkred', 'indigo', 'coral']
for mol_NaX, color in zip(mols_NaX, colors):
    df_NaX = ColumnDataSource(data[data['amount_NaX'] == mol_NaX])
    p.line("pH", "delta_Calcite", legend_label=f'{mol_NaX} mol of NaX', line_width=3, line_cap="round", line_color=color, source=df_NaX)

p.legend.location = 'top_right'
show(p)

Finally, we represent the amount of CaX2 resulting from the reaction NaX + Ca+2 = CaX2 + Na+, where the level of Ca+2 ions grows from the dissolved calcite. We see that at lower pH values, we obtain higher amounts of precipitated CaX2.

hovertool = HoverTool()
hovertool.tooltips = [("pH", "@pH"),
                      ("delta(Calcite)", "@delta_Calcite mol"),
                      ("amount(NaX)", "@amount_NaX mol")]

p = figure(
    title="DEPENDENCE OF CAX2 DISTRIBUTION ON PH",
    x_axis_label=r'PH [-]',
    y_axis_label='CAX2 AMOUNT [MOL]',
    sizing_mode="scale_width",
    plot_height=300)

p.add_tools(hovertool)

colors = ['teal', 'darkred', 'indigo', 'coral']
for mol_NaX, color in zip(mols_NaX, colors):
    df_NaX = ColumnDataSource(data[data['amount_NaX'] == mol_NaX])
    p.line("pH", "amount_CaX2", legend_label=f'{mol_NaX} mol of NaX', line_width=3, line_cap="round", line_color=color, source=df_NaX)

p.legend.location = 'top_right'
show(p)

The numerical experiment performed in this study showed the effect of carbonate mineral solubility (represented by calcite) on the exchangeable cations in the presence of the ion exchanger X- occupied by sodium ions.