How large matrix can LAPACK support

[gyhuang]

Hello
My code is as simple as follows

LaSymmMatDouble A(N, N);
LaVectorDouble v(N);
LaGenMatDouble Eigenvectors(N, N);
LaEigSolve(A, v, Eigenvectors);

However, I find that N can not be more than 40000. I feel puzzled about this. From the Internet a lot of people said LAPACK supports large matrix. Also, it is a little bit slow to calculate eigenvectors.

My problem is simply to calculate eigenvectors and eigenvalue for matrix as large as 40000*40000. Any idea on this?

Thanks.

[admin]

On what kind of machine are you running, how much memory do you have?

[error5772]

Hi!

I have the same question! And want to ask for the integer type of the indizes
in C/C++. Can I use matrices/vectors containing more than 2^32 elements
using "long" indices in these LAPACK-routines?

I have a dual Westmere (2.4GHz) system with 64bit and 48GB of RAM.

Michael

[Julien Langou]

Hi Michael,

As of now, you definitely cannot use one of the dimension of your matrix (so M or N) to be bigger than 2^32, (2^31-1 to be precise).
This is because we use FORTRAN INTEGER in the subroutine interfaces which are mapped by default to 32-bit integers.
This is starting to be a problem with ScaLAPACK and tall and skinny matrices (so typically vectors, which are very tall and very skinny matrices).
First time, someone report some problem with this for LAPACK. Glad to know.

On the other hand, some of our routines have problems when M*N is greater than 2^32 and this is much more of a concern.

Cheers,
Julien.

[rodney]

One solution to the 32-bit integer problem in LAPACK is to recompile the entire library with a compiler option that promotes the fortran INTEGER type to 64-bits. In gfortran, the option is -fdefault-integer-8 and for the Intel compiler the option is -integer-size 64. Note that you must also then recompile and use the included FORTRAN BLAS library with the same option (or find another BLAS library that uses 64-bit integers).

Rodney

[error5772]

Many thanks!

I soon get a double-matrix with more than (2^31 - 1) elements and
more than 16GB RAM. I will try to recompile ATLAS (Goto-BLAS makes
problems) and BLAS 3.3.1 with gfortran using "-fdefault-integer-8" and
write back to that thread.

Michael

[error5772]

True 64bit LAPACK on Mr Goto's true 64bit blas!

Here is what I found out:

The following 3 steps lead to a WORKING LAPACK/Goto-BLAS on a 64bit Dual Intel E5645 ("Westmere"->use NEHALEM)
with 2x6 cores and 2x12 threads using matrices with "long" array-indices - I call it "true 64bit".
1. Install true64bit-Goto-BLAS by editing "Makefile.rule":
TARGET=NEHALEM, CC=gcc, FC=gfortran, BINARY=64, USE_THREAD=1, NUM_THREADS=6, INTERFACE64=1
call the lib "libgoto6t.a" and compile some others "libgotoXt.a" as above with NUM_THREADS=X. Copy them to "/usr/local/lib64".
If you use too many threads (e.g. 12) Goto-BLAS will jump out with "sementation fault"! You can use NUM_THREADS=11
and call that BLAS "libgoto11t.a" - it will be your fastest BLAS, but will sometimes jump out (try again).

2. Install true64bit-LAPACK-3.3.1 by editing "make.inc":
FORTRAN=gfortran -m64 -m128bit-long-double -fdefault-integer-8 -fimplicit-none,
OPTS=-O3 -funroll-all-loops, DRVOPTS=$OPTS, NOOPT=-O0,
LOADER=gfortran, LOADOPTS=-L/usr/local/lib64/ -lgoto6t -lpthread
TIMER=INT_ETIME
BLASLIB=/usr/local/lib64/libgoto6t.a -lpthread

Compile it with: make lapack_install, (make variants), make lapacklib, make tmglib
and make lapack_testing (Should work for "libgoto6t.a, "libgoto11t.a" is not good here.)
make blas_testing does not work (-fdefault-integer-8 is a problem for the testers).
Call it "liblapackgoto.a" and copy to "/usr/local/lib64".

3. Have fun with extremely large matrices! You can link the 6 thread-Goto-BLAS with
"-llapackgoto -lgoto6t -lpthread -lgfortran" or your fastest
"-llapackgoto -lgoto11t -lpthread -lgfortran" that sometimes "jumps out of the threads",
but is unbeatable fast!

No need to compile LAPACKE - just call the the Fortran-function two times:
work space query(resize the workspace arrays) - then run the function. Use your own
LAPACK-header in "lapacke.h"-style for C++ (examples):
#include<complex>
extern "C" {
double dlamch_ ( char* );
void dsyev_( char*,char*,long*,double*,long*,double*,double*,long*,long*);
void zheev_( char*,char*,long*,std::complex<double>*,long*,double*,
std::complex<double>*,long*,double*,long* );
}

These 3 functions can be called in main in the form e.g.:
char sw = 'S';
double tol = 2.0*dlamch_(&sw);
(sorry no example for dsyev_ and zheev_)
and even the complex functions work without problems if you give the
pointer to the first element of your complex<double> array of the
Column-Major-Array of your matrix or the workspace-array.

Had to write the report for other people with problems.
Thanks to the forum for the deciding two tips.

Michael

[MM2011]

I have a question about this too. I am trying to solve a linear least squares problem with 300,000 variables. I want to find the pseudoinverse of a 300,000 X 300,000 matrix. Can LAPACK/ScaLAPACK solve this problem?

[error5772]

This problem can probably be solved by ScaLAPACK with a cluster
consisting of 1000 Cray-Supercomputers . I think its better to find some
approximation or an other way to solve this problem. But of course - there
are problems where you have to solve it only this way.

Good luck!

Michael