{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Parallel Processing in Python\n", "\n", "An obvious way to improve the performance of Python code is to make it run in parallel. Any relatively new computer will have multiple cores, which means that several processors can operate on the same data stored on memory. However, most of the code we have written in the course so far does not take advantage of more than one of them. In addition there is now widespread availability of computing clusters, such as [those offered for use by Amazon](http://aws.amazon.com/ec2/) (Vanderbilt also has [its own cluster](http://www.accre.vanderbilt.edu)). Clusters allow several computers to work together by exchanging data over a network.\n", "\n", "Parallel computing involves breaking a task into several independent sub-tasks, distributing these sub-tasks to available processors or computers, then coordinating the execution of these tasks and combining their outputs in an appropriate way.\n", "\n", "There are several different models for parallel processing, including:\n", "\n", "* **Message passing**: processes or other program components running in parallel communicate by sending and receiving messages, which allows for easy synchronization. \n", "* **Multi-threading**: within a single process, some architectures allow for the existence of several \"threads\", which execute independently, though they share the memory and state of the process in which they reside. Multi-threading is good for increasing *throughput* and reducing *latency*.\n", "* **Task farming**: a master process delegates independent calculations to available processors (task farm), and collects their outputs when complete.\n", "* **Single program, multiple data (SPMD)** Probably the most common type of parallel processing, in which tasks are split up and run simultaneously on multiple processors with different input in order to obtain results faster. All tasks execute their copy of the same program simultaneously.\n", "* **Multiple program, multiple data (MPMD)** Like SPMD, except each task may be executing a different program." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## `multiprocessing`\n", "\n", "The simplest way (though probably not the best) for performing parallel computing in Python is via the built-in process-based library for concurrent computing, called `multiprocessing`. " ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import multiprocessing\n", "import os\n", "import numpy as np" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `multiprocessing` module parallelizes by launching multiple *processes*, each with a seperate interpretor. You may have already heard about *threads*. Processes and threads are not the same:\n", "\n", "* processes are independent of one another, each having their own state, memory and address spaces\n", "* threads share resources, and are therefore interdependent; they are subunits of the same process\n", "\n", "Since processes are independent, they now have independent Global Interpreter Locks (GILs), its best to run multiprocessing on multiple CPUs. You can check how many you have on your machine:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "4" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "multiprocessing.cpu_count()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### `Process` class\n", "\n", "The `Process` class encapsulates a task running in a process. It will usually have a `target` argument that is some callable (function/method) that is executed when the process runs, along with optional arguments that can be passed to the target.\n", "\n", "A `Process` has several methods, with some useful ones being:\n", "\n", "* `is_alive`: Returns `True` if the process is running.\n", "* `join`: Waits for the process to finish its work and terminate. An optional `timeout` argument can be passed.\n", "* `run`: When the process starts, this method is called to invoke the `target`.\n", "* `terminate`: Kills the process forcibly, without any cleanup.\n", "\n", "A `Process` also has several other non-callable attributes, such as `pid`, `name` and `authkey`.\n", "\n", "Here is a trivial example of using the `Process` class, showing that each has its own process ID." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "I am working on job 0 running on PID 18970\n", "I am working on job 1 running on PID 18971\n", "I am working on job 2 running on PID 18972\n", "I am working on job 3 running on PID 18973\n", "I am working on job 4 running on PID 18974\n" ] } ], "source": [ "import os\n", "\n", "def job(n):\n", " print('I am working on job {0} running on PID {1}'.format(n, os.getpid()))\n", "\n", "jobs = []\n", "for i in range(5):\n", " p = multiprocessing.Process(target=job, args=(i,))\n", " jobs.append(p)\n", " p.start()" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[,\n", " ,\n", " ,\n", " ,\n", " ]" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "jobs" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can easily subclass `Process` to our liking:" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "class FibProcess(multiprocessing.Process):\n", " \n", " def __init__(self, n):\n", " self.n = n\n", " multiprocessing.Process.__init__(self)\n", " \n", " def run(self):\n", " a, b = 0, 1\n", " for i in range(self.n):\n", " a, b = b, a + b" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "18975\n", "0\n" ] } ], "source": [ "p = FibProcess(10000)\n", "p.start()\n", "print(p.pid)\n", "p.join()\n", "print(p.exitcode)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### `Queue` class\n", "\n", "Of course, despite being independent, we would like our processes to communicate with one another, for example, to share data. One way to facilitate this in `multiprocessing` is via the `Queue` class, a thread/process safe, first-in-first-out (FIFO) data structure that can store any serializable Python object.\n", "\n", "A `Queue`'s important methods include:\n", "\n", "* `put`: Adds item to `Queue`.\n", "* `get`: Fetches next item from `Queue`.\n", "* `close`: Closes the `Queue`, preventing the addition of more data.\n", "* `empty`: Returns True if the `Queue` is empty.\n", "* `full`: Returns True if full.\n", "* `qsize`: Retuns approximate current number of items in `Queue`.\n", "\n", "A subclass of `Queue` is the `JoinableQueue`, which has additional methods, notably `join`, which waits until all the items have been processed, blocking the addition of new items.\n" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-1\n", "foobar\n" ] } ], "source": [ "from multiprocessing import Queue\n", "\n", "q = Queue()\n", "\n", "q.put(-1)\n", "q.put('foobar')\n", "q.put(5)\n", "\n", "print(q.get())\n", "print(q.get())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here's a toy example of `Queue` usage, where a process is fed items from a function, until it receives a `None` object.\n", "\n", "First, a function to execute a task with items from a `Queue`:" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "def consumer(q):\n", " while True:\n", " thing = q.get()\n", " if thing is None:\n", " break\n", " print('Consuming {}'.format(thing))\n", " print(\"\\nFinished consuming\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Complementing this is another function that provisions the `Queue` with items:" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "def producer(sequence, q):\n", " for thing in sequence:\n", " q.put(thing)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Initialize the `Queue` and start the `consumer` process:" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Consuming 42\n", "Consuming foobar\n", "Consuming True\n", "Consuming range(0, 5)\n", "\n", "Finished consuming\n" ] } ], "source": [ "queue = multiprocessing.Queue()\n", "\n", "consumer_process = multiprocessing.Process(target=consumer, args=[queue])\n", "consumer_process.start()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Feed the `Queue` and process until finished:" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "stuff = [42, 'foobar', True, range(5)]\n", "producer(stuff, queue)" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "queue.put(None)\n", "consumer_process.join()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Two things to be aware of:\n", "\n", "1. if you `terminate` a process that is still accessing a queue, the queue may become corrupted\n", "2. you should make sure that any queue to which a given process has given data is clear before joining the process, or you will get a deadlock condition\n", "\n", "### `Pool` class\n", "\n", "We often have a task that we want to split up among several worker processes in parallel. The `Pool` class creates a number of processes and the methods for passing work to them. A `Pool` has the following key methods:\n", "\n", "* `apply`: Executes a passed function in a process and returns the result.\n", "* `apply_async`: Same as apply, but the result is returned asynchronously via a *callback*\n", "* `map`: A parallel version of `apply`, which splits an iterable of data into chunks and farms chunks out to processes.\n", "* `map_async`: Asynchronous `map`." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Example: parallel bootstrap\n", "\n", "As an example, we will choose a statistical computing task that is [*embarassingly parallel*](http://en.wikipedia.org/wiki/Embarrassingly_parallel). This function generates statistics of bootstrapped samples from a dataset." ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [], "source": [ "def bootstrap(data, nsamples, f):\n", " boot_samples = data[np.random.randint(len(data), size=(nsamples, len(data)))]\n", " return [f(s) for s in boot_samples]" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [], "source": [ "pool = multiprocessing.Pool(processes=4)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [], "source": [ "some_data = np.random.poisson(4, 25)" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [], "source": [ "result = pool.apply_async(bootstrap, (some_data, 1000, np.mean))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The result is an `ApplyResult` object:" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "result" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We may then want to take the result and calculate a confidence interval based on the quantiles." ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [], "source": [ "def bootstrap_ci(boot, alpha=0.05): \n", " \n", " lower_index = int(np.floor((0.5*alpha)*len(boot)))\n", " upper_index = int(np.floor((1.-0.5*alpha)*len(boot)))\n", " return boot[lower_index], boot[upper_index]" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(2.52, 4.0800000000000001)" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bootstrap_ci(np.sort(result.get()))" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [], "source": [ "# Clean up\n", "pool.close()\n", "pool.join()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "But, since we used `Pool.apply`, this is not a parallel task. We need to use `map`." ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [], "source": [ "def mapped_bootstrap(n): \n", " return bootstrap(some_data, n, np.mean)" ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [], "source": [ "pool = multiprocessing.Pool(processes=4)" ] }, { "cell_type": "code", "execution_count": 23, "metadata": {}, "outputs": [], "source": [ "map_result = pool.map_async(mapped_bootstrap, [250]*4)" ] }, { "cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 24, "metadata": {}, "output_type": "execute_result" } ], "source": [ "map_result" ] }, { "cell_type": "code", "execution_count": 25, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[250, 250, 250, 250]" ] }, "execution_count": 25, "metadata": {}, "output_type": "execute_result" } ], "source": [ "parallel_results = map_result.get()\n", "[len(p) for p in parallel_results]" ] }, { "cell_type": "code", "execution_count": 26, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(2.4399999999999999, 4.1200000000000001)" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bootstrap_ci(np.sort(np.ravel(parallel_results)))" ] }, { "cell_type": "code", "execution_count": 27, "metadata": {}, "outputs": [], "source": [ "pool.close()\n", "pool.join()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The multiprocessing package is very useful for highly parallel tasks that do not need to communicate with each other, other than when sending the initial data to the pool of processes and when and collecting the results. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Jupyter parallel\n", "\n", "The IPython architecture consists of four components, which reside in the `ipyparallel` package:\n", "\n", "1. **Engine** The IPython engine is a Python instance that accepts Python commands over a network connection. When multiple engines are started, parallel and distributed computing becomes possible. An important property of an IPython engine is that it blocks while user code is being executed. \n", "\n", "2. **Hub** The hub keeps track of engine connections, schedulers, clients, as well as persist all task requests and results in a database for later use.\n", "\n", "3. **Schedulers** All actions that can be performed on the engine go through a Scheduler. While the engines themselves block when user code is run, the schedulers hide that from the user to provide a fully asynchronous interface to a set of engines.\n", "\n", "4. **Client** The primary object for connecting to a cluster.\n", "\n", "![IPython architecture](images/ipython_architecture.png)\n", "(courtesy Min Ragan-Kelley)\n", "\n", "This architecture is implemented using the ØMQ messaging library and the associated Python bindings in `pyzmq`.\n", "\n", "## Start your engines!\n", "\n", "In order to use IPython for parallel computing, you will need to start the IPython\n", "controller and two or more IPython engines. The simplest way of doing this is\n", "with the [clusters tab](/#tab2), or you can use the `ipcluster` command in a terminal:\n", "\n", " $ ipcluster start --n=4\n", "\n", "This command will start 4 IPython engines on the current host, which is appropriate for many desktop multicore systems. You can also setup IPython clusters that span many nodes in a computing cluster, but this is beyond the scope of this lecture, but you can get more information from \n", "[the IPython.parallel docs](http://ipython.org/ipython-doc/dev/parallel/parallel_process.html)..\n", "\n", "To use the IPython cluster in our Python programs or notebooks, we start by creating an instance of `ipyparallel.Client`:" ] }, { "cell_type": "code", "execution_count": 28, "metadata": {}, "outputs": [], "source": [ "from ipyparallel import Client" ] }, { "cell_type": "code", "execution_count": 32, "metadata": {}, "outputs": [], "source": [ "cli = Client()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This creates a client using the default profile; you can pass an optional `profile=\"my_profile\"` argument if you have a different one running.\n", "\n", "Using the `ids` attribute we can retreive a list of ids for the IPython engines in the cluster:" ] }, { "cell_type": "code", "execution_count": 33, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[0, 1]" ] }, "execution_count": 33, "metadata": {}, "output_type": "execute_result" } ], "source": [ "cli.ids" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can use a `DirectView` object for execution of tasks, which an be accessed simply by indexing the client:" ] }, { "cell_type": "code", "execution_count": 34, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 34, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dv0 = cli[0]\n", "dv0" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The above shows just a single engine, but we want all of them:" ] }, { "cell_type": "code", "execution_count": 35, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 35, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dview = cli[:]\n", "dview" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can get a view on whatever combination of engines we want:" ] }, { "cell_type": "code", "execution_count": 36, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 36, "metadata": {}, "output_type": "execute_result" } ], "source": [ "cli[::2]" ] }, { "cell_type": "code", "execution_count": 37, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 37, "metadata": {}, "output_type": "execute_result" } ], "source": [ "cli[1::2]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `block` flag specifies whether to wait for the result, or return an `AsyncResult` object immediately:" ] }, { "cell_type": "code", "execution_count": 38, "metadata": {}, "outputs": [], "source": [ "dview.block = True" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Finally, since we want to use IPython's parallel magic commands, we set the `DirectView` to be `active`:" ] }, { "cell_type": "code", "execution_count": 39, "metadata": {}, "outputs": [], "source": [ "dview.activate()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Each of these engines are ready to execute tasks. We can selectively run code on individual engines. For example, we can simply use `os.getpid` to return the process ID that the engine is running on. Here is the notebook process:" ] }, { "cell_type": "code", "execution_count": 40, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "18949" ] }, "execution_count": 40, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import os\n", "os.getpid()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here is a single engine's process ID:" ] }, { "cell_type": "code", "execution_count": 41, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "18993" ] }, "execution_count": 41, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dv0.apply_sync(os.getpid)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "And here are all the engines, run simultaneously:" ] }, { "cell_type": "code", "execution_count": 42, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[18993, 18994]" ] }, "execution_count": 42, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dview.apply_sync(os.getpid)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's now consider a useful function that we might want to run in parallel. Here is a version of the approximate Bayesian computing (ABC) algorithm that we have seen in previous lectures." ] }, { "cell_type": "code", "execution_count": 43, "metadata": {}, "outputs": [], "source": [ "import numpy\n", "\n", "def abc(y, N, epsilon=[0.2, 0.8]):\n", " \n", " trace = []\n", "\n", " while len(trace) < N:\n", "\n", " # Simulate from priors\n", " mu = numpy.random.normal(0, 10)\n", " sigma = numpy.random.uniform(0, 20)\n", "\n", " x = numpy.random.normal(mu, sigma, 50)\n", "\n", " #if (np.linalg.norm(y - x) < epsilon):\n", " if ((abs(x.mean() - y.mean()) < epsilon[0]) &\n", " (abs(x.std() - y.std()) < epsilon[1])):\n", " trace.append([mu, sigma])\n", "\n", " return trace" ] }, { "cell_type": "code", "execution_count": 44, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "\n", "y = np.random.normal(4, 2, 50)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's try running this on one of the cluster engines:" ] }, { "cell_type": "code", "execution_count": 45, "metadata": {}, "outputs": [ { "ename": "RemoteError", "evalue": "NameError(name 'numpy' is not defined)", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m\u001b[0;31mNameError\u001b[0m Traceback (most recent call last)\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m()\u001b[0m", "\u001b[0;32m\u001b[0m in \u001b[0;36mabc\u001b[0;34m(y, N, epsilon)\u001b[0m", "\u001b[0;31mNameError\u001b[0m: name 'numpy' is not defined" ] } ], "source": [ "dv0.block = True\n", "dv0.apply(abc, y, 10)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This fails with a `NameError` because NumPy has not been imported on the engine to which we sent the task. Each engine has its own namespace, so we need to import whatever modules we will need prior to running our code:" ] }, { "cell_type": "code", "execution_count": 46, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 46, "metadata": {}, "output_type": "execute_result" } ], "source": [ "cli[0].execute(\"import numpy\")" ] }, { "cell_type": "code", "execution_count": 47, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[[4.729794501294109, 1.955124275924458],\n", " [4.4084785969434535, 3.2710801621330976],\n", " [5.1924225245376245, 2.6639991076064917],\n", " [4.914692595106862, 2.943457251742494],\n", " [5.598432333216517, 2.1019281233409504],\n", " [4.847396141714309, 1.7678739019246792],\n", " [4.775530891472667, 2.3076434075853958],\n", " [5.181119280045877, 2.295797728233673],\n", " [5.11601758179679, 1.8056184659746566],\n", " [5.359761108209431, 2.880568354850497]]" ] }, "execution_count": 47, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dv0.apply(abc, y, 10)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A more efficient way is to simultaneously import modules into the local and the engine namespaces simultaneously, using a context manager:" ] }, { "cell_type": "code", "execution_count": 48, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "importing numpy on engine(s)\n" ] } ], "source": [ "with dview.sync_imports():\n", " import numpy" ] }, { "cell_type": "code", "execution_count": 49, "metadata": {}, "outputs": [], "source": [ "t = dview.apply(abc, y, 10)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Easier yet, you can use the parallel cell magic to import everywhere:" ] }, { "cell_type": "code", "execution_count": 50, "metadata": {}, "outputs": [], "source": [ "%%px\n", "import numpy" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can also use the `require` decorator for functions that you wish to use on engines." ] }, { "cell_type": "code", "execution_count": 51, "metadata": {}, "outputs": [], "source": [ "from ipyparallel import require\n", "\n", "@require(\"numpy\")\n", "def abc(y, N, epsilon=[0.2, 0.8]):\n", "\n", " trace = []\n", "\n", " while len(trace) < N:\n", "\n", " # Simulate from priors\n", " mu = numpy.random.normal(0, 10)\n", " sigma = numpy.random.uniform(0, 20)\n", "\n", " x = numpy.random.normal(mu, sigma, 50)\n", "\n", " #if (np.linalg.norm(y - x) < epsilon):\n", " if ((abs(x.mean() - y.mean()) < epsilon[0]) &\n", " (abs(x.std() - y.std()) < epsilon[1])):\n", " trace.append([mu, sigma])\n", "\n", " return trace" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A simple way to run code on an engine is via the `execute` method:" ] }, { "cell_type": "code", "execution_count": 52, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 52, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dv0.execute('x=3')" ] }, { "cell_type": "code", "execution_count": 53, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3" ] }, "execution_count": 53, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dv0['x']" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Data transfer\n", "\n", "We will often want to send data to our engines, or retrieve objects from them. `DirectView` has `push` and `pull` methods for achieving this.\n", "\n", "Recall that Python namespaces are just dictionaries. So, we can update an engine's namespace by pushing a dictionary:" ] }, { "cell_type": "code", "execution_count": 54, "metadata": {}, "outputs": [], "source": [ "dv0.push({'foo': -3, 'bar': np.arange(10)})" ] }, { "cell_type": "code", "execution_count": 55, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[3, array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])]" ] }, "execution_count": 55, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dv0.pull(('x', 'bar'))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Additionally, `DirectView` objects also have `scatter` and `gather` methods, to distribute data among engines. `scatter` accepts any container or Numpy `array` type, while `gather` assembles the respective return objects in the Client." ] }, { "cell_type": "code", "execution_count": 56, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([ 3.49322828, 5.30494559, 3.6189018 , 6.75997796,\n", " -0.21334747, 2.42591587, 7.25860593, 5.79365828,\n", " 4.3719095 , 4.8857659 , 5.41098905, 2.65143428,\n", " 5.23338454, 7.31421719, 5.97221326, 7.15798182,\n", " 7.17905784, 8.16824998, 6.11871443, 4.17446774,\n", " 9.53597307, 8.11262697, 2.08434239, 3.6807146 ,\n", " 6.60156969, 3.77160205, 1.10850621, 3.0781838 ,\n", " 5.38484312, 5.30269602, 6.78573415, 4.60953635,\n", " 5.43390477, 4.09052607, 7.48890093, 4.50875377,\n", " 4.77165189, 6.6224069 , 2.80725959, 3.68837268,\n", " 6.85987156, 4.05963832, 4.95305724, 3.10880789,\n", " 2.47886851, 2.65170513, 4.8583221 , 3.79538796,\n", " 4.00783596, 10.26270091])" ] }, "execution_count": 56, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Some Gaussian data\n", "y" ] }, { "cell_type": "code", "execution_count": 57, "metadata": { "collapsed": true }, "outputs": [], "source": [ "dview = cli[:2]" ] }, { "cell_type": "code", "execution_count": 58, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[array([ 3.49322828, 5.30494559, 3.6189018 , 6.75997796, -0.21334747,\n", " 2.42591587, 7.25860593, 5.79365828, 4.3719095 , 4.8857659 ,\n", " 5.41098905, 2.65143428, 5.23338454, 7.31421719, 5.97221326,\n", " 7.15798182, 7.17905784, 8.16824998, 6.11871443, 4.17446774,\n", " 9.53597307, 8.11262697, 2.08434239, 3.6807146 , 6.60156969]),\n", " array([ 3.77160205, 1.10850621, 3.0781838 , 5.38484312,\n", " 5.30269602, 6.78573415, 4.60953635, 5.43390477,\n", " 4.09052607, 7.48890093, 4.50875377, 4.77165189,\n", " 6.6224069 , 2.80725959, 3.68837268, 6.85987156,\n", " 4.05963832, 4.95305724, 3.10880789, 2.47886851,\n", " 2.65170513, 4.8583221 , 3.79538796, 4.00783596, 10.26270091])]" ] }, "execution_count": 58, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Send to engines\n", "dview.scatter('y', y)\n", "dview['y']" ] }, { "cell_type": "code", "execution_count": 59, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 59, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Remote execution of function\n", "dview.execute('sum_y = sum(y)')" ] }, { "cell_type": "code", "execution_count": 60, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "249.58457238955828" ] }, "execution_count": 60, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Aggregation on client\n", "sum(dview.gather('sum_y'))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `map` method essentially combines `scatter` and `gather` into a single call:" ] }, { "cell_type": "code", "execution_count": 61, "metadata": {}, "outputs": [ { "data": { "text/plain": [ ">" ] }, "execution_count": 61, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dview.map(lambda x: sum(x**2), np.split(y, 5))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Load balancing\n", "\n", "The `DirectView` objects we have used so far strictly allocate particular tasks to particular engines. This is often inefficient, when tasks take variable amounts of time, leaving some engines idle while some are overworked. We can use a **load balanced** view to distribute memory approximately equally among engines, to minimize idle time." ] }, { "cell_type": "code", "execution_count": 62, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 62, "metadata": {}, "output_type": "execute_result" } ], "source": [ "lb_view = cli.load_balanced_view()\n", "lb_view" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A `LoadBalancedView`, though it works with all the engines (or specified subsets of engines), behaves as if it is working with a single engine.\n", "\n", "If you do not specify the engines when the `LoadBalancedView` is created, it will use all the engines that are available when it assigns tasks." ] }, { "cell_type": "code", "execution_count": 63, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Task 0 ran on process 18993\n", "Task 1 ran on process 18994\n", "Task 2 ran on process 18993\n", "Task 3 ran on process 18994\n", "Task 4 ran on process 18993\n", "Task 5 ran on process 18994\n", "Task 6 ran on process 18993\n", "Task 7 ran on process 18994\n", "Task 8 ran on process 18993\n", "Task 9 ran on process 18994\n" ] } ], "source": [ "for i in range(10):\n", " pid = lb_view.apply_sync(os.getpid)\n", " print('Task {0} ran on process {1}'.format(i, pid))" ] }, { "cell_type": "code", "execution_count": 64, "metadata": {}, "outputs": [], "source": [ "%%px\n", "\n", "def abc(y, N, epsilon=[0.2, 0.8]):\n", "\n", " trace = []\n", "\n", " while len(trace) < N:\n", "\n", " # Simulate from priors\n", " mu = numpy.random.normal(0, 10)\n", " sigma = numpy.random.uniform(0, 20)\n", "\n", " x = numpy.random.normal(mu, sigma, 50)\n", "\n", " #if (np.linalg.norm(y - x) < epsilon):\n", " if ((abs(x.mean() - y.mean()) < epsilon[0]) &\n", " (abs(x.std() - y.std()) < epsilon[1])):\n", " trace.append([mu, sigma])\n", "\n", " return trace" ] }, { "cell_type": "code", "execution_count": 65, "metadata": {}, "outputs": [], "source": [ "tasks = lb_view.map_async(lambda n: abc(y, n), [20]*5)" ] }, { "cell_type": "code", "execution_count": 66, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "['551c609d-e842-4fe3-b25d-f6bf1425991b',\n", " '6c30b36f-47e6-4d72-9de3-56d464ae8756',\n", " '68d30d90-84dd-4b0f-a83d-41454668a660',\n", " 'a3dd0b4a-9497-495c-a85b-cd40860ae881',\n", " 'f7eeb9fd-96dc-4ad5-bcd2-12a22c4efabf']" ] }, "execution_count": 66, "metadata": {}, "output_type": "execute_result" } ], "source": [ "tasks.msg_ids" ] }, { "cell_type": "code", "execution_count": 67, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[ 5.13479778, 1.99074744],\n", " [ 5.83737441, 2.80446401],\n", " [ 4.93086704, 2.78291914],\n", " [ 5.07532547, 2.89403202],\n", " [ 5.57137057, 1.58596582],\n", " [ 5.22463293, 2.65868058],\n", " [ 4.52186877, 2.81636134],\n", " [ 5.49863421, 1.61619842],\n", " [ 5.24600375, 1.74222543],\n", " [ 5.5639975 , 1.67952349]])" ] }, "execution_count": 67, "metadata": {}, "output_type": "execute_result" } ], "source": [ "result = np.concatenate(tasks.get())\n", "result[:10]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Another way that you can dispatch tasks to engines is via the `parallel` decorator. This decorator is a method of the `DirectView` class that controls our engine pool. The decorated function is then disparched to the engines using the `map` method that the decorator adds to the class." ] }, { "cell_type": "code", "execution_count": 68, "metadata": {}, "outputs": [], "source": [ "@lb_view.parallel(block=True)\n", "def abc(y, N, epsilon=[0.2, 0.8]):\n", " \n", " trace = []\n", "\n", " while len(trace) < N:\n", "\n", " # Simulate from priors\n", " mu = numpy.random.normal(0, 10)\n", " sigma = numpy.random.uniform(0, 20)\n", "\n", " x = numpy.random.normal(mu, sigma, 50)\n", "\n", " #if (np.linalg.norm(y - x) < epsilon):\n", " if ((abs(x.mean() - y.mean()) < epsilon[0]) &\n", " (abs(x.std() - y.std()) < epsilon[1])):\n", " trace.append([mu, sigma])\n", "\n", " return trace" ] }, { "cell_type": "code", "execution_count": 69, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[[[5.099682864276421, 1.9600016451620506],\n", " [5.167619695345338, 2.7749820165169115],\n", " [4.65741013915943, 2.744966364750392],\n", " [4.631512967150718, 2.2170096764352754],\n", " [4.912050459991287, 1.9205491182111745],\n", " [4.535967735211762, 2.091694516734597],\n", " [5.300889064301683, 2.9545139934122444],\n", " [4.756953973918762, 1.6123588265695288],\n", " [5.506509682307172, 2.643826970085197],\n", " [5.143915920591466, 1.596264678271515],\n", " [5.058542066042192, 2.2487972003157997],\n", " [4.893977094678111, 1.7820732690578045],\n", " [5.089788438659532, 1.7566732396821627],\n", " [4.828081550671436, 3.4737089991803827],\n", " [5.097180224631751, 2.45113231616515],\n", " [4.830565010283392, 2.4955384226959487],\n", " [5.251808785167854, 1.3646424428460269],\n", " [5.268592126423453, 1.3681561718947721],\n", " [5.103207452428124, 2.378189165650919],\n", " [4.357701820937493, 3.088604471573273],\n", " [5.5757420529349035, 2.3391659934899067],\n", " [4.51170901539765, 1.787305543407789],\n", " [4.805618531567692, 1.9544731262002468],\n", " [4.6593154012349265, 1.496651172318808],\n", " [4.98527699864045, 1.4515522597509323]],\n", " [[5.043037866054284, 2.760736466148954],\n", " [5.501902369652941, 1.8471269647260247],\n", " [4.027288634790987, 2.8583277821183906],\n", " [4.571925195962599, 1.5139056121262096],\n", " [4.655517999519201, 2.0083691902368828],\n", " [5.922888858500221, 1.894627003487559],\n", " [5.09400542931459, 2.0239542911585473],\n", " [4.714090660624263, 2.143535239195755],\n", " [5.202480846756165, 1.5502880999263091],\n", " [4.939644333258536, 1.3568562590203004],\n", " [5.660054555161992, 3.5945100842860356],\n", " [4.796572170783254, 1.3708542397143697],\n", " [4.59672555079491, 2.0861299759442],\n", " [4.463181893885676, 2.9650573934704405],\n", " [5.201618826531734, 1.8309475342930348],\n", " [5.016314691145761, 1.6816660794824845],\n", " [4.690252531051175, 2.0730287337336684],\n", " [3.7847527215534473, 3.5357302516130074],\n", " [5.59678329877705, 2.6486479235522697],\n", " [4.999557808281289, 2.6044297087205526],\n", " [5.122628097911353, 2.2083817217107615],\n", " [4.259974312817068, 2.802690380739563],\n", " [5.36583837194772, 1.711496734943363],\n", " [4.708034862974067, 2.565228125633219],\n", " [4.314408007635806, 1.6187541158897978]],\n", " [[4.809431881211897, 1.5395748158488765],\n", " [5.355231939475894, 2.1214769468543615],\n", " [5.382728115400916, 3.1634158556776137],\n", " [4.6909153950346125, 2.3036456345367573],\n", " [5.129948008754114, 3.2192946990161175],\n", " [4.942738744191026, 1.3527288932880288],\n", " [4.949272644666744, 2.3602313992713952],\n", " [4.894695108329576, 1.2246636397013133],\n", " [5.173654699161258, 2.1567531376704485],\n", " [5.292834233029171, 2.7414480366564797],\n", " [5.010195842122052, 2.5173737425280707],\n", " [4.708587344558397, 2.4279067094023743],\n", " [4.709829846038776, 2.1616290060174115],\n", " [4.721943786965558, 2.8152406502777416],\n", " [5.375790379081728, 1.671267754801502],\n", " [5.1606636252057045, 2.7399390525301337],\n", " [5.137865741913216, 1.6173303592824428],\n", " [5.108503645292948, 1.9547719480972203],\n", " [5.163662951051209, 1.9056116164427928],\n", " [4.505256901736298, 2.311519378292717],\n", " [5.329191548986296, 2.799744816037779],\n", " [4.958438952975342, 1.7562098802836679],\n", " [4.64505225809283, 1.6059599395848534],\n", " [5.330581149600655, 1.9474274174843953],\n", " [5.313382491523115, 1.7845740929617326]],\n", " [[4.842211947920385, 2.627359294119034],\n", " [5.099456335568907, 2.5681590304834967],\n", " [4.718634004127261, 2.48460626196916],\n", " [4.618069615966375, 2.5553027844785325],\n", " [4.767730776300696, 2.0357189308612966],\n", " [5.04003346275795, 2.7533164880175653],\n", " [5.1196081621437814, 1.9081484014271233],\n", " [4.936752863317418, 2.2078852654667425],\n", " [4.9247329122451555, 1.8777748109992132],\n", " [5.110593181374324, 1.494115295816023],\n", " [4.790744834434554, 1.5233470899958856],\n", " [4.261953779448408, 2.6113224156865367],\n", " [5.4863808242316425, 1.722166186742471],\n", " [5.002271190364263, 1.403117428773435],\n", " [5.184409141048819, 1.4056148769259513],\n", " [4.973469075385752, 1.5845828851372623],\n", " [5.09768752043273, 3.0508674312682205],\n", " [4.7866121575582925, 2.491539974940038],\n", " [5.574317182277588, 2.206431566942808],\n", " [4.783901351479342, 3.0675032159658233],\n", " [4.892326071783362, 3.0327085122224484],\n", " [4.365078355384927, 1.4716472281250947],\n", " [4.95156516827712, 2.1233560525511885],\n", " [4.545238260428942, 2.704805297715076],\n", " [5.340427748276136, 1.7598927934344788]]]" ] }, "execution_count": 69, "metadata": {}, "output_type": "execute_result" } ], "source": [ "abc.map([y]*4, [25]*4)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Parallel magics\n", "\n", "The `%px` cell magic is a \"parallel execution\" statement, which will run the code in that cell on all the engines." ] }, { "cell_type": "code", "execution_count": 70, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[stdout:0] 18993\n", "[stdout:1] 18994\n" ] } ], "source": [ "%%px \n", "import os\n", "print(os.getpid())" ] }, { "cell_type": "code", "execution_count": 71, "metadata": {}, "outputs": [], "source": [ "%px b = numpy.random.random()" ] }, { "cell_type": "code", "execution_count": 72, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\u001b[0;31mOut[0:5]: \u001b[0m0.5779141088162374" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\u001b[0;31mOut[1:5]: \u001b[0m0.2809837962562448" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "%px b" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`%pxresult` displays the output of the last request:" ] }, { "cell_type": "code", "execution_count": 73, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\u001b[0;31mOut[0:5]: \u001b[0m0.5779141088162374" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\u001b[0;31mOut[1:5]: \u001b[0m0.2809837962562448" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "%pxresult" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `%pxconfig` magic allows you to configure blocking for the parallel magics." ] }, { "cell_type": "code", "execution_count": 74, "metadata": {}, "outputs": [], "source": [ "# Switch to asynchronous\n", "%pxconfig --block" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Remember that each engine is just another IPython, so anyting you can do in an IPython session you can also do on an engine." ] }, { "cell_type": "code", "execution_count": 75, "metadata": {}, "outputs": [], "source": [ "%px %matplotlib inline" ] }, { "cell_type": "code", "execution_count": 76, "metadata": {}, "outputs": [], "source": [ "%px samples = abc(y, 100)" ] }, { "cell_type": "code", "execution_count": 77, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[output:0]" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgAAAAFyCAYAAACDemKtAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAHNlJREFUeJzt3X+UZGV95/H3FyH0gDps6AAitogTx46urtOuiIoQRU3I\nEQUStbWPItENhk3YibtxOYoaPYmJRseIsMluPCoZbQMKUVcEfyAqJsJCR0VtWlSwFGGcAm0YoPk1\n3/3j3tGaZnqm61dXdz/v1zl1eurWc5/7rTszXZ+697nPjcxEkiSVZa9BFyBJkpaeAUCSpAIZACRJ\nKpABQJKkAhkAJEkqkAFAkqQCGQAkSSqQAUCSpAIZACRJKpABQJKkAhkApAGKiFdFxPaWx90RMRMR\nZ0fEQS3tjqlfP2kP694UEZdExJ9ExEPbqOONEfHJiLil7uvNu2l7XERcFhFbI+LnEXFlREzsot1B\nEfHBiNgSEXdFxDUR8fsL9Pmy+vW7I+JnEfGPEXHgvDZDEfGBiLg2In4REXdExDci4k8jYu/FvldJ\nFf/TSIOXwFnAjcAQ8CzgdcDvRsQTM3Oupd3u1t0HOAQ4Fngv8GcRcUJmXruIGt4O3AxMAS9YqFFE\nnABcBPwr8JZ6+y8BzouIAzPz7+p2DwO+BvxGXcuWut35EfHyzPxYS5+vA84BPg9sBA4D/hswFhFH\nZua9ddM1wCjwmfr9bgeeAWwCngY8KIRI2o3M9OHDx4AewKuAB4AN85b/bb38pfXzY6g+8E7a07r1\na8cCdwI/BPZdRB0j9c8D6+28eYF2lwI/BvZuWfYQ4Hrg31uW/Y+6tmNalgVwJXDTjvWpQsttwGXz\ntvN7dR2nL6L299XbOnjQf58+fKykh6cApOXpMqoPzMd0snJmXk71rf7RLOKbcWY2Ftn1w4GfZ+b9\nLes+ADSBu1vaPQvYmplfbmmXwPlURymOqRc/ETigXt5az2eAbcDLFlHTj+qfaxf5HiThGABpuVpX\n/7y1iz7+iSpEPL/7cn7pcuAJEfG2iHhsRBwREWcBY8A7W9rty86BYIe76prGWtqxQNu7gafMXxgR\n+0TEgRFxWEScCLye6pTA9zt4P1KxDADS8rC2/lB7ZES8lOq8/l3A/+20w8y8CZgFHtujGgHeBlwA\nvJHqsP/3gT8HTs7Mf2lpNwMcFhGPmrf+s6nGDTyyfn59/fyZrY0iYj3V+IE1EfEf5vVxErAVaACf\noDolcUJmbu/urUllcRCgNHgBfLHleVJ9ox3PzJu77Hsb8LAu+2h1L/A9qhBwIdX5//8CfCQijsvM\nq+p2/wicBlwQERupBgG+FHhx/foagMy8NSLOB14VEddRDTA8jOq8/r1UYwTWAD9vqeEy4DiqUwfP\nBZ4MLPqKB0kVA4A0eAn8MdW34fuBLZk506O+H0r14dsr5wBPy8wNOxZExAXAd4C/A44CyMxrI2Ic\n+HvgCqqQczNwRr1sW0uff0R19cO7qAY/JrCZagDji+e1JTO3UoUAgAsj4kzg8xGxLjN/1sP3Kq1q\nngKQlof/l5mXZeZXevXhHxGPpBoY15Nz4xGxD3Aq1WV4v1QPCPws8NS6zY7lFwKHUl2i93SqAYk3\n1C9/r6Xd7Zl5Yv36s4HDM/NVVIMFt2bm7Xso7eNUQedFnb87qTweAZBWr1dSfZu+pEf9HUj1O+Mh\nu3htH6ovFDt9qajDwTU7nkfE8+qavjC/g8z8CfCTut0BVAMFL1hEXWvqn14FILXBIwDSKhQRzwHe\nRHUY/aM96vZnwC+AE1tn3qtnHHwhMJ2Z9+ympnVUh/s/nZl7OirxDqqg8d6W9Q9coO1rqULF1Yt5\nE5IqHgGQBi+6XPf4iBil+v98MPAc4HlUh9tPyF/NpLdwJ9VUvo8G9q8XHRMRb6z/fF5m/jgzt0fE\n31LNL3BlRJxXb/MPqUb1//m8Pr9D9Q2+ARxBNSiwSTXLYWu7N1DNB3Al1RiIE6kG+b0xM69paToR\nEacB/0IVbB5GNWvhccCn6rkPJC2SAUAavF1N8bvYdgn8Rf3ne6lm1bsW+FPgQ5l55yL7/kOq8+87\n+jy2fgB8lepSOzLzryLih1SD+d5MdR3/t3jwZYAA3wBOoQolTeBjwFszszmv3bVUg/1eSPWt/1vA\nH9RjCFpdQTXI8GVU4wPuo7rccCPw/kW+T0m1qCbnkiRJJWlrDEBEnBkRV0XE7fUdvi6KiMfNa3P5\nvDuUPRAR5/a2bEmS1I12BwEeDZwNHEl13m0f4HMRsaalTQL/m+qw3yHAI5h3blCSJA1WW2MAMvP4\n1ucRcQrVyOAxqvNzO9xVT9YhSZKWoW4vAzyA6hv/bfOWvyIitkbEtRHxV/OOEEiSpAHreBBgRATw\naeBhmXlMy/LXUN2e86fAk6juEHZlZv7+Av0cSHUpz43AXEfFSJJUpiHgcODSzGzr7qHdBID/RfXB\n/czd3bAkIn6batavdZl5wy5efznwkY6KkCRJAK/IzLYm/epoHoCIeD9wPHD0Iu5WdiXVZCXr+NU8\n4K1uBNi8eTOjo6OdlFOsjRs3smnTpkGXsaK4zzrjfmuf+6wz7rf2TE9PMzExAfVnaTvaDgD1h/+L\ngGMys7GIVZ5CNU5goaAwBzA6OsqGDRsWaKJdWbt2rfusTe6zzrjf2uc+64z7rWNtn0JvKwDU1/OP\nAycAd0bEwfVLs5k5FxFHAC8HLgZupbpP93uAL2fmt9stTpIk9Ue7RwBOo/o2f/m85a8GzqOaivQ4\nqmlC96eaPvQC4C+7qlKSJPVUu/MA7Paywfp2nsd2U5AkSeo/bwe8go2Pjw+6hBXHfdYZ91v73Ged\ncb8tnYHfDCgiNgDXXHPNNQ78kCSpDVNTU4yNjQGMZeZUO+t6BECSpAIZACRJKpABQJKkAhkAJEkq\nUEdTAUuS2tdoNGg2m22vd88997Dvvvt2tM3h4WFGRkY6WlermwFAkpZAo9Fg/fpR5ubu6mDthwAP\ndLTdoaH9mJmZNgToQQwAkrQEms1m/eG/GWjnxmcXA2d1sB7ANHNzEzSbTQOAHsQAIElLahRoZ86T\n6Q7Xk3bPQYCSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoEMAJIkFcgA\nIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACBJ\nUoEMAJIkFWjvQRegla3RaNBsNrvuZ3h4mJGRkR5UJElaDAOAOtZoNFj/+PXM3T3XdV9Da4aYuW7G\nECBJS8QAoI41m83qw/8kYLibjmDuwjmazaYBQJKWiAFA3RsGDh10EZKkdjgIUJKkAhkAJEkqkAFA\nkqQCGQAkSSqQAUCSpAIZACRJKpABQJKkAhkAJEkqkAFAkqQCGQAkSSqQAUCSpAIZACRJKpABQJKk\nAhkAJEkqkAFAkqQCGQAkSSpQWwEgIs6MiKsi4vaI2BIRF0XE4+a12TcizomIZkTcEREfj4iDelu2\nJEnqRrtHAI4GzgaOBI4D9gE+FxFrWtq8F/g94GTg2cChwCe6L1WSJPXK3u00zszjW59HxCnAz4Ax\n4IqIeDhwKvCyzPxy3ebVwHREPC0zr+pJ1ZIkqSvdjgE4AEjgtvr5GFWo+OKOBpk5AzSAo7rcliRJ\n6pGOA0BEBNXh/isy87v14kOAezPz9nnNt9SvSZKkZaCtUwDznAv8FvCsRbQNqiMFC9q4cSNr167d\nadn4+Djj4+MdFyhJ0moxOTnJ5OTkTstmZ2c77q+jABAR7weOB47OzJ+2vHQL8GsR8fB5RwEOojoK\nsKBNmzaxYcOGTsqRJGnV29WX4qmpKcbGxjrqr+1TAPWH/4uA387MxryXrwHuB57b0v5xwAjwbx1V\nKEmSeq6tIwARcS4wDpwA3BkRB9cvzWbmXGbeHhEfAN4TET8H7gDeB3zNKwAkSVo+2j0FcBrVufzL\n5y1/NXBe/eeNwAPAx4F9gUuA0zsvUZIk9Vq78wDs8ZRBZt4D/En9kCRJy5D3ApAkqUAGAEmSCmQA\nkCSpQAYASZIKZACQJKlABgBJkgrUzb0ApJ6anp7uav3h4WFGRkZ6VI0krW4GAA3eNiBgYmKiq26G\n1gwxc92MIUCSFsEAoMGbo5pf8iRguMM+mjB34RzNZtMAIEmLYADQ8jEMHDroIiSpDA4ClCSpQAYA\nSZIKZACQJKlABgBJkgpkAJAkqUAGAEmSCmQAkCSpQAYASZIK5ERAS6TRaNBsNrvux/nuJbWr0/ts\n+PtmdTMALIFGo8H6x69n7u65rvtyvntJi3czsFfH99kYGtqPmZlpf9+sUgaAJdBsNqsP/27mugfn\nu5fUpl8A24HNwGib604zNzfh75tVzACwlJzrXtJAjAIbBl2ElhkHAUqSVCADgCRJBTIASJJUIAOA\nJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJ\nBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUy\nAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklSgtgNARBwdEZ+KiJsiYntEnDDv9Q/Wy1sfF/euZEmS\n1K1OjgDsD3wDOB3IBdp8FjgYOKR+jHdUnSRJ6ou9210hMy8BLgGIiFig2T2ZubWbwiRJUv/0awzA\nsRGxJSKui4hzI+LX+7QdSZLUgbaPACzCZ4FPADcAjwXeAVwcEUdl5kKnDCRJ0hLqeQDIzPNbnn4n\nIq4FfgAcC3xpofU2btzI2rVrd1o2Pj7O+LjDByRJmpycZHJycqdls7OzHffXjyMAO8nMGyKiCaxj\nNwFg06ZNbNiwod/lSJK0Iu3qS/HU1BRjY2Md9df3eQAi4jDgQODmfm9LkiQtTttHACJif6pv8zuu\nADgiIp4M3FY/3kI1BuCWut3fAN8DLu1FwZIkqXudnAJ4KtWh/Kwf766Xfxj4Y+BJwCuBA4CfUn3w\nvzkz7+u6WkmS1BOdzAPwZXZ/6uB3Oi9HkiQtBe8FIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoEM\nAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACS\nJBXIACBJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQV\nyAAgSVKBDACSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoEMAJIkFcgA\nIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACBJ\nUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKB2g4AEXF0RHwqIm6KiO0RccIu2rwtIn4aEXdFxOcjYl1v\nypUkSb3QyRGA/YFvAKcDOf/FiHgD8F+BPwKeBtwJXBoRv9ZFnZIkqYf2bneFzLwEuAQgImIXTc4A\n3p6Zn67bvBLYArwYOL/zUiVJUq/0dAxARDwGOAT44o5lmXk7cCVwVC+3JUmSOtf2EYA9OITqtMCW\necu31K+pB6anp7vuY3h4mJGRkR5UI0laiXodABYS7GK8QKuNGzeydu3anZaNj48zPj7ez7pWlm1A\nwMTERNddDa0ZYua6GUOAJK0Qk5OTTE5O7rRsdna24/56HQBuofqwP5idjwIcBPz77lbctGkTGzZs\n6HE5q8wcVYw6CRjuop8mzF04R7PZNABI0gqxqy/FU1NTjI2NddRfTwNAZt4QEbcAzwW+BRARDweO\nBM7p5baKNgwcOugiJEkrWdsBICL2B9ZRfdMHOCIingzclpk/Bt4LvCkivg/cCLwd+AnwyZ5ULEmS\nutbJEYCnAl+iOhidwLvr5R8GTs3Md0bEfsA/AAcAXwV+NzPv7UG9kiSpBzqZB+DL7OHywcx8K/DW\nzkqSJEn95r0AJEkqkAFAkqQCGQAkSSqQAUCSpAIZACRJKpABQJKkAi3VvQBWtEajQbPZ7Hj9Xty8\nR9Ly0Onvg9J+D3S6n7xR2dIxAOxBo9Fg/ePXM3f33KBLkTRgjUaD9etHmZu7a9ClLGvd7Kehof2Y\nmZk2BCwBA8AeNJvN6sO/mxvwXE81d6KkFa3ZbNYfapuB0TbXvhg4q/dFLUOd76dp5uYmvFHZEjEA\nLFY3N+Dp/OyBpGVpFGj37qVlnQKodLKftFQcBChJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACS\nJBXIACBJUoEMAJIkFciJgCRJC+rkHgal3fdgpTIASJJ24WZgLyYmJgZdiPrEACBJ2oVfANvxvger\nlwFAkrQb3vdgtXIQoCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUIAOA\nJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJ\nBTIASJJUIAOAJEkF2nvQBfTLBRdcwEUXXdR1P/fff38PqtFSmZ6e7rqP4eFhRkZGelCNJC1fqzYA\nnPFnZ3DLtluIX4/OO9kO22/c3rui1D/bgICJiYmuuxpaM8TMdTOGAEmr2qoNACTkk5J8Tnbex93A\n3/SsIvXTHJDAScBwF/00Ye7COZrNpgFA0qq2egOAyjQMHDroIiRp+XMQoCRJBTIASJJUIAOAJEkF\nMgBIklSgngeAiHhLRGyf9/hur7cjSZI616+rAL4NPBfYcRG+s+lIkrSM9CsA3J+ZW/vUtyRJ6lK/\nxgD8ZkTcFBE/iIjNEfGoPm1HkiR1oB9HAL4OnALMAI8A3gp8JSKemJl39mF76lC38+b3Yt59aRAa\njQbNZrPt9fw3r9Wk5wEgMy9tefrtiLgK+BHwEuCDC623ceNG1q5du9Oy8fFxxsfHe12iejhvvrTS\nNBoN1q8fZW7urkGXIrVlcnKSycnJnZbNzs523F/fpwLOzNmI+B6wbnftNm3axIYNG/pdjqB38+Zf\nD3ypJxVJS6bZbNYf/puB0TbXvhg4q/dFSYuwqy/FU1NTjI2NddRf3wNARDwUeCxwXr+3pTZ1O29+\n+0dQpWVkFGj3S4enALR69GMegHdFxLMj4tER8QzgIqrLACf3sKokSVoi/TgCcBjwUeBAYCtwBfD0\nzLy1D9uSJEkd6McgQEftSZK0zHkvAEmSCmQAkCSpQAYASZIKZACQJKlABgBJkgpkAJAkqUAGAEmS\nCmQAkCSpQAYASZIKZACQJKlABgBJkgpkAJAkqUAGAEmSCmQAkCSpQAYASZIKZACQJKlABgBJkgpk\nAJAkqUB7D7oAabVqNBo0m82u+xkeHmZkZKQHFUnSrxgApD5oNBqsf/x65u6e67qvoTVDzFw3YwiQ\n1FMGAKkPms1m9eF/EjDcTUcwd+EczWbTACCppwwAUj8NA4cOughJejAHAUqSVCADgCRJBTIASJJU\nIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUIGcClHZhenp6oOtLUr8ZAKRW24CAiYmJ\nQVciSX1lAJBazQFJ9zfxuR74Uk8qkqS+MABIu9LtTXyavSpEkvrDQYCSJBXIACBJUoEMAJIkFcgA\nIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoGcCEjSQGQmmzZt4uqrr+5o/cMPP5yTTz6ZiGhrPe/T\nsPx1+nc0PDzMyMhIj6tZvQwAkgbijjvu4PWvfz0RTyDi4LbW3b79JuCfecc73tGf4jQgNwN7dXwv\njqGh/ZiZmTYELJIBQNJAZb6ZzJe0uda7gf8ObAZG21z3YuCsNtfR0vgFsJ3O/l6nmZuboNlsGgAW\nyQAgaQUbBTa0uY6nAJa/Tv5e1S4HAUqSVCADgCRJBTIASJJUIAPASnbtoAtYgdxnHZmcnBx0CSuQ\n+0zLW98CQEScHhE3RMTdEfH1iPjP/dpWsfwwa5/7rCMGgE64z7S89SUARMRLqa7TeQvwFOCbwKUR\nMdyP7UmSpPb06wjARuAfMvO8zLwOOA24Czi1T9uTJElt6HkAiIh9gDHgizuWZWYCXwCO6vX2JElS\n+/oxEdAw8BBgy7zlW4D1u2g/BL2fn/u+e++DnwKdTTNed9Ly5+uBZof9NHrQx676uR341jKqZ5D9\nLLaPPe2z5fSeAH5e/Rj0/PWzs7NMTU31tM9t27bVf/oC1Qxw7fh6/fNi2p/Y52tLtO5PgI/0YLtL\nVe9yWLebbd4ADP7/ylJreb9D7a4b1Zfz3omIRwA3AUdl5pUty98JPCsznzGv/cvZ+X+JJElqzysy\n86PtrNCPIwBN4AFg/t09DuLBRwUALgVeAdwIzPWhHkmSVqsh4HCqz9K29PwIAEBEfB24MjPPqJ8H\n1UHR92Xmu3q+QUmS1JZ+3QzoPcCHI+Ia4CqqqwL2Az7Up+1JkqQ29CUAZOb59TX/b6M6FfAN4AWZ\nubUf25MkSe3pyykASZK0vHkvAEmSCmQAkCSpQAMLABFxdER8KiJuiojtEXHCoGpZKSLizIi4KiJu\nj4gtEXFRRDxu0HUtdxFxWkR8MyJm68e/RsTvDLqulaT+t7c9It4z6FqWs4h4S72fWh/fHXRdy11E\nHBoR/xQRzYi4q/7/umHQdS1n9c325v9b2x4RZy+2j0EeAdifanDg6YADERbnaOBs4EjgOGAf4HMR\nsWagVS1/PwbeQDVF9RhwGfDJiBgdaFUrRH0nz9dS3dRLe/ZtqsHPh9SPZw22nOUtIg6gmgLwHuAF\nwCjwen45D6YW8FR+9W/sEOB5VJ+l5y+2g35dBrhHmXkJcAn8cp4A7UFmHt/6PCJOAX5G9aF2xSBq\nWgky8zPzFr0pIl4HPJ325xstSkQ8FNgMvAY4a8DlrBT3e8VTW/4n0MjM17Qs+9GgilkpMvPW1ucR\n8ULgB5n51cX24RiAle0AqsR326ALWSkiYq+IeBnVvBT/Nuh6VoBzgE9n5mWDLmQF+c361OYPImJz\nRDxq0AUtcy8Ero6I8+tTm1MR8Zo9rqVfqm/C9wrgA+2sZwBYoeqjJu8FrshMzzHuQUQ8MSLuoDrM\neC5wYn2rai2gDkr/CThz0LWsIF8HTqE6lH0a8BjgKxGx/yCLWuaOAF4HzADPB/4eeF9ETAy0qpXl\nRGAt8OF2VhrYKQB17Vzgt4BnDrqQFeI64MlUR01OBs6LiGcbAnYtIg6jCpjPy8z79tRelcxsnY/9\n2xFxFdXh7JcAHxxMVcveXsBVmbnjFNM3I+IJVKFg8+DKWlFOBT6bmbe0s5JHAFagiHg/cDxwbGbe\nPOh6VoLMvD8zf5iZU5n5RqoBbWcMuq5lbAz4DeCaiLgvIu4DjgHOiIh7HbezOJk5C3wPWDfoWpax\nm3nwWJxpYGQAtaw4ETFCNSj8/7S7rkcAVpj6w/9FwDGZ2dhTey1oL2DfQRexjH0B+I/zln2I6hfz\nX6dTiC5KPYjyscB5g65lGfsasH7esvU4EHCxTqW60+7F7a44sABQnxNbB+z4JnFERDwZuC0zfzyo\nupaziDgXGAdOAO6MiB23XJ7NTG+lvICI+Evgs1SXAz6MarDMMVTnG7ULmXknsNPYkoi4E7g1M71y\nYgER8S7g01QfXo8E/gK4H5gcZF3L3CbgaxFxJtUlbEdSXXXy2oFWtQLUR+JOAT6UmdvbXX+QRwCe\nCnyJahR7Au+ul3+YKtHowU6j2leXz1v+avyGsTsHU+2fRwCzwLeA5zuyvW1+69+zw4CPAgcCW6ku\nz336/Eu29CuZeXVEnAj8NdWlpjcAZ2TmxwZb2YpwHPAoOhxf4s2AJEkqkIMAJUkqkAFAkqQCGQAk\nSSqQAUCSpAIZACRJKpABQJKkAhkAJEkqkAFAkqQCGQAkSSqQAUCSpAIZACRJKtD/B9nOYbtXfy0w\nAAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "metadata": { "engine": 0 }, "output_type": "display_data" }, { "data": { "text/plain": [ "[output:1]" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgAAAAFyCAYAAACDemKtAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAG/pJREFUeJzt3X2UZVV55/HvgyANqE1CC6jYKnRsK74wdjkqkVfBl2FG\nQM0YSmsMGjNimAyrk1njYghidJkEjDaDSkyWMwbSWA4aHWFEcHxBAwpoVxTUojXaegGh4XZrtdAU\nb/3MH+e0U11Ud9/Xvrdqfz9r3dVd5+6z71N3VdX5nXP23TsyE0mSVJa9Bl2AJEna8wwAkiQVyAAg\nSVKBDACSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACANUET8fkRsm/V4ICLW\nR8SHIuLgWe2Oq59/3W72vTMiromIP46IJ7RRx7kR8bmIuLvu6127aHtSRHwlIu6NiF9ExE0RMT5P\nu4Mj4uMRsTEitkbEuoj43Z30eXr9/AMRcU9EfCwiDtpNzUfXtT4aEb/Z6vcqqWIAkAYvgT8DxoGz\ngBuAdwDfiIglc9rtat8zgYvrbRcBt0bE81us4b3Ai4DJnbwOABFxCnAtsA9wPvDfgK3AZRFx9qx2\nT6y/j9cCfwP8KbAFuCIiTp/T5zuATwBNYDXwd8DpwJci4vE7qSPq7/W+Fr8/SXPsPegCJAFwTWZO\n1v//nxGxmepgeCrwv+rt0cK+ABdExPHA54HPRcRIZj64m9d/ZmY26rPue3fR7izg58AJmfkIQET8\nHXAbcAbw3+t2ZwKHAy/PzK/V7T4K3Ah8ICI+nZmPRMQ+wPuA6zLzVdtfJCK+CVwF/CHwkXnqeDtw\nGPAx4Ox5npe0G14BkIbTV6gO+M/qZOfMvI7qrP4ZVFcHdte+0WLXTwJ+sf3gX+/7KNXZ+wOz2h0N\n3Lv94F+3S+AK4FDguHrz84AD6+2z6/k81dn9DlcLACLiQKrv7TxgusW6Jc1hAJCG04r6301d9PEP\nVCHild2X82vXAc+NiPdExBERcXhEnAeMAhfOarcvOwaC7bbWNY3OasdO2j4AvHCe7e8D7qK6VSCp\nQ94CkIbD0vry+xKqs+fzqA6W/6fTDjPzzoiYBo7oTYkAvIfqqsS5VGMPAO4HXp+ZV81qtx44MSKe\nnpm3z9p+LNUYg6fVX/+o/vplwKXbG0XESuDJQEbEb2TmL+rtLwD+I/DqzMxqKICkTngFQBq8AL5M\nde/9dqoBcVuA0zLzri77vg94Ypd9zPYQ8EPgU1SX598EfBu4PCJePKvdx4BtwKci4qj6SsE5wGn1\n8/sBZOYmqsv/vx8RfxIRz4qIY4BP1q/167a1i4HPZ+aXe/g9SUXyCoA0eAn8EdXZ8CPAxsxc36O+\nnwBs7FFfUA3Ie3Fmrtq+ISI+BXyfagDgUQCZeWtEjAEfBa6nCjl3UQ3Y+yg7jt5/O9WVj/cDf031\nfqwFfkIVGO6rX+f3gJcCz+3h9yMVywAgDYdvzRnJ37WIeBqwFPiXHvW3D/BW4ILZ2+vR/F8AzoqI\nfTLz4Xr7ZyLiSuBI4HFUHzE8od7th7P23wK8NiIOA54J/Cwzb4+IG6gGEm6pm15IdeXhkYh4Rr3t\nN+p/l0fEvj24YiIVwwAgLV5vpjqbvqZH/R1E9TfjcfM8tw/VLcUdbivWnxZYt/3riHhFXdOX5naQ\nmXcAd9TtDqQaKPipWU2eDryR6rbDXJPAd4BV8zwnaR4GAGkRioiXUw3S+wnVmIJeuAf4JdXZ+rtm\nzQPwBOA1wNSu5huIiBVUl/uvyszdXZX4S6qgcdGsbafN024MeAPwH4A7W/1GJBkApGHQzVD2AE6O\niBGq3+dDgJcDrwA2AKdk5kO72L/qpJrK9xnAAfWm4yLi3Pr/l2Xm7Zm5LSL+muoz+DdFxGX1a/4B\n1aj+/zqnz+9TncE3qCYFOpNqvoB3zGn3Tqr5AG6iGgPxWuAk4NzM/PXVg8y8cp66t39M8JrM3Ly7\n71PS/2cAkAZvp1PvttAugT+v//8QsBm4FfjPwN9n5v0t9v0HVB/R297n8fUD4J+oPp1AZv5FRPyE\najDfu6g+x38L1ccA//ecPr9DNTvgIVQH/k8C787M5px2t1Kd3b+G6qz/FuDfZ+ZnWqxdUgeimpxL\nkiSVpK15ACLinIi4OSK21Ct8fTYinj2nzXVzVih7NCIu6W3ZkiSpG+1OBHQM8CHgJVT36PYBvhgR\nsyfqSKopOg+hmvP7Kcy5NyhJkgarrTEAmXny7K8j4gyqkcGjVJN9bLc1M3e1opgkSRqgbqcCPpDq\njH/u6Ns3RcS9EXFrRPzFnCsEkiRpwDoeBBjVKhxXAU/MzONmbX8b8DOqNcNfQDV7102Z+bs76ecg\n4FXAT4GZjoqRJKlMS6hm0Ly2XlujZd0EgL+hOnC/bFfTb0bECVSzfq3IzA3zPP9G4PKOipAkSQBv\nysy2Jv3qaB6AiPgwcDJwTAtzb99ENVnJCqqJSeb6KcDatWsZGRnppJxirV69mjVr1gy6jAXF96wz\nvm/t8z3rjO9be6amphgfH4f6WNqOtgNAffA/FTguMxst7PJCqnECOwsKMwAjIyOsWuU03u1YunSp\n71mbfM864/vWPt+zzvi+daztW+htBYD68/xjwCnA/RFxSP3UdGbORMThVIt1XA1soloF7IPA1zLz\ne+0WJ0mS+qPdKwBnUp3NXzdn+1uAy6imIj2JaprQA6imD/0U8L6uqpQkST3V7jwAu/zYYL2c5/Hd\nFCRJkvqv23kANEBjY2ODLmHB8T3rjO9b+3zPOuP7tucMfDGgiFgFrFu3bp0DPyRJasPk5CSjo6MA\no5k52c6+XgGQJKlABgBJkgpkAJAkqUAGAEmSCmQAkCSpQAYASZIKZACQJKlABgBJkgpkAJAkqUAG\nAEmSCmQAkCSpQAYASZIKZACQJKlABgBJkgpkAJAkqUAGAEmSCmQAkCSpQAYASZIKZACQJKlABgBJ\nkgpkAJAkqUAGAEmSCmQAkCSpQAYASZIKZACQJKlABgBJkgpkAJAkqUAGAEmSCmQAkCSpQAYASZIK\nZACQJKlABgBJkgpkAJAkqUAGAEmSCmQAkCSpQAYASZIKZACQJKlABgBJkgpkAJAkqUAGAEmSCmQA\nkCSpQHsPugBJ0sLQaDRoNpstt1+2bBnLly/vY0XqhgFAkrRbjUaDlStHmJnZ2vI+S5bsz/r1U4aA\nIWUAkCTtVrPZrA/+a4GRFvaYYmZmnGazaQAYUgYASVIbRoBVgy5CPeAgQEmSCmQAkCSpQAYASZIK\nZACQJKlAbQWAiDgnIm6OiC0RsTEiPhsRz57TZt+I+EhENCPiVxHx6Yg4uLdlS5KkbrR7BeAY4EPA\nS4CTgH2AL0bEfrPaXAT8W+D1wLHAU4F/7L5USZLUK219DDAzT579dUScAdwDjALXR8STgLcCp2fm\n1+o2bwGmIuLFmXlzT6qWJEld6XYMwIFAApvrr0epQsWXtzfIzPVAAziqy9eSJEk90nEAiIigutx/\nfWb+oN58KPBQZm6Z03xj/ZwkSRoC3cwEeAnw28DRLbQNqisFO7V69WqWLl26w7axsTHGxsY6LlCS\npMViYmKCiYmJHbZNT0933F9HASAiPgycDByTmT+f9dTdwOMj4klzrgIcTHUVYKfWrFnDqlVOLylJ\n0nzmOymenJxkdHS0o/7avgVQH/xPBU7IzMacp9cBjwAnzmr/bGA58M2OKpQkST3X1hWAiLgEGANO\nAe6PiEPqp6YzcyYzt0TE/wA+GBG/AH4FXAzc4CcAJEkaHu3eAjiT6l7+dXO2vwW4rP7/auBR4NPA\nvsA1wFmdlyhJknqt3XkAdnvLIDMfBP64fkiSpCHkWgCSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQV\nyAAgSVKBulkLQJK0QDUaDZrNZsvtp6am+liNBsEAIEmFaTQarFw5wszM1kGXogEyAEhSYZrNZn3w\nXwuMtLjX1cB5/StKe5wBQJKKNQK0ugqrtwAWGwcBSpJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCAD\ngCRJBTIASJJUIAOAJEkFciKgArU7B/iuLFu2jOXLl/ekL0mdcV5/dcIAUJhGo8HK56xk5oGZnvS3\nZL8lrL9tvSFAGhDn9VenDACFaTab1cH/dcCybjuDmc/M0Gw2DQDSgDivvzplACjVMuCpgy5CUu84\nr7/a4yBASZIKZACQJKlABgBJkgpkAJAkqUAGAEmSCmQAkCSpQAYASZIKZACQJKlABgBJkgrkTIDq\nWi8WFnFRIUnaswwA6tx9QMD4+HjXXbmokCTtWQYAdW4GSLpfWMhFhSRpjzMAqHsuLCRJC46DACVJ\nKpABQJKkAhkAJEkqkAFAkqQCGQAkSSqQAUCSpAIZACRJKpABQJKkAjkRkBadRqNBs9nsuh/XJ5C6\n1+5aIf7e7TkGAC0qjUaDlc9ZycwDM1335foEUjfuAvZqe62QJUv2Z/36KX/v9gADgBaVZrNZHfxd\nn0AasF8C24C1wEiL+0wxMzPu790eYgDQ4uT6BNKQGAFWDboIzcNBgJIkFcgAIElSgQwAkiQVyAAg\nSVKB2g4AEXFMRFwZEXdGxLaIOGXO8x+vt89+XN27kiVJUrc6uQJwAPAd4Cwgd9LmC8AhwKH1Y6yj\n6iRJUl+0/THAzLwGuAYgImInzR7MzHu7KUySJPVPv8YAHB8RGyPitoi4JCJ+s0+vI0mSOtCPiYC+\nAPwjsAE4AvhL4OqIOCozd3bLQFrUXJ9A0rDpeQDIzCtmffn9iLgV+DFwPPDVne23evVqli5dusO2\nsbExxsYcPqCFzfUJJPXCxMQEExMTO2ybnp7uuL++TwWcmRsiogmsYBcBYM2aNaxa5XSRWnxcn0BS\nL8x3Ujw5Ocno6GhH/fU9AETEYcBBVEtDSeVyfQJJQ6TtABARB1CdzW//BMDhEXEksLl+nE81BuDu\nut0FwA+Ba3tRsCRJ6l4nVwBeRHUpP+vHB+rtlwJ/BLwAeDNwIPBzqgP/uzLz4a6rlSRJPdHJPABf\nY9cfH3x15+VIkqQ9wbUAJEkqkAFAkqQCGQAkSSqQAUCSpAIZACRJKpABQJKkAhkAJEkqkAFAkqQC\nGQAkSSqQAUCSpAIZACRJKpABQJKkAhkAJEkqkAFAkqQCGQAkSSqQAUCSpAIZACRJKpABQJKkAu09\n6AKk7aampoaiD0kqgQFAg3cfEDA+Pj7oSiSpGAYADd4MkMDrgGVd9vUj4KtdVyRJi54BQMNjGfDU\nLvto9qIQSVr8HAQoSVKBDACSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACBJ\nUoGcCVDaBRcokrRYGQCk+bhAkaRFzgAgzccFiiQtcgYAaVdcoEjSIuUgQEmSCmQAkCSpQAYASZIK\nZACQJKlABgBJkgpkAJAkqUAGAEmSCmQAkCSpQAYASZIKZACQJKlABgBJkgpkAJAkqUAGAEmSCmQA\nkCSpQAYASZIKZACQJKlAbQeAiDgmIq6MiDsjYltEnDJPm/dExM8jYmtE/N+IWNGbciVJUi90cgXg\nAOA7wFlAzn0yIt4J/Cfg7cCLgfuBayPi8V3UKUmSemjvdnfIzGuAawAiIuZpcjbw3sy8qm7zZmAj\ncBpwReelSpKkXunpGICIeBZwKPDl7dsycwtwE3BUL19LkiR1ru0rALtxKNVtgY1ztm+sn5Mk7Uaj\n0aDZbLbUdmpqqs/VaLHqdQDYmWCe8QKzrV69mqVLl+6wbWxsjLGxsX7WJUlDpdFosHLlCDMzWwdd\niobMxMQEExMTO2ybnp7uuL9eB4C7qQ72h7DjVYCDgX/e1Y5r1qxh1apVPS5HkhaWZrNZH/zXAiMt\n7HE1cF5/i9JQmO+keHJyktHR0Y7662kAyMwNEXE3cCJwC0BEPAl4CfCRXr6WJC1uI0ArJ0XeAlBn\n2g4AEXEAsILqTB/g8Ig4EticmbcDFwF/FhH/AvwUeC9wB/C5nlQsSZK61skVgBcBX6W6p5/AB+rt\nlwJvzcwLI2J/4G+BA4F/Av5NZj7Ug3olSVIPdDIPwNfYzccHM/PdwLs7K0mSJPWbawFIklQgA4Ak\nSQUyAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUIAOAJEkF\nMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIA\nSJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiS\nVCADgCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklQg\nA4AkSQUyAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVKCeB4CIOD8its15/KDX\nryNJkjq3d5/6/R5wIhD114/06XUkSVIH+hUAHsnMe/vUtyRJ6lK/xgD8VkTcGRE/joi1EfH0Pr2O\nJEnqQD+uANwInAGsB54CvBv4ekQ8LzPv78PrSZIWkampqZbbLlu2jOXLl/exmsWr5wEgM6+d9eX3\nIuJm4GfAG4CP72y/1atXs3Tp0h22jY2NMTY21usSJUlD6S5gL8bHx1veY8mS/Vm/fqqIEDAxMcHE\nxMQO26anpzvur19jAH4tM6cj4ofAil21W7NmDatWrep3OZKkofVLYBuwFhhpof0UMzPjNJvNIgLA\nfCfFk5OTjI6OdtRf3wNARDwBOAK4rN+vJUlaDEYATwj7rR/zALw/Io6NiGdExO8An6X6GODEbnaV\nJEl7SD+uABwGfAI4CLgXuB54aWZu6sNrSZKkDvRjEKCj9iRJGnKuBSBJUoEMAJIkFcgAIElSgQwA\nkiQVyAAgSVKBDACSJBWo7zMBSlLJGo0GzWaz5fbtLIQjdcMAIEl90mg0WLlyhJmZrYMuRXoMA4Ak\n9Umz2awP/q0ubgNwNXBe/4qSagYASeq7dha38RaA9gwHAUqSVCADgCRJBTIASJJUIAOAJEkFMgBI\nklQgA4AkSQUyAEiSVCADgCRJBXIioAXiW9/6Fl//+te77ueOO+7oQTUapF7NFb9s2TKWL1/edT/t\nznW/K72qSWVp93fiwQcfZN99921rn8X4s2kAWCBO/ncns2nzJvbap7uLNtse3tajirTH3QcEjI+P\n96S7JfstYf1t67v6o9ZoNFj5nJXMPDAzNDWpJHcBe3XwO/E44NG29liyZH/Wr59aVD+bBoAF4sEH\nHySPTx49ur0f2se4FvhmT0rSnjYDJPA6YFmXfTVh5jMzNJvNrv6gNZvN6uA/RDWpJL8EttHZWgvt\n7DPFzMz4ovvZNABIC80y4KmDLmKOYaxJBelkrYV29lmcHAQoSVKBDACSJBXIACBJUoEMAJIkFcgA\nIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACBJ\nUoEMAJIkFcgAIElSgQwAkiQVyAAgSVKBDACSJBXIACBJUoH2HnQBkrQQZCYXXHABt9xyS8v7bN68\nuY8VSd0xAEhSC+655x7OOeccIp5PxJNb2mfbttv6XJXUOQOAJLUh831kvqbF1mcBl/SzHKljjgGQ\nJKlABgBJkgpkAJAkqUAGAEmSCmQAWMhuHXQBC5DvWUcmJiYGXcIC5Hum4da3ABARZ0XEhoh4ICJu\njIh/3a/XKpYHs/b5nnXEANAJ3zMNt74EgIj4PeADwPnAC4HvAtdGxLJ+vJ4kSWpPv64ArAb+NjMv\ny8zbgDOBrcBb+/R6kiSpDT0PABGxDzAKfHn7tsxM4EvAUb1+PUmS1L5+zAS4DHgcsHHO9o3Aynna\nLwGYmprqQymLx6OPPgp3At+etXF6ztetuLv+90dAs8uiGj3qq1f9tNLXFqCVqdz3ZE17uh+AX1T/\nXH311S397t1xxx1cfvnlj9m+YcOGntc0rH8LNm3aVP/vGuCuFvb4GdUvKcDVQKvf1w1t7tNu+z21\nz2J5DYDq53wYfzZn1bSk3X2jOjnvnYh4CtWh6qjMvGnW9guBozPzd+a0fyPw2L8skiSpVW/KzE+0\ns0M/rgA0gUeBQ+ZsP5jHXhUAuBZ4E/BTYKYP9UiStFgtAZ5JdSxtS8+vAABExI3ATZl5dv11UF3A\nvDgz39/zF5QkSW3p12qAHwQujYh1wM1UnwrYH/j7Pr2eJElqQ18CQGZeUX/m/z1UtwK+A7wqM+/t\nx+tJkqT29OUWgCRJGm6uBSBJUoEMAJIkFWjgAcBFg9oTEcdExJURcWdEbIuIUwZd07CLiHMi4uaI\n2BIRGyPisxHx7EHXNcwi4syI+G5ETNePb0TEqwdd10JS/9xti4gPDrqWYRYR59fv0+zHDwZd10IQ\nEU+NiH+IiGZEbK1/Z1e1uv9AA4CLBnXkAKpBlWcBDuBozTHAh4CXACcB+wBfjIj9BlrVcLsdeCfV\ntN6jwFeAz0XEyECrWiDqE5k/pPqbpt37HtWA8UPrx9GDLWf4RcSBVNMaPgi8ChgB/pRfz6nZQh+D\nHAS4k/kCbqeaL+DCgRW2QETENuC0zLxy0LUsJHXAvAc4NjOvH3Q9C0VEbAL+S2Z+fNC1DLOIeAKw\nDngHcB7wz5n5J4OtanhFxPnAqZnZ8pmrICL+imrG3eM67WNgVwBcNEgDdCDV1ZPNgy5kIYiIvSLi\ndKq5PL456HoWgI8AV2XmVwZdyALyW/VtzR9HxNqIePqgC1oAXgN8OyKuqG9tTkbE29rpYJC3AHa1\naNChe74claC+ynQRcH1mep9xFyLieRHxK6pLjJcAr62X99ZO1EHpXwHnDLqWBeRG4Ayqy9hnAs8C\nvh4RBwyyqAXgcKqrTOuBVwIfBS6OiPFWO+jXTIDdCLy3rf65BPht4GWDLmQBuA04kuqKyeuByyLi\nWEPA/CLiMKpw+YrMfHjQ9SwUmTl7DvvvRcTNVEspvgHwdtPO7QXcnJnn1V9/NyKeSxUK1rbawaC0\nu2iQ1JWI+DBwMnB8ZraynmvRMvORzPxJZk5m5rlUA9rOHnRdQ2wUeDKwLiIejoiHgeOAsyPiofrq\nk3YjM6eBHwIrBl3LkLuLx65nPAUsb7WDgQWAOiGvA07cvq3+BTkR+Mag6tLiVB/8TwVOyMzGoOtZ\noPYC9h10EUPsS8DzqW4BHFk/vk11NnZkOu1qS+pBlEdQHeC0czcAK+dsW0l19aQlg74F4KJBbarv\ni62gulUCcHhEHAlszszbB1fZ8IqIS4Ax4BTg/ojYftVpOjNdgnoeEfE+4AtUn8p5ItWS3cdR3WvU\nPDLzfmCHcSURcT+wKTPnnqmpFhHvB66iOnA9Dfhz4BFgYpB1LQBrgBsi4hzgCqqPOb+N6uOnLRlo\nAHDRoI68CPgq1TiJpJpHAeBS4K2DKmrInUn1Xl03Z/tbgMv2eDULwyFU781TgGngFuCVjmxvm2f9\nu3cY8AngIOBe4HrgpZm5aaBVDbnM/HZEvBb4K6qPm24Azs7MT7bah4sBSZJUoIFPBSxJkvY8A4Ak\nSQUyAEiSVCADgCRJBTIASJJUIAOAJEkFMgBIklQgA4AkSQUyAEiSVCADgCRJBTIASJJUoP8HW7gI\nXC1lBQEAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": { "engine": 1 }, "output_type": "display_data" } ], "source": [ "%%px \n", "import matplotlib.pyplot as plt\n", "import os\n", "tsamples = numpy.transpose(samples)\n", "plt.hist(tsamples[0])\n", "plt.hist(tsamples[1])\n", "_ = plt.title('PID %i' % os.getpid())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For profiling, we can also use the `%timeit` magic to compare performance on the engines:" ] }, { "cell_type": "code", "execution_count": 78, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[stdout:0] 1 loop, best of 3: 809 ms per loop\n", "[stdout:1] 1 loop, best of 3: 646 ms per loop\n" ] } ], "source": [ "%%px\n", "%%timeit\n", "s = abc(y, 10)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Cython Parallel\n", "\n", "In order to use Cython in parallel on IPython, we need to load and execute `cythonmagic` on all engines." ] }, { "cell_type": "code", "execution_count": 79, "metadata": {}, "outputs": [], "source": [ "%px %load_ext cython" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As an example, let's use the Cythonized Gibbs sampler from the last lecture." ] }, { "cell_type": "code", "execution_count": 80, "metadata": {}, "outputs": [], "source": [ "%%px\n", "%%cython\n", "\n", "from numpy import zeros, random\n", "from numpy cimport *\n", "from libc.math cimport sqrt\n", "gamma = random.gamma\n", "normal = random.normal\n", "\n", "def gibbs(int N=20000, int thin=200):\n", " cdef: \n", " ndarray[float64_t, ndim=2] mat = zeros((N,2))\n", " float64_t x,y = 0\n", " int i,j\n", " for i in range(N):\n", " for j in range(thin):\n", " x = gamma(3, y**2 + 4)\n", " y = normal(1./(x+1), 1./sqrt(2*(x+1)))\n", " mat[i] = x,y\n", "\n", " return mat" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Divide the array by the number of nodes. In this case they divide evenly, a more general partitioning of sizes is easy to do as well." ] }, { "cell_type": "code", "execution_count": 81, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[0, 1]" ] }, "execution_count": 81, "metadata": {}, "output_type": "execute_result" } ], "source": [ "cli.ids" ] }, { "cell_type": "code", "execution_count": 82, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[50000.0, 50000.0]" ] }, "execution_count": 82, "metadata": {}, "output_type": "execute_result" } ], "source": [ "N = 100000\n", "thin = 10\n", "n = N/len(cli.ids)\n", "dv = cli[:]\n", "dv.push(dict(n=n, thin=thin))\n", "# Let's just confirm visually we got what we expect\n", "dv['n']" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can time the execution of the gibbs sampler on the remote nodes" ] }, { "cell_type": "code", "execution_count": 83, "metadata": {}, "outputs": [], "source": [ "%pxconfig --noblock" ] }, { "cell_type": "code", "execution_count": 84, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "100 loops, best of 3: 2.12 ms per loop\n" ] } ], "source": [ "%timeit dv.execute('gibbs(n, thin)')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "But a more realistic (and costly) benchmark must also include the cost of bringing the results back from the cluster engines to our local namespace. For that, we assign the call to the variable `a` on each node and then use the view's `gather` method to pull them back in:" ] }, { "cell_type": "code", "execution_count": 85, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The slowest run took 4.56 times longer than the fastest. This could mean that an intermediate result is being cached.\n", "100 loops, best of 3: 6.19 ms per loop\n" ] } ], "source": [ "%%timeit\n", "dv.execute('a = gibbs(n, thin)')\n", "a = dv.gather('a')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's compare this to the same number of samples executed on a single process:" ] }, { "cell_type": "code", "execution_count": 86, "metadata": {}, "outputs": [], "source": [ "%load_ext cython" ] }, { "cell_type": "code", "execution_count": 87, "metadata": {}, "outputs": [], "source": [ "%%cython\n", "\n", "from libc.math cimport sqrt\n", "from numpy import zeros, random\n", "from numpy cimport *\n", "gamma = random.gamma\n", "normal = random.normal\n", "\n", "def gibbs(int N=20000, int thin=200):\n", " cdef: \n", " ndarray[float64_t, ndim=2] mat = zeros((N,2))\n", " float64_t x,y = 0\n", " int i,j\n", " for i in range(N):\n", " for j in range(thin):\n", " x = gamma(3, y**2 + 4)\n", " y = normal(1./(x+1), 1./sqrt(2*(x+1)))\n", " mat[i] = x,y\n", "\n", " return mat" ] }, { "cell_type": "code", "execution_count": 88, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1 loop, best of 3: 1.07 s per loop\n" ] } ], "source": [ "%%timeit \n", "a = gibbs(N, thin)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Exercise\n", "\n", "Run parallel chains of the `disaster_model` example from PyMC and return the resulting traces to your client, for plotting and summarization." ] }, { "cell_type": "code", "execution_count": 89, "metadata": {}, "outputs": [], "source": [ "# Write your answer here" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## References\n", "\n", "[Scientific Python Lectures](http://github.com/jrjohansson/scientific-python-lectures) by Robert Johansson\n", "\n", "[Using IPython for Parallel Computing](http://ipython.org/ipython-doc/dev/parallel/)" ] } ], "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.5.2" }, "latex_envs": { "bibliofile": "biblio.bib", "cite_by": "apalike", "current_citInitial": 1, "eqLabelWithNumbers": true, "eqNumInitial": 0 } }, "nbformat": 4, "nbformat_minor": 1 }