CLAPACK在Linux下的安装和使用CLAPACK的矩阵求逆C程序

装了3天的CLAPACK，终于算是装好了。具体的步骤我已经忘得差不多了，但是安装过程中最耗时的一个错误我记录下来。是通过CLAPACK的开发小组官方论坛有人解决的。

 

接下来几天准备使用CLAPACK来写个矩阵求逆程序并且测试其在大规模矩阵输入的情况下的性能。甚至可以考虑换一个高性能编译版本的BLAS库。

by » Tue May 10, 2011 11:14 am

Hi, this is my first time using CLAPACK, and the version is 3.2.1, my platform is Linux localhost.localdomain 2.6.35.6-45.fc14.x86_64 #1 SMP Mon Oct 18 23:57:44 UTC 2010 x86_64 x86_64 x86_64 GNU/Linux.

During my installation process, there is a problem:

cd TESTING; make

ZGG: Testing COMPLEX16 Nonsymmetric Generalized Eigenvalue Problem routines
./xeigtstz < zgg.in > zgg.out 2>&1
/bin/sh: line 1: 22252 segmentation fault (core dumped) ./xeigtstz < zgg.in > zgg.out 2>&1
make: *** [zgg.out] fault 139

I don't know whether this problem is the real cause for my following link problem, but maybe this information is helpful.

Although with the above testing problem, I still get blas_LINUX.a lapack_LINUX.a libcblaswr.a tmglib_LINUX.a, and I have correctly put them into /usr/local/lib, and put clapackh and f2c.h and blaswrap.h into /usr/local/include

then I tried to compile a simple program, but it has errors:

g++ clapacktest.cc -lcblas -lclapack -llibf2c

