Inversione matrice grande

Question

Sto osservando l'inverso di una matrice grande, dimensione comune di 1000 x 1000, ma a volte supera 100000 x 100000 (che attualmente non funziona a causa del tempo e della memoria). So che il sentimento normale è "non prendere l'inverso, trovare un altro modo per farlo", ma al momento non è possibile. La ragione di ciò è dovuta all'utilizzo del software già realizzato che si aspetta di ottenere la matrice inversa. (Nota: sto cercando dei modi per cambiare questo, ma ci vorrà molto tempo)

Al momento stiamo usando un metodo di decomposizione LU da recopie numeriche, e attualmente sto testando la libreria di eigen . La libreria eigen sembra essere più stabile e un po 'più veloce, ma sono ancora in fase di test per la precisione. Ho dato una rapida occhiata ad altre librerie come ATLAS e LAPACK ma non ho ancora fatto test sostanziali con questi.

Sembra che la libreria eigen non utilizzi metodi concorrenti per calcolare l'inverso (sebbene per la parte di fattorizzazione LU dell'inverso) e per quanto posso dire ATLAS e LAPACK sono simili in questa limitazione. (Sto attualmente testando la differenza di velocità per eigen con openMP e senza.)

La prima domanda è: qualcuno può spiegare come sarebbe possibile ottimizzare l'inversione della matrice mediante la parallelizzazione. Ho trovato un articolo here che parla di algoritmi paralleli di inversione di matrice, ma non ho capito. Sembra che l'articolo this parli di un altro metodo? Non sono nemmeno sicuro se scaLapack o PETSc sono utili?

Seconda domanda, ho letto l'articolo this di utilizzare le GPU per aumentare le prestazioni, ma non ho mai codificato per GPU e quindi non ho idea di cosa stia cercando di trasmettere, ma i grafici in basso sembravano piuttosto allarmanti. Com'è possibile, e come posso iniziare a implementare qualcosa di simile se è vero?

Ho trovato anche l'articolo this, ho ancora avuto il tempo di leggerlo per capire, ma sembra promettente, poiché la memoria è un problema attuale con il nostro software.

Qualsiasi informazione su questi articoli o sui problemi in generale sarebbe di grande aiuto. E di nuovo mi scuso se questa domanda sembra vaga, cercherò di espandermi di più se necessario.

Answer 1

La prima domanda è: qualcuno può spiegare come sarebbe possibile ottimizzare l'inversione della matrice mediante la parallelizzazione.

Mi azzardo a indovinare che questo e gli argomenti correlati in algebra lineare sono uno degli argomenti più studiati nel calcolo parallelo. Se sei bloccato alla ricerca di un posto dove iniziare a leggere, il buon vecchio Golub and Van Loan ha un capitolo sull'argomento. Se Scalapack e Petsc siano probabilmente utili, sicuramente il primo, probabilmente il secondo. Naturalmente, entrambi dipendono da MPI, ma questo è un dato dato per scontato in questo campo.

Seconda domanda ...

GPU utilizzare se li hai e ti puoi permettere di tradurre il codice nel modello di programmazione supportato dai vostri GPU. Se non hai mai codificato per le GPU e hai accesso a un cluster di CPU di tipo commodity, ti divertirai più rapidamente usando il cluster piuttosto che lottando con una nuova tecnologia.

Per quanto riguarda l'ultimo articolo a cui fai riferimento, ora ha 10 anni in un campo che cambia molto rapidamente (prova a trovare un documento di ricerca di 10 anni sull'uso delle GPU per l'inversione della matrice). Non posso commentare la sua eccellenza o altri attributi, ma le dimensioni dei problemi che hai menzionato mi sembrano ben all'interno delle capacità dei moderni cluster per il calcolo in-core (per usare un termine antico). Se le tue matrici sono molto grandi, sono anche sparse?

Infine, sostengo con forza la tua intenzione evidente di utilizzare codici già esistenti già esistenti piuttosto che provare a sviluppare il tuo.

Answer 2

100000 x 100000 è 80 GB a doppia precisione. È necessaria una libreria che supporti matrici mappate sulla memoria su disco. Non posso raccomandare una libreria particolare e non ho trovato nulla con ricerche di Google veloci. Ma il codice delle ricette numeriche certamente non sarà adeguato.

Answer 3

