Prev: Anybody use GPPTOOL to do programming before
Next: Multiple Key Sort in Fortran
From: David Jones on 11 Jan 2006 06:19
Andrea Campera wrote:
> To be more precise I have to do an operation like (A*B+I)^ -1.

Are you aware that you may be able to speed things up by using a known identity for this type of expression ? See 3.2.2 of the following:
http://www2.imm.dtu.dk/pubdb/views/edoc_download.php/3274/pdf/imm3274.pdf

Some of the libraries suggested may have this type of thing already built in, but you need to be able to recognise that your problem has the right structure.
From: Joost on 11 Jan 2006 06:52
Also here, use a library routine. IIRC the inverse can be computed
using the LAPACK CGESV routine. The usual remark is that often you
might not need the inverse, but you really want to solve A*X=B instead.
CGESV does exactly that, so solve with B=I if you need the inverse.
Also here, CGESV will be part of MKL (typically used on intel chips)
and ACML (used on AMD), and should be relatively fast also with netlib
lapack if combined with atlas or GOTO blas (available for both brands).

And yes, the reference provided by David Jones might be an interesting
additional read...

Joost

From: David Flower on 11 Jan 2006 08:39

andrea campera wrote:
> Hi all,
>
> I have a question regardin matrix multipllication. In my program I call
> a subroutine to multiplies to n*n matrices thousands of times
>
> SUBROUTINE multi(Nc,m1,m2,mtot)
> INTEGER i,j,l,Nc
> complex aa,m1(Nc,Nc),m2(Nc,Nc),mtot(Nc,Nc)
> do i=1,Nc
> do j=1,Nc
> aa=0.
> do l=1,Nc
> aa=aa+m1(i,l)*m2(l,j)
> enddo
> mtot(i,j)=aa
> enddo
> enddo
> return
> end
>
> my question is: is there a faster algorithm to solve this simple
> multiplication? where can I find the code (fortran 77 if possible)?
> Thanks in advance for any help
>
> Regards
>
> Andrea

Well, I can suggest one improvement - lower case L in many fonts can
look identical to 1 (numeral one), so never use it!

Dave Flower

From: Gordon Sande on 11 Jan 2006 08:41
On 2006-01-11 04:14:42 -0400, "Andrea Campera"
<andrea.campera(a)iet.unipi.it> said:

> To be more precise I have to do an operation like (A*B+I)^ -1. up to
> now I have used
>
> call multi(...)
> call summa(...)
> call gaussj(...)
>
> I post also the profile:
>
> Each sample counts as 0.01 seconds.
> % cumulative self self total
> time seconds seconds calls Ks/call Ks/call name

<snip>

> if I increase the dimension Nc of the matrices gaussj and multi become
> the most time consuming routines. I have seen the subroutine cgemm.f of
> Blas for the multiplication, but what about the inversion? is ther
> anything faster than gaussj? I use g77 (3.2 or 3.4) Pentium IV or Xeon
> or Opteron or AMD64 (we have several machines).
>
> Andrea


The fastest way to invert matrices is to not invert matrices. This apparent
piece of self contradictory advice really means that most uses of inversion
are better done some other way.

Solving is usually what is wanted and that can be better done by factoring
and solving with the the factors. The question is "Why to you think you
actually need the inverse?" Other that your simple derivation in matrix
algebra.

Matrix multiply will always be n^3 and inversion and/or solving will
also be n^3. You can fiddle the coefficient of multiply by polishing
code, or having someone else do it for you in a library. You can also
polish inversion but beware the "highly polished" codes that pay no
attention to numerical stability. Many Gauss-Jordan matrix inversion
programs are the trivial algebra coded naively so you have to read and
understand the documentation. Factoring has a bigger reduction in
coefficient over inversion than you are ever going to see with polishing.
And it may be that with a bit of insight into you real problem that
there are other things that can be done. Thee are examples of easy
problems with awkward derivations that are turned into awkward, and
numerically harder, problems with easy derivations. For more details
consult a competent numerical analyst.




From: robin on 11 Jan 2006 17:31
"andrea campera" <andreacampera(a)gmail.com> wrote in message
news:1136934796.383711.298240(a)f14g2000cwb.googlegroups.com...
> I have a question regardin matrix multipllication. In my program I call
> a subroutine to multiplies to n*n matrices thousands of times
>
> SUBROUTINE multi(Nc,m1,m2,mtot)
> INTEGER i,j,l,Nc
> complex aa,m1(Nc,Nc),m2(Nc,Nc),mtot(Nc,Nc)
> do i=1,Nc
> do j=1,Nc
> aa=0.
> do l=1,Nc
> aa=aa+m1(i,l)*m2(l,j)
> enddo
> mtot(i,j)=aa
> enddo
> enddo
> return
> end
>
> my question is: is there a faster algorithm to solve this simple
> multiplication? where can I find the code (fortran 77 if possible)?

You are better off using a Fortran 90 or later compiler because
matrix multiplication is built in. The implementations are,
typically faster than what you have written


First  |  Prev  |  Next  |  Last
Pages: 1 2 3 4
Prev: Anybody use GPPTOOL to do programming before
Next: Multiple Key Sort in Fortran