ulab_ndarray_type is extern now

.github/workflows/build.yml

@@ -1,58 +0,0 @@
-name: Build CI
-
-on:
-  push:
-  pull_request:
-  release:
-    types: [published]
-  check_suite:
-    types: [rerequested]
-
-jobs:
-  test:
-    runs-on: ubuntu-16.04
-    steps:
-    - name: Dump GitHub context
-      env:
-        GITHUB_CONTEXT: ${{ toJson(github) }}
-      run: echo "$GITHUB_CONTEXT"
-    - name: Set up Python 3.5
-      uses: actions/setup-python@v1
-      with:
-        python-version: 3.5
-
-    - name: Versions
-      run: |
-        gcc --version
-        python3 --version
-    - name: Checkout ulab
-      uses: actions/checkout@v1
-
-    - name: Checkout micropython repo
-      uses: actions/checkout@v2
-      with:
-        repository: micropython/micropython
-        path: micropython
-
-    - name: Checkout micropython submodules
-      run: (cd micropython && git submodule update --init)
-
-    - name: Build mpy-cross
-      run: make -C micropython/mpy-cross -j2
-
-    - name: Build micropython unix port
-      run: |
-        make -C micropython/ports/unix -j2 deplibs
-        make -C micropython/ports/unix -j2 USER_C_MODULES=$(readlink -f .)
-
-    - name: Run tests
-      run: env MICROPYTHON_CPYTHON3=python3.5 MICROPY_MICROPYTHON=micropython/ports/unix/micropython micropython/tests/run-tests -d tests
-    - name: Print failure info
-      run: |
-        for exp in *.exp;
-        do testbase=$(basename $exp .exp);
-        echo -e "\nFAILURE $testbase";
-        diff -u $testbase.exp $testbase.out;
-        done
-      if: failure()
-

code/extras.c

@@ -1,33 +0,0 @@
-
-/*
- * This file is part of the micropython-ulab project,
- *
- * https://github.com/v923z/micropython-ulab
- *
- * The MIT License (MIT)
- *
- * Copyright (c) 2020 Zoltán Vörös
-*/
-
-#include <math.h>
-#include <stdlib.h>
-#include <string.h>
-#include "py/obj.h"
-#include "py/runtime.h"
-#include "py/misc.h"
-#include "extras.h"
-
-#if ULAB_EXTRAS_MODULE
-
-STATIC const mp_rom_map_elem_t ulab_filter_globals_table[] = {
-    { MP_OBJ_NEW_QSTR(MP_QSTR___name__), MP_OBJ_NEW_QSTR(MP_QSTR_extras) },
-};
-
-STATIC MP_DEFINE_CONST_DICT(mp_module_ulab_extras_globals, ulab_extras_globals_table);
-
-mp_obj_module_t ulab_filter_module = {
-    .base = { &mp_type_module },
-    .globals = (mp_obj_dict_t*)&mp_module_ulab_extras_globals,
-};
-
-#endif

code/extras.h

@@ -1,23 +0,0 @@
-
-/*
- * This file is part of the micropython-ulab project,
- *
- * https://github.com/v923z/micropython-ulab
- *
- * The MIT License (MIT)
- *
- * Copyright (c) 2020 Zoltán Vörös
-*/
-
-#ifndef _EXTRA_
-#define _EXTRA_
-
-#include "ulab.h"
-#include "ndarray.h"
-
-#if ULAB_EXTRAS_MODULE
-
-mp_obj_module_t ulab_extras_module;
-
-#endif
-#endif

code/fft.c

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,7 +5,7 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #include <math.h>
@@ -14,15 +13,12 @@
 #include <stdlib.h>
 #include <string.h>
 #include "py/runtime.h"
-#include "py/builtin.h"
 #include "py/binary.h"
 #include "py/obj.h"
 #include "py/objarray.h"
 #include "ndarray.h"
 #include "fft.h"
 
-#if ULAB_FFT_MODULE
-
 enum FFT_TYPE {
     FFT_FFT,
     FFT_IFFT,
@@ -80,74 +76,60 @@ void fft_kernel(mp_float_t *real, mp_float_t *imag, int n, int isign) {
 
 mp_obj_t fft_fft_ifft_spectrum(size_t n_args, mp_obj_t arg_re, mp_obj_t arg_im, uint8_t type) {
     if(!MP_OBJ_IS_TYPE(arg_re, &ulab_ndarray_type)) {
-        mp_raise_NotImplementedError(translate("FFT is defined for ndarrays only"));
+        mp_raise_NotImplementedError("FFT is defined for ndarrays only");
     } 
     if(n_args == 2) {
         if(!MP_OBJ_IS_TYPE(arg_im, &ulab_ndarray_type)) {
-            mp_raise_NotImplementedError(translate("FFT is defined for ndarrays only"));
+            mp_raise_NotImplementedError("FFT is defined for ndarrays only");
         }
     }
     // Check if input is of length of power of 2
     ndarray_obj_t *re = MP_OBJ_TO_PTR(arg_re);
     uint16_t len = re->array->len;
     if((len & (len-1)) != 0) {
-        mp_raise_ValueError(translate("input array length must be power of 2"));
+        // TODO: pad the input vector, if the length is not a power of 2
+        mp_raise_ValueError("input array length must be power of 2");
     }
     
-    ndarray_obj_t *out_re = create_new_ndarray(1, len, NDARRAY_FLOAT);
-    mp_float_t *data_re = (mp_float_t *)out_re->array->items;
+    ndarray_obj_t *ndarray_re = ndarray_new_linear_array(len, NDARRAY_FLOAT);
+    mp_float_t *data_re = (mp_float_t *)ndarray_re->array->items;
     
-    if(re->array->typecode == NDARRAY_FLOAT) { 
-        // By treating this case separately, we can save a bit of time.
-        // I don't know if it is worthwhile, though...
-        memcpy((mp_float_t *)out_re->array->items, (mp_float_t *)re->array->items, re->bytes);
-    } else {
     for(size_t i=0; i < len; i++) {
-            *data_re++ = ndarray_get_float_value(re->array->items, re->array->typecode, i);
+        data_re[i] = ndarray_get_float_value(re->array->items, re->array->typecode, re->offset+i*re->strides[0]);
     }
-        data_re -= len;
-    }
-    ndarray_obj_t *out_im = create_new_ndarray(1, len, NDARRAY_FLOAT);
-    mp_float_t *data_im = (mp_float_t *)out_im->array->items;
+    ndarray_obj_t *ndarray_im = ndarray_new_linear_array(len, NDARRAY_FLOAT);
+    mp_float_t *data_im = (mp_float_t *)ndarray_im->array->items;
 
     if(n_args == 2) {
         ndarray_obj_t *im = MP_OBJ_TO_PTR(arg_im);
-        if (re->array->len != im->array->len) {
-            mp_raise_ValueError(translate("real and imaginary parts must be of equal length"));
+        if (re->len != im->len) {
+            mp_raise_ValueError("real and imaginary parts must be of equal length");
         }
-        if(im->array->typecode == NDARRAY_FLOAT) {
-            memcpy((mp_float_t *)out_im->array->items, (mp_float_t *)im->array->items, im->bytes);
-        } else {
         for(size_t i=0; i < len; i++) {
-               *data_im++ = ndarray_get_float_value(im->array->items, im->array->typecode, i);
-            }
-            data_im -= len;
+            data_im[i] = ndarray_get_float_value(im->array->items, im->array->typecode, im->offset+i*im->strides[0]);
         }
     }
-
     if((type == FFT_FFT) || (type == FFT_SPECTRUM)) {
         fft_kernel(data_re, data_im, len, 1);
         if(type == FFT_SPECTRUM) {
             for(size_t i=0; i < len; i++) {
-                *data_re = MICROPY_FLOAT_C_FUN(sqrt)(*data_re * *data_re + *data_im * *data_im);
-                data_re++;
-                data_im++;
+                data_re[i] = MICROPY_FLOAT_C_FUN(sqrt)(data_re[i]*data_re[i] + data_im[i]*data_im[i]);
             }
         }
     } else { // inverse transform
         fft_kernel(data_re, data_im, len, -1);
         // TODO: numpy accepts the norm keyword argument
         for(size_t i=0; i < len; i++) {
-            *data_re++ /= len;
-            *data_im++ /= len;
+            data_re[i] /= len;
+            data_im[i] /= len;
         }
     }
     if(type == FFT_SPECTRUM) {
-        return MP_OBJ_TO_PTR(out_re);
+        return MP_OBJ_TO_PTR(ndarray_re);
     } else {
         mp_obj_t tuple[2];
-        tuple[0] = out_re;
-        tuple[1] = out_im;
+        tuple[0] = ndarray_re;
+        tuple[1] = ndarray_im;
         return mp_obj_new_tuple(2, tuple);
     }
 }
@@ -160,8 +142,6 @@ mp_obj_t fft_fft(size_t n_args, const mp_obj_t *args) {
     }
 }
 
-MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(fft_fft_obj, 1, 2, fft_fft);
-
 mp_obj_t fft_ifft(size_t n_args, const mp_obj_t *args) {
     if(n_args == 2) {
         return fft_fft_ifft_spectrum(n_args, args[0], args[1], FFT_IFFT);
@@ -170,8 +150,6 @@ mp_obj_t fft_ifft(size_t n_args, const mp_obj_t *args) {
     }
 }
 
-MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(fft_ifft_obj, 1, 2, fft_ifft);
-
 mp_obj_t fft_spectrum(size_t n_args, const mp_obj_t *args) {
     if(n_args == 2) {
         return fft_fft_ifft_spectrum(n_args, args[0], args[1], FFT_SPECTRUM);
@@ -179,23 +157,3 @@ mp_obj_t fft_spectrum(size_t n_args, const mp_obj_t *args) {
         return fft_fft_ifft_spectrum(n_args, args[0], mp_const_none, FFT_SPECTRUM);
     }
 }
-
-MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(fft_spectrum_obj, 1, 2, fft_spectrum);
-
-#if !CIRCUITPY
-STATIC const mp_rom_map_elem_t ulab_fft_globals_table[] = {
-    { MP_OBJ_NEW_QSTR(MP_QSTR___name__), MP_OBJ_NEW_QSTR(MP_QSTR_fft) },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_fft), (mp_obj_t)&fft_fft_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_ifft), (mp_obj_t)&fft_ifft_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_spectrum), (mp_obj_t)&fft_spectrum_obj },
-};
-
-STATIC MP_DEFINE_CONST_DICT(mp_module_ulab_fft_globals, ulab_fft_globals_table);
-
-mp_obj_module_t ulab_fft_module = {
-    .base = { &mp_type_module },
-    .globals = (mp_obj_dict_t*)&mp_module_ulab_fft_globals,
-};
-#endif
-
-#endif

code/fft.h

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,12 +5,11 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #ifndef _FFT_
 #define _FFT_
-#include "ulab.h"
 
 #ifndef MP_PI
 #define MP_PI MICROPY_FLOAT_CONST(3.14159265358979323846)
@@ -19,13 +17,8 @@
 
 #define SWAP(t, a, b) { t tmp = a; a = b; b = tmp; }
 
-#if ULAB_FFT_MODULE
-
-extern mp_obj_module_t ulab_fft_module;
-
-MP_DECLARE_CONST_FUN_OBJ_VAR_BETWEEN(fft_fft_obj);
-MP_DECLARE_CONST_FUN_OBJ_VAR_BETWEEN(fft_ifft_obj);
-MP_DECLARE_CONST_FUN_OBJ_VAR_BETWEEN(fft_spectrum_obj);
+mp_obj_t fft_fft(size_t , const mp_obj_t *);
+mp_obj_t fft_ifft(size_t , const mp_obj_t *);
+mp_obj_t fft_spectrum(size_t , const mp_obj_t *);
 
 #endif
-#endif

code/filter.c

@@ -1,101 +0,0 @@
-
-/*
- * This file is part of the micropython-ulab project,
- *
- * https://github.com/v923z/micropython-ulab
- *
- * The MIT License (MIT)
- *
- * Copyright (c) 2020 Jeff Epler for Adafruit Industries
-*/
-
-#include <math.h>
-#include <stdlib.h>
-#include <string.h>
-#include "py/obj.h"
-#include "py/runtime.h"
-#include "py/misc.h"
-#include "filter.h"
-
-#if ULAB_FILTER_MODULE
-mp_obj_t filter_convolve(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_a, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-        { MP_QSTR_v, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-    };
-
-    mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)];
-    mp_arg_parse_all(2, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args);
-
-    if(!MP_OBJ_IS_TYPE(args[0].u_obj, &ulab_ndarray_type) || !MP_OBJ_IS_TYPE(args[1].u_obj, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("convolve arguments must be ndarrays"));
-    }
-
-    ndarray_obj_t *a = MP_OBJ_TO_PTR(args[0].u_obj);
-    ndarray_obj_t *c = MP_OBJ_TO_PTR(args[1].u_obj);
-    int len_a = a->array->len;
-    int len_c = c->array->len;
-    // deal with linear arrays only
-    if(a->m*a->n != len_a || c->m*c->n != len_c) {
-        mp_raise_TypeError(translate("convolve arguments must be linear arrays"));
-    }
-    if(len_a == 0 || len_c == 0) {
-        mp_raise_TypeError(translate("convolve arguments must not be empty"));
-    }
-
-    int len = len_a + len_c - 1; // convolve mode "full"
-    ndarray_obj_t *out = create_new_ndarray(1, len, NDARRAY_FLOAT);
-    mp_float_t *outptr = out->array->items;
-    int off = len_c-1;
-
-    if(a->array->typecode == NDARRAY_FLOAT && c->array->typecode == NDARRAY_FLOAT) {
-        mp_float_t* a_items = (mp_float_t*)a->array->items;
-        mp_float_t* c_items = (mp_float_t*)c->array->items;
-        for(int k=-off; k<len-off; k++) {
-            mp_float_t accum = (mp_float_t)0;
-            int top_n = MIN(len_c, len_a - k);
-            int bot_n = MAX(-k, 0);
-            mp_float_t* a_ptr = a_items + bot_n + k;
-            mp_float_t* a_end = a_ptr + (top_n - bot_n);
-            mp_float_t* c_ptr = c_items + len_c - bot_n - 1;
-            for(; a_ptr != a_end;) {
-                accum += *a_ptr++ * *c_ptr--;
-            }
-            *outptr++ = accum;
-        }
-    } else {
-        for(int k=-off; k<len-off; k++) {
-            mp_float_t accum = (mp_float_t)0;
-            int top_n = MIN(len_c, len_a - k);
-            int bot_n = MAX(-k, 0);
-            for(int n=bot_n; n<top_n; n++) {
-                int idx_c = len_c - n - 1;
-                int idx_a = n+k;
-                mp_float_t ai = ndarray_get_float_value(a->array->items, a->array->typecode, idx_a);
-                mp_float_t ci = ndarray_get_float_value(c->array->items, c->array->typecode, idx_c);
-                accum += ai * ci;
-            }
-            *outptr++ = accum;
-        }
-    }
-
-    return out;
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(filter_convolve_obj, 2, filter_convolve);
-
-#if !CIRCUITPY
-STATIC const mp_rom_map_elem_t ulab_filter_globals_table[] = {
-    { MP_OBJ_NEW_QSTR(MP_QSTR___name__), MP_OBJ_NEW_QSTR(MP_QSTR_filter) },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_convolve), (mp_obj_t)&filter_convolve_obj },
-};
-
-STATIC MP_DEFINE_CONST_DICT(mp_module_ulab_filter_globals, ulab_filter_globals_table);
-
-mp_obj_module_t ulab_filter_module = {
-    .base = { &mp_type_module },
-    .globals = (mp_obj_dict_t*)&mp_module_ulab_filter_globals,
-};
-#endif
-
-#endif

code/filter.h

@@ -1,25 +0,0 @@
-
-/*
- * This file is part of the micropython-ulab project,
- *
- * https://github.com/v923z/micropython-ulab
- *
- * The MIT License (MIT)
- *
- * Copyright (c) 2020 Jeff Epler for Adafruit Industries
-*/
-
-#ifndef _FILTER_
-#define _FILTER_
-
-#include "ulab.h"
-#include "ndarray.h"
-
-#if ULAB_FILTER_MODULE
-
-extern mp_obj_module_t ulab_filter_module;
-
-MP_DECLARE_CONST_FUN_OBJ_KW(filter_convolve_obj);
-
-#endif
-#endif

code/linalg.c

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,7 +5,7 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #include <stdlib.h>
@@ -17,48 +16,6 @@
 #include "py/misc.h"
 #include "linalg.h"
 
