LAPACK and BLAS are originally written in Fortran and meant to be used in Fortran programs. Many vendors supply an optimised version of the LAPACK and BLAS libraries. Naturally, people who are programming in C or C++ and want to make use of the efficient implementation of the LAPACK/BLAS libraries. Here we will demonstrate how this can be done.

At the end you find a complete working example, together with a script to run the same program using various interfaces and LAPACK libraries.

The examples wil be based on the call of dsyev, the fortran code to do this is:

character*1 jobz, uplo integer n,lda,lwork,info double precision, allocatable a(:,:) double precision, allocatable w(:) double precision, allocatable, work(:) double precision wwork n=1000 lda=1000 jobz='V' uplo='U' allocate(a(lda,n)) ... code to fill matrix a ... allocate(w(n)) lwork = -1 ! workspace query call dsyev(jobz, uplo, n, a, lda, w, wwork, lwork, info) lwork = wwork+0.1 ! now the real work: allocate(work(lwork)) call dsyev(jobz, uplo, n, a, lda, w, work, lwork, info) deallocate(work)

You can use the LAPACK library as-is, by mimicking the way Fortran calls Fortran subroutines. Three points are important here:

- use dsyev_ in stead of dsyev
- every parameter must be a pointer
- be aware of strings of varying length: Fortran has it's own way to deal with strings in parameter lists. In general it is safe to use a pointer to a character if the LAPACK subroutine expects only one character.
- matrices in C and Fortran differ: C is row-major, while Fortran is colomn-major.

So, the C-version of above Fortran fragment would become:

void dsyev_(char *jobz, char *uplo, int *n, double *a, int *lda, double *w, double *work, int *lwork, int *info); ... char jobz, uplo; int n,lda,lwork,info; double *a; double *w; double *work; double wwork; n = 1000; lda = 1000; jobz = 'V'; uplo = 'U'; a = malloc(sizeof(double)*n*lda); w = malloc(sizeof(double)*n); lwork = -1; // workspace query dsyev_(&jobz, &uplo, &n, a, &lda, w, &wwork, &lwork, &info); lwork = wwork+0.1; work = malloc(sizeof(double)*lwork); // now the real thing: dsyev_(&jobz, &uplo, &n, a, &lda, w, work, &lwork, &info); free(work);

The newer versions of LAPACK from www./netlib.org/lapack , contain an interface to use LAPACK from C. How to use this interface is described in detail in www.netlib.org/lapack/lapacke , we will give here an example, based on the Fortran code above.

Important things are:

- The name of the subroutine becomes LAPACKE_dsyev
- The function returns the info argument, the info argument is not in the parameter list
- All parameters are used in a C-like fashion: only the arrays are pointers
- In general, there are two versions of the function: for example LAPACKE_dsyev and LAPACKE_dsyev_work. The first one determines the optimal workspace itself, the second one expects that the workspace is provided by the caller. The first one must be called without workspace arrays and dimensions. We will demonstrate the first one here.
- If one or more parameters are 2-dimension matrices, an extra parameter is needed: LAPACK_COL_MAJOR or LAPACK_ROW_MAJOR. The first one for Fortran matrix layout, the second for C matrix layout. This parameter is placed before all other parameters.
- LAPACKE can also be used to call LAPACK subroutines from Intel's MKL library, see below.

Here the example:

#include <lapacke.h> ... char jobz, uplo; int n,lda,info; double *a; double *w; n = 1000; lda = 1000; jobz = 'V'; uplo = 'U'; a = malloc(sizeof(double)*n*lda); w = malloc(sizeof(double)*n); info = LAPACKE_dsyev(LAPACKE_COL_MAJOR,jobz, uplo, n, a, lda, w);

An advantage is, that an header file is provided which prevents trivial mistakes in calling the function.

Here we provide you with working example codes that you can run on the Lisa system. The examples include the usage of MKL from C programs.

Using copy/paste, copy both files (testdsyev.c and doit.txt, see 'attached files' below) to Lisa and give the following commands:

chmod +x doit.txt ./doit.txt

In the output you can see which method is used (LAPACK_PLAIN or LAPACKE), the order of the matrix used; the value of `info`

after the function call and the time it took. In the script `doit.txt`

you can study the linking flags that are used.

Attached files:

The SURFsara Data Archive allows the user to safely archive up to petabytes of valuable research data.

Persistent identifiers (PIDs) ensure the findability of your data. SURFsara offers a PID provisioning service in cooperation with the European Persistent Identifier Consortium (EPIC).

B2SAFE is a robust, secure and accessible data management service. It allows common repositories to reliably implement data management policies, even in multiple administrative domains.

The Data Ingest Service is a service provided by SURFsara for users that want to upload a large amount of data to SURFsara and who not have the sufficient amount...

The Collaboratorium is a visualization and presentation space for science and industry. The facility is of great use for researchers that are faced with...

The ELvis remote visualization cluster provides high-end GPU rendering capabilities for performing (interactive) visualization of large datasets...