{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# 10c. Numerical solution of delay differential equations\n", "\n", "
" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "text/html": [ "\n", "
\n", " \n", " Loading BokehJS ...\n", "
\n" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "application/javascript": [ "(function(root) {\n", " function now() {\n", " return new Date();\n", " }\n", "\n", " const force = true;\n", "\n", " if (typeof root._bokeh_onload_callbacks === \"undefined\" || force === true) {\n", " root._bokeh_onload_callbacks = [];\n", " root._bokeh_is_loading = undefined;\n", " }\n", "\n", "const JS_MIME_TYPE = 'application/javascript';\n", " const HTML_MIME_TYPE = 'text/html';\n", " const EXEC_MIME_TYPE = 'application/vnd.bokehjs_exec.v0+json';\n", " const CLASS_NAME = 'output_bokeh rendered_html';\n", "\n", " /**\n", " * Render data to the DOM node\n", " */\n", " function render(props, node) {\n", " const script = document.createElement(\"script\");\n", " node.appendChild(script);\n", " }\n", "\n", " /**\n", " * Handle when an output is cleared or removed\n", " */\n", " function handleClearOutput(event, handle) {\n", " const cell = handle.cell;\n", "\n", " const id = cell.output_area._bokeh_element_id;\n", " const server_id = cell.output_area._bokeh_server_id;\n", " // Clean up Bokeh references\n", " if (id != null && id in Bokeh.index) {\n", " Bokeh.index[id].model.document.clear();\n", " delete Bokeh.index[id];\n", " }\n", "\n", " if (server_id !== undefined) {\n", " // Clean up Bokeh references\n", " const cmd_clean = \"from bokeh.io.state import curstate; print(curstate().uuid_to_server['\" + server_id + \"'].get_sessions()[0].document.roots[0]._id)\";\n", " cell.notebook.kernel.execute(cmd_clean, {\n", " iopub: {\n", " output: function(msg) {\n", " const id = msg.content.text.trim();\n", " if (id in Bokeh.index) {\n", " Bokeh.index[id].model.document.clear();\n", " delete Bokeh.index[id];\n", " }\n", " }\n", " }\n", " });\n", " // Destroy server and session\n", " const cmd_destroy = \"import bokeh.io.notebook as ion; ion.destroy_server('\" + server_id + \"')\";\n", " cell.notebook.kernel.execute(cmd_destroy);\n", " }\n", " }\n", "\n", " /**\n", " * Handle when a new output is added\n", " */\n", " function handleAddOutput(event, handle) {\n", " const output_area = handle.output_area;\n", " const output = handle.output;\n", "\n", " // limit handleAddOutput to display_data with EXEC_MIME_TYPE content only\n", " if ((output.output_type != \"display_data\") || (!Object.prototype.hasOwnProperty.call(output.data, EXEC_MIME_TYPE))) {\n", " return\n", " }\n", "\n", " const toinsert = output_area.element.find(\".\" + CLASS_NAME.split(' ')[0]);\n", "\n", " if (output.metadata[EXEC_MIME_TYPE][\"id\"] !== undefined) {\n", " toinsert[toinsert.length - 1].firstChild.textContent = output.data[JS_MIME_TYPE];\n", " // store reference to embed id on output_area\n", " output_area._bokeh_element_id = output.metadata[EXEC_MIME_TYPE][\"id\"];\n", " }\n", " if (output.metadata[EXEC_MIME_TYPE][\"server_id\"] !== undefined) {\n", " const bk_div = document.createElement(\"div\");\n", " bk_div.innerHTML = output.data[HTML_MIME_TYPE];\n", " const script_attrs = bk_div.children[0].attributes;\n", " for (let i = 0; i < script_attrs.length; i++) {\n", " toinsert[toinsert.length - 1].firstChild.setAttribute(script_attrs[i].name, script_attrs[i].value);\n", " toinsert[toinsert.length - 1].firstChild.textContent = bk_div.children[0].textContent\n", " }\n", " // store reference to server id on output_area\n", " output_area._bokeh_server_id = output.metadata[EXEC_MIME_TYPE][\"server_id\"];\n", " }\n", " }\n", "\n", " function register_renderer(events, OutputArea) {\n", "\n", " function append_mime(data, metadata, element) {\n", " // create a DOM node to render to\n", " const toinsert = this.create_output_subarea(\n", " metadata,\n", " CLASS_NAME,\n", " EXEC_MIME_TYPE\n", " );\n", " this.keyboard_manager.register_events(toinsert);\n", " // Render to node\n", " const props = {data: data, metadata: metadata[EXEC_MIME_TYPE]};\n", " render(props, toinsert[toinsert.length - 1]);\n", " element.append(toinsert);\n", " return toinsert\n", " }\n", "\n", " /* Handle when an output is cleared or removed */\n", " events.on('clear_output.CodeCell', handleClearOutput);\n", " events.on('delete.Cell', handleClearOutput);\n", "\n", " /* Handle when a new output is added */\n", " events.on('output_added.OutputArea', handleAddOutput);\n", "\n", " /**\n", " * Register the mime type and append_mime function with output_area\n", " */\n", " OutputArea.prototype.register_mime_type(EXEC_MIME_TYPE, append_mime, {\n", " /* Is output safe? */\n", " safe: true,\n", " /* Index of renderer in `output_area.display_order` */\n", " index: 0\n", " });\n", " }\n", "\n", " // register the mime type if in Jupyter Notebook environment and previously unregistered\n", " if (root.Jupyter !== undefined) {\n", " const events = require('base/js/events');\n", " const OutputArea = require('notebook/js/outputarea').OutputArea;\n", "\n", " if (OutputArea.prototype.mime_types().indexOf(EXEC_MIME_TYPE) == -1) {\n", " register_renderer(events, OutputArea);\n", " }\n", " }\n", " if (typeof (root._bokeh_timeout) === \"undefined\" || force === true) {\n", " root._bokeh_timeout = Date.now() + 5000;\n", " root._bokeh_failed_load = false;\n", " }\n", "\n", " const NB_LOAD_WARNING = {'data': {'text/html':\n", " \"
\\n\"+\n", " \"

\\n\"+\n", " \"BokehJS does not appear to have successfully loaded. If loading BokehJS from CDN, this \\n\"+\n", " \"may be due to a slow or bad network connection. Possible fixes:\\n\"+\n", " \"

\\n\"+\n", " \"\\n\"+\n", " \"\\n\"+\n", " \"from bokeh.resources import INLINE\\n\"+\n", " \"output_notebook(resources=INLINE)\\n\"+\n", " \"\\n\"+\n", " \"
\"}};\n", "\n", " function display_loaded() {\n", " const el = document.getElementById(\"d12175cf-5d61-46f6-b21b-7852467e191e\");\n", " if (el != null) {\n", " el.textContent = \"BokehJS is loading...\";\n", " }\n", " if (root.Bokeh !== undefined) {\n", " if (el != null) {\n", " el.textContent = \"BokehJS \" + root.Bokeh.version + \" successfully loaded.\";\n", " }\n", " } else if (Date.now() < root._bokeh_timeout) {\n", " setTimeout(display_loaded, 100)\n", " }\n", " }\n", "\n", " function run_callbacks() {\n", " try {\n", " root._bokeh_onload_callbacks.forEach(function(callback) {\n", " if (callback != null)\n", " callback();\n", " });\n", " } finally {\n", " delete root._bokeh_onload_callbacks\n", " }\n", " console.debug(\"Bokeh: all callbacks have finished\");\n", " }\n", "\n", " function load_libs(css_urls, js_urls, callback) {\n", " if (css_urls == null) css_urls = [];\n", " if (js_urls == null) js_urls = [];\n", "\n", " root._bokeh_onload_callbacks.push(callback);\n", " if (root._bokeh_is_loading > 0) {\n", " console.debug(\"Bokeh: BokehJS is being loaded, scheduling callback at\", now());\n", " return null;\n", " }\n", " if (js_urls == null || js_urls.length === 0) {\n", " run_callbacks();\n", " return null;\n", " }\n", " console.debug(\"Bokeh: BokehJS not loaded, scheduling load and callback at\", now());\n", " root._bokeh_is_loading = css_urls.length + js_urls.length;\n", "\n", " function on_load() {\n", " root._bokeh_is_loading--;\n", " if (root._bokeh_is_loading === 0) {\n", " console.debug(\"Bokeh: all BokehJS libraries/stylesheets loaded\");\n", " run_callbacks()\n", " }\n", " }\n", "\n", " function on_error(url) {\n", " console.error(\"failed to load \" + url);\n", " }\n", "\n", " for (let i = 0; i < css_urls.length; i++) {\n", " const url = css_urls[i];\n", " const element = document.createElement(\"link\");\n", " element.onload = on_load;\n", " element.onerror = on_error.bind(null, url);\n", " element.rel = \"stylesheet\";\n", " element.type = \"text/css\";\n", " element.href = url;\n", " console.debug(\"Bokeh: injecting link tag for BokehJS stylesheet: \", url);\n", " document.body.appendChild(element);\n", " }\n", "\n", " for (let i = 0; i < js_urls.length; i++) {\n", " const url = js_urls[i];\n", " const element = document.createElement('script');\n", " element.onload = on_load;\n", " element.onerror = on_error.bind(null, url);\n", " element.async = false;\n", " element.src = url;\n", " console.debug(\"Bokeh: injecting script tag for BokehJS library: \", url);\n", " document.head.appendChild(element);\n", " }\n", " };\n", "\n", " function inject_raw_css(css) {\n", " const element = document.createElement(\"style\");\n", " element.appendChild(document.createTextNode(css));\n", " document.body.appendChild(element);\n", " }\n", "\n", " const js_urls = [\"https://cdn.bokeh.org/bokeh/release/bokeh-3.2.0.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-gl-3.2.0.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-widgets-3.2.0.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-tables-3.2.0.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-mathjax-3.2.0.min.js\"];\n", " const css_urls = [];\n", "\n", " const inline_js = [ function(Bokeh) {\n", " Bokeh.set_log_level(\"info\");\n", " },\n", "function(Bokeh) {\n", " }\n", " ];\n", "\n", " function run_inline_js() {\n", " if (root.Bokeh !== undefined || force === true) {\n", " for (let i = 0; i < inline_js.length; i++) {\n", " inline_js[i].call(root, root.Bokeh);\n", " }\n", "if (force === true) {\n", " display_loaded();\n", " }} else if (Date.now() < root._bokeh_timeout) {\n", " setTimeout(run_inline_js, 100);\n", " } else if (!root._bokeh_failed_load) {\n", " console.log(\"Bokeh: BokehJS failed to load within specified timeout.\");\n", " root._bokeh_failed_load = true;\n", " } else if (force !== true) {\n", " const cell = $(document.getElementById(\"d12175cf-5d61-46f6-b21b-7852467e191e\")).parents('.cell').data().cell;\n", " cell.output_area.append_execute_result(NB_LOAD_WARNING)\n", " }\n", " }\n", "\n", " if (root._bokeh_is_loading === 0) {\n", " console.debug(\"Bokeh: BokehJS loaded, going straight to plotting\");\n", " run_inline_js();\n", " } else {\n", " load_libs(css_urls, js_urls, function() {\n", " console.debug(\"Bokeh: BokehJS plotting callback run at\", now());\n", " run_inline_js();\n", " });\n", " }\n", "}(window));" ], "application/vnd.bokehjs_load.v0+json": "(function(root) {\n function now() {\n return new Date();\n }\n\n const force = true;\n\n if (typeof root._bokeh_onload_callbacks === \"undefined\" || force === true) {\n root._bokeh_onload_callbacks = [];\n root._bokeh_is_loading = undefined;\n }\n\n\n if (typeof (root._bokeh_timeout) === \"undefined\" || force === true) {\n root._bokeh_timeout = Date.now() + 5000;\n root._bokeh_failed_load = false;\n }\n\n const NB_LOAD_WARNING = {'data': {'text/html':\n \"
\\n\"+\n \"

\\n\"+\n \"BokehJS does not appear to have successfully loaded. If loading BokehJS from CDN, this \\n\"+\n \"may be due to a slow or bad network connection. Possible fixes:\\n\"+\n \"

\\n\"+\n \"\\n\"+\n \"\\n\"+\n \"from bokeh.resources import INLINE\\n\"+\n \"output_notebook(resources=INLINE)\\n\"+\n \"\\n\"+\n \"
\"}};\n\n function display_loaded() {\n const el = document.getElementById(\"d12175cf-5d61-46f6-b21b-7852467e191e\");\n if (el != null) {\n el.textContent = \"BokehJS is loading...\";\n }\n if (root.Bokeh !== undefined) {\n if (el != null) {\n el.textContent = \"BokehJS \" + root.Bokeh.version + \" successfully loaded.\";\n }\n } else if (Date.now() < root._bokeh_timeout) {\n setTimeout(display_loaded, 100)\n }\n }\n\n function run_callbacks() {\n try {\n root._bokeh_onload_callbacks.forEach(function(callback) {\n if (callback != null)\n callback();\n });\n } finally {\n delete root._bokeh_onload_callbacks\n }\n console.debug(\"Bokeh: all callbacks have finished\");\n }\n\n function load_libs(css_urls, js_urls, callback) {\n if (css_urls == null) css_urls = [];\n if (js_urls == null) js_urls = [];\n\n root._bokeh_onload_callbacks.push(callback);\n if (root._bokeh_is_loading > 0) {\n console.debug(\"Bokeh: BokehJS is being loaded, scheduling callback at\", now());\n return null;\n }\n if (js_urls == null || js_urls.length === 0) {\n run_callbacks();\n return null;\n }\n console.debug(\"Bokeh: BokehJS not loaded, scheduling load and callback at\", now());\n root._bokeh_is_loading = css_urls.length + js_urls.length;\n\n function on_load() {\n root._bokeh_is_loading--;\n if (root._bokeh_is_loading === 0) {\n console.debug(\"Bokeh: all BokehJS libraries/stylesheets loaded\");\n run_callbacks()\n }\n }\n\n function on_error(url) {\n console.error(\"failed to load \" + url);\n }\n\n for (let i = 0; i < css_urls.length; i++) {\n const url = css_urls[i];\n const element = document.createElement(\"link\");\n element.onload = on_load;\n element.onerror = on_error.bind(null, url);\n element.rel = \"stylesheet\";\n element.type = \"text/css\";\n element.href = url;\n console.debug(\"Bokeh: injecting link tag for BokehJS stylesheet: \", url);\n document.body.appendChild(element);\n }\n\n for (let i = 0; i < js_urls.length; i++) {\n const url = js_urls[i];\n const element = document.createElement('script');\n element.onload = on_load;\n element.onerror = on_error.bind(null, url);\n element.async = false;\n element.src = url;\n console.debug(\"Bokeh: injecting script tag for BokehJS library: \", url);\n document.head.appendChild(element);\n }\n };\n\n function inject_raw_css(css) {\n const element = document.createElement(\"style\");\n element.appendChild(document.createTextNode(css));\n document.body.appendChild(element);\n }\n\n const js_urls = [\"https://cdn.bokeh.org/bokeh/release/bokeh-3.2.0.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-gl-3.2.0.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-widgets-3.2.0.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-tables-3.2.0.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-mathjax-3.2.0.min.js\"];\n const css_urls = [];\n\n const inline_js = [ function(Bokeh) {\n Bokeh.set_log_level(\"info\");\n },\nfunction(Bokeh) {\n }\n ];\n\n function run_inline_js() {\n if (root.Bokeh !== undefined || force === true) {\n for (let i = 0; i < inline_js.length; i++) {\n inline_js[i].call(root, root.Bokeh);\n }\nif (force === true) {\n display_loaded();\n }} else if (Date.now() < root._bokeh_timeout) {\n setTimeout(run_inline_js, 100);\n } else if (!root._bokeh_failed_load) {\n console.log(\"Bokeh: BokehJS failed to load within specified timeout.\");\n root._bokeh_failed_load = true;\n } else if (force !== true) {\n const cell = $(document.getElementById(\"d12175cf-5d61-46f6-b21b-7852467e191e\")).parents('.cell').data().cell;\n cell.output_area.append_execute_result(NB_LOAD_WARNING)\n }\n }\n\n if (root._bokeh_is_loading === 0) {\n console.debug(\"Bokeh: BokehJS loaded, going straight to plotting\");\n run_inline_js();\n } else {\n load_libs(css_urls, js_urls, function() {\n console.debug(\"Bokeh: BokehJS plotting callback run at\", now());\n run_inline_js();\n });\n }\n}(window));" }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Colab setup ------------------\n", "import os, sys, subprocess\n", "if \"google.colab\" in sys.modules:\n", " cmd = \"pip install --upgrade watermark\"\n", " process = subprocess.Popen(cmd.split(), stdout=subprocess.PIPE, stderr=subprocess.PIPE)\n", " stdout, stderr = process.communicate()\n", "# ------------------------------\n", "\n", "import numpy as np\n", "import scipy.integrate\n", "import scipy.interpolate\n", "import scipy.optimize\n", "\n", "import bokeh.io\n", "import bokeh.plotting\n", "\n", "bokeh.io.output_notebook()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Numerical solution of DDEs\n", "\n", "In general, we may write a system of DDEs as\n", "\n", "\\begin{align}\n", "\\frac{\\mathrm{d}\\mathbf{u}}{\\mathrm{d}t} = &= \\mathbf{f}(t, \\mathbf{u}(t), \\mathbf{u}(t-\\tau_1), \\mathbf{u}(t-\\tau_2), \\ldots, \\tau_n),\n", "\\end{align}\n", "\n", "where $\\tau_1, \\tau_2, \\ldots, \\tau_n$ are the delays with $\\tau_1 < \\tau_2 < \\ldots < \\tau_n$. We define the longest delay to be $\\tau_n = \\tau$ for notational convenience. If we wish to have a solution for time $t \\ge t_0$, we need to specify initial conditions that go sufficiently far back in time,\n", "\n", "\\begin{align}\n", "y(t) = y_0(t) \\text{ for }t_0-\\tau < t \\le t_0.\n", "\\end{align}\n", "\n", "The approach we will take is the **method of steps**. Given $y_0(t)$, we integrate the differential equations from time $t_0$ to time $t_0 + \\tau$. We store that solution and then use it as the past values when we integrate again from time $t_0+\\tau$ to time $t_0 + 2\\tau$. We continue this process until we integrate over the desired time span.\n", "\n", "Because we need to access $y(t-\\tau_i)$ at time $t$, and $t$ can be arbitrary, we often need an **interpolation** of the solution from the previous time steps. We will use [SciPy's functional interface](https://docs.scipy.org/doc/scipy/reference/interpolate.html) for [B-spline interpolation](https://en.wikipedia.org/wiki/B-spline), using the [scipy.interpolate.splrep()](https://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.splrep.html) and [scipy.interpolate.splev()](https://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.splev.html) functions.\n", "\n", "We write the function below to solve DDEs with similar call signature as `scipy.integrate.odeint()`. \n", "\n", " ddeint(func, y0, t, tau, args=(), y0_args=(), n_time_points_per_step=200)\n", "\n", "The input function has a different call signature than `odeint()` because it needs to allow for a function that computes the values of `y` in the past. So, the call signature of the function `func` giving the right hand side of the system of DDEs is\n", " \n", " func(y, t, y_past, *args)\n", " \n", "Furthermore, instead of passing in a single initial condition, we need to pass in a _function_ that computes $y_0(t)$, as defined for initial conditions that go back in time. The function has call signature\n", "\n", " y0(t, *y0_args)\n", " \n", "We also need to specify `tau`, the longest time delay in the system. Finally, in addition to the additional arguments to `func` and to `y0`, we specify how many time points we wish to use for each step to build our interpolants. Note that even though we specify `tau` in the arguments of `ddeint()`, we still need to explicitly pass all time delays as extra `args` to `func`. This is because we could in general have more than one delay, and only the longest is pertinent for the integration." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "def ddeint(func, y0, t, tau, args=(), y0_args=(), n_time_points_per_step=None):\n", " \"\"\"Integrate a system of delay differential equations.\n", " \"\"\"\n", " if np.isscalar(tau):\n", " tau = np.array([tau])\n", " else:\n", " tau = np.array(tau)\n", "\n", " if (tau <= 0).any():\n", " raise RuntimeError(\"All tau's must be greater than zero.\")\n", "\n", " tau_short = np.min(tau)\n", " tau_long = np.max(tau)\n", "\n", " if n_time_points_per_step is None:\n", " n_time_points_per_step = max(\n", " int(1 + len(t) / (t.max() - t.min()) * tau_long), 20\n", " )\n", "\n", " t0 = t[0]\n", "\n", " # Past function for the first step\n", " y_past = lambda time_point: y0(time_point, *y0_args)\n", "\n", " # Integrate first step\n", " t_step = np.linspace(t0, t0 + tau_short, n_time_points_per_step)\n", " y = scipy.integrate.odeint(func, y_past(t0), t_step, args=(y_past,) + args)\n", "\n", " # Store result from integration\n", " y_dense = y.copy()\n", " t_dense = t_step.copy()\n", "\n", " # Get dimension of problem for convenience\n", " n = y.shape[1]\n", "\n", " # Integrate subsequent steps\n", " j = 1\n", " while t_step[-1] < t[-1]:\n", " t_start = max(t0, t_step[-1] - tau_long)\n", " i = np.searchsorted(t_dense, t_start, side=\"left\")\n", " t_interp = t_dense[i:]\n", " y_interp = y_dense[i:, :]\n", "\n", " # Make B-spline\n", " tck = [scipy.interpolate.splrep(t_interp, y_interp[:, i]) for i in range(n)]\n", "\n", " # Interpolant of y from previous step\n", " y_past = (\n", " lambda time_point: np.array(\n", " [scipy.interpolate.splev(time_point, tck[i]) for i in range(n)]\n", " )\n", " if time_point > t0\n", " else y0(time_point, *y0_args)\n", " )\n", "\n", " # Integrate this step\n", " t_step = np.linspace(\n", " t0 + j * tau_short, t0 + (j + 1) * tau_short, n_time_points_per_step\n", " )\n", " y = scipy.integrate.odeint(func, y[-1, :], t_step, args=(y_past,) + args)\n", "\n", " # Store the result\n", " y_dense = np.append(y_dense, y[1:, :], axis=0)\n", " t_dense = np.append(t_dense, t_step[1:])\n", "\n", " j += 1\n", "\n", " # Interpolate solution for returning\n", " y_return = np.empty((len(t), n))\n", " for i in range(n):\n", " tck = scipy.interpolate.splrep(t_dense, y_dense[:, i])\n", " y_return[:, i] = scipy.interpolate.splev(t, tck)\n", "\n", " return y_return" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now that we have a function to solve DDEs, let's solve the simple delay oscillator." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/html": [ "\n", "
\n" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "application/javascript": [ "(function(root) {\n", " function embed_document(root) {\n", " const docs_json = {\"2a536a33-db8e-463b-9534-628a3437465c\":{\"version\":\"3.2.0\",\"title\":\"Bokeh Application\",\"roots\":[{\"type\":\"object\",\"name\":\"Figure\",\"id\":\"p1001\",\"attributes\":{\"x_range\":{\"type\":\"object\",\"name\":\"Range1d\",\"id\":\"p1010\",\"attributes\":{\"end\":200.0}},\"y_range\":{\"type\":\"object\",\"name\":\"DataRange1d\",\"id\":\"p1003\"},\"x_scale\":{\"type\":\"object\",\"name\":\"LinearScale\",\"id\":\"p1011\"},\"y_scale\":{\"type\":\"object\",\"name\":\"LinearScale\",\"id\":\"p1012\"},\"title\":{\"type\":\"object\",\"name\":\"Title\",\"id\":\"p1008\"},\"renderers\":[{\"type\":\"object\",\"name\":\"GlyphRenderer\",\"id\":\"p1036\",\"attributes\":{\"data_source\":{\"type\":\"object\",\"name\":\"ColumnDataSource\",\"id\":\"p1030\",\"attributes\":{\"selected\":{\"type\":\"object\",\"name\":\"Selection\",\"id\":\"p1031\",\"attributes\":{\"indices\":[],\"line_indices\":[]}},\"selection_policy\":{\"type\":\"object\",\"name\":\"UnionRenderers\",\"id\":\"p1032\"},\"data\":{\"type\":\"map\",\"entries\":[[\"x\",{\"type\":\"ndarray\",\"array\":{\"type\":\"bytes\",\"data\":\"AAAAAAAAAACkAJECRArgP6QAkQJECvA/9oDZA2YP+D+kAJECRAoAQM1ANQPVDARA9oDZA2YPCEAfwX0E9xEMQKQAkQJEChBAuCDjgowLEkDNQDUD1QwUQOJgh4MdDhZA9oDZA2YPGEAKoSuErhAaQB/BfQT3ERxANOHPhD8THkCkAJECRAogQK4QukLoCiFAuCDjgowLIkDDMAzDMAwjQM1ANQPVDCRA11BeQ3kNJUDiYIeDHQ4mQOxwsMPBDidA9oDZA2YPKEAAkQJEChApQAqhK4SuECpAFbFUxFIRK0AfwX0E9xEsQCnRpkSbEi1ANOHPhD8TLkA+8fjE4xMvQKQAkQJECjBAqYilIpaKMECuELpC6AoxQLOYzmI6izFAuCDjgowLMkC+qPei3osyQMMwDMMwDDNAyLgg44KMM0DNQDUD1Qw0QNLISSMnjTRA11BeQ3kNNUDc2HJjy401QOJgh4MdDjZA5+ibo2+ONkDscLDDwQ43QPH4xOMTjzdA9oDZA2YPOED7CO4juI84QACRAkQKEDlABRkXZFyQOUAKoSuErhA6QBApQKQAkTpAFbFUxFIRO0AaOWnkpJE7QB/BfQT3ETxAJEmSJEmSPEAp0aZEmxI9QC5Zu2Ttkj1ANOHPhD8TPkA5aeSkkZM+QD7x+MTjEz9AQ3kN5TWUP0CkAJECRApAQKdEmxJtSkBAqYilIpaKQECszK8yv8pAQK4QukLoCkFAsVTEUhFLQUCzmM5iOotBQLbc2HJjy0FAuCDjgowLQkC7ZO2StUtCQL6o96Lei0JAwOwBswfMQkDDMAzDMAxDQMV0FtNZTENAyLgg44KMQ0DK/Crzq8xDQM1ANQPVDERA0IQ/E/5MREDSyEkjJ41EQNUMVDNQzURA11BeQ3kNRUDalGhTok1FQNzYcmPLjUVA3xx9c/TNRUDiYIeDHQ5GQOSkkZNGTkZA5+ibo2+ORkDpLKazmM5GQOxwsMPBDkdA7rS60+pOR0Dx+MTjE49HQPM8z/M8z0dA9oDZA2YPSED5xOMTj09IQPsI7iO4j0hA/kz4M+HPSEAAkQJEChBJQAPVDFQzUElABRkXZFyQSUAIXSF0hdBJQAqhK4SuEEpADeU1lNdQSkAQKUCkAJFKQBJtSrQp0UpAFbFUxFIRS0AX9V7Ue1FLQBo5aeSkkUtAHH1z9M3RS0AfwX0E9xFMQCIFiBQgUkxAJEmSJEmSTEAnjZw0ctJMQCnRpkSbEk1ALBWxVMRSTUAuWbtk7ZJNQDGdxXQW001ANOHPhD8TTkA2JdqUaFNOQDlp5KSRk05AO63utLrTTkA+8fjE4xNPQEA1A9UMVE9AQ3kN5TWUT0BFvRf1XtRPQKQAkQJEClBApSKWilgqUECnRJsSbUpQQKhmoJqBalBAqYilIpaKUECqqqqqqqpQQKzMrzK/ylBAre60utPqUECuELpC6ApRQLAyv8r8KlFAsVTEUhFLUUCydsnaJWtRQLOYzmI6i1FAtbrT6k6rUUC23NhyY8tRQLf+3fp361FAuCDjgowLUkC6QugKoStSQLtk7ZK1S1JAvIbyGsprUkC+qPei3otSQL/K/Crzq1JAwOwBswfMUkDBDgc7HOxSQMMwDMMwDFNAxFIRS0UsU0DFdBbTWUxTQMeWG1tubFNAyLgg44KMU0DJ2iVrl6xTQMr8KvOrzFNAzB4we8DsU0DNQDUD1QxUQM5iOovpLFRA0IQ/E/5MVEDRpkSbEm1UQNLISSMnjVRA0+pOqzutVEDVDFQzUM1UQNYuWbtk7VRA11BeQ3kNVUDZcmPLjS1VQNqUaFOiTVVA27Zt27ZtVUDc2HJjy41VQN76d+vfrVVA3xx9c/TNVUDgPoL7CO5VQOJgh4MdDlZA44KMCzIuVkDkpJGTRk5WQOXGlhtbblZA5+ibo2+OVkDoCqErhK5WQOksprOYzlZA6k6rO63uVkDscLDDwQ5XQO2StUvWLldA7rS60+pOV0Dw1r9b/25XQPH4xOMTj1dA8hrKayivV0DzPM/zPM9XQPVe1HtR71dA9oDZA2YPWED3ot6Lei9YQPnE4xOPT1hA+ubom6NvWED7CO4juI9YQPwq86vMr1hA/kz4M+HPWED/bv279e9YQACRAkQKEFlAArMHzB4wWUAD1QxUM1BZQAT3EdxHcFlABRkXZFyQWUAHOxzscLBZQAhdIXSF0FlACX8m/JnwWUAKoSuErhBaQAzDMAzDMFpADeU1lNdQWkAOBzsc7HBaQBApQKQAkVpAEUtFLBWxWkASbUq0KdFaQBOPTzw+8VpAFbFUxFIRW0AW01lMZzFbQBf1XtR7UVtAGRdkXJBxW0AaOWnkpJFbQBtbbmy5sVtAHH1z9M3RW0Aen3h84vFbQB/BfQT3EVxAIOOCjAsyXEAiBYgUIFJcQCMnjZw0clxAJEmSJEmSXEAla5esXbJcQCeNnDRy0lxAKK+hvIbyXEAp0aZEmxJdQCvzq8yvMl1ALBWxVMRSXUAtN7bc2HJdQC5Zu2Ttkl1AMHvA7AGzXUAxncV0FtNdQDK/yvwq811ANOHPhD8TXkA1A9UMVDNeQDYl2pRoU15AN0ffHH1zXkA5aeSkkZNeQDqL6Syms15AO63utLrTXkA8z/M8z/NeQD7x+MTjE19APxP+TPgzX0BANQPVDFRfQEJXCF0hdF9AQ3kN5TWUX0BEmxJtSrRfQEW9F/Ve1F9AR98cfXP0X0CkAJECRApgQKWRk0ZOGmBApSKWilgqYECms5jOYjpgQKdEmxJtSmBAp9WdVndaYECoZqCagWpgQKj3ot6LemBAqYilIpaKYECqGahmoJpgQKqqqqqqqmBAqzut7rS6YECszK8yv8pgQKxdsnbJ2mBAre60utPqYECuf7f+3fpgQK4QukLoCmFAr6G8hvIaYUCwMr/K/CphQLDDwQ4HO2FAsVTEUhFLYUCx5caWG1thQLJ2ydola2FAswfMHjB7YUCzmM5iOothQLQp0aZEm2FAtbrT6k6rYUC1S9YuWbthQLbc2HJjy2FAt23btm3bYUC3/t36d+thQLiP4D6C+2FAuCDjgowLYkC5seXGlhtiQLpC6AqhK2JAutPqTqs7YkC7ZO2StUtiQLz179a/W2JAvIbyGsprYkC9F/Ve1HtiQL6o96Lei2JAvjn65uibYkC/yvwq86tiQMBb/279u2JAwOwBswfMYkDBfQT3EdxiQMEOBzsc7GJAwp8Jfyb8YkDDMAzDMAxjQMPBDgc7HGNAxFIRS0UsY0DF4xOPTzxjQMV0FtNZTGNAxgUZF2RcY0DHlhtbbmxjQMcnHp94fGNAyLgg44KMY0DJSSMnjZxjQMnaJWuXrGNAymsor6G8Y0DK/Crzq8xjQMuNLTe23GNAzB4we8DsY0DMrzK/yvxjQM1ANQPVDGRAztE3R98cZEDOYjqL6SxkQM/zPM/zPGRA0IQ/E/5MZEDQFUJXCF1kQNGmRJsSbWRA0TdH3xx9ZEDSyEkjJ41kQNNZTGcxnWRA0+pOqzutZEDUe1HvRb1kQNUMVDNQzWRA1Z1Wd1rdZEDWLlm7ZO1kQNe/W/9u/WRA11BeQ3kNZUDY4WCHgx1lQNlyY8uNLWVA2QNmD5g9ZUDalGhTok1lQNola5esXWVA27Zt27ZtZUDcR3AfwX1lQNzYcmPLjWVA3Wl1p9WdZUDe+nfr361lQN6Lei/qvWVA3xx9c/TNZUDgrX+3/t1lQOA+gvsI7mVA4c+EPxP+ZUDiYIeDHQ5mQOLxiccnHmZA44KMCzIuZkDjE49PPD5mQOSkkZNGTmZA5TWU11BeZkDlxpYbW25mQOZXmV9lfmZA5+ibo2+OZkDneZ7neZ5mQOgKoSuErmZA6Zujb46+ZkDpLKazmM5mQOq9qPei3mZA6k6rO63uZkDr361/t/5mQOxwsMPBDmdA7AGzB8weZ0DtkrVL1i5nQO4juI/gPmdA7rS60+pOZ0DvRb0X9V5nQPDWv1v/bmdA8GfCnwl/Z0Dx+MTjE49nQPKJxycen2dA8hrKayivZ0Dzq8yvMr9nQPM8z/M8z2dA9M3RN0ffZ0D1XtR7Ue9nQPXv1r9b/2dA9oDZA2YPaED3EdxHcB9oQPei3ot6L2hA+DPhz4Q/aED5xOMTj09oQPlV5leZX2hA+ubom6NvaED6d+vfrX9oQPsI7iO4j2hA/JnwZ8KfaED8KvOrzK9oQP279e/Wv2hA/kz4M+HPaED+3fp3699oQP9u/bv172hAAAAAAAAAaUA=\"},\"shape\":[400],\"dtype\":\"float64\",\"order\":\"little\"}],[\"y\",{\"type\":\"ndarray\",\"array\":{\"type\":\"bytes\",\"data\":\"2vGgiL3AT7zUwVEJHjv5P5zAol/aQQRA3tOjiPviCEC5NzxR37ALQDsfHnXCYw1AWryvxDJrDkC8BTJ2yAoPQKnKyVx0aw9ADCT89QOmD0DQBatrfckPQNd04sD63g9AIzToH//rD0DyBnF/4vMPQD8ZEt+m+A9Aj8IG1JX7D0CI5fz0Lf0PQGpN45Xf/g9AHctvQQL9D0DBrQ66qAMQQEW3l6dB4A9ACt0uzrk3CkCRtxYSnhcCQJDO6vlHp/g/5/cUkfkD8T+BuxTlwj7oPw2aPNVLB+I/oXuX7ZU33D/go4Nk/EnXP50zOLQcM9Q/orkAhilE0j+oTjYN3A7RP9hNuoTKTdA/OYWGoMuqzz9keAVcyBTPPxww161it84/nA+Q4lR9zj/g60JME1nOP7jFhNYzQ84/C3aavyAzzj+4b9PdPjPOPwdIOqaBW9A/U8Khg7Op1j8x0M5lGObiPwBGNB0OeO8/M5C831Xq9z8gdCraFTsAQMbefv0sCARAhma8b2gQB0BmPy0pqFIJQJ46I7c67wpAQjKYkZQNDEBJUXjKstAMQA8WrEf0Uw1AnJedTGyrDUC41bDsSuUNQJMrrxFcCw5APgupZ0ckDkATBJOMiDQOQGps9WIRPw5A8/Yf+PZFDkA+6qP9zj4OQAaGfa2k2w1Av0evV+RPDEDj5+DevswIQLbZpiQ23QNArCeYy7AE/j+qjDUhmCz2PyNE4MdbavA/CqW8cODF6D/4wdsT6EbjPzyc/O/wKd8/IsDMjrc82j+XiK75SvjWP43G7hZyz9Q/qLEcMN5i0z/7PxUF2nLSP6tA6rYM1dE/bRVmJHFt0T/6JECnbinRP+XmcIwY/dA/xIEqGWni0D/SENRLt+XQP/gsgWRHV9E/kX+x1PcH0z/SV731t37XP3A40HQub+A/SYYOlC6Y6D9MsC9c6DryPwH2DndYcPk/4JKbFAFwAECFJ5uTTs0DQPDhG9RikgZAi53ZpsC1CEDm94ZthkoKQE69VutdbQtAMXte8Uk6DEAUs8d4tMgMQNAusFizKg1ALYmilKRtDUDLFdCbDJsNQPynVAObuQ1AgLRqfk3NDUBcXPoW69MNQOTSxlixvQ1Apk3/dTJdDUBdjc2mIE0MQIDQZ8YO+glAr7hmZeQ2BkC9UyLmXMABQEAhKxswLvs/q02VP+t49D9qgUB+Q/buP+ajxnNv2uc/xlVvD0Hq4j+KI4BvRhHfPzZPP/Xzedo/sxb+cB1g1z+NguNNlUnVP1T+PVpS4tM/2eqQHi/x0j87MPZxc0/SP757BmEt49E/1Fga80Wc0T/zxRuP1XTRP746BaEJd9E/YuV+qd7N0T9H3bT5s+rSP//wPFViuNU/gSh80V242z8ddmVFg2njPyzbw9pURuw/vTzrFQgt9D/9CGOV3i77P93hh/R+GwFAVd09vLQ/BEDC6lhsAdQGQCW29ANY0whAsX6tBvFPCkBXVSsLsmMLQBEdZsjNJwxACIIpTYqxDEDhB5kNXhENQHFIMpl+Uw1A7cU7oHaADUCXWcJiQJ0NQL1KjXVUqg1Az7SqLBKgDUDb6MrW+2UNQHo9fgT+wQxAS3GGa1BHC0DEwE4EcoUIQMJRaFx4mgRAwJ6pyn5VAECO2RaqUwj5PyU6i4tl+vI/lUGO5E/87D/wXcsVM5fmPyv6nXrFIeI/5e4zD5cg3j98xKRP8e/ZP092hiDOFtc/L+0z/HAo1T/Ka14zydnTP4bAarKB99I/NMlDhRVf0j+rv47nOvrRPyvVlmj2vNE/gMSm9kOn0T+TYfjwCc3RP0VXca8Ua9I/HmxA/pYO1D+zhZYMKsHXPyW8dzq7Fd8/MU3zbv/s5T+KCHjxq5TvP1iuoc74DfY/5WFva+kH/T85hre6NekBQDswLQ5f4gRAPLrnfspLB0CZUt8JDicJQGNxy9rghwpAh/gRqWiHC0DARMU8Zj0MQMKeFw5zvQxAwTxgf6oWDUDt2BWUB1QNQEN7yvjUfA1A/8qXuouUDUA/a/nhRZkNQJa172ILfw1A19yJmHAlDUBXDoQXJUYMQFWxK/5FbwpA3+t+OONRB0Dr1/S8wEUDQP089SxZRv4/pCCdZTUm9z+edJSmr5rxP4Eb4lnzEOs/gr1dzzpJ5T//ikNaoULhP+Ujxggz+tw/3hDrbNsv2T+6xT9r9JrWP0ykZWHU2dQ/qi84Uxep0z940TZp2trSPyIuPj2KUNI/6oyOLiP30T+vWotkL8bRP19PicrhwtE/Rsb59scK0j/TJGU/UO3SP+yjUAVlFdU/WPn7IxGv2T+f8MX9qjHhP5LKodJzcOg/YRj4R+t18T/EruyloPb3P4Oh+Y+J7P4/3Y/AuZK/AkC75ZCshpAFQMx7n+si0QdAQJp8+NCICUCZytv0YM0KQH0v0puytwtAmWBsPkNeDEBPefpaTdMMQKq1X8OjJA1A+nZlHR9cDUCNLUEO6n8NQOFFib0Rkg1AEvE4fUOODUC4hAll4GMNQAzgrOSg6QxAiTAieN/MC0DyLbsIW5oJQIPRtylsLAZAsvqQYG8KAkBTPYIZdBb8Px7qb29UbfU/aDpgQzdW8D/GDf3gzUXpP3VdtgVcC+Q/6KXuw7Fp4D+j1pNcw9TbP8l3U92latg/2lFStRkX1j98NbpMK4LUP1cJ16pZb9M/Zw9NmKG10j+lqPRPBjrSP8lnbxGb7NE/lWrbmVvI0T/YLKPrYdjRP1KaveN6RNI/GjJx279z0z8clBdG5CrWPzhtdXrsxNs/i7o8Svb54j8eYia4oR7rP2OuF26SMvM/TQz7GWzm+T/nQBWbC2gAQGfAfk82kwNA2Mb0IfI7BkCPA4ueoVQIQFN1Wbv/6QlAjamJGVwTC0DvYqmxJukLQG/ZBxOvgAxAuBixCOPqDEDLxDVgZDQNQOVX3z/4ZQ1AgF+ljq+EDUC6I6LvCJENQBOvbc8QhA1AXbPnpkZHDUBjIIcbLqgMQDgStvTzRgtAzkbUxli4CEAuO5uW+QMFQMS63pd32gBA86epuFYL+j86rZRfpdbzP6Z9uJCBWO4/MTcfY5Gf5z9+4wRvGObiP2wOIPvUQN8/x/fYjH3C2j8flqcSpLDXP8DvDLFemdU/2ZZEEXct1D/gHhBVozbTP7Xj+6dEkNI/3dmQkski0j9a1zLgCOHRPzGWEDZrytE/U9zVzLXv0T+YW5WFvobSP00WCzAEDdQ/GYExJTFx1z+D5EKDKifeP14Fz7559+Q/7lj+EEgF7j8UVN+phQP1P89m934I3Ps/zrctruRVAUCmVpBahV8EQCbKAnd13wZABRWLs0XRCECK+Sm35UUKQALW1MOTVQtAiYIM5xUYDEAk1asohqEMQCxDGpuMAQ1A8Fw4kK1DDUAbbd5blm8NQPnhNRoTiQ1Ao06+uU2PDUCUWTBir3cNQFL0s/uVJQ1A99iQzG5aDEA6Mh9Pq6wKQMnvYCBkwgdAS2MKYfrWA0AA1oW+wWr/P1O0GZVXIfg/nlrV0vFf8j9TM0m+uTfsP+x2hFLvHuY/rTtEhynb4T+R62yHHdLdPzTuJ+giyNk/OPgCMVwG1z8wZbil2CXVP0enlu1I39M/H41fAgwC0z88gsRmqG3SPw09I7BcDdI/yc5nGmfX0T8ga5uU5s7RPzIQ+SAcDdI/GfXhZMbX0j8sRkyuFcjUP3hG9wQS8tg/5OnxPrJ04D/p3taewDLnP+p6hRU7k/A/rbGhY3vl9j+ih5c7XtH9P+D4grfqPAJArBlvtvghBUCUckjY/ngHQDRip1hQRQlAhI2xpwSbCkBKtlTDwZILQHCzIshhQwxAdZNlHcq/DEAYX8OtcRYNQNLTrC60UQ1AyFsnq0F4DUBY5q6FxYwNQP5yQRYJjA1AUVzUV7pnDUCX5altbPsMQA0UVq4R/AtAdHWEFBr6CUChf88U87gGQP5F/MadpwJA44qUKAQ4/T85ebFrMVj2P8zKJgbTB/E/xtHN0atH6j86oC1UScLkP/7IV3u86eA/uwsm5K6G3D/MoYMk4eXYP7jYARVobNY/g0mQ/0W91D/JIrnedJjTPw2+D89d0tI/sHvjMylO0j+mXSrvWvrRPy7L6IJK0NE/nK9UHP/X0T8bCOkhTDTSP+4q0NUFPtM/zT5m266t1T9s9m/YcL3aPx5HK6kYDeI/yGMcFZCv6T9SR1rflD/yPzzmc9fK0/g/+KlnpErB/z9K7QKqMBsDQJ8GhHp12QVAab3SftoHCEA=\"},\"shape\":[400],\"dtype\":\"float64\",\"order\":\"little\"}]]}}},\"view\":{\"type\":\"object\",\"name\":\"CDSView\",\"id\":\"p1037\",\"attributes\":{\"filter\":{\"type\":\"object\",\"name\":\"AllIndices\",\"id\":\"p1038\"}}},\"glyph\":{\"type\":\"object\",\"name\":\"Line\",\"id\":\"p1033\",\"attributes\":{\"x\":{\"type\":\"field\",\"field\":\"x\"},\"y\":{\"type\":\"field\",\"field\":\"y\"},\"line_color\":\"#1f77b4\",\"line_width\":2}},\"nonselection_glyph\":{\"type\":\"object\",\"name\":\"Line\",\"id\":\"p1034\",\"attributes\":{\"x\":{\"type\":\"field\",\"field\":\"x\"},\"y\":{\"type\":\"field\",\"field\":\"y\"},\"line_color\":\"#1f77b4\",\"line_alpha\":0.1,\"line_width\":2}},\"muted_glyph\":{\"type\":\"object\",\"name\":\"Line\",\"id\":\"p1035\",\"attributes\":{\"x\":{\"type\":\"field\",\"field\":\"x\"},\"y\":{\"type\":\"field\",\"field\":\"y\"},\"line_color\":\"#1f77b4\",\"line_alpha\":0.2,\"line_width\":2}}}}],\"toolbar\":{\"type\":\"object\",\"name\":\"Toolbar\",\"id\":\"p1009\",\"attributes\":{\"tools\":[{\"type\":\"object\",\"name\":\"PanTool\",\"id\":\"p1023\"},{\"type\":\"object\",\"name\":\"WheelZoomTool\",\"id\":\"p1024\"},{\"type\":\"object\",\"name\":\"BoxZoomTool\",\"id\":\"p1025\",\"attributes\":{\"overlay\":{\"type\":\"object\",\"name\":\"BoxAnnotation\",\"id\":\"p1026\",\"attributes\":{\"syncable\":false,\"level\":\"overlay\",\"visible\":false,\"left_units\":\"canvas\",\"right_units\":\"canvas\",\"bottom_units\":\"canvas\",\"top_units\":\"canvas\",\"line_color\":\"black\",\"line_alpha\":1.0,\"line_width\":2,\"line_dash\":[4,4],\"fill_color\":\"lightgrey\",\"fill_alpha\":0.5}}}},{\"type\":\"object\",\"name\":\"SaveTool\",\"id\":\"p1027\"},{\"type\":\"object\",\"name\":\"ResetTool\",\"id\":\"p1028\"},{\"type\":\"object\",\"name\":\"HelpTool\",\"id\":\"p1029\"}]}},\"left\":[{\"type\":\"object\",\"name\":\"LinearAxis\",\"id\":\"p1018\",\"attributes\":{\"ticker\":{\"type\":\"object\",\"name\":\"BasicTicker\",\"id\":\"p1019\",\"attributes\":{\"mantissas\":[1,2,5]}},\"formatter\":{\"type\":\"object\",\"name\":\"BasicTickFormatter\",\"id\":\"p1020\"},\"axis_label\":\"x\",\"major_label_policy\":{\"type\":\"object\",\"name\":\"AllLabels\",\"id\":\"p1021\"}}}],\"below\":[{\"type\":\"object\",\"name\":\"LinearAxis\",\"id\":\"p1013\",\"attributes\":{\"ticker\":{\"type\":\"object\",\"name\":\"BasicTicker\",\"id\":\"p1014\",\"attributes\":{\"mantissas\":[1,2,5]}},\"formatter\":{\"type\":\"object\",\"name\":\"BasicTickFormatter\",\"id\":\"p1015\"},\"axis_label\":\"time\",\"major_label_policy\":{\"type\":\"object\",\"name\":\"AllLabels\",\"id\":\"p1016\"}}}],\"center\":[{\"type\":\"object\",\"name\":\"Grid\",\"id\":\"p1017\",\"attributes\":{\"axis\":{\"id\":\"p1013\"}}},{\"type\":\"object\",\"name\":\"Grid\",\"id\":\"p1022\",\"attributes\":{\"dimension\":1,\"axis\":{\"id\":\"p1018\"}}}],\"frame_width\":450,\"frame_height\":200}}]}};\n", " const render_items = [{\"docid\":\"2a536a33-db8e-463b-9534-628a3437465c\",\"roots\":{\"p1001\":\"ff9658d0-abaf-4dbd-8e53-49d5ea4b9f36\"},\"root_ids\":[\"p1001\"]}];\n", " root.Bokeh.embed.embed_items_notebook(docs_json, render_items);\n", " }\n", " if (root.Bokeh !== undefined) {\n", " embed_document(root);\n", " } else {\n", " let attempts = 0;\n", " const timer = setInterval(function(root) {\n", " if (root.Bokeh !== undefined) {\n", " clearInterval(timer);\n", " embed_document(root);\n", " } else {\n", " attempts++;\n", " if (attempts > 100) {\n", " clearInterval(timer);\n", " console.log(\"Bokeh: ERROR: Unable to run BokehJS code because BokehJS library is missing\");\n", " }\n", " }\n", " }, 10, root)\n", " }\n", "})(window);" ], "application/vnd.bokehjs_exec.v0+json": "" }, "metadata": { "application/vnd.bokehjs_exec.v0+json": { "id": "p1001" } }, "output_type": "display_data" } ], "source": [ "# Right hand side of simple delay oscillator\n", "def delay_rhs(x, t, x_past, beta, n, tau):\n", " return beta / (1 + x_past(t - tau) ** n) - x\n", "\n", "\n", "# Initially, we just have no gene product\n", "def x0(t):\n", " return 0\n", "\n", "\n", "# Specify parameters (All dimensionless)\n", "beta = 4\n", "n = 2\n", "tau = 10\n", "\n", "# Time points we want\n", "t = np.linspace(0, 200, 400)\n", "\n", "# Perform the integration\n", "x = ddeint(delay_rhs, x0, t, tau, args=(beta, n, tau))\n", "\n", "# Plot the result\n", "p = bokeh.plotting.figure(\n", " frame_width=450,\n", " frame_height=200,\n", " x_axis_label=\"time\",\n", " y_axis_label=\"x\",\n", " x_range=[t.min(), t.max()],\n", ")\n", "p.line(t, x[:, 0], line_width=2)\n", "\n", "bokeh.io.show(p)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `ddeint()` function is available as `biocircuits.ddeint()` in the biocircuits package for future use." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Computing environment" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Python implementation: CPython\n", "Python version : 3.11.3\n", "IPython version : 8.12.0\n", "\n", "numpy : 1.24.3\n", "scipy : 1.10.1\n", "bokeh : 3.2.0\n", "jupyterlab: 3.6.3\n", "\n" ] } ], "source": [ "%load_ext watermark\n", "%watermark -v -p numpy,scipy,bokeh,jupyterlab" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.3" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": true, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false }, "varInspector": { "cols": { "lenName": 16, "lenType": 16, "lenVar": 40 }, "kernels_config": { "python": { "delete_cmd_postfix": "", "delete_cmd_prefix": "del ", "library": "var_list.py", "varRefreshCmd": "print(var_dic_list())" }, "r": { "delete_cmd_postfix": ") ", "delete_cmd_prefix": "rm(", "library": "var_list.r", "varRefreshCmd": "cat(var_dic_list()) " } }, "types_to_exclude": [ "module", "function", "builtin_function_or_method", "instance", "_Feature" ], "window_display": false } }, "nbformat": 4, "nbformat_minor": 4 }