Add scipy integration

[h-milz]

As written in feature request #698.

[v923z]

I think this is good to go, the only remaining question is, whether you want to add tests. (If nothing else, we should at least make certain that the arguments and keyword arguments are parsed correctly.) Otherwise, it doesn't make sense to change the float to double.

[h-milz]

Yep, good point. I'll make something up tomorrow and create a new PR I suppose. Btw I found 2 or 3 errors in the docs, these would be in the follow-up PR as well.

[h-milz]

OK, I overhauled the docs, added complex number error handling, and added test cases. As far as these, I am not sure how to handle situations where not all integration algos are built. For example on a fp32 system quadgk will not be built. I can check if i.quadgk is available and just not run the corresponding tests, but the .exp file will contain the expected results just the same. How do you handle conditional compiles?

[v923z]

OK, I overhauled the docs, added complex number error handling, and added test cases. As far as these, I am not sure how to handle situations where not all integration algos are built. For example on a fp32 system quadgk will not be built. I can check if i.quadgk is available and just not run the corresponding tests, but the .exp file will contain the expected results just the same. How do you handle conditional compiles?

Don't worry about that! The unix port should cover all cases, and that's where the tests are going to run.

[v923z]

Can you, please, address the in-line review comments above?

[h-milz]

Are you referring to the commit 6e878c5? I changed a couple of things after my first PR but did not immediately create a new one. The changes are contained in 0a6ad6d, though.

[v923z]

No, I mean the review comments in this thread. There are some that are labelled outdated, though, I don't see that the comment itself was addressed, and as far as I understand, the changes didn't concerned the issue raised. And then are quite a few that are not even labelled outdated, and there was no response to them.

0a6ad6d deals with the documentation, but not the code.

[h-milz]

Sorry Zoltan for sounding dumb but I don't get it. Which inline comments are you referring to? I see none in the diff views, and in the top right of this page I see "No reviews ". Please advise!

[v923z]

Sorry Zoltan for sounding dumb but I don't get it. Which inline comments are you referring to? I see none in the diff views, and in the top right of this page I see "No reviews ". Please advise!

Oh, don't worry, something might very well be just wrong. Do you see the comments after this post? #699 (comment)

You should have a bunch of shorter remarks, like these

Some of those are obsolete, but I believe, some are still open questions.

[h-milz]

The if 0's are out. I had replaced them in my latest commit f446db0 by a check if the value is < 1 and if yes, throw the value error. The if 0's were an intermediate hack because the original checks for negative values were complained by the compiler who rightfully said that the conditions were never met because the variable is a uint16.

[v923z]

I don't think that f446db0 addresses those questions.

Also, can you, please, reply to the comments in place? It's really hard to keep track of what has been fixed, and what needs additional attention. There are some issues that really have to be sorted out. E.g., we can't bring in a feature that depends on the precision, etc.

[h-milz]

