Problem calling LAPACKE from Zig

Help
,
1

Hi,
I try to call lapacke (the C version, or interface, of lapack) to solve a set of linear equations. Here is a little test code:

const std = @import("std");

const lapack = @cImport(@cInclude("lapacke.h"));
const laint = lapack.lapack_int;

const n: laint = 4;
const nrhs: laint = 1;
const lda: laint = n;
const ldb: laint = n;

const matrix_t = [n * n]f32;
const vector_t = [n]f32;

pub fn main() !void {
    var mat: matrix_t = matrix_t{ //
        1.80, 2.88, 2.05, -0.89, //
        5.25, -2.95, -0.95, -3.80, //
        1.58, -2.69, -2.90, -1.04, //
        -1.11, -0.66, -0.59, 0.80, //
    };  
    var vec: vector_t = vector_t{ 9.52, 24.35, 0.77, -6.22 };
    std.debug.print("{any}\n", .{mat});
    std.debug.print("{any}\n", .{vec});

    var ipiv = [n]laint{ 0, 0, 0, 0 };

    const info = lapack.LAPACKE_sgesv( //
        lapack.LAPACK_COL_MAJOR, // matrix orientation
        n, // number of equations in matrix
        nrhs, // number of colums in result vector
        &mat[0], // the matrix to be diagonalized
        lda, // leading dimension of matrix
        &ipiv[0], // pivot indices
        &vec[0], // result vector
        ldb); // leading dimension of result vector

    std.debug.print("{d}\n", .{info});
    std.debug.print("{any}\n", .{mat});
    std.debug.print("{any}\n", .{ipiv});
    std.debug.print("{any}\n", .{vec});
}

When I build it I get the following error:

test.zig:51:9: error: expected type '[*c]f32', found '*f64'
        &mat[0], // the matrix to be diagonalized
        ^~~~~~~
test.zig:51:9: note: pointer type child 'f64' cannot cast into pointer type child 'f32'
/home/fuhl/zig/.zig-cache/o/138acfa2088c11699e6e581cbbe34bc8/cimport.zig:4557:74: note: parameter type declared here
pub extern fn LAPACKE_sgesv(matrix_layout: c_int, n: i32, nrhs: i32, a: [*c]f32, lda: i32, ipiv: [*c]i32, b: [*c]f32, ldb: i32) i32;

This is weird to me. I specify matrix to be [n*n]f32 but it complains that I try to pass f64?
When I change all f32 types to f64 types and change sgesv to dgesv it successfully compiles. However, when running it I get wrong results:
My solution in array vec is:

-1.3674431730055523e0
-9.779755425533203e0
 1.0529064423870595e1
-4.206233505021961e1

But It should be:

 1.0000
-1.0000
 3.0000
-5.0000

The return value of dgesv is 0 which indicates success. I guess I pass something in the wrong way.

Thanks

2

Try to pass &mat (multiple elements pointer) instead of &mat[0] (single element pointer). It this works, do the same for ipiv and vec.

Reasoning: C accepts double (f64) only as function arguments, so mat[0] might be converted to f64 first, even if after that is the address-of (&) operator (of course if that is the case, it is a bug).

3

I tried that as well, but it still gives the same wrong results.

Edit:
An update on the matter. I reimplemented the test code in C and it gives the same numerical results (still does not explain the f64 issue). I thought there might be a bug in my lapacke version, so I compiled one from source, but same result. Still in need of help.

4

Keep the versions without [0] that pass compilation, and change laint to i64 (currently it is translated to i32).

The definitions are:

#ifndef lapack_int 
#if defined(LAPACK_ILP64)
#define lapack_int        int64_t
#else
#define lapack_int        int32_t
#endif
#endif

lapack_int LAPACKE_sgesv( int matrix_layout, lapack_int n, lapack_int nrhs,
                          float* a, lapack_int lda, lapack_int* ipiv, float* b,
                          lapack_int ldb );