diff --git a/.gitignore b/.gitignore index 4d3d00e..d0bbb6d 100644 --- a/.gitignore +++ b/.gitignore @@ -10,3 +10,4 @@ coverage docs/build/ env node_modules +example/*.pdf diff --git a/Project.toml b/Project.toml index 2a25660..935d280 100644 --- a/Project.toml +++ b/Project.toml @@ -1,9 +1,6 @@ name = "InclusiveScans" uuid = "ef74772d-374f-5edb-a43a-0dbb31824c9d" -authors = [ - "Uwe Hernandez Acosta", - " Simeon Ehrig", -] +authors = ["Uwe Hernandez Acosta", "Simeon Ehrig"] version = "0.1.0" [compat] diff --git a/examples/ConstrainedGaussians/Project.toml b/examples/ConstrainedGaussians/Project.toml new file mode 100644 index 0000000..cea9c26 --- /dev/null +++ b/examples/ConstrainedGaussians/Project.toml @@ -0,0 +1,18 @@ +name = "ConstrainedGaussians" +uuid = "72cb5a96-92ad-43ae-9774-2d852a9e5056" +authors = ["Uwe Hernandez Acosta "] +version = "0.1.0" + +[deps] +Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" +DomainSets = "5b8099bc-c8ec-5219-889f-1d9e522a28bf" +QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc" +Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" +StatsPlots = "f3b207a7-027a-5e70-b257-86293d7955fd" + +[compat] +Distributions = "0.25.116" +DomainSets = "0.7.14" +QuadGK = "2.11.1" +Random = "1.11.0" +StatsPlots = "0.15.7" diff --git a/examples/ConstrainedGaussians/src/ConstrainedGaussians.jl b/examples/ConstrainedGaussians/src/ConstrainedGaussians.jl new file mode 100644 index 0000000..f1196a4 --- /dev/null +++ b/examples/ConstrainedGaussians/src/ConstrainedGaussians.jl @@ -0,0 +1,22 @@ +module ConstrainedGaussians + + +using Distributions +using DomainSets +using Random +using QuadGK +using StatsPlots +#using InclusiveScans + +export ConstrainedGaussian1D +export target_dist, filter, norm, normalize +export endpoints +export generate_CPU +export plot_compare + +include("interface.jl") +include("impl.jl") +include("generate.jl") +include("plot_compare.jl") + +end diff --git a/examples/ConstrainedGaussians/src/generate.jl b/examples/ConstrainedGaussians/src/generate.jl new file mode 100644 index 0000000..915ed68 --- /dev/null +++ b/examples/ConstrainedGaussians/src/generate.jl @@ -0,0 +1,36 @@ + +function generate_batch_CPU(rng::AbstractRNG, batch_size::Int, dist::ConstrainedGaussian1D) + rand_args = rand(rng, Uniform(DomainSets.endpoints(dist.domain)...), batch_size) + + # opt this block out to be evaluated on a single element + # and broadcast over everything + # ------- + rand_vals = target_dist.(dist, rand_args) + rel_rand_vals = normalize(dist, rand_vals) + rand_probs = rand(rng, batch_size) + mask = filter.(rel_rand_vals, rand_probs) + # ------- + + accepted_args = rand_args[mask] # this one with InclusiveScan.jl + + return accepted_args +end + +function generate_CPU( + rng::AbstractRNG, + dist::ConstrainedGaussian1D{T}, + N::Int; + batch_size::Int = 1000, +) where {T} + accepted_args = T[] + sizehint!(accepted_args, N + batch_size - 1) + + nrun = 0 + while nrun <= N + accepted_args_batch = generate_batch_CPU(rng, batch_size, dist) + append!(accepted_args, accepted_args_batch) + nrun += length(accepted_args_batch) + end + + return accepted_args +end diff --git a/examples/ConstrainedGaussians/src/impl.jl b/examples/ConstrainedGaussians/src/impl.jl new file mode 100644 index 0000000..d854136 --- /dev/null +++ b/examples/ConstrainedGaussians/src/impl.jl @@ -0,0 +1,39 @@ + +# constrained Gaussian + +struct ConstrainedGaussian1D{T<:Real,DIST,DOM<:AbstractInterval} <: AbstractTestDist{T} + mu::T + sig::T + dist::DIST + domain::DOM + + function ConstrainedGaussian1D( + mu::T, + sig::T, + domain::DOM = Interval(-5 * sig, 5 * sig), + ) where {T<:Real,DOM<:AbstractInterval} + dist = Normal(mu, sig) + return new{T,typeof(dist),DOM}(mu, sig, dist, domain) + end +end + +function ConstrainedGaussian1D(mu::T, sig::T, domain::D) where {T<:Real,D<:Tuple} + return ConstrainedGaussian1D(mu, sig, Interval(domain...)) +end + +endpoints(dist::ConstrainedGaussian1D) = DomainSets.endpoints(dist.domain) + +function norm(d::ConstrainedGaussian1D{T}) where {T} + mu = d.mu + dom = d.domain + + if mu in dom + return target_dist(d, mu) + end + + if mu <= minimum(dom) + return target_dist(d, minimum(dom)) + end + + return target_dist(d, maximum(dom)) +end diff --git a/examples/ConstrainedGaussians/src/interface.jl b/examples/ConstrainedGaussians/src/interface.jl new file mode 100644 index 0000000..0a3ced9 --- /dev/null +++ b/examples/ConstrainedGaussians/src/interface.jl @@ -0,0 +1,16 @@ + + +# general abstract test-distribution + +abstract type AbstractTestDist{T} end + +Base.broadcastable(dist::AbstractTestDist) = Ref(dist) +target_dist(d::AbstractTestDist{T}, x) where {T} = x in d.domain ? pdf(d.dist, x) : zero(T) + +function normalize(d, x) + return x / norm(d) +end + +function filter(rel_rand_val, rand_prob) + return rel_rand_val >= rand_prob +end diff --git a/examples/ConstrainedGaussians/src/plot_compare.jl b/examples/ConstrainedGaussians/src/plot_compare.jl new file mode 100644 index 0000000..41b0197 --- /dev/null +++ b/examples/ConstrainedGaussians/src/plot_compare.jl @@ -0,0 +1,24 @@ +function plot_compare(samples, dist; kwargs...) + P = histogram( + samples; + label = "samples", + xlabel = "x", + ylabel = "normalized event count", + nbins = 100, + normalize = :pdf, + opacity = 0.5, + kwargs..., + ) + + tot_weight, _ = quadgk(x -> target_dist(dist, x), DomainSets.endpoints(dist.domain)...) + plot!( + P, + range(endpoints(dist)...; length = 100), + x -> target_dist(dist, x) / tot_weight; + label = "normalized target dist.", + line = (2, :black, :dash), + alpha = 0.5, + ) + + return P +end diff --git a/examples/Project.toml b/examples/Project.toml new file mode 100644 index 0000000..106f129 --- /dev/null +++ b/examples/Project.toml @@ -0,0 +1,5 @@ +[deps] +BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf" +ConstrainedGaussians = "72cb5a96-92ad-43ae-9774-2d852a9e5056" +InclusiveScans = "ef74772d-374f-5edb-a43a-0dbb31824c9d" +Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" diff --git a/examples/README.md b/examples/README.md new file mode 100644 index 0000000..464b39b --- /dev/null +++ b/examples/README.md @@ -0,0 +1,91 @@ +# Example: Constrained Gaussian + +This repository demonstrates how to generate samples from a **constrained Gaussian distribution**, +a non-trivial extension of the Gaussian distribution confined to a specified interval. +The distribution and sampling functionality (CPU-based) are implemented in the `ConstrainedGaussians.jl` package. + +## Features + +- sampling of constrained Gaussian distributions. +- Comparison of generated samples with the target distribution using histograms. + +--- + +## Setup Instructions + +### 1. Clone the Repository + +Ensure you have cloned the repository containing this example. + +```bash +git clone https://github.com/QEDjl-project/InclusiveScans.jl +cd /example +``` + +### 2. Initialize the Example Environment + +Since `InclusiveScans.jl` and `ConstrainedGaussians.jl` are not registered Julia packages, +you need to initialize the environment manually. Run the following command inside the `example` directory: + +```bash +julia --project=. init.jl +``` + +This script will: + +- Link the source code of `InclusiveScans.jl` and `ConstrainedGaussians.jl` to the example environment. +- Install all other necessary dependencies. +- Instantiate the example environment + +--- + +## Running the Example + +### 1. Execute the Example Script + +To generate samples and visualize the constrained Gaussian distribution, execute the following command inside the `example` directory: + +```bash +julia --project example.jl +``` + +### 2. Output + +- The script generates **1 million samples** from a constrained Gaussian distribution. +- A histogram comparing the generated samples with the target distribution is plotted. +- The histogram is saved as `example_plot_compare.pdf` in the `example` directory. + +### 3. Interactive Exploration (Optional) + +For interactive experimentation, open the Jupyter notebook `example.iypnb`: + +1. Ensure you have the `IJulia` package installed: + + ```julia + using Pkg + Pkg.add("IJulia") + ``` + +2. Start Jupyter and open the notebook: + + ```bash + jupyter notebook + ``` + +This allows you to explore the constrained Gaussian example interactively using Julia's kernel. + +--- + +## Requirements + +Ensure you have the following tools installed: + +- **Julia**: Version 1.6 or later is recommended. +- **Jupyter** (optional): For interactive notebooks. + +--- + +## Troubleshooting + +- **Dependency Issues**: If you encounter problems during initialization, verify that the `init.jl` script executed successfully and installed all required packages. +- **Missing `IJulia`**: Ensure the `IJulia` package is installed to enable Julia kernels in Jupyter. diff --git a/examples/example.ipynb b/examples/example.ipynb new file mode 100644 index 0000000..dad2603 --- /dev/null +++ b/examples/example.ipynb @@ -0,0 +1,355 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "49a6672b-2bef-4997-b8c8-b16bcfaacadd", + "metadata": {}, + "source": [ + "# Example run" + ] + }, + { + "cell_type": "markdown", + "id": "fddf47fa-6939-41c2-9109-e87137b5e91c", + "metadata": {}, + "source": [ + "## Some includes" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "d4c9b18f-bbf5-48d3-ba35-38590f9c4665", + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "\u001b[36m\u001b[1m[ \u001b[22m\u001b[39m\u001b[36m\u001b[1mInfo: \u001b[22m\u001b[39mPrecompiling ConstrainedGaussians [72cb5a96-92ad-43ae-9774-2d852a9e5056] (cache misses: include_dependency fsize change (4))\n" + ] + } + ], + "source": [ + "using Random\n", + "RNG = MersenneTwister(1234)\n", + "using BenchmarkTools\n", + "\n", + "using ConstrainedGaussians" + ] + }, + { + "cell_type": "markdown", + "id": "49115263-017c-4039-9203-9c80fec0906a", + "metadata": {}, + "source": [ + "## Setup" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "id": "f569728a-e50b-4108-8b02-e83f7dc1b186", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "1000000" + ] + }, + "execution_count": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "mu = rand(RNG)\n", + "sig = rand(RNG)\n", + "\n", + "dom = (mu+0.5*sig,mu+3.0*sig)\n", + "d = ConstrainedGaussian1D(mu,sig,dom)\n", + "N = Int(1e6)" + ] + }, + { + "cell_type": "markdown", + "id": "b4fc8c04-987b-408b-abb4-aa32eaf3badb", + "metadata": {}, + "source": [ + "## Generate samples" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "id": "8aaf102e-7a7a-4651-86d5-189b3bf2a612", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "1000288-element Vector{Float64}:\n", + " 1.0742340558435985\n", + " 1.0730781173125767\n", + " 1.9896250113762497\n", + " 3.105113427930066\n", + " 1.8574666488244633\n", + " 1.8159372698439835\n", + " 1.109548243172165\n", + " 1.5770354160050553\n", + " 1.3590575226781938\n", + " 1.1914585945695642\n", + " 1.8759880573744998\n", + " 1.7139305526129174\n", + " 3.381209129317833\n", + " \u22ee\n", + " 1.4397024632246285\n", + " 1.0885387589109017\n", + " 1.5760739229042762\n", + " 1.3555945603290158\n", + " 1.3664612587409264\n", + " 1.1962789622099628\n", + " 2.165842345107823\n", + " 1.2487352157877882\n", + " 1.1257704298484303\n", + " 1.2414722995661573\n", + " 1.1095129273547224\n", + " 1.6673401880835255" + ] + }, + "execution_count": 4, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "samples = generate_CPU(RNG,d,N ; batch_size=1000) # not enough elements????" + ] + }, + { + "cell_type": "markdown", + "id": "80686f69-9c4d-4aa4-b592-b49d09683ce6", + "metadata": {}, + "source": [ + "## Plot samples in histogram and compare to target distribution" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "id": "8054cebb-d5c2-48fc-93f6-9b1383272267", + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdd1gU194H8DNbgF16lSYovVsAFbFhUFEgdlQUS2JJctVc9Y2xJqbnJsYYo4ldVLBrvAawo2JFQUEBRYo06R2WZdn2/rG5G4JLUWGH8v08efLsnpmd/e7I8mNmzpxDSaVSAgAA0FMx6A4AAABAJxRCAADo0VAIAQCgR0MhBACAHg2FEAAAejQUQgAA6NFQCAEAoEdDIQQAgB4NhRAAAHo0FEIAAOjRumohrKur27BhQ9vXF4vFHReme8AuapVEIqE7QmcnlUqxl1qF71qrlLyLumohLCsrO3ToUNvXr6ur67gw3QN2Uav4fD7G5m2ZSCRqaGigO0Vnh+9aq5S8i7pqIQQAAGgXKIQAANCjoRACAECPhkIIAAA9GgohAAD0aCiEAADQo6EQAgBAj9ZzCyGfz09ISBAIBHQHAQAAOvXcQhgXF3f27Nl9+/ZVVlbSnQUAAGjDojsAbfr16/fkyZPi4uK9e/fOnDnT3Nyc7kQAPVFGRkZP+2O0rq6Oy+XSnaJTU7iLmExmv379KIpq97ejuuiQUbm5ud7e3jk5OW1cv6amRlNTs0mjQCA4ffr08+fPWSzWxIkTXV1d2ztmV6JwF0FjPB6Py+V2xPew2xAKhWKxWE1Nre0vsbS01NHRYbPZHZeqs5FKpfgpapnCXZScnJyQkGBvb9/ub9dzjwgJIaqqqjNnzoyKioqLiztz5kxlZeXw4cPpDgXQs4jF4qioKDMzM7qDQGfn7OzcQYNx99xrhDIMBiMgIGDcuHEURV29evXcuXMYOx8AoEfp6YVQxsvLa8aMGWw2++HDh0ePHsXw+QAAPQcK4V/s7e3nz5+vrq6elpa2b9++6upquhMBAIAyoBASSUOlsCxeWBZvpFY4b4qHLodfkJWw59cvCp5fFZbFC8seEYmI7owAANBRUAiJsOyRsDRezMsX8/K12LVzJnqY6KtUlBYcCA3LSX8syL8i5hfSnREAADpKj+41KsfStlcxHil7rEbIQutJp06dSk1NPX61IMCTOFvRmw4AADoQjggVYLPZM2bMGDhwYENDw6nI+y+yculOBADdxJ9//pmUlER3CvgHFELFGAxGYGDgiBEjpFIpn8+nOw4AdBOhoaH37t2jOwX8A06NNouiqNGjR7sbp2lY2tGdBQCUh8/nUxTVeHwcsVgsEAheHfRLKpXW1tbKh2Ti8Xjq6uqvbk1VVZXBaPaog8/nczicxi319fUSiQTDsCkNjghbwWYx6Y4AAEqSk5Pj4eHh7Ozs5OQ0cOBAWaOfn1+fPn1cXV2dnJwePHggaxwzZsyGDRscHR379OkTEBDw9OlTT09PW1vbAQMGFBcXE0Lu3r07ePDgefPmubq6Ghsb7969u8l7SaXSr7/+2sLCwtnZ2cPDIyUlhRBSVlY2atQoOzu7/v37W1tbY3oc5cARYeskglJCKf6LgWKoMtT0lZwHoHsT83KJRKiMd2KqMrn/GNpt+/btPj4+P/74IyEkPz9f1vjzzz87OjoSQk6cOLF48eJHjx4RQmpqau7du/fo0SOKotzd3YODgy9cuGBkZBQUFLRz587PPvtMJBLdv3//gw8+OHjw4PPnzwcNGjRixAgHBwf5ex0+fDgyMvLJkyfa2tpHjx6dM2fOw4cPQ0NDLSwsrl+/TggpKipisfArWhmwl1vBVO8tLLmvcJFUKpEKqzRcP1VyJIBuTdpQdEuqlEJIMTmcvtMbt5iYmBw4cKB///5+fn6mpqayRl1d3R07duTn5wsEgidPnjQ0NKioqBBCFi1aJDul6eXlpa+v36tXL0LIqFGjbt++LXuhnp7e/PnzCSF2dnZ+fn5RUVGNC+HRo0eHDh0qO8Q0NDRMSUkpLS01MTH5/fff9+/fP2HCBGNjYyXsBCAohK1SNRvX3CKpWMBL3qLMMAA9AMWxmkXXey9fvlxLSys8PHzRokVTp04NDQ0tKCgYNGjQBx980L9/f9ncGnV1dbJCqKWlJXuViopK48fyMRp1dXXlUygYGBiUlZU1fq+ioiKJRFJbWyt7Om/ePJFINGvWLAaDcfz48eXLl/v4+Jw4caLJ5UPoCCiEb+L69euFhYWB/s3WSADoiphM5vvvv//++++Xl5e7urreuXMnOzvb3d39s88+I4TEx8e/1tZevnxZWVmpo6NDCElKSpo3b17jpQ4ODra2tl988UWTV82cOXPmzJm1tbVDhgyJjIycNm3a230maB0K4ZvIysrKysoqKy2e5tGgQXcYAGgvv/76q76+vq2tbXZ2dkNDQ9++fVVUVO7duxcVFcVms7///vvX2hqbzV6yZMny5ctv3LiRlpY2ffo/TsOuW7fOx8dHT09vyJAhRUVFMTExmzdvPnDgAEVRjo6OpaWlZWVldnbosq4MKIRvYtq0aYcPHy4qzD/838TFTtXysyIA0KU5OTmdOHHi2LFjxsbGUVFR5ubm5ubmv/zyy969e3V1dbdu3RoeHq6qqkoICQkJsbL6a9CpsWPHamtryx7369dPfjLTwsJi3rx5W7ZsMTQ0vHXrloaGBiEkMDBQ1vXGxcXl3r17O3bsiI6ONjAwGDduHCHEwcEhPDz8zJkz+vr6R44ccXNzU/5O6IF69Az1MoKCaxTFkA+x1kZ8Pv/wodCshNNG9gFz587V09N7rZd3QpihvlWYob5VbzBDvbm5eWxsbPebmPfmzZsffvghBpFpR87OzidPnnRycmr3LeM+wjfE4XBC5gSbG2tXVlYePHiwyWVwAADoKlAI35yamtoM//6WlpZVVVWhoaGlpaV0JwKAzmL48OE4HOwqUAjfiqoKa/bs2VZWVjU1NaGhoSUlJXQnAoCuJDEx8eTJk4QQHo/37bffvvF2JBLJf/7zH/mdG53Ktm3bKioqCCGnT5+WDUfQ2aAQvi0VFZXg4GA7O7va2trQ0NCioiK6EwFAlxEfH3/kyBFCSH19fVhY2BtvRywWr1mzpkkhlEgkHh4eNTU1b5uyGV988UVbMn/55Zeyi0c3btx4/vy5wnWSkpImT57czvnaDIWwHbBYrBkzZtjb2/N4vEOHDqEWAnRpYrFYNl5oY1KpVGFXAIFA0PgrX15ezuPxZI95PF6TIiQWi4uKikQi0avb0dfXlw032oKGhoaqqqomjbW1tfK78l8VHx8vFosbt5SVldXV1b26Zm1tbXl5uexxZWVlfX39q+tUVVUJhX8P+pOdnd3Crzs+n98k7bZt22bMmCF7LJVKCwoKqqurZU95PF5ycnJzm+poKIRvjqIYRCLiJf/CS/6l/tn2Cc6FltzMyqzoPT8syrzxFS/5F0HBVbozAsBrWLhw4ccffzxw4EBvb+++ffumpaXJ2nft2mVsbDx8+PA+ffpcuXJF1ujs7Lxu3To7O7vx48fHxsYOGjRo7ty5Xl5eZmZm4eHhX3/9tbu7e+/evT///HPZ+tu2bbOysvLz8zMxMXn1Pvri4mLZPRgPHjzw+B9bW1sfHx9CSHV1dXBwsL29vaenp4+PT2FhISFELBYvWbLExsZm0KBBGzdufPXjLFiwgBDi4+Pj4eERFxf38OFDW1vbUaNG2dvbT5gwQVaEcnJyLCwsVq5c6eLi8sknn1RXV/v7+7u5uQ0YMGD16tX29vayTSUkJAwYMGDo0KF9+vRZvXo1IeTMmTPnzp3bunWrh4fHq+++fv16CwsLb2/vDz/8UH5vwpw5c/bv308IiY6OtrS09Pf39/T0lJXGBQsWyEY89/DwUFinOxTuI3wLDDbX+d9E+vcfd7Nt5pw6/d/U5+kn70iC37UzqEdXUoDXk5OTc+nSJYlE0vaXsFiswMBAQ0PDVxeJxeLTp09XVlYqfCGTyZwwYYKJiYm8hcfjPXjwICYmRltbe/ny5Zs3b961a1diYuLq1avj4uJsbW0jIiJmzpyZkZGhra1dVVWVnp6ekZHBYrFu37794MGDzz///NChQ1FRUUFBQWvXrn327FlOTo6Dg8OyZcsMDAwmTZq0dOlSBoNRUVExcODAiRMn9u/fX/7W8qlPPT094+LiCCHl5eVDhw6dOXMmIWTdunWamprp6elMJvPzzz9fvXr1oUOHjh49evfu3bS0NE1NzU8++eTVD3jgwIFDhw5du3ZNNrpNeXl5fHy8lpaWRCKZP3/+tm3bNmzYIJFIcnNzTUxMsrKyCCEbNmxgsVgvXrxgMBghISGyQ7r6+vqpU6du27bN39+/rq7Ox8fn7NmzU6ZMiYiIcHZ2XrVqVZP3vXLlyuHDh1NSUgwNDX/66Sf5gaZAIJAdDX/77bebN28OCgoihMiOsw8cOBASEiL74MqHQvhWGOx/3HinokJmBC84efLks2fPwk9dmR3g2LsvXdEAuqT8/Py8vLzXfVVxcbHCQtjQ0JCWltb4bF4ThYWFjQshISQ4OFh2ZObj47Nt2zZCyMWLFydMmGBra0sICQgI0NPTi42NHTt2LCHko48+kk8QYW5u7u/vTwgZNmwYj8dbtGgRIcTCwsLMzCwzM9PAwKBXr15hYWHJycmVlZUMBuPRo0eNC2ETQqFw+vTp06dPX7JkCSHk5MmTX3755bVr1wghBgYGsuOqCxcuzJs3T3b779KlSzdv3tzyXtLT04uOjo6JiSkqKsrLy5OfwmUwGMuWLZM9vnz58oYNG5hMJiFk8eLFssPfuLg42ayKV65cEYvF1tbW0dHRkyZNau6NLly4MGPGDNm/yEcffSQ7gmzMwsJix44dFEWNGTNGX5/+CXxQCNsZk8mcPn368ePHnz2+nZqR39ud7kAAXcrgwYP79u3b5LJWy9hstsIqSAjhcDjLli1r7hIai8UyMjJq0igfVkI+fHZlZaWurq58BT09PVkfSEJI41/i8hey2ewmT2XbmTZtmqam5vTp0zU1NZ89e9bChT2pVPr+++/r6+t/+eWXhBCJRFJeXp6QkPDixQvZCrJhS+UDmRJCGidszq5du3777bdVq1YNGzbsypUrCQkJsnZ1dXX5GAg1NTXyuYXlH6G0tFQqlcrPCVtYWLRQwmXB+vb96yCAw+G8OsDCjh079uzZs2fPnrlz565cufKbb75pNXyHQiFsf0wmc8aMGc/tdMy0m17WBoCWURQlm8+ovWhpab3lIIg2NjZ79+6VPa6rq0tNTX2DIUDFYvGFCxdKSkp0dHQkEon8CEyhL774IjU19dq1a7JhjBgMhr29/TvvvNNkAG4bG5vHjx/LHicmJr66HYqiWCyWvG9ORETE//3f/4WEhBBCLly4oPCtnZ2dY2NjR48eTQi5e/eurNHR0ZHH461du1Y+kpwMm81W+CeLjY2N/DaJtLS0V6/5cTic5cuXL1++PDMz09raes2aNSoqKgr7ECkHCmGHYDKZtjZ9hBW4nRagy5s5c+ZXX331ySefjBkzZufOnYMGDRowYMDrboTJZNrb23/33Xd+fn5HjhyRXzZ71cWLF7///vuff/45MjKSEKKrq+vr6/vtt99+9NFH5eXldnZ2mZmZL1++3Lhx40cffeTl5eXg4GBpaSmbTLgJiqLc3Ny++uorZ2fngIAAFxeX/fv39+7dOzEx8dSpU40nR5Rbs2bN+PHj+Xw+m80+d+6crBLb29tPnz49MDBw1apVHA7n4cOHrq6usj41+/fv19HRsbGxkdVOmffee8/FxeXHH3/s16/fTz/99OoR4apVqzw9Pc3NzW/cuGFjY6OhoWFlZVVeXv7jjz9qa2svWLDg8OHDu3fvvnfv3uvu5zfD3LRpk3LeqX1VV1fv379/xYoVbVy/oaFBNlTuq8S1WRRFMTX6tFs4QgghkvoSSX0xW9e5fTfbcVrYRSAjFArZbDbGGm2BRCKRSqWvNa/6li1bFi1a1ElGrhcKhU5OTubm5oQQiUSipaXl6enJZrNnz5798OHDmJiYQYMG/fzzz7IPWF9fP2zYMNlQ2lKplMPheHl5ybbT0NDg6+sru9ImEAgGDx6so6Mzfvz4mJiYu3fv+vv7+/n5OTg49O7dWyqV9urVy9XVlaIoJpM5atSosrIyHR2dqqqqgoKCgoIC2cvt7e19fHxu3Lhx7dq1hoaGgIAAc3NzfX39MWPGnD9/PjU1ddOmTVpaWvI3lfP39y8sLMzNzXVzcwsICCguLo6KitLT01u9erWRkdGAAQOkUilFUbK+qYQQU1PTiRMnZmdn6+vrT5kyJTY2VnaR8t133+VwOJcvX3706JGRkZGfn5+mpubAgQN1dXUzMzO5XK6Li4v8TdXV1d99990rV648evRo9erVxsbGI0aM4HK5EonE2dnZ1NS0pqbm+vXrMTExmpqaO3bs0NLSUlNT8/Pzy8jIyM/P9/HxkUql+vr6gwYNavxZfvvtt6CgoOZOg78NDLr9hoNut0pUmSKsSOL0DWrfzXYcDLrdKgy63SoMut3V1dXV5ebm2tvb19TUBAcHOzs7v+7kUx0Hg24DAECHEwgEc+bMMTMzc3NzMzMz27BhA92JlAHXCDuSuF5cl9+4oaq6hsvhsNksQghTzYgwsP8BoBPR1dV98OAB3SmUDb+IOwqloiMV1wtyI+QttXWCHQev6mhxZ0/24rL4qqZj2PqvfckdAADaFwphR2FyTbn2ixu3qIhEJg7swsLC0/dYM0ZbqZLXGDsDAAA6SEddIzx37tzatWuDgoLu37+vcAWpVLplyxZvb+9x48bJ79MkhGRmZs6aNWvQoEHLly+Xj8faPbBYrJCQEENDw6KioqN/XBd0yglTAAB6mo4qhIcPHxaLxbdu3WputKTdu3fv3Lnzp59+WrBgwdSpU9PT0wkhYrHYz8+vb9++e/bsefny5eLFixW+tutSV1efP3++gYFBQVH5kRPnO+fkYQAAPUpHnRqVTTUpuyFUoe3bt3/55ZdDhgwZMmTIhQsXdu/e/cMPP1y8eFEgEHzzzTcURW3fvr1Pnz4FBQVNRgLs6tTV1UNCQvZsSczLLzp27FhwcPBr3XQF0P1ERETo6enRnQI6u447R0jPr2ChUJiSkjJ48GDZ08GDB589e5YQkpiY6OnpKbtPy8TExNTUNDk5uZsVQkKItrb27Kk+h//7IDMz8+TJk0FBQU3ugQXoOUJCQq5e7VkTlolEIvz52zKFu+idd94xNTXtiLej5x+jpKREIpE0Hi5WNrtjcXFx46Fj5e2v4vP5hYWF8nFdCSFLly5t4VRqC+PbiurqKIoh6LBJnBVSU2VNCfQ5+mfikydPxGLxu+++S/tt2i3sIpCpq6sTi8W0/0t1ZrIb6luY7eFV69at67g8nVNtba1sPBpoTgu7qOY1f1erqanJhkFvAT2FUDacUl1dnazs1dbWyoqilpZW44mheTxekzFe5TgcjqGhYeM/JI2MjFr+2Wp2ZJlaLkUxVJQ7qEo9h2NpYPr++x4HDx5MT0+PiYnx9/en/TcsRpZpGYPBwMgyLXuDkWV6JnzXWqXMXURPIdTQ0DAwMEhPT5eNq5Senm5paUkI6dOnj3xM9Pr6+ry8PFm7Qkwm08rKSjmB2x+DVZ93QYvBmuxRcyzi4Z2o28K8P8Z4/29Ie4ql7vAhxeLSGhEAoEdQ6hBr8fHxe/bskT0ODg7esWOHVCotLy8/evRocHAwIWTSpEnJycmxsbGEkAMHDtjY2Li6uiozodKomo3VcF6p7rjc3mfjnH/9xDEemphv8qjcQ91xubrjckIxpBJ0KAUAUIaOKoSjR4+mKColJWXq1KkURd2+fZsQEh8fv3v3btkKGzdulB3w2draTpo0STbds66u7s6dOydMmCCbr0ReNbsfimJSLI7sPxs752nTZzJYKlev3YyNS6RYHIrg5BsAgJLQPPtEfn4+l8uV95qREQgEhYWF5ubmLfSl7PyzT7yu+Pj4iIgIQsjUqVP7Ulc4tvMZKjqtvqodYfaJVmH2iVbhGmFb4LvWKiXvIppnnzA1NW1SBQkhqqqqlpaWPe2OAnd3d19fX6lU2gNHvAUAoBHuZelEvL29jYyMdHV1SVEY3VkAAHoKFMLOxdbWlhDCU3zzJAAAtD9MzAsAAD0aCiEAAPRoKIQAANCjoRACAECPhkLY2YlEot27d589e1YiwYz2AADtD71GOyOKwa57tpNQFCGkQSgufnYrK0HYkPvH+JGOFEUxOKZcmxC6MwIAdBMohJ0R136xVCKSPVYnZP7HMw4eDkspEeoVuvkM61+fdZreeAAA3QlOjXZKDLZ8JFKKxTG3tA6aEcxkq966cz82/gnd4QAAuhUUwq7B1tZWNnnvpUvRKWkFdMcBAOg+UAi7jH79+vn6+koJOXc1MT09ne44AADdBAphV+Lt7e012EMilp48ebKgAMeFAADtAIWwixnj6+PqYCYQCMLCwsrKyuiOAwDQ5aEQdjEURU0Y5WJtbc3j8cLDw3k8Ht2JAAC6NhTCrofJZAQFBZmampaXl4eHhzc0NNCdCACgC0Mh7JJUVVXnzJmjr6+fn59/9OhRsVhMdyIAgK4KhbALkjQIK5LZghczAjw4TF7G0/tnwrc3lCcJK5KFFckSQTnd+QAAuhLFhdDe3v7+/ftNGh88eKCnp9fxkaAlFEuDodFHXJUirkrRovKmj7FhiSsSH96Oux0hrkoRFt8SlsTSnREAoCtRPMRaVVWVSCRq0tjQ0FBVVdXxkaAlFFOV02eq/KlVHxJs8M6ff/6pYz1OrY+TsOS+WFBKYzwAgC7nNcYaffDggbGxccdFgTdjY2OzYsUKulMAAHRV/yiE+/fv/+abbwghpaWl06dPV1NTky8qLy+vrKz88MMPlR0QAACgI/2jEFpYWPj6+hJCDh8+7Onp2atXL/kiAwMDFxeX6dOnKzsgAABAR/pHIfT19ZUVQolE8n//93/29vY0pQIAAFASxdcI9+zZo+QcAAAAtGi2s0xeXt758+fz8vIEAkHj9u+//77jU8Gbk9aXCMviJRJJQWGpqYkhRVGNl7J0nCmmWnOvBQDogRQXwv3793/00UcCgUBdXV1FRaXxIhTCzoypYSHmF4l5+XEJzy9ci+/vYh0wZpB8qbgmg2JpsLRxxhsA4G+KC+Hq1auHDx++b98+CwsLJQeCt8HgGKtZBBJCbFQKuUmCpNwGg0zNUaNGyZbyM4/SGQ4AoFNSMLJMcXFxWVnZV199hSrYdRkbG0+fPp3BYNy4cePRo0d0xwEA6LwUFEIdHR11dXU+n6/8NNCObG1t/f39pVJpREQEZrQHAGiOgkKooqKycuXKzZs3Y36frs7d3X348OFisRgz2gMANEfxNUI2m/3w4UMnJycfHx99ff3Gi9BZpmsZPXp0VVXV48ePjxw5MsdX55//mAAA0EwhDAsLEwgEAoHg9OnTTRahEHYtFEVNnDixuro6Kyvr2H+T31/iqaFNdyYAgM5E8TRMqamp5c1Qcj54e0wmc9asWb169Sopqz5xJurVeUUAAHoyTMzbI6iqqs6ePVtTg5OV8/Ls2bNSqZTuRAAAnYXiU6PZ2dlisVjhIisrq47MAx1FS0tr5kTv8Ki0pKQkPT290aNH050IAKBTUFwIBw8eXFRUpHARDia6LiMD7akT/Y7/eTcmJkZbW9vd3Z3uRAAA9Gt20O36+nr509LS0uvXr1+5cuXbb79VVjDoENZWFgEB+ufOnYuMjNTW1raxsaE7EQAAzRQXwsDAwCYtH3744VdffXXy5MklS5Z0fCroKFJhbX8Xq4rSATE37548fnj1/33896jconpCNGlNBwBAg2Znn3jV3LlzP/vssxcvXvTt27fjAkHHYagaNBTdIkW3BptKG6xqhaKq+vTD8qWC2mLNQV8SikljQgAA5XuNQsjj8eT/h65I1WyMqtkY2ePxLk2X8u9vJFIpoZq2AwB0b23tNZqZmblx40YtLS07OzulBAMAAFCG1+g1amRkFBoa2mR6QgAAgC6tTb1GKYoyNjYeOHAgl8tVVjAAAABlaGuvUQAAgG6ppc4ylZWVz549y8vL69Wrl6Ojo4GBgdJiAS0yMjK1dHQNDQ3pDgIAoDyKC6FYLF67du22bdsEAsFf67FY8+fP//XXX9XU1JQYD5SnXiAKCw9XVeMsWLCgV69edMcBAFASxYNub9q0afPmzZMnT/7jjz/u3bt37ty5BQsWhIaGLl26VMn5QGnUVFnOzk719fVHjhypqamhOw4AgJIoOCIUiUTbt29fvXp146kHAwMDXV1dV6xY8eOPP+rq6ioxISjP5EmTeHX1L168OHz48HvvvYejfwDoCRQcERYXF1dWVgYHBzdpDw4OFovFGRkZSgkGNGAQ4bTJAfp62sWF+adOHBU38KQi/l//STCLIQB0TwqOCDU1NSmKys7OdnNza9yenZ1NCNHWxgTn3RNDRYf3bAchZLI7/+CZRyn3Yk9XXZ8wyokQQqQSpnpvjvUcmiMCAHQAxYVw+PDhy5cvNzU1lc/U8/z584ULF9rZ2WG+gu6Kbb1YQ1OTEKJByIK++aGhoSllDaZV3sOGDRPXZgsKoukOCADQIRR3ltm5cyefz/fw8LCzsxsxYoSjo6Ojo2NGRsaBAwf+nqwAui9TU9OpU6cyGIyrV68+fvyY7jgAAB1IcSF0dHRMSkr6+uuvra2t6+rqzM3N161bl5ycPHToUCXnA7rY29uPGTNGKpWeO3cuJzef7jgAAB2l2RvqDQwM1q9fr8wo0Nl4eXlVVlbGxsYeP31u9lhTDK8HAN2S4iPClJSUhw8fNmlMTk5OSEjo+EjQifj5+Tk6OvLr6o+du1tbW0t3HACA9qe4EAYFBZ0/f75JY1xc3NixY0UidKPvQSiKmjJlipmZcWV13ZEjRxoaGuhOBADQzhQUQh6Pl5yc7Ovr26Td19e3pKQkKytLGbmg02Cz2bOmT9TVVs/Pzz916pREIqE7EQBAe1JQCKuqqggh6urqTdplczBVVFQoIRZ0KlwuZ0bgEA6H8/z587S0NHbsveMAACAASURBVLrjAAC0JwWdZQwMDLhc7s2bN11cXBq337x5kxDSu3fvNm5aJBI9ffpUV1fX3Nz81aW1tbVCoVD+lMFgyG7Vr66uFovFf4VjsTQ1Ndv4dtCh9HU15syZ8vDhw7b/AAAAdAkKCqGKisrUqVPXrl3bu3fvgIAAWePNmzf/9a9/+fj4GBsbt2W7mZmZY8eO1dLSys/Pnzx58u+//95khQ8++CAqKkr2mM/nOzo6yrrnDBkyJC8vj8ViEUK8vb3//PPPN/5s0L7MzMzMzMzoTgEA0M4U3z6xZcuWhISEwMBAIyMjMzOzwsLCgoKCvn377t+/v43b3bBhg7+//y+//FJeXu7q6jpz5syRI0c2XiEsLEz+eMSIEe+++6786fnz5729vV//s0CHYbDFvNzaJ/9pbrmq6Ri2/kBlJgIAaC+KC6GBgUFsbOyhQ4cuXbpUUVFhYWExatSo999/v40nKoVC4enTp+Pj4wkhenp6U6dOPX78eJNCKPf8+fN79+6dOHFC3lJXV1dUVIQp8ToPJtdUw3U1kUoVLm0oipEIcWcFAHRVzd5Qz+FwlixZsmTJkjfYaGFhYUNDg5WVleyplZXVlStXmlt5//79/v7+jc+4hoSESCQSFRWVHTt2TJw4sbkXCgSCxpt1dXVF7ew4FLPZKZkoBltxhQQA6AqaLYRvg8fjEUJUVVVlTzkcTnP3YotEosOHD+/cuVPecvHiRVl3jLCwsNmzZ6empiq8LsXn86urq7/99lt5y7Rp00JCQpqL1MLN4KK6OopiCHr8VLRvfL+8uK5OSrEaesAOrKurE4vFGG63BUKhUCwWN+4HB6/C2BStasddpKamxmazW16nQwqh7MisoqLCwMCAEFJWVtZcF5vz589LJBI/Pz95i7xT4pw5c77++ut79+5NnTr11RdyOBxDQ8Po6NeYEqG587qCWi5FMVTQPbX5XdSCgoKCvBf5rs72aj1gBzIYDC6Xi0LYAlkhxJTOrUJ/+FYpcxcpHlnmLenq6lpbW9+5c0f29M6dO/LpnJrYv3///PnzFZZrgUBQVlaG6Q87ucuXL5+7cDs6Jp7uIAAAb6hDjggJIR9//PEnn3zC4XCePHly9+7dAwcOEELS09NHjRqVkpKipaVFCCkqKoqKivruu+/kr0pNTT148KCXl5dUKt25c6eJicnw4cM7KCG0i2HDhmUmXbv74Imh1QNPT0+64wAAvLaOKoRLly5ls9lbt27V09OLjo42NDQkhGhqak6cOFF+/JeRkbF27VoHBwf5q/T09BoaGnbv3s1kMr28vJYtWya/0Aidk5WV1bt+w/44f/v8+fOampqN/zUBALoESqqoT3xQUNCmTZucnJwaN6akpHz88ceXL19WVraW5Obment75+TktHH9mpqaZq8RFlyjKIaKseK7O3qOFnZRyxoKom/cTbr9uJzNZs+bN0/hQELdA4/HwzXCluEaYVu88Xet51DyLlJ8jTAmJqaysrJJY2Vl5dWrVzs+EnQ9w7369e/fXygUHjt27NWfHACAzuw1OstkZWXJznACNEFRVGBgoLW1dW1tbXh4OJ/PpzsRAEBb/eMa4enTp2W39FVUVPz73/9u3GOTx+M9evQoMDBQ2QGhi2AymUFBQQcOHCgsLDx27FhISIhswFgAgE6urb+qjI2NV69e/fHHH3doGuiixDUZAkkDIWTaaNMDRxMzkzJP7H8xxX84RVGEUCq9vFsYmAYAgF7/KIRTp06V3b0+duzYzZs3u7m50ZQKuhKmjhNhqMgea2mrBU+bEHo08mla3rXbSe+M9Gwouc/SdmCqY9oKAOikFB8RXrp0Sck5oOticoyZnL9HDjLrRWa95xweHh6bUmXQV9VZG73jAKBTa/bUaH19fVxcXFZWVl1dXeP2xYsXd3wq6NqsrKwmTZp05syZqKgo1iDKDUeDANCJKS6Ely9fXrBgwcuXL19dhEIIbeHq6lpRUREdHf3fSwl61gWWdiiGANBJKb59YuHChTo6OhcuXMjLyyv/JyXng65rxIgR7u7uQpH41p1YurMAADRLwRFhWVlZTk7O1atXR48erfxA0J34+/tr1d+1HjiU7iAAAM1SUAg5HA6bzcYgSfD2GAyGZ7++qsZGdAcBAGiWglOjXC537ty5svkiAAAAujfFnWUGDx68fv367OzscePGNRn5FJ1lAACgO1FcCDdu3FhSUnL58uVX55pAIQQAgO5EcSFMTU2VSCRKjgLdlVRUI2moULiIYqhSLK6S8wAANKa4EDYebhvgbTBUDQR5Fxu38AXC6Dupbg5mvU10iFSs7rKKrmwAAKSFkWWEQuGlS5eSk5NFItG6desIIc+ePeNyuRYWFkqMB12emuXkJi0F6ekppWEZD1TnzZmoVXmCllQAAHKKb6jPz893d3cPCAj47LPPfvvtN1njnj17Zs6cqcRs0D1ZW1u7ubkJBILwo8cqqzFzIQDQTHEhXLRoEZ/Pj4uLu3DhgrxxxowZsbGxVVVVysoG3RNFURMnTrSysuLV8o5FPGoymC0AgJIpKIS1tbUXL17cunWru7s7RVHydhsbG4lEkpubq8R40D0xmcwZM2YYmxiXV9WFh4c3NDTQnQgAei4FhbC6ulosFltZWTVpF4lEhBD8zoJ2oaqqOid4pq4W5+XLl6dOnUIvZQCgi4LOMoaGhurq6rGxsY6Ojo3bL1++zGKxbGxslJUNujl1dfUZE9yOXKt4lnjjtPDlRH+fxmcgCCEMjgmTa0pXPADoIRQUQjabPXv27DVr1piamqqoqBBCJBJJVFTUypUrp02bpqWlpfSQ0D1RTNVe1sOmqxaEnYpOTHikpSYc6eUqXyptqKBqXjD7TKMxIQD0BIpvn/jpp5+eP38+btw4NTU1kUiko6NTU1PTv3//X3/9Vcn5oDujGKpmftZmZKbO0OPHj9992qBvbeLh4SFbKKxIElc9ozcgAPQEiguhhobG1atX//zzz0uXLpWUlGhoaIwaNWrmzJmyA0SA9mVvbz9p0qQ//vgjMjKSy+U6OTnRnQgAepBmb6hnMBgTJ06cOHGiMtNAj+Xm5lZZWRkdHX3mzBl1dXVLS0u6EwFAT6H4PsIpU6bs27evurpayWmgJxsxYsSQIUNEItHRo0cLCwvpjgMAPYXiQlhbW7tw4UITE5OQkJCrV6+iazsox7hx45ycnOrr648cOVLLw432AKAMigvhpUuXnj179sknn9y6dcvX19fCwmLNmjVpaWlKDgc9DUVRU6ZM6du3b3V1dWUlTkgAgDIoLoSEEHt7+02bNqWnp1+8eHHkyJHbtm2zt7cfMWKEMsNBD8RisebMmbN06VJzM2O6swBAj9BsIZRhMpljx44NDw+/ePGihYXFzZs3lRMLejImk2lgYEB3CgDoKVophMXFxVu3bh0wYMCIESMqKysxPT0AAHQzim+faGhouHjx4uHDh8+ePSsWi728vHbt2jV79mx1dXUl5wMAAOhQiguhi4tLWlqavb39F198ERISYm5uruRYAIQQqVQsFTUzYSFFUUw15cYBgO5JcSGcOnVqYGDg0KFDlZwGQI7B4oprs3hPt8lbhCIxm8X864lEqGY5maXjTE84AOhGFBfC7777Tsk5AJpgalppuH4qf3rjxo1bt25NmfKubFKU+pxzUrGAvnQA0H0021nm0aNH8+fPHzRokPy4cP/+/cePH1dWMIB/4HK5QqHw9OnTL168oDsLAHQrigvhxYsXhwwZcufOHV1d3ZycHFljXV3dp59+KpVKlRgP4C+enp5eXl4ikejYsWP5+fl0xwGA7kNxIVy+fPnEiROTk5PXrl0rbxwzZkx2dnZBQYGysgH8w9ixY93c3AQCQXh4eHlFDd1xAKCbUFAIS0pKnj9/vnr1ajab3XjGcFnfURRCoAtFURMnTrS1teXxeOGnr1VV19KdCAC6AwWFUHbyk8Foukg2IQCHw1FCLACFmExmUFCQpaVlVQ0v7HhEbS1qIQC8LQWF0NDQ0MLCQtYvpvER4e7du/X09Ozs7JSXDuAVbDY7ODjY2EivvKIqLCysvr6e7kQA0LUpuH2CoqiNGzcuXry4oqLC2tpaKBSeP3/+xIkTBw8e/OGHH1isZufyBVAOVVXVWZNHhv03vrCwMDw8PCQkREVFhe5QANBVKa5qCxcu5PP5GzdurKqqIoRMmDBBTU3ts88+W7VqlXLjASimzlULnj7h8Nn43Nzc06dPBwUFMZnM1l8GAPCKZu8jXLZs2cuXL6Ojo0+cOBEVFZWfn79p06bGZ0oB6KWjrTl37lx1dfXU1NSzZ8/ixh4AeDMtnedUV1f38fFRWhSA1yKuy9fWUJ0R4HHo6H8TH0QzGwonjBsp+1uNYrCYWnb4uw0A2qKVaZgAOiemZl8i5ourUow4pUHj7Jii8rjYG5cij4mrUsRVKfysU0SEGw0BoE3Q8wW6JLauK1vXVfbYrg8JNhpz7NixuAzxYL9RhoaG4uSfpVIpjgcBoC1QCKE7sLW1DQ4OzsrK0tPTozsLAHQxKITQTVhbW1tbW9OdAgC6HlwjBACAHu3vI8KGhgYej9fqC3R1dTsyDwAAgFL9XQjDw8Pfe++9Vl+Au7UAAKA7+bsQDh06dNeuXbLHtbW133zzjbGx8ZQpU0xMTEpKSqKiopKSktatW0dTTgAAgA7xdyG0t7e3t7eXPZ41a5a/v//BgwfltyR//vnna9asuXDhwvr162mICfBG0tPTb968OXbsWDMzM7qzAEAnpaCzDI/HO3ny5KpVq5oMzLFq1apbt25lZGQoKxvA2yosLMzOzj58+DDm0QSA5igohDU1NWKx+NWZ3mQtsmG4ATo5MS9HXJM52NXYydqQX11wYPfPeWn3xDWZ4ppMSX0J3ekAoBNRUAh79eplbW29cuVK2Uy8MpWVlUuXLtXV1XV0dFRiPIA3wdSyFZY9EhTdEpbcGT9Iva+RkFeWdmjf1pdPzwsKrtVnnaQ7IAB0IornI9y7d29AQECfPn2GDh1qbGxcWlp69+5dgUBw9OhRzFAPnZ9a74DGT0NsQo4fP/78+fMTMfyQWb6a1VfoCgYAnZDiG+pHjRr15MmTDz/8UCqVPnjwgM/nh4SEPHz4cOrUqUrOB/D2mEzm9OnT+/btW1tbG3bkRGV1Hd2JAKATaXaItb59+/78889vvF2RSLR///6HDx86OjouWbJETU2tyQqRkZFJSUmyx2w2e+XKlbLHPB7v999/z8jIGDx48Ny5cxkMjH0D7YDNZs+aNSssLCw781nY2cTF9pU6Ojp0hwKATqGlMlNSUhITExMREfEG2126dOn+/fsHDx4cGRk5c+bMV1c4efLkpUuXKv5H3h4YGHjjxg1PT89ffvll7dq1b/DWAAqpqKjMmTPHwsK8qoZ/4MCByspKuhMBQKeg+IhQIBB89NFHoaGhEonEzMwsLy+PEBISElJRUdGWulhUVBQaGpqWlta7d+/p06cbGxunpKQ4OTk1WS0gIGDFihWNW+7du5eYmJifn6+qqurl5TVo0KD169draWm96acD+AcVFZVZM6aGFsUWV1UdOnRo/vz5+OkCAMVHhCtWrDh16tS2bdv27dsnb5w7d+6VK1f4fH6rG71//36fPn169+5NCNHQ0Bg0aNCtW7deXe3WrVubNm0KCwurr6+Xt3h7e6uqqhJCHB0dtbS0Hj169AafCqA5qqoqMwM9zczMysvLDx48WFOD+XsBejoFR4QCgSA0NHT79u3vvffejRs35O0uLi4CgSA3N9fOzq7ljRYWFhoaGsqfGhoavno7s42NTVVVFYPB2L59+/fff3/v3j0NDY1XX5ifn6/wLQQCQVlZWePOO9OnT3/33Xebi8Tn85lMpsJFIoGAEEpU19M7ULSwi7oTiYDPoKRTpkw5cuRIYWHhnj175syZo6Gh0ZbX1tXVEUKaDDQBjQmFQrFYLJFI6A7SqfWQ79rbaMddpKKiwmK1MuGggsWlpaV8Pt/Ly6tJO5fLJW27oV5FRUUkEsmfNjQ0qKioNFlnw4YNsgfr1q1zd3c/cODAsmXLXn2h7OjwVWw2m8vlBgUFyVtcXV2bW7nlTVFsNkUx2M2/todoYRd1JxLCFYur2Hm7gwYLj55LLkqrPfTz9eDAgRrqf312lV7D2YZNf/hlRCKRqqoqCmELGAyGWCzuCT9Ib6OHfNfeRjvuorb0uFRQCHV1dVksVmZmZpN75+Pi4gghshOeLZNfVpTJy8ubMmVKswlYrIEDB2ZlZcleeO/ePVm7WCwuKChobohIBoPB4XBmzJjRahgZJpPZ3N8XDAaDohj4A62FXdSdMLkGWv3WEKlUk5BFDnWHwo4UFRY9qXIfO8iXECIsi5eKapvbD7JdhELYAtmxYE/4QXobPeS79jaUvIsUlEoulztu3Lj169cXFBTIv/O5ubkrV6709vY2NjZudaPDhg3j8/k3b94khDx//jw5OdnPz48QkpOTI79YKB/CraKiIjo62s3NjRASEBBw79693NxcQsjFixc1NDTc3d3b4VMCNEIx1SgWh2Jx1LX05y9YOGq0r7unl6yFUK2cQgGA7kfx137btm3Dhw+3tbW1s7OrrKz09fW9e/eumppa40uGLVBTU/v++++nTZv2zjvv3LhxY8OGDbIrf5GRkbt27UpISCCEmJmZDRkyhMPh3LlzZ+TIkXPmzCGEWFpaLl++fNiwYd7e3pcvX/71119bPbcL8DY4HM6oUaPoTgEAdKKam2i3pKRky5YtFy5cKC4u1tHRGTly5Keffmppadn2TWdkZCQmJjo6OspPsZaXl5eVldna2hJCCgsLExMT6+vr7e3tHRwcGr/w8ePH6enp7u7uLbxdbm6ut7d3Tk5OG8PU1NRoamoqXCQouEZRDBXjkW3cVHfVwi7qORqK70mFlapmfgqX8ng8LpeLU6MtkHWWeXUADWgM37VWKXkXNVsIOzkUwnaHLydBIXxrKIRtge9aq5S8ixR3p9m+ffvJk01H6C8tLf3Pf/7T8ZEAAACUR3Eh/Prrr4OCgj7++GOxWCxvLCoqWrNmjbKCAdBJLBY3/uEHgG6s2Rss3nvvvZ07d06dOpXH4ykzEEBnsGvXrm3btpWVldEdBAA6XLOFcPbs2dHR0bdv3/b29m58UyBAT6CtrV1VVRUaGlpcXEx3FgDoWC3dcu/t7X337l0+nz9kyBCM+Qk9hLg2R5B/ZeJwfQs9fmV+4r5fN+QknBDkXxHkXxEVXZMKMTYpQHfTytgzNjY2N2/etLCwGDVq1MWLF5WTCYAuLB17lo4TxVRTUdWYNc3f3s66rl50+FR0fnENxVST1DwX1zUdNRcAurrWB2EzMjK6du1aQEDAqlWrlBAIgEYMFV2VXsNk/3FMR8xauMHV00/IMj4a9fRlnRmlokt3QABof4rHbdmwYYO1tbX8qaqqalhYmJOTU2JiorKCAdCPyWROmzbt7Nmzjx8/Pnr06AQPlrMJ3ZkAoL0pLoRLly5t0kJR1Pr16zs+D0DnwmAwJk+ezGKxHj58eCoyjqXr6uJpT3coAGhPfxdCPp9fUVGhpaWloaFRVFTU3E1UpqamysoG0ClQFBUYGMhms29EJpz847xE1Uw2RjwAdA9/XyM8duyYmZnZr7/+Sghxc3MzawZ9UQFoQ1HU+PHjR3k5SSXSP/74IzY2lu5EANBu/j4iHDFixKFDhwYOHEgI2bFjB5/Ppy8VQGc0ZKCtZi/1K7efXrhwQSQSeXt7050IANrB34XQ2tpa3kFm2rRpNOUB6NQGefRT1bGKjIy8fPlyfX39O++8Q3ciAHhbmO0PoO0Y9bnnHDkqYre6iOgn0WdimcWRnm7/myyMpaZuv5gQzE0B0MX8XQhjY2PDwsJafYHsIiJAD8QyC+CoUhRFedoSbbuM8xcuGTiN5tj+1YmUl/KrVCqhKCa9IQHgdf1dCLOzs//4449WX4BCCD0WxVRjqPw1H6GDs4eDs8c/luJYEKBr+rsQBgUFBQUF0RgFAABA+VofYg0AAKAba6WzDJ/Pr6+vb9yiq4vhFgEAoPtQXAjr6+vXrVt3/Pjx/Pz8JoukUmnHpwIAAFASxYVw5cqV+/fv//jjj2/evKmrqzty5MioqKgHDx5s2LBByfkAupbLly/HxcVNmTLF3h5DkgJ0DYqvER45cuSHH374z3/+Y29v379//9WrV1+/fn3RokVXr15Vcj6ArkRcz2ZKBXzesaNh8Q/uSkV8+X9EIqI7HAAopuCIsKSkpKqqSjZkBovF4vF4svbVq1ebmZnl5ub27t1bqRkBugJKVZf3dLuHIamzLL4V/+L0gXtFj/sO97SiKIpIJUx1C471bLozAoACCgohl8slhAiFQkKIiYnJkydPGreXlJSgEAK8St3xr8nLJrgSk4GPIiIi7ueKG/QdAgMDpXU5wqIYeuMBQHMUnBpVV1e3tLRMSkoihAwdOvTixYsXLlwoLi7esGEDm822srJSekiALmbAgAHBwcGqqqoJCQnh4eECQQPdiQCgWYqvEb733nsPHz4khIwbN87b23v8+PG9evX67bffPvvsMx0dHeUmBOiSrK2t586dq66unpmZefjI6dq6+tZfAwB0oFq9HUIkEsXExGRnZ7u5ubm7uysnVqtyc3O9vb1zcnLauH5NTY2mpqbCRYKCaxTFUDEe2X7puqQWdhHI8Hg8LvevIdbaqLy8PCwsrLQwU5NRMm/5z0ZGRh0XrzMQCoVisVhNTY3uIJ0avmutUvIuan1kGRaLNXr06AULFnSeKgjQVejp6S1atMiyt1lVdd2+ffvS0tLoTgQATbU0skx5eXleXp6s14wcyiHAa+FwOLNnTv7jSE5qieDo0aMTJkzw8PBo/WUAoCyKC2FqauqiRYtu3rz56iKMLAPwulgs5sSx7neyzW/duhUREVFdXe3j4/Nap1gBoOMoLoTTpk2rqKjYunWrra2tioqKkjMBdD8URfn6+urq6kZGRsbExJSXl0+ePJnJxOSFAPRTUAjLysqSkpIiIyMnTJig/EAA3RBFifnF9VknnfWJ2kjdP6LuJdxKcTUutjQ3lC1n6buzNHFjEgA9FBRCdXV1NpuNWSYA2guTa67WO0AqlRBC7NycFpoPzMp+2cfRkcFgEEJEFU8kvJcEhRCAJgoKoZqa2vz58/fu3evl5aX8QADdD8VgsXQc5U+NdYmx9d9LJfWFNGQCgP9RfI1w27Zt8+fP9/b2HjNmDIfDabzo008/VUowAAAAZVBcCM+fPx8REcHj8e7cudNkEQohAAB0JwoKoUgkWrx48YABA7Zs2WJtbY1O3gAA0I0pKISlpaWlpaU//PCDp6en8gMB9FgJCQkPHjzw9/c3NTWlOwtAD6KgEOrr62tra9fV1Sk/DUDPJK4vFlYk56XH5mY83vtbkv+4kW4uf01wT1EMprYdReGOQ4COomCsUTabvWnTps8//7yiokL5gQB6Gia3N0Uk4qoUnwE6/e00G2oL/jh97Py5I8KKZHFVSn3uOQm/mO6MAN2Z4s4yycnJaWlp1tbWHh4eTeZdOnHihFKCAfQULG07lrad7PEU6xmW8fFRUVEPX4jLJLzp06czBEcIwbiGAB1IcSHMzs6WTUNfXl5eXl6u3EgAPZq7u7uRkdGJEyeys7P37t37rofQsjfdmQC6NcWF8NKlS0rOAQByvXv3Xrx48fHjx/Py8g6dTniX1X/gEHSfAegoCq4RFhQUUBQVGRmp/DQAIKOpqTl//nx3d3ehSPzfPy9ERESIRCK6QwF0TwoKoaamJoPBwATKAPRisViBgYHv+g5gsdlxcXH79u2rrKykOxRAN6SgEGpoaIwfP/706dPKTwMATbg6mM8PmaGjo1NQULBnz57MzEy6EwF0N4qvEc6fP//DDz8sLi4ODAw0MjJqvMjX11cpwQDgL6YmvRYvXnzmzJn09PSwsLCRI0eOGDECQz4BtBdK4YzzxsbGRUVFCl/QSWaoz83N9fb2zsnJaeP6NTU1zZ3sFRRcoyiGivHI9kvXJbWwi0CGx+NxuVwlV6C61N0sLRtKRVsqlcbceRhzO54Q0t/VPnD8SEIIxeQ2nteCdkKhUCwWq6mp0R2kU8N3rVVK3kXN9hoVCoVKCwEAzWEbuIt5+UTII4QM629iqutx7uLd+ppiMS+fSEWi6jSNzlQIAboixYXQzc1NyTk6lkQkqX4mFHEVL6wvYXJ6KTkRQBux9d3Z+u7yp04WxGn4ItljqahOVJ1GUy6A7kNxIZQpKipKSkrKy8szMTFxcnIyNzdXWqz2Ja7LFxVcEBs4NLOcYmhYKjUQAAB0GooLoUgk+ve//71r1y75rUsMBmPmzJm7d+9WV1dXYrx2QhFKRV+tz3S6cwAAQKejuBCuX7/+999/nz9//syZM42NjUtLS8+dO/f7778zGIzDhw8rOSIAAEDHUVAIhULhzp07N27cuGnTJnmjj4+Po6PjRx99tHXrVn19feUFBIC2SUxMTE9P9/X11dbWpjsLQFei4Ib6kpKS6urqqVOnNmmfNm2aWCx+8eKFUoIBwOtJSUl58uTJzp07k5OT6c4C0JUoOCLU0tJiMBgZGRmurq6N29PT0wkhTWZlAgDaMFhSSUPNo02yZ+/0EdZnp6a9KA3/9bybvckYbxt1kyFqvf1pjQjQBSgohBoaGj4+PkuXLtXX1x8+fLisMSEh4b333nNycrKxsVFuQgBQjGKoaPbbIH+qScgCLxIXF3fx4sVnVcKia4KAkfnWmMIJoDUKTo0SQnbu3ElR1IgRI0xMTNzd3Xv37j1gwICCgoLQ0FDlxgOA1+Ph4bF48WJjY+OyiqrQE1evX78ukUjoDgXQqSkuhDY2NklJSVu2bBkyZIi2tnb//v2/+eabp0+fenp6KjkfALwuQ0PDhQsXeg3qJ5VKr1+/vn///rKyxGBy4gAAH0ZJREFUMrpDAXReisca7fxea6xRMS+nKiNCz+2jjk7VpWH8w1bRMtboGxNVPElPvnUhTlBZWamiojJ27Fh3d/eODo+xRtsC37VWKXkXKT4ifHuyfqeamprGxsZ79ux5dYXffvvN1dVVVVXVxMRk3bp18rM3Y8eOtf6f4ODgDooH0BNYmukvWTjPzcWxob4u4tzZwwf3V5YVSUV8qYgvFdfTnQ6gs1B8Q31tbe1PP/30xx9/5OfnN5kXu7y8vC3b/eyzz4RCYUlJSWpq6siRI729vZ2cnBqvIJFI9u7dO2DAgMzMTD8/PwsLiw8++IAQkpeX9/3337u7uxNC8HclwBujWFxx7QtSu3usPbFg8S/EPH12PzYr4ZSvt52bg6lU3MDpM5Wl49T6hgC6O8WFcPbs2RERERMmTPD19WWxWhqPVCGJRHLw4MHIyEg1NbV+/fpNmjQpNDT0hx9+aLzO0qVLZQ8cHBz8/f0fPXokX2RqamplZfW6bwoAjTE1rTVcP5U99nAlDr61kZGRT58+vZRMskUW/p5sqaSB3oQAnYSCIicQCKKiojZv3rxixYo322hxcXFlZaX8NkQXF5fbt283t7JAIIiOjl6zZo28JTg4WCwWDxw48LvvvnN2dm7uhVKptKKiQv5UR0enq1y8AVA+DQ2NGTNmPH78+Pz586mpqf3MjG31cCsUACEKC2F1dbVIJBo1atQbb1RWn+TDc2tpaTXXaU0qlS5btszExCQkJETWsmXLFmdnZ6lUunXr1jFjxqSkpCi8hZ/P5xcUFDQ+cFyxYkVzlVvK5wkEgpqamjf+RD1BbW0t3RE6u7q6OrFY3KX/3urbt+/cuXNzc3ONNdN5PB6T3c5fCllnGcxm2jJ811rVjrtITU2NzWa3vI6CQmhgYODo6JiQkDBgwIA3e2MDAwNCSHV1tayGVVRUGBkZKVxz9erVCQkJV65cYTD+6rbj5+cne/DTTz+dOnXq1q1bAQEBr76Qw+GYmpq2tdcoQ12kqopuWq3CLmoZg8HoQr1Gm6OpqWliYlKfnc9UV2e39784eo22Eb5rraJ5hnqKog4cODBv3jx9ff1x48apqqq+7kb19fUNDAweP348YsQIQsjjx4/t7e1fXW3jxo2XL1++evWqlpaWwu1QVFe9uwOgs6OYDS8vNxTcULyQweLYLaSYr/3dB+iKFHeEmTt3bk5OzsSJExkMRpOR7NvSa5TBYLz//vtffvnliRMnnjx5EhERER8fTwjJy8v717/+deTIEXV19W+++Wbnzp1hYWFZWVlZWVk6OjrW1ta5ubn379/38vKSSCTbtm1raGgYNmxYu3xOAGhM1Xyc1LjZL1dd6l4iaSAohNAzKC6EISEhb3mK9vPPP1++fLmTk5O+vv6ePXtkI5SKRKKioiLZQV5qaqqlpeX69etl648YMWLLli0SiWT79u3Lli1jsVgDBw68dOmSrq7u28QAAIUohgqlotKk8fTp02VlZePHj9elOuoOY4BOqKuee8TIMu0Oo120qmuNLPMGDh48+OLFC4qiXAxzJ4R8x9E0fN0t4BphW+C71qpuMrIMAHQ5s2fPHjFiBIPBiE/K3f7b7sePH3fRP5QBXstr3ywPAN0Vi8UaPXq0q6vrmd3J+TVVZ06deBgXO2G8n6GhwT/WY7ApBn51QPeBn2YA+AdDQ8M50955/CQ5+m7y87jYjIcnB/Wz9Hbvq8JmEkKIVMJQM+TaLaQ7JkC7QSEEgKbU7Rd52ZP+/vzo6Oi4uLj4l9K0GrWxY8e6uLhI6ov5L07RHRCgPeEaIQAoxuFw/P39Fy5caG5uXl1dferUqdDQ0KLiErpzAbQzFEIAaImZmdn7778/ceJEdXX1rKys3XsPXYxJrqurozsXQLtBIQSAVlAUNWDAgGXLlg0ZMoQQEv8k69dff3369CnduQDaBwohALSJmpqan5/fkkXz+vY25PP5iYmJdCcCaB/oLAMAr8HI0GBWQP9y9mAdHS1xTWaTpRKRSMrSJbihHroUFEIAeA0US4PBMdKXpJFqIqhuulTEL6c0bIimghljADotFEIAeA0Ui8uxCm5uaV3+TbGgUpl5AN4erhECQLuRSqU37z25d++eSCSiOwtAW+GIEADaTb2g4fb9ZHZqw507d0aNGtW/f3/5nNsAnRZ+RgGg3XDUVGdN9jExMamurj537txvv/2WlJSEkbuhk8MRIQC0G4rJMeNkzfJkPc0ovhmXlZ9UdyTplKGexnCPPnZ9DCiKcO3eZ6r3pjsmwD+gEAJAu2Hp9eNou6ipqQ0eSDymiBMTE2NiYiorKyMeE5MSEy/rGiexgO6MAE2hEAJAh2AymQMHDuzXr9/Dhw9v3rxZUFBwIuWxeUq9z7jpdnZ23Xh+Y+hyUAgBoAMxmUxPT88BAwbExcVd++/z/ILio0ePGhsbjxw50sHBAeUQOgMUQgDocCwWa8iQIU66qUkvNe4lZBcWFh4/ftzIyMjHx8fR0ZHudNDToRACgJKw2SpuRi+cfZmJT2vuPsx8mVwflvznv98bzVVTIYQQpirXfgmOEUH5UAgBQEnULKdIxXWEkGEOxOtdcVLy0/p6gb6ru2wpL3kbiiDQAoUQAJSEYqpSTFXZYwYhAzyG05sHQAaFEAA6CzEvn/zvqFAqlTY+TUqxtRhsDZpyQTeHQggAnQJLs68gL1L2OCuv9Ni5WHMTXU+3vnZWxkTawFDR5VjPoTchdFcohADQKXBs5sof62jmcx8w8uvr/3tfqptW7+Fi7mRSxaExHHRrKIQA0OmYmpquWLEiISEhNja2vLz80tXMqw0F7j69PD09jYyM6E4H3Q0KIQB0RqqqqoMHDx40aNDz58/vXP8zIyX3wYMHDx486NOnj6enp4ODA5PJpDsjdBMohADQeVEUZW9vb2PKfvlYnJDHeZz0PPNpQebTuxoa3P6u9gPcHHS0NQkhTE1rhooO3WGhq0IhBIDOjqFm2MvCeWxv6Uh386SnWXGJaSVlRTE3im7GxFj3Melvo+kw4B01E9yMAW8IhRAAOjuKranaO4AQokaItw3xDiTZ2dnx8fEpKSlZ5aKMq5lacWULl7no6urSnRS6JBRCAOh6LC0tLS0t/fz8EhMTY6+GV9bW8/l8FEJ4MyiEANBVcblcLy+vgRa1/KpsDkkR5Ke8ug5FMdi9hlEMFeXHg64ChRAAujaWrguX1exNhpU5NxlSUwNTB2VGgq4FhRAAujYm14TJNWlu6dHffisXlxmbWbu6urq4uOjooHMpNIVCCADdmYu9Wexz0f+3d+dRTZ15H8CfJGQlCQQCZGEXWUQWFa0s1YBL3WhV6mvFseqxrYy26hxPmTOeMuNpa890anXaTqdHz9se/WPs69LltYooKu4bokhFiQSQhIQEQsiekPX9I23KsOo7JDdwf5+/kuc+Sb7cc8kv9z73PlfZIVF2SKqrfoyJFkydkpqamsRiMhFCBBKNSA3DOiPAGBRCAMBElv+iaFa2rFXa/fCJvLlN2dYobmusOYUIscKwtEn8lDhW5KxdWGcEGINCCACYyGgxS2kIZU1BWYuQzWYTi8WNjY0SiURhdSh+tp+9dDnpkSAtLS01NTUkJATrsAAbUAgBAHhBoVAyMjIyMjL6+vrEYvHDhntNt2pbm2pbm2pPIyTgR6YmJ6YmJ4aH/TKOSKRHkmgR2GYGfgCFEACAO1QqNTMzMzMjXTc7uLlFKpbIW54qO9q6O9oaL1YTNq97KZzDcvVpiTQuKW4F1mGBz0EhBADgFYEUkrQqJwnlvITsdntLS0tTU5PBYAhPW0Wj0ey9DU6dBOuIwB+gEAIAACKTyampqamp/3a5ocuus/c2ep/abHad3sgNDyUQCAQCgcSaRCBR/Z4UjD0ohAAAMAQSLdJJZjp1v81W87+VNx896QhmUBNjo+IiiSmzVnKiZ2OYEIwVKIQAADAEIp1Hi1/VvyU9P7nLVqPVah8r0c/ipso7h/lx9xMTExMTE+Pi4igUmMVtvIJCCAAAzyQrKysrK0utVkskkqbbOpnaoVKpVCrVzZs3SSSSUChMSEiIj4+Pjo4mk8lYhwXPAQohAAA8By6Xy+Vys/mdLgJDaWS0tclan0rlClW7RN0ueYAQIpGChEJ+YsqM2NjY+Ph4EomEdWQwCiiEAADw3EjBMS7tIx4F8VJQboqgzxYplaulcnV7h1ql1rbU329vkyACaXp2+tKX5g54rctGQiyYBDyAQCEEAIDnRubmkLk53qcMhDhTUBZCCCGr1dpc+2370xaFUhMTZu1/ug1CyOWwOo09KAoKYQCBQggAAGOJRqNlvLgx48Whl7osKqv4X3ar4V9H/sfldkULhUKhQCgQhoSwPR0IJCoiEP0XF0AhBAAAfyKQaMhh1j38XNpwzWpzSOp+aQ+mU3gRLEEEMyZ1blzmchaLhWlMfIFCCAAA/kOghJBTtrNYrPKpVplMJpfLOzo6FAqF2WyWWdHTJhl6VE2qbmEymXw+n8fj8fn8qKiosLAwAoGAdfYJCwohAABggEajTZ48efLkyZ6nvb29CoVC+rhG1nRJ1XNH02XXtCLvrDYUMikqnLWgICVh9tYR7kIM/n+gEAIAAPY4HA6Hw0lPn+J2bEIIabW6TmVnZ6dSqVQpVSqD3qByILmGKFTWEIKGPmpKCKJRBQv8m3qCgEIIAACBg0AIoiOEOFw6h8ubMvWXVovFotVquUwbsqkHv8Zssf734e8d2kfCLJXnMsfw8HAulwsDjc8ICiEAAAQ6Op1Op9MRQgjFDV4aZDK5KDd0hnvGBxfF/dqpFHIYJySMww7nsLn8xIiY7LCwMAaD4Z/M4wgUQgAAGN+Cg4N37typepKu7lb1aPQ9vXp1j7anV2+xWDo7LZ2dSrfT6nZWB7ESEUJUKpUTGhIeFlqQNzMqkut5ByJDiOc7aUAhBACAcY9IJPJTl/D//TJ9s9ms0Wh6enrUna2dzRd6dXqtzmQy2k09qKMFMVztRflTEULuPg0lcjY54rc7aSiVSrPZzGazQ0JC8DBvKhRCAACYmBgMBoPBiI6ORllZaNEKT6PJZOrt7TUYDImJiVQqFSHUp6juU92waxo8HZxO11cHqxxO5y9vQqeFx+WyQ8LZbDabzWKz2SwWk8VksdksKpWKiGQCcdzXkXH/BwAAAHh2wcHBwcHB/VsoUQVBoen9WwqXRLa1y3Q6vd5gtPSKpY1nhnwrShCRzWaHJS1lsViZmZneS0HGHSiEAACAawQSncSg928RLVwp+vWx0WjU6/UGg0Gr1RoMBoPBoNPpjEajTqez9ZnU3bd7LRcRQuqnN2LXLB7wzkRqGJES+lgsedwkCQ5mMOg0JjOYTqcxGHQ6jc7iJjGY7ECYKAAKIQAAgGExmUwmkznkIqvV0vP0ismgNxjNQkEEgUTrv9Tt7LP3PkQI3b50vVXaNeC1boeJSGYiIo1Bp9JoFAadSqdR6DQKjUql0yihHNbM+Vv8NueqDwthXV3dvXv30tLSCgoKhuygUCiqq6uZTOaSJUt+PTMYIYQuXrzY2to6a9aszMxM38UDAADwn6DR6MLUl0bttrpsZWtrq9lsNplMJpPJbDabzWajttOk77JYrFYnspqQ1tSHUJ/3JX294pSZ/8UOjfRl/N/4qhDu27dv3759y5cv//jjj0tKSj7++OMBHerr64uKipYvXy6Tyfbs2XP9+nVPLSwrK7t8+fK8efMqKio++OCDN954w0cJAQAA+AGTyRxur8btdlssFrPZbPl3LtkxJjN4yJf4AsHtdo/5m5pMJqFQeP78+ZycHKlUmpKS0tLSIhAI+vd59dVX09LSPvjgA5fLlZeX9+abb27atEkikWRlZbW1tUVGRtbU1JSWlra3t1MolMEfIZPJ8vPzpVLps+RxmqS6llNhmVvG5s+boAwGA8xDMTKTycRgMAJhSCNg2e12p9NJo9FG74pj8L82Ks39T8jINNzSIFYCPWn9GH6cT/YIr127FhoampOTgxCKjY3Nzs4+d+7chg0bvB3cbvfp06ffe+89hBCRSFy+fPnp06c3bdpUWVmZn58fGRmJEBKJRHa7va6uLjc31xchAQAABCZyUtlwvxWc5s4+2cmx/TifFEKFQiEUCr1PhUKhXC7v30Gj0VitVm8fgUCgUCgGvJBAIPD5/AEv9HI4HEaj8aOPPvK2zJs3b/r06UN2djkcTqfTbrf/B3/TxGe322EVjcyzimCPcASePUISiYR1kIAG/2ujGmEVOe325/o+J5FIROIoJ934pBA6nc7+XxZEItH567WZ3g4IIW8fEonkcDie5YVebrfb6XRqNBpvi1qtHrYzKRSFZA23FHg4nU5YRSPzrCIohCNw/grrIAENVtGoRlhF7iA2gZPz7Ctw1CqIfFQI+Xx+V9dvJ8uqVKqFCxf278DlcslkcldXF5fL9XTwjCDy+fzHjx/3f+GAkUUvMpkcEhKyd+/eZ0tEc7inwbjFyOx2O6yikXlGv6AQjoBEIsEY4ajgf21UI64iGgqeNbYf55OrNHJzczs6OiQSCUKot7e3trZWJBIhhCwWi2cfjkgkikSiqqoqT/+qqqqioiKEUGFh4dWrV00mE0Kovr7earUOd7QTAAAAGBM+2SMMCwvbsmXLihUrNmzYcOLEiZUrVyYlJSGEDh06dODAgfr6eoTQn/70pxUrVphMpvb2drFYfPToUYTQtGnT5s6du3Tp0uLi4oMHD+7cuXPAVEAAAADA2PLVdfuffPLJ7t27NRrN22+/ffjwYU9jUVHR+++/73lcWFh46dIll8s1ZcqU2tpaDofjaT9x4sSGDRs0Gs3evXsrKirGJIzRaLx69eqYvNUE5t1BB8O5detWb28v1ikC2tOnTxsbG7FOEdCsVmtNTQ3WKQLd2bNnfXFp33B8ch2hHzzXdYTXrl0rLy+/ceOGr1ONa3Q6vbe3F4YuRrB48eKtW7cuW7YM6yCBa+/evTKZ7LPPPsM6SOBqaGhYu3btzz//jHWQgBYREfHo0aOIiAj/fJyfZnIDAODBOP1hDXAOCiEAAABcg0IIAAAA18brGGFLS0tGRkZ+fv6zdNbpdM3NzZ4p38Bwampq5s6d+ywXn+JWfX19TExMeHg41kECl1QqtVqtycnJWAcJXEaj8eHDh7Nnz8Y6SEC7cuVKbm4umUz+z99qxYoVW7aMMtH0eC2ELpfr8OHDMTExz9LZ4XB0dnY+Y2fcamtrS0hIwDpFQJPL5Vwul0qlYh0kcOn1ervdDr8VRuByuWQyWVxcHNZBAtoYfh0lJCRMmjRp5D7jtRACAAAAYwKOgwEAAMA1KIQAAABwDQohAAAAXINCCAAAANd8Muk2hlwuV3Nzs0wmKygoGG62MLFYfOPGjfj4eJFIhM9b6qjV6sbGxujo6CFPppJKpU+ePPE+zcvLYzAYfkyHPZfLVVtbKxaLuVxuUVHRcBvS1atXJRLJzJkzp06d6ueEgUAmk925c8dut8+ePTs+Pn5wh/r6erVa7XlMpVJffPFFv+YLAM3NzQ0NDQaDIS0t7YUXXhiyj1qtPnfuHI1GW7RoEd7+0RBCWq32zp07crk8IiKisLBw8F0WDAbD7du3vU/T0tL63/V9zLgnEIVCwWazPfc4fPr06ZB9jh49yuVy33rrrfT09HXr1vk5YSB47bXXqFQqi8X685//PGSHffv28fn8+b+SyWR+Toi5ZcuWTZ06df369Xl5eQkJCXK5fHCfrVu3Tp48efPmzVFRUV9//bX/Q2Lr+++/Dw8PX7ly5Zo1a9hs9jfffDO4T3FxcXp6umcreu211/wfEnNxcXGrV6/euHFjbGzs6tWrXS7XgA5isTgiIqK0tHThwoXp6ek6nQ6TnBjatm3bwoUL33jjjcLCQh6P9+TJkwEd7t+/T6PRvF9Hp06d8kWMCVUI+/r6pFKp2WwerhA6nc5JkyZ99913brdbo9FwOJwHDx74PSbG2tvbbTZbaWnpCIVw/fr1/g0VWCQSiffx/Pnzd+3aNaBDS0sLg8Ho7Ox0u90XLlzg8Xg2m82vEbGmUCiMRqPn8ZEjR6Kiogb3KS4uPnTokH9zBSiVSkUikcRi8YD2jRs3vvPOO2632+VyFRYW/v3vf8ciXaBYvnx5eXn5gMb79+/Hxsb6+qMn1BghhUIZ+ar5x48fy+Vyz90DOBxOUVHRqVOn/JUuUMTGxo46X4NGozlz5kx9fb3L5fJPqoDS/4gxj8fr6+sb0KGysjIvL4/H4yGECgsL7XZ7bW2tXyNijc/ne49i8fl8z++Awd1aWlqqqqra2tr8my7gWCyWoKAgNps9oP3UqVOvvvoqQohAIJSUlODw68jL7XabzWbP8bwB7HZ7dXX1zZs3PTs5vjChCuGoFApFZGQkhULxPBUKhXK5HNtIAYhIJCoUiq+++urll1/Oz8/H8x34Ghsbf/rpp9dff31Au1wuj46O9jwmEAh8Ph+3G5LL5frwww83bdo0eLidRqNdvXr1888/z87O/v3vf49JPMzt2LGjqKgoPz//yJEjnl9OXjabTa1We0e8cPt1dPbs2QULFqSmpkZFRW3btm1wBwaD8cUXX5SVlaWkpNy9e9cXGfBVCJ1OZ/9/VxKJ5HA4MMwTmN5+++27d++ePHmyubmZQqHs2bMH60TYUCgUr7zyyp49ezIzMwcsGrAhBQUF4XZD2r59u9Vq9d5wu79vv/22pqamsrLy4cOHx48fr6ys9H88zP3ud7/buXPnkiVLKioqDAZD/0VOp9Ptdns3JNx+HWVkZJSXl2/btu38+fPV1dWDl0okkpMnTz548GDt2rVlZWW+yICvQsjn89Vqtfdwn0ql4vP52EYKQCQSyfOASqWWlJTU19djmwcTKpVq3rx5b7755tatWwcv5fP5XV1d/TsLBAI/pgsU77777q1bt06fPk2n0wcv9W5IMTExBQUF9+/f92+6gJCTk7N06dKDBw8SicQff/yx/yI6nR4aGtrd3e15itutSCAQLFiwYOvWreXl5Z9++umApd6tCCFUWlra0NDgi/EaXBRCk8lkMpkQQmlpaUwm8/r16wghm8126dKlwsJCrNMFBIfD0dPTM7i9rq4Oh5OVd3d3z58/v7S09I9//GP/dq1Wa7PZEEIikejatWueEYuGhgaTyTRjxgxssmLnvffeq66urqqqCgkJ8TZaLBa9Xj+gZ19f38OHD2NjY/0bMIDYbDa9Xu9ZUTabTavVetpFItG5c+c8j8+dOycSibBKGAjUarV3GFWj0QzeP66rqxMKhb64Q85Eu45w+/btRqMRIbRr1y4mk/mPf/yDTCbv2LHD5XJ9/fXXFArl3XffXb9+/TvvvHP+/PlJkybNmTMH68j+dvz48fPnz9+5c6epqUmpVK5Zs0YkEt29ezc3N9dzvsPq1avj4uIiIiJu37594cKFmzdvYh3Z39atW6dUKjs6OjZv3owQ8g5xzZgx4/3331+7du20adMKCgqWLVv28ssvHzx4cMeOHUwmE+vUfnXs2LE9e/asXLly165dnpZ9+/YFBwfv37+/urq6pqZGp9MVFxcXFRVRqdQffviBzWavWrUK28x+dvny5f379+fk5BAIhJ9++onH4y1cuBAhdPz48YqKitbWVoRQeXn5okWLEEJqtfry5cufffYZxqH9bvHixdnZ2eHh4Y8ePTpx4oT3+LlAIDhz5kxhYeGHH34ok8mSkpKkUunhw4f/+c9/+iIGaffu3b54X6x4jnYWFxfHx8cLBIIZM2YQiUQOh5Odne2570leXt7kyZObmppycnI+/fTTMbnf1fii1+tJJFJ+fn5OTo5AIEhNTfWcQJScnOzZrREIBN3d3VqtNjs7+8CBAzj8Ic9gMHJzcwW/SkxM9JxHyuPxZs2aFRYWhhAqKSlxuVwdHR2vv/46Ps8EycrKSkpK8q6l6dOnBwUFsVisjIyM5ORkMpkcGhra1dVlt9uXLVu2f/9+vN2+KiIigkajKZVKl8v1yiuv/O1vf/OsgeDg4ClTpngmYYiJiVm6dGljY2NoaOiXX37pk0vFA1tSUlJXV5dGo0lOTv7yyy+9c1MIhcLZs2ezWCyBQGAwGDzHjf/617/Onz/fFzHgNkwAAABwDRdjhAAAAMBwoBACAADANSiEAAAAcA0KIQAAAFyDQggAAADXoBACAADANSiEAAAAcA0KIQAAAFyDQggAAADXoBACAADANSiEAIxv33zzTVhY2OnTp70t27dvFwqFzc3NGKYCYByBuUYBGN/cbndJScnly5fv3bsXFxd34sSJVatWHThw4K233sI6GgDjAxRCAMY9rVY7bdq0mJiYQ4cOzZw5c968eceOHcM6FADjBhRCACaCW7duzZkzh06nR0ZG1tXVeW9wCgAYFYwRAjARvPDCCzNnztTr9Tt37oQqCMBzgT1CACaCL774Yvv27enp6d3d3fX19TweD+tEAIwbsEcIwLjX0NBQXl7+hz/84eLFi0FBQWvWrHE6nViHAmDcgD1CAMY3o9GYk5PDYrGuX79OoVCuXLlSVFT0l7/8paKiAutoAIwPsEcIwPhWVlamVCqPHj1KoVAQQnPmzKmoqNi9e/fFixexjgbA+AB7hAAAAHAN9ggBAADgGhRCAAAAuAaFEAAAAK5BIQQAAIBrUAgBAADgGhRCAAAAuAaFEAAAAK79Hz/uhfvhPSS6AAAAAElFTkSuQmCC", + "image/svg+xml": [ + "\n", + "\n", + "\n", + " \n", + " \n", + " \n", + "\n", + "\n", + "\n", + " \n", + " \n", + " \n", + "\n", + "\n", + "\n", + " \n", + " \n", + " \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n" + ], + "text/html": [ + "\n", + "\n", + "\n", + " \n", + " \n", + " \n", + "\n", + "\n", + "\n", + " \n", + " \n", + " \n", + "\n", + "\n", + "\n", + " \n", + " \n", + " \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "plot_compare(samples,d)" + ] + }, + { + "cell_type": "markdown", + "id": "bbf15f20-6821-4503-95d9-7a4859be1bbe", + "metadata": {}, + "source": [ + "## Some Benchmarking" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "id": "a282af31-95a0-4c39-904b-de34218532d3", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "BenchmarkTools.Trial: 10000 samples with 1 evaluation per sample.\n", + " Range \u001b[90m(\u001b[39m\u001b[36m\u001b[1mmin\u001b[22m\u001b[39m \u2026 \u001b[35mmax\u001b[39m\u001b[90m): \u001b[39m\u001b[36m\u001b[1m364.083 \u03bcs\u001b[22m\u001b[39m \u2026 \u001b[35m 12.503 ms\u001b[39m \u001b[90m\u250a\u001b[39m GC \u001b[90m(\u001b[39mmin \u2026 max\u001b[90m): \u001b[39m 0.00% \u2026 96.29%\n", + " Time \u001b[90m(\u001b[39m\u001b[34m\u001b[1mmedian\u001b[22m\u001b[39m\u001b[90m): \u001b[39m\u001b[34m\u001b[1m436.125 \u03bcs \u001b[22m\u001b[39m\u001b[90m\u250a\u001b[39m GC \u001b[90m(\u001b[39mmedian\u001b[90m): \u001b[39m 0.00%\n", + " Time \u001b[90m(\u001b[39m\u001b[32m\u001b[1mmean\u001b[22m\u001b[39m \u00b1 \u001b[32m\u03c3\u001b[39m\u001b[90m): \u001b[39m\u001b[32m\u001b[1m486.206 \u03bcs\u001b[22m\u001b[39m \u00b1 \u001b[32m306.643 \u03bcs\u001b[39m \u001b[90m\u250a\u001b[39m GC \u001b[90m(\u001b[39mmean \u00b1 \u03c3\u001b[90m): \u001b[39m10.54% \u00b1 14.32%\n", + "\n", + " \u001b[39m\u2585\u001b[39m\u2584\u001b[39m\u2587\u001b[34m\u2588\u001b[39m\u001b[39m\u2586\u001b[32m\u2584\u001b[39m\u001b[39m\u2582\u001b[39m\u2581\u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m\u2582\n", + " \u001b[39m\u2588\u001b[39m\u2588\u001b[39m\u2588\u001b[34m\u2588\u001b[39m\u001b[39m\u2588\u001b[32m\u2588\u001b[39m\u001b[39m\u2588\u001b[39m\u2588\u001b[39m\u2587\u001b[39m\u2585\u001b[39m\u2583\u001b[39m\u2583\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2581\u001b[39m\u2583\u001b[39m\u2583\u001b[39m\u2581\u001b[39m\u2586\u001b[39m\u2585\u001b[39m\u2585\u001b[39m\u2586\u001b[39m\u2585\u001b[39m\u2585\u001b[39m\u2586\u001b[39m\u2585\u001b[39m\u2587\u001b[39m\u2587\u001b[39m\u2587\u001b[39m\u2588\u001b[39m\u2588\u001b[39m\u2588\u001b[39m\u2587\u001b[39m\u2587\u001b[39m\u2587\u001b[39m \u001b[39m\u2588\n", + " 364 \u03bcs\u001b[90m \u001b[39m\u001b[90mHistogram: \u001b[39m\u001b[90m\u001b[1mlog(\u001b[22m\u001b[39m\u001b[90mfrequency\u001b[39m\u001b[90m\u001b[1m)\u001b[22m\u001b[39m\u001b[90m by time\u001b[39m 1.73 ms \u001b[0m\u001b[1m<\u001b[22m\n", + "\n", + " Memory estimate\u001b[90m: \u001b[39m\u001b[33m1.04 MiB\u001b[39m, allocs estimate\u001b[90m: \u001b[39m\u001b[33m507\u001b[39m." + ] + }, + "execution_count": 6, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "N_bench = 10000\n", + "\n", + "@benchmark generate_CPU($RNG,$d,$N_bench)" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "id": "37a1aecb-27f9-493b-8a62-4e27da048f72", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "BenchmarkTools.Trial: 10000 samples with 1 evaluation per sample.\n", + " Range \u001b[90m(\u001b[39m\u001b[36m\u001b[1mmin\u001b[22m\u001b[39m \u2026 \u001b[35mmax\u001b[39m\u001b[90m): \u001b[39m\u001b[36m\u001b[1m116.000 \u03bcs\u001b[22m\u001b[39m \u2026 \u001b[35m 21.065 ms\u001b[39m \u001b[90m\u250a\u001b[39m GC \u001b[90m(\u001b[39mmin \u2026 max\u001b[90m): \u001b[39m 0.00% \u2026 99.09%\n", + " Time \u001b[90m(\u001b[39m\u001b[34m\u001b[1mmedian\u001b[22m\u001b[39m\u001b[90m): \u001b[39m\u001b[34m\u001b[1m136.666 \u03bcs \u001b[22m\u001b[39m\u001b[90m\u250a\u001b[39m GC \u001b[90m(\u001b[39mmedian\u001b[90m): \u001b[39m 0.00%\n", + " Time \u001b[90m(\u001b[39m\u001b[32m\u001b[1mmean\u001b[22m\u001b[39m \u00b1 \u001b[32m\u03c3\u001b[39m\u001b[90m): \u001b[39m\u001b[32m\u001b[1m174.068 \u03bcs\u001b[22m\u001b[39m \u00b1 \u001b[32m504.927 \u03bcs\u001b[39m \u001b[90m\u250a\u001b[39m GC \u001b[90m(\u001b[39mmean \u00b1 \u03c3\u001b[90m): \u001b[39m20.28% \u00b1 8.61%\n", + "\n", + " \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m\u2581\u001b[39m\u2588\u001b[39m\u2587\u001b[34m\u2584\u001b[39m\u001b[39m\u2581\u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[32m \u001b[39m\u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \n", + " \u001b[39m\u2586\u001b[39m\u2583\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2583\u001b[39m\u2588\u001b[39m\u2588\u001b[39m\u2588\u001b[34m\u2588\u001b[39m\u001b[39m\u2588\u001b[39m\u2587\u001b[39m\u2585\u001b[39m\u2584\u001b[39m\u2584\u001b[39m\u2583\u001b[39m\u2583\u001b[39m\u2583\u001b[39m\u2583\u001b[39m\u2583\u001b[39m\u2583\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[32m\u2582\u001b[39m\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2581\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2581\u001b[39m\u2582\u001b[39m\u2582\u001b[39m\u2581\u001b[39m\u2582\u001b[39m \u001b[39m\u2583\n", + " 116 \u03bcs\u001b[90m Histogram: frequency by time\u001b[39m 244 \u03bcs \u001b[0m\u001b[1m<\u001b[22m\n", + "\n", + " Memory estimate\u001b[90m: \u001b[39m\u001b[33m340.20 KiB\u001b[39m, allocs estimate\u001b[90m: \u001b[39m\u001b[33m18\u001b[39m." + ] + }, + "execution_count": 7, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "@benchmark ConstrainedGaussians.generate_batch_CPU($RNG,$N_bench,$d)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "b9435669-be67-4768-a481-d731325d50fe", + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Julia 1.11.2", + "language": "julia", + "name": "julia-1.11" + }, + "language_info": { + "file_extension": ".jl", + "mimetype": "application/julia", + "name": "julia", + "version": "1.11.2" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/examples/example.jl b/examples/example.jl new file mode 100644 index 0000000..522a3f5 --- /dev/null +++ b/examples/example.jl @@ -0,0 +1,31 @@ + +# some imports + +using BenchmarkTools +using Random +using Plots +RNG = MersenneTwister(1234) + +# import the example +using ConstrainedGaussians + +# setup + +mu = rand(RNG) # central value +sig = rand(RNG) # variance +dom = (mu + 0.5 * sig, mu + 3.0 * sig) # support/domain + +## example distribution to be sampled +d = ConstrainedGaussian1D(mu, sig, dom) + +## number of samples to be generated +N = Int(1e6) + +## generation of samples +samples = generate_CPU(RNG, d, N; batch_size = 1000) + +## plotting: histogram vs target_dist +P = plot_compare(samples, d) + +## save plot +savefig(P, "example_compare.pdf") diff --git a/examples/init.jl b/examples/init.jl new file mode 100644 index 0000000..07e547f --- /dev/null +++ b/examples/init.jl @@ -0,0 +1,7 @@ +import Pkg + +inclusive_scans_jl = Pkg.PackageSpec(path = "..") +constainted_gaussians_jl = Pkg.PackageSpec(path = "./ConstrainedGaussians") + +Pkg.develop([inclusive_scans_jl, constainted_gaussians_jl]) +Pkg.instantiate()