[MRG] Use Scipy cython BLAS API instead of bundled CBLAS #12732

jeremiedbb · 2018-12-06T16:15:33Z

First step to tackle #11638

This PR add a new module in utils, _cython_blas, which contains helpers for the BLAS functions shipped with scipy in scipy.cython_blas.

The scipy functions expect fortran aligned arrays. These helpers allow to specify the memory layout as it is in CBLAS, and can be used as drop in replacement in sklearn.
BLAS functions are not fused (e.g. sgemm and dgemm for float and double matrix matrix multiplication). These helpers are fused, sgemv and dgemv are replaced by _xgemv, which avoids the redundant

if floating is float:
    dot = sdot
else:
    dot = ddot

I still have to fill the docstrings.

For now the module only contains helpers for the BLAS functions which were already used in sklearn (+gemm because I'm using it in another PR :) ). We can add other BLAS functions as we need them.

This module comes with a test suite, which tests the helpers for all type/layout/transpose configurations.

I also added an example of use in sklearn pairwise_fast with the asum function, see sklearn/metrics/pairwise_fast.pyx and sklearn/metrics/setup.py.

I don't think I'll replace all occurences of CBLAS use in sklearn in this PR. I think it would be easier to do it in separate PRs. Removing the bundled CBLAS will only be possible once all replacement are done.

Finally, this PR requires an upgrade of scipy minimal version > 0.16, which seems to be on it's way (#12184)

jeremiedbb · 2018-12-06T16:46:09Z

Note that all the _xxxx_memview functions in the module are just wrappers to be able to tests the module with pytest.

jeremiedbb · 2018-12-06T17:25:22Z

CI is failing because of cython version :(
cython does not support the const keyword for memoryview before version 0.28. Would it be considered to upgrade cython requirement ?

rth · 2018-12-06T19:38:21Z

cython does not support the const keyword for memoryview before version 0.28. Would it be considered to upgrade cython requirement ?

Yes, we did discuss this in #10624 and I think it should be possible. Are you sure this is not affected by the issues from #10624 (comment) though?

jeremiedbb · 2018-12-06T19:41:19Z

Not I don't use const with fused typed memview. But after all it's not necessary, I found a workaround :)

rth · 2018-12-06T19:42:49Z

Maybe it could be worth checking if some discussion about such wrappers did take place at scipy and if not open an issue about it? From a user perspective, we might indeed want 1 cython function for matrix-matrix multiplication with BLAS irrespective of the C/Fortran alignment or dtype, not 4...

jeremiedbb · 2018-12-10T10:18:28Z

There is an old open issue in scipy scipy/scipy#4516 for that. It's been active again this summer, so there may be some progress in that way.

jeremiedbb · 2018-12-17T15:48:06Z

I had to increase rtol up to 1e-3 in tests with float32 to make CI green ! (MKL or OpenBLAS show similar behavior) whereas rtol = 1e-12 is enough for float64.

This is a really high rtol...

jeremiedbb · 2018-12-18T15:01:51Z

So I tested this rtol thing with the bundled CBLAS and it turns out that the issue is already here.
gemv(A,x,y) and A.dot(x) + y are close only up to 1e-3 on float32 data.

However, to mitigate this, I used data such that the result of gemv (or others) was often close to zero, meaning that it was involving differences between close numbers, which can easily lead to numerical precision issues.

Using better conditioned data allows to retrieve a 1e-6 rtol which is fine for float32.

sklearn/utils/_cython_blas.pyx

jeremiedbb · 2018-12-20T10:47:09Z

Since minimal versions requirements have been upgraded in #12746, it's now safe to use scipy cython blas (scipy >= 0.16).

Does that need a what's new entry ?

ogrisel · 2018-12-22T11:25:52Z

Does that need a what's new entry?

I think we should do what's new entries for the individual estimators when there implementation are updated as a consequence of this infrastructure change.

For this PR you can already add an entry for the pairwise_metrics sparse Manhattan distance computation.

ogrisel

LGTM. Thanks @jeremiedbb.

jnothman · 2019-01-29T11:15:27Z

sklearn/utils/_cython_blas.pyx

+        return ddot(&n, x, &incx, y, &incy)
+
+
+cpdef _dot_memview(floating[::1] x, floating[::1] y):


Should we be using const memoryviews to allow read-only input arrays?

fused typed const memoryviews does not work yet, see #10624
However, all the _xxx_memview functions are just python wrappers to be able to test the C functions with pytest. They are not meant to be used in the python code base (if we want to multiply matrices in python we just do numpy dot), we don't want to expose blas functions at the python level.

Ah of course. But is there a reason we should not be using memview interfaces when that would simplify the call?

There might be the small overhead of an additional function call (and I don't really know how memview behave versus pointers performance wise) but I agree it would simplify it.

Currently all blas functions are called within functions where we have access to the pointers, so I'm not sure it's worth making interfaces for what we don't need currently. Maybe we should reconsider doing it when the need comes ?

I don't think there's any function call overhead. Certainly not python functions. There will be cost in accessing members of the memoryview struct, but minimal.

No hurry, you're right, but I still find it strange that we are passed order rather than determining it from the memoryview strides.

You convinced me :) I updated the functions to infer the memory layout.

jnothman · 2019-01-29T11:23:27Z

sklearn/utils/_cython_blas.pyx

+            dgemv(&ta_, &m, &n, &alpha, A, &lda, x, &incx, &beta, y, &incy)
+
+
+cpdef _gemv_memview(BLAS_Order order, BLAS_Trans ta, floating alpha,


shouldn't we determine order from A.strides[0] == A.itemsize rather than having it passed in?

ogrisel · 2019-02-01T16:39:15Z

Merged! Thanks @jeremiedbb!

…n#12732)

…kit-learn#12732)" This reverts commit eb9a022.

…n#12732)

jeremiedbb force-pushed the scipy-blas branch 2 times, most recently from be33b3c to 486cf22 Compare December 6, 2018 16:29

jeremiedbb force-pushed the scipy-blas branch 2 times, most recently from fa78d95 to 2aa4bcb Compare December 14, 2018 14:08

ogrisel reviewed Dec 19, 2018

View reviewed changes

sklearn/utils/_cython_blas.pyx Outdated Show resolved Hide resolved

jeremiedbb changed the title ~~[WIP] Use Scipy cython BLAS API instead of bundled CBLAS~~ [MRG] Use Scipy cython BLAS API instead of bundled CBLAS Dec 20, 2018

jeremiedbb added 12 commits January 14, 2019 14:49

Scipy cython_blas fused helpers

c6cd0a9

cython language_level=3

cc81002

rtol

d823f7c

rtol

5e0656b

rtol

73522a4

alpha beta

d656a3c

enum order & trans

24e0dd0

fix numpy order type

df86d94

change blas functions names

9372276

flake8

1868648

clean up

12b3660

what's new

86556d9

jeremiedbb force-pushed the scipy-blas branch from 461abde to 86556d9 Compare January 14, 2019 14:11

ogrisel approved these changes Jan 29, 2019

View reviewed changes

ogrisel added the Waiting for Reviewer label Jan 29, 2019

jnothman reviewed Jan 29, 2019

View reviewed changes

infer memory layout

0019d25

jnothman approved these changes Jan 31, 2019

View reviewed changes

remove blank line

3564ffa

ogrisel merged commit d0f63a7 into scikit-learn:master Feb 1, 2019

jeremiedbb mentioned this pull request Feb 3, 2019

[MRG] MAINT: Continue moving from CBLAS to scipy cython blas #13084

Merged

thomasjpfan pushed a commit to thomasjpfan/scikit-learn that referenced this pull request Feb 6, 2019

[MRG] Use Scipy cython BLAS API instead of bundled CBLAS (scikit-lear…

f32f6b2

…n#12732)

thomasjpfan pushed a commit to thomasjpfan/scikit-learn that referenced this pull request Feb 7, 2019

[MRG] Use Scipy cython BLAS API instead of bundled CBLAS (scikit-lear…

5014599

…n#12732)

xhluca pushed a commit to xhluca/scikit-learn that referenced this pull request Apr 28, 2019

[MRG] Use Scipy cython BLAS API instead of bundled CBLAS (scikit-lear…

eb9a022

…n#12732)

xhluca pushed a commit to xhluca/scikit-learn that referenced this pull request Apr 28, 2019

Revert "[MRG] Use Scipy cython BLAS API instead of bundled CBLAS (sci…

f366953

…kit-learn#12732)" This reverts commit eb9a022.

xhluca pushed a commit to xhluca/scikit-learn that referenced this pull request Apr 28, 2019

Revert "[MRG] Use Scipy cython BLAS API instead of bundled CBLAS (sci…

a8b1634

…kit-learn#12732)" This reverts commit eb9a022.

jackmitch pushed a commit to jackmitch/scikit-learn that referenced this pull request Jul 2, 2019

[MRG] Use Scipy cython BLAS API instead of bundled CBLAS (scikit-lear…

4cd1358

…n#12732)

koenvandevelde pushed a commit to koenvandevelde/scikit-learn that referenced this pull request Jul 12, 2019

[MRG] Use Scipy cython BLAS API instead of bundled CBLAS (scikit-lear…

bbf1454

…n#12732)

oleksandr-pavlyk mentioned this pull request Feb 14, 2020

[WIP] Make SVC tests independent of SV ordering #12849

Open

jeremiedbb deleted the scipy-blas branch July 20, 2020 14:59

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

[MRG] Use Scipy cython BLAS API instead of bundled CBLAS #12732

[MRG] Use Scipy cython BLAS API instead of bundled CBLAS #12732

jeremiedbb commented Dec 6, 2018 •

edited by ogrisel

Loading

jeremiedbb commented Dec 6, 2018

jeremiedbb commented Dec 6, 2018

rth commented Dec 6, 2018

jeremiedbb commented Dec 6, 2018

rth commented Dec 6, 2018

jeremiedbb commented Dec 10, 2018

jeremiedbb commented Dec 17, 2018 •

edited

Loading

jeremiedbb commented Dec 18, 2018

jeremiedbb commented Dec 20, 2018

ogrisel commented Dec 22, 2018 •

edited

Loading

ogrisel left a comment

jnothman Jan 29, 2019

jeremiedbb Jan 29, 2019

jnothman Jan 29, 2019

jeremiedbb Jan 29, 2019

jnothman Jan 29, 2019

jnothman Jan 29, 2019

jeremiedbb Jan 31, 2019

jnothman Jan 29, 2019

ogrisel commented Feb 1, 2019

		return ddot(&n, x, &incx, y, &incy)


		cpdef _dot_memview(floating[::1] x, floating[::1] y):

		dgemv(&ta_, &m, &n, &alpha, A, &lda, x, &incx, &beta, y, &incy)


		cpdef _gemv_memview(BLAS_Order order, BLAS_Trans ta, floating alpha,

[MRG] Use Scipy cython BLAS API instead of bundled CBLAS #12732

[MRG] Use Scipy cython BLAS API instead of bundled CBLAS #12732

Conversation

jeremiedbb commented Dec 6, 2018 • edited by ogrisel Loading

jeremiedbb commented Dec 6, 2018

jeremiedbb commented Dec 6, 2018

rth commented Dec 6, 2018

jeremiedbb commented Dec 6, 2018

rth commented Dec 6, 2018

jeremiedbb commented Dec 10, 2018

jeremiedbb commented Dec 17, 2018 • edited Loading

jeremiedbb commented Dec 18, 2018

jeremiedbb commented Dec 20, 2018

ogrisel commented Dec 22, 2018 • edited Loading

ogrisel left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

ogrisel commented Feb 1, 2019

jeremiedbb commented Dec 6, 2018 •

edited by ogrisel

Loading

jeremiedbb commented Dec 17, 2018 •

edited

Loading

ogrisel commented Dec 22, 2018 •

edited

Loading