{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from chempy.kinetics.arrhenius import ArrheniusParamWithUnits\n",
    "from chempy.kinetics.rates import MassAction\n",
    "from chempy.kinetics._rates import RTPoly\n",
    "from chempy.reactionsystem import ReactionSystem\n",
    "from chempy.units import patched_numpy, Backend, unit_of, simplified, default_constants as dc, default_units as u\n",
    "from chempy.util._expr import Expr, Constant\n",
    "from chempy.util.parsing import get_parsing_context"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def log10recipoly(*args):\n",
    "    return 10**RTPoly(*args)\n",
    "rpkw = dict(globals_=dict(u=u, log10recipoly=log10recipoly, Constant=Constant, MassAction=MassAction))\n",
    "lgRTp = \"MassAction([Constant(1/u.M/u.s) * log10recipoly([x*u.K**i for i, x in enumerate([13.339, -2.22e3, 7.333e5, -1.065e8])])])\"\n",
    "rs1 = ReactionSystem.from_string(\"OH + OH- -> H2O + O-; %s; name='R27f'\" % lgRTp, rxn_parse_kwargs=rpkw)\n",
    "rs1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rs1.rxns[0].param.rate_coeff({'temperature': 298*u.K})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "be = Backend()\n",
    "variables = {'temperature': 298.15*u.kelvin}\n",
    "rates = rs1['R27f'].rate(dict(variables, **{'OH': 1e-11*u.molar, 'OH-': 1e-7*u.molar}), backend=be)\n",
    "for sk in rs1.substances:\n",
    "    assert sk in rates\n",
    "    assert unit_of(rates[sk]) == u.molar/u.second\n",
    "rates"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class GibbsExpr(Expr):\n",
    "    nargs = 4\n",
    "    parameter_keys = ('temperature',)\n",
    "\n",
    "    def __call__(self, variables, backend=patched_numpy, **kwargs):\n",
    "        dS_over_R, dCp_over_R, dH_over_R, Tref = map(simplified, self.all_args(variables, backend=backend))\n",
    "        T, = self.all_params(variables, backend=backend)\n",
    "        return backend.exp(dS_over_R)*(T/Tref)**dCp_over_R*backend.exp(-dH_over_R/T)\n",
    "\n",
    "GeNH3 = GibbsExpr([\n",
    "    18.8*u.cal/u.K/u.mol/dc.molar_gas_constant,\n",
    "    52*u.cal/u.K/u.mol/dc.molar_gas_constant,\n",
    "    -0.87e3*u.cal/u.mol/dc.molar_gas_constant,\n",
    "    298.15*u.K\n",
    "])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dCp_R = 18.8*u.cal/u.K/u.mol/dc.molar_gas_constant\n",
    "dCp_R, dCp_R.simplified"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "str_rs2 = \"\"\"\n",
    "NH4+ + OH- -> NH3 + H2O; rs1['R27f'].param; name='ammonium_hydroxide'\n",
    "NH3 + H2O -> NH4+ + OH-; MassAction(rs1['R27f'].param.args[0]/GeNH3); name='ammonia_water'\n",
    "\"\"\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "globals_ = get_parsing_context()\n",
    "globals_['rs1'] = rs1\n",
    "globals_['GeNH3'] = GeNH3\n",
    "globals_['Arrh'] = ArrheniusParamWithUnits\n",
    "rs2 = ReactionSystem.from_string(str_rs2, rxn_parse_kwargs=dict(globals_=globals_))\n",
    "rs2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rs2.rxns[-1].rate_expr()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rs2.rxns[-1].param.rate_coeff({'temperature': 298*u.K})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for rk, cd in dict(ammonium_hydroxide={'NH4+': 1e-3*u.molar, 'OH-': 1e-7*u.molar},\n",
    "                   ammonia_water={'NH3': 0.42*u.molar, 'H2O': 55.4*u.molar}).items():\n",
    "    print(rs2[rk].rate(dict(variables, **cd), backend=be))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}