From fbff1452b89a641d729ebde84488a6125be62227 Mon Sep 17 00:00:00 2001 From: Axel Huebl Date: Thu, 20 Oct 2022 21:39:22 -0700 Subject: [PATCH] Helpers: to_numpy/cupy --- MANIFEST.in | 4 ++ src/amrex/Array4.py | 127 ++++++++++++++++++++++++++++++++++ src/amrex/__init__.py | 4 +- src/amrex/space1d/__init__.py | 3 + src/amrex/space2d/__init__.py | 3 + src/amrex/space3d/__init__.py | 3 + tests/test_multifab.py | 25 ++++--- 7 files changed, 155 insertions(+), 14 deletions(-) create mode 100644 src/amrex/Array4.py diff --git a/MANIFEST.in b/MANIFEST.in index 0b0abb03..9a61a6d7 100644 --- a/MANIFEST.in +++ b/MANIFEST.in @@ -6,6 +6,10 @@ recursive-include cmake * recursive-include src * recursive-include tests * +# avoid accidentially copying compiled Python files +global-exclude */__pycache__/* +global-exclude *.pyc + # see .gitignore prune cmake-build* prune .spack-env* diff --git a/src/amrex/Array4.py b/src/amrex/Array4.py new file mode 100644 index 00000000..d4eb98c9 --- /dev/null +++ b/src/amrex/Array4.py @@ -0,0 +1,127 @@ +""" +This file is part of pyAMReX + +Copyright 2022 AMReX community +Authors: Axel Huebl +License: BSD-3-Clause-LBNL +""" + +from .amrex_pybind import ( + Array4_double, + Array4_float, + Array4_int, + Array4_long, + Array4_longdouble, + Array4_longlong, + Array4_short, + Array4_uint, + Array4_ulong, + Array4_ulonglong, + Array4_ushort, +) + + +def array4_to_numpy(self, copy=False, order="F"): + """ + Provide a Numpy view into an Array4. + + Note on the order of indices: + By default, this is as in AMReX in Fortran contiguous order, indexing as + x,y,z. This has performance implications for use in external libraries such + as cupy. + The order="C" option will index as z,y,x and perform better with cupy. + https://github.com/AMReX-Codes/pyamrex/issues/55#issuecomment-1579610074 + + Parameters + ---------- + self : amrex.Array4_* + An Array4 class in pyAMReX + copy : bool, optional + Copy the data if true, otherwise create a view (default). + order : string, optional + F order (default) or C. C is faster with external libraries. + + Returns + ------- + np.array + A numpy n-dimensional array. + """ + import numpy as np + + if order == "F": + return np.array(self, copy=copy).T + elif order == "C": + return np.array(self, copy=copy) + else: + raise ValueError("The order argument must be F or C.") + + +def array4_to_cupy(self, copy=False, order="F"): + """ + Provide a Cupy view into an Array4. + + Note on the order of indices: + By default, this is as in AMReX in Fortran contiguous order, indexing as + x,y,z. This has performance implications for use in external libraries such + as cupy. + The order="C" option will index as z,y,x and perform better with cupy. + https://github.com/AMReX-Codes/pyamrex/issues/55#issuecomment-1579610074 + + Parameters + ---------- + self : amrex.Array4_* + An Array4 class in pyAMReX + copy : bool, optional + Copy the data if true, otherwise create a view (default). + order : string, optional + F order (default) or C. C is faster with external libraries. + + Returns + ------- + cupy.array + A numpy n-dimensional array. + + Raises + ------ + ImportError + Raises an exception if cupy is not installed + """ + import cupy as cp + + if order == "F": + return cp.array(self, copy=copy).T + elif order == "C": + return cp.array(self, copy=copy) + else: + raise ValueError("The order argument must be F or C.") + + +# Array4 helper methods +# +Array4_float.to_numpy = array4_to_numpy +Array4_double.to_numpy = array4_to_numpy +Array4_longdouble.to_numpy = array4_to_numpy + +Array4_short.to_numpy = array4_to_numpy +Array4_int.to_numpy = array4_to_numpy +Array4_long.to_numpy = array4_to_numpy +Array4_longlong.to_numpy = array4_to_numpy + +Array4_ushort.to_numpy = array4_to_numpy +Array4_uint.to_numpy = array4_to_numpy +Array4_ulong.to_numpy = array4_to_numpy +Array4_ulonglong.to_numpy = array4_to_numpy + +Array4_float.to_cupy = array4_to_cupy +Array4_double.to_cupy = array4_to_cupy +Array4_longdouble.to_cupy = array4_to_cupy + +Array4_short.to_cupy = array4_to_cupy +Array4_int.to_cupy = array4_to_cupy +Array4_long.to_cupy = array4_to_cupy +Array4_longlong.to_cupy = array4_to_cupy + +Array4_ushort.to_cupy = array4_to_cupy +Array4_uint.to_cupy = array4_to_cupy +Array4_ulong.to_cupy = array4_to_cupy +Array4_ulonglong.to_cupy = array4_to_cupy diff --git a/src/amrex/__init__.py b/src/amrex/__init__.py index 7313b1da..1b9b2bf4 100644 --- a/src/amrex/__init__.py +++ b/src/amrex/__init__.py @@ -1,7 +1,5 @@ # -*- coding: utf-8 -*- # __version__ is TODO - only in spaceNd submodules -__author__ = ( - "Axel Huebl, Ryan Sandberg, Shreyas Ananthan, Remi Lehe, " "Weiqun Zhang, et al." -) +__author__ = "Axel Huebl, Ryan Sandberg, Shreyas Ananthan, Remi Lehe, Andrew Myers, Reva Jambunathan, Edodardo Zoni, Weiqun Zhang" __license__ = "BSD-3-Clause-LBNL" diff --git a/src/amrex/space1d/__init__.py b/src/amrex/space1d/__init__.py index 6fac7643..db4e7bc9 100644 --- a/src/amrex/space1d/__init__.py +++ b/src/amrex/space1d/__init__.py @@ -32,3 +32,6 @@ # def d_decl(x, y, z): return (x,) + + +from .Array4 import * # noqa diff --git a/src/amrex/space2d/__init__.py b/src/amrex/space2d/__init__.py index c30ed8c5..8bf0fdce 100644 --- a/src/amrex/space2d/__init__.py +++ b/src/amrex/space2d/__init__.py @@ -32,3 +32,6 @@ # def d_decl(x, y, z): return (x, y) + + +from .Array4 import * # noqa diff --git a/src/amrex/space3d/__init__.py b/src/amrex/space3d/__init__.py index d8163ea8..15bd101b 100644 --- a/src/amrex/space3d/__init__.py +++ b/src/amrex/space3d/__init__.py @@ -32,3 +32,6 @@ # def d_decl(x, y, z): return (x, y, z) + + +from .Array4 import * # noqa diff --git a/tests/test_multifab.py b/tests/test_multifab.py index 7ad60c0d..838b7e94 100644 --- a/tests/test_multifab.py +++ b/tests/test_multifab.py @@ -46,17 +46,17 @@ def test_mfab_loop(make_mfab): # numpy representation: non-copying view, including the # guard/ghost region # note: in numpy, indices are in C-order! - marr_np = np.array(marr, copy=False) + marr_np = marr.to_numpy() # check the values at start/end are the same: first component assert marr_np[0, 0, 0, 0] == marr[bx.small_end] - assert marr_np[0, -1, -1, -1] == marr[bx.big_end] + assert marr_np[-1, -1, -1, 0] == marr[bx.big_end] # same check, but for all components for n in range(mfab.num_comp): small_end_comp = list(bx.small_end) + [n] big_end_comp = list(bx.big_end) + [n] - assert marr_np[n, 0, 0, 0] == marr[small_end_comp] - assert marr_np[n, -1, -1, -1] == marr[big_end_comp] + assert marr_np[0, 0, 0, n] == marr[small_end_comp] + assert marr_np[-1, -1, -1, n] == marr[big_end_comp] # now we do some faster assignments, using range based access # this should fail as out-of-bounds, but does not @@ -64,7 +64,7 @@ def test_mfab_loop(make_mfab): # marr_np[24:200, :, :, :] = 42. # all components and all indices set at once to 42 - marr_np[:, :, :, :] = 42.0 + marr_np[()] = 42.0 # values in start & end still match? assert marr_np[0, 0, 0, 0] == marr[bx.small_end] @@ -210,10 +210,11 @@ def test_mfab_ops_cuda_cupy(make_mfab_device): with cupy.profiler.time_range("assign 3 [()]", color_id=0): for mfi in mfab_device: bx = mfi.tilebox().grow(ngv) - marr = mfab_device.array(mfi) - marr_cupy = cp.array(marr, copy=False) + marr_cupy = mfab_device.array(mfi).to_cupy(order="C") # print(marr_cupy.shape) # 1, 32, 32, 32 # print(marr_cupy.dtype) # float64 + # performance: + # https://github.com/AMReX-Codes/pyamrex/issues/55#issuecomment-1579610074 # write and read into the marr_cupy marr_cupy[()] = 3.0 @@ -244,8 +245,11 @@ def set_to_five(mm): for mfi in mfab_device: bx = mfi.tilebox().grow(ngv) - marr = mfab_device.array(mfi) - marr_cupy = cp.array(marr, copy=False) + marr_cupy = mfab_device.array(mfi).to_cupy(order="F") + # print(marr_cupy.shape) # 32, 32, 32, 1 + # print(marr_cupy.dtype) # float64 + # performance: + # https://github.com/AMReX-Codes/pyamrex/issues/55#issuecomment-1579610074 # write and read into the marr_cupy fives_cp = set_to_five(marr_cupy) @@ -266,8 +270,7 @@ def set_to_seven(x): for mfi in mfab_device: bx = mfi.tilebox().grow(ngv) - marr = mfab_device.array(mfi) - marr_cupy = cp.array(marr, copy=False) + marr_cupy = mfab_device.array(mfi).to_cupy(order="C") # write and read into the marr_cupy set_to_seven(marr_cupy)