-#if ULAB_LINALG_MODULE
-
-mp_obj_t linalg_size(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-        { MP_QSTR_axis, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-    };
-
-    mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)];
-    mp_arg_parse_all(1, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args);
-
-    if(!MP_OBJ_IS_TYPE(args[0].u_obj, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("size is defined for ndarrays only"));
-    } else {
-        ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(args[0].u_obj);
-        if(args[1].u_obj == mp_const_none) {
-            return mp_obj_new_int(ndarray->array->len);
-        } else if(MP_OBJ_IS_INT(args[1].u_obj)) {
-            uint8_t ax = mp_obj_get_int(args[1].u_obj);
-            if(ax == 0) {
-                if(ndarray->m == 1) {
-                    return mp_obj_new_int(ndarray->n);
-                } else {
-                    return mp_obj_new_int(ndarray->m);                    
-                }
-            } else if(ax == 1) {
-                if(ndarray->m == 1) {
-                    mp_raise_ValueError(translate("tuple index out of range"));
-                } else {
-                    return mp_obj_new_int(ndarray->n);
-                }
-            } else {
-                    mp_raise_ValueError(translate("tuple index out of range"));
-            }
-        } else {
-            mp_raise_TypeError(translate("wrong argument type"));
-        }
-    }
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(linalg_size_obj, 1, linalg_size);
-
 bool linalg_invert_matrix(mp_float_t *data, size_t N) {
     // returns true, of the inversion was successful, 
     // false, if the matrix is singular
@@ -102,49 +59,46 @@ bool linalg_invert_matrix(mp_float_t *data, size_t N) {
 }
 
 mp_obj_t linalg_inv(mp_obj_t o_in) {
-    // since inv is not a class method, we have to inspect the input argument first
     if(!MP_OBJ_IS_TYPE(o_in, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("only ndarrays can be inverted"));
+        mp_raise_TypeError("only ndarray objects can be inverted");
     }
-    ndarray_obj_t *o = MP_OBJ_TO_PTR(o_in);
-    if(!MP_OBJ_IS_TYPE(o_in, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("only ndarray objects can be inverted"));
+    ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(o_in);
+    if(ndarray->ndim != 2) {
+        mp_raise_ValueError("only two-dimensional tensors can be inverted");
     }
-    if(o->m != o->n) {
-        mp_raise_ValueError(translate("only square matrices can be inverted"));
+    if(ndarray->shape[0] != ndarray->shape[1]) {
+        mp_raise_ValueError("only square matrices can be inverted");
     }
-    ndarray_obj_t *inverted = create_new_ndarray(o->m, o->n, NDARRAY_FLOAT);
+    size_t *shape = m_new(size_t, 2);
+    shape[0] = shape[1] = ndarray->shape[0];
+    ndarray_obj_t *inverted = ndarray_new_dense_ndarray(2, shape, NDARRAY_FLOAT);
     mp_float_t *data = (mp_float_t *)inverted->array->items;
     mp_obj_t elem;
-    for(size_t m=0; m < o->m; m++) { // rows first
-        for(size_t n=0; n < o->n; n++) { // columns next
+    for(size_t m=0; m < ndarray->shape[0]; m++) { // rows first
+        for(size_t n=0; n < ndarray->shape[1]; n++) { // columns next
             // this could, perhaps, be done in single line... 
             // On the other hand, we probably spend little time here
-            elem = mp_binary_get_val_array(o->array->typecode, o->array->items, m*o->n+n);
-            data[m*o->n+n] = (mp_float_t)mp_obj_get_float(elem);
+            elem = mp_binary_get_val_array(ndarray->array->typecode, ndarray->array->items, m*ndarray->shape[1]+n);
+            data[m*ndarray->shape[1]+n] = (mp_float_t)mp_obj_get_float(elem);
         }
     }
     
-    if(!linalg_invert_matrix(data, o->m)) {
-        // TODO: I am not sure this is needed here. Otherwise, 
-        // how should we free up the unused RAM of inverted?
-        m_del(mp_float_t, inverted->array->items, o->n*o->n);
-        mp_raise_ValueError(translate("input matrix is singular"));
+    if(!linalg_invert_matrix(data, ndarray->shape[0])) {
+        // TODO: I am not sure this is needed here. Otherwise, how should we free up the unused RAM of inverted?
+        m_del(mp_float_t, inverted->array->items, ndarray->shape[0]*ndarray->shape[1]);
+        mp_raise_ValueError("input matrix is singular");
     }
     return MP_OBJ_FROM_PTR(inverted);
 }
 
-MP_DEFINE_CONST_FUN_OBJ_1(linalg_inv_obj, linalg_inv);
-
 mp_obj_t linalg_dot(mp_obj_t _m1, mp_obj_t _m2) {
+    return mp_const_none;
+    /*
     // TODO: should the results be upcast?
-    if(!MP_OBJ_IS_TYPE(_m1, &ulab_ndarray_type) || !MP_OBJ_IS_TYPE(_m2, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("arguments must be ndarrays"));
-    }
     ndarray_obj_t *m1 = MP_OBJ_TO_PTR(_m1);
     ndarray_obj_t *m2 = MP_OBJ_TO_PTR(_m2);    
     if(m1->n != m2->m) {
-        mp_raise_ValueError(translate("matrix dimensions do not match"));
+        mp_raise_ValueError("matrix dimensions do not match");
     }
     // TODO: numpy uses upcasting here
     ndarray_obj_t *out = create_new_ndarray(m1->m, m2->n, NDARRAY_FLOAT);
@@ -159,14 +113,12 @@ mp_obj_t linalg_dot(mp_obj_t _m1, mp_obj_t _m2) {
                 v2 = ndarray_get_float_value(m2->array->items, m2->array->typecode, k*m2->n+j);
                 sum += v1 * v2;
             }
-            outdata[j*m1->m+i] = sum;
+            outdata[i*m1->m+j] = sum;
         }
     }
-    return MP_OBJ_FROM_PTR(out);
+    return MP_OBJ_FROM_PTR(out); */
 } 
 
-MP_DEFINE_CONST_FUN_OBJ_2(linalg_dot_obj, linalg_dot);
-
 mp_obj_t linalg_zeros_ones(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args, uint8_t kind) {
     static const mp_arg_t allowed_args[] = {
         { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_obj = MP_OBJ_NULL} } ,
@@ -177,20 +129,22 @@ mp_obj_t linalg_zeros_ones(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw
     mp_arg_parse_all(n_args, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args);
     
     uint8_t dtype = args[1].u_int;
-    if(!MP_OBJ_IS_INT(args[0].u_obj) && !MP_OBJ_IS_TYPE(args[0].u_obj, &mp_type_tuple)) {
-        mp_raise_TypeError(translate("input argument must be an integer or a 2-tuple"));
+    if(!mp_obj_is_int(args[0].u_obj) && !mp_obj_is_type(args[0].u_obj, &mp_type_tuple)) {
+        mp_raise_TypeError("input argument must be an integer or a tuple");
     }
     ndarray_obj_t *ndarray = NULL;
-    if(MP_OBJ_IS_INT(args[0].u_obj)) {
+    if(mp_obj_is_int(args[0].u_obj)) {
         size_t n = mp_obj_get_int(args[0].u_obj);
-        ndarray = create_new_ndarray(1, n, dtype);
-    } else if(MP_OBJ_IS_TYPE(args[0].u_obj, &mp_type_tuple)) {
+        size_t *shape = m_new(size_t, 1);
+        shape[0] = n;
+        ndarray = ndarray_new_dense_ndarray(1, shape, dtype);
+    } else if(mp_obj_is_type(args[0].u_obj, &mp_type_tuple)) {
         mp_obj_tuple_t *tuple = MP_OBJ_TO_PTR(args[0].u_obj);
-        if(tuple->len != 2) {
-            mp_raise_TypeError(translate("input argument must be an integer or a 2-tuple"));
+        size_t *shape = m_new(size_t, tuple->len);
+        for(uint8_t i=0; i < tuple->len; i++) {
+            shape[i] = mp_obj_get_int(tuple->items[i]);
         }
-        ndarray = create_new_ndarray(mp_obj_get_int(tuple->items[0]), 
-                                                  mp_obj_get_int(tuple->items[1]), dtype);
+        ndarray = ndarray_new_dense_ndarray(tuple->len, shape, dtype);
     }
     if(kind == 1) {
         mp_obj_t one = mp_obj_new_int(1);
@@ -205,18 +159,17 @@ mp_obj_t linalg_zeros(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args
     return linalg_zeros_ones(n_args, pos_args, kw_args, 0);
 }
 
-MP_DEFINE_CONST_FUN_OBJ_KW(linalg_zeros_obj, 0, linalg_zeros);
-
 mp_obj_t linalg_ones(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
     return linalg_zeros_ones(n_args, pos_args, kw_args, 1);
 }
 
-MP_DEFINE_CONST_FUN_OBJ_KW(linalg_ones_obj, 0, linalg_ones);
-
 mp_obj_t linalg_eye(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
+    // TODO: this is a bit more problematic in higher dimensions
+    return mp_const_none;
+    /*
     static const mp_arg_t allowed_args[] = {
         { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_INT, {.u_int = 0} },
-        { MP_QSTR_M, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
+        { MP_QSTR_M, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_rom_obj = MP_ROM_PTR(&mp_const_none_obj) } },
         { MP_QSTR_k, MP_ARG_KW_ONLY | MP_ARG_INT, {.u_int = 0} },        
         { MP_QSTR_dtype, MP_ARG_KW_ONLY | MP_ARG_INT, {.u_int = NDARRAY_FLOAT} },
     };
@@ -251,93 +204,98 @@ mp_obj_t linalg_eye(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args)
             i++;
         }
     }
-    return MP_OBJ_FROM_PTR(ndarray);
+    return MP_OBJ_FROM_PTR(ndarray); */
 }
 
-MP_DEFINE_CONST_FUN_OBJ_KW(linalg_eye_obj, 0, linalg_eye);
-
 mp_obj_t linalg_det(mp_obj_t oin) {
-    if(!MP_OBJ_IS_TYPE(oin, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("function defined for ndarrays only"));
+    if(!mp_obj_is_type(oin, &ulab_ndarray_type)) {
+        mp_raise_TypeError("function defined for ndarrays only");
     }
-    ndarray_obj_t *in = MP_OBJ_TO_PTR(oin);
-    if(in->m != in->n) {
-        mp_raise_ValueError(translate("input must be square matrix"));
+    ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(oin);
+    if(ndarray->ndim != 2) {
+        mp_raise_ValueError("only two-dimensional tensors can be inverted");
+    }
+    if(ndarray->shape[0] != ndarray->shape[1]) {
+        mp_raise_ValueError("only square matrices can be inverted");
     }
     
-    mp_float_t *tmp = m_new(mp_float_t, in->n*in->n);
-    for(size_t i=0; i < in->array->len; i++){
-        tmp[i] = ndarray_get_float_value(in->array->items, in->array->typecode, i);
+    mp_float_t *tmp = m_new(mp_float_t, ndarray->shape[0]*ndarray->shape[1]);
+    // TODO: this won't work for sliced arrays
+    for(size_t i=0; i < ndarray->len; i++){
+        tmp[i] = ndarray_get_float_value(ndarray->array->items, ndarray->array->typecode, i);
     }
     mp_float_t c;
-    for(size_t m=0; m < in->m-1; m++){
-        if(MICROPY_FLOAT_C_FUN(fabs)(tmp[m*(in->n+1)]) < epsilon) {
-            m_del(mp_float_t, tmp, in->n*in->n);
+    for(size_t m=0; m < ndarray->shape[0]-1; m++){
+        if(MICROPY_FLOAT_C_FUN(fabs)(tmp[m*(ndarray->shape[1]+1)]) < epsilon) {
+            m_del(mp_float_t, tmp, ndarray->shape[0]*ndarray->shape[1]);
             return mp_obj_new_float(0.0);
         }
-        for(size_t n=0; n < in->n; n++){
+        for(size_t n=0; n < ndarray->shape[1]; n++){
             if(m != n) {
-                c = tmp[in->n*n+m] / tmp[m*(in->n+1)];
-                for(size_t k=0; k < in->n; k++){
-                    tmp[in->n*n+k] -= c * tmp[in->n*m+k];
+                c = tmp[ndarray->shape[0]*n+m] / tmp[m*(ndarray->shape[1]+1)];
+                for(size_t k=0; k < ndarray->shape[1]; k++){
+                    tmp[ndarray->shape[1]*n+k] -= c * tmp[ndarray->shape[1]*m+k];
                 }
             }
         }
     }
     mp_float_t det = 1.0;
                             
-    for(size_t m=0; m < in->m; m++){ 
-        det *= tmp[m*(in->n+1)];
+    for(size_t m=0; m < ndarray->shape[0]; m++){ 
+        det *= tmp[m*(ndarray->shape[1]+1)];
     }
-    m_del(mp_float_t, tmp, in->n*in->n);
+    m_del(mp_float_t, tmp, ndarray->shape[0]*ndarray->shape[1]);
     return mp_obj_new_float(det);
 }
 
-MP_DEFINE_CONST_FUN_OBJ_1(linalg_det_obj, linalg_det);
-
 mp_obj_t linalg_eig(mp_obj_t oin) {
-    if(!MP_OBJ_IS_TYPE(oin, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("function defined for ndarrays only"));
+    if(!mp_obj_is_type(oin, &ulab_ndarray_type)) {
+        mp_raise_TypeError("function defined for ndarrays only");
     }
-    ndarray_obj_t *in = MP_OBJ_TO_PTR(oin);
-    if(in->m != in->n) {
-        mp_raise_ValueError(translate("input must be square matrix"));
+    ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(oin);
+    if(ndarray->ndim != 2) {
+        mp_raise_ValueError("only two-dimensional tensors can be inverted");
     }
-    mp_float_t *array = m_new(mp_float_t, in->array->len);
-    for(size_t i=0; i < in->array->len; i++) {
-        array[i] = ndarray_get_float_value(in->array->items, in->array->typecode, i);
+    if(ndarray->shape[0] != ndarray->shape[1]) {
+        mp_raise_ValueError("only square matrices can be inverted");
+    }
+    mp_float_t *array = m_new(mp_float_t, ndarray->len);
+    // TODO: this won't work for sliced arrays
+    for(size_t i=0; i < ndarray->len; i++) {
+        array[i] = ndarray_get_float_value(ndarray->array->items, ndarray->array->typecode, i);
     }
     // make sure the matrix is symmetric
-    for(size_t m=0; m < in->m; m++) {
-        for(size_t n=m+1; n < in->n; n++) {
+    for(size_t m=0; m < ndarray->shape[0]; m++) {
+        for(size_t n=m+1; n < ndarray->shape[1]; n++) {
             // compare entry (m, n) to (n, m)
             // TODO: this must probably be scaled!
-            if(epsilon < MICROPY_FLOAT_C_FUN(fabs)(array[m*in->n + n] - array[n*in->n + m])) {
-                mp_raise_ValueError(translate("input matrix is asymmetric"));
+            if(epsilon < MICROPY_FLOAT_C_FUN(fabs)(array[m*ndarray->shape[0] + n] - array[n*ndarray->shape[0] + m])) {
+                mp_raise_ValueError("input matrix is asymmetric");
             }
         }
     }
     
     // if we got this far, then the matrix will be symmetric
-    
-    ndarray_obj_t *eigenvectors = create_new_ndarray(in->m, in->n, NDARRAY_FLOAT);
+    size_t *shape = m_new(size_t, 2);
+    shape[0] = shape[1] = ndarray->shape[0];
+    ndarray_obj_t *eigenvectors = ndarray_new_dense_ndarray(2, shape, NDARRAY_FLOAT);
     mp_float_t *eigvectors = (mp_float_t *)eigenvectors->array->items;
     // start out with the unit matrix
-    for(size_t m=0; m < in->m; m++) {
-        eigvectors[m*(in->n+1)] = 1.0;
+    for(size_t m=0; m < ndarray->shape[0]; m++) {
+        eigvectors[m*(ndarray->shape[1]+1)] = 1.0;
     }
     mp_float_t largest, w, t, c, s, tau, aMk, aNk, vm, vn;
     size_t M, N;
-    size_t iterations = JACOBI_MAX*in->n*in->n;
+    size_t iterations = JACOBI_MAX*ndarray->shape[0]*ndarray->shape[0];
     do {
         iterations--;
         // find the pivot here
         M = 0;
         N = 0;
         largest = 0.0;
-        for(size_t m=0; m < in->m-1; m++) { // -1: no need to inspect last row
-            for(size_t n=m+1; n < in->n; n++) {
-                w = MICROPY_FLOAT_C_FUN(fabs)(array[m*in->n + n]);
+        for(size_t m=0; m < ndarray->shape[0]-1; m++) { // -1: no need to inspect last row
+            for(size_t n=m+1; n < ndarray->shape[0]; n++) {
+                w = MICROPY_FLOAT_C_FUN(fabs)(array[m*ndarray->shape[0] + n]);
                 if((largest < w) && (epsilon < w)) {
                     M = m;
                     N = n;
@@ -350,7 +308,7 @@ mp_obj_t linalg_eig(mp_obj_t oin) {
         }
         // at this point, we have the pivot, and it is the entry (M, N)
         // now we have to find the rotation angle
-        w = (array[N*in->n + N] - array[M*in->n + M]) / (2.0*array[M*in->n + N]);
+        w = (array[N*ndarray->shape[0] + N] - array[M*ndarray->shape[0] + M]) / (2.0*array[M*ndarray->shape[0] + N]);
         // The following if/else chooses the smaller absolute value for the tangent 
         // of the rotation angle. Going with the smaller should be numerically stabler.
         if(w > 0) {
@@ -365,26 +323,26 @@ mp_obj_t linalg_eig(mp_obj_t oin) {
         // at this point, we have the rotation angles, so we can transform the matrix
         // first the two diagonal elements
         // a(M, M) = a(M, M) - t*a(M, N)
-        array[M*in->n + M] = array[M*in->n + M] - t * array[M*in->n + N];
+        array[M*ndarray->shape[0] + M] = array[M*ndarray->shape[0] + M] - t * array[M*ndarray->shape[0] + N];
         // a(N, N) = a(N, N) + t*a(M, N)
-        array[N*in->n + N] = array[N*in->n + N] + t * array[M*in->n + N];
+        array[N*ndarray->shape[0] + N] = array[N*ndarray->shape[0] + N] + t * array[M*ndarray->shape[0] + N];
         // after the rotation, the a(M, N), and a(N, M) entries should become zero
-        array[M*in->n + N] = array[N*in->n + M] = 0.0;
+        array[M*ndarray->shape[0] + N] = array[N*ndarray->shape[0] + M] = 0.0;
         // then all other elements in the column
-        for(size_t k=0; k < in->m; k++) {
+        for(size_t k=0; k < ndarray->shape[0]; k++) {
             if((k == M) || (k == N)) {
                 continue;
             }
-            aMk = array[M*in->n + k];
-            aNk = array[N*in->n + k];
+            aMk = array[M*ndarray->shape[0] + k];
+            aNk = array[N*ndarray->shape[0] + k];
             // a(M, k) = a(M, k) - s*(a(N, k) + tau*a(M, k))
-            array[M*in->n + k] -= s*(aNk + tau*aMk);
+            array[M*ndarray->shape[0] + k] -= s*(aNk + tau*aMk);
             // a(N, k) = a(N, k) + s*(a(M, k) - tau*a(N, k))
-            array[N*in->n + k] += s*(aMk - tau*aNk);
+            array[N*ndarray->shape[0] + k] += s*(aMk - tau*aNk);
             // a(k, M) = a(M, k)
-            array[k*in->n + M] = array[M*in->n + k];
+            array[k*ndarray->shape[0] + M] = array[M*ndarray->shape[0] + k];
             // a(k, N) = a(N, k)
-            array[k*in->n + N] = array[N*in->n + k];
+            array[k*ndarray->shape[0] + N] = array[N*ndarray->shape[0] + k];
         }
         // now we have to update the eigenvectors
         // the rotation matrix, R, multiplies from the right
@@ -393,27 +351,29 @@ mp_obj_t linalg_eig(mp_obj_t oin) {
         // R(N, M) = s
         // (M, N) = -s
         // entries. This means that only the Mth, and Nth columns will change
-        for(size_t m=0; m < in->m; m++) {
-            vm = eigvectors[m*in->n+M];
-            vn = eigvectors[m*in->n+N];
+        for(size_t m=0; m < ndarray->shape[0]; m++) {
+            vm = eigvectors[m*ndarray->shape[0]+M];
+            vn = eigvectors[m*ndarray->shape[0]+N];
             // the new value of eigvectors(m, M)
-            eigvectors[m*in->n+M] = c * vm - s * vn;
+            eigvectors[m*ndarray->shape[0]+M] = c * vm - s * vn;
             // the new value of eigvectors(m, N)
-            eigvectors[m*in->n+N] = s * vm + c * vn;
+            eigvectors[m*ndarray->shape[0]+N] = s * vm + c * vn;
         }
     } while(iterations > 0);
     
     if(iterations == 0) { 
         // the computation did not converge; numpy raises LinAlgError
-        m_del(mp_float_t, array, in->array->len);
-        mp_raise_ValueError(translate("iterations did not converge"));
+        m_del(mp_float_t, array, ndarray->len);
+        mp_raise_ValueError("iterations did not converge");
     }
-    ndarray_obj_t *eigenvalues = create_new_ndarray(1, in->n, NDARRAY_FLOAT);
+    size_t *eigen_shape = m_new(size_t, 1);
+    eigen_shape[0] = ndarray->shape[0];
+    ndarray_obj_t *eigenvalues = ndarray_new_dense_ndarray(1, eigen_shape, NDARRAY_FLOAT);
     mp_float_t *eigvalues = (mp_float_t *)eigenvalues->array->items;
-    for(size_t i=0; i < in->n; i++) {
-        eigvalues[i] = array[i*(in->n+1)];
+    for(size_t i=0; i < ndarray->shape[0]; i++) {
+        eigvalues[i] = array[i*(ndarray->shape[0]+1)];
     }
-    m_del(mp_float_t, array, in->array->len);
+    m_del(mp_float_t, array, ndarray->len);
     
     mp_obj_tuple_t *tuple = MP_OBJ_TO_PTR(mp_obj_new_tuple(2, NULL));
     tuple->items[0] = MP_OBJ_FROM_PTR(eigenvalues);
@@ -421,28 +381,3 @@ mp_obj_t linalg_eig(mp_obj_t oin) {
     return tuple;
     return MP_OBJ_FROM_PTR(eigenvalues);
 }
-
-MP_DEFINE_CONST_FUN_OBJ_1(linalg_eig_obj, linalg_eig);
-
-#if !CIRCUITPY
-STATIC const mp_rom_map_elem_t ulab_linalg_globals_table[] = {
-    { MP_OBJ_NEW_QSTR(MP_QSTR___name__), MP_OBJ_NEW_QSTR(MP_QSTR_linalg) },
-    { MP_ROM_QSTR(MP_QSTR_size), (mp_obj_t)&linalg_size_obj },
-    { MP_ROM_QSTR(MP_QSTR_inv), (mp_obj_t)&linalg_inv_obj },
-    { MP_ROM_QSTR(MP_QSTR_dot), (mp_obj_t)&linalg_dot_obj },
-    { MP_ROM_QSTR(MP_QSTR_zeros), (mp_obj_t)&linalg_zeros_obj },
-    { MP_ROM_QSTR(MP_QSTR_ones), (mp_obj_t)&linalg_ones_obj },
-    { MP_ROM_QSTR(MP_QSTR_eye), (mp_obj_t)&linalg_eye_obj },
-    { MP_ROM_QSTR(MP_QSTR_det), (mp_obj_t)&linalg_det_obj },
-    { MP_ROM_QSTR(MP_QSTR_eig), (mp_obj_t)&linalg_eig_obj },    
-};
-
-STATIC MP_DEFINE_CONST_DICT(mp_module_ulab_linalg_globals, ulab_linalg_globals_table);
-
-mp_obj_module_t ulab_linalg_module = {
-    .base = { &mp_type_module },
-    .globals = (mp_obj_dict_t*)&mp_module_ulab_linalg_globals,
-};
-#endif
-
-#endif

code/linalg.h

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,13 +5,12 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #ifndef _LINALG_
 #define _LINALG_
 
-#include "ulab.h"
 #include "ndarray.h"
 
 #if MICROPY_FLOAT_IMPL == MICROPY_FLOAT_IMPL_FLOAT
@@ -23,13 +21,14 @@
 
 #define JACOBI_MAX     20
 
-#if ULAB_LINALG_MODULE || ULAB_POLY_MODULE
 bool linalg_invert_matrix(mp_float_t *, size_t );
-#endif
+mp_obj_t linalg_inv(mp_obj_t );
+mp_obj_t linalg_dot(mp_obj_t , mp_obj_t );
+mp_obj_t linalg_zeros(size_t , const mp_obj_t *, mp_map_t *);
+mp_obj_t linalg_ones(size_t , const mp_obj_t *, mp_map_t *);
+mp_obj_t linalg_eye(size_t , const mp_obj_t *, mp_map_t *);
 
-#if ULAB_LINALG_MODULE
-
-extern mp_obj_module_t ulab_linalg_module;
+mp_obj_t linalg_det(mp_obj_t );
+mp_obj_t linalg_eig(mp_obj_t );
 
 #endif
-#endif

code/micropython.mk

@@ -3,17 +3,11 @@ USERMODULES_DIR := $(USERMOD_DIR)
 
 # Add all C files to SRC_USERMOD.
 SRC_USERMOD += $(USERMODULES_DIR)/ndarray.c
-SRC_USERMOD += $(USERMODULES_DIR)/linalg.c
 SRC_USERMOD += $(USERMODULES_DIR)/vectorise.c
+SRC_USERMOD += $(USERMODULES_DIR)/linalg.c
 SRC_USERMOD += $(USERMODULES_DIR)/poly.c
-SRC_USERMOD += $(USERMODULES_DIR)/fft.c
 SRC_USERMOD += $(USERMODULES_DIR)/numerical.c
-SRC_USERMOD += $(USERMODULES_DIR)/filter.c
-SRC_USERMOD += $(USERMODULES_DIR)/extras.c
+SRC_USERMOD += $(USERMODULES_DIR)/fft.c
 SRC_USERMOD += $(USERMODULES_DIR)/ulab.c
 
-# We can add our module folder to include paths if needed
-# This is not actually needed in this example.
 CFLAGS_USERMOD += -I$(USERMODULES_DIR)
-
-CFLAGS_EXTRA = -DMODULE_ULAB_ENABLED=1

code/ndarray.h

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,7 +5,7 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #ifndef _NDARRAY_
@@ -17,7 +16,10 @@
 #include "py/objstr.h"
 #include "py/objlist.h"
 
-#define PRINT_MAX  10
+#define SWAP(t, a, b) { t tmp = a; a = b; b = tmp; }
+
+#define NDARRAY_NUMERIC   0
+#define NDARRAY_BOOLEAN   1
 
 #if MICROPY_FLOAT_IMPL == MICROPY_FLOAT_IMPL_FLOAT
 #define FLOAT_TYPECODE 'f'
@@ -25,15 +27,10 @@
 #define FLOAT_TYPECODE 'd'
 #endif
 
-#if !CIRCUITPY
-#define translate(x) x
-#endif
-
-#define SWAP(t, a, b) { t tmp = a; a = b; b = tmp; }
-
 extern const mp_obj_type_t ulab_ndarray_type;
 
 enum NDARRAY_TYPE {
+    NDARRAY_BOOL = '?', // this must never be assigned to the typecode!
     NDARRAY_UINT8 = 'B',
     NDARRAY_INT8 = 'b',
     NDARRAY_UINT16 = 'H', 
@@ -43,103 +40,91 @@ enum NDARRAY_TYPE {
 
 typedef struct _ndarray_obj_t {
     mp_obj_base_t base;
-    size_t m, n;
+    uint8_t boolean;
+    uint8_t ndim;
+    size_t *shape;
+    int32_t *strides;
     size_t len;
+    size_t offset;
     mp_obj_array_t *array;
-    size_t bytes;
 } ndarray_obj_t;
 
-mp_obj_t mp_obj_new_ndarray_iterator(mp_obj_t , size_t , mp_obj_iter_buf_t *);
+// this is a helper structure, so that we can return shape AND strides from a function
+typedef struct _ndarray_header_obj_t {
+    size_t *shape;
+    int32_t *strides;
+    int8_t axis;
+    size_t offset;
+} ndarray_header_obj_t;
 
+// various helper functions
+size_t ndarray_index_from_flat(size_t , ndarray_obj_t *, int32_t *);
+size_t ndarray_index_from_contracted(size_t  , ndarray_obj_t * , int32_t * , uint8_t  , uint8_t  );
 mp_float_t ndarray_get_float_value(void *, uint8_t , size_t );
-void fill_array_iterable(mp_float_t *, mp_obj_t );
 
-void ndarray_print_row(const mp_print_t *, mp_obj_array_t *, size_t , size_t );
+// calculates the index (in the original linear array) of an item, if the index in the flat array is given
+// this is the macro equivalent of ndarray_index_from_flat()
+// TODO: This fails, when the last stride is not 1!!!
+#define NDARRAY_INDEX_FROM_FLAT(ndarray, shape_strides, index, _tindex, _nindex) do {\
+    size_t Q;\
+    (_tindex) = (index);\
+    (_nindex) = (ndarray)->offset;\
+    for(size_t _x=0; _x < (ndarray)->ndim; _x++) {\
+        Q = (_tindex) / (shape_strides)[_x];\
+        (_tindex) -= Q * (shape_strides)[_x];\
+        (_nindex) += Q * (ndarray)->strides[_x];\
+    }\
+} while(0)
+
+#define NDARRAY_INDEX_FROM_FLAT2(ndarray, stride_array, shape_strides, index, _tindex, _nindex) do {\
+    size_t Q;\
+    (_tindex) = (index);\
+    (_nindex) = (ndarray)->offset;\
+    for(size_t _x=0; _x < (ndarray)->ndim; _x++) {\
+        Q = (_tindex) / (shape_strides)[_x];\
+        (_tindex) -= Q * (shape_strides)[_x];\
+        (_nindex) += Q * (stride_array)[_x];\
+    }\
+} while(0)
+    
+#define CREATE_SINGLE_ITEM(ndarray, type, typecode, value) do {\
+    (ndarray) = ndarray_new_linear_array(1, (typecode));\
+    type *tmparr = (type *)(ndarray)->array->items;\
+    tmparr[0] = (type)(value);\
+} while(0)
+
+#define RUN_BINARY_LOOP(ndarray, typecode, type_out, type_left, type_right, lhs, rhs, shape, ndim, operator) do {\
+    uint8_t *left = (uint8_t *)(lhs)->array->items;\
+    uint8_t *right = (uint8_t *)(rhs)->array->items;\
+    (ndarray) = ndarray_new_dense_ndarray((ndim), (shape), (typecode));\
+    uint8_t size_left = sizeof(type_left), size_right = sizeof(type_right);\
+    type_out *out = (type_out *)ndarray->array->items;\
+    for(size_t i=0; i < (lhs)->len; i++) {\
+        out[i] = *(type_left *)left + *(type_right *)right;\
+        left += size_left; right += size_right;\
+    }\
+} while(0)
+
+mp_obj_t mp_obj_new_ndarray_iterator(mp_obj_t , size_t , mp_obj_iter_buf_t *);
 void ndarray_print(const mp_print_t *, mp_obj_t , mp_print_kind_t );
-void ndarray_assign_elements(mp_obj_array_t *, mp_obj_t , uint8_t , size_t *);
-ndarray_obj_t *create_new_ndarray(size_t , size_t , uint8_t );
+ndarray_obj_t *ndarray_new_ndarray(uint8_t , size_t *, int32_t *, uint8_t );
+ndarray_obj_t *ndarray_new_dense_ndarray(uint8_t , size_t *, uint8_t );
+ndarray_obj_t *ndarray_new_linear_array(size_t , uint8_t );
+ndarray_obj_t *ndarray_copy_view(ndarray_obj_t *, uint8_t );
 
 mp_obj_t ndarray_copy(mp_obj_t );
-#ifdef CIRCUITPY
-mp_obj_t ndarray_make_new(const mp_obj_type_t *type, size_t n_args, const mp_obj_t *args, mp_map_t *kw_args);
-#else
 mp_obj_t ndarray_make_new(const mp_obj_type_t *, size_t , size_t , const mp_obj_t *);
-#endif
 mp_obj_t ndarray_subscr(mp_obj_t , mp_obj_t , mp_obj_t );
 mp_obj_t ndarray_getiter(mp_obj_t , mp_obj_iter_buf_t *);
 mp_obj_t ndarray_binary_op(mp_binary_op_t , mp_obj_t , mp_obj_t );
 mp_obj_t ndarray_unary_op(mp_unary_op_t , mp_obj_t );
 
 mp_obj_t ndarray_shape(mp_obj_t );
-mp_obj_t ndarray_size(mp_obj_t );
-mp_obj_t ndarray_itemsize(mp_obj_t );
-mp_obj_t ndarray_flatten(size_t , const mp_obj_t *, mp_map_t *);
-
 mp_obj_t ndarray_reshape(mp_obj_t , mp_obj_t );
-MP_DECLARE_CONST_FUN_OBJ_2(ndarray_reshape_obj);
-
 mp_obj_t ndarray_transpose(mp_obj_t );
-MP_DECLARE_CONST_FUN_OBJ_1(ndarray_transpose_obj);
-
-mp_int_t ndarray_get_buffer(mp_obj_t obj, mp_buffer_info_t *bufinfo, mp_uint_t flags);
-//void ndarray_attributes(mp_obj_t , qstr , mp_obj_t *);
-
-
-#define CREATE_SINGLE_ITEM(outarray, type, typecode, value) do {\
-    ndarray_obj_t *tmp = create_new_ndarray(1, 1, (typecode));\
-    type *tmparr = (type *)tmp->array->items;\
-    tmparr[0] = (type)(value);\
-    (outarray) = MP_OBJ_FROM_PTR(tmp);\
-} while(0)
-
-/*  
-    mp_obj_t row = mp_obj_new_list(n, NULL);
-    mp_obj_list_t *row_ptr = MP_OBJ_TO_PTR(row);
-    
-    should work outside the loop, but it doesn't. Go figure! 
-*/
-
-#define RUN_BINARY_LOOP(typecode, type_out, type_left, type_right, ol, or, op) do {\
-    type_left *left = (type_left *)(ol)->array->items;\
-    type_right *right = (type_right *)(or)->array->items;\
-    uint8_t inc = 0;\
-    if((or)->array->len > 1) inc = 1;\
-    if(((op) == MP_BINARY_OP_ADD) || ((op) == MP_BINARY_OP_SUBTRACT) || ((op) == MP_BINARY_OP_MULTIPLY)) {\
-        ndarray_obj_t *out = create_new_ndarray(ol->m, ol->n, typecode);\
-        type_out *(odata) = (type_out *)out->array->items;\
-        if((op) == MP_BINARY_OP_ADD) { for(size_t i=0, j=0; i < (ol)->array->len; i++, j+=inc) odata[i] = left[i] + right[j];}\
-        if((op) == MP_BINARY_OP_SUBTRACT) { for(size_t i=0, j=0; i < (ol)->array->len; i++, j+=inc) odata[i] = left[i] - right[j];}\
-        if((op) == MP_BINARY_OP_MULTIPLY) { for(size_t i=0, j=0; i < (ol)->array->len; i++, j+=inc) odata[i] = left[i] * right[j];}\
-        return MP_OBJ_FROM_PTR(out);\
-    } else if((op) == MP_BINARY_OP_TRUE_DIVIDE) {\
-        ndarray_obj_t *out = create_new_ndarray(ol->m, ol->n, NDARRAY_FLOAT);\
-        mp_float_t *odata = (mp_float_t *)out->array->items;\
-        for(size_t i=0, j=0; i < (ol)->array->len; i++, j+=inc) odata[i] = (mp_float_t)left[i]/(mp_float_t)right[j];\
-        return MP_OBJ_FROM_PTR(out);\
-    } else if(((op) == MP_BINARY_OP_LESS) || ((op) == MP_BINARY_OP_LESS_EQUAL) ||  \
-             ((op) == MP_BINARY_OP_MORE) || ((op) == MP_BINARY_OP_MORE_EQUAL)) {\
-        mp_obj_t out_list = mp_obj_new_list(0, NULL);\
-        size_t m = (ol)->m, n = (ol)->n;\
-        for(size_t i=0, r=0; i < m; i++, r+=inc) {\
-            mp_obj_t row = mp_obj_new_list(n, NULL);\
-            mp_obj_list_t *row_ptr = MP_OBJ_TO_PTR(row);\
-            for(size_t j=0, s=0; j < n; j++, s+=inc) {\
-                row_ptr->items[j] = mp_const_false;\
-                if((op) == MP_BINARY_OP_LESS) {\
-                    if(left[i*n+j] < right[r*n+s]) row_ptr->items[j] = mp_const_true;\
-                } else if((op) == MP_BINARY_OP_LESS_EQUAL) {\
-                    if(left[i*n+j] <= right[r*n+s]) row_ptr->items[j] = mp_const_true;\
-                } else if((op) == MP_BINARY_OP_MORE) {\
-                    if(left[i*n+j] > right[r*n+s]) row_ptr->items[j] = mp_const_true;\
-                } else if((op) == MP_BINARY_OP_MORE_EQUAL) {\
-                    if(left[i*n+j] >= right[r*n+s]) row_ptr->items[j] = mp_const_true;\
-                }\
-            }\
-            if(m == 1) return row;\
-            mp_obj_list_append(out_list, row);\
-        }\
-        return out_list;\
-    }\
-} while(0)
+mp_obj_t ndarray_flatten(size_t , const mp_obj_t *, mp_map_t *);
+mp_obj_t ndarray_itemsize(mp_obj_t );
+mp_obj_t ndarray_strides(mp_obj_t );
+mp_obj_t ndarray_info(mp_obj_t );
 
 #endif

code/ndarray_properties.h

@@ -1,62 +0,0 @@
-
-/*
- * This file is part of the micropython-ulab project, 
- *
- * https://github.com/v923z/micropython-ulab
- *
- * The MIT License (MIT)
- *
- * Copyright (c) 2020 Jeff Epler for Adafruit Industries
- *               2020 Zoltán Vörös   
-*/
-
-#ifndef _NDARRAY_PROPERTIES_
-#define _NDARRAY_PROPERTIES_
-
-#include "py/runtime.h"
-#include "py/binary.h"
-#include "py/obj.h"
-#include "py/objarray.h"
-
-#include "ndarray.h"
-
-typedef struct _mp_obj_property_t {
-    mp_obj_base_t base;
-    mp_obj_t proxy[3]; // getter, setter, deleter
-} mp_obj_property_t;
-
-/* v923z: it is not at all clear to me, why this must be declared; it should already be in obj.h */
-typedef struct _mp_obj_none_t {
-    mp_obj_base_t base;
-} mp_obj_none_t;
-
-const mp_obj_type_t mp_type_NoneType;
-const mp_obj_none_t mp_const_none_obj = {{&mp_type_NoneType}};
-
-MP_DEFINE_CONST_FUN_OBJ_1(ndarray_get_shape_obj, ndarray_shape);
-MP_DEFINE_CONST_FUN_OBJ_1(ndarray_get_size_obj, ndarray_size);
-MP_DEFINE_CONST_FUN_OBJ_1(ndarray_get_itemsize_obj, ndarray_itemsize);
-MP_DEFINE_CONST_FUN_OBJ_KW(ndarray_flatten_obj, 1, ndarray_flatten);
-
-STATIC const mp_obj_property_t ndarray_shape_obj = {
-    .base.type = &mp_type_property,
-    .proxy = {(mp_obj_t)&ndarray_get_shape_obj,
-              (mp_obj_t)&mp_const_none_obj,
-              (mp_obj_t)&mp_const_none_obj},
-};
-
-STATIC const mp_obj_property_t ndarray_size_obj = {
-    .base.type = &mp_type_property,
-    .proxy = {(mp_obj_t)&ndarray_get_size_obj,
-              (mp_obj_t)&mp_const_none_obj,
-              (mp_obj_t)&mp_const_none_obj},
-};
-
-STATIC const mp_obj_property_t ndarray_itemsize_obj = {
-    .base.type = &mp_type_property,
-    .proxy = {(mp_obj_t)&ndarray_get_itemsize_obj,
-              (mp_obj_t)&mp_const_none_obj,
-              (mp_obj_t)&mp_const_none_obj},
-};
-
-#endif

code/numerical.c

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,7 +5,7 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #include <math.h>
@@ -19,8 +18,6 @@
 #include "py/misc.h"
 #include "numerical.h"
 
-#if ULAB_NUMERICAL_MODULE
-
 enum NUMERICAL_FUNCTION_TYPE {
     NUMERICAL_MIN,
     NUMERICAL_MAX,
@@ -31,13 +28,39 @@ enum NUMERICAL_FUNCTION_TYPE {
     NUMERICAL_STD,
 };
 
+// creates the shape and strides arrays of the contracted ndarray
+ndarray_header_obj_t contracted_shape_strides(ndarray_obj_t *ndarray, int8_t axis) {
+    if(axis < 0) axis += ndarray->ndim;
+    if((axis > ndarray->ndim-1) || (axis < 0)) {
+        mp_raise_ValueError("tuple index out of range");
+    }
+    size_t *shape = m_new(size_t, ndarray->ndim-1);
+    int32_t *strides = m_new(int32_t, ndarray->ndim-1);
+    for(size_t i=0, j=0; i < ndarray->ndim; i++) {
+        if(axis != i) {
+            shape[j] = ndarray->shape[j];
+            j++;
+        }
+    }
+    int32_t stride = 1;
+    for(size_t i=0; i < ndarray->ndim-1; i++) {
+        strides[ndarray->ndim-2-i] = stride;
+        stride *= shape[ndarray->ndim-2-i];
+    }
+    ndarray_header_obj_t header;
+    header.shape = shape;
+    header.strides = strides;
+    header.axis = axis;
+    return header;
+}
+
 mp_obj_t numerical_linspace(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
     static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
+        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = MP_ROM_PTR(&mp_const_none_obj) } },
+        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = MP_ROM_PTR(&mp_const_none_obj) } },
         { MP_QSTR_num, MP_ARG_INT, {.u_int = 50} },
-        { MP_QSTR_endpoint, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_rom_obj = mp_const_true} },
-        { MP_QSTR_retstep, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_rom_obj = mp_const_false} },
+        { MP_QSTR_endpoint, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_rom_obj = MP_ROM_PTR(&mp_const_true_obj)} },
+        { MP_QSTR_retstep, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_rom_obj = MP_ROM_PTR(&mp_const_false_obj)} },
         { MP_QSTR_dtype, MP_ARG_KW_ONLY | MP_ARG_INT, {.u_int = NDARRAY_FLOAT} },
     };
 
@@ -46,14 +69,14 @@ mp_obj_t numerical_linspace(size_t n_args, const mp_obj_t *pos_args, mp_map_t *k
 
     uint16_t len = args[2].u_int;
     if(len < 2) {
-        mp_raise_ValueError(translate("number of points must be at least 2"));
+        mp_raise_ValueError("number of points must be at least 2");
     }
     mp_float_t value, step;
     value = mp_obj_get_float(args[0].u_obj);
     uint8_t typecode = args[5].u_int;
     if(args[3].u_obj == mp_const_true) step = (mp_obj_get_float(args[1].u_obj)-value)/(len-1);
     else step = (mp_obj_get_float(args[1].u_obj)-value)/len;
-    ndarray_obj_t *ndarray = create_new_ndarray(1, len, typecode);
+    ndarray_obj_t *ndarray = ndarray_new_linear_array(len, typecode);
     if(typecode == NDARRAY_UINT8) {
         uint8_t *array = (uint8_t *)ndarray->array->items;
         for(size_t i=0; i < len; i++, value += step) array[i] = (uint8_t)value;
@@ -80,126 +103,71 @@ mp_obj_t numerical_linspace(size_t n_args, const mp_obj_t *pos_args, mp_map_t *k
     }
 }
 
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_linspace_obj, 2, numerical_linspace);
-
-void axis_sorter(ndarray_obj_t *ndarray, mp_obj_t axis, size_t *m, size_t *n, size_t *N, 
-                 size_t *increment, size_t *len, size_t *start_inc) {
-    if(axis == mp_const_none) { // flatten the array
-        *m = 1;
-        *n = 1;
-        *len = ndarray->array->len;
-        *N = 1;
-        *increment = 1;
-        *start_inc = ndarray->array->len;
-    } else if((mp_obj_get_int(axis) == 1)) { // along the horizontal axis
-        *m = ndarray->m;
-        *n = 1;
-        *len = ndarray->n;
-        *N = ndarray->m;
-        *increment = 1;
-        *start_inc = ndarray->n;
-    } else { // along vertical axis
-        *m = 1;
-        *n = ndarray->n;
-        *len = ndarray->m;
-        *N = ndarray->n;
-        *increment = ndarray->n;
-        *start_inc = 1;
+mp_obj_t numerical_flat_sum_mean_std(ndarray_obj_t *ndarray, uint8_t optype, size_t ddof) {
+    mp_float_t value;
+    int32_t *shape_strides = m_new(int32_t, ndarray->ndim);
+    shape_strides[ndarray->ndim-1] = ndarray->strides[ndarray->ndim-1];
+    for(uint8_t i=ndarray->ndim-1; i > 0; i--) {
+        shape_strides[i-1] = shape_strides[i] * ndarray->shape[i-1];
     }
+    if(ndarray->array->typecode == NDARRAY_UINT8) {
+        CALCULATE_FLAT_SUM_STD(ndarray, uint8_t, value, shape_strides, optype);
+    } else if(ndarray->array->typecode == NDARRAY_INT8) {
+        CALCULATE_FLAT_SUM_STD(ndarray, int8_t, value, shape_strides, optype);
+    } if(ndarray->array->typecode == NDARRAY_UINT16) {
+        CALCULATE_FLAT_SUM_STD(ndarray, uint16_t, value, shape_strides, optype);
+    } else if(ndarray->array->typecode == NDARRAY_INT16) {
+        CALCULATE_FLAT_SUM_STD(ndarray, int16_t, value, shape_strides, optype);
+    } else {
+        CALCULATE_FLAT_SUM_STD(ndarray, mp_float_t, value, shape_strides, optype);
     }
-
-mp_obj_t numerical_sum_mean_std_iterable(mp_obj_t oin, uint8_t optype, size_t ddof) {
-    mp_float_t value, sum = 0.0, sq_sum = 0.0;
-    mp_obj_iter_buf_t iter_buf;
-    mp_obj_t item, iterable = mp_getiter(oin, &iter_buf);
-    mp_int_t len = mp_obj_get_int(mp_obj_len(oin));
-    while ((item = mp_iternext(iterable)) != MP_OBJ_STOP_ITERATION) {
-        value = mp_obj_get_float(item);
-        sum += value;
-    }
+    m_del(int32_t, shape_strides, ndarray->ndim);
     if(optype == NUMERICAL_SUM) {
-        return mp_obj_new_float(sum);
+        return mp_obj_new_float(value);
     } else if(optype == NUMERICAL_MEAN) {
-        return mp_obj_new_float(sum/len);
-    } else { // this should be the case of the standard deviation
-        // TODO: note that we could get away with a single pass, if we used the Weldorf algorithm
-        // That should save a fair amount of time, because we would have to extract the values only once
-        iterable = mp_getiter(oin, &iter_buf);
-        sum /= len; // this is now the mean!
-        while ((item = mp_iternext(iterable)) != MP_OBJ_STOP_ITERATION) {
-            value = mp_obj_get_float(item) - sum;
-            sq_sum += value * value;
-        }
-        return mp_obj_new_float(MICROPY_FLOAT_C_FUN(sqrt)(sq_sum/(len-ddof)));
+        return mp_obj_new_float(value/ndarray->len);
+    } else {
+        return mp_obj_new_float(MICROPY_FLOAT_C_FUN(sqrt)(value/(ndarray->len-ddof)));
     }
 }   
 
-STATIC mp_obj_t numerical_sum_mean_ndarray(ndarray_obj_t *ndarray, mp_obj_t axis, uint8_t optype) {
-    size_t m, n, increment, start, start_inc, N, len; 
-    axis_sorter(ndarray, axis, &m, &n, &N, &increment, &len, &start_inc);
-    ndarray_obj_t *results = create_new_ndarray(m, n, NDARRAY_FLOAT);
-    mp_float_t sum, sq_sum;
-    mp_float_t *farray = (mp_float_t *)results->array->items;
-    for(size_t j=0; j < N; j++) { // result index
-        start = j * start_inc;
-        sum = sq_sum = 0.0;
+// numerical functions for ndarrays
+mp_obj_t numerical_sum_mean_ndarray(ndarray_obj_t *ndarray, mp_obj_t axis, uint8_t optype) {
+    if(axis == mp_const_none) {
+        return numerical_flat_sum_mean_std(ndarray, optype, 0);
+    } else {
+        int8_t ax = mp_obj_get_int(axis);
+        ndarray_header_obj_t header = contracted_shape_strides(ndarray, ax);
+        ndarray_obj_t *result = ndarray_new_ndarray(ndarray->ndim-1, header.shape, header.strides, NDARRAY_FLOAT);
+        mp_float_t *farray = (mp_float_t *)result->array->items;
+        size_t offset;
+        // iterate along the length of the output array, so as to avoid recursion
+        for(size_t i=0; i < result->array->len; i++) {
+            offset = ndarray_index_from_contracted(i, ndarray, result->strides, result->ndim, header.axis) + ndarray->offset;
             if(ndarray->array->typecode == NDARRAY_UINT8) {
-            RUN_SUM(ndarray, uint8_t, optype, len, start, increment);
+                CALCULATE_SUM(ndarray, uint8_t, farray, ndarray->shape[header.axis], i, ndarray->strides[header.axis], offset, optype);
             } else if(ndarray->array->typecode == NDARRAY_INT8) {
-            RUN_SUM(ndarray, int8_t, optype, len, start, increment);
+                CALCULATE_SUM(ndarray, int8_t, farray, ndarray->shape[header.axis], i, ndarray->strides[header.axis], offset, optype); 
             } else if(ndarray->array->typecode == NDARRAY_UINT16) {
-            RUN_SUM(ndarray, uint16_t, optype, len, start, increment);
+                CALCULATE_SUM(ndarray, uint16_t, farray, ndarray->shape[header.axis], i, ndarray->strides[header.axis], offset, optype);     
             } else if(ndarray->array->typecode == NDARRAY_INT16) {
-            RUN_SUM(ndarray, int16_t, optype, len, start, increment);
-        } else { // this will be mp_float_t, no need to check
-            RUN_SUM(ndarray, mp_float_t, optype, len, start, increment);
+                CALCULATE_SUM(ndarray, int16_t, farray, ndarray->shape[header.axis], i, ndarray->strides[header.axis], offset, optype);
+            } else {
+                CALCULATE_SUM(ndarray, mp_float_t, farray, ndarray->shape[header.axis], i, ndarray->strides[header.axis], offset, optype);       
             }
-        if(optype == NUMERICAL_SUM) {
-            farray[j] = sum;
-        } else { // this is the case of the mean
-            farray[j] = sum / len;
+            if(optype == NUMERICAL_MEAN) farray[i] /= ndarray->shape[header.axis];
         }
+        return MP_OBJ_FROM_PTR(result);
     }
-    if(results->array->len == 1) {
-        return mp_obj_new_float(farray[0]);
-    }
-    return MP_OBJ_FROM_PTR(results);
+    return mp_const_none;
 }
 
-mp_obj_t numerical_std_ndarray(ndarray_obj_t *ndarray, mp_obj_t axis, size_t ddof) {
-    size_t m, n, increment, start, start_inc, N, len; 
-    mp_float_t sum, sum_sq;
-    
-    axis_sorter(ndarray, axis, &m, &n, &N, &increment, &len, &start_inc);
-    if(ddof > len) {
-        mp_raise_ValueError(translate("ddof must be smaller than length of data set"));
-    }
-    ndarray_obj_t *results = create_new_ndarray(m, n, NDARRAY_FLOAT);
-    mp_float_t *farray = (mp_float_t *)results->array->items;
-    for(size_t j=0; j < N; j++) { // result index
-        start = j * start_inc;
-        sum = 0.0;
-        sum_sq = 0.0;
-        if(ndarray->array->typecode == NDARRAY_UINT8) {
-            RUN_STD(ndarray, uint8_t, len, start, increment);
-        } else if(ndarray->array->typecode == NDARRAY_INT8) {
-            RUN_STD(ndarray, int8_t, len, start, increment);
-        } else if(ndarray->array->typecode == NDARRAY_UINT16) {
-            RUN_STD(ndarray, uint16_t, len, start, increment);
-        } else if(ndarray->array->typecode == NDARRAY_INT16) {
-            RUN_STD(ndarray, int16_t, len, start, increment);
-        } else { // this will be mp_float_t, no need to check
-            RUN_STD(ndarray, mp_float_t, len, start, increment);
-        }
-        farray[j] = MICROPY_FLOAT_C_FUN(sqrt)(sum_sq/(len - ddof));
-    }
-    if(results->array->len == 1) {
-        return mp_obj_new_float(farray[0]);
-    }
-    return MP_OBJ_FROM_PTR(results);
+mp_obj_t numerical_argmin_argmax_ndarray(ndarray_obj_t *ndarray, mp_obj_t axis, uint8_t optype) {
+    return mp_const_none;
 }
 
-mp_obj_t numerical_argmin_argmax_iterable(mp_obj_t oin, mp_obj_t axis, uint8_t optype) {
+// numerical function for interables (single axis)
+mp_obj_t numerical_argmin_argmax_iterable(mp_obj_t oin, uint8_t optype) {
     size_t idx = 0, best_idx = 0;
     mp_obj_iter_buf_t iter_buf;
     mp_obj_t iterable = mp_getiter(oin, &iter_buf);
@@ -221,44 +189,34 @@ mp_obj_t numerical_argmin_argmax_iterable(mp_obj_t oin, mp_obj_t axis, uint8_t o
     }    
 }
 
-mp_obj_t numerical_argmin_argmax_ndarray(ndarray_obj_t *ndarray, mp_obj_t axis, uint8_t optype) {
-    size_t m, n, increment, start, start_inc, N, len;
-    axis_sorter(ndarray, axis, &m, &n, &N, &increment, &len, &start_inc);
-    ndarray_obj_t *results;
-    if((optype == NUMERICAL_ARGMIN) || (optype == NUMERICAL_ARGMAX)) {
-        // we could save some RAM by taking NDARRAY_UINT8, if the dimensions 
-        // are smaller than 256, but the code would become more involving 
-        // (we would also need extra flash space)
-        results = create_new_ndarray(m, n, NDARRAY_UINT16);
-    } else {
-        results = create_new_ndarray(m, n, ndarray->array->typecode);
+mp_obj_t numerical_sum_mean_std_iterable(mp_obj_t oin, uint8_t optype, size_t ddof) {
+    mp_float_t value, sum = 0.0, sq_sum = 0.0;
+    mp_obj_iter_buf_t iter_buf;
+    mp_obj_t item, iterable = mp_getiter(oin, &iter_buf);
+    mp_int_t len = mp_obj_get_int(mp_obj_len(oin));
+    while ((item = mp_iternext(iterable)) != MP_OBJ_STOP_ITERATION) {
+        value = mp_obj_get_float(item);
+        sum += value;
     }
-    
-    for(size_t j=0; j < N; j++) { // result index
-        start = j * start_inc;
-        if((ndarray->array->typecode == NDARRAY_UINT8) || (ndarray->array->typecode == NDARRAY_INT8)) {
-            if((optype == NUMERICAL_MAX) || (optype == NUMERICAL_MIN)) {
-                RUN_ARGMIN(ndarray, results, uint8_t, uint8_t, len, start, increment, optype, j);
-            } else {
-                RUN_ARGMIN(ndarray, results, uint8_t, uint16_t, len, start, increment, optype, j);                
+    if(optype ==  NUMERICAL_SUM) {
+        return mp_obj_new_float(sum);
+    } else if(optype == NUMERICAL_MEAN) {
+        return mp_obj_new_float(sum/len);
+    } else { // this should be the case of the standard deviation
+        sum /= len; // this is the mean now
+        iterable = mp_getiter(oin, &iter_buf);
+        while ((item = mp_iternext(iterable)) != MP_OBJ_STOP_ITERATION) {
+            value = mp_obj_get_float(item) - sum;
+            sq_sum += value * value;
         }
-        } else if((ndarray->array->typecode == NDARRAY_UINT16) || (ndarray->array->typecode == NDARRAY_INT16)) {
-            RUN_ARGMIN(ndarray, results, uint16_t, uint16_t, len, start, increment, optype, j);
-        } else {
-            if((optype == NUMERICAL_MAX) || (optype == NUMERICAL_MIN)) {
-                RUN_ARGMIN(ndarray, results, mp_float_t, mp_float_t, len, start, increment, optype, j);
-            } else {
-                RUN_ARGMIN(ndarray, results, mp_float_t, uint16_t, len, start, increment, optype, j);                
+        return mp_obj_new_float(MICROPY_FLOAT_C_FUN(sqrt)(sq_sum/(len-ddof)));
     }
 }
-    }
-    return MP_OBJ_FROM_PTR(results);
-}
 
 STATIC mp_obj_t numerical_function(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args, uint8_t optype) {
     static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none} } ,
-        { MP_QSTR_axis, MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
+        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = MP_ROM_PTR(&mp_const_none_obj)} } ,
+        { MP_QSTR_axis, MP_ARG_OBJ, {.u_rom_obj = MP_ROM_PTR(&mp_const_none_obj)} },
     };
 
     mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)];
@@ -266,10 +224,6 @@ STATIC mp_obj_t numerical_function(size_t n_args, const mp_obj_t *pos_args, mp_m
     
     mp_obj_t oin = args[0].u_obj;
     mp_obj_t axis = args[1].u_obj;
-    if((axis != mp_const_none) && (mp_obj_get_int(axis) != 0) && (mp_obj_get_int(axis) != 1)) {
-        // this seems to pass with False, and True...
-        mp_raise_ValueError(translate("axis must be None, 0, or 1"));
-    }
     
     if(MP_OBJ_IS_TYPE(oin, &mp_type_tuple) || MP_OBJ_IS_TYPE(oin, &mp_type_list) || 
         MP_OBJ_IS_TYPE(oin, &mp_type_range)) {
@@ -278,73 +232,45 @@ STATIC mp_obj_t numerical_function(size_t n_args, const mp_obj_t *pos_args, mp_m
             case NUMERICAL_ARGMIN:
             case NUMERICAL_MAX:
             case NUMERICAL_ARGMAX:
-                return numerical_argmin_argmax_iterable(oin, axis, optype);
+                return numerical_argmin_argmax_iterable(oin, optype);
             case NUMERICAL_SUM:
             case NUMERICAL_MEAN:
                 return numerical_sum_mean_std_iterable(oin, optype, 0);
-            default: // we should never reach this point, but whatever
+            default: // we should never end up here
                 return mp_const_none;
         }
     } else if(MP_OBJ_IS_TYPE(oin, &ulab_ndarray_type)) {
         ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(oin);
         switch(optype) {
             case NUMERICAL_MIN:
-            case NUMERICAL_MAX:
             case NUMERICAL_ARGMIN:
+            case NUMERICAL_MAX:
             case NUMERICAL_ARGMAX:
                 return numerical_argmin_argmax_ndarray(ndarray, axis, optype);
             case NUMERICAL_SUM:
             case NUMERICAL_MEAN:
-                return numerical_sum_mean_ndarray(ndarray, axis, optype);
+                return numerical_sum_mean_ndarray(ndarray, args[1].u_obj, optype);
             default:
-                mp_raise_NotImplementedError(translate("operation is not implemented on ndarrays"));
+                return mp_const_none;
         }
     } else {
-        mp_raise_TypeError(translate("input must be tuple, list, range, or ndarray"));
+        mp_raise_TypeError("input must be tuple, list, range, or ndarray");
     }
     return mp_const_none;
 }
 
-mp_obj_t numerical_min(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    return numerical_function(n_args, pos_args, kw_args, NUMERICAL_MIN);
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_min_obj, 1, numerical_min);
-
-mp_obj_t numerical_max(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    return numerical_function(n_args, pos_args, kw_args, NUMERICAL_MAX);
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_max_obj, 1, numerical_max);
-
-mp_obj_t numerical_argmin(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    return numerical_function(n_args, pos_args, kw_args, NUMERICAL_ARGMIN);
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_argmin_obj, 1, numerical_argmin);
-
-mp_obj_t numerical_argmax(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    return numerical_function(n_args, pos_args, kw_args, NUMERICAL_ARGMAX);
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_argmax_obj, 1, numerical_argmax);
-
 mp_obj_t numerical_sum(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
     return numerical_function(n_args, pos_args, kw_args, NUMERICAL_SUM);
 }
 
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_sum_obj, 1, numerical_sum);
-
 mp_obj_t numerical_mean(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
     return numerical_function(n_args, pos_args, kw_args, NUMERICAL_MEAN);
 }
 
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_mean_obj, 1, numerical_mean);
-
 mp_obj_t numerical_std(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
     static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } } ,
-        { MP_QSTR_axis, MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
+        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = MP_ROM_PTR(&mp_const_none_obj)} } ,
+        { MP_QSTR_axis, MP_ARG_OBJ, {.u_rom_obj = MP_ROM_PTR(&mp_const_none_obj)} },
         { MP_QSTR_ddof, MP_ARG_KW_ONLY | MP_ARG_INT, {.u_int = 0} },
     };
 
@@ -354,405 +280,37 @@ mp_obj_t numerical_std(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_arg
     mp_obj_t oin = args[0].u_obj;
     mp_obj_t axis = args[1].u_obj;
     size_t ddof = args[2].u_int;
-    if((axis != mp_const_none) && (mp_obj_get_int(axis) != 0) && (mp_obj_get_int(axis) != 1)) {
-        // this seems to pass with False, and True...
-        mp_raise_ValueError(translate("axis must be None, 0, or 1"));
-    }
     if(MP_OBJ_IS_TYPE(oin, &mp_type_tuple) || MP_OBJ_IS_TYPE(oin, &mp_type_list) || MP_OBJ_IS_TYPE(oin, &mp_type_range)) {
         return numerical_sum_mean_std_iterable(oin, NUMERICAL_STD, ddof);
     } else if(MP_OBJ_IS_TYPE(oin, &ulab_ndarray_type)) {
         ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(oin);
-        return numerical_std_ndarray(ndarray, axis, ddof);
+        if(axis == mp_const_none) { // calculate for the flat array
+            return numerical_flat_sum_mean_std(ndarray, NUMERICAL_STD, ddof);
         } else {
-        mp_raise_TypeError(translate("input must be tuple, list, range, or ndarray"));
+            int8_t ax = mp_obj_get_int(axis);
+            ndarray_header_obj_t header = contracted_shape_strides(ndarray, ax);
+            ndarray_obj_t *result = ndarray_new_ndarray(ndarray->ndim-1, header.shape, header.strides, NDARRAY_FLOAT);
+            mp_float_t *farray = (mp_float_t *)result->array->items, sum_sq;
+            size_t offset;
+            // iterate along the length of the output array, so as to avoid recursion
+            for(size_t i=0; i < result->array->len; i++) {
+                offset = ndarray_index_from_contracted(i, ndarray, result->strides, result->ndim, header.axis) + ndarray->offset;
+                if(ndarray->array->typecode == NDARRAY_UINT8) {
+                    CALCULATE_STD(ndarray, uint8_t, sum_sq, ndarray->shape[header.axis], ndarray->strides[header.axis], offset);
+                } else if(ndarray->array->typecode == NDARRAY_INT8) {
+                    CALCULATE_STD(ndarray, int8_t, sum_sq, ndarray->shape[header.axis], ndarray->strides[header.axis], offset); 
+                } else if(ndarray->array->typecode == NDARRAY_UINT16) {
+                    CALCULATE_STD(ndarray, uint16_t, sum_sq, ndarray->shape[header.axis], ndarray->strides[header.axis], offset);     
+                } else if(ndarray->array->typecode == NDARRAY_INT16) {
+                    CALCULATE_STD(ndarray, int16_t, sum_sq, ndarray->shape[header.axis], ndarray->strides[header.axis], offset);
+                } else {
+                    CALCULATE_STD(ndarray, mp_float_t, sum_sq, ndarray->shape[header.axis], ndarray->strides[header.axis], offset);
+                }
+                farray[i] = MICROPY_FLOAT_C_FUN(sqrt)(sum_sq / (ndarray->shape[header.axis] - ddof));
+            }
+            return MP_OBJ_FROM_PTR(result);
+        }
+        mp_raise_TypeError("input must be tuple, list, range, or ndarray");
     }
     return mp_const_none;
 }
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_std_obj, 1, numerical_std);
-
-mp_obj_t numerical_roll(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none  } },
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-        { MP_QSTR_axis, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-    };
-
-    mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)];
-    mp_arg_parse_all(2, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args);
-    
-    mp_obj_t oin = args[0].u_obj;
-    int16_t shift = mp_obj_get_int(args[1].u_obj);
-    if((args[2].u_obj != mp_const_none) && 
-           (mp_obj_get_int(args[2].u_obj) != 0) && 
-           (mp_obj_get_int(args[2].u_obj) != 1)) {
-        mp_raise_ValueError(translate("axis must be None, 0, or 1"));
-    }
-
-    ndarray_obj_t *in = MP_OBJ_TO_PTR(oin);
-    uint8_t _sizeof = mp_binary_get_size('@', in->array->typecode, NULL);
-    size_t len;
-    int16_t _shift;
-    uint8_t *array = (uint8_t *)in->array->items;
-    // TODO: transpose the matrix, if axis == 0. Though, that is hard on the RAM...
-    if(shift < 0) {
-        _shift = -shift;
-    } else {
-        _shift = shift;
-    }
-    if((args[2].u_obj == mp_const_none) || (mp_obj_get_int(args[2].u_obj) == 1)) { // shift horizontally
-        uint16_t M;
-        if(args[2].u_obj == mp_const_none) {
-            len = in->array->len;
-            M = 1;
-        } else {
-            len = in->n;
-            M = in->m;
-        }
-        _shift = _shift % len;
-        if(shift < 0) _shift = len - _shift;
-        // TODO: if(shift > len/2), we should move in the opposite direction. That would save RAM
-        _shift *= _sizeof;
-        uint8_t *tmp = m_new(uint8_t, _shift);
-        for(size_t m=0; m < M; m++) {
-            memmove(tmp, &array[m*len*_sizeof], _shift);
-            memmove(&array[m*len*_sizeof], &array[m*len*_sizeof+_shift], len*_sizeof-_shift);
-            memmove(&array[(m+1)*len*_sizeof-_shift], tmp, _shift);
-        }
-        m_del(uint8_t, tmp, _shift);
-        return mp_const_none;
-    } else {
-        len = in->m;
-        // temporary buffer
-        uint8_t *_data = m_new(uint8_t, _sizeof*len);
-        
-        _shift = _shift % len;
-        if(shift < 0) _shift = len - _shift;
-        _shift *= _sizeof;
-        uint8_t *tmp = m_new(uint8_t, _shift);
-
-        for(size_t n=0; n < in->n; n++) {
-            for(size_t m=0; m < len; m++) {
-                // this loop should fill up the temporary buffer
-                memmove(&_data[m*_sizeof], &array[(m*in->n+n)*_sizeof], _sizeof);
-            }
-            // now, the actual shift
-            memmove(tmp, _data, _shift);
-            memmove(_data, &_data[_shift], len*_sizeof-_shift);
-            memmove(&_data[len*_sizeof-_shift], tmp, _shift);
-            for(size_t m=0; m < len; m++) {
-                // this loop should dump the content of the temporary buffer into data
-                memmove(&array[(m*in->n+n)*_sizeof], &_data[m*_sizeof], _sizeof);
-            }            
-        }
-        m_del(uint8_t, tmp, _shift);
-        m_del(uint8_t, _data, _sizeof*len);
-        return mp_const_none;
-    }
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_roll_obj, 2, numerical_roll);
-
-mp_obj_t numerical_flip(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-        { MP_QSTR_axis, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-    };
-
-    mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)];
-    mp_arg_parse_all(1, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args);
-    
-    if(!MP_OBJ_IS_TYPE(args[0].u_obj, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("flip argument must be an ndarray"));
-    }
-    if((args[1].u_obj != mp_const_none) && 
-           (mp_obj_get_int(args[1].u_obj) != 0) && 
-           (mp_obj_get_int(args[1].u_obj) != 1)) {
-        mp_raise_ValueError(translate("axis must be None, 0, or 1"));
-    }
-
-    ndarray_obj_t *in = MP_OBJ_TO_PTR(args[0].u_obj);
-    mp_obj_t oout = ndarray_copy(args[0].u_obj);
-    ndarray_obj_t *out = MP_OBJ_TO_PTR(oout);
-    uint8_t _sizeof = mp_binary_get_size('@', in->array->typecode, NULL);
-    uint8_t *array_in = (uint8_t *)in->array->items;
-    uint8_t *array_out = (uint8_t *)out->array->items;    
-    size_t len;
-    if((args[1].u_obj == mp_const_none) || (mp_obj_get_int(args[1].u_obj) == 1)) { // flip horizontally
-        uint16_t M = in->m;
-        len = in->n;
-        if(args[1].u_obj == mp_const_none) { // flip flattened array
-            len = in->array->len;
-            M = 1;
-        }
-        for(size_t m=0; m < M; m++) {
-            for(size_t n=0; n < len; n++) {
-                memcpy(array_out+_sizeof*(m*len+n), array_in+_sizeof*((m+1)*len-n-1), _sizeof);
-            }
-        }
-    } else { // flip vertically
-        for(size_t m=0; m < in->m; m++) {
-            for(size_t n=0; n < in->n; n++) {
-                memcpy(array_out+_sizeof*(m*in->n+n), array_in+_sizeof*((in->m-m-1)*in->n+n), _sizeof);
-            }
-        }
-    }
-    return out;
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_flip_obj, 1, numerical_flip);
-
-mp_obj_t numerical_diff(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-        { MP_QSTR_n, MP_ARG_KW_ONLY | MP_ARG_INT, {.u_int = 1 } },
-        { MP_QSTR_axis, MP_ARG_KW_ONLY | MP_ARG_INT, {.u_int = -1 } },
-    };
-
-    mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)];
-    mp_arg_parse_all(1, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args);
-    
-    if(!MP_OBJ_IS_TYPE(args[0].u_obj, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("diff argument must be an ndarray"));
-    }
-    
-    ndarray_obj_t *in = MP_OBJ_TO_PTR(args[0].u_obj);
-    size_t increment, N, M;
-    if((args[2].u_int == -1) || (args[2].u_int == 1)) { // differentiate along the horizontal axis
-        increment = 1;
-    } else if(args[2].u_int == 0) { // differtiate along vertical axis
-        increment = in->n;
-    } else {
-        mp_raise_ValueError(translate("axis must be -1, 0, or 1"));
-    }
-    if((args[1].u_int < 0) || (args[1].u_int > 9)) {
-        mp_raise_ValueError(translate("n must be between 0, and 9"));
-    }
-    uint8_t n = args[1].u_int;
-    int8_t *stencil = m_new(int8_t, n+1);
-    stencil[0] = 1;
-    for(uint8_t i=1; i < n+1; i++) {
-        stencil[i] = -stencil[i-1]*(n-i+1)/i;
-    }
-
-    ndarray_obj_t *out;
-    
-    if(increment == 1) { // differentiate along the horizontal axis 
-        if(n >= in->n) {
-            out = create_new_ndarray(in->m, 0, in->array->typecode);
-            m_del(uint8_t, stencil, n);
-            return MP_OBJ_FROM_PTR(out);
-        }
-        N = in->n - n;
-        M = in->m;
-    } else { // differentiate along vertical axis
-        if(n >= in->m) {
-            out = create_new_ndarray(0, in->n, in->array->typecode);
-            m_del(uint8_t, stencil, n);
-            return MP_OBJ_FROM_PTR(out);
-        }
-        M = in->m - n;
-        N = in->n;
-    }
-    out = create_new_ndarray(M, N, in->array->typecode);
-    if(in->array->typecode == NDARRAY_UINT8) {
-        CALCULATE_DIFF(in, out, uint8_t, M, N, in->n, increment);
-    } else if(in->array->typecode == NDARRAY_INT8) {
-        CALCULATE_DIFF(in, out, int8_t, M, N, in->n, increment);
-    }  else if(in->array->typecode == NDARRAY_UINT16) {
-        CALCULATE_DIFF(in, out, uint16_t, M, N, in->n, increment);
-    } else if(in->array->typecode == NDARRAY_INT16) {
-        CALCULATE_DIFF(in, out, int16_t, M, N, in->n, increment);
-    } else {
-        CALCULATE_DIFF(in, out, mp_float_t, M, N, in->n, increment);
-    }
-    m_del(int8_t, stencil, n);
-    return MP_OBJ_FROM_PTR(out);
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_diff_obj, 1, numerical_diff);
-
-mp_obj_t numerical_sort_helper(mp_obj_t oin, mp_obj_t axis, uint8_t inplace) {
-    if(!MP_OBJ_IS_TYPE(oin, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("sort argument must be an ndarray"));
-    }
-
-    ndarray_obj_t *ndarray;
-    mp_obj_t out;
-    if(inplace == 1) {
-        ndarray = MP_OBJ_TO_PTR(oin);
-    } else {
-        out = ndarray_copy(oin);
-        ndarray = MP_OBJ_TO_PTR(out);
-    }
-    size_t increment, start_inc, end, N;
-    if(axis == mp_const_none) { // flatten the array
-        ndarray->m = 1;
-        ndarray->n = ndarray->array->len;
-        increment = 1;
-        start_inc = ndarray->n;
-        end = ndarray->n;
-        N = ndarray->n;
-    } else if((mp_obj_get_int(axis) == -1) || 
-              (mp_obj_get_int(axis) == 1)) { // sort along the horizontal axis
-        increment = 1;
-        start_inc = ndarray->n;
-        end = ndarray->array->len;
-        N = ndarray->n;
-    } else if(mp_obj_get_int(axis) == 0) { // sort along vertical axis
-        increment = ndarray->n;
-        start_inc = 1;
-        end = ndarray->m;
-        N = ndarray->m;
-    } else {
-        mp_raise_ValueError(translate("axis must be -1, 0, None, or 1"));
-    }
-    
-    size_t q, k, p, c;
-
-    for(size_t start=0; start < end; start+=start_inc) {
-        q = N; 
-        k = (q >> 1);
-        if((ndarray->array->typecode == NDARRAY_UINT8) || (ndarray->array->typecode == NDARRAY_INT8)) {
-            HEAPSORT(uint8_t, ndarray);
-        } else if((ndarray->array->typecode == NDARRAY_INT16) || (ndarray->array->typecode == NDARRAY_INT16)) {
-            HEAPSORT(uint16_t, ndarray);
-        } else {
-            HEAPSORT(mp_float_t, ndarray);
-        }
-    }
-    if(inplace == 1) {
-        return mp_const_none;
-    } else {
-        return out;
-    }
-}
-
-// numpy function
-mp_obj_t numerical_sort(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-        { MP_QSTR_axis, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_int = -1 } },
-    };
-
-    mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)];
-    mp_arg_parse_all(1, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args);
-
-    return numerical_sort_helper(args[0].u_obj, args[1].u_obj, 0);
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_sort_obj, 1, numerical_sort);
-
-// method of an ndarray
-mp_obj_t numerical_sort_inplace(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-        { MP_QSTR_axis, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_int = -1 } },
-    };
-
-    mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)];
-    mp_arg_parse_all(1, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args);
-
-    return numerical_sort_helper(args[0].u_obj, args[1].u_obj, 1);
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_sort_inplace_obj, 1, numerical_sort_inplace);
-
-mp_obj_t numerical_argsort(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) {
-    static const mp_arg_t allowed_args[] = {
-        { MP_QSTR_, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = mp_const_none } },
-        { MP_QSTR_axis, MP_ARG_KW_ONLY | MP_ARG_OBJ, {.u_int = -1 } },
-    };
-    mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)];
-    mp_arg_parse_all(1, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args);
-    if(!MP_OBJ_IS_TYPE(args[0].u_obj, &ulab_ndarray_type)) {
-        mp_raise_TypeError(translate("argsort argument must be an ndarray"));
-    }
-
-    ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(args[0].u_obj);
-    size_t increment, start_inc, end, N, m, n;
-    if(args[1].u_obj == mp_const_none) { // flatten the array
-        m = 1;
-        n = ndarray->array->len;
-        ndarray->m = m;
-        ndarray->n = n;
-        increment = 1;
-        start_inc = ndarray->n;
-        end = ndarray->n;
-        N = n;
-    } else if((mp_obj_get_int(args[1].u_obj) == -1) || 
-              (mp_obj_get_int(args[1].u_obj) == 1)) { // sort along the horizontal axis
-        m = ndarray->m;
-        n = ndarray->n;
-        increment = 1;
-        start_inc = n;
-        end = ndarray->array->len;
-        N = n;
-    } else if(mp_obj_get_int(args[1].u_obj) == 0) { // sort along vertical axis
-        m = ndarray->m;
-        n = ndarray->n;
-        increment = n;
-        start_inc = 1;
-        end = m;
-        N = m;
-    } else {
-        mp_raise_ValueError(translate("axis must be -1, 0, None, or 1"));
-    }
-
-    // at the expense of flash, we could save RAM by creating 
-    // an NDARRAY_UINT16 ndarray only, if needed, otherwise, NDARRAY_UINT8
-    ndarray_obj_t *indices = create_new_ndarray(m, n, NDARRAY_UINT16);
-    uint16_t *index_array = (uint16_t *)indices->array->items;
-    // initialise the index array
-    // if array is flat: 0 to indices->n
-    // if sorting vertically, identical indices are arranged row-wise
-    // if sorting horizontally, identical indices are arranged colunn-wise
-    for(uint16_t start=0; start < end; start+=start_inc) {
-        for(uint16_t s=0; s < N; s++) {
-            index_array[start+s*increment] = s;
-        }
-    }
-
-    size_t q, k, p, c;
-    for(size_t start=0; start < end; start+=start_inc) {
-        q = N; 
-        k = (q >> 1);
-        if((ndarray->array->typecode == NDARRAY_UINT8) || (ndarray->array->typecode == NDARRAY_INT8)) {
-            HEAP_ARGSORT(uint8_t, ndarray, index_array);
-        } else if((ndarray->array->typecode == NDARRAY_INT16) || (ndarray->array->typecode == NDARRAY_INT16)) {
-            HEAP_ARGSORT(uint16_t, ndarray, index_array);
-        } else {
-            HEAP_ARGSORT(mp_float_t, ndarray, index_array);
-        }
-    }
-    return MP_OBJ_FROM_PTR(indices);
-}
-
-MP_DEFINE_CONST_FUN_OBJ_KW(numerical_argsort_obj, 1, numerical_argsort);
-
-#if !CIRCUITPY
-STATIC const mp_rom_map_elem_t ulab_numerical_globals_table[] = {
-    { MP_OBJ_NEW_QSTR(MP_QSTR_linspace), (mp_obj_t)&numerical_linspace_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_sum), (mp_obj_t)&numerical_sum_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_mean), (mp_obj_t)&numerical_mean_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_std), (mp_obj_t)&numerical_std_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_min), (mp_obj_t)&numerical_min_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_max), (mp_obj_t)&numerical_max_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_argmin), (mp_obj_t)&numerical_argmin_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_argmax), (mp_obj_t)&numerical_argmax_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_roll), (mp_obj_t)&numerical_roll_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_flip), (mp_obj_t)&numerical_flip_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_diff), (mp_obj_t)&numerical_diff_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_sort), (mp_obj_t)&numerical_sort_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_argsort), (mp_obj_t)&numerical_argsort_obj },    
-};
-
-STATIC MP_DEFINE_CONST_DICT(mp_module_ulab_numerical_globals, ulab_numerical_globals_table);
-
-mp_obj_module_t ulab_numerical_module = {
-    .base = { &mp_type_module },
-    .globals = (mp_obj_dict_t*)&mp_module_ulab_numerical_globals,
-};
-#endif
-
-#endif

