Mentions légales du service

Skip to content
Snippets Groups Projects
Commit 26dfb7d1 authored by hhakim's avatar hhakim
Browse files

Update/Review Faust optimizations jupyter notebook.

parent 6065c635
No related merge requests found
......@@ -13980,7 +13980,7 @@ Indeed, several operations are available in the FAµST API to optimize a Faust o
<h2 id="1.-Faust.optimize_memory">1. Faust.optimize_memory<a class="anchor-link" href="#1.-Faust.optimize_memory">&#182;</a></h2><p>As you probably already know, a Faust can be composed of dense or sparse factors. In the FAµST C++ core, a dense matrix is encoded as a column-major order array in memory while a sparse matrix follows the <a href="https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.csr_matrix.html">CSR format</a>. Although this is not memory efficient, it is totally possible to encode a sparse matrix in the dense format and conversely a dense matrix in the CSR format. There is no reason to do this but in many situations you might end up with a Faust encoded in a non-optimal format in your code.<br>
That's where the <code>Faust.optimize_memory</code> function comes up, it takes a Faust and gives you back a new one for which the format of each factor is ensured to be memory efficient. It means that depending on the number of nonzeros of the considered matrix, the format will be chosen by the function to minimize the memory used for the whole Faust.</p>
<p>Let's show this principle with a very simple example:</p>
<p>I build a Faust creating my own CSR matrix using <a href="https://www.scipy.org/">scipy</a> and a dense matrix using a <a href="https://numpy.org/">numpy</a> array. On purpose I choose to create a CSR matrix of size NxN whose nnz is equal to N^2 (something that obviously must be encoded as a numpy array) and likewise I create a numpy array of size NxN whose nnz is 1 (here again the format is not well chosen because a CSR matrix would be far lighter in memory).</p>
<p>we build a Faust creating my own CSR matrix using <a href="https://www.scipy.org/">scipy</a> and a dense matrix using a <a href="https://numpy.org/">numpy</a> array. On purpose we choose to create a CSR matrix of size NxN whose nnz is equal to N^2 (something that obviously must be encoded as a numpy array) and likewise we create a numpy array of size NxN whose nnz is 1 (here again the format is not well chosen because a CSR matrix would be far lighter in memory).</p>
 
</div>
</div><div class="jp-Cell jp-CodeCell jp-Notebook-cell jp-mod-noOutputs ">
......@@ -14465,10 +14465,10 @@ The following figure illustrates a basic iteration of pruneout on a product AB.
</div>
<div class="jp-Cell-inputWrapper"><div class="jp-InputPrompt jp-InputArea-prompt">
</div><div class="jp-RenderedHTMLCommon jp-RenderedMarkdown jp-MarkdownOutput " data-mime-type="text/markdown">
<p>The example above is trivial but when a Faust contain many factors including very sparse matrices, the benefit can be quite significant.
<p>The example above is trivial but when a Faust contains many factors including very sparse matrices, the benefit can be quite significant.
The reason why is that <code>pruneout</code> repeats the removal operation iteratively until there is no zero column/row eligible to this operation anymore.</p>
<p>Note that of course the resulting Faust gives the same full matrix (i.e. <code>Faust.toarray</code> must give the same result after a <code>Faust.pruneout</code>).</p>
<p>In the next section of code I show with a figure how a sliced DFT Faust can benefit of <code>Faust.pruneout</code> in terms of computation time.</p>
<p>In the next section of code we show with a figure how a sliced DFT Faust can benefit of <code>Faust.pruneout</code> in terms of computation time.</p>
 
