Zengtudor/pybind11
1
0
Fork 0
mirror of https://github.com/pybind/pybind11.git synced 2025-02-16 13:47:53 +00:00
Code Issues Packages Projects Releases Wiki Activity

Making a copy when casting a numpy array with negative strides to Eigen.

Browse Source 
`EigenConformable::stride_compatible` returns false if the strides are
negative. In this case, do not use `EigenConformable::stride`, as it
is {0,0}. We cannot write negative strides in this element, as Eigen
will throw an assertion if we do.

The `type_caster` specialization for regular, dense Eigen matrices now
does a second `array_t::ensure` to copy data in case of negative strides.
I'm not sure that this is the best way to implement this.

I have added "TODO" tags linking these changes to Eigen bug #747, which,
when fixed, will allow Eigen to accept negative strides.
This commit is contained in:
Cris Luengo 2017-04-06 11:34:39 -06:00 committed by Dean Moldovan
parent d400f60c96
commit 627da3f135
4 changed files with 93 additions and 23 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

41
include/pybind11/eigen.h
View File

@ -68,16 +68,22 @@ template <typename T> using is_eigen_other = all_of<
template <bool EigenRowMajor> struct EigenConformable {
bool conformable = false;
EigenIndex rows = 0, cols = 0;
EigenDStride stride{0, 0};
EigenDStride stride{0, 0}; // Only valid if negativestridees is false!
bool negativestrides = false; // If true, do not use stride!
EigenConformable(bool fits = false) : conformable{fits} {}
// Matrix type:
EigenConformable(EigenIndex r, EigenIndex c,
EigenIndex rstride, EigenIndex cstride) :
conformable{true}, rows{r}, cols{c},
stride(EigenRowMajor ? rstride : cstride /* outer stride */,
EigenRowMajor ? cstride : rstride /* inner stride */)
{}
conformable{true}, rows{r}, cols{c} {
// TODO: when Eigen bug #747 is fixed, remove the tests for non-negativity. http://eigen.tuxfamily.org/bz/show_bug.cgi?id=747
if (rstride < 0 || cstride < 0) {
negativestrides = true;
} else {
stride = {EigenRowMajor ? rstride : cstride /* outer stride */,
EigenRowMajor ? cstride : rstride /* inner stride */ };
}
}
// Vector type:
EigenConformable(EigenIndex r, EigenIndex c, EigenIndex stride)
: EigenConformable(r, c, r == 1 ? c*stride : stride, c == 1 ? r : r*stride) {}
@ -86,6 +92,7 @@ template <bool EigenRowMajor> struct EigenConformable {
// To have compatible strides, we need (on both dimensions) one of fully dynamic strides,
// matching strides, or a dimension size of 1 (in which case the stride value is irrelevant)
return
!negativestrides &&
(props::inner_stride == Eigen::Dynamic || props::inner_stride == stride.inner() ||
(EigenRowMajor ? cols : rows) == 1) &&
(props::outer_stride == Eigen::Dynamic || props::outer_stride == stride.outer() ||
@ -138,8 +145,8 @@ template <typename Type_> struct EigenProps {
EigenIndex
np_rows = a.shape(0),
np_cols = a.shape(1),
np_rstride = a.strides(0) / sizeof(Scalar),
np_cstride = a.strides(1) / sizeof(Scalar);
np_rstride = a.strides(0) / static_cast<ssize_t>(sizeof(Scalar)),
np_cstride = a.strides(1) / static_cast<ssize_t>(sizeof(Scalar));
if ((fixed_rows && np_rows != rows) || (fixed_cols && np_cols != cols))
return false;
@ -149,7 +156,7 @@ template <typename Type_> struct EigenProps {
// Otherwise we're storing an n-vector. Only one of the strides will be used, but whichever
// is used, we want the (single) numpy stride value.
const EigenIndex n = a.shape(0),
stride = a.strides(0) / sizeof(Scalar);
stride = a.strides(0) / static_cast<ssize_t>(sizeof(Scalar));
if (vector) { // Eigen type is a compile-time vector
if (fixed && size != n)
@ -255,7 +262,23 @@ struct type_caster<Type, enable_if_t<is_eigen_dense_plain<Type>::value>> {
if (!fits)
return false; // Non-comformable vector/matrix types
value = Eigen::Map<const Type, 0, EigenDStride>(buf.data(), fits.rows, fits.cols, fits.stride);
if (fits.negativestrides) {
// Eigen does not support negative strides, so we need to make a copy here with normal strides.
// TODO: when Eigen bug #747 is fixed, remove this if case, always execute the else part.
// http://eigen.tuxfamily.org/bz/show_bug.cgi?id=747
auto buf2 = array_t<Scalar,array::forcecast || array::f_style>::ensure(src);
if (!buf2)
return false;
// not checking sizes, we already did that
fits = props::conformable(buf2);
value = Eigen::Map<const Type, 0, EigenDStride>(buf2.data(), fits.rows, fits.cols, fits.stride);
} else {
value = Eigen::Map<const Type, 0, EigenDStride>(buf.data(), fits.rows, fits.cols, fits.stride);
}
return true;
}

24
include/pybind11/numpy.h
View File

@ -248,10 +248,10 @@ template <typename T> using is_pod_struct = all_of<
satisfies_none_of<T, std::is_reference, std::is_array, is_std_array, std::is_arithmetic, is_complex, std::is_enum>
>;
template <size_t Dim = 0, typename Strides> size_t byte_offset_unsafe(const Strides &) { return 0; }
template <size_t Dim = 0, typename Strides> ssize_t byte_offset_unsafe(const Strides &) { return 0; }
template <size_t Dim = 0, typename Strides, typename... Ix>
size_t byte_offset_unsafe(const Strides &strides, size_t i, Ix... index) {
return i * strides[Dim] + byte_offset_unsafe<Dim + 1>(strides, index...);
ssize_t byte_offset_unsafe(const Strides &strides, size_t i, Ix... index) {
return static_cast<ssize_t>(i) * strides[Dim] + byte_offset_unsafe<Dim + 1>(strides, index...);
}
/** Proxy class providing unsafe, unchecked const access to array data. This is constructed through
@ -615,18 +615,18 @@ public:
/// Byte offset from beginning of the array to a given index (full or partial).
/// May throw if the index would lead to out of bounds access.
template<typename... Ix> size_t offset_at(Ix... index) const {
template<typename... Ix> ssize_t offset_at(Ix... index) const {
if (sizeof...(index) > ndim())
fail_dim_check(sizeof...(index), "too many indices for an array");
return byte_offset(size_t(index)...);
}
size_t offset_at() const { return 0; }
ssize_t offset_at() const { return 0; }
/// Item count from beginning of the array to a given index (full or partial).
/// May throw if the index would lead to out of bounds access.
template<typename... Ix> size_t index_at(Ix... index) const {
return offset_at(index...) / itemsize();
template<typename... Ix> ssize_t index_at(Ix... index) const {
return offset_at(index...) / static_cast<ssize_t>(itemsize());
}
/** Returns a proxy object that provides access to the array's data without bounds or
@ -692,7 +692,7 @@ protected:
" (ndim = " + std::to_string(ndim()) + ")");
}
template<typename... Ix> size_t byte_offset(Ix... index) const {
template<typename... Ix> ssize_t byte_offset(Ix... index) const {
check_dimensions(index...);
return detail::byte_offset_unsafe(strides(), size_t(index)...);
}
@ -773,8 +773,8 @@ public:
return sizeof(T);
}
template<typename... Ix> size_t index_at(Ix... index) const {
return offset_at(index...) / itemsize();
template<typename... Ix> ssize_t index_at(Ix... index) const {
return offset_at(index...) / static_cast<ssize_t>(itemsize());
}
template<typename... Ix> const T* data(Ix... index) const {
@ -789,14 +789,14 @@ public:
template<typename... Ix> const T& at(Ix... index) const {
if (sizeof...(index) != ndim())
fail_dim_check(sizeof...(index), "index dimension mismatch");
return *(static_cast<const T*>(array::data()) + byte_offset(size_t(index)...) / itemsize());
return *(static_cast<const T*>(array::data()) + byte_offset(size_t(index)...) / static_cast<ssize_t>(itemsize()));
}
// Mutable reference to element at a given index
template<typename... Ix> T& mutable_at(Ix... index) {
if (sizeof...(index) != ndim())
fail_dim_check(sizeof...(index), "index dimension mismatch");
return *(static_cast<T*>(array::mutable_data()) + byte_offset(size_t(index)...) / itemsize());
return *(static_cast<T*>(array::mutable_data()) + byte_offset(size_t(index)...) / static_cast<ssize_t>(itemsize()));
}
/** Returns a proxy object that provides access to the array's data without bounds or

4
include/pybind11/stl_bind.h
View File

@ -350,14 +350,14 @@ vector_buffer(Class_& cl) {
cl.def("__init__", [](Vector& vec, buffer buf) {
auto info = buf.request();
if (info.ndim != 1 || info.strides[0] <= 0 || info.strides[0] % sizeof(T))
if (info.ndim != 1 || info.strides[0] <= 0 || info.strides[0] % static_cast<ssize_t>(sizeof(T)))
throw type_error("Only valid 1D buffers can be copied to a vector");
if (!detail::compare_buffer_info<T>::compare(info) || sizeof(T) != info.itemsize)
throw type_error("Format mismatch (Python: " + info.format + " C++: " + format_descriptor<T>::format() + ")");
new (&vec) Vector();
vec.reserve(info.shape[0]);
T *p = static_cast<T*>(info.ptr);
auto step = info.strides[0] / sizeof(T);
auto step = info.strides[0] / static_cast<ssize_t>(sizeof(T));
T *end = p + info.shape[0] * step;
for (; p < end; p += step)
vec.push_back(*p);

47
tests/test_eigen.py
View File

@ -154,6 +154,53 @@ def test_nonunit_stride_from_python():
np.testing.assert_array_equal(counting_mat, [[0., 2, 2], [6, 16, 10], [6, 14, 8]])
def test_negative_stride_from_python(msg):
from pybind11_tests import (
double_row, double_col, double_complex, double_mat_cm, double_mat_rm,
double_threec, double_threer)
# Eigen doesn't support (as of yet) negative strides. When a function takes an Eigen
# matrix by copy or const reference, we can pass a numpy array that has negative strides.
# Otherwise, an exception will be thrown as Eigen will not be able to map the numpy array.
counting_mat = np.arange(9.0, dtype=np.float32).reshape((3, 3))
counting_mat = counting_mat[::-1, ::-1]
second_row = counting_mat[1, :]
second_col = counting_mat[:, 1]
np.testing.assert_array_equal(double_row(second_row), 2.0 * second_row)
np.testing.assert_array_equal(double_col(second_row), 2.0 * second_row)
np.testing.assert_array_equal(double_complex(second_row), 2.0 * second_row)
np.testing.assert_array_equal(double_row(second_col), 2.0 * second_col)
np.testing.assert_array_equal(double_col(second_col), 2.0 * second_col)
np.testing.assert_array_equal(double_complex(second_col), 2.0 * second_col)
counting_3d = np.arange(27.0, dtype=np.float32).reshape((3, 3, 3))
counting_3d = counting_3d[::-1, ::-1, ::-1]
slices = [counting_3d[0, :, :], counting_3d[:, 0, :], counting_3d[:, :, 0]]
for slice_idx, ref_mat in enumerate(slices):
np.testing.assert_array_equal(double_mat_cm(ref_mat), 2.0 * ref_mat)
np.testing.assert_array_equal(double_mat_rm(ref_mat), 2.0 * ref_mat)
# Mutator:
with pytest.raises(TypeError) as excinfo:
double_threer(second_row)
assert msg(excinfo.value) == """
double_threer(): incompatible function arguments. The following argument types are supported:
1. (numpy.ndarray[float32[1, 3], flags.writeable]) -> arg0: None
Invoked with: array([ 5., 4., 3.], dtype=float32)
"""
with pytest.raises(TypeError) as excinfo:
double_threec(second_col)
assert msg(excinfo.value) == """
double_threec(): incompatible function arguments. The following argument types are supported:
1. (numpy.ndarray[float32[3, 1], flags.writeable]) -> arg0: None
Invoked with: array([ 7., 4., 1.], dtype=float32)
"""
def test_nonunit_stride_to_python():
from pybind11_tests import diagonal, diagonal_1, diagonal_n, block