1012 lines
76 KiB
Plaintext
1012 lines
76 KiB
Plaintext
{
|
|
"cells": [
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": 1,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"%matplotlib notebook\n",
|
|
"import numpy as np\n",
|
|
"from scipy.integrate import odeint\n",
|
|
"import matplotlib.pyplot as plt\n",
|
|
"from matplotlib import animation\n",
|
|
"#from matplotlib import rc\n",
|
|
"#import seaborn as sbs\n",
|
|
"#from scipy.integrate import quad\n",
|
|
"#from scipy.optimize import brentq\n",
|
|
"#rc('font', **{'family': 'serif', 'serif': ['Computer Modern'], 'size': 20})\n",
|
|
"#rc('text', usetex=True)\n",
|
|
"#rc('animation', html='html5')"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": 2,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"# The potential and its first derivative, as callables.\n",
|
|
"V = lambda x: 0.5 * x**2 * (0.5 * x**2 - 1)\n",
|
|
"dVdx = lambda x: x**3 - x"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": 3,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"# The potential energy function on a grid of x-points.\n",
|
|
"xgrid = np.linspace(-1.5, 1.5, 100)\n",
|
|
"Vgrid = V(xgrid)"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": 4,
|
|
"metadata": {
|
|
"scrolled": false
|
|
},
|
|
"outputs": [
|
|
{
|
|
"data": {
|
|
"application/javascript": [
|
|
"/* Put everything inside the global mpl namespace */\n",
|
|
"window.mpl = {};\n",
|
|
"\n",
|
|
"\n",
|
|
"mpl.get_websocket_type = function() {\n",
|
|
" if (typeof(WebSocket) !== 'undefined') {\n",
|
|
" return WebSocket;\n",
|
|
" } else if (typeof(MozWebSocket) !== 'undefined') {\n",
|
|
" return MozWebSocket;\n",
|
|
" } else {\n",
|
|
" alert('Your browser does not have WebSocket support.' +\n",
|
|
" 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n",
|
|
" 'Firefox 4 and 5 are also supported but you ' +\n",
|
|
" 'have to enable WebSockets in about:config.');\n",
|
|
" };\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n",
|
|
" this.id = figure_id;\n",
|
|
"\n",
|
|
" this.ws = websocket;\n",
|
|
"\n",
|
|
" this.supports_binary = (this.ws.binaryType != undefined);\n",
|
|
"\n",
|
|
" if (!this.supports_binary) {\n",
|
|
" var warnings = document.getElementById(\"mpl-warnings\");\n",
|
|
" if (warnings) {\n",
|
|
" warnings.style.display = 'block';\n",
|
|
" warnings.textContent = (\n",
|
|
" \"This browser does not support binary websocket messages. \" +\n",
|
|
" \"Performance may be slow.\");\n",
|
|
" }\n",
|
|
" }\n",
|
|
"\n",
|
|
" this.imageObj = new Image();\n",
|
|
"\n",
|
|
" this.context = undefined;\n",
|
|
" this.message = undefined;\n",
|
|
" this.canvas = undefined;\n",
|
|
" this.rubberband_canvas = undefined;\n",
|
|
" this.rubberband_context = undefined;\n",
|
|
" this.format_dropdown = undefined;\n",
|
|
"\n",
|
|
" this.image_mode = 'full';\n",
|
|
"\n",
|
|
" this.root = $('<div/>');\n",
|
|
" this._root_extra_style(this.root)\n",
|
|
" this.root.attr('style', 'display: inline-block');\n",
|
|
"\n",
|
|
" $(parent_element).append(this.root);\n",
|
|
"\n",
|
|
" this._init_header(this);\n",
|
|
" this._init_canvas(this);\n",
|
|
" this._init_toolbar(this);\n",
|
|
"\n",
|
|
" var fig = this;\n",
|
|
"\n",
|
|
" this.waiting = false;\n",
|
|
"\n",
|
|
" this.ws.onopen = function () {\n",
|
|
" fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n",
|
|
" fig.send_message(\"send_image_mode\", {});\n",
|
|
" if (mpl.ratio != 1) {\n",
|
|
" fig.send_message(\"set_dpi_ratio\", {'dpi_ratio': mpl.ratio});\n",
|
|
" }\n",
|
|
" fig.send_message(\"refresh\", {});\n",
|
|
" }\n",
|
|
"\n",
|
|
" this.imageObj.onload = function() {\n",
|
|
" if (fig.image_mode == 'full') {\n",
|
|
" // Full images could contain transparency (where diff images\n",
|
|
" // almost always do), so we need to clear the canvas so that\n",
|
|
" // there is no ghosting.\n",
|
|
" fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n",
|
|
" }\n",
|
|
" fig.context.drawImage(fig.imageObj, 0, 0);\n",
|
|
" };\n",
|
|
"\n",
|
|
" this.imageObj.onunload = function() {\n",
|
|
" fig.ws.close();\n",
|
|
" }\n",
|
|
"\n",
|
|
" this.ws.onmessage = this._make_on_message_function(this);\n",
|
|
"\n",
|
|
" this.ondownload = ondownload;\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype._init_header = function() {\n",
|
|
" var titlebar = $(\n",
|
|
" '<div class=\"ui-dialog-titlebar ui-widget-header ui-corner-all ' +\n",
|
|
" 'ui-helper-clearfix\"/>');\n",
|
|
" var titletext = $(\n",
|
|
" '<div class=\"ui-dialog-title\" style=\"width: 100%; ' +\n",
|
|
" 'text-align: center; padding: 3px;\"/>');\n",
|
|
" titlebar.append(titletext)\n",
|
|
" this.root.append(titlebar);\n",
|
|
" this.header = titletext[0];\n",
|
|
"}\n",
|
|
"\n",
|
|
"\n",
|
|
"\n",
|
|
"mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n",
|
|
"\n",
|
|
"}\n",
|
|
"\n",
|
|
"\n",
|
|
"mpl.figure.prototype._root_extra_style = function(canvas_div) {\n",
|
|
"\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype._init_canvas = function() {\n",
|
|
" var fig = this;\n",
|
|
"\n",
|
|
" var canvas_div = $('<div/>');\n",
|
|
"\n",
|
|
" canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n",
|
|
"\n",
|
|
" function canvas_keyboard_event(event) {\n",
|
|
" return fig.key_event(event, event['data']);\n",
|
|
" }\n",
|
|
"\n",
|
|
" canvas_div.keydown('key_press', canvas_keyboard_event);\n",
|
|
" canvas_div.keyup('key_release', canvas_keyboard_event);\n",
|
|
" this.canvas_div = canvas_div\n",
|
|
" this._canvas_extra_style(canvas_div)\n",
|
|
" this.root.append(canvas_div);\n",
|
|
"\n",
|
|
" var canvas = $('<canvas/>');\n",
|
|
" canvas.addClass('mpl-canvas');\n",
|
|
" canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n",
|
|
"\n",
|
|
" this.canvas = canvas[0];\n",
|
|
" this.context = canvas[0].getContext(\"2d\");\n",
|
|
"\n",
|
|
" var backingStore = this.context.backingStorePixelRatio ||\n",
|
|
"\tthis.context.webkitBackingStorePixelRatio ||\n",
|
|
"\tthis.context.mozBackingStorePixelRatio ||\n",
|
|
"\tthis.context.msBackingStorePixelRatio ||\n",
|
|
"\tthis.context.oBackingStorePixelRatio ||\n",
|
|
"\tthis.context.backingStorePixelRatio || 1;\n",
|
|
"\n",
|
|
" mpl.ratio = (window.devicePixelRatio || 1) / backingStore;\n",
|
|
"\n",
|
|
" var rubberband = $('<canvas/>');\n",
|
|
" rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n",
|
|
"\n",
|
|
" var pass_mouse_events = true;\n",
|
|
"\n",
|
|
" canvas_div.resizable({\n",
|
|
" start: function(event, ui) {\n",
|
|
" pass_mouse_events = false;\n",
|
|
" },\n",
|
|
" resize: function(event, ui) {\n",
|
|
" fig.request_resize(ui.size.width, ui.size.height);\n",
|
|
" },\n",
|
|
" stop: function(event, ui) {\n",
|
|
" pass_mouse_events = true;\n",
|
|
" fig.request_resize(ui.size.width, ui.size.height);\n",
|
|
" },\n",
|
|
" });\n",
|
|
"\n",
|
|
" function mouse_event_fn(event) {\n",
|
|
" if (pass_mouse_events)\n",
|
|
" return fig.mouse_event(event, event['data']);\n",
|
|
" }\n",
|
|
"\n",
|
|
" rubberband.mousedown('button_press', mouse_event_fn);\n",
|
|
" rubberband.mouseup('button_release', mouse_event_fn);\n",
|
|
" // Throttle sequential mouse events to 1 every 20ms.\n",
|
|
" rubberband.mousemove('motion_notify', mouse_event_fn);\n",
|
|
"\n",
|
|
" rubberband.mouseenter('figure_enter', mouse_event_fn);\n",
|
|
" rubberband.mouseleave('figure_leave', mouse_event_fn);\n",
|
|
"\n",
|
|
" canvas_div.on(\"wheel\", function (event) {\n",
|
|
" event = event.originalEvent;\n",
|
|
" event['data'] = 'scroll'\n",
|
|
" if (event.deltaY < 0) {\n",
|
|
" event.step = 1;\n",
|
|
" } else {\n",
|
|
" event.step = -1;\n",
|
|
" }\n",
|
|
" mouse_event_fn(event);\n",
|
|
" });\n",
|
|
"\n",
|
|
" canvas_div.append(canvas);\n",
|
|
" canvas_div.append(rubberband);\n",
|
|
"\n",
|
|
" this.rubberband = rubberband;\n",
|
|
" this.rubberband_canvas = rubberband[0];\n",
|
|
" this.rubberband_context = rubberband[0].getContext(\"2d\");\n",
|
|
" this.rubberband_context.strokeStyle = \"#000000\";\n",
|
|
"\n",
|
|
" this._resize_canvas = function(width, height) {\n",
|
|
" // Keep the size of the canvas, canvas container, and rubber band\n",
|
|
" // canvas in synch.\n",
|
|
" canvas_div.css('width', width)\n",
|
|
" canvas_div.css('height', height)\n",
|
|
"\n",
|
|
" canvas.attr('width', width * mpl.ratio);\n",
|
|
" canvas.attr('height', height * mpl.ratio);\n",
|
|
" canvas.attr('style', 'width: ' + width + 'px; height: ' + height + 'px;');\n",
|
|
"\n",
|
|
" rubberband.attr('width', width);\n",
|
|
" rubberband.attr('height', height);\n",
|
|
" }\n",
|
|
"\n",
|
|
" // Set the figure to an initial 600x600px, this will subsequently be updated\n",
|
|
" // upon first draw.\n",
|
|
" this._resize_canvas(600, 600);\n",
|
|
"\n",
|
|
" // Disable right mouse context menu.\n",
|
|
" $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n",
|
|
" return false;\n",
|
|
" });\n",
|
|
"\n",
|
|
" function set_focus () {\n",
|
|
" canvas.focus();\n",
|
|
" canvas_div.focus();\n",
|
|
" }\n",
|
|
"\n",
|
|
" window.setTimeout(set_focus, 100);\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype._init_toolbar = function() {\n",
|
|
" var fig = this;\n",
|
|
"\n",
|
|
" var nav_element = $('<div/>')\n",
|
|
" nav_element.attr('style', 'width: 100%');\n",
|
|
" this.root.append(nav_element);\n",
|
|
"\n",
|
|
" // Define a callback function for later on.\n",
|
|
" function toolbar_event(event) {\n",
|
|
" return fig.toolbar_button_onclick(event['data']);\n",
|
|
" }\n",
|
|
" function toolbar_mouse_event(event) {\n",
|
|
" return fig.toolbar_button_onmouseover(event['data']);\n",
|
|
" }\n",
|
|
"\n",
|
|
" for(var toolbar_ind in mpl.toolbar_items) {\n",
|
|
" var name = mpl.toolbar_items[toolbar_ind][0];\n",
|
|
" var tooltip = mpl.toolbar_items[toolbar_ind][1];\n",
|
|
" var image = mpl.toolbar_items[toolbar_ind][2];\n",
|
|
" var method_name = mpl.toolbar_items[toolbar_ind][3];\n",
|
|
"\n",
|
|
" if (!name) {\n",
|
|
" // put a spacer in here.\n",
|
|
" continue;\n",
|
|
" }\n",
|
|
" var button = $('<button/>');\n",
|
|
" button.addClass('ui-button ui-widget ui-state-default ui-corner-all ' +\n",
|
|
" 'ui-button-icon-only');\n",
|
|
" button.attr('role', 'button');\n",
|
|
" button.attr('aria-disabled', 'false');\n",
|
|
" button.click(method_name, toolbar_event);\n",
|
|
" button.mouseover(tooltip, toolbar_mouse_event);\n",
|
|
"\n",
|
|
" var icon_img = $('<span/>');\n",
|
|
" icon_img.addClass('ui-button-icon-primary ui-icon');\n",
|
|
" icon_img.addClass(image);\n",
|
|
" icon_img.addClass('ui-corner-all');\n",
|
|
"\n",
|
|
" var tooltip_span = $('<span/>');\n",
|
|
" tooltip_span.addClass('ui-button-text');\n",
|
|
" tooltip_span.html(tooltip);\n",
|
|
"\n",
|
|
" button.append(icon_img);\n",
|
|
" button.append(tooltip_span);\n",
|
|
"\n",
|
|
" nav_element.append(button);\n",
|
|
" }\n",
|
|
"\n",
|
|
" var fmt_picker_span = $('<span/>');\n",
|
|
"\n",
|
|
" var fmt_picker = $('<select/>');\n",
|
|
" fmt_picker.addClass('mpl-toolbar-option ui-widget ui-widget-content');\n",
|
|
" fmt_picker_span.append(fmt_picker);\n",
|
|
" nav_element.append(fmt_picker_span);\n",
|
|
" this.format_dropdown = fmt_picker[0];\n",
|
|
"\n",
|
|
" for (var ind in mpl.extensions) {\n",
|
|
" var fmt = mpl.extensions[ind];\n",
|
|
" var option = $(\n",
|
|
" '<option/>', {selected: fmt === mpl.default_extension}).html(fmt);\n",
|
|
" fmt_picker.append(option)\n",
|
|
" }\n",
|
|
"\n",
|
|
" // Add hover states to the ui-buttons\n",
|
|
" $( \".ui-button\" ).hover(\n",
|
|
" function() { $(this).addClass(\"ui-state-hover\");},\n",
|
|
" function() { $(this).removeClass(\"ui-state-hover\");}\n",
|
|
" );\n",
|
|
"\n",
|
|
" var status_bar = $('<span class=\"mpl-message\"/>');\n",
|
|
" nav_element.append(status_bar);\n",
|
|
" this.message = status_bar[0];\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.request_resize = function(x_pixels, y_pixels) {\n",
|
|
" // Request matplotlib to resize the figure. Matplotlib will then trigger a resize in the client,\n",
|
|
" // which will in turn request a refresh of the image.\n",
|
|
" this.send_message('resize', {'width': x_pixels, 'height': y_pixels});\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.send_message = function(type, properties) {\n",
|
|
" properties['type'] = type;\n",
|
|
" properties['figure_id'] = this.id;\n",
|
|
" this.ws.send(JSON.stringify(properties));\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.send_draw_message = function() {\n",
|
|
" if (!this.waiting) {\n",
|
|
" this.waiting = true;\n",
|
|
" this.ws.send(JSON.stringify({type: \"draw\", figure_id: this.id}));\n",
|
|
" }\n",
|
|
"}\n",
|
|
"\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_save = function(fig, msg) {\n",
|
|
" var format_dropdown = fig.format_dropdown;\n",
|
|
" var format = format_dropdown.options[format_dropdown.selectedIndex].value;\n",
|
|
" fig.ondownload(fig, format);\n",
|
|
"}\n",
|
|
"\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_resize = function(fig, msg) {\n",
|
|
" var size = msg['size'];\n",
|
|
" if (size[0] != fig.canvas.width || size[1] != fig.canvas.height) {\n",
|
|
" fig._resize_canvas(size[0], size[1]);\n",
|
|
" fig.send_message(\"refresh\", {});\n",
|
|
" };\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_rubberband = function(fig, msg) {\n",
|
|
" var x0 = msg['x0'] / mpl.ratio;\n",
|
|
" var y0 = (fig.canvas.height - msg['y0']) / mpl.ratio;\n",
|
|
" var x1 = msg['x1'] / mpl.ratio;\n",
|
|
" var y1 = (fig.canvas.height - msg['y1']) / mpl.ratio;\n",
|
|
" x0 = Math.floor(x0) + 0.5;\n",
|
|
" y0 = Math.floor(y0) + 0.5;\n",
|
|
" x1 = Math.floor(x1) + 0.5;\n",
|
|
" y1 = Math.floor(y1) + 0.5;\n",
|
|
" var min_x = Math.min(x0, x1);\n",
|
|
" var min_y = Math.min(y0, y1);\n",
|
|
" var width = Math.abs(x1 - x0);\n",
|
|
" var height = Math.abs(y1 - y0);\n",
|
|
"\n",
|
|
" fig.rubberband_context.clearRect(\n",
|
|
" 0, 0, fig.canvas.width, fig.canvas.height);\n",
|
|
"\n",
|
|
" fig.rubberband_context.strokeRect(min_x, min_y, width, height);\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_figure_label = function(fig, msg) {\n",
|
|
" // Updates the figure title.\n",
|
|
" fig.header.textContent = msg['label'];\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_cursor = function(fig, msg) {\n",
|
|
" var cursor = msg['cursor'];\n",
|
|
" switch(cursor)\n",
|
|
" {\n",
|
|
" case 0:\n",
|
|
" cursor = 'pointer';\n",
|
|
" break;\n",
|
|
" case 1:\n",
|
|
" cursor = 'default';\n",
|
|
" break;\n",
|
|
" case 2:\n",
|
|
" cursor = 'crosshair';\n",
|
|
" break;\n",
|
|
" case 3:\n",
|
|
" cursor = 'move';\n",
|
|
" break;\n",
|
|
" }\n",
|
|
" fig.rubberband_canvas.style.cursor = cursor;\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_message = function(fig, msg) {\n",
|
|
" fig.message.textContent = msg['message'];\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_draw = function(fig, msg) {\n",
|
|
" // Request the server to send over a new figure.\n",
|
|
" fig.send_draw_message();\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_image_mode = function(fig, msg) {\n",
|
|
" fig.image_mode = msg['mode'];\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.updated_canvas_event = function() {\n",
|
|
" // Called whenever the canvas gets updated.\n",
|
|
" this.send_message(\"ack\", {});\n",
|
|
"}\n",
|
|
"\n",
|
|
"// A function to construct a web socket function for onmessage handling.\n",
|
|
"// Called in the figure constructor.\n",
|
|
"mpl.figure.prototype._make_on_message_function = function(fig) {\n",
|
|
" return function socket_on_message(evt) {\n",
|
|
" if (evt.data instanceof Blob) {\n",
|
|
" /* FIXME: We get \"Resource interpreted as Image but\n",
|
|
" * transferred with MIME type text/plain:\" errors on\n",
|
|
" * Chrome. But how to set the MIME type? It doesn't seem\n",
|
|
" * to be part of the websocket stream */\n",
|
|
" evt.data.type = \"image/png\";\n",
|
|
"\n",
|
|
" /* Free the memory for the previous frames */\n",
|
|
" if (fig.imageObj.src) {\n",
|
|
" (window.URL || window.webkitURL).revokeObjectURL(\n",
|
|
" fig.imageObj.src);\n",
|
|
" }\n",
|
|
"\n",
|
|
" fig.imageObj.src = (window.URL || window.webkitURL).createObjectURL(\n",
|
|
" evt.data);\n",
|
|
" fig.updated_canvas_event();\n",
|
|
" fig.waiting = false;\n",
|
|
" return;\n",
|
|
" }\n",
|
|
" else if (typeof evt.data === 'string' && evt.data.slice(0, 21) == \"data:image/png;base64\") {\n",
|
|
" fig.imageObj.src = evt.data;\n",
|
|
" fig.updated_canvas_event();\n",
|
|
" fig.waiting = false;\n",
|
|
" return;\n",
|
|
" }\n",
|
|
"\n",
|
|
" var msg = JSON.parse(evt.data);\n",
|
|
" var msg_type = msg['type'];\n",
|
|
"\n",
|
|
" // Call the \"handle_{type}\" callback, which takes\n",
|
|
" // the figure and JSON message as its only arguments.\n",
|
|
" try {\n",
|
|
" var callback = fig[\"handle_\" + msg_type];\n",
|
|
" } catch (e) {\n",
|
|
" console.log(\"No handler for the '\" + msg_type + \"' message type: \", msg);\n",
|
|
" return;\n",
|
|
" }\n",
|
|
"\n",
|
|
" if (callback) {\n",
|
|
" try {\n",
|
|
" // console.log(\"Handling '\" + msg_type + \"' message: \", msg);\n",
|
|
" callback(fig, msg);\n",
|
|
" } catch (e) {\n",
|
|
" console.log(\"Exception inside the 'handler_\" + msg_type + \"' callback:\", e, e.stack, msg);\n",
|
|
" }\n",
|
|
" }\n",
|
|
" };\n",
|
|
"}\n",
|
|
"\n",
|
|
"// from http://stackoverflow.com/questions/1114465/getting-mouse-location-in-canvas\n",
|
|
"mpl.findpos = function(e) {\n",
|
|
" //this section is from http://www.quirksmode.org/js/events_properties.html\n",
|
|
" var targ;\n",
|
|
" if (!e)\n",
|
|
" e = window.event;\n",
|
|
" if (e.target)\n",
|
|
" targ = e.target;\n",
|
|
" else if (e.srcElement)\n",
|
|
" targ = e.srcElement;\n",
|
|
" if (targ.nodeType == 3) // defeat Safari bug\n",
|
|
" targ = targ.parentNode;\n",
|
|
"\n",
|
|
" // jQuery normalizes the pageX and pageY\n",
|
|
" // pageX,Y are the mouse positions relative to the document\n",
|
|
" // offset() returns the position of the element relative to the document\n",
|
|
" var x = e.pageX - $(targ).offset().left;\n",
|
|
" var y = e.pageY - $(targ).offset().top;\n",
|
|
"\n",
|
|
" return {\"x\": x, \"y\": y};\n",
|
|
"};\n",
|
|
"\n",
|
|
"/*\n",
|
|
" * return a copy of an object with only non-object keys\n",
|
|
" * we need this to avoid circular references\n",
|
|
" * http://stackoverflow.com/a/24161582/3208463\n",
|
|
" */\n",
|
|
"function simpleKeys (original) {\n",
|
|
" return Object.keys(original).reduce(function (obj, key) {\n",
|
|
" if (typeof original[key] !== 'object')\n",
|
|
" obj[key] = original[key]\n",
|
|
" return obj;\n",
|
|
" }, {});\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.mouse_event = function(event, name) {\n",
|
|
" var canvas_pos = mpl.findpos(event)\n",
|
|
"\n",
|
|
" if (name === 'button_press')\n",
|
|
" {\n",
|
|
" this.canvas.focus();\n",
|
|
" this.canvas_div.focus();\n",
|
|
" }\n",
|
|
"\n",
|
|
" var x = canvas_pos.x * mpl.ratio;\n",
|
|
" var y = canvas_pos.y * mpl.ratio;\n",
|
|
"\n",
|
|
" this.send_message(name, {x: x, y: y, button: event.button,\n",
|
|
" step: event.step,\n",
|
|
" guiEvent: simpleKeys(event)});\n",
|
|
"\n",
|
|
" /* This prevents the web browser from automatically changing to\n",
|
|
" * the text insertion cursor when the button is pressed. We want\n",
|
|
" * to control all of the cursor setting manually through the\n",
|
|
" * 'cursor' event from matplotlib */\n",
|
|
" event.preventDefault();\n",
|
|
" return false;\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype._key_event_extra = function(event, name) {\n",
|
|
" // Handle any extra behaviour associated with a key event\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.key_event = function(event, name) {\n",
|
|
"\n",
|
|
" // Prevent repeat events\n",
|
|
" if (name == 'key_press')\n",
|
|
" {\n",
|
|
" if (event.which === this._key)\n",
|
|
" return;\n",
|
|
" else\n",
|
|
" this._key = event.which;\n",
|
|
" }\n",
|
|
" if (name == 'key_release')\n",
|
|
" this._key = null;\n",
|
|
"\n",
|
|
" var value = '';\n",
|
|
" if (event.ctrlKey && event.which != 17)\n",
|
|
" value += \"ctrl+\";\n",
|
|
" if (event.altKey && event.which != 18)\n",
|
|
" value += \"alt+\";\n",
|
|
" if (event.shiftKey && event.which != 16)\n",
|
|
" value += \"shift+\";\n",
|
|
"\n",
|
|
" value += 'k';\n",
|
|
" value += event.which.toString();\n",
|
|
"\n",
|
|
" this._key_event_extra(event, name);\n",
|
|
"\n",
|
|
" this.send_message(name, {key: value,\n",
|
|
" guiEvent: simpleKeys(event)});\n",
|
|
" return false;\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.toolbar_button_onclick = function(name) {\n",
|
|
" if (name == 'download') {\n",
|
|
" this.handle_save(this, null);\n",
|
|
" } else {\n",
|
|
" this.send_message(\"toolbar_button\", {name: name});\n",
|
|
" }\n",
|
|
"};\n",
|
|
"\n",
|
|
"mpl.figure.prototype.toolbar_button_onmouseover = function(tooltip) {\n",
|
|
" this.message.textContent = tooltip;\n",
|
|
"};\n",
|
|
"mpl.toolbar_items = [[\"Home\", \"Reset original view\", \"fa fa-home icon-home\", \"home\"], [\"Back\", \"Back to previous view\", \"fa fa-arrow-left icon-arrow-left\", \"back\"], [\"Forward\", \"Forward to next view\", \"fa fa-arrow-right icon-arrow-right\", \"forward\"], [\"\", \"\", \"\", \"\"], [\"Pan\", \"Pan axes with left mouse, zoom with right\", \"fa fa-arrows icon-move\", \"pan\"], [\"Zoom\", \"Zoom to rectangle\", \"fa fa-square-o icon-check-empty\", \"zoom\"], [\"\", \"\", \"\", \"\"], [\"Download\", \"Download plot\", \"fa fa-floppy-o icon-save\", \"download\"]];\n",
|
|
"\n",
|
|
"mpl.extensions = [\"eps\", \"jpeg\", \"pdf\", \"png\", \"ps\", \"raw\", \"svg\", \"tif\"];\n",
|
|
"\n",
|
|
"mpl.default_extension = \"png\";var comm_websocket_adapter = function(comm) {\n",
|
|
" // Create a \"websocket\"-like object which calls the given IPython comm\n",
|
|
" // object with the appropriate methods. Currently this is a non binary\n",
|
|
" // socket, so there is still some room for performance tuning.\n",
|
|
" var ws = {};\n",
|
|
"\n",
|
|
" ws.close = function() {\n",
|
|
" comm.close()\n",
|
|
" };\n",
|
|
" ws.send = function(m) {\n",
|
|
" //console.log('sending', m);\n",
|
|
" comm.send(m);\n",
|
|
" };\n",
|
|
" // Register the callback with on_msg.\n",
|
|
" comm.on_msg(function(msg) {\n",
|
|
" //console.log('receiving', msg['content']['data'], msg);\n",
|
|
" // Pass the mpl event to the overridden (by mpl) onmessage function.\n",
|
|
" ws.onmessage(msg['content']['data'])\n",
|
|
" });\n",
|
|
" return ws;\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.mpl_figure_comm = function(comm, msg) {\n",
|
|
" // This is the function which gets called when the mpl process\n",
|
|
" // starts-up an IPython Comm through the \"matplotlib\" channel.\n",
|
|
"\n",
|
|
" var id = msg.content.data.id;\n",
|
|
" // Get hold of the div created by the display call when the Comm\n",
|
|
" // socket was opened in Python.\n",
|
|
" var element = $(\"#\" + id);\n",
|
|
" var ws_proxy = comm_websocket_adapter(comm)\n",
|
|
"\n",
|
|
" function ondownload(figure, format) {\n",
|
|
" window.open(figure.imageObj.src);\n",
|
|
" }\n",
|
|
"\n",
|
|
" var fig = new mpl.figure(id, ws_proxy,\n",
|
|
" ondownload,\n",
|
|
" element.get(0));\n",
|
|
"\n",
|
|
" // Call onopen now - mpl needs it, as it is assuming we've passed it a real\n",
|
|
" // web socket which is closed, not our websocket->open comm proxy.\n",
|
|
" ws_proxy.onopen();\n",
|
|
"\n",
|
|
" fig.parent_element = element.get(0);\n",
|
|
" fig.cell_info = mpl.find_output_cell(\"<div id='\" + id + \"'></div>\");\n",
|
|
" if (!fig.cell_info) {\n",
|
|
" console.error(\"Failed to find cell for figure\", id, fig);\n",
|
|
" return;\n",
|
|
" }\n",
|
|
"\n",
|
|
" var output_index = fig.cell_info[2]\n",
|
|
" var cell = fig.cell_info[0];\n",
|
|
"\n",
|
|
"};\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_close = function(fig, msg) {\n",
|
|
" var width = fig.canvas.width/mpl.ratio\n",
|
|
" fig.root.unbind('remove')\n",
|
|
"\n",
|
|
" // Update the output cell to use the data from the current canvas.\n",
|
|
" fig.push_to_output();\n",
|
|
" var dataURL = fig.canvas.toDataURL();\n",
|
|
" // Re-enable the keyboard manager in IPython - without this line, in FF,\n",
|
|
" // the notebook keyboard shortcuts fail.\n",
|
|
" IPython.keyboard_manager.enable()\n",
|
|
" $(fig.parent_element).html('<img src=\"' + dataURL + '\" width=\"' + width + '\">');\n",
|
|
" fig.close_ws(fig, msg);\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.close_ws = function(fig, msg){\n",
|
|
" fig.send_message('closing', msg);\n",
|
|
" // fig.ws.close()\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.push_to_output = function(remove_interactive) {\n",
|
|
" // Turn the data on the canvas into data in the output cell.\n",
|
|
" var width = this.canvas.width/mpl.ratio\n",
|
|
" var dataURL = this.canvas.toDataURL();\n",
|
|
" this.cell_info[1]['text/html'] = '<img src=\"' + dataURL + '\" width=\"' + width + '\">';\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.updated_canvas_event = function() {\n",
|
|
" // Tell IPython that the notebook contents must change.\n",
|
|
" IPython.notebook.set_dirty(true);\n",
|
|
" this.send_message(\"ack\", {});\n",
|
|
" var fig = this;\n",
|
|
" // Wait a second, then push the new image to the DOM so\n",
|
|
" // that it is saved nicely (might be nice to debounce this).\n",
|
|
" setTimeout(function () { fig.push_to_output() }, 1000);\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype._init_toolbar = function() {\n",
|
|
" var fig = this;\n",
|
|
"\n",
|
|
" var nav_element = $('<div/>')\n",
|
|
" nav_element.attr('style', 'width: 100%');\n",
|
|
" this.root.append(nav_element);\n",
|
|
"\n",
|
|
" // Define a callback function for later on.\n",
|
|
" function toolbar_event(event) {\n",
|
|
" return fig.toolbar_button_onclick(event['data']);\n",
|
|
" }\n",
|
|
" function toolbar_mouse_event(event) {\n",
|
|
" return fig.toolbar_button_onmouseover(event['data']);\n",
|
|
" }\n",
|
|
"\n",
|
|
" for(var toolbar_ind in mpl.toolbar_items){\n",
|
|
" var name = mpl.toolbar_items[toolbar_ind][0];\n",
|
|
" var tooltip = mpl.toolbar_items[toolbar_ind][1];\n",
|
|
" var image = mpl.toolbar_items[toolbar_ind][2];\n",
|
|
" var method_name = mpl.toolbar_items[toolbar_ind][3];\n",
|
|
"\n",
|
|
" if (!name) { continue; };\n",
|
|
"\n",
|
|
" var button = $('<button class=\"btn btn-default\" href=\"#\" title=\"' + name + '\"><i class=\"fa ' + image + ' fa-lg\"></i></button>');\n",
|
|
" button.click(method_name, toolbar_event);\n",
|
|
" button.mouseover(tooltip, toolbar_mouse_event);\n",
|
|
" nav_element.append(button);\n",
|
|
" }\n",
|
|
"\n",
|
|
" // Add the status bar.\n",
|
|
" var status_bar = $('<span class=\"mpl-message\" style=\"text-align:right; float: right;\"/>');\n",
|
|
" nav_element.append(status_bar);\n",
|
|
" this.message = status_bar[0];\n",
|
|
"\n",
|
|
" // Add the close button to the window.\n",
|
|
" var buttongrp = $('<div class=\"btn-group inline pull-right\"></div>');\n",
|
|
" var button = $('<button class=\"btn btn-mini btn-primary\" href=\"#\" title=\"Stop Interaction\"><i class=\"fa fa-power-off icon-remove icon-large\"></i></button>');\n",
|
|
" button.click(function (evt) { fig.handle_close(fig, {}); } );\n",
|
|
" button.mouseover('Stop Interaction', toolbar_mouse_event);\n",
|
|
" buttongrp.append(button);\n",
|
|
" var titlebar = this.root.find($('.ui-dialog-titlebar'));\n",
|
|
" titlebar.prepend(buttongrp);\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype._root_extra_style = function(el){\n",
|
|
" var fig = this\n",
|
|
" el.on(\"remove\", function(){\n",
|
|
"\tfig.close_ws(fig, {});\n",
|
|
" });\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype._canvas_extra_style = function(el){\n",
|
|
" // this is important to make the div 'focusable\n",
|
|
" el.attr('tabindex', 0)\n",
|
|
" // reach out to IPython and tell the keyboard manager to turn it's self\n",
|
|
" // off when our div gets focus\n",
|
|
"\n",
|
|
" // location in version 3\n",
|
|
" if (IPython.notebook.keyboard_manager) {\n",
|
|
" IPython.notebook.keyboard_manager.register_events(el);\n",
|
|
" }\n",
|
|
" else {\n",
|
|
" // location in version 2\n",
|
|
" IPython.keyboard_manager.register_events(el);\n",
|
|
" }\n",
|
|
"\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype._key_event_extra = function(event, name) {\n",
|
|
" var manager = IPython.notebook.keyboard_manager;\n",
|
|
" if (!manager)\n",
|
|
" manager = IPython.keyboard_manager;\n",
|
|
"\n",
|
|
" // Check for shift+enter\n",
|
|
" if (event.shiftKey && event.which == 13) {\n",
|
|
" this.canvas_div.blur();\n",
|
|
" event.shiftKey = false;\n",
|
|
" // Send a \"J\" for go to next cell\n",
|
|
" event.which = 74;\n",
|
|
" event.keyCode = 74;\n",
|
|
" manager.command_mode();\n",
|
|
" manager.handle_keydown(event);\n",
|
|
" }\n",
|
|
"}\n",
|
|
"\n",
|
|
"mpl.figure.prototype.handle_save = function(fig, msg) {\n",
|
|
" fig.ondownload(fig, null);\n",
|
|
"}\n",
|
|
"\n",
|
|
"\n",
|
|
"mpl.find_output_cell = function(html_output) {\n",
|
|
" // Return the cell and output element which can be found *uniquely* in the notebook.\n",
|
|
" // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n",
|
|
" // IPython event is triggered only after the cells have been serialised, which for\n",
|
|
" // our purposes (turning an active figure into a static one), is too late.\n",
|
|
" var cells = IPython.notebook.get_cells();\n",
|
|
" var ncells = cells.length;\n",
|
|
" for (var i=0; i<ncells; i++) {\n",
|
|
" var cell = cells[i];\n",
|
|
" if (cell.cell_type === 'code'){\n",
|
|
" for (var j=0; j<cell.output_area.outputs.length; j++) {\n",
|
|
" var data = cell.output_area.outputs[j];\n",
|
|
" if (data.data) {\n",
|
|
" // IPython >= 3 moved mimebundle to data attribute of output\n",
|
|
" data = data.data;\n",
|
|
" }\n",
|
|
" if (data['text/html'] == html_output) {\n",
|
|
" return [cell, data, j];\n",
|
|
" }\n",
|
|
" }\n",
|
|
" }\n",
|
|
" }\n",
|
|
"}\n",
|
|
"\n",
|
|
"// Register the function which deals with the matplotlib target/channel.\n",
|
|
"// The kernel may be null if the page has been refreshed.\n",
|
|
"if (IPython.notebook.kernel != null) {\n",
|
|
" IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n",
|
|
"}\n"
|
|
],
|
|
"text/plain": [
|
|
"<IPython.core.display.Javascript object>"
|
|
]
|
|
},
|
|
"metadata": {},
|
|
"output_type": "display_data"
|
|
},
|
|
{
|
|
"data": {
|
|
"text/html": [
|
|
"<img src=\"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAoAAAAHgCAYAAAA10dzkAAAgAElEQVR4nOzdeXyU1dk+8FEKWCtCW7UuPzliBbRacY1iW2mrb6XVYmvrglZRq9Yq1de32MMeAQFFEBVkEVmq7CoonGxsARIgCYGQBBL2kISQsAQSkpB15vr9MZmBQCAJzzxznuX6fj7zR4bJPHftDeeaZ87iARERERG5ikd3AUREREQUXgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAEhERETkMgyARERERC7DAGiA1+tFfn4+SkpKUFpaygcffPDBBx982OBRUlKC/Px8eL1e3VFCGwZAA/Lz8+HxePjggw8++OCDDxs+8vPzdUcJbRgADSgpKQk2kO5PM3zwwQcffPDBR/MegRs4JSUluqOENgyABpSWlsLj8aC0tFR3KURERNRMHL8ZAA1hAxEREdkPx28bBcCJEydCCIG2bdsiIiICycnJZ33t1q1b8dhjj0EIAY/Hg/Hjx5/xmsjIyDPmAnTt2rVFNbGBiIiI7Ifjt00C4Pz589GmTRvMmDED27Ztw8svv4wOHTrg4MGDjb4+JSUF/fr1w7x583DllVeeNQDefPPNKCwsDD4OHz7corrYQERERPbD8dsmATAiIgKvv/568Gev14urr74ao0ePbvJ3hRBnDYDdunUzVBcbiIiIyH44ftsgAFZXV6NVq1ZYvHhxg+efe+459OrVq8nfP1cAvPjii3HVVVehU6dOePrpp5Gbm3vO96qqqmp0FZGbG4iIiMhuGABtEAALCgrg8Xiwfv36Bs+//fbbiIiIaPL3zxYAo6OjsXDhQqSnpyM2Nhbdu3dHx44dcfz48bO+V2PzBt3eQERERHbDAOjiAHi6Y8eO4dJLL8Xnn39+1tfwDiAREZH9MQDaIACa9RVwY+666y7079+/2bWxgYiIiOyH47cNAiDgXwTSt2/f4M9erxfXXHONoUUgpysrK8MPf/hDfPzxx82uiw1ERERkPxy/bRIA58+fj7Zt22LWrFnIysrCK6+8gg4dOqCoqAgA8Oyzzza4c1ddXY20tDSkpaXhqquuQr9+/ZCWloZdu3YFX/Pvf/8bq1evRk5ODtatW4cHH3wQl112GQ4dOtTsuthARERE9sPx2yYBEAAmTJiAjh07ok2bNoiIiEBSUlLwz3r06IE+ffoEf87JyWl0sUaPHj2Cr3nyySdx1VVXoU2bNrjmmmvw5JNPYvfu3S2qiQ1ERERkPxy/bRQArYgNREREZD8cvxkADTGrgaIzDuCNeZvxVWp+SN+XiIjILr5Kzccb8zYjOuNAyN+bAZAB0BCzGmj88h0QUqH/N+khfV8iIiK7kF+nQ0iFj5bvDPl7MwAyABpiVgPNScqFkAovzkwJ6fsSERHZxQszUyCkwtzkc5/SdT4YABkADTGrgZZvK4KQCo98khDS9yUiIrKLhz9ZCyEVVmQVhfy9GQAZAA0xq4HS849BSIWIkctD+r5ERER2cfe7yyGkQkZ+ScjfmwGQAdAQsxqosKQSQipcPyAKdV5fSN+biIjI6uq8PnTqryCkQlFpZcjfnwGQAdAQsxqots6L6+ob/9DxqpC+NxERkdUdPO6/EXJdf4XaOm/I358BkAHQEDMb6M4RyyCkwtaC0N/6JiIisrLM/SUQUuHOEeZMhWIAZAA0xMwG6vmRf/Lrqu0HQ/7eREREVrYq+yCEVPj9R2tNeX8GQAZAQ8xsoOemJ0NIhQUpeSF/byIiIiubn+LfDq3PjGRT3p8BkAHQEDMbqN/CLRBSYcLK0G+ASUREZGWfrNgJIRXe/mqLKe/PAMgAaIiZDTQmNhtCKgz5NjPk701ERGRlgxdnQkiFD2K3m/L+DIAMgIaY2UCz1uVASIV/fJEa8vcmIiKysle+2AghFf67PseU92cAZAA0xMwGis44ACEV/vxpYsjfm4iIyMr+9GkihFSIyTxgyvszADIAGmJmA6XuK4aQCr94b2XI35uIiMjK7hu9EkIqpO47asr7MwAyABpiZgPlFVdASIXOg6Lh8/E0ECIicgefz4fOA6MhpEJecYUp12AAZAA0xMwGqqypg5D+00BKKmpC/v5ERERWdKyiOjj+VdbUmXINBkAGQEPMbqCfR8ZCSIWdRcdNeX8iIiKr2VF0HEIq3PpOnGnXYABkADTE7AZ6YNxqCKmQuOuwKe9PRERkNQk7D0NIhQfHrTbtGgyADICGmN1AvT/bACEVFm3ON+X9iYiIrOabTfkQUuHpaRtMuwYDIAOgIWY30JvzNkNIhSmrd5vy/kRERFYzefVuCKnwv/PTTLsGAyADoCFmN9DIqCwIqTB86TZT3p+IiMhqhi3ZBiEVRkVlmXYNBkAGQEPMbqBpa/dASIW+czeb8v5ERERW8/qcTRBSYdraPaZdgwGQAdAQsxvo27T9EFLhiSnrTXl/IiIiq3l8ynoIqfDdlgLTrsEAyABoiNkNtH73EQip8JsP4k15fyIiIqv59QfxEFJhw54jpl2DAZAB0BCzG2j3oTIIqXDz0FhT3p+IiMhqfjYkBkIq7DlUZto1GAAZAA0xu4GOV9YEd0OvqK415RpERERWUV5VGxz3yqrMG/cYABkADTG7gXw+H24c7P8klHO43JRrEBERWcXew+UQUuGmITGmXocBkAHQkHA00K/eXwUhFZL3Fpt2DSIiIitI2uOf+37/mFWmXocBkAHQkHA00F8mrYOQCkvTzVsNRUREZAVLthRASIW/Tl5n6nUYABkADQlHA/1zdiqEVJiesNe0axAREVnB5wl7IaTCa7M3mXodBkAGQEPC0UCR322FkAqjo7NNuwYREZEVjIr2n4AV+d1WU6/DAGijADhx4kQIIdC2bVtEREQgOTn5rK/dunUrHnvsMQgh4PF4MH78eMPv2ZhwNNDEVbsgpMJbC8w7E5GIiMgK3pqfBiEVPo3fZep1GABtEgDnz5+PNm3aYMaMGdi2bRtefvlldOjQAQcPHmz09SkpKejXrx/mzZuHK6+8stEA2NL3bEw4GmjhxjwIqfC3z5NMuwYREZEVPDMtCUIqfJWab+p1GABtEgAjIiLw+uuvB3/2er24+uqrMXr06CZ/VwjRaAA08p4B4Wig1TsOQUiF3324xrRrEBERWcH/fLgaQiqs2XHI1OswANogAFZXV6NVq1ZYvHhxg+efe+459OrVq8nfbywAnu97VlVVobS0NPjIz883vYGyDpRCSIXbhsWZdg0iIiIr6DYsDkIqZBeaG8wYAG0QAAsKCuDxeLB+/foGz7/99tuIiIho8vcbC4Dn+56RkZHweDxnPMxsoCNlVcFd0atrvaZdh4iISKeq2rrgeFdcXm3qtRgAGQBb9J467gB6vT78dEAUhFTYf+yEadchIiLSKf9oBYRUuGFgFHw+n6nXYgC0QQC00lfApwtXA907agWEVNice9TU6xAREemyKfcohFToPmqF6ddiALRBAAT8Czb69u0b/Nnr9eKaa64xvAjkfN8zIFwN1GtCAoRUiN1aaOp1iIiIdInJLISQCr0mJpp+LQZAmwTA+fPno23btpg1axaysrLwyiuvoEOHDigqKgIAPPvss+jfv3/w9dXV1UhLS0NaWhquuuoq9OvXD2lpadi1a1ez37M5wtVAf5+1EUIqfLFhn6nXISIi0uWL9TkQUuGl/240/VoMgDYJgAAwYcIEdOzYEW3atEFERASSkk7ui9ejRw/06dMn+HNOTk6jizV69OjR7PdsjnA10IBFGRBSYVzcdlOvQ0REpMvYuO0QUmHgogzTr8UAaKMAaEXhaqDxy3dASAX5dbqp1yEiItLlP1+lQ0iFj5bvNP1aDIAMgIaEq4HmJOVCSIUXZqaYeh0iIiJdnp+RDCEV5ibnmn4tBkAGQEPC1UDLtxVBSIWHP1lr6nWIiIh0+cPHayGkwoqs5s/FP18MgAyAhoSrgdLzj0FIhbvfXW7qdYiIiHS5693lEFIhI7/E9GsxADIAGhKuBiosqYSQCp36K9R5zd0ck4iIKNzqvD506u8/BaSotNL06zEAMgAaEq4Gqq3z4rr6vxgHj5v/F4OIiCicDpb6b3Rc11+hts78Y08ZABkADQlnA905wn9rPHO/+bfGiYiIwikjvwRCKtwVpqlODIAMgIaEs4Ee+cR/GsjybeZPjiUiIgqnZfWLHf84ISEs12MAZAA0JJwNxNNAiIjIqcJ5CgjAAAgwABoSzgYavDgTQip8EMvTQIiIyFnGxGZDSIUh32aG5XoMgAyAhoSzgSau2gUhFf5vwRbTr0VERBROby1Ig5AKn8bvCsv1GAAZAA0JZwN9nZoPIRWemday84qJiIis7ulpGyCkwjeb8sNyPQZABkBDwtlA63YdhpAKvx0bb/q1iIiIwuk3Y+MhpMK63YfDcj0GQAZAQ8LZQLsPlUFIhZ8NiTH9WkREROHi8/lw05AYCKmw51BZWK7JAMgAaEg4G6i8qhZC+jeDPl5ZY/r1iIiIwqG0siY4vlVU14bnmgyADIBGhLuBbomMhZAKO4uOh+V6REREZttRdBxCKvw8MjZs12QAZAA0JNwN9D8froaQCmt3HgrL9YiIiMy2ZschCKnwuw/XhO2aDIAMgIaEu4H+9nkShFRYsDEvLNcjIiIy24KUPAip8Oz05LBdkwGQAdCQcDfQ219tgZAKn6zYGZbrERERme3jFTshpMJ/vkoP2zUZABkADQl3A41btgNCKgxYlBGW6xEREZmt/zcZEFLhw2U7wnZNBkAGQEPC3UBzknIhpMKLM1PCcj0iIiKzvTAzBUIqzE3ODds1GQAZAA0JdwOtyj4IIRV+/9HasFyPiIjIbD0/WgshFVZtPxi2azIAMgAaEu4G2lZQCiEV7hi+LCzXIyIiMtvtw5dBSIWsA+ELYwyADICGhLuBjpZXBzfLrKypC8s1iYiIzFJZUxcc145VVIftugyADICGhLuBfD4fOg+KhpAKuUcqwnJNIiIis+w7Ug4hFboMiobP5wvbdRkAGQAN0dFA949ZBSEVkvcWh+2aREREZkjacwRCKvQYsyqs12UAZAA0REcDPTFlPYRU+DZtf9iuSUREZIZv0/ZDSIUnp64P63UZABkADdHRQG/M2wwhFaau2R22axIREZlhyurdEFLhzXmbw3pdBkAGQEN0NNCo6CwIqfDOkq1huyYREZEZIr/bCiEVRkdnh/W6DIAMgIboaKAZiXshpMI/Z6eG7ZpERERmePXLVAipMDNxb1ivywDIAGiIjgaKyTwAIRX+9Gli2K5JRERkhkcnJkJIhZjMwrBelwGQAdAQHQ20OfcohFToPmpF2K5JRERkhntHrYCQCml5x8J6XQZABkBDdDRQYUklhFS4fkAU6rzh2zOJiIgolGrrvLh+QBSEVCgqrQzrtRkAbRYAJ06cCCEE2rZti4iICCQnJ5/z9QsXLkTXrl3Rtm1b3HLLLYiKimrw53369IHH42nweOihh5pdj44Gqq3zolN/peUvDBERUajovKHBAGijADh//ny0adMGM2bMwLZt2/Dyyy+jQ4cOOHiw8cOj161bh1atWmHMmDHIysrC4MGD0bp1a2RmZgZf06dPH/Ts2ROFhYXBx9GjR5tdk64Gumek/5b5ljDfMiciIgqVtLxjEFLhXg1TmhgAbRQAIyIi8Prrrwd/9nq9uPrqqzF69OhGX//EE0/g4YcfbvDcPffcg3/84x/Bn/v06YNHH330vGvS1UC6Js0SERGFis5FjQyANgmA1dXVaNWqFRYvXtzg+eeeew69evVq9HeuvfZajB8/vsFzQ4cOxa233hr8uU+fPmjfvj0uv/xydOnSBa+++iqOHDly1jqqqqpQWloafOTn52tpoH984V82P2tdTlivS0REFCoz67c1e/XL8G9rxgBokwBYUFAAj8eD9esbHhXz9ttvIyIiotHfad26NebOndvguU8//RRXXHFF8Od58+bhu+++Q0ZGBhYvXoybbroJd999N+rq6hp9z8jIyDPmDOpoIF0bZxIREYWKzoMNGABdHgBPt2fPHng8HqxY0fh8BKvcAZxcf3TO/85PC+t1iYiIQuXN+qNNp6wO/9GmDIA2CYBmfQXcmMsuuwxTpkxpVl26GkjX4dlERESh8sSU9RBS4du0/WG/NgOgTQIg4F8E0rdv3+DPXq8X11xzzTkXgTzyyCMNnuvevXuDRSCny8/PxwUXXIDvvvuuWTXpaqANe45ASIVffxAf1usSERGFSo8xqyCkQtKes8+9NwsDoI0C4Pz589G2bVvMmjULWVlZeOWVV9ChQwcUFRUBAJ599ln0798/+Pp169bhe9/7HsaOHYvs7GxERkY22AamrKwM/fr1w4YNG5CTk4MVK1bgjjvuQOfOnVFVVdWsmnQ10L4j5RBSoevgaPh83AyaiIjsxefzocugaAipkHukIuzXZwC0UQAEgAkTJqBjx45o06YNIiIikJSUFPyzHj16oE+fPg1ev3DhQnTp0gVt2rTBzTff3GAj6BMnTuB3v/sdLr/8crRu3RpCCLz88svBQNkcuhqosqYOQvo3gy6pqAnrtYmIiIw6VlEdHMcqaxpfeGkmBkCbBUCr0dlAtw2Lg5AK2YXubV4iIrKnrAOlEFLh9uHLtFyfAZAB0BCdDdTzo7UQUmHV9sZPQiEiIrKqVdkHIaTC7z9aq+X6DIAMgIbobKDnZyRDSIV5yblhvzYREZERc5NzIaTCCzNTtFyfAZAB0BCdDdT/mwwIqfDhsh1hvzYREZER45btgJAKAxZlaLk+AyADoCE6G+ij5TshpIL8Oj3s1yYiIjLiP1+lQ0iFj1fs1HJ9BkAGQEN0NtCClDwIqfDs9OSwX5uIiMiIZ6f7pzEt2Jin5foMgAyAhuhsoISdhyGkwgPjVof92kREREb8dmw8hFRI3HVYy/UZABkADdHZQHsP+zeDvmlIDDeDJiIi2/D5fLhxcAyEVMg5XK6lBgZABkBDdDbQqZtBHy2vDvv1iYiIzkdxud5NoAEGQIAB0BDdDXTniOUQUiFzf4mW6xMREbVU5v4SCKlw17vLtdWge/y2AgZAA3Q3UK8JCRBSIXZroZbrExERtVRMZiGEVOg1MVFbDbrHbytgADRAdwP9c3YqhFSYkbhXy/WJiIhaanrCXgip8NrsTdpq0D1+WwEDoAG6G2jE0m0QUmHE0m1ark9ERNRSw+vHrneVvrFL9/htBQyABuhuoBmJ/k9R/5ydquX6RERELfXql/5vr2Zq/PZK9/htBQyABuhuoNit9fMoJiRouT4REVFL/bF+/nqcxvnrusdvK2AANEB3AwVWUt05Qt9KKiIiopa4c8QyCKmwtUDfDha6x28rYAA0QHcDHbXAXkpERETNdeoetscq9O1hq3v8tgIGQAN0N5DP58NNQ/y7qe/VtJs6ERFRc+05VAYhFX6m+RQr3eO3FTAAGmCFBnpg3GoIqZCwU895ikRERM21duchCKnwoOZz7K0wfuvGAGiAFRrouenJEFJhQUqethqIiIiaY35KLoRU6DMjWWsdVhi/dWMANMAKDdT/mwwIqTBu2Q5tNRARETXHuLjtEFJhwKIMrXVYYfzWjQHQACs00MRVuyCkwr8XbtFWAxERUXP834ItEFJh4qpdWuuwwvitGwOgAVZooEWb8yGkwlNTN2irgYiIqDmenLoeQios3rxfax1WGL91YwA0wAoNlLy3GEIq3D9mlbYaiIiImuNX76+CkAopOcVa67DC+K0bA6ABVmig/KMVEFKh88BoeL36ltQTERGdi9frww0DoyCkwv5jJ7TWYoXxWzcGQAOs0EC1dV506u/fVPNgaaW2OoiIiM6lqLQSQipcPyAKtXVerbVYYfzWjQHQAKs0UPdRKyCkwqbco1rrICIiOpvUfUchpMJ9o1fqLsUy47dODIAGWKWB/jp5HYRUWJpeoLUOIiKis1mypQBCKjw+eb3uUiwzfuvEAGiAVRrojXmbIaTClNW7tdZBRER0NpNX74aQCm/O26y7FMuM3zoxABpglQZ6PyYbQioM/TZTax1ERERnM+TbTAipMCY2W3cplhm/dWIANMAqDfTlhn0QUuHvs1K01kFERHQ2L85MgZAKs5P26S7FMuO3TgyABlilgVZtPwghFXp+tFZrHURERGfz0Pg1EFIhfvtB3aVYZvzWiQHQAKs00I6i4xBS4eeRsVrrICIiOptbImMhpMLOouO6S7HM+K2TrQLgxIkTIYRA27ZtERERgeTk5HO+fuHChejatSvatm2LW265BVFRUQ3+3OfzYciQIbjyyitx0UUX4YEHHsDOnTubXY9VGqisqhZC+vcCPF5Zo7UWIiKi05VW1gTHqfKqWt3lWGb81sk2AXD+/Plo06YNZsyYgW3btuHll19Ghw4dcPBg47eS161bh1atWmHMmDHIysrC4MGD0bp1a2Rmnlwo8d5776F9+/b49ttvkZ6ejl69eqFTp06orGzehspWaqBb34mDkArbC/V/siIiIjpVdmEphFToNixOdykArDV+62KbABgREYHXX389+LPX68XVV1+N0aNHN/r6J554Ag8//HCD5+655x784x//AOC/+3fllVfigw8+CP55SUkJ2rZti3nz5jWrJis10O8/WgshFVZl659bQUREdKqV2UUQUuEPH1tjrrqVxm9dbBEAq6ur0apVKyxevLjB88899xx69erV6O9ce+21GD9+fIPnhg4diltvvRUAsGfPHng8HqSlpTV4zf3334833nij0fesqqpCaWlp8JGfn2+ZBvr7rI0QUuGLDfpXVxHZTZ3Xh72HyxGTWYiJq3ZhbNz2Bo+Jq3YhJrMQew+Xo45nbhO12BfrcyCkwkv/3ai7FAAMgIBNAmBBQQE8Hg/Wr2+4e/jbb7+NiIiIRn+ndevWmDt3boPnPv30U1xxxRUA/F8RezweHDhwoMFrHn/8cTzxxBONvmdkZCQ8Hs8ZDys0UOR3WyGkwnsx+vdXIrI6n8+HtLxjGL50Gx75JAFdB0cH5yc19eg6OBqPfJKA4Uu3IS3vGHw+BkKipoyO9u9XG/ndVt2lAGAABBgAWxQArXwHcOoa/w7r/5qrf4d1IqvaUXQcY2Kz8av3V50R7LoM8ge7txakIfK7rQ0eby1IwyOfJKDLoDOD4q/eX4UxsdnYYYGVjURW1Xeu/8Sqz9bs0V0KAAZAwCYB0CpfAZ/OSg2k0g9ASIW/TFqnuxQiy8ncX4LnZyQ3CG43Do5B37mbsTS9oNlf7Qa+Kl6aXoC+czfjxsExDd7zhZkpyNxfEob/RUT28tgk/5n1URkHmn5xGFhp/NbFFgEQ8C8C6du3b/Bnr9eLa6655pyLQB555JEGz3Xv3v2MRSBjx44N/nlpaaltF4Fszj0KIRXuHbVCdylElrG98Dj+8UVqMKBdPyAKL/13I5ZsKUBFtfGtKCqqa7FkSwFe+u9GXD8gKnidV79M5R1BolPcM3IFhFRIyzumuxQA1hq/dbFNAJw/fz7atm2LWbNmISsrC6+88go6dOiAoqIiAMCzzz6L/v37B1+/bt06fO9738PYsWORnZ2NyMjIRreB6dChA7777jtkZGTg0Ucfte02MAePV0JIhU79FWrqvLrLIdKqpKIG/164Bdf19wey6/or/O/8NOQcLjftmjmHy/HmvM0NrvnvhVtQUsG9Ocndauq8wb8Xh45X6S4HgLXGb11sEwABYMKECejYsSPatGmDiIgIJCUlBf+sR48e6NOnT4PXL1y4EF26dEGbNm1w8803n3Uj6J/85Cdo27YtHnjgAezYsaPZ9VipgbxeHzoP9M9Pyiuu0F0OkTardxwK3m0QUuGfs8N7N25H0XG8+uXJu473jFyBNTsOhe36RFaTV1wBIRU6D4qG1yKr6K00futiqwBoNVZroB5j/BPbN+w5orsUorArr6rFwEUZweD16w/ikbrvqLZ6UvcV49cfxAfrGbQ4wxInIBCF2/rdR4J/J63CauO3DgyABlitgf72eRKEVFi4MU93KURhlZFf0mBlb+R3W3Giuk53WaiorsXQbzODdd0/ZhUy8rlIhNxlwcY8CKnwt8+Tmn5xmFht/NaBAdAAqzXQgPq7H+PitusuhShsYjIPBPfx6z5qBRJ3HdZd0hkSdx1G91ErgvsIxmRaYyUkUTiMjdsOIRUGLsrQXUqQ1cZvHRgADbBaA01e7d8L8M153AuQnM/n8+HT+F3Bu2vPTU9GyQnrLrgoOVGD56af3IpmUvxubiJNrvDGPP8egFNW79ZdSpDVxm8dGAANsFoDRWX49wL886eJukshMlV1rRf9Fm4Jhqmh32ai1gar32vrvA2+En77qy2orrV+3URG/OnTRAipEG2RPQAB643fOjAAGmC1BsrcXwIhFe4csVx3KUSmKa+qxVNTNwS3PZq1Lkd3SS02M3EvOtVvi/HU1A1cHEKOdueIZRBSWWqTdKuN3zowABpgtQYqOVETvLMQik1uiXzEM6wAACAASURBVKymvKoWj09eDyEVbh4ai1XbD+ou6byt2n4QNw+NhZAKj09ezxBIjlReVRscl0orrTNFw2rjtw4MgAZYsYFufScOQipsL+QpBOQs5VW1eHyKP/zdEhmLLRY5UcCItLxjuCUQAqes5wc3cpzswlIIqdBtWJzuUhqw4vgdbgyABlixgf44IQFCKizbVqS7FKKQqaiuxROB8Dc01jLHSYXC5tyjwRD45FSGQHKWuK2FEFKh14QE3aU0YMXxO9wYAA2wYgO9NmcThFT4PGGv7lKIQqKiuhZPTj35te/mXH2bO5tlU+7R4NfBT03dYIk9DIlCYdraPRBS4fU5m3SX0oAVx+9wYwA0wIoN9F5MdnAjXCK7q63zos+M5GD403myh9lS9x3Fz4bEQEiFPjOSbbGqmagpgVXv78dk6y6lASuO3+HGAGiAFRtobnIuhFR4YWaK7lKIDPH5fMGj3boOjsbGnGLdJZluY05xcFPrQYszuE8g2d7z9R/g5iXn6i6lASuO3+HGAGiAFRsocddhCKnwwLjVukshMuSzNf6vjq7rrxCTWai7nLCJyTyA6+q3iJm2do/ucogM+e1Y/3nY6yx2Qo8Vx+9wYwA0wIoNlFdcASEVugyKhtfLuwdkT9EZ7g5Bbg2/5Cxerw+dB/nvaOcVV+gupwErjt/hxgBogBUbqLbOi+sHREFIhaLSSt3lELXY5tyj6FI/aAxenOnKr0F9Ph8GLT759beTVj2TexSWVEJIhZ8OiLLcnFYrjt/hxgBogFUb6Ffvr4KQCikumDNFzlJw7ETw1IDnXb4QorbOG5w/deeIZThQckJ3SUQtkry3GEIq3D9mle5SzmDV8TucGAANsGoDPTMtCUIqfJ2ar7sUomarrvXiz/Vnhvb8aC3KeDIGyqpq0fOjtcEzvmtcHIjJfr5KzYeQCn/7PEl3KWew6vgdTgyABli1gfp/4//q6MNlO3SXQtRsw5ZsC57ykXvEWvOFdNp3pBy3RPr3CBy+dJvucoiabdyyHRBSYcCiDN2lnMGq43c4MQAaYNUGmhS/G0IqvDU/TXcpRM0SlXEgeF4oT7E5U+A0BSEVojMO6C6HqFn+d34ahFSYvHq37lLOYNXxO5wYAA2wagOpdP9g+pdJ63SXQtSkPYfKgqdgjIrO0l2OZY2KygpuiL33cLnucoia9NikdRBSIcqCH1qsOn6HEwOgAVZtoIz8EgipcNe7y3WXQnROJ6rr8ND4NRBS4fEp61296KMpNXVePD7ZfyTeQ+PX8Lg4srw7RyyHkAqZ+0t0l3IGq47f4cQAaIBVG6ikoib4dREHCbKy/3yVHlzlepDbFjWpqLQyuEpafp2uuxyis6qorg2OQyUnanSXcwarjt/hxABogJUb6Of1k8Z3FB3XXQpRowLz2q7rb71TAqxs3a7DwU2y47Zyk2iypu2FxyGkwq3vxOkupVFWHr/DhQHQACs30MOf+LeOWM4J9WRBh8uqcMdw/52sUVGc99dSI+vnA94xfBkOl1XpLofoDMu2FUFIhUc+SdBdSqOsPH6HCwOgAVZuoH/OToWQCtMT9uouhagBn8+Hv8/aGJzLVlXLaQotVVlTh9996J87+dJ/N7rytBSyts8T9kJIhddmb9JdSqOsPH6HCwOgAVZuoFHR/jsEkd9t1V0KUQMLUvIgpELngdHIOmC9vzt2sa2gFDcM9B/7uGBjnu5yiBqI/G4rhFQYHZ2tu5RGWXn8DhcGQAOs3ECzk/ZBSIUXZ6boLoUoKK+4Aj8bEmPZvcHsJrDn581DY5FXzM2zyTpemJkCIRXmJOXqLqVRVh6/w4UB0AArN9DanYcgpMKD41brLoUIAFDn9QW3MXl88nrUefm1pVF1Xh/+Onkd/5uS5TwwbjWEVEjYac0FXlYev8OlRQGwpqYGeXl52L59O4qLi82qyTas3ED7jpRDSIUug6I5P4gsYXr9nKCfDYnh3aoQyj1y8q7qjETO+SX9vF4fugyKhpDKssc6Wnn8DpcmA+Dx48cxadIk3H///bjoootw4YUX4oILLsCFF16Ijh074qWXXkJKiju/ZrRyA9XUedGpfqsI7q9Guu0/dgI31YeU2Un7dJfjOF9u8E/5uGlIDPYfO6G7HHK5otJKCKlw/YAo1Fh0c3crj9/hcs4AOG7cOPzoRz/C3XffjeHDhyM2NhYZGRnYtWsXkpOTMX36dDz//PPo0KEDHnroIezcuTNcdVuC1RvoF++thJAKG3N4t5b08fl8wflAf528Dl5+TRlyXq8Pf6k/duvFmSm8609apeQUQ0iFX76/UncpZ2X18TsczhkAn3rqKWzd2vQq0qqqKkyePBnTp08PWWF2YPUG6v3ZBgip8M2mfN2lkIstTS8IrvrddZAbk5tlZ9Hx4KpglW69s1fJPb5OzYeQCk9P26C7lLOy+vgdDlwEYoDVG0h+7T9ma/zyHbpLIZcqqagJngfKPjTfh8t21B+ttxwlFdY7fovcIdCH/b+x7nGFVh+/w6HZAfD4cT2f3IuLi/H000+jXbt2aN++PV588UWUlZWd83cqKyvx2muv4Uc/+hF+8IMf4LHHHkNRUcMTMTwezxmPefPmtag2qzfQxFW7IKTCWwvSdJdCLhX4EPLbsfHc8DkMqmrr8Jux8ZYffMnZ3pqfBiEVPo3fpbuUs7L6+B0OzQ6A3bp1Q2Fh+M+d7NmzJ7p164akpCQkJCTghhtuQO/evc/5O6+++iquvfZarFy5Eqmpqbj33ntx3333NXiNx+PBzJkzUVhYGHxUVrZssYTVG2jJFv9Xb3+ZtE53KeRCG/YcCR4Gn8J5qGGTvLc4+N99w54jusshFwrMR12aXqC7lLOy+vgdDs0OgM8//zw6duyI7OyGu3qnpaXh97//fcgLA4CsrCx4PB5s3Lgx+FxMTAwuuOACFBQ03lglJSVo3bo1vvrqq+Bz2dnZ8Hg82LDh5HwEj8eDxYsXG6rP6g2UkV8S/DqIKJyqa734bf2dqAGLMnSX4zr9v8kI3nmtrrXmKkxyrjtH+M/5zsgv0V3KWVl9/A6HFs0BHDp0KH784x8jISEBO3bswOOPP44LL7wQjzzyiCnFTZ8+HR06dGjwXG1tLVq1aoVFixY1+jsrV66Ex+PBsWPHGjzfsWNHfPjhh8GfPR4Prr76avz4xz/G3XffjenTpze5cq6qqgqlpaXBR35+vqUb6HhlTfBOQGkl5wNR+Exds7v+w8cylJxg74VbyYma4CD82Zo9usshFyk5cXLcKauq1V3OWTEAnscikJEjR+Kiiy5C69at0bNnTyQnJ5tRV/BaXbp0OeP5yy+/HJMmTWr0d+bMmYM2bdqc8fzdd9+N//znP8Gfhw8fjsTERGzevBnvvfce2rZti48//vic9URGRjY6d9DKDXTXu/4J+On5x5p+MVEIHCytDG5MzDNq9QmcuXzz0FjuBUphsyXvGIRUuPtda3/zxADYggBYVFSEN954A9///vdxxx134OKLL8b8+fPP66JSykaD1KmP7OxsUwPg6YYMGYL/9//+3znrttsdQAB4fIr/6K1v0/brLoVcIjAB/NGJidzzTyOv14deExO5EIzCavHm/RBS4Ykp63WXck4MgC0IgN///vdx2223QSkFwD8X79JLL8WYMWNafNFDhw4hOzv7nI/q6mpTvwI+nVIKHo8HVVVVzf7fYYcGCqzC/HAZt+Ag86XuO7kAYUse7zrrllZ/N0ZIhdR9XIhD5htngy1gAHuM32ZrdgBsbIuUTZs24aqrrsJrr70W0qICAotAUlNTg8/FxcU1axHI119/HXxu+/btZywCOd27776LH/7why2qzw4NNGW1fy7Wv+Zu1l0KOVyd14c/fLwWQiq8/dUW3eVQvX4Lt0BIhYc/WYs63pElk/WduxlCKkxds1t3Kedkh/HbbIY3gs7JycGNN94Yiloa1bNnT9x+++1ITk5GYmIiOnfu3GAbmP3796Nr164N5iK++uqr6NixI1atWoXU1FR0794d3bt3D/75kiVLMG3aNGRmZmLXrl2YNGkSLr74YgwdOrRFtdmhgeK2Fgb/8Scy0+wk/3m0t0TG4nBZ8++kk7kOHa/CLUNjIaTCnKRc3eWQwwU+BC7bVtT0izWyw/httnMGwNzc5v1jcfToUQD+MBZqxcXF6N27Ny655BJceumleOGFFxpsBJ2TkwOPx4P4+Pjgc4GNoH/4wx/i4osvxp///OcGexjGxMTgtttuwyWXXIIf/OAH6NatG6ZMmQKvt2XbJdihgXYdLIOQCj8bEsPzQck0xyqqcduwOAipMD1hr+5y6DSfJ+yFkAq3DYvDsYpq3eWQQ/l8PtxUvwBs18FzH9igmx3Gb7OdMwBeccUVeOWVV5CSknLW15SUlOCzzz7DzTff3OQqWqexQwNV13rRqb9/DhBXApJZhi3ZBiEV/ufD1aip475zVlNT58WD41ZDSIVhS7bpLoccqqi0EkIqXD8gyvL7T9ph/DbbOQPgkSNH8NZbb6F9+/b4yU9+gj/84Q946aWX0LdvXzzzzDO4/fbb0aZNG9x7772IiooKV82WYZcGun/MKp4KQKbJOVyOGwZGQUiFtTsP6S6HzmLNjkMQUuGGgVHYd6RcdznkQOt3+0//6TFmle5SmmSX8dtMTc4BPHr0KE6cOIGvvvoKb775Jv70pz/hoYcewjPPPIOxY8ciMzMzHHVakl0aqM+MZAipMDeZ838o9P45OxVCKjw33bw9QSk0np3u/7fgtdmbdJdCDjQnKRdCKjw/w/r/Fthl/DZTkwHwsssuw+TJkzl/rBF2aaB3lmyFkAojo7J0l0IOE9j2pVN/he2Fx3WXQ03ILiwNTglJ3XdUdznkMO+qbbaZZmCX8dtMTQbAUaNGoV27dujWrRvWruVK0lPZpYG+WJ8DIRX+Pmtj0y8maiafz4c/ferfaFh+be09v+ik/3zl3xv0z58m8oM9hdTfZ6VASIUvNuzTXUqT7DJ+m6lZ28AcOHAAzz//PFq1aoWnnnrKlNW+dmSXBkrYeTh4MDxRqKj0AxBS4cbBMVxgZCNFpZW4cbB/pWZUxgHd5ZCD/GZsPIRUSNx1WHcpTbLL+G2mFu0DmJqail/+8pf4wQ9+gBEjRrTo1AwnsksD7T92Ijj5u5YrNCkEqmu9+NX7q3jKjE0FTmu4f8wqy6/WJHuoqfPipwP8i8EKjp3QXU6T7DJ+m+m8NoKeP38+rrvuOnTq1OmsR7K5gV0ayOv1oevgaAipkHOYq//IuMC+cne9uxzlVbW6y6EWKq+qxZ0jlnPfRgqZvYfLIaRC18HRtjgD3C7jt5nO+ySQqqoqvPvuu2jXrh0efPDBUNZkG3ZqoIfGr4GQCquyD+ouhWyutLIG3eo3fZ7HleW2NTfZv2Kz27A4lFbW6C6HbG5ldhGEVOj5kT3WCthp/DZLiwJgdXU10tPTMWfOHAwYMAB//OMfccUVV+DCCy80qz5Ls1MDBbbq+Jyf9smgsXHbIaTCA+NWc0qBjdXWefHb+jlb4+K26y6HbG7a2j222mLITuO3WZoMgO+88w7++te/4sYbb0Tr1q1xwQUXoEOHDvjlL3+J1157DZMnT8a6devCUavl2KmBxsRmQ0iFQYszdJdCNna4rCp41FNMJhcQ2F10hn8hz01DYnh+MxkycFEGhFT4INYeHybsNH6bpckAePPNN+Opp57CqFGjsHTp0mafD+wGdmqgr1LzIaTC09M26C6FbCxw5NsfJyRwCxEH8Pl8eOSTBAipMHyp9fduI+vq/dkGCKnwdWq+7lKaxU7jt1nOew4g2auBUvcdhZAK3Uet0F0K2dT+YyfQeWA0j3xzmMARcZ0HRdti9SZZ072jVkBIhU259thg3E7jt1kYAA2wUwMdLa+GkP4TAE5U1+kuh2wosIHwk1PX8+6fg/h8PjwxZT039KbzVlFdGxxfjlVU6y6nWew0fpuFAdAAuzVQYOXmtgJ71EvWsftQGa6v3+OLR4g5T+BIv+sHRGHPoTLd5ZDNbC0ogZAKtw2L011Ks9lt/DYDA6ABdmugwLFdKp2T96llXpuzqf44wRTdpZBJXpzpP8br9Tn2WMVJ1rE0vSB4vKBd2G38NgMDoAF2a6C3FqRBSIUJK3fqLoVsJHO//9P9df0Vsg7Yo9ep5bYVlAa/xttaUKK7HLKRT1bshJAK/7dgi+5Sms1u47cZGAANsFsDTVy1C0IqvLUgTXcpZCOBA97/NXez7lLIZP+au7n+Tu9G3aWQjbw1339zYeKqXbpLaTa7jd9mYAA0wG4NFFW/59efbHSbnvRKzz8GIRU69VecG+YCuw+VoVN//13AjHzeBaTmeXSif3pRdIZ9phfZbfw2AwOgAXZroKwDpcGjn4iaIzAv7K35vGvsFv9bfzeH8z2pOXw+H34eGQshFbIL7TEWAvYbv83AAGiA3RroRHVdcI7P0XJ7LNUnfbbknbz7t/dwue5yKEz2nHIXMD3/mO5yyOKKT9lirLLGPluM2W38NgMDoAF2bKDu9Zt1cisPasoLgbt/nDPqOoE5XS/O5F1AOrfAFkL3jV6pu5QWseP4HWoMgAbYsYGenuY/rmfBxjzdpZCFpdXf/bt+QBRyePfPdfYeLg/eBdySx7uAdHYLUvIgpMIz05J0l9Iidhy/Q40B0AA7NlDkd1shpMK7iud+0tk9PyPZdts6UGgFto16gXcB6RxGLPWfDx753VbdpbSIHcfvUGMANMCODTQ7aR+EVHhuerLuUsiiNuce5d0/Qs7h8uDpL2m8C0hn8ex0/4fFOUm5uktpETuO36HGAGiAHRsoJcc/X6P7qBW6SyGL6lN/9+/fC3n3z+3+b8EWCKnw/Ax+YKTG3Vs/r3xjTrHuUlrEjuN3qDEAGmDHBiqpqAmu2DpeWaO7HLKYU+f+7TvCu39ut+/IybuAnAtIpyutPDmelJyw13hix/E71BgADbBrA9397nIIqbAplyuBqaG/z9rIuX/UQGAu4Ev/5ekg1FDqPv90kYiRy3WX0mJ2Hb9DiQHQALs20N8+T4KQCvNT7DVng8wVOAv2uv4Ku3nqB9XbdbAM19WvCOZZ0HSqecm5EFLhb5/bawUwYN/xO5QYAA2wawMNW+JftTVsCVcC00mvzdkEIRVen7NJdylkMewNasw7S/y7Sgxfar+xxK7jdygxABpg1way86c2MsfuQ7zLQ2cXOEaSd4fpVM9Ms++3SXYdv0OJAdAAuzaQnedtkDkCqz05z4vOJjA/lKvDKcDO88ntOn6HEgOgAXZtIDuv3KLQyyuu4EpPatKpK8Tziit0l0Oa2X1HCbuO36Fk+QBYXFyMp59+Gu3atUP79u3x4osvoqzs3F9BTJ06FT169EC7du3g8Xhw7NiZg9r5vO/p7NxAdt27iUJvwKIMCKnwLDcHpyYEFpANXJShuxTSzO57ytp5/A4VywfAnj17olu3bkhKSkJCQgJuuOEG9O7d+5y/M378eIwePRqjR48+awA8n/c9nZ0byK67t1NoHSg5gc4DoyGkQgo/DFATkvf6B/3OA6NRWFKpuxzSyO6nStl5/A4VSwfArKwseDwebNx4cl5STEwMLrjgAhQUFDT5+/Hx8Y0GQKPvG2DnBrLr+Y0UWoFVfE9MWa+7FLKJx6es5y4CZPtz5e08foeKpQPg9OnT0aFDhwbP1dbWolWrVli0aFGTv3+2AHi+71tVVYXS0tLgIz8/37YNtCAlD0IqPD1tg+5SSJPi8mrcODgGQiqs3XlIdzlkE2t2HIKQCjcOjkFxebXuckiT3p9tgJAKCzbm6S7lvDAAWjwAjhw5El26dDnj+csvvxyTJk1q8vfPFgDP930jIyPh8XjOeNixgTbn+lcC3/UuVwK71bhlOyCkwiOfJMDn8+kuh2zC5/Ph4U/WQkiFD5ft0F0OaXLnCP8K4DSbLhxjANQUAKWUjQapUx/Z2dmWC4BOugNYVlUbXMF1lJ/iXaesqha3vhMHIRWiMg7oLodsRqUfgJAKt74Th/KqWt3lUJgVl1cHxw+7/v/PAKgpAB46dAjZ2dnnfFRXV1vuK+DT2b2B7hu9EkIqJO05orsUCrPP1uyBkAq/+SAedV7e/aOWqfP68OsP4iGkwrS1e3SXQ2G2Yc8RCKnwi/dW6i7lvNl9/A4FS38FHFiskZqaGnwuLi4uZItAzvd9A+zeQM/P8K8E/mLDPt2lUBhV1dYhYuRy2+7gT9YQOFEoYuRyVNXW6S6HwuiL9TkQUuGFmSm6Szlvdh+/Q8HSARDwb9dy++23Izk5GYmJiejcuXOD7Vr279+Prl27Ijn55FL0wsJCpKWlYdq0afB4PFi7di3S0tJQXFzc7PdtDrs30KioLAipMOTbTN2lUBgFBu57Rq5Ada1XdzlkU/wg4V6DF2dCSIVR0Vm6Szlvdh+/Q8HyAbC4uBi9e/fGJZdcgksvvRQvvPBCgw2bc3Jy4PF4EB8fH3zubIs1Zs6c2ez3bQ67N9BXqfkQUuHJqdwCxC341R2F0rS1nErgRk/UbwX0dWq+7lLOm93H71CwfAC0Mrs3UHq+/2inO4Yv010KhUlg8n63YZy8T8ZxMZE73T58GYRUyMgv0V3KebP7+B0KDIAG2L2BKqpPrgQ+UlaluxwyGbfvIDNwOyF3OVxWBSEVruuvcKLavnM/7T5+hwIDoAFOaKBfvu9fCbx+N1cCO93anSc38OXWPxQq3FDcXdbtPgwhFX71/irdpRjihPHbKAZAA5zQQC/OTIGQCv9dn6O7FDLZM9OSePwfmSJwLNjfPk/SXQqZbNY6/wrgv8+y7wpgwBnjt1EMgAY4oYFGR2dDSIVBizN0l0ImytxfAiEVrh8QhfyjFbrLIYfJK67A9QOiIKRC5n77zgujpg1clAEhFd6LydZdiiFOGL+NYgA0wAkNtGizfyXw41O4EtjJXp+zCUIqvDlvs+5SyKHemLcZQir0ncsec7LHJ/tXAC/evF93KYY4Yfw2igHQACc0UODO0K3vxHECt0PlHqlAp/7+xT7bCuzbq2RtWwv8/5Z06q+QV8y7zE7k8/mCq763Ftj7Tq8Txm+jGAANcEIDVdXW4af1X93wq0FnGvKtf9PWZ6cnN/1iIgP+9rl/nulQbi7vSHnFFRBS4acDomx/+osTxm+jGAANcEoD9fzIvzVI7NZC3aVQiB0pq0LXwdEQUmHd7sO6yyGHW7fLv0K06+BoFHOluePEZBZCSIXff7RWdymGOWX8NoIB0ACnNNC/F26BkArjuDec4wT2aPvjBO7RRubz+Xx45JME7jXpUOPitkNIhX4Lt+guxTCnjN9GMAAa4JQGmpG41xHL+qmhiupadBvmn6+j0nlKA4XH0vQCCKlw27A4VFTztBknCWwbNjNxr+5SDHPK+G0EA6ABTmmg5L3FEFLh3lErdJdCITSzPtj/6v1VPKeVwqa2zotfvb/KMUGBTrpn5AoIqZCSU6y7FMOcMn4bwQBogFMa6HhlDY+Ec5jaOi9+8Z7/lJcvuMk3hdkX6/2bBf/ivZWorfPqLodC4NQj4MoccI64U8ZvIxgADXBSA/UY4//EvmYHj3JygiVb/F/D3T58GSpr7L1aj+znRHUdbh++DEIqLE0v0F0OhcDqHf6jJH/9QbzuUkLCSeP3+WIANMBJDfTabP9GwZPid+suhQw6dSL++OWciE96fMgFSI7yafwuCKnw2pxNuksJCSeN3+eLAdAAJzVQ4C/36w75y+1mgcPauwziVhykz5GyKnQZ5N+CaP3uI7rLIYNem+OsmwROGr/PFwOgAU5qoMDt/d845Pa+mz0/I5nnO5MlBM6NfWEmdxiwu19/EO+oaUJOGr/PFwOgAU5qoFMn+JY7YIKvW+0oOh78/zHncLnucsjl9h4ux3X1xxDuLDquuxw6T2VVtY5bKOik8ft8MQAa4LQGihi5HEIqbHTAEn+3Cmzq/eqXqbpLIQIA/OOLVMdsHuxWKTn+rcLuGemcrcKcNn6fDwZAA5zWQC/Ub/I5a12O7lLoPBSVVuKGgf5znTfnHtVdDhEAYFPuUQipcMPAKBwsrdRdDp2HwJ6iLzroq3ynjd/ngwHQAKc10FgHHfPjRqOjsyGkwl8nr9NdClEDf5m0DkIqvBeTrbsUOg/B40LjtusuJWScNn6fDwZAA5zWQDGZBxxz0LfblFXV4pbIWAipsGxbke5yiBqI21oIIRVuiYx1xCbCbtPzo7UQUiEms1B3KSHjtPH7fDAAGuC0Bsorrgh+VVNdy9377WTa2j3+Vdxj4+HlsW9kMV6vD78Z619F+nkCj4ezk6raOvx0gH9qSf7RCt3lhIzTxu/zwQBogNMayOfz4dZ34iCkQub+Et3lUDPV1Hlx32j/sW/zknN1l0PUqLnJuRBS4b7RPB7OTjL3l0BIhVvfiXPUht5OG7/PBwOgAU5soKemboCQCgtS8nSXQs30bdp+CKlw5wge+0bWVVlThztH+I+H+24Lj4ezi/kp/uDe+7MNuksJKSeO3y3FAGiAExtoxNJtEFJh6LeZukuhZvD5fPjDx/75OZ+s2Km7HKJz+njFTgip8PAnax11N8nJhnybCSEV3lXbdJcSUk4cv1uKAdAAJzbQos35EFLhsUlcSWoH63b5j327cXAMjvLYN7K44vJqdB3sPx5u3e7DusuhZvjzp4kQUmHx5v26SwkpJ47fLcUAaIATGyhwksRNQ2JQx8UElten/tg33rEluwjcUXp+RrLuUqgJdV4fbhwc48iTXJw4frcUA6ABTmyg2jpv8AD33YfKdJdD57C90B/WO/VX2HeEx76RPew7cvJ4uB0OCxVOs+tgGYRU6Do42nE3BJw4frcUA6ABTm2gXhMTOVHbBvrVb876z9k89o3s5dUv/cfDvf0VN523ssACs0cnJuouJeScOn63BAOgAU5toEGLMyCkXsUsPgAAIABJREFUwoilzpr06yQ89o3sLHA8XOeB0TwezsICiwIHLc7QXUrIOXX8bgkGQAOc2kBfp3IhiNW9F8Nj38jeAsfDvc/j4SwrsADkm035uksJOaeO3y3BAGiAUxtozyH/vI/Og6J5IogFlVXV4uf1x77FbXXO0UzkLrH1x8P9PDIW5TweznKqauvQeaB/Pvjew86bY+zU8bslLB8Ai4uL8fTTT6Ndu3Zo3749XnzxRZSVnXtxwtSpU9GjRw+0a9cOHo8Hx44dO+M1Qgh4PJ4Gj9GjR7eoNqc2kM/nw23D/CeCpOWd+d+O9Po8Ya//2LcPeOwb2ZfX68OvP/AfDzedx8NZzub6r+lvH77MkXs2OnX8bgnLB8CePXuiW7duSEpKQkJCAm644Qb07t37nL8zfvx4jB49GqNHjz5nABw+fDgKCwuDj/Lyln3KcXIDvTAzhf8wW9Cpx77NSeKxb2Rvs5P28Xg4iwp80HxxZoruUkzh5PG7uSwdALOysuDxeLBx48bgczExMbjgggtQUND0CtX4+PhzBsDx48cbqs/JDfRJ/Y79r8/ZpLsUOgWPfSMnqaypwx3DeTycFb02ZxOEVJiw0pknDDl5/G4uSwfA6dOno0OHDg2eq62tRatWrbBo0aImf7+pAPiTn/wEP/rRj3DbbbdhzJgxqK1t2TwUJzdQYv0JE/eNXqm7FKrHY9/IiXg8nDUFvmlYt8uZJ7Y4efxuLksHwJEjR6JLly5nPH/55Zdj0qRJTf7+uQLguHHjEB8fj/T0dEyePBkdOnTAW2+9dc73q6qqQmlpafCRn5/v2AYqq6oNbtbKbRqsIZHHvpEDNTgezqFhw26KSiuDm8yXOXSBDgOgpgAopTxjAcbpj+zsbFMD4OmmT5+O733ve6iqqjrrayIjIxut1akN9ND4NRBSISbzgO5SCMBz03nsGzlT4Hi4PjwezhKiMw5ASIWHxq/RXYppGAA1BcBDhw4hOzv7nI/q6mpTvwI+3datW+HxeLB9+/azvsZNdwABoP83/g2hR0Vl6S7F9bILS4OfyHOPVOguhyik9h0pR6f6bxy2F/J4ON1GRmVBSIUBi5y3AXQAA6DFvwIOLAJJTT151FVcXFxIFoGcbvbs2bjwwgtx9GjzT1VwegMt3JjHzYYt4q0FaRBS4bXZXJRDzvTP2f7j4f5vAY+H0y2wSfdXqc7bADrA6eN3c1g6AAL+bWBuv/12JCcnIzExEZ07d26wDcz+/fvRtWtXJCef/OqgsLAQaWlpmDZtGjweD9auXYu0tDQUFxcDANavX4/x48djy5Yt2LNnD2bPno3LL78czz33XItqc3oDBQ4C78INobU6UHICPx3gP/ZtC/dlJIcK7Dt3w8AoFJZw3rEu1bVedB7kn5O5+9C599y1M6eP381h+QBYXFyM3r1745JLLsGll16KF154ocFG0Dk5OfB4PIiPjw8+d7a5ejNnzgQAbNq0Cffccw/at2+Piy66CDfddBNGjRp1zvl/jXF6A3m9Ptz6jn9D6PR8Bg9dRtV/HfP4lPW6SyEy1eOT13PaiWZb8o5BSIVuw+IcvSrb6eN3c1g+AFqZGxqozwz/woOZidwQWofSyhrcPNR/7NvK7CLd5RCZakVWEYRUuGVoLEora3SX40ozEv0bQD/v8AU5bhi/m8IAaIAbGuij5f49ut6Yt1l3Ka40KX43hFT4nw9X89g3cjyv14cHx62GkAqTV+/WXY4r/WvuZgip8LHD9xp1w/jdFAZAA9zQQGt3HoKQCr98nxtCh1tVbR3uene54ydjE50qsPjs7neXo6qWp92E2y/e828AnbDT2XsyumH8bgoDoAFuaKDSyprghtCHjrdsjiQZMz8lF0Iq3DNyBRfhkGtU13oRMdL/wWdBSp7uclzl4HH/BtDX9Vc47vCv4N0wfjeFAdAAtzTQ/3zo/0ombmuh7lJcw+v14Tdj4yGkwmdr9uguhyispq7xT3347dh4Tn0Io9ithRBS4XcfOncD6AC3jN/nwgBogFsaSH6dDiEVRkdn6y7FNeLq/yG+JTLWsUcxEZ3N8coa3BLpX/y0bBsXP4XLqGj/jgP9v0nXXYrp3DJ+nwsDoAFuaaAFKXnchiTMHqvfiPX9GIZucqf3YrIhpMJfJnEj+nAJbMOzYKPzv3p3y/h9LgyABrilgXYWHYeQCl0HR6OmjnPRzLYxpxhCKnQeGI2Dx7khLrnTwdJKdB7o35B4Y06x7nIcr7rWi66D/f+9dx10/nF8bhm/z4UB0AC3NJDX68Ntw/wbQqfu4z/EZvv7rI2u+RqG6FwC00/+Pmuj7lIcL/DB8/bhy1wx79It4/e5MAAa4KYGCpzT+YnD94bSbUf93dbr+jv7GCai5th9qCy4C8HOIuffldLp4xU7XXXeuJvG77NhADTATQ30xYZ9EFLhyamcB2imt+anQUiFV79M1V0KkSX84wv/h8+3FqTpLsXRnpjin//35YZ9uksJCzeN32fDAGiAmxpoz6Gy4Ly0yhpuzmqGvOIKXD8gimcvE50icDbtTwdEIf9ohe5yHOlEdV1wvuXew+W6ywkLN43fZ8MAaICbGsjn8+GekSsgpELiLmfvEK/LkG8zIaTCM9OSdJdCZClPT9sAIRWGfpupuxRHSth5GEIq3DtqBXw+58//A9w1fp8NA6ABbmugwNeTY2K5NUmoHS6rQpdB/k/g6xiwiRpI3OUPKF0GReNwGU8kCrX367fccdPX7G4bvxvDAGiA2xoocEbnnz5N1F2K44yJ9f8D3Gtioms+gRM1l8/nQ68JCRBS4YPY7brLcZxHJya67sxxt43fjWEANMBtDbT/2AkIqXD9gCjHnxMZTqeeehDL4/aIGhWTefJ0HP77EzqllTXoVL/SuuDYCd3lhI3bxu/GMAAa4MYG6jFmFYRUWJnN45lCZVK8/9zTB8atdsX+W0Tnw+v14bf152NPXr1bdzmOsSKrCEIq/PqDeN2lhJUbx+/TMQAa4MYG6v9NBoRUGL50m+5SHKGypg53jljuuq9fiM5HYBrKXe8u524EITJsyTYIqTBgUYbuUsLKjeP36RgADXBjAy3ZUgAhFXp+tFZ3KY7wZf3+iveNXslj9oiaUF3rRfdRK1y1X53ZHhq/BkIqLE0v0F1KWLlx/D4dA6ABbmygQ8erIKR/vkhxebXucmytps6L+0avhJAKMxP36i6HyBZmJu6FkAq/eI8fmow6Unby33O3ra524/h9OgZAA9zaQL/70P+JMSrjgO5SbG1Biv/rrDtH8OssoubyT5tYBiEVFmzM012Oran0AxBS4aHxa3SXEnZuHb9PxQBogFsbKPK7rRBSYdBid80ZCaXaOm9wQc3UNZzQTtQSU1b7F071GLMKdVw4dd4GLvLP6X5nyVbdpYSdW8fvUzEAGuDWBlq2zb9q7Ddj43WXYluLN++HkAq3DYtDeVWt7nKIbKW8qha3DYuDkArfpu3XXY5t/eYD/6rq5dvct6uDW8fvUzEAGuDWBio5cXLfqMKSSt3l2I7X68MD41ZDSIWJq3bpLofIlias3AkhFR7k9knn5UCJf1/XTv0VSl24r6Jbx+9TMQAa4OYGCuzK/80mbl3SUoF5Nz/nhrZE5620sgY/r99AnfORW+7r1Pzg6UNu5ObxO4AB0AA3N9DoaP/RZW/O26y7FFvx+XzBbRc+XLZDdzlEtjZu2Y7gtlQ8QrFl3py3GUIqvBfjzrPd3Tx+BzAAGuDmBkreWwwhFW59J45bMbRAYP7kz4bE4FgFt9EhMuJYRTV+NiTGtfPYzldNnTd49zQlp1h3OVq4efwOYAA0wM0NVOf14fbh/q0Y1u8+orscWzj1QHu3fuomCrXAtxG9JiTwLmAzrdt9GEIq3DF8mWtXUbt5/A5gADTA7Q3074VbIKTCsCU8Fq45Amdu3jg4BkdctukqkVkOl1XhxsExPKO8Bd5Z4t/Kq9/CLbpL0cbt4zfAAGiI2xsoJrMQQir88v2V/OTdBJ/Phz98vBZCKoyKztJdDpGjjIrKgpAKD3/CuYBN8fl8+OX7/hOIYrcW6i5HG7eP3wADoCFub6CK6lp0HhQNIRW2Fx7XXY6lBcLyz4bE8Ag9ohArLj85F9DNoaY5sgtLIaRCl0HRqKh27x6kbh+/AQZAQ9hAwAszU7ifXRO83pMrfz+I3a67HCJHGhObHTzWjPsCnl1g/8QXZ6boLkUrjt8MgIawgYC5ybmu3kuqOZamF0BIhVsiY1FSwX3/iMxwrKIatwz1r2xV6dwX8GwCC9HmJefqLkUrjt82CIDFxcV4+umn0a5dO7Rv3x4vvvgiysrKzvn6vn37okuXLrjoootw7bXX4l//+hdKSkoavC43Nxd/+MMf8P3vfx+XX345+vXrh9ralt0OZwMBB0srIaT/VJCiUp4Kcrq6U079+Gj5Tt3lEDna+OX+fQEfGLfatatbz6Wo/t/r6/orHDzu7n+vOX7bIAD27NkT3bp1Q1JSEhISEnDDDTegd+/eZ319ZmYmHnvsMSxZsgS7d+/GypUr0blzZ/zlL38Jvqaurg633HILHnzwQaSlpSE6OhqXXXYZBgwY0KLa2EB+j05MhJAKc5Lc/YmyMYEzf299J86Vxy0RhVNpZQ1ufYdnBJ/N7KR9EFLhT5/yGxuO3xYPgFlZWfB4PNi4cWPwuf/f3r1HRXXdewA/hgiaKpJEGpsYp1GB5hqjNA1qeyPp0ipNvLraJmlIa43psulKTPO0xxiV+k7UahOj0XIR05uqIff6PPKID3yAioKoKCOKjjogihmd4TnAzHzvHwNHR0FnGGbOmTnfz1rnD48zzI/t9uzfnLP3b2dkZKBTp04oLy93++ekpaUhNDRUvsOXnp6Oe+65B5cv3ygZ8MUXXyA8PBwNDe5P0GcHcvp81xnoRAmvrs5TOhRVabLZ8WzzZuucI0nkHy3Xo2cXZaOJRepdTFidB50oYXk2r0ccv1WeAKakpCAiIsLlXFNTE0JCQrBhwwa3f05ycjJ69uwp/3nGjBkYNGiQy2vOnTsHQRBw5Ij7W5uxAzmVXK6CTpQQ9VE6aqzaXVV2q68PX4ROlBA7+1u2C5GfVFubMHiW8y5g2uGLSoejGjXWJkRNc1ZtOH2ZVRs4fqs8AZw3bx6io6NvOx8ZGYkVK1a49TOuXr2KPn36YNq0afK5SZMmYdSoUS6vq62thSAISE9Pb/NnWa1WWCwW+TAajZrvQICzrtTwhbugEyVkFHHyNQDUN9rw0wXOWlur9pQqHQ6RpqzcXQqdKOGnC3aivtGmdDiqkH78EnSihPiFu1grEUwAAYUSQFEUIQjCHQ+9Xu91AmixWBAXF4eEhAQ0Nt6Yf9XeBDApKanVWLXcgVrM3noSOlHCe19rt7L8zVbtcQ5AQ+fv4ABE5Gf1jTYMmbcDOlHCP/ecVTocVXj360LoRAlztnLnJoAJIKBQAlhZWQm9Xn/Ho6GhwatHwFVVVRg2bBhGjBiB+nrX1U7tfQTMO4BtO3D2O+hECYNnZWl+3o259sZE9K/5CIpIEV8fuigvwDLXaXsBVpPNjkHNj8UPnuXe7QATQEDlj4BbFoHk5+fL57Kysu66CMRisWDo0KGIj49HbW3tbX/fsgjkypUr8rlVq1YhPDwcVqv7e7SyA93QZLPL8272lFQqHY6i5qc7t6X6xRKWoiBSis3uwC+WOEswLUjXKx2OonaXVPIL+i04fqs8AQScZWBiY2ORl5eHnJwcREVFuZSBKSsrQ0xMDPLynCtQLRYLhgwZgoEDB6K0tBQVFRXyYbM5H8W1lIEZNWoUjh49iszMTERGRrIMjJdmbCqCTpTw1lr3F9IEm/LrdfL2eDuKuTE9kZK2n7wsb3t2yVyndDiKmbz2CHSihJmbipQORTU4fgdAAmgymZCYmIhu3bohPDwcEydOdCkEbTAYIAgCsrOzAQDZ2dltzis0GAzy+86fP49f/vKX6Nq1K3r27In333+fhaC9dMx4Xb7YarXm3QdpR6ETJby4cj8nWhMpzOFw4MUv9kMnSpjyjTbnJ1vqGxHd/KX0uNF89zdoBMfvAEgA1YwdyJXD4cDI5l0vtLjN0KmKKjw21bkrSsGFa0qHQ0QA8s9fg06U8NhUCSUaLH/Ssl3nL5bs5pfSm3D8ZgLoFXag233RXH7hNytylQ7F715LPQSdKOHP/5N/9xcTkd+8/q986EQJf1xzSOlQ/O7XK3KhEyWs3M1yVDfj+M0E0CvsQLe7bKmX74IZrtYoHY7f5J65Cp0ooe+H21Ba2fZe1UTkf6WV1ej74TboRAm5pVeVDsdvzl2tke9+XuFe7S44fjMB9Ao7UOv+kOLcbujvWaeUDsUvmmx2ebXhDE6yJlKl6Rudi9RGLdmjmZWwi7NOQSdKmMBtOm/D8ZsJoFfYgVq3+Wi5XIXfroEyKKk556ATJQyalYXrte7vJU1E/nOtpkGuhbcm16B0OD5ntzvk3Yi2HG27bJpWcfxmAugVdqDW1Tfa8ERSpiYet5hqGjCw+Xf9nwPnlQ6HiO7gXwfOy8WhTTXB/WWtZVrKE0mZ3I2oFRy/mQB6hR2obVP/77gmtob7cIPz90z4x14WfSZSOZvdgdFL90AnSpi24bjS4fhUy9ZvHwb579leHL+ZAHqFHaht+edN0IkSHp+RgRqrZ/UVA0VRmRk/bF7wwu2ViAJDy7aVP5wq4UR5cNbFq7Y24UfTM6ATJeSfZ0mq1nD8ZgLoFXagtjkcDjy7KBs6UcI3+Ualw+lwDocDL3zhLK8wWcM7nxAFojf/XeAs2P5FcBZsTzvs3Af52UXZQfn7dQSO30wAvcIOdGfLdp6GTpQw7vOcoLsIbSosg06UEDM9HeXXtbvFFFEgKr9eh5jpzt0xNhWWKR1Oh3I4HBj3eQ50ooRlO08rHY5qcfxmAugVdqA7q6yyyvviBtMjUnNtI56asx06UcJnO3iBJQpEn+5wfkH9ydztMNcGz9aVLY+4oz5KR2WVVelwVIvjNxNAr7AD3V3LIolXg6gOVct+vz9fnM3VdUQBqr7Rhp8vzg66fYInrM7TxCIXb3H8ZgLoFXaguzNcrZF3BtFXBH477Tt9VZ5AfthgUjocIvLCIYNzsZpOlJBzJvBLVhVfssg7f5z/Tjs7MbUHx28mgF5hB3LPG185J1y/s75Q6VC8UtvQhP/8ZCd3/CAKIi07hDzzyS7UNgR2xYK31x2BTpTwxr8LlA5F9Th+MwH0CjuQe44bzfI+ucZrtUqH026zt56ETpQwbP4OVAdpaRsiramqb8Sw+TugEyXM2XpS6XDa7aKpVt7vuKgsOMvbdCSO30wAvcIO5L5Xkg9AJ0pI2nxC6VDapfDidflR9q5TV5QOh4g60C79FfnR6dGL15UOp12SNp+ATpTwu+SDSocSEDh+MwH0CjuQ+/aeroROlPCj6Rm4FmBbMDU02TFqyZ6geIxNRK1reXw6eukeNDTZlQ7HI6aaBrmszb7TgT+X0R84fjMB9Ao7kPscDgee+3QvdKKEpdtLlA7HI/O2FUMnSoid/W3Q7x9KpFWmmgbEzv4WOlHC/G3FSofjkSXflkAnSnj+s71BV3PVVzh+MwH0CjuQZ7YcLYdOlDB4VlbATLbedeqKvEow80SF0uEQkQ9lnqiQ/79nB8hUj9qGJgyelQWdKGHrsXKlwwkYHL+ZAHqFHcgzTTY7hi/cBZ0oYWGmXulw7uqypV6+IzCTq36JNGHGJueq4B/P/haXLfVKh3NXn2TooRMlDF+4CzY77/65i+M3E0CvsAN5LqPI+Q27/7RtKK2sVjqcNtnsDry8yrlw5Zf/2MuCz0QaUd9oQ8I/nNNVEv95QNVJ1Zkr1eg/bRufULQDx28mgF5hB/Kcw+GQK9X/LvmgauerfNa8TdTjMzJUnagSUccrrazG4zMyVL2frsPhkKsrvLo6T7XXUrXi+M0E0CvsQO1z/rsaeY9gNc5ZyTtnkku+/G++UelwiEgB3+Qb5dIwh1S460/LnOqoj9K560c7cPxmAugVdqD2a1m1Fjdvu6qKKhuv1eInc7dDJ0p4lyVfiDTL4XDgnfWF0IkSfjJ3u6qK2FdbmxA3b3tAVlVQC47fTAC9wg7UfvWNNjzziXNByFxJHdX3zXWNGPn33XItMDUlpkTkf9XWJrkG6Mi/74a5rlHpkAAAc5p3JRq+cBfnJ7cTx28mgF5hB/JOS/X9vh9uw6mKKkVjsTbZ8NtV++W7kpfMdYrGQ0TqUH69Tr7b9vKqA4oXidZXWOQt37grUftx/GYC6BV2IO9N+vIwdKKE/1q2T7Fvsg6HQ94FYMDMTJws578nEd1wotyM/2heFPLO+kLFFlzUN9rwX8v2QSdK+NO/DisSQ7Dg+M0E0CvsQN4ru16HJ//mLGI6ee0RRS6sizJPyXci95RU+v3ziUj9dpdUynfeFmed8vvnOxwOvPnvAuhECU/+LQtl1/mUwhscv5kAeoUdqGPkll5Fv+YLqz8nNDscDizPPiNX/v/60EW/fTYRBZ71hy7I14vl2Wf8+oW1ZeFcvw+3YX/pd3773GDF8ZsJoFfYgTrOurwbF9YtR31fGsZudyBp8wn5M/+xXZ21vohIXZZuL5GvG3/bcgJ2PxSK3txc8kUnSlh/6ILPP08LOH4zAfQKO1DHmis5V7ZFf5SOIxeu+exzrE02vNH8KEUnSkjee9Znn0VEwSd571n5+vHmvwtgbfLd/OWCC9fkuqnzthX77HO0huM3E0CvsAN1LJvdgddSD0EnSnhqznafrAyuqm9E4j8PyNvRbSos6/DPIKLgt6mwTN6G7ZXkA6iq7/gSMacqqvDUHOcK5D+uOazqbekCDcdvJoBeYQfqeNXWJoxe6qy7FTM9Hd904E4c+eevyXX+/mNGBvae5oIPImq/PSWV8pZxI/++G/nnO+7JRdrhi4iZni7XJa1hXdIOxfE7ABJAk8mEV155Bd27d0ePHj3w2muvobq67b1ZTSYTJk+ejOjoaHTp0gWPPvoo3nrrLZjNZpfXCYJw27Fu3TqPYmMH8g1TTQPGp+TJj1imfHMUdQ3tf8RiqW/ERxuP44fN27s9NWc7jhvNd38jEdFdHDNel+/S/XCqhOkbi2Dx4m5gXYMNU745Kl//xqfkwVTT0IERE8DxGwiABDAhIQGDBg3CwYMHsW/fPvTv3x+JiYltvr6oqAi//vWvsWXLFpSWlmLnzp2IiorCb37zG5fXCYKA1NRUVFRUyEd9fb1HsbED+Y7d7sBnO07Le/KOXroH+grP2tnhcCD9+CU83by1m06U8H7aUVzjxZSIOpCppgHvp91I2uLmbUf68UserxLWV1jkJyCPTZXw2Y7TfllkokUcv1WeABYXF0MQBBw+fKPgZUZGBjp16oTycvdXiqalpSE0NBRNTTduoQuCgI0bN3oVHzuQ7+WeuYqn5nwrX1jHfZ6DL/cb2vxG7HA4cMx4HfPTi+Wt5nSihPiFu5B75qqfoyciLck9cxXxC29cd575ZBfmpxfjmPF6m8mgqaYBa3INGPt5jvy+p+Zs5/XKxzh+qzwBTElJQUREhMu5pqYmhISEYMOGDW7/nOTkZPTs2dPlnCAIePjhh/Hggw/i6aefRkpKisff1tiB/OOKpR6TvjwsF2FtqYU1PiUPb609grfWHsFf1h3B5LVH8LOPd8qvaZlHuCjzFPfLJCK/qG+0YVHmKXn+Xsvxs493YnLztarlujU+JU+ugdpSjH7Sl4dxxeLZ0yjyHMdvlSeA8+bNQ3R09G3nIyMjsWLFCrd+xtWrV9GnTx9MmzbN5fzs2bORk5ODI0eO4OOPP0ZYWBg+/fTTO/4sq9UKi8UiH0ajUfMdyJ+uVNUjee9ZPPfpXpcL663Hj6Zn4I2vCiAdu8SJ00SkiBprE7YeK8cbXxXgR9Mz7njNev6zvfjvfedQWWVVOmzNYAKoUAIoimKrizBuPvR6vdcJoMViQVxcHBISEtDYeOdJuTNmzEDv3r3v+JqkpKRWY9VyB1JKyeUqrMk14L/3nZOP5L1nkVFU4dWCESKijlbXYENGUQWS9551uWatyTWg5HLHl7uiu2MCqFACWFlZCb1ef8ejoaHBq0fAVVVVGDZsGEaMGOHW4g5JkiAIAqzWtr+B8Q4gERFR4GMCqPJHwC2LQPLz8+VzWVlZd10EYrFYMHToUMTHx6O2ttatz5o7dy7uv/9+j+JjByIiIgo8HL9VngACzjIwsbGxyMvLQ05ODqKiolzKwJSVlSEmJgZ5eXkAnP+oQ4YMwcCBA1FaWupS5sVmcz4a3LJlC5KTk1FUVIQzZ85gxYoVuO+++zBz5kyPYmMHIiIiCjwcvwMgATSZTEhMTES3bt0QHh6OiRMnuhSCNhgMEAQB2dnZAIDs7Ow25xUaDAYAzlIygwcPRrdu3fC9730PgwYNwsqVK2G32z2KjR2IiIgo8HD8DoAEUM3YgYiIiAIPx28mgF5hByIiIgo8HL+ZAHqFHYiIiCjwcPxmAugVdiAiIqLAw/GbCaBX2IGIiIgCD8dvJoBeYQciIiIKPBy/mQB6hR2IiIgo8HD8ZgLoFXYgIiKiwMPxmwmgV9iBiIiIAg/HbyaAXmEHIiIiCjwcv5kAesVsNkMQBBiNRlgsFh48ePDgwYNHABxGoxGCIMBsNiudSiiGCaAXWjoQDx48ePDgwSPwDqPRqHQqoRgmgF6w2+0wGo0wm80++3bCu4tsK7YV20rtB9uKbRVobWU2m2E0GmG325VOJRTDBFClLBbOT3AX28p9bCv3sa3cx7ZyH9vKfWwr32ICqFLs+O5jW7mPbeU+tpX72FbuY1u5j23lW0wAVYod331sK/exrdzHtnIf28p9bCv3sa18iwmgSlmtViQlJcFqtSodiuqxrdzHtnIf28p9bCv3sa3cx7byLSaARERERBrDBJCIiIhIY5gAEhHo/RIEAAAJZklEQVQREWkME0AiIiIijWECSERERKQxTABVYu7cuRg2bBi6du2KHj16uPWeCRMm3LatzejRo30cqfLa01YOhwMzZsxAr1690KVLF4wYMQKnT5/2caTqYDKZ8Morr6B79+7o0aMHXnvtNVRXV9/xPfHx8bf1rddff91PEfvP559/Dp1Oh7CwMMTFxSEvL++Or09LS0NMTAzCwsLwxBNPYNu2bX6KVHmetFVqaupt/ScsLMyP0Spnz549GDNmDH7wgx9AEARs3Ljxru/Jzs5GbGwsQkND0a9fP6Smpvo+UBXwtK2ys7Nb3c6toqLCTxEHFyaAKjFz5kwsWbIE7733nkcJYEJCAioqKuTj2rVrPo5Uee1pq48//hg9evTApk2bcOzYMYwdOxaPPfYY6uvrfRyt8hISEjBo0CAcPHgQ+/btQ//+/ZGYmHjH98THx2PSpEkufSvYanGtX78eoaGhWL16NU6ePIlJkyYhIiICV65cafX1ubm5CAkJwcKFC1FcXIzp06ejc+fOKCoq8nPk/udpW6WmpiI8PNyl/1y+fNnPUSsjPT0dH330ETZs2OBWUnPu3Dncd999eO+991BcXIxly5YhJCQEmZmZfopYOZ62VUsCWFJS4tK3tLydmzeYAKpMamqqRwnguHHjfByRernbVg6HA7169cKiRYvkc2azGWFhYVi3bp0vQ1RccXExBEHA4cOH5XMZGRno1KkTysvL23xffHw83n77bX+EqJi4uDi8+eab8p/tdjsefvhhLFiwoNXXv/TSS3j++eddzg0ZMiQo74zeytO28uQ6FszcSWr++te/YsCAAS7nfvvb32riac7NPEkAr1+/7qeoghsTQJXxNAHs0aMHIiMjER0djT//+c/47rvvfByherjbVmfPnoUgCCgsLHQ5P3z4cPzlL3/xVXiqkJKSgoiICJdzTU1NCAkJwYYNG9p8X3x8PHr27IkHH3wQAwYMwNSpU1FbW+vrcP2moaEBISEhtw04f/jDHzB27NhW3/Poo49i6dKlLudmzpyJJ5980mdxqkF72io1NRUhISHo06cPevfujbFjx+LEiRP+CFdV3Elqnnnmmdu+bK1evRrh4eG+DE11PEkAdTodevXqhZEjRyInJ8dPEQYfJoAq40kCuG7dOmzevBnHjx/Hxo0b8fjjj+Ppp5+GzWbzcZTq4G5b5ebmQhAEXLp0yeX8iy++iJdeeslX4anCvHnzEB0dfdv5yMhIrFixos33rVq1CpmZmTh+/Di++uorPPLII/jVr37ly1D9qry8HIIgYP/+/S7np0yZgri4uFbf07lzZ6xdu9bl3PLly/H973/fZ3GqQXvaav/+/fjyyy9RWFiI3bt3Y8yYMQgPD4fRaPRHyKrhTlITFRWF+fPnu5zbtm0bBEFAXV2dL8NTFXfa6tSpU1i5ciXy8/ORm5uLiRMn4t5770VBQYGfogwuTAB9SBTFVies3nzo9XqX93jz6KTlTteOHTs6Iny/8mVbBWMC6G57tTcBvNXOnTshCAJKS0s78tdQDBNA97WnrW7V2NiIfv36Yfr06b4IUbWYALrP3QUztxo+fDh+//vf+yCi4McE0IcqKyuh1+vveDQ0NLi8x9u5Mz179sTKlSu9Dd3vfNlWwfgI2N32au8j4FvV1NRAEISgmZjOR8Dua09bteaFF17Ayy+/3NHhqRofAbuvvQngBx98gKFDh/ogouDHBFBlvEkAjUYjOnXqhM2bN3dwVOrk6SKQxYsXy+csFoumFoHk5+fL57Kysu66CORWOTk5EAQBx44d80WYioiLi8PkyZPlP9vtdjzyyCN3XAQyZswYl3PDhg3TzCIQT9rqVjabDTExMXj33Xd9FaIqubsI5IknnnA5l5iYyEUgbho5cmRQTU/xJyaAKnHhwgUUFhZi1qxZ6NatGwoLC1FYWOhSry0mJka+a1NdXY0PPvgABw4cgMFgwI4dO/DjH/8YUVFRsFqtSv0afuFpWwHOMjARERHynMlx48ZpqgxMbGws8vLykJOTg6ioKJcyMGVlZYiJiZHrupWWlmL27NnIz8+HwWDA5s2b0bdvXwwfPlypX8En1q9fj7CwMKxZswbFxcX405/+hIiICLlcyfjx4zF16lT59bm5ubj33nuxePFi6PV6JCUlaaoMjCdtNWvWLGRlZeHs2bMoKCjAyy+/jC5duuDkyZNK/Qp+U11dLV+TBEHAkiVLUFhYiAsXLgAApk6divHjx8uvbykDM2XKFOj1eixfvlwzZWA8baulS5di06ZNOHPmDIqKivD222/jnnvuCchpT2rABFAlWivqLAgCsrOz5dcIgiAXCK2rq8OoUaMQGRmJzp07Q6fTYdKkSZqoteVpWwE3CkE/9NBDCAsLw4gRI1BSUuL/4BVgMpmQmJiIbt26ITw8HBMnTnRJlg0Gg0v7Xbx4EcOHD8cDDzyAsLAw9O/fH1OmTAm6OoAAsGzZMvTp0wehoaGIi4vDwYMH5b+Lj4/HhAkTXF6flpaG6OhohIaGYsCAAZoqBO1JW73zzjvyax966CE899xzOHLkiAJR+19bxYpb2mfChAmIj4+/7T2DBw9GaGgo+vbtq5lC0J621SeffIJ+/fqhS5cueOCBB/Dss89i165dygQfBJgAEhEREWkME0AiIiIijWECSERERKQxTACJiIiINIYJIBEREZHGMAEkIiIi0hgmgEREREQawwSQiIiISGOYABIRERFpDBNAIiIiIo1hAkhERESkMUwAiYhusnbtWnTp0gWXLl2Sz7366qsYOHAgzGazgpEREXUcJoBERDdxOBx48sknMXnyZADAzJkz0bt3b5SVlSkcGRFRx2ECSER0i61btyIsLAxz587F/fffjxMnTigdEhFRh2ICSETUitjYWISGhmL37t1Kh0JE1OGYABIR3SIjIwNdu3ZFSEgI9Hq90uEQEXU4JoBERDcpKChA9+7dsXbtWowaNQovvPCC0iEREXU4JoBERM0MBgN69eqFBQsWAAAOHjyITp06oaCgQOHIiIg6FhNAIiIAJpMJMTExeP31113OP/fccxg9erRCURER+QYTQCIiIiKNYQJIREREpDFMAImIiIg0hgkgERERkcYwASQiIiLSGCaARERERBrDBJCIiIhIY5gAEhEREWkME0AiIiIijWECSERERKQxTACJiIiINIYJIBEREZHGMAEkIiIi0pj/B4DVDxwtMPgoAAAAAElFTkSuQmCC\" width=\"640\">"
|
|
],
|
|
"text/plain": [
|
|
"<IPython.core.display.HTML object>"
|
|
]
|
|
},
|
|
"metadata": {},
|
|
"output_type": "display_data"
|
|
},
|
|
{
|
|
"data": {
|
|
"text/plain": [
|
|
"Text(0, 0.5, '$V(x)$')"
|
|
]
|
|
},
|
|
"execution_count": 4,
|
|
"metadata": {},
|
|
"output_type": "execute_result"
|
|
}
|
|
],
|
|
"source": [
|
|
"plt.plot(xgrid, Vgrid)\n",
|
|
"plt.xlabel('$x$')\n",
|
|
"plt.ylabel('$V(x)$')"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"def deriv(X, t, gamma, delta, omega):\n",
|
|
" \"\"\"Return the derivatives dx/dt and d2x/dt2.\"\"\"\n",
|
|
" \n",
|
|
" x, xdot = X\n",
|
|
" xdotdot = -dVdx(x) -delta * xdot + gamma * np.cos(omega*t)\n",
|
|
" return xdot, xdotdot\n",
|
|
"\n",
|
|
"def solve_duffing(tmax, dt_per_period, t_trans, x0, v0, gamma, delta, omega):\n",
|
|
" \"\"\"Solve the Duffing equation for parameters gamma, delta, omega.\n",
|
|
" \n",
|
|
" Find the numerical solution to the Duffing equation using a suitable\n",
|
|
" time grid: tmax is the maximum time (s) to integrate to; t_trans is\n",
|
|
" the initial time period of transient behaviour until the solution\n",
|
|
" settles down (if it does) to some kind of periodic motion (these data\n",
|
|
" points are dropped) and dt_per_period is the number of time samples\n",
|
|
" (of duration dt) to include per period of the driving motion (frequency\n",
|
|
" omega).\n",
|
|
" \n",
|
|
" Returns the time grid, t (after t_trans), position, x, and velocity,\n",
|
|
" xdot, dt, and step, the number of array points per period of the driving\n",
|
|
" motion.\n",
|
|
" \n",
|
|
" \"\"\"\n",
|
|
" # Time point spacings and the time grid\n",
|
|
"\n",
|
|
" period = 2*np.pi/omega\n",
|
|
" dt = 2*np.pi/omega / dt_per_period\n",
|
|
" step = int(period / dt)\n",
|
|
" t = np.arange(0, tmax, dt)\n",
|
|
" # Initial conditions: x, xdot\n",
|
|
" X0 = [x0, v0]\n",
|
|
" X = odeint(deriv, X0, t, args=(gamma, delta, omega))\n",
|
|
" idx = int(t_trans / dt)\n",
|
|
" return t[idx:], X[idx:], dt, step"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"# Set up the motion for a oscillator with initial position\n",
|
|
"# x0 and initially at rest.\n",
|
|
"x0, v0 = 0, 0\n",
|
|
"tmax, t_trans = 18000, 300\n",
|
|
"omega = 1.4\n",
|
|
"gamma, delta = 0.4, 0.1\n",
|
|
"dt_per_period = 100"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"# Solve the equation of motion.\n",
|
|
"t, X, dt, pstep = solve_duffing(tmax, dt_per_period, t_trans, x0, v0, gamma, delta, omega)\n",
|
|
"x, xdot = X.T"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"# The animation\n",
|
|
"fig, ax = plt.subplots(nrows=2,ncols=2)\n",
|
|
"ax1 = ax[0,0]\n",
|
|
"ax1.plot(xgrid, Vgrid)\n",
|
|
"ln1, = ax1.plot([], [], 'mo')\n",
|
|
"ax1.set_xlabel(r'$x / \\mathrm{m}$')\n",
|
|
"ax1.set_ylabel(r'$V(x) / \\mathrm{J}$')\n",
|
|
"\n",
|
|
"# Position as a function of time\n",
|
|
"ax2 = ax[1,0]\n",
|
|
"ax2.set_xlabel(r'$t / \\mathrm{s}$')\n",
|
|
"ax2.set_ylabel(r'$x / \\mathrm{m}$')\n",
|
|
"ln2, = ax2.plot(t[:100], x[:100])\n",
|
|
"ax2.set_ylim(np.min(x), np.max(x))\n",
|
|
"\n",
|
|
"# Phase space plot\n",
|
|
"ax3 = ax[1,1]\n",
|
|
"ax3.set_xlabel(r'$x / \\mathrm{m}$')\n",
|
|
"ax3.set_ylabel(r'$\\dot{x} / \\mathrm{m\\,s^{-1}}$')\n",
|
|
"ln3, = ax3.plot([], [])\n",
|
|
"ax3.set_xlim(np.min(x), np.max(x))\n",
|
|
"ax3.set_ylim(np.min(xdot), np.max(xdot))\n",
|
|
"\n",
|
|
"# Poincaré section plot\n",
|
|
"ax4 = ax[0,1]\n",
|
|
"ax4.set_xlabel(r'$x / \\mathrm{m}$')\n",
|
|
"ax4.set_ylabel(r'$\\dot{x} / \\mathrm{m\\,s^{-1}}$')\n",
|
|
"ax4.scatter(x[::pstep], xdot[::pstep], s=2, lw=0)\n",
|
|
"scat1 = ax4.scatter([x0], [v0], lw=0, c='m')\n",
|
|
"\n",
|
|
"plt.tight_layout()\n",
|
|
"\n",
|
|
"def animate(i):\n",
|
|
" \"\"\"Update the image for iteration i of the Matplotlib animation.\"\"\"\n",
|
|
" \n",
|
|
" ln1.set_data(x[i], V(x[i]))\n",
|
|
" ln2.set_data(t[:i+1], x[:i+1])\n",
|
|
" ax2.set_xlim(t_trans, t[i])\n",
|
|
" ln3.set_data(x[:i+1], xdot[:i+1])\n",
|
|
" if not i % pstep:\n",
|
|
" scat1.set_offsets(X[i])\n",
|
|
" return\n",
|
|
"\n",
|
|
"anim = animation.FuncAnimation(fig, animate, frames=len(x), interval=1)"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": []
|
|
}
|
|
],
|
|
"metadata": {
|
|
"anaconda-cloud": {},
|
|
"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.5.2"
|
|
}
|
|
},
|
|
"nbformat": 4,
|
|
"nbformat_minor": 1
|
|
}
|