/usr/local/lib/libclapack.a(sgetrf.o): In function `sgetrf_':
sgetrf.c:(.text+0x45c): undefined reference to `f2c_strsm'
sgetrf.c:(.text+0x524): undefined reference to `f2c_sgemm'
/usr/local/lib/libclapack.a(sgetrs.o): In function `sgetrs_':
sgetrs.c:(.text+0x194): undefined reference to `f2c_strsm'
sgetrs.c:(.text+0x1cf): undefined reference to `f2c_strsm'
sgetrs.c:(.text+0x303): undefined reference to `f2c_strsm'
sgetrs.c:(.text+0x33e): undefined reference to `f2c_strsm'
/usr/local/lib/libclapack.a(sgetf2.o): In function `sgetf2_':
sgetf2.c:(.text+0x180): undefined reference to `f2c_isamax'
sgetf2.c:(.text+0x1ca): undefined reference to `f2c_sswap'
sgetf2.c:(.text+0x329): undefined reference to `f2c_sger'
sgetf2.c:(.text+0x3b5): undefined reference to `f2c_sscal'
collect2: ld return 1

what is the problem? thanks in advance ! and my program is pasted below FYI.

#include
using namespace std;
extern "C"
{
#include
#include
}

int
main (void)
{
integer M = 3;
integer N = 1;
real a[9] = { 4, 3, 11, 2, -1, 0, -1, 2, 3 };
real b[3] = { 2, 10, 8 };
integer lda;
integer ldb;
integer INFO;
lda = M;
ldb = M;
integer ipiv[M];

sgesv_ (&M, &N, a, &lda, ipiv, b, &ldb, &INFO);
if (INFO == 0)
{
for (int i = 0; i < M; i++)
{
cout << b[i] << endl;
}
}
else
{
cout << "Failed." << endl;
}
return 0;
}

Posts: 2

Joined: Tue May 10, 2011 11:03 am

---

by » Tue May 10, 2011 4:35 pm

Hi,
Could you try changing the linking sequence and adding a BLAS library
The cblas is just a C wrapper, it needs a BLAS library.

Code:

g++ clapacktest.cc -lclapack -lcblas -lblas -llibf2c


Site Admin

Posts: 244

Joined: Thu Dec 09, 2004 8:07 am

---

by » Tue May 10, 2011 4:50 pm

Thanks !
I put blas_LINUX.a to /usr/local/lib/libblas.a
and "g++ clapacktest.cc -lclapack -lcblas -lblas -llibf2c" works.
maybe this is a good start to CLAPACK, I'd like to try some matrix inversion program later. :)

 

后来在论坛帮助下，写了个矩阵求逆的C程序，似乎这个程序极快，我把输出禁止掉后，使用clock()来统计时间，是620000clocks，大概0.6秒就把一个1000×1000的矩阵求逆了。。。。

程序后面附带有矩阵乘法的验算代码。

 

用这个命令编译（运行的时候把输出重定向到文件，更快。我觉得原因是因为屏幕是字符设备无缓冲吧。）

gcc large_matrix_invert.c -lclapack -lcblas -lblas -llibf2c -lm -Wall -lcblaswr

代码是：

 

 

#include <stdlib.h>
#include <string.h>
#include <stdio.h>
#include <math.h>
#include <time.h>

#include <f2c.h>
#include <clapack.h>
#include <blaswrap.h>

#define N 500

real oldA[N * N];
real A[N * N];
real C[N * N];            //for oldA * A


int
main (void)
{
    integer i, j, n;
    integer lda, arows, acols;

    //integer ipiv[4];

    integer ipiv[N];
    integer lwork;
    integer info;
    float work[N * N];

    arows = acols = N;
    lda = N;
    memset (ipiv, 0, sizeof (int) * N);
    memset (work, 0, sizeof (float) * N * N);
    //memset(C,0,sizeof(float)*N*N);


    srand (time (NULL));
    for (i = 0; i < N * N; i++) {
    A[i] = rand () % 100;
    oldA[i] = A[i];
    }

    printf ("matrix A\r\n\r\n");
    for (i = 0; i < N; i++) {
    for (j = 0; j < N; j++)
     printf ("A[%ld,%ld]:%f\t", i, j, A[j * N + i]);
    printf ("\r\n");
    }

    clock_t st = clock ();

    extern int sgetrf_ (integer * m, integer * n, real * a, integer * lda,
            integer * ipiv, integer * info);

    sgetrf_ (&arows, &acols, A, &lda, ipiv, &info);
    //sgetrf_(&arows,&acols,A,lda,ipiv,&info);

    printf ("\r\n\r\nmatrix A after sgetrf\r\n");
    for (i = 0; i < N; i++) {
    for (j = 0; j < N; j++)
     printf ("A[%ld,%ld]:%f\t", i, j, A[j * N + i]);
    printf ("\r\n");
    }
    printf ("\r\n\r\ninfo for sgetrf_:%ld\r\n\r\n", info);


    //n=acols*arows;

    n = N;
    lda = N;
    lwork = N * N;

    extern int sgetri_ (integer * n, real * a, integer * lda, integer * ipiv,
            real * work, integer * lwork, integer * info);

    sgetri_ (&n, A, &lda, ipiv, work, &lwork, &info);
    //sgetri_(&n,A,lda,ipiv,work,&lwork,&info);


    fprintf (stderr, "%ld\n", (clock () - st));

    printf ("\r\n\r\nMatrix A after sgetri\r\n");
    for (i = 0; i < N; i++) {
    for (j = 0; j < N; j++)
     printf ("A[%ld,%ld]:%f\t", i, j, A[j * N + i]);
    printf ("\r\n");
    }
    printf ("\r\n\r\ninfo for sgetri_:%ld\r\n\r\n", info);


    printf ("\r\noldA * A = \r\n");
    char c = 'n';
    real alpha = 1.0;
    real beta = 0.0;

    extern int sgemm_ (char *transa, char *transb, integer * m, integer *
         n, integer * k, real * alpha, real * a, integer * lda,
         real * b, integer * ldb, real * beta, real * c__,
         integer * ldc);

    sgemm_ (&c, &c, &n, &n, &n, &alpha, oldA, &n, A, &n, &beta, C, &n);
    //f2c_sgemm(&c,&c,&n,&n,&n,&alpha,oldA,&n,A,&n,&beta,C,&n);

    for (i = 0; i < N; i++) {
    for (j = 0; j < N; j++)
     printf ("C[%ld,%ld]:%f\t", i, j, C[j * N + i]);
    printf ("\r\n");
    }

    return (0);
}