[ subroutine @ 13.03.2009. 03:43 ] @
Pozdrav!

Ne mogu da verujem da sam tek sada otkrio ovaj forum!

Ok, imam problem.

Imam komplexnu matricu dimenzija 11x11 koju treba da invertujem. Skinuo sam sa neta kod za inverziju matrice, i on radi u druga dva primera u kojima je koristim. Međutim, u ovom mi kaže da je "matrix is non invertible".
Ako zamenim nekoliko vrsta u matrici (što je dozvoljeno) onda sve bude ok. Na glavnoj dijagonali nemam nultih elemenata, pa me nervira zašto ne mogu da je invertujem u originalnom izgledu.

Da li neko ima ideju zašto se ovo dešava i da li imate neki super-efikasan univerzalni kod za inverziju svih matrica na svetu!
(inače, sedim već 12 sati ispred monitora i pokušavam sve što mi pada na pamet!)

Thanx!
[ Časlav Ilić @ 13.03.2009. 10:51 ] @
To što na glavnoj dijagonali nema nula na početku, ne znači da se neće neka pojaviti tokom proste eliminacije. Tada je potrebno zameniti vrste/kolone, kao što si i učinio, mada bi kôd to trebalo sam da uradi. Štaviše, algoritmi za direktnu eliminaciju često obrću vrste/kolone bez obzira da li se zatekne nula ili ne, tako da na glavnoj dijagonali uvek bude trenutno najveći (po aps. vred.) element matrice, jer je to dobro za numeričku stabilnost. Vidi da li u tom kôdu što koristiš može da se aktivira ovo (na engleskom pivoting).
[ subroutine @ 13.03.2009. 12:14 ] @
Ne, definitivno ih on ne okreće sam. Ali kada ih ja izmenim sve bude ok. Mada, nisam siguran da ih ja izokrećem tako da postavim maximalne elemente na glavnu dijagonalu jer ja ni ne znam koji su mi elementi matrice, znam samo koji su nulti a koji nenulti.
A ovo što si mi rekao za numeričku stabilnost nisam ni znao, a ima PUNO smisla pošto mi rešenje tog sistema jednačina ne konvergira baš najsjajnije.
Ako imaš takav kod kakav si opisao bio bih ti zahvalan za link ili nešto slično. U svakom slučaju hvala puno za objašnjenje.

Pozdrav.
[ Časlav Ilić @ 13.03.2009. 16:49 ] @
Hm, da nemaš možda biblioteku Lapak pri ruci? Možeš onda da upotrebiš njene potprograme zgetrf i zgetri za LU faktorizaciju i inverziju kompleksnih matrica (probni kôd u prilogu, fortran mi nije baš na ruku). zgetrf pritom automatski premeće vrste pri eliminaciji.

Kod mene na Debijanu:
Code:
$ apt-get install liblapack-dev liblapack-doc
$ man zgetrf
$ man zgetri
$ gfortran lapack_test_matinv.f95 -llapack
$ ./a.out
A =
 (  0.3162  0.4957*i) (  0.5710  0.1861*i) (  0.8777  0.8059*i)
 (  0.3292  0.5304*i) (  0.8459  0.8376*i) (  0.0260  0.0697*i)
 (  0.1035  0.6892*i) (  0.1913  0.9955*i) (  0.9950  0.7372*i)
pivoting:
  swap   1   2
  swap   2   2
  swap   3   3
L\U =
 (  0.3292  0.5304*i) (  0.8459  0.8376*i) (  0.0260  0.0697*i)
 (  0.9418 -0.0116*i) ( -0.2354 -0.5930*i) (  0.8524  0.7406*i)
 (  1.0255  0.4414*i) (  0.5223 -0.3093*i) (  0.3248  0.5311*i)
A^-1 =
 ( -1.6646 -5.5931*i) (  4.1295 -2.0745*i) (  0.7524  5.4899*i)
 (  0.1930  3.0004*i) ( -1.6761 -0.0050*i) (  0.2025 -2.8357*i)
 ( -0.0138  0.9749*i) ( -1.4626  0.1170*i) (  0.8380 -1.3704*i)
A * A^-1 =
 (  1.0000 -0.0000*i) (  0.0000  0.0000*i) ( -0.0000  0.0000*i)
 (  0.0000 -0.0000*i) (  1.0000  0.0000*i) ( -0.0000  0.0000*i)
 (  0.0000 -0.0000*i) (  0.0000  0.0000*i) (  1.0000  0.0000*i)

Vidiš gore da su vrste 1 i 2 odmah zamenjene, pošto je abs(A(1, 1)) manje od abs(A(2, 1)).
[ subroutine @ 19.03.2009. 02:46 ] @
Nisam koristio ovaj način niti tu biblioteku, ali mi program ipak radi. Ima slična stvar u Fortranovom Math Library-ju kada se koristi njegov potprogram LINGC i LINGR, ali nemam ni tu biblioteku.

Hvala puno u svakom slučaju, važno mi je da šljaka i ovako.

Pozdrav.