{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "4c6355f0",
   "metadata": {},
   "source": [
    "# Solubility of sodium-chloride in water\n",
    "\n",
    "This tutorial demonstrates how to use Reaktoro to study the solubility of sodium chloride salt in water. In particular,\n",
    "we simulate the precipitation of halite for very concentrated NaCl-brines by preventing excessive complexation of\n",
    "Na<sup>+</sup> and Cl<sup>-</sup> into NaCl(aq).\n",
    "\n",
    "First, we import everything from the `reaktoro` package and load the SUPCRT98 database"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "30ffd2a8",
   "metadata": {},
   "outputs": [],
   "source": [
    "from reaktoro import *\n",
    "\n",
    "db = Database(\"supcrt98.xml\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4a9a0120",
   "metadata": {},
   "source": [
    "To specify aqueous and mineral phases, we use the code below.  We set the Setschenow activity model to observe Halite\n",
    "precipitation. Otherwise, all the NaCl gets split into Na<sup>+</sup>, Cl<sup>-</sup> and NaCl(aq) species. Unlike\n",
    "the default HKF activity model, the Setschenow model uses the ionic strength to compute the activity coefficient\n",
    "of neutral species (e.g., NaCl(aq)). This is similar to the approach used in PHREEQC package."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fbf2303e",
   "metadata": {},
   "outputs": [],
   "source": [
    "editor = ChemicalEditor(db)\n",
    "aqueousphase = editor.addAqueousPhaseWithElements(\"H O Na Cl\")\n",
    "aqueousphase.setActivityModelSetschenow(\"NaCl(aq)\", 1.0)\n",
    "editor.addMineralPhase(\"Halite\")\n",
    "\n",
    "system = ChemicalSystem(editor)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2b7481b1",
   "metadata": {},
   "source": [
    "The temperature and pressure are chosen to be 30 &deg;C and 300 bar (which is realistic for the oil and gas sites).\n",
    "Moreover, we consider a very saline NaCl-brine, the mixture of 1 kg of water and 14.0 mol of NaCl salt."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2adf1470",
   "metadata": {},
   "outputs": [],
   "source": [
    "problem = EquilibriumProblem(system)\n",
    "problem.setTemperature(30.0, \"celsius\")\n",
    "problem.setPressure(300.0, \"bar\")\n",
    "problem.add(\"H2O\", 1.0, \"kg\")\n",
    "problem.add(\"NaCl\", 14.0, \"mol\")\n",
    "\n",
    "state = equilibrate(problem)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3c4dbe11",
   "metadata": {},
   "source": [
    "The full result of the calculation is output below. The obtained ionic strength is about 6.15 molal, and the amount of\n",
    "halite has 7.85 mol."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "01092925",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(state)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "71cc5a28",
   "metadata": {},
   "source": [
    "To access the halite molality only, we use a function below, where molality is calculate as amount of halite (in mol)\n",
    "divided by the mass of water (set to 1 kg the definition of\n",
    "[EquilibriumProblem](https://reaktoro.org/cpp/classReaktoro_1_1EquilibriumProblem.html)):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "756e14fb",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Halite molality = {:4.2f} molal\".format(state.speciesAmount(\"Halite\")))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "16d7e7ca",
   "metadata": {},
   "source": [
    "To calculate the ionic strength directly, we need first to create function `evaluate_I` that accepts the chemical\n",
    "state properties. For that purpose, the code below can be used:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4584eab4",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Access chemical state properties\n",
    "props = state.properties()\n",
    "# Create function for evaluating the ionic strength\n",
    "evaluate_I = ChemicalProperty.ionicStrength(system)\n",
    "# Calculate the ionic strength\n",
    "I = evaluate_I(props)\n",
    "print(\"I = {:4.2f} molal\".format(I.val))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}