diff --git a/CHANGELOG.md b/CHANGELOG.md index 3e4844fc7..192427b81 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -46,7 +46,8 @@ and this project aspires to adhere to [Semantic Versioning](https://semver.org/s #### General - Updated to newer BLT to resolve BLT/FindMPI issues with rpath linking commands when using OpenMPI. - Fixed internal object name string for the Python Iterator object. It used to report `Schema`, which triggered both puzzling and concerned emotions. - +- Fixed a bug with `Node.set` in the Python API that undermined setting NumPy arrays with sliced views and complex striding. General slices should now work with `set`. No changes to the `set_external` case, which requires 1-D effective striding and throws an exception when more complex strides are presented. + #### Relay - Use H5F_ACC_RDONLY in relay::io::is_hdf5_file to avoid errors when checking files that already have open HDF5 handles. diff --git a/src/libs/conduit/python/conduit_python.cpp b/src/libs/conduit/python/conduit_python.cpp index 2bc88badc..055567cfa 100644 --- a/src/libs/conduit/python/conduit_python.cpp +++ b/src/libs/conduit/python/conduit_python.cpp @@ -4007,17 +4007,91 @@ static void PyConduit_Fill_DataArray_From_PyArray(DataArray &conduit_array, PyArrayObject *numpy_array) { - npy_int num_ele = (npy_int)conduit_array.number_of_elements(); + // proper strided iteration, adapted from: + //https://docs.scipy.org/doc/numpy-1.10.0/reference/c-api.iterator.html + + + /* Handle zero-sized arrays specially */ + if (PyArray_SIZE(numpy_array) == 0) + { + // nothing to copy in this case! + return; + } + + /* + * Create and use an iterator to copy out the data to our conduit array + * flag NPY_ITER_READONLY + * - The array is never written to. + * flag NPY_ITER_EXTERNAL_LOOP + * - Inner loop is done outside the iterator for efficiency. + * flag NPY_ITER_NPY_ITER_REFS_OK + * - Reference types are acceptable. + * order NPY_KEEPORDER + * - Visit elements in memory order, regardless of strides. + * This is good for performance when the specific order + * elements are visited is unimportant. + * casting NPY_NO_CASTING + * - No casting is required for this operation. + */ + NpyIter* iter = NpyIter_New(numpy_array, + NPY_ITER_READONLY| + NPY_ITER_EXTERNAL_LOOP| + NPY_ITER_REFS_OK, + NPY_KEEPORDER, + NPY_NO_CASTING, + NULL); + if (iter == NULL) + { + // TODO ERROR! + } + + /* + * The iternext function gets stored in a local variable + * so it can be called repeatedly in an efficient manner. + */ + NpyIter_IterNextFunc *iternext = NpyIter_GetIterNext(iter, NULL); + if (iternext == NULL) + { + NpyIter_Deallocate(iter); + // TODO ERROR! + } + + // The location of the data pointer which the iterator may update + char** dataptr = NpyIter_GetDataPtrArray(iter); + // The location of the stride which the iterator may update + npy_intp* strideptr = NpyIter_GetInnerStrideArray(iter); + // The location of the inner loop size which the iterator may update + npy_intp* innersizeptr = NpyIter_GetInnerLoopSizePtr(iter); + + int idx=0; + do + { + // Get the inner loop data/stride/count values + char* data = *dataptr; + npy_intp stride = *strideptr; + npy_intp count = *innersizeptr; + + // This is a typical inner loop for NPY_ITER_EXTERNAL_LOOP + while (count--) + { + conduit_array[idx] = ((T*)data)[0]; + data += stride; + idx++; + } + // Increment the iterator to the next inner loop + } + while(iternext(iter)); - T *numpy_data = (T*)PyArray_BYTES(numpy_array); + // TODO, add just in case consistency check? + // npy_int num_ele = (npy_int)conduit_array.number_of_elements(); + // if(num_ele != idx) + // { + // // this doesn't match our expectations + // } - for (npy_intp i = 0; i < num_ele; i++) - { - conduit_array[i] = numpy_data[i]; - } + NpyIter_Deallocate(iter); } - //---------------------------------------------------------------------------// // begin Node python special methods //---------------------------------------------------------------------------// @@ -5093,9 +5167,14 @@ PyConduit_Node_set_external(PyConduit_Node* self, index_t stride = (index_t) PyArray_STRIDE(py_arr, 0); int nd = PyArray_NDIM(py_arr); - if (nd > 1) { + if (nd > 1) + { PyErr_SetString(PyExc_TypeError, - "set_external does not handle multidimensional numpy arrays"); + "set_external does not handle multidimensional numpy" + " arrays or multidimensional complex strided views" + " into numpy arrays." + " Views that are effectively 1D-strided are" + "supported."); return (NULL); } diff --git a/src/libs/relay/python/py_src/__init__.py b/src/libs/relay/python/py_src/__init__.py index 284cbb2f2..98f884886 100644 --- a/src/libs/relay/python/py_src/__init__.py +++ b/src/libs/relay/python/py_src/__init__.py @@ -10,7 +10,12 @@ from . import io from . import web -from . import mpi + +# mpi support is optional, so drive on if we can't import +try: + from . import mpi +except: + pass diff --git a/src/tests/conduit/python/t_python_conduit_node.py b/src/tests/conduit/python/t_python_conduit_node.py index 78dc5d7d5..9cc333285 100644 --- a/src/tests/conduit/python/t_python_conduit_node.py +++ b/src/tests/conduit/python/t_python_conduit_node.py @@ -14,7 +14,6 @@ import numpy as np - class Test_Conduit_Node(unittest.TestCase): def test_simple(self): a_val = np.uint32(10) @@ -469,14 +468,73 @@ def test_to_string(self): res_to_yaml = n.to_yaml() res_to_json = n.to_json() - + self.assertEqual(res_to_str_def, res_to_yaml); self.assertEqual(res_to_str_yaml, res_to_yaml); - + self.assertEqual(res_to_str_json, res_to_json); n.print_detailed() + def test_numpy_slice_as_set_input(self): + n = Node() + # slice with non trivial strides + numpy_array = np.array(range(21), dtype='float64') + v = numpy_array.reshape((3, 7)) + print("Input Array") + print(v) + print("Desired Slice") + print(v[:,0]) + n['v'] = v + n['vs'] = v[:,0] + n['vs_expected'] = np.array(v[:,0],np.float64) + print(n) + sdiff = np.setdiff1d(n['vs'], v[:,0]) + print("Set Difference: ",sdiff ) + self.assertEqual(len(sdiff), 0); + # a more complex slice + numpy_array = np.array(range(105), dtype='float64') + v = numpy_array.reshape((3, 7, 5)) + print("Input Array") + print(v) + print("Desired Slice") + print(v[:,0,3:5]) + n['v'] = v + n['vs'] = v[:,0,3:5] + n['vs_expected'] = np.array(v[:,0,3:5],np.float64) + print(n) + sdiff = np.setdiff1d(n['vs'], v[:,0,3:5]) + print("Set Difference: ",sdiff ) + self.assertEqual(len(sdiff), 0); + + + def test_numpy_slice_as_set_external_input(self): + n = Node() + # slice with non trivial strides + numpy_array = np.array(range(21), dtype='float64') + v = numpy_array.reshape((3, 7)) + print("Input Array") + print(v) + print("Desired Slice") + print(v[:,0]) + n['v'] = v + n['vs'].set_external(v[:,0]) + n['vs_expected'] = np.array(v[:,0],np.float64) + print(n) + sdiff = np.setdiff1d(n['vs'], v[:,0]) + print("Set Difference: ",sdiff ) + self.assertEqual(len(sdiff), 0); + # a more complex slice, can't use set external here. + n = Node() + numpy_array = np.array(range(105), dtype='float64') + v = numpy_array.reshape((3, 7, 5)) + with self.assertRaises(TypeError): + n['vs'].set_external(v[:,0,3:5]) + # lets do a 1-d eff slice, this should work since + # it reduces to a 1-D strided case + n['vs'].set_external(v[:,0,0]) + n['vs_expected'] = np.array(v[:,0,0],np.float64) + if __name__ == '__main__': unittest.main()