Repository: rust-ndarray/ndarray
Pull request #1505 (2025-04-20 04:56:12)

Axis permutation to correctly handle cycles using bitmask approach
## Overview
What's up! 😃 Currently, there're two functions to deal with permuting and reversing arrays(`permuted_axes` and `reversed_axes`).
As mentioned from related issue, they're not `in-place`, so I added `permute_axes` and `reverse_axes`. The new implementation for `permute_axes` uses a bitmask-based approach to efficiently track and process cycles in axis permutations!

## Related issue
- #1443 

## Changes
1. **permute_axes (in-place)**: 
   - A new function with a cycle detection algorithm based on bitmasking. 
   - try to maintain the zero-allocation, in-place modification approach
  
2. **reverse_axes (in-place)**:
    - simply use `&mut self`

## Details on permute_axes

### How it works
The new implementation uses a bitmask to track visited axes during permutation cycles:

1. **Cycle Detection**: When permuting axes, we identify cycles of positions that need to be updated.
2. **Bitmask Tracking**: A single `usize` is used as a bitmask to efficiently track which axes have been processed.
4. **In-place Updates**: Each cycle is processed in-place by following the permutation pattern.

### Example

Consider permuting a 2×2 array from `[[1, 2], [3, 4]]` with axes permutation `[1, 0]` (transpose):

```rust
// Original array
let mut arr = array![[1, 2], [3, 4]];
// Original: dims = [2, 2], strides = [2, 1]

// Permute axes [1, 0]
arr.permute_axes([1, 0]);
// Result: dims = [2, 2], strides = [1, 2]
// Array is now [[1, 3], [2, 4]]
```

Here's how the algorithm processes this:

1. Start with axis 1 (new_axis=0):
   - Store initial values: dim[1]=2, stride[1]=1
   - Follow cycle: 1 → 0 → 1
   - Update dim[1] = dim[0] = 2
   - Update stride[1] = stride[0] = 2
   - Mark axes 0 and 1 as visited

2. Process axis 0 (new_axis=1):
   - Already visited in the previous cycle, skip

The bitmask approach ensures each axis is processed exactly once, even in complex permutations with multiple cycles. The key operations are:
- `visited |= 1 << axis`: Mark axis as visited
- `(visited & (1 << axis)) != 0`: Check if axis is visited

@akern40 @nilgoyette added test cases similar to the original functions(`premuted_axes` and `reversed_axes`), but there might be missing cases. please take a look when you guys have some time! 🚀 





