Hello,
I have to invert a large matrix ( N > 10000 x 10000) this matrix is
square and complex. I need to use many times its inverse.
Do you really need the inverse or do you just want to solve several times
a linear system of equations but with the same coefficient matrix?
[1] If you really want the inverse: First call LAPACK ZGETRF, this will
perform an inplace LU factorization of your coefficient matrix. Then call
ZGETRI, this will use the inplace LUfactorization of your initial matrix
and returns you the inverse of the matrix A. Same thing in place as well.
So basically the code is two LAPACK calls. (ZGETRF+ZGETRI). When you want
to use the inverse, for example to solve a linear system of equations, you
just need to multiply by the matrix inverse. This can be done with ZGEMV
(one righthand sides) or ZGEMM (multiple righthand sides). This approach
costs: 2/3*N^3 for ZGETRF, 4/3*N^3 for ZGETRI and 2*N^2 per righthand
sides (ZGEMV or ZGEMM).
[2] However if all you want to do is to solve one linear system of
equations with several righthand sides, this is not the recommended way
to go. A more stable and less expensive approach is the following. First
call ZGETRF. This performs the LU factorization of your coefficient matrix
inplace and cost 2/3*N^3 FLOPS. Keep the NbyN array computed by LAPACK
as long as you want to solve a linear system with the initial matrix. Each
time you want to solve for a given righthand side, you need to call
ZGETRS, ZGETRF uses the LU factorization output form SGETRF. (ZGETRS can
handle one righthand side or multiple, if multiple, the routine is then
more efficient than several call to this routine.) The cost is: 2/3*N^3
for ZGETRF and 2*N^2 per righthand sides. Same thing two lines calls to
LAPACK. As long as you do not touch the output from ZGETRF, you are free
to reuse ZGETRS, i.e. you are free to solve another set of righthand
sides.
So [2] is less expensive and moreover is more stable. (Skip the stability
discussion.)
If you want to see a call to ZGETRF followed to call to ZGETRS look at
the LAPACK routine ZGESV. That's basically all that routine does.
With respect to CLAPACK. Sure you can use CLAPACK. Say that if you have
LAPACK already installed on your system then you better off using it.
CLAPACK is done for arch/OS that do not have a (free/cheap/atall) Fortran
compiler.
For this kind of sizes (N=10,000). Do not forget to use a good BLAS behind
LAPACK, if possible multithreaded.
j
On Sun, 10 Jun 2007, Michel Yvert wrote:
Dear CLAPACK helpers,
This mail is to ask for an advice: I have to invert a large matrix ( N >
10000 x 10000) this matrix is square and complex. I need to use many
times its inverse. I asked Roldan Pozo who suggested to use CLAPACK (I
work in C++ and would like to incorporate this inversion in a C++ job).
Could you indicate me which is the best approach (the routine to use)
for a best efficiency regarding the precision and the calculation time ?
I would appreciate very much a simple example working on a similar task.
Many thanks in advance
Michel Yvert

=======================================================
Michel Yvert
Laboratoire d'AnnecyleVieux de Physique des Particules
B.P. 110
74941 AnnecyleVieux Cedex
France
Tel : +33 4 50 09 16 75 / 04 50 09 16 75
Fax : +33 4 50 27 94 95 / 04 50 27 94 95
http://wwwlapp.in2p3.fr
_______________________________________________
Lapack mailing list
Lapack@Domain.Removed
http://lists.cs.utk.edu/listinfo/lapack
