|
| 1 | +####################################################################### |
| 2 | +# Copyright (c) 2019-present, Blosc Development Team <[email protected]> |
| 3 | +# All rights reserved. |
| 4 | +# |
| 5 | +# This source code is licensed under a BSD-style license (found in the |
| 6 | +# LICENSE file in the root directory of this source tree) |
| 7 | +####################################################################### |
| 8 | + |
| 9 | +# Compute expressions for different array sizes, using the jit decorator. |
| 10 | + |
| 11 | +from time import time |
| 12 | +import blosc2 |
| 13 | +import numpy as np |
| 14 | +import numexpr as ne |
| 15 | + |
| 16 | +niter = 5 |
| 17 | +# Create some data operands |
| 18 | +N = 10_000 # working size of ~1 GB |
| 19 | +dtype = "float32" |
| 20 | +chunks = (100, N) |
| 21 | +blocks = (1, N) |
| 22 | +chunks, blocks= None, None # enforce automatic chunk and block sizes |
| 23 | +cparams = blosc2.CParams(clevel=1, codec=blosc2.Codec.LZ4) |
| 24 | +cparams_out = blosc2.CParams(clevel=1, codec=blosc2.Codec.LZ4) |
| 25 | +print("Using cparams: ", cparams) |
| 26 | +check_result = True |
| 27 | +# Lossy compression |
| 28 | +# filters = [blosc2.Filter.TRUNC_PREC, blosc2.Filter.SHUFFLE] |
| 29 | +# filters_meta = [8, 0] # keep 8 bits of precision in mantissa |
| 30 | +# cparams = blosc2.CParams(clevel=1, codec=blosc2.Codec.LZ4, filters=filters, filters_meta=filters_meta) |
| 31 | +# check_result = False |
| 32 | + |
| 33 | + |
| 34 | +t0 = time() |
| 35 | +na = np.linspace(0, 1, N * N, dtype=dtype).reshape(N, N) |
| 36 | +nb = np.linspace(1, 2, N * N, dtype=dtype).reshape(N, N) |
| 37 | +nc = np.linspace(-10, 10, N, dtype=dtype) # broadcasting is supported |
| 38 | +# nc = np.linspace(-10, 10, N * N, dtype=dtype).reshape(N, N) |
| 39 | +print("Time to create data: ", time() - t0) |
| 40 | + |
| 41 | +def compute_expression_numpy(a, b, c): |
| 42 | + return ((a ** 3 + np.sin(a * 2)) < c) & (b > 0) |
| 43 | + |
| 44 | +t0 = time() |
| 45 | +nout = compute_expression_numpy(na, nb, nc) |
| 46 | +tref = time() - t0 |
| 47 | +print(f"Time to compute with NumPy engine: {tref:.5f}") |
| 48 | + |
| 49 | +nout = ne.evaluate("((na ** 3 + sin(na * 2)) < nc) & (nb > 0)") |
| 50 | +t0 = time() |
| 51 | +for i in range(niter): |
| 52 | + nout = ne.evaluate("((na ** 3 + sin(na * 2)) < nc) & (nb > 0)") |
| 53 | +t1 = (time() - t0) / niter |
| 54 | +print(f"Time to compute with NumExpr: {t1:.5f}") |
| 55 | +print(f"Speedup: {tref / t1:.2f}x") |
| 56 | + |
| 57 | +@blosc2.jit |
| 58 | +def compute_expression_nocompr(a, b, c): |
| 59 | + return ((a ** 3 + np.sin(a * 2)) < c) & (b > 0) |
| 60 | + |
| 61 | +print("\nUsing NumPy operands...") |
| 62 | + |
| 63 | +@blosc2.jit(cparams=cparams_out) |
| 64 | +def compute_expression_compr(a, b, c): |
| 65 | + return ((a ** 3 + np.sin(a * 2)) < c) & (b > 0) |
| 66 | + |
| 67 | +out = compute_expression_compr(na, nb, nc) |
| 68 | +t0 = time() |
| 69 | +for i in range(niter): |
| 70 | + out = compute_expression_compr(na, nb, nc) |
| 71 | +t1 = (time() - t0) / niter |
| 72 | +print(f"Time to compute with NumPy operands and NDArray as result: {t1:.5f}") |
| 73 | +cratio = out.schunk.cratio if isinstance(out, blosc2.NDArray) else 1.0 |
| 74 | +print(f"Speedup: {tref / t1:.2f}x, out cratio: {cratio:.2f}x") |
| 75 | +if check_result: |
| 76 | + np.testing.assert_allclose(out, nout) |
| 77 | + |
| 78 | +out = compute_expression_nocompr(na, nb, nc) |
| 79 | +t0 = time() |
| 80 | +for i in range(niter): |
| 81 | + out = compute_expression_nocompr(na, nb, nc) |
| 82 | +t1 = (time() - t0) / niter |
| 83 | +print(f"Time to compute with NumPy operands and NumPy as result: {t1:.5f}") |
| 84 | +cratio = out.schunk.cratio if isinstance(out, blosc2.NDArray) else 1.0 |
| 85 | +print(f"Speedup: {tref / t1:.2f}x, out cratio: {cratio:.2f}x") |
| 86 | +if check_result: |
| 87 | + np.testing.assert_allclose(out, nout) |
| 88 | + |
| 89 | +print("\nUsing NDArray operands *with* compression...") |
| 90 | +# Create Blosc2 operands |
| 91 | +a = blosc2.asarray(na, cparams=cparams, chunks=chunks, blocks=blocks) |
| 92 | +b = blosc2.asarray(nb, cparams=cparams, chunks=chunks, blocks=blocks) |
| 93 | +c = blosc2.asarray(nc, cparams=cparams) |
| 94 | +# c = blosc2.asarray(nc, cparams=cparams, chunks=chunks, blocks=blocks) |
| 95 | +print(f"{a.chunks=}, {a.blocks=}, {a.schunk.cratio=:.2f}x") |
| 96 | + |
| 97 | +out = compute_expression_compr(a, b, c) |
| 98 | +t0 = time() |
| 99 | +for i in range(niter): |
| 100 | + out = compute_expression_compr(a, b, c) |
| 101 | +t1 = (time() - t0) / niter |
| 102 | +print(f"[COMPR] Time to compute with NDArray operands and NDArray as result: {t1:.5f}") |
| 103 | +cratio = out.schunk.cratio if isinstance(out, blosc2.NDArray) else 1.0 |
| 104 | +print(f"Speedup: {tref / t1:.2f}x, out cratio: {cratio:.2f}x") |
| 105 | +if check_result: |
| 106 | + np.testing.assert_allclose(out, nout) |
| 107 | + |
| 108 | +out = compute_expression_nocompr(a, b, c) |
| 109 | +t0 = time() |
| 110 | +for i in range(niter): |
| 111 | + out = compute_expression_nocompr(a, b, c) |
| 112 | +t1 = (time() - t0) / niter |
| 113 | +print(f"[COMPR] Time to compute with NDArray operands and NumPy as result: {t1:.5f}") |
| 114 | +cratio = out.schunk.cratio if isinstance(out, blosc2.NDArray) else 1.0 |
| 115 | +print(f"Speedup: {tref / t1:.2f}x, out cratio: {cratio:.2f}x") |
| 116 | +if check_result: |
| 117 | + np.testing.assert_allclose(out, nout) |
| 118 | + |
| 119 | +print("\nUsing NDArray operands without compression...") |
| 120 | +# Create NDArray operands without compression |
| 121 | +cparams = cparams_out = blosc2.CParams(clevel=0) |
| 122 | +a = blosc2.asarray(na, cparams=cparams, chunks=chunks, blocks=blocks) |
| 123 | +b = blosc2.asarray(nb, cparams=cparams, chunks=chunks, blocks=blocks) |
| 124 | +c = blosc2.asarray(nc, cparams=cparams) |
| 125 | +# c = blosc2.asarray(nc, cparams=cparams, chunks=chunks, blocks=blocks) |
| 126 | +print(f"{a.chunks=}, {a.blocks=}, {a.schunk.cratio=:.2f}x") |
| 127 | + |
| 128 | +@blosc2.jit(cparams=cparams_out) |
| 129 | +def compute_expression_compr(a, b, c): |
| 130 | + return ((a ** 3 + np.sin(a * 2)) < c) & (b > 0) |
| 131 | + |
| 132 | +out = compute_expression_compr(a, b, c) |
| 133 | +t0 = time() |
| 134 | +for i in range(niter): |
| 135 | + out = compute_expression_compr(a, b, c) |
| 136 | +t1 = (time() - t0) / niter |
| 137 | +print(f"[NOCOMPR] Time to compute with NDArray operands and NDArray as result: {t1:.5f}") |
| 138 | +cratio = out.schunk.cratio if isinstance(out, blosc2.NDArray) else 1.0 |
| 139 | +print(f"Speedup: {tref / t1:.2f}x, out cratio: {cratio:.2f}x") |
| 140 | +if check_result: |
| 141 | + np.testing.assert_allclose(out, nout) |
| 142 | + |
| 143 | +out = compute_expression_nocompr(a, b, c) |
| 144 | +t0 = time() |
| 145 | +for i in range(niter): |
| 146 | + out = compute_expression_nocompr(a, b, c) |
| 147 | +t1 = (time() - t0) / niter |
| 148 | +print(f"[NOCOMPR] Time to compute with NDArray operands and NumPy as result: {t1:.5f}") |
| 149 | +cratio = out.schunk.cratio if isinstance(out, blosc2.NDArray) else 1.0 |
| 150 | +print(f"Speedup: {tref / t1:.2f}x, out cratio: {cratio:.2f}x") |
| 151 | +if check_result: |
| 152 | + np.testing.assert_allclose(out, nout) |
| 153 | + |
| 154 | +print("All results are equal!") |
0 commit comments