Patch:
diff --git a/src/dimension/mod.rs b/src/dimension/mod.rs
index eb07252b..4ada0b12 100644
--- a/src/dimension/mod.rs
+++ b/src/dimension/mod.rs
@@ -269,7 +269,7 @@ fn can_index_slice_impl<D: Dimension>(
 ) -> Result<(), ShapeError>
 {
     // Check condition 3.
-    let is_empty = dim.slice().iter().any(|&d| d == 0);
+    let is_empty = dim.slice().contains(&0);
     if is_empty && max_offset > data_len {
         return Err(from_kind(ErrorKind::OutOfBounds));
     }
diff --git a/src/impl_methods.rs b/src/impl_methods.rs
index ea9c9a0d..9c2c5efc 100644
--- a/src/impl_methods.rs
+++ b/src/impl_methods.rs
@@ -208,7 +208,7 @@ impl<A, D: Dimension> ArrayRef<A, D>
     /// memory layout. Otherwise, the layout of the output array is unspecified.
     /// If you need a particular layout, you can allocate a new array with the
     /// desired memory layout and [`.assign()`](Self::assign) the data.
-    /// Alternatively, you can collectan iterator, like this for a result in
+    /// Alternatively, you can collect an iterator, like this for a result in
     /// standard layout:
     ///
     /// ```
@@ -576,11 +576,7 @@ where
     pub fn slice_move<I>(mut self, info: I) -> ArrayBase<S, I::OutDim>
     where I: SliceArg<D>
     {
-        assert_eq!(
-            info.in_ndim(),
-            self.ndim(),
-            "The input dimension of `info` must match the array to be sliced.",
-        );
+        assert_eq!(info.in_ndim(), self.ndim(), "The input dimension of `info` must match the array to be sliced.",);
         let out_ndim = info.out_ndim();
         let mut new_dim = I::OutDim::zeros(out_ndim);
         let mut new_strides = I::OutDim::zeros(out_ndim);
@@ -648,11 +644,7 @@ impl<A, D: Dimension> LayoutRef<A, D>
     pub fn slice_collapse<I>(&mut self, info: I)
     where I: SliceArg<D>
     {
-        assert_eq!(
-            info.in_ndim(),
-            self.ndim(),
-            "The input dimension of `info` must match the array to be sliced.",
-        );
+        assert_eq!(info.in_ndim(), self.ndim(), "The input dimension of `info` must match the array to be sliced.",);
         let mut axis = 0;
         info.as_ref().iter().for_each(|&ax_info| match ax_info {
             SliceInfoElem::Slice { start, end, step } => {
@@ -1120,8 +1112,7 @@ impl<A, D: Dimension> ArrayRef<A, D>
             // bounds check the indices first
             if let Some(max_index) = indices.iter().cloned().max() {
                 if max_index >= axis_len {
-                    panic!("ndarray: index {} is out of bounds in array of len {}",
-                           max_index, self.len_of(axis));
+                    panic!("ndarray: index {} is out of bounds in array of len {}", max_index, self.len_of(axis));
                 }
             } // else: indices empty is ok
             let view = self.view().into_dimensionality::<Ix1>().unwrap();
@@ -1525,21 +1516,34 @@ impl<A, D: Dimension> ArrayRef<A, D>
     /// }
     /// ```
     pub fn axis_windows(&self, axis: Axis, window_size: usize) -> AxisWindows<'_, A, D>
+    {
+        self.axis_windows_with_stride(axis, window_size, 1)
+    }
+
+    /// Returns a producer which traverses over windows of a given length and
+    /// stride along an axis.
+    ///
+    /// Note that a calling this method with a stride of 1 is equivalent to
+    /// calling [`ArrayRef::axis_windows()`].
+    pub fn axis_windows_with_stride(&self, axis: Axis, window_size: usize, stride_size: usize)
+        -> AxisWindows<'_, A, D>
     {
         let axis_index = axis.index();
 
         ndassert!(
             axis_index < self.ndim(),
-            concat!(
-                "Window axis {} does not match array dimension {} ",
-                "(with array of shape {:?})"
-            ),
+            concat!("Window axis {} does not match array dimension {} ", "(with array of shape {:?})"),
             axis_index,
             self.ndim(),
             self.shape()
         );
 
-        AxisWindows::new(self.view(), axis, window_size)
+        ndassert!(
+            stride_size >0,
+            "Stride size must be greater than zero"
+        );
+
+        AxisWindows::new_with_stride(self.view(), axis, window_size, stride_size)
     }
 
     /// Return a view of the diagonal elements of the array.
@@ -2350,7 +2354,7 @@ impl<A, D: Dimension> ArrayRef<A, D>
     /// The implementation creates a view with strides set to zero for the
     /// axes that are to be repeated.
     ///
-    /// The broadcasting documentation for Numpy has more information.
+    /// The broadcasting documentation for NumPy has more information.
     ///
     /// ```
     /// use ndarray::{aview1, aview2};
@@ -2538,6 +2542,72 @@ where
         unsafe { self.with_strides_dim(new_strides, new_dim) }
     }
 
+    /// Permute the axes in-place.
+    ///
+    /// This does not move any data, it just adjusts the array's dimensions
+    /// and strides.
+    ///
+    /// *i* in the *j*-th place in the axes sequence means `self`'s *i*-th axis
+    /// becomes `self`'s *j*-th axis
+    ///
+    /// **Panics** if any of the axes are out of bounds, if an axis is missing,
+    /// or if an axis is repeated more than once.    
+    ///
+    /// # Example
+    /// ```rust
+    /// use ndarray::{arr2, Array3};
+    ///
+    /// let mut a = arr2(&[[0, 1], [2, 3]]);
+    /// a.permute_axes([1, 0]);
+    /// assert_eq!(a, arr2(&[[0, 2], [1, 3]]));
+    ///
+    /// let mut b = Array3::<u8>::zeros((1, 2, 3));
+    /// b.permute_axes([1, 0, 2]);
+    /// assert_eq!(b.shape(), &[2, 1, 3]);
+    /// ```
+    #[track_caller]
+    pub fn permute_axes<T>(&mut self, axes: T)
+    where T: IntoDimension<Dim = D>
+    {
+        let axes = axes.into_dimension();
+        // Ensure that each axis is used exactly once.
+        let mut usage_counts = D::zeros(self.ndim());
+        for axis in axes.slice() {
+            usage_counts[*axis] += 1;
+        }
+        for count in usage_counts.slice() {
+            assert_eq!(*count, 1, "each axis must be listed exactly once");
+        }
+
+        let dim = self.layout.dim.slice_mut();
+        let strides = self.layout.strides.slice_mut();
+        let axes = axes.slice();
+
+        // The cycle detection is done using a bitmask to track visited positions.
+        // For example, axes from [0,1,2] to [2, 0, 1]
+        // For axis values [1, 0, 2]:
+        // 1 << 1  // 0b0001 << 1 = 0b0010 (decimal 2)
+        // 1 << 0  // 0b0001 << 0 = 0b0001 (decimal 1)
+        // 1 << 2  // 0b0001 << 2 = 0b0100 (decimal 4)
+        //
+        // Each axis gets its own unique bit position in the bitmask:
+        // - Axis 0: bit 0 (rightmost)
+        // - Axis 1: bit 1
+        // - Axis 2: bit 2
+        //
+        let mut visited = 0usize;
+        for (new_axis, &axis) in axes.iter().enumerate() {
+            if (visited & (1 << axis)) != 0 {
+                continue;
+            }
+
+            dim.swap(axis, new_axis);
+            strides.swap(axis, new_axis);
+
+            visited |= (1 << axis) | (1 << new_axis);
+        }
+    }
+
     /// Transpose the array by reversing axes.
     ///
     /// Transposition reverses the order of the axes (dimensions and strides)
@@ -2548,6 +2618,16 @@ where
         self.layout.strides.slice_mut().reverse();
         self
     }
+
+    /// Reverse the axes of the array in-place.
+    ///
+    /// This does not move any data, it just adjusts the array's dimensions
+    /// and strides.
+    pub fn reverse_axes(&mut self)
+    {
+        self.layout.dim.slice_mut().reverse();
+        self.layout.strides.slice_mut().reverse();
+    }
 }
 
 impl<A, D: Dimension> ArrayRef<A, D>
@@ -2702,7 +2782,7 @@ where
 
 impl<A, D: Dimension> ArrayRef<A, D>
 {
-    /// Perform an elementwise assigment to `self` from `rhs`.
+    /// Perform an elementwise assignment to `self` from `rhs`.
     ///
     /// If their shapes disagree, `rhs` is broadcast to the shape of `self`.
     ///
@@ -2714,7 +2794,7 @@ impl<A, D: Dimension> ArrayRef<A, D>
         self.zip_mut_with(rhs, |x, y| x.clone_from(y));
     }
 
-    /// Perform an elementwise assigment of values cloned from `self` into array or producer `to`.
+    /// Perform an elementwise assignment of values cloned from `self` into array or producer `to`.
     ///
     /// The destination `to` can be another array or a producer of assignable elements.
     /// [`AssignElem`] determines how elements are assigned.
@@ -2730,7 +2810,7 @@ impl<A, D: Dimension> ArrayRef<A, D>
         Zip::from(self).map_assign_into(to, A::clone);
     }
 
-    /// Perform an elementwise assigment to `self` from element `x`.
+    /// Perform an elementwise assignment to `self` from element `x`.
     pub fn fill(&mut self, x: A)
     where A: Clone
     {
@@ -3119,8 +3199,7 @@ where
     /// ***Panics*** if not `index < self.len_of(axis)`.
     pub fn remove_index(&mut self, axis: Axis, index: usize)
     {
-        assert!(index < self.len_of(axis), "index {} must be less than length of Axis({})",
-                index, axis.index());
+        assert!(index < self.len_of(axis), "index {} must be less than length of Axis({})", index, axis.index());
         let (_, mut tail) = self.view_mut().split_at(axis, index);
         // shift elements to the front
         Zip::from(tail.lanes_mut(axis)).for_each(|mut lane| lane.rotate1_front());
@@ -3193,15 +3272,16 @@ impl<A, D: Dimension> ArrayRef<A, D>
     /// - All elements equal or greater than the k-th element to its right
     /// - The ordering within each partition is undefined
     ///
+    /// Empty arrays (i.e., those with any zero-length axes) are considered partitioned already,
+    /// and will be returned unchanged.
+    ///
+    /// **Panics** if `k` is out of bounds for a non-zero axis length.
+    ///
     /// # Parameters
     ///
     /// * `kth` - Index to partition by. The k-th element will be in its sorted position.
     /// * `axis` - Axis along which to partition.
     ///
-    /// # Returns
-    ///
-    /// A new array of the same shape and type as the input array, with elements partitioned.
-    ///
     /// # Examples
     ///
     /// ```
@@ -3221,19 +3301,19 @@ impl<A, D: Dimension> ArrayRef<A, D>
         A: Clone + Ord + num_traits::Zero,
         D: Dimension,
     {
-        // Bounds checking
-        let axis_len = self.len_of(axis);
-        if kth >= axis_len {
-            panic!("partition index {} is out of bounds for axis of length {}", kth, axis_len);
-        }
-
         let mut result = self.to_owned();
 
-        // Must guarantee that the array isn't empty before checking for contiguity
-        if result.shape().iter().any(|s| *s == 0) {
+        // Return early if the array has zero-length dimensions
+        if result.shape().contains(&0) {
             return result;
         }
 
+        // Bounds checking. Panics if kth is out of bounds
+        let axis_len = self.len_of(axis);
+        if kth >= axis_len {
+            panic!("Partition index {} is out of bounds for axis {} of length {}", kth, axis.0, axis_len);
+        }
+
         // Check if the first lane is contiguous
         let is_contiguous = result
             .lanes_mut(axis)
@@ -3428,11 +3508,7 @@ mod tests
     fn test_partition_contiguous_or_not()
     {
         // Test contiguous case (C-order)
-        let a = array![
-            [7, 1, 5],
-            [2, 6, 0],
-            [3, 4, 8]
-        ];
+        let a = array![[7, 1, 5], [2, 6, 0], [3, 4, 8]];
 
         // Partition along axis 0 (contiguous)
         let p_axis0 = a.partition(1, Axis(0));
@@ -3442,20 +3518,24 @@ mod tests
         // - Last row should be >= middle row (kth element)
         for col in 0..3 {
             let kth = p_axis0[[1, col]];
-            assert!(p_axis0[[0, col]] <= kth,
+            assert!(
+                p_axis0[[0, col]] <= kth,
                 "Column {}: First row {} should be <= middle row {}",
-                col, p_axis0[[0, col]], kth);
-            assert!(p_axis0[[2, col]] >= kth,
+                col,
+                p_axis0[[0, col]],
+                kth
+            );
+            assert!(
+                p_axis0[[2, col]] >= kth,
                 "Column {}: Last row {} should be >= middle row {}",
-                col, p_axis0[[2, col]], kth);
+                col,
+                p_axis0[[2, col]],
+                kth
+            );
         }
 
         // Test non-contiguous case (F-order)
-        let a = array![
-            [7, 1, 5],
-            [2, 6, 0],
-            [3, 4, 8]
-        ];
+        let a = array![[7, 1, 5], [2, 6, 0], [3, 4, 8]];
 
         // Make array non-contiguous by transposing
         let a = a.t().to_owned();
@@ -3467,12 +3547,69 @@ mod tests
         // - First column should be <= middle column
         // - Last column should be >= middle column
         for row in 0..3 {
-            assert!(p_axis1[[row, 0]] <= p_axis1[[row, 1]],
+            assert!(
+                p_axis1[[row, 0]] <= p_axis1[[row, 1]],
                 "Row {}: First column {} should be <= middle column {}",
-                row, p_axis1[[row, 0]], p_axis1[[row, 1]]);
-            assert!(p_axis1[[row, 2]] >= p_axis1[[row, 1]],
+                row,
+                p_axis1[[row, 0]],
+                p_axis1[[row, 1]]
+            );
+            assert!(
+                p_axis1[[row, 2]] >= p_axis1[[row, 1]],
                 "Row {}: Last column {} should be >= middle column {}",
-                row, p_axis1[[row, 2]], p_axis1[[row, 1]]);
+                row,
+                p_axis1[[row, 2]],
+                p_axis1[[row, 1]]
+            );
         }
     }
+
+    #[test]
+    fn test_partition_empty()
+    {
+        // Test 1D empty array
+        let empty1d = Array1::<i32>::zeros(0);
+        let result1d = empty1d.partition(0, Axis(0));
+        assert_eq!(result1d.len(), 0);
+
+        // Test 1D empty array with kth out of bounds
+        let result1d_out_of_bounds = empty1d.partition(1, Axis(0));
+        assert_eq!(result1d_out_of_bounds.len(), 0);
+
+        // Test 2D empty array
+        let empty2d = Array2::<i32>::zeros((0, 3));
+        let result2d = empty2d.partition(0, Axis(0));
+        assert_eq!(result2d.shape(), &[0, 3]);
+
+        // Test 2D empty array with zero columns
+        let empty2d_cols = Array2::<i32>::zeros((2, 0));
+        let result2d_cols = empty2d_cols.partition(0, Axis(1));
+        assert_eq!(result2d_cols.shape(), &[2, 0]);
+
+        // Test 3D empty array
+        let empty3d = Array3::<i32>::zeros((0, 2, 3));
+        let result3d = empty3d.partition(0, Axis(0));
+        assert_eq!(result3d.shape(), &[0, 2, 3]);
+
+        // Test 3D empty array with zero in middle dimension
+        let empty3d_mid = Array3::<i32>::zeros((2, 0, 3));
+        let result3d_mid = empty3d_mid.partition(0, Axis(1));
+        assert_eq!(result3d_mid.shape(), &[2, 0, 3]);
+
+        // Test 4D empty array
+        let empty4d = Array4::<i32>::zeros((0, 2, 3, 4));
+        let result4d = empty4d.partition(0, Axis(0));
+        assert_eq!(result4d.shape(), &[0, 2, 3, 4]);
+
+        // Test empty array with non-zero dimensions in other axes
+        let empty_mixed = Array2::<i32>::zeros((0, 5));
+        let result_mixed = empty_mixed.partition(0, Axis(0));
+        assert_eq!(result_mixed.shape(), &[0, 5]);
+
+        // Test empty array with negative strides
+        let arr = Array2::<i32>::zeros((3, 3));
+        let empty_slice = arr.slice(s![0..0, ..]);
+        let result_slice = empty_slice.partition(0, Axis(0));
+        assert_eq!(result_slice.shape(), &[0, 3]);
+    }
 }
diff --git a/src/iterators/windows.rs b/src/iterators/windows.rs
index afdaaa89..f3442c0a 100644
--- a/src/iterators/windows.rs
+++ b/src/iterators/windows.rs
@@ -141,7 +141,7 @@ pub struct AxisWindows<'a, A, D>
 
 impl<'a, A, D: Dimension> AxisWindows<'a, A, D>
 {
-    pub(crate) fn new(a: ArrayView<'a, A, D>, axis: Axis, window_size: usize) -> Self
+    pub(crate) fn new_with_stride(a: ArrayView<'a, A, D>, axis: Axis, window_size: usize, stride_size: usize) -> Self
     {
         let window_strides = a.strides.clone();
         let axis_idx = axis.index();
@@ -150,10 +150,11 @@ impl<'a, A, D: Dimension> AxisWindows<'a, A, D>
         window[axis_idx] = window_size;
 
         let ndim = window.ndim();
-        let mut unit_stride = D::zeros(ndim);
-        unit_stride.slice_mut().fill(1);
+        let mut stride = D::zeros(ndim);
+        stride.slice_mut().fill(1);
+        stride[axis_idx] = stride_size;
 
-        let base = build_base(a, window.clone(), unit_stride);
+        let base = build_base(a, window.clone(), stride);
         AxisWindows {
             base,
             axis_idx,


Test patch:
diff --git a/tests/array.rs b/tests/array.rs
index f1426625..3d6fa671 100644
--- a/tests/array.rs
+++ b/tests/array.rs
@@ -2828,3 +2828,81 @@ fn test_slice_assign()
     *a.slice_mut(s![1..3]) += 1;
     assert_eq!(a, array![0, 2, 3, 3, 4]);
 }
+
+#[test]
+fn reverse_axes()
+{
+    let mut a = arr2(&[[1, 2], [3, 4]]);
+    a.reverse_axes();
+    assert_eq!(a, arr2(&[[1, 3], [2, 4]]));
+
+    let mut a = arr2(&[[1, 2, 3], [4, 5, 6]]);
+    a.reverse_axes();
+    assert_eq!(a, arr2(&[[1, 4], [2, 5], [3, 6]]));
+
+    let mut a = Array::from_iter(0..24)
+        .into_shape_with_order((2, 3, 4))
+        .unwrap();
+    let original = a.clone();
+    a.reverse_axes();
+    for ((i0, i1, i2), elem) in original.indexed_iter() {
+        assert_eq!(*elem, a[(i2, i1, i0)]);
+    }
+}
+
+#[test]
+fn permute_axes()
+{
+    let mut a = arr2(&[[1, 2], [3, 4]]);
+    a.permute_axes([1, 0]);
+    assert_eq!(a, arr2(&[[1, 3], [2, 4]]));
+
+    let mut a = Array::from_iter(0..24)
+        .into_shape_with_order((2, 3, 4))
+        .unwrap();
+    let original = a.clone();
+    a.permute_axes([2, 1, 0]);
+    for ((i0, i1, i2), elem) in original.indexed_iter() {
+        assert_eq!(*elem, a[(i2, i1, i0)]);
+    }
+
+    let mut a = Array::from_iter(0..120)
+        .into_shape_with_order((2, 3, 4, 5))
+        .unwrap();
+    let original = a.clone();
+    a.permute_axes([1, 0, 3, 2]);
+    for ((i0, i1, i2, i3), elem) in original.indexed_iter() {
+        assert_eq!(*elem, a[(i1, i0, i3, i2)]);
+    }
+}
+
+#[should_panic]
+#[test]
+fn permute_axes_repeated_axis()
+{
+    let mut a = Array::from_iter(0..24)
+        .into_shape_with_order((2, 3, 4))
+        .unwrap();
+    a.permute_axes([1, 0, 1]);
+}
+
+#[should_panic]
+#[test]
+fn permute_axes_missing_axis()
+{
+    let mut a = Array::from_iter(0..24)
+        .into_shape_with_order((2, 3, 4))
+        .unwrap()
+        .into_dyn();
+    a.permute_axes(&[2, 0][..]);
+}
+
+#[should_panic]
+#[test]
+fn permute_axes_oob()
+{
+    let mut a = Array::from_iter(0..24)
+        .into_shape_with_order((2, 3, 4))
+        .unwrap();
+    a.permute_axes([1, 0, 3]);
+}
diff --git a/tests/windows.rs b/tests/windows.rs
index 6506d830..4d4d0d7d 100644
--- a/tests/windows.rs
+++ b/tests/windows.rs
@@ -294,6 +294,148 @@ fn tests_axis_windows_3d_zips_with_1d()
     assert_eq!(b,arr1(&[207, 261]));
 }
 
+/// Test verifies that non existent Axis results in panic
+#[test]
+#[should_panic]
+fn axis_windows_with_stride_outofbound()
+{
+    let a = Array::from_iter(10..37)
+        .into_shape_with_order((3, 3, 3))
+        .unwrap();
+    a.axis_windows_with_stride(Axis(4), 2, 2);
+}
+
+/// Test verifies that zero sizes results in panic
+#[test]
+#[should_panic]
+fn axis_windows_with_stride_zero_size()
+{
+    let a = Array::from_iter(10..37)
+        .into_shape_with_order((3, 3, 3))
+        .unwrap();
+    a.axis_windows_with_stride(Axis(0), 0, 2);
+}
+
+/// Test verifies that zero stride results in panic
+#[test]
+#[should_panic]
+fn axis_windows_with_stride_zero_stride()
+{
+    let a = Array::from_iter(10..37)
+        .into_shape_with_order((3, 3, 3))
+        .unwrap();
+    a.axis_windows_with_stride(Axis(0), 2, 0);
+}
+
+/// Test verifies that over sized windows yield nothing
+#[test]
+fn axis_windows_with_stride_oversized()
+{
+    let a = Array::from_iter(10..37)
+        .into_shape_with_order((3, 3, 3))
+        .unwrap();
+    let mut iter = a.axis_windows_with_stride(Axis(2), 4, 2).into_iter();
+    assert_eq!(iter.next(), None);
+}
+
+/// Simple test for iterating 1d-arrays via `Axis Windows`.
+#[test]
+fn test_axis_windows_with_stride_1d()
+{
+    let a = Array::from_iter(10..20).into_shape_with_order(10).unwrap();
+
+    itertools::assert_equal(a.axis_windows_with_stride(Axis(0), 5, 2), vec![
+        arr1(&[10, 11, 12, 13, 14]),
+        arr1(&[12, 13, 14, 15, 16]),
+        arr1(&[14, 15, 16, 17, 18]),
+    ]);
+
+    itertools::assert_equal(a.axis_windows_with_stride(Axis(0), 5, 3), vec![
+        arr1(&[10, 11, 12, 13, 14]),
+        arr1(&[13, 14, 15, 16, 17]),
+    ]);
+}
+
+/// Simple test for iterating 2d-arrays via `Axis Windows`.
+#[test]
+fn test_axis_windows_with_stride_2d()
+{
+    let a = Array::from_iter(10..30)
+        .into_shape_with_order((5, 4))
+        .unwrap();
+
+    itertools::assert_equal(a.axis_windows_with_stride(Axis(0), 2, 1), vec![
+        arr2(&[[10, 11, 12, 13], [14, 15, 16, 17]]),
+        arr2(&[[14, 15, 16, 17], [18, 19, 20, 21]]),
+        arr2(&[[18, 19, 20, 21], [22, 23, 24, 25]]),
+        arr2(&[[22, 23, 24, 25], [26, 27, 28, 29]]),
+    ]);
+
+    itertools::assert_equal(a.axis_windows_with_stride(Axis(0), 2, 2), vec![
+        arr2(&[[10, 11, 12, 13], [14, 15, 16, 17]]),
+        arr2(&[[18, 19, 20, 21], [22, 23, 24, 25]]),
+    ]);
+
+    itertools::assert_equal(a.axis_windows_with_stride(Axis(0), 2, 3), vec![
+        arr2(&[[10, 11, 12, 13], [14, 15, 16, 17]]),
+        arr2(&[[22, 23, 24, 25], [26, 27, 28, 29]]),
+    ]);
+}
+
+/// Simple test for iterating 3d-arrays via `Axis Windows`.
+#[test]
+fn test_axis_windows_with_stride_3d()
+{
+    let a = Array::from_iter(0..27)
+        .into_shape_with_order((3, 3, 3))
+        .unwrap();
+
+    itertools::assert_equal(a.axis_windows_with_stride(Axis(1), 2, 1), vec![
+            arr3(&[
+                [[0, 1, 2], [3, 4, 5]],
+                [[9, 10, 11], [12, 13, 14]],
+                [[18, 19, 20], [21, 22, 23]],
+            ]),
+            arr3(&[
+                [[3, 4, 5], [6, 7, 8]],
+                [[12, 13, 14], [15, 16, 17]],
+                [[21, 22, 23], [24, 25, 26]],
+            ]),
+        ]);
+
+    itertools::assert_equal(a.axis_windows_with_stride(Axis(1), 2, 2), vec![
+            arr3(&[
+                [[0, 1, 2], [3, 4, 5]],
+                [[9, 10, 11], [12, 13, 14]],
+                [[18, 19, 20], [21, 22, 23]],
+            ]),
+        ]);
+}
+
+#[test]
+fn tests_axis_windows_with_stride_3d_zips_with_1d()
+{
+    let a = Array::from_iter(0..27)
+        .into_shape_with_order((3, 3, 3))
+        .unwrap();
+    let mut b1 = Array::zeros(2);
+    let mut b2 = Array::zeros(1);
+
+    Zip::from(b1.view_mut())
+        .and(a.axis_windows_with_stride(Axis(1), 2, 1))
+        .for_each(|b, a| {
+            *b = a.sum();
+        });
+    assert_eq!(b1,arr1(&[207, 261]));
+
+    Zip::from(b2.view_mut())
+        .and(a.axis_windows_with_stride(Axis(1), 2, 2))
+        .for_each(|b, a| {
+            *b = a.sum();
+        });
+    assert_eq!(b2,arr1(&[207]));
+}
+
 #[test]
 fn test_window_neg_stride()
 {