{ "cells": [ { "cell_type": "markdown", "id": "a40f6857", "metadata": {}, "source": [ "## Carry out a few single indivual steps of diffusion, \n", "### and directly verify that the values satisfy the diffusion equation\n", "\n", "In this \"PART 2\", we quickly repeat all the steps discussed at length in part1,\n", "but using a much-finer resolution.\n", "This time, we'll start with slightly more complex initial concentrations built from 2 superposed sine waves.\n", "\n", "**We'll also explore the effects of:** \n", "-Spatial resolution (\"delta x\") \n", "-Temporal resolution (\"delta t\") \n", "-Alternate methods of estimating numerical derivatives" ] }, { "cell_type": "markdown", "id": "3d92a505-55d3-4dce-a458-2beafcb4d9c7", "metadata": {}, "source": [ "### TAGS : \"diffusion 1D\", \"under-the-hood\"" ] }, { "cell_type": "code", "execution_count": 1, "id": "e6eea4ef-2b46-46fb-8317-710b5b146496", "metadata": {}, "outputs": [], "source": [ "LAST_REVISED = \"May 3, 2025\"\n", "LIFE123_VERSION = \"1.0.0rc3\" # Library version this experiment is based on" ] }, { "cell_type": "code", "execution_count": 2, "id": "30dffc16-7072-4a62-b606-642138990e5f", "metadata": {}, "outputs": [], "source": [ "#import set_path # Using MyBinder? Uncomment this before running the next cell!" ] }, { "cell_type": "code", "execution_count": 3, "id": "c3d91324-222b-4456-b13c-2759349c7c16", "metadata": {}, "outputs": [], "source": [ "#import sys\n", "#sys.path.append(\"C:/some_path/my_env_or_install\") # CHANGE to the folder containing your venv or libraries installation!\n", "# NOTE: If any of the imports below can't find a module, uncomment the lines above, or try: import set_path \n", "\n", "from life123 import BioSim1D, ChemData, CollectionArray, Numerical, check_version\n", "\n", "import numpy as np" ] }, { "cell_type": "code", "execution_count": 4, "id": "07e1388d-ebf9-4c03-ab0a-878d8b7d3152", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "OK\n" ] } ], "source": [ "check_version(LIFE123_VERSION)" ] }, { "cell_type": "code", "execution_count": null, "id": "73ba0473-38eb-4267-9e33-b2149244bb74", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": 5, "id": "0c291bcf-fd3a-443d-9df1-1b8ce9fefcae", "metadata": {}, "outputs": [], "source": [ "# We'll be considering just 1 chemical species, \"A\"\n", "diffusion_rate = 10.\n", "\n", "chem_data = ChemData(diffusion_rates=diffusion_rate, names=\"A\")" ] }, { "cell_type": "code", "execution_count": null, "id": "617cae4f-a3f0-4f4b-b794-4e3872b3c079", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "82504138-da32-486a-930d-bc9419cb1dc2", "metadata": {}, "source": [ "# BASELINE\n", "This will be our initial system, whose adherence to the diffusion equation we'll test.\n", "Afterwards, we'll tweak individual parameters - and observed their effect on the closeness of the approximation" ] }, { "cell_type": "code", "execution_count": 6, "id": "bca9a060-f62d-433a-9924-1fb2314870f3", "metadata": {}, "outputs": [], "source": [ "# Parameters of the simulation run (the diffusion rate got set earlier, and will never vary)\n", "delta_t = 0.01\n", "n_bins = 300\n", "delta_x = 2 # Note that the number of bins also define the fraction of the sine wave cycle in each bin\n", "algorithm = None # \"Explicit, with 3+1 stencil\"" ] }, { "cell_type": "code", "execution_count": 7, "id": "696122b0-c8c4-4d43-a4fc-319f95221602", "metadata": { "tags": [] }, "outputs": [], "source": [ "# Initialize the system\n", "bio = BioSim1D(n_bins=n_bins, chem_data=chem_data)\n", "\n", "# Initialize the concentrations to 2 superposed sine waves\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=1, amplitude=10, bias=50)\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=2, amplitude=8)" ] }, { "cell_type": "code", "execution_count": 8, "id": "4a4daf15-342b-45b2-8ee1-2324340b79d8", "metadata": {}, "outputs": [ { "data": { "text/html": [ " \n", " " ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "application/vnd.plotly.v1+json": { "config": { "plotlyServerURL": "https://plot.ly" }, "data": [ { "hovertemplate": "Bin number=%{x}
[A]=%{y}", "legendgroup": "", "line": { "color": "darkturquoise", "dash": "solid" }, "marker": { "symbol": "circle" }, "mode": "lines", "name": "", "orientation": "v", "showlegend": false, "type": "scatter", "x": [ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157, 158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169, 170, 171, 172, 173, 174, 175, 176, 177, 178, 179, 180, 181, 182, 183, 184, 185, 186, 187, 188, 189, 190, 191, 192, 193, 194, 195, 196, 197, 198, 199, 200, 201, 202, 203, 204, 205, 206, 207, 208, 209, 210, 211, 212, 213, 214, 215, 216, 217, 218, 219, 220, 221, 222, 223, 224, 225, 226, 227, 228, 229, 230, 231, 232, 233, 234, 235, 236, 237, 238, 239, 240, 241, 242, 243, 244, 245, 246, 247, 248, 249, 250, 251, 252, 253, 254, 255, 256, 257, 258, 259, 260, 261, 262, 263, 264, 265, 266, 267, 268, 269, 270, 271, 272, 273, 274, 275, 276, 277, 278, 279, 280, 281, 282, 283, 284, 285, 286, 287, 288, 289, 290, 291, 292, 293, 294, 295, 296, 297, 298, 299 ], "xaxis": "x", "y": [ 50, 50.544429428667165, 51.088179283950524, 51.630571063807565, 52.170928407051974, 52.70857815921861, 53.242851432961885, 53.77308466117989, 54.29862064106969, 54.81880956733467, 55.333010052784, 55.84059013458597, 56.34092826446227, 56.83341428113769, 57.317450363391465, 57.79245196208926, 58.25784870961231, 58.713085305138996, 59.157622374276286, 59.59093730158345, 60.01252503457715, 60.42189885785704, 60.81859113604284, 61.20215402426825, 61.57216014503327, 61.92820323027551, 62.26989872758102, 62.596884369518115, 62.90882070514144, 63.20539159277945, 63.48630465328596, 63.751291683005014, 64.00010902576838, 64.2325379033164, 64.44838470360496, 64.64748122653477, 64.82968488671305, 64.99487887293242, 65.1429722641268, 65.27390010164028, 65.38762341772012, 65.48412922022146, 65.5634304335874, 65.62556579624368, 65.67059971462322, 65.6986220741107, 65.70974800727217, 65.704117619809, 65.6818956747483, 65.64327123545478, 65.5884572681199, 65.51769020445475, 65.43122946538178, 65.3293569465877, 65.21237646686652, 65.08061318024517, 64.93441295294744, 64.77414170631202, 64.60018472683976, 64.41294594460129, 64.21284718129132, 64.00032736926742, 63.775841742961624, 63.53986100410003, 63.29287046221103, 63.03536915194446, 62.767868928764315, 62.49089354461371, 62.204977705185655, 61.910666110463616, 61.60851248022481, 61.29907856622386, 60.98293315279715, 60.66065104764711, 60.33281206458206, 60.00000000000001, 59.66280160491487, 59.32180555433006, 58.977601415768284, 58.63077861876623, 58.28192542714066, 57.93162791582594, 57.580468954074114, 57.229027196796395, 56.87787608580947, 56.527582862731656, 56.178707595252185, 55.83180221847259, 55.48740959299154, 55.146062581373926, 54.80828314461176, 54.47458146014833, 54.14545506299808, 53.82138801145301, 53.50285007882258, 53.19029597260686, 52.884164582453884, 52.58487825820056, 52.29284211924252, 52.00844339642252, 51.73205080756888, 51.46401396775535, 51.204662835292, 50.954307194392996, 50.71323617540191, 50.481717813388244, 50.2599986458607, 50.04830335027312, 49.846834421928385, 49.65577189281366, 49.47527309182776, 49.305472446787945, 49.146481328529276, 48.9983879373349, 48.86125723186071, 48.735130900642396, 48.62002737619715, 48.51594189165739, 48.42284657979783, 48.34069061424248, 48.269400392563504, 48.20887976090919, 48.159010279725436, 48.11965153006181, 48.09064145988159, 48.071796769724486, 48.062913337001035, 48.063766678128985, 48.07411244765571, 48.093686973444406, 48.12220782693885, 48.1593744274586, 48.20486867941727, 48.25835564129731, 48.31948422515935, 48.387887925409686, 48.46318557549796, 48.544982131167615, 48.63286947883482, 48.72642726762714, 48.82522376357119, 48.92881672437982, 49.03675429325235, 49.148575910068345, 49.263813238324204, 49.38199110613446, 49.502628459594334, 49.625239326778704, 49.74933379063347, 49.874418968999976, 50, 50.125581031000024, 50.25066620936653, 50.374760673221296, 50.497371540405666, 50.618008893865536, 50.736186761675796, 50.851424089931655, 50.96324570674764, 51.07118327562018, 51.17477623642881, 51.27357273237286, 51.36713052116517, 51.455017868832385, 51.53681442450203, 51.61211207459031, 51.68051577484064, 51.74164435870269, 51.79513132058273, 51.8406255725414, 51.87779217306115, 51.90631302655559, 51.92588755234429, 51.936233321871015, 51.937086662998965, 51.92820323027551, 51.90935854011841, 51.88034846993819, 51.840989720274564, 51.791120239090816, 51.7305996074365, 51.659309385757524, 51.57715342020218, 51.484058108342616, 51.37997262380286, 51.26486909935761, 51.13874276813929, 51.00161206266512, 50.85351867147074, 50.69452755321206, 50.52472690817225, 50.344228107186346, 50.15316557807162, 49.951696649726884, 49.74000135413931, 49.518282186611756, 49.2867638245981, 49.04569280560701, 48.795337164708016, 48.53598603224466, 48.26794919243113, 47.9915566035775, 47.70715788075749, 47.41512174179946, 47.11583541754612, 46.80970402739315, 46.497149921177424, 46.17861198854701, 45.85454493700193, 45.525418539851685, 45.19171685538825, 44.85393741862609, 44.51259040700847, 44.16819778152742, 43.82129240474783, 43.47241713726835, 43.122123914190546, 42.77097280320363, 42.41953104592591, 42.068372084174065, 41.71807457285935, 41.36922138123378, 41.02239858423173, 40.678194445669945, 40.337198395085146, 40, 39.66718793541796, 39.3393489523529, 39.017066847202855, 38.70092143377615, 38.3914875197752, 38.0893338895364, 37.79502229481435, 37.509106455386295, 37.23213107123569, 36.964630848055556, 36.70712953778898, 36.460138995899975, 36.22415825703838, 35.99967263073258, 35.787152818708684, 35.58705405539872, 35.39981527316025, 35.22585829368798, 35.065587047052574, 34.919386819754834, 34.7876235331335, 34.670643053412306, 34.568770534618224, 34.48230979554525, 34.41154273188011, 34.35672876454522, 34.31810432525169, 34.295882380191, 34.29025199272783, 34.3013779258893, 34.329400285376764, 34.374434203756316, 34.436569566412594, 34.515870779778524, 34.612376582279865, 34.726099898359706, 34.8570277358732, 35.00512112706758, 35.17031511328693, 35.35251877346523, 35.551615296395035, 35.76746209668359, 35.99989097423161, 36.24870831699498, 36.51369534671404, 36.794608407220544, 37.09117929485854, 37.40311563048187, 37.73010127241896, 38.07179676972448, 38.427839854966706, 38.79784597573175, 39.18140886395714, 39.57810114214295, 39.987474965422834, 40.409062698416534, 40.84237762572369, 41.286914694861004, 41.74215129038767, 42.20754803791071, 42.68254963660851, 43.16658571886228, 43.659071735537715, 44.159409865414, 44.666989947215995, 45.1811904326653, 45.701379358930296, 46.226915338820085, 46.75714856703811, 47.29142184078137, 47.82907159294801, 48.36942893619241, 48.911820716049476, 49.45557057133282 ], "yaxis": "y" } ], "layout": { "autosize": true, "legend": { "tracegroupgap": 0 }, "template": { "data": { "bar": [ { "error_x": { "color": "#2a3f5f" }, "error_y": { "color": "#2a3f5f" }, "marker": { "line": { "color": "#E5ECF6", "width": 0.5 }, "pattern": { "fillmode": "overlay", "size": 10, "solidity": 0.2 } }, "type": "bar" } ], "barpolar": [ { "marker": { "line": { "color": "#E5ECF6", "width": 0.5 }, "pattern": { "fillmode": "overlay", "size": 10, "solidity": 0.2 } }, "type": "barpolar" } ], "carpet": [ { "aaxis": { "endlinecolor": "#2a3f5f", "gridcolor": "white", "linecolor": "white", "minorgridcolor": "white", "startlinecolor": "#2a3f5f" }, "baxis": { "endlinecolor": "#2a3f5f", "gridcolor": "white", "linecolor": "white", "minorgridcolor": "white", "startlinecolor": "#2a3f5f" }, "type": "carpet" } ], "choropleth": [ { "colorbar": { "outlinewidth": 0, "ticks": "" }, "type": "choropleth" } ], "contour": [ { "colorbar": { "outlinewidth": 0, "ticks": "" }, "colorscale": [ [ 0, "#0d0887" ], [ 0.1111111111111111, "#46039f" ], [ 0.2222222222222222, "#7201a8" ], [ 0.3333333333333333, "#9c179e" ], [ 0.4444444444444444, "#bd3786" ], [ 0.5555555555555556, "#d8576b" ], [ 0.6666666666666666, "#ed7953" ], [ 0.7777777777777778, "#fb9f3a" ], [ 0.8888888888888888, "#fdca26" ], [ 1, "#f0f921" ] ], "type": "contour" } ], "contourcarpet": [ { "colorbar": { "outlinewidth": 0, "ticks": "" }, "type": "contourcarpet" } ], "heatmap": [ { "colorbar": { "outlinewidth": 0, "ticks": "" }, "colorscale": [ [ 0, "#0d0887" ], [ 0.1111111111111111, "#46039f" ], [ 0.2222222222222222, "#7201a8" ], [ 0.3333333333333333, "#9c179e" ], [ 0.4444444444444444, "#bd3786" ], [ 0.5555555555555556, "#d8576b" ], [ 0.6666666666666666, "#ed7953" ], [ 0.7777777777777778, "#fb9f3a" ], [ 0.8888888888888888, "#fdca26" ], [ 1, "#f0f921" ] ], "type": "heatmap" } ], "heatmapgl": [ { "colorbar": { "outlinewidth": 0, "ticks": "" }, "colorscale": [ [ 0, "#0d0887" ], [ 0.1111111111111111, "#46039f" ], [ 0.2222222222222222, "#7201a8" ], [ 0.3333333333333333, "#9c179e" ], [ 0.4444444444444444, "#bd3786" ], [ 0.5555555555555556, "#d8576b" ], [ 0.6666666666666666, "#ed7953" ], [ 0.7777777777777778, "#fb9f3a" ], [ 0.8888888888888888, "#fdca26" ], [ 1, "#f0f921" ] ], "type": "heatmapgl" } ], "histogram": [ { "marker": { "pattern": { "fillmode": "overlay", "size": 10, "solidity": 0.2 } }, "type": "histogram" } ], "histogram2d": [ { "colorbar": { "outlinewidth": 0, "ticks": "" }, "colorscale": [ [ 0, "#0d0887" ], [ 0.1111111111111111, "#46039f" ], [ 0.2222222222222222, "#7201a8" ], [ 0.3333333333333333, "#9c179e" ], [ 0.4444444444444444, "#bd3786" ], [ 0.5555555555555556, "#d8576b" ], [ 0.6666666666666666, "#ed7953" ], [ 0.7777777777777778, "#fb9f3a" ], [ 0.8888888888888888, "#fdca26" ], [ 1, "#f0f921" ] ], "type": "histogram2d" } ], "histogram2dcontour": [ { "colorbar": { "outlinewidth": 0, "ticks": "" }, "colorscale": [ [ 0, "#0d0887" ], [ 0.1111111111111111, "#46039f" ], [ 0.2222222222222222, "#7201a8" ], [ 0.3333333333333333, "#9c179e" ], [ 0.4444444444444444, "#bd3786" ], [ 0.5555555555555556, "#d8576b" ], [ 0.6666666666666666, "#ed7953" ], [ 0.7777777777777778, "#fb9f3a" ], [ 0.8888888888888888, "#fdca26" ], [ 1, "#f0f921" ] ], "type": "histogram2dcontour" } ], "mesh3d": [ { "colorbar": { "outlinewidth": 0, "ticks": "" }, "type": "mesh3d" } ], "parcoords": [ { "line": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "type": "parcoords" } ], "pie": [ { "automargin": true, "type": "pie" } ], "scatter": [ { "fillpattern": { "fillmode": "overlay", "size": 10, "solidity": 0.2 }, "type": "scatter" } ], "scatter3d": [ { "line": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "marker": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "type": "scatter3d" } ], "scattercarpet": [ { "marker": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "type": "scattercarpet" } ], "scattergeo": [ { "marker": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "type": "scattergeo" } ], "scattergl": [ { "marker": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "type": "scattergl" } ], "scattermapbox": [ { "marker": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "type": "scattermapbox" } ], "scatterpolar": [ { "marker": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "type": "scatterpolar" } ], "scatterpolargl": [ { "marker": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "type": "scatterpolargl" } ], "scatterternary": [ { "marker": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "type": "scatterternary" } ], "surface": [ { "colorbar": { "outlinewidth": 0, "ticks": "" }, "colorscale": [ [ 0, "#0d0887" ], [ 0.1111111111111111, "#46039f" ], [ 0.2222222222222222, "#7201a8" ], [ 0.3333333333333333, "#9c179e" ], [ 0.4444444444444444, "#bd3786" ], [ 0.5555555555555556, "#d8576b" ], [ 0.6666666666666666, "#ed7953" ], [ 0.7777777777777778, "#fb9f3a" ], [ 0.8888888888888888, "#fdca26" ], [ 1, "#f0f921" ] ], "type": "surface" } ], "table": [ { "cells": { "fill": { "color": "#EBF0F8" }, "line": { "color": "white" } }, "header": { "fill": { "color": "#C8D4E3" }, "line": { "color": "white" } }, "type": "table" } ] }, "layout": { "annotationdefaults": { "arrowcolor": "#2a3f5f", "arrowhead": 0, "arrowwidth": 1 }, "autotypenumbers": "strict", "coloraxis": { "colorbar": { "outlinewidth": 0, "ticks": "" } }, "colorscale": { "diverging": [ [ 0, "#8e0152" ], [ 0.1, "#c51b7d" ], [ 0.2, "#de77ae" ], [ 0.3, "#f1b6da" ], [ 0.4, "#fde0ef" ], [ 0.5, "#f7f7f7" ], [ 0.6, "#e6f5d0" ], [ 0.7, "#b8e186" ], [ 0.8, "#7fbc41" ], [ 0.9, "#4d9221" ], [ 1, "#276419" ] ], "sequential": [ [ 0, "#0d0887" ], [ 0.1111111111111111, "#46039f" ], [ 0.2222222222222222, "#7201a8" ], [ 0.3333333333333333, "#9c179e" ], [ 0.4444444444444444, "#bd3786" ], [ 0.5555555555555556, "#d8576b" ], [ 0.6666666666666666, "#ed7953" ], [ 0.7777777777777778, "#fb9f3a" ], [ 0.8888888888888888, "#fdca26" ], [ 1, "#f0f921" ] ], "sequentialminus": [ [ 0, "#0d0887" ], [ 0.1111111111111111, "#46039f" ], [ 0.2222222222222222, "#7201a8" ], [ 0.3333333333333333, "#9c179e" ], [ 0.4444444444444444, "#bd3786" ], [ 0.5555555555555556, "#d8576b" ], [ 0.6666666666666666, "#ed7953" ], [ 0.7777777777777778, "#fb9f3a" ], [ 0.8888888888888888, "#fdca26" ], [ 1, "#f0f921" ] ] }, "colorway": [ "#636efa", "#EF553B", "#00cc96", "#ab63fa", "#FFA15A", "#19d3f3", "#FF6692", "#B6E880", "#FF97FF", "#FECB52" ], "font": { "color": "#2a3f5f" }, "geo": { "bgcolor": "white", "lakecolor": "white", "landcolor": "#E5ECF6", "showlakes": true, "showland": true, "subunitcolor": "white" }, "hoverlabel": { "align": "left" }, "hovermode": "closest", "mapbox": { "style": "light" }, "paper_bgcolor": "white", "plot_bgcolor": "#E5ECF6", "polar": { "angularaxis": { "gridcolor": "white", "linecolor": "white", "ticks": "" }, "bgcolor": "#E5ECF6", "radialaxis": { "gridcolor": "white", "linecolor": "white", "ticks": "" } }, "scene": { "xaxis": { "backgroundcolor": "#E5ECF6", "gridcolor": "white", "gridwidth": 2, "linecolor": "white", "showbackground": true, "ticks": "", "zerolinecolor": "white" }, "yaxis": { "backgroundcolor": "#E5ECF6", "gridcolor": "white", "gridwidth": 2, "linecolor": "white", "showbackground": true, "ticks": "", "zerolinecolor": "white" }, "zaxis": { "backgroundcolor": "#E5ECF6", "gridcolor": "white", "gridwidth": 2, "linecolor": "white", "showbackground": true, "ticks": "", "zerolinecolor": "white" } }, "shapedefaults": { "line": { "color": "#2a3f5f" } }, "ternary": { "aaxis": { "gridcolor": "white", "linecolor": "white", "ticks": "" }, "baxis": { "gridcolor": "white", "linecolor": "white", "ticks": "" }, "bgcolor": "#E5ECF6", "caxis": { "gridcolor": "white", "linecolor": "white", "ticks": "" } }, "title": { "x": 0.05 }, "xaxis": { "automargin": true, "gridcolor": "white", "linecolor": "white", "ticks": "", "title": { "standoff": 15 }, "zerolinecolor": "white", "zerolinewidth": 2 }, "yaxis": { "automargin": true, "gridcolor": "white", "linecolor": "white", "ticks": "", "title": { "standoff": 15 }, "zerolinecolor": "white", "zerolinewidth": 2 } } }, "title": { "text": "Initial System State
System snapshot at time t=0" }, "xaxis": { "anchor": "y", "autorange": true, "domain": [ 0, 1 ], "range": [ 0, 299 ], "title": { "text": "Bin number" }, "type": "linear" }, "yaxis": { "anchor": "x", "autorange": true, "domain": [ 0, 1 ], "range": [ 32.544724436364255, 67.45527556363574 ], "title": { "text": "[A]" }, "type": "linear" } } }, "image/png": "iVBORw0KGgoAAAANSUhEUgAAA+cAAAFoCAYAAAA1h4qYAAAgAElEQVR4Xu2dCZxdVX34f5l9JpnJLNkggEVFVEAQI9aCgrX611YRcQnIYmpVUNlXBWpMUf8CsipFbKsIIkSFItVC6d8KFqtiWBSoggsKBJJMZs9k9sn/d15yh5ub92bue+e89+6953s+mc9MZu499/y+vzPL951t3jYtQoEABCAAAQhAAAIQgAAEIAABCECgagTmIedVY8+DIQABCEAAAhCAAAQgAAEIQAACOQLIOR0BAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeZUTwOMhAAEIQAACEIAABCAAAQhAAALIOX0AAhCAAAQgAAEIQAACEIAABCBQZQLIeRUTcNl1t8oNa++Wu26+RPZavrSKLeHREIAABCAAAQhAAAIQgAAEIFBNApmS88OPPl029w7sxLMU8Q3qefzeG2bq2u+IVfLWw1fIlWtO2an+W+74oXz2qptk1cq3ybkfO6aoXJYi56Yd0bL/vnvL2utXF/XsuBeX0sa4dbu67v4HHpWTzrt8l+qiOXERi4s6XMVNPRCAAAQgAAEIQAACEIBAdghkSs5NWgJRK0WWg7QmUc6DFwGiLxAE8S7qXCj33X61856ZdBk9c/WX5Z771slFZ5wgxx715pn4g3aHebmIxUUdzpNEhRCAAAQgAAEIQAACEIBA6gkg5zFTWGjkPObteS8rRvTM82cbITeSGh3Vt2lbcG8xbXTxvGLqiPNCTJiLi1hc1FFMjFwLAQhAAAIQgAAEIAABCPhBwBs5D0vVCad+fqfp79Gp7ytPWiMbuntnRqLzTSUPRqoDQQyP3Aaj3NEuVGh0N87U+7gvDgTPjj7LtCWfzOZra3BvwGy2OJ5ev1Heftz5O11y/aVny2GHHDDzuYDn5z754Z2mnwcvNkTbkG/5QL5vx3zsC33bzhVLMAIfvT+cm7nqMPfm4xleHuHHjxWihAAEIAABCEAAAhCAAASKJeCdnBtAYeGKirj5er7PFZLjQnJ+x13377QOPJ80FzMKG0y1j4pvvoSba02JTnOPxhU8P1xnIKkBo9namG+9fb44zXMfe+KpnUb+w1IfnhEQfD7OsoRwHXEEeLZYTNy7L1u0074BURaG6Vx1mCn2s/Es9huU6yEAAQhAAAIQgAAEIAABPwh4J+fRUepAJsNCZSvnhbqOkeaDD9hnZvp5MXKeb4TaPCefxOaLKbg/PKKeL05Tp7n/0Nfun9tBfrY2FnrBwkjtQ4/+dubFgULPMTyWLe7cZTO7Qp/Px7XQLIV8L2IUw9s8Kx+zQnXMNoofzbsfP1qIEgIQgAAEIAABCEAAAhAohgByvmO39XLI+Vw7qxcri0Fig5HocKKjU8Gj4hw8KzzCHIwMz7aWvRQZjT6rkJwX+/m5OnY+3uEXI+biHWfX90J1zFa3idOUcu2oPxcXvg4BCEAAAhCAAAQgAAEIJJ8Acl4GOQ+kNyrM0RHhuWQxbvcJZD38AkO0biOu+UbZ8621Dl9XqI2FRqzDbQ5mKRQr4YWuj8sjuC6Q9eAFidl4B8sGwjKfb4p9oToKrVkP2lKu3fSLZcL1EIAABCAAAQhAAAIQgEAyCSDnZZDzQtO9yyXnhTaBC4R8j90W585ij7PxXHRt+1wj53HWwFdLzuOuny80Jb0YOc83MyGZ3/K0CgIQgAAEIAABCEAAAhBIIgHkPKacF1o3HBW72TY0K1XOTZ1r7/zRTpuVhTvTbKO5Zu23Wde9+7KuXY5ay3f8WjSeQPzziX3cHeTLJeembXvuvmSnneHDXKLPLRRLvjX6pp58uSxURzE7xyfxBwFtggAEIAABCEAAAhCAAASqSwA5jynnwShsdFfwfFKWb7f0YEQ6vL477rT28GZw0ZHqQBbzHT8Wvi+fXAdtyrcOPbh+thcb8u3WHkitOa4u2C2+nHJuZgSYEs1LkK+5pqmHJTzMMMwuPM1/Nh7BM/Pl6IGHf12Wc+ir++ODp0MAAhCAAAQgAAEIQAACrghkSs4D2QzDmetIsLi7tZs6w/XPds559Frzf9MOI6zh3cnjynkQT774grrNzur5ylw7n+fbRC0qutH15fmEN/rssNCWS87DYh19fqE13oViie6Gb+6/6UsX5M5wj67Vn41HMbvHu/omph4IQAACEIAABCAAAQhAIP0EMiXn6U+H2wgKTdd2+xRqgwAEIAABCEAAAhCAAAQgAAFbAsi5LcEE359ven2Cm0vTIAABCEAAAhCAAAQgAAEIeEsAOc9o6gvt4J7RcAkLAhCAAAQgAAEIQAACEIBAqgkg56lOH42HAAQgAAEIQAACEIAABCAAgSwQQM6zkEVigAAEIAABCEAAAhCAAAQgAIFUE0DOU50+Gg8BCEAAAhCAAAQgAAEIQAACWSCAnGchi8QAAQhAAAIQgAAEIAABCEAAAqkmgJynOn00HgIQgAAEIAABCEAAAhCAAASyQAA5z0IWiQECEIAABCAAAQhAAAIQgAAEUk0AOU91+vxp/MqT1siG7l657/ar/QmaSCEAAQhAAAIQgAAEIAABbwhkUs6NyD32xFO7JPHxe28oS2Lvf+BROem8y+WiM06QY496c1me4Xul5ZDzy667VW5Ye7fcdfMlstfypU4Rl9Le/Y5YJW89fIVcueaUndoSnFm/auXb5NyPHeO0nTaVFWqvTZ3m3iAvQT2LOhfyoowtVO6HAAQgAAEIQAACEEg8gUzJ+dPrN8rbjztf8v0xb0TClHKIGHJe/n5eiuzO1SrkfC5Cs3+9HHJ+5uovyz33rdvp+7QcubeLnLshAAEIQAACEIAABCDgnkCm5DwYMS80Qm5kbOWRb3I+Soqcu++Y0RrLIWhpkfPy0y3tCa7lPHhxLTpDIPg8M1NKyxN3QQACEIAABCAAAQikg0Dm5DzOuuTgj/18U4hN2ox07L/v3rL2+tW5LAbXh1Ma3BuIeTTd0boPP/p02dw7MHNZVEDComhG/8PFvNgQbUMxU31na3/wnEC0Dnn1K+SzV9008/hoHPnqMhfnYxm3TnN/MGIajjs8yyGQ85u+dEFudkRQCnGYa2p09OtBfXMJYDSP5r5oGwotq5ht1kYwsyMcf1Bvvhd/Ah6f++SHc0sqghL022AqfPD5Qn09ep25Ps7yj9naW+qPvqAt1196thx2yAE7VWO4L1vcOfM9WeozuA8CEIAABCAAAQhAAAJJJZApOQ8EL87a3HzTZ02SoqOpgRiF6wwkIhC52UbOA5nNJ/thYQrLYliOwjIY/XwcWYnTfhN3IFvhNuWLK4gn3JZCL3bMVmeYZ74ZDyZuU4IN4IJroiKcb/Q2X335PlfKyLlplxHisDxG22raXcpIf6GR6EJybvZVKPQiUr7PR78vgu+BsAwX+r7I9wNstpHzQi9QROsJP3u2fJTCM6k/dGkXBCAAAQhAAAIQgAAE8hHIlJybAPONbIZFJYBQaAqtuf/gA/aZ2ZSrkDCY+3/yi8dyG8DNJueFpCIQ/GA0tdBzAlmKjmYW+nw0yXHaH8h5vtHVuCOWpj0PPfrbnTbuKiRv0bYXus60PdgArRDHaF2z5SL6nFLkPN83UTSXlZLzfLNECuUr+vnZOEW/Bwr96HQ9rX22FwaQc36BQQACEIAABCAAAQhknUDm5DxIWL5py/mmH4cFJxCW8GheeNpvoem+xQhh0L7oswqJYrGfj3bYOO0vRc7j7Ig/1+7jAefgBZXZZjzElfNiRl9LlfNCU+LD/aYUmSx25DyfnBd6bvTzc3EyfSJY1oGcZ/3XAPFBAAIQgAAEIAABCCSBQGblPAo3kJGwAEal2ghMPinJty43PBpfSM4Lrc8Oty2YGl+shBcjlnO1vxg5D+qKvtCRbyQ/rpyb5+eb8RBHdqPPnW1GQXRqezEMTRvD+QyvH8+3Vjrpcp5vjX+4X8bZ04Bp7Un4EU4bIAABCEAAAhCAAASyQsAbOS+0LjoQ8ss+fXJuo7G5NgQziY+ubZ9r5DzOGvhyynm0s+Zbm19ItKLToeOOYM8m/HNJcViCg9kKcZ8714hweLR5rnYUeoEnurFbGuU8iD3O5m+Ffti5ntbOhnBZ+bVCHBCAAAQgAAEIQAACpRDIlJwb6bxyzSl5ORQS6EAIzHrr6JrpQMTz1RkWk9mOeoo7glouOS/EJCpWceW80JrmYkbOw0wMu7V3/mhmbXmQvCiPuHI+1wsl4XX1+daKz/ZNVGhUPp9U5luDP9c3aKG13rPt1h5smBfUHXdau4vj/+KuTZ8r7uDrHKUWlxTXQQACEIAABCAAAQhkkUCm5LzQ2uXgj/5CU3WD+/KNmufb0Tq6W7vpGIWkNd9u7UFHMveYo8H2Wr50l13iC0nqXJ+PdtK47Y8r5/k27QpPkQ6PxEaPpDNtC6Q7mLJeiE++3drzrbHOJ8zB9PXoUWxmd/N8u8zHmdlg2p4v7+ElA/n2Ksh3LFihHySF5L8cch688HTPfesk2kYT0wMP/7rgC11B++NuSljMD858/SvuC1zFPIdrIQABCEAAAhCAAAQgkDQCmZLzsHBEQc82XX0uyci3+VlUaKLry+c659y0L7xuvVwj5+Y5cdofV87zMTZsjcwZ0YvKuXlBJHy+u7k/Oi0831nx0R32446cR1+8CP5f6IWZ6Hr8uZY15Ds/PDgbPtonouu6ZzvnPGhneO19nHPOSx05D56Xbz8C87W4Lyrka6/tD7m5zqi3rZ/7IQABCEAAAhCAAAQgkEQCmZPzYiEXmkpbbD1cvysB12uSYQwBCEAAAhCAAAQgAAEIQCCrBLyX89nOVs5q0isVF3JeKdI8BwIQgAAEIAABCEAAAhBIOwGv5bzQDu5pT2pS2o+cJyUTtAMCEIAABCAAAQhAAAIQSDoBr+U86cmhfRCAAAQgAAEIQAACEIAABCDgBwHk3I88EyUEIAABCEAAAhCAAAQgAAEIJJgAcp7g5NA0CEAAAhCAAAQgAAEIQAACEPCDAHLuR56JEgIQgAAEIAABCEAAAhCAAAQSTAA5T3ByaBoEIAABCEAAAhCAAAQgAAEI+EEAOfcjz0QJAQhAAAIQgAAEIAABCEAAAgkmgJwnODk0DQIQgAAEIAABCEAAAhCAAAT8IICc+5FnooQABCAAAQhAAAIQgAAEIACBBBNAzhOcHJoGAQhAAAIQgAAEIAABCEAAAn4QQM79yDNRQgACEIAABCAAAQhAAAIQgECCCSDnCU4OTYMABCAAAQhAAAIQgAAEIAABPwgg537kmSghAAEIQAACEIAABCAAAQhAIMEEkPMEJ4emQQACEIAABCAAAQhAAAIQgIAfBJBzP/JMlBCAAAQgAAEIQAACEIAABCCQYALIeYKTQ9MgAAEIQAACEIAABCAAAQhAwA8CyLkfeSZKCEAAAhCAAAQgAAEIQAACEEgwAeQ8wcmhaRCAAAQgAAEIQAACEIAABCDgBwHk3I88EyUEIAABCEAAAhCAAAQgAAEIJJgAcp7g5NA0CEAAAhCAAAQgAAEIQAACEPCDAHLuR56JEgIQgAAEIAABCEAAAhCAAAQSTAA5T3ByaBoEIAABCEAAAhCAAAQgAAEI+EEAOfcjz0QJAQhAAAIQgAAEIAABCEAAAgkmgJwnODk0DQIQgAAEIAABCEAAAhCAAAT8IICc+5FnooQABCAAAQhAAAIQgAAEIACBBBNAzhOcHJoGAQhAAAIQgAAEIAABCEAAAn4QQM4d5Pm5nhEHtVBF1gnU186T9gUN0j0wlvVQic8RgWWdzbKpb0SmtzmqkGoyTcD8fBmfmJKtY1OZjpPg3BBoaqiVlsZa6R0ad1MhtWSewO5dzcLfvJlPs3WApp9QSieAnJfObuZOflA5gOhBFci5B0l2HCJy7hhoxqtDzjOeYMfhIeeOgXpQHXLuQZIdhIic20FEzu345e5Gzh1A9KAK5NyDJDsOETl3DDTj1SHnGU+w4/CQc8dAPagOOfcgyQ5CRM7tICLndvyQcwf8fKkCOfcl0+7iRM7dsfShJuTchyy7ixE5d8fSl5qQc18ybRcncm7HDzm344ecO+DnSxXIuS+Zdhcncu6OpQ81Iec+ZNldjMi5O5a+1ISc+5JpuziRczt+yLkdP+TcAT9fqkDOfcm0uziRc3csfagJOfchy+5iRM7dsfSlJuTcl0zbxYmc2/FDzu34IecO+PlSBXLuS6bdxYmcu2PpQ03IuQ9Zdhcjcu6OpS81Iee+ZNouTuTcjh9ybscPOXfAz5cqkHNfMu0uTuTcHUsfakLOfciyuxiRc3csfakJOfcl03ZxIud2/JBzO37IuQN+SapiVLbJb8bHpHtqSgamp2Vg25T0T+rbjo8Hwp+f0q/r58f0nqC0zKuRRTW1srS2Vjr1zXy82LyvrZNlDfXy4vlNUjsynfvcwpqaJIVOWxJIADlPYFIS3CTkPMHJSWDTkPMEJiXhTULOE56ghDQPObdLBHJuxw85d8CvmlWsGx+VB8dG5Bdjo/L42Jg8PTVRseY0yzx5ZWOjvKyuQfZvaJLXNDXJAfWNFXs+D0o+AeQ8+TlKUguR8yRlI/ltQc6Tn6OktRA5T1pGktke5NwuL8i5HT/k3AG/SlXxjIr3Izoq/sDIiDw4PiK/1I/zlZerLO9WVyftOrK9UEe+zQh3e26kuy73cYe+te34vPn/fB0tD4oZSd88PSWbpyalW993T05Kj462b9T/9+oo/KBMy7Nj+jX9/0hoxD24v0mF/dWNKur6tqKxOfdmnkfxkwBy7mfeS40aOS+VnJ/3Ied+5t0mauTchp4/9yLndrlGzu34IecO+JWril+ogP9idExHxrfq22hOlqNlr9p6OdjIsI5av6ahWQ5sKN/IdXTN+UaV9icmxuWxmdH7EelRuY+Wvevq5SBt18FNKuvaxleVsY3lygX1lkYAOS+Nm693Iee+Zr60uJHz0rj5fBdy7nP248eOnMdnle9K5NyOH3LugJ+rKrZs2ya3Dw/KPVu3yM9Vxrdu21l0G3VU+qDciHSTCnmzvFbfd+kIeKVKnA3hzLT6B3V6/UOjhUf3zbr212nb39DcIm9pni8v1pF+SjYJIOfZzGu5okLOy0U2m/Ui59nMazmjQs7LSTc7dSPndrlEzu34IecO+NlW8YCOkH9raEDuHN6y0+ZsnSqxhzW1yGubm8s+Kh4nhjhyHq3HbDb3iL7QYITdzAB4QKW9N/KigxlVP2bBQjmqpU1aa+bFaQrXpIQAcp6SRCWkmch5QhKRkmYg5ylJVIKaiZwnKBkJbgpybpcc5NyOH3LugF8pVQxNb5PvDg/IjVv65cmJFzZxe1l9vbxv/kI5QkeUX1mfrBHlUuQ8H5tf61T4/x4dlv8aGdb3Iztd8s7mBbKydaG8SV+UoKSfAHKe/hxWMgLkvJK00/8s5Dz9Oax0BMh5pYmn83nIuV3ekHM7fsi5A37FVGFGyb+po+TfD42St+losRk1NqPH5VwzXkw7813rSs7DdZvN5769ZUi+Pdwvvw29SLFYp+sfPb9NjlUm++gLFpR0EkDO05m3arUaOa8W+XQ+FzlPZ96q2WrkvJr00/Ns5NwuV8i5HT/k3AG/uaoY1E3Svq1ryW8OjZKb/cvf2Ngi729tk7e3LJAGXU+e9FIOOQ/H/IhuLPftLYNyx9ah3PnrQTHHsxlO79EXMDhbPem9ZOf2Iefpyle1W4ucVzsD6Xo+cp6ufCWhtch5ErKQ/DYg53Y5Qs7t+CHnDvgVquJJnb593WCvivnQzCV/pjuXr9TR4PfPb5VltXVlfLr7qsst5+EW/2BkS24d/r2jW3cKxGwg9wHl91Z9T0k+AeQ8+TlKUguR8yRlI/ltQc6Tn6OktRA5T1pGktke5NwuL8i5HT/k3AG/aBW/mRyXz/d2yw9DYrmypVWOaVsoh+hRYmktlZTzgJE5Pu47Opr+HV2fH16bb6a9f7CtXU5UUa/kjvVpzV212o2cV4t8Op+LnKczb9VqNXJeLfLpfS5ynt7cVbLlyLkdbeTcjh9y7oBfUMWGqUn5Qv9mFcntI+UdNTXyQZXHVW0dYmQy7aUach5m9rBOezeibo6bG9Jj54Jy/II2Oa2tS5bXpWsmQtr7Q5z2I+dxKHFNQAA5py8UQwA5L4YW1xoCyDn9IA4B5DwOpcLXIOd2/JBzB/yMKH5poEeuHezL1WbOI/+IjuqeosKYpaPBqi3n4VT9m057v01F/T91x/egnKAvhJy+sFN2S9lyAQddMLFVIOeJTU0iG4acJzItiW0Ucp7Y1CS2Ych5YlOTqIYh53bpQM7t+CHnFvwm9V5zFNqVfT25s7vNuK1ZD31We1cmRsqjaJIk50HbzGyFqwZ65VYdTZ/QF0nq582TD+gu76eppKdtTb9FV0zsrch5YlOTyIYh54lMS2IbhZwnNjWJbRhyntjUJKphyLldOpBzO37IeYn8/n3rFvmcTmH/4+REbp/1I3VN+fkq5S/SDd+yWpIo5wHr9ZOTcrXOXlirkm5eNDG7339Ap7ufvrBLltSmf0lBWvsUcp7WzFWn3ch5dbin9anIeVozV712I+fVY5+mJyPndtlCzu34IedF8jPrni/s3SS/HB/L3Xl4U7Nc2LFE9qtvKLKm9F2eZDkPS/qVKunfCUn68Tqb4VQdSUfSK9/nkPPKM0/zE5HzNGev8m1HzivPPO1PRM7TnsHKtB85t+OMnNvxQ85j8vuTjpBfrCPld+mIuSkHNjTKpzsWy583pnf39Zihz1yWBjkPGvv01IRcpcsNbtMz081IutkHwEj6ae2dsigDm/MVm7tqXY+cV4t8Op+LnKczb9VqNXJeLfLpfS5ynt7cVbLlyLkdbeTcjh9yPge/vulpuXxgs9ykZ24byXuJTls/v2OR/E3zAgfk01VFmuQ8LOmX9/fIv+oO+lP6ySaV9BNaF8opOpKOpJe//yHn5WecpScg51nKZvljQc7LzzhrT0DOs5bR8sSDnNtxRc7t+CHns/D75pYB+Wxfd+7YrqW6bvnc9kVyrG425mtJo5wHuXpKZz5cppL+PR1JN8WMpJvd3c1IOuekl69HI+flY5vFmpHzLGa1fDEh5+Vjm9WakfOsZtZtXMi5HU/k3I4fcp6H32MTY3Jez8aZdeXnqZSfrEejGaHzuaRZzoO8/V4l/QqV9Dt2SLoZSV+lI+mf0I3jOvVceopbAsi5W55Zrw05z3qG3caHnLvl6UNtyLkPWbaPETm3Y4ic2/FDziP8rtazyi/VteWmHKrryS/rWprpHdiL6T5ZkPOwpF/Wt1nMeemmtMyrUUHvkI+0tst8/ZjihgBy7oajL7Ug575k2k2cyLkbjj7Vgpz7lO3SY0XOS2dn7kTO7fgh5zv4mQ3EPtG9QR7S3dhN+ULnkty0Z8oLBLIk50FUT06MyyX6YszdI8O5T3WqmJs9BczmcRR7Asi5PUOfakDOfcq2fazIuT1D32pAzn3LeGnxIuelcQvuQs7t+CHnSuAburb8H3q7ZVS25Y5Eu27x7rmN3yg7E8iinAcR/mJ8RC7WkfQHx7a/OPOq+kb5gs6aMLvyU0ongJyXzs7HO5FzH7NeeszIeensfL0TOfc188XFjZwXxyt6NXJux89rOe+enpJTNj8v94+O5FaTf7S1Qy5o75K6eX6vLS/UpbIs50HMt+uu7v+gmwCavmEmt5sNAC/SI/PaWI9e0k8a5LwkbN7ehJx7m/qSAkfOS8Lm9U3Iudfpjx08ch4bVd4LkXM7ft7KuTmv/KzeDTI4vU0W67nX1y/eTV7n0ZnlpXQbH+TccBneNi2X6BnpN2zpzx2/1qFifoEKuhF1XrYprucg58Xx8v1q5Nz3HlBc/Mh5cby4WgQ5pxfEIYCcx6FU+Brk3I6fd3JujkX7VM8G+VeVc1Pe0jxfrulaxshojH7ki5wHKJ6cmJBzta+s27EPwQE61f2yRUvFvKfEI4Ccx+PEVW3aPYsAACAASURBVNsJIOf0hGIIIOfF0OJaQwA5px/EIYCcx6GEnM9J6ZY7fiifveqm3HWLOhfKfbdfPXPPmau/LPfcty73//333VvWXr96p/qe6xmZs/4sXPDzsRH5uE5j3zA1Jc06BrpGN307boG/55YXm1Pf5Dzg812d6m7Ouw+mun9A+8yF7Ux1j9N/kPM4lLgmIICc0xeKIYCcF0OLa5Fz+kBcAsh5XFL5r2PkXLlcdt2tcsPau+Xxe2/YhZKR9q/ceOeMrK88aY2sOGhfOfdjx8xcm3U5H9ON3j7Xu1m+ptOUt2nUL69rkH9esrvszaZvRX33+SrnBpKZcXGpbhj3jdBU9wt1qvsxTHWftQ8h50V9i3l/MXLufRcoCgByXhQuLlYCjJzTDeIQQM7jUCp8DXKubPY7YpXcdfMlstfypbuQisp4VNbNDVmW88cmxuTk7uflqcmJ3Hrhk/Uc60+2L2LTtxK+73yW8wDXr/XotXN0qvsj42O5Tx2ku7mbXd2Z6p6/QyHnJXyjeXwLcu5x8ksIHTkvAZrntyDnnneAmOEj5zFBFbjMezm//4FH5cIv/LNs7h2YQfTWw1fIlWtOyf3/8KNPl5NPPFKOPerNuf+b60867/KdRtk39G4/Pipr5ZqBXvm/eoa1KWbTt+uW7CaHNrZkLcyKxVNXO08Wzq+XnsHxij0zqQ+6RY/fM1Pde6enc01cpeeif4pd3XdJ15KOJunu10MKzZQVCgRmIzBvmyxsaZDxySkZGdv+fUWBwGwEGhtqpKWhVvq2TAAKArEILOtskqz+zRsLABfFImD6CaV0At7LebDWPDyl3Yykr1r5ttzUdfPxRWecsIuch0fapzP2l/OgCtMH//i03Dk4lOtZ717YJv+81x7SXltbek/jztzMg3l6zFzW+kupqTX97IL1z8t1Pb25Kjq1f31h+W7yd50dpVaZuftqtL9s058vuHnmUus8oMnJbVJfrwcYamehxzjHm8kK55nfSvrP/IyhQCAOAfM7ib9h4pDy+xrTTyilE0DOI2vKDUqzAZwpZvQ8zsh5lqa1/2FyXE7Y9Jz8Uaexm/J53fTtgzqqSbEnwLT2/Awf16nu5/dslId37Op+oE51v0xPANivvsEeesprYFp7yhNY4eYzrb3CwFP+OKa1pzyBVWg+09qrAD2Fj2Rau13SvJfzfNPUw3Lu05rzn4yOyN91r89t3mXE6PrFbPpm9+21893I+ew0bxkelM/36lR3PSfdFPOi0Pm6v8FCPSfd14Kc+5r50uJGzkvj5utdyLmvmS89buS8dHY+3Ymc22Xbezk3+Mzo+Dve8vrcNPan12+Utx93vlx/6dly2CEHiC+7td+8ZVDO692Y603vnd8qV+vIJcUtAeR8bp4DOtX9Et3n4Bu6Jt2Uznk1ckHnYjlWd3X3sSDnPma99JiR89LZ+Xgncu5j1u1iRs7t+PlyN3Jul2nkfAc/s7Y8KOE15uZzWT/n/DO6Mdc/DfXnwjc7sZ/axppfu2+r/Hcj5/GpmlMCztOp7r/csav7qxuadKr7UnmFZ1PdkfP4fYYrRZBzekExBJDzYmhxrSGAnNMP4hBAzuNQKnwNcm7HL3d3WtecD+v04ZP0mLQfjW7NxXH94t3kHc0LHBChinwEkPPi+4WZ6v45ffGob8eu7hfqi0cf9+jFI+S8+D7j8x3Iuc/ZLz525Lx4Zr7fgZz73gPixY+cx+NU6Crk3I5fauX8ualJOXHjevm1bgC3SI9Ju2nJcnmVbsRFKR8B5Lw0tv0q5kbQv6WibsrrG5vlmkXLZPfautIqTNFdyHmKkpWApiLnCUhCipqAnKcoWQlpKnKekEQkvBnIuV2CkHM7fqmU81/pVOHjNj6b23jrFXUN8s2ly2WZB6LjINVWVSDnVvjkEd3N/bTNG+T3epJAi65FP6e9Uz7c2iFZPuAPObfrM77djZz7lnG7eJFzO34+3o2c+5j14mNGzotnFr4DObfjlzo5v3PrkHxMBceUv2xqka/oVPb5KjqU8hNAzt0wvnawTz6vm8aZYl5c+uKipXKQrknPYkHOs5jV8sWEnJePbRZrRs6zmNXyxoScl5dvVmpHzu0yiZzb8UuVnF812Ctf7O+Rbdrqj7S2y6c7Fgta7qADxKwCOY8JKsZlT09NyNn6ItP/jI3m+vCJeuzaBR2LMvdCE3IeozNwyQwB5JzOUAwB5LwYWlxrCCDn9IM4BJDzOJQKX4Oc2/FLhZxP6Lnlp+vO19/TUXMzBfgLuuv1Bzw9mspBukuuAjkvGV3BG28fHpLP9G2SHl2XvrS2Vi7uXCJ/k6FNDZFz930myzUi51nOrvvYkHP3TLNeI3Ke9Qy7iQ85t+OInNvxS7yc96q0nLhpvTys63Vb582Tr+nGb3+hG2pRKk8AOS8Pc3M2+preTbJWX3wy5YjGFrlMp7pnYcM45Lw8fSartSLnWc1seeJCzsvDNcu1IudZzq672JBzO5bIuR2/RMv5kxMTcvymZ2W97sy+V2293Lx0d3mxrtGlVIcAcl5e7j8fG5GzdYbIU7phXLPMk7P12LWPtrWnesM45Ly8fSZrtSPnWctoeeNBzsvLN4u1I+dZzKr7mJBzO6bIuR2/xMr5vXp2+Uf1DHNzlvmrdbOsm3XEfGENK8wdpLvkKpDzktHFvtEs4bh6oFe+PNQn5uOX1dfrsWu7yQH16TwmEDmPnXouVALIOd2gGALIeTG0uNYQQM7pB3EIIOdxKBW+Bjm345dIOTfnQZ+vI4jT2rp3tbTK1brGvF6ntFOqSwA5rxx/M3pujl17SJdzmJekjs9tGLc4t7QjTQU5T1O2qt9W5Lz6OUhTC5DzNGUrGW1FzpORh6S3Ajm3yxBybscvUXJuZPzTujnW14cGdFKvyHk6rfe0tg4HEVKFCwLIuQuK8eswpxLcqi9UXdzXLWZd+uIa3TCua4m8M0UbxiHn8fPNlYyc0weKI4CcF8eLqxk5pw/EI4Ccx+NU6Crk3I5fYuR8q05fP0mnsf+XTmdvUjW/dtEyeVvLAgfRUYUrAsi5K5LF1bN5ekou7Nkk3x/ZkrvxjU3NclnnMtmjrq64iqpwNXJeBegpfiQj5ylOXhWajpxXAXrKH8nIecoTWKHmI+d2oJFzO36JkPMelY/jN66XX02M5daV37JkDzmwIZ1rbB2kI7FVIOfVTc19+sLVubrcw2yQaMondWbJqQmfWYKcV7fPpO3pyHnaMlbd9iLn1eWfxqcj52nMWuXbjJzbMUfO7fhVXc6fnpqQ927YviP78to6uWXpHvKSunoHUVGFawLIuWuixdc3Ktvkiv4euXawL3fzPrph3OU6iv6axqbiK6vAHch5BSBn6BHIeYaSWYFQkPMKQM7YI5DzjCW0TOEg53ZgkXM7flWV8z9Mjst7VMw36cj5/vUNelTaHrJI19VSkkkAOU9OXn6j3ztn64Zxj4yP5Rp1/II2ubB9sbQl7EQD5Dw5fSYNLUHO05Cl5LQROU9OLtLSEuQ8LZmqbjuRczv+yLkdv6rJ+e91N+r3bHhGulXMzRT2tUv3TN1O1A7Qp6oK5DxZ6TIbxn19S79c0rdZtuixa2bDuDWdi3MnHCSlIOdJyUQ62oGcpyNPSWklcp6UTKSnHch5enJVzZYi53b0kXM7flWR8ycnzIj5M9Krm8Ct0DPMv7V0ucyfxxnmDlJZ1iqQ87LiLbnyTVNT8qnejXL3yHCujkMbm+XyRUtlz9rqLw9BzktOq5c3Iudepr3koJHzktF5eyNy7m3qiwocOS8K1y4XI+d2/Cou54+rmL9v4zO5o6GMmK9dtkdud3ZK8gkg58nO0Q91w7hzdKq7WSZivqfObO+Sk3XDuGru6Y6cJ7vPJK11yHnSMpLs9iDnyc5PEluHnCcxK8lrE3JulxPk3I5fReX8l7o+dqWK+ZBOwX29ju59U0fMEXMHCaxQFch5hUBbPMZMbzfT3G/Q6e7TWo/ZXPEaPZbwIH0hrBoFOa8G9fQ+EzlPb+6q0XLkvBrU0/1M5Dzd+atU65FzO9LIuR2/isn5uvFR+YAelzasU9nNOc3fWLJcGhgxd5C9ylWBnFeOte2TzAthZ+koutk4zsxLOXZ+m3y6c0nF93VAzm0z6df9yLlf+baNFjm3Jejf/ci5fzkvJWLkvBRqL9yDnNvxq4ic/3RsVM8xf1bMMVBvbmqRry3eXermMZXdQeoqWgVyXlHc1g+b0hq+Otgvl/dvlhH93jMnIXxGBf3dLQus645bAXIelxTXGQLIOf2gGALIeTG0uNYQQM7pB3EIIOdxKBW+Bjm341d2Of/x6Ih8cNN6GVc5eFvzfLneiLmDNlNF5Qkg55Vn7uKJz05Oynm9G+Q+/V405XW6pOTqCm0Yh5y7yKA/dSDn/uTaRaTIuQuKftWBnPuV71KjRc5LJbf9PuTcjl9Z5dxsUPUhFfNJfYo53ulLuvaVU8wdJKxKVSDnVQLv6LHf2zokq3u7c8cXNupk9zN0w7iPl3nDOOTcUfI8qQY59yTRjsJEzh2B9Kga5NyjZFuEipxbwEPO7eAFdz/Xs31EzWW5e+sWOWnz8zkxf+/8VrmqaxkrzF0CrkJdyHkVoDt+5ND0NvmHvk1yy/CgzmXZvmHcFfqimTk5oRwFOS8H1ezWiZxnN7fliAw5LwfVbNeJnGc7v66iQ87tSDJybscvd7drOTcjdKfqZlRmzesHdCOqS7uWIuYO8lTtKpDzamfA3fMf0n0gTuvZIE9NTuS+N48xG8Z1LJa2mhp3D9GakHOnODNfGXKe+RQ7DRA5d4rTi8qQcy/SbB0kcm6HEDm34+dczr87PCRn6B/9ZlTub1sXymc7ljhoIVUkgQBynoQsuGvDhB67du1Qn1zd35vbE6JLxXy1fr++R2e6uCrIuSuSftSDnPuRZ1dRIueuSPpTD3LuT65tIkXObeix5tyO3o67XY2cf0unyp7XszEn5ie3tsvf60gcJTsEkPPs5DIciRk9P3vzRvn5uPsN45DzbPaZckWFnJeLbDbrRc6zmddyRoWcl5NudupGzu1yyci5Hb/c3S7k/OtD/XJRX3euPsTcQVISWAVynsCkOGzSWp31crGuR++bnpYGnex+WnunnNLaIfUWxx4i5w4T5EFVyLkHSXYYInLuEKYnVSHnniTaMkzk3A4gcm7Hz4mcf0XF/OIdYn5e+yI5XXeApmSPAHKevZxGI+pVMV/du1Fu1w0dTXlRbb1cs7j0DeOQ8+z3GZcRIucuaWa/LuQ8+zl2HSFy7ppoNutDzu3yipzb8bOW80v7N8vVg325esw0djNqTskmAeQ8m3nNF9VP9Ez0c3SJytNTE7kvH92yQI9C3K1oAMh50ci8vgE59zr9RQePnBeNzPsbkHPvu0AsAMh5LEwFL0LO7fhZyfk/9G2W63VDKVM+pxtJrdIN4CjZJYCcZze3hSL7gn6Pf2nH93irTm//xMIu+UhbuzTFPH8BOfevz9hEjJzb0PPvXuTcv5zbRoyc2xL0437k3C7PyLkdv5Ll/NO93fIvW/pz91/auVSOW9DmoCVUkWQCyHmSs1O+tj05MS6f0bXo9+louil/pmejn7KwU47V49fmKsj5XIT4epgAck5/KIYAcl4MLa41BJBz+kEcAsh5HEqFr0HO7fiVJOfn6nRXszO7KVd3LZP3Ojx6yUE4VFEmAsh5mcCmpFoz1X2NSvrjKuumLK+tk9NU0o9fUHjGDHKekuQmpJnIeUISkZJmIOcpSVSCmomcJygZCW4Kcm6XHOTcjl9Rcm6OSDNnmJuzzGv14y8tWibvanF3JrKDUKiijASQ8zLCTVHV/29kWL7Y3yOPTozlWr2Xbhp3VnuXHK0v0pmfC+GCnKcosQloKnKegCSkqAnIeYqSlZCmIucJSUTCm4Gc2yUIObfjF1vOp/TKUzdvkO9tHZI6/fh63RzqbbpJFMUfAsi5P7mOE+nduqO7kfRfT24fSX+xTnc/S09reJf+XKjZUQFyHock1wQEkHP6QjEEkPNiaHGtIYCc0w/iEEDO41AqfA1ybscvlpxP6lUndT8nd+uImRHzry1ZLm9uanHwZKpIEwHkPE3ZqkxbzWyaH4xskcv11IYnJ7bv7P6y+u2S/o7mBbJbZ7Ns6huRaXMhBQJzEEDO6SLFEEDOi6HFtcg5fSAuAeQ8Lqn81yHndvxiyfnfqpjfo2Juyk0q5n+JmDugnr4qkPP05aySLb5DZ9VcoSPpv5/cLukvr2uQz+65mxw63YCcVzIRKX5W1uR8w9SkdE9PSd/UtGyYnpTNk1PSrccTduvne/Tz3fr/zea9vpmyQE9EaNCTEBprana8Fz0ZQT/Wzzear+143zivRv8v+vbCdS01tdIxr1YW19XKUt0PYh/9/luo9WS5IOdZzm55YmPkvDxcs1Yrcm6X0arL+X5HrCo6gv333VvWXr+66PvKdcNzPdt3YY6WMdkmqzatlx/rRlDm6KRvqJgf1tRcrmZQb8IJIOcJT1BCmneb7klxuUr6n3ackb5/fYOcoyPpb2men5AW0oykEkiDnD+t/Ton2SrURrA3T03JJiPc+rmccO8Q8oHp6apj7lR536ehUV6qS05epu9fVt8oL9GPl9eZOXDpL8h5+nNY6QiQ80oTT+fzkHO7vCVCzh+/94bYUdxyxw/ljrvuT7ycj2yblhNVzP9nbFRa9Bf8N5fsLq9rRMxjJzqDFyLnGUxqGUO6VU90+NJQr/xxfPtI+oEqB+eppB/BzJsyUk931UmS8636O/CR8TF5WH8HrhsbkV+Nj8oGFfFiSpeOXC/SUewltbWySEe2u2r0YxVj87ll+rkG/d1aTBlW+R/ftk3MVoyj2r7cx7n3IiP6tU3aPjMSv1FfIPitLjMxMeQrzfpi+z76otlLdQnKPg1NKu36sY60m/+nqSDnacpWMtqKnCcjD0lvBXJulyHk3I5f7u7oyPmw/kI/fuNz8sD4SG6a3c1L95AV+guc4jcB5Nzv/JcSvdkQ7opnNuamuz+vwmDKwfqz5HyVdGbhlEI02/dUU85/uUPEf6m/9341Nia/2bHRYZT4EpXsxSrWi1WwzXszhbwrJ991uY87zcfm63pdtYsZ3f+TLjP5g779UeP7gwq7GfmfTdxfqaJuRtj30+9T8z36Kn1RLakFOU9qZpLbLuQ8ublJUsuQc7tsIOd2/HaR8yHduemYTc/kRgzaaubJrUv2zI14USCAnNMHiiUQ3q39G1sG5OqBHh3V2z76+Oc6E+f8ji45pIEZOcVyzer1lZJzI6zmKMCHRkflIR0V/4WOiucrZknGwdpPD2hskgP09+ABKq1ZKWaE/Y/K4SkV9j8qCyPuj+v7YM+IIE6zbv2wxhZ5Q3OLHKFve+rRiUkpyHlSMpGediDn6clVNVuKnNvRr7qc2zU/GXcHI+dmjdzKjc/m/mgxv5C/rSPm+2foj5Fk0E5vK5Dz9OauWi3Pd5Ta17b0yzX9vTObYL2rpVWOWbBQ3sh+FtVKU2KeWw45N7/XHlb5fsi86f4pD+rI+GCe4wP20Snd5vfdgSrjBzU2yms9fdHoGR1Zv3dkq74Ny09Gt8qQTp0Pl711zfobdGnKG/XtUH1rq+Kmc8h5Yr51U9MQ5Dw1qapqQ5FzO/ypkPNg07hi1qYXg2XlSWvksSee2umWVSvfJud+7Jjc585c/WW55751uY/zbUZn5NxsZPM+FfMnJsZ1XVyNfGfpnrKvjhpQIBAQQM7pC8USKHTO+ahuNvmNoX65dqBXf/ZsXxdr5Oii9sXyFyrpZp8Lin8EXMj5L1S+zTrxR3RqunlvpnFHi1nv/Sqdtr1CRdy8f7WOjJslXJRdCTyiL2r8t76ocZ/K+k91lkG0vEpf0DCj6m9sml/xpSrIOT22WALIebHE/LweObfLe2Ll3Gz89tmrbpqJ7q6bL5G9li+1i7bA3UbOVxy074yMhy8z7fjKjXfKfbdfnft0vmt/2b1FxfyZ3Do0s2nNbSrmadsYpixgqXQnAsg5HaJYAoXkPKjHSPrturv7dSrpZl2sKY26WdU75i+Q981fqCN0THkvlnmary9Wzp/UF5Mf1iVYj6qQP6QibtaNR0u7vthslma9OifhzTkR70rAevA05sl8v/5MRf1+FXUzuv7rPOvyj9Ap8EcuaJO/aVlQ9hc8kPM09qLqthk5ry7/tDwdObfLVOLkPDxKbUK7/tKz5bBDDrCLco67Z5Pz6Neisr5ehfwNv/mdPKV/GJvRhO+qmJtpaxQIRAkg5/SJYgnMJefh+szo3M06mv5vI1tmPr27brD13vnbp72/iJ9LxeJP3fWzyblZI/2gCviDOnr7sL49Mj6+y27k5sjPA3MSrjKuL+yYNeL0m/J1gz6d9XL/2Fb58dZhuU+nwK/fselj8MR3Ni+QY1oXlu2EBuS8fLnNas3IeVYz6zYu5NyOZ2LkPHze+UVnnCDHHvVmMZ+rlJyHp7WHp7QffvTpcvKJR+baY8r9DzwqJ513uQRT7F/8+G/kKf0jx/wRfNuyPWSvBG32Ytc1uNs1AeTcNdHs11eMnAc0+vUP/u/qMWzf1LXpZjZPUF6r0rVS/9A/Ukfk5jPtPZOdJ5DzzaMTKuA6LV2nVBsZ/6W+DzYSDAduRsQP1pFws7O4kXKz0zilegTM9+v3R4bk9i2DMzNhTGuW698Xx+r37jHz22Q3/dhVQc5dkfSnHuTcn1zbRIqc29ATSYScF1pTXik5DyMM5Dt4UcC0IXixICznwTT7eQ//SvbWtZ53vWhv2V3PX6VAoBCBWl2Saf4YGh4r7qxfiPpLYEFznQyPTOpk2NLKT3Xq7M39/fKdwQEZCW1M9f7WNjmuo0Pe1DK/tIq5K1EEHtFd09fpyKtZI/7Q1hF5XMU8WsxZ3Ac367R0fVuho+Kv1feU5BL4xciI3DrQL9/S793hHftKmNa+RZesHN/eIe9ubbVufJ3+UqqvrZGRcX4nWcP0pIJW/Z00pL+TKBCYjYDpJ5TSCSRGzhd1LpxZ1x2EUw05N88OT2Wfa+T88k3dcmRTqyyrq/6ZrKV3A+6sBIEaPVrPyPnWUX6xVYJ3Fp6xoLlehnUUNLLhc9GhDW+blu8MDMoNA32546+CYl5Q/MDChXJW1yJd38omckWDrcIN5siuX6mE/1ynQq/TXP5cJS5a9tC8vrpJN2zTjcZWqIQfrB8zW6IKyXL0yNuGBuVmFfX/Nzw8U+MSXUZ3bNtCWdXRKS/RAYJSSp2KeX3dPBnhBeNS8Hl5T2tLvQxt3XWTSC9hEHRBAqafUEonkAg5N80PT2sPj1pXYlp7FF9Yzudac27uDY5SKz0N3OkDAaa1+5BltzGWMq19rhb8RjehummwX27Tqe/hY55eo9Obzdr0D+jUWUoyCATHmK3T/QRyo+ITux5jZo7tPFinpB+kb29ob5UDdZ1400Spcy2SETetyE9g09RUbsnK2i0D8rsdG0CaK1+nx9Ydp9Pe36P7SxRTmNZeDC2uNQSY1k4/iEOAae1xKBW+JjFyHjSx0hvCPb1+o1z51e/IlWtOyTUh2CU+mLYeZ7d25NyuE/pyN3LuS6bdxVkOOQ+37l+3bpFbdBO5n0SOeHp5XcP2452a58ufq7RzNJu7nM5WkznG7EGdkm7OE39Ud06PHmPWrBu2HaT5OEjXir9Gd07fXzdu2zO0z0mxu7VXJiqeUg4CZj+Btbo2/Y7Qi2ytepzdu/XFtRNb2+UVMfYPQM7LkZls14mcZzu/ttE9oL/DrtMZev/x8pfYVuX1/YmT8yAb4aPU8k15d5m18Ki9qTc6Wh/nnHOX7aGubBJAzrOZ13JGVW45D9rerTt536ZHsv1A3x7SP/qj5RAdmXtDS4scrsc8mRF2ih2BrbrM4PGJMXlCNxP9X5XwdfriyON6rFm05I4v07cDdFr6QToibtaNz1aQc7u8pPVu8yLbd3U0/V7ddyAor9MXbz7a1iFv0xfYChXkPK0Zr167kfPqsU/6k42YH7thvZgjI7e9+lVJb26i25dYOQ9TK7RhXFLIMnKelEwkux3IebLzk8TWVUrOw7EPTW+Tn4wOy4/1D/3/0WOewju+m+vM2uW/0A3F3qgj64c1zldhZG1Zob5jJPwJle7f6hTk3+qLHr/JHV+2TX6tYm6mrIfLK3S2wgE6En6AvvjxGn0xxOykXmxBzosllq3rzYtst+po+teH+mZ253+JHqF4Zvsiebee0hAtyHm28l+JaJDzSlBO3zPMMZ3HbFqfO57zCH0R/0evfGn6gkhQi6su50a8g2PJ4nAxI+p33HW/rL1+dZzLK3INcl4RzKl/CHKe+hRWPIBqyHk0SHMElxmRu1dH5+7X9736yzdcltTUylt0dO7dul7dTKVt1zXQvhWz4Z4ZBX9S1/M/qe+fUPn+rUp59NzqMJcVZo14U4scquwOamiQZgcb8iHnvvW8wvGa0fRrBjbLkzuOU3ypSvpZ7V3yrpYX1qUj5/SXYgkg58USy/71D+0Qc/N78I36wv0NS5bL3l0t2Q+8jBEi5w7gIucOIHpQBXLuQZIdh5gEOY+G9L8qnf9tRtb1mLaf6dpoM4UtXA7XX85m5LdBJX0v3TX8z3RE+CUq7Qt0PWzai9lA70mdhr5dwnVaukq4kZ/npwqfwLC3SpGZjr6vvr1MhXwfnWmwv05RL0dBzstBNd11fn9ki1ze/4Kkm5kuZ+tI+juaF+ROD2lprJXeoV2XVKQ7alpfLgLIebnIprNes/fFMRuflS36uzEQ80bdG4UN4ezyiZzb8cvdjZw7gOhBFci5B0l2HGIS5Twa4s90vfT/6Nuj+uq5GT3+01T+Y3Y6dWTYSPpeKgcvUTl9UU7c6+Wl+nFSxd2swzfxbJqclHtUcswsgkLFjEy+VON7hcazj46Em7j2i7Epl8sug5y7pJmtur63dUiu6O+ZEtnNcgAAIABJREFU2eXdSPoFi5bIys525DxbqS5rNMh5WfGmqvJf6gvU79/4zC5iboJAzu1SmQg5LzaE/ffdm2ntxULj+qoTQM6rnoLUNSANch6FOqYj6WZKt1mr/nszuqyjzH9QaTdHuBUqHTrKfohuYLXfjhHlF+t662W1dXkvn6f1N+kofJPe06Sv0Oc+VvFvqtH3+v9w2aAj2n26trtP1+Lm3lSuzf97c5/f/n+z9rsn+HpkHXi+Bpid7I2E76vrw/fR9y/V/8fZGbsSnQ85rwTldD/DTHe/QkfS/7DjKLb9tB+f2dopb8+zJj3dkdL6chBAzstBNX11GjFfqWJuZpMdprPlbtSp7GbEPCjIuV1Oqy7nds1Pxt2MnCcjD0lvBXKe9Awlr31plPPZKP5BBf13Ku1G3p/UKeG/2yHxZq1aUooZxe+orZWOebWyTEf3zQsG++uLBWa030xJT3JBzpOcnWS1zcwKuXKwR57asSbdbEh4Xscieessu7snKwJaUw0CyHk1qCfrmY/q7+73bXhBzL+hYh59YRw5t8sZcm7HL3c3cu4AogdVIOceJNlxiFmT80J4zJrtP+lInhnpNqPYvTqNfESXso+qtI/qaLY53C33ce5t+8cjua/J9vf6Zl7Bj5bFulldR22NdNTU6Zt5XyudKt7t5v/m8zri3lmn/1cRN0Jurk9zQc7TnL3Kt92sOb99dEgufn6TPL1jOcr+OhvkHF2TbjZ5pEAgSgA597tPGDE3U9kH9VQXM2KeT8wNIeTcrp8g53b8kHMH/HypAjn3JdPu4vRFzt0R87sm5Nzv/BcbfbAh3CbdEO67w4NydX/vjKSbo/zMxnFv1hMFKBAICCDn/vaF8FT2Q3UZ2k1Ld57KHiaDnNv1E+Tcjh9y7oCfL1Ug575k2l2cyLk7lj7UhJz7kGV3Mebbrf2bWwbk6oFeeW7HCQSv0uUcZ+sRbH/FSLo78CmuCTlPcfIsmh4+x9yI+Y0q5tGp7Mi5BeDIrci5A5ZMa3cA0YMqkHMPkuw4ROTcMdCMV4ecZzzBjsOb7Sg1I+lXqaQHxwS+Xv8gv1DXpL9ajwOk+EsAOfcv9+a4tA9tek426ZIzc1Tq1yObv+Ujwsi5XT9Bzu34MXLugJ8vVSDnvmTaXZzIuTuWPtSEnPuQZXcxxjnn/F+29MuVegSbOeXAlHe1tMpnOhbLEt2jgeIfAeTcr5ybY1JP3LheRvSUlL/SJS5mjXmcgpzHoVT4GuTcjh9y7oCfL1Ug575k2l2cyLk7lj7UhJz7kGV3McaRc/O0LbrZ4nUDPXLVYF/u4S26kaKZ6n5ya7u7xlBTKggg56lIk5NG3jMyLH/b/Vyurnc2L5CvLN4tdr3IeWxUeS9Ezu34IecO+PlSBXLuS6bdxYmcu2PpQ03IuQ9ZdhdjXDkPnviM7uj+yc2b5N6xrblPvaSuXi7rWiqv0ynvFD8IIOd+5Pn7I1vkE93Py6SG++6WBfKlRbuFTjGfmwFyPjej2a5Azu34IecO+PlSBXLuS6bdxYmcu2PpQ03IuQ9ZdhdjsXIePPk/dUTt073dMzu7H2WmuncuTv1RhO7IZrcm5Dy7uQ0iu1VPbjinZ6NOZBc5ccFC+XznkqLE3NSDnNv1E+Tcjh9y7oCfL1Ug575k2l2cyLk7lj7UhJz7kGV3MZYq56YF4/qn+zV69Nq1OtXdfLxg3rzcVPe/a+0QVqO7y1HSakLOk5YRt+3556F+Wd3Xnav0jLYOOVePUyylIOelUHvhHuTcjh9y7oCfL1Ug575k2l2cyLk7lj7UhJz7kGV3MdrIedCK6FT3l9XXyxc6meruLkvJqgk5T1Y+XLbmsv7NM/tKrNFNHz9ssacEcm6XGeTcjh9y7oCfL1Ug575k2l2cyLk7lj7UhJz7kGV3MbqQ86A1/2GmuvdtkmcnzSrV7etUV+t02MU1jKO7y1j1a0LOq58D1y0w09cv6N0kN+rxifP04y/qPhLHzG+zegxyboVPkHM7fsi5A36+VIGc+5Jpd3Ei5+5Y+lATcu5Dlt3F6FLOTavGdHr7VXrs2lcG+2emup+j02I/pCNwKLq7vFWzJuS8mvTdP9uI+ambn5d/3bpF6vTja3VH9nfozuy2BTm3I4ic2/FDzh3w86UK5NyXTLuLEzl3x9KHmpBzH7LsLkbXch607E+TE3Kubij1Ez0j2RQz1f2yrmWyoqHJXeOpqSoEkPOqYC/LQ80cF7Mju9mZ3Yj5DXqG+Zv0LHMXBTm3o4ic2/FDzh3w86UK5NyXTLuLEzl3x9KHmpBzH7LsLsZyyXnQwrt0NM5sLrV+iqnu7rJW3ZqQ8+ryd/V0s4njhzY9Jz8a3SpNOpn9GyrmhzW5OxIRObfLFHJuxw85d8DPlyqQc18y7S5O5NwdSx9qQs59yLK7GMst56alZqr7FTrV/XrdBXpi2zZp1V3dzVT3v2Wqu7tEVrAm5LyCsMv0qK3bpuX4jc/Jz8dHct+P31qyhxzc6HZWC3Julzzk3I4fcu6Any9VIOe+ZNpdnMi5O5Y+1ISc+5BldzFWQs6D1jLV3V3eqlkTcl5N+vbPHpyelpUbn5VfTYxJR02NfGfpnvKK+gb7iiM1IOd2SJFzO37IuQN+vlSBnPuSaXdxIufuWPpQE3LuQ5bdxVhJOQ9a/QNd3/qZ3m55bsdU96N1V/fP6K7uXezq7i6xZawJOS8j3DJX3adifvTGp+XJiQlZWlsrt6mY711XX5anIud2WJFzO37IuQN+vlSBnPuSaXdxIufuWPpQE3LuQ5bdxVgNOTetH9FptVcO9Mq1g325YBbo1NrzO3RX9wXt7oKjprIQQM7LgrXslW6cmpL3bnxG/qCbNe5VWy/fXbqHLK8z28CVpyDndlyRczt+yLkDfr5UgZz7kml3cSLn7lj6UBNy7kOW3cVYLTkPIvi9isIFvRvl/tHtu7rvr9Nrza7ur2podBckNTklgJw7xVmRysySkvfqVHYzW8WcnLB2yZ6yREfOy1mQczu6yLkdP+TcAT9fqkDOfcm0uziRc3csfagJOfchy+5irLacB5GYM5ZXq6T36LRbU8wI+ic7umT+vBp3wVKTEwLIuROMFavkyYlxeb+Keff0lBxQ3yi36oh5u641L3dBzu0II+d2/JBzB/x8qQI59yXT7uJEzt2x9KEm5NyHLLuLMSlybiIamt4mF/dtkm8ND+r+7iLLdGTv4o4l8te6Jp2SHALIeXJyMVdLfjU+Jis3PSOD+r11cEOTfEvF3OzOXomCnNtRRs7t+CHnDvj5UgVy7kum3cWJnLtj6UNNyLkPWXYXY5LkPIjqkfFROW3zBjFT3k05orFFrli0LLeBFaX6BJDz6ucgTgseGhtVMV8v5tg0c365OcfcnGdeqYKc25FGzu34IecO+PlSBXLuS6bdxYmcu2PpQ03IuQ9ZdhdjEuXcRDepb1/VzeIu1/PRR3UcvUWnt5/T3ikfbu0QFN1d/kupCTkvhVpl7zF7OHxQxdx87/yf5vny1UW7SV2FRsyDSJFzu5wj53b8kHMH/HypAjn3JdPu4kTO3bH0oSbk3Icsu4sxqXIeRGg2sDpDR9F/MrZ9wzizmdU1Khpm7SylOgSQ8+pwj/vU/xgZlo92P5d7gevduiTEfL+Uf4X5rq1DzuNmLP91yLkdP+TcAT9fqkDOfcm0uziRc3csfagJOfchy+5iTLqcB5HesXUodza62dTKlA8uWCif6lhcsfWz7oinvybkPLk5/E8V81Uq5qYcN79NLu1aWrXGIud26JFzO37IuQN+vlSBnPuSaXdxIufuWPpQE3LuQ5bdxZgWOTcRmw3j/m9/t3xjy0AOwOKaWvlM52I5qqXVHRBqmpMAcj4noqpc8I+6DORz/Ztzz/54W7tc2L64Ku0IHoqc2+FHzu34IecO+PlSBXLuS6bdxYmcu2PpQ03IuQ9ZdhdjmuQ8iPrRiTHdMO55eXJi+4ZxhzY2y1W6YdzutXXuwFBTQQLIebI6hznZ4ILeTXKjvmhltnv7dMci+ajuzVDtgpzbZQA5t+OHnDvg50sVyLkvmXYXJ3LujqUPNSHnPmTZXYxplPMg+usG++WzOpIeFDNSaEYMKeUlgJyXl28xtQ/rTuwXqph/Z3god9tXF+8mf9OcjKMHkfNiMrnrtci5HT/k3AE/X6pAzn3JtLs4kXN3LH2oCTn3IcvuYkyznBsKZsO483s2yn+Nbs1B2Uc3jLuic5kc3NjkDhI17UQAOU9Gh3hmakJWbXxOfjM5Lm018+RGPSrttQ3NyWictgI5t0sFcm7HDzl3wM+XKpBzXzLtLk7k3B1LH2pCzn3IsrsY0y7nAYl/37pFLurbJBunpnJTe4/XDeMuZMM4dx0lVBNyXhasRVX6Mz294EO68dvA9LS8tK5evrl0uexZW19UHeW+GDm3I4yc2/FDzh3w86UK5NyXTLuLEzl3x9KHmpBzH7LsLsasyLkhYqb4fqGvR27Y0i/T+n+zYdw/6IZxR7JhnLsOozUh505xFl3ZDUMDslpfiDJHpb2hqVn+ZfHuMn9eNQ5Lm73pyHnRqd3pBuTcjh9y7oCfL1Ug575k2l2cyLk7lj7UhJz7kGV3MWZJzgMqZsO4M7s3yK91uq8pbBjnrr+YmpBztzzj1mZk/LzNG2StHitoysmtuiO7zg5JnpZvjwg5j5vZ/Nch53b8kHMH/HypAjn3JdPu4kTO3bH0oSbk3Icsu4sxi3Ju6JiR838Z6pfL+ntyI+pNOtn9rPYuOamtQ9jT3a7/IOd2/Eq5u0+nr6/atF7WjY9K/bx5cmXXMnl3SzI2fisUD3JeSqZfuAc5t+OHnDvg50sVyLkvmXYXJ3LujqUPNSHnPmTZXYxZlfOAkFmDfpaONt47tn3DuJfo+txr9Ni1gxrYMK7UXoScl0qutPuemBiX41XMzeaHHTU1cvOSPeTAhsbSKqvgXci5HWzk3I4fcu6Any9VIOe+ZNpdnMi5O5Y+1ISc+5BldzFmXc4DUnePDMtFeuTU8yo4ZsO4D8xvk7/vXCKtOgpJKY4Acl4cL5urvz+yRc7QJRojsk32rW+Qb+qO7LvXpmPuB3Juk3kR5NyOH3LugJ8vVSDnvmTaXZzIuTuWPtSEnPuQZXcx+iLnhtjIjg3jvq4bxk3p/82GcatV0JM+Pdhdtt3UhJy74ThbLdv0i5f2b5ZrBvtyl/1lU0vuDPPmBG78VigO5NyunyDndvyQcwf8fKkCOfcl0+7iRM7dsfShJuTchyy7i9EnOQ+o/VqnCZ/a/fzMhnGva2yWqxctTdxRVO6y7LYm5Nwtz2htZo+Ek7V//tfo9qUYp+g+CZ9sX5Sb8ZGmgpzbZQs5j/A7/OjTZdniTll7/eqZr5y5+styz33rcv/ff9+9d/qa+dxzPSN2WeBuLwgg516k2WmQyLlTnJmvDDnPfIqdBuijnBuAZsM4c+TaJX2bZcu2bdKo6nOmbhj3MTaMm7N/IedzIir5gmemJuT4jevld5MTuY3fvqz7I7yjOdkbvxUKFjkvuRvkbkTOQ/yMmJsSlvNb7vihfOXGO+W+26/OfW3lSWtkxUH7yrkfO2bmTuTcrhP6cjdy7kum3cWJnLtj6UNNyLkPWXYXo69yHhDsnp6ST+la9Lu2bsl9ig3j5u5byPncjEq54mdjI/Kh7udkQHdmX6RLLm7U9eVp2PgNOS8l23Pfg5zvYGSk+6i3HybPPt8t6x55YmZ0PCrjUVk3tyPnc3c0rhBBzukFxRJAzosl5vf1yLnf+S82et/lPOB1r04hPqtng5jd3c304WN0w7hP6xnSbbo7NmVnAsi5+x7xNZ3Fsaa3W8xZ5vvrxm836o7sS2tr3T+ogjUycm4HGzlXfmEBv+y6W3eSczOafvKJR8qxR705R/r+Bx6Vk867XB6/94YZ8si5XSf05W7k3JdMu4sTOXfH0oeakHMfsuwuRuT8BZZmw7jL9Vz06/R8dFM6dfOti7uWyFEtre6AZ6Am5NxtEv9epdzIuSn/p3m+XKtT2dO08VshGsi5XT/xXs7NenJTrlxzSu59VM73O2KVXHTGCbvI+V03XyJ7LV+au2da1yxRIDAXgWBDD3rLXKTS8/Xyb9JinkCPSU+PqF5LxyampbGekb7qZYAnZ4HAY6Oj8uGnn5UHtm7fS+gvF8yXf95rD/mzhoZEh1ep3xI1uhaav3ntu0L35KS856mn5SfDw7nKPr10iazebbtTZKGYfkIpnYD3cm5Gxjf3DuxCcFHnwtw68zgj58+zIVzpPdCjO+tq54kZ2do8MOZR1NkOtdx/EDFynu3+4zo6Rs5dE812fYycF87vDUMD8jk9zmqrjqibcq5uGHdGW2diO0SlVGi3rmbhb167bvDI+Jis2rReNumeBwtUYq/qWiZ/3ZLOjd8KkTD9hFI6Ae/lPIouOnLOmvPSOxd37kyAae30iGIJIOfFEvP7euTc7/wXGz1yPjsxs2Hc3/dskn8beWHDuMu6loo5fs3XwrR2u8x/dahP/nGgT0zfMv3oSyrmy+vq7CpN4N1Ma7dLCnIe4ReVc3Zrt+tg3P0CAeSc3lAsAeS8WGJ+X4+c+53/YqNHzuMR+5FuGHduz0Z5fsps2WU2jGuV1R1LvNwwDjmP12eiV5lj0k7r3igPjG9fLnGBnl3+CT26L6sFObfLLHI+h5ybL3POuV0n4+7tBJBzekKxBJDzYon5fT1y7nf+i40eOY9PbEz3/viiTnP/6mB/bldts2Hc6s4l8l4VdZ8Kcl58tm/cMiAX921fIvFK3Y392kW7ycv0fZYLcm6XXeTcjl/ubnZrdwDRgyqQcw+S7DhE5Nwx0IxXh5xnPMGOw0POiwf65MSEnNnzvJh1w6aYqclX6lT3F9XVF19ZCu9AzuMnzRzN94nNz8tP9QxzM3H9dN234DTdtyB7k9h3ZYKcx+8n+a5Ezu34IecO+PlSBXLuS6bdxYmcu2PpQ03IuQ9Zdhcjcl4aS7MR6Dd1NPTz/d0yOL1NGvR09FPbO+XU1g6pz/gu1ch5vD6zdnhIVvdulCE9zell9fU6Wr57btTcl4Kc22UaObfjh5w74OdLFci5L5l2Fydy7o6lDzUh5z5k2V2MyLkdyx6zYZyeU/29rUO5il5UWy/XLF4mKxqa7CpO8N3I+ezJMRu9naqj5f89OiK1eunH9QWbc3TEvC7jL9pEqSDndt/EyLkdP+TcAT9fqkDOfcm0uziRc3csfagJOfchy+5iRM7dsPyxitj5umHc07rpV4uuRX97y3w5e2FXJqe6I+eF+8wd+iLNBb2bZGB6Wl6qyxyu0bXlBzY0uulkKasFObdLGHJuxw85d8DPlyqQc18y7S5O5NwdSx9qQs59yLK7GJFzdyxNTZ/Tae7/qBvGBeU03Y37Q/q2uMaMoWajIOe75rFXZfw8ncJ+19btR+59tLVdd/NfnI2ElxgFcl4iuB23Ied2/JBzB/x8qQI59yXT7uJEzt2x9KEm5NyHLLuLETl3xzKoaYMet3ZFf4/cPDw4U/lHVNZOWdgpizIg6cj5zn3mbhVyI+Y9Kuh76HnlZif2LC9riPsdg5zHJZX/OuTcjh9y7oCfL1Ug575k2l2cyLk7lj7UhJz7kGV3MSLn7lhGa3pqckIu16PX/nXHaGqjbhr3wdaFOUnvSrGkI+fbM23Wln+pv1f+Zcv2mRJ/t6BdLuhcJE2aZ4oIcm7XC5BzO37IuQN+vlSBnPuSaXdxIufuWPpQE3LuQ5bdxYicu2NZqKbfq6RfGZJ0I2+rVNI/oWvSO2tqyt8Ax0/wXc7N7utfHuiRf9LlC2OyTXarNaPly3JH6lFeIICc2/UG5NyOH3LugJ8vVSDnvmTaXZzIuTuWPtSEnPuQZXcxIufuWM5V0+/0fPQrVOqCnd2NpP9dW7t8TM+97kiRpPss51/TUfKvqZSbWRGmfFzzd5a+yNKsmwBSdiaAnNv1COTcjh9y7oCfL1Ug575k2l2cyLk7lj7UhJz7kGV3MSLn7ljGrelJlfQv6kj6D0a2bx5mdndf1domp7R1ycIUSLqPcm6WJlyim/09MzmZy9nrGprlC11L9Pxyf84tj9u/g+uQ82KJ7Xw9cm7HDzl3wM+XKpBzXzLtLk7k3B1LH2pCzn3IsrsYkXN3LIutKZD0f1dJ36Y3t+o52B/WM7E/qru7tyVY0n2S8x+NbpXP93XL/06M59L78roG+VTHIvmr5vnFptu765Fzu5Qj53b8kHMH/HypAjn3JdPu4kTO3bH0oSbk3Icsu4sROXfHstSafjM5riPpPWJ2/Q4k/aO6adxHVNSNsCet+CDnD4+PyprebvmFvjdlr9p6Obe9S949v5Xt3mJ2SOQ8JqgClyHndvyQcwf8fKkCOfcl0+7iRM7dsfShJuTchyy7ixE5d8fStiYj6Zf1bZa7R4ZzVZkp7ifpevQP6+Zx8xO0pjnLcv4HzcHnQjkwR9+d3t4pJ85fKHUJfKHEts+V837k3I4ucm7HDzl3wM+XKpBzXzLtLk7k3B1LH2pCzn3IsrsYkXN3LF3V9NjEmB7B1iP37JD0F9fVyyGNTXK0CuKhTdXfETyLcm7OpjezF27ZcTa9eTHkE7nZCwtzewJQiieAnBfPLHwHcm7HDzl3wM+XKpBzXzLtLk7k3B1LH2pCzn3IsrsYkXN3LF3X9KhK+pUqjP+xQ9JN/Qc2NMqbmlpkpUqjmWpdjZIlOTc76H91sFf+3+iwbJyayuH8sJ5XfoZOYU/TDvrV6AdzPRM5n4vQ7F9Hzu34IecO+PlSBXLuS6bdxYmcu2PpQ03IuQ9Zdhcjcu6OZblqMsd2/dNgn6zdMiijuVXp24vZMdxI+jtb5ld0dDcLcm7W9/+HbsT37eGhGZ4rW1rlHN3sbXc9t5xiTwA5t2OInNvxQ84d8POlCuTcl0y7ixM5d8fSh5qQcx+y7C5G5Nwdy3LXNLxtWu5UqVw7NDCzUZl5ZqNuUfbO+Qvk/RWa9p5WOTdr+g2723Tqes/0dC5dZl2/WU/+Id0hf0ltbblT6FX9yLldupFzO37IuQN+vlSBnPuSaXdxIufuWPpQE3LuQ5bdxYicu2NZyZqe1fO2b90yIN/dOjhz9rZ5/h51dfK+ljY5ZsHC3MflKGmS80GV8Nt1dHytsvqVLhMIyusam+WE1nZ5hx6JVs9Gb+XoJoKc22FFzu34IecO+PlSBXLuS6bdxYmcu2PpQ03IuQ9Zdhcjcu6OZTVqMpPcfzY2qvLZL/++dVjM6Lop5gC2Q3Ta+/tb2+TIlgVOp70nXc7NyvEf6/nkt+ooudlUb3zHUgCzhvx989vkgyrlf6ab7FHKSwA5t+OLnNvxQ84d8POlCuTcl0y7ixM5d8fSh5qQcx+y7C5G5Nwdy2rXZNaj/9uwrqNWUf8fFfagmGnvf6tr019c3yhH6GZyyy1H1JMq54/omeTrNO5/1A3egs3dDIO/0J3uj1MhP0rXlFMqRwA5t2ONnNvxQ84d8POlCuTcl0y7ixM5d8fSh5qQcx+y7C5G5NwdyyTVtF6nvX9bp7x/e2hQnp6a2Klpf65TupfppmevaWqSFTq6/irdAb6YkgQ5Ny9EPKIi/tDYmKzXteRmY7etO2YNmFi6dJT8/TpKbqR8b0bJi0mvs2uRczuUyLkdP+TcAT9fqkDOfcm0uziRc3csfagJOfchy+5iRM7dsUxqTet0RPlenfJ+v071/oV+HC1mZP1gHV1e0dAkB+s56q9VeZ/tGLFqyPmfdMf6B42Mj4/IutFRMcfMRUunnkd+qM4M+GvdHO9IRsmr3h2Rc7sUIOd2/JBzB/x8qQI59yXT7uJEzt2x9KEm5NyHLLuLETl3xzINNY3tGHF+UEecHxxTWR8bmdm5PNx+M9r8ah1Rf5+e+d1cI/Ly+iZp3bFxWiXk3Iz236lT9Ndp+x4cHZHe0Kh40M6gja9pbJFD9EWFV9Y3pCEF3rQRObdLNXJuxw85d8DPlyqQc18y7S5O5NwdSx9qQs59yLK7GJFzdyzTWpMR4Yd1VNqMSD+oI9O/HN91VLpZR9cP1NF1U5rqa6Rlep4srqmVJbp+fZG+X6TT5M3/F+lxZOZtvo5iB2Voepv0b5uUAd053eye3j81JQMq2/1T07J12mzftr2Yaenf17PHzU700dKi9R3U0CArdFR/hY6OH6yj/LON7qc1F1lqN3Jul03k3I4fcu6Any9VIOe+ZNpdnMi5O5Y+1ISc+5BldzEi5+5YZqmmB1TSzXrux3Tk+rcT4/KYvlWyvKKuQfZrbNTp9tun2TMqXkn6bp6FnNtxRM7t+CHnDvj5UgVy7kum3cWJnLtj6UNNyLkPWXYXI3LujmXWa+rRUe7NOuotC2rlf/t0OvzUZO5to450b9av9emo+Ab9/3P6Fi5m1HuhbtC20LzXUfV283HurVY6dMS9TafLt+t78zmzUR0ino2ehJzb5RE5t+OHnDvg50sVyLkvmXYXJ3LujqUPNSHnPmTZXYzIuTuWvtRUiTXnvrDMcpzIuV12kXM7fsi5A36+VIGc+5Jpd3Ei5+5Y+lATcu5Dlt3FiJy7Y+lLTci5L5m2ixM5t+OHnNvxQ84d8POlCuTcl0y7ixM5d8fSh5qQcx+y7C5G5NwdS19qQs59ybRdnMi5HT/k3I4fcu6Any9VIOe+ZNpdnMi5O5Y+1ISc+5BldzEi5+5Y+lITcu5Lpu3iRM7t+CHndvyQcwf8fKkCOfcl0+7iRM7dsfShJuTchyy7ixE5d8fSl5qQc18ybRcncm7HDzm344ecO+DnSxXIuS+Zdhcncu6OpQ81Iec+ZNldjMi5O5a+1ISc+5JpuziRczt+yLkdP+TcAT9fqkDOfcm0uziRc3csfagJOfchy+5iRM7dsfSlJuTcl0zbxYmc2/FDzu34IecO+PlSBXLuS6bdxYmcu2PpQ03IuQ9fnAhlAAANNElEQVRZdhcjcu6OpS81Iee+ZNouTuTcjh9ybsePuyEAAQhAAAIQgAAEIAABCEAAAtYEkHNrhFQAAQhAAAIQgAAEIAABCEAAAhCwI4Cc2/HjbghAAAIQgAAEIAABCEAAAhCAgDUB5NwaIRVAAAIQgAAEIAABCEAAAhCAAATsCCDnJfI7c/WX5Z771uXu3n/fvWXt9atLrInbskhg5Ulr5LEnntoptFUr3ybnfuyY3OfoP1nMevEx3f/Ao3LSeZfL4/fesNPNT6/fKG8/7vyZz11/6dly2CEHzPz/ljt+KJ+96qaZ/0fvL74l3JEWAvsdsUqi/eGy626VG9bevVMI4d9Lc/WntMROO+MTOPzo02Vz70DJP0PC94d/d8VvAVemiUD4bxLT7mjO+ZsmTdmkrWkngJyXkEHzh/FXbrxT7rv96tzd5ofWioP2nRGvEqrklowRmK1P0H8yluwSwonKUlSuzR/GJ594pBx71JslKvDBvXfdfInstXypGDFb98gTvEBYQh7SdEtYlvLJ+Wx9YLb+lCYGtDUeAfMz4sqvfkeuXHNK7obgxbzg58xcP0OMqJkS3J/vBaF4LeGqtBAwf7MEg0xB/wj/nOFvmrRkknZmgQByXkIWoz+korJVQpXckjECs/0io/9kLNkW4UT/aDZV5RtND8tVVMajf2hbNIdbE04g3x/NpsmzvUAzV39KeMg0zwGBuWQ8+vWojEdl3UGTqCLhBMK/c0xT+Zsm4QmjeZkigJyXkM7oD61CU1NLqJpbMkIgOgUsPEWM/pORJDsII5+c53uxL/yHUb4/lBnZcpCMFFQxm5yHp7WHp7TP1Z9SEDZNtCQQ/Rtltp8hey1fkltSE8zMmevFH8umcXtCCUR/p/A3TUITRbMySQA5LyGt5ofWRWeckJtyakrwiy/8y6yEarklowSC/hFMEaP/ZDTRJYSVT87NKOj3//OnM8tmglGL3Zd15aaZmj+Sgo+DR0b7VAlN4ZYUECgk59Gmm/4QvCA4V39KQdg00ZKAeUH4HW95/czSu9l+huy5+5LcPhhROY/+TLJsErcnmIB58ea5DT0Fl0rxN02Ck0fTMkEAOS8hjYx8lgDN81vCI5/0H887Qyh8Rs7pC8UQiCvn4WnujJwXQzh715rfNwcfsM/M+nETISPn2cuzq4hM33jo0d/u9OJwvrr5m8YVceqBwK4EkPMSegVrhkuA5vkt4T5D//G8M8wh53OtEWbNub/9pxQ5n6s/+Usz+5HnE3MT9Vw/Q1hznv2+kS/CuGJu7uVvGj/7CFFXhgByXgJndtsuAZpHtxTaKTeYJkj/8agzzBFqvpFzcwu7tdNH8hEoJOf5dloOL71it3b/+lN4aUM0+rk2iGO3dv/6i/kZYkq+Y4H5m8a//kDE1SWAnJfIn3OqSwTnyW3mD6NwiR59RP/xpCMUCDN6lJq5LLxp4FznUnPOuX/9J3pudXjTt9k2azKk5upP/tHMdsTRnw9BtG89fMUux6sFX8t3nGNwTjrnnGe7v+T7fWQiXtS5cGZ6O3/TZLsPEF2yCCDnycoHrYEABCAAAQhAAAIQgAAEIAABDwkg5x4mnZAhAAEIQAACEIAABCAAAQhAIFkEkPNk5YPWQAACEIAABCAAAQhAAAIQgICHBJBzD5NOyBCAAAQgAAEIQAACEIAABCCQLALIebLyQWsgAAEIQAACEIAABCAAAQhAwEMCyLmHSSdkCEAAAhCAAAQgAAEIQAACEEgWAeQ8WfmgNRCAAAQgAAEIQAACEIAABCDgIQHk3MOkEzIEIAABCEAAAhCAAAQgAAEIJIsAcp6sfNAaCEAAAhCAAAQgAAEIQAACEPCQAHLuYdIJGQIQgAAEIAABCEAAAhCAAASSRQA5T1Y+aA0EIAABCEAAAhCAAAQgAAEIeEgAOfcw6YQMAQhAAAIQgAAEIAABCEAAAskigJwnKx+0BgIQgAAEIAABCEAAAhCAAAQ8JICce5h0QoYABCAAAQhAAAIQgAAEIACBZBFAzpOVD1oDAQhAAAIQgAAEIAABCEAAAh4SQM49TDohQwACEIAABCAAAQhAAAIQgECyCCDnycoHrYEABCAAAQhAAAIQgAAEIAABDwkg5x4mnZAhAAEIQAAChQhcdt2t8v3//Kncd/vVQIIABCAAAQhAoIIEkPMKwuZREIAABCCQXgKHH326bO4d2CWAx++9YeZzZ67+sjz06G9TLbbIeXr7KC2HAAQgAIF0E0DO050/Wg8BCEAAAhUiYOT8HW95vZz7sWNmnrjypDWyobs31TIexYecV6hD8RgIQAACEIBAhAByTpeAAAQgAAEIxCCQT86jIlvo/yefeKR89qqbZp4SHm3P92gj/bsv68p96Z771uXe77/v3rL2+tW5j59ev1Heftz5cv2lZ8thhxwwU8V+R6ySi844QY496s25zwX/Dz971cq3ycoj35S7Pyjhe4IYzAsRN6y9u2CbzXXhr+erIxx3+OsxcHMJBCAAAQhAwDsCyLl3KSdgCEAAAhAohUA+OY9+Lp+cG4F96+Er5Mo1p+Qea8TblEC087XFXPPYE0/tItpGrM3IfTFybuoPXgy45Y4f5l4kWNS5cGa0P/hccE0g3cGz8rU5Gme0PUEd4bhLYc49EIAABCAAAZ8IIOc+ZZtYIQABCECgZAKF1pznGzEONlPLN0XcyPBXbrxz1qnwwch5IPSm0WY9uynmc8XIebh9+e4LPnfXzZfIXsuXSr423//Ao3LSeZdLcI0ZkQ8+DoCG28fU+JK7GTdCAAIQgIDHBJBzj5NP6BCAAAQgEJ9AvpHzQGyDUea5prmbp0VHqvO1oJCcP7ehJzfiXg45D6bI5xPr8PNMe42o5yvB1HvkPH6/4koIQAACEIBAQAA5py9AAAIQgAAEYhDIJ+fmNjNiHEizT3IeHTkPI0TOY3QoLoEABCAAAQhECCDndAkIQAACEIBADAJJknPTXDO1PO6GcMEGcbNNa59t5DyY1h6sS49uPBfFh5zH6FBcAgEIQAACEEDO6QMQgAAEIACB4gnkk/NAWoN13ZUaOTetN+05+IB9ZjaaMyP4Zmf3fLu128p5oWeFR88Ni9t+cF+uPch58f2LOyAAAQhAAAKMnNMHIAABCEAAAjEIFNoQLjx6XUk5D0bBg6YbKTc7sbuS8829AzNU8u26Hj1KzVw82+h7DMRcAgEIQAACEPCaAHLudfoJHgIQgAAEIAABCEAAAhCAAASSQAA5T0IWaAMEIAABCEAAAhCAAAQgAAEIeE0AOfc6/QQPAQhAAAIQgAAEIAABCEAAAkkggJwnIQu0AQIQgAAEIAABCEAAAhCAAAS8JoCce51+gocABCAAAQhAAAIQgAAEIACBJBBAzpOQBdoAAQhAAAIQgAAEIAABCEAAAl4TQM69Tj/BQwACEIAABCAAAQhAAAIQgEASCCDnScgCbYAABCAAAQhAAAIQgAAEIAABrwkg516nn+AhAAEIQAACEIAABCAAAQhAIAkEkPMkZIE2QAACEIAABCAAAQhAAAIQgIDXBJBzr9NP8BCAAAQgAAEIQAACEIAABCCQBALIeRKyQBsgAAEIQAACEIAABCAAAQhAwGsCyLnX6Sd4CEAAAhCAAAQgAAEIQAACEEgCAeQ8CVmgDRCAAAQgAAEIQAACEIAABCDgNQHk3Ov0EzwEIAABCEAAAhCAAAQgAAEIJIEAcp6ELNAGCEAAAhCAAAQgAAEIQAACEPCaAHLudfoJHgIQgAAEIAABCEAAAhCAAASSQAA5T0IWaAMEIAABCEAAAhCAAAQgAAEIeE0AOfc6/QQPAQhAAAIQgAAEIAABCEAAAkkggJwnIQu0AQIQgAAEIAABCEAAAhCAAAS8JoCce51+gocABCAAAQhAAAIQgAAEIACBJBBAzpOQBdoAAQhAAAIQgAAEIAABCEAAAl4TQM69Tj/BQwACEIAABCAAAQhAAAIQgEASCCDnScgCbYAABCAAAQhAAAIQgAAEIAABrwkg516nn+AhAAEIQAACEIAABCAAAQhAIAkEkPMkZIE2QAACEIAABCAAAQhAAAIQgIDXBJBzr9NP8BCAAAQgAAEIQAACEIAABCCQBALIeRKyQBsgAAEIQAACEIAABCAAAQhAwGsCyLnX6Sd4CEAAAhCAAAQgAAEIQAACEEgCAeQ8CVmgDRCAAAQgAAEIQAACEIAABCDgNQHk3Ov0EzwEIAABCEAAAhCAAAQgAAEIJIEAcp6ELNAGCEAAAhCAAAQgAAEIQAACEPCaAHLudfoJHgIQgAAEIAABCEAAAhCAAASSQAA5T0IWaAMEIAABCEAAAhCAAAQgAAEIeE0AOfc6/QQPAQhAAAIQgAAEIAABCEAAAkkggJwnIQu0AQIQgAAEIAABCEAAAhCAAAS8JoCce51+gocABCAAAQhAAAIQgAAEIACBJBBAzpOQBdoAAQhAAAIQgAAEIAABCEAAAl4TQM69Tj/BQwACEIAABCAAAQhAAAIQgEASCPx/3t+FIK6cnqQAAAAASUVORK5CYII=", "text/html": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Visualize the initial system state\n", "bio.visualize_system(title_prefix=\"Initial System State\")" ] }, { "cell_type": "code", "execution_count": null, "id": "702d3ab0-20f3-4776-a6aa-697f8adb4629", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "2ba3f95a-2d63-4a11-a150-c5f1ea4a7716", "metadata": {}, "source": [ "#### Now do 4 rounds of single-step diffusion, to collect the system state at a total of 5 time points: t0 (the initial state), plus t1, t2, t3 and t4" ] }, { "cell_type": "code", "execution_count": 9, "id": "28e46cc9-4ddf-4d9c-8ab2-20c08bb3106e", "metadata": {}, "outputs": [], "source": [ "# All the system state will get collected in this object\n", "history = CollectionArray()" ] }, { "cell_type": "code", "execution_count": 10, "id": "122c2273-3f40-4cd0-8ed6-77c03f0f6ed7", "metadata": {}, "outputs": [], "source": [ "# Store the initial state\n", "arr = bio.lookup_species(chem_index=0, copy=True)\n", "history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 11, "id": "45c018bc-5f69-4da2-b8c2-4c412e10bc5e", "metadata": {}, "outputs": [], "source": [ "# Do the 4 rounds of single-step diffusion; accumulate all data in the history object\n", "for _ in range(4):\n", " bio.diffuse(time_step=delta_t, n_steps=1, delta_x=delta_x , algorithm=algorithm)\n", "\n", " arr = bio.lookup_species(chem_index=0, copy=True)\n", " history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 12, "id": "a77b0423-fd2b-4219-8dc0-4f2fa65813be", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(5, 300)" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Now, let's examine the data collected at the 5 time points\n", "all_history = history.get_array()\n", "all_history.shape " ] }, { "cell_type": "code", "execution_count": 13, "id": "128ff22d-0ce2-45ed-8ec0-2a067a613990", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(5, 300)" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compute time derivatives (for each bin)\n", "df_dt_all_bins = np.apply_along_axis(np.gradient, 0, all_history, delta_t)\n", "df_dt_all_bins.shape" ] }, { "cell_type": "code", "execution_count": 14, "id": "f32bd8a9-b062-4db9-9e20-cd230464e6f8", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Let's consider the state at the midpoint in time (t2)\n", "f_at_t2 = all_history[2] # The middle of the 5 time snapshots\n", "f_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 15, "id": "398588e7-2066-4171-b95c-089b83126034", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Computer the second spatial derivative (simple method)\n", "gradient_x_at_t2 = np.gradient(f_at_t2, delta_x)\n", "second_gradient_x_at_t2 = np.gradient(gradient_x_at_t2, delta_x)\n", "second_gradient_x_at_t2.shape" ] }, { "cell_type": "markdown", "id": "8849a5a0-4a36-4e9a-b0d6-34823d147bed", "metadata": {}, "source": [ "#### Now, let's use the computed values to see how the left- and right-hand side of the diffusion equation compare\n", "at all bins (except 2 at each edge), at the middle simulation time t2" ] }, { "cell_type": "code", "execution_count": 16, "id": "8071dbaf-2d87-44e6-b0fc-b361d6094a35", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "lhs = df_dt_all_bins[2] # t2 is the middle point of the 5\n", "lhs.shape" ] }, { "cell_type": "code", "execution_count": 17, "id": "3909d680-2ea2-43ca-9094-bc370ee61c43", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "rhs = diffusion_rate*second_gradient_x_at_t2\n", "rhs.shape" ] }, { "cell_type": "code", "execution_count": 18, "id": "19514cd4-ceb5-450d-8c3e-9b43456f02c0", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.023163289760024783" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Numerical.compare_vectors(lhs, rhs, trim_edges=2) # Euclidean distance, ignoring 2 edge points at each end" ] }, { "cell_type": "markdown", "id": "258e81c8-0720-4a70-a7e0-af06167490d8", "metadata": {}, "source": [ "The above number is a measure of the discrepancy from the perfect match (zero distance) that an ideal solution would provide. \n", "\n", "Notice how vastly closer to zero it is than the counterpart value we got with the very coarse simulation in experiment \"validate_diffusion_1\"" ] }, { "cell_type": "markdown", "id": "5c069ce7-2e2b-4f13-bf2f-32509c865871", "metadata": {}, "source": [ "# VARIATIONS on the BASELINE\n", "We'll now tweak some parameters, and observe the effect on the estimate of the discrepancy from the exact diffusion equation\n", "\n", "The baseline distance was **0.023163289760024783**" ] }, { "cell_type": "markdown", "id": "3002d2cd-5142-483b-86db-126c3f73413f", "metadata": {}, "source": [ "## Variation 1 : reduce spatial resolution\n", "Expectation: greater discrepancy" ] }, { "cell_type": "code", "execution_count": 19, "id": "e3918a7f-10a6-4192-8e00-bfd1b7d6627e", "metadata": {}, "outputs": [], "source": [ "n_bins = 100 # Reducing the spatial resolution" ] }, { "cell_type": "code", "execution_count": 20, "id": "1523c8e7-ac8a-42d9-bc0a-9b9c6ea00f56", "metadata": {}, "outputs": [], "source": [ "# Initialize the system\n", "bio = BioSim1D(n_bins=n_bins, chem_data=chem_data)\n", "\n", "# Initialize the concentrations to 2 superposed sine waves\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=1, amplitude=10, bias=50)\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=2, amplitude=8)" ] }, { "cell_type": "code", "execution_count": 21, "id": "0b2fd9b1-ed93-440d-b67e-ab9a660bef7e", "metadata": {}, "outputs": [], "source": [ "history = CollectionArray() # All the system state will get collected in this object\n", "# Store the initial state\n", "arr = bio.lookup_species(chem_index=0, copy=True)\n", "history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 22, "id": "020e44ca-3dae-47f5-bdb2-377d7efe7ee2", "metadata": {}, "outputs": [], "source": [ "# Do the 4 rounds of single-step diffusion; accumulate all data in the history object\n", "for _ in range(4):\n", " bio.diffuse(time_step=delta_t, n_steps=1, delta_x=delta_x , algorithm=algorithm)\n", "\n", " arr = bio.lookup_species(chem_index=0, copy=True)\n", " history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 23, "id": "0aa5578c-4454-48f4-be62-5fe58ea7c5e1", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(5, 100)" ] }, "execution_count": 23, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Now, let's examine the data collected at the 5 time points\n", "all_history = history.get_array()\n", "all_history.shape " ] }, { "cell_type": "code", "execution_count": 24, "id": "03a64e4c-6b8a-4732-afa8-199028190340", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(100,)" ] }, "execution_count": 24, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compute time derivatives (for each bin) : simple method, using central differences\n", "df_dt_all_bins = np.apply_along_axis(np.gradient, 0, all_history, delta_t)\n", "\n", "# Let's consider the state at the midpoint in time (t2)\n", "f_at_t2 = all_history[2] # The middle of the 5 time snapshots\n", "f_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 25, "id": "4998952c-d051-4d3e-83c2-f443d3e67e47", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(100,)" ] }, "execution_count": 25, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Computer the second spatial derivative (simple method, using central differences)\n", "gradient_x_at_t2 = np.gradient(f_at_t2, delta_x)\n", "second_gradient_x_at_t2 = np.gradient(gradient_x_at_t2, delta_x)\n", "second_gradient_x_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 26, "id": "885a0f51-eaa9-43f8-ae5e-fe8a6bd7bcc9", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.0705630110277808" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compare the left and right hand sides of the diffusion equation\n", "lhs = df_dt_all_bins[2] # t2 is the middle point of the 5\n", "rhs = diffusion_rate*second_gradient_x_at_t2\n", "\n", "Numerical.compare_vectors(lhs, rhs, trim_edges=2) # Euclidean distance, ignoring 2 edge points at each end" ] }, { "cell_type": "markdown", "id": "76a989a9-311c-4e58-9358-5a51b88f705e", "metadata": {}, "source": [ "#### The discrepancy value has indeed gone up from the baseline 0.023163289760024783" ] }, { "cell_type": "code", "execution_count": null, "id": "71eb4699-52f1-4e6b-bc8e-d3f560e65d36", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "976f3b8d-bf0e-4a23-978f-0f2064d654f6", "metadata": {}, "source": [ "## Variation 2 : increase spatial resolution\n", "Expectation: smaller discrepancy" ] }, { "cell_type": "code", "execution_count": 27, "id": "aab88c73-3afa-48f6-959e-1e9043350cc6", "metadata": {}, "outputs": [], "source": [ "n_bins = 900 # Increasing the spatial resolution (from the baseline of 300)" ] }, { "cell_type": "code", "execution_count": 28, "id": "e996032a-63bd-4d41-8cce-b97a158c82b0", "metadata": {}, "outputs": [], "source": [ "# Initialize the system\n", "bio = BioSim1D(n_bins=n_bins, chem_data=chem_data)\n", "\n", "# Initialize the concentrations to 2 superposed sine waves\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=1, amplitude=10, bias=50)\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=2, amplitude=8)" ] }, { "cell_type": "code", "execution_count": 29, "id": "6e26a594-77d6-49d8-9e35-ad2e5d03bdff", "metadata": {}, "outputs": [], "source": [ "history = CollectionArray() # All the system state will get collected in this object\n", "# Store the initial state\n", "arr = bio.lookup_species(chem_index=0, copy=True)\n", "history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 30, "id": "75f95eee-5e93-41a4-8d6a-6663b455faa5", "metadata": {}, "outputs": [], "source": [ "# Do the 4 rounds of single-step diffusion; accumulate all data in the history object\n", "for _ in range(4):\n", " bio.diffuse(time_step=delta_t, n_steps=1, delta_x=delta_x , algorithm=algorithm)\n", "\n", " arr = bio.lookup_species(chem_index=0, copy=True)\n", " history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 31, "id": "2e75b83d-ad36-441e-be9e-e4c7321eb1b4", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(5, 900)" ] }, "execution_count": 31, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Now, let's examine the data collected at the 5 time points\n", "all_history = history.get_array()\n", "all_history.shape " ] }, { "cell_type": "code", "execution_count": 32, "id": "1c5accca-8be9-4410-8cf7-5122774e7306", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(900,)" ] }, "execution_count": 32, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compute time derivatives (for each bin) : simple method, using central differences\n", "df_dt_all_bins = np.apply_along_axis(np.gradient, 0, all_history, delta_t)\n", "\n", "# Let's consider the state at the midpoint in time (t2)\n", "f_at_t2 = all_history[2] # The middle of the 5 time snapshots\n", "f_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 33, "id": "d103fc3f-1578-4d7f-88de-87a7fc896c5a", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(900,)" ] }, "execution_count": 33, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Computer the second spatial derivative (simple method, using central differences)\n", "gradient_x_at_t2 = np.gradient(f_at_t2, delta_x)\n", "second_gradient_x_at_t2 = np.gradient(gradient_x_at_t2, delta_x)\n", "second_gradient_x_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 34, "id": "fab3cfc4-407f-4db0-af21-dda78212933f", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.007721522026370068" ] }, "execution_count": 34, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compare the left and right hand sides of the diffusion equation\n", "lhs = df_dt_all_bins[2] # t2 is the middle point of the 5\n", "rhs = diffusion_rate*second_gradient_x_at_t2\n", "\n", "Numerical.compare_vectors(lhs, rhs, trim_edges=2) # Euclidean distance, ignoring 2 edge points at each end" ] }, { "cell_type": "markdown", "id": "5fbf609b-e96a-4281-bfa2-f2f0f43fae72", "metadata": {}, "source": [ "#### The discrepancy value has indeed gone down from the baseline 0.023163289760024783" ] }, { "cell_type": "code", "execution_count": null, "id": "83f1a02a-13c2-4fda-8fff-37d166a23de5", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "27746f62-229e-4915-9b1d-4dac6190eade", "metadata": {}, "source": [ "## Variation 3 : reduce time resolution\n", "Expectation: greater discrepancy" ] }, { "cell_type": "code", "execution_count": 35, "id": "a143cba9-2e3c-4b0a-a357-7eca88cc454f", "metadata": {}, "outputs": [], "source": [ "n_bins = 300 # Restoring the baseline number of bins\n", "delta_t = 0.02 # Reducing the time resolution (from baselin 0.01)" ] }, { "cell_type": "code", "execution_count": 36, "id": "3fd3d708-dc6f-4387-b49a-090813c52d4f", "metadata": {}, "outputs": [], "source": [ "# Initialize the system\n", "bio = BioSim1D(n_bins=n_bins, chem_data=chem_data)\n", "\n", "# Initialize the concentrations to 2 superposed sine waves\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=1, amplitude=10, bias=50)\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=2, amplitude=8)" ] }, { "cell_type": "code", "execution_count": 37, "id": "cf1d4c86-17be-4c53-b7d6-c14fc80a1604", "metadata": {}, "outputs": [], "source": [ "history = CollectionArray() # All the system state will get collected in this object\n", "# Store the initial state\n", "arr = bio.lookup_species(chem_index=0, copy=True)\n", "history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 38, "id": "b31a59ae-8d90-4cd1-beb6-0a2c518502b2", "metadata": {}, "outputs": [], "source": [ "# Do the 4 rounds of single-step diffusion; accumulate all data in the history object\n", "for _ in range(4):\n", " bio.diffuse(time_step=delta_t, n_steps=1, delta_x=delta_x , algorithm=algorithm)\n", "\n", " arr = bio.lookup_species(chem_index=0, copy=True)\n", " history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 39, "id": "50f13dcb-759d-410d-b7a3-925e475a927c", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(5, 300)" ] }, "execution_count": 39, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Now, let's examine the data collected at the 5 time points\n", "all_history = history.get_array()\n", "all_history.shape " ] }, { "cell_type": "code", "execution_count": 40, "id": "e2f93e1b-4dd9-4163-96a7-30c6552fce91", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 40, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compute time derivatives (for each bin) : simple method, using central differences\n", "df_dt_all_bins = np.apply_along_axis(np.gradient, 0, all_history, delta_t)\n", "\n", "# Let's consider the state at the midpoint in time (t2)\n", "f_at_t2 = all_history[2] # The middle of the 5 time snapshots\n", "f_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 41, "id": "149b6566-22e0-4994-82d4-adc67c4fb08e", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 41, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Computer the second spatial derivative (simple method, using central differences)\n", "gradient_x_at_t2 = np.gradient(f_at_t2, delta_x)\n", "second_gradient_x_at_t2 = np.gradient(gradient_x_at_t2, delta_x)\n", "second_gradient_x_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 42, "id": "9c62a79a-3188-4ffc-9b79-abd2a69e12fa", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.04453006107638132" ] }, "execution_count": 42, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compare the left and right hand sides of the diffusion equation\n", "lhs = df_dt_all_bins[2] # t2 is the middle point of the 5\n", "rhs = diffusion_rate*second_gradient_x_at_t2\n", "\n", "Numerical.compare_vectors(lhs, rhs, trim_edges=2) # Euclidean distance, ignoring 2 edge points at each end" ] }, { "cell_type": "markdown", "id": "c79b51f0-da17-47cf-b4bb-e423951e121a", "metadata": {}, "source": [ "#### The discrepancy value has indeed gone up from the baseline 0.023163289760024783" ] }, { "cell_type": "code", "execution_count": null, "id": "d4675241-060c-42f2-b2ba-8f7811a2f3b6", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "4ff62063-0f8b-4d28-8b47-44bcf0396874", "metadata": {}, "source": [ "## Variation 4 : increase time resolution\n", "Expectation: smaller discrepancy" ] }, { "cell_type": "code", "execution_count": 43, "id": "967b3cb4-7d89-4114-86bb-8dfe76f6cb71", "metadata": {}, "outputs": [], "source": [ "delta_t = 0.005 # Increasing the time resolution (from baselin 0.01)" ] }, { "cell_type": "code", "execution_count": 44, "id": "4ae7056d-6c3d-4db4-93ca-c9688fef52f7", "metadata": {}, "outputs": [], "source": [ "# Initialize the system\n", "bio = BioSim1D(n_bins=n_bins, chem_data=chem_data)\n", "\n", "# Initialize the concentrations to 2 superposed sine waves\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=1, amplitude=10, bias=50)\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=2, amplitude=8)" ] }, { "cell_type": "code", "execution_count": 45, "id": "81a3bd58-8ef4-43fb-9477-b27ab0160d2c", "metadata": {}, "outputs": [], "source": [ "history = CollectionArray() # All the system state will get collected in this object\n", "# Store the initial state\n", "arr = bio.lookup_species(chem_index=0, copy=True)\n", "history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 46, "id": "e32ba5c9-185c-4d84-8c6d-415a55cf98e7", "metadata": {}, "outputs": [], "source": [ "# Do the 4 rounds of single-step diffusion; accumulate all data in the history object\n", "for _ in range(4):\n", " bio.diffuse(time_step=delta_t, n_steps=1, delta_x=delta_x , algorithm=algorithm)\n", "\n", " arr = bio.lookup_species(chem_index=0, copy=True)\n", " history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 47, "id": "b92026f5-82e8-4d4d-b376-a9e9e46d0c00", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(5, 300)" ] }, "execution_count": 47, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Now, let's examine the data collected at the 5 time points\n", "all_history = history.get_array()\n", "all_history.shape " ] }, { "cell_type": "code", "execution_count": 48, "id": "4fdb6d6e-df25-4a1d-a886-0125f82e897c", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 48, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compute time derivatives (for each bin) : simple method, using central differences\n", "df_dt_all_bins = np.apply_along_axis(np.gradient, 0, all_history, delta_t)\n", "\n", "# Let's consider the state at the midpoint in time (t2)\n", "f_at_t2 = all_history[2] # The middle of the 5 time snapshots\n", "f_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 49, "id": "8db28e0f-b63d-43d9-8913-87c62c87739e", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 49, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Computer the second spatial derivative (simple method, using central differences)\n", "gradient_x_at_t2 = np.gradient(f_at_t2, delta_x)\n", "second_gradient_x_at_t2 = np.gradient(gradient_x_at_t2, delta_x)\n", "second_gradient_x_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 50, "id": "01e1ec43-d7b9-4b6b-a57e-2026c400c0ff", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.011808914990091101" ] }, "execution_count": 50, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compare the left and right hand sides of the diffusion equation\n", "lhs = df_dt_all_bins[2] # t2 is the middle point of the 5\n", "rhs = diffusion_rate*second_gradient_x_at_t2\n", "\n", "Numerical.compare_vectors(lhs, rhs, trim_edges=2) # Euclidean distance, ignoring 2 edge points at each end" ] }, { "cell_type": "markdown", "id": "f03489cb-5a2c-490f-a48d-50b3cbe75edd", "metadata": {}, "source": [ "#### The discrepancy value has indeed gone down from the baseline 0.023163289760024783" ] }, { "cell_type": "code", "execution_count": null, "id": "46d94d17-d4e8-4938-b0ac-f43a6f7bf5e8", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "2e42b48c-0c11-4a33-80d7-4c51bbf1fac3", "metadata": {}, "source": [ "## Variation 5 : Use a better (higher-order) method to numerically estimate the SPATIAL derivatives\n", "Expectation: presumably smaller discrepancy (unless dwarfed by errors from approximations in the simulation)" ] }, { "cell_type": "code", "execution_count": 51, "id": "69401e31-a5f6-4faf-b23f-024a73a6e805", "metadata": {}, "outputs": [], "source": [ "n_bins = 300 # Restoring the baseline number of bins\n", "delta_t = 0.01 # Reducing the baseline time resolution" ] }, { "cell_type": "code", "execution_count": 52, "id": "f732c9a0-4b95-4bcf-821e-0221ba7170e9", "metadata": {}, "outputs": [], "source": [ "# Initialize the system\n", "bio = BioSim1D(n_bins=n_bins, chem_data=chem_data)\n", "\n", "# Initialize the concentrations to 2 superposed sine waves\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=1, amplitude=10, bias=50)\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=2, amplitude=8)" ] }, { "cell_type": "code", "execution_count": 53, "id": "98b3b3bf-3ed7-4867-895b-3ae39bfb0424", "metadata": {}, "outputs": [], "source": [ "history = CollectionArray() # All the system state will get collected in this object\n", "# Store the initial state\n", "arr = bio.lookup_species(chem_index=0, copy=True)\n", "history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 54, "id": "38387573-8c4c-4b8f-9e28-f5d3196cc3c6", "metadata": {}, "outputs": [], "source": [ "# Do the 4 rounds of single-step diffusion; accumulate all data in the history object\n", "for _ in range(4):\n", " bio.diffuse(time_step=delta_t, n_steps=1, delta_x=delta_x , algorithm=algorithm)\n", "\n", " arr = bio.lookup_species(chem_index=0, copy=True)\n", " history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 55, "id": "27729fc9-a692-43ee-9ac3-fef64ade8c0e", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(5, 300)" ] }, "execution_count": 55, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Now, let's examine the data collected at the 5 time points\n", "all_history = history.get_array()\n", "all_history.shape " ] }, { "cell_type": "code", "execution_count": 56, "id": "03150f7c-590e-4531-bea4-e18d8df7ffa3", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 56, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compute time derivatives (for each bin) : simple method, using central differences\n", "df_dt_all_bins = np.apply_along_axis(np.gradient, 0, all_history, delta_t)\n", "\n", "# Let's consider the state at the midpoint in time (t2)\n", "f_at_t2 = all_history[2] # The middle of the 5 time snapshots\n", "f_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 57, "id": "dc283a71-7228-486f-8cbc-09d1a2e51c68", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 57, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Computer the second spatial derivative (MORE SOPHISTICATED METHOD, using 5-point stencils)\n", "gradient_x_at_t2 = Numerical.gradient_order4_1d(arr=f_at_t2, dx=delta_x)\n", "second_gradient_x_at_t2 = Numerical.gradient_order4_1d(arr=gradient_x_at_t2, dx=delta_x)\n", "second_gradient_x_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 58, "id": "92b9fc7a-659f-4a62-a0a4-5fe8de345d0e", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.006831207619477847" ] }, "execution_count": 58, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compare the left and right hand sides of the diffusion equation\n", "lhs = df_dt_all_bins[2] # t2 is the middle point of the 5\n", "rhs = diffusion_rate*second_gradient_x_at_t2\n", "\n", "Numerical.compare_vectors(lhs, rhs, trim_edges=2) # Euclidean distance, ignoring 2 edge points at each end" ] }, { "cell_type": "markdown", "id": "abdddfd0-e02b-46fc-a3b4-2b9581740a0b", "metadata": {}, "source": [ "#### Not surprisingly, the discrepancy value (in part caused by poor numeric estimation of spatial derivatives) has indeed gone down from the baseline 0.023163289760024783" ] }, { "cell_type": "code", "execution_count": null, "id": "d54be0a1-f55c-4785-a0fc-22d89e47c0bc", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "38642e67-a039-416d-9945-c0a020d26490", "metadata": {}, "source": [ "## Variation 6 : Use a better (higher-order) method to numerically estimate the TIME derivatives\n", "Expectation: presumably smaller discrepancy (unless dwarfed by errors from approximations in the simulation)\n", "\n", "(Note: we revert to the original method for the *spatial* derivatives)" ] }, { "cell_type": "code", "execution_count": 59, "id": "ba76ad1c-97e7-463c-ac2c-30a283db42dc", "metadata": {}, "outputs": [], "source": [ "# Initialize the system\n", "bio = BioSim1D(n_bins=n_bins, chem_data=chem_data)\n", "\n", "# Initialize the concentrations to 2 superposed sine waves\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=1, amplitude=10, bias=50)\n", "bio.inject_sine_conc(chem_label=\"A\", number_cycles=2, amplitude=8)" ] }, { "cell_type": "code", "execution_count": 60, "id": "9d0a38bf-6c59-4dee-87e8-6aa2e149fb47", "metadata": {}, "outputs": [], "source": [ "history = CollectionArray() # All the system state will get collected in this object\n", "# Store the initial state\n", "arr = bio.lookup_species(chem_index=0, copy=True)\n", "history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 61, "id": "47078710-eb53-44c5-83bc-406dbbd5173b", "metadata": {}, "outputs": [], "source": [ "# Do the 4 rounds of single-step diffusion; accumulate all data in the history object\n", "for _ in range(4):\n", " bio.diffuse(time_step=delta_t, n_steps=1, delta_x=delta_x , algorithm=algorithm)\n", "\n", " arr = bio.lookup_species(chem_index=0, copy=True)\n", " history.store(par=bio.system_time, data_snapshot=arr, caption=f\"State at time {bio.system_time}\")" ] }, { "cell_type": "code", "execution_count": 62, "id": "f8d28717-210e-4adc-8d53-ceeb84b3bb26", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(5, 300)" ] }, "execution_count": 62, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Now, let's examine the data collected at the 5 time points\n", "all_history = history.get_array()\n", "all_history.shape " ] }, { "cell_type": "code", "execution_count": 63, "id": "028d732d-f97b-4103-ae94-85a60cfce4dc", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 63, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compute time derivatives (for each bin) : MORE SOPHISTICATED METHOD, using 5-point stencils\n", "df_dt_all_bins = np.apply_along_axis(Numerical.gradient_order4_1d, 0, all_history, delta_t)\n", "\n", "# Let's consider the state at the midpoint in time (t2)\n", "f_at_t2 = all_history[2] # The middle of the 5 time snapshots\n", "f_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 64, "id": "ca637810-beda-4abf-b232-473b17e1a76d", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(300,)" ] }, "execution_count": 64, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Computer the second spatial derivative (BACK TO SIMPLE METHOD, using central differences)\n", "gradient_x_at_t2 = np.gradient(f_at_t2, delta_x)\n", "second_gradient_x_at_t2 = np.gradient(gradient_x_at_t2, delta_x)\n", "second_gradient_x_at_t2.shape" ] }, { "cell_type": "code", "execution_count": 65, "id": "845f75bd-25fe-4120-8975-e8dcd8e36afc", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.023351269997869635" ] }, "execution_count": 65, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Compare the left and right hand sides of the diffusion equation\n", "lhs = df_dt_all_bins[2] # t2 is the middle point of the 5\n", "rhs = diffusion_rate*second_gradient_x_at_t2\n", "\n", "Numerical.compare_vectors(lhs, rhs, trim_edges=2) # Euclidean distance, ignoring 2 edge points at each end" ] }, { "cell_type": "markdown", "id": "fd1d0295-dfe2-4e7d-9457-cc2fea1a98f2", "metadata": {}, "source": [ "#### Here, the discrepancy value has remained largely the same from the baseline 0.023163289760024783" ] }, { "cell_type": "code", "execution_count": null, "id": "c75683c7", "metadata": {}, "outputs": [], "source": [] } ], "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.9.13" }, "toc-autonumbering": false, "toc-showcode": true, "toc-showmarkdowntxt": false }, "nbformat": 4, "nbformat_minor": 5 }