{ "cells": [ { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Carbonate-rich lakes modeling on the early Earth\n", "\n", "
Written by Svetlana Kyas (ETH Zurich) on Mar 31th, 2022
\n", "\n", "```{attention}\n", "Always make sure you are using the [latest version of Reaktoro](https://anaconda.org/conda-forge/reaktoro). Otherwise, some new features documented on this website will not work on your machine and you may receive unintuitive errors. Follow these [update instructions](updating_reaktoro_via_conda) to get the latest version of Reaktoro!\n", "```\n", "\n", "Reaktoro can be used for various geobiological simulations. One of them is the modeling of carbonate-rich lakes, which were relatively common on the early Earth.\n", "\n", "```{note}\n", "This tutorial is one of two tutorials that follow the paper `Toner2020` and attempts to replicate the geobiological simulations performed in it (see also the second tutorial [**Phosphate accumulation in carbonate-rich brines**](geobiology-phreeqc-fixed-fugacity.ipynb)). This work was done in collaboration with Cara Magnabosco and Laura Murzakhmetov, ETH-Zurich.\n", "```\n", "\n", "## Solubility of phosphate in the hydroxyl-, fluorapatite- and calcite-rich lakes\n", "\n", "Carbonate-rich lakes can be explained by the strong chemical weathering of abundant, fresh volcanic rocks in the CO2-rich atmosphere of the early Earth. Weathering released phosphate from apatites (a group of phosphate minerals) and carbonate alkalinity from other minerals that accumulated in closed basins.\n", "\n", "In Reaktoro, the Earth’s CO2-rich atmosphere can be modelled by fixing the fugacity of the simulated chemical states. In particular, a consequence of early Earth’s CO2-rich atmosphere (corresponding to the partial pressure of CO2 from -2 to 0) is that it would have enhanced the weathering of hydroxyl- and fluorapatite in mafic rocks by lowering the pH of surface waters.\n", "\n", "Below, we load python libraries important to carry out the simulations in this tutorial." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from reaktoro import *\n", "import numpy as np\n", "import pandas as pd\n", "import math as math" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The database is loaded from the extended version of the phreeqc database, including the properties of fluorapatite and hydroxylapatite minerals. The chemical system is composed of an aqueous phase and minerals calcite (CaCO3), fluorapatite (Ca5(F)(PO4)3), and hydroxyapatite (Ca5(PO4)3OH). The corresponding aqueous and chemical properties are also defined below." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "# Fluorapatite\n", "# Ca5(F)(PO4)3 = 5Ca+2 + F- + 3PO4-3\n", "# log_k -59.6\n", "# -analytical_expression -1917.945184 0 87834.57783 631.9611081 0 0\n", "# Hydroxylapatite\n", "# Ca5(OH)(PO4)3 = 5Ca+2 + OH- + 3PO4-3\n", "# log_k -58.517\n", "# -analytical_expression -1.6657 -0.098215 -8219.41 0 0 0\n", "db = PhreeqcDatabase.fromFile('phreeqc-extended.dat')\n", "\n", "# Define the aqueous phase\n", "solution = AqueousPhase(speciate(StringList(\"H O C Na Cl Ca P\")))\n", "solution.set(chain(\n", " ActivityModelHKF(),\n", " ActivityModelDrummond(\"CO2\")\n", "))\n", "\n", "# Define minerals' phases\n", "minerals = MineralPhases(\"Calcite Fluorapatite Hydroxylapatite\")\n", "\n", "# Define the chemical system\n", "system = ChemicalSystem(db, solution, minerals)\n", "\n", "# Define aqueous and chemical properties\n", "aprops = AqueousProps(system)\n", "props = ChemicalProps(system)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To tell the solver that fugacity will be constrained in this chemical system, we need to define equilibrium specifications and define corresponding conditions. The first specifies what will be a constraint and the second by which value (defined below for the range of fugacities)." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "# Define equilibrium specifications\n", "specs = EquilibriumSpecs(system)\n", "specs.temperature()\n", "specs.pressure()\n", "specs.fugacity(\"CO2\")\n", "\n", "# Define conditions to be satisfied at the chemical equilibrium state\n", "conditions = EquilibriumConditions(specs)\n", "conditions.pressure(1.0, \"atm\")\n", "\n", "# Define the equilibrium solver\n", "solver = EquilibriumSolver(specs)\n", "opts = EquilibriumOptions()\n", "opts.epsilon = 1e-13\n", "solver.setOptions(opts)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Below, we perform a series of experiments to determine how much phosphate can be accumulated by abiotic processes in carbonate-rich lakes. In the other words, we calculate the solubility of fluorapatite and hydroxyapatite in the presence of calcite buffer as a function of temperature and CO2 partial pressure.\n", "\n", "The block below defines the array of the partial CO2 pressures and temperatures as well as the data blocks storing results for different temperatures. Finally, we run equilibrium calculations in the loop for different partial CO2 pressure and temperatures and collect the following data:\n", "* pH of the equilibrated state,\n", "* phosphate level in the solution (the amount of element P), and\n", "* calcium level in the solution (the amount of element Ca)." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "# Auxiliary arrays\n", "num_temperatures = 3\n", "num_log10pCO2s = 51\n", "temperatures = np.array([0, 25, 50])\n", "co2pressures = np.linspace(-4.0, 0.0, num=num_log10pCO2s)\n", "\n", "# Output dataframe\n", "data = pd.DataFrame(columns=[\"T\", \"ppCO2\", \"pH\", \"amount_P\", \"amount_Ca\"])\n", "\n", "for ppCO2 in co2pressures:\n", " for T in temperatures:\n", "\n", " conditions.temperature(T, \"celsius\")\n", " conditions.fugacity(\"CO2\", 10 ** ppCO2, 'atm')\n", "\n", " state = ChemicalState(system)\n", " state.set(\"H2O\" , 1.0, \"kg\")\n", " state.set(\"Calcite\" , 10.0, \"mol\")\n", " state.set(\"Fluorapatite\" , 10.0, \"mol\")\n", " state.set(\"Hydroxylapatite\", 10.0, \"mol\")\n", " state.set(\"CO2\" , 100.0, \"mol\")\n", "\n", " # Equilibrate the solution with the given initial chemical state and desired conditions at the equilibrium\n", " solver.solve(state, conditions)\n", "\n", " aprops.update(state)\n", " props.update(state)\n", "\n", " # Collect the value to be added to the dataframe in the following order\n", " # \"T\", \"ppCO2\", \"pH\", \"amount_P\", \"amount_Ca\"\n", " data.loc[len(data)] = [T, ppCO2, float(aprops.pH()),\n", " float(props.elementAmountInPhase(\"P\", \"AqueousPhase\")),\n", " float(props.elementAmountInPhase(\"Ca\", \"AqueousPhase\"))]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Apatites are more soluble at lower pH and weather more rapidly in CO2-acidified stream and rainwater, resulting in potentially high phosphate fluxes to carbonate-rich lakes on the early Earth. Below, we plot the dependency of the phosphates and carbonates solubility on the partial CO2 pressures, which rises with growing log10(pCO2)." ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false, "jupyter": { "outputs_hidden": false }, "tags": [ "hide_input" ] }, "outputs": [ { "data": { "text/html": [ "\n" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "application/javascript": "(function(root) {\n function now() {\n return new Date();\n }\n\n const force = true;\n\n if (typeof root._bokeh_onload_callbacks === \"undefined\" || force === true) {\n root._bokeh_onload_callbacks = [];\n root._bokeh_is_loading = undefined;\n }\n\nconst JS_MIME_TYPE = 'application/javascript';\n const HTML_MIME_TYPE = 'text/html';\n const EXEC_MIME_TYPE = 'application/vnd.bokehjs_exec.v0+json';\n const CLASS_NAME = 'output_bokeh rendered_html';\n\n /**\n * Render data to the DOM node\n */\n function render(props, node) {\n const script = document.createElement(\"script\");\n node.appendChild(script);\n }\n\n /**\n * Handle when an output is cleared or removed\n */\n function handleClearOutput(event, handle) {\n const cell = handle.cell;\n\n const id = cell.output_area._bokeh_element_id;\n const server_id = cell.output_area._bokeh_server_id;\n // Clean up Bokeh references\n if (id != null && id in Bokeh.index) {\n Bokeh.index[id].model.document.clear();\n delete Bokeh.index[id];\n }\n\n if (server_id !== undefined) {\n // Clean up Bokeh references\n const cmd_clean = \"from bokeh.io.state import curstate; print(curstate().uuid_to_server['\" + server_id + \"'].get_sessions()[0].document.roots[0]._id)\";\n cell.notebook.kernel.execute(cmd_clean, {\n iopub: {\n output: function(msg) {\n const id = msg.content.text.trim();\n if (id in Bokeh.index) {\n Bokeh.index[id].model.document.clear();\n delete Bokeh.index[id];\n }\n }\n }\n });\n // Destroy server and session\n const cmd_destroy = \"import bokeh.io.notebook as ion; ion.destroy_server('\" + server_id + \"')\";\n cell.notebook.kernel.execute(cmd_destroy);\n }\n }\n\n /**\n * Handle when a new output is added\n */\n function handleAddOutput(event, handle) {\n const output_area = handle.output_area;\n const output = handle.output;\n\n // limit handleAddOutput to display_data with EXEC_MIME_TYPE content only\n if ((output.output_type != \"display_data\") || (!Object.prototype.hasOwnProperty.call(output.data, EXEC_MIME_TYPE))) {\n return\n }\n\n const toinsert = output_area.element.find(\".\" + CLASS_NAME.split(' ')[0]);\n\n if (output.metadata[EXEC_MIME_TYPE][\"id\"] !== undefined) {\n toinsert[toinsert.length - 1].firstChild.textContent = output.data[JS_MIME_TYPE];\n // store reference to embed id on output_area\n output_area._bokeh_element_id = output.metadata[EXEC_MIME_TYPE][\"id\"];\n }\n if (output.metadata[EXEC_MIME_TYPE][\"server_id\"] !== undefined) {\n const bk_div = document.createElement(\"div\");\n bk_div.innerHTML = output.data[HTML_MIME_TYPE];\n const script_attrs = bk_div.children[0].attributes;\n for (let i = 0; i < script_attrs.length; i++) {\n toinsert[toinsert.length - 1].firstChild.setAttribute(script_attrs[i].name, script_attrs[i].value);\n toinsert[toinsert.length - 1].firstChild.textContent = bk_div.children[0].textContent\n }\n // store reference to server id on output_area\n output_area._bokeh_server_id = output.metadata[EXEC_MIME_TYPE][\"server_id\"];\n }\n }\n\n function register_renderer(events, OutputArea) {\n\n function append_mime(data, metadata, element) {\n // create a DOM node to render to\n const toinsert = this.create_output_subarea(\n metadata,\n CLASS_NAME,\n EXEC_MIME_TYPE\n );\n this.keyboard_manager.register_events(toinsert);\n // Render to node\n const props = {data: data, metadata: metadata[EXEC_MIME_TYPE]};\n render(props, toinsert[toinsert.length - 1]);\n element.append(toinsert);\n return toinsert\n }\n\n /* Handle when an output is cleared or removed */\n events.on('clear_output.CodeCell', handleClearOutput);\n events.on('delete.Cell', handleClearOutput);\n\n /* Handle when a new output is added */\n events.on('output_added.OutputArea', handleAddOutput);\n\n /**\n * Register the mime type and append_mime function with output_area\n */\n OutputArea.prototype.register_mime_type(EXEC_MIME_TYPE, append_mime, {\n /* Is output safe? */\n safe: true,\n /* Index of renderer in `output_area.display_order` */\n index: 0\n });\n }\n\n // register the mime type if in Jupyter Notebook environment and previously unregistered\n if (root.Jupyter !== undefined) {\n const events = require('base/js/events');\n const OutputArea = require('notebook/js/outputarea').OutputArea;\n\n if (OutputArea.prototype.mime_types().indexOf(EXEC_MIME_TYPE) == -1) {\n register_renderer(events, OutputArea);\n }\n }\n if (typeof (root._bokeh_timeout) === \"undefined\" || force === true) {\n root._bokeh_timeout = Date.now() + 5000;\n root._bokeh_failed_load = false;\n }\n\n const NB_LOAD_WARNING = {'data': {'text/html':\n \"\\n\"+\n \"BokehJS does not appear to have successfully loaded. If loading BokehJS from CDN, this \\n\"+\n \"may be due to a slow or bad network connection. Possible fixes:\\n\"+\n \"
\\n\"+\n \"\\n\"+\n \"from bokeh.resources import INLINE\\n\"+\n \"output_notebook(resources=INLINE)\\n\"+\n \"
\\n\"+\n \"\\n\"+\n \"BokehJS does not appear to have successfully loaded. If loading BokehJS from CDN, this \\n\"+\n \"may be due to a slow or bad network connection. Possible fixes:\\n\"+\n \"
\\n\"+\n \"\\n\"+\n \"from bokeh.resources import INLINE\\n\"+\n \"output_notebook(resources=INLINE)\\n\"+\n \"
\\n\"+\n \"