code/numerical.h

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,162 +5,62 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #ifndef _NUMERICAL_
 #define _NUMERICAL_
 
-#include "ulab.h"
 #include "ndarray.h"
 
-#if ULAB_NUMERICAL_MODULE
+mp_obj_t numerical_linspace(size_t , const mp_obj_t *, mp_map_t *);
 
-extern mp_obj_module_t ulab_numerical_module;
+mp_obj_t numerical_sum(size_t , const mp_obj_t *, mp_map_t *);
+mp_obj_t numerical_mean(size_t , const mp_obj_t *, mp_map_t *);
+mp_obj_t numerical_std(size_t , const mp_obj_t *, mp_map_t *);
 
-// TODO: implement minimum/maximum, and cumsum
-//mp_obj_t numerical_minimum(mp_obj_t , mp_obj_t );
-//mp_obj_t numerical_maximum(mp_obj_t , mp_obj_t );
-//mp_obj_t numerical_cumsum(size_t , const mp_obj_t *, mp_map_t *);
-
-#define RUN_ARGMIN(in, out, typein, typeout, len, start, increment, op, pos) do {\
-    typein *array = (typein *)(in)->array->items;\
-    typeout *outarray = (typeout *)(out)->array->items;\
-    size_t best_index = 0;\
-    if(((op) == NUMERICAL_MAX) || ((op) == NUMERICAL_ARGMAX)) {\
-        for(size_t i=1; i < (len); i++) {\
-            if(array[(start)+i*(increment)] > array[(start)+best_index*(increment)]) best_index = i;\
-        }\
-        if((op) == NUMERICAL_MAX) outarray[(pos)] = array[(start)+best_index*(increment)];\
-        else outarray[(pos)] = best_index;\
-    } else{\
-        for(size_t i=1; i < (len); i++) {\
-            if(array[(start)+i*(increment)] < array[(start)+best_index*(increment)]) best_index = i;\
-        }\
-        if((op) == NUMERICAL_MIN) outarray[(pos)] = array[(start)+best_index*(increment)];\
-        else outarray[(pos)] = best_index;\
-    }\
-} while(0)
-
-#define RUN_SUM(ndarray, type, optype, len, start, increment) do {\
+#define CALCULATE_SUM(ndarray, type, farray, shape_ax, index, stride, offset, optype) do {\
     type *array = (type *)(ndarray)->array->items;\
