test/notebooks/Spring.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The autoreload extension is already loaded. To reload it, use:\n",
      "  %reload_ext autoreload\n"
     ]
    }
   ],
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2\n",
    "    \n",
    "import sys\n",
    "import pylab as pl\n",
    "sys.path.insert(0, '..')\n",
    "\n",
    "import pymoca.parser\n",
    "import pymoca.backends.sympy.generator as generator\n",
    "import pylab as pl\n",
    "import sympy\n",
    "import control\n",
    "sympy.init_printing()\n",
    "%load_ext autoreload\n",
    "\n",
    "%autoreload 2\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Modelica Source Code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "modelica_src = '''\n",
    "\n",
    "model System\n",
    "    input Real u;\n",
    "    output Real x(start=1), v_x(start=1);\n",
    "    Spring spring;\n",
    "    Damper damper;\n",
    "equation\n",
    "    spring.x = x;\n",
    "    damper.v = v_x;\n",
    "    der(x) = v_x;\n",
    "    der(v_x) = spring.f + damper.f - u;\n",
    "end System;\n",
    "\n",
    "model Spring\n",
    "    Real x \"displacement\";\n",
    "    Real f \"force\";\n",
    "    parameter Real k = 2.0 \"spring constant\";\n",
    "equation\n",
    "    f = -k*x;\n",
    "end Spring;\n",
    "\n",
    "model Damper\n",
    "    Real v \"velocity\";\n",
    "    Real f \"force\";\n",
    "    parameter Real c = 0.2 \"damping constant\";\n",
    "equation\n",
    "    f = -c*v;\n",
    "end Damper;\n",
    "'''"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Generation of Abstract Syntax Tree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "odict_keys(['u', 'x', 'v_x', 'spring', 'damper'])"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ast = pymoca.parser.parse(modelica_src)\n",
    "ast.classes['System'].symbols.keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{\n",
       "  \"annotation\": [],\n",
       "  \"classes\": {\n",
       "    \"Damper\": {\n",
       "      \"annotation\": [],\n",
       "      \"classes\": {},\n",
       "      \"comment\": \"\",\n",
       "      \"encapsulated\": false,\n",
       "      \"equations\": [\n",
       "        {\n",
       "          \"comment\": \"\",\n",
       "          \"left\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"f\"\n",
       "          },\n",
       "          \"right\": {\n",
       "            \"operands\": [\n",
       "              {\n",
       "                \"operands\": [\n",
       "                  {\n",
       "                    \"child\": [],\n",
       "                    \"indices\": [],\n",
       "                    \"name\": \"c\"\n",
       "                  }\n",
       "                ],\n",
       "                \"operator\": \"-\"\n",
       "              },\n",
       "              {\n",
       "                \"child\": [],\n",
       "                \"indices\": [],\n",
       "                \"name\": \"v\"\n",
       "              }\n",
       "            ],\n",
       "            \"operator\": \"*\"\n",
       "          }\n",
       "        }\n",
       "      ],\n",
       "      \"extends\": [],\n",
       "      \"final\": false,\n",
       "      \"functions\": {},\n",
       "      \"imports\": [],\n",
       "      \"initial_equations\": [],\n",
       "      \"initial_statements\": [],\n",
       "      \"name\": \"Damper\",\n",
       "      \"partial\": false,\n",
       "      \"statements\": [],\n",
       "      \"symbols\": {\n",
       "        \"c\": {\n",
       "          \"class_modification\": {\n",
       "            \"arguments\": [\n",
       "              {\n",
       "                \"redeclare\": false,\n",
       "                \"value\": {\n",
       "                  \"component\": {\n",
       "                    \"child\": [],\n",
       "                    \"indices\": [],\n",
       "                    \"name\": \"value\"\n",
       "                  },\n",
       "                  \"modifications\": [\n",
       "                    {\n",
       "                      \"value\": 0.2\n",
       "                    }\n",
       "                  ]\n",
       "                }\n",
       "              }\n",
       "            ]\n",
       "          },\n",
       "          \"comment\": \"damping constant\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"c\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 10,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [\n",
       "            \"parameter\"\n",
       "          ],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Real\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        },\n",
       "        \"f\": {\n",
       "          \"class_modification\": null,\n",
       "          \"comment\": \"force\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"f\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 9,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Real\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        },\n",
       "        \"v\": {\n",
       "          \"class_modification\": null,\n",
       "          \"comment\": \"velocity\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"v\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 8,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Real\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        }\n",
       "      },\n",
       "      \"type\": \"model\"\n",
       "    },\n",
       "    \"Spring\": {\n",
       "      \"annotation\": [],\n",
       "      \"classes\": {},\n",
       "      \"comment\": \"\",\n",
       "      \"encapsulated\": false,\n",
       "      \"equations\": [\n",
       "        {\n",
       "          \"comment\": \"\",\n",
       "          \"left\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"f\"\n",
       "          },\n",
       "          \"right\": {\n",
       "            \"operands\": [\n",
       "              {\n",
       "                \"operands\": [\n",
       "                  {\n",
       "                    \"child\": [],\n",
       "                    \"indices\": [],\n",
       "                    \"name\": \"k\"\n",
       "                  }\n",
       "                ],\n",
       "                \"operator\": \"-\"\n",
       "              },\n",
       "              {\n",
       "                \"child\": [],\n",
       "                \"indices\": [],\n",
       "                \"name\": \"x\"\n",
       "              }\n",
       "            ],\n",
       "            \"operator\": \"*\"\n",
       "          }\n",
       "        }\n",
       "      ],\n",
       "      \"extends\": [],\n",
       "      \"final\": false,\n",
       "      \"functions\": {},\n",
       "      \"imports\": [],\n",
       "      \"initial_equations\": [],\n",
       "      \"initial_statements\": [],\n",
       "      \"name\": \"Spring\",\n",
       "      \"partial\": false,\n",
       "      \"statements\": [],\n",
       "      \"symbols\": {\n",
       "        \"f\": {\n",
       "          \"class_modification\": null,\n",
       "          \"comment\": \"force\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"f\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 6,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Real\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        },\n",
       "        \"k\": {\n",
       "          \"class_modification\": {\n",
       "            \"arguments\": [\n",
       "              {\n",
       "                \"redeclare\": false,\n",
       "                \"value\": {\n",
       "                  \"component\": {\n",
       "                    \"child\": [],\n",
       "                    \"indices\": [],\n",
       "                    \"name\": \"value\"\n",
       "                  },\n",
       "                  \"modifications\": [\n",
       "                    {\n",
       "                      \"value\": 2.0\n",
       "                    }\n",
       "                  ]\n",
       "                }\n",
       "              }\n",
       "            ]\n",
       "          },\n",
       "          \"comment\": \"spring constant\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"k\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 7,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [\n",
       "            \"parameter\"\n",
       "          ],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Real\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        },\n",
       "        \"x\": {\n",
       "          \"class_modification\": null,\n",
       "          \"comment\": \"displacement\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"x\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 5,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Real\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        }\n",
       "      },\n",
       "      \"type\": \"model\"\n",
       "    },\n",
       "    \"System\": {\n",
       "      \"annotation\": [],\n",
       "      \"classes\": {},\n",
       "      \"comment\": \"\",\n",
       "      \"encapsulated\": false,\n",
       "      \"equations\": [\n",
       "        {\n",
       "          \"comment\": \"\",\n",
       "          \"left\": {\n",
       "            \"child\": [\n",
       "              {\n",
       "                \"child\": [],\n",
       "                \"indices\": [],\n",
       "                \"name\": \"x\"\n",
       "              }\n",
       "            ],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"spring\"\n",
       "          },\n",
       "          \"right\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"x\"\n",
       "          }\n",
       "        },\n",
       "        {\n",
       "          \"comment\": \"\",\n",
       "          \"left\": {\n",
       "            \"child\": [\n",
       "              {\n",
       "                \"child\": [],\n",
       "                \"indices\": [],\n",
       "                \"name\": \"v\"\n",
       "              }\n",
       "            ],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"damper\"\n",
       "          },\n",
       "          \"right\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"v_x\"\n",
       "          }\n",
       "        },\n",
       "        {\n",
       "          \"comment\": \"\",\n",
       "          \"left\": {\n",
       "            \"operands\": [\n",
       "              {\n",
       "                \"child\": [],\n",
       "                \"indices\": [],\n",
       "                \"name\": \"x\"\n",
       "              }\n",
       "            ],\n",
       "            \"operator\": \"der\"\n",
       "          },\n",
       "          \"right\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"v_x\"\n",
       "          }\n",
       "        },\n",
       "        {\n",
       "          \"comment\": \"\",\n",
       "          \"left\": {\n",
       "            \"operands\": [\n",
       "              {\n",
       "                \"child\": [],\n",
       "                \"indices\": [],\n",
       "                \"name\": \"v_x\"\n",
       "              }\n",
       "            ],\n",
       "            \"operator\": \"der\"\n",
       "          },\n",
       "          \"right\": {\n",
       "            \"operands\": [\n",
       "              {\n",
       "                \"operands\": [\n",
       "                  {\n",
       "                    \"child\": [\n",
       "                      {\n",
       "                        \"child\": [],\n",
       "                        \"indices\": [],\n",
       "                        \"name\": \"f\"\n",
       "                      }\n",
       "                    ],\n",
       "                    \"indices\": [],\n",
       "                    \"name\": \"spring\"\n",
       "                  },\n",
       "                  {\n",
       "                    \"child\": [\n",
       "                      {\n",
       "                        \"child\": [],\n",
       "                        \"indices\": [],\n",
       "                        \"name\": \"f\"\n",
       "                      }\n",
       "                    ],\n",
       "                    \"indices\": [],\n",
       "                    \"name\": \"damper\"\n",
       "                  }\n",
       "                ],\n",
       "                \"operator\": \"+\"\n",
       "              },\n",
       "              {\n",
       "                \"child\": [],\n",
       "                \"indices\": [],\n",
       "                \"name\": \"u\"\n",
       "              }\n",
       "            ],\n",
       "            \"operator\": \"-\"\n",
       "          }\n",
       "        }\n",
       "      ],\n",
       "      \"extends\": [],\n",
       "      \"final\": false,\n",
       "      \"functions\": {},\n",
       "      \"imports\": [],\n",
       "      \"initial_equations\": [],\n",
       "      \"initial_statements\": [],\n",
       "      \"name\": \"System\",\n",
       "      \"partial\": false,\n",
       "      \"statements\": [],\n",
       "      \"symbols\": {\n",
       "        \"damper\": {\n",
       "          \"class_modification\": null,\n",
       "          \"comment\": \"\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"damper\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 4,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Damper\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        },\n",
       "        \"spring\": {\n",
       "          \"class_modification\": null,\n",
       "          \"comment\": \"\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"spring\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 3,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Spring\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        },\n",
       "        \"u\": {\n",
       "          \"class_modification\": null,\n",
       "          \"comment\": \"\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"u\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 0,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [\n",
       "            \"input\"\n",
       "          ],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Real\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        },\n",
       "        \"v_x\": {\n",
       "          \"class_modification\": {\n",
       "            \"arguments\": [\n",
       "              {\n",
       "                \"redeclare\": false,\n",
       "                \"value\": {\n",
       "                  \"component\": {\n",
       "                    \"child\": [],\n",
       "                    \"indices\": [],\n",
       "                    \"name\": \"start\"\n",
       "                  },\n",
       "                  \"modifications\": [\n",
       "                    {\n",
       "                      \"value\": 1.0\n",
       "                    }\n",
       "                  ]\n",
       "                }\n",
       "              }\n",
       "            ]\n",
       "          },\n",
       "          \"comment\": \"\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"v_x\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 2,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [\n",
       "            \"output\"\n",
       "          ],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Real\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        },\n",
       "        \"x\": {\n",
       "          \"class_modification\": {\n",
       "            \"arguments\": [\n",
       "              {\n",
       "                \"redeclare\": false,\n",
       "                \"value\": {\n",
       "                  \"component\": {\n",
       "                    \"child\": [],\n",
       "                    \"indices\": [],\n",
       "                    \"name\": \"start\"\n",
       "                  },\n",
       "                  \"modifications\": [\n",
       "                    {\n",
       "                      \"value\": 1.0\n",
       "                    }\n",
       "                  ]\n",
       "                }\n",
       "              }\n",
       "            ]\n",
       "          },\n",
       "          \"comment\": \"\",\n",
       "          \"dimensions\": [\n",
       "            {\n",
       "              \"value\": 1\n",
       "            }\n",
       "          ],\n",
       "          \"final\": false,\n",
       "          \"fixed\": {\n",
       "            \"value\": false\n",
       "          },\n",
       "          \"id\": 0,\n",
       "          \"inner\": false,\n",
       "          \"max\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"min\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"name\": \"x\",\n",
       "          \"nominal\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"order\": 1,\n",
       "          \"outer\": false,\n",
       "          \"prefixes\": [\n",
       "            \"output\"\n",
       "          ],\n",
       "          \"redeclare\": false,\n",
       "          \"start\": {\n",
       "            \"value\": 0\n",
       "          },\n",
       "          \"type\": {\n",
       "            \"child\": [],\n",
       "            \"indices\": [],\n",
       "            \"name\": \"Real\"\n",
       "          },\n",
       "          \"value\": {\n",
       "            \"value\": null\n",
       "          },\n",
       "          \"visibility\": \"private\"\n",
       "        }\n",
       "      },\n",
       "      \"type\": \"model\"\n",
       "    }\n",
       "  },\n",
       "  \"comment\": \"\",\n",
       "  \"encapsulated\": false,\n",
       "  \"equations\": [],\n",
       "  \"extends\": [],\n",
       "  \"final\": false,\n",
       "  \"functions\": {},\n",
       "  \"imports\": [],\n",
       "  \"initial_equations\": [],\n",
       "  \"initial_statements\": [],\n",
       "  \"name\": null,\n",
       "  \"partial\": false,\n",
       "  \"statements\": [],\n",
       "  \"symbols\": {},\n",
       "  \"type\": \"\"\n",
       "}"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ast"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Generation of Sympy Model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "# do not edit, generated by pymoca\n",
      "\n",
      "from __future__ import print_function, division\n",
      "import sympy\n",
      "import sympy.physics.mechanics as mech\n",
      "from pymoca.backends.sympy.runtime import OdeModel\n",
      "from sympy import sin, cos, tan\n",
      "\n",
      "\n",
      "class System(OdeModel):\n",
      "\n",
      "    def __init__(self):\n",
      "\n",
      "        super(System, self).__init__()\n",
      "\n",
      "        # states\n",
      "        x, v_x = mech.dynamicsymbols('x, v_x')\n",
      "        self.x = sympy.Matrix([x, v_x])\n",
      "        self.x0 = {\n",
      "            x : 1.0,\n",
      "            v_x : 1.0,\n",
      "            }\n",
      "\n",
      "        # variables\n",
      "        spring__x, spring__f, damper__v, damper__f = mech.dynamicsymbols('spring.x, spring.f, damper.v, damper.f')\n",
      "        self.v = sympy.Matrix([spring__x, spring__f, damper__v, damper__f])\n",
      "\n",
      "        # constants\n",
      "        self.c = sympy.Matrix([])\n",
      "        self.c0 = {\n",
      "            }\n",
      "\n",
      "        # parameters\n",
      "        spring__k, damper__c = sympy.symbols('spring.k, damper.c')\n",
      "        self.p = sympy.Matrix([spring__k, damper__c])\n",
      "        self.p0 = {\n",
      "            spring__k : 2.0,\n",
      "            damper__c : 0.2,\n",
      "            }\n",
      "\n",
      "        # inputs\n",
      "        u = mech.dynamicsymbols('u')\n",
      "        self.u = sympy.Matrix([u])\n",
      "        self.u0 = {\n",
      "            u : 0,\n",
      "            }\n",
      "\n",
      "        # outputs\n",
      "        x, v_x = mech.dynamicsymbols('x, v_x')\n",
      "        self.y = sympy.Matrix([x, v_x])\n",
      "\n",
      "        # equations\n",
      "        self.eqs = [\n",
      "            spring__f - (- spring__k * spring__x),\n",
      "            damper__f - (- damper__c * damper__v),\n",
      "            spring__x - (x),\n",
      "            damper__v - (v_x),\n",
      "            (x).diff(self.t) - (v_x),\n",
      "            (v_x).diff(self.t) - (spring__f + damper__f - u),\n",
      "            ]\n",
      "\n",
      "        self.compute_fg()\n"
     ]
    }
   ],
   "source": [
    "src_code = generator.generate(ast, 'System')\n",
    "print(src_code)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'t': t, 'x': Matrix([\n",
       "[  x(t)],\n",
       "[v_x(t)]]), 'u': Matrix([[u(t)]]), 'y': Matrix([\n",
       "[  x(t)],\n",
       "[v_x(t)]]), 'p': Matrix([\n",
       "[spring.k],\n",
       "[damper.c]]), 'c': Matrix(0, 0, []), 'v': Matrix([\n",
       "[spring.x(t)],\n",
       "[spring.f(t)],\n",
       "[damper.v(t)],\n",
       "[damper.f(t)]]), 'x0': {x(t): 1.0, v_x(t): 1.0}, 'u0': {u(t): 0}, 'p0': {spring.k: 2.0, damper.c: 0.2}, 'c0': {}, 'eqs': [spring.k*spring.x(t) + spring.f(t), damper.c*damper.v(t) + damper.f(t), spring.x(t) - x(t), damper.v(t) - v_x(t), -v_x(t) + Derivative(x(t), t), -damper.f(t) - spring.f(t) + u(t) + Derivative(v_x(t), t)], 'f': Matrix([\n",
       "[                                 v_x(t)],\n",
       "[-damper.c*v_x(t) - spring.k*x(t) - u(t)]]), 'g': Matrix([\n",
       "[  x(t)],\n",
       "[v_x(t)]])}"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "exec(src_code)\n",
    "model = System()\n",
    "model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAFMAAAA0BAMAAAAEbRnRAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMA74lUMhCZzXbdqyK7ZkQR6/c4AAAACXBIWXMAAA7EAAAOxAGVKw4bAAACxklEQVRIDYWWMWgTURjH/730Lolt0kC7qpnEMZs6ecU6Sq+DpeqQgKIupQEXwcGg0JYuFhUUXDJJncwkCCouLkpoHJxc3MSCNNh2ajV+37289753ubQf3Lvvve/Hy7vcu98dIGMplL1kfrL3xw5VbepmXq9XwuTMRTt6xqZu5s9cKGFKjp2SHTcfddD51TuiPNbhznc94qD+09E9XaDzuzg/1uoPOWhQo66Ne4D/An49DR2p5JqWLNSATBc4nYZWW+PfIsPmW0BAU86loSPh7CVD4i2l+TKQ/YE4nLUWN84u8ujC2ocHDSwDW5ufG/BqMQkHVUPwy4Xn+Tra1D1PR4EWwZGGepG3l2nhCZVv0ZFZZ9BBJ3oq4CP7kUon6HjMSJcaPjt3S40BEw3KCPX/MZKKBrsKng3pTAsodhGhOLgAqvpdZENkokfIRXxZXj0TDrmsTB3VCNuVTbwHaJMFzQVqyjQJhbtW+gvPAVdXFxY7wC9grE1n+p04JOrfXqm8eXZDFehOtVTGd41DolfC2Qh2w/J24UjbLjcxDeyrOrdf4jR1E05iDcWSRdXWHu+vQy7A30f7cq55zbCHPDBTo7udoFE2aCKRl7X1cLlT3AgThOlKdGkoFeMSHWoWNbFEh5plED3ELAzbWV2zqIloF9Am6BvGoAmzaFQYxqAJs2hUGMagrlk0KQ1jUNcsGpWGMahrFo1KwxhUm6Vw/zheswLikIYxqC5ipIms6UjDDKK5PX4G+yEN00eNWegB/4uyJh3DDM6KVy5qDGNRbRZs/27h5453wFNLw1hUmYXK+ZfU3A1Cah3DWFSZhcqFHWrmrjPpGMaiXmwWKvvr1BTUW0kaRqMJs2B+hSd1DKPRhFnwabzGqDSMRhNm+doNDkJmhWE0mjQLYxzCMITyS37ALAoUD0z8ko8/HY4wCwmcPx3iOMIsCvoPQ/zGx12wH7gAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$\\left[\\begin{matrix}\\frac{d}{d t} x{\\left (t \\right )}\\\\\\frac{d}{d t} \\operatorname{v_{x}}{\\left (t \\right )}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡d        ⎤\n",
       "⎢──(x(t)) ⎥\n",
       "⎢dt       ⎥\n",
       "⎢         ⎥\n",
       "⎢d        ⎥\n",
       "⎢──(vₓ(t))⎥\n",
       "⎣dt       ⎦"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t = sympy.Symbol('t')\n",
    "model.x.diff(t)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABkcAAAAyCAMAAADRNPwcAAAAM1BMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADxgEwMAAAAEHRSTlMAEImZRO/dMlQiu6vNZnYgvdTpGQAAAAlwSFlzAAAOxAAADsQBlSsOGwAAEdRJREFUeAHtXemCsyoMdWvdqtf3f9qbhC0qOgh0bOeLPyoihMNJ2JEWhb7KqjZOuQsDwoAwIAwIA6cM9EO3fd88nu3WT56FAWFAGBAGhAE/A9NjKFdvukVGIytC5EEYEAaEAWHglIH++eQNSb8Mp8HlpTAgDPwZBtrhId3GP6NNl5Eb9Novo0u/GBfeqrAX4hQGhIG/x8Cym9j+e3n8F3P0+3qtF7Yesrz+RdIlz8LAv8lAI/3GP6n4G/TaLrOjcjU4cd7iEgaEgT/IQPX8g5mSLBU36LVkK+ultCNihMLAv8PAS5ZH/qSyb9Arbzu4+0/yK5kSBoQBYqCvq7n9/Wl0Yf/NDNylV952cPebsyvihQFh4D4GmkdfFIMsj9yngfekfJteedvB3e/JpkgVBoSB+xkoHxOAqGVbzf2qyIrgPr3ytoO7s2ZPhAkDwsAHMVDTSORZfRAkgZKBgfv0ytsO7s6QKREhDAgDH8nAA783LpfmI8EJqGgG7tMrbzu4OzorElEYEAY+m4FywZFIt3w2SkF3lYEb9crbDu6+mgUJLwwIA1/CQEkbtYZX0cFqu1x/hoEb9crbDu7+M9xKRoQBYWDDwAs+Pp6eYyELJBtivvzxPr3ytoO7gdCybd5z3FaPm0XUNX1yhygS59t4K4pIRIbtDAKcqEhXchYi0z2J9hGQokBERQImynGu2n6cQxdIYtM5If09rz4HaBKS2Mg/6jVWsNXWkQDednA3ROxe2b5TasbZnSTcO2dRjJcakl89yTIWZ0berPqUIxaRFZMswEoiR/tgh+qsX/Enrv3MCEwyI7cp47m+l8PDfwxIXkiBnKyxQQPPMxBaKKIibVMOeE5J5yuK7AEHP1vVe4wqhe+DrCjvZMGHAnjbwd2Ybp+rHWkfZeM2h7xU01E+8V5e3MOeC5Li9fQ3Gmc23rbwohEZQckCjCB1b2m9du23f1ppPzMCk9rohrjGa38/OC4iL6QwTnbgokBERdol/bNHWjrfUGT9HIRY1TuMKo1vf17IN1nwsQB2pla/qRbUss0JrNBXNe9sTfpQH71dpArq0pqUfvEky3ic2Xgzudb3eES5BGwAhT1y7SdnISxJfyh/hXYrJAM0CkRUJJPihXtaOl9RZC+wsQn6BqNK43uDjz8mCz4RwP66anJjBko9W3048j85oWEIyNdf05YPntOf3L94kmU8zmy8bdiIR6QFJQvYAAp75Nq/B4HCeVCh3QnJEhgFIiqSTTLckZbOVxTZcDI2Id9hVGl8bwDyx2TBJwJebhf5czN/nK0+5DVJbxIxX9O+Qhf7kJGDUSQnK5M7AWc23tZZSUCkBCULWOOBic+W/XnN9qV9ZtrPjsAmEuDwV2i5IYVxsoEbBSIq0ibhkMfEdL6iyIbw4A3zBqNK5NsLkzyTBZ8J6BYztVwtbs27rKrZnAc6zTNNPrWvR93O8zgXzTwPUPs3r0c11ZWZRezrqYIVQvCdu0pPpCH+dnw8XF0y07RWMw7LAJLgqtjh1bAV8fBiJ1mCs54fsJmsgQR7BakiSOAz9AgSd5qVkIca8gShygYyxETDmxnrwOm14FFD4/Lc/jPcKU4mae2M4O2Ith20c0SWul08i/BcgA22drRzN82vla6NhvsBJ0WdWWDEZqxqMBg3kbnS/jkCm4U1AvUEGptAY06/TqvlXKNunTFQjK6e587hUH2Qbh5BSAkG0jfKqs8h+ZBAlhM5ARVBqWuX5aXb4VMQnJZU5P78oG+s2XB0W+k5i+zPJbY4ZXGLTT0nWxW2kulGpS24xzotJhv+zKFvsME4PUZawqibj+rhBgblA6vViub+OhzM1FTBPwcor/2C2wUnmo56YJzyqZoG/Kv3Bpl41F1pWyfMDN/3MOhmy35NO7mV9vLh3BiPX+wkS3LSOKoqBoRUEqQOIVXFiE1CBRVIg8B72CIPoV5d8XITaD0Geam0HpSzfbpnODmulTuCt+aYtg20U0Scuk08i/BUgA21dpRUE+MYkuvaaPjFzQJ7ITTMf04Fb7OZ9k8R8CysQcATjloHsDOnX6fVCswOgtTGGBAIWUDP/lwHjKQr5r4BwyzBJqbnXDxR66eQUND+ysAJHYRkWpFzEJyWVOT7vDCfKLPh6JgscmYssiElNkKVyVaVyajqgiqvmeaGIixyy7x9DjaYlR6jLKEYqTF4URdeA1B1rNp31D2wG0gHuw1UFdOKSkutj6pHOlpXedBHTbgT9GWmrkx2WE1SmI6/PWy024Y2sVZ3dpKldkIF13aqLlAbBCbQQg/tBdZt4O6xLSzrlkKBX6P7fvCWoI7YVELrjzlju6LJE36icEbwdkLbBlowok08k6OoLFHHgJoFrmuj4YHaEW0W2EEYUZn0Y1PlvYjgLLjYytVjCWtLpl/UvdIq6ByNyRkDKrojg+WLoOAzl0UJKCvQeAXv6WvuCEgZOIGeTUEDYpW9UBCpyFVqB785zQaTyFlkQ0rsdftOt6o8RgXVD1WHg+q0qYrJLCCDLQfVkF6txhlMpCV0UPnDoAJKl77UOS123287VzU1G2D8cNGurpb+113VLfSfiu1ST3DhJJfyNcLgztsRM+oxyyMwvGchrbN8Pe2F/TZ2kqV24jhEnzWnKg1MBfJAdQdUE+MD5lcm8DChrGwlQD/SyMnzXW8gTisUHRG8ndG2gRaMaBPPQgwWYGOA47kMM3Xwma6thnU7oswCxyA19jSO25FQBFvtA7OPGofLTr//ad0rP2YM2J2ArgIWSdXZwUe4ajs4QUN/UZEFRygkFGGudE6gkjXdThQaCiIZORSP/aWzlWo2W6VlLLJBJfaMxS02leV0q7puVD4koFb6Ekt9jhVqDFpt+rZXK5aDOIOJtISW6vQaxx3qatSftqvxSPkaeiiR+FK1l/SStyMFjlAm1/fbrax52hF32Ki/HTFQzJ2dZKmdhFbNjtXUYCst2JET29Bs59CUNBJgBFOtd9yOXMIZwdsZbRto2r4CEGFFnilLMAkEk/k4cal7BzQaNRrW7QjVyWQWPdiJmtyz/PJeRHgWXHTlaoZF7Sy0+mXnDKp+hPYYsVVraV515h25J3wLi02dupSFgzsGUjonOGDCdlFfV0AkITcJeu85zQYSyFhkg0pshCqTrSqXUZFVqzo1IhtebVrPCIOJsgSzMFJTW4HJw3oD3lQ7MmBhhCyW5VE7oscjtpCaWgZl0MXbET2Eb0xrCcsYJtjJnZ1kyZzQyaSuMLWDDQyZoPmlDRSIiH0MY+fQKAU9ajCpwReSdDpLg5/vV6qygncxOCN4a62WPXvRDDQFNRyRiZchSz32aGhD+Go8Qiti0LFQ81qsHYHF68r2RxRupv3wLKio5hemrmDFkMbEVr9u3I9NHMJUasYFEECMYyhcUDEXNr+9rbt1kYWXEZAycAJgaaOIRncBRBJyQ4b/ntFsIAFWTpnTaOlSkQ0rsZdVmW5V2YyKesJqeeRyNvzKtL4xBhNjCZVtPty+XzXqV+0I1cfd0sL6wsF4RHWtaNUSSjv0XFnxpeywmsQshlXQAHRU4bB1dpv3nUNtpaUTSpkT+jzY7VY9dKhF0DAo8QmW+9XCA03ta6eRqvcA6GlImGbXeOGUM1hO0ZeedriGM4K3E9qgkjPQENQFRDZecpY66h242Uqla6PhXTvSsU62YZJp/0IWTGy6qz0SNGFl9WsUjPqHDgCEIzXj2KTSQxM3RKYVMxOwK2FVDMLDcsklVgkK/GTgBNdo1CwxCQ3lJRW5ycHBPZ/ZQAKsnDJnXJENKrFXCggRkG5VOIeqrS9NNVSQoHrFefhQYzhQ4so7FtV1S9BzMZj6bLtrPQ20R+or05gSGhvI4r4dwb95LlU1raYScJVF7d8CgmmwAIJVPJVBvc0XNoBhkwNX7WZgVpsGVHDzy06yZE7VwVQz4q+KKhMcULQ4RQebVQEb7vt13VAFiTaf4Qq8uuzUUgnbkEx6ZjvyFqfLlA3JHBG8ndC2miyEkZZq7raIVPJr6vJlqaOtiDCzCfvdnK61hrVimVk0r6ZtIRwWL4/2L2SBsQoGhCJnsherX6dVGNviRgnlAQ0/7A8ucTyidyIqJNSaQZcJNndB8zJABaA2V/ghOfQrHOohAyeYFxrGK4l+EOod02wyck9muNcls9mh45LQzcopc14sskoRQSX2vIBs0cFzilXFGJUHgvHCMx+aR0012IExnJqkkbO5hxiMisKsDD2uWwKWKH2xc1HKuprg+xHcANaP9VyV82tux2UZ225A3+q54Mj0NcJa9qj7oCUFLGjOkepjvSUZ47mZaT2N1Q8T7u+By+1+xP638tv/spMsmVMJmylWV9P3BQt83IKNB9Q79azSYDNnKoGqmibViGG40k6sTbTJB/1Aq8p3h/OJ+3QPrwjejmjDNBw0eDhCpMCsqHPxUrPUgB1MFWrY6dpomAyi5GYBoXC9DwqE1eVK+xeyoHKlf+EboWkmm4IRt9Yv0+pIelYe/ThhQDgbsjLfiBE11MsDk4BdF+M8w4Z9bX8HkM7UnMwJlB5oZCf4fkRb4QEIlXun2WTkK1L3D9fMZotuK4+VU+ZUWQ0vsuEl9ocCsoUHNVuKVcUY1R6C8SnHau7akbbwHRnDmUkaOZt7kMGoOM7K8Pm6Jdhu46qDtMFz9KjXXo9ee/3Nt/X6pf1G0hv4sqepOy5HhAhlTWMYHfUIZ47//onhDWpmah5tzkKYy5ilIMwD1OFl38HHGf4rIgtrQVf0y5dH1lL40wGkIDVn4QTBHIDgOPfug0hByM+P4Q0zmz2iKJ8rKt0ksMYZxyKKvAIhyKoOVLOBf/B4FDlEsad6PRJ8gAO91wwfC+Bn/HL3iWj2KqgcsfDoNGd9ae/VZ+aboBGPZuY+ImoxwMydW9o5wumm4WLSUHFieItiLmOWQjCb78cnWnn38HNEqieo3ytQvzhn6Ta2+UVp3wNIQWrOwgniOABxCvwoUhByvT9eHbq9SybMbHbR4jwCVeoTvsZ5RIgv5tovEEK4VUXp02A6ihymWFoU9Ov1SLBJ2HNfM3wsgLcd3O0R6fFaD4Y8AXxeNC1pXlw9N97EO7irFb2Dl+fe0wNm+OplMMsmavrUxLE49aqO8Y+7R/GmJ3RNihaR8djfs2YpBHP5VOzhsMR/JSo/UL/4laHeTOaHwX29kMLUnIcTBOMFwVH63N5IYcj1l8WbHfEqkUCz8SGK8AtUqUfyDmcci9jjOLRXnuwVq/Kqhgs7c/sjhylW7Yv36vW6le0Y9iODvPC2g7vPsmneNbhgsp5rMa/O7qv/QlkBO4sV9K6DieejKZUgATyQH6de1eEBr7ojeYMJXd7Nj2HOLyAkS4GY6Xyxym5h2DPjR7AP5/cJ1W8H63Zsz4RfmPH1QvpNThCIF4RBeHT3RgpBDgLVVM56R/yldI4CX/QPVWmQWC8hP8YMhXDFquKQaKj+yGGKPdWrX/CPBLkAhwJ428HdLmp+FzuF5C/+r25+wqzEZOaSBVgosY77EeyQfwSkKBBRkVT+1VQOfcaxI2TnkZDOTtZbPT4HaBKShMjnek0QrPR2JIC3Hdz9Vm2LcGFAGLiNAXMMLzt0+zYsknA+Bu7UK287uDtf7kSSMCAMfBAD7Bjeg2n0DwIrUIIZuFWvvO3g7mD0ElAYEAa+iAF2DC87WuaLMiBQvQzcq1fednC3F6p4CgPCwJczwI7hLQKXR748x/8G/Hv1ytsO7v43uJdcCgP/GgPsGN7AL2z+NYa+M7/36pW3Hdz9nVwKamFAGDhlgJ+96w7dPo0iL7+AgZv1ytsO7v4C5gSiMCAMXGVAfXNHJ2cX7jDrq1Ik/KcxcLNeV20H++unT6NJ8AgDwkAOBtjZu+7o6ByCRcatDNyrV3Z4Nfw1gDnq+1ZCJHFhQBh4GwPs7F13mPXbUhPBv8XAvXpt2d8Gmj9s+K2cSzrCgDAgDAgD38/AbP8NEfLSy8TW92tUciAMCAPCwK8yUNr/I6dkZ/dv4b8KQxITBoQBYUAY+E4Gyif9baoDPy3j2d/9uYDiEgaEAWFAGBAGiu7x2rYa5Xj417ZCmDAgDAgDwoAwsGKgf9kj7f8HnHLOVVLJJQIAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$\\left \\{ \\operatorname{damper.f}{\\left (t \\right )} : - damper.c \\operatorname{v_{x}}{\\left (t \\right )}, \\quad \\operatorname{damper.v}{\\left (t \\right )} : \\operatorname{v_{x}}{\\left (t \\right )}, \\quad \\operatorname{spring.f}{\\left (t \\right )} : - spring.k x{\\left (t \\right )}, \\quad \\operatorname{spring.x}{\\left (t \\right )} : x{\\left (t \\right )}, \\quad \\frac{d}{d t} \\operatorname{v_{x}}{\\left (t \\right )} : - damper.c \\operatorname{v_{x}}{\\left (t \\right )} - spring.k x{\\left (t \\right )} - u{\\left (t \\right )}, \\quad \\frac{d}{d t} x{\\left (t \\right )} : \\operatorname{v_{x}}{\\left (t \\right )}\\right \\}$$"
      ],
      "text/plain": [
       "⎧                                                                             \n",
       "⎨damper.f(t): -damper.c⋅vₓ(t), damper.v(t): vₓ(t), spring.f(t): -spring.k⋅x(t)\n",
       "⎩                                                                             \n",
       "\n",
       "                     d                                                  d     \n",
       ", spring.x(t): x(t), ──(vₓ(t)): -damper.c⋅vₓ(t) - spring.k⋅x(t) - u(t), ──(x(t\n",
       "                     dt                                                 dt    \n",
       "\n",
       "         ⎫\n",
       ")): vₓ(t)⎬\n",
       "         ⎭"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sol = sympy.solve(model.eqs, list(model.v) + list(model.x.diff(t)))\n",
    "sol"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAAyCAMAAABCrGH/AAAAM1BMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADxgEwMAAAAEHRSTlMA74lUMhBmzavddiJEmbsg4DnoagAAAAlwSFlzAAAOxAAADsQBlSsOGwAABqlJREFUeAHtXNuWqyoQRBQliG7+/2tPdTe3ZGLGRCdL19GHGTQ0FGXT3CpRTcCl1XX9IQOaOG5U03ZdZ/6wnqtoZUBxC6q7s3Ix9hl5P+bkURPdiake54pVe3iuz0z1IOyaG/03Q8X7IZMnprqfhNE28H/tDklwAXViqtmZ0ZJJ/Nk0pVWHTJ2X6vEWCb3FiergD8lwBnV0qvshNJhn2HBrM2ZJOI4f3s5hthw7dAwoD/kOc3t0qpVqLJH1c9Cb40wvhmql+p95DkMzATk+1S5geVXNoBN9yc1jqFaqTREl5TjY/+NTbQLc98m+QRNDcwrVqpOZyMEILnCOT7Wy8NZlqk1Io+FFdXmtH6Z88Lzs9kNolQ5piRgDiCdf5h0cfwWQDxkuZjcb5xaDU1PesInDosZg2PLDa1gslH2a6uHNfJlmKlO+OLebrDKyTpyeRJlPq/wLuxPEamVyZOhpNhKvGC/GuY+zviE7fMpyrP9noDozZiZXzZ3Twjx+nFePOfvBEqeiejZqLkvCtN0UGb22m/Zzrb7BZGMKc44Tc71FfW2i7sf0z5LujgbueP+Z9wBPThVAHvmqluvXgdcjOf/n+1N79ble3EX1197XRfVF9dcY+FpFh/LqA80oCpT35zbFVtXGv1HdzU1Zn/2xA9zNk7dLaGzabl2EbeYmb6/cZ6qhvIuktlWV8W9UK5W21e6x/MXdzhIaG7ehXkEdFvxIoHwm5llqxq9U+2oz7RXo7Z+lPY1vSmgW/ChC+QjJYjN+pVovdLHt1D6WkHbqviihWfKjCOUjJIvN+JXqpS72SNTm+7wJms5lvyChWfCjBOUTJMlW/TB+RfU4addJF0Nycg325b224+gcRC5eu5lOUL2exw5PaNPeaOcm7Lchl/G63tbEB4625HaV0KDQHoUWBKVe4yY6mylwcYOTscm5thotyY9aZwmYQRNG7xgzq3nWi3liXSPxo5aVQC+o9g1ImzlUc3KgA1Ot5hltMMGBZoiz6YmlE1WNFnjqPCPOnZBraNWQVXQj5RhkW38/CQ25Pe1fFwSlXq081T4luLzfShjHUEZCAz9yoyf1gwHq/ubUjUDGY8u1Yp5JiHAsjlhWAi1TbUjAJdrDmATQrhUwI3+jo0fhIygl/RHSYwNCDQ5aKRee+byz3NCxn5Vzwd0kNCM1DV93KAhKvXhGgbbAJSgtu001EuKBM8pQOzU8UuPzlt5JPIxfKebBNJpzztxdom0SbaLWNNgVqs1wyxcd000MjCNOTJI3w58paghqmrkCIoMHTtv06NF4kHIhI19iH28+kNA8QpOSTGgmwlIQ/Ivo5FkFl3wCkjTiogvs4HSPJmYPRyFqYK6QiGqeFG1rMc8TKDBtKFjxn2SbQ3UxLlQj893FwITYmGQnlV418auSwmUYnxsVWF3HheS+x3dsnwvfTULj5yCqkIwALpCqEWeIDyyFhY4ChXLJyZC+WUfRMF0hpYXq9WIerr8LzE98TU+MF6k2HCIYepWMEUXemYeD4JXyMG7gINXXxHLfo2aIfWqQ2klCgxiBsYx9NCMo3ZaiGnmqAAHJBiGO/Lk6nCQ6xiyOwqtgrpBHgsB6MQ87noTqFHyeGL+gmkLrPCB60egRk+gmFHblnaEZ1F6eDvYYfvkj+A0yxCRSdLE7IerIHZwp9ttNEhoZdNmLM4JSL2EmpAyEPFxHB69CNYWVlK81irlC7E7D4noxz0zdGaGagueyEmiRagUeMCZbAlwlxQsk9g2aW0OduKOvibFywNBkL/uKCdR1IYtBo7IwaR8JDR8mOu7zGUGuF12H1avyAKRhUmjIqy27tsDiTRK4Mlro8QZmMC87RaLmWS/mmcCAx14RvdtlJdAy1cY63Y2WJnVVUoQujqNyO/GsNniZTNOc0on+pcjnaCaC+nXfp0CI4nK43CKhwey+d9xREJEiglKvsoxEHoy2p4weLYqLFobFDgjM8EW00GF6HNU70Wq1mMdY7drO8sIhQniiBFqmmihacy2suNaY7iSheQdBFapfIEyLa8mSF4AvLMpH97aqGG+nesPKfScJzUoENOiViUHh5kkqbRnJR/Wq90nmh0f3tqoYb6ZaxsyH6lbd7iWhWYmAFiiKB7AV6GpVybtintq2/jrlVqrbW6BBc5frbk/9DvCr4tciaG3vKiHrqyIxPKbVDLKtRhKLrG1r461Uv0b85qcLJ0VvlrJL9gKlPrNaV3SxfevAa13hV64VDJBXXz9SsYKobVnSj1RcP72yjccV1vLTK/8BCv42bBTkmxkAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$\\left[\\begin{matrix}\\operatorname{v_{x}}{\\left (t \\right )}\\\\- damper.c \\operatorname{v_{x}}{\\left (t \\right )} - spring.k x{\\left (t \\right )} - u{\\left (t \\right )}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡                vₓ(t)                 ⎤\n",
       "⎢                                      ⎥\n",
       "⎣-damper.c⋅vₓ(t) - spring.k⋅x(t) - u(t)⎦"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.x.diff(t).subs(sol)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAe0AAAAyCAMAAABGdXOHAAAAM1BMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADxgEwMAAAAEHRSTlMA74lUMhCZZs3dIrurRHYgTxvb/QAAAAlwSFlzAAAOxAAADsQBlSsOGwAACKpJREFUeAHlXNuCpCoMxHu3t+P/f+1JANvCAQngzI67PkynMVZVEkG0cVS10VYr3sDU3/+qPzXHWSWEFMtGMmBqfssJAMGaVdO2baezAGZCVh7i2lGYTVK1j8R4Q0wGVIn5LScABGtW7ScWMD9tf5XRJlX7SEwwCUmAKiO/5QSAQCZIADMYn7uj6+uhFmTFPeqPfYPI4xrO2fDFmgSIqb6iR6I8ggBCabVfjVLde7zS/i37mspcexLBk3J3rrYv1iRAqPalfiTKIwggFFZ70CNj/UpMeqF7N/fz9uPV9saaVYyIfocoiyCEUFjtaebCNXmpL6n5mkeZlLtT3/bGmgQIfftKv0OURRBCKKz2pqu9bDSe/+x2la0LJUm5O1XbG2sSoLDaDlEWQQihrNrd1nNql229yPC37Pr5avtjzSoGZeRCv0uUQxBEKKv2aJ7KtPbhzLfU1Q96kS3/AaY1KXdu3/bHmgQo69suUQ5BEKG02rpv/yPV9sWaUww+8S7O1tEMmDapOQRBhLJqu0PGVae6e99Ftq6oknLn9m1/rEmAsr7tEuUQBBHKqq3MdKD5h2Zpp1hzisGn49XZ6iQ1iyCEUFjtaWLpQ97tEB+au11l6wIzKXdu31beWJMAZX3bJcoicKQCApkQFJgXKYNd5ja+/+GnKyTgD1TbGyukEtISMo/8Xul3iLIIQgiF1VZvfnJ6BBEK8/b2q2xdkCXl7hyWL9YkQOhYl/qRKI8ggFBa7a5e1/nHfxXp52p7zcNFWQO7knJ3rrYv1iTAT7Uj+pEojyCAUFrtQFZ/a3NS7s7V9gWVBPiptg8p0FZOAAhkQlBgBsgf3gyRxyORZCMJEFMdpzce5QSAQCYEBaZUzbP8IPK4cEk2kgAx1XF641FOAAhkQlBgStU8yw8ijwuXZCMJEFMdpzce5QSAQCYEBaZUzbP8IPK4cEk2kgAx1XF641FOAAhkQlBgStU8yw8ijwuXZCMJEFMdpzce5QSAQCYEBaZUzbP8IPK4cEk2kgAx1XF641FOAAhkQlBgStU8yw8ijwuXZCMJEFMdpzce5QSAQCYEBaZUzbP8IPK4cEk2kgAx1XF641FOAAhkQlBgStU8yw8ijwuXZCMJEFMdpzce5QSAQCYEBaZUzbP8IPK4cEk2kgAx1XF641FOAAhkQlBgStUYv1n/8nl1TDdV76v93n3tVOkFI96dOY0QebfFfrmTZAMAlRwRF/dHwignAAQyISgwIyJOuyW/ULxyCnfzOgmIXL3fkSXp4Wwca/8RUI6Ii/tPmTRf7yQAiTdV2yv51JhTuOXmdRIQOS2Dj7ziEqi2s/bfAZQimp+fg+9c3EwAEn+u2lmFq9NH/9M55n6FyJUybzG7DvgtUG1yOX6ddgCliM7ifmT82HcSgMSfq3ZW4bJG/0/OvhoQuepii+CTqy1FdBb3fxVJLYFqZxFAzNJq06KFgRYtLPU0trR+oSNr7pa65ox1a88f1DKOtM8msenXtYHJGxeuSVj5MPb12prRn8x+5df8DoqlXqdFk+6CSAeJ7Gl0BmnkwRvtWPWSC4hcRRdDJFdbiOguCTUKT38D1c4igJiF1X5RZicqV61mXgRTT2RNNNd48Vt/tVr4o6flc1TpbtOXQ3JRdl2zjqSjwq3jskUF27iXimAmfdnW5mvTTJZiJT36zftDkFr4VdORhmeQpsFGlvzSM3CI3PJcfCRX+wLL7DKI7uJ+70GBant9sdFHADHLqj1yqun/N4xUX37za9jaRr3JWqjHUCMvQ+QWzui4cS8yLwLCxIwa1k510mJ3FXvq5Y3WpLPnoOBL7kCiDkFqrHiJXN9qLyuNvGir2HumvYRwz9v6pcXQdydX71zcSQAxe6vdvd6f7cW167aq53GTBlNdwHr7j9q4xTZyn7YtXATKMI/hrenm/F315g0IbXv+nCl73avfXCdrAoU5lfge/xDUzdVAVxtq2IXsLOZ48w0i33eHP23fPkujA8qKkT+Sh7XaPUaySwAxe6v9FXSZto3rRy/4UT7VRF2k0VXVrvaNXtMycwdv9Yi9woz6Pa88zks3fbaYE8qa+rprKHqNW2m4jyAz9dH4II2+6+MtL0QeV2Kr7XEsq7bgnYs7CSBmUbVpvKSpmO6oematzx1YRE6jO58DpoXq3NEoyxdvvtLbjQs37qPB3hj+7PR7hLpqYO4UussvREGsH0Hw5iFII+Gw43eM5O7rAd4kFFa76O0BPRdTui/rW6KB58f6aqil8gSMh1zdwv28th0fLts8vu9+3vjcRnals+VFzyvAtBSmy1NJ+Sz8CLJ6uMMf0hjVTgx/z3VbOYv7WeKXrbDaDkFq3564o6564OTxvOUhDjrqsqmR5lSmhZ4Q0a2anpjPuoObp8f6QTrN3/iskGwvYhtoHkj+YBoKMzN41RrsI2jlK0jHd2AfaYa653klT99og8j56/X2bSN5/J2Lwmo7BBCzbCSn6c+q+8ayLeaelu53jlTNA0+gTcs4D+y5zGttn6fw3F1N7EH35TziS7aOjm/Hme+0wDQUKwOqptcPAD6CGF3rAGmaWtX1MOhT9aZq49p/SKUgrP38wcX9nsPuJQCJomofgpIeiMFl+0C42UoSdFO1MQRIJTYH7L3agd2+5nICQEistvRJprlts7doviDuapMKsnwQeVyBpDZJgPhzY5zeeJQTAAKZEBSYXjVmxuTd5TTWfJ9mnv077bd/kQraiSHyvSn8GcsGH5kEiKkO07p7ygkAgUwICkyX03xr3htPmeJbMw9rry/0cd8SD7GgnQQi35vCn5Fs6AOTADHVYVp3TzkBIJAJQYHpcv4t3yDyeEiSbCQBYqrj9MajnAAQyISgwJSqeZYfRB4XLslGEiCmOk5vPMoJAIFM+DfKYErVPMcP/nuzTHQsG8mAmGqJhHICQLAm/NN1MCVqnuUD/5ldJjyWjWTAX/Df6P8HEq5B57np+TIAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$\\left [ \\left[\\begin{matrix}0 & 1\\\\- spring.k & - damper.c\\end{matrix}\\right], \\quad \\left[\\begin{matrix}0\\\\-1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}1 & 0\\\\0 & 1\\end{matrix}\\right], \\quad \\left[\\begin{matrix}0\\\\0\\end{matrix}\\right]\\right ]$$"
      ],
      "text/plain": [
       "⎡⎡    0          1    ⎤  ⎡0 ⎤  ⎡1  0⎤  ⎡0⎤⎤\n",
       "⎢⎢                    ⎥, ⎢  ⎥, ⎢    ⎥, ⎢ ⎥⎥\n",
       "⎣⎣-spring.k  -damper.c⎦  ⎣-1⎦  ⎣0  1⎦  ⎣0⎦⎦"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.linearize_symbolic()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "A = [[ 0.   1. ]\n",
       " [-2.  -0.2]]\n",
       "\n",
       "B = [[ 0.]\n",
       " [-1.]]\n",
       "\n",
       "C = [[1. 0.]\n",
       " [0. 1.]]\n",
       "\n",
       "D = [[0.]\n",
       " [0.]]"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ss = control.ss(*model.linearize())\n",
    "ss"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "Input 1 to output 1:\n",
       "      -1\n",
       "---------------\n",
       "s^2 + 0.2 s + 2\n",
       "\n",
       "Input 1 to output 2:\n",
       "      -s\n",
       "---------------\n",
       "s^2 + 0.2 s + 2"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "control.ss2tf(ss)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/jgoppert/anaconda3/envs/phd/lib/python3.6/site-packages/control/freqplot.py:153: MatplotlibDeprecationWarning: pyplot.hold is deprecated.\n",
      "    Future behavior will be consistent with the long-time default:\n",
      "    plot commands add elements without first clearing the\n",
      "    Axes and/or Figure.\n",
      "  plt.hold(True);\n",
      "/home/jgoppert/anaconda3/envs/phd/lib/python3.6/site-packages/matplotlib/__init__.py:805: MatplotlibDeprecationWarning: axes.hold is deprecated. Please remove it from your matplotlibrc and/or style files.\n",
      "  mplDeprecation)\n",
      "/home/jgoppert/anaconda3/envs/phd/lib/python3.6/site-packages/matplotlib/rcsetup.py:155: MatplotlibDeprecationWarning: axes.hold is deprecated, will be removed in 3.0\n",
      "  mplDeprecation)\n",
      "/home/jgoppert/anaconda3/envs/phd/lib/python3.6/site-packages/control/freqplot.py:163: MatplotlibDeprecationWarning: pyplot.hold is deprecated.\n",
      "    Future behavior will be consistent with the long-time default:\n",
      "    plot commands add elements without first clearing the\n",
      "    Axes and/or Figure.\n",
      "  plt.hold(True);\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEOCAYAAACetPCkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAIABJREFUeJzt3XecVdW58PHfc8oUZoahDR1nQBgUUFEQxQpRE1OMsSWahNiNaabc5CWWRBOTa5J7481N7IXYEhW9sQQLNlBBRMBGl5EiKL1Pn3PO8/6x94yHccrp9fl+3J85e+9VnhmWs2bttfdeoqoYY4wxkfKkOwBjjDHZxToOY4wxUbGOwxhjTFSs4zDGGBMV6ziMMcZExToOY4wxUbGOwxhjTFSs4zDGGBMV6ziMMcZExToOY4wxUfGlO4Bk6NWrl44cObLT83V1dZSUlER9LpOlMu5E1hVrWdHmiyZ9d2njOW/tK7V1pap9RZMnknSxtqF4fnZLlizZoaoVESVW1ZzbqqurtStz5syJ6VwmS2Xciawr1rKizRdN+u7SxnPe2ldq60pV+4omTyTpYm1D8fzsgMUa4e9Y0Rx6yaGInAGcMWjQoMv/+c9/dpqutraW0tLSqM9lslTGnci6Yi0r2nzRpO8ubTznrX2ltq5Uta9o8kSSLtY2FM/PburUqUtUdWJEiSPtYbJpsxFH9tRlI474NQeC+vd5a/XPs1dpKBRKePmq1r4SmScXRhw5OcdhTL54o2YH1z+9nDXbagHweISfnFqd5qhMrrO7qozJUjc9u5Jv3rOQxkCQu6ZN4JyjhvKXl9bwzPub0x2ayXE24jAmCz2/bAt3vraWCyYN4/ozxlLk93Ly6ArW7ajlPx57l4P69OCwoeXpDtPkKJscj+JcJrPJy8Slz/TJ8Z0NIX79RgP9iz1ce2wRPo+0ndvbpPxmQQMAN51YTKFXOismKta+EpfHJsczdLPJ8eypyybHoxMIhvTrd7yhY371nK7dXtthmoVrd2rl9Fl6yytr4qornLWvxOXJhclxm+MwJovc8koNC9ft4rdnjmN4v44f9Jo0vA+nHjqAO+Z+yK665hRHaPKBdRzGZIk5q7bxl5c/4Kwjh3D2UUO6TDv99NHUNQe4dU5NiqIz+cQ6DmOywLoddVz1yDscOrAn/3nWYYh0PXcxakAZ500YxoMLNrBxV32KojT5wjoOYzJcbVOAyx9YjM8j3DltAsUF3ojy/eS0UYjA/7z4QZIjNPnG7qqK4lwms7teEpc+k+6qCoSUv73TxNIdQX4+sYgxfSPrNFrNXN3Mc+ta+NXkIkaUR5c3nLWvxOWxu6oydLO7qrKnLrurqnNNLUG9/P5FWjl9lj7wxrqI84Xb19CsE3/3on71lnkaDMb+OhJrX4nLY3dVGWOSojkQ4gf/fJsXVmzl+jPGMG1yVUzllBX5ufqLh/Dexj08vmRTYoM0ecs6DmMyTG1TgO89tIQXV2zlhjPGcPHxw+Mq76wjhzCxsjd/fH4VextaEhSlyWfWcRiTQWq21XLmLfOY+8F2fve1cVwUZ6cBICL85syx7K5vtolykxDWcRiTIZ5ftoWv3TqfPfUtPHjpJL59bGXCyh47uJxvHVPJAwvW897GPQkr1+SnjO84RKRERO4XkbtF5FvpjseYRNtZ28RPHnmHKx9awsEVJfz7Rydw3MH9El7Pz78wmoE9i/jpo+/S0BxMePkmf6Sl4xCRGSKyTUSWtTt+uoisFpEaEfmle/hs4HFVvRz4asqDNSZJQiHl/5Zs4tSbX+WZpZu56pRRzLxyMoN7FSelvvJiP//99SNYu6OOm55bmZQ6TH5I12vV7wNuAR5oPSAiXuBW4DRgE7BIRJ4GhgJL3WT2Z5LJeqrK3NXb+a/Zq1mxeR9HHtSLP55zONUDypJe93EH9+OyE4Zzz7x1fO6Q/kwZ3T/pdZrck5YRh6q+Buxqd3gSUKOqa1W1GXgEOBOnExnqpsn4S2vGdCYYUt7ZFuC8OxZw8X2LqG0K8JdvjOfxK49LSafR6udfGE31gFL+3+Pvs7O2KWX1mtyRtifHRaQKmKWq49z9c4HTVfUyd38acAwwHWd00gjMU9V/dFLeFcAVABUVFRNmzpzZad325Hjm1JUPT47XNitvfBLgpY9a2Fav9CkSvjLCz0lDfQespZFKH+0LcuObjYwo9/CLo4u6jcPaV+Ly2JPjcWxAFbAsbP884J6w/WnA32Ip254cz566cvXJ8caWgL64fIte+eBiHXnNM1o5fZaefdt8/ePDL2pLIBhxLMn0xNubtHL6LP31k0u7TWvtK3F5cuHJ8YhGHOK8ivNbwAhV/a2IHAQMVNW3ou/X2sqs4sARx2TgBlX9grt/tdux3RRFmfauqiyrK5dGHHUtyvIdQRZ+3Mjy3UJjEMoKYPIgH8cP8VHZ05tx7euRVU08vz7AJeMKOGmov9N01r4SlydvRhzA7TgT1yvd/d7Aokh7p07KrOLAEYcPWAsMBwqA94CxsZRtI47sqSubRxz1TQGdX7Nd//zCaj3r1nk6/JeztHL6LB33q1k6/fH39JWVW7W53egi09pXSyCo37r7TR11zbO6cO3OTtNZ+0pcnnwacbytqkeJyDuqeqR77D1VPSKKDi28vIeBKUA/YCtwvareKyJfAv4CeIEZqvr7KMu1EUeW1ZUtIw5VZVu9sm5viNU7Gvmozsv6fSGCCgIML/cwrp+Xcf28DPQ10LMse9pXbbPyu4UN7G1SrjmmmGFln70HxdpX4vLk04hjIc4v87fd/QrgnUh7p1RvNuLInroybcQRCoX0iedf1vk12/WBN9bptU+8r+fcNl/H/vp5rZzujChGXu3MV9z07Ep9eeUW3VPfHHFdmdq+Nu2u12N+/5JO/N2LumFH3WfOW/tKXJ58GnF8C/gGcBRwP3AucJ2qPhZDx5Y0NuKwvwgjSd8cVHY1KjsblO0NIXY0KNvqQ2yrV7bWh2gIfJqv2AfDyjxt2/ByD71ooLxnatbjSKWPa0P858IGSvzCNZOK6FX06cjD2lfi8uTCiCPi23FF5BDgFJyR+cuqmrGPno4ePVpXr17d6fm5c+cyZcqUqM9lslTGnci6Yi2rfb5gSNld38zuumZ21jWzo7aJnbXO1237mli5/hNa/CVs2dvA7voD3xDr9QhDehVT1a+Eqr49CO7ezOnHj2dk/1IG9iz6zDKt3cWcze3r7Y928+17FjKwvIiHLz+WAT2LAGtficwTSbpY21A8PzsRSUzHISJ9usqsqu0f4ksrG3Fk31+E+/fXUlRSQlMAmoJKUxAag0pjABoDSkNAaQjQ9rU+oNS3KPsaAzSpl7oWpbZZqQ90XL4A5YVCmT9En2IffYqkbetb7KFfsdC7UPCGPceQSSsApsMHu4PcvLiR8kJh+qQi+hR5srZ92YgjcgkbcYjIOkBx/v87CNjtfu4FfKSq8b/zOQnydcRx8sknE1Lnr++QOlswpIRCEFQlEAq1fQ6FlEBICYZCBEJKIOikdT6HCIaUFvdzi3uuJRiiJRhi+cpVDD94FC3BEM3BEM0B53hzwN2CIZpaQjQFQzS1BGkKhGhsCdLY4n4NBGlodj7XNwcIRfgMaoHXQ1mRj/JiP7Q0MHRAX3r38NOr2E95jwL6lhTQu6SAPj0K6FdWQL/SQnr3KMDrkaj+XeMZUXR3Plva15INu7loxlv0LingH5cdw4fvv2UjjgTlyYURR5fvqmrtGETkDuBpVX3W3f8icGpM0WWwnz76LivWN3DbqgUozm+z1n5VofVGAfez8xXVsH11vmq79O65kLo3I7iFhMLyhrQ1r3PM+cVP27HW/dZ0rZ1Da95gSNHnn03dD2v58gN2fR7B7/VQ6PdQ4PVQ4PNQ6PNQ6PNS6Hc+9ystoMjvbduK/V62bd7EoaNGUOz30qPAS0mhj5JCL8V+H2VFPkoLfZQUOp+L/J+ume38DzIpdd9vnplQ2ZsHLzuG79y7kHNuf4MfHpaeJ9xNZop0cnyJqk5od2xxpL1TqsR7qequ95vYXteC1+sl/H+T1kvcBxxDcP9zjsun56VdHgk/F5anLY045YWnDU/nkU/TesR5YZcIeOTTPIGWZgoLCtrSfppO8Ah4W4+1nse5tt96PPy8z/NpHq/H+eoTafvc1FBPz7ISfAI+j7N52s0DRCpbbsdN1PlsuFQVbtP+EH9e3EhDIMSPjyrm0L7e7jPFyS5VZf6lqkhvx50NXIfz0F4lcC0wO9Jbt1K92e242VNXpt2OG0naeM5nY/v6eHe9Tr7xGR11zbP65Dubkl5fNravaPLkwu24kb5t9gKcZzeeAJ4E+rvHjDE5bnCvYq49ppjxB/Xix4+8y3/PXk0o0okpk5MiWo9DnbunfpzkWIwxGarELzx06TH86sll3DKnhg+27ud/vjGeksJ0Lelj0inSOY45uHPB4VT1c8kIKlZ2O67dLpmI9DbH8VmtcasqL2wI8MiqZoaUCj88soiBJYldJicb21c0efJpjmNC2HY8cDPwp0ivh6V6szmO7KnL5jiyQ/u4X/tgm47/zWwd9+vn9bmlm5NaVzrKsjmOBMxxqOqSsG2+qv4MZ5ElY0weOnFUBbOuOpERFSVc+dAS/vPZlTQHQukOy6RIRB2HiPQJ2/qJyBeAgUmOzRiTwYb0KmbmlZOZdmwld722lvPuXMBHO+vTHZZJgUjnOMKfIA8A64Dfquq85IYXHZvjsDmORKS3OY7P6i7uRVsCzFjmrF9+4dhCjh0U+6R5NravaPLk0xxHUQfHCiO9HpbqzeY4sqcum+PIDpHE/dHOOj3r1nlaOX2WXvXw25953Xwi60p2WTbHkZjnON7o4NiCCPMaY/LAsD49mPndyfz01Gpmvb+ZL/7lNd74cEe6wzJJ0GXHISIDRWQCUCwiR4rIUe42BeiRkgiNMVnD5/Xw41NH8X/fO45Cv5dv3r2Q659aRn1zJ68vNlmpuwuRXwAuAobi3ILbaj9wTZJiMsZkufHDevHsVSfyp9mr+Pv89cxZvZ0/nXs4x47om+7QTAJ0OeJQ1ftVdSpwkapODdu+qqr/SlGMxpgsVFzg5fozxvLoFccCcP5db3LtE0vZ39jSTU6T6bpbj+PbqvqQiPwHHT85fnMH2dLG7qqyu6oSkd7uqvqseONuCij/WtPMCxsC9CoUvjO2gCP7d3zBIxvbVzR5cv6uKuC77tfrO9oinYFP9WZ3VWVPXXZXVXZIVNzvfLRbP3/zq1o5fZZe+eBi3bynIWl1xVOW3VXV9dbdQk53ul9/E1MXZowxYcYP68W/f3QC98xby/++tIbX1+zg55+vZtrkqgOW7zWZLaKndESkArgcZz2OtjyqeklywjLG5KoCn4fvTxnJlw8bxHVPLuOGf6/gsSWbuPFr4zjqoN7pDs9EINLHO58CXgdeAoLJC8cYky8q+5bwwCWTeGbpZm6ctYKzb3uD848exnFlttZHpou04+ihqtOTGokxJu+ICF85fDBTRvfnry+vYca8dTztUXaVrOPbx1bi8yb2le0mMSL9V5klIl9KaiTGmLxVWujjmi8dyvM/OZHh5R5u+PcKvvK3ebxRY0+eZ6JIO44f43QeDSKyT0T2i8i+ZAbWSkRGiMi9IvJ4KuozxqTPyP5l/HxiEXd8+yj2Nwb45j0L+e6Di+2tuxkm0vU4ylTVo6rFqtrT3e/ZXT4RmSEi20RkWbvjp4vIahGpEZFfdlP3WlW9NJI4jTHZT0Q4fdwgXv6Pk/n556t5fc0OTr35VW56biX77OHBjBDpXVVHdXB4L7BBVbt6Cc19wC3AA2FleYFbgdOATcAiEXka8AI3tct/iapuiyRGY0xuKfJ7+eHnRnHexGH88flV3PnqWh5fvImfnlbN+UcPs/mPNIp0cvw24Chgqbt/GPAe0FdErlTVFzrKpKqviUhVu8OTgBpVXQsgIo8AZ6rqTcBXogvfGJPrBvQs4uavj+fi44Zz46wVXPfkMu57Yz3XfOkQpo7uj4g9/5FqkS7k9Ahwo6oud/fHAL8AbgT+parju8hbBcxS1XHu/rnA6ap6mbs/DThGVX/YSf6+wO9xRij3uB1MR+muAK4AqKiomDBz5sxOvx97JUTm1GWvHMkOmdK+VJW3twWZubqZrfXKoX08fGN0AVXl3qjLijWGePPkwitHIh1xHNLaaQCo6goROVJV18bQ23eUodPeS1V3Ald2V6iq3iUim4EzfD7fhGiDMsZkPhFhwgAfR1R4mbsxwFM1zdywoJFjB3k5Z1QBFT3s8lVKRPJeEuBR4HbgZHe7DZgJFAKLuslbBSwL258MzA7bvxq4OtJ3pESy2buqsqcue1dVdsjU9rW3oVn/9PxKHX3dszrqmmf1N08v1521TTGVFWsM0ebJhXdVRXqpqhj4PnACzohhntt5NOI8HFjbRd4qDrxU5QM+AE4BPgYWAd/UsBFNrOztuOm/lJCqsuxSVWplevva3RjiiZoWXt8UoNALXxrh5/OVfgKNdXapKkIJX3M81g14GNgMtODcQXWpe/xLOJ3Hh8C1ia7XRhzZU5eNOLJDtrSvNVv36eX3L9LK6bN0wo0v6LX3vaBNLcGUxGAjjnZEZBTOrbJjgKKwTmdEFB1a0tmII3P/Ikx0WTbiSK1sa181u4M89kEzq3eH6FcsnDXSz+TBPjwRzsnaiKMbkfQuOJemTgHeByqBG4DfRNo7pXqzEUf21GUjjuyQje0rFArpX2e+qF/+62taOX2Wnvrnufrc0k80FAolJYZ8GnFEegtCsaq+jHP77gZVvQH4XHT9mTHGpI6IcFiFj6d/cAK3fesoQqpc+dDbfPWW+cxdva31j2ITg0gvVc0HTgQeB17BmdT+g6qOTm540bFLVdlzKSHesuxSVWrlQvsKhpQ3Pgnw1Ict7GhQqnt7OHtUAYf0+ewzIHapqhuRDEuAo4FSYCjwd+BfwLGRDmtSvdmlquypyy5VZYdcal9NLUF9YMF6nfT7F7Vy+iz95t0LdPH6nXHHkE+XqiIacWQLG3Fk91+EycxnI4745GL7ag4qczYGeGZtM/ua4bB+Xs4a6WdEL6+NOLrTVa8CPN3VFmnvlOrNRhzZU5eNOLJDLrevuqYWvW1OjY7/zWytnD5LL73vLb3vqZeSUlek6TJ9xNHdK0cmAxtxnsdYSMevCzHGmKzVo8DH96YczLTJldz/xnruem0tL61sYf6exfz41FGMHVye7hAzTnd3VQ0ErgHGAf+L86LBHar6qqq+muzgjDEmVUoLffxg6khenz6Vs0b6WbB2J1/+6zy+++BiVnySknXrskbEcxwiUghcAPwX8FtV/VsyA4uFzXHk1jXoROazOY745GP7ksISXljfwgsbWmgIwIQBXs482M9BPeN7E2/Oz3G4nUohcDbwGM57pX4FDIn0Wlg6NpvjyJ66bI4jO+Rz+9pT16w3v7Bax/36ea2cPksvv3+RLt20J+a6cn6OQ0Tux7lM9RzOk+LLukpvjDG5pryHn5+eVs0lxw9nxvx1zJi/jhdWbOXUQwfw41NGcdjQ/JsD6W5yfBpQB1QDV4WtvSGAagTrjhtjTC5o60BOGM7f569jxrx1nHHLVj53SH9+fMqodIeXUvYcRxTnMlk+XoO2OY7Usfb1WfUtyksftTB7fQt1LTCmt3J2dTEje3c8BxJN2Vk/x5GNm81xZE9dNseRHax9dW5/Y4veOmeNjrtuVtuT6As+3BFX2Vk9x2GMMaZrpYU+vj9lJAcHNvJRQSV3vraW8+96k0lVffjRKSM5YWQ/YlhiO6PZAr3GGJMAhT7h8pNGMG/6VK4/Ywwf7apn2r1vcdZtb/DKqq2td6nmBBtxGGNMAhX5vVx8/HC+ecxBPLZ4E7fP/ZBL7lvM2ME9+dHnRlKQAx2IjTiMMSYJCn1evn1sJXN/MYU/nXs4dU0Brnzoba6b38BT735MIBhKd4gxs47DGGOSyO/18PWJw3jpZyfzv+ePR4AfP/Iup978Ko8u+ojmQPZ1IHY7bhTnMpndLpm49HY77mdZ+0pcnn37a1lTX8S/P2xh/b4QfYqELw33c9JQHwVe6bYsux3XbsdNGLtdMnHp7Xbcz7L2lbg8relCoZDOWbVVz7ltvlZOn6UTbnxRb59bo/sbW+x2XGOMMZ8lIkwZ3Z+TqytYuG4Xt86p4Q/PreL2uR8yZTCMn9RMrx4F6Q6zQ9ZxGGNMGokIx47oy7Ej+vLuxj3cOqeGp1Zs5aU/vMK3j63k0hOG079nUbrDPIBNjhtjTIYYP6wXd39nIjceX8wphw7g7tfXcsKf5nDdk0vZuKs+3eG1sRGHMcZkmGFlHqadcSQ/O62aO179kEcXbeThtzZy5vjBTOyR/ruwMn7EISJfE5G7ReQpEfl8uuMxxphUqepXwh/OOZzX/t9ULpxcxbNLN3PtvAa+99ASlm7am7a4ktpxiMgMEdkmIsvaHT9dRFaLSI2I/LKrMlT1SVW9HLgI+EYSwzXGmIw0qLyYX58xhvnTP8dXRviZV7ODM26Zx7R7F/Lm2p0pf51Jsi9V3QfcAjzQekBEvMCtOOuXbwIWicjTgBe4qV3+S1R1m/v5OjefMcbkpb6lhZxTXcB/fud4HnrzI+6d57xQcUJlb74/5WA8KepAktpxqOprIlLV7vAkoEZV1wKIyCPAmap6E/CV9mWI81rJPwDPqerbyYzXGGOyQVmRn+9NOZiLj6/iscUbuePVtVx6/2KGlgr/mthI/7Lk3oWV9CfH3Y5jlqqOc/fPBU5X1cvc/WnAMar6w07yXwVciLPe+buqekcn6a4ArgCoqKiYMHPmzE5jsid7M6cue3I8O1j7SlyeZCzkFAgpCzcHWLy5iasmlMT0GveMenIcqAKWhe2fB9wTtj8N+Fsi67Qnx7OnLntyPDtY+0pcnlxYyCkdI47JwA2q+gV3/2q3A2s/vxFLXfauqiyry0Yc2cHaV+Ly2NKxsY04fMBaYDhQALwHjE1knTbiyJ66bMSRHax9JS5PLow4kjo5LiIPA1OAfiKyCbheVe8VkR8Cs3HupJqhqssTVN8ZwBlAvYis7CJpOdDZTdD9gB2JiCfFuvqeMrmuWMuKNl806btLG895a1+prStV7SuaPJGk6ypNstpXZcQpI+1hsmkD7or1PFH0upm0dfc9Z2pdsZYVbb5o0sfTfro7b+0rtXWlqn1FkyeSdN20obS3r4x/cjxG/47zfDZK5feUyLpiLSvafNGkj7f9WPvKnLpS1b6iyRNJuq7SpL195dRCTokgIos10gkiY6Jk7cskU6raV66OOOJxV7oDMDnN2pdJppS0LxtxGGOMiYqNOIwxxkTFOg5jjDFRsY7DGGNMVKzjiIKIjBCRe0Xk8XTHYnKDiJSIyP3uYmXfSnc8Jrck63dW3nQcCVpUaq2qXprcSE22i7KtnQ08rs5iZV9NebAm60TTvpL1OytvOg6cRaVODz8QtqjUF4ExwAUiMkZEDhORWe22/qkP2WSp+4iwrQFDgY1usmAKYzTZ6z4ib19JkewVADOGJmBRKWMiEU1bw1kFcyjwLvn1h5yJUZTta0UyYsj3hjqET//aA+d/4iGdJRaRviJyB3Bk6+vgjYlQZ23tX8A5InI7GfAqCZO1OmxfyfqdlTcjjk50tExWp09EqupO4MrkhWNyWIdtTVXrgItTHYzJOZ21r6T8zsr3EccmYFjY/lDgkzTFYnKbtTWTTCltX/necSwCRonIcBEpAM4Hnk5zTCY3WVszyZTS9pU3HYe7qNQCYLSIbBKRS1U1ALQuKrUSmKkJWlTK5C9rayaZMqF92UsOjTHGRCVvRhzGGGMSIyfvqurVq5eOHDmy0/N1dXWUlJREfS6TpTLuRNYVa1nR5osmfXdp4zlv7Su1daWqfUWTJ5J0sbaheH52S5Ys2aGqFRElTsX6tKneqqurtStz5syJ6VwmS2Xciawr1rKizRdN+u7SxnPe2ldq60pV+4omTyTpYm1D8fzsiGK9crtUZYwxJio5eakqVjtrm9jTFGL7/iY8AiLOMzUCeERAQMT5LHz6Gfd8+3Ot+Y0xJpdYxxHmvDsXsHZ7A8x5KWFltnZArV9bO6G2ffl03+l8Pv3cltfTuv/Z9B73fN3+Bv6yfD5ez4HnvB7B4x7zivPZ6x4XAW/Yfmtan/vZK4LX63z1hZ3bsL6ZNZ61eD2Cz+uk9Xs8bfs+97PfK/i8Hvwe56vP66Tz+5w0BV4PuxpD7Khtwu919lvzGGMyl3UcYX5yajWL31vOqFGjCLl3KasqCoTU+QwQUkUV97jzGSAU0rZjIQVUD0jTWoa6aUMKinMuGFIU51go1JpeCboVOZ9b56Sc9K31qCrbm+ooK/IdcC4QCtEUcMtUJRjStnNBt46g+zkUUgJh5wKhT9O37h/gg5WJ+8HPPbCj9ghOR+LzUOhzOpRCv9f96hyr39/Ag+sXUeT3Uuhzzhf7vRQXeCj2eynye+lR4KO4wEOPAh81O4KUbdhNSaGXkgIfZUU+Sgp9+K2TMiZq1nGE+eoRg+m5+wOmTK5KdyhRmzt3LlOmHJPUOkIhpyN7Ze6rHH/CiQSDSkso1NaxBIKhtg6nJRgiEHQ6r5agEnDTBoLOudbzS1es5OCRo2gOKs0B53hzIERz+6+BEE2BIE2tn4OweW8jTYEgjS2htq8NLUGC7Tu5Vovf+MyhIr+H0kI/PYt89Cz207PYT9O+Rl7es4zePfz06lFA39IC+pYU0q+sgIrSQkL27JPJc9ZxmIh5PIIHodArlBYmpun03V8TU0ftdJQndniuOeB0II0tQeqbg9Q3B5i/cDHVYw6jvjlIbVOA2sYAdU0B9jcF2N/Ywr7GAPsaWthb38yWfSFWv/8Jexta6KiP8Aj0X/Ayg3oVMai8iMHlxRzUtwfD+vSgsk+PzjsuY3KEdRwm5xT4nMtc5cX+tmPbe3uZMjqytbicTmkKwZCyt6GFXXVN7KxtZkdtM9v3N7J4+RqKevdj894GVm3ez8srt9EUCLXl9wqMePdVRg0oZezgcsYNKeewIeX0KSlI+PdqTDpYx2FMJ7weoU9JAX1KChgZ1udUtWxgypQj2vZVle0Na85dAAAZ3ElEQVS1TWzcVc+6HfW8sngFzUUlLPt4H88u3dKWbmT/Ug4qaqK+72ZOHNWPsiI/xmQj6ziMiZOI0L+siP5lRUyo7EO//TVMmTIRgL0NLSz/eC/vbdrLm2t3Mv/DWl75x9sU+DycNKqCLx8+kNPHDqK4wJvm78KYyFnHYUwSlRf7OW5kP44b2Y/vTTmYl16ZQ/mII3hu6RaeW7aZl1Zu5bc9VjDt2EqmTa6ioqww3SEb0y3rOIxJIZ9HOLqqD0dX9eG6Lx/KovW7uGfeOv42p4Y7XlvLZScM56pTRlHktxGIyVzWcRiTJh6PcMyIvhwzoi8fbq/l1ldquG3uh8xevoU/nXsEEyp7pztEYzpkTz8ZkwEOrijl5m+M5/5LJtHQHOTcO97g1jk16Q7LmA5Zx2FMBjm5uoLZPz2Jrxw+mP+avZrb5lrnYTJPyjsOEZkhIttEZFnYsRtE5GMRedfdvhR27moRqRGR1SLyhVTHa0yqlRX5+cs3xnPm+MH86fnV3PP62nSHZMwB0jHHcR9wC/BAu+P/o6r/HX5ARMbgLLo+FhgMvCQi1aoaTEWgxqSL1yP8+bwjaAmG+N0zKykp9HHBpIPSHZYxQBpGHKr6GrArwuRnAo+oapOqrgNqgElJC86YDOLzevjf84/khJH9+P0zK9m+vyndIRkDgGgaXtgmIlXALFUd5+7fAFwE7AMWA/+hqrtF5BbgTVV9yE13L/Ccqj7eQZlXAFcAVFRUTJg5c2an9dfW1lJaWhr1uUyWyrgTWVesZUWbL5r03aWN53ws3++WuhDXzmvghCE+Lh6Xnuc8rH0lLk8k6WJtQ/H87KZOnbpEVSdGlDjSpQITuQFVwLKw/QGAF2cE9Htghnv8VuDbYenuBc7prnxbOjZ76rKlYyPz238v16pfztJlH++JKX+8rH0lLo8tHZsgqrpVVYOqGgLu5tPLUZuAYWFJhwKfpDo+Y9Ltqs+NolexnxtnrWhbF8aYdImp4xCRoSLycxF5SkQWichrInKbiHxZRKIuU0QGhe2eBbTecfU0cL6IFIrIcGAU8FYsMRuTzcp7+PnZadW8uXYXL6zYmu5wTJ6L5Zf834EZQDPwR+AC4PvAS8DpwDwROamL/A8DC4DRIrJJRC4F/iQiS0XkfWAq8FMAVV0OzARWAM8DP1C7o8rkqQsmHcTI/qX2YKBJu1hux/2zqi7r4Pgy4F8iUgB0et+gql7QweF7u0j/e5x5D2Pyms/r4awjh/Bfs1ezZW8jA8uL0h2SyVNRjzg66TTCzzerqv1JZEwSnDZmAAAvr7LLVSZ9Yn4AUESWAu1n6fbi3E77O1XdGU9gxpjPGtW/lIP69ODlldv41jGV6Q7H5Kl4nhx/DggC/3T3z3e/7sN5OvyMOMo2xnRARDj10AE8tHAD9c0BehTYC65N6sVzO+7xqnq1qi51t2uBKar6R5znNIwxSXDqof1pDoR4fc2OdIdi8lQ8HUepiBzTuiMik4DWRxYDcUVljOnU0cP7UFbk4yW7LdekSTzj3MuAGSLS2lnsBy4TkRLgprgjM8Z0yO/1MHV0f15ZtY1gSPF6JN0hmTwT84hDVRep6mHAeOBIVT1cVd9S1TpV7fxFUcaYuJ1yaH921jXz7sY96Q7F5KGYOw4RGeC+dPARVd0jImPch/mMMUk2pbo/Po/w0kq7XGVSL545jvuA2TjrZAB8APwk3oCMMd0r7+Hn6Ko+zF29Pd2hmDwUT8fRz70kFQJQ1QDO7bnGmBQYf1Av1mzdT3MglO5QTJ6Jp+OoE5G+uA8BisixOA8AGmNS4JCBZQRCyofba9Mdiskz8dxV9TOct9ceLCLzgQrg3IREZYzp1qGDegKwasu+ts/GpELMHYeqvi0iJwOjAQFWq2pLwiIzxnRpeL8SCrweVm3eD0emOxqTT6LuOETk7E5OVYsIqvqvOGMyxkTA7/UwakApK7fsT3coJs/EMuJofQdVf+A44BV3fyowF7COw5gUOWRgT15fY3dWmdSK5bXqF6vqxTiT4mNU9RxVPQcYm/DojDFdOnRQGdv2N7GztindoZg8Es9dVVWqujlsfytQ3V0mEZkhIttEZFnYsT4i8qKIrHG/9naPi4j8VURqROR9ETkqjniNyTmHDHQmxVfb5SqTQvF0HHNFZLaIXCQiFwLPAHMiyHcfzhKz4X4JvKyqo4CX3X2AL+KsMz4KuAK4PY54jck5hwwqA7B5DpNS8byr6ofAHcAROO+ruktVfxRBvteAXe0Onwnc736+H/ha2PEH1PEm0EtEBsUaszG5pl9pIf1KC1m1eV+6QzF5JJa7qkRVFUBVnwCe6CpNhAa0XvZS1c0i0t89PgTYGJZuk3tsM8YYwJnnWGUjDpNCEt3vdxCRucD/AU+p6kdhxwuAE4ALgTmqel8XZVQBs1R1nLu/R1V7hZ3fraq9ReQZ4CZVnecefxn4f6q6pIMyr8C5nEVFRcWEmTM7f0FvbW0tpaWlUZ/LZKmMO5F1xVpWtPmiSd9d2njOJ+Pf6ZFVTbz0UYA7T+2RtFesW/tKXJ5I0sXahuL52U2dOnWJqk6MKLGqRrUBRcD3gfnAJ8AKYB2wAbgbGB9BGVXAsrD91cAg9/MgnIcJAe4ELugoXVdbdXW1dmXOnDkxnctkqYw7kXXFWla0+aJJ313aeM4n49/p8cUbtXL6LF2zdV/Cy25l7StxeSJJF2sbiudnByzWCPuBWG7HbVTV21T1eKASOAVnPY5KVb1cVd+NtkycV5dc6H6+EHgq7Ph33LurjgX26oF3chmT99omyDfb5SqTGvHcVYWqtqjqZlWNeDUZEXkYWACMFpFN7hoefwBOE5E1wGnuPsCzwFqgBmc08/144jUmF43sX4rXI6zaYhPkJjXieclhTFT1gk5OndJBWgV+kNyIjMluhT4vB1eUOO+sMiYF4hpxGGMywyEDe7LSbsk1KRJXxyEilSJyqvu5WETKEhOWMSYahw7qySd7G9lT35zuUEweiGfN8cuBx3HufAIYCjyZiKCMMdEZM9h59YhNkJtUiGfE8QPgeGAfgKquwXljrjEmxca4CzmtsMtVJgXi6TiaVLVtXCwiPtxlZI0xqVVR5rx6ZMUn1nGY5Iun43hVRK4BikXkNOAx4N+JCcsYE60xg22C3KRGPB3HL4HtwFLguzjPXFyXiKCMMdEbM6gna7btpzkQSncoJsfFs+Z4COehvLtFpA8w1H3uwhiTBmMG96QlqNRsq22bLDcmGeK5q2quiPR0O413gb+LyM2JC80YE40xba8esctVJrniuVRVrqr7gLOBv6vqBODUxIRljInW8H6lFPk9dmeVSbp4Og6fu6jS14FZCYrHGBMjr0cYPbCn3Vllki6ejuO3wGygRlUXicgIYE1iwjLGxGLMoDJWbtmHTTeaZIpn6djHVPVwVf2+u79WVc9JXGjGmGiNGdSTPfUtbN7bmO5QTA6L+a4qESkCLgXG4izuBICqXpKAuIwxMWi9m2rFJ/sY3Ks4zdGYXBXPpaoHgYHAF4BXcd5VZS/KMSaNRg/siYjdWWWSK56OY6Sq/gqoU9X7gS8DhyUmLGNMLEoLfVT26WF3Vpmkimchpxb36x4RGQdswVlLPGYish5n1BIEAqo60X1O5FG37PXA11V1dzz1GJPLxgzuyXK7s8okUTwjjrtEpDfwK5y1wVcAf0pATFNVdbyqTnT3fwm8rKqjgJfdfWNMJ8YP68WGnfV8vKch3aGYHBXPXVX3qOpuVX1VVUeoan9VvSORwbnOBO53P98PfC0JdRiTMz4/ZiAAzy/bkuZITK6K566qQuAcnEtIbeWo6m/jiEeBF0REgTtV9S5ggKpudsveLCK25ocxXajqV8IhA8t4ftlmLj1heLrDMTlIYn1QSESeB/YCS3DmJABQ1T/HHIzIYFX9xO0cXgR+BDytqr3C0uxW1d4d5L0CuAKgoqJiwsyZMzutp7a2ltLS0qjPZbJUxp3IumItK9p80aTvLm0851P17/RkTTNP1bTwP1OL6VUY1wrRgLWvROaJJF2sbSien93UqVOXhE0RdE1VY9qAZbHmjbD8G4CfA6uBQe6xQcDq7vJWV1drV+bMmRPTuUyWyrgTWVesZUWbL5r03aWN53yq/p1Wbd6nldNn6YML1iekPGtficsTSbpY21A8PztgsUb4+zmeP0XeEJGE3X4rIiUiUtb6Gfg8sAxn4v1CN9mFwFOJqtOYXFU9oJQR/UpsnsMkRdRzHCKyFGcuwgdcLCJrgSZAAFXVw2OMZQDwhIi0xvVPVX1eRBYBM0XkUuAj4LwYyzcmb4gIp48byJ2vrWV3XTO9SwrSHZLJIbFMjn8l4VHgvOsKOKKD4zuBU5JRpzG57IvjBnHb3A95ceVWvj5xWLrDMTkklktVW4GzgF8ApwMfq+qG1i2h0RljYjZuSE+G9Cpmtl2uMgkWS8dxPzARZ63xLwIx30VljEkeEeGL4wby+podbN5rDwOaxIml4xijqt9W1TuBc4ETExyTMSZBLjyuChH43ayV6Q7F5JBYOo7Wd1ShqoEExmKMSbBhfXrww6kjeWbpZl5fsz3d4ZgcEUvHcYSI7HO3/cDhrZ9FxN6sZkyGufykEVT17cH1Ty2nKRDsPoMx3Yi641BVr6r2dLcyVfWFfe6ZjCCNMbEr8nu54atjWbujjnteX5fucEwOiP9dBMaYjDdldH9OHzuQv72yhiUbbFUCEx/rOIzJE789cywDexbxnXsX8ta6XekOx2Qx6ziMyRP9exbx6HcnM6C8iAtnvMWCD3emOySTpazjMCaPDOhZxKNXTGZYn2Iu+vtb3PP6WgLBULrDMlnGOg5j8kxFWSEPX34sx4/sx++eWclZt73Bso/3pjssk0Ws4zAmD/UtLeTeCydy6zePYvPeRr56yzx+8I+3WbJhV+uyBsZ0KuYVAI0x2U1E+PLhgzhhZD9uf/VD/rlwA88s3cwRQ8s5+6ihfGHsQAaWF6U7TJOBrOMwJs+V9/Dzyy8ewlWnjOT/lmziwTc3cP3Ty7n+6eUceVAvThxVQdG+IMe2BCnye9MdrskA1nEYYwDoUeBj2uQqpk2uomZbLbOXb+GF5Vu45ZU1hBT+8vYLHDKojLGDyxk3pCeHDCzj4IpSevWwtT7yjXUcxpjPGNm/lJH9R/KDqSPZ19jCjKdfpaFsCMs+3suzSzfz8FsftaXtV1pAVd8ShvXpwbA+PRjaq5iB5UUMKi9iQHkRZYU+3AXaTI7Iio5DRE4H/hfwAveo6h/SHJIxeaNnkZ/x/X1MmXIoAKrKpt0NrNm2n5pttdRsq2XDznreWreLp979mFC7ufVCn4eKskL6lRbSt6SA3iUF9C0poLyHn17FBfTq4adnkZ+exT56FvnZ16Q0BYIU+uyyWKbK+I5DRLzArcBpwCZgkYg8raor0huZMflJRNpGF587ZMAB55oDIbbua2TLvkY+2dPAtn1NbK9tYvv+JnbUNrF5byPLP9nHrvpmmgNdPD8y53kKvB5KCr2UFPooLfRRXOClpMD52sPdivxeit2tyO+lyO+h0O+l0OehyO9l9fYAhR/upNDvocDrocjvocDrpcDnwe8VCnweZ/N6bFQUhYzvOIBJQI27tCwi8ghwJmAdhzEZpsDnaetUutPYEmR3fTN76lvY3xhgb0ML+xpaeGfZSgYdNJx9jS3UNwWpawpQ2xSgocX5vKO2iYaWIA3NztYYCNIS7OIW4iVvRhS73yv4vR58HoFQgJIFL+PzCn6PxznuFXxeD36POMe9Hrwewedx8uza2ciTW97B6+57veJ89bR+9eD1wKaPmnkvsAafV/CI4PWARz5N6/EINRtb2LZ4I14RPO55r0fwirBia4DAiq1tx1vPicDqXUFOCikeT3I7wWzoOIYAG8P2NwHHpCkWY0yCFPm9DCovZlB58QHH++6vYcqUkVGV1RIM0dgSpLHF+doUcD6/+dZixhx2BE3BEM2BEE0B56uzOR1OczBEi3s+EFKaAyHWb9xE/wH92s4HgiECQaUlpG2fawMBgiElEFQCoRD7akNsbdlDIKjO8ZBzPBjStnRBdT7z4Qfdf1PL3+/83DuLOz118VeVgiR3HJLpD/uIyHnAF1T1Mnd/GjBJVX/ULt0VwBUAFRUVE2bOnNlpmbW1tZSWlkZ9LpOlMu5E1hVrWdHmiyZ9d2njOW/tK7V1pap9RZNn//5aSkpLCCkEFULupm37Sm1dPUXFPVAOPB8C6urqKSouJkTYcXerb2jgyCEleGK47DZ16tQlqjoxosSqmtEbMBmYHbZ/NXB1V3mqq6u1K3PmzInpXCZLZdyJrCvWsqLNF0367tLGc97aV2rrSlX7iiZPJOlibUPx/OyAxRrh7+VseOXIImCUiAwXkQLgfODpNMdkjDF5K+PnOFQ1ICI/BGbj3I47Q1WXpzksY4zJWxk/xxELEdkLrOkiSTnQ2etA+wE7Eh5U8nX1PWVyXbGWFW2+aNJ3lzae89a+UltXqtpXNHkiSddVmmS1r0pVrYgoZaTXtLJpA+6K9TxRXOfLpK277zlT64q1rGjzRZM+nvbT3XlrX6mtK1XtK5o8kaTrpg2lvX1lwxxHLP4d5/lslMrvKZF1xVpWtPmiSR9v+7H2lTl1pap9RZMnknRdpUl7+8rJS1XxEJHFGuktacZEydqXSaZUta9cHXHE4650B2BymrUvk0wpaV824jDGGBMVG3EYY4yJinUcxhhjomIdhzHGmKhYxxEFERkhIveKyOPpjsXkBhEpEZH7ReRuEflWuuMxuSVZv7PypuMQkRkisk1ElrU7frqIrBaRGhH5ZVdlqOpaVb00uZGabBdlWzsbeFxVLwe+mvJgTdaJpn0l63dW3nQcwH3A6eEHwlYX/CIwBrhARMaIyGEiMqvd1j/1IZssdR8RtjVgKJ+uNxNMYYwme91H5O0rKTL+JYeJoqqviUhVu8Mdri6oqjcBX0lthCZXRNPWcBYmGwq8S379IWdiFGX7SspKqfneUDtaXXBIZ4lFpK+I3AEcKSJXJzs4k1M6a2v/As4RkdvJgFdJmKzVYftK1u+svBlxdKKjZbI6fSJSVXcCVyYvHJPDOmxrqloHXJzqYEzO6ax9JeV3Vr6PODYBw8L2hwKfpCkWk9usrZlkSmn7yveOw1YXNKlibc0kU0rbV950HCLyMLAAGC0im0TkUlUNAK2rC64EZqqtLmjiZG3NJFMmtC97yaExxpio5M2IwxhjTGJYx2GMMSYq1nEYY4yJinUcxhhjomIdhzHGmKhYx2GMMSYq1nGYrCEiQRF5N2yrSndMiSQiR4rIPXGWcZ+InBu2f4GIXBt/dCAiPxQRez2Kyft3VZns0qCq4zs7KSI+90GobHUN8Lv2B+P8vk4H/hpXVJ+aAcwH/p6g8kyWshGHyWoicpGIPCYi/wZecI/9QkQWicj7IvKbsLTXugvdvCQiD4vIz93jc0Vkovu5n4isdz97ReS/wsr6rnt8ipvncRFZJSL/EBFxzx0tIm+IyHsi8paIlInI6yIyPiyO+SJyeLvvoww4XFXfc/dvEJG7ROQF4AERqXLLedvdjnPTiYjcIiIrROQZoH9YmQKMB94WkZPDRmrvuPV19bP6jnvsPRF5EEBV64H1IjIpEf92JnvZiMNkk2IRedf9vE5Vz3I/T8b5pbtLRD4PjMJZn0CAp0XkJKAO5/09R+K0+7eBJd3UdymwV1WPFpFCYL77ixy3nLE4L5KbDxwvIm8BjwLfUNVFItITaADuAS4CfiIi1UChqr7frq6JwLJ2xyYAJ6hqg4j0AE5T1UYRGQU87OY5CxgNHAYMwFl/YUZYjO+pqrqd5A9Udb6IlAKNXfysdgLXAser6g4R6RMW02LgROCtbn52JodZx2GySWeXql5U1V3u58+72zvufinOL8cy4An3r2ZEJJIXwH0eODxszqDcLasZeEtVN7llvQtUAXuBzaq6CEBV97nnHwN+JSK/AC7BWcGtvUHA9nbHnlbVBvezH7jFHbkEgWr3+EnAw6oaBD4RkVfC8p8OPOd+ng/cLCL/AP6lqpvcjqOjn9UROMvZ7nC/j11hZW4DDun4x2XyhXUcJhfUhX0W4CZVvTM8gYj8hM7XWgnw6WXbonZl/UhVZ7crawrQFHYoiPP/knRUh6rWi8iLOCuyfR1npNBeQ7u64cDv66fAVpxf6h6gMbyKjr4pnE7hHDeGP7iXsr4EvCkip9L5z+qqLsoscmM1eczmOEyumQ1c4l6OQUSGiLNe/GvAWSJS7F7fPyMsz3qcy0IA57Yr63si4nfLqhaRki7qXgUMFpGj3fRlItL6x9k9OJPUi9r9Bd9qJTCyi7LLcUYzIWAa4HWPvwac787HDAKmunWXAz53IR9E5GBVXaqqf8S53HQInf+sXga+LiJ93ePhl6qq+ewlNZNnbMRhcoqqviAihwIL3PnqWuDbqvq2iDyKs7b3BuD1sGz/DcwUkWlA+KWee3AuQb3tTjRvB77WRd3NIvIN4G8iUozzl/mpQK2qLhGRfXRyR5KqrhKRchEpU9X9HSS5Dfg/ETkPmMOno5EngM8BS4EPgFfd46cBL4Xl/4mITMUZHa0AnlPVpk5+VstF5PfAqyISxLmUdZFbzvHAbzB5zV6rbvKSiNyA8wv9v1NU32BgLnCIO2roKM1Pgf2qGtezHG5Z9wD3qOqb8ZYVVuaRwM9UdVqiyjTZyS5VGZNkIvIdYCFwbWedhut2Dpw7iZmqXpbITsPVD/hVgss0WchGHMYYY6JiIw5jjDFRsY7DGGNMVKzjMMYYExXrOIwxxkTFOg5jjDFRsY7DGGNMVP4/3n1SHQYSofYAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f2b1b870518>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "control.bode(ss[0,0], omega=pl.logspace(-1,1,100));"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEKCAYAAAAFJbKyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAFoJJREFUeJzt3X2wXXV97/H318SAJDwIBw4PCQZjLhYskuZIJXrpiahFdEB78Qq5trbUZjI00t5bH/AyIwXbEWt779XGXorAiFaMDggEiEapOWgngiQEIk+pCSLEcIUATTigPOV7/9grckjOw1o5Z5+118n7NbPnrLX2b+/1/Q3M/mSt31q/FZmJJEllvaLuAiRJzWJwSJIqMTgkSZUYHJKkSgwOSVIlBockqRKDQ5JUicEhSarE4JAkVTK5zp1HxBXAe4BHM/MNg7zfC1wP/KzY9K3MvGik7+3q6sqZM2eWquHpp59m6tSpZUtuDPvVHBOxT2C/mmbNmjVbMvPgMm1rDQ7gy8AS4CvDtPlhZr6nypfOnDmT1atXl2rb19dHb29vla9vBPvVHBOxT2C/miYifl62ba2nqjLzB8ATddYgSaqmCWMcJ0bEXRHx7Yg4tu5iJGlPF3XPjhsRM4Ebhxjj2A/Ynpn9EXEq8PnMnD3E9ywEFgJ0d3fPXbp0aan99/f3M23atN2svnPZr+aYiH0C+9U08+fPX5OZPaUaZ2atL2AmcHfJtg8CXSO1mzt3bpa1cuXK0m2bxH41x0TsU6b9ahpgdZb83e7oU1URcWhERLF8Aq1Ta4/XW5Uk7dlqDY6I+DrwI+DoiNgUEX8aEYsiYlHR5Azg7oi4C/gCcGaRjFJjXHLLRlZt3PKybas2buGSWzbWVJE0OrVejpuZZ43w/hJal+tKjXXc9P1ZfNValiyYw7xZXazauOU361IT1X0fhzThzZvVxZIFczjn6mW8+w2H8+07Jv0mRKQmMjikcTBvVhcHzriJbz2Y/NnvfsbQUKN19OC4NFGs2riFLf3PceSBr+JfbntolzEPqUkMDqnNdoxpHNU1ldccNJUlC+aw+Kq1hocay+CQ2mzdpq0sWTCHffdunRneMeaxbtPWmiuTdo9jHFKbLfq9WQD80/qXts2b1eU4hxrLIw5JUiUGhySpEoNDklSJwSFJqsTgkCRVYnBIkioxOCRJlRgckqRKDA5JUiUGhySpEoNDklSJwSFJqsTgkCRVYnBIkioxOCRJlRgckqRKDA5JUiUGhySpEoNDklSJwSFJqqTW4IiIKyLi0Yi4e4j3IyK+EBEbImJdRPzOeNcoSXq5uo84vgycMsz77wJmF6+FwP8dh5okScOoNTgy8wfAE8M0OR34SrbcChwQEYeNT3WSpMHUfcQxkiOAhwesbyq2SZJqMrnuAkYQg2zLQRtGLKR1Oovu7m76+vpK7aC/v7902yaxX51n27ZtPBvP7lJ/k/s0HPs1cXV6cGwCZgxYnw5sHqxhZl4KXArQ09OTvb29pXbQ19dH2bZNYr86z2XLL2OfyfvsUn+T+zQc+zVxdfqpqmXAHxVXV70Z2JqZj9RdlCTtyWo94oiIrwO9QFdEbAIuAF4JkJmXAMuBU4ENwDPAn9RTqSRph1qDIzPPGuH9BP58nMqRJJXQ6aeqJEkdxuCQJFVicEiSKjE4JEmVGBySpEoMDklSJQaHJKkSg0OSVInBIUmqxOCQJFVicEiSKjE4JEmVGBySpEoMDklSJQaHJKkSg0OSVInBIUmqxOCQJFVicEiSKjE4JEmVGBySpEoMDklSJQaHJKkSg0OSVInBIUmqxOCQJFVSa3BExCkRsT4iNkTEeYO8/8cR8VhE3Fm8PlxHnZKkl0yua8cRMQn4IvAOYBNwe0Qsy8x7d2r6jcxcPO4FSpIGVecRxwnAhsx8IDOfA5YCp9dYjySphDqD4wjg4QHrm4ptO/svEbEuIq6OiBnjU5okaSi1naoCYpBtudP6DcDXM/PZiFgEXAm8bdAvi1gILATo7u6mr6+vVBH9/f2l2zaJ/eo827Zt49l4dpf6m9yn4diviavO4NgEDDyCmA5sHtggMx8fsPol4LNDfVlmXgpcCtDT05O9vb2liujr66Ns2yaxX53nsuWXsc/kfXapv8l9Go79mrjqPFV1OzA7Io6KiCnAmcCygQ0i4rABq6cB941jfZKkQdR2xJGZL0TEYmAFMAm4IjPviYiLgNWZuQw4NyJOA14AngD+uK56JUktdZ6qIjOXA8t32vapAcufBD453nVJkobmneOSpEoMDklSJQaHJKkSg0OSVInBIUmqxOCQJFVicEiSKikVHMUU6JIklT7i2BARn4uIY9pajSSp45UNjuOAfwcui4hbI2JhROzXxrokSR2qVHBk5lOZ+aXMnAd8HLgAeCQiroyI17W1QklSRyk9xhERp0XEtcDngX8AXkvreRnLh/2wJGlCKTvJ4U+BlcDnMnPVgO1XR8RJY1+WJKlTjRgcxRVVX87MiwZ7PzPPHfOqJEkda8RTVZn5IjB/HGqRJDVA2VNVqyJiCfAN4OkdGzPzjrZUJUnqWGWDY17xd+DpqgTeNrblSJI6XangyExPVUmSgAqPjo2IdwPHAnvv2DbUgLkkaeIqex/HJcAHgI8AAbwfeE0b65IkdaiyU47My8w/Ap7MzAuBE4EZ7StLktSpygbHr4q/z0TE4cDzwFHtKUmS1MnKjnHcGBEHAJ8D7qB1RdVlbatKktSxyl5V9eli8ZqIuBHYOzO3tq8saWJKsu4SpFGrclXVPGDmjs9EBJn5lTbVJU04R+57JCseXMH9T9zP6w98fd3lSLut7FVVXwX+Hngr8Kbi1dPGuqQJ56Nv+igH7HUAH7vlYzz9/NMjf0DqUGUHx3uAt2TmOZn5keI16skNI+KUiFgfERsi4rxB3t8rIr5RvH9bRMwc7T6luhy494FcfNLFPPTUQ3z61k+T6WkrNVPZ4LgbOHQsd1zMuvtF4F3AMcBZgzya9k9pXQL8OuB/A58dyxqk8famQ9/Eojcu4qYHbuK6DdfVXY60W8qOcXQB90bEj4Fnd2zMzNNGse8TgA2Z+QBARCwFTgfuHdDmdOCvi+WrgSUREek/1dRgC397Iav/32o+8+PPcNzBx9VdjlRZ2eD46zbs+wjg4QHrm4DfHapNZr4QEVuBg4AtbagHrvkwvPDsyO0a4NjHHoNfXl53GWOuUf165T7Qex4c+PJbnia9YhIX/+eLOeOGMzjn5nM4JA9hWd+ymopsn8cee8x+jbN9p+zLhfMubPt+yl6Oe0sb9h2D7Wo32rQaRiwEFgJ0d3fT19dXqoj+/v7ftP2dB9cy6cXnSn2u0+21fTtPP/OLussYc03q196/fpRnHriNtXP+ju2Tpuzy/gf3/yDXPnktj7z4CL/c/MsaKmyv7du3269xNvUVU+l7rq/9O8rMIV/AvxV/nwK2DXg9BWwb7rMjvWhNW7JiwPongU/u1GYFcGKxPJnWkUaM9N1z587NslauXFm6bZPYrw6w/juZF+yXeeNfDdusUX2qwH41C7A6S/5+Dzs4nplvLf7um5n7DXjtm5n7jTKzbgdmR8RRETEFOBPY+fhvGfChYvkM4PtFB6XO959+H05cDLd/Ce69vu5qpDFT6lRVRBw4yOanMvP53d1xtsYsFtM6qpgEXJGZ90TERbSSbxlwOfDViNgAPEErXKTmOPkCeOhHcP1H4LA3wqtn1l2RNGplB8fvoDUb7pO0xh0OAB6JiEeBP8vMNbuz88xcDizfadunBiz/mtYU7lIzTZ4CZ1wBl5wEV58Nf/Kd1japwcrex/Ed4NTM7MrMg2jde/FN4Bzgn9pVnDQhvHomnPYF+MUa+OHf112NNGql7xzPzBU7VjLzu8BJmXkrsFdbKpMmkmPfC4cdDw//uO5KpFEre6rqiYj4BLC0WP8A8GRx9/f2tlQmTTSDXJIrNVHZI44FwHTgOuB64Mhi2yTgv7anNElSJyp7A+AWWs8bH8yGsStHktTpyl6OezDwceBYYO8d2zPzbW2qS5LUocqeqvoacD+t54xfCDxI6wY+SdIepmxwHJSZlwPPZ+YtmXk28OY21iVJ6lBlr6racYf4IxHxbmAzrcFySdIepmxw/E1E7A/8FfCPwH7Af29bVZKkjlX2qqobi8WtwPz2lSNJ6nRlr6o6itbluDMHfiZH9wRASVIDlT1VdR2tmWpvwDvFJWmPVjY4fp2ZX2hrJZKkRigbHJ+PiAuA7wK/eSh3Zt7RlqokSR2rbHD8NvCHwNt46VRVFuuSpD1I2eB4H/DazHyuncVIkjpf2TvH76L11D9J0h6u7BFHN3B/RNzOy8c4vBxXkvYwZYPjgrZWIUlqjLJ3jt/S7kIkSc0wbHBExFO0rp7a5S0gM3O/tlQlSepYwwZHZu47XoVIkpqh7FVVkiQBBockqSKDQ5JUSS3BEREHRsT3IuKnxd9XD9HuxYi4s3gtG+86JUm7quuI4zzgXzNzNvCvxfpgfpWZxxcvbzaUpA5QV3CcDlxZLF8JvLemOiRJFdUVHN2Z+QhA8feQIdrtHRGrI+LWiDBcJKkDROZg9/eNwRdH3AwcOshb5wNXZuYBA9o+mZm7jHNExOGZuTkiXgt8Hzg5MzcOsb+FwEKA7u7uuUuXLi1VZ39/P9OmTSvVtknsV+eZc8fHeXHSq1j3xgtftr3JfRqO/WqW+fPnr8nMnlKNM3PcX8B64LBi+TBgfYnPfBk4o8z3z507N8tauXJl6bZNYr860Jfennnl6btsbnSfhmG/mgVYnSV/w+s6VbUM+FCx/CHg+p0bRMSrI2KvYrkLeAtw77hVKEkaVF3BcTHwjoj4KfCOYp2I6ImIy4o2vwWsjoi7gJXAxZlpcEhSzcpOqz6mMvNx4ORBtq8GPlwsr6L1yFpJUgfxznFJUiUGhySpEoNDklSJwSFJqsTgkCRVYnBIkioxOCRJlRgckqRKDA5JUiUGhySpEoNDklSJwSFJqsTgkCRVYnBIkioxOCRJlRgckqRKDA5JUiUGhySpEoNDklSJwSFJqsTgkCRVYnBIkioxOCRJlRgckqRKDA5JUiUGhySpklqCIyLeHxH3RMT2iOgZpt0pEbE+IjZExHnjWaMkaXB1HXHcDfwB8IOhGkTEJOCLwLuAY4CzIuKY8SlPkjSUyXXsNDPvA4iI4ZqdAGzIzAeKtkuB04F7216gJGlInTzGcQTw8ID1TcU2SVKN2nbEERE3A4cO8tb5mXl9ma8YZFsOs7+FwEKA7u5u+vr6ypRJf39/6bZNYr86z5xtW3lx0nOs26n+JvdpOPZr4mpbcGTm20f5FZuAGQPWpwObh9nfpcClAD09Pdnb21tqJ319fZRt2yT2qwNt2B+mTN2l/kb3aRj2a+Lq5FNVtwOzI+KoiJgCnAksq7kmSdrj1XU57vsiYhNwInBTRKwoth8eEcsBMvMFYDGwArgP+GZm3lNHvZKkl9R1VdW1wLWDbN8MnDpgfTmwfBxLkySNoJNPVUmSOpDBIUmqxOCQJFVicEiSKjE4JEmVGBySpEoMDklSJQaHJKkSg0OSVInBIUmqxOCQJFVicEiSKjE4JEmVGBySpEoMDklSJQaHJKkSg0OSVInBIUmqxOCQJFVicEiSKjE4JEmVGBxSm11yy0ZWbdzysm2rNm7hkls21lSRNDoGh9Rmx03fn8VXrWXbs88DrdBYfNVajpu+f82VSbvH4JDabN6sLpYsmMPPHnuaBx9/hsVXrWXJgjnMm9VVd2nSbplcdwHSnmDerC5uOPpsrv/JY3yw90hDQ43mEYc0DlZt3MIFG17HMb0f4F9ue2iXMQ+pSQwOqc12jGksWTCH//HOo1myYA6Lr1preKixagmOiHh/RNwTEdsjomeYdg9GxE8i4s6IWD2eNUpjZd2mrS8b09gx5rFu09aaK5N2T11jHHcDfwD8c4m28zPTf5qpsRb93qxdts2b1eU4hxqrluDIzPsAIqKO3UuSRiEys76dR/QBH83MQU9DRcTPgCeBBP45My8d5rsWAgsBuru75y5durRUDf39/UybNq1i5Z3PfjXHROwT2K+mmT9//prMHHLoYKC2HXFExM3AoYO8dX5mXl/ya96SmZsj4hDgexFxf2b+YLCGRahcCtDT05O9vb2ldtDX10fZtk1iv5pjIvYJ7NdE1rbgyMy3j8F3bC7+PhoR1wInAIMGhyRpfHTs5bgRMTUi9t2xDLyT1qC6JKlGtYxxRMT7gH8EDgb+A7gzM38/Ig4HLsvMUyPitcC1xUcmA1dl5t+W/P7HgJ+XLKcLmIhXbdmv5piIfQL71TRHZ+a+ZRrWOjjeCSJiddkBoSaxX80xEfsE9qtpqvSrY09VSZI6k8EhSarE4Cgu4Z2A7FdzTMQ+gf1qmtL92uPHOCRJ1XjEIUmqxOAAIuLTEbGumIX3u8VlwY0XEZ+LiPuLvl0bEQfUXdNolZ1ZuSki4pSIWB8RGyLivLrrGQsRcUVEPBoRE+a+q4iYERErI+K+4v+/v6i7prEQEXtHxI8j4q6iXxeW+pynqiAi9svMbcXyucAxmbmo5rJGLSLeCXw/M1+IiM8CZOYnai5rVCLit4DttGZWHnKesyaIiEnAvwPvADYBtwNnZea9tRY2ShFxEtAPfCUz31B3PWMhIg4DDsvMO4obk9cA750A/60CmJqZ/RHxSuDfgL/IzFuH+5xHHMCO0ChMpTWpYuNl5ncz84Vi9VZgep31jIXMvC8z19ddxxg5AdiQmQ9k5nPAUuD0mmsatWI+uSfqrmMsZeYjmXlHsfwUcB9wRL1VjV629BerryxeI/7+GRyFiPjbiHgY+G/Ap+qupw3OBr5ddxF6mSOAhwesb2IC/BhNdBExE5gD3FZvJWMjIiZFxJ3Ao8D3MnPEfu0xwRERN0fE3YO8TgfIzPMzcwbwNWBxvdWWN1K/ijbnAy/Q6lvHK9OnCWKwB9JMiKPdiSoipgHXAH+505mKxsrMFzPzeFpnJE6IiBFPL9b1BMBxV2G23quAm4AL2ljOmBmpXxHxIeA9wMnZkAGtsZhZuSE2ATMGrE8HNtdUi0ZQjAFcA3wtM79Vdz1jLTP/o3hG0imMMKHsHnPEMZyImD1g9TTg/rpqGUsRcQrwCeC0zHym7nq0i9uB2RFxVERMAc4EltVckwZRDCJfDtyXmf+r7nrGSkQcvONqy4h4FfB2Svz+eVUVEBHXAEfTulrn58CizPxFvVWNXkRsAPYCHi823dr0q8WGmlm53qp2X0ScCvwfYBJwRdkZoDtZRHwd6KU1i+wvgQsy8/JaixqliHgr8EPgJ7R+JwD+Z2Yur6+q0YuI44Araf3/9wrgm5l50YifMzgkSVV4qkqSVInBIUmqxOCQJFVicEiSKjE4JEmVGBzSGIiIF4vZle+OiBtGMxNxRDwYEV1jWZ80lgwOaWz8KjOPL2aDfQL487oLktrF4JDG3o8YMFlhRHwsIm4vnoty4YDt10XEmuI5CAtrqVTaDQaHNIaKZ2ycTDF1SPFMlNm0plA/HphbPK8C4OzMnAv0AOdGxEE1lCxVZnBIY+NVxdTUjwMHAt8rtr+zeK0F7gBeTytIoBUWd9F6VsqMAduljmZwSGPjV8XU1K8BpvDSGEcAnynGP47PzNdl5uUR0UtrQrkTM/ONtIJl7zoKl6oyOKQxlJlbgXOBjxbTcK8Azi6e40BEHBERhwD7A09m5jMR8XrgzbUVLVW0xzyPQxovmbm2OAV1ZmZ+tXhO+o9aM3PTD3wQ+A6wKCLWAetpna6SGsHZcSVJlXiqSpJUicEhSarE4JAkVWJwSJIqMTgkSZUYHJKkSgwOSVIlBockqZL/DwPqZ+IoWU+vAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f2b1876ee10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "control.rlocus(ss[0,0], klist=pl.logspace(-2,1,1000));\n",
    "pl.grid()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simulation of Sympy Model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAAyCAMAAABCrGH/AAAAM1BMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADxgEwMAAAAEHRSTlMA74lUMhBmzavddiJEmbsg4DnoagAAAAlwSFlzAAAOxAAADsQBlSsOGwAABqlJREFUeAHtXNuWqyoQRBQliG7+/2tPdTe3ZGLGRCdL19GHGTQ0FGXT3CpRTcCl1XX9IQOaOG5U03ZdZ/6wnqtoZUBxC6q7s3Ix9hl5P+bkURPdiake54pVe3iuz0z1IOyaG/03Q8X7IZMnprqfhNE28H/tDklwAXViqtmZ0ZJJ/Nk0pVWHTJ2X6vEWCb3FiergD8lwBnV0qvshNJhn2HBrM2ZJOI4f3s5hthw7dAwoD/kOc3t0qpVqLJH1c9Cb40wvhmql+p95DkMzATk+1S5geVXNoBN9yc1jqFaqTREl5TjY/+NTbQLc98m+QRNDcwrVqpOZyMEILnCOT7Wy8NZlqk1Io+FFdXmtH6Z88Lzs9kNolQ5piRgDiCdf5h0cfwWQDxkuZjcb5xaDU1PesInDosZg2PLDa1gslH2a6uHNfJlmKlO+OLebrDKyTpyeRJlPq/wLuxPEamVyZOhpNhKvGC/GuY+zviE7fMpyrP9noDozZiZXzZ3Twjx+nFePOfvBEqeiejZqLkvCtN0UGb22m/Zzrb7BZGMKc44Tc71FfW2i7sf0z5LujgbueP+Z9wBPThVAHvmqluvXgdcjOf/n+1N79ble3EX1197XRfVF9dcY+FpFh/LqA80oCpT35zbFVtXGv1HdzU1Zn/2xA9zNk7dLaGzabl2EbeYmb6/cZ6qhvIuktlWV8W9UK5W21e6x/MXdzhIaG7ehXkEdFvxIoHwm5llqxq9U+2oz7RXo7Z+lPY1vSmgW/ChC+QjJYjN+pVovdLHt1D6WkHbqviihWfKjCOUjJIvN+JXqpS72SNTm+7wJms5lvyChWfCjBOUTJMlW/TB+RfU4addJF0Nycg325b224+gcRC5eu5lOUL2exw5PaNPeaOcm7Lchl/G63tbEB4625HaV0KDQHoUWBKVe4yY6mylwcYOTscm5thotyY9aZwmYQRNG7xgzq3nWi3liXSPxo5aVQC+o9g1ImzlUc3KgA1Ot5hltMMGBZoiz6YmlE1WNFnjqPCPOnZBraNWQVXQj5RhkW38/CQ25Pe1fFwSlXq081T4luLzfShjHUEZCAz9yoyf1gwHq/ubUjUDGY8u1Yp5JiHAsjlhWAi1TbUjAJdrDmATQrhUwI3+jo0fhIygl/RHSYwNCDQ5aKRee+byz3NCxn5Vzwd0kNCM1DV93KAhKvXhGgbbAJSgtu001EuKBM8pQOzU8UuPzlt5JPIxfKebBNJpzztxdom0SbaLWNNgVqs1wyxcd000MjCNOTJI3w58paghqmrkCIoMHTtv06NF4kHIhI19iH28+kNA8QpOSTGgmwlIQ/Ivo5FkFl3wCkjTiogvs4HSPJmYPRyFqYK6QiGqeFG1rMc8TKDBtKFjxn2SbQ3UxLlQj893FwITYmGQnlV418auSwmUYnxsVWF3HheS+x3dsnwvfTULj5yCqkIwALpCqEWeIDyyFhY4ChXLJyZC+WUfRMF0hpYXq9WIerr8LzE98TU+MF6k2HCIYepWMEUXemYeD4JXyMG7gINXXxHLfo2aIfWqQ2klCgxiBsYx9NCMo3ZaiGnmqAAHJBiGO/Lk6nCQ6xiyOwqtgrpBHgsB6MQ87noTqFHyeGL+gmkLrPCB60egRk+gmFHblnaEZ1F6eDvYYfvkj+A0yxCRSdLE7IerIHZwp9ttNEhoZdNmLM4JSL2EmpAyEPFxHB69CNYWVlK81irlC7E7D4noxz0zdGaGagueyEmiRagUeMCZbAlwlxQsk9g2aW0OduKOvibFywNBkL/uKCdR1IYtBo7IwaR8JDR8mOu7zGUGuF12H1avyAKRhUmjIqy27tsDiTRK4Mlro8QZmMC87RaLmWS/mmcCAx14RvdtlJdAy1cY63Y2WJnVVUoQujqNyO/GsNniZTNOc0on+pcjnaCaC+nXfp0CI4nK43CKhwey+d9xREJEiglKvsoxEHoy2p4weLYqLFobFDgjM8EW00GF6HNU70Wq1mMdY7drO8sIhQniiBFqmmihacy2suNaY7iSheQdBFapfIEyLa8mSF4AvLMpH97aqGG+nesPKfScJzUoENOiViUHh5kkqbRnJR/Wq90nmh0f3tqoYb6ZaxsyH6lbd7iWhWYmAFiiKB7AV6GpVybtintq2/jrlVqrbW6BBc5frbk/9DvCr4tciaG3vKiHrqyIxPKbVDLKtRhKLrG1r461Uv0b85qcLJ0VvlrJL9gKlPrNaV3SxfevAa13hV64VDJBXXz9SsYKobVnSj1RcP72yjccV1vLTK/8BCv42bBTkmxkAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$\\left[\\begin{matrix}\\operatorname{v_{x}}{\\left (t \\right )}\\\\- damper.c \\operatorname{v_{x}}{\\left (t \\right )} - spring.k x{\\left (t \\right )} - u{\\left (t \\right )}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡                vₓ(t)                 ⎤\n",
       "⎢                                      ⎥\n",
       "⎣-damper.c⋅vₓ(t) - spring.k⋅x(t) - u(t)⎦"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAOoAAAAVBAMAAAC6fm+xAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAZpkQ3Ynvq81UMrtEdiLw+n06AAAACXBIWXMAAA7EAAAOxAGVKw4bAAADgklEQVRIDb1WPWgUQRT+NnNr7rzcGYOKlV6KoHaLaLQJXhEhSgwhaSwMxlYEg9iJ8RAl4O+BksLCLLYWiVikEPQEURu5w0S0Ca5gp8VpEI1R8Hszk3PuNhc7X/HezPe9d9/OzJvdA+B14n/avTlRUwcbNT8I4EWNMDB/RWNb++OUQbLtErfF6UxfKKB3qQQ8l6FXpXMtO6Nn4y5mxscstBJjGZ800hrFCL8E1bU5BCNaCqT9UTrXkhGwJ8J+FzNjUZsMgKaqN7h1A1CFWCnVkkVSournScdUOwhOEa8w1tsxqDMX1lD1ckCqCozUl3FGtS3AUaM6KoA41wY42cmHHmWsN1nj+zVUE5HZv+76Ms6oOsu1BHqtowKIcyyb50Q67JoDmuG/VD8yLTEIJEuNpVT7DpSLrqp38wQuHLKZmRzazi2fB4Y0kBiT0PP285uKPs/V12r418CLJ7MV+DmpMZVqYfzVCNeqvlK1UlNtLQDrxpCURLF0uzkcbNXT7gKDGvSOJArNVS3/gKnsCHgFOuhK7MWmSh9Vs9+A6UFR1ZdmMgQyi9CXVJITRdvct2VmzA/8xVTUXNXysl3sCKTypkz8MI4Hd0WVazWquAjcP0lKfQGPw9i6kjkcTFuAQSE5w9D0XC1/ijnSEeurdNYCXOXI3WG0HrGvirP1qgdKzHRUgQ0VIk1VLU9V9Yt5rirwmwh3lt00pbsJEwTkXFF+GuHlkr/MoexwmavDM05qNh1yuIaq5rnDbVUEaMvX6vgIPzmh6iNge8BozlW/KhJnSF1vCemlm47iDmw3EeFBBe+RCZqrWp7d5BdSoe0mXYlbmSrmRJVviZ0SzVXV99VbYk73LkmUtt+huKFdMjOdWC4+kXb7u1b1Q5MrzvJD7MOxHrrBWmXbUrrKbwLVEkV1skFV5ZnnLUqybM+njRyYD4y+dfsWejq5BVTdcXqiXWJvILkrZnmeSfYB85KhELpSnZ97zR+jqhp5Gzao6vLd73QY0B7ZvImOl7WKMfIexS1hwY8xSnZWTKLe3L9vxIfpnDAd4khWTHS8o+o7cG3o5cxw1be/pmqqte/rlmrLckguHemEJl864ajO01vFHmtMFWKUu9YWTV9uyLFfdfMLdZyzVp7eKma+6vaxXd5VfV4SJtUp3jH9cvQiB7HD+cN6MN8fpwzS/B9MbygZmd4S7g0DfwAp5fV1Wcy6qgAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$$\\left \\{ \\operatorname{v_{x}}{\\left (t \\right )} : x[1], \\quad x{\\left (t \\right )} : x[0]\\right \\}$$"
      ],
      "text/plain": [
       "{vₓ(t): x[1], x(t): x[0]}"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import sympy\n",
    "ss_sub = {}\n",
    "ss_sub.update({model.x[i]: sympy.DeferredVector('x')[i] for i in range(len(model.x))})\n",
    "#ss_sub.update({model.y[i]: sympy.DeferredVector('y')[i] for i in range(len(model.y))})\n",
    "ss_sub"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "res = model.simulate(t0=10, tf=20, dt=0.001)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0,0.5,'x')"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEKCAYAAAAFJbKyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAIABJREFUeJzsnXdYVFfawH9nqNJBYECKIAgCduw9tmia6Ympm83GJFuyu9mSsiW72exustl8W7LZ9Gx6NcYYS4yJPXajCIoNlI4oSO/M+f64YIwiMDBz7wXP73l4Bph7z/teZpj3nrcKKSUKhUKhUHQVi9EKKBQKhaJ3oQyHQqFQKOxCGQ6FQqFQ2IUyHAqFQqGwC2U4FAqFQmEXynAoFAqFwi6U4VAoFAqFXSjDoVAoFAq7UIZDoVAoFHbharQCziA4OFjGxMR069yamhq8vb0dq5DJUdfc97nYrhfUNdvL7t27T0kpQ7pybJ80HDExMezatatb565fv54ZM2Y4ViGTo66573OxXS+oa7YXIUROV49VriqFQqFQ2IUyHAqFQqGwC2U4FAqFQmEXynAoFAqFwi6U4VAoFAqFXSjDoVAoFAq7UIZDoVAoFHbRJ+s4FBch1Sch6ysozwMPH4gYAxGpYFH3RgqFo1GGQ9G7qT4JX/0B9r4HsuW7z4UkwZzHIWGuIaopFH0VZTi6QUZBBW9vy2FLViknqxoI8HJjbEwQt4yPZnxsEEIIo1W8OMjfBe/fCrWnYNw9MPJWCE6A+gpt97H5H/DuDTBuEcx7EiwuRmusUPQJlOGwg/qmFp5YcYB3tufi5ebCtIQQIgL6UVLVwKYjJ1mWVsi8lDD+fM1Q+vt4GK1u3yZ3G7x1DXgHw6L1EDbs2+fcPGHkLTD0Ovjyj7DtOagqhuv/By7qLa9Q9BT1X9RFymsb+d7/drI3r5zvT47lZ3MG4+fpdub5+qYWXvv6GP9cc4Sr//s1b9w1jkEhPgZq3IcpOQjv3AC+4fD9z8EntP3jXD1g3l/APwJWPworHoQr/wVqR6hQ9AgVOewCVfVN3Pm/nRworOSF20bz+yuTv2M0ADzdXPjhjHg+um8itQ0t3PjiNvLKag3SuA/TUA0f3qEZhTs+vbDROJuJP4Kpv4Rv3oAdLzlfR4Wij6MMRydIKXnwwzQyCip47tbRzBsa3uHxI6ICeH/RBJpabNz52g4qapt00vQi4fOH4dRhuO4VCIjq+nkzfwsJ8+CL30FxuvP0UyguApTh6ISXN2Wz5sAJHr0siTnJ1i6dM9jqyyt3jiG3rJZff5yGlNLJWl4kHNsEe96CyQ/AoBn2nSsELHgO+gXCkkUIW7MzNFQoLgqU4eiAwyeqeHr1IealhPH9yTF2nTs2JoiH5g1h9f4TvL091zkKXkw0N8Dyn0PAQJj+cPfW8A6GK/4BJQeIzF/mWP0UiosIZTgugM0meWRJOt4ervz5mqHdSrG9e0osUwcH8+TKTIoq6pyg5UXEzleh9Ahc/gy4e3V/nSGXQeLlxBx/H8qVQVcouoMyHBfgo9157M45zW8vT+52aq3FIvjLNcNokZI/LjvgYA0vIhqqYNPfIXY6DJ7T8/XmP4mQNlj3l56vpVBchCjD0Q51jS3835rDjI4O4LrRET1aKyrIiwdmDebz/cVsOHzSQRpeZGz9L9SWwqzHHLNeQDT5kZdD2vtwQhl0hcJelOFoh9e3HOdEZQMPz09ySBX4D6YMIjrIiydXHcRmU4Fyu6ivhK3/gSFXQGSqw5bNjb4OPHxh7RMOW1OhuFhQhuMcapokz68/yiWJIYyLDXLImu6uFn4xN4HMokqWpRU6ZM2Lhm/ehIZKmPoLhy7b7OYHkx6AQyugKM2haysUfR1DDYcQ4jUhRIkQIuMCzwshxL+FEEeFEPuEEKOdrdP6vCYq65t5cE6iQ9e9cvgAUgb48cyaQzS12By6dp+lpQm2PQ8Dp0CEE176cfeAuy98/S/Hr61Q9GGM3nG8Dszr4Pn5wODWr0XA885UpqG5hTU5zUyO78+wSH+Hrm2xCB6ck0BeWR2fqV1H19i/FCrzYdJPnLN+vwAYcxfs/wTKsp0jQ6HogxhqOKSUG4GyDg5ZALwpNbYBAUKIjku3e8Cnewspb5DcOy3OKetfkhhKotWXFzZkqVhHV9j+PPQfDIOd2BZ9wg/B4gpb/uM8GQpFH8PoHUdnRAB5Z/2c3/o7h2NramToquv4hdcqpg4OdoYILBbB/TPiOHyimrUHS5wio89QnA4Fu2Hs3c4dxuQXDsNvgr3vQN1p58lRKPoQZu+O215KU7u36kKIRWjuLKxWK+vXr7dLUF2zZLQs5xLPw2zYsMFePbuMr00S3E/w1LJvcC3p5zQ59lBdXW3338vZDD78IuHCjS1VkTQ7Qbezr9nHMpoxzW9xdPGfyI+6yuGyzIAZX2Nno67ZiUgpDf0CYoCMCzz3IrDwrJ8PAeGdrZmamiq7xfu3yponh3TvXDt4ZVO2HPjQcpmeX+50WV1h3bp1RqvwXRpqpPxLlJSL73aaiPOu+eXZUv5rlJQtLU6TaSSme411QF2zfQC7ZBc/t83uqloG3NGaXTUBqJBSFjlNWthw+tUVaa27ncj1qZF4ull4a2uOU+X0Wg58Cg0VMPpO/WSO/QGUZcEx5+02FYq+gtHpuO8BW4FEIUS+EOJuIcR9Qoj7Wg9ZCWQDR4GXgR86VaGwYQgklDi3mti/nxtXj4zg07QC1Xa9Pfa+A0GDIGaKfjKTF4BXf9j5in4yFYpeitFZVQullOFSSjcpZaSU8lUp5QtSyhdan5dSyh9JKeOklMOklLucqlDb+NHifU4VA3D7xIHUN9n4aHde5wdfTFQWwvHNWsBaz0l9bp7azPJDq6BaJS4oFB1hdleVvvhF0OTqC0XONxwpA/xJHRjIO9tz1byOs8lYAkgYer3+skfeCrIF0j/SX7ZC0YtQhuNshKDaJ1a3CXELx0Vz7FQNu3NUGugZMhZD+AgIjtdfdugQGDAK9r6nv2yFohehDMc5VPsM0mIcLc6fEDd/aBhe7i58tCvf6bJ6BaVZULgHht1gnA4jboET6Wq8rELRAcpwnEO1Tyw010PpUafL8vZw5fJh4axIL6K2UY0yJX0xICDlWuN0GHY9WNzUrkOh6ABlOM6h2idW+0anO87rUyOpbmhm9f5iXeSZmv1LYOAk8HdKc4Cu4RUEifMg/UOtyaJCoTgPZTjOodYrAlw8dMmsAhgXG0R0kJdyV506CicPQpIJKrdHLISak5CtajoUivZQhuMcpMUVQpN023EIIbg+NZItWaUUlF/Ec8kPrdAeh1xmrB4A8bPBw0/rmqtQKM5DGY72CBumGQ6d0mSvHqm5Zlbsu4jbrR9cAWHDISDaaE3A1QOGXA4HP4PmRqO1UShMhzIc7RE2HGpPQZU+cYfo/l6MiPTnszTndVMxNVUnIG+HNh7WLKRcA/UVkL3OaE0UCtNh9u64xnB2Bbmf08Z/fIcrRwzgiRWZHD9VQ0ywty4yTcPhVYDU7vLNwqBLwDNAK0hMuNRobRTdpbIQ0t7XuhGUZUNTLXgFa7VCCXMh8XJwdTday16HMhztYU3RHov36fahcdmwcJ5YkcnyfYX8eOZgXWSahoMrIWDgt393M+DqDklXwP5Poalea0nSBTIKKliWVsiu42Xkna6jxSYJ9fVgWIQ/c5KtzBwSiquL2ug7neqTJBz6D2xYq3UDCE2BiFRw99Jayhz+HNLeBZ8wuOQRGHU7WFyM1rrXoAxHe3j6aU32dCwCGxDQj7ExgXyWVnRxGY6Gasher3Wn1bM3VVdIuQb2vA1ZX3W6GzpQWMlfV2Wy6cgp3FwEI6MCmJkYiouLoKi8jtX7i/lodz5RQf34xZxEFowcgDDb9fYVjqyBJfcQVl+pzZUffx8ExX73GJsNstbCpr/DZz+Fve/C9a+Bf6QxOvcylOG4EG0Bch25YvgAHlu2n0PFVSSG+eoq2zCObYSWBq12wmzETod+QVp21QUMh80meXFjNs98cQhfT1cevWwIN46JIsDru+6PphYbX2WW8N/1R/nZB3v5+Jt8/nHTSIJ9PPS4kouHr/8Fax4D61B2DV3EuPl3tH+cxQKDZ0P8LNj3Aaz4BbwwFRa+B9ET9NW5F6L2zBfCOhTKjkFjjW4i5w8LwyJg+cWUXXX0S3D3gSgT/rO6uEHiZXD4i3aLARuaW/jJ+3t46vODzE2xsu6XM1g0Le48owHg5mJh3tAwlv5wMo8vSGHn8TKufHYz6fkVelxJ30dKWPcXWPN7baf4gzXUenchQ08IGHEzLNqgFX++eTUc+dL5+vZylOG4EKHJgISSg/qJ9PVkfGx/VmVcJFXkUsLRNdqdvVkDlEMu14ZKHd/8nV/XN7Vw9+u7WLGviEfmD+G5W0a3azDOxWIR3DExhsX3TcIiBAtf3sau42XO0v7iYdt/YcNTMPI2uO4VcLNzLHNwPNz1ufb4/sLzXm/Fd1GG40K0BWpPZOgq9tIUK0dLqsk66dwphKag9CiU52ruArMyaAa49tPqTFppsUl++v4eNh89xdPXD+fe6XF2xyuGRvjz8f2TCPX14I7XdqgOyT0hczms/o3WdeCqZ7sf5PYJgTuWQWAsvLcQivX93+9NKMNxIQIGai4UJ08DPJe5KWEAfLH/hK5yDeFoq0sgfraxenSEu5dm2A6uOFMQ+sfP9rN6/wkeuzKZG8ZEdXvpMH9P3r93AqG+Htzz5i6On9LPLdpnOHUEliyCiNFwzYta7KIneAXB7Uu0//33FkKt2g22hzIcF8Ji0VqPnNivq9gBAf0YHul/cTQ9PPolBCdA4ECjNemYIZdDVSEU7uGTPfm8uTWHH0yJ5a7JsZ2f2wmhvp68ftc4pJR873871Chhe2hugMXf19ycN72tGXlH4B8JN78N1cWw5B4tA0vxHZTh6Ahriuaq0nlC39xkK3vzyimuqNdVrq401Wl+ZDPvNtpImAfCQuk3S3l0SQbjYoN4eP4Qhy0fE+zNy3eMIf90Hb9cnKYmQnaVrx7Xaq0W/Bf8Bjh27YhUmP+UdnOz+RnHrt0HUIajI0JToO60bq1H2ri01V21JrMPu6uOf63NPTFzfKMNryBs0ZOo2rMUbw8X/rNwlMOL+MbEaMZozYETvLr5mEPX7pPk7YStz8GYu53XGDP1Lhh6Hax/Updx0r0JZTg64kyAXF93VXyoD4OCvfmiL7urjn4Jrp4wcLLRmnSJzS7jiLHl8H+z/Qj161oVub3cPSWWuclWnlx1kP2FKk33grQ0wfKfgW84zP6D8+QIAZf9Hbz6wyf3aa4xBaAMR8dYk7XHEn0NhxCCuSlhbM0q7bs+7+z12tAme9MmDeBQcRWPHdLiMNNsO5wmRwjBU9cNJ8DLnV99tI+mFuVbb5etz2ku5Mue1ro8OBOvIC1Tq2Q/bHzaubJ6EcpwdES/QPCL0H3HATA3xUqzTbL2UB90V1WXwMlMrX7D5Egp+d3SDMo9BtAcnASHVztVXqC3O3++ZigHiip5fn2WU2X1SqqKYcPftOaESTp1U064VBvutfmfWhaXQhmOTrGmGGI4RkYGEOrrwZcHSnSX7XSObdQeB5nfcKxML2bH8TJ+dekQXBPnQu5WqK90qsxLU8JYMHIAz649wpETVU6V1etY9xdoaYRLn9BX7pzHwc0LVv5S92QZM6IMR2eEJsPJQ7rPn7ZYBDOHhLLx8Mm+57I4tgE8/bW5JyamvqmFv6zMJCncj5vGRsHgS8HWrMuMjseuTMHL3ZXHlu1XWVZtlGTCnrdg7N1aE1I98QmFmb/VXKwHluor24Qow9EZ1qFgazJki3rJkFCqGprZ2ddaUhzbCDFTTd/G+uWN2RSU1/H7K5JxsQiIGg8e/nDkC6fLDvJ255dzE9iSVcrK9D6cJGEPax7TCvOm/doY+WPv1m52Vv9GSye/iFGGozPaAuQGuKumxAfj7mJh3cE+5K46nQOnj0PsNKM16ZDy2kZe2pjN3GQrE+P6a790cYX4mVrbbh2Kwm4ZP5DkcD+eWHGA2sZmp8szNXk74MhqmPJz8O5vjA4WF5j3JFQWwLbnjdHBJBhqOIQQ84QQh4QQR4UQD7fz/PeEECeFEHtbv36gu5L9B4PFTffMKgBvD1fGDwriq75kOI5v0h5Nbjhe2phNdWMzD85N+O4Tgy+F6hNQnOZ0HVwsgj9dnUJRRT0vXOyB8g1PaWmx4+81Vo+YyVpB6OZ/XtTtSAwzHEIIF+A5YD6QDCwUQiS3c+gHUsqRrV+v6KokaO0MQhIN2XEAzBoSSvbJmr7Tx+jYRvAOgRDHVV47mlPVDby+5ThXDB/AkLBz0j3jZwNC23XoQOrAIC4fHs4rm49RUtWHOwl0RMFure5n4o/B3QRjlWf/ARqrYOPfjdbEMIzccYwDjkops6WUjcD7wAID9bkwoclwQt9mh23MHGIFYG1f2HVIqRmO2Gnmm/Z3Fi+sz6K+qYWfzW5nEqNPiNZQz8lpuWfzy7mJNDbbeParo7rJNBUbntbmv4+7x2hNNEKTYOStsPNlze3aDVpskvqmFhqaWxyrm04YOQEwAsg76+d8YHw7x10nhJgGHAZ+LqXMa+cY52JNgfQPtfYj/QJ1FR3d34v4UB/WHSrh+1N63lTPUEqPQlWRqd1Up6obeGtbDlePiiAuxKf9gwbP1dpQ1JwC72Cn6xQb7M3CcdG8tyOX70+JJTbYBHfdelG0Dw6vghmPgoeJpmJe8iikf6QVBS54rtPDy2oa+TyjmE1HTnKgqJL81nn0AIFebgwK8WFcbBCzk6yMjg4w/VhhIw1He3+Zc/MOPwPek1I2CCHuA94AZra7mBCLgEUAVquV9evXd0up6urq884NKm1hOLBn9btUBKR0a92eEO/VyJqj1az6ch39XB3/hmrvmp3BgIJVJADbSjyo10FeR1zompccaaSx2cZYr7IL/k18K/uTiiTzs2c5EXaJcxVtJbWfjQ+F5KG3N/HDkfa3PNHrNXY0SQeeob9LP7Y1pdBsp/7OvuZ46xwi9rzLdvep1PcLa/eYomoby7Ib2VnUQrOEIE9BfICFYTGueLqCTUJZvSSvvIKXNpzm+fVZhHkJ5sa4MS3SFVeLff/vur3OUkpDvoCJwOqzfn4EeKSD412Aiq6snZqaKrvLunXrzv9lRYGUj/lJuf2lbq/bE7ZmnZIDH1ouV6UXOWX9dq/ZGXxwh5TPJEtps+kjrwPau+bq+iY5/A+r5T1v7Oz45JYWKf8WL+VHdzlHuQvwzOqDcuBDy2VGQbnd5+r2GjuSigIp/xgk5aqHu3W606+5olDKx0OkXPqj85+qa5S//SRdDnpkhUz+3Sr5+6Xp8kBhhbR18N6vqm+SH+zMlVf9Z7Mc+NByOemvX8nPM+z7n+/JNQO7ZBc/v42McewEBgshYoUQ7sDNwLKzDxBChJ/141VApo76fYtvuOai0nkaYBupAwPx9XRl7cFe3H5ESsjZomWlmHQb/uGuPCrqmrh3elzHB1osMHiOFrC16eejvnvqIHw9XS+eWMeOl0DajM+kuhB+4TDmLkh7D8q+7Wi8LbuU+f/cxDvbc7h1fDQbfn0Jf1wwlKRwvw5dUD4ertw4JoqlP5zE63eNxdfTlXvf2s29b+2itNpcDRYNMxxSymbgx8BqNIPwoZRyvxDicSHEVa2HPSCE2C+ESAMeAL5niLJCaC3WDcqscnOxMC0hhHWHTmKz9dIq4rJsqCmB6IlGa9IuzS02Xtl0jDEDA0kd2IU4VtxMqK+Agm+cr1wr/v3cuGtyLJ/vLyazyLltTwynsQZ2/Q+GXAGBMUZrc2Em/wyEC2x6Biklr2zK5paXt+HmIlh8/yQeXzCUYB8Pu5YUQjAjMZTPfjKFh+YNYd2hk1z+783szjFP+q+hdRxSypVSygQpZZyU8s+tv/u9lHJZ6/ePSClTpJQjpJSXSCkPGqasNUVreWDQNLCZiaGcrGrgQG/9wMjZoj0OnGSsHhdgRXoRBeV1ne822oibCQjI+sqpep3L9yfH4OPhyn/W9vFdR9p7UF8OE39ktCYd4xcOqd9Dpr3H3z5YwxMrMpmbHMaKB6YyOrpniTRuLhbunxHHkvsn4e5q4aYXt/HhLv1zg9pDVY53FWsyNFZDeY4h4qcmaNk7Gw6fNER+j8ndqhVwBSd0fqwBvPb1cQaFeDNrSGjXTvAKggGj4Ki+hiPAy53vTYphZUYRh/tqA0SbTavMHjBaa/Nicpon/oQWG4RmvMIPZ8Tx31tH4+3huLyjoRH+fPaTKUyM68+vF+/j2a+OGN6/TBmOrmIdqj0a5K4K9fUkZYBf7zUcOVs0N5UJ4xv78stJyyvnjgkDsdiTxRI/Cwp2aWnaOnL3lFi83Fx4tq/uOo6t11K3J/zQlO+Xs2lqsfHAyhI+aZ7Ebe4b+PXUEPveQ13Ev58br945lmtHRfDMmsP8aXmmocZDGY6u0lbpXGJMISDA9IQQvsk5TWV9LxvuVFUMp4+ZNr7x5tYcvNxduDY10r4T42Zpwdu2NvE6Eejtzh2TYli+r5BjfaWjwNnsek3bnSZf1fmxBiKl5NEl6axML0ZMfgA3W71WFOgk3F0tPHPjCO6aHMNrXx/jLyuNMx7KcHQVDx8IjDVsxwGa4Wi2SbYcLTVMh25xJr5hPsNxuqaRz9IKuWZUBH6ebvadHDkG3H11d1cBfH9yLG4uFl7elK27bKdSVQwHV8LIW8DVvqCy3vzrqyN8tDufB2YN5vr5cyBhPmx/ERprnSZTCMHvr0jmjokDeXnTMZ5efchpsjpCGQ57MGioUxujBwbi4+Ha+9xVOVvAzRvCRhityXl8uCuPhmYbd0yMsf9kFzdtGFXWWt2H+4T4enB9aiSLd+dzsspcqZo9Ys9bIFsg9S6jNemQj3fn888vj3B9aiQ/b2tNM/mnUFcGe99xqmwhBH+4MoWF46L47/osXv/6WOcnORhlOOwhNBnKsgzrxe/mYmFyfH82Hj5peHDMLnK3QtRYrS25iWixSd7ensO42CASw7rZziJuJlTkGTKv5Z6pg2hqsfHGluO6y3YKthbY/abWkqZ/F7PbDCCjoIJHPklnUlx//nrtsG9rMwZO1IL5W56FFue2wbdYBE9cPYzZSVb+uPwAX+zXd2aLMhz2YE3RfNonjcsKnp4QSkF5HVkne4lvu65c26VFmy8Nd8PhEvLK6rhj4sDuLxI/S3vMWusYpewgNtibeSlhvLn1ODUNfWBeR9ZaqMg19W6jvLaR+97eTbC3O88uHIWbyzkfoZN/qmVe6jAl0MUieHbhKIZHBvDA+3vYk6tfkoYyHPZgbe1TZVCnXIBpvS0tN287IE0Z33h3ex4hvh5cmtJ+n6EuERgDQXG613O0sWjaICrrm3l/pzny+3vErv9pLfeHXGG0Ju1is0l+/sFeTlTW89yto+nfXmFfwnxths/W53RxX/Zzd+HVO8cQ6uvJord2c7penzozZTjsIWgQuHoamlkVGehFXIh37zEcOVu0QVgRY4zW5DuU19tYd6iE60ZHnn/XaC/xs+D4ZmjWP9YwKjqQ8bFBvLopu3fPpq8sgsOfa+3KXd2N1qZdXvv6GOsOneR3VyQz6kLFfRaL1iKl8BvI36mLXsE+Hrxy5xhqG5r5954G6puc3wZHGQ57sLhoabkG9axqY3pCKNuzS3V5g/SY3K1aoZy7l9GafIevC5tpsUluHGNnCm57xM2EplrtWg3gvulxFFbUs2xvoSHyHcK+D7Sg+KjbjdakXQ4VV/G3zw8xO8nK7RM6cW2OWAie/rDtv/ooByRYffnHTSMZHGCxu6Nud1CGw16sKYa6qgCmJ4bQ0GxjW7bJ03Kb6rReTiZzU0kp2ZTfzNiYQAZdaOaGPcRM1XZVBsQ5AGYkhpBg9eHVzcd6V9JEG1JqLUaixkNwvNHanEdDcws/+2Avvp6uPHndsM5nZXj4wOg74cAyqMjXR0lgbkoYtyR54NrTHXQXUIbDXqwpWrO+auNcReNjg/BwtZjfXVXwDdiaTFf4tyvnNMW1khvHRDlmQQ8fiJ4AR40xHEIIvj85lgNFlWw/Zp5GeF2mcI+WcDJiodGatMs/vzxCZlElT103vOsNC8fdA0jY4byCQCNRhsNeQlvHopcYV8/h6ebChEFaWq6pafPxRo4zVo9z+GBnHp4ucPnw8M4P7ipxM+FEOlQZ0/r+6lERBHm789pm/XP6e8zed8HFA1KuMVqT80jLK+fFDVncPDaK2cnWrp8YEA1JV8Lu17VOv30MZTjsxQSZVaBVkWedrCGvzHlVqj0mf6eWUODd32hNzlDd0MyKfUWMC3fFy92BdSUGpuWCdjNx6/ho1mSeILfUxO+Jc2lugIzFkHQF9AswWpvv0NRi4+El6YT4evDo5Un2LzDhh1qH37T3Ha+cwSjDYS8+oVrKoIEV5KDFOQA2HjHprkNKyNthut3G8rRC6ppamBbp4GJE6zDwCjbMcADcNmEgrhbB/7b0ol3H4dVak8gRtxityXm8tvkYmUWV/PGqofa3owEtZhM+Era/YNg4BmehDEd3CE021FUFMCjYm8jAfmw4ZFLDUZ6rxYKixhqtyXf4YFce8aE+xPk7+K1vsWjuqqy1hn1IWP08uWL4AD7alU9Vb2mEmfYe+IRBnD6z27tKbmkt//jyMHOTrcwb2s06HyG0Xcepw5Bt3A2FM1CGoztYh0LJQV3Hhp6LEILpCSFsySqlsdmEdzNn4hvmMRzZJ6vZk1vODamRnWfGdIf4WVB7Cor3OX7tLvL9ybFUNzTz4S79snm6TfVJOPIFDL9RS3U3CVJKfrM0HVeLhT8uSOnZYinXgI9Vmy/Sh1CGoztYk6G57jtzho1gekII1Q3N7M7Rdx5El8jboTU2DO3hP54DWbq3ECG0QLJTiJupPRpURQ4wLNKfsTGBvL7lGC1mHzOc/hHYmrVOuCbi84xiNh05xa8uTSTcv1/PFnN1h7E/0ObTnzSmk60zUIajO5ggswpgUnwwrhZhzjhH/g6IGG2axoZSSpbuKWBSXH+sfp7OEeITCmHDDEvLbeP7k2PJK6vjy0xjMry6zL4PtBhAaDemJUOjAAAgAElEQVQCz06isdnGk58fJNHqy22dFfp1ldS7wMUddrzkmPVMgDIc3SFkCAiL4QFyHw9XxsQEmi/O0VQHxenavAqT8E1uOblltVw90km7jTbiZkLeNmgwbqzr3JQwIgP7mTs199RRKNqrualMxFvbcsgpreWRy4bg4qgKbJ8QGHYD7H1Pa/rZB1CGozu4e2lppgYbDtDajxwoqqSkst5oVb6lcK/mgjBRRtXSPQV4uFq6H+jsKnGztGs/vtm5cjrAxSL43qQYth8rI6OgwjA9OiRjMSAg5VqjNTlDRW0T//7qCFMHBzMjsYuz57vK+HuhqUabN9IHUIaju4QmG9rssI3pCW1puacM1uQs8ndojyYJjDc221i+r5A5yVZ8u5NWaQ/RE8DNy5CpgGdz49govN1deM2AIT+dIiWkL4aYKeDnwCLMHvKfdUeorG/i0cuc4DoLH6GNFtjxkqFJNY5CGY7uYh2qBccNrgpNCvclxNeD9YdKDNXjO+Tv1Mbs+oQYrQkAGw+f5HRtE9c4Kyh+Nq4eWu8qAwPkAH6ebtwwJorlaUWUVJloNwpa1lnpERh6ndGanCG3tJY3tuRwQ2okSeF+zhEy4T4tTf3QSuesryPKcHQXazIgtbRcA2lLy9105JQ5smikhLydptltAHyyt4BALzemJehkyOJnQVm24Vl335sUQ5PNxttbcwzV4zzSF4PFFZIXGK3JGZ5afRAXi+AXcxOdJyTxcvCP0uaS93KU4eguba1HDM6sAs1dVVHXRFq+CQJvFXlQXQxR5ohvVNY38eWBE1w5YkDP5250lbi29iPG7jpigr2ZNcTK29tzaWwxwU0FaMWRGUu0v5FXkNHaALA75zQr9hWxaNog52XcgZZhOO4eOL4Jio0dzdBTlOHoLgExWp2CCQLkU+KDsQjMkV1lssK/zzOKaWi2Oa92oz36x4F/tOFpuQB3T4mlrKaRrUUmGS2btw0q82HY9UZrAmhp2n9ecYBQXw/unT7I+QJH36HFwLb37oJAQw2HEGKeEOKQEOKoEOLhdp73EEJ80Pr8diFEjP5aXgCLBUKHmMJwBHq7MyIqwBxt1vN2gmu/b3dkBrN0TwED+3sxKkrHBnpCQPxMOLYRWoxt/TFhUBBJ4X58cbzJHLM60hdr74/Ey4zWBIBVGcV8k1vOL+YmOLbp5YXoFwgjboZ9H0GNiRJa7MQwwyGEcAGeA+YDycBCIUTyOYfdDZyWUsYD/wCe0lfLTrCmaIbDBP+Q0xNCSMsv53RNo7GKnCn8c3L2Uhcoqqhja3YpV4+McE6LkY6ImwWNVbqND70QQgjunhJLQbVk81GDP6hamuDAUkicr80wMZiG5haeXHWQIWG+XJ/qoNksXWHcvdDSALv/p59MB2PkjmMccFRKmS2lbATeB86Nli0A3mj9fjEwS+j+CdABoSlQVwbVxlfoTk8IQUrYZOSHQ1M9FO0zjZtq2d5CpHRii5GOGDQdhIvhabkAV44Ix89dGF8QmL0BaktN46Z6a2sOuWW1PHpZkuOK/bpC6BAYdAnsfNXwHWl3MdJwRAB5Z/2c3/q7do+RUjYDFYB5hjtYWzdIJnBXDY8MINDLzdi03KI0beKfSQzH0r2FjIgKIDbYW3/hnv7a38HgADmAh6sLs6JdWXfoJEdLqo1TJGOx9neJn22cDq2U1zby7NqjTEsI0S/b7mwm3A9VRXDgU/1lOwAjGwm1Z+LP9fl05RjtQCEWAYsArFYr69ev75ZS1dXVXT7XrbGSyUDWlmXk5Rvf3TPR38aXGYWsDTmNxY6NmT3X3BGReUuJB77Oa6bpRM/X6wmF1TYyi+q4ZYh7u9fmqGvuiIEug4jJe48tXyyjyd1JtQFdZFxQI59ZBH/+6GvuTOni+FMHYmlpYFLGUk6GTObQ5q26yOzoNX4vs4HKumbmhjj/fdAu0o1x/QbQtOZv7CkNdtiyeryvAS2rwIgvYCKw+qyfHwEeOeeY1cDE1u9dgVOA6Gzt1NRU2V3WrVtn3wlPJ0i55N5uy3Mki3flyYEPLZfp+eV2nWf3NV+I92+T8h/DHLNWD3lm9UEZ+/ByeaKirt3nHXbNHZG3U8rH/KTc95HzZXXCunXr5K8/SpOJv10pT9c06K9Axifa3+LoWt1EXug1Pn6qWsY/ukI+tDhNN13aZdsL2t8kb5fDluzJ+xrYJbv4+W2kq2onMFgIESuEcAduBpadc8wy4M7W768H1rZeoHloC5CbgKkJ2p2LYdlV+btM4aaSUrIsrZCJcf0JdWZefmcMGKVl0Rg4FfBs7poSQ32TjXd35OovPGMxeIdC7DT9ZZ/DU58fxM3FwoNzEoxVZMRCcPftlam5hhkOqcUsfoy2q8gEPpRS7hdCPC6EuKr1sFeB/kKIo8CDwHkpu4ZjTdb67LcYnycf6utJygA/YwxHRT5UFZqi8G9ffgXHS2u5asQAYxWxuMCgGZrhMMH9zpAwP6bEB/PmlhyaWnQc/lVfAYe/gKHXGj6waXdOGSvTi7l3WpyxNxUAnn4w6jbY/wlUFhmri50YWschpVwppUyQUsZJKf/c+rvfSymXtX5fL6W8QUoZL6UcJ6XMNlLfdglN0VLryrKM1gTQsqu+yTlNpd6jQ/PM09hwWVohbi6CeSkmaKAXN0sLgpqgISZoBYHFlfWsTNfxg+rgCu1/ZKix2VRSSp5YkYnVz4N7psUaqssZxi/Smh7ues1oTexCVY73lLZCtxPmaCEwPSGEZptky9FSfQXn7wRXT22QkYG02CTL9xUyPSEUfy/ja0nOTAU0QVouaO+PQSHevLb5mH4FgRkfQ0C04fNZVqQXsSe3nF/MTdSn2K8rBA2ChHma4WgyWTPKDlCGo6eEJGr5+ifMcUc5emAgPh6u+rur8ndqPn2DC/92HCvjRGUDC0Ya7KZqwz8CQpJMkZYLYLEI7pocS1p+hT4jh2tKIXu9NnfDwBKshuYWnvpcK/a7bnSkYXq0y4T7tFn1GR8brUmXUYajp7h6QP9407gi3FwsTI7vz4ZDJfrdUTY3aDUcJnFTebm7MDvJarQq3xI3E3K2QmOt0ZoAcN3oCPz7ufGqHgWBmcu0wVYGt1B/c0sOeWV1/PbyZH2L/bpC7HTt5mL786aIhXUFZTgcgTXFNK4q0KYCFlbU61fsVZQGLY2GB8Ybm22sTC9iTrKVfu7G19WcIX6m5uPP2WK0JgB4ubuycFw0q/cXk1fmZGOW8TH0H2yoC/N0TSPPrj3CjMQQpgx2XM2EwxBCmxBYnG6a90hndGo42ukfhRBihlO06a1Yk7UBLfWVRmsCwPRErRJWN3eVSTribjpykoq6JuOzqc5l4GQt/mMSdxVoszpcLIKXNzkx36SqWBuhO/Q6Q91Uz649SnVDs3Mm+zmK4TdpqdvbXzBaky7RlR3Hh0KIh4RGPyHEs8Bfna1Yr8I6VHssyTRWj1YiAvoxONRHP8ORt0NrI+7r5HnenbAsrRD/fm5MHWyOyYNncOsHAyeZJkAOEObvyTWjIvhgZx6nqhucI+TAp4DU0nAN4vipGt7adpybxkaTYPU1TI9OcfeC0XfCweXaTajJ6YrhGA9EAVvQivYKgcnOVKrXEdq6KTPBUKc2pieEsD27jNpGHepL8ndClLG7jbrGFtYcOMFlw8JxdzWhBzZuFpw6pNW7mIR7p8fR2GLjf86aS57xsXZTFeLEqXqd8NTnB3F3sfDzOYMN06HLjLsHhAW2mb8gsCv/YU1AHdAP8ASOSSl1rB7qBQREaxWgJsmsAs1d1dhiY3t2mXMFVRRAZYHhbqovM09Q29hiPjdVG/GtUwFNtOuIC/Fh/tAw3tyaQ5Wj637K8yBvu6G7jcOnW1iVUcx90+MI9TW42K8r+EfCsBtg9+taNpqJ6Yrh2IlmOMYCU9DmZix2qla9DSG0OIdJWo8AjI0JwtPNwjpnd8s9E98wNjD+6d5CrH4ejIs1xzjS8wgZAr4DTBXnALh/ejxV9c28vc3B7pH9n2iPKcYYDptN8v7BRsL8PPnBVB0m+zmKyT+DplrYYe655F0xHHe3VnM3SSmLpZQLgN7ZC9iZhCZrriqTpNN5urkwJT6YrzKdnJabvxNcPAzNmqmobWLD4RKuGD7AfKmWbbRNBcxab6oZDMMi/Zk6OJhXNx+jvqnFcQtnfAwDRkOQMRXay9IKya6w8atLE82VYdcZoUMg8XLY/iI0VBmtzQXp1HBIKXe187u3nKNOL8aaovXkqSw0WpMzzE6yUlBex8FiJ74B83bAgJHg6u48GZ3w+f4imlqked1UbSTMg4YKzYVjIu6fEcep6gYW73ZQ/KU0C4r2Gla7UdeoFfvF+Fm4xoghXj1l6oNQX665rEyKCaOIvZQzrUfM466aOSQUgK8ynTSh0CSFf8vSConp78XwSH9D9eiUQTPAxR0Of260Jt9h4qD+jIwK4MWNWY5pfpixRHtMuabna3WDlzdlU1RRz8Ih7ljMugPtiMgxEDMVtj6n/Y+ZEGU4HEVoa464iTKrQv08GRHpz5eZTopzFKdrhW0GFv6VVNazNauUq0YM0H+uuL14+ELMFDi82mhNvoMQgp/MjCevrI6PHbHr2L8Eoidq7VZ05kRlPc+vz2L+0DASg3qRi+pcpj6oNcdMe99oTdpFGQ5H0S8Q/CJNteMAmJVkJS2/nJNVTrhzOdMR1zjDsXxfETYJV5mlN1VnJMyDU4c1d46JmDkklJFRATy79igNzT2IdZw4oLXfMchN9fTqQ7TYJA/PH2KIfIcx6BIIHwmb/2GKkQ3nogyHI7EmmyolF2BWUihSwrqDTth15O/UjKWfce3Ll6UVkhTuR3yoiYu7zmbwXO3xyBfG6nEOQggenJNAQXkdH+7M6/5C+5dotQjJCxynXBfJKKjg42/y+d7kGAb2N2DOvCMRAqY/BKePQdp7RmtzHspwOJLQZO1usrnRaE3OkBzuxwB/T750RpzD4MK/3NJa9uaVmz8ofjZBsVpqrsniHABTBwczNiaQ/6w72r0MKym1+EbsNPAJdbyCHYqWPP7ZAQK93PnRJfG6ynYaifO1zLQNfzPVZwoow+FYrEPB1gSlR4zW5AxCCGYlWdl05JRj0y0ri6Aiz9DA+LK0AgCuHGGCgU32kHApHP/aNL3N2tB2HYmcqGzgne3dqOsoStMGmhlQu7HkmwJ2HC/j15cm4t/PBHNYHIEQcMlvoCIX9pgrkVUZDkdibW09YkJ3VV1TC1uzHFiNmm9sfENKyZI9BYyLDSIy0MsQHbpNwjztBiN7ndGanMfEuP5MiuvP8+uPUtNgp28942OwuELSlc5R7gJU1DXx11WZjIoO4MYxUbrKdjrxsyBqAmz8u6kGPSnD4Uj6DwaLm6larANMGNQfL3cXx7qr8nZoqaXhwx23ph3sy68g+2QN1/bGPP3IceAZYLrsqjZ+dWkip6obeXGDHQF8m00zHHGzwEvf6v3/++IQZTWN/GnB0N6ZftsRQsDM30BVIez+n9HanEEZDkfi6q41dDNZZpWnmwtTBwez9qADq8jzd2lZH64ejlnPTj7ZU4C7q4X5w3qZmwrAxRUGz9EMh82B7kMHMSo6kKtGDOClTdkUltd17aScr7WeZcNvdK5y55BRUMFb23K4fcJAhkaYvI6nu8RO0742/l0rMjYBynA4mrBhULzPaC3OY3aSlaKKevYXOsCv3twIhXsMq99oarHxWVohs5NCe68/O2GeNi604BujNWmXX89LxCa19NYukf4huPtA4mXOVewsWmyS3y7NIMjbnQfnGteBVxfmPA61pbDpGaM1AZThcDxhw6H6BFQ5qVq7m8xOsuJiEazKKOr5Ym2FfwYFxjcdOUlpTSPXjDLZ7Gh7iJ+lzao/tMJoTdolMtCLu6fE8smeAtLyyjs+uKke9n+qxTbc9Ys3vb7lOHvzyvndFcm99waiqwwYBSMWai3XTx83WhtlOBxOm8/fZLuOQG93JgwKYlV6cc/dVW2BcYN2HEu+KSDQy43pCSYb2GQP/QK1KvLMz0zTGPNcfjgjjv7e7vzxs/3YbB3oeGS11oNLRzfV8VM1PL36ILOTQntXOnZPmPU7LflgzWNGa6IMh8Np6xJblGasHu0wf2g42adqOHyih7PI83a0Fv7p/w9bWd/EmgMnuHLEAHMObLKHpCuh9Cic7KI7SGd8Pd14eP4Qvskt5/2OigL3fQg+VoidroteNpvkoY/34eZi4Ymrh5m/1Yyj8BsAk38KB5YaPpu8l//nmRBPfwiMMd2OA2BuihUh6Lm7Kn+n1ojNAD5PL6ah2cbVvTGb6lySrgSEtuswKdenRjJhUBBPrsqkpKqddNC601oV/NDrwaJPb6h3tuew/VgZv7s8mTD/XjCgyZFM+ol207b8QUOLApXhcAZhw6HIfIYj1NeTsQOD+DyjuPuLtBX+GeSm+mRPAbHB3oyKCjBEvkPxDdP+jpnLjNbkgggh+PM1w6hvsvGn5ZnnH3DgU2hp1M1NdfhEFU+syGRaQgg3jOnFMa7u4u4Nl/8dTmbCln8bpoYyHM4gfLjWY8YkqXNnM29oGAeLq8g+2U13lYGFf4XldWw7VsrVIyP6jnsi6Uptd2qCgOeFiAvx4YeXxPFZWiFf7D/npmPfhxCcAOEjnK5HfVMLD7y3B19PV/5+w/C+8x6wl8T5kHSV1orEoGaZhhgOIUSQEGKNEOJI62PgBY5rEULsbf0y723ZuYS1/hMVm6sQEDTDAbCqu7uO/J2GFf4t3p2PlPTO4TwXYsgV2mPmcmP16IQfzognOdyPh5ekf+uyKs/V6jeG36gVqjmZv6zM5GBxFX+/YUTvmCHuTOb/Tauh+uynWvGlzhi143gY+EpKORj4qvXn9qiTUo5s/bpKP/V6iEkzqwAGBPRjZFRA9+MceTsNKfyz2SQf7c5jUlx/ovv3shYjHREUqyVUmDjOAeDuauHfC0dS09DMrz7ap2Xm7ftQe3LYDU6X/+neAt7cmsPdU2KZkahvA0VT4hcOl/4Zjm+Crf/RXbxRhmMB8Ebr928AVxukh3PwDQPvUFPGOQDmDw0jo6CS3NJa+040sPBvW3YpeWV13DS2j/UiAs3tkLcdqnoQe9KB+FBffnN5EhsOn+TljVmw521tUl1gjFPlpudX8OvF+xgXG8RD83r5nA1HMup2bcf61eNQuFdX0UYZDquUsgig9fFCtxCeQohdQohtQojeZVzCh5tyxwFw+XCtTcdn++ycj25g4d8Hu/Lw83Tl0pQw3WU7naQrAWn6XQfA7RMGMi8ljK9WL9XieKNuc6q8kqp6Fr21i2AfD/576+jen4LtSISAq54F72D4+Ae6dlt2ddbCQogvgfb+y39jxzLRUspCIcQgYK0QIl1K2W40SAixCFgEYLVaWb9+vb0qA1BdXd3tc88mtjGAqJJ1bFq7BmkxX1VrQqCFdzYfJpk8ampqunTNEfmfMRjYmtdMw8nOj3cUNU2SFftqmR7pyravNzlkTUe9zg5BSsZ6RdP09WvsrR3sFBGOvN6rwyW27A1U2/qxOtef/qcds+651DRJntxRT2mtjd+M9yRj11a7zjfVa+xEAgb9iBFpv6f05euojvmJPtcspdT9CzgEhLd+Hw4c6sI5rwPXd2X91NRU2V3WrVvX7XO/Q8YSKR/zk7LgG8es52De2npcDnxouUzPL+/6NX/4PSmfSXKqXu3xxpZjZ3R1FA57nR3Fhr9p75fTuU5Z3qHXW1chW/5klR89dq2c+tRaWVxR57i1W6ltaJbXP/+1jH90hdxwqKRba5juNXYmW5+X8jE/eey1H3R7CWCX7OJnuFH7vmXAna3f3wl8eu4BQohAIYRH6/fBwGTAXIMuOiKsNUBu0jjH5cPCcbUIPt1b0PWT8nca46bamUfKAL++2/0UtAI60FqTm539n2BprmPYFT+itLqB217ZTlmN44rRKuubuPO1HezKOc0/bhrJtN7cWkYvxt8LI2/DemIjNPSwM0QXMMpwPAnMEUIcAea0/owQYowQ4pXWY5KAXUKINGAd8KSUsvcYjsBYcPc1bZwj0NudGYkhLEsrxNaVXkkGFf5lFFSwv7CybwbFzyYoVjPK6YuN1qRz9rwNwYkkpl7CK3eOJbeslutf2EJemZ3JFu1QUlnPzS9u45vc0/z75lFcMfwi6UPVU4SAK/6P3alPg4eP08UZYjiklKVSyllSysGtj2Wtv98lpfxB6/dbpJTDpJQjWh9fNULXbmOxaGmWJt1xAFw9KoITlQ0cLOtCHrhBhX/vbM/B083CghF9qHbjQgy7AU6kQ0k7Fdpm4eQh7b0w6jYQgolx/Xnr7vGcqmrgmv9uYcexsm4vvT27lMv+vZnjpTW8cucYrrxYmhc6ClcPmt38dBGlUhScSfhwbRqgCYf1gNZq3cfDla2FXRgRakDhX0VtE0v3FHL1yAj8vcyXYOBwUq4BYTH3ruObN7UOrcNvOvOrcbFBfHz/JLw9XLj5pa0888Uhu+bb1zQ088TyAyx8eRt+nq4s/dFkVathcpThcCZhw6GpVuuAakI83VyYNzSMncXNnc+Xztuhe+HfR7vzqGtq4faJA3WTaSg+oTBoBmQsNmer9cZazU2VdCX4Wr/z1GCrLysemMo1oyJ5du1RZj2zgfd25Hb4vqqoa+KVTdnM+Pt6Xtl8jJvHRfPpjyeTYPV19pUoeojT0nEVaMNXQCuaCzHnhLKbx0axeHc+y/cVctPY6PYPaqrXrmH8fbrpZbNJ3t6WQ+rAQFIG9OGg+LkMuxGW3ge5W2HgJKO1+S77l0B9OYy9p92nfTxceebGEVyfGslfVmbyyJJ0nlh+gIlx/RkeGUCIrwctNklxRT1p+eVsyy6lqUUyOb4/L9w2mtSB+s4qV3QfZTicSUgiuHlr40FH3Gy0Nu2SOjCQAd6C93bkXdhwFO7ROqBGT9RNr01HT3G8tJafz0nQTaYpSL4KVv5KcwmZzXDsfAVCkjrVa2Jcf5b9eDLf5J5m8e4CtmWX8mVmyZnnXSyC+BAfvjcphqtGRDAs8iK6MegjKMPhTCwuWtfQQnPOlQatbfb0KDfeO1jOweJKhoS1E1zLbS28ihqvm15vbT1OsI/7maaMFw3u3jDsekh7H+Y/pc13MQP5u7UbiMv+3qWGhkIIUgcGndlF1De1UFrTiJtF4O/lhoerPrM7FM5BxTicTcTo1lYdTUZrckEmDXDF3cXC+zsuMOUtdxsEJ4J3f130OXaqhq8OlnDz2OiL8wNm9O3QXGeumo6dr4C7z3eC4vbg6eZCREA/Qv08L87XtI+hDIezGTAKmuuhxLwlKL7ugkuHhvHJnoLzs2FsNsjbBtETdNPn5U3ZuFks3DHpIgmKn8uA0WAdqrmrzEDNKc2IDb8RPPVJ91SYG2U4nE3EaO2xwLzuKoCF46KoqGtiWdo5jQ9PHtQGUukU3zhV3cDi3flclxpx8c5cEELrfFq4R9utGs2Ol7XmluPvN1oThUlQhsPZBMZCv0BTxzkAJg7qz5AwX17bfKytN5hGW3xDpx3Hm1uO09Ri4wdTB+kiz7QMvxFcPWGnwXWvjbWw82VImA8hF1miguKCKMPhbITQ3FUFe4zWpEOEEHx/SiwHi6vYklX67RO528AnzOkzFwBqG5t5c1sOs5OsxIU4v22CqfEK0irJ096H2u5XY/eYtPegthQm/cQ4HRSmQxkOPRgwWotxNPa8l48zuWrEAIJ93Hl187Fvf5nbGt/QYTTou9tzKa9t4r7pF/luo40J92tB8m/e6PxYZ2Brga3PQUSq+VKDFYaiDIceRIwG2WIOf3UHeLq5cOv4gaw9WMLRkmqoyIeKXF3iG7WNzTy/Posp8cGqEKwNawrEToMdr0BLF9rCOJqDy6EsS9tt6HDjoOg9KMOhBwNaA+Qmj3MA3D5xIB6uFp5fn6XtNkCX+MabW3MorWm8+Ar+OmP8/VCZD5nL9JVrs8H6p6D/YG20rUJxFspw6IFfOPiGmz6zCiDYx4PbJgxk6d4CKg9v0nL3rUOdKrO6oZkXN2QxIzGE1IGBTpXV60i4FIIGwdf/1Ld/VeanULIfpj+kFbIqFGehDIdeDBjdK3YcAPdOH4SrRVB1eLM2I8LFuQ0GXtmUzenaJn4+W+02zsPiAlN/AUVpcOQLfWS27TaCE2HotfrIVPQqlOHQi4hRWpfcunKjNemUUF9P7h7Tn/D6LE4Hj3aqrMLyOl7YkMXlw8MZERXgVFm9luE3QcBAWP+kPruO/UvgZCbMULsNRfsow6EXbQOQCnYZq0cXWRRTjEVI/pfv3AFKT31+EJuEh+cNcaqcXo2Lm7brKPwGjn7pXFlNdfDlH7UhZMnXOFeWoteiDIdeRKRqQ3rydhitSZcIOLGdZuHOi9lBbMk65RQZu3PK+HRvIYumDiIqyMspMvoMIxZCQLT2oe7MwWBbn9My6S79qzbFUqFoB/XO0AuP1iBz3najNekaxzchoscRHODP458doMXmWBdJQ3MLD32czgB/T+6fEefQtfskru4w+w/aaNm97zhHRtUJ2PwPGHIFxE51jgxFn0AZDj2JGg/5u0w7SvYMdeVQtA+X2Kn85vIkDhZX8drZRYEO4NmvjnK0pJq/XDsMbw/V3b9LpFyrvYe++hM0VDl+/c8f0uauzHnc8Wsr+hTKcOhJ1DhorDZ1p1ygtT+VhJgpzB8axpxkK09/cUgrCnQAe3JP8/yGLK4bHalmS9uDEJoLqaYE1v3FsWtnLof9n2jpt/3VDlDRMcpw6ElUa4Dc7O6qY5u0BnsRYxBC8OdrhuLl7sIvPtxLY7OtR0ufrmnkx+/uIdzfk99fkewghS8iIlNh7A9g2/PfFmj2lNoyWPELsA6DyT91zJqKPo0yHHoSMBB8rOYPkB/fpNVvuGltzUN9PfnrNcNIy6/gD5/t7/ayTS02fvrBXk5WNaG5Cv8AAA8QSURBVPDfW0fj7+XmKI0vLmb/Efyj4NMfQWNNz9ay2eCT+6CuDBb8R8vgUig6QRkOPRFC23WYecdRd1rrqRXz3eDo/GHh3Dc9jne35/LGluN2Lyul5NEl6Ww8fJI/LkhheKSq2eg2Hj7ah3xpFix7oGe1HV//E46shrl/hgEjHaejok+jDIfeRI2H08e1DBYzkrMFkO1m1fzq0kRmJ1l5bNl+3tuR2+UlW2yS3y7N4KPd+TwwazALx0U7UOGLlEHTYdbvIGMxbHm2e2tkfAxfPa4F3cfd41j9FH0aZTj0Jmq89phvUnfV8c2t8Y3U855ysQieu3UUMxJDeGRJOk+vPthpmm55bSP3vrWbd7bnct/0OH4+e7CzNL/4mPIgJC+ANb+zf8zsoc81F1X0BLj6v6r7rcIulOHQm/AR4OLuuMCmozm2SXOnuXq0+7SHqwsv3p7KTWOieG5dFtc+v4Vdx8u+OzUQLZ7x0a48Lv3nRjYcLuEPVybz8PwhCPUB5TiEgGtfhvjZmsvq63917raSEna/Ae/fAqHJsPA9cOunj76KPoMhCfRCiBuAPwBJwDgpZbt9OIQQ84B/AS7AK1LKJ3VT0lm4emh3820jWc1E9UmtwGzmbzs8zMPVhSevG8bEuP48sSKT61/YSnyoD6nRgfh4ulJcWc+Wo6c4XdvEsAh/XrljLMMi/XW6iIsMVw+46W1t97Dm95CzFeY/BYEDzz+2PE/bnez/BAZdAje9BR6++uus6PUYVXmVAVwLvHihA4QQLsBzwBwgH9gphFgmpTR5EUQXiJkCm/4P6ivB089obb4le732GDer00OFEFw9KoLZyVaW7S1kVUYRXx08QU1DCyG+HkwdHMLVowZwSWKo2mU4G7d+cMPrWruQdX+GZ0dD4nyInQ7eIVBzkuT9S2BT6/3ZzN/BlJ+rBoaKbmOI4ZBSZgKdfaCMA45KKbNbj30fWAD0DcOx8Wktu2rwHKO1+ZastdAvSHOndREfD1duGR/NLeNVwNtQhIBJP9ZiHttfgH0fQuZnZ54OcAuA1Lu0aX4BUQYqqugLmLnXQwSQd9bP+cD4Cx0shFgELAKwWq2sX7++W0Krq6u7fW5XsbQ0MEW4kr/xXbILjM+br66uZv26dUzMXEWFfzIHNm4yWiWno8frbBgec2DMbDwaTuHaXEOzqw+lTR74ePnC3iwgy2gNdaFPv8YXQK9rdprhEEJ8CYS189RvpJSfdmWJdn53wciflPIl4CWAMWPGyBkzZnRFzfNYv3493T3XLnLGEt2SQ7Qesjph/fr1zEi2wobThE68mdDRxuvkbHR7nU3CxXa9oK7ZmTjNcEgpZ/dwiXzg7D11JFDYwzXNQ8xU2PSMeeIcWWu1x0GXGKuHQqEwPWZOx90JDBZCxAoh3IGbgWUG6+Q4YqaAbDFPFXnWWghOUP5vhULRKYYYDiHENUKIfGAisEIIsbr19wOEECsBpJTNwI+B1UAm8KGUsvuNksxG5FitnuO48fEES0sjHP9a7TYUCkWXMCqr6hPgk3Z+XwhcdtbPK4GVOqqmH+5eEDFGq9Q2GP+KA9BcB3EzjVZFoVD0Aszsqur7xE6Fwj1aY0EDCSrbBS4emvtMoVAoOkEZDiOJmwXS9m3hnRFISfCpnRA7Teu6qlAoFJ2gDIeRRKSCpz8c/dI4HU4doV99MSRcapwOCoWiV6EMh5G4uMKgGXB0bc9mKvSEw6u0x4R5xshXKBS9DmU4jCZ+NlQVQkmmMfIPr6baO1al4SoUii6jDIfRtDUUzPpKf9m1ZZC7jVPBY/WXrVAoei3KcBiNfwSEJBkT5zj6JcgWSvuP0V+2QqHotSjDYQbiZ2kjWxuq9ZWb+Rn4WKnyVVP5FApF11GGwwwkzIOWRn13HQ3VcOQLSLoKhHobKBSKrqM+McxA9ERtDsbB5frJPLIamush5Wr9ZCoUij6BMhxmwMUVEi+Dw19Ac6M+MvcvBe9QzWgpFAqFHSjDYRaSroCGCji+0fmy2txUyQvU+FCFQmE3ynCYhUGXgJs3ZOrgrlJuKoVC0QOU4TALbp4weDYcXAG2FufK2vch+IQpN5VCoegWynCYieSroaYEjjnRXVVdAkfWwIiblZtKoVB0C2U4zETifPDw03YEzmLfh9rkwZG3OE+GQqHo0yjDYSbc+mkB68xl0Fjr+PWlhL3vaF15QxIdv75CobgoUIbDbIy4GRqr4ZATBh8W7YWSA2q3oVAoeoQyHGYjehL4R8Hedx2/9s5XwM0Lhl7n+LUVCsVFgzIcZsNigZG3QtZaKM1y3Lq1ZZC+GIbfBP0CHbeuQqG46FCGw4ykfk/LeNr1muPW/OYNrXZj3D2OW/P/27v/WKvrOo7jzxcQLiIEBK78MpTBMkhRFJnTEZSCzCAL+sc2Vy6KGFOWqY1VZluRaCZO10zadWnSD0EryyACRZs/gLj82BVhgUCAkBrCTAfcd398v7Q77r1yv+fe7/keznk9trPzvd/zPee+3zvwfe37434+ZlaTHByVqNfAZPDBf/yycy6SHz8KryyGYVdC3aiOf56Z1TQHR6UaNwveOwQNj3f8sxqWwKHdcPncjn+WmdU8B0elOmc8DL4Env9pxwY+PH4M1twDA8fAiKs7rz4zq1kOjkolwYTb4NCujh11bPw1vL0DJtyafKaZWQc5OCrZiKtg0EXw3N1w9L3s73//CPztBzDoYhh5TefXZ2Y1qZDgkDRT0hZJTZLanPBa0k5JmyRtkLS2nDVWBAk+c0dy1PH3Rdnf//y9cHgfXPPj5DZfM7NOUNTeZDPweaA9o/lNjIgxEdFmwFS18z6VDH645h54a0f737dvI7xwX/J3G0PH5VWdmdWgQoIjIhojYmsRv/u0NPmH0LU7LP1qcmvtqRz9Lyz7GvQ4C6YsyL8+M6splX7+IoDlktZJmlV0MYU5czB89j7Y8wr8ZX4yWGFbmpqS0DjQCNMfgB59y1enmdUExQfthDrywdJfgbNbeWl+RDyVbrMauCUiWr1+IWlQROyVNABYAcyNiFZPb6XBMgugrq5u7JIlS0qq+8iRI/Ts2bOk9+Zt+PaHGbrnD7x+zgx2nPulFndJqek4I197gIH7V7J9+JfZM7R9M/xVcs95qbWea61fcM9ZTZw4cV27LwlERGEPYDVwSTu3vYMkZE657dixY6NUq1atKvm9uTt+POLJORHf6xXx6MyINxqT9U1NEbtejvj5p5PXVv0oWddOFd1zTmqt51rrN8I9ZwWsjXbuu7uVFE1lIOkjQJeIOJwuXw3cWXBZxerSBabdnwwbsvJOePCyZArYpqPw7pvw4b7whcXwyRlFV2pmVayQ4JB0HXA/0B94WtKGiJgsaRDwcERMBeqAZUpOx3QDfhURzxRRb0WRYPxsGD0DNv8O9m8CdYEhlybDpZ9RW4fmZlZ+hQRHRCwDlrWyfi8wNV3+J3BhmUs7ffTsnwSImVmZVfpdVWZmVmEcHGZmlomDw8zMMnFwmJlZJg4OMzPLxMFhZmaZODjMzCwTB4eZmWWS2yCHRZJ0EHi9xLf3A/7dieWcDtxz9au1fsE9Z/WxiOjfng2rMjg6QtLaqLFJo9xz9au1fsE958mnqszMLBMHh5mZZeLgaOmhogsogHuufrXWL7jn3Pgah5mZZeIjDjMzy6Smg0PSLyQdkLS52bq+klZI2pY+9ymyxs7URr8LJb0qaaOkZZJ6F1ljZ2ut52av3SIpJPUrora8tNWzpLmStkraIumuourLQxv/tsdIelHSBklrJY0rssbOJmmopFWSGtPv9KZ0fe77sJoODqAemHLSutuBlRExAliZ/lwt6mnZ7wpgdERcALwGfLvcReWsnpY9I2kocBWwq9wFlUE9J/UsaSIwHbggIkYBdxdQV57qafk93wV8PyLGAN9Nf64mx4BvRsT5wHhgjqRPUIZ9WE0HR0Q8B7x10urpwCPp8iPA58paVI5a6zcilkfEsfTHF4EhZS8sR218xwD3ArcCVXeRr42eZwMLIuL9dJsDZS8sR230HECvdPlMYG9Zi8pZROyLiPXp8mGgERhMGfZhNR0cbaiLiH2QfDHAgILrKaevAH8uuoi8SZoG/CsiGoqupYxGAldKeknSs5IuLbqgMrgZWChpN8kRVrUdTf+fpGHARcBLlGEf5uAwACTNJzn0fazoWvIkqQcwn+TURS3pBvQhOaXxLeA3klRsSbmbDcyLiKHAPGBxwfXkQlJP4Ang5oh4pxy/08HR0huSBgKkz1V1SN8aSTcA1wLXR/Xfnz0cOBdokLST5NTceklnF1pV/vYASyPxMtBEMq5RNbsBWJou/xaoqovjAJI+RBIaj0XEiV5z34c5OFr6Pck/ONLnpwqsJXeSpgC3AdMi4t2i68lbRGyKiAERMSwihpHsUC+OiP0Fl5a3J4FJAJJGAt2p/gEA9wIT0uVJwLYCa+l06RHjYqAxIn7S7KX892ERUbMP4HFgH3CUZAdyI3AWyZ0I29LnvkXXmXO/24HdwIb08bOi68y755Ne3wn0K7rOMnzP3YFHgc3AemBS0XWWoecrgHVAA8m5/7FF19nJPV9BcgPAxmb/f6eWYx/mvxw3M7NMfKrKzMwycXCYmVkmDg4zM8vEwWFmZpk4OMzMLBMHh1kJJPWW9I2i6zArgoPDrDS9AQeH1SQHh1lpFgDD07keFra1kaSukuolbZa0SdK8dP1wSc9IWidpjaSPp+vr0nlRGtLH5WXqx6zduhVdgNlp6naSeUzGnGK7McDgiBgNySmudP1DwNcjYpuky4AHSYbFWAQ8GxHXSeoK9MynfLPS+S/HzUqQDmP9xxOB8AHb9QHWAn8CngaWAz2Ag8DWZpueERHnSzoIDIl03gyzSuQjDrMcRcTbki4EJgNzgC+SzBPxn3YcrZhVJF/jMCvNYeCjzVdIevXkjdL5zLtExBPAd0hG4n0H2CFpZrqN0nCBZFC62en6rpJ6nfyZZkVzcJiVICLeBF5IL3ovTAOitYmRBgOrJW0gmRf7xCx01wM3SmoAtpBM9wlwEzBR0iaSkV1H5diGWUl8jcOsE0i6FjgvIhYVXYtZ3hwcZmaWiU9VmZlZJg4OMzPLxMFhZmaZODjMzCwTB4eZmWXi4DAzs0wcHGZmlsn/AFXzdj1Qh0vTAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f2b185fb710>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "pl.plot(res['t'], res['x'])\n",
    "pl.grid()\n",
    "pl.xlabel('t, sec')\n",
    "pl.ylabel('x')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAAyCAMAAABCrGH/AAAAM1BMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADxgEwMAAAAEHRSTlMA74lUMhBmzavddiJEmbsg4DnoagAAAAlwSFlzAAAOxAAADsQBlSsOGwAABqlJREFUeAHtXNuWqyoQRBQliG7+/2tPdTe3ZGLGRCdL19GHGTQ0FGXT3CpRTcCl1XX9IQOaOG5U03ZdZ/6wnqtoZUBxC6q7s3Ix9hl5P+bkURPdiake54pVe3iuz0z1IOyaG/03Q8X7IZMnprqfhNE28H/tDklwAXViqtmZ0ZJJ/Nk0pVWHTJ2X6vEWCb3FiergD8lwBnV0qvshNJhn2HBrM2ZJOI4f3s5hthw7dAwoD/kOc3t0qpVqLJH1c9Cb40wvhmql+p95DkMzATk+1S5geVXNoBN9yc1jqFaqTREl5TjY/+NTbQLc98m+QRNDcwrVqpOZyMEILnCOT7Wy8NZlqk1Io+FFdXmtH6Z88Lzs9kNolQ5piRgDiCdf5h0cfwWQDxkuZjcb5xaDU1PesInDosZg2PLDa1gslH2a6uHNfJlmKlO+OLebrDKyTpyeRJlPq/wLuxPEamVyZOhpNhKvGC/GuY+zviE7fMpyrP9noDozZiZXzZ3Twjx+nFePOfvBEqeiejZqLkvCtN0UGb22m/Zzrb7BZGMKc44Tc71FfW2i7sf0z5LujgbueP+Z9wBPThVAHvmqluvXgdcjOf/n+1N79ble3EX1197XRfVF9dcY+FpFh/LqA80oCpT35zbFVtXGv1HdzU1Zn/2xA9zNk7dLaGzabl2EbeYmb6/cZ6qhvIuktlWV8W9UK5W21e6x/MXdzhIaG7ehXkEdFvxIoHwm5llqxq9U+2oz7RXo7Z+lPY1vSmgW/ChC+QjJYjN+pVovdLHt1D6WkHbqviihWfKjCOUjJIvN+JXqpS72SNTm+7wJms5lvyChWfCjBOUTJMlW/TB+RfU4addJF0Nycg325b224+gcRC5eu5lOUL2exw5PaNPeaOcm7Lchl/G63tbEB4625HaV0KDQHoUWBKVe4yY6mylwcYOTscm5thotyY9aZwmYQRNG7xgzq3nWi3liXSPxo5aVQC+o9g1ImzlUc3KgA1Ot5hltMMGBZoiz6YmlE1WNFnjqPCPOnZBraNWQVXQj5RhkW38/CQ25Pe1fFwSlXq081T4luLzfShjHUEZCAz9yoyf1gwHq/ubUjUDGY8u1Yp5JiHAsjlhWAi1TbUjAJdrDmATQrhUwI3+jo0fhIygl/RHSYwNCDQ5aKRee+byz3NCxn5Vzwd0kNCM1DV93KAhKvXhGgbbAJSgtu001EuKBM8pQOzU8UuPzlt5JPIxfKebBNJpzztxdom0SbaLWNNgVqs1wyxcd000MjCNOTJI3w58paghqmrkCIoMHTtv06NF4kHIhI19iH28+kNA8QpOSTGgmwlIQ/Ivo5FkFl3wCkjTiogvs4HSPJmYPRyFqYK6QiGqeFG1rMc8TKDBtKFjxn2SbQ3UxLlQj893FwITYmGQnlV418auSwmUYnxsVWF3HheS+x3dsnwvfTULj5yCqkIwALpCqEWeIDyyFhY4ChXLJyZC+WUfRMF0hpYXq9WIerr8LzE98TU+MF6k2HCIYepWMEUXemYeD4JXyMG7gINXXxHLfo2aIfWqQ2klCgxiBsYx9NCMo3ZaiGnmqAAHJBiGO/Lk6nCQ6xiyOwqtgrpBHgsB6MQ87noTqFHyeGL+gmkLrPCB60egRk+gmFHblnaEZ1F6eDvYYfvkj+A0yxCRSdLE7IerIHZwp9ttNEhoZdNmLM4JSL2EmpAyEPFxHB69CNYWVlK81irlC7E7D4noxz0zdGaGagueyEmiRagUeMCZbAlwlxQsk9g2aW0OduKOvibFywNBkL/uKCdR1IYtBo7IwaR8JDR8mOu7zGUGuF12H1avyAKRhUmjIqy27tsDiTRK4Mlro8QZmMC87RaLmWS/mmcCAx14RvdtlJdAy1cY63Y2WJnVVUoQujqNyO/GsNniZTNOc0on+pcjnaCaC+nXfp0CI4nK43CKhwey+d9xREJEiglKvsoxEHoy2p4weLYqLFobFDgjM8EW00GF6HNU70Wq1mMdY7drO8sIhQniiBFqmmihacy2suNaY7iSheQdBFapfIEyLa8mSF4AvLMpH97aqGG+nesPKfScJzUoENOiViUHh5kkqbRnJR/Wq90nmh0f3tqoYb6ZaxsyH6lbd7iWhWYmAFiiKB7AV6GpVybtintq2/jrlVqrbW6BBc5frbk/9DvCr4tciaG3vKiHrqyIxPKbVDLKtRhKLrG1r461Uv0b85qcLJ0VvlrJL9gKlPrNaV3SxfevAa13hV64VDJBXXz9SsYKobVnSj1RcP72yjccV1vLTK/8BCv42bBTkmxkAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$\\left[\\begin{matrix}\\operatorname{v_{x}}{\\left (t \\right )}\\\\- damper.c \\operatorname{v_{x}}{\\left (t \\right )} - spring.k x{\\left (t \\right )} - u{\\left (t \\right )}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡                vₓ(t)                 ⎤\n",
       "⎢                                      ⎥\n",
       "⎣-damper.c⋅vₓ(t) - spring.k⋅x(t) - u(t)⎦"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sympy.init_printing()\n",
    "model.f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAOIAAAAyCAMAAAB28PfKAAAAM1BMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADxgEwMAAAAEHRSTlMA74lUMhCZZs3dIrurRHYgTxvb/QAAAAlwSFlzAAAOxAAADsQBlSsOGwAABCxJREFUaAXtWtuWoyAQROSiCCz//7Vb1XjbORtjjGPOZMJDQlqgq7qbBiGqKShavWXR5NaopjXG2LdkqCyotaBo3pPexMocoGid9vrn2OUIxa5VyvZhstJl321zaDodoOgbktLdZdREkY0ulqsoDpE622P6nrFLPqbygBeLUEwF4XptuYyiLY7MUsnXElTqMoqhbhPM9buFCymKF9+Z4i8IVFXTTfvG6UYNA/OMP5bBn0lRl81FVZd+d/HSD9tcR1H13MC9YPN+IUWrc46Xb8NdbEoX/eORfmB387iS1/b4UHyt/c/R/vHiOXZ87SgfL77W/udo/3jxHDu+dpSPF19r/3O0P+HFKC9VWzDs0PRbz//7zAyNHCv89+ER4TMUd2yJuyNoT37ZfoLiHoseQZtOftn+XoqH0OrHg3vT2t9L8RDaQ8G9QXKTIl59PV59kx6CwVuwRS3apDUOiW12PCqGIAQ8Go+NW5dzu8pCRNs+8PocnM6mBjeqLvOeZlGRdB6SKJ3wAAcwOtwgLcgmsniQ5b19i2KH8QaA1CryLEMPqA24lupwb6NV4u2Nw2EVCRe5p0ILFephuSiyQJtDKjvykrRPDYYZZCpKtSsQQ2dVkQEH16FrPCrxgizgintBJiOpQMSdnC9tUAxUgOvjAFK8xvDFtKpHLRnKeD5FAYcJhfaqNzmrDANBtsruZWgbtpRzr7EKky0qeFfvgWnBo0LDYyRnpFVFhkYsDVtHPMUI8xWq7fq5dERsS+MYFwgWga3LH8goqTJ6bxRQM4ZliJrqUP5WbuVREfz78VWjE//1BDdWVyqq/bgWL3hsbDzmEgQTkElB7V9/LRSnp8t3GkohalzRYBQ1IEha4ULZePlWBZGuNBKQeZUP+5gZY3uLmKgacazKXKoqnIzbyHAznnpqLeMvyPhT+o96b1NEPCCpiE8kMcpB/3J6itgl7yoAOYsg4oTk5B0L0YbR7ZNs49vKRZBAXVUnFeLcBBXQOuNZXR0tyKCi9h913abIpIIQ5Yekcc/0JiFOETMJI0oE9KgePbyaigzfqR373ClsChN1qg2r6qiiOhc8aPkZzwiHrp2RiZYxw32Zi18BDPRJlsBguBqeDS8+SUUFJIcqwKU4FhdJq1FcaQsjVzaxSEQ0xZ7SQZlHPkP7VbWqqLO90zLYjCdTjeWiMSOrqh3zI/MQym0vYrnzWcyQSqqLD3I0+0iJnumvCkL0bJhi1uNiz8SrBrbA+smA3lMs+psQuTisqlVF5oCqdbJOz3g4uuBYIRPVSmvvxT1bFBdQj+xRVlNxGeDk2iN4oPq2FxdcO3dUdaGpi8rS+RtqO/FMmndQrHN/6nDzW3Nlqf/muNnmlAc78cy67lNs+8LJf7e00Wcnc/du06ca7MUzK7lPcW76Uyu/hOIv+LPm2//l9i/zEyHhFdd0mQAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$$\\left[\\begin{matrix}0 & 1\\\\- spring.k & - damper.c\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡    0          1    ⎤\n",
       "⎢                    ⎥\n",
       "⎣-spring.k  -damper.c⎦"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.f.jacobian(model.x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Lyapunov Stability"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using the energy function as a candidate function, we can show that the system has global exponential stability (GES)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQ4AAAArCAMAAAC+eZamAAAAM1BMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADxgEwMAAAAEHRSTlMAdt3NVDLvmRCrIolEZrsgN5RvdgAAAAlwSFlzAAAOxAAADsQBlSsOGwAABVNJREFUaAXdWuuagygMBUW5yvL+T7tJuChqC7ZO/Xb5YREiJ4TkEJhh7OaifWtA5eehJfO/6ddLeyqjbMt8ICFm9cFXz38SxP06KC11+E+aw/2R2vyPxr1/7aoR7R9RxzPmsJwvE2POejNxrhXUtHLWcsYUl/DEBmOgC+tQhORcbBgWqUNwDWPcWx4xx+gY8zAhy/QAHGA91Pwo2Dhjm4Mnk8x7NE4wOF+QYCas7KmAOrhxoYNzr5nrCXOYADpOihkwgIbqEibBBqi5CdvkCL3QAD9gA3QAQRG9YU9o4Iqp263BnjCHCrME/2BA4jRFG/4BN8CW2IZekRoWtBybMUym6Cj4zuTGU6jhrscT5mDOh0BEEDcID9EhaN44q+gKqUGji0wUFHzDnoPmGEonxUovv9grO8zxJcJBZ4gHIExaa9ogFC41RkgsED84n9gAhlAQTEggSDapoCuZ5E65Lf5akCQT1s3db21zfItwUAUJE8IEH5RbLpgJzhbfoSBLYp0a0FNscpQNdWAIZTn6qDwGMMc2qkpHZ6Vtjm8RDop4XGtOvo4hM83Alutau8AMUGRssCODDZm2F00uogJ6ER1qgGSzCVeIARhoHWtt7611mONLhIMmkE4snHJsFxznEq3jVmLQC24YscHoBQWd5jYlXrgBMY8SkJ+ck0QiH7DLuvMsCNIsUs9h7DgxFQQc8YDChxAga4AyDbBpYKWvXMktN9TRGHxIJGuIrZOw7rJHY+TcnRHw/QQlujB2Cly87tJ5LMXUI2/E7bFtXo8xWkBhtDNVqLo9QkuiIKDgGcpM7Aidl6yBbNhRLCZhvnfkJVsjV9IWjgeAuohDS93/6i0PTP35pUIZU9bArxwixBDGHo2EXrjsMhwo6IBeIbmFQk4Bv2kLV3nFsI/Kh+ZYEXCUU5RI+8wcqT4i/+w5jUIIjxRrMjkPSSk8KFXlM3OsCDjYOYqNOVGvQ1daff+i4ORrHMdwmAMWHJJT0Djtg4/H4iriUaJljmUMM+xNOuAxs5QVAZvOUZaAXr8LFTtuSo7oMux9FQXLvwDZ05mvDEubMbyV7H/Zk2nLHGBcWt/9dwUCKucoDvPsp0LFQmwg3YpqL81LWrJ/kaMnT6dtDkrKNnlF/nL9PUdRmF4+FCp0whm3OQZpOyeqyNQBuXuZRfTbYSbvfaO2wrPjWz58gRIGtuypqqDXFYru2x5pbArWCicpuuYrG3NEwb13bHTKI2nwqB5z7FGGYA7sYIdNebMGGfuL34luiaoBkhvDEQhyMOzZnAOi4N4csbV6wimC0ns34sEyXkpUAi9QgL8Js5L90YtQcIkFWPURJpEcHAKZoDzkAyqFtEKnVYYMSZ4kVS9Q5Ontw0/s4WDhPJhjR3lpX5VwexJzPbn3+g7vgKuWtMmq+TQTfIGy/G0wvLOr0pzD9TudgzdyKTaMzx3jfnF7zKHKdrTQfe0GgKptlP0XT73n9Dnhl/yx6NNjjiKsJD/NP5ooZYQvKkZK7/ereXG8fLhKnx2PcKrKU1DsDSxQ4+kVQxPlotZn4gpDLyX+Z/19bXTrlkV7DvivYZcZyEmGsxW6jJI16v+FP8tBKbcG/R9WktXFTKV1Jba+fAR7GWXF663Ff97Qx2Sid4Aot9lrui4HP4O9inJtDig9Eot/bY6rwA/Bdqo5nG5tnR9/LvYQbEthF/Z5U+uLW/ofgm3qTn9CaErdLvAQbGse8nAkbH1xS/9DsC3d461bS+r2/odgW/MQ6BvmkDa2Pvu2/yHYltoTRQr/tTkegm1ZwwzWWolXUT8tD8E25zjGe7pfm+Me2H8BIHsqiQ8jahkAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$V = \\frac{spring.k}{2} x^{2}{\\left (t \\right )} + \\frac{1}{2} \\operatorname{v_{x}}^{2}{\\left (t \\right )}$$"
      ],
      "text/plain": [
       "              2        2   \n",
       "    spring.k⋅x (t)   vₓ (t)\n",
       "V = ────────────── + ──────\n",
       "          2            2   "
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "v_x, spring_k, t, V, x = sympy.symbols('v_x, spring.k, t, V, x')\n",
    "V_eq = (v_x(t)**2 + spring_k*x(t)**2)/2\n",
    "sympy.Eq(V, V_eq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXEAAAArBAMAAAB8y6WyAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMrvvmVREEIndzSJ2Zqtw62vAAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAFNElEQVRYCdVYXWgcVRT+9ncmO9lN8E0UM7Y+tdRGDbb4O30LabAj+lJ8cKQvKYg7oA+2FHcFH0wJdF8ENRC2WEWj2KX6YMCfiVRRadMVRWhjyaIPKg24tKSJVazn3Hs32bmZTGifnAP3zj3nnvPtd+f+nR2AJbXVEc8EVn1uAkkLymU/qczHk0ocJ5LJ3BqeXUom8wON9HIimRsPI1tPJPPsMnqqiWSesVHxEsm8UsXTiSSOiottuUYSuReC4pSZROIwFl/+/Oswc6vJ+nzYqGmfafoNqJvAC+QirQIuukTGfnv9beDWKw5wSviXIgJXgSybm7cMrRpuoBEPL5FxiQC5aBIZa1wjr1yLqvcBYxJGnZobSUEO66Yu33h4Rj7UwGH6ZS6aRMfeR14lKkWHRtAGJrSobvVHoeRv5vLV4XNBNzAYmQ47k6xcwqLHqt6T9Jyl0kOj5mv1eWWPekwKo1mL6tvEpsNrzBl5hCZ9UJYwmB6rev9swKpSe4YKX6sFjxrRYknKBSe6O9aqw4eZC+TXCeAdVUJYeqzqPOaJxYLTwHM77w5gOqGobqWXuz6eGHBhDY/uwtilT+eCsS89jG2/fc7Hi2f2Yuy3MbnaJs7048KVZukvVwHo8B3mR0Zh7QIh53+5vAi8Qu5cYmMV5EIgFgt2kM7XarHOHUenWU5wc03STSD/IHb6ONAotFHdEhj/eLlBVP9wTTs/iU9Q3ea+xf4zbqEOPAJQLSUMT1sqkHbnCfS2wci8yXBBldhYGYmyIxYLXiOdVhpyNdWx/tHj0jZokhdnmnbKfRPWCopty/2Q6o9cY0/KncIzHLcbJfI958uDmC06vGJuNUeQtsHIvMmwX5XYWPJhqdgl8byfal5pqbZQ16rDb7AcJ0PGA/b55MWZZmAYV5FbRsk2sARz+defz/psYqF+lkpz7ZoOwZvT0+89MD3dou2Ih5AJkPHkJsMAIMpGsdc7Qiz7aj+RN0DQ4nfXMRe9osp4wL3sxZlmP09vtsa/S8dkevBfdhEzTkO0WUGpFmIegu+sFood6EfGo3FSWWMeF0t+LD1/O+JJ05lvw0e+xuq4eNFTome1ojk1VpidzDRpeukwKjd8s4aKs8JuYsaJRFXE5K/OI3/SPWaTFoZfW+c09nNEglZLmSYAB8mVC+JiuZ8kfdkXT9pCZj3Xr3aoMOkV7SNjCen6UzLTzFT5PY1YbsHBD/5uoBdkYulrAsQFUw6dssEeNunwnXeebmFL0ecdehSPd+3QuFjGI8nKpclnUbZ1hKoqGaOFD8wPcL5VlZnmFx4WgHdN9Nmpe/AdjL0gEyifyLbwQj9hLBB9Y4jHsA6+wzxr54dMcRQ/ZgTAMLlyiY0V/XQKtmSD5sjaQUMv8E9GCy+kl765bdGTmSZFzAFPfo/ys6cb6B29qCaa1tirFwOGOMTVgs+1Dt9hbvz+1R2zYomeongcVyU2lvG6pKchlZkum97kOzpCxiNsqyZjtM5tHb7DXPkpZL5L1U3NHZGxqa2OClKPotI3z7jCcaRpF1a4P+3fxQYd3lAvSjnLXA5mEMq4omP7XE4ru2S7aMdmuelQQCeWdu3GUnwU++7k7nh4hRzOcjeILfvaf2iZvkeTU8wsWzVCj/3XgpAercTDK2QenhxiCCQcS4tzJNQtb+r5kE1X/hf/5mhx8mWfNOHvoTKtTBpz8T2UrvHEifwemq0njrjIUquJ/CTKWaon08qEvfaK+B4q0sqkMRffQ0VamTDmMksVaWXCmMssVaSVyWL+H3K8qUT9INXvAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$$\\frac{d}{d t} V{\\left (t \\right )} = - \\left(damper.c \\operatorname{v_{x}}{\\left (t \\right )} + u{\\left (t \\right )}\\right) \\operatorname{v_{x}}{\\left (t \\right )}$$"
      ],
      "text/plain": [
       "d                                        \n",
       "──(V(t)) = -(damper.c⋅vₓ(t) + u(t))⋅vₓ(t)\n",
       "dt                                       "
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dV_dt_eq = (sympy.Matrix([V_eq]).jacobian(model.x).dot(model.f)).simplify()\n",
    "sympy.Eq(V(t).diff(t), dV_dt_eq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEKCAYAAAD9xUlFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAIABJREFUeJzt3Xd8lfXd//HX55yTQRYBQgIkQBAIQ3YiQ0SC4i24uLVq3dpKqavVau1th/a+2/vnfbeO6t26cOEqlDqpexEXGxVkg+y9ZCSQeb6/PxJTigwNuXKd8X4+HnnkjOvkvL8GzzvX+l7mnENERAQg4HcAERGJHCoFERGpp1IQEZF6KgUREamnUhARkXoqBRERqadSEBGReioFERGpp1IQEZF6Ib8DfFdZWVkuPz+/Qa8tKysjNTW1cQNFOI05PmjM8eFYxjx37tztzrnWR1su6kohPz+fOXPmNOi1JSUlFBcXN26gCKcxxweNOT4cy5jNbM23WU6bj0REpJ5KQURE6qkURESknkpBRETqqRRERKSeSkFEROqpFEREpF7clMKO0gqeW1xBeVWN31FERCJW3JTCtC938M6aaq58YhZ7yqv8jiMiEpHiphTO7tuOH/dJYu6arxjzl09YsGG335FERCJO3JQCwJB2IZ4dO4j9lTWMeeATfvPyF2zYtd/vWCIiESPq5j46VoOPa8WbNw3j3neW8deZa3lu5loGd2rF0C6tOD63OXmZzchMSSQxFCAxGCAxFCAYML9ji4g0ibgrBYDMlER+N6YX404+jufnruf1LzZx99vLDrt8MGD1BdEqNZH8rFSOb5fBsK6tGdAhk1Awrla4RCSGxWUpfC2vRQo3jSzgppEFfFVWyZfbStmwaz979ldRWeOorA5TVROmorqGyuowldVhtu6tYNX2Mj5Yto0/v7+Cts2TuXRQB64a2om0pLj+zykiMUCfYnVapCZSlNqSom+5/O79VXy8fDuTZq/l7reX8dT0NfzmzB6M6ZfraU4RES9pu0cDNW+WwJl92vLM1YN46boTyc1sxo2TPufmyZ+zv1LnQohIdFIpNIL+HVrwwrUnctPIrrz02QYue3wmu/ZV+h1LROQ7Uyk0kmDAuGlkAQ9eMoAv1u/movEz2L1PJ8mJSHRRKTSy0b3b8vhVRazcVsbYp2drWg0RiSoqBQ8M69qaP32/H3PWfMWtz8/HOed3JBGRb8WzUjCzJ8xsq5ktOMpyJ5hZjZmd71UWP5zZpy0//7du/GPeRp6e/q2uly0i4jsv1xQmAKOOtICZBYE/AG95mMM31w7vzMge2fz+1UV8sV5zLYlI5POsFJxzHwI7j7LYT4AXgK1e5fBTIGDcc0E/WqUl8vO/z6OiWvsXRCSy+bZPwcxygXOBh/3K0BSapyTwP+f1ZumWvfz5vRV+xxEROSI/z2i+D/gP51yN2ZEnnDOzccA4gJycHEpKShr0hqWlpQ1+7bEIAEPbhXiwZAVtKjeQl950XezXmP2kMccHjdkjzjnPvoB8YMFhnlsFrK77KqV2E9K/H+1nFhYWuoaaOnVqg197rHaUVrg+//mWu+TR6S4cDjfZ+/o5Zr9ozPFBY/5ugDnuW3xu+7b5yDnXyTmX75zLB54HrnPOvexXHq+1TE3k5tMK+GTFDt5auNnvOCIih+TlIakTgelANzNbb2ZXm9k1ZnaNV+8Z6S4d1IFuOen892uLtdNZRCKSZ/sUnHMXf4dlr/IqRyQJBQP85qweXP74LCbOXMtVQzv5HUlE5F/ojOYmdlKXLAYf15K/TP2SfZXVfscREfkXKoUmZmbceno3tpdWMGHaar/jiIj8C5WCDwo7tuSU7tk88sFKdu/XTKoiEjlUCj65+bQCdu+v4tkZmhdJRCKHSsEnvXKbM7ygNU9+skrTa4tIxFAp+Oi64s5sL61k8px1fkcREQFUCr4a2KklAzpk8sgHK6mqCfsdR0REpeAnM+O64i5s2LWfV+dv9DuOiIhKwW+ndM+mW046j3ywUldoExHfqRR8FggYPzwpnyWb9zJj5dEuPyEi4i2VQgQY0y+XFikJPPHJKr+jiEicUylEgOSEIJcO6si7i7ewdsc+v+OISBxTKUSIy4d0JGjGU9NX+x1FROKYSiFC5GQkc2aftkyevY7SCk2UJyL+UClEkB8M7cTeimqe18lsIuITlUIE6dc+k/4dMnlq+hrCYR2eKiJNT6UQYa4cks+q7WVM+3KH31FEJA6pFCLM6N5taJmayDMzVvsdRUTikEohwiSFglxY1J53Fm1h0+79fscRkTijUohAlw7qgAMmztIOZxFpWiqFCNS+ZQojumUzcdZazZ4qIk1KpRChLhvcgW17K3h74Ra/o4hIHPGsFMzsCTPbamYLDvP8pWY2v+5rmpn19SpLNBpekE1ei2ba4SwiTcrLNYUJwKgjPL8KGO6c6wP8HhjvYZaoEwwYlw7qyIyVO1mxda/fcUQkTnhWCs65D4HDzgXtnJvmnPuq7u4MIM+rLNHqwqI8EoMBnp2x1u8oIhInImWfwtXAG36HiDSt0pI4o3cbXpi7njLNhyQiTcC8vNqXmeUDrzrneh1hmRHAg8BJzrlDnsZrZuOAcQA5OTmFkyZNalCe0tJS0tLSGvRavyz/qob/N7Ocq45PpLh9wnd+fTSO+VhpzPFBY/5uRowYMdc5V3TUBZ1znn0B+cCCIzzfB/gSKPi2P7OwsNA11NSpUxv8Wr+Ew2E36r4P3aj7PnThcPg7vz4ax3ysNOb4oDF/N8Ac9y0+Y33bfGRmHYAXgcudc8v8yhHpzIzLB3dk8aY9fLp2l99xRCTGeXlI6kRgOtDNzNab2dVmdo2ZXVO3yB1AK+BBM/vczOZ4lSXajenXjrSkEM9MX+13FBGJcSGvfrBz7uKjPD8WGOvV+8eS1KQQ5xfm8dzMNfz6zJ60Tk/yO5KIxKhIOfpIjuKKIR2pqnFMnKXDU0XEOyqFKHFc6zSGF7Tm2RlrqKzWfEgi4g2VQhS5amg+W/dW8MaCTX5HEZEYpVKIIsO7tqZTVioTpq32O4qIxCiVQhQJBIwrh3Tks7W7+HydDk8VkcanUogy3yvMIy0pxFNaWxARD6gUokx6cgLnF+bx6vyNbN1b7nccEYkxKoUodOWJ+VTVOJ7T7Kki0shUClGoU1Yqp3bP5pkZa9hfWeN3HBGJISqFKHVNcWd2llUyec46v6OISAxRKUSpE/JbUtixBeM/XElVjU5mE5HGoVKIYtcO78yGXft5df5Gv6OISIxQKUSxU7pnU5CTxkMlXxIOe3exJBGJHyqFKBYIGNcM78yyLaW8v2Sr33FEJAaoFKLc2X3bkdeiGfe9t+zrq9mJiDSYSiHKJQQD3HhqVxZs2MObCzb7HUdEopxKIQac2z+X41qncs87y6jRvgUROQYqhRgQCga45bRurNhayiufb/A7johEMZVCjBjdqw0922Zw7zvLKK/SWc4i0jAqhRgRCBi/PrMH67/az2MfrfQ7johEKZVCDBnaJYvRvdrwwNQv2bhrv99xRCQKqRRizK/O6EHYOe58fbHfUUQkCqkUYkz7lilcW9yZV+dvYqpOaBOR78izUjCzJ8xsq5ktOMzzZmb/Z2YrzGy+mQ3wKku8uba4M91y0rntxfmUVekQVRH59rxcU5gAjDrC86OBrnVf44CHPMwSV5JCQe6+oC/bSyt5bnGl33FEJIp4VgrOuQ+BnUdYZAzwtKs1A8g0s7Ze5Yk3vfOac/2ILkzbWM3zc9f7HUdEokTIx/fOBQ68Qsz6usc2HbygmY2jdm2CnJwcSkpKGvSGpaWlDX5tNOobdBQ0d/zyhXmUbVhKx4yg35GaRLz9nkFjjhdNMWY/S8EO8dghN4A758YD4wGKiopccXFxg96wpKSEhr42WpVVTeXOuWHGLzJevG4wORnJfkfyXDz+njXm+NAUY/bz6KP1QPsD7ucBulpMI8tIMh69oohd+yq54vFZ7NqnfQwicnh+lsIU4Iq6o5AGA7udc9/YdCTHrndecx69oohV28u49LGZbNtb4XckEYlQXh6SOhGYDnQzs/VmdrWZXWNm19Qt8jqwElgBPApc51UWgRO7ZPHIFYWs3FbG9x6axvIte/2OJCIRyLN9Cs65i4/yvAOu9+r95ZtGdMtm4rjBXD1hNmf/5WNuP6snF5/QgUDgULt3RCQe6YzmONOvfSZv3DiMoo4t+fVLCzj3oWlM+3K7rtomIoC/Rx+JT7Izknn6hwN58bMN3PXWEi55dCbd26Qzpl8uxd1aU5CTTlBrDyJxSaUQpwIB4/zCPM7s3ZYp8zbw3My1/OHNJfzhzSUkhQJ0yU4jOz2JFimJpCaFCBiYGcFA7VdyQpCUxNqvNhnJdGyVSn5WCkmh+DgXQiRWqRTiXLPEIN8/oQPfP6EDm3eX88mK7SzZvIdlW0rZVlrBsi2l7KusJuwg7BzOQXU4THlV+Bs/KzEUoH/7TE4uaM2/988lN7OZDyMSkWOhUpB6bZon873CvG+1bDjsKK+uoayiho279rN6RxlfrN/NjFU7uOutpdz99lKGF7Tmp6d2ZUCHFh4nF5HGolKQBgkEjJTEECmJIVqnJ9G3fSZj+uUCsG7nPp6fu56np6/mvAencfrxOfzXOb1o0zz2z6YWiXY6+kgaXfuWKfzstAI+/o9TuPX0bnywbBsj7/2AFz/VxHwikU6lIJ5JTQpx/YguvH3TcI5vl8HNk+fxq5e+oKK6xu9oInIYKgXxXIdWKTw3dhDXDO/MX2eu5eoJcyitqPY7logcgkpBmkQoGOC20d2554K+TF+5g0sencFXZZqcTyTSqBSkSX2vMI/xlxeyZPNernxyFnvLq/yOJCIHUClIkzu1Rw4PXjKARRv3MPapOZRXaR+DSKRQKYgvRvbM4Z4L+zJr9U5u+fs8zb0kEiEOWwpm9hczO7Epw0h8GdMvl9tGdee1+Zv48/sr/I4jIhx5TWE5cI+ZrTazP5hZv6YKJfFj3MnHcV7/XO59ZxlvfKFrLIn47bCl4Jy73zk3BBgO7ASeNLPFZnaHmRU0WUKJaWbGnef1pl/7TG59fj6rt5f5HUkkrh11n4Jzbo1z7g/Ouf7AJcC5wGLPk0ncSE4I8sClAwgGjBsmfqqT20R8dNRSMLMEMzvbzJ4D3gCWAd/zPJnEldzMZtx9QV8WbNjD/7y+xO84InHrSDuaTzOzJ4D1wDhqr6nc2Tn3fefcy00VUOLHaT1z+OHQTkyYtpp3F23xO45IXDrSmsKvgOlAD+fc2c6555xz2uArnrptdHd6tM3gthe/YKfOeBZpckfa0TzCOfeoc25nUwaS+JYYCnDvhX3Zvb+S219Z4Hcckbijk9ck4vRom8FNIwt4bf4mpszb6HcckbjiaSmY2SgzW2pmK8zstkM838HMpprZZ2Y238zO8DKPRI8fn3wc/dpncvvLC9i6p9zvOCJxw7NSMLMg8AAwGugJXGxmPQ9a7DfA5LrDXS8CHvQqj0SXUDDAPRf2pbyqhl+++IWmwRBpIl6uKQwEVjjnVjrnKoFJwJiDlnFARt3t5oC2FUi9zq3TuPX0bry3ZCsvfbbB7zgiccG8+gvMzM4HRjnnxtbdvxwY5Jy74YBl2gJvAy2AVGCkc27uIX7WOGoPiyUnJ6dw0qRJDcpUWlpKWlpag14braJ9zGHnuHNmORtLw9x5UjMyk4/+d0y0j7khNOb4cCxjHjFixFznXNHRlgs16Kd/O3aIxw5uoIuBCc65e8xsCPCMmfVyzoX/5UXOjQfGAxQVFbni4uIGBSopKaGhr41WsTDmjr1KOeP+j3h1SzqPXlGE2aH+af1TLIz5u9KY40NTjNnLzUfrgfYH3M/jm5uHrgYmAzjnpgPJQJaHmSQKdW6dxs//rRvvLt7KK59rC6OIl7wshdlAVzPrZGaJ1O5InnLQMmuBUwHMrAe1pbDNw0wSpX54UicGdMjkt1MW6mgkEQ95VgrOuWrgBuAtaifQm+ycW2hmvzOzc+oWuwX4kZnNAyYCVzkdZiKHEAwYd13Ql/1VNfzqpQU6GknEI17uU8A59zq1cyYd+NgdB9xeBAz1MoPEjtrNSAXc+foSpszbyJh+uX5HEok5OqNZosrVJx1H/683I+3VZiSRxqZSkKgSDBh3nd+XfZU1/EabkUQanUpBok6X7DRuOa2Atxdt0dxIIo1MpSBRaeyw2rmRfjtlIdv2VvgdRyRmqBQkKgUDxt0X9KndjPSy5kYSaSwqBYlaXbLTufm0At5auIVX52/yO45ITFApSFQbe1In+rbP5I5XFmgzkkgjUClIVAsFA9xzQR/KKmu4/WUdjSRyrFQKEvW6ZKfzs5EFvLlwM699oc1IIsdCpSAx4UfDvt6MtJA9FVpbEGkolYLEhFAwwN3n96G0vJpHv6igJqxiEGkIlYLEjK456fz2nJ58sb2G+99d5ncckaikUpCYcsnADgzLDfF/76/gnUVb/I4jEnVUChJTzIzLeybSJ685N//tc5Zu3ut3JJGoolKQmJMYNB6+rJCUpCBXPjGLjbv2+x1JJGqoFCQmtctsxoQfDKSsopqrnpzF7n1VfkcSiQoqBYlZPdpm8MjlhazaXsblT8xk175KvyOJRDyVgsS0E7tk8fBlhSzZtJeLH53JjlJNhSFyJCoFiXmn9sjhsSuLWLmtlO89NI0VW7XzWeRwVAoSF04uaM1ffzSI0opqzn1gGu/qcFWRQ1IpSNwo7NiSV244iQ6tUhj79Bx++eJ8Siuq/Y4lElFUChJXcjOb8cK1J/Lj4ccxafY6Tr2nhMmz12laDJE6npaCmY0ys6VmtsLMbjvMMhea2SIzW2hmf/UyjwhAckKQX47uwQvXnki7zGb84oX5nHpPCY99tFJHKEncC3n1g80sCDwAnAasB2ab2RTn3KIDlukK/BIY6pz7ysyyvcojcrABHVrw4rUn8uaCzTz28Sr++7XF/M8bSxiY35JhBVn0zm1Oj7YZtEpNxMz8jivSJDwrBWAgsMI5txLAzCYBY4BFByzzI+AB59xXAM65rR7mEfkGM2N077aM7t2WhRt389r8Tby7eAt/fHNp/TKJoQA5GUlkJCeQFAqQFAqSGAqQEAyQELS67wESQ0YoECApFKBN82TyWqTQKSuVLtlpBAMqFYkOXpZCLrDugPvrgUEHLVMAYGafAEHgP51zb3qYSeSwjm/XnOPbNecXo7qza18lCzfuYenmvWzZU87mPeWUVVRTUR2moirMrn2VVNU4qmrCVIcdldXh+tv7KqsprwrX/9zUxCB922dySvdsTj++De1bpvg4SpEjM68uX2hmFwCnO+fG1t2/HBjonPvJAcu8ClQBFwJ5wEdAL+fcroN+1jhgHEBOTk7hpEmTGpSptLSUtLS0Br02WmnMTc85R1kVbN8fZkNpmJW7wyzdWcP60tr/17pmBji1YwJFOUFCjbQG4feY/aAxfzcjRoyY65wrOtpyXq4prAfaH3A/D9h4iGVmOOeqgFVmthToCsw+cCHn3HhgPEBRUZErLi5uUKCSkhIa+tpopTFHjrU79vHGgk1MnLWWh+ftIzezGTeO7MJ5/XMJBY/tmI9IHbOXNGZveHn00Wygq5l1MrNE4CJgykHLvAyMADCzLGo3J630MJOIbzq0SuHHwzvz/i3FPHnVCWSlJfKL5+cz+v6PmL16p9/xRAAPS8E5Vw3cALwFLAYmO+cWmtnvzOycusXeAnaY2SJgKnCrc26HV5lEIkEgYIzons3L1w/l4csGsK+yhgsens6vX/pCJ9OJ77zcfIRz7nXg9YMeu+OA2w64ue5LJK6YGaN6tWVY19bc+84ynvxkFdO+3MFfLunP8e2a+x1P4pTOaBbxWWpSiNvP6smkcUPYX1nDuQ9O49kZa/yOJXFKpSASIQZ2asnrNw7jxM6t+M3LC7jjlQVU14SP/kKRRqRSEIkgLVMTefzKExh38nE8PX0NP3xqDnvKddU4aToqBZEIEwwYvzqjB/97Xm+mrdjOhQ9PZ+vecr9jSZxQKYhEqIsGduDJH5zAmh37uPDh6azbuc/vSBIHVAoiEWxY19Y8O3YQO8squeDh6bpqnHhOpSAS4Qo7tuBvPx5Cddhx4SMzWLBht9+RJIapFESiQI+2GTx/zRCaJQS5+NEZfLr2K78jSYxSKYhEifysVP7248G0TE3k8sdmMmOlTv6XxqdSEIkieS1SmPzjIbTNbMaVT8zig2Xb/I4kMUalIBJlcjKS+du4wXRuncaPnprD2ws3+x1JYohKQSQKtUpLYuKPBtOzXQbXPvcpMzZpIj1pHCoFkSjVPCWBZ8cOorBjCx6ZV8HkOeuO/iKRo1ApiESxtKQQT/1gIMe3CvKL5+fzzPTVfkeSKKdSEIlyzRKD3FiYxMgeOdz+ykLGf/il35EkiqkURGJAQsB46LIBnNWnLXe+voT7312OV9dfl9jm6UV2RKTpJAQD3H9Rf5ITgvzp3WXsq6rmtlHdMTO/o0kUUSmIxJBgwPjj9/qQnBDgkQ9WUl5Zw2/PPp5AQMUg345KQSTGBALG78f0ollCkEc/WsXe8mruPK83yQlBv6NJFFApiMQgs9prMmQkJ3DPO8tYvrWUhy8vJDezmd/RJMJpR7NIjDIzfnJqVx69oojV28s4+88f86GmxZCjUCmIxLjTeubw8g1DaZWayBVPzOL2lxewr1JnQMuheVoKZjbKzJaa2Qozu+0Iy51vZs7MirzMIxKvOrdO4x8/OYmrT+rEszPXcMb9HzF16Va/Y0kE8qwUzCwIPACMBnoCF5tZz0Mslw78FJjpVRYRgeSEILef1ZO/jh2MmfGDJ2dz1ZOzWLpZV3OTf/JyTWEgsMI5t9I5VwlMAsYcYrnfA38EdGVykSYwpHMr3rrpZH59Rg/mrv6K0+/7kLFPzWH26p064U08PfooFzhwhq71wKADFzCz/kB759yrZvZzD7OIyAESQwF+dPJxnF+Yx1PTVzNh2mrefXgLnbJSObd/LiN75NC9TbrOb4hD5tVfBmZ2AXC6c25s3f3LgYHOuZ/U3Q8A7wNXOedWm1kJ8HPn3JxD/KxxwDiAnJycwkmTJjUoU2lpKWlpaQ16bbTSmOPDsY65otoxa3M10zZWs3hnGID0RChoESQvLUBeeoCsZkbzJCMj0QhFQFno9/zdjBgxYq5z7qj7bb1cU1gPtD/gfh6w8YD76UAvoKTuNPw2wBQzO+fgYnDOjQfGAxQVFbni4uIGBSopKaGhr41WGnN8aIwxn173fcuecj5avp2Pl29j3vrdfLqyjIP/dkwMBkgKBUhKCJAUChIK1hZFQjBAQjBAKPj1bSMUCJAYCpCdnkS7zGa0bZ5MQU46XXPSSAo1/IQ6/Z694WUpzAa6mlknYANwEXDJ108653YDWV/fP9Kagog0nZyMZM4vzOP8wjwA9lVW8+XWMjbt3s+20gq2761kf1UNFdU1VFSHqagKUx0OU13jqKoJU1UTpjrsqKwOU14VprqmmvKqMHNW7+SrfVX175MQNApy0hlyXCtOLmjNwE4tddZ1BPCsFJxz1WZ2A/AWEASecM4tNLPfAXOcc1O8em8RaTwpiSF65zWnd17zY/5Z+ytr2LBrP0s272Hhxj18vnYXT09fw2MfryIlMcjpx7fh3P65DO2SRTACNlHFI0+nuXDOvQ68ftBjdxxm2WIvs4iI/5olBumSnUaX7DTO6tMOqC2KGat28PbCzbw6fxMvfbaBds2TuWpoPt8/oQPNmyX4nDq+aO4jEfFVs8QgI7plM6JbNr89+3imLtnK09PX1F8X4tLBHblmeGdapib6HTUuqBREJGIkJwQZ3bsto3u3ZcGG3Tz20Uoe+2glf525lrHDOnH1SZ1IT9aag5c095GIRKReuc2576L+vHXTyQzrmsV97y6n+K4SJs9eRzisk+y8olIQkYjWNSedhy4rZMoNQ8nPSuUXL8zn3IemsXJ3jd/RYpJKQUSiQp+8TJ6/Zgj3XtiXjbv28/vp5fzyxfl8VVbpd7SYolIQkahhZpw3II/3bxnO6fkhJs9Zzyn3aJNSY1IpiEjUSU9O4KLuSbz205Pokp3GL16YzwWPTGfxpj1+R4t6KgURiVrd22Twt3FDuOv8PqzaXsZZf/6Y/351EaUVuohQQ6kURCSqBQLGBUXtef+W4Xz/hPY8/skqTr2nhNfmb9JU4A2gUhCRmJCZksid5/bmhWtPJCstiev/+ilXPDGLVdvL/I4WVVQKIhJTBnRowSvXD+U/z+7J52t3cfqfPuTed5ZRXqVDWL8NlYKIxJxQMMBVQzvx3i3DGd27Df/33nJOv+9DSnRd6qNSKYhIzMrOSOb+i/rz3NhBBAPGVU/O5tpn5/LltlK/o0UslYKIxLyhXbJ448Zh3Hp6N6Yu3crIez/gJxM/Y8lmHcJ6ME2IJyJxISkU5PoRXbiwqD2Pf7yKZ6av5h/zNjKsaxaXDurIqT2ySQjq72SVgojEldbpSdw2ujvXDD+Op6evYeKstVzz7Fyy05M4u287RvVqw4AOLeL2Ij8qBRGJS5kpifz01K5cV9yZkqXbmDR7Hc9MX8PjH68iKy2RIZ2zGJjfgqL8lnTJToubtQiVgojEtVAwwMieOYzsmcPe8ipKlm7jnUVbmLVqJ/+Yt7F2mYCRn5VK59aptMlIpnV6EllpSaQnJ5AUCpCcECQ5IUBSKEgoaCQEAyQGA/W3E+q/1942i9y1EJWCiEid9OQEzu7bjrP7tsM5x7qd+/l07Vcs37qX5VtK+XJbGdO/3MGe8mObRiMUMEJBIykUJDs9ibaZzWibkUyHVin0aJtO9zYZtG2e7Et5qBRERA7BzOjQKoUOrVK+8Vx5VQ3bSysoq6ihvKruqzpMVXWYqpowVWFHVXWY6nCYypp/3q6qcbXP14SprnHsr6phy55yNu0uZ9HGPWwvrah/j+bNEijq2IKBnVoy6LhW9GqX0STjVimIiHxHyQlB8lp8syyO1Z7yKpZt3svizXtZuGE3s1bt5L0ltSfcZSSHGN0xQHFxo7/tv1ApiIhEiIzkBIryW1KU37L+sa17ypmxaicfLdtGi+ptnmfwdHe6mY2Bu8vaAAAFX0lEQVQys6VmtsLMbjvE8zeb2SIzm29m75lZRy/ziIhEm+yMZM7p2467LujL4Lbe/x3vWSmYWRB4ABgN9AQuNrOeBy32GVDknOsDPA/80as8IiJydF6uKQwEVjjnVjrnKoFJwJgDF3DOTXXO7au7OwPI8zCPiIgchZelkAusO+D++rrHDudq4A0P84iIyFF4uYHqUAfYHvIySGZ2GVAEDD/M8+OAcQA5OTmUlJQ0KFBpaWmDXxutNOb4oDHHh6YYs5elsB5of8D9PGDjwQuZ2Ujg18Bw51zFwc8DOOfGA+MBioqKXHEDj8kqKSmhoa+NVhpzfNCY40NTjNnLzUezga5m1snMEoGLgCkHLmBm/YFHgHOcc7r6hYiIzzwrBedcNXAD8BawGJjsnFtoZr8zs3PqFrsLSAP+bmafm9mUw/w4ERFpAp4e9Oqcex14/aDH7jjg9kgv319ERL4bc+6Q+34jlpltA9Y08OVZwPZGjBMNNOb4oDHHh2MZc0fnXOujLRR1pXAszGyOc67I7xxNSWOODxpzfGiKMcfHVSNERORbUSmIiEi9eCuF8X4H8IHGHB805vjg+Zjjap+CiIgcWbytKYiIyBHEbCmY2RNmttXMFhzwWEsze8fMltd9b+FnxsZ2mDHfZWZL6q5Z8ZKZZfqZsbEdaswHPPdzM3NmluVHNq8cbsxm9pO665csNLOYmob+MP+2+5nZjLoTX+eY2UA/MzYmM2tvZlPNbHHd7/PGusc9/wyL2VIAJgCjDnrsNuA951xX4L26+7FkAt8c8ztAr7prViwDftnUoTw2gW+OGTNrD5wGrG3qQE1gAgeN2cxGUDs1fR/n3PHA3T7k8tIEvvl7/iPwX865fsAdxNb1WKqBW5xzPYDBwPV116Px/DMsZkvBOfchsPOgh8cAT9Xdfgr49yYN5bFDjdk593bdlCMQg9esOMzvGeBPwC84zMy80ewwY74W+N+vJ5WMtbnEDjNmB3x9NfvmHGLCzWjlnNvknPu07vZeaqcKyqUJPsNithQOI8c5twlq/6MD2T7naWo/JA6uWVE3t9YG59w8v7M0oQJgmJnNNLMPzOwEvwM1gZuAu8xsHbVrRrG2FgyAmeUD/YGZNMFnWLyVQtwys19Tu0r6nN9ZvGRmKdROxX7H0ZaNMSGgBbWbGm4FJpvZoa5pEkuuBX7mnGsP/Ax43Oc8jc7M0oAXgJucc3ua4j3jrRS2mFlbgLrvMbWKfThmdiVwFnCpi/1jkDsDnYB5Zraa2s1ln5pZG19TeW898KKrNQsIUztPTiy7Enix7vbfqb0EcMwwswRqC+E559zX4/T8MyzeSmEKtf+QqPv+io9ZmoSZjQL+g9prVuw72vLRzjn3hXMu2zmX75zLp/bDcoBzbrPP0bz2MnAKgJkVAInE/mRxG/nn1RpPAZb7mKVR1a3lPQ4sds7de8BT3n+GOedi8guYCGwCqqj9YLgaaEXtHvvldd9b+p2zCca8gtprZX9e9/Ww3zm9HvNBz68GsvzO2QS/50TgWWAB8Clwit85m2DMJwFzgXnUbm8v9DtnI473JGp3pM8/4P/dM5riM0xnNIuISL1423wkIiJHoFIQEZF6KgUREamnUhARkXoqBRERqadSEGkEZpZpZtf5nUPkWKkURBpHJqBSkKinUhBpHP8LdK6b2/8uv8OINJROXhNpBHUzWb7qnOvlcxSRY6I1BRERqadSEBGReioFkcaxF0j3O4TIsVIpiDQC59wO4BMzW6AdzRLNtKNZRETqaU1BRETqqRRERKSeSkFEROqpFEREpJ5KQURE6qkURESknkpBRETqqRRERKTe/wfxuZF/OX1AjgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f2b18636240>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "V_lam = sympy.lambdify(model.x, V_eq.subs(model.p0))\n",
    "V_vals = [V_lam(*x) for x in res['x']]\n",
    "pl.plot(res['t'], V_vals)\n",
    "pl.xlabel('t')\n",
    "pl.ylabel('V')\n",
    "pl.grid()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Parameter Sensitivity"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAJgAAAAyBAMAAABSaw+aAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMA74lUMhCZZs3dIrurRHbLQ9+lAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACmElEQVRIDe2Xz2vUQBTHvyQds79iBcWT0j2IemrjX7Ap1HNTEaEg2yCCnnT1JL004EkvFr3YIrI3ERF7EAta6B78A7wKQoMHT6XdKhUEa32zP2rezGyTQg4e+iAJ8+aTD5OZl8kuRna/I5cY293C8YmLubhweeICTuSjkpZxRXbl3JlMcnHtrq+Bquw+ppsaZEgUPLGopRWZcwxWpEGGxB3gvJZWZOUW7G0NMiRWgFqg5hXZcAvuD5UxtX8Da57aochmQrh/VMbQFj9JFqsdiuxNFUNZarhE0GQ1TRZmlNHIUmW5PiYtgJ11AWppC1BegpOpND4Ac2mlQUV7JFLn1dSmor2p5ZXVxBNM+RpkSBQ9oe8QquzSxkfDrXpKjK77WlaVacBBEoeyg8xWlz2cs/9xzkoNOarZTENLssbVfNXRVJpZbEnWKHsI0CdWRFlkSdYkc0LAbgOjGWSMNcmKTXS2yKkMMsaaZC9JUqwChU/pNsZyWX399dcY34AvqysxrNAou7fVqPzyel2MZTJRdZ4XIywTWKPDiXp3KJdTwF4PY5nMCqxtu0nfge7Xwp5XLL3mZtAtLtlkbEcmHixQPPUFCksEnKbjMR1DbTr1o89Qe6Zh9bOcZSMDhmPCSCZ26Mpke7cDlXkm+8cqskmfbqKhu20EcAc8prszC/eRNyafgrFMZgebOBrIBbAi2x+4AHgWUt3E18nFWSZb81bxAjhLNduq06kqcUPMeTQRJ+mksEx2daN+owF8BkrLdC34BpFMTctT96cGY5lMMhTyFZEhq3tgiFuR7GOsSSZfXhn7vujl4J1kGGuS4b3E9t+CnBHcfiupJGuUdTe8clPCaZFkjbLkVpwmS7JGWZpgUP94vn8R8/zz+hehwq9PsjM3hQAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$$\\left[\\begin{matrix}0 & 0\\\\- x{\\left (t \\right )} & - \\operatorname{v_{x}}{\\left (t \\right )}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡  0      0   ⎤\n",
       "⎢             ⎥\n",
       "⎣-x(t)  -vₓ(t)⎦"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.f.jacobian(model.p)"
   ]
  }
 ],
 "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.6.3"
  },
  "name": ""
 },
 "nbformat": 4,
 "nbformat_minor": 1
}