Handle different versions of BLAS and LAPACK

Description

Beside being implementations, BLAS and LAPACK are also API standard for basic linear algebra operations (such as vector and matrix multiplication).

Many implementations of these API exist. The reference implementation of BLAS and LAPACK is very stable but is not as fast as optimized ones such as ATLAS and OpenBLAS.

Implementations of BLAS (both the Fortran and the C interface):

Implementations of LAPACK (Fortran interface):

How to switch from an implementation to the other

It is just trivial. It can be used like any other software using update-alternatives. If you don't like update-alternatives' CLI, you can use any of it's front-end programs instead, e.g. rover (text-based), galternatives (GTK), kalternatives (Qt).

BLAS

In jessie and stretch

update-alternatives --config libblas.so.3

In (upcoming) buster and unstable/sid

update-alternatives --config libblas.so.3-<multiarch>

where <multiarch> is the multiarch path for you architecture (e.g. x86_64-linux-gnu for amd64).

Example

There are 3 choices for the alternative libblas.so.3 (providing /usr/lib/libblas.so.3).

  Selection    Path                                     Priority   Status
------------------------------------------------------------
* 0            /usr/lib/openblas-base/libopenblas.so.0   40        auto mode
  1            /usr/lib/atlas-base/atlas/libblas.so.3    35        manual mode
  2            /usr/lib/libblas/libblas.so.3             10        manual mode
  3            /usr/lib/openblas-base/libopenblas.so.0   40        manual mode

/!\ libgslcblas0 is not in that list, because it does not implement the Fortran interface.

LAPACK

In jessie and stretch

update-alternatives --config liblapack.so.3

In (upcoming) buster and unstable/sid

update-alternatives --config liblapack.so.3-<multiarch>

where <multiarch> is the multiarch path for you architecture (e.g. x86_64-linux-gnu for amd64).

Example

There are 3 choices for the alternative liblapack.so.3 (providing /usr/lib/liblapack.so.3).

  Selection    Path                                      Priority   Status
------------------------------------------------------------
* 0            /usr/lib/openblas-base/liblapack.so.3      40        auto mode
  1            /usr/lib/atlas-base/atlas/liblapack.so.3   35        manual mode
  2            /usr/lib/lapack/liblapack.so.3             10        manual mode
  3            /usr/lib/openblas-base/liblapack.so.3      40        manual mode

Getting the best performance out of ATLAS and OpenBLAS

The binary packages of ATLAS and OpenBLAS distributed by Debian are generic packages, which are not optimized for your specific machine.

For ATLAS, the optimal way of using it is to recompile it locally on the machine on which it is to be used. See the README.Debian file for details on how to achieve that easily.

For OpenBLAS, if you are on amd64 or i386, there is no need to recompile it, since the binary includes kernels optimized for several CPU microarchitectures, and the selection is done at runtime. For non-x86 architectures however, OpenBLAS should be recompiled locally for optimal performance, see the README.Debian file for instructions.

Using the GSL with an optimized BLAS

The GSL comes with its own BLAS implementation (libgslcblas0), which is not optimized for speed.

It is however possible to use an optimized BLAS implementation (OpenBLAS or ATLAS) as the backend for GSL matrix linear algebra when compiling your program.

Provided that you have installed libblas-dev, libatlas-base-dev or libopenblas-dev, you can use -lgsl -lblas at link time (instead of -lgsl -lgslcblas); or, if you are using pkg-config, you can get the right libraries with:

pkg-config --define-variable=GSL_CBLAS_LIB=-lblas --libs gsl

At runtime, the resulting binary will use the BLAS implementation selected by the alternatives system described above.

(WIP, experimental) BLAS64 / LAPACK64, aka. BLAS/LAPACK with 64-bit indexing

This is an experimental on-going work.

Background

Netlib BLAS (Fortran) uses INTEGER as the type of index.

SUBROUTINE SCOPY(N,SX,INCX,SY,INCY)
INTEGER INCX,INCY,N

By default, INTEGER is int32_t, but when the index of your vector of matrix goes beyond the range of int32_t, one can use the -i8 fortran compiler option to change INTEGER into int64_t.

intel-mkl calls the 32-bit index version as LP64 interface (long and pointer are 64bit), and the 64-bit index version as ILP64 interface (int, long and pointer are 64bit). However this term may cause confusion in the C world.

We call BLAS and LAPACK interfaces that support 64-bit index as "BLAS64" and "LAPACK64" respectively.

Packages that provide BLAS64 API/ABI

Packages that provide LAPACK64 API/ABI

https://lists.debian.org/debian-devel/2018/10/msg00339.html