diff --git a/README.md b/README.md index 1400f8888..597c16e9f 100644 --- a/README.md +++ b/README.md @@ -58,7 +58,11 @@ python -m pip install dace ``` However, you may want to install the latest version from the [GitHub repository](https://github.com/spcl/dace). To run NPBench with DaCe, you have to select as framework (see details below) -either `dace_cpu` or `dace_gpu`. +either `dace_cpu` or `dace_gpu` (which build the SDFG and `auto_optimize` it), or +`dace_canonicalize_cpu` / `dace_canonicalize_gpu` (which instead run DaCe's +`canonicalize` pipeline, leaving `sdfg.openmp_array_reductions` on so whole-buffer +WCR accumulators lower to OpenMP array-section reductions rather than per-element +atomics). ### DPNP diff --git a/framework_info/dace_canonicalize_cpu.json b/framework_info/dace_canonicalize_cpu.json new file mode 100644 index 000000000..a90fceb04 --- /dev/null +++ b/framework_info/dace_canonicalize_cpu.json @@ -0,0 +1,10 @@ +{ + "framework": { + "simple_name": "dace_canonicalize_cpu", + "full_name": "DaCe Canonicalize CPU", + "prefix": "dc", + "postfix": "dace", + "class": "DaceCanonicalizeFramework", + "arch": "cpu" + } +} diff --git a/framework_info/dace_canonicalize_gpu.json b/framework_info/dace_canonicalize_gpu.json new file mode 100644 index 000000000..af2e998ab --- /dev/null +++ b/framework_info/dace_canonicalize_gpu.json @@ -0,0 +1,10 @@ +{ + "framework": { + "simple_name": "dace_canonicalize_gpu", + "full_name": "DaCe Canonicalize GPU", + "prefix": "dc", + "postfix": "dace", + "class": "DaceCanonicalizeFramework", + "arch": "gpu" + } +} diff --git a/npbench.db b/npbench.db new file mode 100644 index 000000000..3e58a3653 Binary files /dev/null and b/npbench.db differ diff --git a/npbench/benchmarks/cavity_flow/cavity_flow_cupy.py b/npbench/benchmarks/cavity_flow/cavity_flow_cupy.py index dc2cdd7e2..3ad5ec48f 100644 --- a/npbench/benchmarks/cavity_flow/cavity_flow_cupy.py +++ b/npbench/benchmarks/cavity_flow/cavity_flow_cupy.py @@ -22,7 +22,7 @@ def build_up_b(b, rho, dt, u, v, dx, dy): def pressure_poisson(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) pn = p.copy() for q in range(nit): @@ -39,8 +39,8 @@ def pressure_poisson(nit, p, dx, dy, b): def cavity_flow(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu): - un = np.empty_like(u) - vn = np.empty_like(v) + un = np.zeros_like(u) + vn = np.zeros_like(v) b = np.zeros((ny, nx)) for n in range(nt): diff --git a/npbench/benchmarks/cavity_flow/cavity_flow_dace.py b/npbench/benchmarks/cavity_flow/cavity_flow_dace.py index 6ff03beb8..ce11d5e1d 100644 --- a/npbench/benchmarks/cavity_flow/cavity_flow_dace.py +++ b/npbench/benchmarks/cavity_flow/cavity_flow_dace.py @@ -30,7 +30,7 @@ def build_up_b(b: dc.float64[ny, nx], rho: dc.float64, dt: dc.float64, @dc.program def pressure_poisson(p: dc.float64[ny, nx], dx: dc.float64, dy: dc.float64, b: dc.float64[ny, nx]): - pn = np.empty_like(p) + pn = np.zeros_like(p) pn[:] = p.copy() for q in range(nit): @@ -51,8 +51,8 @@ def cavity_flow(nt: dc.int64, nit: dc.int64, u: dc.float64[ny, nx], v: dc.float64[ny, nx], dt: dc.float64, dx: dc.float64, dy: dc.float64, p: dc.float64[ny, nx], rho: dc.float64, nu: dc.float64): - un = np.empty_like(u) - vn = np.empty_like(v) + un = np.zeros_like(u) + vn = np.zeros_like(v) b = np.zeros((ny, nx)) for n in range(nt): diff --git a/npbench/benchmarks/cavity_flow/cavity_flow_dpnp.py b/npbench/benchmarks/cavity_flow/cavity_flow_dpnp.py index 70cb508b8..58753fd13 100644 --- a/npbench/benchmarks/cavity_flow/cavity_flow_dpnp.py +++ b/npbench/benchmarks/cavity_flow/cavity_flow_dpnp.py @@ -1,74 +1,74 @@ -import dpnp as np - -def build_up_b(b, rho, dt, u, v, dx, dy): - b[1:-1, 1:-1] = ( - rho * ( - (1 / dt * ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx) + - (v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy))) - - ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx))**2 - - 2 * ((u[2:, 1:-1] - u[0:-2, 1:-1]) / (2 * dy) * - (v[1:-1, 2:] - v[1:-1, 0:-2]) / (2 * dx)) - - ((v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy))**2 - ) - ) - -def pressure_poisson(nit, p, dx, dy, b): - pn = np.empty_like(p) - for q in range(nit): - pn[:] = p.copy() - p[1:-1, 1:-1] = ( - ((pn[1:-1, 2:] + pn[1:-1, 0:-2]) * dy**2 + - (pn[2:, 1:-1] + pn[0:-2, 1:-1]) * dx**2) / - (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / - (2 * (dx**2 + dy**2)) * b[1:-1, 1:-1] - ) - - p[:, -1] = p[:, -2] # dp/dx = 0 at x = 2 - p[0, :] = p[1, :] # dp/dy = 0 at y = 0 - p[:, 0] = p[:, 1] # dp/dx = 0 at x = 0 - p[-1, :] = 0 # p = 0 at y = 2 - -def cavity_flow(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu): - un = np.empty_like(u) - vn = np.empty_like(v) - b = np.zeros((ny, nx)) - - for n in range(nt): - un[:] = u.copy() - vn[:] = v.copy() - - build_up_b(b, rho, dt, u, v, dx, dy) - pressure_poisson(nit, p, dx, dy, b) - - u[1:-1, 1:-1] = ( - un[1:-1, 1:-1] - un[1:-1, 1:-1] * dt / dx * - (un[1:-1, 1:-1] - un[1:-1, 0:-2]) - - vn[1:-1, 1:-1] * dt / dy * - (un[1:-1, 1:-1] - un[0:-2, 1:-1]) - dt / (2 * rho * dx) * - (p[1:-1, 2:] - p[1:-1, 0:-2]) + nu * - (dt / dx**2 * - (un[1:-1, 2:] - 2 * un[1:-1, 1:-1] + un[1:-1, 0:-2]) + - dt / dy**2 * - (un[2:, 1:-1] - 2 * un[1:-1, 1:-1] + un[0:-2, 1:-1])) - ) - - v[1:-1, 1:-1] = ( - vn[1:-1, 1:-1] - un[1:-1, 1:-1] * dt / dx * - (vn[1:-1, 1:-1] - vn[1:-1, 0:-2]) - - vn[1:-1, 1:-1] * dt / dy * - (vn[1:-1, 1:-1] - vn[0:-2, 1:-1]) - dt / (2 * rho * dy) * - (p[2:, 1:-1] - p[0:-2, 1:-1]) + nu * - (dt / dx**2 * - (vn[1:-1, 2:] - 2 * vn[1:-1, 1:-1] + vn[1:-1, 0:-2]) + - dt / dy**2 * - (vn[2:, 1:-1] - 2 * vn[1:-1, 1:-1] + vn[0:-2, 1:-1])) - ) - - u[0, :] = 0 - u[:, 0] = 0 - u[:, -1] = 0 - u[-1, :] = 1 # set velocity on cavity lid equal to 1 - v[0, :] = 0 - v[-1, :] = 0 - v[:, 0] = 0 - v[:, -1] = 0 +import dpnp as np + +def build_up_b(b, rho, dt, u, v, dx, dy): + b[1:-1, 1:-1] = ( + rho * ( + (1 / dt * ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx) + + (v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy))) - + ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx))**2 - + 2 * ((u[2:, 1:-1] - u[0:-2, 1:-1]) / (2 * dy) * + (v[1:-1, 2:] - v[1:-1, 0:-2]) / (2 * dx)) - + ((v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy))**2 + ) + ) + +def pressure_poisson(nit, p, dx, dy, b): + pn = np.zeros_like(p) + for q in range(nit): + pn[:] = p.copy() + p[1:-1, 1:-1] = ( + ((pn[1:-1, 2:] + pn[1:-1, 0:-2]) * dy**2 + + (pn[2:, 1:-1] + pn[0:-2, 1:-1]) * dx**2) / + (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / + (2 * (dx**2 + dy**2)) * b[1:-1, 1:-1] + ) + + p[:, -1] = p[:, -2] # dp/dx = 0 at x = 2 + p[0, :] = p[1, :] # dp/dy = 0 at y = 0 + p[:, 0] = p[:, 1] # dp/dx = 0 at x = 0 + p[-1, :] = 0 # p = 0 at y = 2 + +def cavity_flow(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu): + un = np.zeros_like(u) + vn = np.zeros_like(v) + b = np.zeros((ny, nx)) + + for n in range(nt): + un[:] = u.copy() + vn[:] = v.copy() + + build_up_b(b, rho, dt, u, v, dx, dy) + pressure_poisson(nit, p, dx, dy, b) + + u[1:-1, 1:-1] = ( + un[1:-1, 1:-1] - un[1:-1, 1:-1] * dt / dx * + (un[1:-1, 1:-1] - un[1:-1, 0:-2]) - + vn[1:-1, 1:-1] * dt / dy * + (un[1:-1, 1:-1] - un[0:-2, 1:-1]) - dt / (2 * rho * dx) * + (p[1:-1, 2:] - p[1:-1, 0:-2]) + nu * + (dt / dx**2 * + (un[1:-1, 2:] - 2 * un[1:-1, 1:-1] + un[1:-1, 0:-2]) + + dt / dy**2 * + (un[2:, 1:-1] - 2 * un[1:-1, 1:-1] + un[0:-2, 1:-1])) + ) + + v[1:-1, 1:-1] = ( + vn[1:-1, 1:-1] - un[1:-1, 1:-1] * dt / dx * + (vn[1:-1, 1:-1] - vn[1:-1, 0:-2]) - + vn[1:-1, 1:-1] * dt / dy * + (vn[1:-1, 1:-1] - vn[0:-2, 1:-1]) - dt / (2 * rho * dy) * + (p[2:, 1:-1] - p[0:-2, 1:-1]) + nu * + (dt / dx**2 * + (vn[1:-1, 2:] - 2 * vn[1:-1, 1:-1] + vn[1:-1, 0:-2]) + + dt / dy**2 * + (vn[2:, 1:-1] - 2 * vn[1:-1, 1:-1] + vn[0:-2, 1:-1])) + ) + + u[0, :] = 0 + u[:, 0] = 0 + u[:, -1] = 0 + u[-1, :] = 1 # set velocity on cavity lid equal to 1 + v[0, :] = 0 + v[-1, :] = 0 + v[:, 0] = 0 + v[:, -1] = 0 diff --git a/npbench/benchmarks/cavity_flow/cavity_flow_legate.py b/npbench/benchmarks/cavity_flow/cavity_flow_legate.py index 31860bb3a..7acb865a9 100644 --- a/npbench/benchmarks/cavity_flow/cavity_flow_legate.py +++ b/npbench/benchmarks/cavity_flow/cavity_flow_legate.py @@ -22,7 +22,7 @@ def build_up_b(b, rho, dt, u, v, dx, dy): def pressure_poisson(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) pn = p.copy() for q in range(nit): @@ -39,8 +39,8 @@ def pressure_poisson(nit, p, dx, dy, b): def cavity_flow(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu): - un = np.empty_like(u) - vn = np.empty_like(v) + un = np.zeros_like(u) + vn = np.zeros_like(v) b = np.zeros((ny, nx)) for n in range(nt): diff --git a/npbench/benchmarks/cavity_flow/cavity_flow_numba_n.py b/npbench/benchmarks/cavity_flow/cavity_flow_numba_n.py index ecd1a7eee..3255d7cf3 100644 --- a/npbench/benchmarks/cavity_flow/cavity_flow_numba_n.py +++ b/npbench/benchmarks/cavity_flow/cavity_flow_numba_n.py @@ -16,7 +16,7 @@ def build_up_b(b, rho, dt, u, v, dx, dy): @nb.jit(nopython=True, parallel=False, fastmath=True) def pressure_poisson(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) pn = p.copy() for q in range(nit): @@ -34,8 +34,8 @@ def pressure_poisson(nit, p, dx, dy, b): @nb.jit(nopython=True, parallel=False, fastmath=True) def cavity_flow(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu): - un = np.empty_like(u) - vn = np.empty_like(v) + un = np.zeros_like(u) + vn = np.zeros_like(v) b = np.zeros((ny, nx)) for n in range(nt): diff --git a/npbench/benchmarks/cavity_flow/cavity_flow_numba_np.py b/npbench/benchmarks/cavity_flow/cavity_flow_numba_np.py index 477f017c6..e8a3cce77 100644 --- a/npbench/benchmarks/cavity_flow/cavity_flow_numba_np.py +++ b/npbench/benchmarks/cavity_flow/cavity_flow_numba_np.py @@ -16,7 +16,7 @@ def build_up_b(b, rho, dt, u, v, dx, dy): @nb.jit(nopython=True, parallel=True, fastmath=True) def pressure_poisson(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) pn = p.copy() for q in range(nit): @@ -33,9 +33,9 @@ def pressure_poisson(nit, p, dx, dy, b): @nb.jit(nopython=True, parallel=True, fastmath=True) -def nopython_mode(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu): - un = np.empty_like(u) - vn = np.empty_like(v) +def cavity_flow(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu): + un = np.zeros_like(u) + vn = np.zeros_like(v) b = np.zeros((ny, nx)) for n in range(nt): diff --git a/npbench/benchmarks/cavity_flow/cavity_flow_numpy.py b/npbench/benchmarks/cavity_flow/cavity_flow_numpy.py index 9890acb39..23235d8d6 100644 --- a/npbench/benchmarks/cavity_flow/cavity_flow_numpy.py +++ b/npbench/benchmarks/cavity_flow/cavity_flow_numpy.py @@ -22,7 +22,7 @@ def build_up_b(b, rho, dt, u, v, dx, dy): def pressure_poisson(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) pn = p.copy() for q in range(nit): @@ -39,8 +39,8 @@ def pressure_poisson(nit, p, dx, dy, b): def cavity_flow(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu): - un = np.empty_like(u) - vn = np.empty_like(v) + un = np.zeros_like(u) + vn = np.zeros_like(v) b = np.zeros((ny, nx)) for n in range(nt): diff --git a/npbench/benchmarks/cavity_flow/cavity_flow_pythran.py b/npbench/benchmarks/cavity_flow/cavity_flow_pythran.py index 78c5008ad..5eb983e33 100644 --- a/npbench/benchmarks/cavity_flow/cavity_flow_pythran.py +++ b/npbench/benchmarks/cavity_flow/cavity_flow_pythran.py @@ -26,7 +26,7 @@ def build_up_b(b, rho, dt, u, v, dx, dy): # pythran export pressure_poisson(int64, float64[:,:], float64, float64, # float64[:,:]) def pressure_poisson(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) pn = p.copy() for q in range(nit): @@ -46,8 +46,8 @@ def pressure_poisson(nit, p, dx, dy, b): # float64[:,:], float64, float64, float64, # float64[:,:], float64, float64) def cavity_flow(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu): - un = np.empty_like(u) - vn = np.empty_like(v) + un = np.zeros_like(u) + vn = np.zeros_like(v) b = np.zeros((ny, nx)) for n in range(nt): diff --git a/npbench/benchmarks/channel_flow/channel_flow_cupy.py b/npbench/benchmarks/channel_flow/channel_flow_cupy.py index 4cbef94dd..5142982b1 100644 --- a/npbench/benchmarks/channel_flow/channel_flow_cupy.py +++ b/npbench/benchmarks/channel_flow/channel_flow_cupy.py @@ -40,7 +40,7 @@ def build_up_b(rho, dt, dx, dy, u, v): def pressure_poisson_periodic(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) for q in range(nit): pn = p.copy() diff --git a/npbench/benchmarks/channel_flow/channel_flow_dace.py b/npbench/benchmarks/channel_flow/channel_flow_dace.py index 7654d7c01..c519332d7 100644 --- a/npbench/benchmarks/channel_flow/channel_flow_dace.py +++ b/npbench/benchmarks/channel_flow/channel_flow_dace.py @@ -47,7 +47,7 @@ def build_up_b(rho: dc.float64, dt: dc.float64, dx: dc.float64, dy: dc.float64, @dc.program def pressure_poisson_periodic(p: dc.float64[ny, nx], dx: dc.float64, dy: dc.float64, b: dc.float64[ny, nx]): - pn = np.empty_like(p) + pn = np.zeros_like(p) for q in range(nit): pn[:] = p.copy() diff --git a/npbench/benchmarks/channel_flow/channel_flow_dpnp.py b/npbench/benchmarks/channel_flow/channel_flow_dpnp.py index 4e5f12052..d128ce43a 100644 --- a/npbench/benchmarks/channel_flow/channel_flow_dpnp.py +++ b/npbench/benchmarks/channel_flow/channel_flow_dpnp.py @@ -1,142 +1,142 @@ -import dpnp as np - -def build_up_b(rho, dt, dx, dy, u, v): - b = np.zeros_like(u) - b[1:-1, - 1:-1] = (rho * (1 / dt * ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx) + - (v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy)) - - ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx))**2 - 2 * - ((u[2:, 1:-1] - u[0:-2, 1:-1]) / (2 * dy) * - (v[1:-1, 2:] - v[1:-1, 0:-2]) / (2 * dx)) - - ((v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy))**2)) - - # Periodic BC Pressure @ x = 2 - b[1:-1, -1] = (rho * (1 / dt * ((u[1:-1, 0] - u[1:-1, -2]) / (2 * dx) + - (v[2:, -1] - v[0:-2, -1]) / (2 * dy)) - - ((u[1:-1, 0] - u[1:-1, -2]) / (2 * dx))**2 - 2 * - ((u[2:, -1] - u[0:-2, -1]) / (2 * dy) * - (v[1:-1, 0] - v[1:-1, -2]) / (2 * dx)) - - ((v[2:, -1] - v[0:-2, -1]) / (2 * dy))**2)) - - # Periodic BC Pressure @ x = 0 - b[1:-1, 0] = (rho * (1 / dt * ((u[1:-1, 1] - u[1:-1, -1]) / (2 * dx) + - (v[2:, 0] - v[0:-2, 0]) / (2 * dy)) - - ((u[1:-1, 1] - u[1:-1, -1]) / (2 * dx))**2 - 2 * - ((u[2:, 0] - u[0:-2, 0]) / (2 * dy) * - (v[1:-1, 1] - v[1:-1, -1]) / - (2 * dx)) - ((v[2:, 0] - v[0:-2, 0]) / (2 * dy))**2)) - - return b - - -def pressure_poisson_periodic(nit, p, dx, dy, b): - pn = np.empty_like(p) - - for q in range(nit): - pn = p.copy() - p[1:-1, 1:-1] = (((pn[1:-1, 2:] + pn[1:-1, 0:-2]) * dy**2 + - (pn[2:, 1:-1] + pn[0:-2, 1:-1]) * dx**2) / - (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / - (2 * (dx**2 + dy**2)) * b[1:-1, 1:-1]) - - # Periodic BC Pressure @ x = 2 - p[1:-1, -1] = (((pn[1:-1, 0] + pn[1:-1, -2]) * dy**2 + - (pn[2:, -1] + pn[0:-2, -1]) * dx**2) / - (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / - (2 * (dx**2 + dy**2)) * b[1:-1, -1]) - - # Periodic BC Pressure @ x = 0 - p[1:-1, - 0] = (((pn[1:-1, 1] + pn[1:-1, -1]) * dy**2 + - (pn[2:, 0] + pn[0:-2, 0]) * dx**2) / (2 * (dx**2 + dy**2)) - - dx**2 * dy**2 / (2 * (dx**2 + dy**2)) * b[1:-1, 0]) - - # Wall boundary conditions, pressure - p[-1, :] = p[-2, :] # dp/dy = 0 at y = 2 - p[0, :] = p[1, :] # dp/dy = 0 at y = 0 - - -def channel_flow(nit, u, v, dt, dx, dy, p, rho, nu, F): - udiff = 1 - stepcount = 0 - - while udiff > .001: - un = u.copy() - vn = v.copy() - - b = build_up_b(rho, dt, dx, dy, u, v) - pressure_poisson_periodic(nit, p, dx, dy, b) - - u[1:-1, - 1:-1] = (un[1:-1, 1:-1] - un[1:-1, 1:-1] * dt / dx * - (un[1:-1, 1:-1] - un[1:-1, 0:-2]) - - vn[1:-1, 1:-1] * dt / dy * - (un[1:-1, 1:-1] - un[0:-2, 1:-1]) - dt / (2 * rho * dx) * - (p[1:-1, 2:] - p[1:-1, 0:-2]) + nu * - (dt / dx**2 * - (un[1:-1, 2:] - 2 * un[1:-1, 1:-1] + un[1:-1, 0:-2]) + - dt / dy**2 * - (un[2:, 1:-1] - 2 * un[1:-1, 1:-1] + un[0:-2, 1:-1])) + - F * dt) - - v[1:-1, - 1:-1] = (vn[1:-1, 1:-1] - un[1:-1, 1:-1] * dt / dx * - (vn[1:-1, 1:-1] - vn[1:-1, 0:-2]) - - vn[1:-1, 1:-1] * dt / dy * - (vn[1:-1, 1:-1] - vn[0:-2, 1:-1]) - dt / (2 * rho * dy) * - (p[2:, 1:-1] - p[0:-2, 1:-1]) + nu * - (dt / dx**2 * - (vn[1:-1, 2:] - 2 * vn[1:-1, 1:-1] + vn[1:-1, 0:-2]) + - dt / dy**2 * - (vn[2:, 1:-1] - 2 * vn[1:-1, 1:-1] + vn[0:-2, 1:-1]))) - - # Periodic BC u @ x = 2 - u[1:-1, -1] = ( - un[1:-1, -1] - un[1:-1, -1] * dt / dx * - (un[1:-1, -1] - un[1:-1, -2]) - vn[1:-1, -1] * dt / dy * - (un[1:-1, -1] - un[0:-2, -1]) - dt / (2 * rho * dx) * - (p[1:-1, 0] - p[1:-1, -2]) + nu * - (dt / dx**2 * - (un[1:-1, 0] - 2 * un[1:-1, -1] + un[1:-1, -2]) + dt / dy**2 * - (un[2:, -1] - 2 * un[1:-1, -1] + un[0:-2, -1])) + F * dt) - - # Periodic BC u @ x = 0 - u[1:-1, - 0] = (un[1:-1, 0] - un[1:-1, 0] * dt / dx * - (un[1:-1, 0] - un[1:-1, -1]) - vn[1:-1, 0] * dt / dy * - (un[1:-1, 0] - un[0:-2, 0]) - dt / (2 * rho * dx) * - (p[1:-1, 1] - p[1:-1, -1]) + nu * - (dt / dx**2 * - (un[1:-1, 1] - 2 * un[1:-1, 0] + un[1:-1, -1]) + dt / dy**2 * - (un[2:, 0] - 2 * un[1:-1, 0] + un[0:-2, 0])) + F * dt) - - # Periodic BC v @ x = 2 - v[1:-1, -1] = ( - vn[1:-1, -1] - un[1:-1, -1] * dt / dx * - (vn[1:-1, -1] - vn[1:-1, -2]) - vn[1:-1, -1] * dt / dy * - (vn[1:-1, -1] - vn[0:-2, -1]) - dt / (2 * rho * dy) * - (p[2:, -1] - p[0:-2, -1]) + nu * - (dt / dx**2 * - (vn[1:-1, 0] - 2 * vn[1:-1, -1] + vn[1:-1, -2]) + dt / dy**2 * - (vn[2:, -1] - 2 * vn[1:-1, -1] + vn[0:-2, -1]))) - - # Periodic BC v @ x = 0 - v[1:-1, - 0] = (vn[1:-1, 0] - un[1:-1, 0] * dt / dx * - (vn[1:-1, 0] - vn[1:-1, -1]) - vn[1:-1, 0] * dt / dy * - (vn[1:-1, 0] - vn[0:-2, 0]) - dt / (2 * rho * dy) * - (p[2:, 0] - p[0:-2, 0]) + nu * - (dt / dx**2 * - (vn[1:-1, 1] - 2 * vn[1:-1, 0] + vn[1:-1, -1]) + dt / dy**2 * - (vn[2:, 0] - 2 * vn[1:-1, 0] + vn[0:-2, 0]))) - - # Wall BC: u,v = 0 @ y = 0,2 - u[0, :] = 0 - u[-1, :] = 0 - v[0, :] = 0 - v[-1, :] = 0 - - udiff = (np.sum(u) - np.sum(un)) / np.sum(u) - stepcount += 1 - - return stepcount +import dpnp as np + +def build_up_b(rho, dt, dx, dy, u, v): + b = np.zeros_like(u) + b[1:-1, + 1:-1] = (rho * (1 / dt * ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx) + + (v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy)) - + ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx))**2 - 2 * + ((u[2:, 1:-1] - u[0:-2, 1:-1]) / (2 * dy) * + (v[1:-1, 2:] - v[1:-1, 0:-2]) / (2 * dx)) - + ((v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy))**2)) + + # Periodic BC Pressure @ x = 2 + b[1:-1, -1] = (rho * (1 / dt * ((u[1:-1, 0] - u[1:-1, -2]) / (2 * dx) + + (v[2:, -1] - v[0:-2, -1]) / (2 * dy)) - + ((u[1:-1, 0] - u[1:-1, -2]) / (2 * dx))**2 - 2 * + ((u[2:, -1] - u[0:-2, -1]) / (2 * dy) * + (v[1:-1, 0] - v[1:-1, -2]) / (2 * dx)) - + ((v[2:, -1] - v[0:-2, -1]) / (2 * dy))**2)) + + # Periodic BC Pressure @ x = 0 + b[1:-1, 0] = (rho * (1 / dt * ((u[1:-1, 1] - u[1:-1, -1]) / (2 * dx) + + (v[2:, 0] - v[0:-2, 0]) / (2 * dy)) - + ((u[1:-1, 1] - u[1:-1, -1]) / (2 * dx))**2 - 2 * + ((u[2:, 0] - u[0:-2, 0]) / (2 * dy) * + (v[1:-1, 1] - v[1:-1, -1]) / + (2 * dx)) - ((v[2:, 0] - v[0:-2, 0]) / (2 * dy))**2)) + + return b + + +def pressure_poisson_periodic(nit, p, dx, dy, b): + pn = np.zeros_like(p) + + for q in range(nit): + pn = p.copy() + p[1:-1, 1:-1] = (((pn[1:-1, 2:] + pn[1:-1, 0:-2]) * dy**2 + + (pn[2:, 1:-1] + pn[0:-2, 1:-1]) * dx**2) / + (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / + (2 * (dx**2 + dy**2)) * b[1:-1, 1:-1]) + + # Periodic BC Pressure @ x = 2 + p[1:-1, -1] = (((pn[1:-1, 0] + pn[1:-1, -2]) * dy**2 + + (pn[2:, -1] + pn[0:-2, -1]) * dx**2) / + (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / + (2 * (dx**2 + dy**2)) * b[1:-1, -1]) + + # Periodic BC Pressure @ x = 0 + p[1:-1, + 0] = (((pn[1:-1, 1] + pn[1:-1, -1]) * dy**2 + + (pn[2:, 0] + pn[0:-2, 0]) * dx**2) / (2 * (dx**2 + dy**2)) - + dx**2 * dy**2 / (2 * (dx**2 + dy**2)) * b[1:-1, 0]) + + # Wall boundary conditions, pressure + p[-1, :] = p[-2, :] # dp/dy = 0 at y = 2 + p[0, :] = p[1, :] # dp/dy = 0 at y = 0 + + +def channel_flow(nit, u, v, dt, dx, dy, p, rho, nu, F): + udiff = 1 + stepcount = 0 + + while udiff > .001: + un = u.copy() + vn = v.copy() + + b = build_up_b(rho, dt, dx, dy, u, v) + pressure_poisson_periodic(nit, p, dx, dy, b) + + u[1:-1, + 1:-1] = (un[1:-1, 1:-1] - un[1:-1, 1:-1] * dt / dx * + (un[1:-1, 1:-1] - un[1:-1, 0:-2]) - + vn[1:-1, 1:-1] * dt / dy * + (un[1:-1, 1:-1] - un[0:-2, 1:-1]) - dt / (2 * rho * dx) * + (p[1:-1, 2:] - p[1:-1, 0:-2]) + nu * + (dt / dx**2 * + (un[1:-1, 2:] - 2 * un[1:-1, 1:-1] + un[1:-1, 0:-2]) + + dt / dy**2 * + (un[2:, 1:-1] - 2 * un[1:-1, 1:-1] + un[0:-2, 1:-1])) + + F * dt) + + v[1:-1, + 1:-1] = (vn[1:-1, 1:-1] - un[1:-1, 1:-1] * dt / dx * + (vn[1:-1, 1:-1] - vn[1:-1, 0:-2]) - + vn[1:-1, 1:-1] * dt / dy * + (vn[1:-1, 1:-1] - vn[0:-2, 1:-1]) - dt / (2 * rho * dy) * + (p[2:, 1:-1] - p[0:-2, 1:-1]) + nu * + (dt / dx**2 * + (vn[1:-1, 2:] - 2 * vn[1:-1, 1:-1] + vn[1:-1, 0:-2]) + + dt / dy**2 * + (vn[2:, 1:-1] - 2 * vn[1:-1, 1:-1] + vn[0:-2, 1:-1]))) + + # Periodic BC u @ x = 2 + u[1:-1, -1] = ( + un[1:-1, -1] - un[1:-1, -1] * dt / dx * + (un[1:-1, -1] - un[1:-1, -2]) - vn[1:-1, -1] * dt / dy * + (un[1:-1, -1] - un[0:-2, -1]) - dt / (2 * rho * dx) * + (p[1:-1, 0] - p[1:-1, -2]) + nu * + (dt / dx**2 * + (un[1:-1, 0] - 2 * un[1:-1, -1] + un[1:-1, -2]) + dt / dy**2 * + (un[2:, -1] - 2 * un[1:-1, -1] + un[0:-2, -1])) + F * dt) + + # Periodic BC u @ x = 0 + u[1:-1, + 0] = (un[1:-1, 0] - un[1:-1, 0] * dt / dx * + (un[1:-1, 0] - un[1:-1, -1]) - vn[1:-1, 0] * dt / dy * + (un[1:-1, 0] - un[0:-2, 0]) - dt / (2 * rho * dx) * + (p[1:-1, 1] - p[1:-1, -1]) + nu * + (dt / dx**2 * + (un[1:-1, 1] - 2 * un[1:-1, 0] + un[1:-1, -1]) + dt / dy**2 * + (un[2:, 0] - 2 * un[1:-1, 0] + un[0:-2, 0])) + F * dt) + + # Periodic BC v @ x = 2 + v[1:-1, -1] = ( + vn[1:-1, -1] - un[1:-1, -1] * dt / dx * + (vn[1:-1, -1] - vn[1:-1, -2]) - vn[1:-1, -1] * dt / dy * + (vn[1:-1, -1] - vn[0:-2, -1]) - dt / (2 * rho * dy) * + (p[2:, -1] - p[0:-2, -1]) + nu * + (dt / dx**2 * + (vn[1:-1, 0] - 2 * vn[1:-1, -1] + vn[1:-1, -2]) + dt / dy**2 * + (vn[2:, -1] - 2 * vn[1:-1, -1] + vn[0:-2, -1]))) + + # Periodic BC v @ x = 0 + v[1:-1, + 0] = (vn[1:-1, 0] - un[1:-1, 0] * dt / dx * + (vn[1:-1, 0] - vn[1:-1, -1]) - vn[1:-1, 0] * dt / dy * + (vn[1:-1, 0] - vn[0:-2, 0]) - dt / (2 * rho * dy) * + (p[2:, 0] - p[0:-2, 0]) + nu * + (dt / dx**2 * + (vn[1:-1, 1] - 2 * vn[1:-1, 0] + vn[1:-1, -1]) + dt / dy**2 * + (vn[2:, 0] - 2 * vn[1:-1, 0] + vn[0:-2, 0]))) + + # Wall BC: u,v = 0 @ y = 0,2 + u[0, :] = 0 + u[-1, :] = 0 + v[0, :] = 0 + v[-1, :] = 0 + + udiff = (np.sum(u) - np.sum(un)) / np.sum(u) + stepcount += 1 + + return stepcount diff --git a/npbench/benchmarks/channel_flow/channel_flow_numba_n.py b/npbench/benchmarks/channel_flow/channel_flow_numba_n.py index 9c80f27ae..96347011d 100644 --- a/npbench/benchmarks/channel_flow/channel_flow_numba_n.py +++ b/npbench/benchmarks/channel_flow/channel_flow_numba_n.py @@ -43,7 +43,7 @@ def build_up_b(rho, dt, dx, dy, u, v): @nb.jit(nopython=True, parallel=False, fastmath=True) def pressure_poisson_periodic(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) for q in range(nit): pn = p.copy() diff --git a/npbench/benchmarks/channel_flow/channel_flow_numba_np.py b/npbench/benchmarks/channel_flow/channel_flow_numba_np.py index ccd6ac43d..54bb35dd2 100644 --- a/npbench/benchmarks/channel_flow/channel_flow_numba_np.py +++ b/npbench/benchmarks/channel_flow/channel_flow_numba_np.py @@ -43,7 +43,7 @@ def build_up_b(rho, dt, dx, dy, u, v): @nb.jit(nopython=True, parallel=True, fastmath=True) def pressure_poisson_periodic(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) for q in range(nit): pn = p.copy() diff --git a/npbench/benchmarks/channel_flow/channel_flow_numpy.py b/npbench/benchmarks/channel_flow/channel_flow_numpy.py index 4cbef94dd..5142982b1 100644 --- a/npbench/benchmarks/channel_flow/channel_flow_numpy.py +++ b/npbench/benchmarks/channel_flow/channel_flow_numpy.py @@ -40,7 +40,7 @@ def build_up_b(rho, dt, dx, dy, u, v): def pressure_poisson_periodic(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) for q in range(nit): pn = p.copy() diff --git a/npbench/benchmarks/channel_flow/channel_flow_pythran.py b/npbench/benchmarks/channel_flow/channel_flow_pythran.py index 23df35070..73937d824 100644 --- a/npbench/benchmarks/channel_flow/channel_flow_pythran.py +++ b/npbench/benchmarks/channel_flow/channel_flow_pythran.py @@ -44,7 +44,7 @@ def build_up_b(rho, dt, dx, dy, u, v): # pythran export pressure_poisson_periodic(int64, float64[:,:], float64, # float64, float64[:,:]) def pressure_poisson_periodic(nit, p, dx, dy, b): - pn = np.empty_like(p) + pn = np.zeros_like(p) for q in range(nit): pn = p.copy() diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_cupy.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_cupy.py index c1dc6870e..b5c22b2ac 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_cupy.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_cupy.py @@ -8,7 +8,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_dpnp.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_dpnp.py index fd0c640c2..7a916a1ed 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_dpnp.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_dpnp.py @@ -7,7 +7,7 @@ def conv2d(input, weights): C_out = weights.shape[3] # Output channels H_out = H - K + 1 W_out = W - K + 1 - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Perform convolution manually by iterating over the kernel dimensions for i in range(K): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_jax.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_jax.py index 6a5055170..78db5cd3e 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_jax.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_jax.py @@ -11,7 +11,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = jnp.empty((N, H_out, W_out, C_out), dtype=jnp.float32) + output = jnp.zeros((N, H_out, W_out, C_out), dtype=jnp.float32) def row_update(output, i): def col_update(output, j): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_n.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_n.py index 5d1d322f8..ad5a6e2ce 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_n.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_n.py @@ -11,7 +11,7 @@ def conv2d(input, weights): W_out = input.shape[2] - K + 1 C_in = input.shape[3] C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_np.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_np.py index bdaabae46..37f735b74 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_np.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_np.py @@ -11,7 +11,7 @@ def conv2d(input, weights): W_out = input.shape[2] - K + 1 C_in = input.shape[3] C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_npr.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_npr.py index 774659564..405251504 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_npr.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_npr.py @@ -11,7 +11,7 @@ def conv2d(input, weights): W_out = input.shape[2] - K + 1 C_in = input.shape[3] C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in nb.prange(H_out): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_o.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_o.py index 55d8722bf..b57e375fb 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_o.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_o.py @@ -10,7 +10,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_op.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_op.py index eb81880c3..f0e644757 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_op.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_op.py @@ -10,7 +10,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_opr.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_opr.py index dae6cc52f..7d53c057c 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_opr.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numba_opr.py @@ -10,7 +10,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in nb.prange(H_out): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numpy.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numpy.py index 7ef66e95b..bc8cd3fc6 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numpy.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_numpy.py @@ -8,7 +8,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_pythran.py b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_pythran.py index b84f835d3..0a26dade2 100644 --- a/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_pythran.py +++ b/npbench/benchmarks/deep_learning/conv2d_bias/conv2d_pythran.py @@ -9,7 +9,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_cupy.py b/npbench/benchmarks/deep_learning/lenet/lenet_cupy.py index dbe621c5e..1382f1b56 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_cupy.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_cupy.py @@ -12,7 +12,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): @@ -28,7 +28,7 @@ def conv2d(input, weights): # 2x2 maxpool operator, as used in LeNet-5 def maxpool2d(x): - output = np.empty( + output = np.zeros( [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], dtype=x.dtype) for i in range(x.shape[1] // 2): diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_jax.py b/npbench/benchmarks/deep_learning/lenet/lenet_jax.py index f1dbbd4e1..7d31f383c 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_jax.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_jax.py @@ -16,7 +16,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = jnp.empty((N, H_out, W_out, C_out), dtype=jnp.float32) + output = jnp.zeros((N, H_out, W_out, C_out), dtype=jnp.float32) def row_update(output, i): def col_update(output, j): @@ -47,7 +47,7 @@ def col_update(output, j): # 2x2 maxpool operator, as used in LeNet-5 @jax.jit def maxpool2d(x): - output = jnp.empty( + output = jnp.zeros( [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], dtype=x.dtype) diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_numba_n.py b/npbench/benchmarks/deep_learning/lenet/lenet_numba_n.py index 5a9d5d224..259885e10 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_numba_n.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_numba_n.py @@ -16,7 +16,7 @@ def conv2d(input, weights): W_out = input.shape[2] - K + 1 C_in = input.shape[3] C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): @@ -45,7 +45,7 @@ def maxpool2d(x): # output = np.empty( # [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], # dtype=x.dtype) - output = np.empty( + output = np.zeros( (x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]), dtype=x.dtype) for i in range(x.shape[1] // 2): diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_numba_np.py b/npbench/benchmarks/deep_learning/lenet/lenet_numba_np.py index 84066995e..dc41e139b 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_numba_np.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_numba_np.py @@ -16,7 +16,7 @@ def conv2d(input, weights): W_out = input.shape[2] - K + 1 C_in = input.shape[3] C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): @@ -45,7 +45,7 @@ def maxpool2d(x): # output = np.empty( # [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], # dtype=x.dtype) - output = np.empty( + output = np.zeros( (x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]), dtype=x.dtype) for i in range(x.shape[1] // 2): diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_numba_npr.py b/npbench/benchmarks/deep_learning/lenet/lenet_numba_npr.py index 70c5e2a00..bfbc8f482 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_numba_npr.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_numba_npr.py @@ -16,7 +16,7 @@ def conv2d(input, weights): W_out = input.shape[2] - K + 1 C_in = input.shape[3] C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in nb.prange(H_out): @@ -45,7 +45,7 @@ def maxpool2d(x): # output = np.empty( # [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], # dtype=x.dtype) - output = np.empty( + output = np.zeros( (x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]), dtype=x.dtype) for i in nb.prange(x.shape[1] // 2): diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_numba_o.py b/npbench/benchmarks/deep_learning/lenet/lenet_numba_o.py index 776a07329..6673c45ac 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_numba_o.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_numba_o.py @@ -15,7 +15,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): @@ -32,7 +32,7 @@ def conv2d(input, weights): # 2x2 maxpool operator, as used in LeNet-5 @nb.jit(nopython=False, forceobj=True, parallel=False, fastmath=True) def maxpool2d(x): - output = np.empty( + output = np.zeros( [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], dtype=x.dtype) for i in range(x.shape[1] // 2): diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_numba_op.py b/npbench/benchmarks/deep_learning/lenet/lenet_numba_op.py index 8afd5bad9..ca4401605 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_numba_op.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_numba_op.py @@ -15,7 +15,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): @@ -32,7 +32,7 @@ def conv2d(input, weights): # 2x2 maxpool operator, as used in LeNet-5 @nb.jit(nopython=False, forceobj=True, parallel=True, fastmath=True) def maxpool2d(x): - output = np.empty( + output = np.zeros( [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], dtype=x.dtype) for i in range(x.shape[1] // 2): diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_numba_opr.py b/npbench/benchmarks/deep_learning/lenet/lenet_numba_opr.py index 79f430bff..74cb3f19a 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_numba_opr.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_numba_opr.py @@ -15,7 +15,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in nb.prange(H_out): @@ -32,7 +32,7 @@ def conv2d(input, weights): # 2x2 maxpool operator, as used in LeNet-5 @nb.jit(nopython=False, forceobj=True, parallel=True, fastmath=True) def maxpool2d(x): - output = np.empty( + output = np.zeros( [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], dtype=x.dtype) for i in nb.prange(x.shape[1] // 2): diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_numpy.py b/npbench/benchmarks/deep_learning/lenet/lenet_numpy.py index 5c8f4afca..51053edf8 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_numpy.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_numpy.py @@ -12,7 +12,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): @@ -28,7 +28,7 @@ def conv2d(input, weights): # 2x2 maxpool operator, as used in LeNet-5 def maxpool2d(x): - output = np.empty( + output = np.zeros( [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], dtype=x.dtype) for i in range(x.shape[1] // 2): diff --git a/npbench/benchmarks/deep_learning/lenet/lenet_pythran.py b/npbench/benchmarks/deep_learning/lenet/lenet_pythran.py index 51a750219..601b65e25 100644 --- a/npbench/benchmarks/deep_learning/lenet/lenet_pythran.py +++ b/npbench/benchmarks/deep_learning/lenet/lenet_pythran.py @@ -13,7 +13,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): @@ -36,7 +36,7 @@ def conv2d(input, weights): # 2x2 maxpool operator, as used in LeNet-5 # pythran export maxpool2d(float32[:,:,:,:]) def maxpool2d(x): - output = np.empty( + output = np.zeros( [x.shape[0], x.shape[1] // 2, x.shape[2] // 2, x.shape[3]], dtype=x.dtype) for i in range(x.shape[1] // 2): diff --git a/npbench/benchmarks/deep_learning/mlp/mlp.py b/npbench/benchmarks/deep_learning/mlp/mlp.py index 27966a0ad..7fe4f4589 100644 --- a/npbench/benchmarks/deep_learning/mlp/mlp.py +++ b/npbench/benchmarks/deep_learning/mlp/mlp.py @@ -9,7 +9,7 @@ def initialize(C_in, N, S0, S1, S2): mlp_sizes = [S0, S1, S2] # [300, 100, 10] # Inputs - input = np.random.rand(N, C_in).astype(np.float32) + input = rng.random((N, C_in), dtype=np.float32) # Weights w1 = rng.random((C_in, mlp_sizes[0]), dtype=np.float32) b1 = rng.random((mlp_sizes[0], ), dtype=np.float32) diff --git a/npbench/benchmarks/deep_learning/mlp/mlp_numba_n.py b/npbench/benchmarks/deep_learning/mlp/mlp_numba_n.py index 0388915ee..6b9b74e18 100644 --- a/npbench/benchmarks/deep_learning/mlp/mlp_numba_n.py +++ b/npbench/benchmarks/deep_learning/mlp/mlp_numba_n.py @@ -12,7 +12,7 @@ def relu(x): def softmax(x): new_shape = (x.shape[0], 1) # tmp_max = np.max(x, axis=-1, keepdims=True) - tmp_max = np.empty(new_shape, dtype=x.dtype) + tmp_max = np.zeros(new_shape, dtype=x.dtype) for i in range(x.shape[1]): tmp_max[:, 0] = np.max(x[:, i]) tmp_out = np.exp(x - tmp_max) diff --git a/npbench/benchmarks/deep_learning/mlp/mlp_numba_np.py b/npbench/benchmarks/deep_learning/mlp/mlp_numba_np.py index 56f1b15d1..da9dc2d46 100644 --- a/npbench/benchmarks/deep_learning/mlp/mlp_numba_np.py +++ b/npbench/benchmarks/deep_learning/mlp/mlp_numba_np.py @@ -12,7 +12,7 @@ def relu(x): def softmax(x): new_shape = (x.shape[0], 1) # tmp_max = np.max(x, axis=-1, keepdims=True) - tmp_max = np.empty(new_shape, dtype=x.dtype) + tmp_max = np.zeros(new_shape, dtype=x.dtype) for i in range(x.shape[1]): tmp_max[:, 0] = np.max(x[:, i]) tmp_out = np.exp(x - tmp_max) diff --git a/npbench/benchmarks/deep_learning/mlp/mlp_numba_npr.py b/npbench/benchmarks/deep_learning/mlp/mlp_numba_npr.py index 42255eb2c..2b9172832 100644 --- a/npbench/benchmarks/deep_learning/mlp/mlp_numba_npr.py +++ b/npbench/benchmarks/deep_learning/mlp/mlp_numba_npr.py @@ -12,7 +12,7 @@ def relu(x): def softmax(x): new_shape = (x.shape[0], 1) # tmp_max = np.max(x, axis=-1, keepdims=True) - tmp_max = np.empty(new_shape, dtype=x.dtype) + tmp_max = np.zeros(new_shape, dtype=x.dtype) for i in nb.prange(x.shape[1]): tmp_max[:, 0] = np.max(x[:, i]) tmp_out = np.exp(x - tmp_max) diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_cupy.py b/npbench/benchmarks/deep_learning/resnet/resnet_cupy.py index bbf3112e9..b934366a7 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_cupy.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_cupy.py @@ -12,7 +12,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_jax.py b/npbench/benchmarks/deep_learning/resnet/resnet_jax.py index e573da9d9..3a34a6ed3 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_jax.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_jax.py @@ -15,7 +15,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = jnp.empty((N, H_out, W_out, C_out), dtype=jnp.float32) + output = jnp.zeros((N, H_out, W_out, C_out), dtype=jnp.float32) def row_update(output, i): def col_update(output, j): diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_numba_n.py b/npbench/benchmarks/deep_learning/resnet/resnet_numba_n.py index 0bd9b8f75..dd5bd7243 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_numba_n.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_numba_n.py @@ -16,7 +16,7 @@ def conv2d(input, weights): W_out = input.shape[2] - K + 1 C_in = input.shape[3] C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): @@ -43,10 +43,10 @@ def conv2d(input, weights): @nb.jit(nopython=True, parallel=False, fastmath=True) def batchnorm2d(x, eps=1e-5): # mean = np.mean(x, axis=0, keepdims=True) - mean = np.empty(x.shape, dtype=x.dtype) + mean = np.zeros(x.shape, dtype=x.dtype) mean[:] = np.sum(x, axis=0) / x.shape[0] # std = np.std(x, axis=0, keepdims=True) - std = np.empty(x.shape, dtype=x.dtype) + std = np.zeros(x.shape, dtype=x.dtype) std[:] = np.sqrt(np.sum((x - mean)**2, axis=0) / x.shape[0]) return (x - mean) / np.sqrt(std + eps) diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_numba_np.py b/npbench/benchmarks/deep_learning/resnet/resnet_numba_np.py index abe88bbe1..427c31a4a 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_numba_np.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_numba_np.py @@ -16,7 +16,7 @@ def conv2d(input, weights): W_out = input.shape[2] - K + 1 C_in = input.shape[3] C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): @@ -43,10 +43,10 @@ def conv2d(input, weights): @nb.jit(nopython=True, parallel=True, fastmath=True) def batchnorm2d(x, eps=1e-5): # mean = np.mean(x, axis=0, keepdims=True) - mean = np.empty(x.shape, dtype=x.dtype) + mean = np.zeros(x.shape, dtype=x.dtype) mean[:] = np.sum(x, axis=0) / x.shape[0] # std = np.std(x, axis=0, keepdims=True) - std = np.empty(x.shape, dtype=x.dtype) + std = np.zeros(x.shape, dtype=x.dtype) std[:] = np.sqrt(np.sum((x - mean)**2, axis=0) / x.shape[0]) return (x - mean) / np.sqrt(std + eps) diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_numba_npr.py b/npbench/benchmarks/deep_learning/resnet/resnet_numba_npr.py index 672100759..f4266aa06 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_numba_npr.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_numba_npr.py @@ -16,7 +16,7 @@ def conv2d(input, weights): W_out = input.shape[2] - K + 1 C_in = input.shape[3] C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in nb.prange(H_out): @@ -43,10 +43,10 @@ def conv2d(input, weights): @nb.jit(nopython=True, parallel=True, fastmath=True) def batchnorm2d(x, eps=1e-5): # mean = np.mean(x, axis=0, keepdims=True) - mean = np.empty(x.shape, dtype=x.dtype) + mean = np.zeros(x.shape, dtype=x.dtype) mean[:] = np.sum(x, axis=0) / x.shape[0] # std = np.std(x, axis=0, keepdims=True) - std = np.empty(x.shape, dtype=x.dtype) + std = np.zeros(x.shape, dtype=x.dtype) std[:] = np.sqrt(np.sum((x - mean)**2, axis=0) / x.shape[0]) return (x - mean) / np.sqrt(std + eps) diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_numba_o.py b/npbench/benchmarks/deep_learning/resnet/resnet_numba_o.py index 12e1dde5b..26d336f8d 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_numba_o.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_numba_o.py @@ -15,7 +15,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_numba_op.py b/npbench/benchmarks/deep_learning/resnet/resnet_numba_op.py index 02b6e7e04..6beed94de 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_numba_op.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_numba_op.py @@ -15,7 +15,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_numba_opr.py b/npbench/benchmarks/deep_learning/resnet/resnet_numba_opr.py index 9352b6662..63143f140 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_numba_opr.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_numba_opr.py @@ -15,7 +15,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in nb.prange(H_out): diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_numpy.py b/npbench/benchmarks/deep_learning/resnet/resnet_numpy.py index 313607b50..69113d51c 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_numpy.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_numpy.py @@ -12,7 +12,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/resnet/resnet_pythran.py b/npbench/benchmarks/deep_learning/resnet/resnet_pythran.py index 2d1529d4a..e5e1ec77d 100644 --- a/npbench/benchmarks/deep_learning/resnet/resnet_pythran.py +++ b/npbench/benchmarks/deep_learning/resnet/resnet_pythran.py @@ -14,7 +14,7 @@ def conv2d(input, weights): H_out = input.shape[1] - K + 1 W_out = input.shape[2] - K + 1 C_out = weights.shape[3] - output = np.empty((N, H_out, W_out, C_out), dtype=np.float32) + output = np.zeros((N, H_out, W_out, C_out), dtype=np.float32) # Loop structure adapted from https://github.com/SkalskiP/ILearnDeepLearning.py/blob/ba0b5ba589d4e656141995e8d1a06d44db6ce58d/01_mysteries_of_neural_networks/06_numpy_convolutional_neural_net/src/layers/convolutional.py#L88 for i in range(H_out): diff --git a/npbench/benchmarks/deep_learning/softmax/softmax_numba_n.py b/npbench/benchmarks/deep_learning/softmax/softmax_numba_n.py index 50b67edc1..ffb9c6db6 100644 --- a/npbench/benchmarks/deep_learning/softmax/softmax_numba_n.py +++ b/npbench/benchmarks/deep_learning/softmax/softmax_numba_n.py @@ -7,7 +7,7 @@ def softmax(x): new_shape = (x.shape[0], x.shape[1], x.shape[2], 1) # tmp_max = np.max(x, axis=-1, keepdims=True) - tmp_max = np.empty(new_shape, dtype=x.dtype) + tmp_max = np.zeros(new_shape, dtype=x.dtype) for i in range(x.shape[3]): tmp_max[:, :, :, 0] = np.max(x[:, :, :, i]) tmp_out = np.exp(x - tmp_max) diff --git a/npbench/benchmarks/deep_learning/softmax/softmax_numba_np.py b/npbench/benchmarks/deep_learning/softmax/softmax_numba_np.py index 366adb76d..7e8050f75 100644 --- a/npbench/benchmarks/deep_learning/softmax/softmax_numba_np.py +++ b/npbench/benchmarks/deep_learning/softmax/softmax_numba_np.py @@ -7,7 +7,7 @@ def softmax(x): new_shape = (x.shape[0], x.shape[1], x.shape[2], 1) # tmp_max = np.max(x, axis=-1, keepdims=True) - tmp_max = np.empty(new_shape, dtype=x.dtype) + tmp_max = np.zeros(new_shape, dtype=x.dtype) for i in range(x.shape[3]): tmp_max[:, :, :, 0] = np.max(x[:, :, :, i]) tmp_out = np.exp(x - tmp_max) diff --git a/npbench/benchmarks/deep_learning/softmax/softmax_numba_npr.py b/npbench/benchmarks/deep_learning/softmax/softmax_numba_npr.py index 31dd65567..36310a442 100644 --- a/npbench/benchmarks/deep_learning/softmax/softmax_numba_npr.py +++ b/npbench/benchmarks/deep_learning/softmax/softmax_numba_npr.py @@ -7,7 +7,7 @@ def softmax(x): new_shape = (x.shape[0], x.shape[1], x.shape[2], 1) # tmp_max = np.max(x, axis=-1, keepdims=True) - tmp_max = np.empty(new_shape, dtype=x.dtype) + tmp_max = np.zeros(new_shape, dtype=x.dtype) for i in nb.prange(x.shape[3]): tmp_max[:, :, :, 0] = np.max(x[:, :, :, i]) tmp_out = np.exp(x - tmp_max) diff --git a/npbench/benchmarks/mandelbrot1/mandelbrot1_legate.py b/npbench/benchmarks/mandelbrot1/mandelbrot1_legate.py index 7fab1e6e3..8b7c70cb8 100644 --- a/npbench/benchmarks/mandelbrot1/mandelbrot1_legate.py +++ b/npbench/benchmarks/mandelbrot1/mandelbrot1_legate.py @@ -8,7 +8,7 @@ def linspace(start, stop, num, dtype): - X = np.empty((num, ), dtype=dtype) + X = np.zeros((num, ), dtype=dtype) dist = (stop - start) / (num - 1) for i in range(num): X[i] = start + i * dist diff --git a/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_n.py b/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_n.py index 4c2e4bcf2..dd698c46b 100644 --- a/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_n.py +++ b/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_n.py @@ -10,7 +10,7 @@ @nb.jit(nopython=True, parallel=False, fastmath=True) def linspace(start, stop, num, dtype): - X = np.empty((num, ), dtype=dtype) + X = np.zeros((num, ), dtype=dtype) dist = (stop - start) / (num - 1) for i in range(num): X[i] = start + i * dist diff --git a/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_np.py b/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_np.py index 89844a117..f47416170 100644 --- a/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_np.py +++ b/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_np.py @@ -10,7 +10,7 @@ @nb.jit(nopython=True, parallel=True, fastmath=True) def linspace(start, stop, num, dtype): - X = np.empty((num, ), dtype=dtype) + X = np.zeros((num, ), dtype=dtype) dist = (stop - start) / (num - 1) for i in range(num): X[i] = start + i * dist diff --git a/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_npr.py b/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_npr.py index a38baa314..70cb889a9 100644 --- a/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_npr.py +++ b/npbench/benchmarks/mandelbrot1/mandelbrot1_numba_npr.py @@ -10,7 +10,7 @@ @nb.jit(nopython=True, parallel=True, fastmath=True) def linspace(start, stop, num, dtype): - X = np.empty((num, ), dtype=dtype) + X = np.zeros((num, ), dtype=dtype) dist = (stop - start) / (num - 1) for i in nb.prange(num): X[i] = start + i * dist diff --git a/npbench/benchmarks/mandelbrot2/mandelbrot2_legate.py b/npbench/benchmarks/mandelbrot2/mandelbrot2_legate.py index 7a4ea6098..4be3b0945 100644 --- a/npbench/benchmarks/mandelbrot2/mandelbrot2_legate.py +++ b/npbench/benchmarks/mandelbrot2/mandelbrot2_legate.py @@ -8,7 +8,7 @@ def linspace(start, stop, num, dtype): - X = np.empty((num, ), dtype=dtype) + X = np.zeros((num, ), dtype=dtype) dist = (stop - start) / (num - 1) for i in range(num): X[i] = start + i * dist diff --git a/npbench/benchmarks/mandelbrot2/mandelbrot2_numba_n.py b/npbench/benchmarks/mandelbrot2/mandelbrot2_numba_n.py index f0c00f7df..7f5f06a7f 100644 --- a/npbench/benchmarks/mandelbrot2/mandelbrot2_numba_n.py +++ b/npbench/benchmarks/mandelbrot2/mandelbrot2_numba_n.py @@ -10,8 +10,8 @@ @nb.jit(nopython=True, parallel=False, fastmath=True) def mgrid(xn, yn): - Xi = np.empty((xn, yn), dtype=np.int64) - Yi = np.empty((xn, yn), dtype=np.int64) + Xi = np.zeros((xn, yn), dtype=np.int64) + Yi = np.zeros((xn, yn), dtype=np.int64) for i in range(xn): Xi[i, :] = i for j in range(yn): @@ -21,7 +21,7 @@ def mgrid(xn, yn): @nb.jit(nopython=True, parallel=False, fastmath=True) def linspace(start, stop, num, dtype): - X = np.empty((num, ), dtype=dtype) + X = np.zeros((num, ), dtype=dtype) dist = (stop - start) / (num - 1) for i in range(num): X[i] = start + i * dist diff --git a/npbench/benchmarks/mandelbrot2/mandelbrot2_pythran.py b/npbench/benchmarks/mandelbrot2/mandelbrot2_pythran.py index 455c3f2ce..107a60e0a 100644 --- a/npbench/benchmarks/mandelbrot2/mandelbrot2_pythran.py +++ b/npbench/benchmarks/mandelbrot2/mandelbrot2_pythran.py @@ -8,8 +8,8 @@ def mgrid(xn, yn): - Xi = np.empty((xn, yn), dtype=np.int64) - Yi = np.empty((xn, yn), dtype=np.int64) + Xi = np.zeros((xn, yn), dtype=np.int64) + Yi = np.zeros((xn, yn), dtype=np.int64) for i in range(xn): Xi[i, :] = i for j in range(yn): diff --git a/npbench/benchmarks/nbody/nbody_jax.py b/npbench/benchmarks/nbody/nbody_jax.py index ca6a635db..d34dfbab8 100644 --- a/npbench/benchmarks/nbody/nbody_jax.py +++ b/npbench/benchmarks/nbody/nbody_jax.py @@ -92,8 +92,8 @@ def nbody(mass, pos, vel, N, Nt, dt, G, softening): acc = getAcc(pos, mass, G, softening) # calculate initial energy of system - KE = jnp.empty(Nt + 1, dtype=jnp.float64) - PE = jnp.empty(Nt + 1, dtype=jnp.float64) + KE = jnp.zeros(Nt + 1, dtype=jnp.float64) + PE = jnp.zeros(Nt + 1, dtype=jnp.float64) ke, pe = getEnergy(pos, vel, mass, G) KE = KE.at[0].set(ke) PE = PE.at[0].set(pe) diff --git a/npbench/benchmarks/nbody/nbody_numba_n.py b/npbench/benchmarks/nbody/nbody_numba_n.py index c256c38fb..f20c31b78 100644 --- a/npbench/benchmarks/nbody/nbody_numba_n.py +++ b/npbench/benchmarks/nbody/nbody_numba_n.py @@ -106,8 +106,8 @@ def nbody(mass, pos, vel, N, Nt, dt, G, softening): # calculate initial energy of system # KE = np.ndarray(Nt+1, dtype=np.float64) # PE = np.ndarray(Nt+1, dtype=np.float64) - KE = np.empty(Nt + 1, dtype=np.float64) - PE = np.empty(Nt + 1, dtype=np.float64) + KE = np.zeros(Nt + 1, dtype=np.float64) + PE = np.zeros(Nt + 1, dtype=np.float64) KE[0], PE[0] = getEnergy(pos, vel, mass, G) t = 0.0 diff --git a/npbench/benchmarks/nbody/nbody_pythran.py b/npbench/benchmarks/nbody/nbody_pythran.py index abb4e876c..580982628 100644 --- a/npbench/benchmarks/nbody/nbody_pythran.py +++ b/npbench/benchmarks/nbody/nbody_pythran.py @@ -95,8 +95,8 @@ def nbody(mass, pos, vel, N, Nt, dt, G, softening): # calculate initial energy of system # KE = np.ndarray(Nt+1, dtype=np.float64) # PE = np.ndarray(Nt+1, dtype=np.float64) - KE = np.empty(Nt + 1, dtype=np.float64) - PE = np.empty(Nt + 1, dtype=np.float64) + KE = np.zeros(Nt + 1, dtype=np.float64) + PE = np.zeros(Nt + 1, dtype=np.float64) KE[0], PE[0] = getEnergy(pos, vel, mass, G) t = 0.0 diff --git a/npbench/benchmarks/polybench/adi/adi_cupy.py b/npbench/benchmarks/polybench/adi/adi_cupy.py index 7371f9070..48ebc3909 100644 --- a/npbench/benchmarks/polybench/adi/adi_cupy.py +++ b/npbench/benchmarks/polybench/adi/adi_cupy.py @@ -5,9 +5,9 @@ def kernel(TSTEPS, N, u): - v = np.empty(u.shape, dtype=u.dtype) - p = np.empty(u.shape, dtype=u.dtype) - q = np.empty(u.shape, dtype=u.dtype) + v = np.zeros(u.shape, dtype=u.dtype) + p = np.zeros(u.shape, dtype=u.dtype) + q = np.zeros(u.shape, dtype=u.dtype) DX = 1.0 / N DY = 1.0 / N diff --git a/npbench/benchmarks/polybench/adi/adi_dace.py b/npbench/benchmarks/polybench/adi/adi_dace.py index 18c87a768..346a7359f 100644 --- a/npbench/benchmarks/polybench/adi/adi_dace.py +++ b/npbench/benchmarks/polybench/adi/adi_dace.py @@ -9,9 +9,9 @@ @dc.program def kernel(TSTEPS: dc.int64, u: dc.float64[N, N]): - v = np.empty(u.shape, dtype=u.dtype) - p = np.empty(u.shape, dtype=u.dtype) - q = np.empty(u.shape, dtype=u.dtype) + v = np.zeros(u.shape, dtype=u.dtype) + p = np.zeros(u.shape, dtype=u.dtype) + q = np.zeros(u.shape, dtype=u.dtype) DX = 1.0 / np.float64(N) DY = 1.0 / np.float64(N) diff --git a/npbench/benchmarks/polybench/adi/adi_dask.py b/npbench/benchmarks/polybench/adi/adi_dask.py index d57f60a1b..7354b6ba6 100644 --- a/npbench/benchmarks/polybench/adi/adi_dask.py +++ b/npbench/benchmarks/polybench/adi/adi_dask.py @@ -7,9 +7,9 @@ def kernel(TSTEPS, N, np_u): u = np.from_array(np_u, chunks='auto') - v = np.empty(u.shape, dtype=u.dtype) - p = np.empty(u.shape, dtype=u.dtype) - q = np.empty(u.shape, dtype=u.dtype) + v = np.zeros(u.shape, dtype=u.dtype) + p = np.zeros(u.shape, dtype=u.dtype) + q = np.zeros(u.shape, dtype=u.dtype) DX = 1.0 / N DY = 1.0 / N diff --git a/npbench/benchmarks/polybench/adi/adi_dpnp.py b/npbench/benchmarks/polybench/adi/adi_dpnp.py index a8ba990a1..f2657589e 100644 --- a/npbench/benchmarks/polybench/adi/adi_dpnp.py +++ b/npbench/benchmarks/polybench/adi/adi_dpnp.py @@ -1,9 +1,9 @@ import dpnp as np def kernel(TSTEPS, N, u): - v = np.empty(u.shape, dtype=u.dtype) - p = np.empty(u.shape, dtype=u.dtype) - q = np.empty(u.shape, dtype=u.dtype) + v = np.zeros(u.shape, dtype=u.dtype) + p = np.zeros(u.shape, dtype=u.dtype) + q = np.zeros(u.shape, dtype=u.dtype) DX = 1.0 / N DY = 1.0 / N diff --git a/npbench/benchmarks/polybench/adi/adi_legate.py b/npbench/benchmarks/polybench/adi/adi_legate.py index f4a4eeb9b..212c96f41 100644 --- a/npbench/benchmarks/polybench/adi/adi_legate.py +++ b/npbench/benchmarks/polybench/adi/adi_legate.py @@ -5,9 +5,9 @@ def kernel(TSTEPS, N, u): - v = np.empty(u.shape, dtype=u.dtype) - p = np.empty(u.shape, dtype=u.dtype) - q = np.empty(u.shape, dtype=u.dtype) + v = np.zeros(u.shape, dtype=u.dtype) + p = np.zeros(u.shape, dtype=u.dtype) + q = np.zeros(u.shape, dtype=u.dtype) DX = 1.0 / N DY = 1.0 / N diff --git a/npbench/benchmarks/polybench/adi/adi_numba_n.py b/npbench/benchmarks/polybench/adi/adi_numba_n.py index ca551cda2..17bbb07fb 100644 --- a/npbench/benchmarks/polybench/adi/adi_numba_n.py +++ b/npbench/benchmarks/polybench/adi/adi_numba_n.py @@ -7,9 +7,9 @@ @nb.jit(nopython=True, parallel=False, fastmath=True) def kernel(TSTEPS, N, u): - v = np.empty(u.shape, dtype=u.dtype) - p = np.empty(u.shape, dtype=u.dtype) - q = np.empty(u.shape, dtype=u.dtype) + v = np.zeros(u.shape, dtype=u.dtype) + p = np.zeros(u.shape, dtype=u.dtype) + q = np.zeros(u.shape, dtype=u.dtype) DX = 1.0 / N DY = 1.0 / N diff --git a/npbench/benchmarks/polybench/adi/adi_numba_np.py b/npbench/benchmarks/polybench/adi/adi_numba_np.py index a31294d69..4402e78ea 100644 --- a/npbench/benchmarks/polybench/adi/adi_numba_np.py +++ b/npbench/benchmarks/polybench/adi/adi_numba_np.py @@ -7,9 +7,9 @@ @nb.jit(nopython=True, parallel=True, fastmath=True) def kernel(TSTEPS, N, u): - v = np.empty(u.shape, dtype=u.dtype) - p = np.empty(u.shape, dtype=u.dtype) - q = np.empty(u.shape, dtype=u.dtype) + v = np.zeros(u.shape, dtype=u.dtype) + p = np.zeros(u.shape, dtype=u.dtype) + q = np.zeros(u.shape, dtype=u.dtype) DX = 1.0 / N DY = 1.0 / N diff --git a/npbench/benchmarks/polybench/adi/adi_numpy.py b/npbench/benchmarks/polybench/adi/adi_numpy.py index c5df9a73e..74b3e49dd 100644 --- a/npbench/benchmarks/polybench/adi/adi_numpy.py +++ b/npbench/benchmarks/polybench/adi/adi_numpy.py @@ -5,9 +5,9 @@ def kernel(TSTEPS, N, u): - v = np.empty(u.shape, dtype=u.dtype) - p = np.empty(u.shape, dtype=u.dtype) - q = np.empty(u.shape, dtype=u.dtype) + v = np.zeros(u.shape, dtype=u.dtype) + p = np.zeros(u.shape, dtype=u.dtype) + q = np.zeros(u.shape, dtype=u.dtype) DX = 1.0 / N DY = 1.0 / N diff --git a/npbench/benchmarks/polybench/adi/adi_pythran.py b/npbench/benchmarks/polybench/adi/adi_pythran.py index 9aa421cdc..045bfcde7 100644 --- a/npbench/benchmarks/polybench/adi/adi_pythran.py +++ b/npbench/benchmarks/polybench/adi/adi_pythran.py @@ -6,9 +6,9 @@ # pythran export kernel(int, int, float64[:,:]) def kernel(TSTEPS, N, u): - v = np.empty(u.shape, dtype=u.dtype) - p = np.empty(u.shape, dtype=u.dtype) - q = np.empty(u.shape, dtype=u.dtype) + v = np.zeros(u.shape, dtype=u.dtype) + p = np.zeros(u.shape, dtype=u.dtype) + q = np.zeros(u.shape, dtype=u.dtype) DX = 1.0 / N DY = 1.0 / N diff --git a/npbench/benchmarks/polybench/atax/atax_legate.py b/npbench/benchmarks/polybench/atax/atax_legate.py index 2cb73da1f..37e4f7119 100644 --- a/npbench/benchmarks/polybench/atax/atax_legate.py +++ b/npbench/benchmarks/polybench/atax/atax_legate.py @@ -10,10 +10,10 @@ def kernel(A, x): def init_data(M, N, datatype): fn = datatype(N) - A = np.empty((M, N), dtype=datatype) - x = np.empty((N, ), dtype=datatype) - y = np.empty((N, ), dtype=datatype) - tmp = np.empty((M, ), dtype=datatype) + A = np.zeros((M, N), dtype=datatype) + x = np.zeros((N, ), dtype=datatype) + y = np.zeros((N, ), dtype=datatype) + tmp = np.zeros((M, ), dtype=datatype) # for i in range(N): # x[i] = 1 + (i / fn) # for i in range(M): diff --git a/npbench/benchmarks/polybench/bicg/bicg_legate.py b/npbench/benchmarks/polybench/bicg/bicg_legate.py index 041020ebd..b4206085c 100644 --- a/npbench/benchmarks/polybench/bicg/bicg_legate.py +++ b/npbench/benchmarks/polybench/bicg/bicg_legate.py @@ -9,11 +9,11 @@ def kernel(A, p, r): def init_data(M, N, datatype): - A = np.empty((N, M), dtype=datatype) - s = np.empty((M, ), dtype=datatype) - q = np.empty((N, ), dtype=datatype) - p = np.empty((M, ), dtype=datatype) - r = np.empty((N, ), dtype=datatype) + A = np.zeros((N, M), dtype=datatype) + s = np.zeros((M, ), dtype=datatype) + q = np.zeros((N, ), dtype=datatype) + p = np.zeros((M, ), dtype=datatype) + r = np.zeros((N, ), dtype=datatype) # for i in range(M): # p[i] = (i % M) / M # for i in range(N): diff --git a/npbench/benchmarks/polybench/cholesky/cholesky.py b/npbench/benchmarks/polybench/cholesky/cholesky.py index 1fe67bd77..969ae8706 100644 --- a/npbench/benchmarks/polybench/cholesky/cholesky.py +++ b/npbench/benchmarks/polybench/cholesky/cholesky.py @@ -4,7 +4,7 @@ def initialize(N, datatype=np.float64): - A = np.empty((N, N), dtype=datatype) + A = np.zeros((N, N), dtype=datatype) for i in range(N): A[i, :i + 1] = np.fromfunction(lambda j: (-j % N) / N + 1, (i + 1, ), dtype=datatype) diff --git a/npbench/benchmarks/polybench/cholesky/cholesky_legate.py b/npbench/benchmarks/polybench/cholesky/cholesky_legate.py index 48886107d..140a04818 100644 --- a/npbench/benchmarks/polybench/cholesky/cholesky_legate.py +++ b/npbench/benchmarks/polybench/cholesky/cholesky_legate.py @@ -19,7 +19,7 @@ def kernel2(A): def init_data(N, datatype): - A = np.empty((N, N), dtype=datatype) + A = np.zeros((N, N), dtype=datatype) # for i in range(N): # for j in range(i + 1): # A[i, j] = (-j % N) / N + 1 diff --git a/npbench/benchmarks/polybench/cholesky2/cholesky2.py b/npbench/benchmarks/polybench/cholesky2/cholesky2.py index 1fe67bd77..969ae8706 100644 --- a/npbench/benchmarks/polybench/cholesky2/cholesky2.py +++ b/npbench/benchmarks/polybench/cholesky2/cholesky2.py @@ -4,7 +4,7 @@ def initialize(N, datatype=np.float64): - A = np.empty((N, N), dtype=datatype) + A = np.zeros((N, N), dtype=datatype) for i in range(N): A[i, :i + 1] = np.fromfunction(lambda j: (-j % N) / N + 1, (i + 1, ), dtype=datatype) diff --git a/npbench/benchmarks/polybench/cholesky2/cholesky2_legate.py b/npbench/benchmarks/polybench/cholesky2/cholesky2_legate.py index 48886107d..140a04818 100644 --- a/npbench/benchmarks/polybench/cholesky2/cholesky2_legate.py +++ b/npbench/benchmarks/polybench/cholesky2/cholesky2_legate.py @@ -19,7 +19,7 @@ def kernel2(A): def init_data(N, datatype): - A = np.empty((N, N), dtype=datatype) + A = np.zeros((N, N), dtype=datatype) # for i in range(N): # for j in range(i + 1): # A[i, j] = (-j % N) / N + 1 diff --git a/npbench/benchmarks/polybench/cholesky2/cholesky2_pythran.py b/npbench/benchmarks/polybench/cholesky2/cholesky2_pythran.py index 5b99b9255..3daa37350 100644 --- a/npbench/benchmarks/polybench/cholesky2/cholesky2_pythran.py +++ b/npbench/benchmarks/polybench/cholesky2/cholesky2_pythran.py @@ -2,5 +2,5 @@ # pythran export kernel(float64[:,:]) -def kernel2(A): +def kernel(A): A[:] = np.linalg.cholesky(A) + np.triu(A, k=1) diff --git a/npbench/benchmarks/polybench/correlation/correlation_legate.py b/npbench/benchmarks/polybench/correlation/correlation_legate.py index a486046eb..318bc7ec7 100644 --- a/npbench/benchmarks/polybench/correlation/correlation_legate.py +++ b/npbench/benchmarks/polybench/correlation/correlation_legate.py @@ -29,7 +29,7 @@ def kernel(M, float_n, data): def init_data(M, N, datatype): float_n = datatype(N) - data = np.empty((N, M), dtype=datatype) + data = np.zeros((N, M), dtype=datatype) # for i in range(N): # for j in range(M): # data[i, j] = (i * j) / M + i diff --git a/npbench/benchmarks/polybench/covariance/covariance_legate.py b/npbench/benchmarks/polybench/covariance/covariance_legate.py index 0d5433a6b..049189b5e 100644 --- a/npbench/benchmarks/polybench/covariance/covariance_legate.py +++ b/npbench/benchmarks/polybench/covariance/covariance_legate.py @@ -20,7 +20,7 @@ def kernel(M, float_n, data): def init_data(M, N, datatype): float_n = datatype(N) - data = np.empty((N, M), dtype=datatype) + data = np.zeros((N, M), dtype=datatype) # for i in range(N): # for j in range(M): # data[i, j] = (i * j) / M + i diff --git a/npbench/benchmarks/polybench/deriche/deriche_cupy.py b/npbench/benchmarks/polybench/deriche/deriche_cupy.py index 6e71a2129..8fbb41073 100644 --- a/npbench/benchmarks/polybench/deriche/deriche_cupy.py +++ b/npbench/benchmarks/polybench/deriche/deriche_cupy.py @@ -13,14 +13,14 @@ def kernel(alpha, imgIn): b2 = -np.exp(-2.0 * alpha) c1 = c2 = 1 - y1 = np.empty_like(imgIn) + y1 = np.zeros_like(imgIn) y1[:, 0] = a1 * imgIn[:, 0] y1[:, 1] = a1 * imgIn[:, 1] + a2 * imgIn[:, 0] + b1 * y1[:, 0] for j in range(2, imgIn.shape[1]): y1[:, j] = (a1 * imgIn[:, j] + a2 * imgIn[:, j - 1] + b1 * y1[:, j - 1] + b2 * y1[:, j - 2]) - y2 = np.empty_like(imgIn) + y2 = np.zeros_like(imgIn) y2[:, -1] = 0.0 y2[:, -2] = a3 * imgIn[:, -1] for j in range(imgIn.shape[1] - 3, -1, -1): diff --git a/npbench/benchmarks/polybench/deriche/deriche_dace.py b/npbench/benchmarks/polybench/deriche/deriche_dace.py index 2272eca41..626fe515e 100644 --- a/npbench/benchmarks/polybench/deriche/deriche_dace.py +++ b/npbench/benchmarks/polybench/deriche/deriche_dace.py @@ -22,14 +22,14 @@ def kernel(alpha: dc.float64, imgIn: dc.float64[W, H]): c1 = 1 c2 = 1 - y1 = np.empty_like(imgIn) + y1 = np.zeros_like(imgIn) y1[:, 0] = a1 * imgIn[:, 0] y1[:, 1] = a1 * imgIn[:, 1] + a2 * imgIn[:, 0] + b1 * y1[:, 0] for j in range(2, H): y1[:, j] = (a1 * imgIn[:, j] + a2 * imgIn[:, j - 1] + b1 * y1[:, j - 1] + b2 * y1[:, j - 2]) - y2 = np.empty_like(imgIn) + y2 = np.zeros_like(imgIn) y2[:, -1] = 0.0 y2[:, -2] = a3 * imgIn[:, -1] for j in range(H - 3, -1, -1): diff --git a/npbench/benchmarks/polybench/deriche/deriche_jax.py b/npbench/benchmarks/polybench/deriche/deriche_jax.py index 978c64f78..7b799564a 100644 --- a/npbench/benchmarks/polybench/deriche/deriche_jax.py +++ b/npbench/benchmarks/polybench/deriche/deriche_jax.py @@ -15,7 +15,7 @@ def kernel(alpha, imgIn): b2 = -jnp.exp(-2.0 * alpha) c1 = c2 = 1 - y1 = jnp.empty_like(imgIn) + y1 = jnp.zeros_like(imgIn) y1 = y1.at[:, 0].set(a1 * imgIn[:, 0]) y1 = y1.at[:, 1].set(a1 * imgIn[:, 1] + a2 * imgIn[:, 0] + b1 * y1[:, 0]) @@ -27,7 +27,7 @@ def horizontal_forward(j, y1): y1 = lax.fori_loop(2, imgIn.shape[1], horizontal_forward, y1) - y2 = jnp.empty_like(imgIn) + y2 = jnp.zeros_like(imgIn) y2 = y2.at[:, -1].set(0.0) y2 = y2.at[:, -2].set(a3 * imgIn[:, -1]) @@ -42,7 +42,7 @@ def horizontal_backward(j, y2): imgOut = c1 * (y1 + y2) - y1 = jnp.empty_like(imgOut) + y1 = jnp.zeros_like(imgOut) y1 = y1.at[0, :].set(a5 * imgOut[0, :]) y1 = y1.at[1, :].set(a5 * imgOut[1, :] + a6 * imgOut[0, :] + b1 * y1[0, :]) @@ -54,7 +54,7 @@ def vertical_forward(i, y1): y1 = lax.fori_loop(2, imgIn.shape[0], vertical_forward, y1) - y2 = jnp.empty_like(imgOut) + y2 = jnp.zeros_like(imgOut) y2 = y2.at[-1, :].set(0.0) y2 = y2.at[-2, :].set(a7 * imgOut[-1, :]) diff --git a/npbench/benchmarks/polybench/deriche/deriche_legate.py b/npbench/benchmarks/polybench/deriche/deriche_legate.py index 5ef7d9559..c5454a9e5 100644 --- a/npbench/benchmarks/polybench/deriche/deriche_legate.py +++ b/npbench/benchmarks/polybench/deriche/deriche_legate.py @@ -17,7 +17,7 @@ def kernel(alpha, imgIn): b2 = -np.exp(-2.0 * alpha) c1 = c2 = 1 - y1 = np.empty_like(imgIn) + y1 = np.zeros_like(imgIn) y1[:, 0] = a1 * imgIn[:, 0] # y1[:, 1] = a1 * imgIn[:, 1] + a2 * imgIn[:, 0] + b1 * y1[:, 0] y1[:, 1] = b1 * y1[:, 0] @@ -26,7 +26,7 @@ def kernel(alpha, imgIn): y1[:, j] = (a1 * imgIn[:, j] + a2 * imgIn[:, j - 1] + b1 * y1[:, j - 1] + b2 * y1[:, j - 2]) - y2 = np.empty_like(imgIn) + y2 = np.zeros_like(imgIn) y2[:, -1] = 0.0 y2[:, -2] = a3 * imgIn[:, -1] for j in range(imgIn.shape[1] - 3, -1, -1): @@ -55,10 +55,10 @@ def kernel(alpha, imgIn): def init_data(W, H, datatype): alpha = datatype(0.25) - imgIn = onp.empty((W, H), dtype=datatype) - imgOut = onp.empty((W, H), dtype=datatype) - y1 = onp.empty((W, H), dtype=datatype) - y2 = onp.empty((W, H), dtype=datatype) + imgIn = onp.zeros((W, H), dtype=datatype) + imgOut = onp.zeros((W, H), dtype=datatype) + y1 = onp.zeros((W, H), dtype=datatype) + y2 = onp.zeros((W, H), dtype=datatype) for i in range(W): for j in range(H): imgIn[i, j] = ((313 * i + 991 * j) % 65536) / 65535.0 diff --git a/npbench/benchmarks/polybench/deriche/deriche_numba_n.py b/npbench/benchmarks/polybench/deriche/deriche_numba_n.py index ae0d185d5..4e759f806 100644 --- a/npbench/benchmarks/polybench/deriche/deriche_numba_n.py +++ b/npbench/benchmarks/polybench/deriche/deriche_numba_n.py @@ -15,14 +15,14 @@ def kernel(alpha, imgIn): b2 = -np.exp(-2.0 * alpha) c1 = c2 = 1 - y1 = np.empty_like(imgIn) + y1 = np.zeros_like(imgIn) y1[:, 0] = a1 * imgIn[:, 0] y1[:, 1] = a1 * imgIn[:, 1] + a2 * imgIn[:, 0] + b1 * y1[:, 0] for j in range(2, imgIn.shape[1]): y1[:, j] = (a1 * imgIn[:, j] + a2 * imgIn[:, j - 1] + b1 * y1[:, j - 1] + b2 * y1[:, j - 2]) - y2 = np.empty_like(imgIn) + y2 = np.zeros_like(imgIn) y2[:, -1] = 0.0 y2[:, -2] = a3 * imgIn[:, -1] for j in range(imgIn.shape[1] - 3, -1, -1): diff --git a/npbench/benchmarks/polybench/deriche/deriche_numba_np.py b/npbench/benchmarks/polybench/deriche/deriche_numba_np.py index 10ffceae4..2cd900e2f 100644 --- a/npbench/benchmarks/polybench/deriche/deriche_numba_np.py +++ b/npbench/benchmarks/polybench/deriche/deriche_numba_np.py @@ -15,14 +15,14 @@ def kernel(alpha, imgIn): b2 = -np.exp(-2.0 * alpha) c1 = c2 = 1 - y1 = np.empty_like(imgIn) + y1 = np.zeros_like(imgIn) y1[:, 0] = a1 * imgIn[:, 0] y1[:, 1] = a1 * imgIn[:, 1] + a2 * imgIn[:, 0] + b1 * y1[:, 0] for j in range(2, imgIn.shape[1]): y1[:, j] = (a1 * imgIn[:, j] + a2 * imgIn[:, j - 1] + b1 * y1[:, j - 1] + b2 * y1[:, j - 2]) - y2 = np.empty_like(imgIn) + y2 = np.zeros_like(imgIn) y2[:, -1] = 0.0 y2[:, -2] = a3 * imgIn[:, -1] for j in range(imgIn.shape[1] - 3, -1, -1): diff --git a/npbench/benchmarks/polybench/deriche/deriche_numpy.py b/npbench/benchmarks/polybench/deriche/deriche_numpy.py index b0511a583..831b4c92e 100644 --- a/npbench/benchmarks/polybench/deriche/deriche_numpy.py +++ b/npbench/benchmarks/polybench/deriche/deriche_numpy.py @@ -13,14 +13,14 @@ def kernel(alpha, imgIn): b2 = -np.exp(-2.0 * alpha) c1 = c2 = 1 - y1 = np.empty_like(imgIn) + y1 = np.zeros_like(imgIn) y1[:, 0] = a1 * imgIn[:, 0] y1[:, 1] = a1 * imgIn[:, 1] + a2 * imgIn[:, 0] + b1 * y1[:, 0] for j in range(2, imgIn.shape[1]): y1[:, j] = (a1 * imgIn[:, j] + a2 * imgIn[:, j - 1] + b1 * y1[:, j - 1] + b2 * y1[:, j - 2]) - y2 = np.empty_like(imgIn) + y2 = np.zeros_like(imgIn) y2[:, -1] = 0.0 y2[:, -2] = a3 * imgIn[:, -1] for j in range(imgIn.shape[1] - 3, -1, -1): diff --git a/npbench/benchmarks/polybench/deriche/deriche_pythran.py b/npbench/benchmarks/polybench/deriche/deriche_pythran.py index 96e157e95..411b74b61 100644 --- a/npbench/benchmarks/polybench/deriche/deriche_pythran.py +++ b/npbench/benchmarks/polybench/deriche/deriche_pythran.py @@ -14,14 +14,14 @@ def kernel(alpha, imgIn): b2 = -np.exp(-2.0 * alpha) c1 = c2 = 1 - y1 = np.empty_like(imgIn) + y1 = np.zeros_like(imgIn) y1[:, 0] = a1 * imgIn[:, 0] y1[:, 1] = a1 * imgIn[:, 1] + a2 * imgIn[:, 0] + b1 * y1[:, 0] for j in range(2, imgIn.shape[1]): y1[:, j] = (a1 * imgIn[:, j] + a2 * imgIn[:, j - 1] + b1 * y1[:, j - 1] + b2 * y1[:, j - 2]) - y2 = np.empty_like(imgIn) + y2 = np.zeros_like(imgIn) y2[:, -1] = 0.0 y2[:, -2] = a3 * imgIn[:, -1] for j in range(imgIn.shape[1] - 3, -1, -1): diff --git a/npbench/benchmarks/polybench/doitgen/doitgen_legate.py b/npbench/benchmarks/polybench/doitgen/doitgen_legate.py index 5b5573780..102c810af 100644 --- a/npbench/benchmarks/polybench/doitgen/doitgen_legate.py +++ b/npbench/benchmarks/polybench/doitgen/doitgen_legate.py @@ -16,12 +16,12 @@ def kernel(NR, NQ, NP, A, C4): def init_data(NR, NQ, NP, datatype): - A = onp.empty((NR, NQ, NP), dtype=datatype) - C4 = onp.empty(( + A = onp.zeros((NR, NQ, NP), dtype=datatype) + C4 = onp.zeros(( NP, NP, ), dtype=datatype) - sum = onp.empty((NP, ), dtype=datatype) + sum = onp.zeros((NP, ), dtype=datatype) for i in range(NR): for j in range(NQ): for k in range(NP): diff --git a/npbench/benchmarks/polybench/durbin/durbin_cupy.py b/npbench/benchmarks/polybench/durbin/durbin_cupy.py index 1a84a9ed2..a2d256bb4 100644 --- a/npbench/benchmarks/polybench/durbin/durbin_cupy.py +++ b/npbench/benchmarks/polybench/durbin/durbin_cupy.py @@ -3,7 +3,7 @@ def kernel(r): - y = np.empty_like(r) + y = np.zeros_like(r) alpha = -r[0] beta = 1.0 y[0] = -r[0] diff --git a/npbench/benchmarks/polybench/durbin/durbin_dace.py b/npbench/benchmarks/polybench/durbin/durbin_dace.py index 213723b3d..bb42da42d 100644 --- a/npbench/benchmarks/polybench/durbin/durbin_dace.py +++ b/npbench/benchmarks/polybench/durbin/durbin_dace.py @@ -15,7 +15,7 @@ def flip(A: dc.float64[M]): @dc.program def kernel(r: dc.float64[N]): - y = np.empty_like(r) + y = np.zeros_like(r) alpha = -r[0] beta = 1.0 y[0] = -r[0] diff --git a/npbench/benchmarks/polybench/durbin/durbin_jax.py b/npbench/benchmarks/polybench/durbin/durbin_jax.py index 491ccacd4..8359ed7b0 100644 --- a/npbench/benchmarks/polybench/durbin/durbin_jax.py +++ b/npbench/benchmarks/polybench/durbin/durbin_jax.py @@ -6,7 +6,7 @@ @jax.jit def kernel(r): - y = jnp.empty_like(r) + y = jnp.zeros_like(r) alpha = -r[0] beta = 1.0 y = y.at[0].set(-r[0]) diff --git a/npbench/benchmarks/polybench/durbin/durbin_legate.py b/npbench/benchmarks/polybench/durbin/durbin_legate.py index 7a71cb8d7..710a9eaef 100644 --- a/npbench/benchmarks/polybench/durbin/durbin_legate.py +++ b/npbench/benchmarks/polybench/durbin/durbin_legate.py @@ -7,7 +7,7 @@ def kernel(r): - y = np.empty_like(r) + y = np.zeros_like(r) alpha = -r[0] beta = 1.0 y[0] = -r[0] @@ -25,8 +25,8 @@ def kernel(r): def init_data(N, datatype): - r = onp.empty((N, ), dtype=datatype) - y = onp.empty((N, ), dtype=datatype) + r = onp.zeros((N, ), dtype=datatype) + y = onp.zeros((N, ), dtype=datatype) for i in range(N): r[i] = N + 1 - i diff --git a/npbench/benchmarks/polybench/durbin/durbin_numba_n.py b/npbench/benchmarks/polybench/durbin/durbin_numba_n.py index b90702fb5..2a703cf0f 100644 --- a/npbench/benchmarks/polybench/durbin/durbin_numba_n.py +++ b/npbench/benchmarks/polybench/durbin/durbin_numba_n.py @@ -5,7 +5,7 @@ @nb.jit(nopython=True, parallel=False, fastmath=True) def kernel(r): - y = np.empty_like(r) + y = np.zeros_like(r) alpha = -r[0] beta = 1.0 y[0] = -r[0] diff --git a/npbench/benchmarks/polybench/durbin/durbin_numba_np.py b/npbench/benchmarks/polybench/durbin/durbin_numba_np.py index 0ab3d053b..97618dd16 100644 --- a/npbench/benchmarks/polybench/durbin/durbin_numba_np.py +++ b/npbench/benchmarks/polybench/durbin/durbin_numba_np.py @@ -5,7 +5,7 @@ @nb.jit(nopython=True, parallel=True, fastmath=True) def kernel(r): - y = np.empty_like(r) + y = np.zeros_like(r) alpha = -r[0] beta = 1.0 y[0] = -r[0] diff --git a/npbench/benchmarks/polybench/durbin/durbin_numpy.py b/npbench/benchmarks/polybench/durbin/durbin_numpy.py index cff58112d..4dc6253c0 100644 --- a/npbench/benchmarks/polybench/durbin/durbin_numpy.py +++ b/npbench/benchmarks/polybench/durbin/durbin_numpy.py @@ -3,7 +3,7 @@ def kernel(r): - y = np.empty_like(r) + y = np.zeros_like(r) alpha = -r[0] beta = 1.0 y[0] = -r[0] diff --git a/npbench/benchmarks/polybench/durbin/durbin_pythran.py b/npbench/benchmarks/polybench/durbin/durbin_pythran.py index e0bea63ea..cf793d63b 100644 --- a/npbench/benchmarks/polybench/durbin/durbin_pythran.py +++ b/npbench/benchmarks/polybench/durbin/durbin_pythran.py @@ -2,7 +2,7 @@ def flip(A): - B = np.empty_like(A) + B = np.zeros_like(A) for i in range(B.shape[0]): B[i] = A[-1 - i] return B @@ -11,7 +11,7 @@ def flip(A): # pythran export kernel(float64[:]) def kernel(r): - y = np.empty_like(r) + y = np.zeros_like(r) alpha = -r[0] beta = 1.0 y[0] = -r[0] diff --git a/npbench/benchmarks/polybench/gemm/gemm_legate.py b/npbench/benchmarks/polybench/gemm/gemm_legate.py index 0798f3bed..edece63f9 100644 --- a/npbench/benchmarks/polybench/gemm/gemm_legate.py +++ b/npbench/benchmarks/polybench/gemm/gemm_legate.py @@ -13,17 +13,17 @@ def init_data(NI, NJ, NK, datatype): alpha = datatype(1.5) beta = datatype(1.2) - C = np.empty((NI, NJ), dtype=datatype) + C = np.zeros((NI, NJ), dtype=datatype) # for i in range(NI): # for j in range(NJ): # C[i, j] = ((i * j + 1) % NI) / NI C[:] = np.random.randn(NI, NJ) - A = np.empty((NI, NK), dtype=datatype) + A = np.zeros((NI, NK), dtype=datatype) # for i in range(NI): # for k in range(NK): # A[i, k] = (i * (k + 1) % NK) / NK A[:] = np.random.randn(NI, NK) - B = np.empty((NK, NJ), dtype=datatype) + B = np.zeros((NK, NJ), dtype=datatype) # for k in range(NK): # for j in range(NJ): # C[i, j] = (k * (j + 2) % NJ) / NJ diff --git a/npbench/benchmarks/polybench/gramschmidt/gramschmidt.py b/npbench/benchmarks/polybench/gramschmidt/gramschmidt.py index c7423ca62..a01df9747 100644 --- a/npbench/benchmarks/polybench/gramschmidt/gramschmidt.py +++ b/npbench/benchmarks/polybench/gramschmidt/gramschmidt.py @@ -8,7 +8,11 @@ def initialize(M, N, datatype=np.float64): rng = default_rng(42) A = rng.random((M, N), dtype=datatype) - while np.linalg.matrix_rank(A) < N: - A = rng.random((M, N), dtype=datatype) + # Add a diagonal-dominance term so A is deterministically full column rank and + # well-conditioned. A plain random matrix is only full-rank in expectation (hence the + # former reject-sampling ``while matrix_rank(A) < N`` loop) and can still be poorly + # conditioned, which makes the Gram-Schmidt QR numerically unstable; the added diagonal + # gives cond(A) ~= 1.5 without the nondeterministic resampling. + A[:N, :N] += N * np.eye(N, dtype=datatype) return A diff --git a/npbench/benchmarks/polybench/k2mm/k2mm_legate.py b/npbench/benchmarks/polybench/k2mm/k2mm_legate.py index 1829872ab..0dbb206fa 100644 --- a/npbench/benchmarks/polybench/k2mm/k2mm_legate.py +++ b/npbench/benchmarks/polybench/k2mm/k2mm_legate.py @@ -11,11 +11,11 @@ def init_data(NI, NJ, NK, NL, datatype): alpha = datatype(1.5) beta = datatype(1.2) - tmp = np.empty((NI, NJ), dtype=datatype) - A = np.empty((NI, NK), dtype=datatype) - B = np.empty((NK, NJ), dtype=datatype) - C = np.empty((NJ, NL), dtype=datatype) - D = np.empty((NI, NL), dtype=datatype) + tmp = np.zeros((NI, NJ), dtype=datatype) + A = np.zeros((NI, NK), dtype=datatype) + B = np.zeros((NK, NJ), dtype=datatype) + C = np.zeros((NJ, NL), dtype=datatype) + D = np.zeros((NI, NL), dtype=datatype) # for i in range(NI): # for j in range(NK): # A[i, j] = ((i * j + 1) % NI) / NI diff --git a/npbench/benchmarks/polybench/k3mm/k3mm_legate.py b/npbench/benchmarks/polybench/k3mm/k3mm_legate.py index 5d7709fbb..659c6d39d 100644 --- a/npbench/benchmarks/polybench/k3mm/k3mm_legate.py +++ b/npbench/benchmarks/polybench/k3mm/k3mm_legate.py @@ -47,13 +47,13 @@ def kernel(A, B, C, D): def init_data(NI, NJ, NK, NL, NM, datatype): - E = np.empty((NI, NJ), dtype=datatype) - A = np.empty((NI, NK), dtype=datatype) - B = np.empty((NK, NJ), dtype=datatype) - F = np.empty((NJ, NL), dtype=datatype) - C = np.empty((NJ, NM), dtype=datatype) - D = np.empty((NM, NL), dtype=datatype) - G = np.empty((NI, NL), dtype=datatype) + E = np.zeros((NI, NJ), dtype=datatype) + A = np.zeros((NI, NK), dtype=datatype) + B = np.zeros((NK, NJ), dtype=datatype) + F = np.zeros((NJ, NL), dtype=datatype) + C = np.zeros((NJ, NM), dtype=datatype) + D = np.zeros((NM, NL), dtype=datatype) + G = np.zeros((NI, NL), dtype=datatype) # for i in range(NI): # for j in range(NK): # A[i, j] = ((i * j + 1) % NI) / (5 * NI) diff --git a/npbench/benchmarks/polybench/lu/lu.py b/npbench/benchmarks/polybench/lu/lu.py index 1fe67bd77..969ae8706 100644 --- a/npbench/benchmarks/polybench/lu/lu.py +++ b/npbench/benchmarks/polybench/lu/lu.py @@ -4,7 +4,7 @@ def initialize(N, datatype=np.float64): - A = np.empty((N, N), dtype=datatype) + A = np.zeros((N, N), dtype=datatype) for i in range(N): A[i, :i + 1] = np.fromfunction(lambda j: (-j % N) / N + 1, (i + 1, ), dtype=datatype) diff --git a/npbench/benchmarks/polybench/ludcmp/ludcmp.py b/npbench/benchmarks/polybench/ludcmp/ludcmp.py index ef9b96791..1df052c5f 100644 --- a/npbench/benchmarks/polybench/ludcmp/ludcmp.py +++ b/npbench/benchmarks/polybench/ludcmp/ludcmp.py @@ -4,7 +4,7 @@ def initialize(N, datatype=np.float64): - A = np.empty((N, N), dtype=datatype) + A = np.zeros((N, N), dtype=datatype) for i in range(N): A[i, :i + 1] = np.fromfunction(lambda j: (-j % N) / N + 1, (i + 1, ), dtype=datatype) diff --git a/npbench/benchmarks/polybench/symm/symm.py b/npbench/benchmarks/polybench/symm/symm.py index a2de70551..4d4dfbaf3 100644 --- a/npbench/benchmarks/polybench/symm/symm.py +++ b/npbench/benchmarks/polybench/symm/symm.py @@ -10,7 +10,7 @@ def initialize(M, N, datatype=np.float64): dtype=datatype) B = np.fromfunction(lambda i, j: ((N + i - j) % 100) / M, (M, N), dtype=datatype) - A = np.empty((M, M), dtype=datatype) + A = np.zeros((M, M), dtype=datatype) for i in range(M): A[i, :i + 1] = np.fromfunction(lambda j: ((i + j) % 100) / M, (i + 1, ), diff --git a/npbench/benchmarks/polybench/symm/symm_cupy.py b/npbench/benchmarks/polybench/symm/symm_cupy.py index f3cc3c3e0..a9067d52a 100644 --- a/npbench/benchmarks/polybench/symm/symm_cupy.py +++ b/npbench/benchmarks/polybench/symm/symm_cupy.py @@ -3,7 +3,7 @@ def kernel(alpha, beta, C, A, B): - temp2 = np.empty((C.shape[1], ), dtype=C.dtype) + temp2 = np.zeros((C.shape[1], ), dtype=C.dtype) C *= beta for i in range(C.shape[0]): for j in range(C.shape[1]): diff --git a/npbench/benchmarks/polybench/symm/symm_dace.py b/npbench/benchmarks/polybench/symm/symm_dace.py index fa970f9da..fca79850e 100644 --- a/npbench/benchmarks/polybench/symm/symm_dace.py +++ b/npbench/benchmarks/polybench/symm/symm_dace.py @@ -8,7 +8,7 @@ def kernel(alpha: dc.float64, beta: dc.float64, C: dc.float64[M, N], A: dc.float64[M, M], B: dc.float64[M, N]): - temp2 = np.empty((N, ), dtype=C.dtype) + temp2 = np.zeros((N, ), dtype=C.dtype) C *= beta for i in range(M): for j in range(N): diff --git a/npbench/benchmarks/polybench/symm/symm_jax.py b/npbench/benchmarks/polybench/symm/symm_jax.py index abeafceea..a8c75f003 100644 --- a/npbench/benchmarks/polybench/symm/symm_jax.py +++ b/npbench/benchmarks/polybench/symm/symm_jax.py @@ -5,7 +5,7 @@ @jax.jit def kernel(alpha, beta, C: jax.Array, A: jax.Array, B: jax.Array): - temp2 = jnp.empty((C.shape[1], ), dtype=C.dtype) + temp2 = jnp.zeros((C.shape[1], ), dtype=C.dtype) C *= beta def row_update(i, arrays): diff --git a/npbench/benchmarks/polybench/symm/symm_numba_n.py b/npbench/benchmarks/polybench/symm/symm_numba_n.py index 0ad1265be..84cde5609 100644 --- a/npbench/benchmarks/polybench/symm/symm_numba_n.py +++ b/npbench/benchmarks/polybench/symm/symm_numba_n.py @@ -5,7 +5,7 @@ @nb.jit(nopython=True, parallel=False, fastmath=True) def kernel(alpha, beta, C, A, B): - temp2 = np.empty((C.shape[1], ), dtype=C.dtype) + temp2 = np.zeros((C.shape[1], ), dtype=C.dtype) C *= beta for i in range(C.shape[0]): for j in range(C.shape[1]): diff --git a/npbench/benchmarks/polybench/symm/symm_numba_npr.py b/npbench/benchmarks/polybench/symm/symm_numba_npr.py index 79fd006b6..2aabfb006 100644 --- a/npbench/benchmarks/polybench/symm/symm_numba_npr.py +++ b/npbench/benchmarks/polybench/symm/symm_numba_npr.py @@ -5,7 +5,7 @@ @nb.jit(nopython=True, parallel=True, fastmath=True) def kernel(alpha, beta, C, A, B): - temp2 = np.empty((C.shape[1], ), dtype=C.dtype) + temp2 = np.zeros((C.shape[1], ), dtype=C.dtype) C *= beta for i in range(C.shape[0]): for j in nb.prange(C.shape[1]): diff --git a/npbench/benchmarks/polybench/symm/symm_numpy.py b/npbench/benchmarks/polybench/symm/symm_numpy.py index bfb653b27..9ac391386 100644 --- a/npbench/benchmarks/polybench/symm/symm_numpy.py +++ b/npbench/benchmarks/polybench/symm/symm_numpy.py @@ -3,7 +3,7 @@ def kernel(alpha, beta, C, A, B): - temp2 = np.empty((C.shape[1], ), dtype=C.dtype) + temp2 = np.zeros((C.shape[1], ), dtype=C.dtype) C *= beta for i in range(C.shape[0]): for j in range(C.shape[1]): diff --git a/npbench/benchmarks/polybench/symm/symm_pythran.py b/npbench/benchmarks/polybench/symm/symm_pythran.py index 6272ebdf0..b3e3e50d6 100644 --- a/npbench/benchmarks/polybench/symm/symm_pythran.py +++ b/npbench/benchmarks/polybench/symm/symm_pythran.py @@ -4,7 +4,7 @@ # pythran export kernel(float64, float64, float64[:,:], float64[:,:], float64[:,:]) def kernel(alpha, beta, C, A, B): - temp2 = np.empty((C.shape[1], ), dtype=C.dtype) + temp2 = np.zeros((C.shape[1], ), dtype=C.dtype) C *= beta for i in range(C.shape[0]): for j in range(C.shape[1]): diff --git a/npbench/benchmarks/spmv/spmv_cupy.py b/npbench/benchmarks/spmv/spmv_cupy.py index 39666b47a..1167ef0a5 100644 --- a/npbench/benchmarks/spmv/spmv_cupy.py +++ b/npbench/benchmarks/spmv/spmv_cupy.py @@ -5,7 +5,7 @@ # Matrix-Vector Multiplication with the matrix given in Compressed Sparse Row # (CSR) format def spmv(A_row, A_col, A_val, x): - y = np.empty(A_row.size - 1, A_val.dtype) + y = np.zeros(A_row.size - 1, A_val.dtype) for i in range(A_row.size - 1): cols = A_col[A_row[i]:A_row[i + 1]] diff --git a/npbench/benchmarks/spmv/spmv_dace.py b/npbench/benchmarks/spmv/spmv_dace.py index 53f8e0334..47cc6e81e 100644 --- a/npbench/benchmarks/spmv/spmv_dace.py +++ b/npbench/benchmarks/spmv/spmv_dace.py @@ -11,7 +11,7 @@ def spmv(A_row: dc.uint32[M + 1], A_col: dc.uint32[nnz], A_val: dc.float64[nnz], x: dc.float64[N]): # y = np.empty(A_row.size - 1, A_val.dtype) - y = np.empty(M, A_val.dtype) + y = np.zeros(M, A_val.dtype) # for i in range(A_row.size - 1): for i in range(M): diff --git a/npbench/benchmarks/spmv/spmv_jax.py b/npbench/benchmarks/spmv/spmv_jax.py index 25cc00874..fcdd712ac 100644 --- a/npbench/benchmarks/spmv/spmv_jax.py +++ b/npbench/benchmarks/spmv/spmv_jax.py @@ -7,7 +7,7 @@ # (CSR) format @jax.jit def spmv(A_row, A_col, A_val, x): - y = jnp.empty(A_row.size - 1, dtype=A_val.dtype) + y = jnp.zeros(A_row.size - 1, dtype=A_val.dtype) def row_update(i, y): diff --git a/npbench/benchmarks/spmv/spmv_numba_n.py b/npbench/benchmarks/spmv/spmv_numba_n.py index fcaf4b902..d6b6ccb21 100644 --- a/npbench/benchmarks/spmv/spmv_numba_n.py +++ b/npbench/benchmarks/spmv/spmv_numba_n.py @@ -7,7 +7,7 @@ # (CSR) format @nb.jit(nopython=True, parallel=False, fastmath=True) def spmv(A_row, A_col, A_val, x): - y = np.empty(A_row.size - 1, A_val.dtype) + y = np.zeros(A_row.size - 1, A_val.dtype) for i in range(A_row.size - 1): cols = A_col[A_row[i]:A_row[i + 1]] diff --git a/npbench/benchmarks/spmv/spmv_numba_np.py b/npbench/benchmarks/spmv/spmv_numba_np.py index 50e45b2f0..ba93ebbeb 100644 --- a/npbench/benchmarks/spmv/spmv_numba_np.py +++ b/npbench/benchmarks/spmv/spmv_numba_np.py @@ -7,7 +7,7 @@ # (CSR) format @nb.jit(nopython=True, parallel=True, fastmath=True) def spmv(A_row, A_col, A_val, x): - y = np.empty(A_row.size - 1, A_val.dtype) + y = np.zeros(A_row.size - 1, A_val.dtype) for i in range(A_row.size - 1): cols = A_col[A_row[i]:A_row[i + 1]] diff --git a/npbench/benchmarks/spmv/spmv_numba_npr.py b/npbench/benchmarks/spmv/spmv_numba_npr.py index 3ed06204d..2784bc569 100644 --- a/npbench/benchmarks/spmv/spmv_numba_npr.py +++ b/npbench/benchmarks/spmv/spmv_numba_npr.py @@ -7,7 +7,7 @@ # (CSR) format @nb.jit(nopython=True, parallel=True, fastmath=True) def spmv(A_row, A_col, A_val, x): - y = np.empty(A_row.size - 1, A_val.dtype) + y = np.zeros(A_row.size - 1, A_val.dtype) for i in nb.prange(A_row.size - 1): cols = A_col[A_row[i]:A_row[i + 1]] diff --git a/npbench/benchmarks/spmv/spmv_numpy.py b/npbench/benchmarks/spmv/spmv_numpy.py index f51a2836b..ed4518e61 100644 --- a/npbench/benchmarks/spmv/spmv_numpy.py +++ b/npbench/benchmarks/spmv/spmv_numpy.py @@ -5,7 +5,7 @@ # Matrix-Vector Multiplication with the matrix given in Compressed Sparse Row # (CSR) format def spmv(A_row, A_col, A_val, x): - y = np.empty(A_row.size - 1, A_val.dtype) + y = np.zeros(A_row.size - 1, A_val.dtype) for i in range(A_row.size - 1): cols = A_col[A_row[i]:A_row[i + 1]] diff --git a/npbench/benchmarks/spmv/spmv_pythran.py b/npbench/benchmarks/spmv/spmv_pythran.py index d27417d83..0a17ce1e5 100644 --- a/npbench/benchmarks/spmv/spmv_pythran.py +++ b/npbench/benchmarks/spmv/spmv_pythran.py @@ -6,7 +6,7 @@ # (CSR) format # pythran export spmv(uint32[:], uint32[:], float64[:], float64[:]) def spmv(A_row, A_col, A_val, x): - y = np.empty(A_row.size - 1, A_val.dtype) + y = np.zeros(A_row.size - 1, A_val.dtype) for i in range(A_row.size - 1): cols = A_col[A_row[i]:A_row[i + 1]] diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_cupy.py b/npbench/benchmarks/stockham_fft/stockham_fft_cupy.py index f4978dc1c..94716acfd 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_cupy.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_cupy.py @@ -6,7 +6,7 @@ def stockham_fft(N, R, K, x, y): # Generate DFT matrix for radix R. # Define transient variable for matrix. i_coord, j_coord = np.mgrid[0:R, 0:R] - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y # to avoid overwriting the input. @@ -21,7 +21,7 @@ def stockham_fft(N, R, K, x, y): yv = np.reshape(y, (R**i, R, R**(K - i - 1))) tmp_perm = np.transpose(yv, axes=(1, 0, 2)) # Twiddle Factor multiplication - D = np.empty((R, R**i, R**(K - i - 1)), dtype=np.complex128) + D = np.zeros((R, R**i, R**(K - i - 1)), dtype=np.complex128) tmp = np.exp(-2.0j * np.pi * ii_coord[:, :R**i] * jj_coord[:, :R**i] / R**(i + 1)) D[:] = np.repeat(np.reshape(tmp, (R, R**i, 1)), R**(K - i - 1), axis=2) diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_dace.py b/npbench/benchmarks/stockham_fft/stockham_fft_dace.py index 2acc20e98..f33e3a3e1 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_dace.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_dace.py @@ -35,7 +35,7 @@ def stockham_fft(x: dc.complex128[R**K], y: dc.complex128[R**K]): i_coord = np.ndarray((R, R), dtype=np.uint32) j_coord = np.ndarray((R, R), dtype=np.uint32) mgrid1(i_coord, j_coord) - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat[:] = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y # to avoid overwriting the input. @@ -46,9 +46,9 @@ def stockham_fft(x: dc.complex128[R**K], y: dc.complex128[R**K]): jj_coord = np.ndarray((R, N), dtype=np.uint32) mgrid2(ii_coord, jj_coord) - tmp_perm = np.empty_like(y) - D = np.empty_like(y) - tmp = np.empty_like(y) + tmp_perm = np.zeros_like(y) + D = np.zeros_like(y) + tmp = np.zeros_like(y) # Main Stockham loop for i in range(K): diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_dpnp.py b/npbench/benchmarks/stockham_fft/stockham_fft_dpnp.py index 502c7dfac..086669412 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_dpnp.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_dpnp.py @@ -5,7 +5,7 @@ def stockham_fft(N, R, K, x, y): # Generate DFT matrix for radix R. # Define transient variable for matrix. i_coord, j_coord = np.mgrid[0:R, 0:R] - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y to avoid overwriting the input. # y[:] = x[:] @@ -19,7 +19,7 @@ def stockham_fft(N, R, K, x, y): yv = np.reshape(y, (R**i, R, R**(K - i - 1))) tmp_perm = np.transpose(yv, axes=(1, 0, 2)) # Twiddle Factor multiplication - D = np.empty((R, R**i, R**(K - i - 1)), dtype=np.complex128) + D = np.zeros((R, R**i, R**(K - i - 1)), dtype=np.complex128) tmp = np.exp(-2.0j * np.pi * ii_coord[:, :R**i] * jj_coord[:, :R**i] / R**(i + 1)) np.copyto(D, np.repeat(np.reshape(tmp, (R, R**i, 1)), R**(K - i - 1), axis=2)) diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_numba_n.py b/npbench/benchmarks/stockham_fft/stockham_fft_numba_n.py index 261ca2c97..09fc63b9b 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_numba_n.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_numba_n.py @@ -4,8 +4,8 @@ @nb.jit(nopython=True, parallel=False, fastmath=True) def mgrid(xn, yn): - Xi = np.empty((xn, yn), dtype=np.uint32) - Yi = np.empty((xn, yn), dtype=np.uint32) + Xi = np.zeros((xn, yn), dtype=np.uint32) + Yi = np.zeros((xn, yn), dtype=np.uint32) for i in range(xn): Xi[i, :] = i for j in range(yn): @@ -19,7 +19,7 @@ def stockham_fft(N, R, K, x, y): # Generate DFT matrix for radix R. # Define transient variable for matrix. i_coord, j_coord = mgrid(R, R) - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y # to avoid overwriting the input. @@ -36,7 +36,7 @@ def stockham_fft(N, R, K, x, y): # tmp_perm = np.transpose(yv, axes=(1, 0, 2)) tmp_perm = np.transpose(yv, axes=(1, 0, 2)).copy() # Twiddle Factor multiplication - D = np.empty((R, R**i, R**(K - i - 1)), dtype=np.complex128) + D = np.zeros((R, R**i, R**(K - i - 1)), dtype=np.complex128) tmp = np.exp(-2.0j * np.pi * ii_coord[:, :R**i] * jj_coord[:, :R**i] / R**(i + 1)) # D[:] = np.repeat(np.reshape(tmp, (R, R**i, 1)), R ** (K-i-1), axis=2) diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_numba_np.py b/npbench/benchmarks/stockham_fft/stockham_fft_numba_np.py index bdbe067ca..be7af6747 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_numba_np.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_numba_np.py @@ -4,8 +4,8 @@ @nb.jit(nopython=True, parallel=True, fastmath=True) def mgrid(xn, yn): - Xi = np.empty((xn, yn), dtype=np.uint32) - Yi = np.empty((xn, yn), dtype=np.uint32) + Xi = np.zeros((xn, yn), dtype=np.uint32) + Yi = np.zeros((xn, yn), dtype=np.uint32) for i in range(xn): Xi[i, :] = i for j in range(yn): @@ -19,7 +19,7 @@ def stockham_fft(N, R, K, x, y): # Generate DFT matrix for radix R. # Define transient variable for matrix. i_coord, j_coord = mgrid(R, R) - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y # to avoid overwriting the input. @@ -36,7 +36,7 @@ def stockham_fft(N, R, K, x, y): # tmp_perm = np.transpose(yv, axes=(1, 0, 2)) tmp_perm = np.transpose(yv, axes=(1, 0, 2)).copy() # Twiddle Factor multiplication - D = np.empty((R, R**i, R**(K - i - 1)), dtype=np.complex128) + D = np.zeros((R, R**i, R**(K - i - 1)), dtype=np.complex128) tmp = np.exp(-2.0j * np.pi * ii_coord[:, :R**i] * jj_coord[:, :R**i] / R**(i + 1)) # D[:] = np.repeat(np.reshape(tmp, (R, R**i, 1)), R ** (K-i-1), axis=2) diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_numba_npr.py b/npbench/benchmarks/stockham_fft/stockham_fft_numba_npr.py index c367593e2..24cc581a4 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_numba_npr.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_numba_npr.py @@ -4,8 +4,8 @@ @nb.jit(nopython=True, parallel=True, fastmath=True) def mgrid(xn, yn): - Xi = np.empty((xn, yn), dtype=np.uint32) - Yi = np.empty((xn, yn), dtype=np.uint32) + Xi = np.zeros((xn, yn), dtype=np.uint32) + Yi = np.zeros((xn, yn), dtype=np.uint32) for i in nb.prange(xn): Xi[i, :] = i for j in nb.prange(yn): @@ -19,7 +19,7 @@ def stockham_fft(N, R, K, x, y): # Generate DFT matrix for radix R. # Define transient variable for matrix. i_coord, j_coord = mgrid(R, R) - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y # to avoid overwriting the input. @@ -36,7 +36,7 @@ def stockham_fft(N, R, K, x, y): # tmp_perm = np.transpose(yv, axes=(1, 0, 2)) tmp_perm = np.transpose(yv, axes=(1, 0, 2)).copy() # Twiddle Factor multiplication - D = np.empty((R, R**i, R**(K - i - 1)), dtype=np.complex128) + D = np.zeros((R, R**i, R**(K - i - 1)), dtype=np.complex128) tmp = np.exp(-2.0j * np.pi * ii_coord[:, :R**i] * jj_coord[:, :R**i] / R**(i + 1)) # D[:] = np.repeat(np.reshape(tmp, (R, R**i, 1)), R ** (K-i-1), axis=2) diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_numba_o.py b/npbench/benchmarks/stockham_fft/stockham_fft_numba_o.py index 6b41deeb7..7f36085a4 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_numba_o.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_numba_o.py @@ -8,7 +8,7 @@ def stockham_fft(N, R, K, x, y): # Generate DFT matrix for radix R. # Define transient variable for matrix. i_coord, j_coord = np.mgrid[0:R, 0:R] - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y # to avoid overwriting the input. @@ -23,7 +23,7 @@ def stockham_fft(N, R, K, x, y): yv = np.reshape(y, (R**i, R, R**(K - i - 1))) tmp_perm = np.transpose(yv, axes=(1, 0, 2)) # Twiddle Factor multiplication - D = np.empty((R, R**i, R**(K - i - 1)), dtype=np.complex128) + D = np.zeros((R, R**i, R**(K - i - 1)), dtype=np.complex128) tmp = np.exp(-2.0j * np.pi * ii_coord[:, :R**i] * jj_coord[:, :R**i] / R**(i + 1)) D[:] = np.repeat(np.reshape(tmp, (R, R**i, 1)), R**(K - i - 1), axis=2) diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_numba_op.py b/npbench/benchmarks/stockham_fft/stockham_fft_numba_op.py index 248c6d044..518dc15dc 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_numba_op.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_numba_op.py @@ -8,7 +8,7 @@ def stockham_fft(N, R, K, x, y): # Generate DFT matrix for radix R. # Define transient variable for matrix. i_coord, j_coord = np.mgrid[0:R, 0:R] - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y # to avoid overwriting the input. @@ -23,7 +23,7 @@ def stockham_fft(N, R, K, x, y): yv = np.reshape(y, (R**i, R, R**(K - i - 1))) tmp_perm = np.transpose(yv, axes=(1, 0, 2)) # Twiddle Factor multiplication - D = np.empty((R, R**i, R**(K - i - 1)), dtype=np.complex128) + D = np.zeros((R, R**i, R**(K - i - 1)), dtype=np.complex128) tmp = np.exp(-2.0j * np.pi * ii_coord[:, :R**i] * jj_coord[:, :R**i] / R**(i + 1)) D[:] = np.repeat(np.reshape(tmp, (R, R**i, 1)), R**(K - i - 1), axis=2) diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_numpy.py b/npbench/benchmarks/stockham_fft/stockham_fft_numpy.py index d122e4029..6244bea3c 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_numpy.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_numpy.py @@ -6,7 +6,7 @@ def stockham_fft(N, R, K, x, y): # Generate DFT matrix for radix R. # Define transient variable for matrix. i_coord, j_coord = np.mgrid[0:R, 0:R] - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y # to avoid overwriting the input. @@ -21,7 +21,7 @@ def stockham_fft(N, R, K, x, y): yv = np.reshape(y, (R**i, R, R**(K - i - 1))) tmp_perm = np.transpose(yv, axes=(1, 0, 2)) # Twiddle Factor multiplication - D = np.empty((R, R**i, R**(K - i - 1)), dtype=np.complex128) + D = np.zeros((R, R**i, R**(K - i - 1)), dtype=np.complex128) tmp = np.exp(-2.0j * np.pi * ii_coord[:, :R**i] * jj_coord[:, :R**i] / R**(i + 1)) D[:] = np.repeat(np.reshape(tmp, (R, R**i, 1)), R**(K - i - 1), axis=2) diff --git a/npbench/benchmarks/stockham_fft/stockham_fft_pythran.py b/npbench/benchmarks/stockham_fft/stockham_fft_pythran.py index 164088714..e822ffb00 100644 --- a/npbench/benchmarks/stockham_fft/stockham_fft_pythran.py +++ b/npbench/benchmarks/stockham_fft/stockham_fft_pythran.py @@ -3,8 +3,8 @@ # pythran export mgrid(int, int) def mgrid(xn, yn): - Xi = np.empty((xn, yn), dtype=np.uint32) - Yi = np.empty((xn, yn), dtype=np.uint32) + Xi = np.zeros((xn, yn), dtype=np.uint32) + Yi = np.zeros((xn, yn), dtype=np.uint32) for i in range(xn): Xi[i, :] = i for j in range(yn): @@ -19,7 +19,7 @@ def stockham_fft(N, R, K, x, y): # Define transient variable for matrix. # i_coord, j_coord = np.mgrid[0:R, 0:R] i_coord, j_coord = mgrid(R, R) - dft_mat = np.empty((R, R), dtype=np.complex128) + dft_mat = np.zeros((R, R), dtype=np.complex128) dft_mat = np.exp(-2.0j * np.pi * i_coord * j_coord / R) # Move input x to output y # to avoid overwriting the input. @@ -36,7 +36,7 @@ def stockham_fft(N, R, K, x, y): yv = y.reshape(R**i, R, R**(K - i - 1)) tmp_perm = np.transpose(yv, axes=(1, 0, 2)) # Twiddle Factor multiplication - D = np.empty((R, R**i, R**(K - i - 1)), dtype=np.complex128) + D = np.zeros((R, R**i, R**(K - i - 1)), dtype=np.complex128) tmp = np.exp(-2.0j * np.pi * ii_coord[:, :R**i] * jj_coord[:, :R**i] / R**(i + 1)) # D[:] = np.repeat(np.reshape(tmp, (R, R**i, 1)), R ** (K-i-1), axis=2) diff --git a/npbench/benchmarks/weather_stencils/vadv/vadv_dpnp.py b/npbench/benchmarks/weather_stencils/vadv/vadv_dpnp.py index 7feacce8d..073653101 100644 --- a/npbench/benchmarks/weather_stencils/vadv/vadv_dpnp.py +++ b/npbench/benchmarks/weather_stencils/vadv/vadv_dpnp.py @@ -7,9 +7,9 @@ # Adapted and optimized for runtime efficiency def vadv(utens_stage, u_stage, wcon, u_pos, utens, dtr_stage): I, J, K = utens_stage.shape[0], utens_stage.shape[1], utens_stage.shape[2] - ccol = np.empty((I, J, K), dtype=utens_stage.dtype) - dcol = np.empty((I, J, K), dtype=utens_stage.dtype) - data_col = np.empty((I, J), dtype=utens_stage.dtype) + ccol = np.zeros((I, J, K), dtype=utens_stage.dtype) + dcol = np.zeros((I, J, K), dtype=utens_stage.dtype) + data_col = np.zeros((I, J), dtype=utens_stage.dtype) # Use np.sum for slicing and cumulative operations for k in range(1): diff --git a/npbench/benchmarks/weather_stencils/vadv/vadv_jax.py b/npbench/benchmarks/weather_stencils/vadv/vadv_jax.py index a6d6a19cf..1f8108a77 100644 --- a/npbench/benchmarks/weather_stencils/vadv/vadv_jax.py +++ b/npbench/benchmarks/weather_stencils/vadv/vadv_jax.py @@ -10,9 +10,9 @@ @jax.jit def vadv(utens_stage, u_stage, wcon, u_pos, utens, dtr_stage): I, J, K = utens_stage.shape[0], utens_stage.shape[1], utens_stage.shape[2] - ccol = jnp.empty((I, J, K), dtype=utens_stage.dtype) - dcol = jnp.empty((I, J, K), dtype=utens_stage.dtype) - data_col = jnp.empty((I, J), dtype=utens_stage.dtype) + ccol = jnp.zeros((I, J, K), dtype=utens_stage.dtype) + dcol = jnp.zeros((I, J, K), dtype=utens_stage.dtype) + data_col = jnp.zeros((I, J), dtype=utens_stage.dtype) def loop1(k, loop_vars): ccol, dcol = loop_vars diff --git a/npbench/benchmarks/weather_stencils/vadv/vadv_numba_n.py b/npbench/benchmarks/weather_stencils/vadv/vadv_numba_n.py index fca378939..7aa3615d1 100644 --- a/npbench/benchmarks/weather_stencils/vadv/vadv_numba_n.py +++ b/npbench/benchmarks/weather_stencils/vadv/vadv_numba_n.py @@ -13,9 +13,9 @@ def vadv(utens_stage, u_stage, wcon, u_pos, utens, dtr_stage): # ccol = np.ndarray((I, J, K), dtype=utens_stage.dtype) # dcol = np.ndarray((I, J, K), dtype=utens_stage.dtype) # data_col = np.ndarray((I, J), dtype=utens_stage.dtype) - ccol = np.empty((I, J, K), dtype=utens_stage.dtype) - dcol = np.empty((I, J, K), dtype=utens_stage.dtype) - data_col = np.empty((I, J), dtype=utens_stage.dtype) + ccol = np.zeros((I, J, K), dtype=utens_stage.dtype) + dcol = np.zeros((I, J, K), dtype=utens_stage.dtype) + data_col = np.zeros((I, J), dtype=utens_stage.dtype) for k in range(1): gcv = 0.25 * (wcon[1:, :, k + 1] + wcon[:-1, :, k + 1]) diff --git a/npbench/benchmarks/weather_stencils/vadv/vadv_numba_np.py b/npbench/benchmarks/weather_stencils/vadv/vadv_numba_np.py index 1774e6f5e..f589d01e8 100644 --- a/npbench/benchmarks/weather_stencils/vadv/vadv_numba_np.py +++ b/npbench/benchmarks/weather_stencils/vadv/vadv_numba_np.py @@ -13,9 +13,9 @@ def vadv(utens_stage, u_stage, wcon, u_pos, utens, dtr_stage): # ccol = np.ndarray((I, J, K), dtype=utens_stage.dtype) # dcol = np.ndarray((I, J, K), dtype=utens_stage.dtype) # data_col = np.ndarray((I, J), dtype=utens_stage.dtype) - ccol = np.empty((I, J, K), dtype=utens_stage.dtype) - dcol = np.empty((I, J, K), dtype=utens_stage.dtype) - data_col = np.empty((I, J), dtype=utens_stage.dtype) + ccol = np.zeros((I, J, K), dtype=utens_stage.dtype) + dcol = np.zeros((I, J, K), dtype=utens_stage.dtype) + data_col = np.zeros((I, J), dtype=utens_stage.dtype) for k in range(1): gcv = 0.25 * (wcon[1:, :, k + 1] + wcon[:-1, :, k + 1]) diff --git a/npbench/benchmarks/weather_stencils/vadv/vadv_pythran.py b/npbench/benchmarks/weather_stencils/vadv/vadv_pythran.py index 709bc3aa5..28985e53c 100644 --- a/npbench/benchmarks/weather_stencils/vadv/vadv_pythran.py +++ b/npbench/benchmarks/weather_stencils/vadv/vadv_pythran.py @@ -13,9 +13,9 @@ def vadv(utens_stage, u_stage, wcon, u_pos, utens, dtr_stage): # ccol = np.ndarray((I, J, K), dtype=utens_stage.dtype) # dcol = np.ndarray((I, J, K), dtype=utens_stage.dtype) # data_col = np.ndarray((I, J), dtype=utens_stage.dtype) - ccol = np.empty((I, J, K), dtype=utens_stage.dtype) - dcol = np.empty((I, J, K), dtype=utens_stage.dtype) - data_col = np.empty((I, J), dtype=utens_stage.dtype) + ccol = np.zeros((I, J, K), dtype=utens_stage.dtype) + dcol = np.zeros((I, J, K), dtype=utens_stage.dtype) + data_col = np.zeros((I, J), dtype=utens_stage.dtype) for k in range(1): gcv = 0.25 * (wcon[1:, :, k + 1] + wcon[:-1, :, k + 1]) diff --git a/npbench/infrastructure/__init__.py b/npbench/infrastructure/__init__.py index 65083a2e3..42ea568b8 100644 --- a/npbench/infrastructure/__init__.py +++ b/npbench/infrastructure/__init__.py @@ -7,6 +7,7 @@ from .cupy_framework import * from .dace_framework import * +from .dace_canonicalize_framework import * from .legate_framework import * from .numba_framework import * from .pythran_framework import * diff --git a/npbench/infrastructure/dace_canonicalize_framework.py b/npbench/infrastructure/dace_canonicalize_framework.py new file mode 100644 index 000000000..89b59efb2 --- /dev/null +++ b/npbench/infrastructure/dace_canonicalize_framework.py @@ -0,0 +1,87 @@ +# Copyright 2021 ETH Zurich and the NPBench authors. All rights reserved. +import importlib +import traceback + +import numpy as np + +from npbench.infrastructure import Benchmark +from npbench.infrastructure.dace_framework import DaceFramework +from typing import Callable, Sequence, Tuple + + +class DaceCanonicalizeFramework(DaceFramework): + """Runs DaCe's ``canonicalize`` pipeline instead of ``auto_optimize``. + + Reuses the ``*_dace.py`` benchmark implementations (``postfix`` ``dace``) and every + calling convention of :class:`DaceFramework`; only the SDFG-optimization step differs. + The CPU and GPU variants share this class and differ only by the ``arch`` field in + their framework JSON, exactly as ``dace_cpu`` / ``dace_gpu`` share :class:`DaceFramework`. + + Contrast with ``dace_cpu`` / ``dace_gpu`` (WCR-config OFF): the canonicalize pipeline + leaves ``sdfg.openmp_array_reductions = True`` (WCR-config ON), so a whole-buffer WCR + accumulator of a parallel map lowers to an OpenMP ``reduction(op:A[0:n])`` array-section + clause instead of per-element atomics. + """ + + def copy_func(self) -> Callable: + """Copy inputs to the device: cupy for the GPU variant, ``np.copy`` otherwise.""" + if self.info["arch"] == "gpu": + import cupy + + def cp_copy_func(arr): + darr = cupy.asarray(arr) + cupy.cuda.stream.get_current_stream().synchronize() + return darr + + return cp_copy_func + return np.copy + + def implementations(self, bench: Benchmark) -> Sequence[Tuple[Callable, str]]: + """Build the SDFG for ``bench``, run the canonicalize pipeline + target + finalization, turn the WCR-config ON, then compile.""" + import dace # noqa: F401 + from dace.transformation.passes.canonicalize import canonicalize + from dace.transformation.passes.canonicalize.finalize import finalize_for_target + + module_pypath = "npbench.benchmarks.{r}.{m}".format(r=bench.info["relative_path"].replace('/', '.'), + m=bench.info["module_name"]) + postfix = self.info["postfix"] if "postfix" in self.info.keys() else self.fname + module_str = "{m}_{p}".format(m=module_pypath, p=postfix) + func_str = bench.info["func_name"] + + try: + module = importlib.import_module(module_str) + ct_impl = vars(module)[func_str] + except Exception as e: + print("Failed to load the DaCe implementation.") + raise e + + target = "gpu" if self.info["arch"] == "gpu" else "cpu" + + sdfg = ct_impl.to_sdfg(simplify=True) + sdfg._name = "canonicalize" + + canonicalize(sdfg, + validate=True, + target=target, + peel_limit=4, + break_anti_dependence=True, + interchange_carry_with_map=True, + scatter_to_guarded_maps=True) + finalize_for_target(sdfg, target) + + # WCR-config ON. ``canonicalize`` already sets this on every nested SDFG at the end of + # the pipeline; set it again explicitly right before compile for clarity. + for nested in sdfg.all_sdfgs_recursive(): + nested.openmp_array_reductions = True + + implementations = [] + try: + dc_exec = sdfg.compile() + implementations.append((dc_exec, sdfg._name)) + except Exception as e: + print("Failed to compile DaCe {a} canonicalize implementation.".format(a=self.info["arch"])) + print(e) + traceback.print_exc() + + return implementations diff --git a/npbench/infrastructure/dace_framework.py b/npbench/infrastructure/dace_framework.py index 6cbf141f5..274ca0639 100644 --- a/npbench/infrastructure/dace_framework.py +++ b/npbench/infrastructure/dace_framework.py @@ -274,6 +274,11 @@ def copy_to_gpu(sdfg): else: gpu_time1 = [0] fe_time += gpu_time1[0] + # WCR-config OFF: dace_cpu / dace_gpu keep the default per-element atomic WCR + # lowering (contrast dace_canonicalize_cpu / dace_canonicalize_gpu, which turn ON + # OpenMP array-section reductions). Set it explicitly for clarity before compile. + for nested in sdfg.all_sdfgs_recursive(): + nested.openmp_array_reductions = False try: dc_exec, compile_time = util.benchmark("__npb_result = sdfg.compile()", out_text="DaCe compilation time",