-    type value;\
-    for(size_t j=0; j < (len); j++) {\
-        value = array[(start)+j*(increment)];\
-        sum += value;\
+    (farray)[(index)] = 0.0;\
+    for(size_t j=0; j < (shape_ax); j++, (offset) += (stride)) {\
+        (farray)[(index)] += array[(offset)];\
     }\
 } while(0)
 
-#define RUN_STD(ndarray, type, len, start, increment) do {\
+// TODO: this can be done without the NDARRAY_INDEX_FROM_FLAT macro
+// Welford algorithm for the standard deviation
+#define CALCULATE_FLAT_SUM_STD(ndarray, type, value, shape_strides, optype) do {\
     type *array = (type *)(ndarray)->array->items;\
-    mp_float_t value;\
-    for(size_t j=0; j < (len); j++) {\
-        sum += array[(start)+j*(increment)];\
-    }\
-    sum /= (len);\
-    for(size_t j=0; j < (len); j++) {\
-        value = (array[(start)+j*(increment)] - sum);\
-        sum_sq += value * value;\
-    }\
-} while(0)
-
-#define CALCULATE_DIFF(in, out, type, M, N, inn, increment) do {\
-    type *source = (type *)(in)->array->items;\
-    type *target = (type *)(out)->array->items;\
-    for(size_t i=0; i < (M); i++) {\
-        for(size_t j=0; j < (N); j++) {\
-            for(uint8_t k=0; k < n+1; k++) {\
-                target[i*(N)+j] -= stencil[k]*source[i*(inn)+j+k*(increment)];\
-            }\
+    (value) = 0.0;\
+    mp_float_t m = 0.0, mtmp;\
+    size_t index, nindex;\
+    for(size_t j=0; j < (ndarray)->len; j++) {\
+        NDARRAY_INDEX_FROM_FLAT((ndarray), (shape_strides), j, index, nindex);\
+        if((optype) == NUMERICAL_STD) {\
+            mtmp = m;\
+            m = mtmp + (array[nindex] - mtmp) / (j+1);\
+            (value) += (array[nindex] - mtmp) * (array[nindex] - m);\
+        } else {\
+            (value) += array[nindex];\
         }\
     }\
 } while(0)
 
-#define HEAPSORT(type, ndarray) do {\
+// we calculate the standard deviation in two passes, in order to avoid negative values through truncation errors
+// We could do in a single pass, if we resorted to the Welford algorithm above
+#define CALCULATE_STD(ndarray, type, sq_sum, shape_ax, stride, offset) do {\
     type *array = (type *)(ndarray)->array->items;\
-    type tmp;\
-    for (;;) {\
-        if (k > 0) {\
-            tmp = array[start+(--k)*increment];\
-        } else {\
-            q--;\
-            if(q == 0) {\
-                break;\
+    mp_float_t x, ave = 0.0;\
+    (sq_sum) = 0.0;\
+    size_t j, _offset = (offset);\
+    for(j=0; j < (shape_ax); j++, _offset += (stride)) {\
+        ave += array[_offset];\
     }\
-            tmp = array[start+q*increment];\
-            array[start+q*increment] = array[start];\
-        }\
-        p = k;\
-        c = k + k + 1;\
-        while (c < q) {\
-            if((c + 1 < q)  &&  (array[start+(c+1)*increment] > array[start+c*increment])) {\
-                c++;\
-            }\
-            if(array[start+c*increment] > tmp) {\
-                array[start+p*increment] = array[start+c*increment];\
-                p = c;\
-                c = p + p + 1;\
-            } else {\
-                break;\
-            }\
-        }\
-        array[start+p*increment] = tmp;\
+    ave /= j;\
+    for(j=0; j < (shape_ax); j++, (offset) += (stride)) {\
+        x = array[(offset)] - ave;\
+        (sq_sum) += x * x;\
     }\
 } while(0)
 
-// This is pretty similar to HEAPSORT above; perhaps, the two could be combined somehow
-// On the other hand, since this is a macro, it doesn't really matter
-// Keep in mind that initially, index_array[start+s*increment] = s
-#define HEAP_ARGSORT(type, ndarray, index_array) do {\
-    type *array = (type *)(ndarray)->array->items;\
-    type tmp;\
-    uint16_t itmp;\
-    for (;;) {\
-        if (k > 0) {\
-            k--;\
-            tmp = array[start+index_array[start+k*increment]*increment];\
-            itmp = index_array[start+k*increment];\
-        } else {\
-            q--;\
-            if(q == 0) {\
-                break;\
-            }\
-            tmp = array[start+index_array[start+q*increment]*increment];\
-            itmp = index_array[start+q*increment];\
-            index_array[start+q*increment] = index_array[start];\
-        }\
-        p = k;\
-        c = k + k + 1;\
-        while (c < q) {\
-            if((c + 1 < q)  &&  (array[start+index_array[start+(c+1)*increment]*increment] > array[start+index_array[start+c*increment]*increment])) {\
-                c++;\
-            }\
-            if(array[start+index_array[start+c*increment]*increment] > tmp) {\
-                index_array[start+p*increment] = index_array[start+c*increment];\
-                p = c;\
-                c = p + p + 1;\
-            } else {\
-                break;\
-            }\
-        }\
-        index_array[start+p*increment] = itmp;\
-    }\
-} while(0)
-
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_linspace_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_min_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_max_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_argmin_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_argmax_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_sum_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_mean_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_std_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_roll_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_flip_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_diff_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_sort_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_sort_inplace_obj);
-MP_DECLARE_CONST_FUN_OBJ_KW(numerical_argsort_obj);
-
-#endif
 #endif