</div>
</div><div class="jp-Cell jp-CodeCell jp-Notebook-cell ">
......
%% Cell type:markdown id:f4a6da37 tags:
# Optimizing a Faust
When a Faust object is instantiated it might be possible to optimize its structure or configuration before using it elsewhere. It applies in all cases, whether you gave explicitly the matrices to define the Faust or you obtained it through a FAµST algorithm or an operation like Faust-slicing.
Indeed, several operations are available in the FAµST API to optimize a Faust object, in memory space or in computation time.
In this notebook, we'll present three functions that aim at fulfilling these tasks:
[1. Faust.optimize_memory](#1.-Faust.optimize_memory)
[2. Faust.optimize_time](#2.-Faust.optimize_time)
[3. Faust.pruneout](#3.-Faust.pruneout)
%% Cell type:markdown id:9b61e964 tags:
## 1. Faust.optimize_memory
As you probably already know, a Faust can be composed of dense or sparse factors. In the FAµST C++ core, a dense matrix is encoded as a column-major order array in memory while a sparse matrix follows the [CSR format](https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.csr_matrix.html). Although this is not memory efficient, it is totally possible to encode a sparse matrix in the dense format and conversely a dense matrix in the CSR format. There is no reason to do this but in many situations you might end up with a Faust encoded in a non-optimal format in your code.
That's where the ``Faust.optimize_memory`` function comes up, it takes a Faust and gives you back a new one for which the format of each factor is ensured to be memory efficient. It means that depending on the number of nonzeros of the considered matrix, the format will be chosen by the function to minimize the memory used for the whole Faust.
Let's show this principle with a very simple example:
I build a Faust creating my own CSR matrix using [scipy](https://www.scipy.org/) and a dense matrix using a [numpy](https://numpy.org/) array. On purpose I choose to create a CSR matrix of size NxN whose nnz is equal to N^2 (something that obviously must be encoded as a numpy array) and likewise I create a numpy array of size NxN whose nnz is 1 (here again the format is not well chosen because a CSR matrix would be far lighter in memory).
we build a Faust creating my own CSR matrix using [scipy](https://www.scipy.org/) and a dense matrix using a [numpy](https://numpy.org/) array. On purpose we choose to create a CSR matrix of size NxN whose nnz is equal to N^2 (something that obviously must be encoded as a numpy array) and likewise we create a numpy array of size NxN whose nnz is 1 (here again the format is not well chosen because a CSR matrix would be far lighter in memory).
%% Cell type:code id:27b0be31 tags:
``` python
import numpy as np
import pyfaust as pf
from scipy.sparse import random
N = 64
S = random(N, N, 1, format='csr') # density = 1 <=> nnz = N^2
M = np.zeros((N,N))
M[0,0] = 1
# now create the Faust
F = pf.Faust([S, M]) # S is in CSR format, M is in column-major order array
```
%% Cell type:markdown id:56689697 tags:
Let's check the size (in bytes) of the Faust F we've just created.
%% Cell type:code id:2570500a tags:
``` python
F.nbytes
```
%% Cell type:markdown id:e7066f2b tags:
Now use the ``optimize_memory`` function to get a memory-efficient Faust equal to F and check its new size.
%% Cell type:code id:17d4fd86 tags:
``` python
G = F.optimize_memory()
G.nbytes
```
%% Cell type:markdown id:918552e4 tags:
Let us verify that F and G are really equal.
%% Cell type:code id:f195506b tags:
``` python
np.allclose(F.toarray(), G.toarray())
```
%% Cell type:markdown id:d58d51d4 tags:
As you can see the memory consumed by F is uselessly greater than the memory used for G. What ``optimize_memory`` did is simply to use a sparse CSR format for M which is formally a sparse matrix (it contains a few nonzeros for a lot of zeros), and convert S to an array because a dense matrix consumes more in memory if it is encoded as a sparse CSR matrix. The change of format can be verified printing the Faust-s:
%% Cell type:code id:5e825711 tags:
``` python
F
```
%% Cell type:code id:82ee7a15 tags:
``` python
G
```
%% Cell type:markdown id:7942efa4 tags:
You'll find the API documentation of ``Faust.optimize_memory`` [here](https://faustgrp.gitlabpages.inria.fr/faust/last-doc/html/classpyfaust_1_1Faust.html#a86e18a39e5b6a1091e232f3249a20191).
%% Cell type:markdown id:950a24f2 tags:
## 2. Faust.optimize_time
%% Cell type:markdown id:c9d2dcd2 tags:
The goal of ``Faust.optimize_time`` is not to be memory efficient like ``optimize_storage`` does but rather time efficient. The problem to consider here is the [matrix chain problem](https://en.wikipedia.org/wiki/Matrix_chain_multiplication). If you want to multiply an arbitrary number n of matrices, because of the associative property of the matrix product, you have many ways to do it. In fact, as many ways as number of possible parenthesis pair positions in the product. This number of possible positions is given by the [catalan number](https://en.wikipedia.org/wiki/Catalan_number) $C_{n-1}$. From an optimization point of view the question is what set of positions would be the most efficient to compute the matrix product? Or equivalently which order of matrix products will allow to compute the whole product of n matrices in the minimum time?
That is the problem ``optimize_time`` tackles with several solutions:
- The ``DYNPROG`` method which is the [classic dynamic programming approach](https://en.wikipedia.org/wiki/Matrix_chain_multiplication#A_dynamic_programming_algorithm) to solve the problem. A good reference to look at is:
**Cormen, Thomas H; Leiserson, Charles E; Rivest, Ronald L; Stein, Clifford (2001). "15.2: Matrix-chain multiplication". Introduction to Algorithms. Second Edition. MIT Press and McGraw-Hill. pp. 331–338.**
Note that neverthless the method implemented in FAµST is an extended version of the one discussed in the given reference because FAµST takes also into account the complexity of the sparse matrices that have to be computed differently from the complexity of a dense matrix product.
- FAµST provides other methods in exprimental packages, basically these methods use the Torch library as a backend to compute the product with some variations.
Anyway, ``optimize_time`` is able to select the best method available for your Faust. It's capable to benchmark these methods automatically on ``Faust.toarray`` or on a Faust-matrix product. Note that the GPU backend is not yet included in the benchmark.
That's what the following examples show:
%% Cell type:code id:2a65e466 tags:
``` python
import pyfaust as pf
from pyfaust.demo import get_data_dirpath
from os.path import join
from scipy.io import loadmat
from scipy.sparse import csr_matrix
F = pf.Faust(join(get_data_dirpath(), 'F_DYNPROG.mat'))
M = csr_matrix(loadmat(join(get_data_dirpath(), 'M_DYNPROG.mat'))['M'])
G = F.optimize_time(mat=M)
```
%% Cell type:markdown id:a41520b8 tags:
The code above is creating a Faust ``F`` based on the file ``F_DYNPROG.mat`` (provided by the FAµST data archive). Likewise a matrix ``M`` based on the file ``M_DYNPROG.mat`` is created. As you can guess from the file names the best method found by ``optimize_time`` is the ``DYNPROG`` method. Hence ``G`` is a version of the Faust ``F`` for which the ``DYNPROG`` method will be used to compute ``G.toarray()`` or ``G@M`` (for any matrix ``M``). By default, ``F`` is configured to compute the both of these products by simply multiplying the matrices from the right to the left (aka ``DEFAULT_L2R``).
Note that M is passed as ``mat`` argument of ``optimize_time`` which in consequence runs the benchmark on ``F@M`` (without ``mat`` argument, ``optimize_time`` runs the benchmark on ``F.toarray()``).
Now let us compare the computation times respectively for F and G.
%% Cell type:code id:eb43f5b5 tags:
``` python
timeit F.toarray()
```
%% Cell type:code id:ae123712 tags:
``` python
timeit G.toarray()
```
%% Cell type:code id:ee46450a tags:
``` python
timeit F@M
```
%% Cell type:code id:10d41d81 tags:
``` python
timeit G@M
```
%% Cell type:markdown id:0c7e8834 tags:
G is definitely faster than F!
%% Cell type:markdown id:064153b1 tags:
For further details about ``Faust.optimize_time``, you can look at the API reference [here](https://faustgrp.gitlabpages.inria.fr/faust/last-doc/html/classpyfaust_1_1Faust.html#a0348f4cfe22b920fa741860773240174) and maybe the [FaustMulMode](https://faustgrp.gitlabpages.inria.fr/faust/last-doc/html/classpyfaust_1_1FaustMulMode.html) too for a tour of all methods available.
%% Cell type:markdown id:87433857 tags:
## 3. Faust.pruneout
As ``Faust.optimize_time``, the member function ``Faust.pruneout`` is mainly dedicated to optimize the Faust product time (for ``Faust.toarray()`` or Faust-matrix product). Basically, this function removes zero rows (resp. columns) and corresponding columns (resp. rows) in consecutive pairs of factors.
The following figure illustrates a basic iteration of pruneout on a product AB. The submatrices covered by the blue grid are removed by the pruneout operation because the row block of B is composed only of zeros. Hence the resulting product A'B' is less costly.
%% Cell type:markdown id:b8d84ce6 tags:
<img align="center" src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAc4AAAHBCAYAAADgsFtlAAABhWlDQ1BJQ0MgcHJvZmlsZQAAKJF9kT1Iw0AcxV9TpVoqHewgIpKhOlkQFdFNq1CECqFWaNXB5NIPoUlD0uLiKLgWHPxYrDq4OOvq4CoIgh8gbm5Oii5S4v+SQosYD4778e7e4+4dINRLTLM6RgFNr5ipRFzMZFfEwCu6EUQYg5iWmWXMSlISnuPrHj6+3sV4lve5P0ePmrMY4BOJZ5hhVojXiSc3KwbnfeIIK8oq8TnxiEkXJH7kuuLyG+eCwwLPjJjp1BxxhFgstLHSxqxoasQTxFFV0ylfyLisct7irJWqrHlP/sJQTl9e4jrNASSwgEVIEKGgig2UUEGMVp0UCynaj3v4+x2/RC6FXBtg5JhHGRpkxw/+B7+7tfLjY25SKA50vtj2xxAQ2AUaNdv+Prbtxgngfwau9Ja/XAemPkmvtbToERDeBi6uW5qyB1zuAH1PhmzKjuSnKeTzwPsZfVMW6L0Fgqtub819nD4AaeoqeQMcHALDBcpe83h3V3tv/55p9vcDmrdyt+LwhsoAAAAGYktHRAD/AJcAL2Cq798AAAAJcEhZcwAALiMAAC4jAXilP3YAAAAHdElNRQflCQoPHzgG9FO5AAAgAElEQVR42uydeXhU1fnHv9kh+2QBAiaEYZFdcCQWRIQ6EQVcAOOPVWvVgJWqta0JFJcWiwlFapVWCbTUBbCEpQqkYkZZDTIQ2RISJRlCQgiQZCaZZEIm2/39Ee9w586ZzD6ZCe/neXwk98y98825N3PmfM/7vseH4zgOBEEQBEFYhS91AUEQBEHQwEkQBEEQNHASBEEQBA2cBEEQBEEDJ0EQBEHQwEkQBEEQNHBSFxAEQRAEDZwEQRAEQQMnQRAEQdDASRAEQRA0cBIEQRAEDZyEm8jOzoaPj4/Rf1lZWdQxBEEQNHAS5gZOMQqFgjqGIAjCA/Ch3VE8C41Gg6ioKGabWq2GRCKhTiIIgqAZJ8GabaakpCAlJcXwM9m1hCeTnp5ussTA+i85ORnp6elQqVTUaQQNnITjCAdHuVwOuVzOHFQJwltRKBTIzMzE4MGD6csg4ZWQVetBqFQqDB482PCzWq0GACPrtrS0FFKplDqL8MgZZ2ZmptEXP/HzzZpl0jNN0IyTcMpsMyUlBRKJBBKJhOxawivJzc01+q+0tBSlpaVITU01+9wTBA2chE0IrVjht3Wya4meglQqRUZGBnUEQQMn4Tj5+flGNpZwlin8t0qlotQUokdBkeIEDZyEw7NN3qYVfrAIB0+adRLe/AVxyZIlZp9tgvAGKDjIQ4iKioJGo7H6GzofOEQQnoI4OMgSUqkU27dvh0wmo84jaMZJ2D7btHbQBDqLJNCsk/BmpFIp0tLSKJqW8Er8qQu6H+GapUQiMfsNPD8/3zDAZmdnk8VFeDTidBThM6xSqbBkyRIsWbIEGzZsMIm0JQhPhqzabkZcYi8tLc1s1GFmZibS09MNP1MJPsKTEFu15j5aVCoVMjMzjdJQKJeT8CbIqu1mxJZrV7NIcRvZtYQ3wkpJoVxOggZOwq6BUyaTdRkoIZVKjdrpw4YgCIIGzlsKcU6mNWuWwteIcz8JwhvQaDRGSw78l0KC8BYoOKgbkUgkyM3NNZpxWiI1NdXodbTGSXgqycnJzOPiAh6Uy0nQwEnYNHCyIg+dfQ5BdAfWVLjivzzSF0CCBk6CIIguBkuZTAa5XI7U1FQaNAmvg9JRCIIgCMIGKDiIIAiCIGjgJAiCIAgaOAmCIAiCBk6CIAiCoIGTIAiCIGjgJKxB++1R6MvLvVJHxZoMaL89SjeRIAjCmwfOtvp6VKzJgFI60Cs+1IsXL8T1z7Z6pY6qDz9APQ2cBEEQXeLRBRB0BedQtmI5ghISvKZDk1SXSAdBEATNOLsHdc4+RM2YiSHr/+E1HVq8eCGub9vq0TqKFy9ExZoMq2f76px93f77VKzJgL68HNpvj6JiTYbH6CIIgmacHkX8q+le16Hab48iZMyYHqGj7A/LoS8vR9yS57v996n68AOjL1Pt9fUoWfYrxC193iufE4IgaOAkehj8oHn7J1vgHxHhMQ7EHQePGH72i4iAOmcfDZwEQdDASTif4sULDf9uKjgHfUU5dOfOAQBCxowxGnzK/rAcunPnPGrQBICoGTNNjnlCBDNBEDRwEj0QoW2rryhHUHyC4VhQ/M3gK+3Ro9AVnEOf+Qs8atAkCIKggZNwK8IZpe7cOZNZppC4pc+j6sMPEH7PZOYsjyAIggZOgviJ8MmTDQNq2R+WI2T0GK9KByIIgnA1Hp2Oov32qOE/fqYk/Nkc17dtNVrTs/Rzd2KNNmfrDRkzxsieNTdDDYpPQMmyX3n870MQBEEzzp8Qf7gKcw+7SvDXV5QbDa6Wfu5OrNHmbL3WRqEOWf8PFD46CxVrMqw+pzt+H4IgCHfiw3EcR93gPJTSgR6RW+gpOgiCIHoaVOTdifCzqIh7JpMOgiAIGjiJrqhYk4HixQsRMnoMgkePueV1EARB9FTIqnUSbfX1aK+v75YIVG1TC4ovqTt/aGzAjVo1eg8caPf1fqioxe3x0UbHwoIDMWJgFN1ogiBo4HT3wLn7cAkGxIZ6ZWeVX2vAjkPFiO8TgUB/48l6S1sHVFc0CO4VgIiQIESEBJmcX1xei5jIYDQ1tyKhT7jp9a9r0dHBobmlDcMTok3a63V6VNU2orW9HbfHxxhpuD0hCs0tbThythxBgX64/bZo9I0y7me+PSq8N/z9fHHH4L4m73H47CXEhIegRqvDlLE3B98fytV4dPIQj703ldWNmD3FffqURVfp04NwGUkj+lEneDKcmzl+vorzVo6fr+IOfF/OvfHvw1x9Y7PheH1js9Gxd7Yf5wouVhudKzy2/cB5bvuB80btwmMFF6u5d7YfN2oXHhO/nzUaxO0sDcJzxBre3/m9x98beo6JngA9W56P22ecn35xGP2CW7zyS0bJtU7d/SX+2HNKjYfHd1qX/L+DA2/OAPecUmPcwBDERwUZ/Zsnr0QLAJg0JNzo3zwVaj1OX9Lh4fFRRv/maWrpsFpDdGgAs134viyNwvf98mwjHhzruU7B1aZALHpkitve72vlBXBa2vOUcIENGD4Q9ycNpY4gq9bY4vJWG4K355JG9INWp0fmtmPgfID0eRMRzrBm12UrUVPfhIXy0RiVGGPSnn2wCPk/XoVsWD+kTB1h0l5YVoMtigLERATjlZQkk3ZrNGRszUOVWodVT09htmcfLEJufhlemjuBqbGwrAbv7VQiwC8Q61+We/S9cedz5c3PMeH5nzP0bHk2FFXriM3tA/hwnqtBq9OjubUdcVEhqKhuYL6moroBoxJjcL6smtl+vqwaIxNj0dDUQjecIAiiO2acZNV24mqrVtgeHOjL1CA8xtIgPEZWrTFk1RIu+1Amq9YLZk208N3jgoNY7WINLI1CDWKNFBxEARwEPVtENwUHefsa58dfncPqZ03XC7U6PdbtUCIowB+PTBrKXC9cl62En68P+keHMtc0sw8W4UptI9o7OOaaZmFZDfbk/Yjm1na88ngSU4Or1l3X7zqFZXPGe/S9oTVOoidgzbPVVl+Pqg0foOrDDzD8ky0IpyphbsXtRd5/vPAjtJUFXtlZhZXNuHZdi4OHDhnZsjzXq2qhaWrDqMgGVJWYtjfWanG2QoeZd0RBoag0bVfrceS0BmMTgqFQaE3am1o6cKa4GpIQfyiPHWG2X6mqBQCzGi+VqaFtbodSqUNVienAWlSiRdX1ZhR1aKBou6mxrKwRCkWtR1u17hzIGhoaoFAU0CcI4RKrFjD/LOsKzqFsxXLa7o+sWu9h7WcnurRBzdmkQuuTZZMKrViWlSu8LsvKJauWrFri1rBqyzPf5q588I/O1w5K4OqPHqFOczMUVWsjQQF+eOXxJKzboYRWpzfYn9MnSDEqMQbhIUEm7dkHiwDAYH2+kpKE/SdUKCyrMViw+0+oDPYs/zr+PN4G5u3ZUYkxmD5BinXZSmY7S4NQI0uDWKNYA0EQnkH8q+mIW/o8dUR3ugIUVWsbfHQpH7Ua4OeDCdJQo2hV3jbdc0qNvhEBCPT3MYpW5dlzSo24yEBU1bUYRczy5JVo0dLG4Vp9q0nULtAZeXviYiNa2zhmuyWNvAZfX0AS4s/UmFeihUbXhqt1HJ65L9qjrVqKqiV6ilVrbVStUjqQ1ji7AbevcQ4bOsyrgyqKtacgl3cGyRTWHcNVdSPmzmQXF6jzPwdF/kWsf3E6sz1uSA3e+vQoVi76OTNQJ2miHsve+wryCcPwyAOmO51odXqc3nwYCVEheGTGJKZeSxo1/kVQ5F/Esw8mMTXEDanB33aeQHR0AORyzy6A4E7CwsKQlCSnTxDC659lwnbIqrWTddlKPDJpKFY9PcXIEuXJPliE3oH+WP/idGY7b89uW/mYkWUqHBTX7VBi/YsPoHegv4llyrevenoKHp40zGDb2qoRADb8dgZTA68x67cPYUBMGN10giAIkFVrM1+ebURre4tRMQFxsQFxMQFxO6uggbAYgfj1gHExAla7+JriggeWNIrPEV+PCiAYQ1Yt4bIPZbJqaeBk2RDebNU+tToHry4wtTX5GaA0LhK9A/1N8jT59hl3D8HRcxXMPM112UpMHn0bcpSlzDzN7INFuNHSBlVVHbOdnyECMAoEslYjr2HobVG4cFltpJHyOHvWc0ygRzzLNHB2D2TV2kiAvx9zLTA8JAjSuEgo8i9i+gQps33G3UOQuS0Pz864g3ntZ2fcgczPjmFG0mDmeuT0CVIo8ssgjYtkto9KjEFNfRNq6pvs0si/x9avC8y2EwTRvWi/PWr4DwB0584Z/Uz0wBmnt1u125V16B3QYRIFy1uf4xJCmbVreevz/lGRzHahlfp1YZ1JXVlh++nyRgAwiYLlrVYAJlawLRofHh9lYvWSVWsMWbWEyz6ULVi1SulAs21JKnom3QIVQLCN93d+b1KAQFwsQFyAQPx6cTurYIGwQAGrXfye4oIGjmoUX5MKIFABBIKeLaITvzfffPNNdw7UlTWNGBAb6rVfNJRFVzFrYqeV+mluAS5Xd87ihOuFQYH+uHNoP6zboURwr0AcOHXJaL1Q2D50QBT+/nm+yZrlxFED8GluAfx8ffBxboFJ+6jEWJwvq8H5shrsP3HRZE2zT2SwQxrFGv53/CIevmeIx94Xdz9X3v4cE6BnmSCr1l0ILct9Z9SoUOuxdFoc87UXrt3A7vxaLJP3Z9aNrW1sw4YDV7FkWj9Eh5qm1Da1dGC9ogqzZdEY2rcX8z0+/OYq4qMDMfOOKGa7oxp5DfGSEMyfGElWLVm1RDdbtQRZtV5p1QqtT5bFKbQ+zdWN5Y9fvq7tsr2+sZlZV1ZopbJq2zpTY31jM/fK+gNkb5GdRtCzRVCtWvsQ1nUV140FjGvPsurGCmvLDogN67I9PCTIpK4sYFx7llVX1pkaw0OCMCgugm48QRAEWbW2889DtRgxIMAkopWPSB03MIQZ0cpHxcpHRUJRWGc2qtZcO3AzapZ/H3HtWT5qlseZGimqlqxagqxagqxau3jpvW/Mtn2Rd4Gbt2q32fbL17XcrBX/4S5f19rVznEcN2/Vbu6LvAtm29M2fMOlbXC+RoqqJTuNoGeLIKvWLtra25lbbRWW1eDCZTVWLprMrBur1emxMec0Pnz5IWzMOW1SN9ZSO2/Prlw0GRcuq03qyvL2rGxYP8iG9XOJRoIgCIKsWpv58mwjwoM7jKxQcV1X8c/iOrG2/iy0aXl7VvyzuPas+GdHNZJVS1YtQVYtQVYtRdVSVC3ZaQRZtQRZte4iZeoIXKltxN92nmAWbB+VGIPJo2/D02v2Mguyh4cE4bkZ47D03f/huRnjmO2vPJ6Ep9fsxeQx8czas6+kJOG9nUpcqW1kFmx3hkZeQ69Af7rpBEEQZNXaZ9U+ODbUYHVGh/kbWaI8vNV5Z2Ioiq802Rxly9uxw/v3xvdlOmY7b8fWNrQxo2wd1SjU8PV5HZZOiyGrlqxagqxagmwI261aqlVL9hbZaQQ9W7cuNOO0EdodhWacNOMkaMZJM076NmUDL/w112zbR/vPcotXf2ESZHPzd6/k5ry+w2x7fWMzN+f1Hdzxoitm2xev/pz7aP9ZsxqWbzzALd94wG6NBReruXmrdpsEI9GMk2YFhOc9Wx03bnAt169zbVotx3W0U+dRcJBn0tDUwsyh1Or0UFXVQS5LxP4TKmZ7jrIUafMmYlPOGea1N+WcQdr8Scg5XsLModx/QgW5bBBUVXXM9sKyGsREBCMmItgujfx7LLh/tNl2giA8h6ai87j273+h/uABdOh01CHucgXIqrWNL882orW9xcjGFOc8inMoxe3iHEqhzRofFcTM4xRek9UuvqbYarWkUXyO+Hpk1ZJVS3ieVXv9449Q9ubriJ71MBJWrERAv37UgWTVemZwkDBwxlwOJB+MY65dGIzDytMUnsfK0xS2m8vTtFajOQ3C677xr2/JqiWrlvCwZ6v28/9yZ6bdx1364xtca3U1dV5PDQ5SFl1F0gjv/Va0ftcpLJszHgCQsTUPVWodVj09xSQHEgA+/uocFPllWP/iA8z2wrIavPXpUaxcNJmZp6nV6bHsvf2QywbhyQfGMNtf23wY/aJCsXzBRKZeSxqzDxYhN78ML82dwNRQWFaD93YqEeAXiPUvyz32vrj7ufL255gAPcsEWbXutGofHBtqsD4D/H0wYVCoSQ4l3943IgCB/j4mEbC8NRoXGYiquhZmDmVeiRYtbRyu1bcy8zgr1HqcUDWitZ1jtlvSyGvw9QUkIf5MjXklWmh0bbh4vR3LkmPJqiWrlvAgq5Ygq9ZrrFqx9Sm2OcXtLKtVeA7LahWew7Jaheew2i1pFB9jaRQeW/vZCbJqyaol6NkiKKrWdvSt7UYbPAMw2mhavAE0AJONpoWbUAMw2WhauAk1AJONpoWbULPaWRrEm2GLNYg1ijUEBfjRzScIT5v4tLSgtboa7fV1QEcHdQhZtZ7J+txq/PK+aBNbFAB25ddCo2vD/J/FMtvzSrQ4W96EmeMkTNu0Qq3HvjNqjI0PYdqmTS0d2PZdNSTB/phzVzSzfXd+LQBgtoytcc8pNbTN7ZhyezhTQ16JFqXXmzG4Ty8jDRRVS1Yt4XlWbd03X+PHZ3+J8KS7Megv7yAoPp460A24vXL3sKHDvHrh+xNlDqbed59JoI1Wp8fRSiUSAvwweOQwZqDNWY0S98b6IDQ6FHJGUfbsg0W4984+aO/gIJebFmUvLKvB2LoL0Le2IWliElPDd1XH4MMBU++byAwGOqtRoqa+CUlJo5kaNf5FaPK9ihHD+hlpLNaeglw+3mPvi7LoqlvfLywsDElJcvoEIbr1WfYLDUWvxEEI7NsXvv60EYO7IKvWRsYO7mOwRIUDFm+Npi+YZGSJ8vDW6EtzJxgGSfGgCQAvzZ1gZNsKB839J1RYvmCikS0r1pA+byLS5k80aRdqWP3sVKZGXkNG6jSmRoIgPIteAwei/69eQPTsufCLCKcOcZcrQFatbXx5thFThgcbigkAYNZ9FRYTEBcjAIwLGrCKEQgLELAKJggLGrA0iAsesDQIj7E0CI/tzq/HbFmEx94XsmqJHvOhTFG1NHCybIiekMep1emRue0YOB8gfR7bFl2X3WmLLpSzbdHsg0XI//EqZMP6MffTLCyrwRZFAWIigpn7aVrSwM9CgwL88cikoUwN67KV8PP1Qf/oUKaG7INFuFLbiO9/qMZHK2Z47H2hPE6ip0DPludDprgDcD6AD+e5GsJDghAU4I+r6kbEx4YxXxMfGwZF/kXIZYOY7SMTY5GbX4aw4EC64W7g+x+uoaCs1m3vV9+oR0RokNk2+HRGVLM2Mm9uaYO+tR3gYPYa19SdO/n0jTIOLDtwqgIAcEPfAj9fXwQGsD+KWlrb0NLWjtDe7Ou3t3fgRksbegf6w8+PvfLUeEOPQH8/s+9xQ98CfUs7IsN6M9vrGm6gg+MQFR5s9nxtUwtiwoOZGtrbO1CjbUJ4cCB6BwWa7afwkM72uOgQGjhpxmkMWbWduNqqFWqIDg1gtgvfl6VR+L7blXV4IinSY++Lu63ajE+P4cDpSgQF+CE6wvgDtba+CX6+Pp0fmh2cSfsNfQs0DXrERPRGTf0N9I+5+aUmPDgQowZF48jZciTfJUXuSRXuHZtgMnAdPnsJd0j74YzqKm6/LdpkYDpTeg29Av3R3NIGALhjcF+TD+ofLtfi9tui8cPlWkwZO5DZnjR8AI6cLTfR0NzSZjiuLK5kajh89hJuvy0aVzU6Ew33yxKQV1Bh9Hqx48Gv60+fIDVKvxI7Ks/NGIeNOaeN0q+Ejgp/vjD9SuiodKWBbx+ZGMvUwGt8dsYdJilgQo2vPJ6ETTlnmBrEGv+0+Rj+88eHrXsQ29vRqlbDx88P/pERgC+ljblnxkLJvTax9rMTLi82wCqI4GjBA0eLMtC2YuLn4CT30ZeFJv1o6Wdxvzq6yTnrZ0c1dIcmT+yn7tD0zvbj3GMr/mv1c3ij5AJ3fFACd/ree7gbJSVUmYBq1XomqX/Zj7W/murS9cT2Do65pllYVoM9eT+iubWd+e3akga+XRoXid6B/kwN67KVGHpbFC5cVhtpENbopXUh4J3/5CM2sjeenD7So2Yt58uqHdJwq8zuLPVTd86C5/1xHz57Y6ZVz2FL5WUUL16EoH59MSjjLwhMSKDZIFm1nsfHRzV4crKE2bbrZC00TRYKIFToMPOOKPMFEE5rMDYhuOsCCCH+mCOLtkvDweJ6nK3Q4dn7+jHbK9R67DhRg8cnxBhppAIIplZt/9gIPDl9JAAgPetA5/GfUnnE/G3nCRSW1SDrtw+ZGfivIPOzY9j86ixmoJlWp8fTa/Yibf4kJA2PY15jyTs5GJkYa0h5slXDnmMl2Pp1AbatfIzZXlndgKXv/g8fvvwQBphZM5//1n+x4P7ReHjiEGa7N/STJY2u7qf/e2OP1VZth04H9f9y4BcSgoh774VvaBiNamTVeh5v/Otbl2/h5eptxOzZ6oysWrZVK+xX1n0T9qu5dv5+HD9fybyv/P04XnSFed+EGljPljUa+HZzzxav8fJ1rVkNwvNZGryhn6zV6Mp++r8395IX6uHQwGkj7+/83uSPpruLvNujQdzO0iC8Jg2c9q1xivvV0fvAehbE1xQ/C5Y0eOK6a3f0k6esu9qyxkncImucPWVbMT7iFIBJNKow6tXV24o5ooFvvzMxBMVXbjA17DmlxvD+wdh7SovfPNiHrFqBVVutbcbPRvazea2Ob58+QcpcqxMW8RcX4xev5e0/oWKuxfHnAehSgyevu7qznxxdG3ZmP+04UGJ9VC1BVq23zDh5lm88wC3feMDsaz/af5ZbvPpzpl3Df2Odt2o3067hv7EuXv0F99H+s2bfw1ENx4uucHNe32G2vb6xmZvz+g5uRdZhmnEK+P0Hh7jn//qV2fZ3dyi559bmWLgvXzj4bHzu0LNhSePx85VWPRvHi66YvUbq2n3cuzuUHt1Pjmp0dj+RVev5UK1aO1mXrcRC+WgslI82qSvLf1PtHeiP9S9OZ9aN5b+pblv5GLNuLP9Ndf2LD6B3oD+zbqwzNBw9V4HNr85itvMaNr86C5erG+imC4gICUJLa5vJfeP7tb2DQ7IskXnftDo9VFV1kMsSDTMhMftPqLDg/tFdtstlg6CqqjO5b7yGmIhgxEQE260xR1mKtHkTsSnnDFPDppwzSJs/CTnHS5gasg8WQS4bhPYOzmP7yRkaXd1PBEXV9girtrW9xcgaFRcoEBc0ENeNZRU0EBYgEL+edU1xwQJHNYjbxT9TVK0xr286jCHxMajRNhrZb+K9UsX7mor3ShW3C627UYkxJtcTX5O196r4HLGNaatGlgbhNVkaxNd0VIMr+snZGp3VT4dPXcF/Vz9GoxNZtT3Hqn3yz/vsilo0RAQWXek6IvB8pcsiJy1p6CoikIKDjHnr4+8MUbX8/TAXaWkpmtmdEdf2auwq4lqswdzz50n95CqNzuinlNe/IC+UgoN61oyzq11C9p1Ro0Ktx9Jp7PyxC9duYHd+LZbJ+zNzKJtaOrBeUYXZsmgM7duLeY0Pv7mK+OhAzLwjyiUaahvbsOHAVSyZ1g/Rof4mQVE04zSecfJ5nCs2HQQArH52KvP1H391Dor8Mqx/8QFm/mH2wSLk5pfhpbkTmMUzCstq8N5OJeSyQczCFVqdHq9tPox+UaFYvmAiU4OjGgvLavDWp0exctFkpkatTo9l7+2HXDYITz4wxiUanNFPrtboaD8t+NM+bH19JgjPhTaytpH/Fiig8R/AjL6bJBvQZfTdWY0KO1Y90GX03Y5V87Ap5wzihrCj7/6YOg7ny6qhAVym4ZOVczurniTfbKeNrI0JCgoyWU/k/82q2CRcqxPeN41GA5VKhYrqdoxKjMH5smrmh+35smqMTIxFhWCtmT9XJpOhoroB/aJCoW9tg1anN/lAd0Qjaz2RpVG8nugKDRXVDTb3kzM08n09dPhol/eTr48PjUxk1fa8qNqeUI/Ukgbxz2TVsq1ae/NpT548ycnlcg4AB4ALj4jk0tLSrKpj/JuMj0zO/fmjT5p9lnpC3rH4HHvqPdujcfsXCqO+7hUSxr38yu9c2k+z0neRF0pWbc8LDnpwbKgh+IZHnCPJB9+MGxhiEggkDMaRj4qEorCOmafJBwDx1xHnabpag7B92zENliXHklUrsGrDQnrBzx8mM3s+ICQ+Now5s1+/7SukL0nB8GFDMf0Xf8CCWfdi73/+ifT0dKSlpUH24FOGgBFx8IhKpcL4O2UIj+6HQ7n7cMM3HGkr38S+LX9HWloaMjIyjIJQKqobmO6DJY3C92XlSAoDYViBO+7QIGxnaXBUo0qlwuix4zBIKsW2rZ9hV341fMsP443XVxr62hX9dLKo2mqrVvvtUQTFJyCIatTSjNPTZ5w8aRu+4dI2fGP2tV/kXeDmrdpttv3ydS03a8V/uMvXtWZfM2/Vbu6LvAtm212tgW9f/cl3NOMU8Po/v+WeWWM+t+/tLd9yL77/FTPIKyMjgwPAzVv+L6NZV1paGieRSAwzk7QN35jMaPhzvzx0klu+8YBhRiM8l5/NvPj+V9zqLXl2aeQ1PLc2x2yOZMHFai517T5mgIu7NJjrJ2do5Pv693/NNtIo7mtn99PCP+2z/rkflMCVZ75NU0DK4/QOsg8WQTasH2TD+jHzuwrLanDhshorF01m5lhqdXpszDmND19+CBtzTjPzu9ZlK7Fy0WRcuKxm5ne5WoOwvbi8lm66CH1LO/O+aXV6NLe2Iy4qhLnWplAocJt0OH4mG22oFgMAcrkcGo0G2dnZZt9ToVBAJpPhtgTj/TPF54rXPG3VKF5PZGFpPdEdGizhiEaFQoFx4+9EsKSvkUZWX7uyn7oiSXUJ8a+m0x+jmyGr1g6rNjy4w8gaFedMinMkbc2ZFNq0vD0r/ln8ns7WQHmclq3a2x7bFhUAACAASURBVPpG4oq63shiFOfqsey74NBwTLlvGr7c97mJfefj44NH5z+LhamvMC3IqKgoJE28Fw/84g94JSXJyAr08fFBWloaFi/9neEYK3fQGo3CYyyLUXiMlb8oPOYODa7QKJFEof/QcTh24H8mGvm+7iOb4/R+Kr5Yj8/enGXVc1i8eCGiZsxEn/kLaDQjq9ZzeeGvuUzLhV/0N5f/xR83l/8lPG4uT5M/bi7/y1kaKI/T+uAgYV+a61fh/Xxn+3EOAJeWlsYMGAHAPTr/WbMBJQAMgUDi+wqAe2bpiyb33h6N4udPqIH1/AmfOdbz5w4NztYIwBAIJNbI3wdX9NOjK3aTVevh0MBpI12t9XlD7U9LGsytedLAyR44+fu2YuMBbvmmA2b79Z3tx7nlGw8YPnSFA6dwrY7Vxre/u0Nptl34YW7u2bJFo73ricJ11+7Q4EyN5vpa+CXFFRpojdPzcXse548XfoS2ssBrZ+jHz1Vj8/YakwjXCrUeqks63Bbmjz99uMvsTiTyYaF47R+fm92JRD6sNzbuPISrjCjbvBItBoQAKtVFbN5+3eka+PZFEyLxZtY+I+vY/0YLPtxShSF9Az3WqnVnfrBerwdwc9NgzgfwcXDRo7nJfD3gkYmxyM0vMx9dWV9n2V1ygkaHHS4P0OCoxhs619Zt7uA8vIMIsmpt5f2d39u836En7LloSYOlPM7j56s8+t55k1U78Z77OLlcbnKvcnNzOQDca2v/adbeGyObyI2RTTRpX/La+xwAbuOnO8iqdZJG6ai7uKnTfm6iceOnOzgA3JLX3ierlqxaGjitHTiFf9jeXPvTlnXX3394gDvwfTkNnIKBc8MXZ0z6TdyXrHs7c+ELHACutLTU6N7yaQ5dJfLzKRLvb91v1C5MkRCew7q31mh0tNiAuzW4QuMf//QWB4A7c67ISI+wr13RT/Ns2FaMBs7ugaJq7Yiq5aNLtxzrDD9fOJFdGOBgcT3OVujw7H39mHVhK9R67DhRg8cnxJjYrrx1uunQVYyND8HU4ez6uI5quHCt+afatXFd1M+9grulERggCfJoq9adUbXLPzyIy7WN+PvLySYl0/joyKAAfzwyaahJ2TWNRoMB8QMR03cAXn3rb1g2/wFkZmYaCiDwifVXahvR3sEZRWFqNBoMHjwYktj+WPjSW4iM7Y9W1QGjc/lozT15P6K5td2ksIA1GvmoVj9fH/SPDmWWlss+WISVL/0CP55VQvwx4uPjg4n33IdHfrXaoEF57AiSk5Mhl8uRm5vrNA35P16FbFg/s+2sfhRGtZrrJ41GA6l0MHpH9sHCl/6MXzx2n1GhCr6vrdFoi4al7yisLoCglA5E3NLnKSXFzVCtWhvha7YWltVgjKZzD0BWXVmtTo+jlUo8fn8k9IH+eITxB7UuW4mlc27HhctqPC1PYv7BPX5/HFRVdUiamMSsq+mIBq1Oj6M7lHjruTtxtOAyXpmRxNT4VuoYfLy/APePHI6p4+M98r64u1ZtYGAgggL9mEW+w0OCEBTgj6vqRsTHhpm0SyQSvLHuX1i94tf49YLp+PWCzmPCD2N+TTNZlmhybm5uLh6bMxerXphtOCY8FwDiY8NQpdYhLirELo38NRT5FyGXDTK77qq70Wq2j4IC/VyuwRLm+tGafpJIJFAocvHAjEew9jePY+1vzPe1pX6yRQPVqqU1zh5p1faE2p+2rrse+L6crFqRVbtm6/Eut/CyZnus32R8xOWdvmDWvutqnezkyZPc5i++c2ibMU/Y6sxeDbasJzqq8dCRPC7t7/vc0k8UVUtWbY+zarcr69A7oMMkIpUvUBAdGmBSzAAwLlAgLmbA27Z8gQJxMQPeMuWvW9vYyqw9a62GcQmhzHahBrHGkmud94ys2k74bcVkt/fpctNoZ2/I7IqNrT1hc21bNdhTbMBb+intw8PYmzGHrFoPxu0Dp7LoqldbtW9uzsObT09itmVszUOVWodVT0/ptj0XLWmwdy9B3gr11Hvn7ufqz58cR3yfMDw5faThAxKASfUb4QeuNC4SvQP9mfd2XbYSQ2+LwoXLauY6mKV1suyDRbjR0gZVVR1zTdNdGt2hwZH1RG/op0WrcvDpazMsPoPab4+iePFCDP9kC8LvmUyjmRuhWrU20tDU4rW1P8X7HbIQ7iVIWMeoxBjU1Dehpr6J+YUoPCQI0rhIKPIvYvoEKfMa0ydIsfXrArPtIxNjUVhWY3YtcPoEKRT5ZZDGRTK/ELlDozs0xMeGGf5Oemo/WUPFmgwUL16IkNFjEDx6DP0R9vQZp7dbtXtONaCpRW9kc4rrurKsWOExlhUrPCauKyu2UVm1bS1pELezNAjPEWsgq5Zt1T45faTBnuO/eIhnGbw9N32ClLl9ltDSc3RrKv4Lj3gm5GqN7tDgzC28PLmfUl7fg+w/Pdzl89dWX4/2+nraTqy7oOAg24ODvC2h3BkBIH/f/T0FBzEKIDhaiGLvgRPcqn9/49TAMVcXyxBf0x0aXBFgZ4/G97Z9w127ds2l/TRn5ecUfUPBQT1rjXP9rlNYNmc8tDo9MrcdA+cDpM+byLRc1mUrUVPfhIXy0Uzbx1IOWmFZDbYoChATEcxcK7FGg7PWXYcOiMGUOxJojVOwxnnxqhq/STHtN37WwduFKYw0oHU7lJhx9xAsz/gAzzwwHAsWLHDqOhk/M6qobmCu1dmi8ei5CqaGddlKTB59G3KUpS7V4Mr1RFs17vnX23jwwQcN98sV/fT2J0qr1jgBgGtuRltDA3x79YJfaAjgQ6tv7oBq1dpIWVkjFIpaNLV04EpV5x6VBw8dYhYPuFSmhra5HUqlDlUlpoNWUYkWVdebUdShgaKt0nQdQ61HVZUWOo0fFAqtSbslDU0tHSi9qEaAvw927lMwiywcPaVGiC+w88vDqBLVtuUt5GC0oeCHBoS210BbSbVqAUBTr0NY7yDml41RiTFQ5F+EIv8iNvx2BnOdbMbdQ5C5LQ/6i9+itjaWuQ7GB2mZWydb+MfP0Hr9B4Q/da9Je8rUEUh9538YlRjjsMbNr7K3uHp2xh14es1epJn50uYMDfx64voXp5vtB0v9tOy9ryCXJTpF43uvXURtba1L+8mWWrVNReehUeQiePgIRE6bBt/QMBrVyKolq9Zeq9bZJc3+800xt/dYCVm1Aqv2lb9/4/D2bksztnNnCovN3vuuyjFu/ngrJ3voqS7LMbp6CzpXboPnzrKV1mo8evQoV1pa6tJ+siWP89pH/+aOD0rgSn79AtdS5d2frd4EDZw2svazEx5f+9MVgzkVeTfm06/Oc69t+takn8Q/i++lpeL6rlhPdLZGlgbxezqqoTvWXT2ln2al77L6Oaz9/L/cmWn3cZf++AbXWl1NIxoVQPBMPj6qweNJEUxbdM8pNQL8fDBBGsq0RfecUsPXF5CE+Jts+cXbohpdGzo6wNx2rEKtx4mLjWht40yKF/Aadud32kizZdFM+3jPqU77eMrt4UyNeSValF5vxuA+vYw0UlStMV8rL+DTb1TYnD7dsEYltN5Y62TTJ0hNIjUvXLgA+PfClsNlLl1PdESjcC3vuRnjsDHnNFMDH0W688AZzJ12h90aunPd1ZLGL498jxOldZgzdazD/cT3t1jDvDf34rM3Z5EdSlZtz5lxpm84bHYT3tVb8rgX3//KbPv2A+e51LX7zG7SW3CxmntubY7ZDXDrG5u5F9//int7y7dm25dvOsCt2Oj8jYJpxmn6fr94+0vDz2kbvuHSNnxj9vVf5F3g5q0y3S7ql7/8JZeVlcUdL7pi3Sbn5ysNx/Ly8rjnnnvO8PNza3O63OTcXo085jY5FzLz9x9zs578LbdmzRq7NLy7Q8mlrjVvVdrTT2Ic6adf/vKX3Ct/znK4n+at2s19kXeB2WaLVUvQjNMrYOVx8jO57iy5J2zn9ZjTGB8VxNQgfF+xBppxumbGec8992DC1JlIGCfHszPuYOYGCnP/NuWcMcxSrl+/juLiYkyZMsUwi+GLa4g1uGvGGeunQVrmJsy+bzT+vuYNmzTw7SMTY5kaeI229hNrtmdvP/3s/kcwYMxU/OnlJ10241zwp31W745C0IzTK4ODWOsU7i7y7o7goL3HSiiPkzHjdHRd7NPdudwfPsxx63qiO9fuSkpKuLNnz3qUJkf66TfrdnLFxcUes8bZoddzLdevc211Go5rb6epYE+dcfakPM6u9hJ0V+1Pd+y5eKW2EZU1Ojw+ZTjlcQre762Pj2Hx9OF2z5wcnd1NHdsfnx8pwhPy8Tav1bljdpe59Sjaq4vx0Y4v8cgjj+Avv3+y29Zdzc3uPHEW/PmRi1bncdZ98zV+fPaXCE+6G4P+8g6C4uNBuB6yam1EuJH1rpO10DS1Yf7PYpmBOK7eyLqppQPbvquGJMQfc2TRTL2WNOaVaHG2QoeZd0QxNVSo9dh3WoM4SS+M6h9MVq3Aqn131znsMbOLxd92nkBhWQ2yfvsQ2/I/VoKtXxfgzwvGonfv3oiLizNqr6xuwNJ3/4cPX34IA8zUXX1sxTb41xZix4a3mO3pWQcAABmp0+zSqCy6gszPjmHzq7OY+Ydanb4zP3H+JCQNj2Ne45FXN6Pqwims/MX9ePTRR23WyPfTtpWPMdut6af5b/0XC+4fjYcnDnG4n1Qqlcn9ckY/LXknByMTY/HS3Ak2WbUNyuO4mJ6G0DFjEJ++AgFxcTSqkVXrmVat0GKxVM7OVTlo7t5z8S/bjlMep+j9HvhtNjPIiu83c7mBfHvBxWrugWf+xP31r39l3tvL17XM+8bfr6OnLnDLMrcxNfDvbUmDuXZew/Hzlcxni9dwvOiKWY38tYXP1t69e7n9+/fbpNHc821tPwmv42g/LVmyxOh+uaKfbAkOarlaxVVnb+fqDhzg2nWN5KG6CRo47Rg4Xb1ptKUcNFdsbG1p3ZWiatlrnLbWYBW333nvg9yDz73l1rW77sqZ/Orb09wzzzzD9enTx2tzS8c9sIh74YUXXNpPVKuW1jh7nFX78VENpo0MNrE1+ajWvhEBCPT3McnT5NvvTAxB8ZUbzDzNPafUGN4/GN+XNTLzNPNKtGhp43CtvpXZzkfmAjCJmLVGI68hLjIQVXUtRhopqtZ8VK21O2qw1slyc3OhbeuFS41BbokWNbdW567dURoufY/zl2rw+Ny5Hh1Va07DHz/YhZCQEKTOneqyfpqz8nPseutRskPJqu15Vi2Lj/af5Rav/txsjpmjOWj1jc3c4tVfcB/tP2tWw/KNB7jlGw/YrbHgYjU3b9VuE1uLZpzsGae1/f7uDiX33NqcLq5Xad2zUXTF7DVS1+7rMj/RkobOZ+MLm58N4+fz8y6fz9S1+7iMLUeon7roJ8rj9HyolL6N1NTfwLpsJTP6tHegP9a/OB3rdihNNpouLKvB0XMV2PzqLGY7/01086uzcLTgMgrLapjt6198AL0D/U2iAflv0wvlo7FQPtpujftPqLBt5WPYf0JlooFgU1hWg5iIYMREBDP7rLCsBu0dHJJlicz7ptXpkaMsRdq8idiUc4b5HptyziBt/iTkHC9hbmKefbAIctkgtHdwdmtwdJPz/SdUkMsGQVVV16XGgIBA6qcu+skm2tvRWl2NNrUa6GinP0aKqvXcqNpRtwUYFSAQFwsQbxot3hRa3M7amFpYoIDVLn5PcUED8XvaqlF8zcLKZgDAqAG9yKoVWLW/myczsvTE9pzQuhNagbw9V155FX/etBd/+d1ihIcEmbxefE2h1dfe0gSVSgVVQ7DRNW3VILxmeEgQcxNo4TVZGoXniK/Hek9HNdrST87QoNFocK6wGN+Utrm8nzbnnEfOmrlWPYfNpSU4m3w/gm6Lx+2bP0KvwYNpVCOr1nOtWqujEouuOBQRePx8pU1Ri3ZFTprRKNTwwrtfUQEE0fv935tfdBnN3NWuHxkbsrmp037OAeAAcBKJhEtLSzMJMmHd20NH8jjpqLsM54aEhRvOtUVDV1HftkRcW4r6trQziSfsjmJOQ8aGbG6MbKKhryMjJUZ97Yp+ssWq1V+u4M5Mu48rnv8Ep790iTxUN0EDpwNrnI7W1bRU09IZdTedVfvzoy8LaI1T9H4Pd1Hhpas1z9LSUi4kNJyLHjCEO3OuqPMDOiODA2A0eLLWyUpLSzmJRMKNG38nN/t3G7iP9p81Ode29e4vuqyt/NzanC5rK6eu3ddlbWVLa56OanTGuqs5DXxfJw4dyc3+XRZX39jM7Gtn95MtA2d7YyNXnb2dU+fs49obtDSiuQnayNpG+I2seevzZ/FBeOkv202iZHnr86m7I/Hbd3ebRMHyVumiCZF4M2uf2fan7o7Clr1HUXjONEp2zyk1Jg4MQUV1Bf70YYVJlKyzND51dxR2fVeAqxW9oa0M9lir1p2VgxoaGtCkb4VWp2dWdIqJCDb8WxzZ+emWbdA1avHH9Z/hhxoOYwGkpaVBo9EgKysLGRkZRutkwvOzs7Oh0Wjwp79/hvDouM6qUsteNjrXGg3itTpW1aiK6gaMSozB+bJqZtWp82XVGJkYi4rqBqvW8mztJ2s0musnZ2jg+3rxig8gGzsS+0+oTO6TK/rJlo2sfUNCEPN4ik3PbsWaDETcMxnh90wmy5WsWvfNOJ2dc9YdOWm2aqSoWtP3W/xWjl35tNJRd3Hjxt9pcp9yc3M5ANyTv1ltNp9WLpdz0mGjTHJ6P9+TwwHgtm/f7vScXnv2j+3uvGNnaJg67edcXOLtJhr4+7R9+3aX9NO8N/e69tkdlMCVZ75N00ayat3HG//61iVVTtxZBcUejTRwmr7fL97+0q4KTpGREi4lJYX54QmAe2bpi2bXyULCwrmJU6cz7xt/riuqSHW1ybmta3nuqHTlDA29QsK4x2bPZWoAwP380Sdd0k+2FHmngbN7oKhaO6Jq+Vq1Yk5ebMSRH+vxm+kDmO21jW3YcOAqlkzrh+hQf5vbAeCv+ytx77AI3DWIreGTvOsAgMWT+jhVIxVAYEfVbk6fDq1Oj9c2H0a/qFAsXzCR+fqMrXmoUuuw6ukpiAjthbS0NIPVx0dT5uaXYePvZpi08TbiezuVyPrdTGa7VqdHRGgv3PfwIhz84hOLGlg1VXkNL82dwLQcsw8WIf/Hq5AN68e0TXmNctkgZrut/WSPRmdqeO/XDzD7OvtgEZ6YNhLPLH0Rmz74m9P7yZZatcWLFyJqxkz0mb+A2RYyZgziX003Oq6UDkTc0ucNx9vq61G14QOEjB6DqBm0nZk1uH2Nc9jQYV69O8rRijxo/CXMnRwCNSqseWFylzs5fLJybudODsnsqifm2vlQ9zUvjMP+EyrEDWHv9DBrSueAqAGcqjG86CoAePTuKO4kLCzMaC2wX1Qo9K1tzHU0rU6P5tZ2xEWFmF3n4tfJzMGvk5ldt/rpum0dnFUaRjEGJUtrdZawtJZnaz/Zo9GZGiz1dU19k13PjiWNvj4+Vl9L++1RhIwZ41jcxh+WQ19ejrglz9OIaO3aMnWBbcRE9DYMUsIBiR+IRiXGdJYtExQgEOZrDYgNwyuPJxkVILDUzg+afH7YKylJJgUKhPlh/IDpTI0EG75fly+YyOwzYb+mL5jETIrn7+2TP+/cvaO4vNbkCxEAvDR3AgDgRHGViYb/HuwsCDDh9jirNIiT/3kN/HuIk//5n/kdRMy1vzR3gsmzZW8/2avRWRo6B8cbJhomDosyeU9n9pO68Ybbnl9+0Lz9ky3wj4igP2grIavWTquWj1iNjwoyKRwA3IxYvX9UpEnxAuBmxKp8VCQUhXVm2x8eH4WvC+uYtWf5AgUV6s4/fnNRtY5qfHh8FHad0GDkgN64M9Fzo2rdbdVu/F8xkkb2NZq5Cz98ATDrmd4+9m5Eh/dG3tGDRl+IFAoFkpOT8draf2KMbCJSpo4wSbJPTk6G7kYrHn9pDV5JSTIMBmMlWiQnJyM3NxdJE++1qEH4vqyaqsL3ZSX6W2oXflFjFQKwpp8c1eiohuTkZJRf02LXF/uMNFSVnDb0tcZ/gNP76WRRtdVWrdh2LV688Obfb8E5+EVEICg+AQAMti1/Tnt9PXTnztGgaQ8UHORYHqcjdTWtzuN0Yd1Na2t/Zh/8gcs96bkJ1t21rZi5Plux8QC3fNMBZr/yuYDPr95qFFySlpbGSSQSQ+BI2oZvTIJH+HO/PHSSW77xgCHQRHiuNRr4QJjlGw+YzT80p8Ha9oKL1UYabe0nZ2h0RAPf1yve32WkQdzXzu4nW/I4xYE+5ZlvG/47fd9krmjRAsPP17ZuMZxT8PBM7vigBO7iinSK9LEDsmodsOgcqaup1emxMec0Pnz5IWzMOc20RN1Rd9Pa2p8FZddwpuQ63Xhrv5D6AD5mvJzU1FT0DgnD7g1v4nL5JQBAZmYmMjMzkZqa2uV1U1NTIZFI8MqyZ1BXXdXluV1pcBb/eff3eGLaSFMry8cHzy163OiYQqGAj48PkpOT3arR3nvF9/XHf11hsa89hfhX0w3/BcUnGGaZ8a+mmwQQxS19Hte3bYU6Zx/9wZJV63qr9lasVdvU0oH95xowW+aZlo43WbXrspWI79WA9JdToVJ1rnlKJBKkpqYiIyPDsk2652s8/9xTqL1WCQDoHRKGF5f9yhD96Q4blD+W9ecXoVAosP3AeaN2Hx8fDBuThB/OHjexk+VyOXb+d6/HW7UA8Ermx8j+x+uGLzjhEZF4fukSQ1+7wtJ2xKoVYk1UbcWaDFRv24pRn+9FUEICjYhk1bouj9ORuprm2oXHXVV309Han11tqXar5nEK+4zVr+Jj4n49efIk99ZmhV1J9CdPnuTUarXDGmwteGBPor+7NTpbA9/X7uinh17d4ZScTKE929U5BQ/P5Aoenkn+qyfncSqLrnp1Osr6XaewbM54ZtvHX52DIr8M6198gJmDpiy6gszPjmHzq7OY7VqdHk+v2Yu0+ZOQNDyO2b7svf2QywbhyQfYIegrNh0EAKx+dqpdGgvLavDWp0exctFkk3D/rn737sbdz5Wy6Co++PwMNqdPR2FZDfbk/Yjm1naTGYtwZhMU4I9HJg1lplGsy1bCz9cH/aNDmbl92QeLcKW2Ee0dnEkaEX/f3KGhq/xET9Ho7f20aFUOPn1thsMzTmvP0ZeXo/DRWYidv8Cm69zKUK1aGzlbXIfN26+YRLg2tXTgyCk1EiIC8O4n+0wiXA1RtMNC8do/PjeJcOWtUfmw3ti48xCuiiJcebs1IYzDkROFiGyrMmmvUOuh0+gAAJu377NZI6/hjr6ByMr+xkQjX6fXU61ad9eqNawrxYahSq1DXFQI88tIeEgQggL8cVXdiPjYMPbaVGwYFPkXIZcNYraPTIxFbn4ZkmWJZs93tQZLeILGntBPHS6eyySpLhn9HJSQgDtPnaXR0JPXOHvCjLOlvdVovcXSfobO3kuQ1S6+pnhNyBl7LtKM0/j9/vLZCTw8aTByT6pw79gE1Dc244fLtZgydqDRaw+fvYTbb4tGRGgvHDlbjnvHJqBX4M3vrGdKrwEA7hjc1/DavlE3K0NdUzcarit8LU9zS5vhumINJZfrUFHdgCs1DZCEBSHQ3x9XahsQ38d4rbqu4QbaOzhERwTjSk0DYiJ6IzDAn9leW98EP18fRIb1NrS3tLahpv4G+seEMdsB4EpNA/pKQtDS1oam5jZERxinNtXWNyG4lz8C/f1xTaND/xjjgeuGvsVwHn8tP7+bXx7b2zsM59XWNyEsONDodxBr0Le0m2gUaqjRNqGvJNSshtr6JpPfge+ryLDehv931R4WHGT0OwCdBRD+mfaAxWdQ++1RFC9eiOGfbKGC7T194OwpeZx84Ex0aAAzB5IPxhmXEMpsFwbjiAN7hDPUh8dH4XR5IwDjPE1he21jKzNP01qNk4aEMzUINeaVaHFZ3YEnkiI9dsbp7uCg786q4Ovri5a2DgT6d/ZrY3MbQnsZf1gLj7W0daBF3wI/nw636Azr7YvxA3vTJ503fSiHD8T9SUO7fE3FmgxUfdhZJo/yMGnG6RUzTn7WZamupiPrifwssas1TXfX/rymbqYZZw95jgnPxZpnq62+Hu319RQJ201QHqedWKo/Kt5LkIVwL0Fz7cJ9+ky+dYrqbtqqkb8GX/uThaW6mgRBuB//iAgaNG+lGWdPsGqnDA82sj7FNqc4R1KcQwkY50iyciiF57DyOIXnsNrFx1hWrPAYS6PwWIVaj2MlN8iq/YmvlRfAaS/RJwjh/A9lK6xaopuhPE7bWPvZCWYOZHfvJShsN6dB+L727LlIeZw95zkmOHqWCcrjdBdvf3ocL8wex1wv7O69BLU6PTK3HQPnA6TPm8jUsC5biZr6JiyUj7Z5z0WKqu05zzHhudjybHHNzWhraIBvr17wCw0BfGj1zS1Wubvf0NvzOC9X1uPgoXqTHMqmlg6UXlQjwN8HO/cpTHIoAeDoKTVCfIGdXx5GFSOHMq9Ei2C04eiJc5C0VZq0V6j1aG1shKq+Fl/kNDA1XKnqzLM8eOiQSTsAXCpTQ9vcDqVSh6oSU41FJVpUXW9GUYcGCpEGyuO8SUNDAxSKAvoEIVxi1QLWPctNReehUeQiePgIRE6bBt/QMOpAsmq9w6r1tHJhZNXSc0zcGlbttY/+zR0flMCV/PoFrqWKnkl3QfN6GwkK8DPaBJdVjEC80bS4GIF4o2lxMQLxRtPiYgThIUFdahC3szSINYo1iDUWltXgYlU9PQAE4UmWYWQkeiUOQkBMDHz8/alD3OUKUFStbfAFEJpaOrA7v9O2nC2LZtqie0512qJTbg9nWrd5JVqUXm/G4D69mOXvKtR6HP5Bi/BefswSfZY08JG1AX4+mCANZWrYc0oNX19AEuLP1JBXooVG14aODiDALxAPjg310eANzQAAIABJREFUWKuWomqJnmLVUlSth39hcfcbDhs6zKuDKhoCSzBiRD+EBQfgu6pj8OGAqfexA3HOajoDcZKS2IE4Gv8iNPlexYhh/SBnBAMVltXggrYAMRHBkMtNi0FrdXqLGgrrjuGquhFzZ7IDljT+RVDkX8SzDyYxNcYNqcHfdp5A8oREXFM3Qy733OAgdxIWFoakJDl9ghBe/ywTtkNWrY0MiA3FhcsarNuhRPq8iUibP9HIEuXhrdHVz041skR5eAs0I3Wa0c/CQXP/CRVWPzvVyLYVDprWaHhk0lCsenoKs51/zw2/ncHUyGvI+u1DAICa+hv0ABAEQa4AWbW2UVjZjGMlDVgwMdqlxQbEBREcLXjgjKIMvE3tiZBVS/SYD2UbrFqupQVt9fXwDQyAX1g44EtzoR45cHp7/tvBUxUAgKnj45kzQE/eS5Bvl8ZFonegP1PDumwlht4WhQuX1SYaKI+z5zzHhOdiy7NV983X+PHZXyI86W4M+ss7CIqPpw50A/T1xA6CewWYHLN2Hz++9iuLkYmxKCyr6fL8KrUOvQL87NpLMDwkCNK4SCjyL2L6BCnzPaZPkGLr1wVm2wmC8Bz8QkPRK3EQAvv2hS9F1fbcGae3W7XflzXhfOUNLLon2ui4s7fw4l/L48xtxOzd6oys2puQVUu47EPZBqu29dpV1B85goCYWIQlTYBvcIhbZ8bm6OluDFm1dugvu1qPyzX1Lt002tUbV9u6uXZhWQ2yD/yIN5+e5LH3haxaoifgic+WsugqlEVVRseGD4xGeHAgBsSEYkBsKLRNLSi+pP7p/7Umrx2REIUBsaE94h7RwGnnt6yQ3v6G7cCEAxaPo+uJ/OB5o6UNqqo65pomP9jZq4Fvn3H3EBw9V8HUINQY6BdAa5w0cBK34LO1+tPjWLHobod+p6JLalTWNGBATBjksgSvHkTJqrWRkmud2of0DcSWY517WC6cyF6zPFhcj7MVOjx7Xz9mgYQKtR47TtTg8QkxzOIETS0d2HToKsbGh2DqcPYO745quHDtBnbn12KZvL9FjYWXW8mq/QmyagmXfSh7YAEEZwYGVlY3YveRCwCApBFx3vkFlGrV2q7/+PkqQ61XVl1ZYa1Xc3Vjheex6sYKz2PVjRWe56gGazW+8a9vqVYt1aolPOnZamvjWq5f51prazmuvc1lmlxVp3rXoQvcnz/5zuv+niiq1g52HC42WKPiurK8xQp0riey6saK1xPFdWPF64/iurG8heosDXz7ppwzXWqMiehNN58gPIjmsos4dfddKHzsETRfLPM6/bOnDMGKRXejsroRqz89jqJLarJqe6JV+31ZE/RtHCYOCTGxNE9f0iE6rDMkXFz3lY+KvTMxBMVXbjBrz+45pcbw/sH4vqzRJOIVuBkVW9vQZhLx6iwNivN1GBgdxGynqFqyagnPsmpbKi+jePEiBPXri0EZf0FgQoLHW7WW3gcAnnxwFMKDA8mq7WlWLYt3dyi51LX7zJ9bdIWb8/oOE0tUaM/OeX0Hd/x8pdlrPLc2h3t3h9Jsu6MaCi5Wc0+9vYfZTtuKkVVLeNaz1d7YyFVnb+fUOfu49gatyzTtOnSBO19W65bf//L1Bu7Pn3zntvcjq9YNVNfrcPhMucnx7INF6B8dihfnJpnUleWtz6PnKrD51VnMurG8Pbv51Vk4WnDZpG4sb8++NHcC+keHmtS2dYYG3p79d/osZB8qMmknCMKz8A0JQczjKZA8NMOlm1gPiA1FQ5N7nMIBsaFYsehufJ1/yTADJau2B0TVFl5pQmyYr8EKFRcsENd5Ff8srhvLqj0rLkAg/ln8no5qYNWmVZyvw6Qh4QaN+881YLYsgqxasmoJD7Fq3UVldSMU+eV46sGRbn3foktqfJ1/yeOsW8rjtEM/AFy6pjE6Ls6R5Gdv0ydITQoLCGeYz80Yh405p5l5mnwAEH8dcZ6meNZpr4au8jj/mXMaD04YjI05pxHWqzd++393eex9oTxOoifgqc9Wd9WqrqxuxEf7C7FszniPGTxpxmnHjBPozOP8JO86AGDxpD7M15682IgjP9bjN9MHMNtrG9uw4cBVLJnWD9Gh7DqTf91fiXuHReCuQeygHEc1WMrjbGrpwIffVOGpyX1xQtVMwUE04yRuwRlndw6cAKBtasH6Xacw+96hGDEwqtv7gjaytpFwwYxz1pTOwUhjZrYXqFFhzQuTu5ztfbJybueMM5k941zzwjjsP6FC3BD2jNMRDYVlNTirUWHHqgdMSvoJNX72x3nIPlSEuAFtkMs9d8bpTmgja6KnPMte8bkbHIgVi+7G+l2nOsv8dXPVISqnbweHz5RjYL8wo/zK7INFhp/N1X3lfxbnafI5lsKBS5ynyQf78IOnuPasrRpM8jRFGsQan5kxDm9uzqObTxC3MJXVjaisaTSqWxsWHIQRA6N+amswHHdFfdplc8Zj/a5TuF82sFtnnmTV2sj3ZU2o0LTg0fGRRsf54Jz4qCDmziV88M39oyKZO5MIg3W+Lqxj5mnyAUIV6s5oV3GeprUaxg0MYbbzGuSjIqEorDPRSHmcNyGrlnDZh7KHWrV8vdkRA6MwfGCUxfVG/vUNTXrDQCqXOSfPlLdtu23mSXmczs3jfG5tThfnVlqXx1l0xew1Utfus5jH2ZWGL/IucPNW7e4ih0rLzVrxH+7ydS3lcVIeJ3ELPMvuIvfkJe79nd9z7+/8nrt8vcGha9Xr9NyfP/mOq9fpKY/TWyi7Wm9yrLCsBu0dHJJlicwcS61OjxxlKdLmTcSmnDPM627KOYO0+ZOQc7yEmUOZfbAIctkgtHdwzDxPSxoKy2pw4bIaKxdNZuZ5anV6bMw5jQ9ffggbc05THidBeDgVazKgLy+H9tujqFiTgYo1GVDn7PNIrXJZApbNGY8nHxyF3UcuYP2uU6isbrTrWuHBgQbblqxaL6DkWgu+K2nAfSNCDVaqOAdSnFNpTc6kME+TldcpvqY4r9OSBmfklm5X1uGJpEiyasmqJTzEqlVKByJu6fNQ5+xD1IyZaK+vx/VtWxG39HnEv5ru0b+ntqkFH39ZCAB2R+vyeZ5uj/YlG8I+q9abdkcx125pdxThcdodhaxawvOereODErjT9002Olae+bbJMU/mfFmtQyX2ck9e4nYdukBWrTfwSkoStigKsEVRwCwckDJ1BG60tGHZe18xixvwO5rMf+u/zOIGfLTtsvf240ZLG3MjbGs0XKltxN92nmC2j0qMweTRt+HpNXuZGnkNT6/Zi8lj4ml3FILwQKJmzDQ5pi8v9xr9IwZGGUrsffTlebss4MqaBrfurEJWrR1WLQAEBXA4fUkHAMwIWN7q7BsRgEB/H5MIWN5ujYsMRFVdC3OnkrwSLVraOFyrb2XulsLbreY0WL9bSiiKrzR1sWNLb3xfpkNwYBAeHh/mkfeFrFriVrZqhbZsxZoMVH34AZJU3vd8KouuQvGT9WprlaDVnx7HikV3k1XrqTbKf74pNrI+xVar2PpkWa3Cc1hWqvAclpUqPkesQdwu1iC+JkuD8Jr1jc3cc2u+JKuWrFrCw6za8sy3Taza44MSvPb3tzdi9nxZrdsi/8mqtZGm5lac/LHKyPoUbkQtLhzAW6bAzaIFwuIGQtuWj3QVFzcQb0TNKm4g1MBqF2pgaRRrEGsMDwnC2MF96AEgCMKlCCNmbbFfRwyMQlhwkFsqL5FVa6dVO6SvqY2wK78WGl0b5v8slln3Na9Ei7PlTZg5TmJiq/LW6r4zaoyND2Fau00tHdj2XTUkwf6Yc1c0U9+WY9UAgIUTY5ntB4vrcbZCh2fv68fUWKHWY8eJGjw+IcZEIxVAIKuWIKvWndhTJcgdli3VqrW1w05VQKdvgfxng42Oa3V6HK1UIiHAD4NHDjMJ9gGAsxol7o31QWh0KOSMYJ/sg0W4984+aO/gIJebBvMUltVgbN0F6FvbkDTRNJinsKwGYzQqAGDWtuU1Pn5/JPSB/niEoWFdthJL59yOC5fVeFqkoVh7CnL5eI+8L1SrlugpUK3amyybMx6rPz2OAbHWr3nOvneoywvSk1VrI8G9AnDoTIVRAQKh9Zm+YJLBMhUPSNMnSPHS3AmGQVI8aALAS3MnGFmmwkFx/wkVli+YaGTbittfSUkysm1ZGp98YAxTA6/x4YlDmBoIgvAsklSXTPI1419N7xGzTeHguX7XKWit3Eibn53aW1zBKlegO6zaxoZGr72JNQ3tuFjdgso6PQZEds74rta3IDYsAH6+PobXCds1TW0I8vdBcKDfzes0tiKslx+C/H3R3sFB29wOSfBNA6Cppb1zoP7pnJrGVsSEBhjaxedomtqMzhcfa2huR3Cgr5FGYXt7B4e2Dg5B/je/S+nbOuDv62M4J8Cfw7iEYI+8L6FhoW61as+XVuLwd2fok5twOlN+dgdGDh5AHSGA31bMWguWL67gqlmn2wdOgiAIgrAVW6sEffTleYwYGOWSpUGyagmCIAiPZ8TAKAyICYMi37riDk89OBKKfNdY1jRwEgRBEF7B7ClDUHyp1ur1yxEJ0S4JtqKBkyAIgvAannxwFD7aX2j1QOuKWScNnARBEITXEB4cCLlsoNV1bUckRDu9ji0NnARBEIRXkTSiHxqa9FZZtrOnDMHXTp510sBJEARBeB38htjW4sy8Tho4CYIgCK8jPDjQ6tq0tg6yNHASBEEQPRJrU05s3aKMBk6CIAiix2JtyknSiDjsPlxCAydBEARxa2NtyknSiH6orGmggZMgCIIgbCl04IwgIRo4CYIgiFti1jn73qFWl+yjgZMgCILo8bNOS4UOBsSGOsWupYGTIAiC8HruvyvBqkIHA2LCHLZraeAkCIIgvB4+5cTShtdyWYLDdi0NnARBEETPmHXKBuLrk10PigNiQ9HQpKeBkyAIgiBGDIxCUXmtVa+1NDOlgZMgCIK4JbBmDdOamSkNnARBEMQtgTVrmCMGRjkUXUsDJ0EQBNFzZpxOSjmhgZMgCIK4dQZPK+zasOAgu9NSaOAkCIIgehRJI/pZLMFnzWto4CQIgiBuCaxZw3RknZMGToIgCIKggZMgCIK4lRkQE2axdi0NnARBEATxE9bUrrV3cKWBkyAIguhx8LVru2L4wCgU08BJEARBENZhb4AQDZwEQRBEj8WRmrQ0cBIEQRC3FEkj4uyyYmngJAiCIG5Jhg+MgrKoigZOgiAIgrAGawKE7Cm9509d27NQqVRQqVSGn+VyOXUKQRCEGToDhBoxIDaUZpw9hezsbPj4+Bj9l5WVZfb1WVlZSE5ONvxHEARBdI2tAUQ0cHrBwClGoVBQxxAEQTiBpBH9UHyplgbOnoJGo2EOnNnZ2dBoNNRBBEEQFnBk+zAaOL18tpmSkoKUlBTDz+bs2oyMDHAcR51HEASBm2uYNHDeIggHR7lcbhTow5qJitsoMIggCML50MDpoahUKuTn55udcebn5xtFz7IGXOHrCYIgbkVckctJA6cXzDZTUlIgkUggkUgs2rUKhQIKhQJSqRSpqanUkQRB3NJYk8tJA2cPQWjFCi1XS3atXC4Hx3EoLS2lTiQIwmt44okn4OPj4xVaaeD0QMQ2rHCWKfy3SqWi1BSCILwelUrVZdyGp0GVgzx8tsnbtDy8Xcu/Jjs7m4KACMKD4b/cymQySCQSZGdnQ6VSQS6XQyaTGb02PT0dAJCamgqNRoOsrCykpKRALpcbCpqkpKQYLcOwjguvA3Qu66hUKshkMqSmphp9pgi/sPOfKzKZDHK5nPm6rKwsKBQKQ0ocf02pVNqlpq60CicAycnJkMlkyMjI8NybyhEeh0Qi4QBY9Z9EIqEOIwgPhv9b3bBhAyeTyYz+fuVyOfO1GRkZhs+BtLQ0ozb+Z/E5wuPC9xR/nojfk+M4LjU1lfnZcvLkSaPXpaSkmP0s2rBhQ5eazB2Xy+Um12JpdIT3d37vULsYsmo9cLZpS3EDc0USCILwLJYsWQKpVIrt27cbllwUCgUzyC8zMxNyuRxpaWkOOUpLlixBamoqcnNzjd5TGLGflZVl0HDy5Emo1WrDjJefufKv4z9r+Hzx0tJSw6w5PT3drsIsubm5SEtLE07mkJubS1YtYbutw9uyYitHaKvwD2l2djalnhCEhyOXy7F9+3bDv3m7Mzs72yQCXiKRGF7r6HvylqdUKjVa4uE/W/hjKSkphmMZGRkGS1alUhmdyw/o/DU3bNiAu+66y+zv0hOhgdOD4Nc0eFJTU836/JmZmYZvg/wslbUeQRCEZyD8Esx/KebTx8Q464uw8D2Fa5CsL+v8hhKsL+lSqdRordbce5jLLaeBk3AZYsu1qz+elJQUIxvlVvmmRxA9aSD1pKh4VrASP8gTNHB6xcApk8nM2rT8t0eZTGZYq8jKyqKBkyC8CGfMzoRrlfYikUig0WgsRrJKpVKoVCqTdUyhBnOzWmdp9RQoOMiD/oiE3z6tsWqsLcFHEN74t2DOxvRWxAMH/7OtwT/C62RmZjqsi/8cycrKMgyKKpUKycnJeOKJJwzHeJ382idLg/h3sUerN3yO0YzTQ5BIJEaRZF3NNnlSU1NN1k0IortIT0+36sORtwTFuX88WVlZRtfpKbv9KBQKpKenQy6XG3I5rf2SzPcb/2UiKioKUqnUKX/zqamphjiJu+66y/A+fK4p/x4ZGRnIz89Hfn6+IddSo9EYfo8NGzYY7qetWoXPwV133QWpVIqTJ0967s2kLCuCIJxBWlqa1fnHYOT+mbtOT8njTE1NNcmplMlkFvMcedRqtVG+ZUpKiuGYXC63O4+S4zju5MmTRjmaUqmUS0tL49RqNVODVCo1yrncvn273Vr51wvfX9wvnpbH6cPR5o0EQbhgxim27VQqFdOGKy0tNZpxiK/j7R9RfKRqWloa0tLSDJaoTCajNDI3UHRJjeJLasyeMsTsa9bvOoVlc8aTVUsQRPfCSmJXqVTIzMw0SrvKysoyCkrJyMhARkaG1xT8tgWJRGKU7E+4noamFgyIDXXqNSk4iCAItyGVSq2qQUqbsROeDA2cBEF0+yxMDG3GTjgLZVEVkkb0c+o1yaolCMJt5OfnG61f/n979xvT1n3vcfwbIAGFgBOHJvFS4jQjKizJqgaNPgipNpXqRtC7imqTlmRppaurgrM8vKJQHq9DXN37CMXk0VWaZnsG0jRQq0uuqsSZVDQSbSYzvUEZxqJml8SdzR8FLhF7QOyBY8w59vkdH/u8X9KkQV1DTk79Od/POed3kh/OLlJ4D2PnMhJriy2tSMXuUoITQO5td44yvuB58i0p8YexA2aYCEakzu3U9e9Q1QIw3bFjx+Sjjz5Ku9IMYNbEqRfBCUCJpqamF/4XP5/56NEjaWtrE6fTmfKxWoB5E+cT3edAqWoBKLHVMxWTb0lpa2uTpqYmpk8YLhCMyOGqCsPfl4nThsL9XvnjD8/KajSq/GeFenskdten9GdMdXfJ5JXL/MXmiVS3pDB1Qs00GZFanecvCU4TrEajEurtkdFjbuUBYYTYXZ+Eenukpu+qlDgcpoR01IDtkm47v9zRKYvjfgn3e9khASTMPJ7XfeEPwanY4rhfvr50UZanp/Pmd/5Ld5ccOH9Byk+eKpjtXOJwSHVHp4SvefPq78Kuvv32203Pko1PoYDpg0QGFwaJcI4zK5HhIXE2t4ir3SOjx9x58fsuT0+L61PPpu+HenvkwM8uyHJoOjEdlp88Jc7mFiUTeviaV9f7a9nOzuYWCfd7JXzNK0d/+St2Tgt4++23U34/+VFhqe7lBMwwEYxIQ52L4DRTdUdn3gV9+clTUnrkyKbvxyvOeEA9i0Zl8splcbV7DP8zTnV3rYd3m8fw7VzZ2Cgxn48d0yK0PEsz/jg9HokHo40GZqXWvX+b14Tl/XMnCE6kD05Xu2fLf/bal3cSXxc7HBIZHjI0OOOh+eqNm0rOr5ZWH5HFcb+sRqOmnL9FZvbt2yf19fXS1NQkH374IaEJRcGpLRQrd+8iOJGZVLWpkecLp7q7ZNHvVxaaIiJlzyfppXG/VJ5p5C/VZPGnmgBWkUkoEpxQdzSXdK6x4VFwy9fGfD5ZHPfLgfMXtg1NPe+bbHceXfAEQB0tF/3MzC3oXqOW4ERWam/c1PV6V7tHwv1eqTzTmPaiIL3vu9EzE+5LBZAHB/aB2W0v+ll/TWZPTSE4kRE9VWhlY2PiXOlUd1fKC5Qyed9ky6H1armY85uArU0En8iV915PP3E+npfWN2sITtPrgKQb8Rf9fkMC4P9+82uJDA9tmr6Sv5fqNVrCLpf3OVZ3dErM55PJK5flxG9/Z/h2XvT7pcThyKt7VAHkH4Izm6OaSxc3fR3q/ceFEXrOzaWanJLDIvl7qV6zHWdzi0x1d+X0qtOavqvy4N13JNTbo/mKXa3bOTI8xEVBgM2NjE1vextKNuc3RUR2rPHgO1u59/r35aXzF/LuHlQtU+nEpYvy2pd3tqyBARS+voH729a0g7cnpdbtzHg5PpbcsxlXm0fmfvNrUxZ4N9M3/V5xtXsITcDGtC6hl+0atgSn3YKz3SPO5hb5Oqn+zGeh3h55Fo3qWo0IQOG59Ydpeate/fKnVLUAgILwyWdfycc/fyPta0bG1i+QbKrPvJ3i4qA8MxqYVfje4bT//P1zJ5StxAEA2QgEI1J3ZP+2r9NyqwrBabC+gfs5/fm17v3KwitdMI4GZp8/SeAQOwEAy7k1FsxowXaC0yTZHq0AAIwzM7cgItuvTavlVhWCU4G//S2q6XFJhWbyryty+vRpdgAAljN456GmadOImpbgzMDevQ5parLfxFmp8NwqAKieNrXeqqIFt6MAAAp+2hy8PSmtZ48b8jOZOHWiqgUAawgEI5qmTZH1RQ8Ov7SH4MwFqloAsM60ud19myLrFwVt95gxPahqAQB5R08YjgbCWS14wMSZJapaAHY2M7cgg3ce5vS2vNjSiowGwpqmTa0LIxCcClHVArCr2NKKXP/iQc5XENPyBJQ4FQsjUNUCADSFpp7AUiVe0Wq6IEjjrSpMnIpR1QKwo3hoTgQjhqy+k4mZuQVdixhovVWF4FSMqhaAHUOz9ezxxOSWq6r2+hcPNIfmaGBWKnaXKvldqWoBAGlD8616d+LBzzNzC1KRg+BMDu/tjIwF5YNz31PyuzBx6kRVC8AuBm9PSq17fyI0RdYXEmh9sybnv0f60DT2vk2CM0tUtQDsEpoiYuj9j2b8HnpuVckUVS0A4IWJTURMnyy3Ck09v8ennz+QD/5J7XM5mTh1oqoFUMgCwYhhj98yOzTja9catSYtwWkQqloAhRyat8aCeRmaItrXriU4mTiZOAFkLba0YlrwpNM3cF9q3ft1n1uNX3VrBoKTiZOJEyA0c74qUPx3aD17XPPVs3Hxezb1/nsEJwAg4ynvynuvb3uP5GhgVsmqQSNj0zIaCGv6HVIF7shY0NRJmeDUiaoWQKGFpp6FBYxciWdmbkGuf/FAGupcGQdfLiZlglMnqloAhRSaG1cF0hJ0tQbUofFHk4lIRlNmpqFPcAIATAtNkexXDYpftSsi8v65E1kFnt7VhAjOHKKqBZDvBm9PyuGqClNCZzQwK6OBsIiIHK6qyDow47+/SO5WNSI4daKqBZDvoSmS2apAh6sqpG/g/qavNy42MDO3IDOP5zf9O7Xu/YaE5capdebxfE6vACY4AcAmsg2d5LANBCMyv7SyISSdSpfpGxmbtsSqRgSnTlS1APLRzNyC4asCmXl+MT4p5zo0Cc4MUNUCgLmyqZdV4OkoAGAD8XORsQ3Vaj7oG7gvFeW7LBOaTJwZoKoFkK9azx6XTz9/YIm6czvZLMFHcFoMVS2AQpg6zV40QI9sluAzA1UtANhw6rTqlPnJZ1+JiMjHP3/DsuHOxKkTVS0Apk7jXf/8zzK/tGzZKZPgzAJVLYBCmTqtcK5z8PakBKafSFO9WxrqDuXF9iM4AYCp01SxpZVEXdxQ57LUFbMEpwJUtQCYOjMLy1t/mE4sydd69vim5foIzgJGVQvAjlNn38B9OVxVIbVup6bbQ9aX91uQieCTxM94q96dd9MlwQkAyHjqrHU7ZWLDo8HSBvPzBeDz4Z5RglMxqloAdps646sN1WmcNgsdwakTVS2AQps6B29PygfnvrflayaCEWmoc7GxnmMBBACw+dSZ/AxNMHEaiqoWQKFpqnfL4O3JLS/cCQQj0lR/hA1FcGaGqhZAoWmoOySfjH0lrZI6OOeXlvP21hGCk4mTiRNATqZOEJxMnEychpqZW5DBOw9t+Wevde+3/Nqi+b19naZs3+2mThCcgLHB+XhBGupcebPWppFGOahSyszVfZg6CU4lqGqBFycVqDwwCZv6d8nUSXAajqoWqczPz8vDhw8lNkNlCWM9fWruRTlMnQQnYIqKigo5ffo00xcMNzFw3/QGYePUGV81CARnxqhqARS6jVMnqwYRnFmjqkUqVLVQd7Ce/mN6NRqV8DWvhPu9UnvjplSeaTR86gTBCRiOqhaqpKtqF8f9MvVxl5QeMX5Vn/jUGVtaYdUggjPboz+qWgDWEBkeEmdzi7jaPTJ6zG3oe8enzsrdu1g1iODMDlUtUqGqhbqD9a0/pqs7OpX+7KZ6t3zy2VcF+UxNghPIMapaqGL2VbXJU2fr2eP8JRCc2R79UdUCsI90z+kkOKEJVS22cu/ePapaKDhY52Oa4AQKFFUtVMhlVQuC06CjP6paACA4oRlVLbZCVQs1B+tbf0zH7vo2fb3o9//jv1kDFkIAwQkoRVULFdJVtROXLm76OtTbk/j/DY+CbDyC0ypHf1S1YOKENSZOwpHgzAtUtWDihFUmTuRGEZsAAAAmTmWoarEVqlqo+czhY5rgzHNUtdgKVS1UoKq1HqpaAACYONXiq8T5AAAGq0lEQVShqkUqFbt3yef/83uJzZSxMWDwZw4f0wRnnqOqRSp1bqfcchy25b4BtahqCU4mTiZOAAqsPX0qq/PzUlRWJsV7ykV2cCaO4GTiZOK0uKmpoIyMPGFDwOCDdW0f00uBP8u3I/8tu2vrZO+PfiRFeyrYeAQnYG1Hj7qpamE4rVXtot8v33ivyv53/lkqftBAcBKcVjr6o6oFYMEP8717pezoK7Kzqkp2lPDRTnBaCFUttkJVCzUH69o+pp0/flecP36XDUZwAvmDqhYqcFUtwZn3Zmf/SlULAAQntDp06CBVLVKiqoUKWqvatZUVWY1GpWjXTimuqBQp4nYUghOwOKpaqKC1qo367sj//uu/SGXDG/LKv/+HlFZXs/EITmugqgVgRcV79kjZ0Vdk18GDUsRVtQSnlVDVYitUtVBBa1Vb5nbLdy7/QnZWvSTFjko2HMEJWB9VLVTQWtXuPHhIqn7yUzYYwWk9VLUAQHBCB6pacFAFc/erVTYCwQlwUAVoNRHTuADCs2fy/5GI7CgulpK9DpGiYjYewclUkUtUtYC1PZ36i/zp7bek9OVqefW/rkvZd7/LRiE4mSpyiaqWgyrkar/SVtUWlZVJ2dFXpPTQQSnauZMNR3ACHFTBnrRWtSV798l3Lv9CisvLpcS5jw1HcDJV5BpVLfsGLD5xlpdzOwrByVRhJVS17Buw9sSpQqi3RxxnGqXyTKOpP3equ0tWo1Gp6btqyb8TVgEGAKQU7vdK9K5PyXuvRqMS6u2R0WNuiSX9jJc7OmVx3C/hfi8TZyGgqgX7Bszdr9JXtX/84VlxnGmUo7/8Vd78mRbH/TL1cZeUHjmSOpgcDqnu6JSp7i5xNrds+TqCM09Q1YJ9A2ZKV9VGhodkeXpaXJ96Et8L9fbIgZ9dkOXQdGJaLD95SpzNLYZMieFr3qzfLzI8JM7mFnG1e2T0mDvla5zNLRLu90r4mtdyBwUEJwDkqcjwkJSfPLVpIovXm/FwehaNyuSVy+Jq90h1R2dWP2+qu2s9qNs8Wb2P1t+jsrFRYj6f5bY7wakTVS3YN2DufrWaNjhd7Z6U33/tyzuJr4sdDokMD2UVnPHQfPXGTSlxOEz5s5dWH5HFcb+sRqOm/UyCUwGqWrBvwEyZXFWbqkZdnp7OKjQX/X5TQ1NEpOz5JL007jf9yl6Ck6mCiRNAQvJ5xYZHwS1fG/P5ZHHcLwfOX9g2NPW8rxa7T56y5PYjOJkqmDg5qIKl9yvjn45Se+Omrte72j0S7vdK5ZnGtBcF6X3f7TyLRglOgIMqQB8VCyDoqT0rGxsT50anurteuBgp0/fVYjm0Xi8XW+j8JsHJVKEZVS1gPZVnGrM6d6lHdUenxHw+mbxyWU789ndZvVfyggeLfn/K8F30+6XE4ZByi1W2BCdThbb/QKlqOahCjvarrataZ3NLYnk6My7aqem7Kg/efUdCvT1ZXaE7cenipq9DvT2J/7/xvGhkeMhSFwXF7VhbW1tj19Sub+C+XHnPfsE5+jw4G+oOsROwb8BC+9W9178vL52/kPU9mlYTu+uTiUsX5bUv77ByEFNFfqKqBazJ1eaR8DWvuNo8lrrXMVvf9HvF1e6xXGgSnBmgqgUHVTB3v0p/Va2r3SPLoWn5+tLFrM89WkWot0eeRaNZr1BEcAIcVMGGtFxVm08LvGtR3dEp0mHd34/gZKrQhKoWAAhOpgodqGo5qEKu9qtVNgLBCXBQBWilYgEEEJxMFSagqt3e4aoK6RvgQw4gOMFUIVS1WrS+WcNGAGygiE0AAAATpzJUtQBAcEIHqloAsDeqWgAAmDjVoaoFAIITOlDVAoC9UdUCAMDEqQ5VLQAQnNCBqhYA7I2qFgAAJk51qGoBgOCEDlS1AGBvVLUAADBxqkNVCwAEJ3SgqgUAghNMnEycAEBwMnEycQKA8bg4CAAAJk51qGoBgOCEDlS1AGBvVLUAADBxqkNVCwAEJ3SgqgUAe6OqBQCAiVMdqloAIDihA1UtANgbVS0AAEyc6lDVAgDBCR2oagHA3qhqAQBg4lSHqhYACE7oQFULAPZGVQsAABOnOlS1AEBwQgeqWgCwN6paAACYONWhqgUAghM6UNUCAMEJvCC2tCITwUji60AwInVuJxsGAMHJJtDnUfAb+bf/DNriz1pzcNemr2vfrGEHAGB7O9bW1tbYDAAAaMNVtQAAEJwAABCcAAAQnAAAEJwAABCcAACA4AQAgOAEAIDgBACA4AQAgOAEAIDgBAAABCcAADr8HfpNUax6lBayAAAAAElFTkSuQmCC"/>
%% Cell type:markdown id:827d9797 tags:
The example above is trivial but when a Faust contain many factors including very sparse matrices, the benefit can be quite significant.
The example above is trivial but when a Faust contains many factors including very sparse matrices, the benefit can be quite significant.
The reason why is that ``pruneout`` repeats the removal operation iteratively until there is no zero column/row eligible to this operation anymore.
Note that of course the resulting Faust gives the same full matrix (i.e. ``Faust.toarray`` must give the same result after a ``Faust.pruneout``).
In the next section of code I show with a figure how a sliced DFT Faust can benefit of ``Faust.pruneout`` in terms of computation time.
In the next section of code we show with a figure how a sliced DFT Faust can benefit of ``Faust.pruneout`` in terms of computation time.
%% Cell type:code id:270523bf tags:
``` python
%matplotlib inline
import pyfaust
from timeit import timeit
from numpy.random import rand
from numpy.random import permutation
import numpy
from numpy import median
import matplotlib.pyplot as plt
times_Fsub = []
times_Fpr = []
n_exps = numpy.arange(5,15)
N = 3
for n_exp in n_exps:
n = 2**n_exp
F = pyfaust.dft(n)
k = 10
I = permutation(n)[:k]
I.sort()
I = list(I)
Fsubset = F[:,I]
Fpruned = Fsubset.pruneout()
x = rand(n)
xsubset = x[I]
time = median(timeit(lambda: Fsubset@xsubset, number=N))
times_Fsub += [time]
time = median(timeit(lambda: Fpruned@xsubset, number=N))
times_Fpr += [time]
plt.rc("font", family="serif", size=18)
plt.semilogy(n_exps, times_Fsub, marker='+', label="Fsubset*xsubset")
plt.semilogy(n_exps, times_Fpr, marker='o', label="Fpruned*xsubset")
plt.grid(True)
plt.legend()
plt.xlabel("n_exps")
plt.ylabel("time (s)")
plt.tight_layout()
plt.show()
```
%% Cell type:markdown id:e94876cc tags:
``Faust.pruneout`` is referenced [here](https://faustgrp.gitlabpages.inria.fr/faust/last-doc/html/classpyfaust_1_1Faust.html#a2a2940af5acbf0531affd10cf03076f1) in the API documentation.
**Note**: this notebook was executed using the following pyfaust version:
%% Cell type:code id:969e808c tags:
``` python
import pyfaust
pyfaust.version()
```
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment