{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# 8. More about 1. Brief Tour of E-Cell4 Simulations\n", "\n", "Once you read through [1. Brief Tour of E-Cell4 Simulations](tutorial1.ipynb), it is NOT difficult to use `World` and `Simulator`.\n", "`volume` and `{'C': 60}` is equivalent of the `World` and solver is the `Simulator` below." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3deXxU9b3/8ddnspMEwhLWCAiigOxEwLWtG1ot0Fvbq/Vaqt5627q07ty2Lr2/tnq91WoXrdalWK24tlq1pYpr3QEBZadsQZYECCQhkGXm8/tjJhgphCFk5pCZ97OPeZw5J2d5D9JPDt/5nu/X3B0REUkfoaADiIhIcqnwi4ikGRV+EZE0o8IvIpJmVPhFRNJMZtAB4tGtWzfv379/0DFERNqVOXPmbHb34j23t4vC379/f2bPnh10DBGRdsXM1uxtu5p6RETSjAq/iEiaUeEXEUkzKvwiImlGhV9EJM0ktPCbWZGZPWVmS8xssZkda2ZdzOwlM1seW3ZOZAYREfmsRN/x3wX8zd0HAyOBxcA0YJa7DwJmxdZFRCRJElb4zawjcBLwAIC717v7NmAyMD2223RgSqIyvLzmZZ5Y+kSiTi8i0i4l8o5/AFABPGRmH5rZ/WaWD/Rw9w0AsWX3vR1sZpeY2Wwzm11RUdGqAC+uepHbZ99O5a7KVn4EEZHUk8jCnwmMAe5x99HADg6gWcfd73P3UncvLS7+lyeO43LpqEvZ2biTBz56oFXHi4ikokQW/nXAOnd/L7b+FNFfBJvMrBdAbFmeqAADiwbypYFfYsbSGZTXJuwyIiLtSsIKv7tvBMrM7KjYplOARcBzwNTYtqnAs4nKAPCdkd8h7GHuW3BfIi8jItJuJLpXz+XAo2a2ABgF/Ay4FTjNzJYDp8XWE6aksISvDPoKTy97mrLqskReSkSkXUho4Xf3ebF2+hHuPsXdK919i7uf4u6DYsuticwAcMmIS8gIZXDPvHsSfSkRkUNeWjy5271Dd74++Os8v/J5llUuCzqOiEig0qLwA1w8/GIKsgu4a+5dQUcREQlU2hT+TjmduHjYxbyx7g0+2PhB0HFERAKTNoUf4Pwh59O9Q3funHMn7h50HBGRQKRV4c/NzOXSUZeyYPMCZq2dFXQcEZFApFXhB5g0cBIDOg3grrl30RBpCDqOiEjSpV3hzwxlcuXYK1ldtZonlz4ZdBwRkaRLu8IP8LmSzzGu5zjumX8PVfVVQccREUmqtCz8ZsY1pdewvW479y+4P+g4IiJJlZaFH2BI1yFMGjiJRxY/wrrqdUHHERFJmrQt/ACXj76czFAmd869M+goIiJJk9aFv0d+D7559DeZuXomczbNCTqOiEhSpHXhB7hw2IX0zO/Jre/fSjgSDjqOiEjCpX3hz8vM4+qxV7Nk6xKeWfFM0HFERBIu7Qs/wMT+ExnbYyy/mvsrttdtDzqOiEhCqfAT7d45bdw0ttdv5575GrNfRFKbCn/M4C6DOWfQOcxYMoOlW5cGHUdEJGFU+Ju5fPTlFGYX8tP3fkrEI0HHERFJCBX+Zopyi7hq7FV8WP4hz/3zuaDjiIgkhAr/HiYfMZlRxaO4Y/Yd+qJXRFKSCv8eQhbiRxN+xPb67fxy7i+DjiMi0uZU+PfiqC5H8fXBX+fJZU8yr3xe0HFERNqUCv8+XDb6Mnrk9+DH7/yYhrAmbBGR1KHCvw/5WfncMOEGVmxbwQMfPxB0HBGRNpPQwm9mq83sIzObZ2azY9u6mNlLZrY8tuycyAwH46SSk5jYfyL3LbiPldtXBh1HRKRNJOOO/wvuPsrdS2Pr04BZ7j4ImBVbP2RNGzeN3Mxcfvz2j9W3X0RSQhBNPZOB6bH304EpAWSIW7e8blxbei1zy+cyY8mMoOOIiBy0RBd+B/5uZnPM7JLYth7uvgEgtuy+twPN7BIzm21msysqKhIcs2VTjpjC8X2O5865d1JWVRZoFhGRg5Xown+8u48BzgQuNbOT4j3Q3e9z91J3Ly0uLk5cwjiYGTcfezOZlskNb9+gJh8RadcSWvjdfX1sWQ78CRgHbDKzXgCxZXkiM7SVnvk9ufaYa5mzaQ6PLXks6DgiIq2WsMJvZvlmVtj0Hjgd+Bh4Dpga220q8GyiMrS1KUdM4cQ+J3LnnDtZtX1V0HFERFolkXf8PYB/mNl84H3gBXf/G3ArcJqZLQdOi623C2bGzcfdTE5mDtPenEZDRA92iUj7s9/Cb2ZfbXbn/iMze8bMxuzvOHdf6e4jY6+j3f2nse1b3P0Udx8UW249+I+RPN07dOfmY29m0ZZF3DNPk7aISPsTzx3/De5ebWYnABOJdsFM64p3ar9T+fIRX+b+j+5nzqY5QccRETkg8RT+cGx5FnCPuz8LZCcuUvtw/bjrKSks4Qdv/oDq+uqg44iIxC2ewv+Jmd0LfA140cxy4jwupeVn5XPLibewqXYTN719E+4edCQRkbjEU8C/BswEznD3bUAX4NqEpmonRhaP5IoxV/DSmpd4fOnjQccREYnLfgu/u9cS7Wt/QmxTI7A8kaHak28e/U1O6HMCt31wG0u2Lgk6jojIfsXTq+cm4Hrgv2ObsoBHEhmqPQlZiJ+e8FM653Tmmtevoaa+JuhIIiItiqep58vAJGAH7H4atzCRodqbLrlduO1zt7Gueh0//McPNaSDiBzS4in89R795tJh91O4soexPcZy1direKXsFR78+MGg44iI7FM8hf+JWK+eIjP7FvAy8LvExmqfLhh6AWf2P5Nfzv0lb3/ydtBxRET2Kp4vd38OPAU8DRwF3Ojuv0p0sPaoaUiHgUUDue7N6yir1hDOInLoias/vru/5O7Xuvs17v5SokO1Zx2yOnDXF+7C3bl81uV6uEtEDjn7LPxmVm1mVc2WVc3XkxmyvenbsS+/+PwvWFO1hmvfuJbGSGPQkUREdttn4Xf3Qnfv2GzZsfl6MkO2R+N6jeMHE37AW5+8xe2zbw86jojIbpn728HM+u5tu7uvbfs4qeWrR36VldtW8sjiR+jbsS/nDT4v6EgiIvsv/MALzd7nAocDS4GjE5IoxVxTeg3ratZxy3u30D2vO6f0OyXoSCKS5uLp1TO82WsQ0ekT/5H4aKkhI5TBbSfdxvDi4Vz3xnXM3TQ36EgikuYOeJRNd58LHJOALCkrLzOPX5/8a3oX9ObyVy5nReWKoCOJSBqLZ6yeq5q9rjGzPwIVSciWUjrnduaeU+8hJyOHS166hLIq9fEXkWDEc8df2OyVQ7TNf3IiQ6WqksISfnf672iINPCff/9PNu7YGHQkEUlD1h4mECktLfXZs2cHHaPNLNqyiItnXky3vG48dMZDdMvrFnQkEUlBZjbH3Uv33B5PU89LZlbUbL2zmc1s64DpZGjXodx96t1sqt3ERTMvoqJWLWcikjzxNPUUx2beAsDdK4HuiYuUHkZ3H83dp9zNxh0buWjmRWzasSnoSCKSJuKabL35Q1xm1o/YEM1ycEp7lnLvafdSXlvORTMvUpu/iCRFPIX/h8A/zOwPZvYH4A0+nY1LDtLo7qO597R72bprKxf89QJWbV8VdCQRSXHxPMD1N2AM8DjwBDDW3eNu4zezDDP70Myej60fbmbvmdlyM3vczLJbGz5VjOo+igcnPkh9uJ6pf53Kwi0Lg44kIiksni93DTgDGOPufwE6mNm4A7jG94DFzdb/F/hF7CngSuDiAzhXyhrSdQgPn/kweZl5XPS3i3hn/TtBRxKRFBVPU8/dwLFA0whj1cBv4jm5mZUAZwH3x9YNOJnoxC4A04EpB5A3pfXr2I+Hz3yY3gW9+e7L3+VPy/8UdCQRSUHxFP7x7n4psAt29+qJt3nmTuA6oGn28a7ANndvGqB+HdBnbwea2SVmNtvMZldUpE93xx75PXj4zIcZ12scN759I3fNvUuTt4tIm4qn8DeYWQafTrZezKeFfJ/M7Gyg3N3nNN+8l1332kPI3e9z91J3Ly0uLo4jZuoozC7k16f8mnOOPIf7P7qfq1+7mh0NO4KOJSIpIp7C/0vgT0B3M/sp0ZE5fxbHcccDk8xsNTCDaBPPnUQnbW8aDroEWH+godNBViiLGyfcyLWl1/Jq2auc/8L5rK3SFAgicvDi6dXzKNHmmluADcAUd38yjuP+291L3L0/cC7wirufD7wKnBPbbSrwbCuzpzwz4xtHf4N7T7uXLbu2cO4L5/J62etBxxKRdq6lOXe7NL2AcuAx4I/Apti21roeuMrMVhBt83/gIM6VFsb3Gs+Ms2dQUlDCZa9cxu2zb6ch0hB0LBFpp/Y5SJuZrSLa/r7Xdnl3H5DIYM2l2iBtrVUXruP/Pvg/Hl/6OCO6jeC2z91Gn4K9fjcuIrLPQdo0Omc7NHP1TG5++2YcZ9q4aUweOJloT1kRkU+1enTO2MGTzOznsdfZbR9PDsTE/hN5atJTDO4ymBveuoErX7uSrbu2Bh1LRNqJeJ7cvZXo07eLYq/vmdktiQ4mLetT0IcHTn+Aq8dezRvr3mDKn6fw/MrnaQ//ghORYO23qcfMFgCj3KNPEcX69H/o7iOSkA9QU8/+rKhcwU3v3MSCigUc3+d4bphwg9r+ReTgmnqAombvO7VNJGkrR3Q+gofPeJhp46Yxd9Ncpvx5Cr+d/1vqwnVBRxORQ1A8hf8W4EMz+72ZTQfmEN8DXJJEGaEMzh9yPs9NeY6TSk7iN/N+w5Q/T+GVta+o+UdEPiOuXj1m1gs4hmjXzvfcPakzhqip58C9u+FdbnnvFlZuX8nYHmO5pvQahnUbFnQsEUmiA+7OaWZjWjqhu89to2z7pcLfOg2RBp5Z9gx3z7+brbu2MrH/RL478rsMKEraIxgiEqDWFP4IsBBoGhqzeUdxd/eT2zzlPqjwH5wdDTt46OOH+MOiP7CzcSdfHPBFvj3i2/Tv1D/oaCKSQK0p/FcCXwG2Ex1k7U/uXpPQlPugwt82KndV8tDCh5ixZAa7GndxWr/TuHj4xQztOjToaCKSAK1+ctfMDic6CctkYA3wM3efl5CU+6DC37Y279zMo4sfZcaSGdQ01DCh1wQuGHoBJ/Q5gZDF29FLRA51BzVkg5kdTXSEzQuA69z9ibaPuG8q/IlRXV/NE0uf4I+L/0j5znL6dezHeYPP4+wBZ9MpR712Rdq71jT1DCBa7CcDZUSbe553912JDLo3KvyJ1RBp4KXVL/Ho4kdZsHkBORk5TOw/kXOOPIdRxaM0DpBIO9XaL3cXEB0vv4o9Zspy9zsSkHOvVPiTZ9GWRTy97GleWPUCOxp20LewL18a+CXOHnA2JYUlQccTkQPQmsJ/M/uYFhHA3X/cZun2Q4U/+Wobavn7mr/zl3/+hfc3vg/A8G7Dmdh/IhP7T6Rnfs+AE4rI/mhYZmm19TXr+euqvzJz9UwWb10MwLCuwzi578mc3PdkBnQaoOYgkUOQCr+0idXbV/Py2pd5de2rLNi8AIDe+b05seRETuhzAsf0PIb8rPyAU4oIqPBLApTXlvP6utd5c92bvLvhXXY27iTTMhnWbRjje42ntGcpI4tHkpeZF3RUkbTUmjb+77n7XWZ2vLu/lfCELVDhP/TVh+v5sPxD3tvwHu9ueJeFWxYS8QiZlsnQbkMZVTyKkcUjGVk8kh75PYKOK5IWWlP457n7KDOb6+4tjtuTaCr87U91fTXzyucxZ9Mc5pbPZeHmhdRH6gHontedod2GcnTXoxnSZQhHdTmKHh166HsCkTa2r8Kf2cIxi81sNVAcm4xl97mIjtWTtIlYpP0pzC7kxJITObHkRAAawg0s3rqYBRULWLhlIR9v/pjXyl7bvX/nnM4M6jyII4qO4IjORzCg0wAGdBpA59zOAX0CkdS1z8Lv7ueZWU9gJjApeZEkFWVlZDGieAQjij+9X6ipr2H5tuUs2bqEpVuXsrxyOX9e8WdqG2t371OUU0S/jv12vw4rPIySghJKCksoyinSvxJEWiHeIRuygSNjq0vdvSGhqfagpp70EfEIG3ZsYOW2lazavoqV21eytnota6rWUF5b/pl98zLz6J3fm14FveiV34seHXrQM78n3Tt03/0qyCrQLwdJW61p6mk68HPAw8Bqos08h5nZVHd/o81TStoLWYg+BX3oU9BndzNRk9qGWtbXrGddzTrKqstYX7Oe9TXr2bBjAws3L6SyrvJfzpebkUvXvK50y+tG19yudM7tTJfcLnTJ7UJRbhGdczpTlFNEp5xOdMrppF8Ukhb2W/iBO4DT3X0pgJkdCTwGjG3pIDPLBd4AcmLXecrdb4qN9jkD6ALMBS5w9/rWfwRJFx2yOnBE5+h3AHtTF65j045NlNeWU7GzgvLacjbv3MzmnZup2FlBWU0Z8yvmU1lXScQjez1HhmVQmF1Ix+yOdMzuSEF2AYXZhRRkFZCflU9BdgEFWQXkZeaRn5VPflY+HTI70CGrA3mZeZ955WTk6JeIHJLiKfxZTUUfwN2XmVlWHMfVASe7e01s/3+Y2V+Bq4BfuPsMM/stcDFwT2vCizSXk5FD34596duxb4v7RTxCdX01lbsqqayrZNuubWyv3872uuirqr6KqroqqhqqqKmvoby2nJr6GnY07mBHw44DytT0CyAnI+cz77MzssnJzCE7lL17PSuURXZGNtmhbLIysnYvMy3zs8tQJpmWGV02f8W2ZYQyyLToMsNir1AGIQuRaZmELLR7vennIQvtXjez3UtJTfEU/tlm9gDwh9j6+UQnXG+RR788aJq4JSv2cuBk4Oux7dOBm1HhlyQKWWh3005/+h/QseFImNrGWmobatnRuIPahlp2Nu6ktqGW2sbo+6bXrsZd1IXr2Nm4k7pwHXWNdewK76I+XE9duI6quirqw/XUR+qjy9j2hkgDDeEGGr0xMX8AByCE7f5f9H3zbTRbsvdtbrHzRDX/Of7ptH62xwvAmn39+Jntu98b4JjHzofHzt7Em63F3se+02z+K832WO7e3sLXn3s7fm/naPHYfZx/z3Ncf8qvOHrgMXGcOX7xFP7vAJcCV8QyvQHcHc/JzSyD6C+JI4DfAP8Etrnv/hu9Duizj2MvAS4B6Nu35Ts4kWTJCEWbggqzC4mEw+ysrWZnYxW76qqoq82iYWcGDXVG4y4nXBcmUt9IpD4EDeCNDg2OhcNYYyOhcCOhcAOhSAMZu19hMr2RDA8T8kbcGsEbwcJErBHzMG4RCEVoxGg0aDSjESNs0BBbRohuDwNhMxqBiFlsO0T49GcRIGIQxnAgbNGaFMaINHvvsf0iRMto03FObJtF94nuH+XNztH8FWl2PqfpOIudO/oewO3T9d3HW1PN/HQbzc695/bYpmb7Nf/V8Oky0rTjHsf6p5v/5Wd7W99T8181+9q3pXPsqqtt4aets9/C7+51RNv5D3gYZncPA6PMrAj4EzBkb7vt49j7gPsg2qvnQK8tsj+RcJiqygqqKzexY1sFddWbqa/eSri2Et+5Ddu1nVB9NVkNVWQ11pAdriUvXEOu76SD7yTfdpEPHMjIRPWeyS7Lpp5sGsiiIZRNg2UTtiwaLZuGjFzqQoVEQtmxVyYeytr9IiMbz8iCUBZkZGKxJaEsLCPW5JORRSgjE0IZhDKyIJRJKJSBZWQSysjEQhlYKBPLyCAUyiCUmYVZCAtlxH4eiu4fChEKNe0fIhTb38yi+4ZCmIWi1zLDzAiFQru37z6P2afHmWEhzfIWtHju+A+au28zs9eACUCRmWXG7vpLgPXJyCDpIRIOU7l5A9s2raVm81rqKjcSrt5IqKacrF2byanfSkHjNjpGttPRqykyp2gf59rhudRYPjtD+ewK5bMrsxPVub0JZ+UTySqA7ALIzsdyCwnlFJCZm09GTgFZeQVk5eaTnVdAdl4BObn55OR1IDevgOyMDLKT+ici8q8SVvjNrBhoiBX9POBU4H+BV4FziPbsmUp0oheRuNTX7WJT2XIq1y1jZ8UqIpVryKpZT4edGylqLKdbZAtdLUzXPY6rIp9toSJqMjqzJe9wNuZ2JpLXFcvvRmZ+F7IKu5HXqZi8jl0p6NSNwqKu5GdlH9DdvEh70WLhj7XR3+ru17bi3L2A6bFzhIAn3P15M1sEzDCznwAfAg+04tySwiLhMBvXLqdi9UfsXL8Y2/pP8mtW063uE7r7Zg4z57DYvg2eQUWoG5VZ3VlfOJI1Bb0IdepDVlEfOnTtTVH3vnTu3oeOuR3oGOinEjl0tFj43T1sZmPNzPwAx2929wXA6L1sXwmMO7CYkqoqKzbwyZL3qVnzIRkViymqWUFJ41p6Wz29Y/tsJ5+NmSWs6ziKVZ36k9n1cPJ7DqRrySC69exH78zM3fuKyP7F09TzIfCsmT0J7O7E7O7PJCyVpKTtlZtZM/91dqz+gNzyBfSuXUIPttA0DFsFndmYczjzu36ZUPfBFB52ND0HjKBzcS86BZpcJLXEU/i7AFuI9r9v4oAKv+yTRyJsWLOMT+a/jK95m+7bF9A/UkbTEG1l1puywlGs6jGc/L5jKBkyjuLiXhQHmlokPcTTnfPCZASR9q9i/WpWf/AirHyNw7bPoTeb6U20qWZ13tG80/0sCgYeR9/hx3NY52672+lFJLniGaTtSKJP1vZw92FmNgKY5O4/SXg6OaQ1NtSz9IOXqfror/Tc9AaHR1ZTDFRSyKqCMaw57Di6Dz+FfkeNYWRGRtBxRSQmnqae3wHXAvdC9EtbM/sjoMKfhnbuqGbJW8/SuPA5Bm1/i6OpocEzWJYzjHdLrqB41JkcfvR4xqjQixyy4in8Hdz9/T0GbAp+EBFJmrpdtSx64xnCHz3D0Kp/MNrq2E4+yzodT8aQsxh07CSO7tQl6JgiEqd4Cv9mMxtI05AVZucAGxKaSgLnkQhL577K9nceZvCWlxjNDiop5KOuE+kw6hwGTziDY7Jzgo4pIq0QT+G/lOiYOYPN7BNgFdEROiUFbd+yicUzf0fPFTMYHCljp2ezsNNJZI0+l6HHT2K8ir1IuxdPr56VwKlmlg+E3L068bEk2VbMf4vKV37JiG2zmGANLMs8kvePvpkhp06lVM04Iiklnl49XYGbgBMAN7N/AP/j7lsSHU4SKxIOM3/WY+TM/i1D6z+i1nOY1+0sun3uvzhyxHFBxxORBImnqWcG0TH4vxJbPx94nOiga9IONdTXMe/F+ylecA+jI2VsoJh3j7iSIWddxvjO3YKOJyIJFteTu+7+/5qt/8TMpiQqkCROY0M9c5//LSULfs0xvolVof7MHnsbo864kF5ZGixYJF3EU/hfNbNzgSdi6+cALyQukrS1SDjM3Bfvp+fcXzDON7A8cxDzjv0fRn7haxyuSTFE0s4+C7+ZVdM0iWV0gvRHYj8KEZ1L96aEp5OD9vFbfyHnlZspDa/gnxmHM2/CPYw85VzNgiSSxvZZ+N29MJlBpG19snIhFU9dzajad9hIMbPH3MqYsy4hpCdqRdJeXDNwxcbn6d98fw3LfGjauaOa+Y/dxOiyhykig3cGXsHor06jZ57mkhKRqHi6cz4IjAAW0jQRvYZlPiR99PozdH3teiZ4ObM7nUq/c2/n2N79g44lIoeYeO74J7j70IQnkVbbtnkjy/9wBcdsn8naUB8WnvYYpcd9MehYInKIiqfwv2NmQ919UcLTyAGbN2sGJW9ezyiv5p3DLmL0+T+hr5p1RKQF8RT+6USL/0agjmgvH3f3ES0fJolUU1XJoocuY1zl86wM9Wf7lx/n2OETgo4lIu1APIX/QeAC4CM+beOXAC2f9yZ5z/4npZFNvNP7G4yZ+r/k5HYIOpaItBPxFP617v5cwpPIfnkkwnuP38qYJbdTaZ1YcuYMjp1wRtCxRKSdiafwL4nNuPUXok09gLpzJltNVSXL7v0GE3a8wbwOE+h/8XSGdusZdCwRaYfiKfx5RAv+6c22qTtnEpUtn0/ksfMZEf6Edwddyfiv36gnb0Wk1eIZj//C1pzYzA4DHgZ6Ev1u4D53v8vMuhAd3bM/sBr4mrtXtuYa6WDerBkMfOP7NFoWS05/mAnHfynoSCLSzsXzANdDxKZdbM7dL9rPoY3A1e4+18wKgTlm9hLwTWCWu99qZtOAacD1B5w8xXkkwnuP/YRxy+5gZeYA8r/xGMP6HRV0LBFJAfE09Tzf7H0u8GVg/f4OcvcNxObmdfdqM1sM9AEmA5+P7TYdeA0V/s9obKhnzr3/xYTNzzC34ESGfPcx8vI1dJKItI14mnqebr5uZo8BLx/IRcysPzAaeA/oEfulgLtvMLPu+zjmEuASgL59+x7I5dq12prtLP/NOYzf+T7v9jyfcd/6lQZWE5E21ZpvCAcBcVdiMysAnga+7+5V8R7n7ve5e6m7lxYXF7ciZvuzfcsm1t45kWG1H/De0B8x4dt3q+iLSJuLp42/aVz+JhuJs2nGzLKIFv1Hm3X/3GRmvWJ3+72A8gPMnJLKP1lF7QOTGBBez/zjfsX4iRcEHUlEUlQ8TT2talw2MwMeABa7+x3NfvQcMBW4NbZ8tjXnTyXrVy2Bh79EcaSKZac9xJgTJgUdSURSWEszcLXYnOPua/dz7uOJDfVgZvNi235AtOA/YWYXA2uBr8YfN/V8snIhGQ9PIo+drJ/yJMNGnxR0JBFJcS3d8b/Ap1MvNnGgGOgOtNj47O7/2OPY5k45gIwpq2zFR2Q/Mpls6tn8b08xaMRxQUcSkTTQ0tSLw5uvx3rmXA+cCvwsoanSwLoVH5PzyCSyaKTynKcZOGx80JFEJE3st1ePmQ0ys98DfwXmAEPd/VeJDpbKNpatIOORKWTRwLavPcMAFX0RSaKW2viHAT8EjgZuAy5293CygqWqzRvLaHjwSxR5DZv+7SmOGHpM0JFEJM201MY/Hygj2tY/DhgX7agT5e5XJDZa6tleuZmq+86mZ2QLa7/4CINHnhB0JBFJQy0V/v2NxSMHYNfOHay7ZwqDwmUsPeUhho8/ff8HiYgkQEtf7k5PZpBUFgmHWfSb8xhT/xGzx/2c0pMmBx1JRNKYBnVPMI9EeP/ebzOm5nXePeJKSs/6VtCRRCTNqfAn2HuP38qE8id4t/u/M0N2jNEAAAsRSURBVP7rNwYdR0REhT+RFrz6FMcsuY0POxzHuP+6R7NmicghIZ5+/Eea2Swz+zi2PsLMfpT4aO3bmsVzOPy1y1id2Z8jv/OYRtkUkUNGPLegvwP+G2gAcPcFwLmJDNXebdu8kcwnzqPOcsif+iT5hUVBRxIR2S2ewt/B3d/fY1tjIsKkgnBjI2W/O4/iyBY2n/0gPfsOCjqSiMhnxFP4N5vZQGJj8pvZOcSmVJR/9f6DVzG8bi7zRtzA4FKNRScih5545ty9FLgPGGxmnwCrgP9IaKp26sOZ0zl2/XTe6zKJ8V/5ftBxRET2Kp6JWFYCp5pZPhBy9+rEx2p/1iydx5FvX8eyrCMZdcm9QccREdmneKZezAG+AvQHMpvG63H3/0losnZkV20NkcenUm/ZdPrmDHJyOwQdSURkn+Jp438WmEz0C90dzV4SM//+73B4ZDVln/sFPUoGBh1HRKRF8bTxl7j7GQlP0k7NfuF3jN/6HO/0+gbHfuGcoOOIiOxXPHf8b5vZ8P3vln4+WbmQwe/fwJKsoZRe+POg44iIxKWliVg+BiKxfS40s5VAHdF5dN3dRyQn4qGpob6Omj9eSKGFKLrgYbKyc4KOJCISl5aaevoAo5IVpL2Z/YcfcGzjUuaMu4OxekhLRNqRlgr/Kndfk7Qk7ciSD15m3NoH+KBoIsd88eKg44iIHJCWCn93M7tqXz909zsSkOeQV1NVSeGL32VTqJjBF/026DgiIgespcKfARQQbdOXmEW/v5yxkXKWffFxenfqEnQcEZED1lLh33AwD2mZ2YPA2UC5uw+LbesCPE70YbDVwNfcvbK110i2Ba89zbitf+HdXuczYfzEoOOIiLRKS905D/ZO//fAnv3/pwGz3H0QMCu23i5UbdtCj9euZU2ohFFT/y/oOCIirdZS4T+ooSXd/Q1g6x6bJwNNk7hPB6YczDWSacn0y+nmW6k7+zfk5uUHHUdEpNX2Wfjdfc+i3RZ6uPuG2Pk3AN33taOZXWJms81sdkVFRQKixG/Bq08xrvIF3u9zAUeO+XygWUREDtYhOwmsu9/n7qXuXlpcXBxYjh3V2yh+fRprQiWMvuDWwHKIiLSVZBf+TWbWCyC2LE/y9Q/YR49cTy8qqJ14h5p4RCQlJLvwPwdMjb2fSnTkz0PWsrmvc8zGx3mv6xSGqBePiKSIhBV+M3sMeAc4yszWmdnFwK3AaWa2HDgttn5IaqivI+OF77HFOjPkgrR8Vk1EUlQ8wzK3iruft48ftYuJaOfM+AkTwquYe+yvGVPUNeg4IiJt5pD9cjdIG9cuZ8Q/7+XDDscxZuIFQccREWlTKvx7sf7xKzGcnv9+Z9BRRETanAr/Hua/8gRjdrzJvAHfole/o4KOIyLS5lT4m9lVW0O3N3/EmlAJY8+9Meg4IiIJocLfzIczfkwf30T1ybeQnZMbdBwRkYRQ4Y/ZsGYpo9c8xJyCzzPshElBxxERSRgV/pgNT16LY/T599uDjiIiklAq/MDCt15gTM3rzOt3IT0POyLoOCIiCZX2hb+xoZ68WT9gA8WM1he6IpIG0r7wz/nTnQyIrGbD+B+S26Eg6DgiIgmX1oW/atsWjlz0KxZmD2f0xKn7P0BEJAWkdeFf+PiNdPJqcs66FQul9R+FiKSRtK1261cvZez6GcwpOp0jRp4QdBwRkaRJ28K/4elphAnR96u3BB1FRCSp0rLwL539CmOrX2HeYf9Bj5KBQccREUmqtCv8HokQmfkjNlPE8K+p+6aIpJ+0K/zzX3mcIQ0L+efQyyjo2DnoOCIiSZdWhT/c2EjR2z+lzHozZsoVQccREQlEWhX+OX+5m/6RMirGXUdWdk7QcUREApE2hX9XbQ395t/Jsswj9bCWiKS1tCn88575P3qwhYYv3KyHtUQkraVFBazatoXBK+5nQe4xHH38WUHHEREJVFoU/oVP/4wiauhw5s1BRxERCVzKF/7Kig2MWPsIc/NP0tAMIiIEVPjN7AwzW2pmK8xsWiKvtfSpH5NLHV3PvjmRlxERaTeSXvjNLAP4DXAmMBQ4z8yGJuJa5Z+sYtTGp5hbdDr9hoxNxCVERNqdIO74xwEr3H2lu9cDM4DJibjQqmduJoMIvSf/OBGnFxFpl4Io/H2Asmbr62LbPsPMLjGz2WY2u6KiolUX8qJ+zO59Pn0GDGldUhGRFJQZwDVtL9v8Xza43wfcB1BaWvovP4/HhAv+pzWHiYiktCDu+NcBhzVbLwHWB5BDRCQtBVH4PwAGmdnhZpYNnAs8F0AOEZG0lPSmHndvNLPLgJlABvCguy9Mdg4RkXQVRBs/7v4i8GIQ1xYRSXcp/+SuiIh8lgq/iEiaUeEXEUkzKvwiImnG3Fv1bFRSmVkFsKaVh3cDNrdhnPZAnzk96DOnvoP9vP3cvXjPje2i8B8MM5vt7qVB50gmfeb0oM+c+hL1edXUIyKSZlT4RUTSTDoU/vuCDhAAfeb0oM+c+hLyeVO+jV9ERD4rHe74RUSkGRV+EZE0k9KFP5mTugfNzA4zs1fNbLGZLTSz7wWdKVnMLMPMPjSz54POkgxmVmRmT5nZkth/72ODzpRoZnZl7O/1x2b2mJnlBp2prZnZg2ZWbmYfN9vWxcxeMrPlsWXntrhWyhb+ZE7qfohoBK529yHABODSFP+8zX0PWBx0iCS6C/ibuw8GRpLin93M+gBXAKXuPozocO7nBpsqIX4PnLHHtmnALHcfBMyKrR+0lC38JHFS90OBu29w97mx99VEi8G/zGWcasysBDgLuD/oLMlgZh2Bk4AHANy93t23BZsqKTKBPDPLBDqQgrP2ufsbwNY9Nk8GpsfeTwemtMW1UrnwxzWpeyoys/7AaOC9YJMkxZ3AdUAk6CBJMgCoAB6KNW/db2b5QYdKJHf/BPg5sBbYAGx3978Hmypperj7Boje3AHd2+KkqVz445rUPdWYWQHwNPB9d68KOk8imdnZQLm7zwk6SxJlAmOAe9x9NLCDNvrn/6Eq1q49GTgc6A3km9l/BJuqfUvlwp92k7qbWRbRov+ouz8TdJ4kOB6YZGariTblnWxmjwQbKeHWAevcvelfc08R/UWQyk4FVrl7hbs3AM8AxwWcKVk2mVkvgNiyvC1OmsqFP60mdTczI9ruu9jd7wg6TzK4+3+7e4m79yf63/cVd0/pO0F33wiUmdlRsU2nAIsCjJQMa4EJZtYh9vf8FFL8C+1mngOmxt5PBZ5ti5MGMuduMqThpO7HAxcAH5nZvNi2H8TmN5bUcjnwaOyGZiVwYcB5Esrd3zOzp4C5RHuvfUgKDt1gZo8Bnwe6mdk64CbgVuAJM7uY6C/Ar7bJtTRkg4hIeknlph4REdkLFX4RkTSjwi8ikmZU+EVE0owKv4hImknZ7pwirWFmXYkOhgXQEwgTHSIBoNbd0+XBIUlh6s4psg9mdjNQ4+4/DzqLSFtSU49InMysJrb8vJm9bmZPmNkyM7vVzM43s/fN7CMzGxjbr9jMnjazD2Kv44P9BCJRKvwirTOS6DwAw4k+MX2ku48jOjz05bF97gJ+4e7HAF8hTYaOlkOf2vhFWueDpuFyzeyfQNMwwR8BX4i9PxUYGh1eBoCOZlYYmy9BJDAq/CKtU9fsfaTZeoRP/38VAo51953JDCayP2rqEUmcvwOXNa2Y2agAs4jspsIvkjhXAKVmtsDMFgHfDjqQCKg7p4hI2tEdv4hImlHhFxFJMyr8IiJpRoVfRCTNqPCLiKQZFX4RkTSjwi8ikmb+P0QFscEdUFPxAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "%matplotlib inline\n", "from ecell4.prelude import *\n", "\n", "with reaction_rules():\n", " A + B == C | (0.01, 0.3)\n", "\n", "run_simulation(10.0, {'C': 60}, volume=1.0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here we give you a breakdown for `run_simulation`.\n", "`run_simulation` use ODE simulator by default, so we create `ode.World` step by step." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 8.1. Creating ODEWorld\n", "\n", "You can create `World` like this." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "from ecell4_base.core import *\n", "from ecell4_base import *" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "w = ode.World(Real3(1, 1, 1))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`Real3` is a coordinate vector.\n", "In this example, the first argument for `ode.World` constructor is a cube.\n", "Note that you can NOT use volume for `ode.World` argument, like `run_simulation` argument.\n", "\n", "Now you created a cube box for simulation, next let's throw molecules into the cube." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.0 60\n" ] } ], "source": [ "w = ode.World(Real3(1, 1, 1))\n", "w.add_molecules(Species('C'), 60)\n", "print(w.t(), w.num_molecules(Species('C'))) # must return (0.0, 60)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "Use `add_molecules` to add molecules, `remove_molecules` to remove molecules, `num_molecules` to know the number of molecules.\n", "First argument for each method is the `Species` you want to know.\n", "You can get current time by `t` method.\n", "However the number of molecules in ODE solver is real number, in these `_molecules` functions work only for integer number.\n", "When you handle real numbers in ODE, use `set_value` and `get_value`." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 8.2. How to Use Real3\n", "\n", "Before the detail of `Simulator`, we explaing more about `Real3`." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "(1.0, 2.0, 3.0)\n" ] } ], "source": [ "pos = Real3(1, 2, 3)\n", "print(pos) # must print like \n", "print(tuple(pos)) # must print (1.0, 2.0, 3.0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can not print the contents in `Real3` object directly.\n", "You need to convert `Real3` to Python tuple or list once." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.7320508075688772\n", "9.0\n" ] } ], "source": [ "pos1 = Real3(1, 1, 1)\n", "x, y, z = pos[0], pos[1], pos[2]\n", "pos2 = pos1 + pos1\n", "pos3 = pos1 * 3\n", "pos4 = pos1 / 5\n", "print(length(pos1)) # must print 1.73205080757\n", "print(dot_product(pos1, pos3)) # must print 9.0" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can use basic function like `dot_product`.\n", "Of course, you can convert `Real3` to numpy array too." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[1. 2. 3.]\n" ] } ], "source": [ "import numpy\n", "a = numpy.asarray(tuple(Real3(1, 2, 3)))\n", "print(a) # must print [ 1. 2. 3.]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`Integer3` represents a triplet of integers." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(1, 2, 3)\n" ] } ], "source": [ "g = Integer3(1, 2, 3)\n", "print(tuple(g))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Of course, you can also apply simple arithmetics to `Integer3`." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(5, 7, 9)\n", "(3, 3, 3)\n", "(2, 4, 6)\n", "32\n", "3.7416573867739413\n" ] } ], "source": [ "print(tuple(Integer3(1, 2, 3) + Integer3(4, 5, 6))) # => (5, 7, 9)\n", "print(tuple(Integer3(4, 5, 6) - Integer3(1, 2, 3))) # => (3, 3, 3)\n", "print(tuple(Integer3(1, 2, 3) * 2)) # => (2, 4, 6)\n", "print(dot_product(Integer3(1, 2, 3), Integer3(4, 5, 6))) # => 32\n", "print(length(Integer3(1, 2, 3))) # => 3.74165738677" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 8.3. Creating and Running ODE Simulator\n", "\n", "You can create a `Simulator` with `Model` and `World` like " ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "with reaction_rules():\n", " A + B > C | 0.01 # equivalent to create_binding_reaction_rule\n", " C > A + B | 0.3 # equivalent to create_unbinding_reaction_rule\n", "\n", "m = get_model()\n", "\n", "sim = ode.Simulator(w, m)\n", "sim.run(10.0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "then call `run` method, the simulation will run.\n", "In this example the simulation runs for 10 seconds.\n", "\n", "You can check the state of the `World` like this." ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "10.0 30\n" ] } ], "source": [ "print(w.t(), w.num_molecules(Species('C'))) # must return (10.0, 30)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can see that the number of the `Species` `C` decreases from 60 to 30.\n", "\n", "`World` describes the state at a timepoint, so you can NOT tack the transition during the simulation with the `World`.\n", "To obtain the time-series result, use `Observer`." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[0.0, 0.0, 60.0], [0.1, 1.7722206143224584, 58.22777938567754], [0.2, 3.4860124975248006, 56.51398750247521], [0.30000000000000004, 5.137633294323578, 54.862366705676436], [0.4, 6.7240908315276045, 53.27590916847241], [0.5, 8.2431297778128, 51.75687022218721], [0.6000000000000001, 9.693203786964157, 50.30679621303585], [0.7000000000000001, 11.073435610343822, 48.926564389656185], [0.8, 12.383567710238625, 47.61643228976138], [0.9, 13.623905934263576, 46.37609406573643], [1.0, 14.795258697983998, 45.204741302016004], [1.1, 15.898873915953542, 44.101126084046456], [1.2000000000000002, 16.936375649258018, 43.06362435074198], [1.3, 17.909702127696086, 42.090297872303914], [1.4000000000000001, 18.821046480898254, 41.17895351910174], [1.5, 19.672801194656238, 40.32719880534375], [1.6, 20.467507011344757, 39.532492988655235], [1.7000000000000002, 21.207806726765288, 38.7921932732347], [1.8, 21.896404106006866, 38.10359589399312], [1.9000000000000001, 22.53602795042438, 37.4639720495756], [2.0, 23.129401196247105, 36.87059880375288], [2.1, 23.679214810305954, 36.320785189694035], [2.2, 24.188106166226664, 35.81189383377333], [2.3000000000000003, 24.6586415307847, 35.341358469215294], [2.4000000000000004, 25.093302260285835, 34.906697739714154], [2.5, 25.494474296223075, 34.505525703776904], [2.6, 25.864440553798364, 34.13555944620162], [2.7, 26.205375812342282, 33.7946241876577], [2.8000000000000003, 26.519343739941075, 33.480656260058915], [2.9000000000000004, 26.808295712922448, 33.19170428707754], [3.0, 27.074071122047478, 32.925928877952515], [3.1, 27.318398889554317, 32.68160111044567], [3.2, 27.54289995329749, 32.45710004670249], [3.3000000000000003, 27.749090505157007, 32.250909494842965], [3.4000000000000004, 27.93838580001201, 32.06161419998797], [3.5, 28.11210437846467, 31.887895621535307], [3.6, 28.27147257094361, 31.728527429056367], [3.7, 28.41762917272733, 31.582370827272648], [3.8000000000000003, 28.551630198838193, 31.448369801161785], [3.9000000000000004, 28.674453644762753, 31.325546355237226], [4.0, 28.78700419370472, 31.212995806295257], [4.1000000000000005, 28.890117823751574, 31.109882176248405], [4.2, 28.98456627912797, 31.015433720872007], [4.3, 29.071061378816175, 30.928938621183804], [4.4, 29.150259143439115, 30.849740856560864], [4.5, 29.222763727609138, 30.77723627239084], [4.6000000000000005, 29.289131150117587, 30.71086884988239], [4.7, 29.349872818534465, 30.650127181465514], [4.800000000000001, 29.40545884814348, 30.5945411518565], [4.9, 29.456321177785775, 30.543678822214204], [5.0, 29.502856487234418, 30.49714351276556], [5.1000000000000005, 29.545428922270702, 30.454571077729277], [5.2, 29.584372634768023, 30.415627365231956], [5.300000000000001, 29.619994145881115, 30.380005854118863], [5.4, 29.652574540951854, 30.347425459048125], [5.5, 29.68237150503112, 30.317628494968858], [5.6000000000000005, 29.709621208023663, 30.290378791976316], [5.7, 29.73454004842796, 30.26545995157202], [5.800000000000001, 29.757326264497603, 30.24267373550238], [5.9, 29.77816142142091, 30.22183857857907], [6.0, 29.797211782822988, 30.202788217176995], [6.1000000000000005, 29.814629574557085, 30.185370425442898], [6.2, 29.830554148384813, 30.16944585161517], [6.300000000000001, 29.84511305275851, 30.154886947241472], [6.4, 29.858423017523837, 30.141576982476145], [6.5, 29.870590858963563, 30.12940914103642], [6.6000000000000005, 29.88171431121031, 30.118285688789673], [6.7, 29.891882789671147, 30.108117210328835], [6.800000000000001, 29.901178091733776, 30.098821908266206], [6.9, 29.909675039664805, 30.090324960335177], [7.0, 29.917442070267278, 30.082557929732705], [7.1000000000000005, 29.92454177553783, 30.075458224462153], [7.2, 29.931031398254603, 30.06896860174538], [7.300000000000001, 29.9369632861354, 30.063036713864584], [7.4, 29.942385307931307, 30.057614692068675], [7.5, 29.94734123456418, 30.0526587654358], [7.6000000000000005, 29.951871088176286, 30.048128911823696], [7.7, 29.95601146173637, 30.043988538263612], [7.800000000000001, 29.9597958116382, 30.04020418836178], [7.9, 29.963254725533886, 30.036745274466096], [8.0, 29.966416167464857, 30.033583832535125], [8.1, 29.969305702187025, 30.030694297812957], [8.200000000000001, 29.971946700432824, 30.02805329956716], [8.3, 29.974360526710818, 30.025639473289164], [8.4, 29.976566711112177, 30.023433288887805], [8.5, 29.978583106472595, 30.021416893527388], [8.6, 29.9804260321266, 30.01957396787338], [8.700000000000001, 29.98211040538864, 30.017889594611344], [8.8, 29.98364986180096, 30.016350138199023], [8.9, 29.985056865101367, 30.014943134898616], [9.0, 29.98634280778428, 30.013657192215703], [9.1, 29.987518103055105, 30.012481896944877], [9.200000000000001, 29.988592268910782, 30.0114077310892], [9.3, 29.98957400501743, 30.01042599498255], [9.4, 29.990471262999563, 30.00952873700042], [9.5, 29.9912913107032, 30.008708689296782], [9.600000000000001, 29.9920407909477, 30.00795920905228], [9.700000000000001, 29.992725775237364, 30.007274224762618], [9.8, 29.993351812863917, 30.006648187136065], [9.9, 29.993923975794377, 30.006076024205605], [10.0, 29.994446899704982, 30.005553100295]]\n" ] } ], "source": [ "w = ode.World(Real3(1, 1, 1))\n", "w.add_molecules(Species('C'), 60)\n", "sim = ode.Simulator(w, m)\n", "\n", "obs = FixedIntervalNumberObserver(0.1, ('A', 'C'))\n", "sim.run(10.0, obs)\n", "print(obs.data()) # must return [[0.0, 0.0, 60.0], ..., [10.0, 29.994446899691276, 30.005553100308752]]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "There are several types of `Observer`s for E-Cell4.\n", "`FixedIntervalNumberObserver` is the simplest `Observer` to obtain the time-series result.\n", "As its name suggests, this `Observer` records the number of molecules for each time-step.\n", "The 1st argument is the time-step, the 2nd argument is the molecule types.\n", "You can check the result with `data` method, but there is a shortcut for this." ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dd5xU9b3/8ddntrILCyy9ShEFVIoiYklMVBSDURONPbElmHuNGpOfVxOTa4wxMbmJxhT1YkWjWPFaoyJKDHYEBKUISq9LWcoWtn1+f5yzsCAss8vOnN2Z9/PxOI9T5szMex7oZ85+53u+X3N3REQkfcSiDiAiIsmlwi8ikmZU+EVE0owKv4hImlHhFxFJM5lRB4hHx44dvU+fPlHHEBFpUT766KP17t5p9+MtovD36dOH6dOnRx1DRKRFMbOlezquph4RkTSjwi8ikmZU+EVE0owKv4hImlHhFxFJMwkt/GbWzsyeNrP5ZjbPzI42s0Izm2xmC8N1+0RmEBGRXSX6iv9O4BV3HwgMBeYBNwBT3H0AMCXcFxGRJElY4TezAuCrwP0A7l7h7sXAGcCE8LQJwJmJysDc52H6Awl7eRGRliiRV/z9gCLgQTObaWb3mVk+0MXdVwOE6857erKZjTOz6WY2vaioqHEJ5jwFr/0SSjY07vkiIikokYU/EzgcuNvdhwMlNKBZx93Hu/sIdx/RqdOX7jiOz9dvhIoSmHZ7454vIpKCEln4VwAr3P39cP9pgi+CtWbWDSBcr0tYgs4DYeh58OF9sGV1wt5GRKQlSVjhd/c1wHIzOzg8dCIwF3geuDg8djHwXKIyAPC1G6CmGt76n4S+jYhIS5HoQdquAh41s2zgC+BSgi+bJ83scmAZ8J2EJmjfBw7/HsyYAMdcBYV9E/p2IiLNXUK7c7r7rLCdfoi7n+num9x9g7uf6O4DwvXGRGYA4KvXQSwTpt6W8LcSEWnu0uPO3YJuMHIczH4C1n4adRoRkUilR+EHOO5ayC2A12+OOomISKTSp/DnFQbFf+GrsGRa1GlERCKTPoUf4KgfQpvuMPkmcI86jYhIJNKr8Ge1gq//DFZOh3kvRJ1GRCQS6VX4AYZeAB0Phik3Q3Vl1GlERJIu/Qp/RiaMvhk2LILpD0adRkQk6dKv8AMcNAb6fAWm/g7KiqNOIyKSVOlZ+M3glFuhbBP8+09RpxERSar0LPwA3YbCsAvg/Xtg05Ko04iIJE36Fn6AE34RDOXw+q+iTiIikjTpXfgLusMxV8Onz8LSd6JOIyKSFOld+AGOvQYKesI//ysYvllEJMWp8Gfnwcm/hjVzYMbDUacREUk4FX6AQ74NBxwLb9wS9PQREUlhKvwQdO8cc1tQ9DVmv4ikOBX+Wt2GwBGXwAf3wppPok4jIpIwKvx1nfBLyG0LL/0UamqiTiMikhAq/HXlFcLoX8Py9+DjiVGnERFJCBX+3Q27EHodBZN/CaWJnw5YRCTZVPh3F4vB2D8FP/S+cUvUaUREmpwK/550PSyYrWv6g7Ds/ajTiIg0KRX+vfn6z6GgB7xwDVRVRJ1GRKTJqPDvTU4bOO12KJoHb/856jQiIk0moYXfzJaY2Rwzm2Vm08NjhWY22cwWhuv2icywXw46BQ75Frz1P1D0WdRpRESaRDKu+L/u7sPcfUS4fwMwxd0HAFPC/eZrzO+DSdpfuEZ9+0UkJUTR1HMGMCHcngCcGUGG+LXpAiffCsvegQ/vizqNiMh+S3Thd+A1M/vIzMaFx7q4+2qAcN15T080s3FmNt3MphcVFSU45j4MvwgOPAlevwk2fhFtFhGR/ZTown+sux8OnApcaWZfjfeJ7j7e3Ue4+4hOnTolLmE8zOCbf4FYFjz3IzX5iEiLltDC7+6rwvU64FlgJLDWzLoBhOt1iczQZNr2gDG/haVvwwfjo04jItJoCSv8ZpZvZm1qt4GTgU+A54GLw9MuBp5LVIYmN+xCGHByMEfv+oVRpxERaZREXvF3AaaZ2cfAB8BL7v4KcBsw2swWAqPD/ZahtsknKxee+T5UV0adSESkwfZZ+M3sO3Wu3H9hZpPM7PB9Pc/dv3D3oeFyiLvfGh7f4O4nuvuAcN2yRkIr6BYU/9WzYOrvok4jItJg8Vzx/9Ldt5rZccApBF0w705srGZu8OlBT59/3w5L34k6jYhIg8RT+KvD9Vjgbnd/DshOXKQWYsxt0L4PTLoCyjdHnUZEJG7xFP6VZva/wDnAy2aWE+fzUltOG/j2vbBlJTx/FbhHnUhEJC7xFPBzgFeBMe5eDBQC1yU0VUvR60g48b9h7nO6q1dEWox9Fn53LyXoa39ceKgKUF/GWsdcDQeOhld/DqtnR51GRGSf4unVcxNwPfCz8FAW8I9EhmpRYjH41j2Q1wGeugTKt0SdSESkXvE09XwLOB0ogR1347ZJZKgWJ78jnP0AbFoC//cfGtJBRJq1eAp/hbs7wYBrtXfhyu4OOAZOvgXmvwhv3xF1GhGRvYqn8D8Z9uppZ2Y/AF4H7k1srBZq1H/CoWfBlFtg0ZSo04iI7FE8P+7+EXgaeAY4GPhvd/9rooO1SGZw+l+h8yB45nLYuDjqRCIiXxJXf3x3n+zu17n7/3P3yYkO1aJl58O5/wj69U88Tzd3iUizs9fCb2ZbzWxLnfWWuvvJDNnidOgP5z4CGxbB05dBdVXUiUREdthr4Xf3Nu5eUGddUHc/mSFbpL5fhW/8ERa9DpN/GXUaEZEdMvd1gpn13tNxd1/W9HFSzIhLYf1n8N5dUNgPRv4g6kQiIvsu/MBLdbZzgb7AAuCQhCRKNSf/Jujf//J10KYrDPpm1IlEJM3F06vnsDrLAILpE6clPlqKiGXAWfdDzxHw9OWw9N2oE4lImmvwKJvuPgM4MgFZUld2Hpz/BLTrFfT0WTcv6kQiksbiaeP/SZ3dGHA4UJSwRKkqvwNc9Azcfwo8fCZc9s+g3V9EJMniueJvU2fJIWjzPyORoVJW+z7wveegugImnAGbV0SdSETSkHkLmEBkxIgRPn369KhjNJ1Vs2DCN6F1Z7jkZWjTJepEIpKCzOwjdx+x+/F4hmWebGbt6uy3N7NXmzpgWuk+DC58CrasggmnwdY1UScSkTQST1NPp3DmLQDcfRPQOXGR0kTvUXDh07B5JTw0NvgSEBFJgrgmW697E5eZHUA4RLPspz7HwncnBVf8D41Vm7+IJEU8hf9GYJqZPWJmjwBvsXM2LtlfvUfBd5+FkvVBj5/1mtVSRBIrnhu4XiHowvkE8CRwhLvH3cZvZhlmNtPMXgz3+5rZ+2a20MyeMLPsxoZPGb1GwiUvQvV2eOAUWDUz6kQiksLi+XHXgDHA4e7+ApBnZiMb8B7XAHXvWPo9cEd4F/Am4PIGvFbq6jYULnsVsvLhodPg8zejTiQiKSqepp67gKOB88P9rcDf43lxM+sJjAXuC/cNOIFgYheACcCZDcib2jr0h8tfhXa94dGzYcYjUScSkRQUT+E/yt2vBMphR6+eeJtn/gz8F1A7+3gHoNjdaweoXwH02NMTzWycmU03s+lFRWl0o3BBd7jslWBY5+d/BK/frMnbRaRJxVP4K80sg52TrXdiZyHfKzM7DVjn7h/VPbyHU/fYQ8jdx7v7CHcf0alTpzhippDctnDBk3DEJTDtdnjqe7B9a9SpRCRFxFP4/wI8C3Q2s1sJRub8bRzPOxY43cyWAI8TNPH8mWDS9toxgnoC6sC+JxlZcNqf4ZTfwvyX4b6TYMPnUacSkRQQT6+eRwmaa34HrAbOdPen4njez9y9p7v3Ac4D3nD3C4E3gbPD0y4Gnmtk9tRnBkdfGXT33LYO7v06LHgl6lQi0sLVN+duYe0CrAMmAo8Ba8NjjXU98BMzW0TQ5n//frxWeuh3PIybCu0OgInnwmu/gOrKqFOJSAu110HazGwxQfv7Htvl3T1pYwqn3CBtjVVZDq/dCB/eBz1GwNkPQPsDok4lIs3U3gZp2+t4/O7eN7GRpMGycmHsn6DPcfD81XD3sXDq72HYBUGzkIhIHOKagcvMTjezP4bLaYkOJftwyLfgh9Og2xB47j/hiYuCIR9EROIQz527txHcfTs3XK4xs98lOpjsQ/sD4OIXYPQtsPA1+PtImP0ktID5FUQkWvFc8X8DGO3uD7j7AwTDN4xNbCyJSywDjr0arngrmMZx0g+CO343LY06mYg0Y/FOtt6uznbbRASR/dB5UDDOz6l/gKXvwt+Pgn/9IfgxWERkN/EU/t8BM83sITObAHxEfDdwSTLFMuCoK+BHH8BBp8Cbt8JdR8H8l9T8IyK7iOcGronAKGBSuBzt7o8nOpg0UtuecM6EYFL3jBx4/IJgkpeVH+37uSKSFuq7gevw2gXoRjCg2nKge3hMmrN+X4P/eDvo/rn+M7j3BHjqEihaEHEwEYnaXvvxA9OBT4HaoTHrdhR3grF3pDnLyIIjvw9DzoW3/wLv3QWf/h8c9h04/nroeGDUCUUkAvXduXstcBawmWCQtWfdfVsSs+2gO3ebSMkGeOdO+OBeqCyDwWfAcddC92FRJxORBNjbnbt7Lfx1ntiXYBKWM4ClwG/dfVZCUu6FCn8T27YO3rs7GPph+5agWWjUf8KBoyEWb0cvEWnu9lb44/lxdzHBCJqvASOBg5o+niRV685w0k1w7Sdw0q+Cdv/HzoG/jYD37oGyTVEnFJEEqq+ppx/BcMpnEPyo+zjworsnvXO4rvgTrLoS5j4X/BWwcjpk5gbDQhxxCfQ6SuMAibRQDW7qMbMaYDbB1f4Wdpspy91vT0DOPVLhT6JVs2DGBJj9FFRsDe4IHno+DDkH2veJOp2INEBjCv+v2Mu0iADufnOTpdsHFf4IbN8W/BXw8URY8u/gWI8j4JBvB38NtN3jVMki0ow0+sfd5kCFP2LFy+CTZ+CTSbBmdnCs++EwcCwMPA06HazmIJFmSIVfmsb6RTDvuWAe4JXhv0nb3jBgdLD0OQ5y2kSbUUQAFX5JhC2r4bNXYOFk+GIqVJZALDNoEup7PPQ5FnqOhOy8qJOKpKXGtPFf4+53mtmx7v52whPWQ4W/BajaDsveg8X/Cr4EVs0Erwm+CLoPD3oH9TwSeo2Egu5RpxVJC40p/LPcfZiZzXD3SMfmUeFvgco3w/IPYOnbwVDRq2ZC9fbgsTbdgi+D7sOh6xDoeljwZaDfCUSaVIPn3AXmmdkSoJOZza77WgSTrQ9p4oySSnLb7mz3B6iqCH4YXvFh8CWwcgYseHnn+XkdoPPgcBkU/GDc8WDI7xBNfpEUVt9k6+ebWVfgVeD05EWSlJSZDT1HBEut8i2wbi6smRN8KaydC7MehYo6Q0K1KoQOB4ZLP2jfN1z6QF6h/koQaYT6rvhx9zXAUDPLZudQDQvcvTLhyST15RZA71HBUqumBjYvD4aSXv9ZMJzExi+C3w0+fmzX52flQ7te0LZXMA9BQY/g/oI23YKloBvkFOjLQWQ39RZ+ADM7HngYWELQzNPLzC5297cSnE3SUSwWTCTf/oCdzUS1KkqCewo2LYGNi4MviOJlwbJqBpRu+PLrZbYKxiZq3SVY53WA/E6Q3zHYzisM/qrIK4RW7fVFIWlhn4UfuB042d0XAJjZQcBE4Ij6nmRmucBbQE74Pk+7+03haJ+PA4XADOC77l7R+I8gaSM7P2j/7zxoz49XlsOWlbB1DWxdHSzb1gajkW5dE3xZLP8AStcHPY72xDKC3ydatYPcdsFfJTkFO9c5bYIlOx+y20BO63A7P/gLJDtv5zozV18i0izFU/izaos+gLt/ZmZZcTxvO3CCu28Lz59mZv8EfgLc4e6Pm9k9wOXA3Y0JL7KLrFzo0D9Y6lNTA+XFULox+BIo3RiMSFpWuy4OHi/fHPwOsWV1MHz19m3B+EUNyhR+AWS1CpbMVpCZEx7LDabHzAyXjOyd6x3bWRDLCo9lButYVng8s87jmcF+LDweywjXdbYttnPfMsLj4bbF9rCtL61UFU/hn25m9wOPhPsXEky4Xi8P+onW/kqXFS61M3ddEB6fAPwKFX5JplgsaNrJKwQaOAtZTXXw43NFyc4vgorSYL9iG1SWBpPcVJQE66qy4PGq7cF2ZTlUlQf7ZZuC3k7V28NjFcG6uiJYaqoS8vEb9HEtAzDcYrhl4LXbGFgMJ4abBcexneeG+27ByO9O7THCY+H2jufUfQ122f7yuTu3d57Lju2dx6nzGuz1uNc5zh6O73ysrr08p/bLcrdu8r7Hc3Z/Sdvj4a7n/Y3OPft9KcH+iKfw/wdwJXA1wSd8C7grnhc3swyCL4kDgb8DnwPF7l77X/QKYI+jfZnZOGAcQO/eveN5O5HEi4VNQbltqalxyiqrKamoonR7sC6rqKassprSimrKK4OlrKKa8qqacL+G7VU719uraqioqgnX1VRU1VBRXUNllVNVVYVXV+DVlcHQ2TUVUF2JVVdhXkkm1WRSTRbVZFJFptXsOBajhqxwnUk1GdSQQQ0x27kfPFaDhY8FZT08Dydm4brOYzGcDGpqy/qO7TqlHgMyrGZHua59LFbn8Rg1QYk2djwvBsE5tmvJN3aeU/f8nefU7Cj3hoOBhYV35x8tvus5dV6fPR6vy7+0X3tg16+VvdwTtdfX3fdzASormr4lfJ+F3923E7TzN3gYZnevBoaZWTvgWWBPjbN7/MTuPh4YD8ENXA19b5F9qalxNpdVsqm0gk2llWwuq6C4tJLNZcGypayKreWVbCmvZNv2KraVV7G1vIpt26so2V5FSUV1g98zOyNGTlaMnMwMcjJ3bmdnxsjJiJGXnUnbDCM7M0Z2ZgZZMSMrI0ZWppEZi5GdGSMrI9jOyjAyM2JkxozMmJFRZzszw4hZ8NyYGRm154RL7bGMGGTEYsSMHceCNZgZGRbsm7HjMQvPjdnOczB2Pc7O83asCQu9mpAiF88V/35z92IzmwqMAtqZWWZ41d8TWJWMDJIeamqcjaUVrN1Sztot5RRt3b5jWb+tgvXbtrOhpIKNJRUUl1ZQU88lRX52BgWtsmiTm0nrnEza5WXTszCP1tmZ5Odk0jong7ycYDs/O4O87AxaZWcG66wMWoXr3KwMcrNi5GZmEIup6En0Elb4zawTUBkW/VbAScDvgTeBswl69lxMMNGLSFwqqmpYWVzGso2lrNhUyopNZawqLmN1cTmrNpexdks5ldVfruYFuZl0bJNDx/wcBnRuTfv8bDrkZ1OYn027vCza5WXTPi+btq2yaNsqi4LcTDIzNP+wpKZ6C3/YRn+bu1/XiNfuBkwIXyMGPOnuL5rZXOBxM/sNMBO4vxGvLSmspsZZWVzGoqJtfL5uG4vXl7BkQwlL1peyanPZLr+bZcaMbu1y6d62FUf2KaRr21y6tc2lc5tcOhfk0KUgl46ts8nJzIjuA4k0M/u6c7fazI4wM/MGjt/s7rOB4Xs4/gXBpO0ibCypYN7qLcxdtYX5a7by2dqtLFy3lfLKnf3s27bKom/HfEb2LaR3YR69C/PoVZhHr8JWdG6TS4aaT0QaJJ6mnpnAc2b2FFBSe9DdJyUslaSkzWWVzFpezJwVxcxesZk5KzezenP5jsc7t8nh4K5tuGDkAQzo0poDO7emf6fWFOZnR5haJPXEU/gLgQ0E/e9rOaDCL3vl7qzYVMb7izfy4eKNzFi2iYXrdg6+1q9jPkf2KeTQHgUM7taWwd0LVOBFkiSe7pyXJiOItHxrt5Tz9qL1TFu0nvc+38Cq8Gq+bassDu/djtOHdmd47/Yc1rMtbVvFc/O3iCRCPIO0HURwZ20Xdz/UzIYAp7v7bxKeTpq1quoapi/dxNQFRUxdsI75a4LhDNrnZXFM/478sF8hR/XtwIDOrdWNUaQZiaep517gOuB/IfjR1sweA1T401BZRTVvLSzi1U/X8Mb8dRSXVpIZM47sU8gNpw7kKwM6MqhrgQq9SDMWT+HPc/cPdrvbLvpBRCRptldV868FRbw4ezWvz1tLaUU1bVtlceLAzowe3IXjBnSkTa6abkRaingK/3oz6084tIKZnQ2sTmgqiZy7M3N5MZNmrOCFj1ezuayS9nlZnDGsB2MP68ZR/QrJ0g1OIi1SPIX/SoIxcwaa2UpgMcEInZKCiksrmDRjJRM/WMbCddvIzYpxyiFdOXN4D447sKOKvUgKiKdXzxfASWaWD8TcvYEDkktL8MnKzTz49hJemL2KiqoahvZqx23fPoyxQ7qpGUckxcTTq6cDcBNwHOBmNg34tbvvYZ47aUlqapzX563lvmmL+WDxRlplZfCdI3pywVG9OaR726jjiUiCxNPU8zjBGPxnhfsXAk8QDLomLVBldQ3Pz1rFPf/6nIXrttGjXStu/MYgzjmyl/rXi6SBuO7cdfdb6uz/xszOTFQgSZyq6homzVzJX99YyPKNZQzs2oY7zxvG2MO6aSRKkTQST+F/08zOA54M988GXkpcJGlqNTXOC7NXccfkz1iyoZQhPdvyq28ewgkDO2tSDJE0tNfCb2ZbCbpwGsEE6f8IH4oRzKV7U8LTyX575/P1/O7l+cxZuZlB3Qq493sjOGmQCr5IOttr4Xf3NskMIk1r6YYSbnlxLq/PW0ePdq2449yhnDG0h+6oFZH4ZuAKx+fpU/d8DcvcPJVVVHP31EXc89YXZMWM68cM5NJj+5CbpYlIRCQQT3fOB4AhwKdA7ewYGpa5GXrrsyJ+/uwcVmwq48xh3fnZNwbRpSA36lgi0szEc8U/yt0HJzyJNNqmkgpueWkuk2aspF+nfB4fN4pR/TpEHUtEmql4Cv+7ZjbY3ecmPI002JR5a7n+mTkUl1Zw1QkHcuXXD1SzjojUK57CP4Gg+K8BthP08nF3H5LQZFKvbdur+M2Lc3n8w+UM7NqGhy8byeDuBVHHEpEWIJ7C/wDwXWAOO9v4JUJzVmzmRxNnsGxjKT88vj/Xjh5ATqau8kUkPvEU/mXu/nzCk8g+uTsT3lnCb1+eT4fW2Twx7mhG9i2MOpaItDDxFP754YxbLxA09QDqzpls27ZXcd1TH/PPT9Zw4sDO/PE7Q2mvyclFpBHiKfytCAr+yXWOqTtnEn1RtI1xj3zE4vUl3PiNQXz/K311562INFo84/Ff2pgXNrNewMNAV4LfBsa7+51mVkgwumcfYAlwjrtvasx7pIMp89by48dnkZUZ45HLR3JM/45RRxKRFi6eG7geJJx2sS53v2wfT60CfuruM8ysDfCRmU0GLgGmuPttZnYDcANwfYOTpzh35/5pi7n15Xkc0r2Aey46gp7t86KOJSIpIJ6mnhfrbOcC3wJW7etJ7r6acG5ed99qZvOAHsAZwNfC0yYAU1Hh30VVdQ03vzCXR95byphDunLHucNola1eOyLSNOJp6nmm7r6ZTQReb8ibmFkfYDjwPtAl/FLA3VebWee9PGccMA6gd+/eDXm7Fq20ooorH53BmwuKuOKr/bh+zEANrCYiTSquQdp2MwCIuxKbWWvgGeDH7r4l3h8l3X08wSTvjBgx4ktNTamouLSCyx76kFnLi/nNmYdy0agDoo4kIikonjb+2nH5a60hzqYZM8siKPqP1un+udbMuoVX+92AdQ3MnJLWbC7new+8z5L1pdx14RGMObRr1JFEJEXF09TTqHH5Lbi0vx+Y5+6313noeeBi4LZw/VxjXj+VLN9Yyvn3vsemkgoeuvRIjjlQPXdEJHHqm4Gr3uYcd1+2j9c+lnCoBzObFR77OUHBf9LMLgeWAd+JP27qWbqhhPPHv0dJRTUTx41iSM92UUcSkRRX3xX/S+ycerGWA52AzkC93Uzcfdpuz63rxAZkTFmL15dwwb3vUV5ZzWM/OIpDureNOpKIpIH6pl48rO5+2DPneuAk4LcJTZUGlqwv4bzx71JZ7Tz2g1EM6qaRNUUkOWL7OsHMBpjZQ8A/gY+Awe7+10QHS2Wrisu48L73qaiqYaKKvogkWX1t/IcCNwKHAH8ALnf36mQFS1VFW7dz0X3vs6WskonjRnFwV81pLyLJVV8b/8fAcoK2/pHAyLp98N396sRGSz2byyr57v3vs3pzOY9cPpJDe6hNX0SSr77Cv6+xeKQByiurGffwdD4v2saDl4xkRB+Noy8i0ajvx90JyQySympqnJ8++THvL97IX84fznED1E9fRKKzzx93Zf+4O7e8NJeX5qzmxm8M4vSh3aOOJCJpToU/wSa8s4QH317CZcf25ftf6Rt1HBERFf5EmrpgHb9+cS6jB3fhF2MHadYsEWkW4unHf5CZTTGzT8L9IWb2i8RHa9kWrt3KVY/N5OCuBfz53GEaWllEmo14rvjvBX4GVAK4+2zgvESGauk2lVRw+YTp5GRlcN/FI8jPaczo1yIiiRFP4c9z9w92O1aViDCpoLrGufrxmazZXM747x1Bj3atoo4kIrKLeC5F15tZf8Ix+c3sbMIpFeXL/vTaAv69cD2/P+swDu/dPuo4IiJfEk/hv5JgJqyBZrYSWAxclNBULdQrn6zmrqmfc/7I3px7ZPpMFykiLUs8E7F8AZxkZvlAzN23Jj5Wy7No3TZ++uTHDO3Vjl+dPjjqOCIiexXP1Is5wFlAHyCztkuiu/86oclakPLKan702AxysjK456LDycmsd6oCEZFIxdPU8xywmWBI5u2JjdMy3fzCXOav2cpDlx5Jt7b6MVdEmrd4Cn9Pdx+T8CQt1PMfr2LiB8v44fH9+drBnaOOIyKyT/F053zHzA7b92npZ+mGEn4+aQ5HHNCen558UNRxRETiUt9ELJ8ANeE5l5rZFwRNPQa4uw9JTsTmqbK6hmsen0XM4C/nDycrQ6NfiEjLUF9TTw9gWLKCtDR/fWMRs5YX87cLhusmLRFpUeor/IvdfWnSkrQgHy3dxN/eWMi3D+/BaUM0zLKItCz1Ff7OZvaTvT3o7rcnIE+zt217Fdc+MYvu7Vpx8+mHRB1HRKTB6iv8GUBrgjZ9Cd360lxWbCrliSuOpk1uVtRxREQarL7Cv3p/btIysweA04B17n5oeKwQeILgZrAlwDnuvqmx75Fs//qsiIkfLOeKr/bjSM2ZK+o/MY0AAAh6SURBVCItVH1dUfb3Sv8hYPf+/zcAU9x9ADAl3G8RtpRXcsMzs+nfKZ9rR6vrpoi0XPUV/hP354Xd/S1g426HzwBqJ3GfAJy5P++RTLe+OI+1W8r50znDyM3SkAwi0nLttfC7++5Fuyl0cffV4euvBvZ6q6uZjTOz6WY2vaioKAFR4jd1wTqemL6cK47vz7Be7SLNIiKyv5rtXUfuPt7dR7j7iE6dOkWWo2R7FTc++wn9O+VzzYkDIsshItJUkl3415pZN4BwvS7J799gt0/+jJXFZdx21hA18YhISkh24X8euDjcvphg5M9m6+PlxTz49mIuPKq3evGISMpIWOE3s4nAu8DBZrbCzC4HbgNGm9lCYHS43yxVVtdww6Q5dGqTw/WnDow6johIk4lnWOZGcffz9/LQfvUWSpb7py1m3uot3HPRERToRi0RSSHN9sfdKK0sLuPO1xdy0qAujDm0a9RxRESalAr/Htzywlwc19y5IpKSVPh38+b8dbzy6RquOmEAPdvnRR1HRKTJqfDXUV5ZzU3Pf0r/Tvn84Cv9oo4jIpIQCftxtyW651+fs2xjKY99/yiyM/WdKCKpSdUttGJTKXdP/ZyxQ7pxzIEdo44jIpIwKvyh3708HzO48RuDoo4iIpJQKvzAu59v4KU5q/mP4w+ku+bPFZEUl/aFv6q6hptf+JQe7VpxxfH6QVdEUl/aF/6JHy5n/pqt3Dh2kAZhE5G0kNaFf0t5JXdM/oyj+hZyqu7QFZE0kdaF/+9vLmJTaQW/PG0wZppTXkTSQ9oW/uUbS3lw2hK+NbwHh/ZoG3UcEZGkSdvC/4dXFxCLwXWnHBx1FBGRpErLwj9z2SZe+HgVP/hKP7q1VfdNEUkvaVf43Z3fvjyPjq1zuOL4/lHHERFJurQr/FPmrePDJZv48UkDaJ2joYpEJP2kVeGvrnF+/8p8+nbM59wje0UdR0QkEmlV+J+ZsYKF67Zx3SkHk5WRVh9dRGSHtKl+5ZXV3DH5M4b2aqebtUQkraVN4X/43SWs3lzODWMG6mYtEUlraVH4t5RXctfUzzn+oE4c3b9D1HFERCKVFoX/vn8vpri0UjdriYiQBoV/Y0kF9//7C049tKuGZhARIaLCb2ZjzGyBmS0ysxsS+V53T11EWWU1Pxl9UCLfRkSkxUh64TezDODvwKnAYOB8MxuciPdas7mch99dypnDezCgS5tEvIWISIsTxRX/SGCRu3/h7hXA48AZiXijv76xkOoa58cn6mpfRKRWFIW/B7C8zv6K8NguzGycmU03s+lFRUWNeqPehXl8/yv96N0hr3FJRURSUBSD1eypE71/6YD7eGA8wIgRI770eDw0CJuIyJdFccW/Aqg7UE5PYFUEOURE0lIUhf9DYICZ9TWzbOA84PkIcoiIpKWkN/W4e5WZ/Qh4FcgAHnD3T5OdQ0QkXUUyIL27vwy8HMV7i4iku5S/c1dERHalwi8ikmZU+EVE0owKv4hImjH3Rt0blVRmVgQsbeTTOwLrmzBOS6DPnB70mVPf/n7eA9y90+4HW0Th3x9mNt3dR0SdI5n0mdODPnPqS9TnVVOPiEiaUeEXEUkz6VD4x0cdIAL6zOlBnzn1JeTzpnwbv4iI7CodrvhFRKQOFX4RkTST0oU/mZO6R83MepnZm2Y2z8w+NbNros6ULGaWYWYzzezFqLMkg5m1M7OnzWx++O99dNSZEs3Mrg3/u/7EzCaaWW7UmZqamT1gZuvM7JM6xwrNbLKZLQzX7ZvivVK28CdzUvdmogr4qbsPAkYBV6b4563rGmBe1CGS6E7gFXcfCAwlxT+7mfUArgZGuPuhBMO5nxdtqoR4CBiz27EbgCnuPgCYEu7vt5Qt/CRxUvfmwN1Xu/uMcHsrQTH40lzGqcbMegJjgfuizpIMZlYAfBW4H8DdK9y9ONpUSZEJtDKzTCCPFJy1z93fAjbudvgMYEK4PQE4syneK5ULf1yTuqciM+sDDAfejzZJUvwZ+C+gJuogSdIPKAIeDJu37jOz/KhDJZK7rwT+CCwDVgOb3f21aFMlTRd3Xw3BxR3QuSleNJULf1yTuqcaM2sNPAP82N23RJ0nkczsNGCdu38UdZYkygQOB+529+FACU30539zFbZrnwH0BboD+WZ2UbSpWrZULvxpN6m7mWURFP1H3X1S1HmS4FjgdDNbQtCUd4KZ/SPaSAm3Aljh7rV/zT1N8EWQyk4CFrt7kbtXApOAYyLOlCxrzawbQLhe1xQvmsqFP60mdTczI2j3nefut0edJxnc/Wfu3tPd+xD8+77h7il9Jejua4DlZnZweOhEYG6EkZJhGTDKzPLC/85PJMV/0K7jeeDicPti4LmmeNFI5txNhjSc1P1Y4LvAHDObFR77eTi/saSWq4BHwwuaL4BLI86TUO7+vpk9Dcwg6L02kxQcusHMJgJfAzqa2QrgJuA24Ekzu5zgC/A7TfJeGrJBRCS9pHJTj4iI7IEKv4hImlHhFxFJMyr8IiJpRoVfRCTNpGx3TpHGMLMOBINhAXQFqgmGSAAodfd0uXFIUpi6c4rshZn9Ctjm7n+MOotIU1JTj0iczGxbuP6amf3LzJ40s8/M7DYzu9DMPjCzOWbWPzyvk5k9Y2Yfhsux0X4CkYAKv0jjDCWYB+AwgjumD3L3kQTDQ18VnnMncIe7HwmcRZoMHS3Nn9r4RRrnw9rhcs3sc6B2mOA5wNfD7ZOAwcHwMgAUmFmbcL4Ekcio8Is0zvY62zV19mvY+f9VDDja3cuSGUxkX9TUI5I4rwE/qt0xs2ERZhHZQYVfJHGuBkaY2Wwzmwv8MOpAIqDunCIiaUdX/CIiaUaFX0Qkzajwi4ikGRV+EZE0o8IvIpJmVPhFRNKMCr+ISJr5/+kzSnXvxQwJAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "show(obs)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This plots the time-series result easily.\n", "\n", "We explained the internal of `run_simulation` function.\n", "When you change the `World` after creating the `Simulator`, you need to indicate it to `Simulator`.\n", "So do NOT forget to call `sim.initialize()` after that." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 8.4. Switching the Solver\n", "\n", "It is NOT difficult to switch the solver to stochastic method, as we showed `run_simulation`." ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAcw0lEQVR4nO3df7RcdXnv8feHACI/UkQIpGAMuEANhF895Ye4vChci4LELsGi1pVyWTftvVZo9baEtveiXLRplxdK7w9qFGisSETQFfxRFFNSr9WFngCXSKLVhhACCTlUkKCVkPDcP/aeZOcwM2efmdl7z+z9ea111pm9z8zsZyY5z9nz7Of7/SoiMDOz5tir6gDMzKxcTvxmZg3jxG9m1jBO/GZmDePEb2bWMHtXHUAehx56aMydO7fqMMzMRsrq1aufiojDJu8ficQ/d+5cxsfHqw7DzGykSHq03X6XeszMGsaJ38ysYZz4zcwaxonfzKxhnPjNzBqm0MQv6WBJd0j6oaR1ks6UdIikeyT9OP3+iiJjMDOzPRV9xn8DcHdEvA44CVgHLAZWRsSxwMp028zMSlJYH7+kmcCbgN8BiIjtwHZJC4Cz07stA1YBVxYSxN8vhi1rdm/PvwjGLi3kUGZmo6LIM/5jgAngFkkPSPq0pAOAwyNiM0D6fVa7B0taJGlc0vjExET/0WxZA2vu6P95zMxGXJEjd/cGTgU+GBH3SbqBaZR1ImIpsBRgbGyst9Vi3rZk9+1bzu/pKczM6qbIM/5NwKaIuC/dvoPkD8GTkmYDpN+3FhiDmZlNUljij4gtwGOSXpvuOgdYC9wFLEz3LQRWFBWDmZm9VNGTtH0QuFXSvsB64FKSPza3S7oM2AhcXHAMZmaWUWjij4gHgbE2PzqnyOOamVlnIzEt88BsWbP7Iq9bO82soZqT+OdftPt2q7ffid/MGqg5iX/s0t2J3q2dZtZgnqTNzKxhnPjNzBrGid/MrGGaU+OfzB0+ZtZQzUz87vAxswZrZuJ3h4+ZNZhr/GZmDdPMM/7JsvX+LNf+zayGnPiz9f4s1/7NrKac+LP1/izX/s2splzjNzNrGCd+M7OGceI3M2sYJ34zs4Zx4jczaxgnfjOzhnHiNzNrGPfxd+MZPM2shpz4O/EMnmZWU078nXgGTzOrKSf+vDyRm5nVRKGJX9IGYBuwE9gREWOSDgE+D8wFNgDvjoini4yjb57IzcxqpIwz/jdHxFOZ7cXAyohYImlxun1lCXH0zhO5mVmNVNHOuQBYlt5eBryzghgGp1UCuuV8GL+l6mjMzKZUdOIP4BuSVktalO47PCI2A6TfZ7V7oKRFksYljU9MTBQcZo/mXwRHzE9ub1kDa+6oNh4zsxyKLvWcFRFPSJoF3CPph3kfGBFLgaUAY2NjUVSAfXHnj5mNoEITf0Q8kX7fKulLwGnAk5JmR8RmSbOBrUXGUCoP+DKzEVBYqUfSAZIOat0G3gr8ALgLWJjebSGwoqgYSuWyj5mNiCLP+A8HviSpdZzPRcTdkr4P3C7pMmAjcHGBMZTHZR8zGxFTJn5JFwN3R8Q2SX8GnApcGxH3d3tcRKwHTmqz/1+Bc3qM18zM+pTnjP+/RsQXJL0R+A3gE8CNwOmFRjbqJo/0dc3fzIZEnhr/zvT7+cCNEbEC2Le4kGogW+8H1/zNbKjkOeN/XNIngXOBv5D0MjyPf3eTR/q65m9mQyRPAn838HXgvIh4BjgE+KNCozIzs8JMmfgj4hckvfZvTHftAH5cZFBmZlacKRO/pKtJJlG7Kt21D/DZIoMyM7Pi5Cn1/CZwIfBz2DUa96AigzIzs+LkSfzbIyJIJlxrjcI1M7MRlSfx35529Rws6T8C3wQ+VWxYZmZWlCnbOSPiE5L+PfAs8Frgv0XEPYVHZmZmhcg1V0+a6J3s++GZO81sSHRM/JK2kdT1lX7f9SMgImJmwbHVR3bNXq/Ta2YV65j4I8KdO4PimTvNbIjkmZ1zTrv9EbFx8OE0hMs+ZlahPDX+r2Zu7wccDfwIOL6QiOrOZR8zq1ierp752W1JpwK/W1hEdeeyj5lVbNorcEXE/ZJ+vYhgGsllHzMrWZ4a/4cym3uRrMA1UVhETeKyj5lVIM8Zf7a7ZwdJzf/OYsJpGJd9zKwCeWr8Hy0jEDMzK0eeUs89wMXpIixIegWwPCJ+o+jgGmv8ls5LNfo6gJn1Kc8kbYe1kj5ARDwNzCouJGPNHbtr/lleu9fMBiBPjX+npDmtAVuSXs2eUzhYEY6YD5d+dc99vg5gZgOQJ/H/KfBtSf+Ybr8JWFRcSA3Wau3csiZJ/N3uM1mrBDS5TOTSkJlNkufi7t3poK0zSCZo+8OIeCrvASTNAMaBxyPiAklHA8tJFm2/H3h/RGzvKfo6ybZ2HjF/z+1298nKtoK2ykRHzHeLqJm1lefiroDzgGMi4hpJcySdFhHfy3mMK4B1QGs2z78Aro+I5ZL+BrgMuLGH2Osl29o53ftM/gTQKhO5NGRmbeQp9fwf4EXgLcA1wDaSPv4pR+9KOgo4H/gY8KH0j8hbgPemd1kGfAQn/v51KhN1Kg1luRxkVq5unXtZBf1u5kn8p0fEqZIegKSrR9K+OZ//r4A/ZvcgsFcCz0TEjnR7E3BkuwdKWkR6LWHOnLYThFpLpzJRp9JQlstBZuXLlmQ7KfB3M0/ifyGt07cWWz+M5BNAV5IuALZGxGpJZ7d2t7lr2w6hiFgKLAUYGxtzF1E3nUpAecpHLgeZVaNd515Wgb+beRL/XwNfAmZJ+hhwEfBnOR53FnChpLeTTOc8k+QTwMGS9k7P+o8Cnugpchuc6XYKdeKSkdlImHIAV0TcSlKu+XNgM/DOiPhCjsddFRFHRcRc4BLgHyLifcC9JH88ABYCK3qM3QZh/kXtP25mB4t1GlDW6f5mNtS6rbl7SGZzK3Bb9mcR8dMej3klsFzStcADwE09Po8NwnQ7hTpxychsZHQr9axm92LrkwVwTN6DRMQqYFV6ez1wWu4IzcxsoLottn50mYHYEMozktjMOut0fazi36lcK3BJupBkqgaAVRHxleJCsqGQZySxmXXXqW0zz+9UgX8Y8ozcXUIyWOvWdNcVks6KiKsKi8qql6cV1MymNtX1sU7etmTwsaTynPG/HTg5Il4EkLSM5KKsE7/taaq2ULMmyJZ3hrRMmmc+foCDM7d/pYhAbMTlaQs1a4Js+/OQlknznPH/OfCApHtJOnzehM/2bbK8baFmTdBreackeaZlvk3SKpI6v4ArI2JL0YFZAxQxIngU1yNo9z7U/TXXwZB27OTRsdQj6dTWFzCbZEK1x4BfTfeZ9aeIEcHZ5xyVMtPk96EJr7kOOv3/HdLyTla3M/5x4GFgIt3ODuQKkumVzfpTxIjgUVyPIPs+NOU118GQl3Q66Zb4Pwy8C/g3khWzvhQRz5USlZmZFabbyN3rgevTpRLfA6yU9Cjw8Yh4sKwArQaybZ691J+nevwItM+9xFQx9/Oe9ft+255GuJbfSZ7ZOR8hmUHzGyRz7BxXdFBWI9k2z17qz3kePwLtcy/RLeZ+3rN+3297qRGu5XfSbXbOY0imU15AclF3OfCxiPhlSbFZHWTbPHupP+d9/CjWWjvF3M971u/7be2N4v+vLrrV+H8CPERytv8sMAf4z8myuRAR1xUendVPtgzRz0flMss73dpOe2277CXmXkoOwzqa2i2olepW6rmGZOWtF4EDSdbNzX6ZTc/k0b39fFQus7zT6aN+P22XvcQ83ZLDMI+mdgtqpbpd3P1IiXFYEwx64rcyP363O1Y/bZeDjKOTYR9N7RbUyuSaltlsaAx6jYA8o4fzllKK6jjyuggAfO6+jax48PFd2wtOPpL3nj6nmIONYqfYNOSdpM2setnSxaDKO3lGD+cppRTVcVTEax5RKx58nLWbnwVg7eZn9/gjMHCj2Ck2Dd26eq6IiBvSuff/qcygzNoqao2AXssvZXQceV2EPcybPZPP/+6Z/NYnv1v8wWrWyZPVrdRzKXAD8D8Bz81jo6lTV0v254P6GO+BUy+VtxOpzXs3ubSzdvOzzJs9c4/t1h+AdmWfUktDI6Zb4l8naQNwmKSHMvsFREScWGhkZv3K8/F8UB/js8/RKhE48edberDDe9cq7bSS/bzZM1lw8pEAu74Du8o/k5N69vGd7tNU3bp63iPpCODrwIXlhWQ2IGWWSTxwqrOpSiZd3rtWaWey954+Z1cS71b2KbU0NEK6dvWk8+6fJGlfdk/V8KOIeKHwyMzMrBB5Flv/d8BngA0kZZ5XSVoYEd8qODaz0TXsLZhFXo/ooxXyyW2/5KnnnueaT373JTX9fk11TaDuLZxZefr4rwPeGhE/ApB0HHAb8GvdHiRpP+BbwMvS49wREVens30uBw4B7gfeHxHbe38JZkMmW7MexlbAoq9HZOv603z9Tz33PL/YvhPYs6bfrzzXBPqJe9TkSfz7tJI+QET8s6R9cjzueeAtEfFcev9vS/p74EPA9RGxXNLfAJcBN/YSvNlQGvYWzDKuR/TRCrn/vjPa1vX7kfeaQJ1bOLPyJP5xSTcBf5duvw9YPdWDIiKA1sIt+6RfrZW73pvuXwZ8BCd+s0o8ue2XHPj0OjZ8/I0De865L6znuVe8nsMH9ozdnfOLr3HKz77Jwx+fscf+/7J9Jw/8yrnAmbnuX3bcVcozcvc/kSzBeDlwBbAW+L08Ty5phqQHga3APcC/AM9ExI70LpuAtp/lJC2SNC5pfGJiot1dzKxPK3a+gbXx6oE+59p4NSt2vmGgz9nNghnfYZ4efcn+eXqUBTO+k/v+ZcddpSnP+CPieZI6/7SnYY6IncDJkg4mmenz9e3u1uGxS4GlAGNjY23vY2b9Wbn/21m5/9sHWlpplVIWDewZuzv8oP3goFM4vs0kegdM4/5lx12lUiZpi4hnJK0CzgAOlrR3etZ/FPBEGTGYlWXYR4xm4xt050zLlB00Xcx9Yf3gO45yduz0E/coKWySNkmHpWf6SHo5cC6wDrgXaF0uX0iy0ItZbZQ6mVgPsvENsnOmZcHJR+76YzLd1/9PL38zG/Y5JtkY5Dz9OSZd6yfuUdP1jF/SDGBJRPxRD889G1iWPsdewO0R8RVJa4Hlkq4FHgBu6uG5zYbasI8Y7TQidhByd9C0sav0dOmZg+84mqJjp5+4R81UI3d3Svo1SUq7dHKLiIeAU9rsX0+yaLuZmVUgT43/AWCFpC8AP2/tjIgvFhaV2QBMrrV3Mohabhl181GVrZtntd73ru/d5NlVWzX/vKNsh30EdUXy1PgPAf6VpP/+HenXBUUGZTYI2Vp2J4Oq5RZdNx9V2bp5VvZ97/jeTV4zOFvzz7NQihex6ShPO+cQD0E0626qWvYga7lF1s1HVbZunjX5fW/73k0eAT255j+dWT9tD3kmaTuOZGTt4RFxgqQTgQsj4trCozPLKKpNslMpImuqhT66lXf6aREc9tbQ0pVYuqlza2eeUs+ngKuAF2DXRdtLigzKrJ0i2iQ7lSKyOh0rT3mn3xbBYW8NLVWJpZu6t3bmubi7f0R8T1J2345OdzYr0qDbJDuVIrLyLPSR5/l7jXnYW0NLU2Lppu6tnXkS/1OSXkM6tYKki4DNhUZltdSuy6afj9BlfhRvVw7qpXun35in6pDJ6tTVNCxdR63XMizxNEmexP8BkjlzXifpceARkhk6zaZl8hqq/ayDmmt+9QHp1KEz3e6dfmPudKw8a85mDUPXUfb4wxBP0+Tp6lkPnCvpAGCviNhWfFhWV9nSSD8focv8KJ6nHDTd5+kl5rwdMlnD2mk0qPfUejPlxV1Jr5T018D/BVZJukHSK4sPzczMipCn1LOcZAnFd6Xb7wM+TzLpmllX/Y5oLaNlsiqdYu7nPfMI4mJM59rKKMg1cjci/ntEPJJ+XQscXHRgVg/9jmgto2WyCt1i7uc98wjiwcsz+njU5Dnjv1fSJcDt6fZFQP0XpbSB6bfOXEbLZNmmirmf92xY6/qjqpdrK8OuY+KXtI2khVMkC6R/Nv3RXiRr6V5deHRWe1ONmu2nZXKUSh3Z96Epr7kORrUE1DHxR8RBZQZizZOnDNFPy+SolDomx9iE11wH022vHSbKM81+Oj/PXDJ/KMqclnlsbCzGx8fLOpwNUOtsyKUHa4ph+j8vaXVEjE3en2eStpuBE4GHgRfT3QF4Pn7bZdhHiZqVadi7zPJc3D0jIuYVHomNtGEeJWpWpjJHlfcqT+L/rqR5EbG28GhspLmbxGw0uszyJP5lJMl/C/A8SZdPRMSJhUZmlfI88GaD0Wvnz0e//DAAV7/j+IHHlCfx3wy8H1jD7hq/1Vy2dDOsH1fNhl0/nT9rn+i+bGg/8iT+jRFxV2ER2NDyPPBm/RnWwV95Ev8PJX0O+DJJqQcot53TzMwGJ0/ifzlJwn9rZp/bOWsibxtmESNszawaeebj72mtM0mvAj4DHEFybWBpRNwg6RCS2T3nAhuAd0fE070cw/qXpw2ziBG2ZladPAO4biFddjErIv7DFA/dAXw4Iu6XdBCwWtI9wO8AKyNiiaTFwGLgymlHbgMznUnQzGz05Sn1fCVzez/gN4EnpnpQRGwmXZs3IrZJWgccCSwAzk7vtgxYhRN/Xzq1XnYq42S5RGNWjSrLp3lKPXdmtyXdBnxzOgeRNBc4BbgPODz9o0BEbJY0q8NjFgGLAObM8dlmN51aLzuVcbJcojErX9Xl0zxn/JMdC+TOxJIOBO4E/iAinpWU63ERsZRkkXfGxsamnkmu4Tq1Xno0rdnwqbp8mqfG35qXv2ULOUszkvYhSfq3Zto/n5Q0Oz3bnw1snWbMtZRnpOx0O3BcxjGzdqZcejEiDoqImZmv4yaXf9pRcmp/E7AuIq7L/OguYGF6eyGwopfA6ya7ZF6nJd2y98ma3IHTSvYu45hZO91W4Or6OSQiNk7x3GeRTvUg6cF0358AS4DbJV0GbAQuzh9uveUZKesOHDPrV7dSz1fZvfRiSwCHAbOAGd2eOCK+PemxWedMI0YzMxugbksvzs9up505VwLnAh8vNCpr2+rlmr2ZDUKei7vHAn8KnA78D+DyiHih6MCarFNd3jV7MxuEbjX+E0gS/vHAXwKXRcTOsgJrMtfpzaxI3c74/x/wGEmt/zTgtGwPfkRcXmxo9Zdtz3QZx8zK0i3xTzUXj/UpO7LWZRwzK0u3i7vLygykqTyy1szK1suUDdYHl3fMrGpTjty1wcqOvnV5x8yq4DP+Cri8Y2ZVytPHfxxwI8l0yidIOhG4MCKuLTy6mnB5x8yGSZ5Sz6eAq4AXACLiIeCSIoOqG5d3zGyY5Cn17B8R35s0j/6OguKpLZd3zGxY5Dnjf0rSa0jn5Jd0EemSimZmNnrynPF/gGQlrNdJehx4BPjtQqMyM7PC5Flzdz1wrqQDgL0iYlvxYZmZWVHydPW8DHgXMBfYu1Xrj4hrCo3MzMwKkafUswL4GbAaeL7YcMzMrGh5Ev9REXFe4ZGYmVkp8nT1fEfS/KnvZmZmo6DbQiw/AF5M73OppPUkpR4BEREnlhPi6MmO1AWP1jWz4dKt1HMkcHJZgdRJdp598GhdMxsu3RL/IxHxaGmR1IxH6prZsOqW+GdJ+lCnH0bEdQXEY2ZmBeuW+GcAB5LU9G0KnoHTzEZFt8S/uZ9BWpJuBi4AtkbECem+Q4DPkwwG2wC8OyKe7vUYw8Tr55rZqOiW+Ps90/9b4H8Bn8nsWwysjIglkhan21f2eZyh4bq+mY2Cbon/nH6eOCK+JWnupN0LgLPT28uAVYxw4nd5x8xGUccBXBHx0wKOd3hEbE6ffzMwq9MdJS2SNC5pfGJiooBQ+ucFVsxsFA3tmrsRsZRkOmjGxsai4nA6cnnHzEZN2Yn/SUmzI2KzpNnA1pKP35PJI3FbXN4xs1GUZ66eQboLWJjeXkgy8+fQy5Z0slzeMbNRVNgZv6TbSC7kHippE3A1sAS4XdJlwEbg4qKOP2gu6ZhZXRSW+CPiPR1+1Fe3kJmZ9WdoL+5Wza2aZlZXZdf4R4ZbNc2srnzG34Xr+mZWRz7jNzNrGCd+M7OGceI3M2sYJ34zs4Zx4jczaxgnfjOzhml8O6cnYDOzpmn8Gb8nYDOzpmn8GT94oJaZNUvjz/jNzJqmkWf8noDNzJqskWf8noDNzJqskWf84Lq+mTVXI8/4zcyazInfzKxhnPjNzBqmMTV+d/KYmSUac8bvTh4zs0RjzvjBnTxmZtCgM34zM0s48ZuZNUwliV/SeZJ+JOknkhZXEYOZWVOVnvglzQD+N/A2YB7wHknzyo7DzKypqri4exrwk4hYDyBpObAAWDvoA330yw+z9omkk8ctnGZmiSpKPUcCj2W2N6X79iBpkaRxSeMTExN9H9QtnGZmiSrO+NVmX7xkR8RSYCnA2NjYS36ex9XvOL6Xh5mZ1VoVZ/ybgFdlto8CnqggDjOzRqoi8X8fOFbS0ZL2BS4B7qogDjOzRiq91BMROyT9PvB1YAZwc0Q8XHYcZmZNVcmUDRHxNeBrVRzbzKzpPHLXzKxhnPjNzBrGid/MrGGc+M3MGkYRPY2NKpWkCeDRHh9+KPDUAMMZBX7NzeDXXH/9vt5XR8Rhk3eOROLvh6TxiBirOo4y+TU3g19z/RX1el3qMTNrGCd+M7OGaULiX1p1ABXwa24Gv+b6K+T11r7Gb2Zme2rCGb+ZmWU48ZuZNUytE3+TFnWX9CpJ90paJ+lhSVdUHVNZJM2Q9ICkr1QdSxkkHSzpDkk/TP+9z6w6pqJJ+sP0//UPJN0mab+qYxo0STdL2irpB5l9h0i6R9KP0++vGMSxapv4G7io+w7gwxHxeuAM4AM1f71ZVwDrqg6iRDcAd0fE64CTqPlrl3QkcDkwFhEnkEznfkm1URXib4HzJu1bDKyMiGOBlel232qb+Mks6h4R24HWou61FBGbI+L+9PY2kmRQ+0WGJR0FnA98uupYyiBpJvAm4CaAiNgeEc9UG1Up9gZeLmlvYH9quGpfRHwL+Omk3QuAZentZcA7B3GsOif+XIu615GkucApwH3VRlKKvwL+GHix6kBKcgwwAdySlrc+LemAqoMqUkQ8DnwC2AhsBn4WEd+oNqrSHB4RmyE5uQNmDeJJ65z4cy3qXjeSDgTuBP4gIp6tOp4iSboA2BoRq6uOpUR7A6cCN0bEKcDPGdDH/2GV1rUXAEcDvwocIOm3q41qtNU58TduUXdJ+5Ak/Vsj4otVx1OCs4ALJW0gKeW9RdJnqw2pcJuATRHR+jR3B8kfgjo7F3gkIiYi4gXgi8AbKo6pLE9Kmg2Qft86iCetc+Jv1KLukkRS910XEddVHU8ZIuKqiDgqIuaS/Pv+Q0TU+kwwIrYAj0l6bbrrHGBthSGVYSNwhqT90//n51DzC9oZdwEL09sLgRWDeNJK1twtQwMXdT8LeD+wRtKD6b4/Sdc3tnr5IHBrekKzHri04ngKFRH3SboDuJ+ke+0Bajh1g6TbgLOBQyVtAq4GlgC3S7qM5A/gxQM5lqdsMDNrljqXeszMrA0nfjOzhnHiNzNrGCd+M7OGceI3M2uY2rZzmvVC0itJJsMCOALYSTJFAsAvIqIpA4esxtzOadaBpI8Az0XEJ6qOxWyQXOoxy0nSc+n3syX9o6TbJf2zpCWS3ifpe5LWSHpNer/DJN0p6fvp11nVvgKzhBO/WW9OIlkHYD7JiOnjIuI0kumhP5je5wbg+oj4deBdNGTqaBt+rvGb9eb7relyJf0L0JomeA3w5vT2ucC8ZHoZAGZKOihdL8GsMk78Zr15PnP7xcz2i+z+vdoLODMi/q3MwMym4lKPWXG+Afx+a0PSyRXGYraLE79ZcS4HxiQ9JGkt8HtVB2QGbuc0M2scn/GbmTWME7+ZWcM48ZuZNYwTv5lZwzjxm5k1jBO/mVnDOPGbmTXM/wfYjBu0++ZBYwAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "from ecell4 import *\n", "\n", "with reaction_rules():\n", " A + B == C | (0.01, 0.3)\n", "\n", "m = get_model()\n", "\n", "# ode.World -> gillespie.World\n", "w = gillespie.World(Real3(1, 1, 1))\n", "w.add_molecules(Species('C'), 60)\n", "\n", "# ode.Simulator -> gillespie.Simulator\n", "sim = gillespie.Simulator(w, m)\n", "obs = FixedIntervalNumberObserver(0.1, ('A', 'C'))\n", "sim.run(10.0, obs)\n", "\n", "show(obs, step=True)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "`World` and `Simulator` never change the `Model` itself, so you can switch several `Simulator`s for one `Model`." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.1" } }, "nbformat": 4, "nbformat_minor": 1 }