code/poly.c

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,7 +5,7 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #include "py/obj.h"
@@ -16,116 +15,73 @@
 #include "linalg.h"
 #include "poly.h"
 
-#if ULAB_POLY_MODULE
-bool object_is_nditerable(mp_obj_t o_in) {
-    if(MP_OBJ_IS_TYPE(o_in, &ulab_ndarray_type) || 
-      MP_OBJ_IS_TYPE(o_in, &mp_type_tuple) || 
-      MP_OBJ_IS_TYPE(o_in, &mp_type_list) || 
-      MP_OBJ_IS_TYPE(o_in, &mp_type_range)) {
-        return true;
-    }
-    return false;
-}
-
-size_t get_nditerable_len(mp_obj_t o_in) {
-    if(MP_OBJ_IS_TYPE(o_in, &ulab_ndarray_type)) {
-        ndarray_obj_t *in = MP_OBJ_TO_PTR(o_in);
-        return in->array->len;
-    } else {
-        return (size_t)mp_obj_get_int(mp_obj_len_maybe(o_in));
+void fill_array_iterable(mp_float_t *array, mp_obj_t oin) {
+    mp_obj_iter_buf_t buf;
+    mp_obj_t item, iterable = mp_getiter(oin, &buf);
+    while((item = mp_iternext(iterable)) != MP_OBJ_STOP_ITERATION) {
+        *array++ = mp_obj_get_float(item);
     }
 }
 
 mp_obj_t poly_polyval(mp_obj_t o_p, mp_obj_t o_x) {
-    // TODO: return immediately, if o_p is not an iterable
-    // TODO: there is a bug here: matrices won't work, 
-    // because there is a single iteration loop
-    size_t m, n;
-    if(MP_OBJ_IS_TYPE(o_x, &ulab_ndarray_type)) {
-        ndarray_obj_t *ndx = MP_OBJ_TO_PTR(o_x);
-        m = ndx->m;
-        n = ndx->n;
-    } else {
-        mp_obj_array_t *ix = MP_OBJ_TO_PTR(o_x);
-        m = 1;
-        n = ix->len;
-    }
-    // polynomials are going to be of type float, except, when both 
-    // the coefficients and the independent variable are integers
-    ndarray_obj_t *out = create_new_ndarray(m, n, NDARRAY_FLOAT);
-    mp_obj_iter_buf_t x_buf;
-    mp_obj_t x_item, x_iterable = mp_getiter(o_x, &x_buf);
-
-    mp_obj_iter_buf_t p_buf;
-    mp_obj_t p_item, p_iterable;
-
-    mp_float_t x, y;
-    mp_float_t *outf = (mp_float_t *)out->array->items;
+    // we always return floats: polynomials are going to be of type float, except, 
+    // when both the coefficients and the independent variable are integers; 
     uint8_t plen = mp_obj_get_int(mp_obj_len_maybe(o_p));
     mp_float_t *p = m_new(mp_float_t, plen);
-    p_iterable = mp_getiter(o_p, &p_buf);
-    uint16_t i = 0;    
-    while((p_item = mp_iternext(p_iterable)) != MP_OBJ_STOP_ITERATION) {
-        p[i] = mp_obj_get_float(p_item);
-        i++;
+    fill_array_iterable(p, o_p);
+    ndarray_obj_t *ndarray;
+    mp_float_t *array;
+    if(MP_OBJ_IS_TYPE(o_x, &ulab_ndarray_type)) {
+        ndarray_obj_t *input = MP_OBJ_TO_PTR(o_x);
+        ndarray = ndarray_copy_view(input, NDARRAY_FLOAT);
+        array = (mp_float_t *)ndarray->array->items;
+    } else { // at this point, we should have a 1-D iterable
+        size_t len = mp_obj_get_int(mp_obj_len_maybe(o_x));
+        ndarray = ndarray_new_linear_array(len, NDARRAY_FLOAT);
+        array = (mp_float_t *)ndarray->array->items;
+        fill_array_iterable(array, o_x);
     }
-    i = 0;
-    while ((x_item = mp_iternext(x_iterable)) != MP_OBJ_STOP_ITERATION) {
-        x = mp_obj_get_float(x_item);
+    mp_float_t x, y;
+    for(size_t i=0; i < ndarray->len; i++) {
+        x = array[i];
         y = p[0];
         for(uint8_t j=0; j < plen-1; j++) {
             y *= x;
             y += p[j+1];
         }
-        outf[i++] = y;
+        array[i] = y;
     }
     m_del(mp_float_t, p, plen);
-    return MP_OBJ_FROM_PTR(out);
+    return MP_OBJ_FROM_PTR(ndarray);
 }
 
-MP_DEFINE_CONST_FUN_OBJ_2(poly_polyval_obj, poly_polyval);
-
 mp_obj_t poly_polyfit(size_t n_args, const mp_obj_t *args) {
-    if((n_args != 2) && (n_args != 3)) {
-        mp_raise_ValueError(translate("number of arguments must be 2, or 3"));
+    if(n_args != 3) {
+        mp_raise_ValueError("number of arguments must be 3");
     }
-    if(!object_is_nditerable(args[0])) {
-        mp_raise_ValueError(translate("input data must be an iterable"));
+    if(!MP_OBJ_IS_TYPE(args[0], &ulab_ndarray_type) && !MP_OBJ_IS_TYPE(args[0], &mp_type_tuple) &&
+      !MP_OBJ_IS_TYPE(args[0], &mp_type_list) && !MP_OBJ_IS_TYPE(args[0], &mp_type_range) &&
+      !MP_OBJ_IS_TYPE(args[1], &ulab_ndarray_type) && !MP_OBJ_IS_TYPE(args[1], &mp_type_tuple) &&
+      !MP_OBJ_IS_TYPE(args[1], &mp_type_list) && !MP_OBJ_IS_TYPE(args[1], &mp_type_range)) {
+        mp_raise_ValueError("input data must be a 1D iterable");
     }
-    uint16_t lenx = 0, leny = 0;
-    uint8_t deg = 0;
+    uint16_t lenx, leny;
+    uint8_t deg;
     mp_float_t *x, *XT, *y, *prod;
 
-    if(n_args == 2) { // only the y values are supplied
-        // TODO: this is actually not enough: the first argument can very well be a matrix, 
-        // in which case we are between the rock and a hard place
-        leny = (uint16_t)mp_obj_get_int(mp_obj_len_maybe(args[0]));
-        deg = (uint8_t)mp_obj_get_int(args[1]);
-        if(leny < deg) {
-            mp_raise_ValueError(translate("more degrees of freedom than data points"));
-        }
-        lenx = leny;
-        x = m_new(mp_float_t, lenx); // assume uniformly spaced data points
-        for(size_t i=0; i < lenx; i++) {
-            x[i] = i;
-        }
-        y = m_new(mp_float_t, leny);
-        fill_array_iterable(y, args[0]);
-    } else if(n_args == 3) {
     lenx = (uint16_t)mp_obj_get_int(mp_obj_len_maybe(args[0]));
-        leny = (uint16_t)mp_obj_get_int(mp_obj_len_maybe(args[0]));
+    leny = (uint16_t)mp_obj_get_int(mp_obj_len_maybe(args[1]));
     if(lenx != leny) {
-            mp_raise_ValueError(translate("input vectors must be of equal length"));
+        mp_raise_ValueError("input vectors must be of equal length");
     }
     deg = (uint8_t)mp_obj_get_int(args[2]);
     if(leny < deg) {
-            mp_raise_ValueError(translate("more degrees of freedom than data points"));
+        mp_raise_ValueError("more degrees of freedom than data points");
     }
     x = m_new(mp_float_t, lenx);
     fill_array_iterable(x, args[0]);
     y = m_new(mp_float_t, leny);
     fill_array_iterable(y, args[1]);
-    }
     
     // one could probably express X as a function of XT, 
     // and thereby save RAM, because X is used only in the product
@@ -159,7 +115,7 @@ mp_obj_t poly_polyfit(size_t  n_args, const mp_obj_t *args) {
         m_del(mp_float_t, x, lenx);
         m_del(mp_float_t, y, lenx);
         m_del(mp_float_t, prod, (deg+1)*(deg+1));
-        mp_raise_ValueError(translate("could not invert Vandermonde matrix"));
+        mp_raise_ValueError("could not invert Vandermonde matrix");
     } 
     // at this point, we have the inverse of X^T * X
     // y is a column vector; x is free now, we can use it for storing intermediate values
@@ -173,7 +129,7 @@ mp_obj_t poly_polyfit(size_t  n_args, const mp_obj_t *args) {
     // XT is no longer needed
     m_del(mp_float_t, XT, (deg+1)*leny);
     
-    ndarray_obj_t *beta = create_new_ndarray(deg+1, 1, NDARRAY_FLOAT);
+    ndarray_obj_t *beta = ndarray_new_linear_array(deg+1, NDARRAY_FLOAT);
     mp_float_t *betav = (mp_float_t *)beta->array->items;
     // x[0..(deg+1)] contains now the product X^T * y; we can get rid of y
     m_del(float, y, leny);
@@ -194,22 +150,3 @@ mp_obj_t poly_polyfit(size_t  n_args, const mp_obj_t *args) {
     }
     return MP_OBJ_FROM_PTR(beta);
 }
-
-MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(poly_polyfit_obj, 2, 3, poly_polyfit);
-
-#if !CIRCUITPY
-STATIC const mp_rom_map_elem_t ulab_poly_globals_table[] = {
-    { MP_OBJ_NEW_QSTR(MP_QSTR___name__), MP_OBJ_NEW_QSTR(MP_QSTR_poly) },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_polyval), (mp_obj_t)&poly_polyval_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_polyfit), (mp_obj_t)&poly_polyfit_obj },
-};
-
-STATIC MP_DEFINE_CONST_DICT(mp_module_ulab_poly_globals, ulab_poly_globals_table);
-
-mp_obj_module_t ulab_poly_module = {
-    .base = { &mp_type_module },
-    .globals = (mp_obj_dict_t*)&mp_module_ulab_poly_globals,
-};
-#endif
-
-#endif

code/poly.h

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,20 +5,13 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #ifndef _POLY_
 #define _POLY_
 
-#include "ulab.h"
-
-#if ULAB_POLY_MODULE
-
-extern mp_obj_module_t ulab_poly_module;
-
-MP_DECLARE_CONST_FUN_OBJ_2(poly_polyval_obj);
-MP_DECLARE_CONST_FUN_OBJ_VAR_BETWEEN(poly_polyfit_obj);
+mp_obj_t poly_polyval(mp_obj_t , mp_obj_t );
+mp_obj_t poly_polyfit(size_t  , const mp_obj_t *);
 
 #endif
-#endif

code/ulab.c

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,7 +5,7 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #include <math.h>
@@ -18,27 +17,82 @@
 #include "py/obj.h"
 #include "py/objarray.h"
 
-#include "ulab.h"
 #include "ndarray.h"
-#include "ndarray_properties.h"
-#include "linalg.h"
 #include "vectorise.h"
+#include "linalg.h"
 #include "poly.h"
-#include "fft.h"
-#include "filter.h"
 #include "numerical.h"
-#include "extras.h"
+#include "fft.h"
 
-STATIC MP_DEFINE_STR_OBJ(ulab_version_obj, "0.34.0");
+#define ULAB_VERSION 1.0
+
+typedef struct _mp_obj_float_t {
+    mp_obj_base_t base;
+    mp_float_t value;
+} mp_obj_float_t;
+
+mp_obj_float_t ulab_version = {{&mp_type_float}, ULAB_VERSION};
+
+MP_DEFINE_CONST_FUN_OBJ_1(ndarray_shape_obj, ndarray_shape);
+MP_DEFINE_CONST_FUN_OBJ_1(ndarray_strides_obj, ndarray_strides);
+MP_DEFINE_CONST_FUN_OBJ_2(ndarray_reshape_obj, ndarray_reshape);
+MP_DEFINE_CONST_FUN_OBJ_1(ndarray_transpose_obj, ndarray_transpose);
+MP_DEFINE_CONST_FUN_OBJ_KW(ndarray_flatten_obj, 1, ndarray_flatten);
+MP_DEFINE_CONST_FUN_OBJ_1(ndarray_itemsize_obj, ndarray_itemsize);
+MP_DEFINE_CONST_FUN_OBJ_1(ndarray_info_obj, ndarray_info);
+
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_acos_obj, vectorise_acos);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_acosh_obj, vectorise_acosh);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_asin_obj, vectorise_asin);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_asinh_obj, vectorise_asinh);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_atan_obj, vectorise_atan);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_atanh_obj, vectorise_atanh);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_ceil_obj, vectorise_ceil);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_cos_obj, vectorise_cos);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_erf_obj, vectorise_erf);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_erfc_obj, vectorise_erfc);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_exp_obj, vectorise_exp);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_expm1_obj, vectorise_expm1);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_floor_obj, vectorise_floor);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_gamma_obj, vectorise_gamma);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_lgamma_obj, vectorise_lgamma);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_log_obj, vectorise_log);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_log10_obj, vectorise_log10);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_log2_obj, vectorise_log2);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_sin_obj, vectorise_sin);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_sinh_obj, vectorise_sinh);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_sqrt_obj, vectorise_sqrt);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_tan_obj, vectorise_tan);
+MP_DEFINE_CONST_FUN_OBJ_1(vectorise_tanh_obj, vectorise_tanh);
+
+MP_DEFINE_CONST_FUN_OBJ_1(linalg_inv_obj, linalg_inv);
+MP_DEFINE_CONST_FUN_OBJ_2(linalg_dot_obj, linalg_dot);
+MP_DEFINE_CONST_FUN_OBJ_KW(linalg_zeros_obj, 0, linalg_zeros);
+MP_DEFINE_CONST_FUN_OBJ_KW(linalg_ones_obj, 0, linalg_ones);
+MP_DEFINE_CONST_FUN_OBJ_KW(linalg_eye_obj, 0, linalg_eye);
+MP_DEFINE_CONST_FUN_OBJ_1(linalg_det_obj, linalg_det);
+MP_DEFINE_CONST_FUN_OBJ_1(linalg_eig_obj, linalg_eig);
+
+MP_DEFINE_CONST_FUN_OBJ_2(poly_polyval_obj, poly_polyval);
+MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(poly_polyfit_obj, 2, 3, poly_polyfit);
+
+MP_DEFINE_CONST_FUN_OBJ_KW(numerical_linspace_obj, 2, numerical_linspace);
+MP_DEFINE_CONST_FUN_OBJ_KW(numerical_sum_obj, 1, numerical_sum);
+MP_DEFINE_CONST_FUN_OBJ_KW(numerical_mean_obj, 1, numerical_mean);
+MP_DEFINE_CONST_FUN_OBJ_KW(numerical_std_obj, 1, numerical_std);
+
+
+MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(fft_fft_obj, 1, 2, fft_fft);
+MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(fft_ifft_obj, 1, 2, fft_ifft);
+MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(fft_spectrum_obj, 1, 2, fft_spectrum);
 
 STATIC const mp_rom_map_elem_t ulab_ndarray_locals_dict_table[] = {
-    { MP_ROM_QSTR(MP_QSTR_flatten), MP_ROM_PTR(&ndarray_flatten_obj) },
+    { MP_ROM_QSTR(MP_QSTR_shape), MP_ROM_PTR(&ndarray_shape_obj) },
+    { MP_ROM_QSTR(MP_QSTR_strides), MP_ROM_PTR(&ndarray_strides_obj) },
     { MP_ROM_QSTR(MP_QSTR_reshape), MP_ROM_PTR(&ndarray_reshape_obj) },
     { MP_ROM_QSTR(MP_QSTR_transpose), MP_ROM_PTR(&ndarray_transpose_obj) },
-    { MP_ROM_QSTR(MP_QSTR_shape), MP_ROM_PTR(&ndarray_shape_obj) },
-    { MP_ROM_QSTR(MP_QSTR_size), MP_ROM_PTR(&ndarray_size_obj) },
+    { MP_ROM_QSTR(MP_QSTR_flatten), MP_ROM_PTR(&ndarray_flatten_obj) },
     { MP_ROM_QSTR(MP_QSTR_itemsize), MP_ROM_PTR(&ndarray_itemsize_obj) },
-//    { MP_ROM_QSTR(MP_QSTR_sort), MP_ROM_PTR(&numerical_sort_inplace_obj) },
 };
 
 STATIC MP_DEFINE_CONST_DICT(ulab_ndarray_locals_dict, ulab_ndarray_locals_dict_table);
@@ -52,37 +106,55 @@ const mp_obj_type_t ulab_ndarray_type = {
     .getiter = ndarray_getiter,
     .unary_op = ndarray_unary_op,
     .binary_op = ndarray_binary_op,
-    .buffer_p = { .get_buffer = ndarray_get_buffer, },
     .locals_dict = (mp_obj_dict_t*)&ulab_ndarray_locals_dict,
 };
 
-#if !CIRCUITPY
 STATIC const mp_map_elem_t ulab_globals_table[] = {
     { MP_OBJ_NEW_QSTR(MP_QSTR___name__), MP_OBJ_NEW_QSTR(MP_QSTR_ulab) },
-    { MP_ROM_QSTR(MP_QSTR___version__), MP_ROM_PTR(&ulab_version_obj) },
+    { MP_ROM_QSTR(MP_QSTR___version__), MP_ROM_PTR(&ulab_version) },
     { MP_OBJ_NEW_QSTR(MP_QSTR_array), (mp_obj_t)&ulab_ndarray_type },
-    #if ULAB_LINALG_MODULE
-    { MP_ROM_QSTR(MP_QSTR_linalg), MP_ROM_PTR(&ulab_linalg_module) },
-    #endif
-    #if ULAB_VECTORISE_MODULE
-    { MP_ROM_QSTR(MP_QSTR_vector), MP_ROM_PTR(&ulab_vectorise_module) },
-    #endif
-    #if ULAB_NUMERICAL_MODULE
-    { MP_ROM_QSTR(MP_QSTR_numerical), MP_ROM_PTR(&ulab_numerical_module) },
-    #endif
-    #if ULAB_POLY_MODULE
-    { MP_ROM_QSTR(MP_QSTR_poly), MP_ROM_PTR(&ulab_poly_module) },
-    #endif
-    #if ULAB_FFT_MODULE
-    { MP_ROM_QSTR(MP_QSTR_fft), MP_ROM_PTR(&ulab_fft_module) },
-    #endif
-    #if ULAB_FILTER_MODULE
-    { MP_ROM_QSTR(MP_QSTR_filter), MP_ROM_PTR(&ulab_filter_module) },
-    #endif
-    #if ULAB_EXTRAS_MODULE
-    { MP_ROM_QSTR(MP_QSTR_extras), MP_ROM_PTR(&ulab_extras_module) },
-    #endif
+    { MP_OBJ_NEW_QSTR(MP_QSTR_ndinfo), (mp_obj_t)&ndarray_info_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_acos), (mp_obj_t)&vectorise_acos_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_acosh), (mp_obj_t)&vectorise_acosh_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_asin), (mp_obj_t)&vectorise_asin_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_asinh), (mp_obj_t)&vectorise_asinh_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_atan), (mp_obj_t)&vectorise_atan_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_atanh), (mp_obj_t)&vectorise_atanh_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_ceil), (mp_obj_t)&vectorise_ceil_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_cos), (mp_obj_t)&vectorise_cos_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_erf), (mp_obj_t)&vectorise_erf_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_erfc), (mp_obj_t)&vectorise_erfc_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_exp), (mp_obj_t)&vectorise_exp_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_expm1), (mp_obj_t)&vectorise_expm1_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_floor), (mp_obj_t)&vectorise_floor_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_gamma), (mp_obj_t)&vectorise_gamma_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_lgamma), (mp_obj_t)&vectorise_lgamma_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_log), (mp_obj_t)&vectorise_log_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_log10), (mp_obj_t)&vectorise_log10_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_log2), (mp_obj_t)&vectorise_log2_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_sin), (mp_obj_t)&vectorise_sin_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_sinh), (mp_obj_t)&vectorise_sinh_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_sqrt), (mp_obj_t)&vectorise_sqrt_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_tan), (mp_obj_t)&vectorise_tan_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_tanh), (mp_obj_t)&vectorise_tanh_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_polyval), (mp_obj_t)&poly_polyval_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_polyfit), (mp_obj_t)&poly_polyfit_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_inv), (mp_obj_t)&linalg_inv_obj },
+    { MP_ROM_QSTR(MP_QSTR_dot), (mp_obj_t)&linalg_dot_obj },
+    { MP_ROM_QSTR(MP_QSTR_zeros), (mp_obj_t)&linalg_zeros_obj },
+    { MP_ROM_QSTR(MP_QSTR_ones), (mp_obj_t)&linalg_ones_obj },
+    { MP_ROM_QSTR(MP_QSTR_eye), (mp_obj_t)&linalg_eye_obj },
+    { MP_ROM_QSTR(MP_QSTR_det), (mp_obj_t)&linalg_det_obj },
+    { MP_ROM_QSTR(MP_QSTR_eig), (mp_obj_t)&linalg_eig_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_linspace), (mp_obj_t)&numerical_linspace_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_sum), (mp_obj_t)&numerical_sum_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_mean), (mp_obj_t)&numerical_mean_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_std), (mp_obj_t)&numerical_std_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_fft), (mp_obj_t)&fft_fft_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_ifft), (mp_obj_t)&fft_ifft_obj },
+    { MP_OBJ_NEW_QSTR(MP_QSTR_spectrum), (mp_obj_t)&fft_spectrum_obj },
     // class constants
+    { MP_ROM_QSTR(MP_QSTR_bool), MP_ROM_INT(NDARRAY_BOOL) },
     { MP_ROM_QSTR(MP_QSTR_uint8), MP_ROM_INT(NDARRAY_UINT8) },
     { MP_ROM_QSTR(MP_QSTR_int8), MP_ROM_INT(NDARRAY_INT8) },
     { MP_ROM_QSTR(MP_QSTR_uint16), MP_ROM_INT(NDARRAY_UINT16) },
@@ -95,10 +167,9 @@ STATIC MP_DEFINE_CONST_DICT (
     ulab_globals_table
 );
 
-mp_obj_module_t ulab_user_cmodule = {
+const mp_obj_module_t ulab_user_cmodule = {
     .base = { &mp_type_module },
     .globals = (mp_obj_dict_t*)&mp_module_ulab_globals,
 };
 
 MP_REGISTER_MODULE(MP_QSTR_ulab, ulab_user_cmodule, MODULE_ULAB_ENABLED);
-#endif

code/ulab.h

@@ -1,36 +0,0 @@
-
-/*
- * This file is part of the micropython-ulab project, 
- *
- * https://github.com/v923z/micropython-ulab
- *
- * The MIT License (MIT)
- *
- * Copyright (c) 2019-2020 Zoltán Vörös
-*/
-
-#ifndef __ULAB__
-#define __ULAB__
-
-// vectorise (all functions) takes approx. 3 kB of flash space
-#define ULAB_VECTORISE_MODULE (1)
-
-// linalg adds around 6 kB
-#define ULAB_LINALG_MODULE (1)
-
-// poly is approx. 2.5 kB
-#define ULAB_POLY_MODULE (1)
-
-// numerical is about 12 kB
-#define ULAB_NUMERICAL_MODULE (1)
-
-// FFT costs about 2 kB of flash space
-#define ULAB_FFT_MODULE (1)
-
-// the filter module takes about 1 kB of flash space
-#define ULAB_FILTER_MODULE (1)
-
-// user-defined modules
-#define ULAB_EXTRAS_MODULE (0)
-
-#endif

code/vectorise.c

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,7 +5,7 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #include <math.h>
@@ -22,153 +21,72 @@
 #define MP_PI MICROPY_FLOAT_CONST(3.14159265358979323846)
 #endif
     
-#if ULAB_VECTORISE_MODULE
 mp_obj_t vectorise_generic_vector(mp_obj_t o_in, mp_float_t (*f)(mp_float_t)) {
     // Return a single value, if o_in is not iterable
-    if(mp_obj_is_float(o_in) || MP_OBJ_IS_INT(o_in)) {
+    if(mp_obj_is_float(o_in) || mp_obj_is_integer(o_in)) {
         return mp_obj_new_float(f(mp_obj_get_float(o_in)));
     }
-    mp_float_t x;
     if(MP_OBJ_IS_TYPE(o_in, &ulab_ndarray_type)) {
         ndarray_obj_t *source = MP_OBJ_TO_PTR(o_in);
-        ndarray_obj_t *ndarray = create_new_ndarray(source->m, source->n, NDARRAY_FLOAT);
-        mp_float_t *dataout = (mp_float_t *)ndarray->array->items;
+        ndarray_obj_t *ndarray;
+        if(source->len == source->array->len) { // this is a dense array
+            ndarray = ndarray_new_dense_ndarray(source->ndim, source->shape, NDARRAY_FLOAT);
+            mp_float_t *array = (mp_float_t *)ndarray->array->items;
             if(source->array->typecode == NDARRAY_UINT8) {
-            ITERATE_VECTOR(uint8_t, source, dataout);
+                ITERATE_VECTOR(uint8_t, source, array);
             } else if(source->array->typecode == NDARRAY_INT8) {
-            ITERATE_VECTOR(int8_t, source, dataout);
+                ITERATE_VECTOR(int8_t, source, array);
             } else if(source->array->typecode == NDARRAY_UINT16) {
-            ITERATE_VECTOR(uint16_t, source, dataout);
+                ITERATE_VECTOR(uint16_t, source, array);
             } else if(source->array->typecode == NDARRAY_INT16) {
-            ITERATE_VECTOR(int16_t, source, dataout);
+                ITERATE_VECTOR(int16_t, source, array);
             } else {
-            ITERATE_VECTOR(mp_float_t, source, dataout);
+                ITERATE_VECTOR(mp_float_t, source, array);
+            } 
+        } else {
+            ndarray = ndarray_copy_view(source, NDARRAY_FLOAT);
+            mp_float_t *dataout = (mp_float_t *)ndarray->array->items;
+            ITERATE_VECTOR(mp_float_t, ndarray, dataout);
         }
         return MP_OBJ_FROM_PTR(ndarray);
     } else if(MP_OBJ_IS_TYPE(o_in, &mp_type_tuple) || MP_OBJ_IS_TYPE(o_in, &mp_type_list) || 
-        MP_OBJ_IS_TYPE(o_in, &mp_type_range)) { // i.e., the input is a generic iterable
+        MP_OBJ_IS_TYPE(o_in, &mp_type_range)) { // i.e., the input is a generic, one-dimensional iterable
             mp_obj_array_t *o = MP_OBJ_TO_PTR(o_in);
-            ndarray_obj_t *out = create_new_ndarray(1, o->len, NDARRAY_FLOAT);
+            ndarray_obj_t *out = ndarray_new_linear_array(o->len, NDARRAY_FLOAT);
             mp_float_t *dataout = (mp_float_t *)out->array->items;
+            
             mp_obj_iter_buf_t iter_buf;
             mp_obj_t item, iterable = mp_getiter(o_in, &iter_buf);
-            size_t i=0;
+            mp_float_t x;
             while ((item = mp_iternext(iterable)) != MP_OBJ_STOP_ITERATION) {
                 x = mp_obj_get_float(item);
-                dataout[i++] = f(x);
+                *dataout++ = f(x);
             }
         return MP_OBJ_FROM_PTR(out);
     }
     return mp_const_none;
 }
 
-
 MATH_FUN_1(acos, acos);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_acos_obj, vectorise_acos);
-
 MATH_FUN_1(acosh, acosh);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_acosh_obj, vectorise_acosh);