Per quanto riguarda la prima domanda (come parallellize calcolando l'inverso):

Presumo che si sta calcolando l'inversa facendo una decomposizione LU della matrice e quindi utilizzando la decomposizione di risolvere A * B = I dove A è la tua matrice originale, B è la matrice per cui risolvi e io sono la matrice di identità. Quindi B è l'inverso.

L'ultimo passaggio è facile da eseguire in parallelo. Dividi la tua matrice identità lungo le colonne. Se hai CPU p e la tua matrice è n-by-n, allora ogni parte ha colonne n/p e n righe. Chiamiamo le parti I1, I2, ecc. Su ogni CPU, risolve un sistema del modulo A * B1 = I1, questo ti dà le parti B1, B2, ecc., E puoi combinarle per formare B che è l'inverso .

Answer 4

Una decompressione LU su una GPU può essere ~ 10 volte più veloce rispetto a una CPU. Sebbene ciò stia cambiando, le GPU sono state tradizionalmente progettate attorno all'aritmetica a precisione singola, e così su hardware più vecchio l'aritmetica a precisione singola è generalmente molto più veloce dell'aritmetica a doppia precisione. Inoltre, i requisiti di archiviazione e le prestazioni saranno fortemente influenzati dalla struttura delle vostre matrici. Una decompilazione LU di matrice 100.000 x 100.000 sparse è un problema ragionevole da risolvere e non richiederà molta memoria.

A meno che non si voglia diventare specialisti e spendere un molto di tempo per l'aggiornamento dell'hardware, consiglio vivamente l'uso di una libreria commerciale. Suggerirei CULA tools. Dispongono di librerie GPU sparse e dense e in effetti il ​​loro free library offre SGETRF - una routine di decompilazione LU a singola precisione (densa). Dovrai pagare per le loro librerie a doppia precisione.

Answer 5

So che è il vecchio post - ma in realtà - OpenCL (si scarica quello pertinente basato sulla scheda grafica) + OpenMP + Vectorization (non in questo ordine) è la strada da percorrere.

In ogni caso, per me la mia esperienza con Matrix ha a che fare con overheads dalla duplicazione di double double array dentro e fuori il sistema e anche per riempire o inizializzare le matrici con 0 prima di iniziare qualsiasi computazione - specialmente quando sto lavorando con la creazione di .xll per l'utilizzo di Excel.

Se dovessi reprioritize cima -

1. cercare di vettorizzare il codice (Visual Studio 2012 e Intel C++ ha autovettorizazione - io non sono sicuro di MinGW o GCC, ma penso che ci sono bandiere per il compilatore per analizzare i loop for per generare i giusti codici assembly da utilizzare al posto dei normali registri per contenere i dati, per popolare i registri vettoriali del processore. Penso che Excel stia facendo ciò perché quando ho monitorato i thread di Excel mentre eseguivo MINVERSE() , Noto solo 1 thread. Non conosco molto linguaggio assembly - quindi non so come vettorializzare manualmente ... (non ho avuto il tempo di imparare questo ancora, ma mi piacerebbe farlo!)
2. Parallelamente con OpenMP (omp pragma) o libreria MPI o pthreads (parallel_for) - molto semplice - ma ...ecco il trucco - mi rendo conto che se la tua classe matrix è completamente single threaded in primo luogo - quindi parallelizzare l'operazione come mat multiply o inverse è scrappable - il parallelismo cuz peggiorerà la velocità a causa dell'inizializzazione o della copia o dell'accesso solo al non- classe di matrice parallelizzata. Ma ... dove la parallelizzazione aiuta è - se stai progettando la tua classe matrix e parallelizzi la sua operazione di costruzione (padding con 0s ecc.), Allora il tuo calcolo di LU (A^-1) = Sarò anche più veloce. È anche matematicamente semplice ottimizzare anche la scomposizione LU e ottimizzare la sostituzione inversa all'indietro per il caso speciale dell'identità. (Cioè non perdere tempo creando una matrice di identità - analizza il tuo per (riga = colonna) e valuta di essere una funzione con 1 e il resto con 0.)
3. Una volta che è stato parallelizzato (sugli strati esterni) - le operazioni di matrice che richiedono elemento per elemento possono essere mappate per essere calcolate da GPU (SSSSSS) - centinaia di processori per calcolare gli elementi - battere !. C'è ora codice Monte Carlo di esempio disponibile sul sito Web di ATI - utilizzando OpenCL di ATI - non preoccuparti di portare il codice a qualcosa che usa GeForce - tutto ciò che devi fare è ricompilare lì.

Per 2 e 3 anche se - ricorda che le spese generali sono sostenuti in modo che nessun punto a meno che non si sta movimentazione F * K * G matrici enormi - ma vedo 100k^2? wow ...

Gene