{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "### NMODL SympySolver - sparse\n", "\n", "This notebook describes the implementation of the `sparse` part of the `SympySolverVisitor`, which solves the systems of ODEs defined in `DERIVATIVE` blocks when these ODEs are *linear* and *coupled*.\n", "\n", "For a higher level overview of the approach to solving ODEs in NMODL, please see the [nmodl-odes-overview](nmodl-odes-overview.ipynb) notebook. \n", "\n", "For a more general tutorial on using the NMODL python interface, please see the [tutorial notebook](nmodl-python-tutorial.ipynb).\n", "\n", "***\n", "\n", "#### Implementation\n", "The `SympySolverVisitor` for solver method `sparse` does the following:\n", "\n", "* Construct the implicit Euler equations to convert the sysytem of ODEs to a `LINEAR` block\n", "* This `LINEAR` block is then solved as described in [nmodl-linear-solver](nmodl-linear-solver.ipynb)\n", "\n", "***\n", "\n", "#### Implementation Tests\n", "The unit tests may be helpful to understand what these functions are doing\n", " - `SympySolverVisitor` tests are located in [test/visitor/sympy_solver.cpp](https://github.com/BlueBrain/nmodl/blob/master/test/visitor/sympy_solver.cpp), and tests using the `sparse` solver method have the tag \"`[sparse]`\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Examples" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%%capture\n", "! pip install nmodl" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "import nmodl.dsl as nmodl\n", "\n", "\n", "def run_sympy_solver(mod_string, cse=False):\n", " # parse NMDOL file (supplied as a string) into AST\n", " driver = nmodl.NmodlDriver()\n", " AST = driver.parse_string(mod_string)\n", " # run SymtabVisitor to generate Symbol Table\n", " nmodl.symtab.SymtabVisitor().visit_program(AST)\n", " # constant folding, inlining & local variable renaming passes\n", " nmodl.visitor.ConstantFolderVisitor().visit_program(AST)\n", " nmodl.visitor.InlineVisitor().visit_program(AST)\n", " nmodl.visitor.LocalVarRenameVisitor().visit_program(AST)\n", " # run SympySolver visitor\n", " nmodl.visitor.SympySolverVisitor().visit_program(AST)\n", " # return contents of new DERIVATIVE block as list of strings\n", " return nmodl.to_nmodl(\n", " nmodl.visitor.AstLookupVisitor().lookup(\n", " AST, nmodl.ast.AstNodeType.DERIVATIVE_BLOCK\n", " )[0]\n", " )" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "##### Ex. 1\n", "2 coupled linear ODEs" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "DERIVATIVE states {\n", " LOCAL old_mc, old_m, tmp0, tmp1, tmp2\n", " old_mc = mc\n", " old_m = m\n", " tmp0 = a*dt\n", " tmp1 = b*dt\n", " tmp2 = 1/(tmp0+tmp1+1)\n", " mc = tmp2*(old_m*tmp1+old_mc*tmp1+old_mc)\n", " m = tmp2*(old_m*tmp0+old_m+old_mc*tmp0)\n", "}\n" ] } ], "source": [ "ex1 = \"\"\"\n", "BREAKPOINT {\n", " SOLVE states METHOD sparse\n", "}\n", "STATE {\n", " mc m\n", "}\n", "DERIVATIVE states {\n", " mc' = -a*mc + b*m\n", " m' = a*mc - b*m\n", "}\n", "\"\"\"\n", "print(run_sympy_solver(ex1, cse=True))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "##### Ex. 2\n", "5 coupled linear ODEs" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "DERIVATIVE ihkin {\n", " EIGEN_LINEAR_SOLVE[5]{\n", " LOCAL alpha, beta, k3p, k4, k1ca, k2, old_c1, old_o1, old_o2, old_p0, old_p1\n", " }{\n", " evaluate_fct(v, cai)\n", " old_c1 = c1\n", " old_o1 = o1\n", " old_o2 = o2\n", " old_p0 = p0\n", " old_p1 = p1\n", " }{\n", " X[0] = c1\n", " X[1] = o1\n", " X[2] = o2\n", " X[3] = p0\n", " X[4] = p1\n", " F[0] = -old_c1\n", " F[1] = -old_o1\n", " F[2] = -old_o2\n", " F[3] = -old_p0\n", " F[4] = -old_p1\n", " J[0] = -alpha*dt-1\n", " J[5] = beta*dt\n", " J[10] = 0\n", " J[15] = 0\n", " J[20] = 0\n", " J[1] = alpha*dt\n", " J[6] = -beta*dt-dt*k3p-1\n", " J[11] = dt*k4\n", " J[16] = 0\n", " J[21] = 0\n", " J[2] = 0\n", " J[7] = dt*k3p\n", " J[12] = -dt*k4-1\n", " J[17] = 0\n", " J[22] = 0\n", " J[3] = 0\n", " J[8] = 0\n", " J[13] = 0\n", " J[18] = -dt*k1ca-1\n", " J[23] = dt*k2\n", " J[4] = 0\n", " J[9] = 0\n", " J[14] = 0\n", " J[19] = dt*k1ca\n", " J[24] = -dt*k2-1\n", " }{\n", " c1 = X[0]\n", " o1 = X[1]\n", " o2 = X[2]\n", " p0 = X[3]\n", " p1 = X[4]\n", " }{\n", " }\n", "}\n" ] } ], "source": [ "ex2 = \"\"\"\n", "STATE {\n", " c1 o1 o2 p0 p1\n", "}\n", "BREAKPOINT {\n", " SOLVE ihkin METHOD sparse\n", "}\n", "DERIVATIVE ihkin {\n", " LOCAL alpha, beta, k3p, k4, k1ca, k2\n", " evaluate_fct(v, cai)\n", " c1' = (-1*(alpha*c1-beta*o1))\n", " o1' = (1*(alpha*c1-beta*o1))+(-1*(k3p*o1-k4*o2))\n", " o2' = (1*(k3p*o1-k4*o2))\n", " p0' = (-1*(k1ca*p0-k2*p1))\n", " p1' = (1*(k1ca*p0-k2*p1))\n", "}\n", "\"\"\"\n", "print(run_sympy_solver(ex2))" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.7.7" } }, "nbformat": 4, "nbformat_minor": 2 }