Zoltan, I am truly sorry but I see zero comments anywhere, and I also received no notifications from the platform (although I get everything else from github from other activities). In particular, Notifications in this thread are set to "subscribed" ("You’re receiving notifications because you authored the thread. "), and my other notification settings are default. I have no idea what open questions there might be, but I am more than willing to resolve any. At the moment I see two commits with a red cross (6e878c5 and f446db0), and if I click them they are free of any comments although I believe this is exactly where I should see them. Please see the screenshots. I also see no review comments here in this thread, as I know them from elsewhere (and in other pull requests in this repo). There is also no "Changes requested" box as outlined in the github documentation (yes, I reviewed the doc to make sure I'm not completely off track). In "Files changed", there is no open conversation (see screenshot below). I also checked that is not a Firefox problem. It's the same in Brave, and the github Android app tells me "Reviews: None requested".

I'm pretty desperated. Please advise!

The list of commits. No conversation bubbles.

The "Files changed" tab. No conversations.

[v923z]

I'm pretty desperated. Please advise!

Well, my advice would be not to despair, because we can solve this. I will copy the comments into this thread, and link the relevant lines of the code.

[v923z]

On the unix port, the float implementation defaults to double, so you don't actually need this. I would really like to take this out, because here we override a previous definition for no reason.
https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/build.sh#L62

[v923z]

If the float implementation is single-precision, then the tolerance must surely be different, right?

https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/code/scipy/integrate/integrate.c#L42

[v923z]

This constant is hard-coded, and this might fail with single precision.

https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/code/scipy/integrate/integrate.c#L568

[v923z]

Use math.isclose here, so that the result is independent of the precision!

https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/tests/2d/scipy/integrate.py#L26

[v923z]

You really don't need this. The tests are for the most general case, when all modules, methods and constants are compiled into the firmware, so it's a safe assumption that the method that you want to test is available. There is no point in not keeping the test scripts lean.

https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/tests/2d/scipy/integrate.py#L8-L10

[v923z]

Instead of letting pass the exception, you could just do import scipy. See here: https://github.com/v923z/micropython-ulab/?tab=readme-ov-file#usage

https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/tests/2d/scipy/integrate.py#L13-L19

[v923z]

Uh, this is a bit misleading. You could just do from ulab.scipy import integrate, right? Or is there something that I overlooked?

https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/tests/2d/scipy/integrate.py#L12

[h-milz]

Strange. Of my last 4 commits only 2 seem to have made it into the PR.

[v923z]

I don't think you need an #if #endif in the makefile: that should be taken care of during the linking.

[v923z]

Use math.isclose here, so that the result is independent of the precision!
https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/tests/2d/scipy/integrate.py#L26

As before. Won't fix for 32 bit for now.

The tests don't run for 32 bits, anyway. However, you still have to do this with math.isclose, because here you're effectively comparing floating point numbers.

[v923z]

This will work on the Unix build as well as on platforms that support fp64 just fine, like the ESP32 with PSRAM (which is my focus at the moment). Would that be something you could live with?

Otherwise I would have to change quite a few things in constants handling. We could do this in a later step, though.

I get that, but we can't include a feature that works only with certain settings. I don't quite see why it's so complicated to allow single-precision floats. I get that the constants have to be defined differently, but this should work irrespective of what the precision is. Something similar happens in https://github.com/v923z/micropython-ulab/blob/master/code/numpy/linalg/linalg_tools.h, so these two cases just have to be treated separately.

[v923z]

This constant is hard-coded, and this might fail with single precision.
https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/code/scipy/integrate/integrate.c#L568

Correct. If we restrict the build to 64 bit platforms, won't fix for now.

This is actually not safe: https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/code/ulab.h#L420. The user won't necessarily know what the precision is, and they will be wondering, why that particular method is missing.

[v923z]

Instead of letting pass the exception, you could just do import scipy. See here: https://github.com/v923z/micropython-ulab/?tab=readme-ov-file#usage
https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/tests/2d/scipy/integrate.py#L13-L19

resolved.

I think you took out the wrong try - except clause: https://github.com/h-milz/micropython-ulab/blob/af290c155b8e6fb3dbc7c1843228a92f23e8e1d1/tests/2d/scipy/integrate.py#L5. You have to do that on the import (otherwise, you won't be able to run the code in numpy.scipy), but not on the actual function call.

[v923z]

Compilation fails here: https://github.com/v923z/micropython-ulab/actions/runs/12313775838/job/34392672220?pr=699#step:8:710

[h-milz]

Use math.isclose here, so that the result is independent of the precision!
https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/tests/2d/scipy/integrate.py#L26

As before. Won't fix for 32 bit for now.

The tests don't run for 32 bits, anyway. However, you still have to do this with math.isclose, because here you're effectively comparing floating point numbers.

OK, this is in now. Please check.

[h-milz]

This will work on the Unix build as well as on platforms that support fp64 just fine, like the ESP32 with PSRAM (which is my focus at the moment). Would that be something you could live with?
Otherwise I would have to change quite a few things in constants handling. We could do this in a later step, though.

I get that, but we can't include a feature that works only with certain settings. I don't quite see why it's so complicated to allow single-precision floats. I get that the constants have to be defined differently, but this should work irrespective of what the precision is. Something similar happens in https://github.com/v923z/micropython-ulab/blob/master/code/numpy/linalg/linalg_tools.h, so these two cases just have to be treated separately.

This is now resolved as well. The code builds cleanly e.g. for fp32 nanbox and appears to work properly:

MicroPython v1.25.0-preview.114.gbdda91fe7 on 2024-12-14; linux [GCC 13.3.0] version
Use Ctrl-D to exit, Ctrl-E for paste mode
>>> import math as m
>>> m.pi
3.141593
>>> from ulab import scipy
>>> f = lambda x: x**2 +2*x +4
>>> scipy.integrate.quad(f, 0, 2)
(14.66669, 2.048229e-06)

[h-milz]

Compilation fails here: https://github.com/v923z/micropython-ulab/actions/runs/12313775838/job/34392672220?pr=699#step:8:710

Seems resolved - I did not change anything, my build worked cleanly. Strange.

[h-milz]

Instead of letting pass the exception, you could just do import scipy. See here: https://github.com/v923z/micropython-ulab/?tab=readme-ov-file#usage
https://github.com/h-milz/micropython-ulab/blob/f446db0ea7925b82d4948b089ae52f2eb3392c4e/tests/2d/scipy/integrate.py#L13-L19

resolved.

I think you took out the wrong try - except clause: https://github.com/h-milz/micropython-ulab/blob/af290c155b8e6fb3dbc7c1843228a92f23e8e1d1/tests/2d/scipy/integrate.py#L5. You have to do that on the import (otherwise, you won't be able to run the code in numpy.scipy), but not on the actual function call.

I get what you mean. On the other hand, this test is less than likely to work on CPython scipy anyway because the function names are different, and they mean something different. For example, scipy.integrate has tanhsinh as a .py module and Gauss-Kronrod is found in _quadpack.c . Romberg may be similar, and while scipy.integrate has a simpson method, I use adaptive simpson which tunes the integration steps according to the integrand's behaviour.

Bummer. So the approximate mapping would be

ulab.scipy.integrate         CPython scipy.integrate
quad                         tanhsinh
romberg                      romb
simpson                      simpson    (not quite, I use adaptive)
quadgk                       no direct match, it's hidden in the old Fortran derived Quadpack code, and quadpack.py is marked deprecated. 

In order for the tests to work with CPython, I would have to rename the methods, and make sure they have the same user interface.

[h-milz]

I don't think you need an #if #endif in the makefile: that should be taken care of during the linking.

Can you specify please? Which makefile are you referring to?

[v923z]

In order for the tests to work with CPython, I would have to rename the methods, and make sure they have the same user interface.

Well, the title of this very PR is "Add scipy integration", so I assumed that the methods' names are compliant. I really don't understand what the intention with this was, if it would so spectacularly break compatibility. Any function or method that is derived from/inspired by numpy/scipy should work with the syntax given here: https://github.com/v923z/micropython-ulab/?tab=readme-ov-file#usage. There are at least two reasons for that: one is simply compatibility at the user level (people develop code on a computer and then want to run it on a microcontroller), and the second is compatibility in testing. All tests should pass with numpy/scipy, as well as ulab, with the above-mentioned header in the test script.

Using scripts that can be run in CPython and micropython is the only way to ensure that the tests are correct. This is one more reason to use math.isclose, whenever floats are involved: there is no guarantee that you're going to get the exact same results in CPython and micropython, but a test is OK, if its results are close enough.

[v923z]

I don't think you need an #if #endif in the makefile: that should be taken care of during the linking.

Can you specify please? Which makefile are you referring to?

This one here: https://github.com/h-milz/micropython-ulab/blob/cb7166385451e2128702a8526087193e5139c69d/code/micropython.mk#L5

You should just add the path to the source list, and then the file will be compiled, but if the functions are not used, then they won't be linked, so they won't take flash space.

[v923z]

Compilation fails here: https://github.com/v923z/micropython-ulab/actions/runs/12313775838/job/34392672220?pr=699#step:8:710

Seems resolved - I did not change anything, my build worked cleanly. Strange.

Yes, the compilation succeeds.

At this point, the only outstanding issues are the makefile, and the function names/compatibility with scipy. When you change those, you'll have to update the documentation, too.

[h-milz]

I think it is now as close to CPython's scipy.integrate as feasible. In fact, we're faster than that because tanhsinh will only be available in v1.15.0. :-)

I also hope that the expectation management now works better.

The test when run with python3 now says "AttributeError: module 'scipy.integrate' has no attribute 'tanhsinh'" ... my Ubuntu 24.04LTS still has scipy 1.11.something.

[v923z]

I think it is now as close to CPython's scipy.integrate as feasible. In fact, we're faster than that because tanhsinh will only be available in v1.15.0. :-)

I also hope that the expectation management now works better.

Yes, this is OK now. Thanks for sticking with it till the very end!

[h-milz]

My pleasure :-)

Btw it's not in the README.md yet. ;-)