detection of the ILP mode of BLAS

Question

the blas mathematical library is distributed either in the i32lp64 mode (that is, with integer*4 integers) or in the ilp64 mode (with integer*8 or 8-bytes integers).

The question is how to distinguish between these two BLAS modes (i32lp64 vs ilp64) in a short Fortran routine, and without giving the segfault crash.

Answer 1

You can perform a test at runtime by being rather sneaky; it requires understanding the calling conventions on x86_64 and determining what the values passed to certain functions will look like when interpreted as either 32-bit or 64-bit numbers. libblastrampoline does this whenever it loads a new BLAS library in order to determine whether a library is ILP64 or just LP64.

The detection code for BLAS is here, and is well-commented. Summarizing it here, we will call a known BLAS function isamax(blas_int n, float * X, blas_int incx) with specially-crafted arguments:

int64_t n = 0xffffffff00000003;
float X[3] = {1.0f, 2.0f, 1.0f};
int64_t incx = 1;

The most important here is n, which appears to be a large, negative number when viewed as a 32-bit number, but appears to be a positive number when viewed as a 64-bit number. Here we are using the FORTRAN-style API but from C, so as FORTRAN expects scalars to be passed as pointers we'll actually do something like:

int64_t max_idx = isamax(&n, X, &incx);

Now when an LP64 (e.g. 32-bit) BLAS library dereferences n, because this is a little-endian system, it will see the lower 32 bits of that value, which is just 3. This will then return a value of 2 in max_idx, which is the correct return value. However, if N is seen to be negative, this function returns 0. This error handling behavior allows you to disambiguate how the library is interpreting your arguments.

LBT contains a similar method for looking at an LAPACK library as well.

The programming language Julia uses LBT internally, and it is quite easy to load multiple BLAS/LAPACK libraries and learn all sorts of details about their internal ABIs:

julia> using LinearAlgebra, OpenBLAS_jll
       LinearAlgebra.BLAS.lbt_forward(OpenBLAS_jll.libopenblas_path; clear=true, verbose=true);
Generating forwards to /Users/sabae/.julia/juliaup/julia-1.9.2+0.aarch64.apple.darwin14/lib/julia/libopenblas64_.dylib (clear: 1, verbose: 1, suffix_hint: '(null)')
 -> Autodetected symbol suffix "64_"
 -> Autodetected interface ILP64 (64-bit)
 -> Autodetected normal complex return style
 -> Autodetected gfortran calling convention
Processed 4949 symbols; forwarded 4856 symbols with 64-bit interface and mangling to a suffix of "64_"

So you can see it has dynamically determined the symbol suffix (e.g. dgemm_ is actually called dgemm_64_ in this library), the interface type (ILP64), the complex number return style (some libraries/compilers return complex values in an implicit argument), etc... Compare with a 32-bit build of OpenBLAS:

julia> using LinearAlgebra, OpenBLAS32_jll
       LinearAlgebra.BLAS.lbt_forward(OpenBLAS32_jll.libopenblas_path; clear=true, verbose=true);
Generating forwards to /Users/sabae/.julia/artifacts/1f789e077dba9b9c9f09f7ef54f80239bbdf04f3/lib/libopenblas.0.3.21.dylib (clear: 1, verbose: 1, suffix_hint: '(null)')
 -> Autodetected symbol suffix ""
 -> Autodetected interface LP64 (32-bit)
 -> Autodetected normal complex return style
 -> Autodetected gfortran calling convention
Processed 4949 symbols; forwarded 4856 symbols with 32-bit interface and mangling to a suffix of ""