{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Multi-fluid EOS\n",
    "\n",
    "Peering into the innards of teqp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:02.562349Z",
     "iopub.status.busy": "2025-10-15T23:12:02.562226Z",
     "iopub.status.idle": "2025-10-15T23:12:02.768729Z",
     "shell.execute_reply": "2025-10-15T23:12:02.768395Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'0.23.1'"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import timeit, json\n",
    "import pandas\n",
    "import numpy as np\n",
    "import teqp\n",
    "teqp.__version__"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Ancillary Equations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ancillary equations are provided along with multiparameter equations of state. The give a good *approximation* to the phase equilibrium densities.  There are routines in teqp to use the ancillary equations provided with the EOS. First a class containing the ancillary equations is obtained, then methods on that class are called"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:02.793544Z",
     "iopub.status.busy": "2025-10-15T23:12:02.793264Z",
     "iopub.status.idle": "2025-10-15T23:12:02.822384Z",
     "shell.execute_reply": "2025-10-15T23:12:02.821841Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Densities are: 27357.335621492966 42.04100696197727 mol/m^3\n"
     ]
    }
   ],
   "source": [
    "model = teqp.build_multifluid_model([\"Methane\"], teqp.get_datapath())\n",
    "anc = model.build_ancillaries()\n",
    "T = 100.0 # [K]\n",
    "rhoL, rhoV = anc.rhoL(T), anc.rhoV(T)\n",
    "print('Densities are:', rhoL, rhoV, 'mol/m^3')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "But those densities do not correspond to the *true* phase equilibrium solution, so we need to polish the solution:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:02.823680Z",
     "iopub.status.busy": "2025-10-15T23:12:02.823538Z",
     "iopub.status.idle": "2025-10-15T23:12:02.827193Z",
     "shell.execute_reply": "2025-10-15T23:12:02.826519Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "VLE densities are: 27357.147019094475 42.04798227835163 mol/m^3\n"
     ]
    }
   ],
   "source": [
    "Niter = 10\n",
    "rhoLtrue, rhoVtrue = model.pure_VLE_T(T, rhoL, rhoV, Niter)\n",
    "print('VLE densities are:', rhoLtrue, rhoVtrue, 'mol/m^3')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And looking the densities, they are slightly different after the phase equilibrium calculation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Ammonia-Water\n",
    "\n",
    "Tillner-Roth and Friend provided a hard-coded model that is in a form not compatible with the other multi-fluid models. It is available via the high-level factory function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:02.828461Z",
     "iopub.status.busy": "2025-10-15T23:12:02.828321Z",
     "iopub.status.idle": "2025-10-15T23:12:02.832448Z",
     "shell.execute_reply": "2025-10-15T23:12:02.831348Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-0.09731055757504622"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "AW = teqp.AmmoniaWaterTillnerRoth()\n",
    "AW.get_Ar01(300, 300, np.array([0.9, 0.0]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Pure fluid loading"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:02.833660Z",
     "iopub.status.busy": "2025-10-15T23:12:02.833513Z",
     "iopub.status.idle": "2025-10-15T23:12:04.863875Z",
     "shell.execute_reply": "2025-10-15T23:12:04.863315Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "25 ms ± 93 μs per loop (mean ± std. dev. of 7 runs, 10 loops each)\n"
     ]
    }
   ],
   "source": [
    "# By default teqp looks for fluids relative to the set of fluids in ROOT/dev/fluids\n",
    "# The name (case-sensitive) should match the .json file, without the json extension.\n",
    "%timeit model = teqp.build_multifluid_model([\"Methane\"], teqp.get_datapath())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:04.865455Z",
     "iopub.status.busy": "2025-10-15T23:12:04.865298Z",
     "iopub.status.idle": "2025-10-15T23:12:06.883062Z",
     "shell.execute_reply": "2025-10-15T23:12:06.882431Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "24.8 ms ± 92.3 μs per loop (mean ± std. dev. of 7 runs, 10 loops each)\n"
     ]
    }
   ],
   "source": [
    "# And if you provide valid aliases, alias lookup will be used to resolve the name\n",
    "# But beware, this is rather a lot slower than the above because all fluid files need to be read\n",
    "# in to build the alias map\n",
    "%timeit model = teqp.build_multifluid_model([\"n-C1H4\"], teqp.get_datapath())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So, how to make it faster? Only do it once and cache"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:06.884810Z",
     "iopub.status.busy": "2025-10-15T23:12:06.884625Z",
     "iopub.status.idle": "2025-10-15T23:12:08.893688Z",
     "shell.execute_reply": "2025-10-15T23:12:08.893251Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "24.4 ms ± 60.8 μs per loop (mean ± std. dev. of 7 runs, 10 loops each)\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "['1,2-DICHLOROETHANE',\n",
       " '1,2-dichloroethane',\n",
       " '1-BUTENE',\n",
       " '1-Butene',\n",
       " '100-41-4',\n",
       " '10024-97-2',\n",
       " '102687-65-0',\n",
       " '106-42-3',\n",
       " '106-97-8',\n",
       " '106-98-9']"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Here is the set of possible aliases to absolute paths of files\n",
    "# Building this map takes a little while (somewhat faster in C++) due to all the file reads\n",
    "# If you know your files will not change, good idea to build this alias map yourself.\n",
    "%timeit aliasmap = teqp.build_alias_map(teqp.get_datapath())\n",
    "aliasmap = teqp.build_alias_map(teqp.get_datapath())\n",
    "list(aliasmap.keys())[0:10] # the first 10 aliases in the dict"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:08.895092Z",
     "iopub.status.busy": "2025-10-15T23:12:08.894960Z",
     "iopub.status.idle": "2025-10-15T23:12:12.365499Z",
     "shell.execute_reply": "2025-10-15T23:12:12.364988Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "415 μs ± 5.39 μs per loop (mean ± std. dev. of 7 runs, 1,000 loops each)\n"
     ]
    }
   ],
   "source": [
    "# Then load the absolute paths from the alias map, \n",
    "# which will guarantee that you hit exactly what you were looking for,\n",
    "# resolving aliases as needed\n",
    "identifiers = [aliasmap[n] for n in [\"n-C1H4\"]]\n",
    "%timeit model = teqp.build_multifluid_model(identifiers, teqp.get_datapath())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "At some point soon teqp will support in-memory loading of JSON data for the pure components, without requiring reads from the operating system"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:12.367115Z",
     "iopub.status.busy": "2025-10-15T23:12:12.366976Z",
     "iopub.status.idle": "2025-10-15T23:12:12.370223Z",
     "shell.execute_reply": "2025-10-15T23:12:12.369729Z"
    }
   },
   "outputs": [],
   "source": [
    "# And you can also load the JSON that teqp is loading for the pure fluids\n",
    "pureJSON = teqp.collect_component_json(['Neon','Hydrogen'], teqp.get_datapath())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Mixture model loading"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:12.371381Z",
     "iopub.status.busy": "2025-10-15T23:12:12.371267Z",
     "iopub.status.idle": "2025-10-15T23:12:12.374543Z",
     "shell.execute_reply": "2025-10-15T23:12:12.373982Z"
    }
   },
   "outputs": [],
   "source": [
    "# Load the default JSON for the binary interaction parameters\n",
    "BIP = json.load(open(teqp.get_datapath()+'/dev/mixtures/mixture_binary_pairs.json'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:12.375530Z",
     "iopub.status.busy": "2025-10-15T23:12:12.375422Z",
     "iopub.status.idle": "2025-10-15T23:12:12.380684Z",
     "shell.execute_reply": "2025-10-15T23:12:12.380319Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'BibTeX': 'Kunz-JCED-2012',\n",
       " 'CAS1': '74-82-8',\n",
       " 'CAS2': '74-84-0',\n",
       " 'F': 1.0,\n",
       " 'Name1': 'Methane',\n",
       " 'Name2': 'Ethane',\n",
       " 'betaT': 0.996336508,\n",
       " 'betaV': 0.997547866,\n",
       " 'function': 'Methane-Ethane',\n",
       " 'gammaT': 1.049707697,\n",
       " 'gammaV': 1.006617867}"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# You can obtain interaction parameters either by pairs of names, where name is the name that teqp uses, the [\"INFO\"][\"NAME\"] field\n",
    "params, swap_needed = teqp.get_BIPdep(BIP, ['Methane','Ethane'])\n",
    "params"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:12.381939Z",
     "iopub.status.busy": "2025-10-15T23:12:12.381823Z",
     "iopub.status.idle": "2025-10-15T23:12:12.387053Z",
     "shell.execute_reply": "2025-10-15T23:12:12.386674Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'BibTeX': 'Kunz-JCED-2012',\n",
       " 'CAS1': '74-82-8',\n",
       " 'CAS2': '74-84-0',\n",
       " 'F': 1.0,\n",
       " 'Name1': 'Methane',\n",
       " 'Name2': 'Ethane',\n",
       " 'betaT': 0.996336508,\n",
       " 'betaV': 0.997547866,\n",
       " 'function': 'Methane-Ethane',\n",
       " 'gammaT': 1.049707697,\n",
       " 'gammaV': 1.006617867}"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Or also by CAS#\n",
    "params, swap_needed = teqp.get_BIPdep(BIP, ['74-82-8','74-84-0'])\n",
    "params"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:12.388153Z",
     "iopub.status.busy": "2025-10-15T23:12:12.388040Z",
     "iopub.status.idle": "2025-10-15T23:12:12.566431Z",
     "shell.execute_reply": "2025-10-15T23:12:12.565990Z"
    },
    "tags": [
     "raises-exception"
    ]
   },
   "outputs": [
    {
     "ename": "ValueError",
     "evalue": "Can't match the binary pair for: 74-82-8/Ethane",
     "output_type": "error",
     "traceback": [
      "\u001b[31m---------------------------------------------------------------------------\u001b[39m",
      "\u001b[31mValueError\u001b[39m                                Traceback (most recent call last)",
      "\u001b[36mCell\u001b[39m\u001b[36m \u001b[39m\u001b[32mIn[13]\u001b[39m\u001b[32m, line 2\u001b[39m\n\u001b[32m      1\u001b[39m \u001b[38;5;66;03m# But mixing is not allowed\u001b[39;00m\n\u001b[32m----> \u001b[39m\u001b[32m2\u001b[39m params, swap_needed = \u001b[43mteqp\u001b[49m\u001b[43m.\u001b[49m\u001b[43mget_BIPdep\u001b[49m\u001b[43m(\u001b[49m\u001b[43mBIP\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43m[\u001b[49m\u001b[33;43m'\u001b[39;49m\u001b[33;43m74-82-8\u001b[39;49m\u001b[33;43m'\u001b[39;49m\u001b[43m,\u001b[49m\u001b[33;43m'\u001b[39;49m\u001b[33;43mEthane\u001b[39;49m\u001b[33;43m'\u001b[39;49m\u001b[43m]\u001b[49m\u001b[43m)\u001b[49m\n\u001b[32m      3\u001b[39m params\n",
      "\u001b[31mValueError\u001b[39m: Can't match the binary pair for: 74-82-8/Ethane"
     ]
    }
   ],
   "source": [
    "# But mixing is not allowed\n",
    "params, swap_needed = teqp.get_BIPdep(BIP, ['74-82-8','Ethane'])\n",
    "params"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Estimation of interaction parameters\n",
    "\n",
    "Estimation of interaction parameters can be used when no mixture model is present.  The ``flags`` keyword argument allows the user to control how estimation is applied. The ``flags`` keyword argument should be a dictionary, with keys of ``\"estimate\"`` to provide the desired estimation scheme as-needed. For now, the only allowed estimation scheme is ``Lorentz-Berthelot``. \n",
    "\n",
    "If it is desired to force the estimation, the ``\"force-estimate\"`` to force the use of the provided estimation scheme for all binaries, even when a proper mixture model is available. The value associated with ``\"force-estimate\"`` is ignored."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:12.568047Z",
     "iopub.status.busy": "2025-10-15T23:12:12.567923Z",
     "iopub.status.idle": "2025-10-15T23:12:12.574408Z",
     "shell.execute_reply": "2025-10-15T23:12:12.573900Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'F': 0.0, 'betaT': 1.0, 'betaV': 1.0, 'gammaT': 1.0, 'gammaV': 1.0}"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "params, swap_needed = teqp.get_BIPdep(BIP, ['74-82-8','74-84-0'], flags={'force-estimate':'yes', 'estimate': 'Lorentz-Berthelot'})\n",
    "params"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:12.575588Z",
     "iopub.status.busy": "2025-10-15T23:12:12.575445Z",
     "iopub.status.idle": "2025-10-15T23:12:12.581520Z",
     "shell.execute_reply": "2025-10-15T23:12:12.580972Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'BibTeX': 'Kunz-JCED-2012',\n",
       " 'CAS1': '74-82-8',\n",
       " 'CAS2': '74-84-0',\n",
       " 'F': 1.0,\n",
       " 'Name1': 'Methane',\n",
       " 'Name2': 'Ethane',\n",
       " 'betaT': 0.996336508,\n",
       " 'betaV': 0.997547866,\n",
       " 'function': 'Methane-Ethane',\n",
       " 'gammaT': 1.049707697,\n",
       " 'gammaV': 1.006617867}"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# And without the force, the forcing is ignored\n",
    "params, swap_needed = teqp.get_BIPdep(BIP, ['74-82-8','74-84-0'], flags={'estimate': 'Lorentz-Berthelot'})\n",
    "params"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-10-15T23:12:12.582528Z",
     "iopub.status.busy": "2025-10-15T23:12:12.582416Z",
     "iopub.status.idle": "2025-10-15T23:12:12.613459Z",
     "shell.execute_reply": "2025-10-15T23:12:12.612998Z"
    }
   },
   "outputs": [],
   "source": [
    "# And the same flags can be passed to the multifluid model constructor\n",
    "model = teqp.build_multifluid_model(\n",
    "    ['74-82-8','74-84-0'], \n",
    "    teqp.get_datapath(), \n",
    "    flags={'force-estimate':'yes', 'estimate': 'Lorentz-Berthelot'})"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Edit Metadata",
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}