From 7ae0e901a5f3b370c3cec6eb27640fedaf0c2a8d Mon Sep 17 00:00:00 2001 From: Francesc Alted Date: Sun, 19 Jul 2026 14:32:42 +0200 Subject: [PATCH 1/8] Add blosc2.random: chunk-parallel NumPy-quality random constructors Spawns one SeedSequence stream per chunk and generates chunks concurrently in a thread pool (NumPy Generator fill loops release the GIL), giving full PCG64 quality with ~3.2x measured speedup over asarray(np.random...) while staying seed-reproducible regardless of thread count. Covers default_rng().{random,integers,normal,uniform} for now. --- RELEASE_NOTES.md | 9 +++ doc/reference/array_operations.rst | 1 + doc/reference/random.rst | 14 ++++ src/blosc2/__init__.py | 1 + src/blosc2/random.py | 118 +++++++++++++++++++++++++++++ tests/test_random.py | 96 +++++++++++++++++++++++ 6 files changed, 239 insertions(+) create mode 100644 doc/reference/random.rst create mode 100644 src/blosc2/random.py create mode 100644 tests/test_random.py diff --git a/RELEASE_NOTES.md b/RELEASE_NOTES.md index fcd5d826e..3578a1349 100644 --- a/RELEASE_NOTES.md +++ b/RELEASE_NOTES.md @@ -4,6 +4,15 @@ XXX version-specific blurb XXX +### New features + +- New `blosc2.random` module: seedable, NumPy-quality random `NDArray` + constructors (`default_rng`, `Generator.random/integers/normal/uniform`). + Each chunk gets its own independent `SeedSequence`-spawned stream and is + generated concurrently in a thread pool, giving full `PCG64` quality with + genuinely parallel generation (measured ~3.2x faster than + `asarray(np.random.default_rng(...).random(...))` on a 100M-element array). + ## Changes from 4.9.0 to 4.9.1 A small hot-fix release for the Arrow interop work in 4.9.0: a real diff --git a/doc/reference/array_operations.rst b/doc/reference/array_operations.rst index 2c2adb15b..b40ed4310 100644 --- a/doc/reference/array_operations.rst +++ b/doc/reference/array_operations.rst @@ -7,5 +7,6 @@ Operations with arrays ufuncs reduction_functions linalg + random additional_funcs index_funcs diff --git a/doc/reference/random.rst b/doc/reference/random.rst new file mode 100644 index 000000000..413c00b1a --- /dev/null +++ b/doc/reference/random.rst @@ -0,0 +1,14 @@ +Random +------ +Chunk-parallel, NumPy-quality random :ref:`NDArray ` constructors. + +.. currentmodule:: blosc2.random + +.. autosummary:: + + default_rng + Generator + +.. autofunction:: blosc2.random.default_rng +.. autoclass:: blosc2.random.Generator + :members: diff --git a/src/blosc2/__init__.py b/src/blosc2/__init__.py index 3050887bb..cea161a22 100644 --- a/src/blosc2/__init__.py +++ b/src/blosc2/__init__.py @@ -595,6 +595,7 @@ def _raise(exc): from .linalg import tensordot, vecdot, permute_dims, matrix_transpose, matmul, transpose, diagonal, outer from .utils import linalg_funcs as linalg_funcs_list from . import fft +from . import random # Registry for postfilters postfilter_funcs = {} diff --git a/src/blosc2/random.py b/src/blosc2/random.py new file mode 100644 index 000000000..5282f2656 --- /dev/null +++ b/src/blosc2/random.py @@ -0,0 +1,118 @@ +####################################################################### +# Copyright (c) 2019-present, Blosc Development Team +# All rights reserved. +# +# SPDX-License-Identifier: BSD-3-Clause +####################################################################### +"""Chunk-parallel NumPy-quality random generators producing :class:`blosc2.NDArray`. + +Each chunk gets its own independent stream, spawned from a root +:class:`numpy.random.SeedSequence` via :meth:`~numpy.random.SeedSequence.spawn`, and +generation runs in a thread pool (NumPy ``Generator`` fill loops release the GIL). +Results are reproducible for a given ``(seed, call order, shape, chunks)``; changing +the chunk layout changes the values, since the chunk grid determines the stream +assignment. +""" + +from __future__ import annotations + +from concurrent.futures import FIRST_COMPLETED, ThreadPoolExecutor, wait +from typing import Any + +import numpy as np + +import blosc2 + + +def _chunk_slices(a: blosc2.NDArray): + for info in a.iterchunks_info(): + yield tuple( + slice(c * s, min((c + 1) * s, sh)) + for c, s, sh in zip(info.coords, a.chunks, a.shape, strict=False) + ) + + +def _bounded_map(fn, n: int, max_workers: int): + # ponytail: window bounds peak memory to O(workers * chunk); a plain + # executor.map would submit all n chunks upfront and could buffer O(n). + window = 2 * max_workers + with ThreadPoolExecutor(max_workers=max_workers) as ex: + pending = {ex.submit(fn, i): i for i in range(min(window, n))} + next_i = len(pending) + while pending: + done, _ = wait(pending, return_when=FIRST_COMPLETED) + for fut in done: + del pending[fut] + yield fut.result() + if next_i < n: + pending[ex.submit(fn, next_i)] = next_i + next_i += 1 + + +class Generator: + """Chunk-parallel counterpart of :class:`numpy.random.Generator`. + + Create via :func:`default_rng`. Each method call draws a full + :class:`~blosc2.NDArray` in one shot, generating its chunks concurrently. + """ + + def __init__(self, seed=None): + self._root = seed if isinstance(seed, np.random.SeedSequence) else np.random.SeedSequence(seed) + + def _fill(self, method: str, args: tuple, kwargs: dict, shape, **b2_kwargs: Any) -> blosc2.NDArray: + if shape is None: + raise TypeError("shape is required") + dtype = getattr(np.random.default_rng(0), method)(*args, size=1, **kwargs).dtype + dst = blosc2.empty(shape, dtype=dtype, **b2_kwargs) + + slices = list(_chunk_slices(dst)) + call_ss = self._root.spawn(1)[0] + chunk_streams = call_ss.spawn(len(slices)) + + def gen(i): + rng = np.random.default_rng(chunk_streams[i]) + chunk_shape = tuple(s.stop - s.start for s in slices[i]) + return i, getattr(rng, method)(*args, size=chunk_shape, **kwargs) + + if len(slices) <= 1 or blosc2.nthreads == 1: + results = (gen(i) for i in range(len(slices))) + else: + results = _bounded_map(gen, len(slices), blosc2.nthreads) + + for i, buf in results: + dst[slices[i]] = buf + return dst + + def random(self, shape, dtype=np.float64, **kwargs: Any) -> blosc2.NDArray: + """Draw uniform floats in [0, 1). Mirrors :meth:`numpy.random.Generator.random`.""" + return self._fill("random", (), {"dtype": dtype}, shape, **kwargs) + + def integers( + self, low, high=None, *, shape, dtype=np.int64, endpoint: bool = False, **kwargs: Any + ) -> blosc2.NDArray: + """Draw random integers. Mirrors :meth:`numpy.random.Generator.integers`.""" + return self._fill("integers", (low, high), {"dtype": dtype, "endpoint": endpoint}, shape, **kwargs) + + def normal(self, loc=0.0, scale=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a normal distribution. Mirrors :meth:`numpy.random.Generator.normal`.""" + return self._fill("normal", (loc, scale), {}, shape, **kwargs) + + def uniform(self, low=0.0, high=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a uniform distribution. Mirrors :meth:`numpy.random.Generator.uniform`.""" + return self._fill("uniform", (low, high), {}, shape, **kwargs) + + +def default_rng(seed=None) -> Generator: + """Construct a chunk-parallel :class:`Generator`, mirroring :func:`numpy.random.default_rng`. + + Parameters + ---------- + seed: int, array-like, SeedSequence or None + Seed for the root :class:`numpy.random.SeedSequence`. ``None`` (default) draws + entropy from the OS. + + Returns + ------- + out: :class:`Generator` + """ + return Generator(seed) diff --git a/tests/test_random.py b/tests/test_random.py new file mode 100644 index 000000000..f80e7dea9 --- /dev/null +++ b/tests/test_random.py @@ -0,0 +1,96 @@ +####################################################################### +# Copyright (c) 2019-present, Blosc Development Team +# All rights reserved. +# +# SPDX-License-Identifier: BSD-3-Clause +####################################################################### + +import numpy as np +import pytest + +import blosc2 + + +def test_reproducible_same_seed(): + a = blosc2.random.default_rng(42).random((100, 100), chunks=(17, 17)) + b = blosc2.random.default_rng(42).random((100, 100), chunks=(17, 17)) + np.testing.assert_array_equal(a[:], b[:]) + + +def test_successive_calls_differ(): + rng = blosc2.random.default_rng(42) + a = rng.random((1000,)) + b = rng.random((1000,)) + assert not np.array_equal(a[:], b[:]) + + +def test_different_seeds_differ(): + a = blosc2.random.default_rng(1).random((1000,)) + b = blosc2.random.default_rng(2).random((1000,)) + assert not np.array_equal(a[:], b[:]) + + +def test_seed_none_works(): + a = blosc2.random.default_rng(None).random((10,)) + assert a.shape == (10,) + + +@pytest.mark.parametrize("nthreads", [1, 2]) +def test_random_distribution_sanity(nthreads, monkeypatch): + monkeypatch.setattr(blosc2, "nthreads", nthreads) + a = blosc2.random.default_rng(0).random((200, 200), chunks=(31, 31))[:] + assert a.dtype == np.float64 + assert a.min() >= 0.0 + assert a.max() < 1.0 + assert abs(a.mean() - 0.5) < 0.05 + + +def test_normal_distribution_sanity(): + a = blosc2.random.default_rng(0).normal(loc=5.0, scale=2.0, shape=(500, 500))[:] + assert abs(a.mean() - 5.0) < 0.1 + assert abs(a.std() - 2.0) < 0.1 + + +def test_integers_bounds_dtype_endpoint(): + a = blosc2.random.default_rng(0).integers( + 0, 10, shape=(1000,), dtype=np.int32, endpoint=True, chunks=(100,) + )[:] + assert a.dtype == np.int32 + assert a.min() >= 0 + assert a.max() <= 10 + + +def test_uniform_bounds(): + a = blosc2.random.default_rng(0).uniform(-3.0, 3.0, shape=(1000,))[:] + assert a.min() >= -3.0 + assert a.max() < 3.0 + + +def test_blosc2_kwargs_passthrough(tmp_path): + urlpath = tmp_path / "r.b2nd" + a = blosc2.random.default_rng(0).random( + (100, 100), chunks=(20, 20), cparams={"clevel": 5}, urlpath=str(urlpath) + ) + assert a.chunks == (20, 20) + assert a.schunk.cparams.clevel == 5 + assert urlpath.exists() + b = blosc2.open(str(urlpath)) + np.testing.assert_array_equal(a[:], b[:]) + + +def test_nd_shape(): + a = blosc2.random.default_rng(0).random((4, 5, 6), chunks=(2, 3, 3)) + assert a.shape == (4, 5, 6) + + +def test_shape_required(): + with pytest.raises(TypeError): + blosc2.random.default_rng(0).normal() + + +def test_scheduling_independent_of_nthreads(monkeypatch): + monkeypatch.setattr(blosc2, "nthreads", 1) + serial = blosc2.random.default_rng(7).random((300, 300), chunks=(17, 23))[:] + monkeypatch.setattr(blosc2, "nthreads", 8) + parallel = blosc2.random.default_rng(7).random((300, 300), chunks=(17, 23))[:] + np.testing.assert_array_equal(serial, parallel) From 6381ffb2937ae46faa2f3d0f4d8adc762bf371a7 Mon Sep 17 00:00:00 2001 From: Francesc Alted Date: Sun, 19 Jul 2026 14:37:24 +0200 Subject: [PATCH 2/8] blosc2.random: add 30 more Generator distributions Fills out the easy-win half of numpy.random.Generator's surface: every distribution whose numpy signature is (params..., size=None) drops straight into the existing _fill/_bounded_map machinery as a one-line wrapper. Still out: choice/permutation/shuffle/permuted (sequential, not per-chunk independent), the multivariate draws and dirichlet (vector-valued per draw), and bytes (not an NDArray). gamma/standard_gamma rename numpy's distribution-shape parameter to shape_param (positional-only) since it collides with blosc2's array-shape `shape` keyword. --- RELEASE_NOTES.md | 4 +- src/blosc2/random.py | 137 +++++++++++++++++++++++++++++++++++++++++++ tests/test_random.py | 55 +++++++++++++++++ 3 files changed, 195 insertions(+), 1 deletion(-) diff --git a/RELEASE_NOTES.md b/RELEASE_NOTES.md index 3578a1349..2ac21ec78 100644 --- a/RELEASE_NOTES.md +++ b/RELEASE_NOTES.md @@ -7,7 +7,9 @@ XXX version-specific blurb XXX ### New features - New `blosc2.random` module: seedable, NumPy-quality random `NDArray` - constructors (`default_rng`, `Generator.random/integers/normal/uniform`). + constructors, mirroring most of `numpy.random.Generator` (`random`, + `integers`, `normal`, `uniform`, and 30 further distributions such as + `poisson`, `gamma`, `beta`, `standard_normal`, `binomial`, `exponential`...). Each chunk gets its own independent `SeedSequence`-spawned stream and is generated concurrently in a thread pool, giving full `PCG64` quality with genuinely parallel generation (measured ~3.2x faster than diff --git a/src/blosc2/random.py b/src/blosc2/random.py index 5282f2656..f0f7c4ed3 100644 --- a/src/blosc2/random.py +++ b/src/blosc2/random.py @@ -101,6 +101,143 @@ def uniform(self, low=0.0, high=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: """Draw samples from a uniform distribution. Mirrors :meth:`numpy.random.Generator.uniform`.""" return self._fill("uniform", (low, high), {}, shape, **kwargs) + def beta(self, a, b, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a beta distribution. Mirrors :meth:`numpy.random.Generator.beta`.""" + return self._fill("beta", (a, b), {}, shape, **kwargs) + + def binomial(self, n, p, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a binomial distribution. Mirrors :meth:`numpy.random.Generator.binomial`.""" + return self._fill("binomial", (n, p), {}, shape, **kwargs) + + def chisquare(self, df, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a chi-square distribution. Mirrors :meth:`numpy.random.Generator.chisquare`.""" + return self._fill("chisquare", (df,), {}, shape, **kwargs) + + def exponential(self, scale=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from an exponential distribution. Mirrors :meth:`numpy.random.Generator.exponential`.""" + return self._fill("exponential", (scale,), {}, shape, **kwargs) + + def f(self, dfnum, dfden, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from an F distribution. Mirrors :meth:`numpy.random.Generator.f`.""" + return self._fill("f", (dfnum, dfden), {}, shape, **kwargs) + + def gamma(self, shape_param, /, scale=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a gamma distribution. Mirrors :meth:`numpy.random.Generator.gamma`. + + ``shape_param`` is numpy's ``shape`` (the gamma distribution's shape parameter, + usually called *k*); renamed here, and made positional-only, to avoid colliding + with the array-shape ``shape`` keyword. + """ + return self._fill("gamma", (shape_param, scale), {}, shape, **kwargs) + + def geometric(self, p, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a geometric distribution. Mirrors :meth:`numpy.random.Generator.geometric`.""" + return self._fill("geometric", (p,), {}, shape, **kwargs) + + def gumbel(self, loc=0.0, scale=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a Gumbel distribution. Mirrors :meth:`numpy.random.Generator.gumbel`.""" + return self._fill("gumbel", (loc, scale), {}, shape, **kwargs) + + def hypergeometric(self, ngood, nbad, nsample, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a hypergeometric distribution. Mirrors :meth:`numpy.random.Generator.hypergeometric`.""" + return self._fill("hypergeometric", (ngood, nbad, nsample), {}, shape, **kwargs) + + def laplace(self, loc=0.0, scale=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a Laplace distribution. Mirrors :meth:`numpy.random.Generator.laplace`.""" + return self._fill("laplace", (loc, scale), {}, shape, **kwargs) + + def logistic(self, loc=0.0, scale=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a logistic distribution. Mirrors :meth:`numpy.random.Generator.logistic`.""" + return self._fill("logistic", (loc, scale), {}, shape, **kwargs) + + def lognormal(self, mean=0.0, sigma=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a log-normal distribution. Mirrors :meth:`numpy.random.Generator.lognormal`.""" + return self._fill("lognormal", (mean, sigma), {}, shape, **kwargs) + + def logseries(self, p, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a logarithmic series distribution. Mirrors :meth:`numpy.random.Generator.logseries`.""" + return self._fill("logseries", (p,), {}, shape, **kwargs) + + def negative_binomial(self, n, p, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a negative binomial distribution. Mirrors :meth:`numpy.random.Generator.negative_binomial`.""" + return self._fill("negative_binomial", (n, p), {}, shape, **kwargs) + + def noncentral_chisquare(self, df, nonc, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a noncentral chi-square distribution. Mirrors :meth:`numpy.random.Generator.noncentral_chisquare`.""" + return self._fill("noncentral_chisquare", (df, nonc), {}, shape, **kwargs) + + def noncentral_f(self, dfnum, dfden, nonc, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a noncentral F distribution. Mirrors :meth:`numpy.random.Generator.noncentral_f`.""" + return self._fill("noncentral_f", (dfnum, dfden, nonc), {}, shape, **kwargs) + + def pareto(self, a, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a Pareto distribution. Mirrors :meth:`numpy.random.Generator.pareto`.""" + return self._fill("pareto", (a,), {}, shape, **kwargs) + + def poisson(self, lam=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a Poisson distribution. Mirrors :meth:`numpy.random.Generator.poisson`.""" + return self._fill("poisson", (lam,), {}, shape, **kwargs) + + def power(self, a, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a power distribution. Mirrors :meth:`numpy.random.Generator.power`.""" + return self._fill("power", (a,), {}, shape, **kwargs) + + def rayleigh(self, scale=1.0, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a Rayleigh distribution. Mirrors :meth:`numpy.random.Generator.rayleigh`.""" + return self._fill("rayleigh", (scale,), {}, shape, **kwargs) + + def standard_cauchy(self, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a standard Cauchy distribution. Mirrors :meth:`numpy.random.Generator.standard_cauchy`.""" + return self._fill("standard_cauchy", (), {}, shape, **kwargs) + + def standard_exponential( + self, *, shape, dtype=np.float64, method="zig", **kwargs: Any + ) -> blosc2.NDArray: + """Draw samples from a standard exponential distribution. + + Mirrors :meth:`numpy.random.Generator.standard_exponential`. + """ + return self._fill("standard_exponential", (), {"dtype": dtype, "method": method}, shape, **kwargs) + + def standard_gamma(self, shape_param, /, *, shape, dtype=np.float64, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a standard gamma distribution. + + Mirrors :meth:`numpy.random.Generator.standard_gamma`; see :meth:`gamma` for why the + distribution's ``shape`` parameter is renamed ``shape_param`` here. + """ + return self._fill("standard_gamma", (shape_param,), {"dtype": dtype}, shape, **kwargs) + + def standard_normal(self, *, shape, dtype=np.float64, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a standard normal distribution. + + Mirrors :meth:`numpy.random.Generator.standard_normal`. + """ + return self._fill("standard_normal", (), {"dtype": dtype}, shape, **kwargs) + + def standard_t(self, df, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a standard Student's t distribution. Mirrors :meth:`numpy.random.Generator.standard_t`.""" + return self._fill("standard_t", (df,), {}, shape, **kwargs) + + def triangular(self, left, mode, right, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a triangular distribution. Mirrors :meth:`numpy.random.Generator.triangular`.""" + return self._fill("triangular", (left, mode, right), {}, shape, **kwargs) + + def vonmises(self, mu, kappa, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a von Mises distribution. Mirrors :meth:`numpy.random.Generator.vonmises`.""" + return self._fill("vonmises", (mu, kappa), {}, shape, **kwargs) + + def wald(self, mean, scale, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a Wald distribution. Mirrors :meth:`numpy.random.Generator.wald`.""" + return self._fill("wald", (mean, scale), {}, shape, **kwargs) + + def weibull(self, a, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a Weibull distribution. Mirrors :meth:`numpy.random.Generator.weibull`.""" + return self._fill("weibull", (a,), {}, shape, **kwargs) + + def zipf(self, a, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a Zipf distribution. Mirrors :meth:`numpy.random.Generator.zipf`.""" + return self._fill("zipf", (a,), {}, shape, **kwargs) + def default_rng(seed=None) -> Generator: """Construct a chunk-parallel :class:`Generator`, mirroring :func:`numpy.random.default_rng`. diff --git a/tests/test_random.py b/tests/test_random.py index f80e7dea9..467858970 100644 --- a/tests/test_random.py +++ b/tests/test_random.py @@ -94,3 +94,58 @@ def test_scheduling_independent_of_nthreads(monkeypatch): monkeypatch.setattr(blosc2, "nthreads", 8) parallel = blosc2.random.default_rng(7).random((300, 300), chunks=(17, 23))[:] np.testing.assert_array_equal(serial, parallel) + + +_DIST_CASES = { + "beta": (2.0, 5.0), + "binomial": (10, 0.3), + "chisquare": (3.0,), + "exponential": (2.0,), + "f": (5.0, 2.0), + "gamma": (2.0, 3.0), + "geometric": (0.3,), + "gumbel": (0.0, 1.0), + "hypergeometric": (10, 10, 5), + "laplace": (0.0, 1.0), + "logistic": (0.0, 1.0), + "lognormal": (0.0, 1.0), + "logseries": (0.5,), + "negative_binomial": (5, 0.5), + "noncentral_chisquare": (3.0, 2.0), + "noncentral_f": (5.0, 2.0, 1.0), + "pareto": (3.0,), + "poisson": (4.0,), + "power": (2.0,), + "rayleigh": (1.0,), + "standard_cauchy": (), + "standard_exponential": (), + "standard_gamma": (2.0,), + "standard_normal": (), + "standard_t": (5.0,), + "triangular": (0.0, 0.5, 1.0), + "vonmises": (0.0, 4.0), + "wald": (1.0, 1.0), + "weibull": (1.5,), + "zipf": (2.0,), +} + + +@pytest.mark.parametrize(("method", "args"), _DIST_CASES.items(), ids=_DIST_CASES.keys()) +def test_distribution_reproducible_and_finite(method, args): + def draw(): + rng = getattr(blosc2.random.default_rng(0), method) + return rng(*args, shape=(200,), chunks=(37,))[:] + + a, b = draw(), draw() + np.testing.assert_array_equal(a, b) + assert np.all(np.isfinite(a)) + + +def test_poisson_mean_sanity(): + a = blosc2.random.default_rng(0).poisson(4.0, shape=(5000,))[:] + assert abs(a.mean() - 4.0) < 0.2 + + +def test_standard_normal_dtype(): + a = blosc2.random.default_rng(0).standard_normal(shape=(10,), dtype=np.float32) + assert a.dtype == np.float32 From 1d0c05c36beec20e6f22af9aa97a36c9083a8bd8 Mon Sep 17 00:00:00 2001 From: Francesc Alted Date: Sun, 19 Jul 2026 14:49:00 +0200 Subject: [PATCH 3/8] blosc2.random: add choice and the vector-output distributions choice (replace=True, 1-D/scalar `a`) drops straight into the existing per-chunk-independent machinery, same as any scalar distribution; replace=False and multi-dim `a` raise NotImplementedError since they'd need cross-chunk coordination this module doesn't do. dirichlet/multinomial/multivariate_hypergeometric/multivariate_normal each draw a whole length-k vector per element, so output shape is shape + (k,) with the trailing dim pinned to k and never split across chunks. Factored the shared chunk-loop out of _fill into _fill_chunks so _fill_vector can reuse it with a different chunk-shape function. --- RELEASE_NOTES.md | 13 ++++--- src/blosc2/random.py | 92 ++++++++++++++++++++++++++++++++++++++++---- tests/test_random.py | 64 ++++++++++++++++++++++++++++++ 3 files changed, 156 insertions(+), 13 deletions(-) diff --git a/RELEASE_NOTES.md b/RELEASE_NOTES.md index 2ac21ec78..e2645245c 100644 --- a/RELEASE_NOTES.md +++ b/RELEASE_NOTES.md @@ -8,11 +8,14 @@ XXX version-specific blurb XXX - New `blosc2.random` module: seedable, NumPy-quality random `NDArray` constructors, mirroring most of `numpy.random.Generator` (`random`, - `integers`, `normal`, `uniform`, and 30 further distributions such as - `poisson`, `gamma`, `beta`, `standard_normal`, `binomial`, `exponential`...). - Each chunk gets its own independent `SeedSequence`-spawned stream and is - generated concurrently in a thread pool, giving full `PCG64` quality with - genuinely parallel generation (measured ~3.2x faster than + `integers`, `normal`, `uniform`, `choice`, 30 further scalar distributions + such as `poisson`/`gamma`/`beta`/`standard_normal`/`binomial`/`exponential`, + and the vector-valued `dirichlet`/`multinomial`/`multivariate_normal`/ + `multivariate_hypergeometric`, whose trailing per-draw dimension is kept + whole within a chunk). Each chunk gets its own independent + `SeedSequence`-spawned stream and is generated concurrently in a thread + pool, giving full `PCG64` quality with genuinely parallel generation + (measured ~3.2x faster than `asarray(np.random.default_rng(...).random(...))` on a 100M-element array). ## Changes from 4.9.0 to 4.9.1 diff --git a/src/blosc2/random.py b/src/blosc2/random.py index f0f7c4ed3..2eff993da 100644 --- a/src/blosc2/random.py +++ b/src/blosc2/random.py @@ -59,20 +59,14 @@ class Generator: def __init__(self, seed=None): self._root = seed if isinstance(seed, np.random.SeedSequence) else np.random.SeedSequence(seed) - def _fill(self, method: str, args: tuple, kwargs: dict, shape, **b2_kwargs: Any) -> blosc2.NDArray: - if shape is None: - raise TypeError("shape is required") - dtype = getattr(np.random.default_rng(0), method)(*args, size=1, **kwargs).dtype - dst = blosc2.empty(shape, dtype=dtype, **b2_kwargs) - + def _fill_chunks(self, dst, method, args, kwargs, chunk_shape_of) -> blosc2.NDArray: slices = list(_chunk_slices(dst)) call_ss = self._root.spawn(1)[0] chunk_streams = call_ss.spawn(len(slices)) def gen(i): rng = np.random.default_rng(chunk_streams[i]) - chunk_shape = tuple(s.stop - s.start for s in slices[i]) - return i, getattr(rng, method)(*args, size=chunk_shape, **kwargs) + return i, getattr(rng, method)(*args, size=chunk_shape_of(slices[i]), **kwargs) if len(slices) <= 1 or blosc2.nthreads == 1: results = (gen(i) for i in range(len(slices))) @@ -83,6 +77,38 @@ def gen(i): dst[slices[i]] = buf return dst + def _fill(self, method: str, args: tuple, kwargs: dict, shape, **b2_kwargs: Any) -> blosc2.NDArray: + if shape is None: + raise TypeError("shape is required") + dtype = getattr(np.random.default_rng(0), method)(*args, size=1, **kwargs).dtype + dst = blosc2.empty(shape, dtype=dtype, **b2_kwargs) + return self._fill_chunks(dst, method, args, kwargs, lambda sl: tuple(s.stop - s.start for s in sl)) + + def _fill_vector( + self, method: str, args: tuple, kwargs: dict, shape, k: int, **b2_kwargs: Any + ) -> blosc2.NDArray: + # Each draw produces a whole length-k vector, so the trailing dimension must + # never be split across chunks (a chunk can only hold complete draws). + if shape is None: + raise TypeError("shape is required") + shape = (shape,) if isinstance(shape, int) else tuple(shape) + full_shape = shape + (k,) + dtype = getattr(np.random.default_rng(0), method)(*args, size=1, **kwargs).dtype + + chunks = b2_kwargs.get("chunks") + if chunks is not None: + if tuple(chunks)[-1] != k: + raise ValueError(f"chunks[-1] must equal the trailing vector length {k}") + else: + # ponytail: reuse empty()'s own chunk-sizing heuristic instead of + # reimplementing it, then pin the non-splittable trailing dim to its full length. + b2_kwargs["chunks"] = blosc2.empty(full_shape, dtype=dtype).chunks[:-1] + (k,) + + dst = blosc2.empty(full_shape, dtype=dtype, **b2_kwargs) + return self._fill_chunks( + dst, method, args, kwargs, lambda sl: tuple(s.stop - s.start for s in sl[:-1]) + ) + def random(self, shape, dtype=np.float64, **kwargs: Any) -> blosc2.NDArray: """Draw uniform floats in [0, 1). Mirrors :meth:`numpy.random.Generator.random`.""" return self._fill("random", (), {"dtype": dtype}, shape, **kwargs) @@ -238,6 +264,56 @@ def zipf(self, a, *, shape, **kwargs: Any) -> blosc2.NDArray: """Draw samples from a Zipf distribution. Mirrors :meth:`numpy.random.Generator.zipf`.""" return self._fill("zipf", (a,), {}, shape, **kwargs) + def choice(self, a, *, shape, p=None, replace: bool = True, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from ``a`` with replacement. Mirrors :meth:`numpy.random.Generator.choice`. + + Only ``replace=True`` (numpy's default) and 1-D or scalar-int ``a`` are + supported: sampling without replacement, or along an axis of a multi-dimensional + ``a``, needs whole-array coordination across chunks that this module doesn't do. + """ + if not replace: + raise NotImplementedError( + "blosc2.random.Generator.choice only supports replace=True: sampling " + "without replacement needs whole-array coordination across chunks" + ) + if np.ndim(a) > 1: + raise NotImplementedError("blosc2.random.Generator.choice only supports 1-D or scalar `a`") + return self._fill("choice", (a,), {"p": p, "replace": True}, shape, **kwargs) + + def dirichlet(self, alpha, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a Dirichlet distribution. Mirrors :meth:`numpy.random.Generator.dirichlet`. + + Output shape is ``shape + (len(alpha),)``; see :meth:`_fill_vector` — the + trailing dimension holds one full draw and is never split across chunks. + """ + return self._fill_vector("dirichlet", (alpha,), {}, shape, len(alpha), **kwargs) + + def multinomial(self, n, pvals, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a multinomial distribution. Mirrors :meth:`numpy.random.Generator.multinomial`. + + Output shape is ``shape + (len(pvals),)``; see :meth:`dirichlet` for the + trailing-dimension note. + """ + return self._fill_vector("multinomial", (n, pvals), {}, shape, len(pvals), **kwargs) + + def multivariate_hypergeometric(self, colors, nsample, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a multivariate hypergeometric distribution. + + Mirrors :meth:`numpy.random.Generator.multivariate_hypergeometric`. Output shape is + ``shape + (len(colors),)``; see :meth:`dirichlet` for the trailing-dimension note. + """ + return self._fill_vector( + "multivariate_hypergeometric", (colors, nsample), {}, shape, len(colors), **kwargs + ) + + def multivariate_normal(self, mean, cov, *, shape, **kwargs: Any) -> blosc2.NDArray: + """Draw samples from a multivariate normal distribution. + + Mirrors :meth:`numpy.random.Generator.multivariate_normal`. Output shape is + ``shape + (len(mean),)``; see :meth:`dirichlet` for the trailing-dimension note. + """ + return self._fill_vector("multivariate_normal", (mean, cov), {}, shape, len(mean), **kwargs) + def default_rng(seed=None) -> Generator: """Construct a chunk-parallel :class:`Generator`, mirroring :func:`numpy.random.default_rng`. diff --git a/tests/test_random.py b/tests/test_random.py index 467858970..2b9833ece 100644 --- a/tests/test_random.py +++ b/tests/test_random.py @@ -149,3 +149,67 @@ def test_poisson_mean_sanity(): def test_standard_normal_dtype(): a = blosc2.random.default_rng(0).standard_normal(shape=(10,), dtype=np.float32) assert a.dtype == np.float32 + + +def test_choice_int_a_bounds_and_reproducible(): + a = blosc2.random.default_rng(0).choice(5, shape=(1000,), chunks=(137,)) + b = blosc2.random.default_rng(0).choice(5, shape=(1000,), chunks=(137,)) + np.testing.assert_array_equal(a[:], b[:]) + assert a[:].min() >= 0 + assert a[:].max() < 5 + + +def test_choice_array_a(): + a = blosc2.random.default_rng(0).choice(np.array(["x", "y", "z"]), shape=(20,)) + assert set(a[:].tolist()) <= {"x", "y", "z"} + + +def test_choice_replace_false_not_implemented(): + with pytest.raises(NotImplementedError): + blosc2.random.default_rng(0).choice(5, shape=(10,), replace=False) + + +def test_choice_2d_a_not_implemented(): + with pytest.raises(NotImplementedError): + blosc2.random.default_rng(0).choice(np.zeros((3, 3)), shape=(10,)) + + +_VECTOR_DIST_CASES = { + "dirichlet": (([1.0, 2.0, 3.0],), 3), + "multinomial": ((10, [0.2, 0.3, 0.5]), 3), + "multivariate_hypergeometric": (([3, 4, 5], 6), 3), + "multivariate_normal": (([0.0, 0.0], [[1.0, 0.0], [0.0, 1.0]]), 2), +} + + +@pytest.mark.parametrize( + ("method", "args", "k"), + [(m, a, k) for m, (a, k) in _VECTOR_DIST_CASES.items()], + ids=_VECTOR_DIST_CASES.keys(), +) +def test_vector_distribution_shape_reproducible_and_finite(method, args, k): + def draw(): + rng = getattr(blosc2.random.default_rng(0), method) + return rng(*args, shape=(40,)) + + a = draw() + assert a.shape == (40, k) + assert a.chunks[-1] == k + b = draw() + np.testing.assert_array_equal(a[:], b[:]) + assert np.all(np.isfinite(a[:])) + + +def test_dirichlet_sums_to_one(): + a = blosc2.random.default_rng(0).dirichlet([1.0, 2.0, 3.0], shape=(50,), chunks=(7, 3))[:] + np.testing.assert_allclose(a.sum(axis=-1), 1.0) + + +def test_vector_dist_chunks_mismatch_raises(): + with pytest.raises(ValueError): + blosc2.random.default_rng(0).dirichlet([1.0, 2.0, 3.0], shape=(50,), chunks=(7, 2)) + + +def test_vector_dist_scalar_leading_shape(): + a = blosc2.random.default_rng(0).dirichlet([1.0, 2.0, 3.0], shape=()) + assert a.shape == (3,) From e83249f68429e6adffe68b50454315a1d3f6d2fc Mon Sep 17 00:00:00 2001 From: Francesc Alted Date: Sun, 19 Jul 2026 14:59:48 +0200 Subject: [PATCH 4/8] docs: add NumPy-compatibility table and usage example for blosc2.random doc/reference/random.rst now documents the module's consistent deviations from numpy.random.Generator (shape vs size, no out=, gamma's shape_param rename, choice's replace=True restriction, vector-output trailing dim) plus a full coverage table for all 43 Generator methods (39 implemented). examples/ndarray/random-constructor.py is a runnable walkthrough: basic draws, reproducibility semantics, blosc2 kwargs passthrough, a few distributions including a vector-valued one, and a timing comparison against asarray(np.random...). RELEASE_NOTES.md entry updated to match the current 39/43 coverage. --- RELEASE_NOTES.md | 29 +++-- doc/reference/random.rst | 172 ++++++++++++++++++++++++- examples/ndarray/random-constructor.py | 73 +++++++++++ 3 files changed, 264 insertions(+), 10 deletions(-) create mode 100644 examples/ndarray/random-constructor.py diff --git a/RELEASE_NOTES.md b/RELEASE_NOTES.md index e2645245c..cbf80cd2d 100644 --- a/RELEASE_NOTES.md +++ b/RELEASE_NOTES.md @@ -7,16 +7,27 @@ XXX version-specific blurb XXX ### New features - New `blosc2.random` module: seedable, NumPy-quality random `NDArray` - constructors, mirroring most of `numpy.random.Generator` (`random`, - `integers`, `normal`, `uniform`, `choice`, 30 further scalar distributions - such as `poisson`/`gamma`/`beta`/`standard_normal`/`binomial`/`exponential`, - and the vector-valued `dirichlet`/`multinomial`/`multivariate_normal`/ - `multivariate_hypergeometric`, whose trailing per-draw dimension is kept - whole within a chunk). Each chunk gets its own independent - `SeedSequence`-spawned stream and is generated concurrently in a thread - pool, giving full `PCG64` quality with genuinely parallel generation - (measured ~3.2x faster than + constructors. Each chunk gets its own independent `SeedSequence`-spawned + stream and is generated concurrently in a thread pool, giving full `PCG64` + quality with genuinely parallel generation (measured ~3x faster than `asarray(np.random.default_rng(...).random(...))` on a 100M-element array). + Covers 39 of `numpy.random.Generator`'s 43 public methods (full + compatibility table in `doc/reference/random.rst`): + - Core: `random`, `integers`, `normal`, `uniform`, `choice` + (`replace=True` only). + - 30 scalar distributions: `beta`, `binomial`, `chisquare`, + `exponential`, `f`, `gamma`, `geometric`, `gumbel`, `hypergeometric`, + `laplace`, `logistic`, `lognormal`, `logseries`, `negative_binomial`, + `noncentral_chisquare`, `noncentral_f`, `pareto`, `poisson`, `power`, + `rayleigh`, `standard_cauchy`, `standard_exponential`, + `standard_gamma`, `standard_normal`, `standard_t`, `triangular`, + `vonmises`, `wald`, `weibull`, `zipf`. + - 4 vector-valued distributions (output shape `shape + (k,)`, one draw + per trailing vector): `dirichlet`, `multinomial`, + `multivariate_hypergeometric`, `multivariate_normal`. + - Not implemented: `permutation`/`permuted`/`shuffle` (whole-array + shuffling needs cross-chunk coordination this module doesn't do) and + `bytes` (returns raw `bytes`, not an `NDArray`). ## Changes from 4.9.0 to 4.9.1 diff --git a/doc/reference/random.rst b/doc/reference/random.rst index 413c00b1a..996e4fea5 100644 --- a/doc/reference/random.rst +++ b/doc/reference/random.rst @@ -1,6 +1,9 @@ Random ------ -Chunk-parallel, NumPy-quality random :ref:`NDArray ` constructors. +Chunk-parallel, NumPy-quality random :ref:`NDArray ` constructors. Each chunk +gets its own independent :class:`~numpy.random.SeedSequence`-spawned stream and is +generated concurrently in a thread pool, so generation itself parallelizes (not just +compression). .. currentmodule:: blosc2.random @@ -12,3 +15,170 @@ Chunk-parallel, NumPy-quality random :ref:`NDArray ` constructors. .. autofunction:: blosc2.random.default_rng .. autoclass:: blosc2.random.Generator :members: + +NumPy compatibility +~~~~~~~~~~~~~~~~~~~ +:class:`Generator` mirrors :class:`numpy.random.Generator` method-for-method: same +method names, same distribution parameters, same math. A few differences apply +consistently across the module: + +- ``shape`` replaces numpy's ``size`` and is **required** (no implicit scalar draws), + and keyword-only on every method except :meth:`~Generator.random`. +- There is no ``out=`` parameter — each call allocates its own destination + :class:`~blosc2.NDArray`, filled chunk by chunk. +- Extra keyword arguments (``chunks``, ``cparams``, ``urlpath``, ...) are forwarded to + :func:`blosc2.empty`. +- Results are reproducible for a given ``(seed, call order, shape, chunks)``; changing + the chunk layout changes the values, since chunks map one-to-one to independent + ``SeedSequence`` streams. +- :meth:`~Generator.gamma` and :meth:`~Generator.standard_gamma` rename numpy's + distribution-shape parameter (itself called ``shape``) to ``shape_param``, and make it + positional-only, to avoid colliding with the array-shape ``shape`` keyword. +- :meth:`~Generator.choice` only supports ``replace=True`` (numpy's default) and a 1-D + or scalar-int ``a``; sampling without replacement, or along an axis of a + multi-dimensional ``a``, needs whole-array coordination across chunks that this + module doesn't do. +- The vector-valued distributions (:meth:`~Generator.dirichlet`, + :meth:`~Generator.multinomial`, :meth:`~Generator.multivariate_hypergeometric`, + :meth:`~Generator.multivariate_normal`) draw one length-``k`` vector per element: + output shape is ``shape + (k,)``, and that trailing dimension is always kept whole + within a chunk. + +Coverage of :class:`numpy.random.Generator`'s public methods: + +.. list-table:: + :header-rows: 1 + :widths: 25 10 45 + + * - Method + - Status + - Notes + * - ``beta`` + - ✅ + - + * - ``binomial`` + - ✅ + - + * - ``bytes`` + - ❌ + - returns raw ``bytes``, not an ``NDArray`` — outside this module's contract + * - ``chisquare`` + - ✅ + - + * - ``choice`` + - ✅ + - ``replace=True`` and 1-D/scalar ``a`` only + * - ``dirichlet`` + - ✅ + - vector output, see above + * - ``exponential`` + - ✅ + - + * - ``f`` + - ✅ + - + * - ``gamma`` + - ✅ + - distribution-shape parameter renamed ``shape_param`` + * - ``geometric`` + - ✅ + - + * - ``gumbel`` + - ✅ + - + * - ``hypergeometric`` + - ✅ + - + * - ``integers`` + - ✅ + - + * - ``laplace`` + - ✅ + - + * - ``logistic`` + - ✅ + - + * - ``lognormal`` + - ✅ + - + * - ``logseries`` + - ✅ + - + * - ``multinomial`` + - ✅ + - vector output, see above + * - ``multivariate_hypergeometric`` + - ✅ + - vector output, see above + * - ``multivariate_normal`` + - ✅ + - vector output, see above + * - ``negative_binomial`` + - ✅ + - + * - ``noncentral_chisquare`` + - ✅ + - + * - ``noncentral_f`` + - ✅ + - + * - ``normal`` + - ✅ + - + * - ``pareto`` + - ✅ + - + * - ``permutation`` + - ❌ + - whole-array shuffle is inherently sequential, not per-chunk independent + * - ``permuted`` + - ❌ + - same as ``permutation`` + * - ``poisson`` + - ✅ + - + * - ``power`` + - ✅ + - + * - ``random`` + - ✅ + - + * - ``rayleigh`` + - ✅ + - + * - ``shuffle`` + - ❌ + - same as ``permutation`` + * - ``standard_cauchy`` + - ✅ + - + * - ``standard_exponential`` + - ✅ + - + * - ``standard_gamma`` + - ✅ + - distribution-shape parameter renamed ``shape_param`` + * - ``standard_normal`` + - ✅ + - + * - ``standard_t`` + - ✅ + - + * - ``triangular`` + - ✅ + - + * - ``uniform`` + - ✅ + - + * - ``vonmises`` + - ✅ + - + * - ``wald`` + - ✅ + - + * - ``weibull`` + - ✅ + - + * - ``zipf`` + - ✅ + - diff --git a/examples/ndarray/random-constructor.py b/examples/ndarray/random-constructor.py new file mode 100644 index 000000000..8a12f3c62 --- /dev/null +++ b/examples/ndarray/random-constructor.py @@ -0,0 +1,73 @@ +####################################################################### +# Copyright (c) 2019-present, Blosc Development Team +# All rights reserved. +# +# SPDX-License-Identifier: BSD-3-Clause +####################################################################### + +# This example shows how to use `blosc2.random` to build blosc2 arrays with the same +# API shape as `numpy.random.Generator`, but generated chunk-by-chunk in parallel. + +from time import time + +import numpy as np + +import blosc2 + +# `default_rng` mirrors `numpy.random.default_rng`: same seeding, same distributions. +rng = blosc2.random.default_rng(42) + +# Every method returns a blosc2.NDArray, not a numpy array. `shape` replaces numpy's +# `size` and is required (no implicit scalar draws) and keyword-only for every +# distribution except `random`, where it stays positional to mirror numpy's own order. +a = rng.random((4, 4)) +print("*** rng.random((4, 4)) ***") +print(a[:]) + +b = rng.integers(0, 10, shape=(10,), dtype=np.int32) +print("\n*** rng.integers(0, 10, shape=(10,), dtype=np.int32) ***") +print(b[:]) + +# Reproducibility: same seed -> same array. Successive calls on the *same* generator +# draw different arrays (the "generator state advances" semantics numpy users expect). +c1 = blosc2.random.default_rng(7).normal(shape=(1000,)) +c2 = blosc2.random.default_rng(7).normal(shape=(1000,)) +same_seed = blosc2.random.default_rng(7) +d1 = same_seed.normal(shape=(1000,)) +d2 = same_seed.normal(shape=(1000,)) +print("\n*** Reproducibility ***") +print(f"Two fresh generators, same seed, first call: equal = {np.array_equal(c1[:], d1[:])}") +print(f"One generator, two successive calls: equal = {np.array_equal(d1[:], d2[:])}") + +# Extra kwargs are forwarded to blosc2.empty(): chunks, cparams, urlpath, ... +e = rng.uniform(-1, 1, shape=(1_000, 1_000), chunks=(100, 1_000), cparams={"clevel": 5}) +print(f"\n*** rng.uniform(..., chunks=(100, 1000), cparams=clevel5) -> cratio {e.schunk.cratio:.1f}x ***") + +# Most of numpy.random.Generator's other distributions are covered too. +poisson = rng.poisson(4.0, shape=(10,)) +gamma = rng.gamma(2.0, scale=1.5, shape=(10,)) # shape_param renamed, see docs +choice = rng.choice(["red", "green", "blue"], shape=(10,)) +print("\n*** A few more distributions ***") +print("poisson:", poisson[:]) +print("gamma: ", gamma[:]) +print("choice: ", choice[:]) + +# Vector-valued distributions draw a whole length-k vector per element: output shape +# is `shape + (k,)`, and that trailing dimension is always kept whole within a chunk. +dirichlet = rng.dirichlet([1.0, 2.0, 3.0], shape=(5,)) +print("\n*** rng.dirichlet([1, 2, 3], shape=(5,)) -> shape", dirichlet.shape, "***") +print(dirichlet[:]) +print("rows sum to 1:", np.allclose(dirichlet[:].sum(axis=-1), 1.0)) + +# Parallel generation vs. building the array via numpy + asarray(). +N = 100_000_000 +print(f"\n*** Timing: {N:_} float64 uniform draws ***") +t0 = time() +blosc2.asarray(np.random.default_rng(42).random(N), cparams={"clevel": 1}) +t_numpy = time() - t0 +print(f"asarray(np.random...): {t_numpy:.3f} s") + +t0 = time() +blosc2.random.default_rng(42).random((N,), cparams={"clevel": 1}) +t_blosc2 = time() - t0 +print(f"blosc2.random.default_rng...: {t_blosc2:.3f} s ({t_numpy / t_blosc2:.1f}x faster)") From 2dfb86bb739bc70cd83017c24bdb3f3529e8f497 Mon Sep 17 00:00:00 2001 From: Francesc Alted Date: Sun, 19 Jul 2026 17:52:48 +0200 Subject: [PATCH 5/8] blosc2.random: add permutation, permuted, and shuffle Unlike the rest of the module these are not chunk-parallel: whole-array shuffling is inherently sequential, so each call materializes the full array into memory and shuffles it single-threaded with one spawned SeedSequence stream, then writes the result back out. shuffle additionally requires its argument to already be a blosc2.NDArray, since it mutates in place and returns None to match numpy's contract. Brings coverage of numpy.random.Generator's public methods to 42/43 (only bytes, which returns raw bytes rather than an NDArray, remains out of scope). Docs, coverage table, and RELEASE_NOTES updated to match. --- RELEASE_NOTES.md | 12 ++++++--- doc/reference/random.rst | 18 +++++++++----- src/blosc2/random.py | 40 +++++++++++++++++++++++++++++ tests/test_random.py | 54 ++++++++++++++++++++++++++++++++++++++++ 4 files changed, 114 insertions(+), 10 deletions(-) diff --git a/RELEASE_NOTES.md b/RELEASE_NOTES.md index cbf80cd2d..dc98b3e53 100644 --- a/RELEASE_NOTES.md +++ b/RELEASE_NOTES.md @@ -11,7 +11,7 @@ XXX version-specific blurb XXX stream and is generated concurrently in a thread pool, giving full `PCG64` quality with genuinely parallel generation (measured ~3x faster than `asarray(np.random.default_rng(...).random(...))` on a 100M-element array). - Covers 39 of `numpy.random.Generator`'s 43 public methods (full + Covers 42 of `numpy.random.Generator`'s 43 public methods (full compatibility table in `doc/reference/random.rst`): - Core: `random`, `integers`, `normal`, `uniform`, `choice` (`replace=True` only). @@ -25,9 +25,13 @@ XXX version-specific blurb XXX - 4 vector-valued distributions (output shape `shape + (k,)`, one draw per trailing vector): `dirichlet`, `multinomial`, `multivariate_hypergeometric`, `multivariate_normal`. - - Not implemented: `permutation`/`permuted`/`shuffle` (whole-array - shuffling needs cross-chunk coordination this module doesn't do) and - `bytes` (returns raw `bytes`, not an `NDArray`). + - `permutation`, `permuted`, `shuffle`: unlike the rest of the module, + these are not chunk-parallel — whole-array shuffling is inherently + sequential, so they materialize the full array and shuffle it + single-threaded. `shuffle` additionally requires its argument to + already be an `NDArray`, since it mutates in place and returns `None`, + matching numpy. + - Not implemented: `bytes` (returns raw `bytes`, not an `NDArray`). ## Changes from 4.9.0 to 4.9.1 diff --git a/doc/reference/random.rst b/doc/reference/random.rst index 996e4fea5..4fbcc32fc 100644 --- a/doc/reference/random.rst +++ b/doc/reference/random.rst @@ -43,6 +43,12 @@ consistently across the module: :meth:`~Generator.multivariate_normal`) draw one length-``k`` vector per element: output shape is ``shape + (k,)``, and that trailing dimension is always kept whole within a chunk. +- :meth:`~Generator.permutation`, :meth:`~Generator.permuted`, and + :meth:`~Generator.shuffle` are **not** chunk-parallel: whole-array shuffling is + inherently sequential, so these load the full array into memory and shuffle it + single-threaded rather than generating chunks independently. :meth:`~Generator.shuffle` + additionally requires its argument to already be a :class:`~blosc2.NDArray`, since it + mutates it in place and returns ``None``, matching numpy. Coverage of :class:`numpy.random.Generator`'s public methods: @@ -129,11 +135,11 @@ Coverage of :class:`numpy.random.Generator`'s public methods: - ✅ - * - ``permutation`` - - ❌ - - whole-array shuffle is inherently sequential, not per-chunk independent + - ✅ + - single-threaded, full-materialization — see above * - ``permuted`` - - ❌ - - same as ``permutation`` + - ✅ + - single-threaded, full-materialization — see above * - ``poisson`` - ✅ - @@ -147,8 +153,8 @@ Coverage of :class:`numpy.random.Generator`'s public methods: - ✅ - * - ``shuffle`` - - ❌ - - same as ``permutation`` + - ✅ + - single-threaded, full-materialization; requires an ``NDArray`` — see above * - ``standard_cauchy`` - ✅ - diff --git a/src/blosc2/random.py b/src/blosc2/random.py index 2eff993da..d16388b40 100644 --- a/src/blosc2/random.py +++ b/src/blosc2/random.py @@ -32,6 +32,10 @@ def _chunk_slices(a: blosc2.NDArray): ) +def _materialize(x) -> np.ndarray: + return x[:] if isinstance(x, blosc2.NDArray) else np.asarray(x) + + def _bounded_map(fn, n: int, max_workers: int): # ponytail: window bounds peak memory to O(workers * chunk); a plain # executor.map would submit all n chunks upfront and could buffer O(n). @@ -314,6 +318,42 @@ def multivariate_normal(self, mean, cov, *, shape, **kwargs: Any) -> blosc2.NDAr """ return self._fill_vector("multivariate_normal", (mean, cov), {}, shape, len(mean), **kwargs) + def permutation(self, x, *, axis: int = 0, **kwargs: Any) -> blosc2.NDArray: + """Return a shuffled copy of ``x`` (or of ``arange(x)`` if ``x`` is an int). + + Mirrors :meth:`numpy.random.Generator.permutation`. + + Unlike the rest of this module, permuting is inherently sequential: it loads + ``x`` fully into memory and shuffles it single-threaded, rather than generating + chunks independently in parallel. + """ + rng = np.random.default_rng(self._root.spawn(1)[0]) + arr = np.arange(x) if isinstance(x, int | np.integer) else _materialize(x) + return blosc2.asarray(rng.permutation(arr, axis=axis), **kwargs) + + def permuted(self, x, *, axis: int | None = None, **kwargs: Any) -> blosc2.NDArray: + """Return a copy of ``x`` with elements permuted along ``axis`` (flattened if ``None``). + + Mirrors :meth:`numpy.random.Generator.permuted`. See :meth:`permutation` for the + single-threaded, full-materialization caveat. + """ + rng = np.random.default_rng(self._root.spawn(1)[0]) + return blosc2.asarray(rng.permuted(_materialize(x), axis=axis), **kwargs) + + def shuffle(self, x: blosc2.NDArray, *, axis: int = 0) -> None: + """Shuffle ``x`` in place along ``axis``. Mirrors :meth:`numpy.random.Generator.shuffle`. + + ``x`` must be a :class:`~blosc2.NDArray`: its full contents are read into memory, + shuffled single-threaded (see :meth:`permutation`), and written back. Returns + ``None``, matching numpy. + """ + if not isinstance(x, blosc2.NDArray): + raise TypeError("blosc2.random.Generator.shuffle requires a blosc2.NDArray to shuffle in place") + rng = np.random.default_rng(self._root.spawn(1)[0]) + arr = x[:] + rng.shuffle(arr, axis=axis) + x[:] = arr + def default_rng(seed=None) -> Generator: """Construct a chunk-parallel :class:`Generator`, mirroring :func:`numpy.random.default_rng`. diff --git a/tests/test_random.py b/tests/test_random.py index 2b9833ece..03488c688 100644 --- a/tests/test_random.py +++ b/tests/test_random.py @@ -213,3 +213,57 @@ def test_vector_dist_chunks_mismatch_raises(): def test_vector_dist_scalar_leading_shape(): a = blosc2.random.default_rng(0).dirichlet([1.0, 2.0, 3.0], shape=()) assert a.shape == (3,) + + +def test_permutation_int_is_a_permutation_and_reproducible(): + a = blosc2.random.default_rng(0).permutation(10) + b = blosc2.random.default_rng(0).permutation(10) + np.testing.assert_array_equal(a[:], b[:]) + assert sorted(a[:].tolist()) == list(range(10)) + + +def test_permutation_successive_calls_differ(): + rng = blosc2.random.default_rng(0) + a = rng.permutation(1000) + b = rng.permutation(1000) + assert not np.array_equal(a[:], b[:]) + + +def test_permutation_array_axis(): + arr = blosc2.arange(12).reshape((3, 4)) + p = blosc2.random.default_rng(1).permutation(arr, axis=1) + assert p.shape == (3, 4) + for row_before, row_after in zip(arr[:], p[:], strict=True): + assert sorted(row_after.tolist()) == sorted(row_before.tolist()) + + +def test_permutation_kwargs_passthrough(): + p = blosc2.random.default_rng(0).permutation(1000, chunks=(100,), cparams={"clevel": 5}) + assert p.chunks == (100,) + assert p.schunk.cparams.clevel == 5 + + +def test_permuted_is_a_permutation(): + a = blosc2.random.default_rng(2).permuted(np.arange(10)) + assert sorted(a[:].tolist()) == list(range(10)) + + +def test_permuted_flattens_when_axis_none(): + a = blosc2.random.default_rng(2).permuted(blosc2.arange(12).reshape((3, 4))) + assert a.shape == (3, 4) + assert sorted(a[:].flatten().tolist()) == list(range(12)) + + +def test_shuffle_in_place_preserves_multiset(): + a = blosc2.arange(10) + before = a[:].copy() + result = blosc2.random.default_rng(3).shuffle(a) + after = a[:] + assert result is None + assert sorted(before.tolist()) == sorted(after.tolist()) + assert not np.array_equal(before, after) + + +def test_shuffle_requires_ndarray(): + with pytest.raises(TypeError): + blosc2.random.default_rng(3).shuffle(np.arange(5)) From 930b500d8af2f69323d41d8eefd6f491247e7c67 Mon Sep 17 00:00:00 2001 From: Francesc Alted Date: Sun, 19 Jul 2026 18:13:16 +0200 Subject: [PATCH 6/8] Fix CI: fall back to serial generation when OS threads are unavailable The WASM/Pyodide CI jobs (PR #678) reported a real nthreads count but can't actually spawn OS threads, so ThreadPoolExecutor.submit() raised "RuntimeError: can't start new thread" on the very first chunk, crashing test_random_distribution_sanity[2] and test_scheduling_independent_of_nthreads. _bounded_map now catches that RuntimeError and falls back to serial generation instead of propagating it. Since the failure happens on the very first submit (thread creation is refused outright, not exhausted), nothing has been yielded yet, so restarting fully serial is safe and correct. --- src/blosc2/random.py | 31 ++++++++++++++++++++----------- tests/test_random.py | 29 +++++++++++++++++++++++++++++ 2 files changed, 49 insertions(+), 11 deletions(-) diff --git a/src/blosc2/random.py b/src/blosc2/random.py index d16388b40..008fa2b62 100644 --- a/src/blosc2/random.py +++ b/src/blosc2/random.py @@ -40,17 +40,26 @@ def _bounded_map(fn, n: int, max_workers: int): # ponytail: window bounds peak memory to O(workers * chunk); a plain # executor.map would submit all n chunks upfront and could buffer O(n). window = 2 * max_workers - with ThreadPoolExecutor(max_workers=max_workers) as ex: - pending = {ex.submit(fn, i): i for i in range(min(window, n))} - next_i = len(pending) - while pending: - done, _ = wait(pending, return_when=FIRST_COMPLETED) - for fut in done: - del pending[fut] - yield fut.result() - if next_i < n: - pending[ex.submit(fn, next_i)] = next_i - next_i += 1 + try: + with ThreadPoolExecutor(max_workers=max_workers) as ex: + pending = {ex.submit(fn, i): i for i in range(min(window, n))} + next_i = len(pending) + while pending: + done, _ = wait(pending, return_when=FIRST_COMPLETED) + for fut in done: + del pending[fut] + yield fut.result() + if next_i < n: + pending[ex.submit(fn, next_i)] = next_i + next_i += 1 + except RuntimeError: + # ponytail: some sandboxes (e.g. wasm32/Pyodide without pthreads) report a + # real nthreads count but can't actually start OS threads at all — submit() + # then raises "can't start new thread" on the very first call, before + # anything has been yielded. Fall back to serial generation rather than + # crashing; upgrade path is none needed, this is the correct behavior. + for i in range(n): + yield fn(i) class Generator: diff --git a/tests/test_random.py b/tests/test_random.py index 03488c688..263b87170 100644 --- a/tests/test_random.py +++ b/tests/test_random.py @@ -96,6 +96,35 @@ def test_scheduling_independent_of_nthreads(monkeypatch): np.testing.assert_array_equal(serial, parallel) +def test_falls_back_to_serial_when_threads_unavailable(monkeypatch): + # Some sandboxes (e.g. wasm32/Pyodide without pthreads) report a real nthreads + # count but can't actually start OS threads: ThreadPoolExecutor.submit() raises + # RuntimeError. _bounded_map must catch that and fall back to serial generation + # instead of propagating it (regression for the CI failure this triggered). + class _NoThreadsExecutor: + def __init__(self, max_workers=None): + pass + + def __enter__(self): + return self + + def __exit__(self, *exc): + return False + + def submit(self, fn, *args): + raise RuntimeError("can't start new thread") + + monkeypatch.setattr(blosc2.random, "ThreadPoolExecutor", _NoThreadsExecutor) + monkeypatch.setattr(blosc2, "nthreads", 8) + fallback = blosc2.random.default_rng(7).random((300, 300), chunks=(17, 23))[:] + + # nthreads=1 takes the serial fast path directly, never touching + # ThreadPoolExecutor, so this is a trustworthy baseline despite the patch above. + monkeypatch.setattr(blosc2, "nthreads", 1) + serial = blosc2.random.default_rng(7).random((300, 300), chunks=(17, 23))[:] + np.testing.assert_array_equal(fallback, serial) + + _DIST_CASES = { "beta": (2.0, 5.0), "binomial": (10, 0.3), From 95ab1b8af31ad43b12485c343bf81147b4eb89aa Mon Sep 17 00:00:00 2001 From: Francesc Alted Date: Sun, 19 Jul 2026 18:13:43 +0200 Subject: [PATCH 7/8] Fix _fill_vector shape normalization for numpy integer scalars Copilot review on PR #678 (src/blosc2/random.py:98): _fill_vector normalized shape with isinstance(shape, int), but blosc2.empty() and the rest of this module accept np.integer too. shape=np.int64(5) fell into tuple(shape) and raised "'numpy.int64' object is not iterable" since numpy scalars aren't iterable. Widened the check to int | np.integer, matching how permutation() already handles its own int-vs-array-like input, and added the regression test Copilot also flagged as missing. --- src/blosc2/random.py | 2 +- tests/test_random.py | 5 +++++ 2 files changed, 6 insertions(+), 1 deletion(-) diff --git a/src/blosc2/random.py b/src/blosc2/random.py index 008fa2b62..fd36d708a 100644 --- a/src/blosc2/random.py +++ b/src/blosc2/random.py @@ -104,7 +104,7 @@ def _fill_vector( # never be split across chunks (a chunk can only hold complete draws). if shape is None: raise TypeError("shape is required") - shape = (shape,) if isinstance(shape, int) else tuple(shape) + shape = (shape,) if isinstance(shape, int | np.integer) else tuple(shape) full_shape = shape + (k,) dtype = getattr(np.random.default_rng(0), method)(*args, size=1, **kwargs).dtype diff --git a/tests/test_random.py b/tests/test_random.py index 263b87170..b5fdf9c0d 100644 --- a/tests/test_random.py +++ b/tests/test_random.py @@ -244,6 +244,11 @@ def test_vector_dist_scalar_leading_shape(): assert a.shape == (3,) +def test_vector_dist_numpy_integer_shape(): + a = blosc2.random.default_rng(0).dirichlet([1.0, 2.0, 3.0], shape=np.int64(5)) + assert a.shape == (5, 3) + + def test_permutation_int_is_a_permutation_and_reproducible(): a = blosc2.random.default_rng(0).permutation(10) b = blosc2.random.default_rng(0).permutation(10) From 460030f3429fc27070ad37cb0ff5b0605e29b1dd Mon Sep 17 00:00:00 2001 From: Francesc Alted Date: Wed, 22 Jul 2026 17:46:40 +0200 Subject: [PATCH 8/8] docs: rewrite the pandas engine guide around the readable-UDF pitch MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit The guide opened with `add_one(col) -> col + 1`, which measures 0.53x on a 1M x 8 frame — i.e. it led with a case where `engine=blosc2.jit` is twice as slow as a plain `apply`. The actual payoff only appeared in a benchmark block at the bottom. Restructure so the first screenful explains why and when the engine helps: - Replace the running example with Yeo-Johnson (the power transform behind sklearn's PowerTransformer). Its parameter is fitted per column, so applying it across every column is the correct semantics rather than an artifact, and it is multi-statement, which is what the readability argument rests on. - Add a two-panel plot measuring the two conditions separately: speedup vs row count (break-even between 100k and 1M rows) and vs the number of fused operations (2.1x for one, 3.8x for five). - Keep `col + 1` as an explicit counter-example under "when not to use it". - Move the `axis=` discussion and Limitations to the end as reference material. Add an honest comparison against `pd.eval`/numexpr, which a reader arriving from pandas' Enhancing performance page will ask about: numexpr is somewhat faster on in-memory frames (3.14x vs 2.13x), and the argument for the Blosc2 engine is that you write a Python function rather than a quoted string. Blosc2's compressed/out-of-core strength does not apply on this path, since pandas materialises each column as a NumPy array first. Document three gotchas found while measuring: - Plain `apply` passes a Series (`std()` defaults to ddof=1) while the engine passes an ndarray (ddof=0), so identical source computes different numbers. - `np.where` evaluates both arms, so each must be clamped to its own domain; without this the Yeo-Johnson example raises a negative base to a fractional power for ~159k of 1M elements and floods the run with NaN warnings. - `np.sign` is unsupported on traced expressions. Also say "evaluates" rather than "compiles": with the default `jit=None`, plain expressions go to miniexpr's bytecode interpreter, not a JIT-compiled kernel. The win here is fusion, not compilation. bench/bench_pandas_engine.py gains the two sweeps, a numexpr reference series and the plot, so every figure quoted in the guide is reproducible. Co-Authored-By: Claude Opus 4.8 --- bench/bench_pandas_engine.py | 164 +++++++++++++++++++--- doc/guides/pandas_engine.md | 202 +++++++++++++++++++++------ doc/guides/pandas_engine/speedup.png | Bin 0 -> 51870 bytes 3 files changed, 305 insertions(+), 61 deletions(-) create mode 100644 doc/guides/pandas_engine/speedup.png diff --git a/bench/bench_pandas_engine.py b/bench/bench_pandas_engine.py index b73a1b083..898be4dce 100644 --- a/bench/bench_pandas_engine.py +++ b/bench/bench_pandas_engine.py @@ -8,10 +8,22 @@ # Benchmark: DataFrame.apply(f, engine=blosc2.jit) vs plain DataFrame.apply(f) # # engine=blosc2.jit calls the vectorized function once per column (the -# default axis=0), so the win comes from the Blosc2/numexpr compute engine -# (operator fusion, multi-threading) beating plain NumPy on a -# multi-operation elementwise expression over a large 1D array. This script -# measures that on a 1,000,000-row, 8-column frame. +# default axis=0), so the win comes from the Blosc2 compute engine (operator +# fusion, multi-threading) beating plain NumPy, which evaluates the function +# body one operation at a time and allocates a full-size temporary at each +# step. +# +# Two conditions must both hold for the engine to pay off, and this script +# measures each one separately: +# +# * enough rows -- below ~50k the per-call setup dominates (rows sweep) +# * enough ops -- a single operation has nothing to fuse (ops sweep) +# +# numexpr is measured alongside as a reference point: on in-memory frames it +# is somewhat faster than the Blosc2 engine, so the argument for +# engine=blosc2.jit is that you write a readable Python function instead of a +# quoted expression string, not that it wins a raw speed race. See +# doc/guides/pandas_engine.md. # # Note: axis=1 (row-wise) is NOT a good fit for this engine. It still calls # the function once per row in a Python loop either way, and for a handful @@ -22,27 +34,76 @@ # # Each measurement is the minimum of NRUNS repetitions to reduce noise. +from pathlib import Path from time import perf_counter +import numexpr import numpy as np import pandas as pd import blosc2 -NRUNS = 3 +NRUNS = 5 NROWS = 1_000_000 NCOLS = 8 +ROW_SWEEP = (1_000, 10_000, 100_000, 1_000_000, 5_000_000) + +OUT_DIR = Path(__file__).resolve().parent.parent / "doc" / "guides" / "pandas_engine" + +# dataviz reference palette, same values as bench/optim_tips/common.py +COLOR_TIP = "#1baf7a" +INK = "#0b0b0b" +MUTED = "#898781" +GRID = "#e1e0d9" -def make_df(): + +def make_df(nrows=NROWS, ncols=NCOLS): rng = np.random.default_rng(0) return pd.DataFrame( - {f"c{i}": rng.random(NROWS) for i in range(NCOLS)}, + {f"c{i}": rng.normal(size=nrows) for i in range(ncols)}, ) -def transform(col): - return np.sin(col) * np.cos(col) + col**2 - np.sqrt(np.abs(col)) + np.exp(-col) +def yeo_johnson(col, lam=0.5): + """Power transform behind sklearn's PowerTransformer, applied per column. + + np.where evaluates both arms over the whole column, so each is clamped to + its own domain: without the np.maximum calls, the unselected arm raises a + negative base to a fractional power and floods the run with NaNs and + "invalid value encountered in power" warnings. + """ + pos = (np.power(np.maximum(col, 0.0) + 1.0, lam) - 1.0) / lam + neg = -(np.power(np.maximum(-col, 0.0) + 1.0, 2.0 - lam) - 1.0) / (2.0 - lam) + return np.where(col >= 0, pos, neg) + + +# The same transform as a single numexpr expression: legal, but this is what +# the readability argument is about. +YEO_JOHNSON_NX = ( + "where(c >= 0, " + "((maximum(c, 0.0) + 1.0) ** lam - 1.0) / lam, " + "-((maximum(-c, 0.0) + 1.0) ** (2.0 - lam) - 1.0) / (2.0 - lam))" +) + + +def numexpr_apply(df, expr=YEO_JOHNSON_NX, lam=0.5): + return pd.DataFrame( + {col: numexpr.evaluate(expr, local_dict={"c": df[col].values, "lam": lam}) for col in df}, + index=df.index, + ) + + +# Prefixes of one expression, so the number of fused operations is the only +# variable: mixing in different *kinds* of operation would measure the cost of +# the operations rather than the benefit of fusing them. +OPS_SWEEP = ( + ("1", lambda col: np.sin(col)), + ("2", lambda col: np.sin(col) * np.cos(col)), + ("3", lambda col: np.sin(col) * np.cos(col) + col**2), + ("4", lambda col: np.sin(col) * np.cos(col) + col**2 - np.sqrt(np.abs(col))), + ("5", lambda col: np.sin(col) * np.cos(col) + col**2 - np.sqrt(np.abs(col)) + np.exp(-col)), +) def timeit(fn): @@ -55,18 +116,89 @@ def timeit(fn): return best, result +def speedup(df, func): + """Plain apply vs engine=blosc2.jit, returning (speedup, t_plain, t_engine).""" + t_plain, _ = timeit(lambda: df.apply(func)) + t_engine, _ = timeit(lambda: df.apply(func, engine=blosc2.jit)) + return t_plain / t_engine, t_plain, t_engine + + +def save_plot(row_speedups, ops_speedups, out_path): + import matplotlib + + matplotlib.use("Agg") + import matplotlib.pyplot as plt + + fig, (ax_rows, ax_ops) = plt.subplots(1, 2, figsize=(8, 3.2)) + + ax_rows.semilogx(ROW_SWEEP, row_speedups, "o-", color=COLOR_TIP, linewidth=2) + ax_rows.set_xlabel("rows (log scale)", color=INK, fontsize=9) + ax_rows.set_ylabel("speedup vs plain apply", color=INK, fontsize=9) + ax_rows.set_title(f"{NCOLS} columns, Yeo-Johnson", color=MUTED, fontsize=9) + + labels = [name for name, _ in OPS_SWEEP] + ax_ops.plot(range(len(labels)), ops_speedups, "o-", color=COLOR_TIP, linewidth=2) + ax_ops.set_xticks(range(len(labels))) + ax_ops.set_xticklabels(labels) + ax_ops.set_xlabel("operations fused into one pass", color=INK, fontsize=9) + ax_ops.set_title(f"{NROWS:,} rows x {NCOLS} columns", color=MUTED, fontsize=9) + + for ax, values in ((ax_rows, row_speedups), (ax_ops, ops_speedups)): + # Break-even: below this line the engine is a net loss. + ax.axhline(1.0, color=MUTED, linestyle="--", linewidth=1) + ax.set_ylim(0, max(values) * 1.25) + ax.yaxis.set_major_formatter(lambda v, _pos: f"{v:g}x") + ax.yaxis.grid(True, color=GRID, linewidth=0.8) + ax.set_axisbelow(True) + ax.spines["top"].set_visible(False) + ax.spines["right"].set_visible(False) + ax.spines["left"].set_color(GRID) + ax.spines["bottom"].set_color(GRID) + ax.tick_params(labelsize=9, colors=MUTED) + + fig.suptitle( + "df.apply(f, engine=blosc2.jit): when it pays off", + fontsize=11, + color=INK, + ) + fig.tight_layout(rect=[0, 0, 1, 0.90]) + out_path.parent.mkdir(parents=True, exist_ok=True) + fig.savefig(out_path, dpi=150) + plt.close(fig) + + def main(): df = make_df() - t_plain, result_plain = timeit(lambda: df.apply(transform)) - t_engine, result_engine = timeit(lambda: df.apply(transform, engine=blosc2.jit)) + t_plain, result_plain = timeit(lambda: df.apply(yeo_johnson)) + t_engine, result_engine = timeit(lambda: df.apply(yeo_johnson, engine=blosc2.jit)) + t_numexpr, result_numexpr = timeit(lambda: numexpr_apply(df)) pd.testing.assert_frame_equal(result_engine, result_plain) - - print(f"rows={NROWS}, cols={NCOLS}") - print(f"plain df.apply(f): {t_plain:.4f} s") - print(f"df.apply(f, engine=blosc2.jit): {t_engine:.4f} s") - print(f"speedup: {t_plain / t_engine:.1f}x") + pd.testing.assert_frame_equal(result_numexpr, result_plain) + + print(f"rows={NROWS}, cols={NCOLS}, transform=Yeo-Johnson") + print(f"plain df.apply(f): {t_plain:.4f} s") + print(f"df.apply(f, engine=blosc2.jit): {t_engine:.4f} s {t_plain / t_engine:.2f}x") + print(f"numexpr per column: {t_numexpr:.4f} s {t_plain / t_numexpr:.2f}x") + + print("\nrows sweep (speedup vs plain apply):") + row_speedups = [] + for nrows in ROW_SWEEP: + sp, tp, te = speedup(make_df(nrows=nrows), yeo_johnson) + row_speedups.append(sp) + print(f" {nrows:>9,} rows: plain {tp:.4f} s engine {te:.4f} s {sp:.2f}x") + + print("\nops sweep (speedup vs plain apply):") + ops_speedups = [] + for name, func in OPS_SWEEP: + sp, tp, te = speedup(df, func) + ops_speedups.append(sp) + print(f" {name.replace(chr(10), ' '):>20}: plain {tp:.4f} s engine {te:.4f} s {sp:.2f}x") + + out_path = OUT_DIR / "speedup.png" + save_plot(row_speedups, ops_speedups, out_path) + print(f"\nplot saved to {out_path}") if __name__ == "__main__": diff --git a/doc/guides/pandas_engine.md b/doc/guides/pandas_engine.md index fde3c366b..d80dd5fd5 100644 --- a/doc/guides/pandas_engine.md +++ b/doc/guides/pandas_engine.md @@ -1,50 +1,169 @@ # Using Blosc2 as a pandas Engine -pandas' `DataFrame.apply` and `Series.map` accept an `engine=` argument: a -callable exposing a `__pandas_udf__` attribute that pandas dispatches to -instead of running the Python-level per-row/per-element loop. `blosc2.jit` -is such an engine. - -The contract is different from a plain `apply`/`map` callback: the function -passed to `engine=blosc2.jit` must be **vectorized** — it is called *once* -with a full NumPy array (a column, a row, or the whole array, depending on -`axis`), not once per element. This is the same contract `@blosc2.jit` -already has everywhere else in the library; using it as a pandas engine -just changes who supplies the array. +pandas' `DataFrame.apply` and `Series.map` accept an `engine=` argument, and +`blosc2.jit` is one such engine. Instead of running your function once per +element in a Python loop, pandas hands it the **whole column at once**, and +Blosc2 evaluates the entire function body in a single multi-threaded pass +over the data. + +The result is typically **2-3x faster than a plain `apply`**, with the +function itself left exactly as you wrote it. + +## An example + +Yeo-Johnson is the power transform behind scikit-learn's `PowerTransformer`, +used to make skewed features more normally distributed. Its parameter is +fitted per column, so applying it to every column of a feature matrix is +precisely what you want: ```python import numpy as np import pandas as pd + import blosc2 -df = pd.DataFrame( - { - "a": np.arange(1_000_000, dtype=np.float64), - "b": np.arange(1_000_000, dtype=np.float64), - } -) +rng = np.random.default_rng(0) +df = pd.DataFrame({f"c{i}": rng.normal(size=1_000_000) for i in range(8)}) @blosc2.jit -def add_one(col): - return col + 1 +def yeo_johnson(col, lam=0.5): + # np.where evaluates both arms, so clamp each to its own domain + pos = (np.power(np.maximum(col, 0.0) + 1.0, lam) - 1.0) / lam + neg = -(np.power(np.maximum(-col, 0.0) + 1.0, 2.0 - lam) - 1.0) / (2.0 - lam) + return np.where(col >= 0, pos, neg) + + +result = df.apply(yeo_johnson, engine=blosc2.jit) +``` + +`result` is a DataFrame of the same shape and column names as `df`, with the +transform applied to every column — the same thing plain `df.apply(yeo_johnson)` +returns, only computed differently. On the machine below it takes 0.056 s +instead of 0.119 s, a **2.1x** speedup. + +## Why it is faster + +Evaluated by NumPy, that function body is a sequence of separate steps: raise +to a power, subtract, divide, do it again for the negative branch, then select. +Each step walks the full column and allocates a new full-size temporary array +to hold its result. + +Blosc2 does not execute the steps one at a time. It captures the whole +expression first, then makes a single pass over the data, computing every step +on one small piece while that piece is still in cache, and spreading the pieces +across cores. The intermediate arrays are never created. +## When it pays off -result = df.apply(add_one, engine=blosc2.jit) +Two things have to be true, and the plot below measures each one on its own. + +![Speedup vs rows and vs number of fused operations](pandas_engine/speedup.png) + +**Enough rows.** Setting up the compute engine costs a fixed amount per call. +At a few hundred thousand rows that setup is still larger than anything it +saves, and the engine is a net loss — at 100,000 rows it runs at 0.70x. Break-even +falls between 100,000 and 1,000,000 rows. + +**Enough arithmetic.** The more operations there are to fuse into one pass, the +more temporaries are avoided and the bigger the win — from 2.1x for a single +operation up to 3.8x for five. + +Beyond a few million rows the speedup flattens and then eases off (1.9x at 5M +above): the arrays no longer fit in cache and the whole computation becomes +limited by memory bandwidth, which fusion can reduce but not eliminate. + +## When not to use it + +**Trivial expressions.** Arithmetic intensity matters more than the raw number +of operations. A single cheap operation over a large array is limited by memory +bandwidth, not by computation, so there is nothing for the engine to win back: + +```python +result = df.apply(lambda col: col + 1, engine=blosc2.jit) # 0.53x — slower! ``` -`axis=0` (the default) calls the function once per column; `axis=1` calls it -once per row. **Use `axis=0`** (or restructure the computation so it works -column-wise): the win comes from the Blosc2/numexpr compute engine (operator -fusion, multi-threading) processing one large 1D array per call, and that -only happens for columns. `axis=1` still calls the function once per row — -same as plain pandas — and for a handful of columns, the overhead of -wrapping each tiny row array for the compute engine outweighs any benefit, -so `engine=blosc2.jit` with `axis=1` is typically *slower* than plain -`apply(axis=1)`. See the benchmark below. +That runs at roughly **half** the speed of a plain `apply`. Reach for the engine +when the function does real work per element, such as transcendental functions +or several chained operations. + +**Small frames.** See the plot above: below a few hundred thousand rows, use a +plain `apply`. + +## Compared to `pd.eval` and numexpr + +pandas' own [Enhancing performance](https://pandas.pydata.org/docs/user_guide/enhancingperf.html) +guide describes `pd.eval(..., engine="numexpr")`, which fuses expressions using +the same underlying idea. It is worth being straightforward about how the two +compare on the example above: -`Series.map(func, engine=blosc2.jit)` works the same way: `func` is called -once with the Series' full underlying array. +| approach | time | vs plain apply | +| --- | --- | --- | +| plain `df.apply(f)` | 0.1194 s | 1.00x | +| `df.apply(f, engine=blosc2.jit)` | 0.0561 s | 2.13x | +| `numexpr.evaluate(...)` per column | 0.0380 s | 3.14x | + +On an in-memory DataFrame, numexpr is somewhat faster. The reason to prefer +`engine=blosc2.jit` is not raw speed but that **you write a Python function +rather than a quoted string**. The numexpr equivalent of `yeo_johnson` has to +become a single expression: + +```python +"where(c >= 0, " +"((maximum(c, 0.0) + 1.0) ** lam - 1.0) / lam, " +"-((maximum(-c, 0.0) + 1.0) ** (2.0 - lam) - 1.0) / (2.0 - lam))" +``` + +The Blosc2 version keeps its intermediate variables and its name, can call +helper functions, can be unit-tested and reused elsewhere under `@blosc2.jit`, +and is checked by your editor and linters. That is the trade being offered. + +Note also that Blosc2's characteristic strength — computing directly over +compressed, potentially larger-than-memory arrays — does not come into play on +this path, because pandas materialises each column as a plain NumPy array +before the engine ever sees it. + +## Gotchas + +**`std()` silently changes meaning.** A plain `apply` passes your function a +pandas Series, whose `.std()` defaults to `ddof=1`. The engine passes a NumPy +array, whose `.std()` defaults to `ddof=0`. The same source code therefore +computes slightly different numbers depending on the engine. Pass `ddof` +explicitly if it matters: + +```python +z = (col - col.mean()) / col.std(ddof=0) +``` + +**No Python `if` on values.** The function is traced rather than executed +statement by statement, so branching on array contents raises +`ValueError: The truth value of an array ... is ambiguous`. Use `np.where`, +nesting it where you would have used `elif`. (This is the same restriction +numexpr has.) Branching on a scalar *parameter* is fine. + +**`np.where` evaluates both arms.** Unlike a real `if`, both branches are +computed over the whole column and only then selected between, so each one runs +on values it was never meant to see. In `yeo_johnson` above, `np.power(col + 1.0, +0.5)` would hit a negative base wherever `col < -1` — around 159,000 elements in +a million-row standard normal — producing NaNs and a +`RuntimeWarning: invalid value encountered in power`. The final answer is still +correct, because `np.where` discards them, but the warning is noise and the work +is wasted. That is why each arm is clamped to its own domain with `np.maximum`. +The same caveat applies to numexpr's `where()`. + +**`np.sign` is not supported** on the traced expressions and raises a +`TypeError`. Express it with `np.where` instead. + +## Row-wise (`axis=1`) + +`axis=0`, the default, calls the function once per column and is where the +benefit lies. `axis=1` calls it once per row — the same Python-level loop plain +pandas would run, except each call now also pays the cost of wrapping a tiny +array for the compute engine. For a handful of columns that overhead outweighs +any gain, making `engine=blosc2.jit` with `axis=1` typically *slower* than a +plain `apply(axis=1)`. + +Use `axis=0`, or restructure the computation so it works column-wise. ## Limitations @@ -61,23 +180,16 @@ once with the Series' full underlying array. the pandas-side API surface, not of the Blosc2 engine. The two entry points that do reach the engine are `DataFrame.apply` and `Series.map`. -## Benchmark +`Series.map(func, engine=blosc2.jit)` works the same way as `DataFrame.apply`: +`func` is called once with the Series' full underlying array. -`bench/bench_pandas_engine.py` compares `df.apply(f, engine=blosc2.jit)` -against plain `df.apply(f)` (`axis=0`, the default) on a 1,000,000-row, -8-column frame, for a multi-operation elementwise expression -(`sin(x)*cos(x) + x**2 - sqrt(|x|) + exp(-x)`). Measured on the development -machine (Apple M4, conda env with pandas 3.0.3): +## Reproducing these numbers -``` -rows=1000000, cols=8 -plain df.apply(f): 0.1114 s -df.apply(f, engine=blosc2.jit): 0.0260 s -speedup: 4.3x -``` - -Run the script for the numbers on your machine: +All figures on this page come from `bench/bench_pandas_engine.py`, measured on +an Apple M4 with pandas 3.0.3 and 8 threads. Run it yourself with: ``` python bench/bench_pandas_engine.py ``` + +It prints the comparison table, both sweeps, and regenerates the plot above. diff --git a/doc/guides/pandas_engine/speedup.png b/doc/guides/pandas_engine/speedup.png new file mode 100644 index 0000000000000000000000000000000000000000..b26d37a3e8eda63a649fef80fe73d455559e944f GIT binary patch literal 51870 zcmeFZXH-*L)HWJK;0OXL0)n)of^?LsbVNitN)KHi5F#KDx(Of$6;V2&3(|Yz@I{&2nD+N z`9poY<)m-llD;b`GX%rJ$D=! zYk*nrzLt4#_Qo8Ob&A#g!5?j11}!>HhC1@R?SwUlfxJ$^?U>vBBGX9=Tabn;#qoPM~M9VYqzI|RMC2V(i}Q2HH+>%SvAUOLJDj#}lvr z)6e4gUru!Y9`zJm>^bxAQG~(u2!?-;Kxq$v3yP05jHpW>k;q+qNiJ8r=bC@nm9siyX696$10 zhkcve$qeMjw%BvYkqsnzcjR#+1IRp2gN8R*>6@FgSB9?% z(b3UCoMQ74Ffo4ddX_EpUx;8$FCxy8^@Q}=D(Ub zRI)YKk+vSF&^poVWn1RhjUEo9Iks@0SD$P3{d9IFmf5MTt*zfcrt3z$+J1WHpRZSU zx5H05QsTLainhHF=uAZ~Po)S>4vyRCCZ8UYr#KtK7 zR`1yvq7#kWebSX79z9SHSwW+c8)7S_GUVKEsHn`kYR|A3a|Qvng_lM$fGmWhq@+%^ zl>^py_)~yH{1ywV(N; zJz&tA`r=C=DSXhxH#(Yxo)6nr^Uo_K<-CZB5|)=Yf601PDd5pIPbTCwq9p zKmVLP#G8bBvi+mrlJCFz9pY*F_;J#wV=(gz7~{m}7!F}oRXbh<&l}etyiJfnALM_g zuFnN6D6}Xbj}ETQLL+Ux$cL-BvwVL^I--B@4Z&E2U|EV@^Y83w z130T9VvK+%kwQ^;@9&p5IfYHsNGWB)1=}hJfq5j=>(zE-{u<1YC4v3Qt)7ilVAW;e z@{Dr~jEw`W3)Z>252E|y#;LgTEC>@kjG7ZPpExvR*g@TnVpqRaq0Anc zEJHlToCI{BNBFu#h8Rf{eeu<+<+b>$`-< zL)k>&rZiml@#NXCMy|xJ$X2%Ffc5z0a2B;->4EIK{zF8*;alT$9RhaiTFC-eHs(5T z=g|R_bgtMY-C56DzmE?0R(fpMci{R-T#gPMfmLF0TWdWbf0n+i2pQEBFr5npb;EGd zdj%K@lH9St+0hLjNS=SmHctZeI7B+D^_2PYu6=YWgl;KpCIps?N^3WlR&!U25($cDc~3SmeqFn+LgCkj#V6m?n|~! zp2FPHb_v?A*!tPS5486>*COFsGxLu%ZRt}f4elcX{-%XmrS;bkQW28%D9iM|P?Q#jQg*3j?x=v?4A13rQ_dl~`$-ZMehQNqd*i z0ve(W(V`q~%DeC0?GkX?)cuaXIv=1YR>W&ot0;;XX)hd~q4Z!EZ z;RD&o__d%Ev=5P%Pc?%zLx*({-Kk_1t!?d-c{HNq(6Jk)uvj=!(maV6t#RyrT(wd0 z1*>s<;Ih|(#DdkiHK-eMhZI81aDPFWSxM#~$(Mk}^7?|)j2$(wW@kxZ<)wpFx!fR4 zvEOGv15ZD{zSu@upB&^@hI?#9vvBXWY|nS8v;a}REiy84WjADRW!g?gTzuduE_ZyY zJpi@^1T{Z7S=pz}1cb)X&Zm|AV<5tDx6jLssJcjFj*07&5@aGCbC|}j5-8aMSb1gF zrUAm)UR{(MI+qKC!!{@s3JEhdHr`?qe4aqq-u97pZ1K(+CILA`V3!ZfME$UE`k1?xOG}^BZ_c)? zzI^p6YtS)lQ=wy9DU7P`qbqaumXoG85J`^%5%GZ4`nCx8FWipuvo{UXo!*DT(xmMN zPWQi@Fkn|b@|vEWKFsA*3)>%>YjBSr4)~s{vAsic8bKqdyN*^&>;gkltb?$;Mfaov zG7PFw(H^wzJ%glfwy-&-`E>qy5JSz%&d$Di;eyk9e?*nd=g*(hXFX?|s}5);_!v80hj_)<1Rw_Lc#!iPFTj=QIk*;+V%YCEo#qqILTI`ud|@lt4Aa z2faYQufh{f2sPb`R7FY2p_LFxLeR$~VS}M&V8d7zxa;@1m{w)m;&q|I3+jn_}Mn-zK=m0maJDL`*G!VL!-w2d+ounlUkl&cE)pu;8XU}FIAC82h%11)GUTE!_XwW8}Qi0qVa^^ba(6uW3dVpH$SB^F zig5i5H3Ws&pBv}-!7MOvW6Jp$dGx1m)4{~5L_z}5HD&wdQHiZs&Na8~Kd?W~!_2Vo zTx-v+%axn?7;oaAlzrWjiM#n`ui$S?w3NGa9R-yU!}TZny8;3Vp5{m+F&LPmf439d zEh8ymw@#RJsH%9;n3E~%rr&1OHvv59JM$7V=$wj`OGx-Sao~mG({Z{na88-Q3K}-}2JKcuE_~P<6{gY_X!^BjVuatH4=n=9?uCS#l{Dyab?rg;W72G(EGfG*wAN` z%`UR?X+y>-NCbl^cC$V{B>3v=*ER}?M@L}CN$#-B%7f6kIQj5LDs+4$v%2E=)J~L zbb+!nD)S8mQLS*c>H9TM&*5o(#$&^#2q)dv)wh*z5#<#M1?oDZ9epN3X3U2uPc@F} zNVX#S5UKMUq7sES5CQND%WFVsgjJ{z$CLP~rYsY(6i#`h3nE zTrRp{6TEEk9#^Mekd;4Q(tgdkcNRm=fw>nM^LowD$0-4B^ zU!nYVU-`lLXADJ%T>J~BcA?kEmxS^}-q+)XVzkb&RlO*)%%85?)H1-Np;^4DSn@xR zz$$Rlz>zg&zaf`4K)<_|B;Nhmxh(c zImUj(z2CFxY$;I|ft2I{qZUv- zmH|iqnGX6Oe>^YA%`?O|!)p1^cr-7GF;GFmSnfBRx5HQs`&1Rty%wO5V3C+h&mbrV z>S;SITHRdFtymW`O%J9k`>z;AZ#o*hV^gzkP*O8x&yHpByt7_WZhd2#m}V4P5p4|Da6Qg}K6Bc+?7FMW6gZ*C=u!NjUq)A1r2axI@54g??D1Ns zoAV!0Tj|np|DJu1YC?ajNFi+h=_xz`m9?HnknKZudYQgtY)p3f-=}9pNp`(Z*ccd8 z9vrxoSAVuvRZYzZg_>_+l!c4KACYRxKet0>rs8Bb2gg9RCEk~H*@gDE7aTo0jN8WO+64}*ae{UG$JFiFxoVyERIGrnyv|nv zso+Y5=Z$&W-Nz&ujy%{{e)OlZM_|W?{CA%9@QJDyCvWnbDpHJ@8o5O$fED3y-6Hh{9WJTG zK_H`z1I#k(EHZhW!TuX18nm*z={QM_9LVxhtimdv(Ypq=;5E43Fj~BP-oeKDg#~Zt z_=rRoJ7LBRHP~dFz9mzoRNRzVPH@-MQFULg(zXOD(&f)@y0%3T;Wotw785})uA{*I z!nET_q3rp13Q(21v<#nSx<@*;Adtl?01!+`6u5=&Er-~uk|xX)nTo*`%cmPV(yM%Y zw2VE}$?{q1-7c9fwc(}Fd*8V-iwc^wiB2u8^om`_KrVkhF)>k*i3`PvlJEE=n;`SB z;i^2x4K=m7qZ+1|O#UYsil<0l%DV5mjlY5O{0aHwIQOuMX)UgviP-a?!3oZbD$4VN zLSIG*z->&F#Lo?%AEg3wpO+y3={iY4Gs2mt0j%v8v&55QVA1|K!d6~25P7kDh&)_Q z8*Wr^A1hO%gyQ;&d(&T3yvD<{=woHNh)Egi;A!G9T;UFZ-VNdCK?g0ur0L5J>o$bS z`|Hx@Q}}F|@|zOS-WclmVp)J?Om%(<(;o&H8hYDf(0>vDD~>jQ^~EcZ3pN8zIL}{$ zjrrsRx*5so5q&)tPQi*j9<%w8g!)m7#y=m&y;+4_Fr|{hPsrg#Tdx@P zB?D?3zdQEgBT*+(66%(}iM{2N^P&=i<*i@T964#~4V}8DNmW))8YFgIU81@+6pH4M z?W z1@dK0D;h$uw*2VV__hPmyw>ADtzha?o_oMMz(o@&KTDhk5&oJPX!j%u@K%sfz5(=Rbt z@9l!HKo>tH<>W^WzR4RLf%UquI^iSpxpc6rw)7IbpTmB=Q!*5)kgxDcO%+(g!(bF5 zenrT>VY87dYhSlPRo?95rJ*q=BR5RDB4r)}wWky`eUy}lfi4gg9zW(ikL4Ve7cCj7 zK+6zmdD|?A@HJxo3$w)ijcJz`E$!P1))%f^@nFj)b(P{KN9YJMDkxeB=@Stt-+Jld zMQ2Mj3^|eUFJdUFCt|s&TI=F)=fn({;u4QM_Hb&uk^8sX2L!IC**Inb*Ii9ZqCBRku74CPI?R^dq#8zTIU)~>&W)~zX)dVng2oB*UNaMCUl9E(uc0{_dpeMc zp0&+DX}r?#kXo{>)n$S(@uYPnjcevGmJtzEH47)Sg%bq~LjZFlvU^?LqDw>9EQSvB zU5gP|+fmZLs@&TYFj`h~<6`p4UH92915k|32QzC)_h!;Yy#U{uRym^_G2$_MoI?5| zKKuS{>UOshJzQ-O+uO!?O|Tu#뗌(rlB&{3VrfHuP#y&Qc*AQ%-v!Q*`LF%6LciVEL=}OY+v#C+KJ2 z{KaeSvGUG%ZUW(bE@3UQy?VlJ+DBL_+tGnS!A3n%%bSdHVg+~81**78?(fl4SjA_9~}SnI9z#;uw+#K^xc$oNNHU0uYwj5?VKnGsVH zl_|Hh*1=Bu|E_A8m$h-H2#of;N3P)GNtYE6Z1Syz)T(+1daPF23p#OD5%88(ls#qI ze5`y)E>BW^R*za*W5UYo3CUJ>j4RHyngHh5oWhF}$~C@;d8FG?x z@_1rq*jZgYpBXoF5)l<6iP(?hxcmSvS3M8QwD);q)Jdm8#w?QW8uRO!@VOdLggA-? zRQl(A$V$c!fHhb~1L%Aay4t}cJU>5QyvQetkyWgv*Gg{*#WFPLi~f7loZJnuJ58(V!J*end4m) z5Tb7;X#+ctZ@kp?-dL`T4gRI!g~14{N~}t_pu*8}STJY)Ov*pTe`ClI#%=2}O|a-q z`6^fuDPF(+lDO&E!a1PQiMqB00+mpR7!#xY1g!sfyEyo)4XiKaCDuvYYb;pJL#Di^ zO9BiaNjCfjurz3|N7c=O0Jp}6b=BXBE83uukUUX}112H@9oV?oS}krWk1DUEM0NL0xajCPFI?*J}(1V?y#aC(>JxtTNG_Ybenb3gU-+FrmOhK4pBe^eM znp$N|v$K?ZGRthha9b-|>S7f86!aO-bHEPz14b;)sjUkCv1*I1uK%~B+rMZ2f63zg zp9}uKXF~sv4gddQ!#sVdYz?bp>Xt?@muE)H!ck7h7QkaW5q3V%Ln<$G@PDbMs|*I8 z(C0x01_qMospNsg>rBgENwB4r6_LnD)*@Joz<5=d$)BB_(sk*3_b+j9xNdEbg@tM( zhn_DnE#Ykru2YT%zBcBj@&&x2C&&TI@OkEuKGrUT8rG zl!Z_IRR&O~ADRIStRA+i$RYnwzy0q3!esLPrb>T8$IIB3mY^7FdSc?=v(4Tq9DVWq z!&zcW+m4Ry@J(qN)A?o z=@Nf{rL3Hs@SWB~;--eFa^QyDlKMwfoO|Hbl(<#39jbA{d@mahGsta-h|mn=UVQSc z=Zt6Tj~{l2>x0qt&NSXcyW3h7#yK4sk>f2etATo9M~eV%Ut3et+VcXkhL<>&+-9*x z(fhG%mR2=(A`T7?ej}|h7MD8Dfbt5lHMZCL;sg|Ljsp)q*!c@6%dQ)s5AP+JB^Bw8 z%#SzQhsCJ+>#q#c9#%VcTb9zkZsYoA6eMm)sMv;*wj_qHg9j2nF}WU>Gz=e<DRh--WO~L+tj=n;jLF|fe zk$Yj5?~i%44*gCU-@dWyuZHy6V!5Rs4x+pYt+tUXZJ5Q?AJc;A5r2blAItpnUaAOk zNR2&O1?IPTcTpR&8zo>96rW_Garn?UN7-36Tt4_kIj&`)XVyKVMXzI@zeq26^k@)% zq0&7`H;$LE^HvT7s$TjahAGnXkShat2IRgsbgn(X485>ZHBqs8Q!Y1HL?f!)f7QF* zAHkD=qPFRi$pm38r#M z+4gIj30w>!ckWSbTb{~mFzs*O?k!s4Q*JepwtJmeJr_FA^-E;M$Fa$;tbaopK`CtX znpF>{#%*R9S!O~)Sb0_DH8>2ToOo3U<3(BWog>Xc$gY`f`w!$A>17VO@=}roH?`YY72dU1txkhl0M9JYUh!~j4r4d zfZG(Z-`AM|KvHx6)EE(=cD+mtQ`X`gbCxsQR=IYj+b906OKSQR8SkyFOIJoUXBxfS z8J(~;j@11vlr|H2un88k^z*0KaAtIW-rw3qP)zNUZC$3eneA015pypM+9hd(f3NOR zY!t?S`Q(UOCad?`pI^E+uB+MK%Ti1pE7(rWTfDG$hxB&i3(M;8Vxvn*wd|AMSF%!& z{7&ndd-J~Jn3|fJr*rN7uG)C>7w^@m-y$NlnceQ=h4!ZEt|!5NwNDEoNk@et_gUJ) zKBpgP`9nPxw)n$ypUbY&+cZ5XQ$7BSDUSE8;~qvtM@MJL{i9!UcTMx?$Mi$F(n(}9 z4EN~ec^(mUbhCRqX`_Yl3(>H3%htz75MzRv7Z(rCc#jeRkg6oJn=;0&X1^eI1Mf_$w_$08nqi`rC*}7Zmu*eq=MXdP25z< z`=$;K&)wg`j2P;^eb*%#@SO0iLE}ge^B~GI^rNH*CA+Ys#DGSj8!z*>^4CC>itc=_?vlI*=#$4;mrtq!ZDfdk_&00(y_l=nE ztFn!5An}d$)xr7{_iV{oapyjM%*w)5srzvG`f|GwP(~2G_aZuo5gJ-wG>3XsOh%5LBLfhiV^18$w}UGq3ON9<2G%d&Qv=z4K-q9 zR^E3}G7Ya8Y!>^wk34jsqYoAbSBS_ z+~lT4q34uAV-dq@zqaJ)`S>q*Bp7TPFuhS0I?u*-8)4VF!?AY+W*;ZEVS3-CwX|el zN9!$Mdv^kqKQ0fZ6?Qgmz1sRU*AbWwa>+&07eYPcqOs2&e=CZUSe`OAxz~c^W$t{fe08ilbz0?+isC_d1#BA_h$>2JAtTzdAeI3@D zfnWrs#6uLsstn8%yfGsVth_0ui2DY502kfMLD;9vmEOEt3Eh#Xv}rw*mfljar+=Cz z{Uvi-I6QooK17Wt|JON#Zb_|cYJUnBj2H-!3SL|givN5$;^_v199u&-J9t`kD`_-2 zmW695^+&cHeykv^(R1$W)A(*b8(_mUHO^7>_O|L{3B6}1H-^KUO4==J!vRSrWQMl&FbaUO`AG8r~DeoLRno=|JFdLO*YOuwBD z`C+pbEWf|Rz>zJs^-oz`5icXd;N<6j?wK(JNS0+PjSzf%gI3a1W~iED_PRFK%os?kXM z_DyaaZLbTRE6xew$ESbYwM%n552sA%9qz_PL$DD7^lQhjJ?6(d(E~}*=d7>Ovy}cM zve=Da@o9wY>Fungg|@qQMcmZoB+NvNkfo;kZ`61%-|GDkmoq%gO=w8z?+;=lO1ZNt z{96Bsk6&dU3MFsBUI2?{mRBaX8CK3eFRpZqxJKZr(zqXXUT`_{n742=FCX&F2HlN75PrLp zvXujVY*49!&!-^ZB7LJSi>D_6KNOvpTn}BhUakC-8s|);+Q=nj%0n|Xgd6Svop2~0 zu5jnCeVZz9`H^?;*Ty;Wc;RD14NGfRf``|aSK2cJAT&Pwxo<>b!KJ=BUbKm2z+W*XC!$v@a9e7_`C1Zb;V$+08}D)yZcTpBE|Gl~ z2_F+h#1R9l7G6=N0A42rfJQW>Zh0TwNFMWF-WBGUn(7Ej%Xm;scb=9nbB@1pEiDWC zu9nR~j8JCT7svNqLj}VBc}$rtM>!xh<+p20Oj)*YdT_*f*MQOVm^IcZQfvld5 z-gDUV;+0W8T*^*$0;QfPaF*icFlf?`3oAR?q-11D|6-YS6i-rEw#7Gs)Z{9o#M~_5 z?U!L^?qYi-il!e>f}7$jdski8b`llo;rutCnsVPg14`?aLx0!uXtMLVs(4Homh~{Q z8YSp4RbU%wi+?wpNj0v#s!!uAZgJtcBY}`a^%kL<0oSwZ#g!fVs5OoA{OdfZv!?sT z&4&lSt1fut{9P0dyB!!N+2RcKIovGma_bQO>k&$QrXxpXbF=2z+OJJ7fby`B`a;6= z&6b#FEM#g;lK0a;-fKuOyl_oYRU_E#uzGUxDv3N+@M^ckwm$iL-)}c(Xy7}>oANLN z`q!DHO8M({CBx8X1x^X0^H<7^b2`iNzKl7g?;Fq?gRgrX*82>Rj}BszSr8{G`M>hY z*fl)A#)0>rtj2kX!1PPT0yZYv?#C{BEb0js#5J%mUI~zNKey;{!NjB)ALB&0UEJIT zTbm1M36S^cv>YSeZx#+Nf^U@m;zJDI;Z?b5?be*?cD`^VvZG_u+L~O#NQu}L`Zof- z9eK07kp5PJ)wRhRV+mI+RJab$%;gHkHOT5P*ugIp(Zz{;Fp8DB5Gt(h^Tx8F_G5_dv_9HRz`F3{FlKV;mm9f z3jPCQy)U{}{}@wM=!4W;O3l>W!-HzJC9C0KXvE>e@}0gKr}!|Z7g&vg>#Cj!KwlTiBGH0il^4rJp9q`pf8ismH7Vsqcy7I zNud1q+U$+m#+IYNiDTp96rUKXb<%7bM*R~~Lr&I^E7`5`v7JVfBGf5@wjC0-VN}!j zBXF)A$=1h`*;TU;!25Iz8pk_Pf3{pcp-90(0cIcdh%#3ZDKWT=jTB;it+u-+Io*C3 zmfTzDVo?;yW)E$1qVv^-pLpDoRHkH~i-T|Zxo6w*NgdS$ln~KdJW9929SJc~IWq^x z;U{81M_2++{CDx)Kz-#0&7F$t)6?DU8l)zBn_n&nsk^*5WnPuptAc6p^M!J%@bL9Y z_F+$qAb%J%Z)YpfWe9dwRN{=aCTj^PJEcqRApkYcf^t)QUXq%o5e zLzH#-%!i&zeu+lO$%Ko@GpFc4-j z=UuZ<68dq8@{Ry$cj%RLGzbgh8485K=DZyffX3KB1JHNBFx0Uv2I-OajiN3^YYL)L z5Je>ZDzcXEkHTlgt;^Ic%ZtAlSTEVQdqeZ|OE^8}#B;8dOq&DkIVXB7C69N~InH?* z&J4LhQlD$64Nf6v!X&ly8j^f3M)DZseVY$9jM4xL8gZWy*FWyJQgLy$;uVyrS}oGx z6`jZ@f#WUQon2iam&;|CQy3LamlV?l{j{%A-|B}y6|BXdp^Z=C-7+l^d&RU!qh5Rv zwxZhF+OncfNoQo`_&H)&v&`v1EZdr0-9Z7HD4ao4ic#Hni5`+wU1??m7cUv73)NfW z8OstAJJfNh!uQnQ7}4{sVOQg*U>L)&sA`1h*pFn<`Iwr%zLUjvuW3I@S~Rce>i6us zo4dfa^y6Y+$Bjsos_4TOe9*H7$WLLIAmaPkM54-UdlzE_GqS|fY5C)qJ;)mi-`F{i zS{>M%b>`mCrO%ltp3!7G6jF@e4EoHd8&+w>YRNk^cD%+6?7gfG-N_rHyi`CPJUgFyLzz-}?lRspE^hrXRK?!jcWBCY zya>h)6qoowow03Csl>#bRADm#<}^jTbn?Mcf>ES6&#wbLl)>73V@f}_4Xu@oQ9JwHUD%NLWauVQP1a8=Rw4G)3 z#EYYO|Ea$NU)Rx&$1z0 za@xd=^THL6bq$Bs3oh%$;I&pQg8sVCVDQI}k4N<6va3#TJ@L}&1(8w+HfAhhy!{aa zx89?&z^|rpYrn9dyz17XYN_Pt+IF!HHCxz8|SYzBa|hP^yQYjQkiv}!g8n1__>uJp#`upQ04P_Y2e`8=JwP=yEvTI1_xbhhB}vEdrf-!` z3f&{<&zu>ed+|PTQY+t0Sl5am2b#Fgd%rz1y=9=v%{Otf0=0}A=23srXtblpsoFYe zD4au~ZO)|127RYOKCo5s2OLc~@7dIt9aiCk{1%Y(Mou2zH^@a^vw9naS>&@+6PAm?tJu^la-wq1#Y(Q>vd7gf4!$qp-@XMqarltZ{Q z(^&TT%hycJENS96*3(Tx!D-Rnk3tmVSRCKZ%gYmrd86!JY_J|=(QEOt$u1Jv4VxOUc-`_?-h^O@P%TMc3?iknU{bBb>YR`n_Bn;`M6h7Fwo z-L^GCrM-fKdd)A~BzM;xHsxmN>94*iDR~9dUcF}viAZMk{TIX=IoPT3+v&@P6I3kN z(h*T-ime2Tu&c^feVZ*6ydcy?B_;zPdFCHGpgr8swua*U;BAijB z9GI;PB_Q=Y5p;uj+TH0?lAw9M+~1Th)Sf7e>%#n;AMu-(zY6bb=D(hrdfX~?9DT>D zbUP}!D=aWW&Mef@XJ^s5_zuLE7<=;bm0&7>X09;UbRveho7v{!&lZ{Clnixm)=E5b z^4fbKJ=Lmx&)Vp~(3ezJ6#pUrEW@jBLpCFsatX@zk77$gLqY&ZkTmSnku{lNi$DrK z_h^R#3bu>>#pX;9<#Yi*+%tdMqQtO{{TBh6#y&ojw>J&xg88lS(@Ok?p0~`lwVCX{ z160n#n768Kj30U5?pB%g6)uxEmHI2;zuzDjBQhVFyvLK16@~@sZZS2q&|8Pr_g1dM zRy!V>dXq*Uv;K@;OvCNsb{vh|% zd=GZCI7UbFY;w^ww{P}A+{*akUdBl@f+=frX6r+loa-$1o!`;HBe?+MsH+r2{!K&X z@k(*23)Ry`iSEXk9I?Dk7@Z*EpP8#cz%(QmcWvP+D-HK(YJbp7133scwPB;ax4h;U z>S&!2ZT1&G|Ll~?#ga&I{tC_^%lflbk{F39s8`|8(D^+zFPpAgdC0twd4Qy`aup0v zj4Kb@Qrp&FXpQJ4L7F4-+;{r7GEa4^~P^S9e8_kH!&-Yijg)e-Bh zcRy58_W@KVV|CEtj5W{^7E#Lx33yxis|tXHjc%FH7ifmGE^`fS(}T=K?Zv@Nx+AnT zcZDbCezI!>pjnNqaib#0B1-!=*h2k%t9l%l4Dw!j2h~7!!38pV%nP04Gp5S+S_BBt7Vt>g?v79K{WiZHeJI2tSDkjw0Z?%Q1tQ2)(#RlQ z)qke&AJ~V`Dygwzqd!fUQ_?u^+G+omsbp1E9!99C8Ta~;3-C|CAD)TJ{%45)! z{#D)nd-VUL0lB?khgYFK`|EXu{`k{4uuXbxFxypuMGKF=h*p66eX@t_E(??a{9@ z*3|vYq#TuIJ%ICT^u`ziT;A{TBI~=)4m}R?9FX_(6=fL6nvD7O#Vz+J8pS z)SGBE&kx*GxyyPcdzCp+$EdxffLL zy;;uOW%BQk5i^Bf8EW4wt!z1z5PTm6Blth|L zKUV7MUW;Q_rszUAxk~;#Gs-zAzPz^Y-wT0EXD3?4{)x+04|ps)jw$l$U=v1v|Kp-N z6dxTQU-P!jy%WAMHC=4qjieXe9TMYM^PFq1`bW`AL;@wH$Yb;CrlDenqwt@xS0{BH zh{c}crPWm#-?~S_TPh$>srtRFW`KVNMrv&YU+ifQk<-Ce*$N8?2nd%$Er9%w2TAq- zP|~W28bP75vX&b_{XE_18zt#TX*!@D8dwMK9qoI2dyAMVP=AqKzwG)LCP)ILG$T8Y z$cGO`5{F~|sk(@31HV9X)PYpU{$rfq!}Vw{95pTb80#}IrW!;r6ub59%&UI+}(WboDo22Ts3J%|7sHvUzvj3HpavwD4SP^}9}H zPby=r_&tEjwaG9gJmG{pojBaR8u9=DdEM8Y=Q=~ofUL~nZ}NX)c+;%OR*1BA`i#kT~(ZwaYUP!RzWPGn)*K7+L0E2DU_1tUN(*2619HG zCI~I&4J6DP80RYg{feExv`3fe-8ZWCpB8)YI~b30skLw|L=*<%c+L zwckX1axu~p>{;F#7W~lklM{ZrqbB2;huB2*Da&8SE%Pm}S6f)V(dIbN$H#|9$vH_9 zXrZu_cD#lW2(bfb8u1O^|6Twv!S%1RSj?Tf{BF6J4y0++d%BQQblz3`Pl$@hjQ-#6 z7s$k%vDToB#k&F$|L^U9Q=fh&E*P^o4QDBtoGM!Zby^D0VO?{Hn|rmg9#{f)Y7zjP)9lu1(3-p!;kLOxSyAsXb}d6A9L=Zx zeCjE{S;wApHOf*n0_0lO7WV#P=ZT>I^?>Fbz{G1RASiH-TF)P?bsmyN^fsq0#*ID; z^GC);MpXmS|5+&SO^D@}oR#0d)A3QcNJ?3yLl;<9R`!OZq#=?c`Y48sn< z{|5na|IK*S-Pe`YbssbA&;Dy#O)pj_)~kw!KpTZm@wvdQA2trbKck-z>;TcG0l@pO zRQ|7i<|3w+u`HAbtnlg3w})8*8_w(EmF^ z#LOVH;sQ$Xf0ixC+>;TTdrT~mwKrjDBY6Jv?@6GD817|H3zGig_O{~+$nx>y$I;{x z<=EJFvfbLJL9H4%ErsPYS}Z>1u?F;VRh3!)Gxt)Dp>&li!&z2XBVfe5Y_4!pMz%WQ zF)&*L8$b+S9M5iMcGilCiAl=Nep(C=Q4AIqO29QebAU#efeg8l&zOdRmoI>C zs@s@9v?FQloq$eb9C*(Gbh^F>PQVTyd1H#-9jz#J-%dLnX1B0-J}F?K8s*SM(VuH4 zO-%w{f}hPmTc)P0c~k-`@9XI30Zo;zYoi}7UAlA)(2x49jsAVOs&OI)h*nb-VAzL1 zxso#}1?;h4H(I;G(og)~L9ADui!XpNO$e@`-0P>y>}?-J4^!-Xm?nJlX8%U~v4q3? zyey%A^)p~*=1_Lv2_#C~oM?bZP_efH=NN7ebcP5IrpcrMQ<1emIf2!50Rq?az$hT| z&Tn!hzo*7N$wG$i{Mp}c?%*XJoJ;$-8NKwZ@zbQL?0UiZswVcz}{XMp8tpAFE)M z55)tqy|sF@SgJgnwUUir`Ep1Ipm`hOQ9r%?%I$9IY3Gsa#9G)JuUVkq3W{m&;+yKW z3f!Lg`)aqNj8^a8vQG~lz6ah`V%i{A%g8S9RJQU>z5l9l8=<-*ds6t+)_;pA>$3?= z-6a1|d9*&cCovYV)fV((D6Y$?zV^r5oZAayDakM8(I#p%?T%Ssr%DY71gXgt^x z3P@dr0kJv~P$l&Pbg-8+<|Leg$sGX1x~DC_ewJycKGm6>%_qlyOI7t=dm!P_5&_Q~ zjk>l4#EImuM)xG@PS4Tt9}F1J?gC`4{Y^03@3pall?{79ohf|hjy2#eH$+4p08R3R zr*CkGzL&+;-q!zyQ=(Aumm%;7q7qV_$0Oa{_D@v=1qBcBBaQ-8vwZ*a#_(LKsWb+*lC_};HeIQga0qaz5*)hy=!+Afuo2V zQBhPdC~0YFL;-1}8-Q;7cujnU+n$tCrGGT-}6MEP@bxVJX6K@+NR!oJP3`BHTdPoDL@f@`MbOi-nsI( z$H++~6cqyOgvQ4+i2IDWMEv3uOXQsCkZ!H-q2}B`j`e&E6&-m$>EeY8LLJ9svJ}f;ue=Qca z$}(Vhr~pu4+PM_J@BoF%dewWG*E}Z_1?Oam>Cxc<#$iHwUGE6)W`=Z&?=JKStPjB% zWB13ffifk7)23Cj*U_ku<)!m%^)G^p$`$?Mu;tV^#n0H9=-Qke>;vIFoX2Ss(c);epR7hzI1#~J%Z~r-3Fkr%?~eGDTH}%R5{3thq+8^6Y8eSPg#s4rgYRk&jMjXz1ibuy`c1 zTH5*YeV}|tStPiBGkxN{*{9Vp7kRHmR@U4sp}D!ok&o-4?Pxjk=G>0#8>(E>@$m}! z+XZ(uT+gCZa~Adgp-h@t+gdmtwqapH88-e{BY~JQ?amLvG<+} zqU>5E$?xsWRKxX39yk#H|8)X=qbE?klI05yItvPKL&vfS;Y;LSQ#>E>tXIec+_oGy zDEYbvD*!c98XiNxmmeSHJ~1Ex3?E=2w{IERw7>pZ81#kqaYna-+^f0ebyQsD3@WDhBX0N{p8L z4Wh&coe84rocI#E#C`A}3Wr0ir=zc|tl|xRv4+pY=Iiqo?xUU|3;?CkKzhq$PtIoJ zFv1S*i^6E#Ukl}TiYMoxp8MYKayk6xvQXCaunAme)Xerw{)3)RJ`?+Xk>)nV>l9~| zU6nfEwW&YPo)xgO&@}jPpw-#-%=#idGiz>p$wf^j6^^vA+a?vE-CTsCmtv~fMB4R1 z6aqtX;6ud4N)8x{EVvU*!Sx;*9c?+xYuSAtJR{7VM9EP!VsWn0gBSk0yAr0>Rg4X6 z0al&q9u|ktVrRu471#`y8d6gJG;3)7LV5nGYAWFWy*Zaj#fR5*GV`||+`s0I#r-X4 ze~^5&(&g?M!~5gwVr)lQ!W8)?IqDv13GYT;<9dDovGI78Eu-Xv#ivD zLuTcot`1@Ek13|iEc!(pKG`AlWIF4MNBc*=D5YD~|574+JtK7=JMZS8VOy4W{)VG* zjOf--vODEP&xhZ}gzti>ghqzqV!tLPSOG?f1Lm!f@i8`X)Nzzm`ee>i z_oaqYg3hep%{{#-?8+j-7<62$Zn%|~wKa^pti5^p?B|3b#!1ssP{(msGx@6gmlMU5 z3bi>R0f8?Jom13CAI(m>GG+ADW-AHWewofxJjT{`vTS(c<1tjRysCLj^r6L^CcGyG z^AWc7MZVSC#6+FcyQ~KG>3J3{7i;Bl9pM8M`5!tzZ%5~Nc=nL+t;+=f_jx@WrFrt# z`KJ78?OR-254nNN*PVZKy<(A<>UT>V#njr=o!Q>9Y=w_*zJVgL&@Cjn$hRm&$zXB( ztT=6J9nJ^U<@2Ozzjl*5#OZ8lYrX&eGgh}I>85IJ7*EgS%BD8YVi)(p%6_i( zRJmeq%Og#(N&~fQ>qlyt%Tszw)n7Lel%&9CyIaesI@xB<{m8Ncx5f4r_@|cSS<>8U)(|GG&1_GQaW@*P?#BwS537B2b=kSYwikQc|{R32Y zF?aJwvT~M(*J=)xDq-BteQm; zedVZVij#y|m_~$F1&co1akXFh`;0Gq1+y7Fu;Zhv0TGxP{7-kKyOaA=gw|wcdO|4p z&BQJuEw-T3^h?%3taI2&!@f*4EiZ3x+?+qLj^@n~sD>Jvj3R+F{t_C_A}LU>|wThhxe8ueC2#>ABF~ z*J#r{i$c{p2tF6iKAs>THQ@X`4)Z!=U`1&(i{!`j^Vr~^9!ZX3{~qoFkF5pf8#it+ z63OrE?Zu+{+c;}%4hNV6iAYX=NKf?9J>2#^(#jfcosw5r$-GswyLlA0KJ(e^yT3BM zh{)I8N(?PqdRkg*)3pN`M~u^COL`(0^R@DKO+mHZOS@aV-~G(a^zH_2?Bb?Km04Yq zPIrZ|Q|2o*hG2%qmK*1~T4R0OSOlxv8@z>Z!&Grl&elH@?_oi zyKw$IeW%!74sx6G9+93n(F`(fCwI{Bt>Z6kN zTzO4NoF9($JXks9p4~L5DkcfN;aOYLlo+mO<{FC?q8C=IzU6Wg44=gNRhDa%5NVDc zkCmfQwoAwvK3{9qmlXkx{Ez~(evMBR$Zc>Z3X7Wep!NVCy%8OQG{I~mtnERLdTUN-sZYheKewETzqD9=^Rv;8cX*7+`p7>#gDfHxT&(=$Q#!CL^ zQzU7DbTBQu@6jHCfq_8^n&(Lw@!b%Sm0F^-C5;-NgGt18a%n`+eH|wdk958IUu8n@ zlQYYj_%tEPanwVn6z}AcXi7@@Uj-*dP-oMD@0U3qTx>Y(sx~9R*Hd0w`=d_C3b+os zO1oxVu#$qZrTDg0@cvruPf0vsiA#sp@_UKguKGDGXOYjW*L3GB4QjOWbT4h4Cby&` zBk9a4z-PBy8h~E5ucmi-ZeMbKpd%_v-LHz+pvW-7OhbbGf(oH6D&FXwAiG~6#sp4M zzAz;{PoK_bL&otvmou0OqO zZ0AZbOWDO+(_BuH`*oiGZvlZ4?onE+o%65^?VIKXP1bsJN2D2lJm%xFORtScJqmSaXgf-|EqIHKTMeHh zO0d1;9G<0OplnBdW%!H2v7BM796huer1`XX{*q=5cfJfhajTO1j7-gycX<@=K27Dz zdO!BEB3>urVu6tFyl(U!<@gk&WvW&+FeE0JQ&Bl&0w~zUHN5S)NgKJxuR261{<5<44~Fi=0_v;K zj}B^--KvvtVw!a%vB()v*v*b5=^XM(9E_$Y9h|DU)~k=a)C|+~u>9tY+`=ho^B9~< zZplW8Aw|qn_tG=DKRy+U7>ZX)iUp5Zxc;5rLn5zRN`>;v(sdBHagauzlK%1L+_Y)4dO*-@2=CAQ1b7Z5FB2M`f0%jOt_%Uz76J za7OinqOlZqy_#ZkM;X7}>^WK@X8n90ZqN=T*ynwXJd?i&+eTu@lKS{OsyE(G{f9)~ zyLR>TB`bF9p?+*ddzN2$Yko~szS{CNzQi6PZ9vh1m|$7vF!93F0?H=!g zVF4$Ch=}|i51_n!nI$1XAr=;$$N;x(%Le*d+few4*u`?Tjml>{jh$`jt-F9KwNi!4 zwJH?q$z*l$R&Rt$vYydm$StN`P*e8z-fvfaxoKyFKZ4aB(XbJVz&7$4%?8%@<~BT% zxCPn$H-;<{WlL;FC&3coA}8nTcxYJ|X-4L(0+>QM@HTc67aw=e5mn8r>kwcaQ{nh^ z|5dV-*^u8Qmow;8V0UayFpDK|?smd7XgsOC)tsI%dOdt+F?PEmBsh`Uvoz{)-AC!_ zL+(^}Wu3pP(i;Z82P*1Mxl7>ygx7R*#3;n_zZ<`RO8QnF5phF#;rD&$VxJU>CfnZA zx#Jc|{opt%%e9M0>Q5VUhCPb&@z%xG@uBm+(qn>OKB;Os;9788taApk@*IqMES$F@ zB4Q`Fl{14YW)|5pe-f!ZGjvwiiEsU9UnDO)5Ye*H%U#UuWKuDCT^dX)s*706uFrJ3 zejDw(bzYblHIS}GY`*xBDRzPmjScEQ`5`}OU*hXLUk0Z7Sl)97a=9xT1*Wm;=3PV= zwnfCSXxmo-wCg0ea3Kaau|h}86x`AKwQ6;`s{Pbzuc%At@qzxTyT4B<4YA&N*GU<^ zw!M+ro20Ij(D_eSwpx9E@zf(U9;?Wx7(-7+#df=nE&M8xD9RokEjYB+o=3d4n0k&M z@2s-hXCVGT)2B|OU!b<{cHXl>e_H!$=+5=CEG@Q!oTDs;Pglp&I2W$)R^xTlE?2f* zNWY4Wd}&THY%bO3q^*~gV6mC+nORpS_1)9d^Q#~?-9zEy+6XFay1g_C&ASk_m{%#{ z?7V+)K<3(vQ$&V3!k^{)NRXFf;1RjaX^!&10#ZrOF2 zNWZCG(+qA|cp$7cgK=`gh2Lo-?9xP-V9jEW$m~jroQp;=xK;GL-KY=UL}aD=zFuLj z(5aO*cKyi_C2TlJ+Fwb>`9%E8IJM1!NeTe*W2qAFJBj>pMWL6z4o8o7M9dR$F7d}# zhCZaH=shZAnZGigb!nhBrk|N*(aIvhc4=*AJe|I_;=P8-w315WT8VQ~YUyI6#M=Y> zt$MWld()ZorQ|Jp!fhA>8>0K>73Z&+I-XaMG$yZ{7s|My(P{4aAvp23+nrnyUO#Tf z;R>)}jp!(c@d8VH!jb!C*BX?0qO5x@wb3W?$5mC;)Xe6mml^n%lIcfFIc8^Pktz4P z*ypW&&jLW=nx-)-)vGGRo(VT*>a$osXcjRV8wz66jH= z)Q#;hLW{4G9ri%d*(qv!S260wjohI`dqMKcyy1z53tB&ymTU_}xw!65cla3bd)FkO zsN$BSrL*p*H5ATdf2=5d6%<=hT2b~Yy`sYH!5E<&SC3O186Cm7#-%FOEk>d(LMIpi zC+pqz)p}s9bLV%9oVyLmDpHQ@5qb#IxWCks9{hWEclyYP^4^QnnK`ylQF1KDkKObY z{ZidEenjw2Yt3M!8LY6EB_`6Mx_h?{nM(i6*XNH8gv0$mfk?_rgV$>oBif=&LBlL- zZ%^mGvC7m~ACM_{D|THVza=OHZ~CyEwk$g0;S#giWrZ)~d8NXGnr@zksBcx>c67cK zrEh{fpJ9L69+%-SaS+Z3FBzTVUqDCU>*%+SZ%=4G+FP90wHrX|29)01U!Qbc+gf!_ ztsYAn*Sz3(gtDrJ_cRI)_S5zDU>OUOIQ;(6n`9!U*QG=#{k}>k;i+vItqa441BMjM z&&>_8Qj8XRD>+Aq$V!WjiVZ5GMX6h& z=EaY8f~hBAi8SmV5snT^FiWRUesbGp?LO{1 zYZLD+BJ~S(Z0l3Ax3`Ss36HqYcmj0?d)R%_-D<|G(59mJqR`t9@So|xC*4-Sd*%o7 zVCE=1+5Cfu+p4{)yU^lV1@V~kZW7RDM0R?eej3&Ilr6_PDB@fg*mb7kMsFm>;h6_` z)M&8T^FP>o=ysdca!KW0=GTb)=~zXy78aqHV@^C1!(2^kcwxjzJ-l~8!5#X<9TAVj z(b?xhX#0OEYms1-zV8J|n!kYJAt(LfAo<;oA3v0q0bwUbzY7kxaA9XVU+Dpu@s9i6 z*@fHcdNCTQYR8so1Rkn}hMCrv^L3{w3a@V0oJ+DAy#2u6Uvd&~K;Yq^15f0uIK@9G zXS)_Rbn^}8)$q4hoA^G&{O*rp(xV8cKI_d;TcLEgW&dsewPn5Ck4UKfm zKltE(PwXtM3pfs9eZmJ?YPhh1+KjSiH&2bptAeb2S9Xbjd4t zI-kK)_9J0w$<-eN@GAdi{$+ua5iMcq%_D-#u^J~WL;)E7{afAG?#C%^M~rDp;JvlO z=v+sUS=A7<46m&PMw67dgyp+My7{yS!KZZR7%EBbxHi&x_GW7-V{L|R*45Q9`tgC6 zy*pwIH-$JfPHhetMMW8UwESuCmv}$QU^^lKa*ka*7vmBlzFnEffPcu&&|nloLLSiG zd)NAMh7d5<++Z?+P^fmbStGqo{K?NTDn$C#N>KCd8tPQre91X-DMqlo*R;f0!<%o_ zcX+L$+QljiIMM|<%9j?GZcsem{$Qa?5#q9WMX)#+rcIE*Jk zl~xxrpXBzQZwdrToz;lp;gH2^w{kq~K2cv|UBKSd;P-poY&m5oPMnAt4}Txh5ie@O zQ8P34$~-G+_&;zwA)PF>nLKC$WCa~3lZq_yW9C^mtg16&U#%JYLrBJ) zr&c2mrd~x7WjLHbW!NjdsZB{sllv^~+EEMl6IGUTU$zz}#4IS4En5TkVcHdx^yYn6 zbMY#gKYpC091}_Ia&dLtcHBIR8rS6*DQE1#`mn9ey7y$Vy?9;}d9)R4Xo{n%?JKq; zJcM~S1E`E&y%%E=92!#1Pe^?FYfuNE2F_24vC%Fwoic4p2Eq!t*<;?j_>pfwMVR^; zcuA+wql}k~jAuXwmrwkAEIJp%xb$)Iym>i^m)J^IYH{?1A{om4mNFDlocOHrIvR@)fJm1qw3U6W+k zIFm;&EFh!VL>snESJFC`%_W6hyw>}*x8!nPifB*vF%%z@UthZNwV@I_qr;<#Z|}mu z=3y8$yd zKKXg@Wc*wj3Y@DL^Dg$mQHvX6bVf;0<3g(3 zJyTy6ou+=hE0rwPr#cg3_WCT5grf7pHwO2kb}AjmA0JGo&4q0Sij@3B`m8=E*nJO5 z5`)7?!%|8Y7n)u$(u0Ai`f3q`FfIe!@&`@{p zh3u0jJ(Qy}|6>LytuZ)3S;#j;Ok0mz@^w?_=$G;6($bLo1z(&?=@He=tEibMlxZ7YK*Bn41&GB7D8gr}s? ztDOfK*z0)nbmIPe#_3#I=S?gBohpw$%s%Jb4?5`?f{~mf#%rUlf|8c@Mvp~%d(Ary z*2o7yB8G#QBDV3GN>bA}Pb6G~WROz=@>A1en zm>%r%b&Dhij^?nVB&K%TR-M<3C;4O1fK9F$A5o`0;6ka(Fh|c~Zg~XncvXDx7&R7F z7I2Lo%*ht%F)w&ZZZnOXD=KIrZ8%2OVhQ`wnz_1WnCLa(AG6HJ*$iF6b~nashk~&q z{#`f>j?h>dVUIgU%dl@>oj-TG=WQnR2674t;VG$UuPcIrPVk+MOsLw6VfG)kcmpJP zIj|4nezKx+lDvxN_6kr~9L_N!%DO*o&0+p)rF-0aE0zukAjciaKg>o?V5PnPm0NWfBa; z7?$IQSn!k`R6|&p)FO-x>Vcv{w{p)rxq#31mZJJ0;cBXm^jguCN>wDeLyc^c6y=B6 zo&I`vr|+5FDYP$w{HI76xv-(4q9Ws&jpD_xUz4_t=xJF!gSXV2{m4i&3CB>1k;;iw z78HsP;U?=+u!vOtib%M5^QHo@lUpAvj&SYFj&YK)UprE5`fe^qjw;0fD7e;A7SBDi za7MPhbhaQqKPXB8J*JjmE2^!#Lo5u%9jf;?svlT?W*`?@gE$p@mA<^g?36AO8OpDF z_i)W}=f{#%s6i5c`}*RGmi5J;!h6R`bDAgJcCI^Sdz~#b>ULUUd#bMTGB{D>yOp0v ze{v%4M&`_)rD*-Sl!>0c8KL{aCl;n<=q?2l&S>8ui1C<2?%#&S#+L3>c?QT^qc+Im zzvXfYC4~>MSE3yi8&IhDvZbDVcglwlW zDEPG*gPaA1wBW5s13lXr04^XcgNj3b+7}_6&y7Qg$H2q0?DttVj4?3n?{dfngLz{r zKQe%MzgXj9n1(*kwly~!v`M7>E;F{o6+M1o{IS(uDI*V)o@N6r*R$+<114{5JOS!Z zBUrh8#(G5_Ls^NYE8b{tX#op4f9C?;CuX^iL~ydJd?*RNrwh_?BIRf$SS#kt2zF2E zh0>P>HiyoQ&Ci8+de{r~F0iq=nR8-WO6Si1J@=<1XpN()1sA)Ih0w~Yki{ug-|Q*r zh`bVQ$&YV(;Syf6>7gPrGLoHRhkE3&+V8w`%Q2Siai1sGb33@0aMI$Y3Y z@}v0+C80r#ApW<$3*^UG1`un_c3biKM5xH62+8&r=K)5P?EgsqlxkO(xu)3V1 zlv+J4u+1F}mAXdnH?01)ozok|7oDU23`A+vXD|2#XA>iD_F!?5>m_+ZHp@fT^`}p} zDa^SIbDTJi>i7@?3%+q73aKyw0S{V)vEb6lMfy%GHJi8pO*cp2nGz#@t5!dBzp+PG zFHT632E`LenvK;)|D0*;y@Kk9_yA0cj-?OE9gE@F)FoH3FsNjp5LjkcE}^DppCd*J zL%8{n2}s}LO%4O}%KdJ$*hTcgJyV#0Q~okQu=XSuj+>u0bF0OCB>?TdinU&t=PcAE z52t3|l_Od%|MNij$l-zb$yfYOEc;yLskKZkGjAdl&&lcOR_yb_wCcVWZRxLV@Sl84 zkRzciUq@hW4M4EFNHU?mrmjnxNb$PbP{bv4p+VF~Bi`S#kujtyv)oP?v)RUa_+;1} zH>CAL&_R4w*2lG61B<929WenS5v^2F^J`OW80UFSpmK+N|19&SM1a`s7|NaWo%l%A z1O6tN;f42qUlZT4SqY;w%s*CfQX;i1B!p%Cvuj&Z2!2XcH?-bmY>C#xZtM-l zb;y=fUq|%=qNE%x0P|*p3eC&B&BDJ1A%{sJnp-L~bd7` zcT~7dwzJiJGIkMQc+(%_OsNq@_7R)p5;&K`X}+D?usUB4(FqiCU&Z!xes0p+ii*zN zl{X>m{&|+gA#jQrd2dFXASRw7&F=-S0~=5f*J_UTY8)XQ%YC_Od(!A=H>wiO-;Oyd zRBhwrBIF0mWELts#uU!!Q+ujA(hnF58%V8?=90FAPgkYn zG4TL(HpAgqPDaKpq_$l|ilv%i0$X{g7IT0(J#QabZ! z6hu15j~`zB;GTRD5mfb!Y*y%hyr-0AJn(gfvyD()kJp-%3J}?)cnvJIVE+gr6 zu!e`pN~EBLYco=LZnmojm1qsr>%kZ^Bwd#xe0uY;(CkB?NJjnjN1H*#q5ot~AGkF# z)jQKykeow zxh)oF7W=caS(u5xx%Z_ZH(hdmyJ*QhDx*8Ms0Fc$IY4vcRL0Ui876fbaQF%*5sG0d zGTVtiPQu#o-BrQveNAA=V2tZYK=K#3BMt;Ns1=wHV5jn63*YfewfOACD`q}H{5ZT8 zN+qmaxJ1oRthA_W3y#pP^Rgr=h3dVf(g9k*=vY=^%Fc)iWjIi9ymhZ|dH?=|p(j6H(bDI2wOtYTt@NW@%A#>3ADkzeA%9=CbTD8TvUiu7*8 zNe^I#JT~6<^4p?<@z&kCTm`fy3N>CH7x^p@1BPD`ibeZB&+dF*k_wE>J|c8%H#?!1 zvoJo#><=8?^<(C0raJ@E$RnR~oL-f1qxe(&=lOi*Mkv0$U8mV{T{lw?U(N^y5sLPFA~tR5#rfTJ1+48kCTuBEol zRrb-N0IOocMX6`)s5=7Lii)2`%)RWeP(L{p*|Sd8%Q96 zjw=Pq{;b{v-Bh#z8B$D{S13&&NF%UJk{E%m8|78XqBTLaW25o(BueF#61_p>AI$Ig5XC=`pV7!@b1iY!>(~lUAF9fd zk`jfFw{<{_f)wG8b>^A=nJgT~onus@qM!(7@%~D3A}SZPUppRoVGJno<5x9&8vdrD zVnjG|P@29;Pq#RffWO8QrzrbjR!=-~J~|Vk!1!!Aj}~~OL`Ji>Brav(VEV*NJv-Yi z@>A(ix@F;Q>f81&TUv*UonnWuV2LXG;APaOdTBWo`i+H#*J(hM8s!||KN zl-^g1i3|&)BlsRcirrLH8g2i(D+2A*RX6UvXbxvq%Sx7q;f4K=e!gfdAXLiAUxvi! zQ;PRX2v1jC@sFA4`> zEy9i9z460$!S>pnf*LSsmDjf!I`u{GErl&Bp&{{>!`{g{G{d}=e}k-K*%CaC^Xx^I zL)@uvCVJ^ebL~>4MlOO%4@AF=_$+WSPJ6ir*@?Fx0@uOZ6@`lU^?w7=fbi0INSgJs zzkGK$yOSPNk_demd~1jF^X_J^o?~ho8Cq_8v1#4d)TH(8UE-u_@tn_DE|tgKw$0U( zrH}DjYZQEoph(__IF1wWno4$e+oTXm^9zEg98j9gN(|~>+@(B!{$>k)PI(iiZ;%^xwh=85A_k{3-~A?)Z`#^@nHl>zL)oWGsuur^tKZKx|>osc`F;J!pPWp?QzgU zr1`Z`zQY71uX*uZ_Q5>9e7MH4%7I6bbi8y*m2Wh!FeLpSli@r>3EGhuoQFA`^IAmh zBr1ZLfoJqY{+5BTU13#(QGu7>q^0HLeCIioddR|pjSl8F8Qf*zC-p72BIpnUX}~JO z$3IAOyRV&(7M}42salARvM3WoElkak{!t9GYu$LNk*iw*!p+d$4{k?QzHJ18yBat? zn>!7`%yKf{_b>qXKS@TW?5kpdnVjfqQNQsQgL8&AHM1 z=JAO4Iv|4k)>{#C+Zc1=1W9Yg_iY3`&^cPPO5eHyOj3W8rXUOzW9#Dk!?er1US+PD z$$lSrmQ1~`7-#eP6a1YBVE&5_3KA3uU-}0JjHt-nRb%?ebAaOT_Ym-75x(zSZzgMX zb+x&3oxeC$tN{AM@84|*(s}aK87{n#^AZbq5QYLrwly=1!Dxp9o^P^BNab-}G!-Ke z$(fKvJb$dmSrGX~Dqb0omz1o0dr6SrcSVc#@}(XdPEU}pq0{(x0kfyxl#i7R$kKIuP&5exM;`T@0%{77Gz|NN< zUec34c)lL)6=LW|oAt)C{n6m4nEfMPr@d^= zJdCG1QYSPthwL12Y8OagXs%`m*p6s}{ch~h5jc1=0&~Cx4^q#T($$IeqodvI{7I0( z!9M{N);o-G>Zs;v3uT#uXMw&-mB;4fGtz7AN4`_9FNA*lSki5BSJ%R1uAXA%sKy0% zxy)+^Gm526NLXv83!XfxF+T7pgM#+X+TcX+%n|7Lq+h*y?=t8zyfM2&N&u(rh{BMN zGkB$Fp}x3LO<~uiaPU|bvhh43qIj8TG0Q|P_@njh?{qRq9dri&o|Z>C267z_#oT{BUl=ZKGsaiQy!}v6aC_Z#GMWd2 z#`6cEgU-;foMm5|!%E3WwWKL-QGD!8x*L$<9Lqm#XOO#5L4M^bc5xs-v1qxgZr^*h zn_P^W|c&{o)NrX;0ZZ&dOXG?tf_uJBdKdWcZN4hsnD5|#>RMBR2=f9px9=bHa zq{LsxLT(6=bQ9l{SsYeQp|Y;v0=;;w-u%n*Lc;7OQT7SB@>8Fbc6LlVg3Mbib1!fi z|2(+hlBQ}Wd?zJ(v@KfFja2Qg_1l5BDb)UAUFMM#4?BG7mFsn>n9euW;{D6uigd4D z-xnKX7#yQ<%qoq-4Yqy3sZ(k}2F5WCXe15?b?XfaHCMc#ll1$mf;U@QT7E7sbH)q1 zMgl8A&N+QEpVw^)gqQEH3Eu>NV_rUb7Ju4(9a4+es@*QBXlrj3J(}UjX}{VF!@Ad2 zdW*0gewy6_F|0-GxINs3VHS&$gY_?VXn(C%T>ZU6;YTk{TaM%n;-`3SQY4_oI}r4LZnKq*>-ELF)fSMR9&!a&Spze$w}We2_5hbx!xM0@&)e}_Pa`+qtOHTO11Nfi()bWr zIx}d)!uip7`HC&~?6hMN&Ip0H(#1{c)1uwxNBZp`*Zfd1IAl;6Gzi!OF<}bvLVdO? za^&zFX+3`3*ja1iY`g%4jpeBlJzi`gd;C2(2i7TFiH);S4UyNKmUm{+#+XfW|DpY{5pL8e+V4>%ko42jo54 zeH+bSbrF#38n0Q6iYTPD^Vtu($Wm&78rQe}`MyPY+lRR<)VlB>_5J6m+5P}GA^2q| z_43Lcd!^g%YoAz6+eN!QZ-n z20fgdRdGBia?p3NQUSmFWuzU0i;4di$Tze!!-on#@kh$v%~wD|e!1Zo1xlT$$+n>d z#Jb?nqIhPW!*1(igM}dcF~68AwYpUv#ecv#>Hp->Zh1-yw-7un^w21bf#?^a^->y$(XTd7-=q|6frs-v|>* z7+U~r7)(Xz{4QL$U{R__EH4PfAZH{K~c04d1@C0*;E%A~O0yxqI6Gi#R%e?IsxsE>vP! z2-u3qkJ^`A)&^sYD5#&n!2dOZ*=UT?K+Y40l$y^8i-iOXa1Pr9XW+X|A8>EdF|1s}=ANJo!#n{XOP49>4Nsv5UN0hZa{9W1a*{<%V6s?ivVGiL20KD%+0YnTu zcJ&007f^A@!))jG-d-CFm}hf1j!C==gQ_l{Mwi2=VIiu0#3r^4q_NDsV}u?Vs?!M! zG;;zjk2#S&SSv-OZES1|Q7xKJD&Pwe#N`pi4-!a=@VJ2F;dl9BrH2h(e<;K)@BX6q zsOVJKGoDmnLdq>e|(M)D!2k+%>ME^9P9FU;T821G^r6-}?C^J_p3D zT5i6mB*vaTlD5adAa$<>5)v#q_<}ziS(Us9jzeJ}D1IFj6vS&crY|We*?#IdOxZyw z&yHL}u!+=R*Ct7q+t0@UfUTM~pJtkh{8so%8a`2QiV7Y4{(U5u$+u=udvhRi^c0xZ z??g4iQoCYTdE*Md&4*6$qe@Kq2^8NqT+UM+YR}%@s^Lg`YLmv6%%?q2G+JklPcJlD7F1lf{9Uh>nBb&PwG;4;5sDp)_tW2{b9}^dv!# z0Udq#9>3AtofzilsJPdNvfV!ObdEp7Ewe$C*wYX+ybt#m(^_cP@G24Z*YPjS>EfK8l8vi%X}k2&p%fm->pB~A4EIfo+xf)-v067FC^kx zJ(9zKWqZssA9*Pp3_!==HxM-O*UX*($0sGtkVVmv5i^)bAPj5(+X8-htu$eA5j3x) z64}69TesD&;(dPY81h+5kY_w|-dfP12Ga%bb&$)}%uTT)RHuz;fG7qGyYc#YNBMds z9s2KuGcFwKwP6b5h%MIyHW(snYQ!$9)O>;AkvuYO(b4KhEEG2DStr+jd4q)l+I2X| z+|Y8R6(rx1*GT_;t{!7^Ui~SBfE$Ww8Xmjxsd#TMFVpT6*)_TS&epwBO14F47s_iBn-95y3WPv+)*ezQ@ ztx8;&D$Q?chsj;C9Z{+r5smBDLS?Wa;x<3|fm%#jT^5Hyxq71kk)PVEp#Nuqp9)=a zJb9Enwp5<+OV4nVpQc64R10mpyV4WXc26dm<0~xEVW`!bTB=>-2Y=U3C;wmm!s@}7;94z=$SV4@sFej&rVH5LJ zslHRA3?L8(*4xDBXrI|1`snE)l*%VaA2RV70QW30v1-u6exrr=-}hRLNE<->&MBZ+ za#NpdEr0Rm@L)1b`On5eHWd%@t8cOktMwZQ8Dny*YaEAjShUXr)1Zv0D69or&YDu} z#x3=?1cUpN)~A+~!wPf#xvXR}oe?lC9vXlxtm~kQnjdu~L{NzDtAN6l4Y?+~wuUSa z($Az3I0SGU>}^LuNaWet<0|MU5>rOj7$$h|M~M=|3(Tb`#GAUtm5 zskvD9m7`zHhj&CDRF9NL-&7Xa$pH_y=9&c%I&v?E1=B;q6XNUl^T%BSFkn$3BGPu6 zZV!PeeG5~f*JDiX!Ugoe<0b3<$<$IJKh`&jRX2X@l@A0(2*d@3UgqE5cwMG z#q4aYe9h*A#``iY{kd1Bx>N5sxH^TU;l6^q&e`r#*(a~^Dj>h;3}^Z;W)9&T-=Y4M5%+0W|O8%L4R0WrjDa`pVbtpRQfWydU<1;eeO z=-@8182+BMsxQ_vxKozd2gPq&ShLYCknfWD61(fs zH4vCx)aRYF_w_7JsmSzJxul-=yIJepo5Q*_;_!_Q%VuJ9OXL5l^LMBp<`QY^I{g5H zCn@Z!*RQia9A5&T#g+k3WCQjA6@hEv;_ol@(>BQovZu+>n#p*?66}vXo?dk&s@H4x zMimcNYj3gcShOwCb{v+2Sk{N9E z8grPSzjKGvbLUqG9mb9z7KD}wrlte9k;^wx#6i+!ad@g=agh~h8^g)}hd zK`bDe-8)-yOMq`;tfu6G(9YBH+q9zX1dr2Ln_1N$w6OCWJ$wjUzUe6`%|G*NTn4cQ zds}G>Z{+0V!#Kr`MvU>A;NS(aP8`9b7S_6C5IR=oesh6w{n}CD_Nxig<6?AcE7T68 z&@gUOpB?61UdIfR6hL=%9qt=$qeCHRKEy|ycgLZZYjqYEx~^Mcp89GR=e8hw@>J^{ zX1=c?|4*t^g-2?8+-3mfXe3GgixIv`{expFqEdiy&gptP`_wdylxy1Tz4RbAlXFvc zzg?OBe4s*&hXj>HN3Q!W@s7?PVMC6QwI;_^0UyGAgG|rC8W}dr!j( z05~rzxep-lYamEggZ2&G@fImv9S`^23k^O$n`)`vookcou9C5%_Uw!=w6WUEnsWR# zE*8DhK=C?G?n$&{hQ+lvWRDEW{sh-_rSynz{Z`}WyRo&5rz7MUFaw+u0eMyI>n5)e z3o%yK*C($l%bAz61&i8Xa(U;o?tmi)_ZJ;4wGRelg-(x1lY|O#B(?QrPte@-;hwoy z<|U|vpmFCI$Bb@l=TBU&S{(m-k9t*sqC8;suKP#Y$(WvikCM&2CXKukUlY~xZad#; zA{L{|awdsP3z=nkF22*EJP85eoxUq*{B-W?&zUXY`QuaOuQi7ZE`*b%+9acjpM38= zQBb5PYh`QzNq&fnm7w5|`Go{3cD>kMJL~^)m3&9Y4a?QzZ|yRsYu^Hr@6_IMb;OIi zb=(U}52rj7P%O$Ww(o7;iCe$^lK)zm(Q@d$_=kv!(iVev!ylRbP3{TH(gX_!TtsEa z7X6=CZ)2Q-Q&Cc|C<9t+7~Aihoplf2ef95BNxgZZgqVI!iY~cbU8g{7s(VeDq~6pi zw80PWPK)RBmCHMCZ{a=5V3GCSqX0wjU{3u~m93>w>uKFNk)5$B&ow>8edrba+~jCQcP2CO)08+Wp1EBG%<#mz+1+2(g;z&&A6<=`qJZIucL3 z17yH49F8+*X+BR(^z~t+@>_c&SH;0fsnpim(BsHZ??ZmxntlF6>+^9JM>~}(Cs*?< zum&|=4`1gkkfInStaSQgjPo}0H=YwFThvUK;{*g1{~O zeaN`_FRM(1j0v8EtpBUL_l#<4@4|h9f{F;JfG7yKEg*eM5d|rVqJs1qsx;|EARr}x z1yB*NP^8yDLa(70MS4eSLJ^Q&1wt=(#&^8u+;Pr0_kOycZaxT`u?Z{5TI>Ix&wQTW z1lizAXJ;zY!Fh^MMK!ImOvMi(FEg%)v%Uimcj#1SmhI_1!N z^U<#z71Ruz3NlXgBYM`Zrx2YRcVMSL`xG`;m1brW8{-qI*Eqhwp2uauiL}xDLOX}^>u82HKN$8ZJlJX< zez?iC(gSGCt5OH~phDO8Y-#8ood<5D)9=Q+WH4KyO?2S9Z_l9TGz1^!-vF_#$)RP}hDY$BlsMTP4%M>MqO9F!` zp8`SY6Q+owP)B!jsd6>2gu=1PR7{{L;>Io+Fm*OQPR(|=og(oW)3>Iv;cEk24%Cha zEI0V=MErH<5eY!6 z{vo_1yhEbjc=;dYSO}6%VJ2-keVZe`djF)%6mdlTmebE$e2t!K+BcSba(6?CqE2uD z_K&^m7lajh1OON3z1^MRDwv_37QDB`C8*)p0XI*5K3@`X^gye-yZbkIU%xC3m$z^n zZas`fglmyFZaF5j;qy<>*a%3NN8f)K`vstRgZk%~Py*G8NXtW%`Mp>Hoz9p0o|Cq} zFU$^?#}`@-Ha`*V5e_ArNI!&7$~-pvZPv7-z%|+ObwB@L^>I^k8+U$Xae2T8`|d8v zEPqjJUc;qpdhXz~fIT415gk@B`%nnBj#B1aglVRw&>YOx-(&y?@YZgBgy}>2$AFUt zYuSLU@{=bYCmQG!lIs;)-;0M|#zZMjPR?E<9X+z@FtIwR&|Q+Z`2KAa5<`Il*3_&D zcnWE)(a=gJLvs^rA4X{&FGnq33!q#2QX4==#*#LZz9I&?xf=trjYLpc;re&3f_;~% zJkUer#1|awZht&8z`PWT{Pg!wV{h^Mf_TYKP2I#xIrmmcu;=K{)lL?TBrx#$5gWDZ}W_>eCXXGYi-3_5u&~$VDuo zL?DeSSGUmepWV%6!TdJBDyBec6^#2p)^F;u3TejkN0UQOsXAT-xsEVt=nr)a!<8WnZxzQ|K5$w=v$O{ft}~%F{FvH2er#J$XyZg1yC&k&$PSvye<^ zhMm0`;YvufH{yCMq@9qgZbOt=P}yF&t>*HP*k7vv7BI8?tYZO%XG4g*wttT>T@a3t z-D~CUHnbkE#(S;4yK}HL;U{o7Z*<55y6Y$Dy)k+wC!aHC_}i3Kj~qEVsa<{G4!mR`}l#XjrY=-p^n7Zv41ep4ZilL~s2M6bhl#-~*}F z2D<29<)nnOcxbhKkoTWhe>X5CHG_Mx_5Ew{Nu|TM+s>)mXz1S0j#gJ^>K5_%u&}TG zsFRBvg{)hyf=aSM2}-d_+^2@o4A3}v0~=!~R#*t>MVGnBwloZ6YgUc_f)T}i&8bZ) zVHCQM+mb%S+T{*hrKP36gpt6K8GS9vIPfLAb$&E*HlJtXs~3BLGX7c&KOQZq>=4;m z7=Pa|VmIbXokujA2WYPoeX^rDJNut+Pv{VZs@ z(e!>G8K`v0jT`pk&*qS`jNhIsN+Tt|eTvpxhZdlI~*pF=~h!qN}{ zA*>cyqOTg1xjc7uLulluw2B~I|K#FqigvAC@)qYkcbd`01^`Yl?9RA*E5nWaWQZ6s zKb7SvOoZWQkU=`0{8JbaMjCSNXaZ(j-A^?HZUbu%QnT4w*rV<7Jt-uc=M(Vny|-N! zuotR$uE$WwGMw^y%Om&WQACf%e=2*6r~CWMJ`sm|)6(5A-_89`9dq0J*I+3+%fTTi z;kxj#dT&y2Y8C8O+c0i^3c{;cnY7Q@UlV%!r|BCf3xC7Q&KVtea;8GMqFiJntHb&4 z`H3Z*;)`XUnVd9#r3qz{pMWy#nPE2NBPc5)pS_j1rB1Ew3cwh?QpZOC(%}L$6(1fg zyIdv&ga?Rfs+UoUn7_Y-=n?wwVZs~r*n2&hW02O6+sMd>+!xu~tUg46#1^XhObPGh zwjE0L;6ge7Ph{!f1NDShz6+XRQ`t#r;xc%U$&t*+<%|3KrvwR z?qdHBSkN4!#EscOBsoA&uCce5nym`bpQ?I+wa(V^Avrv8fu&pS&afoUpDyy>M_99_ zF+Jbh+!U*e3M`kfGflK39fHRn*LL0YeZD*yrAEguHmJTar{P}ePIjiP60;{1PtVnm zXh6gGOcMay4C}n&QJ^Ljy#EDZ-)o)abjpl~1fBnWVyZE0j@~PCd3Mg`565o@*opk5 zme`g#)l(LOq9>~A&7e1ba#ONxj8~fi^osfZ z{VoW9T$SX~6LL6AnKp|Yy-@69q@)nuVX_hCDMkljQ}yUgSSzFjImm=sI6d0vp28qM zd6)`nd$Pnvk`F^g9DOlCi!61*K&tVp=E*e0Xdp26ienI7jQ`0!-u%C4}dY~5=Z*i_9n3@l1BrUe=YM1RJaayxy)cs|8#bF?C zx*_`+6l>J%OM3s!kKi)YBJBTtMZcUL3|X|taO>f1KY?5ySnfHpo=mol17_RjU3`2z z^)vHl&&Vwa^R1OJ@cCu&y6$A`2b$kqk2_l&-<6b8(QBHoA%bp-M0}!d9inc(bwCvV zsGh61o$U}?X}0ecd~4Ef%&F}7%`cO+-j)Sl;S^_jcX+*=*5%ut2VaGD2 z+GBaiLKd)fTfnJgw({#6XTk$S#yA)zZvq^})VPNoLjKyI8QLC^<2>5LgD{g;$2zJs z?%Z^;vGQ%~buTGpsQ7yY^^73M*dee;4pJ79Ij6^^w_ZZ{Fe)Cl+FxhiU|;-KlyAGY z(~k?f!`A!n2@KNSyy@~!lY_r+Ob-DOehdpdg0;iZlIO`y>$&atlOS$@93`HiQ1m|g z=G*sUvs3o@-GET5BM#S&jfNkR8XAUwE{(bortHd{C5goAZ6r}&QU{O$hzVo{c6?s@4A>OmAwxR#RB~^n z7CaQ!$=Hp<>J5SoXM%Yi7^j8ND{-Dn&$S|~!{ASusqaOq87+flxMgh-$&@d?V^>V= zbn}`jnYN2x1ezH6HW1q_NPb#`d<(yJhGcTRpjcMM=?$Z`X*FeTqUeUU*`$UkZZL9n zXN+kTI9lZ$l~N5+?l3vfh(3o06vqTt+*@*|0gEPELDTC(^9i`%fC%ylTQOi9ihSn; z9T#*N_5GRTqx zCu?k{Z46|yn?IaGxAolcbSXDf`cFe0QF2l!*LpUx8E($0p|DD&fdv(O`BW5{*1o`6 zgWBXjYy8ub|B(M2?DAQ_;{NwhA_A`azYj;NIL%=SZDt0xYan=qDnyL`7;aM}vjY-k zA3mDkOw=3}g+o1S`UKhx>B2&}?@zl=g?iyF<9-;G^1&h&0b>lwvl$`}XKmjFmd{jJ zP{V3)@jOfsfw3Nl@S`yBBHy)5{?Gb)=!FnT?&O9c8}D@oaJ(Zy1#*suhN-X6Q}VN0 zA335!cRw~T^!M+pUr+#+i>Wp^TpGY-=vQ6jiV&=O{J$D8OM}RxFhDy$9x1VVRQLr{ z4~{qMMc~hKUhF=IT^1b;Jev^!2fuc7?gS|Icb}64ylX;phn3 z_pxTI&qbfRDaGcqe3fEMDesQRYw|tnIdS=1QazUX24@uZ5c}@N=kucVPmj2Y97+1& z??Ux}`=z6DOYH)Do)*))?s4d2OGBm3OOkF=aqbQwSHh{Ic#CSh_z?unx_z@;#-mF_RhkE}cTFw&ocptwuGpdyPBeZS`1z&s~f{_nv zUe#~+1)c(LNo+mhf+xxU zrd63rvwNj654kh`=y-9t&?>3OMhSDtd-@+AY)jJXGCpzmV(M05C@6VfbSdb@$XjXY zaTtjeK3Q5=7`s#o=Irx$*pRS|j2`R&PW`8H@XE&Gu->JIZqNW?Nt zx>xr;O?(9WG9S4Lt>8*kl|k^bK9+8!KI3k2w|F7I{qi2lD~k7?8y8$6ANY;%KVWu) z+nmiQ1mY6CCHc}w@`X=CO>M5Ufi&Ou%-w#NzGV+(w!2|xpZ%6a)MhY9YNY62d01hYEQo8UC)q5 z9#?QW!-VDG`85JH*$j4bXS9PbRy*rjt5~e-*dz)5I^$yO8utB%z9d9u=p9b-3yj zEw^=wk*&wP>#o5$@|QpNV(KYMXpCKtr{1TY|Cv**=w%aADK+)LW-vWVuL>?ZMK6c_ z9fPRSC9VPsLNFUK(;}xTvN*!Nw=4^0b9IOFAuZufRQ#+>FD7^gF76Yu-&RAOCtsf z-m%~cj?XnN9t6vjzkZALW@dr1g^w8#k16x$-`okwX8-XIEDSQ2 z&eYs{qFa2sMci%48Y3TphmJjg{PsyKKgngaukgE~^M>D`Y(PjTn9=oyV-yU?>Dw{7 z={pgNxM{EB0K0(po%KZrc$;#^>jH@zb~rx+@|ti-k23f4m*2nBjSu@nVu#DSZ+Scz z0M7H`7o@P=3QCMneMov@%zV6tTg0aCl*I1BMqe3<*lJ#6-9L8ziU4}M?hdgG6^OAP zHYxYEjNY5;-_(0|R>XOt=>WZJTlXG=f7ZrB3=4k-t&u{A0O`Fx7PL^YR*sx(e(Yy( z|HpO6AWV6!=GhIF4Z`JZ3MGHu>LjmDP*hQ7i~F=59PsV%=l0`z6Sm+k0-OKs^0nzr zM%U-3ldJC>1lY!x?=_9Z>Y_%6V?YE`$2!F!O3UC#c%!N3M&G~U7j z85gQiH%>v3NoiM;h1}0#l}e(KkG?q?U(*u7kMB3~Y&u3?+GX9_OkfyIN9k?4{i>nR ztC+E2IrE^TDuU0@@*y1u9%fgaTQR(J9W_!V`*8lbsvxf(J_fTmWc-D*(&h_SrOnY} zhO*a~@^<1cG48R~M^FVJFxU(T%jv%`yj zKC$y%;OuImRu&_x(|(&47sh^Im1=F8@K897&>);PDB4dq1}}d9H2Xn|r$K-l7Vi8T z%?io38x6uIp4Q>k*JW>_31&8~|28aOrc@&uyaxg5C1oUy4HCIuQ1 z7u`-3&-L1vm5=jUY#RB!of~w8s(Pe*TgP}vG%{yl>J||X?r`@PFk)n^!Afp*tAP9X6?WPgcMI( z&%WCN5(*3@Hn)s1&5!XL<)TzoMm4xT>~nSIihEu*eH9O&)Br-rU^^_J#6m0C#84Xd zv>Hvb75PvD?Q5m!i^Js&sFX`m68-JhpZlZcy9#mxUr|0drUfsh&SaP%iuXB`lGu%L zuARbnSL@-eZ354En-kTFpGmYng8a{z_kgv82uq+AI?aCAb~b=+aNgAL6tSJ?5K@4* z>&P>#+>ITVyn9-Df>JQh9zFl;Vw3`>kd$Hxi;M~ruVKtG_Ogs@pJ+gKwsN8Hcu0O_ zSS7l*_qC8j?P(u!sc7W0@5iR3?B@Kl*P*l1e76BaS?b2+Z{1BSxVx|Yml}fm{`f3l z-HckIiY+JonXNvz&Hj{f;qvHB?^W{LR?Q_1=c+R7&MsD;&|5D3q)e=f&|m6%*p=Z7 zn3{JPyHMwqZ-+}o-EIR606^na*XPOs*L>6GU>eSI8-;A$wr~ODZ;u}n#;Zz*cjmpt zyu*7LiqS>a7m_Z8g!WdiL<<_=!@Q*}vlU5k>vPUG8id*{FMFGN&Tl46VP*FKP86*e zy;{BCWv;;Uv#7zZa^iCC69g~CBDi5B<;tTKCR0*!QTMr4qk1!J+tUj4`@*J79^>>n zyzk6-h^Oh}sd+&=LeFD?r&%&6FbFK}(Rip_ zdr(*<;j%E!>C7R0fRc#%G$wu}yX@J&Vs?nr*CAb&JNr}GWwGr=T*0+Z%a;3tF@uqZ z{n#Li2>QYu((LIVF%BP*wqW+Q+!0ES3Ei~zY2@0utSosoCW*1wTyWo&Zc|Nn><;a> z&pVA!opYhgJ^Rrk?n^+ILEMY2&hA3(46RE9rJ-X*rM`+W7}!GUY+~`#57o; zPnqvkOB!dAZ;=8O*y8do%JB&A(>5hL*Y-u(9$bmOjfQa3nQv+CAzVINalKi3qr-bt z50w4hGiBYamCKgZa63sWE|hy#rHR91Ry6#!8hQI=v{M*CL-l$`=}wQkNn#ZRh<-uiP~V{(2pcNQbyN)vE~VQi~4z&CUKZJ|D{3 zjo1%gw;bW{-bANWdfgIj4&6Rbs$$_}A-|m}-So~$GlAG=KMWc-_Wo6VUrjg8L5h8W zUHfqlee^nyNYBFXp^_s79aZ-EDgAW-4cADX6J+0OO8SG|&?goU9!N6sYVnCK|Dbkj zFty)d7*$g*E@^-pv|i_#8hIJm57@?`pUwAY6+wf?PZH3bgjSS_h|1`@0C8o$9x>k} zA&ptsXc6XmM|F3dB}x-tzg32QWpN5m6LP8>g7wk-hI>IYt~35BC!~7A{JrH0yocRf z_qOtE;~(ZO43*Y?{+w>Mt#5r~vnNZB?;0js6RYy=RWL&9Ncw2PC-{%9aDIc|_K}e> z`n~-9U;8Ip&gqd@L~HGfgSX_+TE72xj0K7_PAm&r9R69(*7P>Q~~Bi-g25p zWv&H>+omvO^paMwQ=O1lCl0WhE+S8-i z_urlEOEdC`-s-#MiEa%PSX5`M=zcepGm7ymVSdZ>KR%N@{S-M4vlDxruU!_I3d*10MJ#OQP?F=*|kOUBW^#3ioDG8gV zE$|a~`1+FR#ect(x$^(bcdGhcf9m?uH8V=YZz_cIcHQF@%^3)ty7!fvLf0GLO%{bM{hu1yC~t>y4V7dmSI!Smlrb&+IbdC zs6KOi)yLZ+jhU%rKKy6zEW-#}QSx@QzAX<~?wXYarcLkFFYyDThX-!F=_5XeXoBBZOcz%>5c zIRP>dQ|w{xnxQx)4Kc)78hB4{JzlO@3+C6!ty!wvGz>PUuq5KM44j&-Z9sZv0M`H8 z^$3+mkfksQkv=Ga5ee;|zkJC8-!7*SbXd9X8YOmtr?>K|IU+vQA1J&e%fYK4zup`O*D5Im?lVLX}(G0yFI z4eXw+u;O0C$2DUYakk)o%acOy=oDG^>LaGys&`?Wj)1{5StUh+MK=p^my9fUmhBCKXS6vOaR2zjC)i@$ zmpR*kSH7AG*le3(VKQ?Bdu|q{?3Z8 zD0rB|$}KwLmWEch8%zMyvLskIg5vW6(Wt0=sjMaHrrq3}&a(R`dUFu80aJnLcMMz$ zrPFZ*Ape!@D?Wal$FV>2Fnjnsci(D2x5aE~w7t{xmuBo(Qo1Qx%QR7Za~=})fBJCh zT$^z_`ggW}f63khUw}6N=!rmVS>AZCo6-*mvctywZ0rz(JBox%-r1k4N2erlA z{z;X)UYY40KH|M=$Ncf@9e1bj(r`eG44FZB#Ox-W4{r3CzGdhM0p*03ufXd28)|`0 zO*nhzp)KgyBYyI{7qm8118QM>g--#5oxf-0a8$#zLzfnktp4&(K)%jHBZLi>3#px% zE6kK?{@pp=B;MkdK`x?fk_Ki7$O;zyS?wK($YT3gR-v?s zWQ3}#23O!Yey49GGy)|xi`CaVziN9_-+g_o^b*0B;;25_7GSQ65cFL7Gy8()BQ~KH z2)5I&967vJUP{~)gp)$Ye+TlrLx(M4ycVyI)8kJv@symhq4m6@$*7~BXX3|OJa+&Y z50}FOpy&)8Xrya2QPFEICpp(lbMR8-Fq-|m?LCW$(1oK|00e{!g(mvgdRqRp97QgC ze7#iPTn-io3mRhQ#eH1jz||wsF}aCeZQ{yS)X3Di;kLx%h^_hhDneI#IeFUT4-AqT zgh3f|QWHB~=?6?$wOHvtKNuqCA>7^8YXv-LeH9&w0hmS6%GLo?m;khr<}s>CmlzEx zs?VhDXqcx-^K@1?y5E_$GlCOwDgGRv4tbJgJL_#Hv+!x?q)r#v4hoHhEaGm_!Klkk zpiBQY@ZO%J_PIlm@NWeeOO`yufGoKT0sY99v?r0{;hQsES8UPfWt;7qmq&-)N*snt zp4)Fj9!@9%lbY%QN3jXaml+HOV=f>K@E1?%EV!{$v`&Ud?O^pF4b0{Mg&(rH3~YYO z&`dRMDdZ2s$OW;CygYhbD1vO*wwuQj?E|-2`pVrLHlygJwk|nG1~);1W+ps@MU=;R z8_e-TLi^~gv4F@HdJd_DSkI+W`)Xi^LsX+5IC8nBxfs2>k&7ydPei4lyKR75*_Uft1;gH8MZ_DD2^TF{;+8;T(1t45 z0d9|x^XXC=o40#qufqDl8G4XoH*;fMhuq<|-v;M@oH^@Ykus!_rVQ7hnd`^z1LN3+ zog}B8O|?ea4;xl``1SrtdJs62X9uvuh)WcSg0Rr=_xaf#U4;f<`A;I&QV>fI+i$t@PLl z*%K5ou8Fhjb7wRiQapn0sec3iO{l4~B#cmO&+g|HkwjWG5rWlXWQUzQB<)K*JY_54 z)^G1_NpQ-#!kqO zaojL@qj0=`7pIL}zTjcOnU&0blf=HmJujcp?-*2((z^^nO40z1kKm6d6~ynil04+p zmMcAEik3@h6QKj{VcHOk8obET&KA|1XiIBua=rFZfAOQuQJ?Bvx=rrgRQI*1CGm`f zvp&S;XWd1L{YT^TL)JZ@Uu+U?2f>=Qr5!dL*y6R4?r+I3*eNU?sZ1lB-^R=@0p-db zZhvttd#kfw??Ai6-l$0WQ)u2NdoXx4K=i=Ha%tf|CgD}Q#(wYzVM-Tnjc+Ht` zwlvHFA7)7`-f15k%uR1menkqcOU_NZBn)F&7Ev(^0*gi5`vW#{xiq-5U9%E-)pM9{iQ z%g{RipHG;P1#>auTXOxfUGZ?Z`){~+{;Q(v6SD$vA8Da8QZW?4fbF=y^gZTC74do)4Y9< zS$n1j+y7?^Fqir&im0T7sOC&#x|1*7vBX=aWSlxVgrkG}#@s9YboE{%OIdyvHjCa} zw^8&oyXcSm=H0o+MO}Im)DirludHENNRpVcZlypdaiBZq8#2T_E?bpERTZOt($a>l z&X#4{psu3itw#9vtqypos_%~;mK_UNoqqfN+{K?L{_#fz7Qj#8)UQ)S!MN2edfVxd$=#oxHkR>o z^54$Urm5>)6?&UmmZ6iEWW#lTHGpadHO(PmC;WNzrUcXFkfPvN!i&U)+gv<3uXL>4 z?z6iUT@YQ@ushSXd~Q&S8ndr?T8`~)mz1z+{fFn%@2BWmY_|6=WfE+=rbN{{uM;dd zCA#K*7)|v@E6Z~9vj21n<1X(AC%m2W>I({Xp_^5WXQe`{IcHPcL(s@1s{WE&J7Ou+ z(kYkbvEXu}=S{G6;k)`V-o}&ngQ6K{TW{bEvYKj5LSAeDiGM1m7?gWp7hG ziX`wZmE@5QsB~yeNhVhAAQ~0wx8y8b7yWRB>rC6be4I;Ozvr%GF+rPjIU_)U_DFTH z@&(hAhOoE2%6C*ghYgCy)@5A6M{k*#fV-9<*{{?SFW8-w$NeMXwY9<*^bcuuJCrGH z`_W}qJMYp$4jyNv&}(UWXq(`(q@8UM3x5_tG7)t#(jt!9Gss6^=IiO2BfdWs`tSIB zUwd@GwhZHrIq>>yDxV_nj~huBDt1-KoVAPNjTiUMCqnL{!O(l#T#V}r;Z*Z)t+cP4 z3w@%h&d!q&)ggnG#le30g?p~=1PNY`z3bIxXyrw6B3?a>($KTxP~GmmsV4U%)0rif z&-^iNiP5)n$e`I(dnzcqI?A-kK>@eaGnQ&z+8gV$P$>P$=_KW3*b(clU{t2qt3Mqr z3roMxdsIzp*h}U-RV-AJ^kF^loGjWqu#lT}cU={_;jM%^t6LBAEOG>B)Wvi=-%?Q( z^|EUdwVdDkDY54yq=-Kwv7n@Eggh4m#NYP5%zDyP5o6=~@0-&yiCrnviC{fwOgWaC z|6$M9FO&}{RPYaT6h2EPOdvC=P$QHg#kbSJXP203%Xg?GtSc>zYwNPaMvw6J$!Qw~ zuwB2DV&0r^k}1-IBXG$7%Jq+F_4yvJ!^+sXMQWem4!Smww6vl zA}_A7aAKcEt8Et8ELh#sMAXb;IC^&ug6oagFQllZsTE3m;dZVYP2F@I!IwfhM{FnV{lK#;*h;zR882Cw}%Id-4_E|L{ z9`^Ut*Y+6?I*nh?{n2c$^f0^GP6Pa(%!Qd)Gn6Ka*lTPOMRX zfLLA5A+lqf{3NGQnsk1yw6Q*P>gQ?6~h?H1fR=% z>D2eLv{%-101q4cU{1f=albJ*8kXT_$WbFM-2{FXX;b3pw94Fd5J}a0#$IhFZVQU` z7{yARo!^V;+N#AT86-PLv&s+NK6z(G301ZZGsE5tV~HbxzdPB6TIP zUUZ>H|Gzq@z1cfKwt~2<6nqrHq7?OX?oQ|P_jrdv^fQx5`ZiL)Xk4ZxYFuzlSWr~_ z?l+~P=qnOF$!PpN?h~nJXjs@Y{&{5`IJbLMaGK>3K;;;tESbHlng7)#;5~$yVook{ zFAU@*+n8;=m0s4&un{(dD-$y&a~UW z%U5FFrjO?zD$KT9v36PBdX?YLsNbFZ%5#oiQ@(ig`pVDR{q0*$N$wmXZ`5VaGFfUY zf8{*W8_&|%<-sDmZM8gSrEdM=cf7y3oB~=|S-In;DIG#3*CpE1Ksc0|ap$!Y%&cGX zNXVS)+|Kcenz$kIC}F__-C=Y7Q~ZM@ zf@NI_R~L=U93S=@Q#poSo3g03TD4y+AROuP&*WAa6VzsX5{o*xeOmm7ifprz80rJ} zb%xt18V$(1jIqv93jCqXxQR5piM`0RS6F9(b;+JYTo+}teLwB#RoU82#{0W1Mh9m_ zf67lknYqy(Bvg7&S~wljS#yIO*=o@mu`xS*e%{9|v$HHoX{5^`S@|VB6q&ygr8n&| z;UjB(NtMufzx>kn>_AsHzE?(VLve;Jzy-#GWC?hm?7TiX;&6+{7jqlagHSyd?%{Fe z=)JiI-SJ0m5$YkcvOun8XpwDQVGKTcc_Z$EIkj|!-i4t{G5AkE9Vz#5iQ>*vPsEg8 zdcV2;E2W9BzrJGWhWL=OP2J<1Eu%GKR%P9KXMg=?Y1PQ#x1^_0&5XR&-(n3^=n*Q8 z{&Z|}{s?a2dl$r-lT<814svWB{jTs!;ut8)o_SctnXYw7-kF2PaRsPc&$83n8E7AN zU;J2>%Ihf~*;H0#o#bZz>iWloU#j+uj4F;#%r*T3@M$*j`I%u9Nn;9$E?IW{*|FuR zvq~i%^hL@p6JKMpU-Kj%VT@XpxhqT?uf6hkH!8xN=Za5Zh9Zz+6(TU%w_mwbYFk8S zmmjogtSaoLxfk|Gb$>*&8L`o$Y-j5QJ~?a7Rfa`om8n?79&v5{wHiRV^ODg{o0Nv* z(C}^xOmPjJ>eV>?VC~b+B5M&tk14N0;tuWL zHForR7hgwPoFq}wt%bkX%-y~{=Jrx)w0#-3?e>+BJLBQxa&zhV*mpnsZIHui!=c)` zle%i)F^jR^@9F8e*e~&X{Eapgbx^n2tHv)w+25aU-qQLvcxzXz9503JGk^(y+Ic1q zazILnY}wgZoF{9vlFT@eE9~_M&%A93g-w+QFeusOu%oeY#kvPOM%G>+EzB;KH4H`9*ND%b$AJjf#`^D5-(n~Ru*81UZ_pD~ z$VISbw*H_bl;=dxIb2!==w)WvX&bMzk6SC>FRmy=f768v$=Lw3J5Rn76 zhtzsTlU+qd{dXmEsu!}2tDx5QfL93;RK|3(M<3Fn)aiktth4nR_G030qD>xw4@`NR z&<$D$eh$o?3Q#KC$?1OPz?c!`I3QyeOA}mSkas0$2cmP&FGzb=h`m^OnghqZbiOag zXtc-(QlF<$`5L|n>HrIVX?^x*C1_LSB#j18?n2yFe%;mf zC)$_ur|$P`Bo~R@gKCn5io;go#5KyVuEwZ&PD1k41(hgn57jT(AMR1^y0~hB+S_{W z@VoJoAsFd~4*IeUh((ofMnsqoSA-jBOvWuPcfm80L30rViODj-qsQr7m{%GX&H86FW%bP?LHyFWv;gIf zYs+9U8xpHRdYcK^I%7 zX|*apyrEvV_Jt6`P&o?0m|@CT`EMPHhW>2V74?dDwYJ0Rvz)94{M8CY8D-ak@0Uxu zdk0*Z5l|%15C3^o%T**D<*UuT&*DF;d@qOE92MhYD|)7CYN2z8nnmbz%SGL!O%=4Y zURlOxe>UtANlMuDQZQ5bFEkM&d=X;|42UR&kv#SwlEUGCa3hBglI?uI44=vk)qnZ zQ>ad_V|kr#{+v9}e_qFQ6W9lvhEe?r3<)-dA3KTmHimCNQ7yX^!r7?tRCz66rgG1{ zW^X$sKiM=HCR}fIT`ez*EjeOza%xt!fsp2Y-DODgoFX=Tn-AiU^HZ0%CIx*ixN6se zDlpP!oJ9ScNx#SI%E|Y2;t{W38Jv7%6?~@cc%u8oJ(<{G9Oo)6#{G#1vmfV<<}Y#= z1DXV>BAMR>C=ovjz!gqZG|75?X`F#h9+6Lu@o{Gwo$^4j7a#pj3J%f@ zE3DDw0dUW<%Hd4OPM%9QKE0$uP^6Y7u7iJwg}G5)NeK`6d;-uu9pG~0TNo$*Z#7nT w|A)z^-rxlJN8kpsB-4cd|MT6`yBmjfmq<+C++^#?uln%=rTf|Up1%6O0H&Sas{jB1 literal 0 HcmV?d00001