Modulo (#734)

* add modulo operator

* fix module loops

* add in-place modulo operator

* update readme

* add test files, update documentation

.github/workflows/build.yml

@@ -20,7 +20,7 @@ jobs:
     strategy:
         matrix:
             os:
-                - ubuntu-20.04
+                - ubuntu-24.04
                 - macOS-latest
             dims: [1, 2, 3, 4]
     runs-on: ${{ matrix.os }}
@@ -61,7 +61,7 @@ jobs:
     strategy:
         matrix:
             os:
-                - ubuntu-20.04
+                - ubuntu-24.04
                 - macOS-latest
             dims: [1, 2, 3, 4]
     runs-on: ${{ matrix.os }}

README.md

@@ -36,7 +36,7 @@ detected and handled.
 ## ndarray methods
 
 `ulab` implements `numpy`'s `ndarray` with the `==`, `!=`, `<`, `<=`, `>`, `>=`, `+`, `-`, `/`, `*`, `**`,
-`+=`, `-=`, `*=`, `/=`, `**=` binary operators, and the `len`, `~`, `-`, `+`, `abs` unary operators that
+`%`, `+=`, `-=`, `*=`, `/=`, `**=`, `%=` binary operators, and the `len`, `~`, `-`, `+`, `abs` unary operators that
 operate element-wise. Type-aware `ndarray`s can be initialised from any `micropython` iterable, lists of
 iterables via the `array` constructor, or by means of the `arange`, `concatenate`, `diag`, `eye`,
 `frombuffer`, `full`, `linspace`, `logspace`, `ones`, or `zeros`  functions.
@@ -112,6 +112,7 @@ of the user manual.
 1. `MaixPy` https://github.com/sipeed/MaixPy
 1. `OpenMV` https://github.com/openmv/openmv
 1. `pimoroni-pico` https://github.com/pimoroni/pimoroni-pico
+1. `Tulip Creative Computer` https://github.com/shorepine/tulipcc 
 
 ## Compiling
 

code/ndarray.c

@@ -563,6 +563,9 @@ ndarray_obj_t *ndarray_new_dense_ndarray(uint8_t ndim, size_t *shape, uint8_t dt
 ndarray_obj_t *ndarray_new_ndarray_from_tuple(mp_obj_tuple_t *_shape, uint8_t dtype) {
     // creates a dense array from a tuple
     // the function should work in the general n-dimensional case
+    if(_shape->len > ULAB_MAX_DIMS) {
+        mp_raise_ValueError(MP_ERROR_TEXT("maximum number of dimensions is " MP_STRINGIFY(ULAB_MAX_DIMS)));
+    }
     size_t *shape = m_new0(size_t, ULAB_MAX_DIMS);
     for(size_t i = 0; i < _shape->len; i++) {
         shape[ULAB_MAX_DIMS - 1 - i] = mp_obj_get_int(_shape->items[_shape->len - 1 - i]);
@@ -583,43 +586,10 @@ void ndarray_copy_array(ndarray_obj_t *source, ndarray_obj_t *target, uint8_t sh
     }
     #endif
 
-    #if ULAB_MAX_DIMS > 3
+    ITERATOR_HEAD();
-    size_t i = 0;
-    do {
-    #endif
-        #if ULAB_MAX_DIMS > 2
-        size_t j = 0;
-        do {
-        #endif
-            #if ULAB_MAX_DIMS > 1
-            size_t k = 0;
-            do {
-            #endif
-                size_t l = 0;
-                do {
         memcpy(tarray, sarray, target->itemsize);
         tarray += target->itemsize;
-                    sarray += source->strides[ULAB_MAX_DIMS - 1];
+    ITERATOR_TAIL(source, sarray);
-                    l++;
-                } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-            #if ULAB_MAX_DIMS > 1
-                sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS-1];
-                sarray += source->strides[ULAB_MAX_DIMS - 2];
-                k++;
-            } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-            #endif
-        #if ULAB_MAX_DIMS > 2
-            sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS-2];
-            sarray += source->strides[ULAB_MAX_DIMS - 3];
-            j++;
-        } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-        #endif
-    #if ULAB_MAX_DIMS > 3
-        sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS-3];
-        sarray += source->strides[ULAB_MAX_DIMS - 4];
-        i++;
-    } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-    #endif
 }
 
 ndarray_obj_t *ndarray_new_view(ndarray_obj_t *source, uint8_t ndim, size_t *shape, int32_t *strides, int32_t offset) {
@@ -673,20 +643,7 @@ ndarray_obj_t *ndarray_copy_view_convert_type(ndarray_obj_t *source, uint8_t dty
     uint8_t complex_size = 2 * sizeof(mp_float_t);
     #endif
 
-    #if ULAB_MAX_DIMS > 3
+    ITERATOR_HEAD()
-    size_t i = 0;
-    do {
-    #endif
-        #if ULAB_MAX_DIMS > 2
-        size_t j = 0;
-        do {
-        #endif
-            #if ULAB_MAX_DIMS > 1
-            size_t k = 0;
-            do {
-            #endif
-                size_t l = 0;
-                do {
         mp_obj_t item;
         #if ULAB_SUPPORTS_COMPLEX
         if(source->dtype == NDARRAY_COMPLEX) {
@@ -715,27 +672,7 @@ ndarray_obj_t *ndarray_copy_view_convert_type(ndarray_obj_t *source, uint8_t dty
         ndarray_set_value(dtype, array, 0, item);
         #endif
         array += ndarray->itemsize;
-                    sarray += source->strides[ULAB_MAX_DIMS - 1];
+    ITERATOR_TAIL(source, sarray);
-                    l++;
-                } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-            #if ULAB_MAX_DIMS > 1
-                sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS-1];
-                sarray += source->strides[ULAB_MAX_DIMS - 2];
-                k++;
-            } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-            #endif
-        #if ULAB_MAX_DIMS > 2
-            sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS-2];
-            sarray += source->strides[ULAB_MAX_DIMS - 3];
-            j++;
-        } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-        #endif
-    #if ULAB_MAX_DIMS > 3
-        sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS-3];
-        sarray += source->strides[ULAB_MAX_DIMS - 4];
-        i++;
-    } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-    #endif
     return ndarray;
 }
 
@@ -762,20 +699,7 @@ mp_obj_t ndarray_byteswap(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_
         return MP_OBJ_FROM_PTR(ndarray);
     } else {
         uint8_t *array = (uint8_t *)ndarray->array;
-        #if ULAB_MAX_DIMS > 3
+        ITERATOR_HEAD();
-        size_t i = 0;
-        do {
-        #endif
-            #if ULAB_MAX_DIMS > 2
-            size_t j = 0;
-            do {
-            #endif
-                #if ULAB_MAX_DIMS > 1
-                size_t k = 0;
-                do {
-                #endif
-                    size_t l = 0;
-                    do {
             if(self->dtype == NDARRAY_FLOAT) {
                 #if MICROPY_FLOAT_IMPL == MICROPY_FLOAT_IMPL_FLOAT
                 SWAP(uint8_t, array[0], array[3]);
@@ -789,27 +713,7 @@ mp_obj_t ndarray_byteswap(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_
             } else {
                 SWAP(uint8_t, array[0], array[1]);
             }
-                        array += ndarray->strides[ULAB_MAX_DIMS - 1];
+        ITERATOR_TAIL(ndarray, array);
-                        l++;
-                    } while(l < ndarray->shape[ULAB_MAX_DIMS - 1]);
-                #if ULAB_MAX_DIMS > 1
-                    array -= ndarray->strides[ULAB_MAX_DIMS - 1] * ndarray->shape[ULAB_MAX_DIMS-1];
-                    array += ndarray->strides[ULAB_MAX_DIMS - 2];
-                    k++;
-                } while(k < ndarray->shape[ULAB_MAX_DIMS - 2]);
-                #endif
-            #if ULAB_MAX_DIMS > 2
-                array -= ndarray->strides[ULAB_MAX_DIMS - 2] * ndarray->shape[ULAB_MAX_DIMS-2];
-                array += ndarray->strides[ULAB_MAX_DIMS - 3];
-                j++;
-            } while(j < ndarray->shape[ULAB_MAX_DIMS - 3]);
-            #endif
-        #if ULAB_MAX_DIMS > 3
-            array -= ndarray->strides[ULAB_MAX_DIMS - 3] * ndarray->shape[ULAB_MAX_DIMS-3];
-            array += ndarray->strides[ULAB_MAX_DIMS - 4];
-            i++;
-        } while(i < ndarray->shape[ULAB_MAX_DIMS - 4]);
-        #endif
     }
     return MP_OBJ_FROM_PTR(ndarray);
 }