-
 MATH_FUN_1(asin, asin);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_asin_obj, vectorise_asin);
-
 MATH_FUN_1(asinh, asinh);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_asinh_obj, vectorise_asinh);
-
 MATH_FUN_1(atan, atan);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_atan_obj, vectorise_atan);
-
 MATH_FUN_1(atanh, atanh);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_atanh_obj, vectorise_atanh);
-
 MATH_FUN_1(ceil, ceil);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_ceil_obj, vectorise_ceil);
-
 MATH_FUN_1(cos, cos);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_cos_obj, vectorise_cos);
-
-MATH_FUN_1(cosh, cosh);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_cosh_obj, vectorise_cosh);
-
 MATH_FUN_1(erf, erf);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_erf_obj, vectorise_erf);
-
 MATH_FUN_1(erfc, erfc);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_erfc_obj, vectorise_erfc);
-
 MATH_FUN_1(exp, exp);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_exp_obj, vectorise_exp);
-
 MATH_FUN_1(expm1, expm1);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_expm1_obj, vectorise_expm1);
-
 MATH_FUN_1(floor, floor);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_floor_obj, vectorise_floor);
-
 MATH_FUN_1(gamma, tgamma);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_gamma_obj, vectorise_gamma);
-
 MATH_FUN_1(lgamma, lgamma);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_lgamma_obj, vectorise_lgamma);
-
 MATH_FUN_1(log, log);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_log_obj, vectorise_log);
-
 MATH_FUN_1(log10, log10);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_log10_obj, vectorise_log10);
-
 MATH_FUN_1(log2, log2);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_log2_obj, vectorise_log2);
-
 MATH_FUN_1(sin, sin);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_sin_obj, vectorise_sin);
-
 MATH_FUN_1(sinh, sinh);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_sinh_obj, vectorise_sinh);
-
 MATH_FUN_1(sqrt, sqrt);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_sqrt_obj, vectorise_sqrt);
-
 MATH_FUN_1(tan, tan);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_tan_obj, vectorise_tan);
-
 MATH_FUN_1(tanh, tanh);
-MP_DEFINE_CONST_FUN_OBJ_1(vectorise_tanh_obj, vectorise_tanh);
-
-#if !CIRCUITPY
-STATIC const mp_rom_map_elem_t ulab_vectorise_globals_table[] = {
-    { MP_OBJ_NEW_QSTR(MP_QSTR___name__), MP_OBJ_NEW_QSTR(MP_QSTR_vector) },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_acos), (mp_obj_t)&vectorise_acos_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_acosh), (mp_obj_t)&vectorise_acosh_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_asin), (mp_obj_t)&vectorise_asin_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_asinh), (mp_obj_t)&vectorise_asinh_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_atan), (mp_obj_t)&vectorise_atan_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_atanh), (mp_obj_t)&vectorise_atanh_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_ceil), (mp_obj_t)&vectorise_ceil_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_cos), (mp_obj_t)&vectorise_cos_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_erf), (mp_obj_t)&vectorise_erf_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_erfc), (mp_obj_t)&vectorise_erfc_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_exp), (mp_obj_t)&vectorise_exp_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_expm1), (mp_obj_t)&vectorise_expm1_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_floor), (mp_obj_t)&vectorise_floor_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_gamma), (mp_obj_t)&vectorise_gamma_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_lgamma), (mp_obj_t)&vectorise_lgamma_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_log), (mp_obj_t)&vectorise_log_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_log10), (mp_obj_t)&vectorise_log10_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_log2), (mp_obj_t)&vectorise_log2_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_sin), (mp_obj_t)&vectorise_sin_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_sinh), (mp_obj_t)&vectorise_sinh_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_sqrt), (mp_obj_t)&vectorise_sqrt_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_tan), (mp_obj_t)&vectorise_tan_obj },
-    { MP_OBJ_NEW_QSTR(MP_QSTR_tanh), (mp_obj_t)&vectorise_tanh_obj },
-};
-
-STATIC MP_DEFINE_CONST_DICT(mp_module_ulab_vectorise_globals, ulab_vectorise_globals_table);
-
-mp_obj_module_t ulab_vectorise_module = {
-    .base = { &mp_type_module },
-    .globals = (mp_obj_dict_t*)&mp_module_ulab_vectorise_globals,
-};
-#endif
-
-#endif

code/vectorise.h

@@ -1,4 +1,3 @@
-
 /*
  * This file is part of the micropython-ulab project, 
  *
@@ -6,23 +5,42 @@
  *
  * The MIT License (MIT)
  *
- * Copyright (c) 2019-2020 Zoltán Vörös
+ * Copyright (c) 2019 Zoltán Vörös
 */
     
 #ifndef _VECTORISE_
 #define _VECTORISE_
 
-#include "ulab.h"
 #include "ndarray.h"
 
-#if ULAB_VECTORISE_MODULE
-
-mp_obj_module_t ulab_vectorise_module;
+mp_obj_t vectorise_acos(mp_obj_t );
+mp_obj_t vectorise_acosh(mp_obj_t );
+mp_obj_t vectorise_asin(mp_obj_t );
+mp_obj_t vectorise_asinh(mp_obj_t );
+mp_obj_t vectorise_atan(mp_obj_t );
+mp_obj_t vectorise_atanh(mp_obj_t );
+mp_obj_t vectorise_ceil(mp_obj_t );
+mp_obj_t vectorise_cos(mp_obj_t );
+mp_obj_t vectorise_erf(mp_obj_t );
+mp_obj_t vectorise_erfc(mp_obj_t );
+mp_obj_t vectorise_exp(mp_obj_t );
+mp_obj_t vectorise_expm1(mp_obj_t );
+mp_obj_t vectorise_floor(mp_obj_t );
+mp_obj_t vectorise_gamma(mp_obj_t );
+mp_obj_t vectorise_lgamma(mp_obj_t );
+mp_obj_t vectorise_log(mp_obj_t );
+mp_obj_t vectorise_log10(mp_obj_t );
+mp_obj_t vectorise_log2(mp_obj_t );
+mp_obj_t vectorise_sin(mp_obj_t );
+mp_obj_t vectorise_sinh(mp_obj_t );
+mp_obj_t vectorise_sqrt(mp_obj_t );
+mp_obj_t vectorise_tan(mp_obj_t );
+mp_obj_t vectorise_tanh(mp_obj_t );
 
 #define ITERATE_VECTOR(type, source, out) do {\
     type *input = (type *)(source)->array->items;\
     for(size_t i=0; i < (source)->array->len; i++) {\
-                (out)[i] = f(input[i]);\
+        *(out)++ = f(*input++);\
     }\
 } while(0)
 
@@ -32,4 +50,3 @@ mp_obj_module_t ulab_vectorise_module;
     }
     
 #endif
-#endif

docs/manual/source/conf.py

@@ -22,7 +22,7 @@ copyright = '2019, Zoltán Vörös'
 author = 'Zoltán Vörös'
 
 # The full version, including alpha/beta/rc tags
-release = '0.32'
+release = '0.26'
 
 
 # -- General configuration ---------------------------------------------------

docs/manual/source/ulab.rst

@@ -1,10 +1,10 @@
 Introduction
 ============
 
-In the `last
-chapter <https://micropython-usermod.readthedocs.io/en/latest/usermods_15.html>`__
-of the usermod documentation, I mentioned that I have another story, for
-another day. The day has come, so here is my story.
+In
+https://micropython-usermod.readthedocs.io/en/latest/usermods_14.html, I
+mentioned that I have another story, for another day. The day has come,
+so here is my story.
 
 Enter ulab
 ----------
@@ -68,11 +68,9 @@ The main points of ``ulab`` are
 -  polynomial fits to numerical data
 -  fast Fourier transforms
 
-At the time of writing this manual (for version 0.32), the library adds
+At the time of writing this manual (for version 0.26), the library adds
 approximately 30 kB of extra compiled code to the micropython
-(pyboard.v.11) firmware. However, if you are tight with flash space, you
-can easily shave off a couple of kB. See the section on `customising
-ulab <#Custom_builds>`__.
+(pyboard.v.11) firmware.
 
 Resources and legal matters
 ---------------------------
@@ -145,35 +143,6 @@ can always be queried as
 
 If you find a bug, please, include this number in your report!
 
-Customising ``ulab``
---------------------
-
-``ulab`` implements a great number of functions, and it is quite
-possible that you do not need all of them in a particular application.
-If you want to save some flash space, you can easily exclude arbitrary
-functions from the firmware. The
-`https://github.com/v923z/micropython-ulab/blob/master/code/ulab.h <ulab.h>`__
-header file contains a pre-processor flag for all functions in ``ulab``.
-The default setting is 1 for each of them, but if you change that to 0,
-the corresponding function will not be part of the compiled firmware.
-
-The first couple of lines of the file look like this
-
-.. code:: c
-
-   // vectorise (all functions) takes approx. 3 kB of flash space
-   #define ULAB_VECTORISE_ACOS (1)
-   #define ULAB_VECTORISE_ACOSH (1)
-   #define ULAB_VECTORISE_ASIN (1)
-   #define ULAB_VECTORISE_ASINH (1)
-   #define ULAB_VECTORISE_ATAN (1)
-   #define ULAB_VECTORISE_ATANH (1)
-
-In order to simplify navigation in the file, each flag begins with
-``ULAB_``, continues with the sub-module, where the function itself is
-implemented, and ends with the function’s name. Each section displays a
-hint as to how much space you can save by un-setting the flag.
-
 Basic ndarray operations
 ------------------------
 
@@ -197,10 +166,14 @@ Methods of ndarrays
 
 `.reshape <#.reshape>`__
 
+`.rawsize\*\* <#.rawsize>`__
+
 `.transpose <#.transpose>`__
 
 `.flatten\*\* <#.flatten>`__
 
+`.asbytearray <#.asbytearray>`__
+
 Matrix methods
 --------------
 
@@ -266,11 +239,6 @@ FFT routines
 
 `spectrum\*\* <#spectrum>`__
 
-Filter functions
-----------------
-
-`convolve <#convolve>`__
-
 ndarray, the basic container
 ============================
 
@@ -424,8 +392,8 @@ Methods of ndarrays
 .shape
 ~~~~~~
 
-The ``.shape`` method (property) returns a 2-tuple with the number of
-rows, and columns.
+The ``.shape`` method returns a 2-tuple with the number of rows, and
+columns.
 
 .. code::
         
@@ -435,11 +403,11 @@ rows, and columns.
     
     a = np.array([1, 2, 3, 4], dtype=np.int8)
     print("a:\n", a)
-    print("shape of a:", a.shape)
+    print("shape of a:", a.shape())
     
     b= np.array([[1, 2], [3, 4]], dtype=np.int8)
     print("\nb:\n", b)
-    print("shape of b:", b.shape)
+    print("shape of b:", b.shape())
 
 .. parsed-literal::
 
@@ -455,74 +423,6 @@ rows, and columns.
     
 
 
-.size
-~~~~~
-
-The ``.size`` method (property) returns an integer with the number of
-elements in the array.
-
-.. code::
-        
-    # code to be run in micropython
-    
-    import ulab as np
-    
-    a = np.array([1, 2, 3], dtype=np.int8)
-    print("a:\n", a)
-    print("size of a:", a.size)
-    
-    b= np.array([[1, 2], [3, 4]], dtype=np.int8)
-    print("\nb:\n", b)
-    print("size of b:", b.size)
-
-.. parsed-literal::
-
-    a:
-     array([1, 2, 3], dtype=int8)
-    size of a: 3
-    
-    b:
-     array([[1, 2],
-    	 [3, 4]], dtype=int8)
-    size of b: 4
-    
-    
-
-
-.itemsize
-~~~~~~~~~
-
-The ``.itemsize`` method (property) returns an integer with the siz
-enumber of elements in the array.
-
-.. code::
-        
-    # code to be run in micropython
-    
-    import ulab as np
-    
-    a = np.array([1, 2, 3], dtype=np.int8)
-    print("a:\n", a)
-    print("itemsize of a:", a.itemsize)
-    
-    b= np.array([[1, 2], [3, 4]], dtype=np.float)
-    print("\nb:\n", b)
-    print("itemsize of b:", b.itemsize)
-
-.. parsed-literal::
-
-    a:
-     array([1, 2, 3], dtype=int8)
-    itemsize of a: 1
-    
-    b:
-     array([[1.0, 2.0],
-    	 [3.0, 4.0]], dtype=float)
-    itemsize of b: 8
-    
-    
-
-
 .reshape
 ~~~~~~~~
 
@@ -559,6 +459,41 @@ consistent with the old, a ``ValueError`` exception will be raised.
     
 
 
+.rawsize
+~~~~~~~~
+
+The ``rawsize`` method of the ``ndarray`` returns a 5-tuple with the
+following data
+
+1. number of rows
+2. number of columns
+3. length of the storage (should be equal to the product of 1. and 2.)
+4. length of the data storage in bytes
+5. datum size in bytes (1 for ``uint8``/``int8``, 2 for
+   ``uint16``/``int16``, and 4, or 8 for ``floats``, see `ndarray, the
+   basic container <#ndarray,-the-basic-container>`__)
+
+**WARNING:** ``rawsize`` is a ``ulab``-only method; it has no equivalent
+in ``numpy``.
+
+.. code::
+        
+    # code to be run in micropython
+    
+    import ulab as np
+    
+    a = np.array([1, 2, 3, 4], dtype=np.float)
+    print("a: \t\t", a)
+    print("rawsize of a: \t", a.rawsize())
+
+.. parsed-literal::
+
+    a: 		 array([1.0, 2.0, 3.0, 4.0], dtype=float)
+    rawsize of a: 	 (1, 4, 4, 16, 4)
+    
+    
+
+
 .flatten
 ~~~~~~~~
 
@@ -600,6 +535,83 @@ verbatim copy of the contents.
     
 
 
+.asbytearray
+~~~~~~~~~~~~
+
+The contents of an ``ndarray`` can be accessed directly by calling the
+``.asbytearray`` method. This will simply return a pointer to the
+underlying flat ``array`` object, which can then be manipulated
+directly.
+
+**WARNING:** ``asbytearray`` is a ``ulab``-only method; it has no
+equivalent in ``numpy``.
+
+In the example below, note the difference between ``a``, and ``buffer``:
+while both are designated as an array, you recognise the micropython
+array from the fact that it prints the typecode (``b`` in this
+particular case). The ``ndarray``, on the other hand, prints out the
+``dtype`` (``int8`` here).
+
+.. code::
+        
+    # code to be run in micropython
+    
+    import ulab as np
+    
+    a = np.array([1, 2, 3, 4], dtype=np.int8)
+    print('a: ', a)
+    buffer = a.asbytearray()
+    print("array content:", buffer)
+    buffer[1] = 123
+    print("array content:", buffer)
+
+.. parsed-literal::
+
+    a:  array([1, 2, 3, 4], dtype=int8)
+    array content: array('b', [1, 2, 3, 4])
+    array content: array('b', [1, 123, 3, 4])
+    
+    
+
+
+This in itself wouldn’t be very interesting, but since ``buffer`` is a
+proper micropython ``array``, we can pass it to functions that can
+employ the buffer protocol. E.g., all the ``ndarray`` facilities can be
+applied to the results of timed ADC conversions.
+
+.. code::
+        
+    # code to be run in micropython
+    
+    import pyb
+    import ulab as np
+    
+    n = 100
+    
+    adc = pyb.ADC(pyb.Pin.board.X19)
+    tim = pyb.Timer(6, freq=10)
+    
+    a = np.array([0]*n, dtype=np.uint8)
+    buffer = a.asbytearray()
+    adc.read_timed(buffer, tim)
+    
+    print("ADC results:\t", a)
+    print("mean of results:\t", np.mean(a))
+    print("std of results:\t", np.std(a))
+
+.. parsed-literal::
+
+    ADC results:	 array([48, 2, 2, ..., 0, 0, 0], dtype=uint8)
+    mean of results:	 1.22
+    std of results:	 4.744639
+    
+
+
+Likewise, data can be read directly into ``ndarray``\ s from other
+interfaces, e.g., SPI, I2C etc, and also, by laying bare the
+``ndarray``, we can pass results of ``ulab`` computations to anything
+that can read from a buffer.
+
 .transpose
 ~~~~~~~~~~
 
@@ -2977,40 +2989,6 @@ Y_2(k) &=& -\frac{i}{2}\left(Y(k) - Y^*(N-k)\right)
 :math:`Y_1, Y_2`, and :math:`Y`, respectively, are the Fourier
 transforms of :math:`y_1, y_2`, and :math:`y = y_1 + iy_2`.
 
-Filter routines
-===============
-
-numpy:
-https://docs.scipy.org/doc/numpy/reference/generated/numpy.convolve.html
-
-convolve
---------
-
-Returns the discrete, linear convolution of two one-dimensional
-sequences.
-
-Only the ``full`` mode is supported, and the ``mode`` named parameter is
-not accepted. Note that all other modes can be had by slicing a ``full``
-result.
-
-.. code::
-        
-    # code to be run in micropython
-    
-    import ulab as np
-    
-    x = np.array((1,2,3))
-    y = np.array((1,10,100,1000))
-    
-    print(np.convolve(x, y))
-
-.. parsed-literal::
-
-    array([1.0, 12.0, 123.0, 1230.0, 2300.0, 3000.0], dtype=float)
-    
-    
-
-
 Extending ulab
 ==============
 

docs/ulab-change-log.md

@@ -1,100 +1,4 @@
 
-Tue, 18 Feb 2020
-
-version 0.34.0
-
-    split ulab into multiple modules
-
-Sun, 16 Feb 2020
-
-version 0.33.2
-
-    moved properties into ndarray_properties.h, implemented pointer arithmetic in fft.c to save some time
-
-Fri, 14 Feb 2020
-
-version 0.33.1
-
-    added the __name__attribute to all sub-modules
-
-Thu, 13 Feb 2020
-
-version 0.33.0
-
-    sub-modules are now proper sub-modules of ulab
-
-Mon, 17 Feb 2020
-
-version 0.32.1
-
-    temporary fix for issue #40
-    
-Tue, 11 Feb 2020
-
-version 0.32.0
-
-    added itemsize, size and shape attributes to ndarrays, and removed rawsize
-
-Mon, 10 Feb 2020
-
-version 0.31.0
-
-    removed asbytearray, and added buffer protocol to ndarrays, fixed bad error in filter.c
-    
-Sun, 09 Feb 2020
-
-version 0.30.2
-
-    fixed slice_length in ndarray.c
-
-Sat, 08 Feb 2020
-
-version 0.30.1
-
-    fixed typecode error, added variable inspection, and replaced ternary operators in filter.c
-    
-Fri, 07 Feb 2020
-
-version 0.30.0
-
-    ulab functions can arbitrarily be excluded from the firmware via the ulab.h configuration file
-    
-Thu, 06 Feb 2020
-
-version 0.27.0
-
-    add convolve, the start of a 'filter' functionality group
-
-Wed, 29 Jan 2020
-
-version 0.26.7
-
-    fixed indexing error in linalg.dot
-
-Mon, 20 Jan 2020
-
-version 0.26.6
-
-    replaced MP_ROM_PTR(&mp_const_none_obj), so that module can be compiled for the nucleo board
-
-Tue, 7 Jan 2020
-
-version 0.26.5
-
-    fixed glitch in numerical.c, numerical.h
-
-Mon, 6 Jan 2020
-
-version 0.26.4
-
-    switched version constant to string
-
-Tue, 31 Dec 2019
-
-version 0.263
-
-    changed declaration of ulab_ndarray_type to extern
-
 Fri, 29 Nov 2019
 
 version 0.262

docs/ulab-manual.ipynb

@@ -24,11 +24,11 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 14,
+   "execution_count": 23,
    "metadata": {
     "ExecuteTime": {
-     "end_time": "2020-02-11T19:06:35.427133Z",
-     "start_time": "2020-02-11T19:06:35.418598Z"
+     "end_time": "2019-11-05T16:01:28.817558Z",
+     "start_time": "2019-11-05T16:01:28.810704Z"
     }
    },
    "outputs": [
@@ -66,7 +66,7 @@
     "author = 'Zoltán Vörös'\n",
     "\n",
     "# The full version, including alpha/beta/rc tags\n",
-    "release = '0.32'\n",
+    "release = '0.26'\n",
     "\n",
     "\n",
     "# -- General configuration ---------------------------------------------------\n",
@@ -120,11 +120,10 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 16,
+   "execution_count": null,
    "metadata": {
     "ExecuteTime": {
-     "end_time": "2020-02-11T20:30:43.287360Z",
-     "start_time": "2020-02-11T20:30:40.308932Z"
+     "start_time": "2019-11-06T17:37:29.723Z"
     }
    },
    "outputs": [],
@@ -293,11 +292,11 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 17,
+   "execution_count": 3,
    "metadata": {
     "ExecuteTime": {
-     "end_time": "2020-02-16T14:53:49.098172Z",
-     "start_time": "2020-02-16T14:53:49.093201Z"
+     "end_time": "2019-11-04T20:57:36.730820Z",
+     "start_time": "2019-11-04T20:57:36.723446Z"
     }
    },
    "outputs": [],
@@ -311,11 +310,11 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 18,
+   "execution_count": 4,
    "metadata": {
     "ExecuteTime": {
-     "end_time": "2020-02-16T14:53:53.396267Z",
-     "start_time": "2020-02-16T14:53:53.375754Z"
+     "end_time": "2019-11-04T20:57:38.576560Z",
+     "start_time": "2019-11-04T20:57:38.521927Z"
     }
    },
    "outputs": [],
@@ -384,20 +383,13 @@
     "ip.register_magics(PyboardMagic)"
    ]
   },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "## pyboard"
