Invertire una matrice
Innuendo
Io finora usavo il metodo, data una matrice A, di moltiplicare 1/det(A) per
la trasposta dei suoi complementi algebrici, ma non mi sembra un
procedimento molto rapido specialmente in caso di matrici che non siano le
semplici 2x2....
grazie in anticipo
Innuendo
--
My soul is painted like the wings of butterflies
Fairy tales of yesterday will grow but never die
I can fly - my friends
=====================
ICQ UIN : 29447269
EMAIL : tetrag...@yahoo.it
Giovanni Anzani
> Esiste un metodo rapido per invertire una matrice?
> Io finora usavo il metodo, data una matrice A, di moltiplicare 1/det(A) per
> la trasposta dei suoi complementi algebrici, ma non mi sembra un
> procedimento molto rapido specialmente in caso di matrici che non siano le
> semplici 2x2....
> grazie in anticipo
> Innuendo
Farebbe piacere anche a me conoscerlo, rispondeteci!
per me il problema e' ridurre i tempi: 36 sec per
invertire 9 matrici 8x8 e con Cramer siamo la' 34 sec.
Che se po' fa'?
Grazie Giovanni
Giovanni Anzani
Ciao Innuendo,
prima t'ho considerato un po' pochino, sai l'emozione di trovare finalmente qualcuno con i miei stessi problemi....
Continuo adesso, ma off-line.
Ti premetto che passavo di qua da it.comp.software.cad, che non sono un matematico ma un quasi architetto e temo che, per quello che qua dirotti in seguito, mi prendero' qualche calcio nei denti.
In attesa che qualche "it.scienza.matematico" ci indirizzi verso metodi alternativi ti diro' cosa ho fatto sino ad ora io [vedrai che calci...], se puo' esserti utile rispondimi che lo e' e nel prossimo messaggio taglio-incollo un po' di codice che nel mio caso e' in autolisp, ma con un po' di pazienza te lo puoi tradurre in quello che vuoi.
Cio' premesso.
Se la tua non e' una pura curiosita' matematica, ma vuoi farti il tuo bel programmino di inversione, dovresti specificare il "linguaggio" da te usato, perche' la riduzione del tempo di calcolo dipende molto da quello; ad esempio le mie belle matriciotte 8x8 excel me le inverte in tempo reale, mentre quel "dinosauro ignorantone" dell'autolisp ha prima voluto che gli spiegassi che cavolo e' una matrice, la sua trasposta, un suo minore ... risultato ci mette 4 sec.
Non sarebbe neanche tanto se non fosse che 1°) ne devo fare 9 di fila (e siamo a 36 sec) e con gli spiccioli della procedura saliamo a 40 sec. e che 2°) che gli racconto in sede di discussione di tesi mentre attendiamo la fine della procedura?
-----
Non avendo questo dato (la tua lingua) partiro' dal presupposto che sei nella condizione peggiore, ovvero la mia.
premessa: per me (in lisp) una matrice A ad es. e' fatta così:
A = ((11 12 13)(21 22 23)(31 32 33))
Ed eccoti qualche truce suggerimento [calci liberi!]:
----
1) ti sei fatto una procedura generale per la matrici nxn? Se ti e' possibile passa ad una procedura (n<k)x(n<k) dove k e il rango massimo delle matrici con cui hai a che fare; lo so piange il cuore a fare qualcosa di limitato, ma cosi' facendo forse puoi semplificare un po'. Questo vale sia per la procedura di inversione sia per quella del calcolo del determinante.
----
2) l'estrazione dei minori e' pure un passaggio che se non curato porta via tempo; ti spiego:
per il minore (ad es. "di centro 23") si potrebbe:
[MINORE (A)
tolgo-riga (2) = ((11 12 13)(31 32 33))
trasposta = ((11 31)(12 32)(13 33))
tolgo-riga (3) = ((11 31)(12 32))
trasposta = ((11 12)(31 32))]
ma se consideri che di tale minore devi solo calcolare il det allora siccome det(A) = det(A)t:
[MINORE TRASPOSTO (A)
tolgo-riga (2) = ((11 12 13)(31 32 33))
trasposta = ((11 31)(12 32)(13 33))
tolgo-riga (3) = ((11 31)(12 32))]
e cosi togliamo una trasposizione di torno per ogni minore che hai da calcolare.
3) altra nota dolente il det. In breve: sali di livello!; mi spiego:
l'ideale sarebbe che tu arrivassi a definire il det di livello (k-1) cosi' eviteresti di dover "scavare in basso" fino al det da te definito.
esempio, un conto e' [limitiamoci al caso di M6 6x6]:
------
[DET3X3 ...procedura...]
[DET6X6 (M6)
ripeti (6) minore (1*(M6)) = (6*(M5))
ripeti (5) minore (6*(M5)) = (5*6*(M4))
ripeti (4) minore (5*6*(M4)) = (4*5*6*(M3))
DET3x3 (4*5*6*(M3))] ;;; e qua ti calcoli 120 det3x3
------
[DET5X5 ...procedura...]
[DET6X6 (M6)
ripeti (6) minore (1*(M6)) = (6*(M5))
DET5x5 (6*(M5))] ;;; e qua ti calcoli 6 det5x5
------
E' vero che quanto a velocità: [DET3X3 ...] < [DET5X5 ...]
ma (almeno nel mio caso): 120*[DET3X3 ...] >> [DET5X5 ...]
Se ti serve, fino al 5x5 ho la formula diretta del det.
------
Questo e' quanto ho da dire per il momento
P.S. Qualcuno ha e puo' passarmi/ci lo sviluppo di det [A]6x6 o [A]7x7, o in alternativa sa indicarmi con quale programma posso fare calcolo simbolico con simboli del tipo:
|a11 ... a16| |a11 ... a17|
det | ... ... ... | o magari det | ... ... ... |
|a61 ... a66| |a71 ... a77|
Ciao Giovanni
Subatomic
Innuendo <tetragr...@SPAMMAREyahoo.it> wrote in message
8dnr6t$ne$1...@fe1.cs.interbusiness.it...
> Esiste un metodo rapido per invertire una matrice?
> Io finora usavo il metodo, data una matrice A, di moltiplicare 1/det(A)
per
> la trasposta dei suoi complementi algebrici, ma non mi sembra un
> procedimento molto rapido specialmente in caso di matrici che non siano le
> semplici 2x2....
>
Se la tua necessità è quella di invertire *manualmente* una matrice (i.e.
devi fare un esame di Geometria)
purtroppo devi fare proprio in questo modo.
Alternativamente detta A la tua matrice da invertire puoi *costruire* la
matrice AI
affiancando alla matrice A la matrice unità con tante righe e colonne quante
sono quelle di A.
Fatto questo mediante *trasformazioni elementari* (inversione di riga
esclusa...mi pare)
trasformi la tua matrice
"n x 2n" di modo da ottenere
al posto di A la matrice unità......è chiaro che al posto di I otterrai
A^-1.
( (AI)^-1 = I A^-1 )
Nel caso tu possa fare questo con un pc trovi una marea infinita di
procedure o programmi che invertono matrici,
ma a riguardo dovresti essere più preciso.
Ciao
Giovanni Anzani
Grazie intanto per la risposta e spiegherotti che debbo fa':
N.B. il mio messaggio:
"Re: ... continua Invertire una matrice" - ven. 21/4/2000 ore 01.53
ignora il tuo semplicemente perche', contrariamente a quanto potrebbe sembrare dall'orario, lo ignorava.
E' un po' lunghetto, ma se ti/vi va di leggerlo, forse fa capire qualcosa dei miei problemi, quanto meno per la parte "informatica".
Venendo al tuo messaggio:
-----
> Vedi "Numerical Methods that (usually) Work" di Forman Acton.
-----
Grazie provero' a cercarlo e, compatibilmente al fatto che "je parle francais mais non anglais", vedro' di capirci qualcosa [oltre ai Newsgroup ho i miei angeli custodi in carne ed ossa].
Chiedere comunque qualcosa di italiano e' pretenzioso vero? [lo dico sinceramente, non c'e' retorica, risentimento o ironia]
------
> Comunque, spesso chi chiede di invertire le matrici non ha davvero
> bisogno di farlo. Cosa vuoi fare _davvero_? Magari puoi farlo senza
> invertire le matrici.
------
Ti/vi esporro' la parte matematica in notazione "matriciale sintetica". [perdonate la terminologia so' architetto che ce volete fa', e quanto a Tex e LaTex conosco solo i fumetti]
Ho un sistema omogeneo di (n) equazioni in (9) incognite.
K9 · X9 = N9
K9 (n · 9) = ((# # # # # # # # #) ... (# # # # # # # # #))
X9 (9 · 1) = ((x1)(x2)(x3)(x4)(x5)(x6)(x7)(x8)(x9))
N9 (n · 1) = ((0)...(0))
A questo punto per determinare "approssimativamente" le X pongo una ad una tali incognite come note (x1 = 1, . , x9 = 1); le tre matrici del sistema omogeneo di (n) equazioni diverranno in (8) incognite [ne avro' 9 ma trascuriamolo d'ora in poi]:
K8 · X8 = N8
- K8 (n · 8) = ((# # # # # # # #) ... (# # # # # # # #))
- X8 (8 · 1) = ((x1)(x2)(x3)(x4)(x5)(x6)(x7)(x8))
- N8 (n · 1) = ((#)...(#))
Posta:
- K8t (8 · n) = la trasposta di K8
E quindi poste:
- A8 (8 · 8) = K8t (8 · n) · K8 (n · 8)
- B8 (8 · 1) = K8t (8 · n) · N8 (n · 1)
Premoltiplicando la mia equazione " K8 · X8 = N8 " per " K8t " otterro':
A8 · X8 = B8
e quindi posta:
- A8i (8 · 8) = l'inversa di A8
Avro' due alternative:
X8 = A8i · B8
X8 = Cramer (A8 B8)
E prescindendo da cio' che ho esposto in sunto per hora ho scritto questi programmi:
X8 = Non_so_come (A8 B8)
X8 = So_ancora_meno_come_fare (K8 N8)
Vorei poter scriverne altri:
X8 = Grazie_Adam_o_chi_per_lui_1 (A8 B8)
X8 = Grazie_Adam_o_chi_per_lui_2 (K8 N8)
------
Tornando alla parte informatica, cio' che vorrei fare e': come minimo velocizzare cio' che ho gia' fatto e (come ho scritto nell'altro mio messaggio citato all'inizio) a tal proposito mi farebbe comodo lo sviluppo "simbolico" di detA (6 · 6) o se non e' utopico di detA (7 · 7).
Con Derive [da un mio amico "licenziato" che aveva molto tempo libero da dedicarmi] ho fatto il detA (5 · 5) e con santa pazienza l'ho messo in lisp, il problema è che non potendo o non riuscendo a settare come simboli "a11" "a12" ... ed essendoci dovuti limitare ad "a" "b" ... ho penato non poco a "raccogliere" ed "ordinare".
Quindi o voi matematici pazienti e gentili [si, un po' mi vergogno di leccare] ditemi, con quale programma posso fare (o voi per me se ...issimi) tale sviluppo, con "simboli" appropriati e magari ben "raccolto" ed "ordinato"?
------
Questo e' tutto
Grazie ad Adam,
Grazie da e Ciao da Giovanni
Adam Atkinson
>Chiedere comunque qualcosa di italiano e'
>pretenzioso vero?
No. E' solo che non ho quasi niente in italiano.
>> Comunque, spesso chi chiede di invertire le matrici non ha davvero
>> bisogno di farlo. Cosa vuoi fare _davvero_? Magari puoi farlo senza
>> invertire le matrici.
[snip]
Se vuoi risolvere delle equazioni simultanee, non e' obbligatorio
invertire una matrice.
Ti avverto che non so quasi niente dell'analisi numerica, ma ho
sofgliato cose come il libro di Acton in cui dicono spesso cose del
tipo "la regola di Cramer dovrebbe essere tolto dai programmi
liceali!"
Prova anche uno dei libri della serie "Numerical Recipes in
(linguaggio)"
--
Adam Atkinson (gh...@mistral.co.uk)
VOLCANO MISSING FEARED DEAD
Giovanni Rana
serve un algoritmo per risolvere numericamente un sistema di equazioni
lineari di n equazioni in n incognite, con il determinante diverso da zero
(se no non avrebbe senso parlare di regola di Cramer).
In tal caso, inizio dalla letteratura italiana sull'argomento, che è
vastissima:
- V. Comincioli, " Analisi Numerica", McGraw-Hill;
- F.Fontanella, A. Pasquali, "Calcolo Numerico", Pitagora Editrice;
- N.S. Bachalov, "Metodi Numerici", Editori Riuniti;
Ce ne sono molti altri, ma direi che per cominicare ti bastano questi: il
più semplice dovrebbe essere il Comincioli, che per inciso è l'unico che
ho letto ( mi baso più che altro su libri americani e inglesi ).
Volendo essere molto sintetici, l' idea è che per calcolare un determinante
devi fare o(n!) moltiplicazioni. Per cui è chiaro che la regola di Cramer,
per n che sale, diviene di nessuna utilità pratica ( l' utilità teorica è
invece indiscussa, perchè permette di dimostrare il teorema della
caratterizzazione del rango di una matrice ). I metodi alternativi si
dividono in metodi diretti e metodi iterativi: esempio di metodo diretto
è il metodo di Gauss, in cui riduci a scalini la matrice A dei coefficienti
tramite operazioni elementari avendo cura di trasfomare anche il
termine noto secondo quelle operazioni), dopodichè hai un sistema
di equazioni che risolvi per sostituizone partendo dall'ultima riga.
Tali metodi son molto più veloci o(n^3)), e si possono specializzare
se A è di forma particolare (es.tridiagonale), arrivando a o(3*n)
moltiplicazioni. Però hanno il grave difetto che , se A è qualunque,
allora gli errori di macchina possono accumularsi, essere amplificati,
e il computer ti restituisce valori tipo 10^5 per un'incognita che doveva
valere 1. Devo dire che , da quanto ho capito del tuo post, tu lavori
con matrici piccolissime, delle 9x9, mentre io trovo problemi
con tali metodi per delle 200x200: per cui magari tu potresti usare un
metodo diretto. Credo che tu sappia come ridurre a scalini una matrice,
dunque sai pure implementare l' algoritmo di Gauss, comunque se così
non fosse, puoi scrivere al NG e ti spedisco l'algoritmo.
L'altra classe di metodi sono i metodi iterativi, cioè metodi che, partendo
da una soluzione di tentativo, la modificano a mano a mano, invertendo in
pratica una matrice che non è A ma è "simile" ad A, anche se molto più
semplice da invertire, ed arrivano in tot passi alla soluzione. Anche
questi non convergono qualunque sia A, però certo funzionano molto
più spesso dei metodi diretti, e nel tuo caso (ordine "piccolo"), dovresti
avere ottimi risultati. Se hai il Fortran o il C, cerca nelle "Numerical
Recipes" il codice GAUSSJ se ti interessa un metodo diretto, o il codice
RELAX se vuoi usare un metodo iterativo. Se hai MATLAB, infine, ti
posso spedire il codice del metodo iterativo di Gauss-Siedel, che non è
il migliore dell'universo, però ha il notevole (secondo me) pregio di avere
una logica facilmente comprensibile anche da chi non ha studiato
Analisi Numerica.
Ciao
Giovanni Anzani
Grazie anche per il tuo secondo messaggio, da cui traggo:
> Se vuoi risolvere delle equazioni simultanee, non e' obbligatorio
> invertire una matrice.
dovresti se non ti secca dirmi cosa sono le eq. simult. perche' proprio non lo so, sapendolo potrei o rassegnarmi a dover invertire a sperare in soluzioni alternative.
Grazie Ciao Giovanni.
Giovanni Anzani
Innanzitutto grazie per il tempo dedicatomi e per la dovizia di particolari, come da te consentitomi torno a scriverti riferendomi al tuo messaggio precedente:
1) Hai scritto:
> o(n!)
> o(n^3)
Scusa la mia profonda ignoranza (ho gia' detto che sono un semi-architetto?)
cosa vuole dire la "o"?
2) Hai scritto:
> se non ho capito male ti serve un algoritmo per risolvere numericamente un
> sistema di equazioni lineari di n equazioni in n incognite, con il determinante
> diverso da zero (se no non avrebbe senso parlare di regola di Cramer).
Io ho un sistema omogeneo di (2*n) equazioni in (9) incognite [dove n<20 circa].
K9 · X9 = N9 dove :
K9 (2*n · 9) = ( (# # 1 0 0 0 # # #) (0 0 0 # # 1 # # #) ... (# # 1 0 0 0 # # #) (0 0 0 # # 1 # # #))
X9 (9 · 1) = ((x1)(x2)(x3)(x4)(x5)(x6)(x7)(x8)(x9))
N9 (2*n · 1) = ((0)(0)...(0)(0))
A questo punto per determinare "approssimativamente" le x pongo una ad una tali incognite come note (x1 = 1, . , x9 = 1) e contestualmente porto la relativa colonna di K9 ai termini noti N9; le tre matrici del sistema omogeneo di (2*n) equazioni diverranno 9 sistemi non omogenei in (8) incognite :
K8 · X8 = N8 e posti: K8t (la trasposta di K8), A8 = K8t · K8, B8 = K8t · N8
Premoltiplicando la mia equazione " K8 · X8 = N8 " per " K8t " otterro':
A8 · X8 = B8
A quanto sapevo sino ad ora avevo due alternative:
X8 = Inversa (A8) · B8
X8 = Cramer (A8 B8)
3) Hai scritto:
> Credo che tu sappia come ridurre a scalini una matrice,
> dunque sai pure implementare l'algoritmo di Gauss, comunque se così
> non fosse, puoi scrivere al NG e ti spedisco l'algoritmo.
La risposta, mi vergogno un po' a dirlo, e' doppiamente negativa, ti prego pertanto di insegnami ad usare un metodo diretto spedendomi tale l'algoritmo!
SE non ti spiace mi prendero' poi la briga di "tradurlo" in autolisp e di farmene una procedura per il mio programma che, compatibilmente con altre paternita' che tu dovessi segnalarmi si chiamera' a scelta:
MA-GAUSS-GR
MA-GAUSS-come_vuoi_tu'
4) Hai scritto:
> I metodi alternativi si dividono in metodi diretti e metodi iterativi...
> ... tu potresti usare un metodo diretto.
Concordo, non credo di dover passare alle iterazioni, avro' tra un po' a che fare con matrici 12x12, ma sono sempre piccole; inolre dopo il metodo di Cross a tecnica delle costruzioni vorrei nel mio caso evitare la convergenza per iterazione.
5) Hai scritto:
> Tali metodi [quelli diretti] si possono specializzare
> se A è di forma particolare (es.tridiagonale), arrivando a o(3*n)
> moltiplicazioni.
L'unica particolarita' delle mie 9 matrici A8 (8x8) e quella di derivare tutte
da una matrice A9 (nx9) formata da coppie di n/2 righe del tipo che riporto sotto:
r. n1 = (# # 1 0 0 0 # # #)
r. n2 = (0 0 0 # # 1 # # #)
a cui tolgo volta volta una colonna.
Non so se questo puo' dar luogo ad una particolarizzazione, ma sarebbe gia' qualcosa per me credo passare a Gauss se, da quello che ho capito (a parte la "o" di cui ti ho gia' chiesto), mi abbatterebbe un po' i tempi di calcolo anche in vista del passaggio ad un problema perfettamente analogo a quello che ti ho posto, ma con 12->11 incongnite in luogo delle attuali 9->8.
Tu che ti diletti con 200x200 sorriderai sicuramente a questa mia fissa dell'abbatimento di tempi dell'ordine dei 9x(4 sec.)+ spiccioli, ma tieni presente che le procedure che vado sviluppando devono essere integrate in programmi CAD di uso frequente e nel disegnare al Cad sarebbe un po' scomodo attendere.
(Ti dico questo perche', scusa la digressione, quando mi disegnai per Progettazione I tutta una piazza 3D con 486 16Mb ram furono dolori)
Per il momento mi fermo qui, mi collego, spedisco, sincronizzo, e metto me e il computer a nanna.
Ciao Giovanni
Adam Atkinson
puo' essere utile anche a terzi)
On 21-Apr-00 21:01:50, Giovanni Anzani said:
>> Se vuoi risolvere delle equazioni simultanee, non e' obbligatorio
>> invertire una matrice.
>dovresti se non ti secca dirmi cosa sono le eq. simult.
??!!!
O magari non si chiamano "equazioni simultanee" in italiano.
comunque....
qualcosa come
3x+4y+5z = 17
4x+2y+3x = 12
x-15y-z = 3
ma magari con piu' equazioni e piu' variabili.
> perche' proprio non
>lo so, sapendolo potrei o rassegnarmi a dover invertire a sperare in
>soluzioni alternative.
Si _puo'_ farlo invertendo una matrice, o usando la (lentissima)
regola di Cramer.
Pero', si puo' anche fare con l'eliminazione Gaussiana. Riscrivi una
delle equazioni nella forma x = (formula con y e z dentro). Poi
sostituisci questo per x nelle altre 2 equazioni. Poi ottieni da
un'altra equazione una formula per y in termini di z e lo sostituisci
nell'ultima equazione. Trovi z, e poi y, e poi x. Puoi generalizzare
questo procedimento per qualsiasi numero di variabili.
Probabilmente c'e' gia' un module Perl per fare questo cose. E
immagino che le facciano anche Mathematica e molte altre cose.
--
Adam Atkinson (gh...@mistral.co.uk)
It is a sobering thought, for example, that when Mozart was my age, he had
been dead for two years. (T. Lehrer)
Giovanni Anzani
Provero' ad "eliminare" le procedure sinora fatte buttandomi sul Gauss, ho capito spero tutto, tra l'altro con un po' di fortuna dovrei ricevere da Giovanni Rana un uovo di Pasqua con dentro tale algoritmo.
[ho visto altrove proteste sul virgolettare, e te virgolettare _cosi'_, faccio qualcosa di male a virgolettare "così"? se si scusami/scusatemi]
Hai scritto:
> puo' essere utile anche a terzi)
A proposito di terzi mi domando che fine abbia fatto Innuendo, lancia il sasso [buon per me che ho trovato lo spunto] e nasconde la mano, dico questo bonariamente, spero anzi che le votre esaurienti risposte siano state utili anche a lui.
Con questo ti/vi saluto e se non vi spiace tornerei al mio ruolo di lurker, in quanto su i.s.m ho molto da imparare ma poco da insegnare.
Ciao a Tutti Giovanni A.
Adam Atkinson
>[ho visto altrove proteste sul virgolettare, e te virgolettare _cosi'_,
>faccio qualcosa di male a virgolettare "così"? se si scusami/scusatemi]
Le tue virgolette sembrano normali. Usi un prodotto Microsoft? Forse
in passato hai usato "smart quotes" che sostituisce le virgolette con
dei caratteri non-standard che appaiono come quadretti su altri
sistemi. _Questo_ non e' 'virgolettare', e' il corsivo o grassetto.
--
Adam Atkinson (gh...@mistral.co.uk)
This .sig intentionally left blank.
Piero
se il sistema assume la forma
Ax=b
allora la soluzioneè x= inv(A)*b
e in matlab c'è anche il comando \.
x=A\ b. Che invece di usare l'inversa usa un metodo di eliminazione.
Poi una volta trovati i risultati che fai scrivere un file di testo formattaro e te lo
leggi dal tuo programma.
Potresti anche automatizzare il tuo tramite il DDE, ad esempio tramite una
funzione simile ad sendkeys di VBapplication.
Ti crei prima un file di testo con i dati in input, poi lo fai leggere a Matbal, calcoli i risultati
produci un file di testo con i risultati e lo rileggi con il programma client.
diamond
sui miei appunti di algebra lineare ;-)
Data una matrice Aij = {aij} con i=j e 1 <= i <=n.
Si calcoli il seguente determinante det(A-xI) e lo si ponga uguale a zero; x
è uno scalare e I è la matrice identica; dal calcolo risulta che
det(A-xI) = (-1)^n * x^n + (-1)^n-1 * Bn-1 * x^n-1 + ... + (-1) * B1 * x +
B0 = 0 (1)
dove Bn-1, ....., B1, B0 sono scalari che si ricavano dai conti.
A questo punto si applica il Teorema di Cayley-Hamilton: sostituendo nella
(1) la matrice A risulta ancora
(-1)^n * A^n + (-1)^n-1 * Bn-1 * A^n-1 + ... + (-1) * B1 * A + B0 * I = 0
(2)
dove I è ancora la matrice identica.
A questo punto si moltiplicano entrrambi i membri della (2) per A^-1
(matrice inversa di A) e facendo alcuni calcoli si trova la matrice inversa
di A
(sarà del tipo A^-1 = C1 * A^n-1 + C2 * A^n-2 + .... + Cn * I)
senza calcolare i minori complementari, ma facendo soltanto delle
moltiplicazioni della matrice A stessa.
Diamond
-----------
Remember when you were young, you shone like the sun.
Shine on you crazy Diamond.
Now there's a look in your eyes, like black holes in the sky.
Shine on you crazy Diamond.
....
However that may be, I'm linving in the dark.
-----------
"Innuendo" <tetragr...@SPAMMAREyahoo.it> ha scritto nel messaggio
news:8dnr6t$ne$1...@fe1.cs.interbusiness.it...
> Esiste un metodo rapido per invertire una matrice?
> Io finora usavo il metodo, data una matrice A, di moltiplicare 1/det(A)
per
> la trasposta dei suoi complementi algebrici, ma non mi sembra un
> procedimento molto rapido specialmente in caso di matrici che non siano le
> semplici 2x2....
>
> grazie in anticipo
> Innuendo
> --
Roberto Zanasi
wrote:
> O magari non si chiamano "equazioni simultanee" in italiano.
Noi di solito li chiamiamo "sistemi di equazioni".
Giovanni Rana
Giovanni Anzani <giova...@tiscalinet.it> wrote in message
8dqhkg$urt$2...@pegasus.tiscalinet.it...
>Caro Giovanni Rana,
>Innanzitutto grazie per il tempo dedicatomi e per la dovizia di
>particolari, come da te consentitomi torno a scriverti
>riferendomi al tuo messaggio precedente:
Ciao, vedo che hai molte domande, per cui il mio post sarà lunghetto,
come avverto nel subject.
>1) Hai scritto:
>> o(n!)
>> o(n^3)
>Scusa la mia profonda ignoranza (ho gia' detto che sono un
>semi-architetto?) cosa vuole dire la "o"?
Prego figurati, siamo in due ad essere "ignoranti". Difatti io, sbagliando,
ho usato il simbolo o( ) (o piccolo di Landau) che indica (nel nostro
caso in cui si ha a che fare con successioni divergenti) che una data
successione an ha ordine d'infinito rispettivamente SUPERIORE a n! o a
n^3, mentre avrei dovuto usare il simbolo O( ) ( O grande di Landau) che
invece indica che an ha ordine d'infinito non maggiore( nel nostro caso
uguale) a n! o a n^3.
In altre parole, per n-->+inf il metodo di Cramer richiede di fare circa n!
fra moltiplicazioni e divisioni, mentre il metodo di Gauss ne richiede
n^3/3+n^2- n/3, che ha lo stesso ordine d'infinito di n^3.
>2) Hai scritto:
>> se non ho capito male ti serve un algoritmo per risolvere numericamente
>>un sistema di equazioni lineari di n equazioni in n incognite, con il
>>determinante diverso da zero (se no non avrebbe senso parlare di regola
>>di Cramer).
>Io ho un sistema omogeneo di (2*n) equazioni in (9) incognite [dove n<20
[snip]
Ma...credo che ci eravamo già capiti. Tu parti da un sistema omogeneo
sovradeterminato( cioè in cui il numero delle equazioni è superiore a
quello delle incognite), e giungi a 9 sistemi di 8 equazioni ed 8 incognite
ciascuno, usando tra l'altro la prima trasformata di Gauss ( il fatto che
moltiplichi K8 e N8 per K8t, così da avere un sistema di 8 equazioni
in 8 incognite). Mi sbagliavo sull' ordine dei tuoi 9 sistemi ( ordine 8
invece di 9), ma il ragionamento non cambia, io mi riferivo ai
sistemi finali e non a K9 · X9 = N9.
Solo un' osservazione: il tuo sistema di partenza, essendo omogeneo e
sovradeterminato, ammetterà una o infinite soluzioni, ma delle quali puoi
determinare esattamente una base: perchè parli di soluzioni approssimate?
Io ti capirei se tu stessi cercando di risolvere un sistema non omogeneo
sovradeterminato, che spesso non ammette soluzione, e del quale dunque
ti cerchi una soluzione "approssimata" nel senso della minima norma
euclidea del residuo. Comunque resta il fatto che per risolvere 9 sistemi
quadrati Gauss ti aiuta molto, sempre che det(A8)<>0 ( condizione
necessaria anche per l'uso di Cramer o dell'inversione ): det(A8)<>0
se rango(K8)=8 per ogni K8.
>3) Hai scritto:
>> Credo che tu sappia come ridurre a scalini una matrice,
>> dunque sai pure implementare l'algoritmo di Gauss, comunque se così
>> non fosse, puoi scrivere al NG e ti spedisco l'algoritmo.
>La risposta, mi vergogno un po' a dirlo, e' doppiamente negativa, ti prego
>pertanto di insegnami ad usare un metodo diretto spedendomi tale
>l'algoritmo! SE non ti spiace mi prendero' poi la briga di "tradurlo" in
autolisp e di
>farmene una procedura per il mio programma che, compatibilmente con altre
>paternita' che tu dovessi segnalarmi si chiamera' a scelta:
>MA-GAUSS-GR
>MA-GAUSS-come_vuoi_tu'
Non c'è problema: eccoti l' algoritmo di Gauss, però ti rinnovo l'invito a
"spulciare" le Numerical Recipes, che anche se son carenti per i metodi
iterativi (non ci stanno nè BiGCstab nè GMRES, o io non riesco a trovarli,
ma ne dubito, ho recentemente imparato a leggere :-) ), vanno forte per i
metodi diretti. Per il nome, ti ringrazio per la correttezza, ma non c'è
bisogno che mi citi, l' algoritmo è molto noto e dovuto interamente a Gauss,
io faccio solo lo scriba. Ipotesi: det(A)<>0
Dette x_1, x_2, ...x_n le n incognite del tuo sistema numerate nell'ordine
in cui appaiono in tutte le equazioni del sistema, scambia fra loro la
equazione1 e la prima equazione n cui x_1 appare con coefficiente non nullo:
esempio,
0*x_1+ 2*x_2+ 3*x_3=0
3*x_1 + 6*x_2+0*x_3=5
...
A prescindere dall 'altra equazione, la prima equazione in cui x_1 compare
con coefficiente non nullo è la seconda, per cui scambi eq.1 con eq.2. In
generale, avrai fatto lo scambio e1<->ej. Adesso dividi (e1)' (cioè la
vecchia ej, che ora sta in prima posizione) per a(j,1)=a'(1,1)<>0: ottieni
in prima posizione un'equazione (e1)'' che ha come primo coefficiente 1.
Ciò fatto, per ogni equazione ei seguente alla prima:
- se a(i,1)<>0, ad ei sottrai (e1)''*a(i,1);
- se a(i,1)=0, ovviamente passi avanti;
Esempio: sempre il sistema di prima, dopo aver scambiato e1 con e2, ed aver
diviso (e1)' per a'(1,1), hai
1*x_1+2*x_2+ 0*x_3=5/3
0*x_1+ 2*x_2+ 3*x_3=0
2*x_1+0*x_2+1*x_3=1
Esamini ei=e2, vedi che a(i,1)=a(2,1)=0 e passi avanti. Esamini ei=e3, vedi
che a(i,1)=a(3,1)=2, per cui ad e3 sottrai (e1)''*2, e ottieni
(e3)'=0*x_1-4*x_2+1*x_3= -7/3. Insomma, come è chiaro, otttieni una A' in
cui la prima colonna ha il primo termine, che è il termine più in alto,
unitario (pivot) e tutti gli altri nulli, cioè x_1 appare (con coeff. non
nullo) solo nella prima equazione. E' poi ovvio come proseguire: cerchi la
prima equazione e'j fra e'2,....e'n in cui x_2 compare con coefficiente non
nullo, la scambi con e'2 e così via.
Alla fine ottieni una matrice Af dei coefficienti a scalini senza righe
nulle ( per ipotesi era det(A) <>0 ), e anche il termine noto è stato
trasformato. Nell'ultima equazione avrai
0*x_1+...+0*x_(n-1)+1*x_n=m =>x_n=m. Poi nella penultima
equazione appaiono solo x_(n-1) ed x_n, che però hai già
calcolato, ergo hai solo un' incognita, che ricavi, e così via.
>4) Hai scritto:
>> I metodi alternativi si dividono in metodi diretti e metodi iterativi...
>> ... tu potresti usare un metodo diretto.
>Concordo, non credo di dover passare alle iterazioni, avro' tra un po' a
>che fare con matrici 12x12, ma sono sempre piccole; inolre dopo il metodo
>di Cross a tecnica delle costruzioni vorrei nel mio caso evitare la
>convergenza per iterazione.
Metodo di Cross? Che metodo è? Mi puoi dare un riferimento bibliografico?
Son sempre alla ricerca di nuovi metodi iterativi per la soluzione di
sistemi lineari, per motivi di studio. E' un metodo stazionario o
instazionario? O usa il multigrid?
[snip]
>L'unica particolarita' delle mie 9 matrici A8 (8x8) e quella di derivare
>tutte da una matrice A9 (nx9) formata da coppie di n/2 righe del tipo
>che riporto sotto:
>r. n1 = (# # 1 0 0 0 # # #)
>r. n2 = (0 0 0 # # 1 # # #)
> cui tolgo volta volta una colonna.
Se le tue A8 non conservano questa struttura (ovviamente con 8
equazioni e non n), allora quest' informazione non aiuta.
>Tu che ti diletti con 200x200 sorriderai sicuramente a questa mia fissa
>dell'abbatimento di tempi dell'ordine dei 9x(4 sec.)+ spiccioli, ma tieni
>presente che le procedure che vado sviluppando devono essere integrate in
>programmi CAD di uso frequente e nel disegnare al Cad sarebbe un po'
scomodo attendere.
Non sorrido, anzi ti capisco. Il tuo è un bel casino, perchè per numeri
bassi i vari metodi "seri" non son tanto differenti, ergo non è facile
trovare quello ottimale: sarebbe ancora semplice se dovessi risolvere
un sol gruppo di 9 sistemi una volta per tutte, ma credo che tu
debba risolverne più gruppi durante uno stesso disegno. Purtroppo
non posso dedicare molto tempo al problema, comunque già così
passi dalle circa 11!=39916800 moltiplicazioni di Cramer alle circa
560 di Gauss.
>Ciao Giovanni
Ciao
Giovanni Rana
diamond <diam...@libero.it> wrote in message
NVhM4.19260$xt2.2...@news.infostrada.it...
> Propongo un altro modo per calcolare l'inversa di una matrice; l'ho
trovato
> sui miei appunti di algebra lineare ;-)
>
> Data una matrice Aij = {aij} con i=j e 1 <= i <=n.
>
> Si calcoli il seguente determinante det(A-xI) e lo si ponga uguale a zero;
x
> è uno scalare e I è la matrice identica; dal calcolo risulta che
>
> det(A-xI) = (-1)^n * x^n + (-1)^n-1 * Bn-1 * x^n-1 + ... + (-1) * B1 * x +
> B0 = 0 (1)
>
> dove Bn-1, ....., B1, B0 sono scalari che si ricavano dai conti.
>
> A questo punto si applica il Teorema di Cayley-Hamilton: sostituendo nella
> (1) la matrice A risulta ancora
>
> (-1)^n * A^n + (-1)^n-1 * Bn-1 * A^n-1 + ... + (-1) * B1 * A + B0 * I = 0
> (2)
>
> dove I è ancora la matrice identica.
>
> A questo punto si moltiplicano entrrambi i membri della (2) per A^-1
> (matrice inversa di A) e facendo alcuni calcoli si trova la matrice
inversa
> di A
>
> (sarà del tipo A^-1 = C1 * A^n-1 + C2 * A^n-2 + .... + Cn * I)
>
> senza calcolare i minori complementari, ma facendo soltanto delle
> moltiplicazioni della matrice A stessa.
>
>
>
> Diamond
Sì, se le cose fossero così semplici, effettivamente si farebbe come hai
detto tu: peccato che a norma di un noto teorema di algebra lineare, i
(-1)^i*Bi sono le somme dei minori principali di ordine i di A, che non
sono di calcolo tanto immediato. Si tratta difatti di calcolare determinanti
di matrici i x i: già se A è una 8x8 il calcolo diviene pesante anche per un
computer, ed assolutamente fuori dal mondo per un uomo. Se però A
è molto sparsa o è a banda stretta, allora è facile trovare i Bi: ancora,
il teorema di Cayley-Hamilton permette di invertire le matrici in forma
compagna, perchè in tal caso i Bi son gli elementi di una certa riga o
di una certa colonna, ma per il resto non è utile per matrici grosse.
Ciao
diamond
-----------
Remember when you were young, you shone like the sun.
Shine on you crazy Diamond.
Now there's a look in your eyes, like black holes in the sky.
Shine on you crazy Diamond.
....
However that may be, I'm linving in the dark.
-----------
"Giovanni Rana" <pan...@studenti.unina.it> ha scritto nel messaggio
news:8dsmvg$age$1...@news.unina.it...
>
> Sì, se le cose fossero così semplici, effettivamente si farebbe come hai
> detto tu: peccato che a norma di un noto teorema di algebra lineare, i
> (-1)^i*Bi sono le somme dei minori principali di ordine i di A, che non
> sono di calcolo tanto immediato. Si tratta difatti di calcolare
determinanti
> di matrici i x i: già se A è una 8x8 il calcolo diviene pesante anche per
un
> computer, ed assolutamente fuori dal mondo per un uomo. Se però A
> è molto sparsa o è a banda stretta, allora è facile trovare i Bi: ancora,
> il teorema di Cayley-Hamilton permette di invertire le matrici in forma
> compagna, perchè in tal caso i Bi son gli elementi di una certa riga o
> di una certa colonna, ma per il resto non è utile per matrici grosse.
>
> Ciao
>
>
Effettivamente invertire una matrice non è cosa semplice soprattutto se si
tratta di matrici molto "grosse"!!!!! Io volevo soltanto proporre uno modo
diverso per calcolarle. Nientaltro!
Buona Pasqua per chi ci crede. Ciao.
Innuendo
sasso
>[buon per me che ho trovato lo spunto] e nasconde la mano, dico questo
bonariamente, spero anzi che le votre esaurienti >risposte siano state utili
anche a lui.
Si ho "lurkato" la vostra interessantissima discussione, anche se ad un
certo pto è andata un po' al di fuori di quello che chiedevo dato che non
cercavo (cmq potrei avere un interesse in futuro) un algoritmo per PC. La
mia esigenza nasceva dal fatto che, dovendo fare prossimamente un esame di
Programmazione Lineare (quindi simplessi a go-go), ho la necessità di
invertire matrici velocemente. Il metodo che uso adesso può andare benissimo
per semplici matrici 2x2, ma son dolori a partire dalle 3x3 in su. Per
esempio se mi capita una base 4x4, per determinare le soluzioni del primale
e del duale (x=AB^-1*Bb,y=c*AB^-1) devo faticare non poco.
>Con questo ti/vi saluto e se non vi spiace tornerei al mio ruolo di lurker,
in
>quanto su i.s.m ho molto da imparare ma poco da insegnare.
Beh cmq spero che in questo ng rimanga spazio anche per chi ha da porre
questioni più "terra-terra" e non sia solo una "elite" per superlaureati in
matematica!!! :)
>Ciao a Tutti Giovanni A.
Ciao a tutti anche da parte mia.
I.
Giovanni Anzani
grazie per il tuo consiglio di e su come usare Matlab, quello che mi serve sarebbe proprio un bel "x=A\ b", ma aime' me lo devo fa' da me coll'autolisp 'sto metodo di eliminazione e ti dico perche':
A livello concettuale hai ragione da vendere, in sintesi (me lo dico da solo) e' una ca##ata usare un linguaggio di programmazione che non si presta come l'autolisp per lavorare con le matrici, senza scomodare Matlab, lo stesso excel "si comporta meglio", ma:
Cio' che sto facendo con queste matrici e' preparare la mia tesi di laurea in architettura e sarebbe meglio che lasciassi fuori da questo Matlab o altri programmi; gia' nella commissione potrebbero storcere il naso per il fatto che la mia tesi e' un programma e forse me la cavo perche' il "linguaggio" che uso e' quello di un CAD.
Un altra possibile contestazione potra' essere che rispetto ad alcune soluzioni gia' esistenti sul mercato, con la mia tesi devo "portarmi nello zaino" Autocad (e gia e' un bel "peso economico") e non vorrei alimentare tale critica infilandoci anche Matlab o chi per lui.
L'unica soluzione potrebbe essere un "*.exe" in c++ o altro, ma vorrei provare a restare in autolisp, perche' credo di far prima a svilupparmi l'eliminazione in lisp che non ad impararmi l'"interfacciamento verso l'esterno" dell'autolisp.
Quanto alle tue, per me, arcane parole:
> Potresti anche automatizzare il tuo tramite il DDE, ad esempio tramite una
> funzione simile ad sendkeys di VBapplication.
Spero di averti sufficientemente spiegato le mie ragioni.
Ti ringrazio, ti saluto
Ciao Giovanni
Giovanni Anzani
Nuovamente grazie per la tua risposta a cui tagliuzzando qua e la direttamente rispondo:
-----
> Prego figurati, siamo in due ad essere "ignoranti".
-----
-----
> Tu parti da un sistema omogeneo sovradeterminato [...],
> e giungi a 9 sistemi di 8 equazioni ed 8 incognite ciascuno,
> [...]
> sovradeterminato, ammetterà una o infinite soluzioni, ma delle quali puoi
> determinare esattamente una base: perchè parli di soluzioni approssimate?
> Io ti capirei se tu stessi cercando di risolvere un sistema non omogeneo
> sovradeterminato, che spesso non ammette soluzione, e del quale dunque
> ti cerchi una soluzione "approssimata" nel senso della minima norma
> euclidea del residuo.
non so cosa vuol dire: "ma delle quali puoi determinare esattamente una base"
non so cosa vuol dire: "minima norma euclidea del residuo"
non sapevo che da qualche parte facevo una cosa chiamata "prima trasformata di Gauss".
Alla luce di quanto detto sopra riguardo al capirsi ti rispondero' che di quello che faccio ho capito questo [scusa se quello che ti dico lo hai gia' capito]:
Premesse:
1) escluse in partenza combinazioni lineari.
2) le soluzioni che mi servono possono differire per un fattore di proporzionalita'.
3) le equazioni sono di numero pari e almeno 8.
Cio' premesso potro' avere 2 possibili casi:
- a) 2x4 equazioni 9 incognite sistema omogeneo
alla luce di 2) ho una soluzione per me "univoca" si potrebbe dire "banale" e quindi lasciamola perdere.
- b) 2x[n>4] eqazioni 9 incogite sistema omogeneo
Io avevo pensato di prendere le eq. 8 a 8 con ogni possibile permutazione, procedere come se fossi nel caso a) e poi "mediare il risultato" salvo stabilire come.
Mi e' stata invece indicata (molto meglio di come qua la riporto) questa strada:
poni le incognite una ad una =1 sostituendole nella rispettiva colonna della matrice (2n)x9 ottenendo cosi' dei termini noti "da portare di la'".
Passo cosi' ad un sistema:
- b2) 2x[n>4] eqazioni 8 incogite sistema non omogeneo.
Premoltiplico tutto per la trasposta non so a quale scopo (prima trasformata di Gauss ?).
Questa operazione mi pare di capire che "ottimizzi" qualcosa, ma non ho capito cosa, comunque ottengo un sistema:
- b3) 8 eqazioni 8 incogite sistema non omogeneo
che avrà soluzione univoca e ben determinabile e questo lo so fare pure io.
Il punto e che cosi' operando passo da questa fase della moltiplicazione per la trasposta (di cui non ho capito il significato) e che inoltre ottengo 9 "situazioni" b3) e quindi da cio' deducevo che in quanto non esisteva un'unica soluzione nessuna di esse era esatta, per questo parlavo, forse con un termine inadatto, di soluzioni approssimate.
----
> [...] eccoti l' algoritmo di Gauss, però ti rinnovo l'invito a "spulciare" le Numerical Recipes
> [...].
grazie per l'algoritmo, me lo studiacchio subito; quanto alle Numerical Recipes, mi ero dimenticato di dirti che non ho ne il Fortran ne il C.
N.B. Per chi fosse interessato Giovanni R. in un messaggio precedente scriveva:
> [...] Se hai il Fortran o il C, cerca nelle "Numerical Recipes" il codice GAUSSJ [...]
Quanto alla spiegazione dell'algoritmo, preferisco leggermelo e provarci un po' prima di chiedertene qualche delucidazione ulteriore, ma eventualmente sappi fin d'ora che tale eventualita' sara da inputarsi escusivamende alla mia ignoranza e non alla chiarezza della tua spiegazione (che dici esagero?). Salamelecchi a parte, nei prossimi 4-5 giorni, in cui mio malgrado, non credo che potro' leggere altro che il newsgroup it.all-newsgroup.off-line (causa matrimonio fuori regione di amici) vedro' di approfittarene per rimettere il autolisp tale procedura almeno concettualmente. La tua spiegazione c'e, la stampante c'e, 1+1=3 mi portero' una stampina e carta e penna si travano pure nelle toilette degli eurostar (questa e' proprio brutta).
----
> Metodo di Cross? Che metodo è? Mi puoi dare un riferimento bibliografico?
> Son sempre alla ricerca di nuovi metodi iterativi per la soluzione di
> sistemi lineari, per motivi di studio. E' un metodo stazionario o
> instazionario? O usa il multigrid?
----
> così passi dalle circa 11!=39916800 moltiplicazioni di Cramer alle circa
> 560 di Gauss.
Ti faro' sapere l'abbattimento dei tempi ottenuto.
----
Ciao Giovanni A.
Ancora grazie