@@ -1438,43 +1342,10 @@ mp_obj_t ndarray_flatten(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_a
     uint8_t *array = (uint8_t *)ndarray->array;
 
     if(memcmp(order, "C", 1) == 0) { // C-type ordering
-        #if ULAB_MAX_DIMS > 3
+        ITERATOR_HEAD();
-        size_t i = 0;
-        do {
-        #endif
-            #if ULAB_MAX_DIMS > 2
-            size_t j = 0;
-            do {
-            #endif
-                #if ULAB_MAX_DIMS > 1
-                size_t k = 0;
-                do {
-                #endif
-                    size_t l = 0;
-                    do {
             memcpy(array, sarray, self->itemsize);
             array += ndarray->strides[ULAB_MAX_DIMS - 1];
-                        sarray += self->strides[ULAB_MAX_DIMS - 1];
+        ITERATOR_TAIL(self, sarray);
-                        l++;
-                    } while(l <  self->shape[ULAB_MAX_DIMS - 1]);
-                #if ULAB_MAX_DIMS > 1
-                    sarray -= self->strides[ULAB_MAX_DIMS - 1] * self->shape[ULAB_MAX_DIMS-1];
-                    sarray += self->strides[ULAB_MAX_DIMS - 2];
-                    k++;
-                } while(k <  self->shape[ULAB_MAX_DIMS - 2]);
-                #endif
-            #if ULAB_MAX_DIMS > 2
-                sarray -= self->strides[ULAB_MAX_DIMS - 2] * self->shape[ULAB_MAX_DIMS-2];
-                sarray += self->strides[ULAB_MAX_DIMS - 3];
-                j++;
-            } while(j <  self->shape[ULAB_MAX_DIMS - 3]);
-            #endif
-        #if ULAB_MAX_DIMS > 3
-            sarray -= self->strides[ULAB_MAX_DIMS - 3] * self->shape[ULAB_MAX_DIMS-3];
-            sarray += self->strides[ULAB_MAX_DIMS - 4];
-            i++;
-        } while(i <  self->shape[ULAB_MAX_DIMS - 4]);
-        #endif
     } else { // 'F', Fortran-type ordering
         #if ULAB_MAX_DIMS > 3
         size_t i = 0;
@@ -1527,6 +1398,13 @@ mp_obj_t ndarray_itemsize(mp_obj_t self_in) {
 }
 #endif
 
+#if NDARRAY_HAS_NDIM
+mp_obj_t ndarray_ndim(mp_obj_t self_in) {
+    ndarray_obj_t *self = MP_OBJ_TO_PTR(self_in);
+    return MP_OBJ_NEW_SMALL_INT(self->ndim);
+}
+#endif
+
 #if NDARRAY_HAS_SHAPE
 mp_obj_t ndarray_shape(mp_obj_t self_in) {
     ndarray_obj_t *self = MP_OBJ_TO_PTR(self_in);
@@ -1611,7 +1489,7 @@ ndarray_obj_t *ndarray_from_mp_obj(mp_obj_t obj, uint8_t other_type) {
 
     if(mp_obj_is_int(obj)) {
         int32_t ivalue = mp_obj_get_int(obj);
-        if((ivalue < -32767) || (ivalue > 32767)) {
+        if((ivalue < -32768) || (ivalue > 65535)) {
             // the integer value clearly does not fit the ulab integer types, so move on to float
             ndarray = ndarray_new_linear_array(1, NDARRAY_FLOAT);
             mp_float_t *array = (mp_float_t *)ndarray->array;
@@ -1619,7 +1497,7 @@ ndarray_obj_t *ndarray_from_mp_obj(mp_obj_t obj, uint8_t other_type) {
         } else {
             uint8_t dtype;
             if(ivalue < 0) {
-                if(ivalue > -128) {
+                if(ivalue >= -128) {
                     dtype = NDARRAY_INT8;
                 } else {
                     dtype = NDARRAY_INT16;
@@ -1770,6 +1648,12 @@ mp_obj_t ndarray_binary_op(mp_binary_op_t _op, mp_obj_t lobj, mp_obj_t robj) {
             return ndarray_inplace_ams(lhs, rhs, rstrides, op);
             break;
         #endif
+        #if NDARRAY_HAS_INPLACE_MODULO
+        case MP_BINARY_OP_INPLACE_MODULO:
+            COMPLEX_DTYPE_NOT_IMPLEMENTED(lhs->dtype);
+            return ndarray_inplace_modulo(lhs, rhs, rstrides);
+            break;
+        #endif
         #if NDARRAY_HAS_INPLACE_MULTIPLY
         case MP_BINARY_OP_INPLACE_MULTIPLY:
             COMPLEX_DTYPE_NOT_IMPLEMENTED(lhs->dtype);
@@ -1825,6 +1709,12 @@ mp_obj_t ndarray_binary_op(mp_binary_op_t _op, mp_obj_t lobj, mp_obj_t robj) {
             return ndarray_binary_add(lhs, rhs, ndim, shape, lstrides, rstrides);
             break;
         #endif
+        #if NDARRAY_HAS_BINARY_OP_MODULO
+        case MP_BINARY_OP_MODULO:
+            COMPLEX_DTYPE_NOT_IMPLEMENTED(lhs->dtype);
+            return ndarray_binary_modulo(lhs, rhs, ndim, shape, lstrides, rstrides);
+            break;
+        #endif
         #if NDARRAY_HAS_BINARY_OP_MULTIPLY
         case MP_BINARY_OP_MULTIPLY:
             return ndarray_binary_multiply(lhs, rhs, ndim, shape, lstrides, rstrides);

code/ndarray.h

@@ -232,6 +232,10 @@ mp_obj_t ndarray_dtype(mp_obj_t );
 mp_obj_t ndarray_itemsize(mp_obj_t );
 #endif
 
+#if NDARRAY_HAS_NDIM
+mp_obj_t ndarray_ndim(mp_obj_t );
+#endif
+
 #if NDARRAY_HAS_SIZE
 mp_obj_t ndarray_size(mp_obj_t );
 #endif
@@ -709,4 +713,89 @@ ndarray_obj_t *ndarray_from_mp_obj(mp_obj_t , uint8_t );
 #endif /* ULAB_MAX_DIMS == 4 */
 #endif /* ULAB_HAS_FUNCTION_ITERATOR */
 
+
+// iterator macro for traversing arrays over all dimensions
+#if ULAB_MAX_DIMS == 1
+#define ITERATOR_HEAD()\
+    size_t _l_ = 0;\
+    do {
+
+#define ITERATOR_TAIL(_source_, _source_array_)\
+    (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 1];\
+    _l_++;\
+    } while(_l_ < (_source_)->shape[ULAB_MAX_DIMS - 1]);
+
+#endif /* ULAB_MAX_DIMS == 1 */
+
+#if ULAB_MAX_DIMS == 2
+#define ITERATOR_HEAD()\
+    size_t _k_ = 0;\
+    do {\
+        size_t _l_ = 0;\
+        do {
+
+#define ITERATOR_TAIL(_source_, _source_array_)\
+            (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 1];\
+            _l_++;\
+        } while(_l_ < (_source_)->shape[ULAB_MAX_DIMS - 1]);\
+        (_source_array_) -= (_source_)->strides[ULAB_MAX_DIMS - 1] * (_source_)->shape[ULAB_MAX_DIMS - 1];\
+        (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 2];\
+        _k_++;\
+    } while(_k_ < (_source_)->shape[ULAB_MAX_DIMS - 2]);
+#endif /* ULAB_MAX_DIMS == 2 */
+
+#if ULAB_MAX_DIMS == 3
+#define ITERATOR_HEAD()\
+    size_t _j_ = 0;\
+    do {\
+        size_t _k_ = 0;\
+        do {\
+            size_t _l_ = 0;\
+            do {
+
+#define ITERATOR_TAIL(_source_, _source_array_)\
+                (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 1];\
+                _l_++;\
+            } while(_l_ < (_source_)->shape[ULAB_MAX_DIMS - 1]);\
+            (_source_array_) -= (_source_)->strides[ULAB_MAX_DIMS - 1] * (_source_)->shape[ULAB_MAX_DIMS - 1];\
+            (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 2];\
+            _k_++;\
+        } while(_k_ < (_source_)->shape[ULAB_MAX_DIMS - 2]);\
+        (_source_array_) -= (_source_)->strides[ULAB_MAX_DIMS - 2] * (_source_)->shape[ULAB_MAX_DIMS - 2];\
+        (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 3];\
+        _j_++;\
+    } while(_j_ < (_source_)->shape[ULAB_MAX_DIMS - 3]);
+
+#endif /* ULAB_MAX_DIMS == 3 */
+
+#if ULAB_MAX_DIMS == 4
+#define ITERATOR_HEAD()\
+    size_t _i_ = 0;\
+    do {\
+        size_t _j_ = 0;\
+        do {\
+            size_t _k_ = 0;\
+            do {\
+                size_t _l_ = 0;\
+                do {
+
+#define ITERATOR_TAIL(_source_, _source_array_)\
+                    (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 1];\
+                    _l_++;\
+                } while(_l_ < (_source_)->shape[ULAB_MAX_DIMS - 1]);\
+                (_source_array_) -= (_source_)->strides[ULAB_MAX_DIMS - 1] * (_source_)->shape[ULAB_MAX_DIMS - 1];\
+                (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 2];\
+                _k_++;\
+            } while(_k_ < (_source_)->shape[ULAB_MAX_DIMS - 2]);\
+            (_source_array_) -= (_source_)->strides[ULAB_MAX_DIMS - 2] * (_source_)->shape[ULAB_MAX_DIMS - 2];\
+            (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 3];\
+            _j_++;\
+        } while(_j_ < (_source_)->shape[ULAB_MAX_DIMS - 3]);\
+        (_source_array_) -= (_source_)->strides[ULAB_MAX_DIMS - 3] * (_source_)->shape[ULAB_MAX_DIMS - 3];\
+        (_source_array_) += (_source_)->strides[ULAB_MAX_DIMS - 4];\
+        _i_++;\
+    } while(_i_ < (_source_)->shape[ULAB_MAX_DIMS - 4]);
+#endif /* ULAB_MAX_DIMS == 4 */
+
+
 #endif

code/ndarray_operators.c

@@ -248,6 +248,105 @@ mp_obj_t ndarray_binary_add(ndarray_obj_t *lhs, ndarray_obj_t *rhs,
 }
 #endif /* NDARRAY_HAS_BINARY_OP_ADD */
 
+#if NDARRAY_HAS_BINARY_OP_MODULO
+mp_obj_t ndarray_binary_modulo(ndarray_obj_t *lhs, ndarray_obj_t *rhs,
+                                        uint8_t ndim, size_t *shape, int32_t *lstrides, int32_t *rstrides) {
+
+    ndarray_obj_t *results = NULL;
+    uint8_t *larray = (uint8_t *)lhs->array;
+    uint8_t *rarray = (uint8_t *)rhs->array;
+
+    if(lhs->dtype == NDARRAY_UINT8) {
+        if(rhs->dtype == NDARRAY_UINT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_UINT8);
+            BINARY_LOOP(results, uint8_t, uint8_t, uint8_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_INT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_INT16);
+            BINARY_LOOP(results, int16_t, uint8_t, int8_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_UINT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_UINT16);
+            BINARY_LOOP(results, uint16_t, uint8_t, uint16_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_INT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_INT16);
+            BINARY_LOOP(results, int16_t, uint8_t, int16_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_FLOAT) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            MODULO_FLOAT_LOOP(results, mp_float_t, uint8_t, mp_float_t, larray, lstrides, rarray, rstrides);
+        }
+    } else if(lhs->dtype == NDARRAY_INT8) {
+        if(rhs->dtype == NDARRAY_UINT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_INT16);
+            BINARY_LOOP(results, int16_t, int8_t, uint8_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_INT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_INT8);
+            BINARY_LOOP(results, int8_t, int8_t, int8_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_UINT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_INT16);
+            BINARY_LOOP(results, int16_t, int8_t, int16_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_INT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_INT16);
+            BINARY_LOOP(results, int16_t, int8_t, int16_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_FLOAT) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            MODULO_FLOAT_LOOP(results, mp_float_t, int8_t, mp_float_t, larray, lstrides, rarray, rstrides);
+        }
+    } else if(lhs->dtype == NDARRAY_UINT16) {
+        if(rhs->dtype == NDARRAY_UINT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_UINT8);
+            BINARY_LOOP(results, uint16_t, uint16_t, uint8_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_INT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            BINARY_LOOP(results, mp_float_t, uint16_t, int8_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_UINT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_UINT16);
+            BINARY_LOOP(results, uint16_t, uint16_t, uint16_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_INT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            BINARY_LOOP(results, mp_float_t, uint16_t, int16_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_FLOAT) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            MODULO_FLOAT_LOOP(results, mp_float_t, uint16_t, mp_float_t, larray, lstrides, rarray, rstrides);
+        }
+    } else if(lhs->dtype == NDARRAY_INT16) {
+        if(rhs->dtype == NDARRAY_UINT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_INT16);
+            BINARY_LOOP(results, int16_t, int16_t, uint8_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_INT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_INT16);
+            BINARY_LOOP(results, int16_t, int16_t, int8_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_UINT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            BINARY_LOOP(results, mp_float_t, int16_t, uint16_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_INT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_INT16);
+            BINARY_LOOP(results, int16_t, int16_t, int16_t, larray, lstrides, rarray, rstrides, %);
+        } else if(rhs->dtype == NDARRAY_FLOAT) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            MODULO_FLOAT_LOOP(results, mp_float_t, int16_t, mp_float_t, larray, lstrides, rarray, rstrides);
+        }
+    } else if(lhs->dtype == NDARRAY_FLOAT) {
+        if(rhs->dtype == NDARRAY_UINT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            MODULO_FLOAT_LOOP(results, mp_float_t, mp_float_t, uint8_t, larray, lstrides, rarray, rstrides);
+        } else if(rhs->dtype == NDARRAY_INT8) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            MODULO_FLOAT_LOOP(results, mp_float_t, mp_float_t, int8_t, larray, lstrides, rarray, rstrides);
+        } else if(rhs->dtype == NDARRAY_UINT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            MODULO_FLOAT_LOOP(results, mp_float_t, mp_float_t, uint16_t, larray, lstrides, rarray, rstrides);
+        } else if(rhs->dtype == NDARRAY_INT16) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            MODULO_FLOAT_LOOP(results, mp_float_t, mp_float_t, int16_t, larray, lstrides, rarray, rstrides);
+        } else if(rhs->dtype == NDARRAY_FLOAT) {
+            results = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            MODULO_FLOAT_LOOP(results, mp_float_t, mp_float_t, mp_float_t, larray, lstrides, rarray, rstrides);
+        }
+    }
+
+    return MP_OBJ_FROM_PTR(results);
+}
+#endif /* NDARRAY_HAS_BINARY_OP_MODULO */
+
 #if NDARRAY_HAS_BINARY_OP_MULTIPLY
 mp_obj_t ndarray_binary_multiply(ndarray_obj_t *lhs, ndarray_obj_t *rhs,
                                             uint8_t ndim, size_t *shape, int32_t *lstrides, int32_t *rstrides) {
@@ -1074,6 +1173,29 @@ mp_obj_t ndarray_inplace_ams(ndarray_obj_t *lhs, ndarray_obj_t *rhs, int32_t *rs
 }
 #endif /* NDARRAY_HAS_INPLACE_ADD || NDARRAY_HAS_INPLACE_MULTIPLY || NDARRAY_HAS_INPLACE_SUBTRACT */
 
+
+#if NDARRAY_HAS_INPLACE_MODULO
+mp_obj_t ndarray_inplace_modulo(ndarray_obj_t *lhs, ndarray_obj_t *rhs, int32_t *rstrides) {
+    if((lhs->dtype != NDARRAY_FLOAT) && (rhs->dtype == NDARRAY_FLOAT)) {
+        mp_raise_TypeError(MP_ERROR_TEXT("results cannot be cast to specified type"));
+    }
+    if(lhs->dtype == NDARRAY_FLOAT) {
+        if(rhs->dtype == NDARRAY_UINT8) {
+            INLINE_MODULO_FLOAT_LOOP(lhs, uint8_t, larray, rarray, rstrides);
+        } else if(rhs->dtype == NDARRAY_UINT8) {
+            INLINE_MODULO_FLOAT_LOOP(lhs, int8_t, larray, rarray, rstrides);
+        } else if(rhs->dtype == NDARRAY_UINT16) {
+            INLINE_MODULO_FLOAT_LOOP(lhs, uint16_t, larray, rarray, rstrides);
+        } else if(rhs->dtype == NDARRAY_INT16) {
+            INLINE_MODULO_FLOAT_LOOP(lhs, int16_t, larray, rarray, rstrides);
+        } else {
+            INLINE_MODULO_FLOAT_LOOP(lhs, mp_float_t, larray, rarray, rstrides);
+        }
+    }
+    return MP_OBJ_FROM_PTR(lhs);
+}
+#endif /* NDARRAY_HAS_INPLACE_MODULO */
+
 #if NDARRAY_HAS_INPLACE_TRUE_DIVIDE
 mp_obj_t ndarray_inplace_divide(ndarray_obj_t *lhs, ndarray_obj_t *rhs, int32_t *rstrides) {
 

code/ndarray_operators.h

@@ -12,6 +12,7 @@
 
 mp_obj_t ndarray_binary_equality(ndarray_obj_t *, ndarray_obj_t *, uint8_t , size_t *,  int32_t *, int32_t *, mp_binary_op_t );
 mp_obj_t ndarray_binary_add(ndarray_obj_t *, ndarray_obj_t *, uint8_t , size_t *, int32_t *, int32_t *);
+mp_obj_t ndarray_binary_modulo(ndarray_obj_t *, ndarray_obj_t *, uint8_t , size_t *, int32_t *, int32_t *);
 mp_obj_t ndarray_binary_multiply(ndarray_obj_t *, ndarray_obj_t *, uint8_t , size_t *, int32_t *, int32_t *);
 mp_obj_t ndarray_binary_more(ndarray_obj_t *, ndarray_obj_t *, uint8_t , size_t *, int32_t *, int32_t *, mp_binary_op_t );
 mp_obj_t ndarray_binary_power(ndarray_obj_t *, ndarray_obj_t *, uint8_t , size_t *, int32_t *, int32_t *);
@@ -21,6 +22,7 @@ mp_obj_t ndarray_binary_logical(ndarray_obj_t *, ndarray_obj_t *, uint8_t , size
 mp_obj_t ndarray_binary_floor_divide(ndarray_obj_t *, ndarray_obj_t *, uint8_t , size_t *, int32_t *, int32_t *);
 
 mp_obj_t ndarray_inplace_ams(ndarray_obj_t *, ndarray_obj_t *, int32_t *, uint8_t );
+mp_obj_t ndarray_inplace_modulo(ndarray_obj_t *, ndarray_obj_t *, int32_t *);
 mp_obj_t ndarray_inplace_power(ndarray_obj_t *, ndarray_obj_t *, int32_t *);
 mp_obj_t ndarray_inplace_divide(ndarray_obj_t *, ndarray_obj_t *, int32_t *);
 
@@ -537,3 +539,176 @@ mp_obj_t ndarray_inplace_divide(ndarray_obj_t *, ndarray_obj_t *, int32_t *);
 } while(0)
 
 #endif /* ULAB_MAX_DIMS == 4 */
+
+#define MODULO_FLOAT1(results, array, type_left, type_right, larray, lstrides, rarray, rstrides)\
+({\
+    size_t l = 0;\
+    do {\
+        *(array)++ = MICROPY_FLOAT_C_FUN(fmod)(*((type_left *)(larray)), *((type_right *)(rarray)));\
+        (larray) += (lstrides)[ULAB_MAX_DIMS - 1];\
+        (rarray) += (rstrides)[ULAB_MAX_DIMS - 1];\
+        l++;\
+    } while(l < (results)->shape[ULAB_MAX_DIMS - 1]);\
+})
+
+#if ULAB_MAX_DIMS == 1
+#define MODULO_FLOAT_LOOP(results, type_out, type_left, type_right, larray, lstrides, rarray, rstrides) do {\
+    type_out *array = (type_out *)(results)->array;\
+    MODULO_FLOAT1((results), (array), type_left, type_right, (larray), (lstrides), (rarray), (rstrides));\
+} while(0)
+#endif /* ULAB_MAX_DIMS == 1 */
+
+#if ULAB_MAX_DIMS == 2
+#define MODULO_FLOAT_LOOP(results, type_out, type_left, type_right, larray, lstrides, rarray, rstrides) do {\
+    type_out *array = (type_out *)(results)->array;\
+    size_t l = 0;\
+    do {\
+        MODULO_FLOAT1((results), (array), type_left, type_right, (larray), (lstrides), (rarray), (rstrides));\
+        (larray) -= (lstrides)[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+        (larray) += (lstrides)[ULAB_MAX_DIMS - 2];\
+        (rarray) -= (rstrides)[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+        (rarray) += (rstrides)[ULAB_MAX_DIMS - 2];\
+        l++;\
+    } while(l < (results)->shape[ULAB_MAX_DIMS - 2]);\
+} while(0)
+#endif /* ULAB_MAX_DIMS == 2 */
+
+#if ULAB_MAX_DIMS == 3
+#define MODULO_FLOAT_LOOP(results, type_out, type_left, type_right, larray, lstrides, rarray, rstrides) do {\
+    type_out *array = (type_out *)(results)->array;\
+    size_t k = 0;\
+    do {\
+        size_t l = 0;\
+        do {\
+            MODULO_FLOAT1((results), (array), type_left, type_right, (larray), (lstrides), (rarray), (rstrides));\
+            (larray) -= (lstrides)[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+            (larray) += (lstrides)[ULAB_MAX_DIMS - 2];\
+            (rarray) -= (rstrides)[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+            (rarray) += (rstrides)[ULAB_MAX_DIMS - 2];\
+            l++;\
+        } while(l < (results)->shape[ULAB_MAX_DIMS - 2]);\
+        (larray) -= (lstrides)[ULAB_MAX_DIMS - 2] * (results)->shape[ULAB_MAX_DIMS - 2];\
+        (larray) += (lstrides)[ULAB_MAX_DIMS - 3];\
+        (rarray) -= (rstrides)[ULAB_MAX_DIMS - 2] * (results)->shape[ULAB_MAX_DIMS - 2];\
+        (rarray) += (rstrides)[ULAB_MAX_DIMS - 3];\
+        k++;\
+    } while(k < (results)->shape[ULAB_MAX_DIMS - 3]);\
+} while(0)
+#endif /* ULAB_MAX_DIMS == 3 */
+
+#if ULAB_MAX_DIMS == 4
+#define MODULO_FLOAT_LOOP(results, type_out, type_left, type_right, larray, lstrides, rarray, rstrides) do {\
+    type_out *array = (type_out *)(results)->array;\
+    size_t j = 0;\
+    do {\
+        size_t k = 0;\
+        do {\
+            size_t l = 0;\
+            do {\
+                MODULO_FLOAT1((results), (array), type_left, type_right, (larray), (lstrides), (rarray), (rstrides));\
+                (larray) -= (lstrides)[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+                (larray) += (lstrides)[ULAB_MAX_DIMS - 2];\
+                (rarray) -= (rstrides)[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+                (rarray) += (rstrides)[ULAB_MAX_DIMS - 2];\
+                l++;\
+            } while(l < (results)->shape[ULAB_MAX_DIMS - 2]);\
+            (larray) -= (lstrides)[ULAB_MAX_DIMS - 2] * (results)->shape[ULAB_MAX_DIMS - 2];\
+            (larray) += (lstrides)[ULAB_MAX_DIMS - 3];\
+            (rarray) -= (rstrides)[ULAB_MAX_DIMS - 2] * (results)->shape[ULAB_MAX_DIMS - 2];\
+            (rarray) += (rstrides)[ULAB_MAX_DIMS - 3];\
+            k++;\
+        } while(k < (results)->shape[ULAB_MAX_DIMS - 3]);\
+            (larray) -= (lstrides)[ULAB_MAX_DIMS - 3] * (results)->shape[ULAB_MAX_DIMS - 3];\
+            (larray) += (lstrides)[ULAB_MAX_DIMS - 4];\
+            (rarray) -= (rstrides)[ULAB_MAX_DIMS - 3] * (results)->shape[ULAB_MAX_DIMS - 3];\
+            (rarray) += (rstrides)[ULAB_MAX_DIMS - 4];\
+        j++;\
+    } while(j < (results)->shape[ULAB_MAX_DIMS - 4]);\
+} while(0)
+#endif /* ULAB_MAX_DIMS == 4 */
+
+
+#define INPLACE_MODULO_FLOAT1(results, type_right, larray, rarray, rstrides)\
+({\
+    size_t l = 0;\
+    do {\
+        *((mp_float_t *)larray) = MICROPY_FLOAT_C_FUN(fmod)(*((mp_float_t *)(larray)), *((type_right *)(rarray)));\
+        (larray) += (results)->strides[ULAB_MAX_DIMS - 1];\
+        (rarray) += (rstrides)[ULAB_MAX_DIMS - 1];\
+        l++;\
+    } while(l < (results)->shape[ULAB_MAX_DIMS - 1]);\
+})
+
+
+#if ULAB_MAX_DIMS == 1
+#define INPLACE_MODULO_FLOAT_LOOP(results, type_right, larray, rarray, rstrides) do {\
+    INPLACE_MODULO_FLOAT1((results), type_right, (larray), (rarray), (rstrides));\
+} while(0)
+#endif /* ULAB_MAX_DIMS == 1 */
+
+
+#if ULAB_MAX_DIMS == 2
+#define INLINE_MODULO_FLOAT_LOOP(results, type_right, larray, rarray, rstrides) do {\
+    size_t l = 0;\
+    do {\
+        INPLACE_MODULO_FLOAT1((results), type_right, (larray), (rarray), (rstrides));\
+        (larray) -= (results)->strides[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+        (larray) += (results)->strides[ULAB_MAX_DIMS - 2];\
+        (rarray) -= (rstrides)[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+        (rarray) += (rstrides)[ULAB_MAX_DIMS - 2];\
+        l++;\
+    } while(l < (results)->shape[ULAB_MAX_DIMS - 2]);\
+} while(0)
+#endif /* ULAB_MAX_DIMS == 2 */
+
+#if ULAB_MAX_DIMS == 3
+#define INLINE_MODULO_FLOAT_LOOP(results, type_right, larray, rarray, rstrides) do {\
+    size_t k = 0;\
+    do {\
+        size_t l = 0;\
+        do {\
+            INPLACE_MODULO_FLOAT1((results), type_right, (larray), (rarray), (rstrides));\
+            (larray) -= (results)->strides[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+            (larray) += (results)->strides[ULAB_MAX_DIMS - 2];\
+            (rarray) -= (rstrides)[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+            (rarray) += (rstrides)[ULAB_MAX_DIMS - 2];\
+            l++;\
+        } while(l < (results)->shape[ULAB_MAX_DIMS - 2]);\
+        (larray) -= (results)->strides[ULAB_MAX_DIMS - 2] * (results)->shape[ULAB_MAX_DIMS - 2];\
+        (larray) += (results)->strides[ULAB_MAX_DIMS - 3];\
+        (rarray) -= (rstrides)[ULAB_MAX_DIMS - 2] * (results)->shape[ULAB_MAX_DIMS - 2];\
+        (rarray) += (rstrides)[ULAB_MAX_DIMS - 3];\
+        k++;\
+    } while(k < (results)->shape[ULAB_MAX_DIMS - 3]);\
+} while(0)
+#endif /* ULAB_MAX_DIMS == 3 */
+
+#if ULAB_MAX_DIMS == 4
+#define INLINE_MODULO_FLOAT_LOOP(results, type_right, larray, rarray, rstrides) do {\
+    size_t j = 0;\
+    do {\
+        size_t k = 0;\
+        do {\
+            size_t l = 0;\
+            do {\
+                INPLACE_MODULO_FLOAT1((results), type_right, (larray), (rarray), (rstrides));\
+                (larray) -= (results)->strides[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+                (larray) += (results)->strides[ULAB_MAX_DIMS - 2];\
+                (rarray) -= (rstrides)[ULAB_MAX_DIMS - 1] * (results)->shape[ULAB_MAX_DIMS - 1];\
+                (rarray) += (rstrides)[ULAB_MAX_DIMS - 2];\
+                l++;\
+            } while(l < (results)->shape[ULAB_MAX_DIMS - 2]);\
+            (larray) -= (results)->strides[ULAB_MAX_DIMS - 2] * (results)->shape[ULAB_MAX_DIMS - 2];\
+            (larray) += (results)->strides[ULAB_MAX_DIMS - 3];\
+            (rarray) -= (rstrides)[ULAB_MAX_DIMS - 2] * (results)->shape[ULAB_MAX_DIMS - 2];\
+            (rarray) += (rstrides)[ULAB_MAX_DIMS - 3];\
+            k++;\
+        } while(k < (results)->shape[ULAB_MAX_DIMS - 3]);\
+            (larray) -= (results)->strides[ULAB_MAX_DIMS - 3] * (results)->shape[ULAB_MAX_DIMS - 3];\
+            (larray) += (results)->strides[ULAB_MAX_DIMS - 4];\
+            (rarray) -= (rstrides)[ULAB_MAX_DIMS - 3] * (results)->shape[ULAB_MAX_DIMS - 3];\
+            (rarray) += (rstrides)[ULAB_MAX_DIMS - 4];\
+        j++;\
+    } while(j < (results)->shape[ULAB_MAX_DIMS - 4]);\
+} while(0)
+#endif /* ULAB_MAX_DIMS == 4 */

code/ndarray_properties.c

@@ -6,7 +6,7 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2021 Zoltán Vörös
+ * Copyright (c) 2021-2025 Zoltán Vörös
  *
 */
 
@@ -42,6 +42,11 @@ void ndarray_properties_attr(mp_obj_t self_in, qstr attr, mp_obj_t *dest) {
                 dest[0] = ndarray_itemsize(self_in);
                 break;
             #endif
+            #if NDARRAY_HAS_NDIM
+            case MP_QSTR_ndim:
+                dest[0] = ndarray_ndim(self_in);
+                break;
+            #endif
             #if NDARRAY_HAS_SHAPE
             case MP_QSTR_shape:
                 dest[0] = ndarray_shape(self_in);

code/ndarray_properties.h

@@ -7,7 +7,7 @@
  * The MIT License (MIT)
  *
  * Copyright (c) 2020 Jeff Epler for Adafruit Industries
- *               2020-2021 Zoltán Vörös
+ *               2020-2025 Zoltán Vörös
 */
 
 #ifndef _NDARRAY_PROPERTIES_
@@ -36,6 +36,10 @@ MP_DEFINE_CONST_FUN_OBJ_1(ndarray_flatiter_make_new_obj, ndarray_flatiter_make_n
 MP_DEFINE_CONST_FUN_OBJ_1(ndarray_itemsize_obj, ndarray_itemsize);
 #endif
 
+#if NDARRAY_HAS_NDIM
+MP_DEFINE_CONST_FUN_OBJ_1(ndarray_ndim_obj, ndarray_ndim);
+#endif
+
 #if NDARRAY_HAS_SHAPE
 MP_DEFINE_CONST_FUN_OBJ_1(ndarray_shape_obj, ndarray_shape);
 #endif

code/numpy/compare.c

@@ -260,48 +260,14 @@ static mp_obj_t compare_isinf_isfinite(mp_obj_t _x, uint8_t mask) {
         }
         uint8_t *xarray = (uint8_t *)x->array;
 
-        #if ULAB_MAX_DIMS > 3
+        ITERATOR_HEAD();
-        size_t i = 0;
-        do {
-        #endif
-            #if ULAB_MAX_DIMS > 2
-            size_t j = 0;
-            do {
-            #endif
-                #if ULAB_MAX_DIMS > 1
-                size_t k = 0;
-                do {
-                #endif
-                    size_t l = 0;
-                    do {
             mp_float_t value = *(mp_float_t *)xarray;
             if(isnan(value)) {
                 *rarray++ = 0;
             } else {
                 *rarray++ = isinf(value) ? mask : 1 - mask;
             }
-                        xarray += x->strides[ULAB_MAX_DIMS - 1];
+        ITERATOR_TAIL(x, xarray);
-                        l++;
-                    } while(l < x->shape[ULAB_MAX_DIMS - 1]);
-                #if ULAB_MAX_DIMS > 1
-                    xarray -= x->strides[ULAB_MAX_DIMS - 1] * x->shape[ULAB_MAX_DIMS-1];
-                    xarray += x->strides[ULAB_MAX_DIMS - 2];
-                    k++;
-                } while(k < x->shape[ULAB_MAX_DIMS - 2]);
-                #endif
-            #if ULAB_MAX_DIMS > 2
-                xarray -= x->strides[ULAB_MAX_DIMS - 2] * x->shape[ULAB_MAX_DIMS-2];
-                xarray += x->strides[ULAB_MAX_DIMS - 3];
-                j++;
-            } while(j < x->shape[ULAB_MAX_DIMS - 3]);
-            #endif
-        #if ULAB_MAX_DIMS > 3
-            xarray -= x->strides[ULAB_MAX_DIMS - 3] * x->shape[ULAB_MAX_DIMS-3];
-            xarray += x->strides[ULAB_MAX_DIMS - 4];
-            i++;
-        } while(i < x->shape[ULAB_MAX_DIMS - 4]);
-        #endif
-
         return MP_OBJ_FROM_PTR(results);
     } else {
         mp_raise_TypeError(MP_ERROR_TEXT("wrong input type"));

code/numpy/io/io.c

@@ -338,7 +338,7 @@ static mp_obj_t io_loadtxt(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw
         buffer[read] = '\0';
         offset = buffer;
         while(*offset != '\0') {
-            if(*offset == comment_char) {
+            while(*offset == comment_char) {
                 // clear the line till the end, or the buffer's end
                 while((*offset != '\0')) {
                     offset++;
@@ -425,7 +425,7 @@ static mp_obj_t io_loadtxt(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw
         offset = buffer;
 
         while(*offset != '\0') {
-            if(*offset == comment_char) {
+            while(*offset == comment_char) {
                 // clear the line till the end, or the buffer's end
                 while((*offset != '\0')) {
                     offset++;
@@ -619,48 +619,14 @@ static mp_obj_t io_save(mp_obj_t file, mp_obj_t ndarray_) {
 
     uint8_t *array = (uint8_t *)ndarray->array;
 
-    #if ULAB_MAX_DIMS > 3
+    ITERATOR_HEAD();
-    size_t i = 0;
-    do {
-    #endif
-        #if ULAB_MAX_DIMS > 2
-        size_t j = 0;
-        do {
-        #endif
-            #if ULAB_MAX_DIMS > 1
-            size_t k = 0;
-            do {
-            #endif
-                size_t l = 0;
-                do {
         memcpy(buffer+offset, array, sz);
         offset += sz;
         if(offset == ULAB_IO_BUFFER_SIZE) {
             stream_p->write(stream, buffer, offset, &error);
             offset = 0;
         }
-                    array += ndarray->strides[ULAB_MAX_DIMS - 1];
+    ITERATOR_TAIL(ndarray, array);
-                    l++;
-                } while(l <  ndarray->shape[ULAB_MAX_DIMS - 1]);
-            #if ULAB_MAX_DIMS > 1
-                array -= ndarray->strides[ULAB_MAX_DIMS - 1] * ndarray->shape[ULAB_MAX_DIMS-1];
-                array += ndarray->strides[ULAB_MAX_DIMS - 2];
-                k++;
-            } while(k <  ndarray->shape[ULAB_MAX_DIMS - 2]);
-            #endif
-        #if ULAB_MAX_DIMS > 2
-            array -= ndarray->strides[ULAB_MAX_DIMS - 2] * ndarray->shape[ULAB_MAX_DIMS-2];
-            array += ndarray->strides[ULAB_MAX_DIMS - 3];
-            j++;
-        } while(j <  ndarray->shape[ULAB_MAX_DIMS - 3]);
-        #endif
-    #if ULAB_MAX_DIMS > 3
-        array -= ndarray->strides[ULAB_MAX_DIMS - 3] * ndarray->shape[ULAB_MAX_DIMS-3];
-        array += ndarray->strides[ULAB_MAX_DIMS - 4];
-        i++;
-    } while(i <  ndarray->shape[ULAB_MAX_DIMS - 4]);
-    #endif
-
     stream_p->write(stream, buffer, offset, &error);
     stream_p->ioctl(stream, MP_STREAM_CLOSE, 0, &error);
 
@@ -751,16 +717,32 @@ static mp_obj_t io_savetxt(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw
     char *buffer = m_new(char, ULAB_IO_BUFFER_SIZE);
     int error;
 
+    size_t len_comment;
+    char *comments;
+
+    if(mp_obj_is_str(args[5].u_obj)) {
+        const char *_comments = mp_obj_str_get_data(args[5].u_obj, &len_comment);
+        comments = (char *)_comments;
+    } else {
+        len_comment = 2;
+        comments = m_new(char, len_comment);
+        comments[0] = '#';
+        comments[1] = ' ';
+    }
+
     if(mp_obj_is_str(args[3].u_obj)) {
         size_t _len;
-        if(mp_obj_is_str(args[5].u_obj)) {
-            const char *comments = mp_obj_str_get_data(args[5].u_obj, &_len);
-            stream_p->write(stream, comments, _len, &error);
-        } else {
-            stream_p->write(stream, "# ", 2, &error);
-        }
         const char *header = mp_obj_str_get_data(args[3].u_obj, &_len);
-        stream_p->write(stream, header, _len, &error);
+
+        stream_p->write(stream, comments, len_comment, &error);
+
+        // We can't write the header in the single chunk, for it might contain line breaks
+        for(size_t i = 0; i < _len; header++, i++) {
+            stream_p->write(stream, header, 1, &error);
+            if((*header == '\n') && (i < _len)) {
+                stream_p->write(stream, comments, len_comment, &error);
+            }
+        }
         stream_p->write(stream, "\n", 1, &error);
     }
 
@@ -799,16 +781,19 @@ static mp_obj_t io_savetxt(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw
     } while(k < ndarray->shape[ULAB_MAX_DIMS - 2]);
     #endif
 
-    if(mp_obj_is_str(args[4].u_obj)) {
+    if(mp_obj_is_str(args[4].u_obj)) { // footer string
         size_t _len;
-        if(mp_obj_is_str(args[5].u_obj)) {
-            const char *comments = mp_obj_str_get_data(args[5].u_obj, &_len);
-            stream_p->write(stream, comments, _len, &error);
-        } else {
-            stream_p->write(stream, "# ", 2, &error);
-        }
         const char *footer = mp_obj_str_get_data(args[4].u_obj, &_len);
-        stream_p->write(stream, footer, _len, &error);
+
+        stream_p->write(stream, comments, len_comment, &error);
+
+        // We can't write the header in the single chunk, for it might contain line breaks
+        for(size_t i = 0; i < _len; footer++, i++) {
+            stream_p->write(stream, footer, 1, &error);
+            if((*footer == '\n') && (i < _len)) {
+                stream_p->write(stream, comments, len_comment, &error);
+            }
+        }
         stream_p->write(stream, "\n", 1, &error);
     }
 

code/numpy/numerical.c

@@ -546,10 +546,6 @@ static mp_obj_t numerical_argmin_argmax_ndarray(ndarray_obj_t *ndarray, mp_obj_t
         }
 
         m_del(int32_t, strides, ULAB_MAX_DIMS);
-
-        if(results->len == 1) {
-            return mp_binary_get_val_array(results->dtype, results->array, 0);
-        }
         return ulab_tools_restore_dims(ndarray, results, keepdims, _shape_strides);
     }
     // we should never get to this point

code/numpy/poly.c

@@ -197,42 +197,9 @@ mp_obj_t poly_polyval(mp_obj_t o_p, mp_obj_t o_x) {
 
         // TODO: these loops are really nothing, but the re-impplementation of
         // ITERATE_VECTOR from vectorise.c. We could pass a function pointer here
-        #if ULAB_MAX_DIMS > 3
+        ITERATOR_HEAD();
-        size_t i = 0;
-        do {
-        #endif
-            #if ULAB_MAX_DIMS > 2
-            size_t j = 0;
-            do {
-            #endif
-                #if ULAB_MAX_DIMS > 1
-                size_t k = 0;
-                do {
-                #endif
-                    size_t l = 0;
-                    do {
             *array++ = poly_eval(func(sarray), p, plen);
-                        sarray += source->strides[ULAB_MAX_DIMS - 1];
+        ITERATOR_TAIL(source, sarray);
-                        l++;
-                    } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-                #if ULAB_MAX_DIMS > 1
-                    sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS-1];
-                    sarray += source->strides[ULAB_MAX_DIMS - 2];
-                    k++;
-                } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-                #endif
-            #if ULAB_MAX_DIMS > 2
-                sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS-2];
-                sarray += source->strides[ULAB_MAX_DIMS - 3];
-                j++;
-            } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-            #endif
-        #if ULAB_MAX_DIMS > 3
-            sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS-3];
-            sarray += source->strides[ULAB_MAX_DIMS - 4];
-            i++;
-        } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-        #endif
     } else {
         // o_x had better be a one-dimensional standard iterable
         ndarray = ndarray_new_linear_array(mp_obj_get_int(mp_obj_len_maybe(o_x)), NDARRAY_FLOAT);

code/numpy/random/random.c

@@ -57,7 +57,7 @@ const mp_obj_type_t random_generator_type = {
 void random_generator_print(const mp_print_t *print, mp_obj_t self_in, mp_print_kind_t kind) {
     (void)kind;
     random_generator_obj_t *self = MP_OBJ_TO_PTR(self_in);
-    mp_printf(MP_PYTHON_PRINTER, "Gnerator() at 0x%p", self);
+    mp_printf(MP_PYTHON_PRINTER, "Generator() at 0x%p", self);
 }
 
 mp_obj_t random_generator_make_new(const mp_obj_type_t *type, size_t n_args, size_t n_kw, const mp_obj_t *args) {
@@ -149,12 +149,9 @@ static mp_obj_t random_normal(size_t n_args, const mp_obj_t *pos_args, mp_map_t
             ndarray = ndarray_new_linear_array((size_t)mp_obj_get_int(size), NDARRAY_FLOAT);
         } else if(mp_obj_is_type(size, &mp_type_tuple)) {
             mp_obj_tuple_t *_shape = MP_OBJ_TO_PTR(size);
-            if(_shape->len > ULAB_MAX_DIMS) {
-                mp_raise_ValueError(MP_ERROR_TEXT("maximum number of dimensions is " MP_STRINGIFY(ULAB_MAX_DIMS)));
-            }
             ndarray = ndarray_new_ndarray_from_tuple(_shape, NDARRAY_FLOAT);
         } else { // input type not supported
-            mp_raise_TypeError(MP_ERROR_TEXT("shape must be None, and integer or a tuple of integers"));
+            mp_raise_TypeError(MP_ERROR_TEXT("shape must be None, an integer or a tuple of integers"));
         }
     } else {
         // return single value
@@ -221,27 +218,16 @@ static mp_obj_t random_random(size_t n_args, const mp_obj_t *pos_args, mp_map_t
     mp_obj_t out = args[2].u_obj;
 
     ndarray_obj_t *ndarray = NULL;
-    size_t *shape = m_new(size_t, ULAB_MAX_DIMS);
+    size_t *shape = m_new0(size_t, ULAB_MAX_DIMS);
-    uint8_t ndim = 1;
 
     if(size != mp_const_none) {
         if(mp_obj_is_int(size)) {
             shape[ULAB_MAX_DIMS - 1] = (size_t)mp_obj_get_int(size);
         } else if(mp_obj_is_type(size, &mp_type_tuple)) {
             mp_obj_tuple_t *_shape = MP_OBJ_TO_PTR(size);
-            if(_shape->len > ULAB_MAX_DIMS) {
+            ndarray = ndarray_new_ndarray_from_tuple(_shape, NDARRAY_FLOAT);
-                mp_raise_ValueError(MP_ERROR_TEXT("maximum number of dimensions is " MP_STRINGIFY(ULAB_MAX_DIMS)));
-            }
-            ndim = _shape->len;
-            for(size_t i = 0; i < ULAB_MAX_DIMS; i++) {
-                if(i >= ndim) {
-                    shape[ULAB_MAX_DIMS - 1 - i] = 0;
-                } else {
-                    shape[ULAB_MAX_DIMS - 1 - i] = mp_obj_get_int(_shape->items[i]);
-                }
-            }
         } else { // input type not supported
-            mp_raise_TypeError(MP_ERROR_TEXT("shape must be None, and integer or a tuple of integers"));
+            mp_raise_TypeError(MP_ERROR_TEXT("shape must be None, an integer or a tuple of integers"));
         }
     }
 
@@ -267,7 +253,8 @@ static mp_obj_t random_random(size_t n_args, const mp_obj_t *pos_args, mp_map_t
         }
     } else { // out == None
         if(size != mp_const_none) {
-            ndarray = ndarray_new_dense_ndarray(ndim, shape, NDARRAY_FLOAT);
+            mp_obj_tuple_t *_shape = MP_OBJ_TO_PTR(size);
+            ndarray = ndarray_new_ndarray_from_tuple(_shape, NDARRAY_FLOAT);
         } else {
             // return single value
             mp_float_t value;
@@ -336,13 +323,9 @@ static mp_obj_t random_uniform(size_t n_args, const mp_obj_t *pos_args, mp_map_t
         return mp_obj_new_float(value);
     } else if(mp_obj_is_type(size, &mp_type_tuple)) {
         mp_obj_tuple_t *_shape = MP_OBJ_TO_PTR(size);
-        // TODO: this could be reduced, if the inspection was in the ndarray_new_ndarray_from_tuple function
-        if(_shape->len > ULAB_MAX_DIMS) {
-            mp_raise_ValueError(MP_ERROR_TEXT("maximum number of dimensions is " MP_STRINGIFY(ULAB_MAX_DIMS)));
-        }
         ndarray = ndarray_new_ndarray_from_tuple(_shape, NDARRAY_FLOAT);
     } else { // input type not supported
-        mp_raise_TypeError(MP_ERROR_TEXT("shape must be None, and integer or a tuple of integers"));
+        mp_raise_TypeError(MP_ERROR_TEXT("shape must be None, an integer or a tuple of integers"));
     }
 
     mp_float_t *array = (mp_float_t *)ndarray->array;

code/numpy/vector.c

@@ -91,43 +91,10 @@ static mp_obj_t vector_generic_vector(size_t n_args, const mp_obj_t *pos_args, m
 
             mp_float_t (*func)(void *) = ndarray_get_float_function(source->dtype);
 
-            #if ULAB_MAX_DIMS > 3
+            ITERATOR_HEAD();
-            size_t i = 0;
-            do {
-            #endif
-                #if ULAB_MAX_DIMS > 2
-                size_t j = 0;
-                do {
-                #endif
-                    #if ULAB_MAX_DIMS > 1
-                    size_t k = 0;
-                    do {
-                    #endif
-                        size_t l = 0;
-                        do {
                 mp_float_t value = func(sarray);
                 *tarray++ = f(value);
-                            sarray += source->strides[ULAB_MAX_DIMS - 1];
+            ITERATOR_TAIL(source, sarray);
-                            l++;
-                        } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-                    #if ULAB_MAX_DIMS > 1
-                        sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS-1];
-                        sarray += source->strides[ULAB_MAX_DIMS - 2];
-                        k++;
-                    } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-                    #endif /* ULAB_MAX_DIMS > 1 */
-                #if ULAB_MAX_DIMS > 2
-                    sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS-2];
-                    sarray += source->strides[ULAB_MAX_DIMS - 3];
-                    j++;
-                } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-                #endif /* ULAB_MAX_DIMS > 2 */
-            #if ULAB_MAX_DIMS > 3
-                sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS-3];
-                sarray += source->strides[ULAB_MAX_DIMS - 4];
-                i++;
-            } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-            #endif /* ULAB_MAX_DIMS > 3 */
         #else
         if(source->dtype == NDARRAY_UINT8) {
             ITERATE_VECTOR(uint8_t, target, tarray, tstrides, source, sarray);
@@ -171,43 +138,10 @@ static mp_obj_t vector_generic_vector(mp_obj_t o_in, mp_float_t (*f)(mp_float_t)
 
             mp_float_t (*func)(void *) = ndarray_get_float_function(source->dtype);
 
-            #if ULAB_MAX_DIMS > 3
+            ITERATOR_HEAD();
-            size_t i = 0;
-            do {
-            #endif
-                #if ULAB_MAX_DIMS > 2
-                size_t j = 0;
-                do {
-                #endif
-                    #if ULAB_MAX_DIMS > 1
-                    size_t k = 0;
-                    do {
-                    #endif
-                        size_t l = 0;
-                        do {
                 mp_float_t value = func(sarray);
                 *array++ = f(value);
-                            sarray += source->strides[ULAB_MAX_DIMS - 1];
+            ITERATOR_TAIL(source, sarray);
-                            l++;
-                        } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-                    #if ULAB_MAX_DIMS > 1
-                        sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS-1];
-                        sarray += source->strides[ULAB_MAX_DIMS - 2];
-                        k++;
-                    } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-                    #endif /* ULAB_MAX_DIMS > 1 */
-                #if ULAB_MAX_DIMS > 2
-                    sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS-2];
-                    sarray += source->strides[ULAB_MAX_DIMS - 3];
-                    j++;
-                } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-                #endif /* ULAB_MAX_DIMS > 2 */
-            #if ULAB_MAX_DIMS > 3
-                sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS-3];
-                sarray += source->strides[ULAB_MAX_DIMS - 4];
-                i++;
-            } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-            #endif /* ULAB_MAX_DIMS > 3 */
         #else
         if(source->dtype == NDARRAY_UINT8) {
             ITERATE_VECTOR(uint8_t, array, source, sarray);
@@ -327,43 +261,11 @@ mp_obj_t vector_around(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_arg
 
     mp_float_t (*func)(void *) = ndarray_get_float_function(source->dtype);
 
-    #if ULAB_MAX_DIMS > 3
+    ITERATOR_HEAD();
-    size_t i = 0;
-    do {
-    #endif
-        #if ULAB_MAX_DIMS > 2
-        size_t j = 0;
-        do {
-        #endif
-            #if ULAB_MAX_DIMS > 1
-            size_t k = 0;
-            do {
-            #endif
-                size_t l = 0;
-                do {
         mp_float_t f = func(sarray);
         *narray++ = MICROPY_FLOAT_C_FUN(round)(f * mul) / mul;
-                    sarray += source->strides[ULAB_MAX_DIMS - 1];
+    ITERATOR_TAIL(source, sarray);
-                    l++;
+
-                } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-            #if ULAB_MAX_DIMS > 1
-                sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS-1];
-                sarray += source->strides[ULAB_MAX_DIMS - 2];
-                k++;
-            } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-            #endif
-        #if ULAB_MAX_DIMS > 2
-            sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS-2];
-            sarray += source->strides[ULAB_MAX_DIMS - 3];
-            j++;
-        } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-        #endif
-    #if ULAB_MAX_DIMS > 3
-        sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS-3];
-        sarray += source->strides[ULAB_MAX_DIMS - 4];
-        i++;
-    } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-    #endif
     return MP_OBJ_FROM_PTR(ndarray);
 }
 
@@ -631,46 +533,13 @@ static mp_obj_t vector_exp(mp_obj_t o_in) {
             mp_float_t *array = (mp_float_t *)ndarray->array;
             uint8_t itemsize = sizeof(mp_float_t);
 
-            #if ULAB_MAX_DIMS > 3
+            ITERATOR_HEAD();
-            size_t i = 0;
-            do {
-            #endif
-                #if ULAB_MAX_DIMS > 2
-                size_t j = 0;
-                do {
-                #endif
-                    #if ULAB_MAX_DIMS > 1
-                    size_t k = 0;
-                    do {
-                    #endif
-                        size_t l = 0;
-                        do {
                 mp_float_t real = *(mp_float_t *)sarray;
                 mp_float_t imag = *(mp_float_t *)(sarray + itemsize);
                 mp_float_t exp_real = MICROPY_FLOAT_C_FUN(exp)(real);
                 *array++ = exp_real * MICROPY_FLOAT_C_FUN(cos)(imag);
                 *array++ = exp_real * MICROPY_FLOAT_C_FUN(sin)(imag);
-                            sarray += source->strides[ULAB_MAX_DIMS - 1];
+            ITERATOR_TAIL(source, sarray);
-                            l++;
-                        } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-                    #if ULAB_MAX_DIMS > 1
-                        sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS-1];
-                        sarray += source->strides[ULAB_MAX_DIMS - 2];
-                        k++;
-                    } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-                    #endif /* ULAB_MAX_DIMS > 1 */
-                #if ULAB_MAX_DIMS > 2
-                    sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS-2];
-                    sarray += source->strides[ULAB_MAX_DIMS - 3];
-                    j++;
-                } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-                #endif /* ULAB_MAX_DIMS > 2 */
-            #if ULAB_MAX_DIMS > 3
-                sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS-3];
-                sarray += source->strides[ULAB_MAX_DIMS - 4];
-                i++;
-            } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-            #endif /* ULAB_MAX_DIMS > 3 */
             return MP_OBJ_FROM_PTR(ndarray);
         }
     }
@@ -921,20 +790,7 @@ mp_obj_t vector_sqrt(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args)
                 mp_float_t *array = (mp_float_t *)ndarray->array;
                 uint8_t itemsize = sizeof(mp_float_t);
 
-                #if ULAB_MAX_DIMS > 3
+                ITERATOR_HEAD();
-                size_t i = 0;
-                do {
-                #endif
-                    #if ULAB_MAX_DIMS > 2
-                    size_t j = 0;
-                    do {
-                    #endif
-                        #if ULAB_MAX_DIMS > 1
-                        size_t k = 0;
-                        do {
-                        #endif
-                            size_t l = 0;
-                            do {
                     mp_float_t real = *(mp_float_t *)sarray;
                     mp_float_t imag = *(mp_float_t *)(sarray + itemsize);
                     mp_float_t sqrt_abs = MICROPY_FLOAT_C_FUN(sqrt)(real * real + imag * imag);
@@ -942,47 +798,15 @@ mp_obj_t vector_sqrt(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args)
                     mp_float_t theta = MICROPY_FLOAT_CONST(0.5) * MICROPY_FLOAT_C_FUN(atan2)(imag, real);
                     *array++ = sqrt_abs * MICROPY_FLOAT_C_FUN(cos)(theta);
                     *array++ = sqrt_abs * MICROPY_FLOAT_C_FUN(sin)(theta);
-                                sarray += source->strides[ULAB_MAX_DIMS - 1];
+                ITERATOR_TAIL(source, sarray);
-                                l++;
+
-                            } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-                        #if ULAB_MAX_DIMS > 1
-                            sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS-1];
-                            sarray += source->strides[ULAB_MAX_DIMS - 2];
-                            k++;
-                        } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-                        #endif /* ULAB_MAX_DIMS > 1 */
-                    #if ULAB_MAX_DIMS > 2
-                        sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS-2];
-                        sarray += source->strides[ULAB_MAX_DIMS - 3];
-                        j++;
-                    } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-                    #endif /* ULAB_MAX_DIMS > 2 */
-                #if ULAB_MAX_DIMS > 3
-                    sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS-3];
-                    sarray += source->strides[ULAB_MAX_DIMS - 4];
-                    i++;
-                } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-                #endif /* ULAB_MAX_DIMS > 3 */
                 return MP_OBJ_FROM_PTR(ndarray);
             } else if(source->dtype == NDARRAY_FLOAT) {
                 uint8_t *sarray = (uint8_t *)source->array;
                 ndarray_obj_t *ndarray = ndarray_new_dense_ndarray(source->ndim, source->shape, NDARRAY_COMPLEX);
                 mp_float_t *array = (mp_float_t *)ndarray->array;
 
-                #if ULAB_MAX_DIMS > 3
+                ITERATOR_HEAD();
-                size_t i = 0;
-                do {
-                #endif
-                    #if ULAB_MAX_DIMS > 2
-                    size_t j = 0;
-                    do {
-                    #endif
-                        #if ULAB_MAX_DIMS > 1
-                        size_t k = 0;
-                        do {
-                        #endif
-                            size_t l = 0;
-                            do {
                     mp_float_t value = *(mp_float_t *)sarray;
                     if(value >= MICROPY_FLOAT_CONST(0.0)) {
                         *array++ = MICROPY_FLOAT_C_FUN(sqrt)(value);
@@ -991,27 +815,8 @@ mp_obj_t vector_sqrt(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args)
                         array++;
                         *array++ = MICROPY_FLOAT_C_FUN(sqrt)(-value);
                     }
-                                sarray += source->strides[ULAB_MAX_DIMS - 1];
+                ITERATOR_TAIL(source, sarray);
-                                l++;
+
-                            } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-                        #if ULAB_MAX_DIMS > 1
-                            sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS-1];
-                            sarray += source->strides[ULAB_MAX_DIMS - 2];
-                            k++;
-                        } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-                        #endif /* ULAB_MAX_DIMS > 1 */
-                    #if ULAB_MAX_DIMS > 2
-                        sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS-2];
-                        sarray += source->strides[ULAB_MAX_DIMS - 3];
-                        j++;
-                    } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-                    #endif /* ULAB_MAX_DIMS > 2 */
-                #if ULAB_MAX_DIMS > 3
-                    sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS-3];
-                    sarray += source->strides[ULAB_MAX_DIMS - 4];
-                    i++;
-                } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-                #endif /* ULAB_MAX_DIMS > 3 */
                 return MP_OBJ_FROM_PTR(ndarray);
             } else {
                 mp_raise_TypeError(MP_ERROR_TEXT("input dtype must be float or complex"));
@@ -1071,45 +876,12 @@ static mp_obj_t vector_vectorized_function_call(mp_obj_t self_in, size_t n_args,
         uint8_t *sarray = (uint8_t *)source->array;
         uint8_t *narray = (uint8_t *)ndarray->array;
         
-        #if ULAB_MAX_DIMS > 3
+        ITERATOR_HEAD();
-        size_t i = 0;
-        do {
-        #endif
-            #if ULAB_MAX_DIMS > 2
-            size_t j = 0;
-            do {
-            #endif
-                #if ULAB_MAX_DIMS > 1
-                size_t k = 0;
-                do {
-                #endif
-                    size_t l = 0;
-                    do {
             avalue[0] = mp_binary_get_val_array(source->dtype, sarray, 0);
             fvalue = MP_OBJ_TYPE_GET_SLOT(self->type, call)(self->fun, 1, 0, avalue);
             ndarray_set_value(self->otypes, narray, 0, fvalue);
-                        sarray += source->strides[ULAB_MAX_DIMS - 1];
             narray += ndarray->itemsize;
-                        l++;
+        ITERATOR_TAIL(source, sarray);
-                    } while(l < source->shape[ULAB_MAX_DIMS - 1]);
-                #if ULAB_MAX_DIMS > 1
-                    sarray -= source->strides[ULAB_MAX_DIMS - 1] * source->shape[ULAB_MAX_DIMS - 1];
-                    sarray += source->strides[ULAB_MAX_DIMS - 2];
-                    k++;
-                } while(k < source->shape[ULAB_MAX_DIMS - 2]);
-                #endif /* ULAB_MAX_DIMS > 1 */
-            #if ULAB_MAX_DIMS > 2
-                sarray -= source->strides[ULAB_MAX_DIMS - 2] * source->shape[ULAB_MAX_DIMS - 2];
-                sarray += source->strides[ULAB_MAX_DIMS - 3];
-                j++;
-            } while(j < source->shape[ULAB_MAX_DIMS - 3]);
-            #endif /* ULAB_MAX_DIMS > 2 */
-        #if ULAB_MAX_DIMS > 3
-            sarray -= source->strides[ULAB_MAX_DIMS - 3] * source->shape[ULAB_MAX_DIMS - 3];
-            sarray += source->strides[ULAB_MAX_DIMS - 4];
-            i++;
-        } while(i < source->shape[ULAB_MAX_DIMS - 4]);
-        #endif /* ULAB_MAX_DIMS > 3 */
         
         return MP_OBJ_FROM_PTR(ndarray);
     } else if(mp_obj_is_type(args[0], &mp_type_tuple) || mp_obj_is_type(args[0], &mp_type_list) ||

code/scipy/integrate/integrate.c

@@ -41,21 +41,21 @@
 
 #if MICROPY_FLOAT_IMPL == MICROPY_FLOAT_IMPL_DOUBLE
 ULAB_DEFINE_FLOAT_CONST(etolerance, MICROPY_FLOAT_CONST(1e-14), 0x283424dcUL, 0x3e901b2b29a4692bULL);
-#define MACHEPS        MICROPY_FLOAT_CONST(1e-17)
+#define ULAB_MACHEPS        MICROPY_FLOAT_CONST(1e-17)
 #else
 ULAB_DEFINE_FLOAT_CONST(etolerance, MICROPY_FLOAT_CONST(1e-8), 0x358637cfUL, 0x3e7010c6f7d42d18ULL);
-#define MACHEPS        MICROPY_FLOAT_CONST(1e-8)
+#define ULAB_MACHEPS        MICROPY_FLOAT_CONST(1e-8)
 #endif
 
-#define ZERO           MICROPY_FLOAT_CONST(0.0)
+#define ULAB_ZERO           MICROPY_FLOAT_CONST(0.0)
-#define POINT_TWO_FIVE MICROPY_FLOAT_CONST(0.25)
+#define ULAB_POINT_TWO_FIVE MICROPY_FLOAT_CONST(0.25)
-#define ONE            MICROPY_FLOAT_CONST(1.0)
+#define ULAB_ONE            MICROPY_FLOAT_CONST(1.0)
-#define TWO            MICROPY_FLOAT_CONST(2.0)
+#define ULAB_TWO            MICROPY_FLOAT_CONST(2.0)
-#define FOUR           MICROPY_FLOAT_CONST(4.0)
+#define ULAB_FOUR           MICROPY_FLOAT_CONST(4.0)
-#define SIX            MICROPY_FLOAT_CONST(6.0)
+#define ULAB_SIX            MICROPY_FLOAT_CONST(6.0)
-#define TEN            MICROPY_FLOAT_CONST(10.0)
+#define ULAB_TEN            MICROPY_FLOAT_CONST(10.0)
-#define FIFTEEN        MICROPY_FLOAT_CONST(15.0)
+#define ULAB_FIFTEEN        MICROPY_FLOAT_CONST(15.0)
-#define EPS_5          MICROPY_FLOAT_CONST(1e-5)
+#define ULAB_EPSILON_5      MICROPY_FLOAT_CONST(1e-5)
 
 
 static mp_float_t integrate_python_call(const mp_obj_type_t *type, mp_obj_t fun, mp_float_t x, mp_obj_t *fargs, uint8_t nparams) {
@@ -68,7 +68,7 @@ static mp_float_t integrate_python_call(const mp_obj_type_t *type, mp_obj_t fun,
 
 // sign helper function
 int sign(mp_float_t x) {
-    if (x >= ZERO) 
+    if (x >= ULAB_ZERO) 
         return 1;
     else
         return -1;
@@ -85,7 +85,7 @@ mp_float_t exp_sinh_opt_d(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_
     mp_obj_t fargs[1];
     mp_float_t h2 = integrate_python_call(type, fun, a + d/2, fargs, 0) - integrate_python_call(type, fun, (a + d*2)*4, fargs, 0);
     int i = 1, j = 32;                   // j=32 is optimal to find r
-    if (isfinite(h2) && MICROPY_FLOAT_C_FUN(fabs)(h2) > EPS_5) {    // if |h2| > 2^-16
+    if (isfinite(h2) && MICROPY_FLOAT_C_FUN(fabs)(h2) > ULAB_EPSILON_5) {    // if |h2| > 2^-16
         mp_float_t r, fl, fr, h, s = 0, lfl, lfr, lr = 2;
         do {                                  // find max j such that fl and fr are finite
             j /= 2;
@@ -118,8 +118,8 @@ mp_float_t exp_sinh_opt_d(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_
             if (s > eps) {                    // if sum of |h| > eps
                 h = lfl - lfr;                // use last fl and fr before the sign change
                 r = lr;                       // use last r before the sign change
-                if (h != ZERO)                   // if last diff != 0, back up r by one step
+                if (h != ULAB_ZERO)           // if last diff != 0, back up r by one step
-                    r /= TWO;
+                    r /= ULAB_TWO;
                 if (MICROPY_FLOAT_C_FUN(fabs)(lfl) < MICROPY_FLOAT_C_FUN(fabs)(lfr))
                     d /= r;                   // move d closer to the finite endpoint
                 else
@@ -135,8 +135,8 @@ mp_float_t exp_sinh_opt_d(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_
 mp_float_t tanhsinh(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, uint16_t n, mp_float_t eps, mp_float_t *e) {
     const mp_obj_type_t *type = mp_obj_get_type(fun);
     mp_obj_t fargs[1];
-    const mp_float_t tol = TEN*eps;
+    const mp_float_t tol = ULAB_TEN * eps;
-    mp_float_t c = ZERO, d = ONE, s, sign = ONE, v, h = TWO;
+    mp_float_t c = ULAB_ZERO, d = ULAB_ONE, s, sign = ULAB_ONE, v, h = ULAB_TWO;
     int k = 0, mode = 0;                   // Tanh-Sinh = 0, Exp-Sinh = 1, Sinh-Sinh = 2
     if (b < a) {                                // swap bounds
         v = b;
@@ -145,8 +145,8 @@ mp_float_t tanhsinh(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, u
         sign = -1;
     }
     if (isfinite(a) && isfinite(b)) {
-        c = (a+b)/TWO;
+        c = (a+b) / ULAB_TWO;
-        d = (b-a)/TWO;
+        d = (b-a) / ULAB_TWO;
         v = c;
     }
     else if (isfinite(a)) {
@@ -165,20 +165,20 @@ mp_float_t tanhsinh(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, u
     }
     else {
         mode = 2;                               // Sinh-Sinh
-        v = ZERO;
+        v = ULAB_ZERO;
     }
     s = integrate_python_call(type, fun, v, fargs, 0);
     do {
-        mp_float_t p = ZERO, q, fp = ZERO, fm = ZERO, t, eh;
+        mp_float_t p = ULAB_ZERO, q, fp = ULAB_ZERO, fm = ULAB_ZERO, t, eh;
-        h /= TWO;
+        h /= ULAB_TWO;
         t = eh = MICROPY_FLOAT_C_FUN(exp)(h);
-        if (k > ZERO)
+        if (k > ULAB_ZERO)
            eh *= eh;
         if (mode == 0) {                        // Tanh-Sinh
             do {
-                mp_float_t u = MICROPY_FLOAT_C_FUN(exp)(ONE / t - t); // = exp(-2*sinh(j*h)) = 1/exp(sinh(j*h))^2
+                mp_float_t u = MICROPY_FLOAT_C_FUN(exp)(ULAB_ONE / t - t); // = exp(-2*sinh(j*h)) = 1/exp(sinh(j*h))^2
-                mp_float_t r = TWO * u / (ONE + u);                   // = 1 - tanh(sinh(j*h))
+                mp_float_t r = ULAB_TWO * u / (ULAB_ONE + u);                   // = 1 - tanh(sinh(j*h))
-                mp_float_t w = (t + ONE / t) * r / (ONE + u);         // = cosh(j*h)/cosh(sinh(j*h))^2
+                mp_float_t w = (t + ULAB_ONE / t) * r / (ULAB_ONE + u);         // = cosh(j*h)/cosh(sinh(j*h))^2
                 mp_float_t x = d*r;
                 if (a+x > a) {                  // if too close to a then reuse previous fp
                     mp_float_t y = integrate_python_call(type, fun, a+x, fargs, 0);
@@ -196,11 +196,11 @@ mp_float_t tanhsinh(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, u
             } while (MICROPY_FLOAT_C_FUN(fabs)(q) > eps*MICROPY_FLOAT_C_FUN(fabs)(p));
         }
         else {
-            t /= TWO;
+            t /= ULAB_TWO;
             do {
-                mp_float_t r = MICROPY_FLOAT_C_FUN(exp)(t - POINT_TWO_FIVE / t);    // = exp(sinh(j*h))
+                mp_float_t r = MICROPY_FLOAT_C_FUN(exp)(t - ULAB_POINT_TWO_FIVE / t);    // = exp(sinh(j*h))
                 mp_float_t x, y, w = r;
-                q = ZERO;
+                q = ULAB_ZERO;
                 if (mode == 1) {                // Exp-Sinh
                     x = c + d/r;
                     if (x == c)                 // if x hit the finite endpoint then break
@@ -210,8 +210,8 @@ mp_float_t tanhsinh(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, u
                         q += y/w;
                 }
                 else {                          // Sinh-Sinh
-                    r = (r - ONE / r) / TWO;              // = sinh(sinh(j*h))
+                    r = (r - ULAB_ONE / r) / ULAB_TWO;              // = sinh(sinh(j*h))
-                    w = (w + ONE / w) / TWO;              // = cosh(sinh(j*h))
+                    w = (w + ULAB_ONE / w) / ULAB_TWO;              // = cosh(sinh(j*h))
                     x = c - d*r;
                     y = integrate_python_call(type, fun, x, fargs, 0);
                     if (isfinite(y))    // if f(x) is finite, add to local sum
@@ -221,7 +221,7 @@ mp_float_t tanhsinh(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, u
                 y = integrate_python_call(type, fun, x, fargs, 0);
                 if (isfinite(y))        // if f(x) is finite, add to local sum
                     q += y*w;
-                q *= t + POINT_TWO_FIVE / t;                   // q *= cosh(j*h)
+                q *= t + ULAB_POINT_TWO_FIVE / t;                   // q *= cosh(j*h)
                 p += q;
                 t *= eh;
             } while (MICROPY_FLOAT_C_FUN(fabs)(q) > eps*MICROPY_FLOAT_C_FUN(fabs)(p));
@@ -319,12 +319,12 @@ mp_float_t qromb(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, uint
     for (i = 1; i < n; ++i) {
         unsigned long long k = 1UL << i;
         unsigned long long s = 1;
-        mp_float_t sum = ZERO;
+        mp_float_t sum = ULAB_ZERO;
         mp_float_t *Rt;
-        h /= TWO;
+        h /= ULAB_TWO;
         for (j = 1; j < k; j += 2)
             sum += integrate_python_call(type, fun, a+j*h, fargs, 0);
-        Ru[0] = h*sum + Ro[0] / TWO;
+        Ru[0] = h*sum + Ro[0] / ULAB_TWO;
         for (j = 1; j <= i; ++j) {
             s <<= 2;
             Ru[j] = (s*Ru[j-1] - Ro[j-1])/(s-1);
@@ -408,17 +408,17 @@ mp_float_t as(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, mp_floa
               mp_float_t fb, mp_float_t v, mp_float_t eps, int n, mp_float_t t) {
     const mp_obj_type_t *type = mp_obj_get_type(fun);
     mp_obj_t fargs[1];
-    mp_float_t h = (b-a) / TWO;
+    mp_float_t h = (b-a) / ULAB_TWO;
-    mp_float_t f1 = integrate_python_call(type, fun, a + h / TWO, fargs, 0);
+    mp_float_t f1 = integrate_python_call(type, fun, a + h / ULAB_TWO, fargs, 0);
-    mp_float_t f2 = integrate_python_call(type, fun, b - h / TWO, fargs, 0);
+    mp_float_t f2 = integrate_python_call(type, fun, b - h / ULAB_TWO, fargs, 0);
-    mp_float_t sl = h*(fa + FOUR * f1 + fm) / SIX;
+    mp_float_t sl = h*(fa + ULAB_FOUR * f1 + fm) / ULAB_SIX;
-    mp_float_t sr = h*(fm + FOUR * f2 + fb) / SIX;
+    mp_float_t sr = h*(fm + ULAB_FOUR * f2 + fb) / ULAB_SIX;
     mp_float_t s = sl+sr;
-    mp_float_t d = (s-v) / FIFTEEN;
+    mp_float_t d = (s-v) / ULAB_FIFTEEN;
     mp_float_t m = a+h;
     if (n <= 0 || MICROPY_FLOAT_C_FUN(fabs)(d) < eps)
         return t + s + d; // note: fabs(d) can be used as error estimate
-    eps /= TWO;
+    eps /= ULAB_TWO;
     --n;
     t = as(fun, a, m, fa, f1, fm, sl, eps, n, t);
     return as(fun, m, b, fm, f2, fb, sr, eps, n, t);
@@ -430,7 +430,7 @@ mp_float_t qasi(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, int n
     mp_float_t fa = integrate_python_call(type, fun, a, fargs, 0);
     mp_float_t fm = integrate_python_call(type, fun, (a+b)/2, fargs, 0);
     mp_float_t fb = integrate_python_call(type, fun, b, fargs, 0);
-    mp_float_t v = (fa + FOUR * fm + fb) * (b-a) / SIX;
+    mp_float_t v = (fa + ULAB_FOUR * fm + fb) * (b-a) / ULAB_SIX;
     return as(fun, a, b, fa, fm, fb, v, eps, n, 0);
 }
 
@@ -562,8 +562,8 @@ mp_float_t gk(mp_float_t (*fun)(mp_float_t), mp_float_t c, mp_float_t d, mp_floa
     };
     const mp_obj_type_t *type = mp_obj_get_type(fun);
     mp_obj_t fargs[1];
-    mp_float_t p = ZERO; // kronrod quadrature sum
+    mp_float_t p = ULAB_ZERO; // kronrod quadrature sum
-    mp_float_t q = ZERO; // gauss quadrature sum
+    mp_float_t q = ULAB_ZERO; // gauss quadrature sum
     mp_float_t fp, fm;
     mp_float_t e;
     int i;
@@ -581,24 +581,24 @@ mp_float_t gk(mp_float_t (*fun)(mp_float_t), mp_float_t c, mp_float_t d, mp_floa
         p += (fp + fm) * weights[i];
     }
     *err = MICROPY_FLOAT_C_FUN(fabs)(p - q);
-    e = MICROPY_FLOAT_C_FUN(fabs)(2 * p * MACHEPS); // optional, to take 1e-17 MachEps prec. into account
+    e = MICROPY_FLOAT_C_FUN(fabs)(2 * p * ULAB_MACHEPS); // optional, to take 1e-17 MachEps prec. into account
     if (*err < e)
         *err = e;
     return p;
 }
 
 mp_float_t qakro(mp_float_t (*fun)(mp_float_t), mp_float_t a, mp_float_t b, int n, mp_float_t tol, mp_float_t eps, mp_float_t *err) {
-    mp_float_t c = (a+b) / TWO;
+    mp_float_t c = (a+b) / ULAB_TWO;
-    mp_float_t d = (b-a) / TWO;
+    mp_float_t d = (b-a) / ULAB_TWO;
     mp_float_t e;
     mp_float_t r = gk(fun, c, d, &e);
     mp_float_t s = d*r;
     mp_float_t t = MICROPY_FLOAT_C_FUN(fabs)(s*tol);
-    if (tol == ZERO)
+    if (tol == ULAB_ZERO)
         tol = t;
     if (n > 0 && t < e && tol < e) {
-        s = qakro(fun, a, c, n-1, t / TWO, eps, err);
+        s = qakro(fun, a, c, n-1, t / ULAB_TWO, eps, err);
-        s += qakro(fun, c, b, n-1, t / TWO, eps, &e);
+        s += qakro(fun, c, b, n-1, t / ULAB_TWO, eps, &e);
         *err += e;
         return s;
     }

code/ulab.c

@@ -33,7 +33,7 @@
 #include "user/user.h"
 #include "utils/utils.h"
 
-#define ULAB_VERSION 6.7.2
+#define ULAB_VERSION 6.9.0
 #define xstr(s) str(s)
 #define str(s) #s
 

code/ulab.h

@@ -117,6 +117,10 @@
 #define NDARRAY_HAS_BINARY_OP_LESS_EQUAL    (1)
 #endif
 
+#ifndef NDARRAY_HAS_BINARY_OP_MODULO
+#define NDARRAY_HAS_BINARY_OP_MODULO        (1)
+#endif
+
 #ifndef NDARRAY_HAS_BINARY_OP_MORE
 #define NDARRAY_HAS_BINARY_OP_MORE          (1)
 #endif
@@ -161,6 +165,10 @@
 #define NDARRAY_HAS_INPLACE_ADD             (1)
 #endif
 
+#ifndef NDARRAY_HAS_INPLACE_MODULO
+#define NDARRAY_HAS_INPLACE_MODU            (1)
+#endif
+
 #ifndef NDARRAY_HAS_INPLACE_MULTIPLY
 #define NDARRAY_HAS_INPLACE_MULTIPLY        (1)
 #endif
@@ -245,6 +253,10 @@
 #define NDARRAY_HAS_ITEMSIZE            (1)
 #endif
 
+#ifndef NDARRAY_HAS_NDIM
+#define NDARRAY_HAS_NDIM                (1)
+#endif
+
 #ifndef NDARRAY_HAS_RESHAPE
 #define NDARRAY_HAS_RESHAPE             (1)
 #endif

docs/manual/source/conf.py

@@ -27,7 +27,7 @@ copyright = '2019-2025, Zoltán Vörös and contributors'
 author = 'Zoltán Vörös'
 
 # The full version, including alpha/beta/rc tags
-release = '6.7.2'
+release = '6.9.0'
 
 
 # -- General configuration ---------------------------------------------------

docs/manual/source/ulab-ndarray.rst

@@ -1814,12 +1814,12 @@ array.
 Binary operators
 ================
 
-``ulab`` implements the ``+``, ``-``, ``*``, ``/``, ``**``, ``<``,
+``ulab`` implements the ``+``, ``-``, ``*``, ``/``, ``**``, ``%``,
-``>``, ``<=``, ``>=``, ``==``, ``!=``, ``+=``, ``-=``, ``*=``, ``/=``,
+``<``, ``>``, ``<=``, ``>=``, ``==``, ``!=``, ``+=``, ``-=``, ``*=``,
-``**=`` binary operators, as well as the ``AND``, ``OR``, ``XOR``
+``/=``, ``**=``, ``%=`` binary operators, as well as the ``AND``,
-bit-wise operators that work element-wise. Note that the bit-wise
+``OR``, ``XOR`` bit-wise operators that work element-wise. Note that the
-operators will raise an exception, if either of the operands is of
+bit-wise operators will raise an exception, if either of the operands is
-``float`` or ``complex`` type.
+of ``float`` or ``complex`` type.
 
 Broadcasting is available, meaning that the two operands do not even
 have to have the same shape. If the lengths along the respective axes

docs/ulab-change-log.md

@@ -1,3 +1,39 @@
+Fri, 06 Jun 2025
+
+version 6.8.0
+
+    expose ndim property
+
+Fri, 06 Jun 2025
+
+version 6.7.7
+
+    fix ndarray type inference for micropython objects
+
+Thu, 29 May 2025
+
+version 6.7.6
+
+    loadtxt can deal with multi-line comments
+
+Thu, 29 May 2025
+
+version 6.7.5
+
+    fix typo and shape in radnom module
+
+Sun, 16 Mar 2025
+
+version 6.7.4
+
+    re-name integration constants to avoid name clash with EPS ports
+
+Sun, 26 Jan 2025
+
+version 6.7.3
+
+    fix keepdims for min, max, argmin, argmax
+    
 Sun, 19 Jan 2025
 
 version 6.7.2

docs/ulab-convert.ipynb

@@ -14,7 +14,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 1,
+   "execution_count": 2,
    "metadata": {
     "ExecuteTime": {
      "end_time": "2022-02-09T06:27:15.118699Z",
@@ -61,7 +61,7 @@
     "author = 'Zoltán Vörös'\n",
     "\n",
     "# The full version, including alpha/beta/rc tags\n",
-    "release = '6.7.2'\n",
+    "release = '6.9.0'\n",
     "\n",
     "\n",
     "# -- General configuration ---------------------------------------------------\n",
@@ -217,7 +217,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 4,
+   "execution_count": 3,
    "metadata": {
     "ExecuteTime": {
      "end_time": "2022-02-09T06:27:21.647179Z",
@@ -258,7 +258,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 5,
+   "execution_count": 4,
    "metadata": {
     "ExecuteTime": {
      "end_time": "2022-02-09T06:27:42.024028Z",

docs/ulab-ndarray.ipynb

@@ -2599,7 +2599,7 @@
    "source": [
     "# Binary operators\n",
     "\n",
-    "`ulab` implements the `+`, `-`, `*`, `/`, `**`, `<`, `>`, `<=`, `>=`, `==`, `!=`, `+=`, `-=`, `*=`, `/=`, `**=`  binary operators, as well as the `AND`, `OR`, `XOR` bit-wise operators that work element-wise. Note that the bit-wise operators will raise an exception, if either of the operands is of `float` or `complex` type.\n",
+    "`ulab` implements the `+`, `-`, `*`, `/`, `**`, `%`, `<`, `>`, `<=`, `>=`, `==`, `!=`, `+=`, `-=`, `*=`, `/=`, `**=`, `%=` binary operators, as well as the `AND`, `OR`, `XOR` bit-wise operators that work element-wise. Note that the bit-wise operators will raise an exception, if either of the operands is of `float` or `complex` type.\n",
     "\n",
     "Broadcasting is available, meaning that the two operands do not even have to have the same shape. If the lengths along the respective axes are equal, or one of them is 1, or the axis is missing, the element-wise operation can still be carried out. \n",
     "A thorough explanation of broadcasting can be found under https://numpy.org/doc/stable/user/basics.broadcasting.html. \n",

docs/ulab-numerical.ipynb

@@ -50,7 +50,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 5,
+   "execution_count": 3,
    "metadata": {
     "ExecuteTime": {
      "end_time": "2022-01-07T19:16:37.453883Z",
@@ -245,115 +245,11 @@
     "\n",
     "**WARNING:** Difference to `numpy`: the `out` keyword argument is not implemented.\n",
     "\n",
-    "These functions follow the same pattern, and work with generic iterables, and `ndarray`s. `min`, and `max` return the minimum or maximum of a sequence. If the input array is two-dimensional, the `axis` keyword argument can be supplied, in which case the minimum/maximum along the given axis will be returned. If `axis=None` (this is also the default value), the minimum/maximum of the flattened array will be determined.\n",
+    "These functions follow the same pattern, and work with generic iterables, and `ndarray`s. `min`, and `max` return the minimum or maximum of a sequence. If the input array is two-dimensional, the `axis` keyword argument can be supplied, in which case the minimum/maximum along the given axis will be returned. If `axis=None` (this is also the default value), the minimum/maximum of the flattened array will be determined. The functions also accept the `keepdims=True` or `keepdims=False` keyword argument. The latter case is the default, while the former keeps the dimensions (the number of axes) of the supplied array. \n",
     "\n",
     "`argmin/argmax` return the position (index) of the minimum/maximum in the sequence."
    ]
   },
-  {
-   "cell_type": "code",
-   "execution_count": 8,
-   "metadata": {
-    "ExecuteTime": {
-     "end_time": "2020-10-17T21:26:22.507996Z",
-     "start_time": "2020-10-17T21:26:22.492543Z"
-    }
-   },
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "array([1.0, 2.0, 3.0], dtype=float64)\n",
-      "array([], dtype=float64)\n",
-      "[] 0\n",
-      "array([1.0, 2.0, 3.0], dtype=float64)\n",
-      "array([], dtype=float64)\n",
-      "\n",
-      "\n"
-     ]
-    }
-   ],
-   "source": [
-    "%%micropython -unix 1\n",
-    "\n",
-    "from ulab import numpy as np\n",
-    "\n",
-    "a = np.array([1, 2, 3])\n",
-    "print(a)\n",
-    "print(a[-1:-1:-3])\n",
-    "try:\n",
-    "    sa = list(a[-1:-1:-3])\n",
-    "    la = len(sa)\n",
-    "except IndexError as e:\n",
-    "    sa = str(e)\n",
-    "    la = -1\n",
-    "    \n",
-    "print(sa, la)\n",
-    "\n",
-    "a[-1:-1:-3] = np.ones(0)\n",
-    "print(a)\n",
-    "\n",
-    "b = np.ones(0) + 1\n",
-    "print(b)"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 9,
-   "metadata": {
-    "ExecuteTime": {
-     "end_time": "2020-10-17T21:54:49.123748Z",
-     "start_time": "2020-10-17T21:54:49.093819Z"
-    }
-   },
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "array([], dtype=float64)\n",
-      "\n",
-      "\n"
-     ]
-    }
-   ],
-   "source": [
-    "%%micropython -unix 1\n",
-    "\n",
-    "from ulab import numpy as np\n",
-    "\n",
-    "a = np.array([1, 2, 3])\n",
-    "print(a[0:1:-3])"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 81,
-   "metadata": {
-    "ExecuteTime": {
-     "end_time": "2020-10-17T20:59:58.285134Z",
-     "start_time": "2020-10-17T20:59:58.263605Z"
-    }
-   },
-   "outputs": [
-    {
-     "data": {
-      "text/plain": [
-       "(0,)"
-      ]
-     },
-     "execution_count": 81,
-     "metadata": {},
-     "output_type": "execute_result"
-    }
-   ],
-   "source": [
-    "a = np.array([1, 2, 3])\n",
-    "np.ones(0, dtype=uint8) / np.zeros(0, dtype=uint16)\n",
-    "np.ones(0).shape"
-   ]
-  },
   {
    "cell_type": "code",
    "execution_count": 10,
@@ -400,6 +296,43 @@
     "print('min of b (axis=1):', np.min(b, axis=1))"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "a: array([[0.0, 1.0, 2.0, 3.0],\n",
+      "       [4.0, 5.0, 6.0, 7.0],\n",
+      "       [8.0, 9.0, 10.0, 11.0]], dtype=float64)\n",
+      "\n",
+      "min of a (axis=1):\n",
+      " array([[0.0],\n",
+      "       [4.0],\n",
+      "       [8.0]], dtype=float64)\n",
+      "\n",
+      "min of a (axis=0):\n",
+      " array([[0.0, 1.0, 2.0, 3.0]], dtype=float64)\n",
+      "\n",
+      "\n"
+     ]
+    }
+   ],
+   "source": [
+    "%%micropython -unix 1\n",
+    "\n",
+    "from ulab import numpy as np\n",
+    "\n",
+    "a = np.array(range(12)).reshape((3, 4))\n",
+    "\n",
+    "print('a:', a)\n",
+    "print('\\nmin of a (axis=1):\\n', np.min(a, axis=1, keepdims=True))\n",
+    "print('\\nmin of a (axis=0):\\n', np.min(a, axis=0, keepdims=True))"
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {},

tests/2d/numpy/modulo.py

@@ -0,0 +1,26 @@
+try:
+    from ulab import numpy as np
+except:
+    import numpy as np
+
+
+dtypes = (np.uint8, np.int8, np.uint16, np.int16, np.float)
+
+for dtype1 in dtypes:
+    x1 = np.array(range(6), dtype=dtype1).reshape((2, 3))
+    for dtype2 in dtypes:
+        x2 = np.array(range(1, 4), dtype=dtype2)
+        print(x1 % x2)
+
+print()
+print('=' * 30)
+print('inplace modulo')
+print('=' * 30)
+print()
+
+for dtype1 in dtypes:
+    x1 = np.array(range(6), dtype=dtype1).reshape((2, 3))
+    for dtype2 in dtypes:
+        x2 = np.array(range(1, 4), dtype=dtype2)
+        x1 %= x2
+        print(x1)

tests/2d/numpy/modulo.py.exp

@@ -0,0 +1,105 @@
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=uint8)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=uint16)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int8)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0, 0, 1],
+       [0, 2, 0]], dtype=uint8)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=uint16)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+
+==============================
+inplace modulo
+==============================
+
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=uint8)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0, 0, 1],
+       [0, 2, 0]], dtype=uint8)
+array([[0, 0, 1],
+       [0, 0, 0]], dtype=int16)
+array([[0.0, 0.0, 1.0],
+       [0.0, 0.0, 0.0]], dtype=float64)
+array([[0.0, 0.0, 1.0],
+       [0.0, 0.0, 0.0]], dtype=float64)
+array([[0.0, 0.0, 1.0],
+       [0.0, 0.0, 0.0]], dtype=float64)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0, 1, 2],
+       [0, 0, 2]], dtype=int16)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)
+array([[0.0, 1.0, 2.0],
+       [0.0, 0.0, 2.0]], dtype=float64)

tests/2d/numpy/random.py

@@ -0,0 +1,10 @@
+try:
+    from ulab import numpy as np
+except ImportError:
+    import numpy as np
+
+rng = np.random.Generator(1234)
+
+for generator in (rng.normal, rng.random, rng.uniform):
+    random_array = generator(size=(1, 2))
+    print("array shape:", random_array.shape)

tests/2d/numpy/random.py.exp

@@ -0,0 +1,3 @@
+array shape: (1, 2)
+array shape: (1, 2)
+array shape: (1, 2)