-   ]
-  },
   {
    "cell_type": "code",
-   "execution_count": 111,
+   "execution_count": 520,
    "metadata": {
     "ExecuteTime": {
-     "end_time": "2020-02-16T18:36:59.172039Z",
-     "start_time": "2020-02-16T18:36:59.144651Z"
+     "end_time": "2019-10-20T06:48:01.610725Z",
+     "start_time": "2019-10-20T06:48:00.856261Z"
     }
    },
    "outputs": [],
@@ -409,11 +401,11 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 110,
+   "execution_count": 501,
    "metadata": {
     "ExecuteTime": {
-     "end_time": "2020-02-16T17:34:35.250747Z",
-     "start_time": "2020-02-16T17:34:35.241871Z"
+     "end_time": "2019-10-19T13:36:42.010602Z",
+     "start_time": "2019-10-19T13:36:42.003900Z"
     }
    },
    "outputs": [],
@@ -485,7 +477,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "In the [last chapter](https://micropython-usermod.readthedocs.io/en/latest/usermods_15.html) of the usermod documentation, I mentioned that I have another story, for another day. The day has come, so here is my story.\n",
+    "In https://micropython-usermod.readthedocs.io/en/latest/usermods_14.html, I mentioned that I have another story, for another day. The day has come, so here is my story.\n",
     "\n",
     "## Enter ulab\n",
     "\n",
@@ -513,7 +505,7 @@
     "- polynomial fits to numerical data\n",
     "- fast Fourier transforms\n",
     "\n",
-    "At the time of writing this manual (for version 0.32), the library adds approximately 30 kB of extra compiled code to the micropython (pyboard.v.11) firmware. However, if you are tight with flash space, you can easily shave off a couple of kB. See the section on [customising ulab](#Custom_builds).\n",
+    "At the time of writing this manual (for version 0.26), the library adds approximately 30 kB of extra compiled code to the micropython (pyboard.v.11) firmware. \n",
     "\n",
     "## Resources and legal matters\n",
     "\n",
@@ -576,42 +568,6 @@
     "If you find a bug, please, include this number in your report!"
    ]
   },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "## Customising `ulab`\n",
-    "\n",
-    "`ulab` implements a great number of functions, which are organised in sub-modules. E.g., functions related to Fourier transforms are located in the `ulab.fft` sub-module, so you would import `fft` as\n",
-    "\n",
-    "```python\n",
-    "import ulab\n",
-    "from ulab import fft\n",
-    "```\n",
-    "by which point you can get the FFT of your data by calling `fft.fft(...)`. \n",
-    "\n",
-    "The idea of such grouping of functions and methods is to provide a means for granularity: It is quite possible that you do not need all functions in a particular application. If you want to save some flash space, you can easily exclude arbitrary sub-modules from the firmware. The [https://github.com/v923z/micropython-ulab/blob/master/code/ulab.h](ulab.h) header file contains a pre-processor flag for each sub-module. The default setting is 1 for each of them, but if you change that to 0, the corresponding sub-module will not be part of the compiled firmware. \n",
-    "\n",
-    "The first couple of lines of the file look like this\n",
-    "\n",
-    "```c\n",
-    "// vectorise (all functions) takes approx. 3 kB of flash space\n",
-    "#define ULAB_VECTORISE_MODULE (1)\n",
-    "\n",
-    "// linalg adds around 6 kB\n",
-    "#define ULAB_LINALG_MODULE (1)\n",
-    "\n",
-    "// poly is approx. 2.5 kB\n",
-    "#define ULAB_POLY_MODULE (1)\n",
-    "```\n",
-    "\n",
-    "In order to simplify navigation in the header, each flag begins with `ULAB_`, and continues with the name of the sub-module. This name is also the `.c` file, where the sub-module is implemented. So, e.g., the linear algebra routines can be found in `linalg.c`, and the corresponding compiler flag is `ULAB_LINALG_MODULE`. Each section displays a hint as to how much space you can save by un-setting the flag.\n",
-    "\n",
-    "At first, having to import everything in this way might appear to be overly complicated, but there is a very good reason behind all this: you can find out at the time of importing, whether a function is part of your `ulab` firmware, or not. The alternative, namely, that you do not have to import anything beyond `ulab`, could be catastrophic: you would learn only at run time that a particular function is not in the firmware, and that is most probably too late.\n",
-    "\n",
-    "The standard sub-modules, `vector`, `linalg`, `numerical`, `poly`, and `fft` are all `numpy`-compatible. User-defined functions that accept `ndarray`s as their argument should be implemented in the `extra` sub-module, or its sub-modules."
-   ]
-  },
   {
    "cell_type": "markdown",
    "metadata": {},
@@ -637,10 +593,14 @@
     "\n",
     "[.reshape](#.reshape)\n",
     "\n",
+    "[.rawsize<sup>**</sup>](#.rawsize)\n",
+    "\n",
     "[.transpose](#.transpose)\n",
     "\n",
     "[.flatten<sup>**</sup>](#.flatten)\n",
     "\n",
+    "[.asbytearray](#.asbytearray)\n",
+    "\n",
     "## Matrix methods\n",
     "\n",
     "[size](#size)\n",
@@ -699,11 +659,7 @@
     "\n",
     "[ifft<sup>**</sup>](#ifft)\n",
     "\n",
-    "[spectrum<sup>**</sup>](#spectrum)\n",
-    "\n",
-    "## Filter functions\n",
-    "\n",
-    "[convolve](#convolve)"
+    "[spectrum<sup>**</sup>](#spectrum)"
    ]
   },
   {
@@ -712,13 +668,13 @@
    "source": [
     "# ndarray, the basic container\n",
     "\n",
-    "The `ndarray` is the underlying container of numerical data. It is derived from micropython's own `array` object, but has a great number of extra features starting with how it can be initialised, which operations can be done on it, and which functions can accept it as an argument. One important property of an `ndarray` is that it is also a proper `micropython` iterable.\n",
+    "The `ndarray` is the underlying container of numerical data. It is derived from micropython's own `array` object, but has a great number of extra features starting with how it can be initialised, how operations can be done on it, and which functions can accept it as an argument.\n",
     "\n",
-    "Since the `ndarray` is a binary container, it is also compact, meaning that it takes only a couple of bytes of extra RAM in addition to what is required for storing the numbers themselves. `ndarray`s are also type-aware, i.e., one can save RAM by specifying a data type, and using the smallest reasonable one. Five such types are defined, namely `uint8`, `int8`, which occupy a single byte of memory per datum, `uint16`, and `int16`, which occupy two bytes per datum, and `float`, which occupies four or eight bytes per datum. The precision/size of the `float` type depends on the definition of `mp_float_t`. Some platforms, e.g., the PYBD, implement `double`s, but some, e.g., the pyboard.v.11, don't. You can find out, what type of float your particular platform implements by looking at the output of the [.itemsize](#.itemsize) class property.\n",
+    "Since the `ndarray` is a binary container, it is also compact, meaning that it takes only a couple of bytes of extra RAM in addition to what is required for storing the numbers themselves. `ndarray`s are also type-aware, i.e., one can save RAM by specifying a data type, and using the smallest reasonable one. Five such types are defined, namely `uint8`, `int8`, which occupy a single byte of memory per datum, `uint16`, and `int16`, which occupy two bytes per datum, and `float`, which occupies four or eight bytes per datum. The precision/size of the `float` type depends on the definition of `mp_float_t`. Some platforms, e.g., the PYBD, implement `double`s, but some, e.g., the pyboard.v.11, don't. You can find out, what type of float your particular platform implements by looking at the output of the [.rawsize](#.rawsize) class method.\n",
     "\n",
     "On the following pages, we will see how one can work with `ndarray`s. Those familiar with `numpy` should find that the nomenclature and naming conventions of `numpy` are adhered to as closely as possible. I will point out the few differences, where necessary.\n",
     "\n",
-    "For the sake of comparison, in addition to the `ulab` code snippets, sometimes the equivalent `numpy` code is also presented. You can find out, where the snippet is supposed to run by looking at its first line, the header.\n",
+    "For the sake of comparison, in addition to `ulab` code snippets, sometimes the equivalent `numpy` code is also presented. You can find out, where the snippet is supposed to run by looking at its first line, the header.\n",
     "\n",
     "Hint:  you can easily port existing `numpy` code, if you `import ulab as np`."
    ]
@@ -886,16 +842,16 @@
    "source": [
     "### .shape\n",
     "\n",
-    "The `.shape` method (property) returns a 2-tuple with the number of rows, and columns."
+    "The `.shape` method returns a 2-tuple with the number of rows, and columns."
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 12,
+   "execution_count": 283,
    "metadata": {
     "ExecuteTime": {
-     "end_time": "2020-02-11T19:01:40.377272Z",
-     "start_time": "2020-02-11T19:01:40.364448Z"
+     "end_time": "2019-10-16T15:30:33.810628Z",
+     "start_time": "2019-10-16T15:30:33.796088Z"
     }
    },
    "outputs": [
@@ -923,111 +879,11 @@
     "\n",
     "a = np.array([1, 2, 3, 4], dtype=np.int8)\n",
     "print(\"a:\\n\", a)\n",
-    "print(\"shape of a:\", a.shape)\n",
+    "print(\"shape of a:\", a.shape())\n",
     "\n",
     "b= np.array([[1, 2], [3, 4]], dtype=np.int8)\n",
     "print(\"\\nb:\\n\", b)\n",
-    "print(\"shape of b:\", b.shape)"
-   ]
-  },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "### .size\n",
-    "\n",
-    "The `.size` method (property) returns an integer with the number of elements in the array."
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 10,
-   "metadata": {
-    "ExecuteTime": {
-     "end_time": "2020-02-11T06:32:22.721112Z",
-     "start_time": "2020-02-11T06:32:22.713111Z"
-    }
-   },
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "a:\n",
-      " array([1, 2, 3], dtype=int8)\n",
-      "size of a: 3\n",
-      "\n",
-      "b:\n",
-      " array([[1, 2],\n",
-      "\t [3, 4]], dtype=int8)\n",
-      "size of b: 4\n",
-      "\n",
-      "\n"
-     ]
-    }
-   ],
-   "source": [
-    "%%micropython -unix 1\n",
-    "\n",
-    "import ulab as np\n",
-    "\n",
-    "a = np.array([1, 2, 3], dtype=np.int8)\n",
-    "print(\"a:\\n\", a)\n",
-    "print(\"size of a:\", a.size)\n",
-    "\n",
-    "b= np.array([[1, 2], [3, 4]], dtype=np.int8)\n",
-    "print(\"\\nb:\\n\", b)\n",
-    "print(\"size of b:\", b.size)"
-   ]
-  },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "### .itemsize\n",
-    "\n",
-    "The `.itemsize` method (property) returns an integer with the siz enumber of elements in the array."
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 13,
-   "metadata": {
-    "ExecuteTime": {
-     "end_time": "2020-02-11T19:05:04.296601Z",
-     "start_time": "2020-02-11T19:05:04.280669Z"
-    }
-   },
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "a:\n",
-      " array([1, 2, 3], dtype=int8)\n",
-      "itemsize of a: 1\n",
-      "\n",
-      "b:\n",
-      " array([[1.0, 2.0],\n",
-      "\t [3.0, 4.0]], dtype=float)\n",
-      "itemsize of b: 8\n",
-      "\n",
-      "\n"
-     ]
-    }
-   ],
-   "source": [
-    "%%micropython -unix 1\n",
-    "\n",
-    "import ulab as np\n",
-    "\n",
-    "a = np.array([1, 2, 3], dtype=np.int8)\n",
-    "print(\"a:\\n\", a)\n",
-    "print(\"itemsize of a:\", a.itemsize)\n",
-    "\n",
-    "b= np.array([[1, 2], [3, 4]], dtype=np.float)\n",
-    "print(\"\\nb:\\n\", b)\n",
-    "print(\"itemsize of b:\", b.itemsize)"
+    "print(\"shape of b:\", b.shape())"
    ]
   },
   {
@@ -1078,6 +934,54 @@
     "print('a (1 by 16):', a.reshape((1, 16)))"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### .rawsize\n",
+    "\n",
+    "The `rawsize` method of the `ndarray` returns a 5-tuple with the following data\n",
+    "\n",
+    "1. number of rows\n",
+    "2. number of columns\n",
+    "3. length of the storage (should be equal to the product of 1. and 2.)\n",
+    "4. length of the data storage in bytes \n",
+    "5. datum size in bytes (1 for `uint8`/`int8`, 2 for `uint16`/`int16`, and 4, or 8 for `floats`, see [ndarray, the basic container](#ndarray,-the-basic-container))\n",
+    "\n",
+    "**WARNING:** `rawsize` is a `ulab`-only method; it has no equivalent in `numpy`."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 510,
+   "metadata": {
+    "ExecuteTime": {
+     "end_time": "2019-10-19T17:44:26.983908Z",
+     "start_time": "2019-10-19T17:44:26.764912Z"
+    }
+   },
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "a: \t\t array([1.0, 2.0, 3.0, 4.0], dtype=float)\n",
+      "rawsize of a: \t (1, 4, 4, 16, 4)\n",
+      "\n",
+      "\n"
+     ]
+    }
+   ],
+   "source": [
+    "%%micropython -unix 1\n",
+    "\n",
+    "import ulab as np\n",
+    "\n",
+    "a = np.array([1, 2, 3, 4], dtype=np.float)\n",
+    "print(\"a: \\t\\t\", a)\n",
+    "print(\"rawsize of a: \\t\", a.rawsize())"
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {},
@@ -1131,6 +1035,109 @@
     "print(\"b flattened (F): \\t\", b.flatten(order='F'))"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### .asbytearray\n",
+    "\n",
+    "The contents of an `ndarray` can be accessed directly by calling the `.asbytearray` method. This will simply return a pointer to the underlying flat `array` object, which can then be manipulated directly.\n",
+    "\n",
+    "**WARNING:** `asbytearray` is a `ulab`-only method; it has no equivalent in `numpy`.\n",
+    "\n",
+    "In the example below, note the difference between `a`, and `buffer`: while both are designated as an array, you recognise the micropython array from the fact that it prints the typecode (`b` in this particular case). The `ndarray`, on the other hand, prints out the `dtype` (`int8` here)."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 211,
+   "metadata": {
+    "ExecuteTime": {
+     "end_time": "2019-11-01T14:25:39.008074Z",
+     "start_time": "2019-11-01T14:25:38.876546Z"
+    }
+   },
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "a:  array([1, 2, 3, 4], dtype=int8)\n",
+      "array content: array('b', [1, 2, 3, 4])\n",
+      "array content: array('b', [1, 123, 3, 4])\n",
+      "\n",
+      "\n"
+     ]
+    }
+   ],
+   "source": [
+    "%%micropython -unix 1\n",
+    "\n",
+    "import ulab as np\n",
+    "\n",
+    "a = np.array([1, 2, 3, 4], dtype=np.int8)\n",
+    "print('a: ', a)\n",
+    "buffer = a.asbytearray()\n",
+    "print(\"array content:\", buffer)\n",
+    "buffer[1] = 123\n",
+    "print(\"array content:\", buffer)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "This in itself wouldn't be very interesting, but since `buffer` is a proper micropython `array`, we can pass it to functions that can employ the buffer protocol. E.g., all the `ndarray` facilities can be applied to the results of timed ADC conversions."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 563,
+   "metadata": {
+    "ExecuteTime": {
+     "end_time": "2019-10-20T07:29:00.153589Z",
+     "start_time": "2019-10-20T07:28:50.210383Z"
+    }
+   },
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "ADC results:\t array([48, 2, 2, ..., 0, 0, 0], dtype=uint8)\n",
+      "mean of results:\t 1.22\n",
+      "std of results:\t 4.744639\n",
+      "\n"
+     ]
+    }
+   ],
+   "source": [
+    "%%micropython -pyboard 1\n",
+    "\n",
+    "import pyb\n",
+    "import ulab as np\n",
+    "\n",
+    "n = 100\n",
+    "\n",
+    "adc = pyb.ADC(pyb.Pin.board.X19)\n",
+    "tim = pyb.Timer(6, freq=10)\n",
+    "\n",
+    "a = np.array([0]*n, dtype=np.uint8)\n",
+    "buffer = a.asbytearray()\n",
+    "adc.read_timed(buffer, tim)\n",
+    "\n",
+    "print(\"ADC results:\\t\", a)\n",
+    "print(\"mean of results:\\t\", np.mean(a))\n",
+    "print(\"std of results:\\t\", np.std(a))"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Likewise, data can be read directly into `ndarray`s from other interfaces, e.g., SPI, I2C etc, and also, by laying bare the `ndarray`, we can pass results of `ulab` computations to anything that can read from a buffer."
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {},
@@ -3981,11 +3988,11 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 114,
+   "execution_count": 458,
    "metadata": {
     "ExecuteTime": {
-     "end_time": "2020-02-16T18:38:07.294862Z",
-     "start_time": "2020-02-16T18:38:07.233842Z"
+     "end_time": "2019-10-19T13:07:43.168629Z",
+     "start_time": "2019-10-19T13:07:43.130341Z"
     }
    },
    "outputs": [
@@ -3993,13 +4000,13 @@
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "real part:\t array([5119.996, -5.004663, -5.004798, ..., -5.005482, -5.005643, -5.006577], dtype=float)\r\n",
-      "\r\n",
-      "imaginary part:\t array([0.0, 1631.333, 815.659, ..., -543.764, -815.6588, -1631.333], dtype=float)\r\n",
-      "\r\n",
-      "real part:\t array([5119.996, -5.004663, -5.004798, ..., -5.005482, -5.005643, -5.006577], dtype=float)\r\n",
-      "\r\n",
-      "imaginary part:\t array([0.0, 1631.333, 815.659, ..., -543.764, -815.6588, -1631.333], dtype=float)\r\n",
+      "real part:\t array([5119.996, -5.004663, -5.004798, ..., -5.005482, -5.005643, -5.006577], dtype=float)\n",
+      "\n",
+      "imaginary part:\t array([0.0, 1631.333, 815.659, ..., -543.764, -815.6588, -1631.333], dtype=float)\n",
+      "\n",
+      "real part:\t array([5119.996, -5.004663, -5.004798, ..., -5.005482, -5.005643, -5.006577], dtype=float)\n",
+      "\n",
+      "imaginary part:\t array([0.0, 1631.333, 815.659, ..., -543.764, -815.6588, -1631.333], dtype=float)\n",
       "\n"
      ]
     }
@@ -4008,20 +4015,16 @@
     "%%micropython -pyboard 1\n",
     "\n",
     "import ulab as np\n",
-    "from ulab import numerical\n",
-    "from ulab import vector\n",
-    "from ulab import fft\n",
-    "from ulab import linalg\n",
     "\n",
-    "x = numerical.linspace(0, 10, num=1024)\n",
-    "y = vector.sin(x)\n",
-    "z = linalg.zeros(len(x))\n",
+    "x = np.linspace(0, 10, num=1024)\n",
+    "y = np.sin(x)\n",
+    "z = np.zeros(len(x))\n",
     "\n",
-    "a, b = fft.fft(x)\n",
+    "a, b = np.fft(x)\n",
     "print('real part:\\t', a)\n",
     "print('\\nimaginary part:\\t', b)\n",
     "\n",
-    "c, d = fft.fft(x, z)\n",
+    "c, d = np.fft(x, z)\n",
     "print('\\nreal part:\\t', c)\n",
     "print('\\nimaginary part:\\t', d)"
    ]
@@ -4265,56 +4268,6 @@
     "where $N$ is the length of $y_1$, and $Y_1, Y_2$, and $Y$, respectively, are the Fourier transforms of $y_1, y_2$, and $y = y_1 + iy_2$."
    ]
   },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "# Filter routines"
-   ]
-  },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.convolve.html\n",
-    "\n",
-    "## convolve\n",
-    "Returns the discrete, linear convolution of two one-dimensional sequences.\n",
-    "\n",
-    "Only the ``full`` mode is supported, and the ``mode`` named parameter is not accepted. Note that all other modes can be had by slicing a ``full`` result."
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 5,
-   "metadata": {
-    "ExecuteTime": {
-     "end_time": "2020-02-10T18:46:06.538207Z",
-     "start_time": "2020-02-10T18:46:06.525851Z"
-    }
-   },
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "array([1.0, 12.0, 123.0, 1230.0, 2300.0, 3000.0], dtype=float)\n",
-      "\n",
-      "\n"
-     ]
-    }
-   ],
-   "source": [
-    "%%micropython -unix 1\n",
-    "\n",
-    "import ulab as np\n",
-    "\n",
-    "x = np.array((1,2,3))\n",
-    "y = np.array((1,10,100,1000))\n",
-    "\n",
-    "print(np.convolve(x, y))"
-   ]
-  },
   {
    "cell_type": "markdown",
    "metadata": {},
@@ -4360,7 +4313,7 @@
     "    NDARRAY_FLOAT = 'd',\n",
     "};\n",
     "```\n",
-    "The ambiguity is caused by the fact that not all platforms implement `double`, and there one has to take `float`s. But you haven't actually got to be concerned by this, because at the very beginning of `ndarray.h`, this is already taken care of: the pre-processor figures out, what the `float` implementation of the hardware platform is, and defines the  `NDARRAY_FLOAT` typecode accordingly. All you have to keep in mind is that wherever you would use `float` or `double`, you have to use `mp_float_t`. That type is defined in `py/mpconfig.h` of the micropython code base.\n",
+    "The ambiguity is caused by the fact that not all platforms implement `double`, and there one has to take `float`s. But you haven't actually got to be concerned by this, because at the very beginning of `ndarray.h`, this is already taken care of: the preprocessor figures out, what the `float` implementation of the hardware platform is, and defines the  `NDARRAY_FLOAT` typecode accordingly. All you have to keep in mind is that wherever you would use `float` or `double`, you have to use `mp_float_t`. That type is defined in `py/mpconfig.h` of the micropython code base.\n",
     "\n",
     "Therefore, a 4-by-5 matrix of type float can be created as\n",
     "\n",

tests/00smoke.py

@@ -1,2 +0,0 @@
-from ulab import linalg
-print(linalg.eye(3))

tests/00smoke.py.exp

@@ -1,3 +0,0 @@
-array([[1.0, 0.0, 0.0],
-	 [0.0, 1.0, 0.0],
-	 [0.0, 0.0, 1.0]], dtype=float)