View Full Version : [C/C++] Consigli ottimizzazione velocità
Come da oggetto vorrei dei consigli per ottimizzare la velocità del mio programma, principalmente il "collo di bottiglia" del mio programma è la valutazione del seno e del coseno di una certa quantità, (quantità che varia volta per volta ovviamente),funzioni che occupano oltre l'80% del tempo processore. (Ovviamente le funzioni seno e coseno sono richiamate moooooooolte volte :muro: )
Il programma è già multithreaded, quindi sfrutta più core contemporaneamente.
Inoltre se possibile sono graditi dei consigli sull'allocazione e deallocazione rapida di vettori di dimensioni non particolarmente elevate qualche decina di elementi per intenderci, ma ciò apporterebbe solo marginali incrementi prestazionali.
Non so se è possibile utilizzare qualche libreria in maniera facile adatta allo scopo, per effettuare le operazioni trigonometriche.
SO Windows 7 RC Build 7100
IDE Visual Studio Team System 2010 Beta 1 (mi da un incremento di prestazioni rispetto a Visual Studio 2008 di circa il 20%)
Sicuro che non puoi usare qualcosaltro?
Per esempio il coseno dell'angolo tra 2 vettori si trova col prodotto scalare, che è decisamente più veloce di qualsiasi funzione trigonometrica.
Altro non saprei proprio, senza sapere codice e algoritmo ci vorrebbe la palla di vetro :asd:
Sicuro che non puoi usare qualcosaltro?
Per esempio il coseno dell'angolo tra 2 vettori si trova col prodotto scalare, che è decisamente più veloce di qualsiasi funzione trigonometrica.
Altro non saprei proprio, senza sapere codice e algoritmo ci vorrebbe la palla di vetro :asd:
Innanzitutto grazie della risposta, comunque praticamente devo scindere un'esponente con argomento complesso nella parte reale ed immaginaria, quindi valuto il seno e coseno. Ad ogni modo posto la parte di codice di interesse che così è tutto più chiaro e non serve la palla di cristallo ;)
void derivsU(double x ,double y[],double dydx[],double lam)
{
/* La parte iniziale del vettore contiene i dati legati alla parte
* reale dell'ampiezza, mentre la seconda parte contiene la parte
* immaginaria.
*/
int i,j;
int n=N;
double q1,q2,esp;
double K;double Kmod=(1+cos(2*pi/lam*x));
for(i=1;i<=n;i++){
dydx[i]=0;
dydx[i+n]=0;
for(j=1;j<=n;j++){
K=Kmod*intePerParam[i][j];
esp=x*beta[i][j];
q1=K*sin(esp);
q2=K*cos(esp);
dydx[i]-=q1*y[j]+q2*y[j+n];
dydx[i+n]+=q2*y[j]-q1*y[j+n];
}
}
}
Vorrei sottolineare come abbia provato anche altre soluzioni, ad esempio q1*=sin(esp) e q2*=cos(esp) ma i risultati sono analoghi, quindi il problema è proprio della valutazione del seno e del coseno e non del prodotto (che dovrebbe tra le funzioni più veloci per un processore dopo somma e sottrazione)
La funzione che ho riportato è richiamata anche oltre 3.000.000 di volte (è un dato indicativo), mentre N può assumere valore che spaziano fra qualche decina a 100-200
Ciò che mi lascia perplesso inoltre è che le funzioni che occupano la maggior parte del tempo processore sono sin(esp) e cos(esp) mentre cos(2*pi/lam*x) no, risulta trascurabile (almeno stando all'analizzatore di performance di Visual Studio)
Dunque, mi sembra già ottimizzato a bestia :asd:
poi vista la complicatezza dei calcoli che svolgi, scommetto che non puoi scambiare imprecisione per velocità no?
Tipo che so, sostituire seno e coseno con una table di valori, calcolare meno valori e poi sfocarli ecc ecc.
L'unica cosa da fare sarebbe riorganizzare l'algoritmo in modo di svolgere meno volte la parte critica, comunque non capisco nemmeno che stai facendo (dannati numeri complessi :asd:) quindi non saprei come...
Comunque stai facendo un filtro su un'immagine praticamente, quindi ti consiglierei CUDA, se ben usato queste robe se le mangia (a precisione singola però).
Ciò che mi lascia perplesso inoltre è che le funzioni che occupano la maggior parte del tempo processore sono sin(esp) e cos(esp) mentre cos(2*pi/lam*x) no, risulta trascurabile (almeno stando all'analizzatore di performance di Visual Studio)
Ciao,
io ti consiglio di introdurre nel tuo codice le liberie ACML, utilizzando la funzione vettoriale vrda_sincos. In pratica, raggruppi in un vettore tutti gli elementi su cui vuoi fare seni e coseni e con UNA SOLA chiamata li calcoli tutti di un botto, sfruttando anche le proprietà trigonometriche. ;)
Nel tuo caso dei raggruppare i vari esp
Link utili:
http://developer.amd.com/cpu/Libraries/acml/Pages/default.aspx
http://developer.amd.com/cpu/Libraries/acml/onlinehelp/Documents/vrda_005fsincos.html#vrda_005fsincos
PS: Perchè nei cicli for parti da 1 e non da 0 ?
Dunque, mi sembra già ottimizzato a bestia :asd:
poi vista la complicatezza dei calcoli che svolgi, scommetto che non puoi scambiare imprecisione per velocità no?
Tipo che so, sostituire seno e coseno con una table di valori, calcolare meno valori e poi sfocarli ecc ecc.
L'unica cosa da fare sarebbe riorganizzare l'algoritmo in modo di svolgere meno volte la parte critica, comunque non capisco nemmeno che stai facendo (dannati numeri complessi :asd:) quindi non saprei come...
Comunque stai facendo un filtro su un'immagine praticamente, quindi ti consiglierei CUDA, se ben usato queste robe se le mangia (a precisione singola però).
Si diciamo che l'ho ottimizzato per benino, quando ho avuto il programma ci impiegava 16 volte di più per fare le stesse cose, poi ho fatto qualche ottimizzazione... :sofico:
Comunque non si tratta di un filtro su un'immagine, è l'applicazione di una teoria dei campi elettromagnetici.
La possibilità di utilizzare una look-up table non credo faccia al caso mio, comunque grazie per il suggerimento
Ciao,
io ti consiglio di introdurre nel tuo codice le liberie ACML, utilizzando la funzione vettoriale vrda_sincos. In pratica, raggruppi in un vettore tutti gli elementi su cui vuoi fare seni e coseni e con UNA SOLA chiamata li calcoli tutti di un botto, sfruttando anche le proprietà trigonometriche. ;)
Nel tuo caso dei raggruppare i vari esp
Link utili:
http://developer.amd.com/cpu/Libraries/acml/Pages/default.aspx
http://developer.amd.com/cpu/Libraries/acml/onlinehelp/Documents/vrda_005fsincos.html#vrda_005fsincos
PS: Perchè nei cicli for parti da 1 e non da 0 ?
Ma il mio processore è Intel, non so se vadano bene ugualmente le librerie della AMD.
La domanda dell'indice di partenza è una bella domanda, praticamente il programma mi è stato dato e poi io l'ho ottimizzato (pesantemente, ma togliendo il collo di bottiglia del seno e coseno potrei aumentare la velocità di altre 4-5 volte immagino) e nel programma di partenza tutti gli array sono utilizzati con indice di partenza 1 (sprecando quindi la memoria per il primo elemento)
Ma il mio processore è Intel, non so se vadano bene ugualmente le librerie della AMD
Si, vanno benissimo ;) Anche se perdi un pò di performance rispetto a processori AMD, guadagni comunque molto rispetto al codice originale. RIguardo l'indice, va bene anche se parti da 1,alla fine sprechi un elemento, va bè. Solo che è molto facile confondersi poi..
Si, vanno benissimo ;) Anche se perdi un pò di performance rispetto a processori AMD, guadagni comunque molto rispetto al codice originale. RIguardo l'indice, va bene anche se parti da 1,alla fine sprechi un elemento, va bè. Solo che è molto facile confondersi poi..
Tu per caso le hai già sfruttate in qualche ambito? che miglioramenti hai ottenuto?
(Grazie comunque del suggerimento vedrò di utilizzare la libreria :) )
Tu per caso le hai già sfruttate in qualche ambito? che miglioramenti hai ottenuto?
(Grazie comunque del suggerimento vedrò di utilizzare la libreria :) )
Le sto utilizzando proprio in questo periodo ed i miglioramenti sono di un buon 25%. Ma tutto dipende dal numero di seni e coseni che devi calcolare. Nel tuo caso ad esempio, se n è molto piccolo, non avrai un grande vantaggio.
Le sto utilizzando proprio in questo periodo ed i miglioramenti sono di un buon 25%. Ma tutto dipende dal numero di seni e coseni che devi calcolare. Nel tuo caso ad esempio, se n è molto piccolo, non avrai un grande vantaggio.
Un 25% non si butta di certo via :D
Comunque N varia da qualche decina a 100-200, per molto piccolo su che ordine di grandezza siamo?
Un 25% non si butta di certo via :D
Comunque N varia da qualche decina a 100-200, per molto piccolo su che ordine di grandezza siamo?
Si va bene direi. Per i piccolo intendo intorno a 10. Comunque nel tuo caso, siccome dici che la funzione derivsU la chiami molte volte, devi fare N*M, dove M è il numero di volte che chiami tale funzione.
Un altro mio consiglio è di fare un unrolling del loop più interno con passo 2. Proprio in questi giorni ho a che fare con una funzione chiamata moltissime volte. Un unrolling con passo 2 mi ha fatto raddoppiare la velocità di esecuzione.
Esempio:
for (i = 0; i < 100; i++)
a += vect[i];
//Unroll passo 2:
for (i = 0; i < 100; i += 2)
{
a += vect[i];
a += vect[i+1];
}
Ovviamente devi gestire il resto nel caso che il numero di cicli non sia divisibile per 2.
Fatto con il passo2, prova con 4, ed 8.
Ultimo consiglio, dichiara double dydx[] come restrict ( double* restrict dydx) e metti -restrict nei flags di compilazione.
anche io ultimamente mi sono sbattuto parecchio per l'ottimizzazione di qualche programmino numerico in ambito scientifico e ho notato che la miglior ottimizzazione si ottiene solamente utilizzando l'algoritmo migliore, ma immagino che tu sia già arrivato a questo punto ;)
prima di tutto immagino che esp possa avere qualsiasi valore giusto? perchè se (come nel mio caso ad esempio) sai che può avere solo valori molto bassi puoi tranquillamente utilizzare gli sviluppi in serie di sin e cos. Se, come credo, questo non è il tuo caso puoi provare a calcolare solo il sin (o il cos) e poi usare sqrt(1- sin^2), io, prima di rendermi conto che effettivamente lo sviluppo era la cosa migliore da fare, avevo ottenuto un discreto miglioramento (con un test "sintetico" veniva fuori un 30% se non sbaglio).
_Claudio
21-07-2009, 01:37
Hai provato a cercare qualche libreria che calcola seno e coseno in tempi più rapidi? Immagino usi la math.h...
Altrimenti l'unica soluzione è rivedere le chiamate perchè chiamare 3000000 di volte una funzione è veramente pesante (e nelle ottimizzazioni si parte sempre dal numero maggiore, dannate derivate parziali)... che so magari con qualche stratagemma riesci a dimezzare, in fin dei conti i calcoli nel campo complesso derivano dalla circonferenza goniometrica che è simmetrica rispetto l'origine quindi in sostanza cambiano i segni ma le "forme" nei 4 pezzetti di piano sono sempre le stesse...
Se beta viene aggiornato meno spesso di quanto tu calcoli derivsU, puoi salvarti il seno e coseno di beta dentro un'altro array (caching), avendo quindi un costo costante per ogni chiamata di derivsU.
Inoltre mi pare ci siano istruzioni vettoriali (SIMD) che calcolano 4 seni e coseni in una sola volta, potrebbero interessarti.
Ringrazio tutti delle risposte.
Ho provato ad effettuare l'unroll del ciclo più esterno (un unroll completo, ho fatto giusto una prova su un caso particolare, un unroll di 64) ma le prestazioni non miglioravano (anzi peggioravano leggermente, anche se ciò non ha molto senso)
per l'opzione restrict, in visual studio la parola chiave è __restrict ma se la uso da errore, inoltre su MSDN per utilizzare restrict è indicata la modifica da fare al compilatore ma tale opzione non è presente nel mio.... (eppure ho consultato la sezione della versione che utilizzo di Visual studio )
per quanto riguarda la possibilità di utilizzare cos(x)=sqrt(1-sin(x)^2) in questo caso perderei l'informazione legata al segno, oppure dovrei fare degli if, ma a quanto ho letto gli if rallentano l'esecuzione del codice.
Per quanto riguarda le librerie, in realtà avevo provato ad utilizzare Intel Compiler, che aveva Math kernel Library, avevo utilizzato su macchina virtuale. Il problema è che non è gratuito il compilatore dell'Intel nè le librerie, sarebbe un po' scocciante ogni tot giorni formattare e reinstallare oppure ricopiare un'immagine della partizione. Si potrebbe pensare di compilare su macchina virtuale ed eseguire sull'host ma il problema è che su macchina virtuale non posso fare l'analisi delle performance (non sono emulate delle cose che non ricordo esattamente).
Per quanto riguarda le librerie dell'AMD, le ho scaricate installate, l'istruzione che è riportata sulla guida
cl driver.c -Ic:/acml4.3.0/win64_mp/include
c:/acml4.3.0/win64_mp/lib/libacml_mp_dll.lib
non mi funziona da il seguente errore
driver.c
c1: fatal error C1083: Cannot open source file: 'driver.c': No such file or directory
Ho provato a inserire gli header manualmente e le librerie ma purtroppo non sono molto pratico, e la compilazione termina (se uso una delle funzioni del file di intestazione) con un classico unresolved external symbol
:help:
Ho provato a seguire qualche guida per utilizzare le librerie ma niente da fare
Se qualcuno può indicarmi come utilizzare le librerie ACML, o indicarmi qualche link o guida :)
_Claudio
22-07-2009, 20:39
È chiaro che comunque c'è una discrepanza computazionale, la richiesta di 3000000 di chiamate a funzione non è eseguibile a bassa latenza su un normale computer.
O prevedi il fatto che quando si eseguono i calcoli vi è un tempo di attesa oppure non rimane che adottare tecnologie hardware più potenti o addirittura cluster di computer, se le finanze lo permettono.
Una domanda: i calcoli devono essere svolti per forza sequenzialmente?
Non è possibile dividere quei 3000000 di chiamate su più thread (core)?
Ringrazio tutti delle risposte.
Ho provato ad effettuare l'unroll del ciclo più esterno (un unroll completo, ho fatto giusto una prova su un caso particolare, un unroll di 64) ma le prestazioni non miglioravano (anzi peggioravano leggermente, anche se ciò non ha molto senso)
Ti avevo suggerito quello più interno :) Quando fai l'unroll le prestazioni possono peggiorare, tutto dipende dalla cache. Però in generale, siccome togli un sacco di if alla fine di ogni ciclo, si migliora. Inoltre un unroll di 64 è esagerato,( gestire i resti da 1 a 63 è pazzesco !) infatti ti avevo suggerito di partire con 2,4,8 ecc...
l'opzione restrict, in visual studio la parola chiave è __restrict ma se la uso da errore, inoltre su MSDN per utilizzare restrict è indicata la modifica da fare al compilatore ma tale opzione non è presente nel mio.... (eppure ho consultato la sezione della versione che utilizzo di Visual studio )
Purtroppo non conosco Visual Studio,ma comunque restrict è presente nello standard C99, prova a guardare se è già configurato per tale standard.
Per quanto riguarda le librerie, in realtà avevo provato ad utilizzare Intel Compiler, che aveva Math kernel Library, avevo utilizzato su macchina virtuale. Il problema è che non è gratuito il compilatore dell'Intel nè le librerie, sarebbe un po' scocciante ogni tot giorni formattare e reinstallare oppure ricopiare un'immagine della partizione. Si potrebbe pensare di compilare su macchina virtuale ed eseguire sull'host ma il problema è che su macchina virtuale non posso fare l'analisi delle performance (non sono emulate delle cose che non ricordo esattamente).
Si, infatti Tra ACML e MKL ti ho consigliato le prime proprio perché gratuite.
Per quanto riguarda le librerie dell'AMD, le ho scaricate installate, l'istruzione che è riportata sulla guida
cl driver.c -Ic:/acml4.3.0/win64_mp/include
c:/acml4.3.0/win64_mp/lib/libacml_mp_dll.lib
non mi funziona da il seguente errore
driver.c
c1: fatal error C1083: Cannot open source file: 'driver.c': No such file or directory
Ho provato a inserire gli header manualmente e le librerie ma purtroppo non sono molto pratico, e la compilazione termina (se uso una delle funzioni del file di intestazione) con un classico unresolved external symbol
L'errore è chiaro, avrai qualche path sbagliato. Ma anche su questo non posso aiutarti non essendo pratico di Visual Studio.
È chiaro che comunque c'è una discrepanza computazionale, la richiesta di 3000000 di chiamate a funzione non è eseguibile a bassa latenza su un normale computer.
O prevedi il fatto che quando si eseguono i calcoli vi è un tempo di attesa oppure non rimane che adottare tecnologie hardware più potenti o addirittura cluster di computer, se le finanze lo permettono.
Una domanda: i calcoli devono essere svolti per forza sequenzialmente?
Non è possibile dividere quei 3000000 di chiamate su più thread (core)?
Il problema non è tanto la latenza fra una chiamata e l'altra quanto l'esecuzione del seno e del coseno, stando all'analizzatore di Perforance di Visual Studio oltre l'80% del tempo processore è occupato nell'esecuzione di queste due funzioni. Ho già previsto il programma della possibilità di sfruttare multithreading, che ovviamente dimezza esattamente il tempo di esecuzione (il mio processore è Dual Core) ma basta cambiare un parametro e può lanciare anche 3,4,5 .... thread tutti in parallelo e sincronizzati. Il mio obiettivo è quindi cercare di utilizzare dei metodi per ridurre il costo computazionale del calcolo del seno e del coseno, utilizzando magari qualche libreria.
_Claudio
23-07-2009, 18:41
Il problema non è tanto la latenza fra una chiamata e l'altra quanto l'esecuzione del seno e del coseno, stando all'analizzatore di Perforance di Visual Studio oltre l'80% del tempo processore è occupato nell'esecuzione di queste due funzioni. Ho già previsto il programma della possibilità di sfruttare multithreading, che ovviamente dimezza esattamente il tempo di esecuzione (il mio processore è Dual Core) ma basta cambiare un parametro e può lanciare anche 3,4,5 .... thread tutti in parallelo e sincronizzati. Il mio obiettivo è quindi cercare di utilizzare dei metodi per ridurre il costo computazionale del calcolo del seno e del coseno, utilizzando magari qualche libreria.
Visto che è il calcolo del seno e del coseno "pesante"... hai provato a vedere se puoi costruirtela da te la libreria?
Se usi angoli entro una certa approssimazione puoi ad esempio tabularti... a costo di usare un buon 5-6MB di memoria però hai un accesso a costo costante senza calcoli, meglio di così non saprei proprio.
Visto che è il calcolo del seno e del coseno "pesante"... hai provato a vedere se puoi costruirtela da te la libreria?
Se usi angoli entro una certa approssimazione puoi ad esempio tabularti... a costo di usare un buon 5-6MB di memoria però hai un accesso a costo costante senza calcoli, meglio di così non saprei proprio.
Purtroppo non ci sono limiti per per gli angoli quindi tabulare i risultati non può essere effettuato, ciò che mi servirebbe è riuscire ad utilizzare una libreria che sfrutti le nuove potenzialità dei processori come la KML o (meglio come suggerisce Unrue ;) ) la ACML della AMD che è gratuita ma nell'utilizzare le librerie esterne e linkarle al programma non sono capace.
Se qualcuno può postarmi qualche guida e o suggerimento per utilizzare librerie gliene sarei grato :help:
Attualmente non sono pratico nella creazione di una libreria, menche meno nella creazione di una libreria altamente ottimizzata per il calcolo del seno e coseno, non conosco così affondo l'assembly da poter sfruttare a basso livello le istruzioni processore ed ottenere un risultato migliore di quanto non faccia la funzione standard. Ho provato ad utilizzare qualche altra funzione trovata in rete, dei metodi di approssimazione del seno e del coseno, fra quei pochi che sono riusciti a far funzionare risultavano tutti più lenti.
Ti avevo suggerito quello più interno :) Quando fai l'unroll le prestazioni possono peggiorare, tutto dipende dalla cache. Però in generale, siccome togli un sacco di if alla fine di ogni ciclo, si migliora. Inoltre un unroll di 64 è esagerato,( gestire i resti da 1 a 63 è pazzesco !) infatti ti avevo suggerito di partire con 2,4,8 ecc...
Purtroppo non conosco Visual Studio,ma comunque restrict è presente nello standard C99, prova a guardare se è già configurato per tale standard.
Si, infatti Tra ACML e MKL ti ho consigliato le prime proprio perché gratuite.
L'errore è chiaro, avrai qualche path sbagliato. Ma anche su questo non posso aiutarti non essendo pratico di Visual Studio.
A quanto dice MSDN _restrict è analogo della keyword restrict definita dallo standard C99 ma è utilizzabile sia a codice C che C++, sono riuscito ad utilizzare la keyword, praticamente non può essere utilizzata per double dydx[] ma per double *dydx, ho provato ad utilizzarla ma nessun miglioramento probabilmente è legata alla mancata modificata effettuata al compilatore,stando a MSDN dovrei cambiare le runtime library ma la runtime library che mi serve non c'è nell'elenco.
Per l'unroll ho provato a fare un unroll di passo 2 per un'altra funzione che anch'essa richiamata svariate milioni di volte ma non prevede la necessità di utilizzare seno e coseno e quindi risulta ben più veloce ma anche in questo caso non ho riscontrato miglioramenti. Nel corpo di questa funzione vi sono vari prodotti e somme. Ho provato a fare un unroll di 64 semplicemente per fare una cosa più estrema e verificare se c'erano dei miglioramenti, non pensando ad eventuali riscontri negativi per quanto riguarda la cache.
Devo cercare di utilizzare quelle librerie che mi hai consigliato, ma sono un po' in alto mare, non capisco perchè non ci sia una guida facile facile su come utilizzare delle librerie non fornite direttamente con l'IDE :muro:
A quanto dice MSDN _restrict è analogo della keyword restrict definita dallo standard C99 ma è utilizzabile sia a codice C che C++, sono riuscito ad utilizzare la keyword, praticamente non può essere utilizzata per double dydx[] ma per double *dydx, ho provato ad utilizzarla ma nessun miglioramento probabilmente è legata alla mancata modificata effettuata al compilatore,stando a MSDN dovrei cambiare le runtime library ma la runtime library che mi serve non c'è nell'elenco.
Per l'unroll ho provato a fare un unroll di passo 2 per un'altra funzione che anch'essa richiamata svariate milioni di volte ma non prevede la necessità di utilizzare seno e coseno e quindi risulta ben più veloce ma anche in questo caso non ho riscontrato miglioramenti. Nel corpo di questa funzione vi sono vari prodotti e somme. Ho provato a fare un unroll di 64 semplicemente per fare una cosa più estrema e verificare se c'erano dei miglioramenti, non pensando ad eventuali riscontri negativi per quanto riguarda la cache.
Devo cercare di utilizzare quelle librerie che mi hai consigliato, ma sono un po' in alto mare, non capisco perchè non ci sia una guida facile facile su come utilizzare delle librerie non fornite direttamente con l'IDE :muro:
Scusa un attimo,ma questo errore:
cl driver.c -Ic:/acml4.3.0/win64_mp/include
c:/acml4.3.0/win64_mp/lib/libacml_mp_dll.lib
driver.c
c1: fatal error C1083: Cannot open source file: 'driver.c': No such file or directory
Il cui comando hai preso da qua presumo:
http://biowulf.nih.gov/doc/acml/Linking_002fWindows.html
Il file driver.c magari non va sostituito con il nome del tuo sorgente?:fagiano: Quello è solo un esempio presumo..
Riguardo il restrict, se non hai miglioramenti vuol dire che il compilatore non fa ottimizzazioni su quel loop, oppure che le fa ed è riuscito a capire che tale puntatore ha un utilizzo ristretto, quindi aggiungendola tu non cambia nulla. Io comunque la lascerei, non si sa mai..
Purtroppo non ci sono limiti per per gli angoli quindi tabulare i risultati non può essere effettuato.
Beh, non vedo il motivo: dato che seno e coseno sono funzioni periodiche la tabella deve coprire solo l'intervallo 0...pigreco.
Quindi basta riportare l'angolo in questo intervallo (una divisione) e regolare il segno di conseguenza.
Tutto sta a decidere la risoluzione dei valori della tabella, ma come dice _Claudio in 5-6 mb ci stanno milioni di valori... difficile notare la differenza in precisione.
Mentre quella in velocità sarebbe decisamente evidente :read:
Scusa un attimo,ma questo errore:
cl driver.c -Ic:/acml4.3.0/win64_mp/include
c:/acml4.3.0/win64_mp/lib/libacml_mp_dll.lib
driver.c
c1: fatal error C1083: Cannot open source file: 'driver.c': No such file or directory
Il cui comando hai preso da qua presumo:
http://biowulf.nih.gov/doc/acml/Linking_002fWindows.html
Il file driver.c magari non va sostituito con il nome del tuo sorgente?:fagiano: Quello è solo un esempio presumo..
Riguardo il restrict, se non hai miglioramenti vuol dire che il compilatore non fa ottimizzazioni su quel loop, oppure che le fa ed è riuscito a capire che tale puntatore ha un utilizzo ristretto, quindi aggiungendola tu non cambia nulla. Io comunque la lascerei, non si sa mai..
Effettivamente..... :mc: :mc: si era solo un esempio :doh: comunque ho provato a quindi a mettere tutti i file del mio programma e mettere i percorsi corretti, non da alcun errore del tipo non riesco a trovare o aprire quel file o quella libreria, ma da lo stesso identico errore che mi da compilando direttamente dall'IDE ossia unresolved external linker, eppure ho provato sia con le librerie nella cartella win64 che in quela win64_mp, ho provato persino a mettere tutte le librerie ma mi da errore dicendo che alcune librerie vanno in conflitto chissà come mai :D
Comunque ho creato una stupida libreria, una libreria che come unica funzione ha questa
int SempliceSomma(int a,int b);
una libreria statica (.lib)
poi ho creato un semplice programma
#include <stdio.h>
#include <windows.h>
#include <IntestazioneLibreria.h>
int main()
{
printf("Prova Compilazione\n");
printf("%d \n",SempliceSomma(3,6));
system("PAUSE");
return 0;
}
e il tutto era compilato correttamente, ho ovviamente fatto le modifiche necessarie per far vedere file di intestaizone e libreira.
Poi ho fatto gli stessi identici passi nel programma vero e proprio ma l'unresolved link rimane :muro: :muro:
riporto l'errore esatto comunque:
error LNK2019: unresolved external symbol _fastsincos referenced in function _derivsU
Ma tu usi la stessa identica versione delle librerie? non vorrei che fosse una versione buggata.
Beh, non vedo il motivo: dato che seno e coseno sono funzioni periodiche la tabella deve coprire solo l'intervallo 0...pigreco.
Quindi basta riportare l'angolo in questo intervallo (una divisione) e regolare il segno di conseguenza.
Tutto sta a decidere la risoluzione dei valori della tabella, ma come dice _Claudio in 5-6 mb ci stanno milioni di valori... difficile notare la differenza in precisione.
Mentre quella in velocità sarebbe decisamente evidente :read:
Creare i valori tabellati non sarebbe un grosso problema
int i;
double *ValoriTabellati;
ValoriTabellati=(double *)malloc(sizeof(double)*40000001);
for (i=0;i<=40000000;i++)
{
ValoriTabellati[i]=sin(2*3.14*i/40000000);
}
Anche il tempo di calcolo non è eccessivo (istantaneo) il tempo di premere start... (inizio ad avere dei dubbi sull'analizzatore di performance di Visual Studio....:muro: )
con fmod dovrei ottenere il resto della divisione fra esp e 2*pi, poi dovrei ottenere l'indice corrispondente all'angolo ed arrotondarlo con floor, e poi accedere a quell'indice.
Implemento il tutto e faccio sapere i risultati :)
Ho implementato il codice ed ho fatto qualche prova, innanzitutto riporto il codice:
int i;
ValoriTabellatiCos=(double *)malloc(sizeof(double)*(NTABELLATI+1));
ValoriTabellatiSin=(double *)malloc(sizeof(double)*(NTABELLATI+1));
for (i=0;i<=NTABELLATI;i++)
{
ValoriTabellatiCos[i]=cos(DUEPERPISUNTABELLATI*i);
ValoriTabellatiSin[i]=sin(DUEPERPISUNTABELLATI*i);
}
Nella funzione derivsU invece ho sostituito al posto di
q1=K*sin(esp);
q2=K*cos(esp);
le seguenti righe di codice
if (esp==0)
{
q1=0;
q2=K;
}
else
{
ValoreAssoluto=fabs(esp);
IndiceRicerca=(int) floor((fmod(ValoreAssoluto,DUEPERPI)*NTABELLATISUDUEPERPI)+0.5);
q1=K*esp/ValoreAssoluto*ValoriTabellatiSin[IndiceRicerca];
q2=K*ValoriTabellatiCos[IndiceRicerca];
}
Ho fatto qualche prova per valutare tempi di esecuzione e correttezza dell'algoritmo, ho preso in esame un caso in cui la funzione era richiamata circa 3.700.000 volte, ed N=26 per un totale di circa 190.000.000 chiamate alle due funzioni trigonometriche incriminate cos(esp) e sin(esp).
Senza Tabella 1.04
Con Tabella da 700000 valori 1.02
Con Tabella da 400000 valori 0.52
Con Tabella da 40000 valori 0.37
Con Tabella da 4000 valori 0.39
Quindi al di sotto di 40000 valori il tempo è praticamente identivo (il tempo lo rilevo a mano e poi avevo decine di pagine internet aperte e visual studio aperto quindi qualche secondo di differenza ci può stare comunque anche per un 4000 valori l'errore massimo relativo è dell'ordine dello 0.01% quindi del tutto trascurabile (farò alter prove per verificare tempi ed accuratezza dei risultati).
Nel caso questa dovesse risultare un'ottimo approccio penso di utilizzare una tabella da 40000 valori che mi da un errore relativo di due ordini di grandezza inferiore.
Se qualcuno ha suggerimenti per migliorare quella porzione di codice, in particolare per eliminare l'if che dovrebbe far perdere un bel po' di tempo o per ovviare a qualche chiamata a funzioni fmod e floor non esiti a postare :D
Qualsiasi commento e soluzione è ben accetta.
(Ovviamente i ringraziamenti a tutti coloro che stanno cercando e che cercheranno di aiutarmi sono scontati ;) )
Comunque inizio ad avere dei forti dubbi sull'analizzatore di performace di Visual Studio 2010 beta 1, perchè il tempo per creare le due tabelle da 40.000.000 valori quindi 80.000.000 chiamate alle funzioni seno e coseno (e altrettante divisioni) risulta essere pressochè nullo 1-2 secondi :mbe:
Vedi che era più veloce :D
Cmq l'if lo puoi proprio togliere, dato che 0.0/0.0 coi float è definito come #inf e cos0 e' 1 (valore che puoi includere nella tabella).
Potresti anche togliere floor() che dato il minimo intervallo di angolo tra il valore approssimato per difetto e per eccesso non porta ad alcuna imprecisione.
Per il resto se vuoi misurare il tempo è meglio che usi GetQueryPerformanceCounter & co su windows che dà la massima precisione...
Inoltre dovresti considerare di ridurre le chiamate a funzione, che hanno si un peso irrisorio ma quando n fai 3.700.000 iniziano a farsi sentire come overhead :asd:
Prova a mettere quella che è derivsU all'interno di un ciclo while di una funzione più ampia, per esempio.
Vedi che era più veloce :D
Cmq l'if lo puoi proprio togliere, dato che 0.0/0.0 coi float è definito come #inf e cos0 e' 1 (valore che puoi includere nella tabella).
Potresti anche togliere floor() che dato il minimo intervallo di angolo tra il valore approssimato per difetto e per eccesso non porta ad alcuna imprecisione.
Per il resto se vuoi misurare il tempo è meglio che usi GetQueryPerformanceCounter & co su windows che dà la massima precisione...
Inoltre dovresti considerare di ridurre le chiamate a funzione, che hanno si un peso irrisorio ma quando n fai 3.700.000 iniziano a farsi sentire come overhead :asd:
Prova a mettere quella che è derivsU all'interno di un ciclo while di una funzione più ampia, per esempio.
il problema è che non è l'unica funzione che può essere chiamata, era solo un esempio, infatti è nominata come derivsU, ma ce ne sono anche altre derivsCR, derivsCH, derivsL, derivsUL ecc... per le quali ho comunque verificato che il collo di bottiglia era sempre cos(esp) e sin(esp), per evitare di fare degli if o degli switch utilizzo dei puntatori a funzione, ed è per questo che ho estrapolato la parte di codice in delle funzioni a sè stanti.
per quanto riguarda 0/0 si evvero da #inf, però poi quel valore lo utilizzo nelle successive due istruzioni di codice
dydx[i]-=q1*y[j]+q2*y[j+n];
dydx[i+n]+=q2*y[j]-q1*y[j+n];
e quindi dydx[i] anzichè essere pari a q2*y[j+n] diventa pari a #inf, compromettendo tutti i risultati perchè poi l'#inf è propagato a y[] e quindi per tutta l'esecuzione del programma
Utilizzerò GetQueryPerformanceCounter come mi hai suggerito per meglio misurare i tempi di esecuzione ;)
Uhm già allora l'if serve... comunque quanto pesa dipende dall'incidenza di casi true sul totale: cioè se esp è 0 pochissime volte il valutare l'if potrebbe diventare più pesante di fare il calcolo e basta.
Comunque non saprei come renderlo più leggero, probabilmente l'if non pesa poi tanto.
Cmq i puntatori a funzione pesano molto di più delle chiamate a funzione normali, secondo me potresti farci qualcosa.
Certo, perderesti tantissimo in mantenibilità e gestibilità...
Ma tu usi la stessa identica versione delle librerie? non vorrei che fosse una versione buggata.
Io utilizzo la versione Linux delle ACML. In pratica a te manca l'equivalente in Linux del flag -L e -l .. Non basta dargli l'header
Io utilizzo la versione Linux delle ACML. In pratica a te manca l'equivalente in Linux del flag -L e -l .. Non basta dargli l'header
Io gli do l'header e dico al linker dove andare a pescare le librerie, ho provato a fare la stessa cosa con le librerie SDL, (delle librerie multipiattaforma per creare applicazioni grafiche in C) e funziona. Tu che compilatore usi?? gcc immagino?
L'idea di utilizzare delle tabelle con i valori già precalcolati è interessante, però riuscire ad utilizzare delle librerie specificatamente pensate per sfruttare le nuove potenzialità dei processori, poter utilizzare in maniera più veloce i vettori mi piaceva ancor di più. :eek:
Io gli do l'header e dico al linker dove andare a pescare le librerie, ho provato a fare la stessa cosa con le librerie SDL, (delle librerie multipiattaforma per creare applicazioni grafiche in C) e funziona. Tu che compilatore usi?? gcc immagino?
Sto usando vari compilatori tra i quali Intel, PGI e GCC. Ma devi per forza usare Visual Studio? Prova con Gcc per Windows a vedere se ti funziona.
A distanza di tempo ed usando Visual studio 2010 beta 2 ho riprovato ad utilizzare le librerie ACML e (probabilmente sbagliavo qualcosa oppure Visual Studio 2010 beta 1 aveva qualche problema) la procedura per chi volesse utilizzare tali librerie con Visual studio è la seguente:
1) installare ACML
2) nel programma inserire gli header necessari, ad esempio
#include <acml.h>
3) poi inserire nel progetto i percorsi per gli header
Project->NomeProgetto Properties->Configuration Properties->VC++ Directories->Include Directories ed inserire il percorso opportuno
4) indicare i file .lib necessari, ad esempio
Project->NomeProgetto Properties->Configuration Properties->VC++ Directories->Linker->Input->Additional Dependencies ed inserire C:\AMD\acml4.3.0\win64_mp\lib\libacml_mp_dll.lib o i file .lib opportuni
5) copiare tutti i file .dll nella cartella di lavoro del progetto
(Testato con Windows 7)
Ciao,
visto che lavori su processore Intel, forse sono più indicate le MKL ( Math Kernel Library), sviluppate da Intel stessa. A suo tempo ti indicai le ACML in quanto le MKL ancora non le conoscevo :)
Inoltre ti consiglio anche di utilizzare il compilatore Intel ( icc), versione free. In genere è più performante di gcc, sopratutto su architetture Intel.
HyperText
06-12-2009, 10:18
Dico una stupidaggine: hai provato ad utilizzare funzioni inline?
@HyperText
No non ho provato ad usare funzioni inline, devo approfondire dato che se definisco la funzione con la parola chiave __inline il compilatore mi da un unresolved link.
Utilizzando ACML sono riuscito ad ottenere oltre il 50% di riduzione dei tempi di calcolo (rispetto a prima dell'utilizzo di ACML, quindi complessivamente un codice 30 volte più veloce :sofico:), ho trasformato i calcoli in calcoli fra vettori e matrici. Il punto debole del programma rimane sempre il calcolo del seno e del coseno (che effettuo tramite la funzione fastsincos()) che stando a Visual Studio mi occupa il 70% del tempo.
@Unrue
Sai per caso di qualche buona guida per le funzioni BLAS?
Finora sono riuscito a trovare tutte le funzioni che mi servivano tranne una, il prodotto elemento per elemento (element-wise) fra due vettori (o matrici) tu conosci quale funzione devo usare? (ho letto che si parla anche di Hadamard product, ho trovato anche una funzione ma non c'è in ACML)
Attualmente faccio il prodotto elemento per elemento tramite il classico ciclo for e ciò mi comporta oltre il 15% di tempo, se potessi farlo tramite BLAS otterrei una ulteriore riduzione di tempo :D (tutte le funzioni BLAS che ho utilizzato comportano meno dello 0.1% del tempo :eek: ).
Comunque ho letto che le GPU consentono di effettuare il calcolo del seno e del coseno con 1 ciclo :eek: quindi se riuscissi ad utilizzare la GPU avrei un incremento prestazionale enorme. Purtroppo ACML-GPU non è thread-safe (quindi non potrei usare più thread) ed inoltre consente solo di effettuare prodotti fra matrici e non di usare sin e cos quindi non molto utile nel mio caso.
Per quanto riguarda MKL e il compilatore Intel non sono gratuiti per Windows mi sembra, anche se in MKL si ha la possibilità di fare operazioni direttamente con i numeri complessi, che mi eviterebbe di passare per le funzioni sin e cos, quindi potrebbe essere molto interessante.
@HyperText
Utilizzando ACML sono riuscito ad ottenere oltre il 50% di riduzione dei tempi di calcolo (rispetto a prima dell'utilizzo di ACML, quindi complessivamente un codice 30 volte più veloce :sofico:), ho trasformato i calcoli in calcoli fra vettori e matrici. Il punto debole del programma rimane sempre il calcolo del seno e del coseno (che effettuo tramite la funzione fastsincos()) che stando a Visual Studio mi occupa il 70% del tempo.
E' un ottimo risultato! Considera che quei calcoli li devi fare, quindi non puoi azzerare quel tempo :) E con MKL farai sicuramente meglio.
@Unrue
Sai per caso di qualche buona guida per le funzioni BLAS?
La migliore guida è la guida ufficiale:
http://www.netlib.org/blas/
Finora sono riuscito a trovare tutte le funzioni che mi servivano tranne una, il prodotto elemento per elemento (element-wise) fra due vettori (o matrici) tu conosci quale funzione devo usare? (ho letto che si parla anche di Hadamard product, ho trovato anche una funzione ma non c'è in ACML)
Attualmente faccio il prodotto elemento per elemento tramite il classico ciclo for e ciò mi comporta oltre il 15% di tempo, se potessi farlo tramite BLAS otterrei una ulteriore riduzione di tempo :D (tutte le funzioni BLAS che ho utilizzato comportano meno dello 0.1% del tempo :eek: ).
Non sono esperto di BLAS, ma sicuramente c'è la funzione di cui hai bisogno, Quelle che iniziano con DGEM mi pare.
Comunque ho letto che le GPU consentono di effettuare il calcolo del seno e del coseno con 1 ciclo :eek: quindi se riuscissi ad utilizzare la GPU avrei un incremento prestazionale enorme. Purtroppo ACML-GPU non è thread-safe (quindi non potrei usare più thread) ed inoltre consente solo di effettuare prodotti fra matrici e non di usare sin e cos quindi non molto utile nel mio caso.
Non sempre le GPU ti danno benefici, dipende dai casi. E magari in questo caso potresti non "threadizzare" la parte di codice che passi alla GPU. Tanto poi ci pensa la GPU stessa a parallelizzarla.
Per quanto riguarda MKL e il compilatore Intel non sono gratuiti per Windows mi sembra, anche se in MKL si ha la possibilità di fare operazioni direttamente con i numeri complessi, che mi eviterebbe di passare per le funzioni sin e cos, quindi potrebbe essere molto interessante.
Devi usare per forza Visual Studio e Windows? Passare al compilatore Intel ( e a Linux per la versione free) secondo me ti darebbe un grosso boost. C'è la versione gratuita che comprende anche le MKL.:
http://software.intel.com/en-us/articles/non-commercial-software-download/
PS: Hai lasciato l'unroll dei loop come ti avevo suggerito a suo tempo?
Considera anche che non puoi ottimizzare all'infinito, arrivi ad un certo punto che più di lì non vai. 30 volte più veloce mi sembra già più che ottimo.
The3DProgrammer
07-12-2009, 00:18
in genere buona pratica è eliminare tutti i rapporti, calcolare prima il reciproco e poi fare il prodotto col divisore.
double rec = 1/x;
double z=y*x;
invece di
double z = y/x;
ora non ricordo bene l'assembly che viene generato nei due casi, ma se non mi sbaglio il primo caso è abbastanza + veloce (mi pare di ricordare che per il calcolo del reciproco c'è un'istruzione ad hoc o qualcosa del genere)
bye
Purtroppo non conosco Visual Studio,ma comunque restrict è presente nello standard C99, prova a guardare se è già configurato per tale standard.
Il compilatore MS non è compatibile con C99.
Ci puoi fare vedere il codice ora ? Potremmo magari vedere se si può fare qualche trasformazione sul calcolo matriciale.
Un ebook interessante: http://www.agner.org/optimize/optimizing_cpp.pdf
@Tommo
Grazie del link, leggerò sicuramente il libro ;)
@Unrue
Sicuramente ho già ottenuto un più che ottimo risultato, però più ottimizzo più posso ottenere risultati accurati, insomma finchè si può... ottimizzo :D
Purtroppo il codice è parallelizzato invarie parti (elaborazione e salvataggio dei dati) facendo uso di funzioni specefiche per Windows, dovrei quindi riscriverlo, inoltre il codice dovrà girare anche su windows. Mi costerebbe un po' di fatica, comunque vedo che si può fare, magari riesco a sfruttare le MKL.
L'unroll del codice l'ho tolto in quanto mi da tempi analoghi (se non leggermente peggiori, forse il compilatore fa un unroll più efficiente, o il processore ha degli algoritmi predettivi efficaci per i salti nei loop).
@cionci
Ecco il codice:
void derivsU(double x , double y[], double dydx[],double lam)
{
int i;
int n=N;
int nquadro=n*n;
double Kmod=(1+fastcos(2*pi/lam*x));
double **xdbeta;
double **mcos, **msin;
xdbeta=matrixd(1,n,1,n);
mcos=matrixd(1,n,1,n);
msin=matrixd(1,n,1,n);
dcopy(nquadro,&beta[1][1],1,&xdbeta[1][1],1);
dscal(nquadro,x,&xdbeta[1][1],1);
vrda_sincos(nquadro,&xdbeta[1][1],&msin[1][1],&mcos[1][1]);
for (i=1;i<=nquadro;i++){
msin[1][i]*=intePerParam[1][i];
mcos[1][i]*=intePerParam[1][i];
}
dgemv('t',n,n,-Kmod,&mcos[1][1],n,&y[n+1],1,0,&dydx[1],1);
dgemv('t',n,n,-Kmod,&msin[1][1],n,&y[1],1,1,&dydx[1],1);
dgemv('t',n,n,Kmod,&mcos[1][1],n,&y[1],1,0,&dydx[1+n],1);
dgemv('t',n,n,-Kmod,&msin[1][1],n,&y[n+1],1,1,&dydx[n+1],1);
free_matrixd(xdbeta,1,n,1,n);
free_matrixd(mcos,1,n,1,n);
free_matrixd(msin,1,n,1,n);
}
Le matrici sono salvate in memoria in modo che tutte le righe siano memorizzate in spazi consecutivi di memoria, in questo modo posso trattarle indifferentemente come matrici o come vettori.
Essenzialmente il codice prevede la copia della matrice beta nella matrice xdbeta (dcopy) che poi moltiplico per x (dscal). poi calcolo il seno e il coseno di questa matrice (vrda_sincos), moltiplico elemento per elemento le due matrici ottenute per la matrice intePerParam (tramite il ciclo for).
Poi effettuo dei prodotti matriciali (o meglio matrice-vettore) e salvo il risultato in dydx.
Sicuramente può essere ottimizzato, magari evitando di creare e liberare le matrici ad ogni chiamata alla funzione, ma il punto lento è il calcolo del seno e del coseno. Ho utilizzato AMD CodeAnalyst e i risultati sono:
vrda_sincos 65% del tempo (funzione per il calcolo del seno e del coseno)
zupmtr_ circa 20% del tempo (dovrebbe essere una funzione che consente di salvare le matrici, io non faccio esplicita richiesta a questa funzione, sarà richiamata da un'altra funzione BLAS ma non so quale sia)
il mio programma occupa meno del 10% del tempo di cui il 13% nella creazione delle matrici, e il 70% nel fare il prodotto elemento per elemento. Qundi in parole pavore c'è da concentrarsi sul calcolo dei seni e dei coseni e su quella funzione il resto darebbe risultati pressochè trascurabili.
Sicuramente ho già ottenuto risultati ottimi, vedo solo se è possibile migliorare ulteriormente :D , devo informarmi meglio sulla possibilità di utilizzare le schede video che dovrebbero essere delle schegge nelle funzioni trigonometriche.
Per chi avesse bisogno di utilizzare le funzioni BLAS ottime guide sono anche le seguenti:
http://www.intel.com/software/products/mkl/docs/mklqref/
e una volta individuata la funzione altre informazioni (molto più dettagliate) sulla funzione e sui parametri possono essere trovati in
http://techpubs.sgi.com/library/tpl/cgi-bin/browse.cgi?db=man&coll=0650&pth=/cat3
Non conosco le grandezze in gioco. Ti ci mette tanto tempo perché devi eseguire tante volte quella funzione ?
Sicuramente ho già ottenuto risultati ottimi, vedo solo se è possibile migliorare ulteriormente :D , devo informarmi meglio sulla possibilità di utilizzare le schede video che dovrebbero essere delle schegge nelle funzioni trigonometriche.
Secondo me ti conviene calcolare i flop/s che ottieni con il tuo processore. In modo da valutare se sei arrivato al limite con il processore che hai.
Per fare questo, devi calcolare il numero di operazioni floating point che effettui e dividere per il tempo che impieghi. Poi lo confronti con flop/s teorici del tuo processore. Se raggiungi il 30/40% di questo valore, vuol dire che hai già raggiunto il massimo possibile con questo processore. Questo perché comunque intervengono altri fattori limitanti quali ad esempio la memoria.
Adesso però, avendo chiamate BLAS, è più complicato, dovresti vedere l'assembly generato e calcolarle da lì. Non è banalissimo, ma secondo me conviene farlo. Altrimenti rischi di continuare a cercare di ottimizzare quando poi scopri che più di lì è impossibile andare. Basta farlo facendo un ciclo di for. Altrimenti è un bordello.
Sono riuscito ad usare un trucco per eliminare il calcolo continuo dei seni e dei coseni ottenendo quindi un notevole boost prestazionale (precalcolo i seni e i coseni prima di dare in pasto al solutore che richiamerà milioni di volte la funzione, i seni e i coseni sono calcolati in modo esatto).
Il codice risulta quindi (ho inserito sia il codice con le librerie BLAS che quello senza)
void derivsU(double x , double y[], double dydx[],double lam)
{
#if !defined NOT_USE_ACML
int n=N;
int Indice=(int) floor(x*2/h+0.5)+1;
double Kmod=(1+fastcos(2*pi/lam*x));
dgemv('t',n,n,-Kmod,&ValoriPreCalcolatiCos[Indice][1][1],n,&y[n+1],1,0,&dydx[1],1);
dgemv('t',n,n,-Kmod,&ValoriPreCalcolatiSin[Indice][1][1],n,&y[1],1,1,&dydx[1],1);
dgemv('t',n,n,Kmod,&ValoriPreCalcolatiCos[Indice][1][1],n,&y[1],1,0,&dydx[1+n],1);
dgemv('t',n,n,-Kmod,&ValoriPreCalcolatiSin[Indice][1][1],n,&y[n+1],1,1,&dydx[n+1],1);
#else
int i,j;
int n=N;
int Indice=(int) floor(x*2/h+0.5)+1;
double Kmod=(1+cos(2*pi/lam*x));
double q1,q2;
for(i=1;i<=n;i++){
dydx[i]=0;
dydx[i+n]=0;
for(j=1;j<=n;j++)
{
q1=Kmod*ValoriPreCalcolatiSin[Indice][i][j];
q2=Kmod*ValoriPreCalcolatiCos[Indice][i][j];
dydx[i]-=q1*y[j]+q2*y[j+n];
dydx[i+n]+=q2*y[j]-q1*y[j+n];
}
}
#endif
}
Ora per ottimizzare ulteriormente non mi resta che evitare anche di fare somme e moltiplicazioni ("dydx[i]-=q1*y[j]+q2*y[j+n];") e poi il programma sarà diventerà ancora più veloce :p .
Scherzo ovviamente, meglio di così non credo si possa fare per quella funzione.
Scrivo questo messaggio semplicemente per suggerire a chiunque avesse problemi di questo tipo di pensare in generale a tutto il programma nella sua interezza senza focalizzarsi eccessivamente sulla singola funzione. Nel mio caso ad esempio i sin e cos dipendevano sì dalla x e quindi variavano ad ogni chiamata ma gli stessi identici valori di x erano usati migliaia e migliaia di volte.
Per la cronaca il tutto è inserito in un ciclo che serve per risolvere un sistema di equazioni differenziali del prim'ordine, il solutore è un solutore Runge Kutta del 4° ordine RK4 a step costante, utilizzare un Runge kutta del 4° 5° RK45 a step variabile non risolve molto dato che il problema mi sembra si presenti moderatamente stiff, il passo di soluzione quindi resta piccolo (circa il doppio) e in più in questo caso dovrei calcolare sempre seni e coseni.
Se qualcuno avesse da suggerirmi un solutore più efficiente, mi sembra di aver capito che i solutori impliciti dovrebbero adattarsi meglio a quelli espliciti (come il RK) in queste situazioni, ben venga :D
Inoltre in determinati contesti solo una parte dell'elaborazione risulta essere stiff mentre la restante no, in questi casi quindi è opportuno cambiare solutore, nel mio caso passo dal RK4 al RK45 che consente di aumentare enormemente le prestazioni per la parte non stiff.
In definitiva quindi se non si può ottimizzare ulteriormente la funzione più onerosa occorre ingegnarsi per capire come è possibile ridurre il numero di chiamate a tale funzione. :sofico:
Ad ogni modo mi posso ritenere più che soddisfatto.
Grazie a tutti.
Ottimo!
Quanto è il boost rispetto al codice originale?
Ottimo!
Quanto è il boost rispetto al codice originale?
Con un test sintetico in un tipico scenario di utilizzo
senza BLAS e senza precalcolo 8'33''
senza BLAS e con precalcolo 45'' :cool: :eek:
con BLAS e senza precalcolo 3'05''
con BLAS e con precalcolo 53''
il test l'ho fatto con tanti altri programmi avviati, che non facevano niente ma possono aver modificato leggermente i valori. Le librerie BLAS risultano quindi estremamente valide per velocizzare il calcolo dei sin e cos, per i prodotti matriciali penso che il compilatore di Visual Studio riesca ad ottimizzare maggiormente, anche se di poco.
Dando uno sguardo al codice assembly, da quel poco che ne capisco, usa funzioni tipo xmms o qualcosa del genere quindi sfrutta le istruzioni SSE che velocizzano, anche le librerie BLAS (per l'esattezza l'implementazione dell'AMD ossia ACML) sicuramente le sfruttano ma il dover richiamare ogni volta delle DLL esterne probabilmente causa quel leggero degrado prestazionale. Almeno credo.
Con un test sintetico in un tipico scenario di utilizzo
senza BLAS e senza precalcolo 8'33''
senza BLAS e con precalcolo 45'' :cool: :eek:
con BLAS e senza precalcolo 3'05''
con BLAS e con precalcolo 53''
il test l'ho fatto con tanti altri programmi avviati, che non facevano niente ma possono aver modificato leggermente i valori. Le librerie BLAS risultano quindi estremamente valide per velocizzare il calcolo dei sin e cos, per i prodotti matriciali penso che il compilatore di Visual Studio riesca ad ottimizzare maggiormente, anche se di poco.
Dando uno sguardo al codice assembly, da quel poco che ne capisco, usa funzioni tipo xmms o qualcosa del genere quindi sfrutta le istruzioni SSE che velocizzano, anche le librerie BLAS (per l'esattezza l'implementazione dell'AMD ossia ACML) sicuramente le sfruttano ma il dover richiamare ogni volta delle DLL esterne probabilmente causa quel leggero degrado prestazionale. Almeno credo.
In realtà le BLAS sono molto famose per i prodotti matriciali, strano che vadano leggermente peggio. Forse aumentando la dimensionalità del test case andranno sicuramente meglio. Non ricordo se giò te lo consigliai, ma se sei su processori Intel, forse vale la pena provare anche limplementazione BLAS delle MKL, che sono ancora più spinte. ( le trovi insieme al compilatore Intel e la versione non pro è gratuita)
Hai controllato anche l'accuratezza dei risultati? Hai differenze tra gli output originali e quelli ottimizzati? Se si, di quanto?
PS:in realtà gli xmms sono registri vettoriali, non funzioni :), utilizzati dalle istruzioni SSE.
i sin e cos precalcolati sono calcolati in maniera esatta, cioè con le normali funzioni sin e cos (o con vrda_sincos se uso ACML), quindi i risultati sono pressoché identici. La massima variazione che ho riscontrato nel test è di 1.2488e-008%, sarà dovuto a qualche errore di troncamento sulle ultime cifre.
Le MKL mi sembrano che siano gratis solo per Linux :( , e il mio programma fa uso di funzioni di parallelizzazioni scritte per windows.
Comunque il codice Assembly che Visual Studio genera per quella funzione è il seguente:
249: void derivsU(double x , double y[], double dydx[],double lam)
250: {
000000013F3F1300 40 53 push rbx
000000013F3F1302 57 push rdi
000000013F3F1303 41 55 push r13
000000013F3F1305 41 57 push r15
000000013F3F1307 48 83 EC 68 sub rsp,68h
000000013F3F130B 0F 29 7C 24 30 movaps xmmword ptr [rsp+30h],xmm7
000000013F3F1310 0F 29 74 24 40 movaps xmmword ptr [rsp+40h],xmm6
000000013F3F1315 66 0F 28 F0 movapd xmm6,xmm0
251: #if !defined NOT_USE_ACML
252: /* La parte iniziale del vettore contiene i dati legati alla parte
253: * reale dell'ampiezza, mentre la seconda parte contiene la parte
254: * immaginaria.
255: */
256: int n=N;
257: int Indice=(int) floor(x*2/h+0.5)+1;
258: double Kmod=(1+fastcos(2*pi/lam*x));
259:
260: dgemv('t',n,n,-Kmod,&ValoriPreCalcolatiCos[Indice][1][1],n,&y[n+1],1,0,&dydx[1],1);
261: dgemv('t',n,n,-Kmod,&ValoriPreCalcolatiSin[Indice][1][1],n,&y[1],1,1,&dydx[1],1);
262:
263: dgemv('t',n,n,Kmod,&ValoriPreCalcolatiCos[Indice][1][1],n,&y[1],1,0,&dydx[1+n],1);
264: dgemv('t',n,n,-Kmod,&ValoriPreCalcolatiSin[Indice][1][1],n,&y[n+1],1,1,&dydx[n+1],1);
265: #else
266: int i,j;
267: int n=N;
268: int Indice=(int) floor(x*2/h+0.5)+1;
269: double Kmod=(1+cos(2*pi/lam*x));
000000013F3F1319 F2 0F 10 05 CF 79 03 00 movsd xmm0,mmword ptr [__real@401921fb54442d18 (13F428CF0h)]
000000013F3F1321 4D 8B F8 mov r15,r8
000000013F3F1324 48 8B FA mov rdi,rdx
000000013F3F1327 F2 0F 5E C3 divsd xmm0,xmm3
000000013F3F132B F2 0F 59 C6 mulsd xmm0,xmm6
000000013F3F132F E8 1C 3D 04 00 call cos (13F435050h)
000000013F3F1334 F2 0F 59 35 44 79 03 00 mulsd xmm6,mmword ptr [__real@4000000000000000 (13F428C80h)]
000000013F3F133C 48 63 1D 3D F1 03 00 movsxd rbx,dword ptr [N (13F430480h)]
000000013F3F1343 F2 0F 5E 35 7D 05 04 00 divsd xmm6,mmword ptr [h (13F4318C8h)]
000000013F3F134B 66 0F 28 F8 movapd xmm7,xmm0
000000013F3F134F F2 0F 58 3D B9 D0 03 00 addsd xmm7,mmword ptr [VersioneFileSalvataggio+20h (13F42E410h)]
000000013F3F1357 F2 0F 58 35 81 79 03 00 addsd xmm6,mmword ptr [__real@3fe0000000000000 (13F428CE0h)]
000000013F3F135F 66 0F 28 C6 movapd xmm0,xmm6
000000013F3F1363 E8 F8 86 02 00 call floor (13F419A60h)
000000013F3F1368 0F 28 74 24 40 movaps xmm6,xmmword ptr [rsp+40h]
000000013F3F136D F2 0F 2C C0 cvttsd2si eax,xmm0
000000013F3F1371 FF C0 inc eax
000000013F3F1373 4C 63 E8 movsxd r13,eax
270: double q1,q2;
271:
272: for(i=1;i<=n;i++){
000000013F3F1376 48 83 FB 01 cmp rbx,1
000000013F3F137A 0F 8C 34 03 00 00 jl derivsU+3B4h (13F3F16B4h)
251: #if !defined NOT_USE_ACML
252: /* La parte iniziale del vettore contiene i dati legati alla parte
253: * reale dell'ampiezza, mentre la seconda parte contiene la parte
254: * immaginaria.
255: */
256: int n=N;
257: int Indice=(int) floor(x*2/h+0.5)+1;
258: double Kmod=(1+fastcos(2*pi/lam*x));
259:
260: dgemv('t',n,n,-Kmod,&ValoriPreCalcolatiCos[Indice][1][1],n,&y[n+1],1,0,&dydx[1],1);
261: dgemv('t',n,n,-Kmod,&ValoriPreCalcolatiSin[Indice][1][1],n,&y[1],1,1,&dydx[1],1);
262:
263: dgemv('t',n,n,Kmod,&ValoriPreCalcolatiCos[Indice][1][1],n,&y[1],1,0,&dydx[1+n],1);
264: dgemv('t',n,n,-Kmod,&ValoriPreCalcolatiSin[Indice][1][1],n,&y[n+1],1,1,&dydx[n+1],1);
265: #else
266: int i,j;
267: int n=N;
268: int Indice=(int) floor(x*2/h+0.5)+1;
269: double Kmod=(1+cos(2*pi/lam*x));
000000013F3F1380 48 89 AC 24 98 00 00 00 mov qword ptr [y],rbp
000000013F3F1388 48 89 74 24 60 mov qword ptr [rsp+60h],rsi
000000013F3F138D 4D 8D 4F 08 lea r9,[r15+8]
000000013F3F1391 4C 89 64 24 58 mov qword ptr [rsp+58h],r12
273: dydx[i]=0;
000000013F3F1396 49 F7 DF neg r15
000000013F3F1399 4C 8D 14 DD 00 00 00 00 lea r10,[rbx*8]
000000013F3F13A1 4C 89 74 24 50 mov qword ptr [rsp+50h],r14
000000013F3F13A6 48 8B EB mov rbp,rbx
000000013F3F13A9 48 89 9C 24 90 00 00 00 mov qword ptr [x],rbx
000000013F3F13B1 33 C0 xor eax,eax
000000013F3F13B3 66 66 66 66 66 0F 1F 84 00 00 00 00 00 nop word ptr [rax+rax]
000000013F3F13C0 49 89 01 mov qword ptr [r9],rax
274: dydx[i+n]=0;
275: for(j=1;j<=n;j++)
000000013F3F13C3 BE 01 00 00 00 mov esi,1
000000013F3F13C8 4B 89 04 0A mov qword ptr [r10+r9],rax
000000013F3F13CC 48 83 FB 04 cmp rbx,4
000000013F3F13D0 0F 8C 00 02 00 00 jl derivsU+2D6h (13F3F15D6h)
270: double q1,q2;
271:
272: for(i=1;i<=n;i++){
000000013F3F13D6 48 8B 05 33 06 04 00 mov rax,qword ptr [ValoriPreCalcolatiSin (13F431A10h)]
000000013F3F13DD 4B 8D 0C 0F lea rcx,[r15+r9]
000000013F3F13E1 4C 8B F3 mov r14,rbx
000000013F3F13E4 4A 8B 04 E8 mov rax,qword ptr [rax+r13*8]
000000013F3F13E8 4D 8B E2 mov r12,r10
000000013F3F13EB 49 F7 DE neg r14
000000013F3F13EE 4C 8B 04 01 mov r8,qword ptr [rcx+rax]
000000013F3F13F2 48 8B 05 1F 06 04 00 mov rax,qword ptr [ValoriPreCalcolatiCos (13F431A18h)]
000000013F3F13F9 49 F7 DC neg r12
000000013F3F13FC 4A 8B 04 E8 mov rax,qword ptr [rax+r13*8]
273: dydx[i]=0;
000000013F3F1400 48 8D 6C DF 10 lea rbp,[rdi+rbx*8+10h]
000000013F3F1405 48 8B 14 01 mov rdx,qword ptr [rcx+rax]
000000013F3F1409 49 8B C8 mov rcx,r8
000000013F3F140C B8 02 00 00 00 mov eax,2
000000013F3F1411 48 2B C3 sub rax,rbx
000000013F3F1414 4C 8B DA mov r11,rdx
000000013F3F1417 49 2B CA sub rcx,r10
000000013F3F141A 48 C1 E0 03 shl rax,3
000000013F3F141E 4D 2B DA sub r11,r10
000000013F3F1421 48 2B CF sub rcx,rdi
000000013F3F1424 48 89 84 24 A8 00 00 00 mov qword ptr [lam],rax
000000013F3F142C 48 2B C7 sub rax,rdi
000000013F3F142F 4C 2B DF sub r11,rdi
000000013F3F1432 4C 03 C0 add r8,rax
000000013F3F1435 48 03 C2 add rax,rdx
000000013F3F1438 48 8B 94 24 A8 00 00 00 mov rdx,qword ptr [lam]
000000013F3F1440 48 89 44 24 20 mov qword ptr [rsp+20h],rax
000000013F3F1445 48 8B C3 mov rax,rbx
000000013F3F1448 48 8B 5C 24 20 mov rbx,qword ptr [rsp+20h]
000000013F3F144D 48 C1 E8 02 shr rax,2
000000013F3F1451 48 8D 34 85 01 00 00 00 lea rsi,[rax*4+1]
000000013F3F1459 0F 1F 80 00 00 00 00 nop dword ptr [rax]
276: {
277: q1=Kmod*ValoriPreCalcolatiSin[Indice][i][j];
000000013F3F1460 66 0F 28 D7 movapd xmm2,xmm7
278: q2=Kmod*ValoriPreCalcolatiCos[Indice][i][j];
000000013F3F1464 66 0F 28 DF movapd xmm3,xmm7
281: dydx[i+n]+=q2*y[j]-q1*y[j+n];
000000013F3F1468 48 83 C5 20 add rbp,20h
000000013F3F146C 48 FF C8 dec rax
000000013F3F146F F2 41 0F 59 5C 2B D8 mulsd xmm3,mmword ptr [r11+rbp-28h]
000000013F3F1476 F2 0F 59 54 29 D8 mulsd xmm2,mmword ptr [rcx+rbp-28h]
000000013F3F147C 66 0F 28 CB movapd xmm1,xmm3
000000013F3F1480 66 0F 28 C2 movapd xmm0,xmm2
000000013F3F1484 F2 41 0F 59 44 2C D8 mulsd xmm0,mmword ptr [r12+rbp-28h]
000000013F3F148B F2 0F 59 4D D8 mulsd xmm1,mmword ptr [rbp-28h]
000000013F3F1490 F2 0F 58 C8 addsd xmm1,xmm0
000000013F3F1494 F2 41 0F 10 01 movsd xmm0,mmword ptr [r9]
000000013F3F1499 F2 0F 5C C1 subsd xmm0,xmm1
000000013F3F149D F2 41 0F 11 01 movsd mmword ptr [r9],xmm0
000000013F3F14A2 F2 41 0F 59 5C 2C D8 mulsd xmm3,mmword ptr [r12+rbp-28h]
000000013F3F14A9 F2 0F 59 55 D8 mulsd xmm2,mmword ptr [rbp-28h]
000000013F3F14AE F2 0F 5C DA subsd xmm3,xmm2
000000013F3F14B2 66 0F 28 D7 movapd xmm2,xmm7
000000013F3F14B6 F2 43 0F 58 1C 0A addsd xmm3,mmword ptr [r10+r9]
000000013F3F14BC F2 43 0F 11 1C 0A movsd mmword ptr [r10+r9],xmm3
000000013F3F14C2 66 0F 28 DF movapd xmm3,xmm7
000000013F3F14C6 F2 0F 59 54 29 E0 mulsd xmm2,mmword ptr [rcx+rbp-20h]
000000013F3F14CC F2 41 0F 59 5C 2B E0 mulsd xmm3,mmword ptr [r11+rbp-20h]
000000013F3F14D3 66 0F 28 C2 movapd xmm0,xmm2
000000013F3F14D7 F2 42 0F 59 44 F5 E0 mulsd xmm0,mmword ptr [rbp+r14*8-20h]
000000013F3F14DE 66 0F 28 CB movapd xmm1,xmm3
000000013F3F14E2 F2 0F 59 4D E0 mulsd xmm1,mmword ptr [rbp-20h]
000000013F3F14E7 F2 0F 58 C8 addsd xmm1,xmm0
000000013F3F14EB F2 41 0F 10 01 movsd xmm0,mmword ptr [r9]
000000013F3F14F0 F2 0F 5C C1 subsd xmm0,xmm1
000000013F3F14F4 F2 41 0F 11 01 movsd mmword ptr [r9],xmm0
000000013F3F14F9 F2 42 0F 59 5C F5 E0 mulsd xmm3,mmword ptr [rbp+r14*8-20h]
000000013F3F1500 F2 0F 59 55 E0 mulsd xmm2,mmword ptr [rbp-20h]
000000013F3F1505 F2 0F 5C DA subsd xmm3,xmm2
000000013F3F1509 66 0F 28 D7 movapd xmm2,xmm7
000000013F3F150D F2 43 0F 58 1C 0A addsd xmm3,mmword ptr [r10+r9]
000000013F3F1513 F2 43 0F 11 1C 0A movsd mmword ptr [r10+r9],xmm3
000000013F3F1519 66 0F 28 DF movapd xmm3,xmm7
000000013F3F151D F2 0F 59 54 29 E8 mulsd xmm2,mmword ptr [rcx+rbp-18h]
000000013F3F1523 F2 41 0F 59 5C 2B E8 mulsd xmm3,mmword ptr [r11+rbp-18h]
000000013F3F152A 66 0F 28 C2 movapd xmm0,xmm2
000000013F3F152E F2 41 0F 59 44 2C E8 mulsd xmm0,mmword ptr [r12+rbp-18h]
000000013F3F1535 66 0F 28 CB movapd xmm1,xmm3
000000013F3F1539 F2 0F 59 4D E8 mulsd xmm1,mmword ptr [rbp-18h]
000000013F3F153E F2 0F 58 C8 addsd xmm1,xmm0
000000013F3F1542 F2 41 0F 10 01 movsd xmm0,mmword ptr [r9]
000000013F3F1547 F2 0F 5C C1 subsd xmm0,xmm1
000000013F3F154B F2 41 0F 11 01 movsd mmword ptr [r9],xmm0
000000013F3F1550 F2 41 0F 59 5C 2C E8 mulsd xmm3,mmword ptr [r12+rbp-18h]
000000013F3F1557 F2 0F 59 55 E8 mulsd xmm2,mmword ptr [rbp-18h]
000000013F3F155C F2 0F 5C DA subsd xmm3,xmm2
000000013F3F1560 66 0F 28 D7 movapd xmm2,xmm7
000000013F3F1564 F2 43 0F 58 1C 0A addsd xmm3,mmword ptr [r10+r9]
000000013F3F156A F2 43 0F 11 1C 0A movsd mmword ptr [r10+r9],xmm3
000000013F3F1570 66 0F 28 DF movapd xmm3,xmm7
000000013F3F1574 F2 41 0F 59 54 28 E0 mulsd xmm2,mmword ptr [r8+rbp-20h]
000000013F3F157B F2 0F 59 5C 2B E0 mulsd xmm3,mmword ptr [rbx+rbp-20h]
000000013F3F1581 66 0F 28 C2 movapd xmm0,xmm2
000000013F3F1585 F2 0F 59 44 2A E0 mulsd xmm0,mmword ptr [rdx+rbp-20h]
000000013F3F158B 66 0F 28 CB movapd xmm1,xmm3
000000013F3F158F F2 0F 59 4D F0 mulsd xmm1,mmword ptr [rbp-10h]
000000013F3F1594 F2 0F 58 C8 addsd xmm1,xmm0
000000013F3F1598 F2 41 0F 10 01 movsd xmm0,mmword ptr [r9]
000000013F3F159D F2 0F 5C C1 subsd xmm0,xmm1
000000013F3F15A1 F2 41 0F 11 01 movsd mmword ptr [r9],xmm0
000000013F3F15A6 F2 0F 59 5C 2A E0 mulsd xmm3,mmword ptr [rdx+rbp-20h]
000000013F3F15AC F2 0F 59 55 F0 mulsd xmm2,mmword ptr [rbp-10h]
000000013F3F15B1 F2 0F 5C DA subsd xmm3,xmm2
000000013F3F15B5 F2 43 0F 58 1C 0A addsd xmm3,mmword ptr [r10+r9]
000000013F3F15BB F2 43 0F 11 1C 0A movsd mmword ptr [r10+r9],xmm3
000000013F3F15C1 0F 85 99 FE FF FF jne derivsU+160h (13F3F1460h)
000000013F3F15C7 48 63 1D B2 EE 03 00 movsxd rbx,dword ptr [N (13F430480h)]
000000013F3F15CE 48 8B AC 24 90 00 00 00 mov rbp,qword ptr [x]
274: dydx[i+n]=0;
275: for(j=1;j<=n;j++)
000000013F3F15D6 48 3B F3 cmp rsi,rbx
000000013F3F15D9 0F 8F A4 00 00 00 jg derivsU+383h (13F3F1683h)
000000013F3F15DF 4B 8D 0C 0F lea rcx,[r15+r9]
000000013F3F15E3 48 8D 04 33 lea rax,[rbx+rsi]
000000013F3F15E7 48 8D 14 C7 lea rdx,[rdi+rax*8]
000000013F3F15EB 48 8B 05 1E 04 04 00 mov rax,qword ptr [ValoriPreCalcolatiSin (13F431A10h)]
000000013F3F15F2 4A 8B 04 E8 mov rax,qword ptr [rax+r13*8]
000000013F3F15F6 4C 8B 04 08 mov r8,qword ptr [rax+rcx]
000000013F3F15FA 48 8B 05 17 04 04 00 mov rax,qword ptr [ValoriPreCalcolatiCos (13F431A18h)]
000000013F3F1601 4A 8B 04 E8 mov rax,qword ptr [rax+r13*8]
000000013F3F1605 4D 2B C2 sub r8,r10
000000013F3F1608 4C 8B 1C 08 mov r11,qword ptr [rax+rcx]
000000013F3F160C 48 8B C3 mov rax,rbx
000000013F3F160F 48 8B CB mov rcx,rbx
000000013F3F1612 4D 2B DA sub r11,r10
000000013F3F1615 48 2B C6 sub rax,rsi
000000013F3F1618 4C 2B C7 sub r8,rdi
000000013F3F161B 4C 2B DF sub r11,rdi
000000013F3F161E 48 F7 D9 neg rcx
000000013F3F1621 48 FF C0 inc rax
276: {
277: q1=Kmod*ValoriPreCalcolatiSin[Indice][i][j];
000000013F3F1624 66 0F 28 D7 movapd xmm2,xmm7
278: q2=Kmod*ValoriPreCalcolatiCos[Indice][i][j];
000000013F3F1628 66 0F 28 DF movapd xmm3,xmm7
000000013F3F162C 48 83 C2 08 add rdx,8
000000013F3F1630 48 FF C8 dec rax
000000013F3F1633 F2 41 0F 59 5C 13 F8 mulsd xmm3,mmword ptr [r11+rdx-8]
000000013F3F163A F2 41 0F 59 54 10 F8 mulsd xmm2,mmword ptr [r8+rdx-8]
279:
280: dydx[i]-=q1*y[j]+q2*y[j+n];
000000013F3F1641 66 0F 28 C3 movapd xmm0,xmm3
000000013F3F1645 66 0F 28 CA movapd xmm1,xmm2
000000013F3F1649 F2 0F 59 42 F8 mulsd xmm0,mmword ptr [rdx-8]
000000013F3F164E F2 0F 59 4C CA F8 mulsd xmm1,mmword ptr [rdx+rcx*8-8]
000000013F3F1654 F2 0F 58 C8 addsd xmm1,xmm0
000000013F3F1658 F2 41 0F 10 01 movsd xmm0,mmword ptr [r9]
000000013F3F165D F2 0F 5C C1 subsd xmm0,xmm1
000000013F3F1661 F2 41 0F 11 01 movsd mmword ptr [r9],xmm0
281: dydx[i+n]+=q2*y[j]-q1*y[j+n];
000000013F3F1666 F2 0F 59 5C CA F8 mulsd xmm3,mmword ptr [rdx+rcx*8-8]
000000013F3F166C F2 0F 59 52 F8 mulsd xmm2,mmword ptr [rdx-8]
000000013F3F1671 F2 0F 5C DA subsd xmm3,xmm2
000000013F3F1675 F2 43 0F 58 1C 0A addsd xmm3,mmword ptr [r10+r9]
000000013F3F167B F2 43 0F 11 1C 0A movsd mmword ptr [r10+r9],xmm3
000000013F3F1681 75 A1 jne derivsU+324h (13F3F1624h)
270: double q1,q2;
271:
272: for(i=1;i<=n;i++){
000000013F3F1683 49 83 C1 08 add r9,8
000000013F3F1687 48 FF CD dec rbp
000000013F3F168A B8 00 00 00 00 mov eax,0
000000013F3F168F 48 89 AC 24 90 00 00 00 mov qword ptr [x],rbp
000000013F3F1697 0F 85 23 FD FF FF jne derivsU+0C0h (13F3F13C0h)
000000013F3F169D 4C 8B 74 24 50 mov r14,qword ptr [rsp+50h]
000000013F3F16A2 4C 8B 64 24 58 mov r12,qword ptr [rsp+58h]
000000013F3F16A7 48 8B 74 24 60 mov rsi,qword ptr [rsp+60h]
000000013F3F16AC 48 8B AC 24 98 00 00 00 mov rbp,qword ptr [y]
282: }
283: }
284: #endif
285: }
000000013F3F16B4 0F 28 7C 24 30 movaps xmm7,xmmword ptr [rsp+30h]
000000013F3F16B9 48 83 C4 68 add rsp,68h
000000013F3F16BD 41 5F pop r15
000000013F3F16BF 41 5D pop r13
000000013F3F16C1 5F pop rdi
000000013F3F16C2 5B pop rbx
000000013F3F16C3 C3 ret
Come si può notare usa un sacco di registri xmm, in particolare usa questi registri all'interno dei cicli, quindi immagino che faccia eseguire alla CPU più operazioni in simultanea per ogni thread (dovrebbero essere 2 trattandosi di dati di tipo double).
Comunque non pensavo che ci fossero tutte queste istruzioni per fare due prodotti e due somme..., vabbè che ottimizza e sfrutta SSE però... :D
Come si può notare usa un sacco di registri xmm, in particolare usa questi registri all'interno dei cicli, quindi immagino che faccia eseguire alla CPU più operazioni in simultanea per ogni thread (dovrebbero essere 2 trattandosi di dati di tipo double).
Comunque non pensavo che ci fossero tutte queste istruzioni per fare due prodotti e due somme..., vabbè che ottimizza e sfrutta SSE però... :D
Si, sono molte, ma considera che per ogni xmm ti risparmi 4 operazioni singole se sei in singola precisione o due in doppia. Insomma, senza le SSE, avresti almeno il doppio delle istruzioni, o il doppio delle iterazioni.
E con le prossime AVX arriveranno a 8 e 4 rispettivamente.. ;)
Si, sono molte, ma considera che per ogni xmm ti risparmi 4 operazioni singole se sei in singola precisione o due in doppia. Insomma, senza le SSE, avresti almeno il doppio delle istruzioni, o il doppio delle iterazioni.
E con le prossime AVX arriveranno a 8 e 4 rispettivamente.. ;)
Bè allora mi comprerò un bel Sandy Bridge o un Bulldozer da 8 core :sofico: , così passando da 2 o 8 core e con le istruzioni AVX riduco il tempo di elaborazione ad 1/8 (o anche meno considerando l'aumento di IPC che c'è con due generazioni di mezzo di CPU) :D .
Comunque chiaramente ora il tempo di elaborazione è molto più che accettabile, sono pienamente soddisfatto dei progressi che ho fatto. Magari un giorno imparerò l'assembly per futuri programmi CPU-limited, anche se dubito che riuscierei ad ottenere un codice assembly più performante di quello generato dai compilatori odierni.
Prima di precalcolarti tutti i seni e coseni, avrei provato anche ad usare il buon vecchio sviluppo di Taylor (limitandoti al terzo ordine che da risultati più che accettabili). Certo..se non hai problemi di memoria, una LUT è sempre la scelta ideale..in caso contrario c'è Taylor! ;)
Prima di precalcolarti tutti i seni e coseni, avrei provato anche ad usare il buon vecchio sviluppo di Taylor (limitandoti al terzo ordine che da risultati più che accettabili). Certo..se non hai problemi di memoria, una LUT è sempre la scelta ideale..in caso contrario c'è Taylor! ;)
Ma taylor presuppone uno sviluppo intorno ad un certo punto iniziale e man mano che ci si allontana dal punto iniziale l'approssimazione peggiora, nel mio caso non c'è un limite al valore all'argomento del seno coseno, o comunque se c'è è molto grande (anche svariate centinaia). Quindi dovrei riportare il valore in un range tramite la funzione fmod, ma in tal caso dato converrebbe forse creare un vettore con i valori tabellati nel range con la funzione sin e cos dato che il numero di valori tabellati è molto più esiguo.
O forse mi sbaglio?
Vorrei sottolineare, forse è sfuggito, che i valori precalcolati sono esattamente i valori che mi servirebbero e sono calcolati in maniera esatta tramite il semplice sin e cos. Inoltre si ottiene un notevole risparmio di tempo anche se i vettori hanno una dimensione molto ragguardevole (ho provato anche con vettori da 500MB, probabilmente ciò è dovuto alla cache da 6MB e al fatto che la funzione utilizzerà i valori della tabella in maniera sequenziale, quindi saranno caricati chessò 4-5MB della tabella nella cache e sfruttati tutti, poi altri 4-5MB e così via).
Aggiungo altre informazioni che potrebbero risultare interessanti ad altri utenti:
Ho provato la versione trial delle librerie MKL.
1) A differenza delle librerie ACML, è possibile con il compilatore della Microsoft sia utilizzare le librerie dinamiche (permesso anche dalle ACML) che quelle statiche. O forse sono io incapace, eppure nel manuale delle ACML non si fa menzione di utilizzo delle librerie statiche con il compilatore Microsoft.
2) Le librerie MKL risultano più efficienti, mi consentono di ottenere una ulteriore riduzione del tempo. Nello stesso scenario del precedente post, nel quale riuscivo ad ottenere come migliore risultato 45'', con le MKL ottengo 35'' (un'altra ragguardevole riduzione del tempo)
3) Il compilatore Microsoft secondo me fa comunque un ottimo lavoro! , considerando che deve interpretare il codice che ho scritto e cercare di sfruttare SSE nel miglior modo possibile senza sapere che alla fine il codice si tratta di prodotti vettore-matrice, mentre utilizzando le MKL ho compattato il tutto a 4 prodotti vettore-matrice.
Cmq anche a me è successo di dover ottimizzare in qualche modo il calcolo di seni e coseni e me la sono cavata con una look up table moooooolto piccola (100-500 valori per il seno e altrettanti per il coseno) + un'interpolazione cubica. La precisione, per quel che dovevo fare io, era ottimo (in valore assoluto la differenza tra il risultato "vero" e quello trovato in questa maniera era ~10^-6 - 10^-9)
rеpne scasb
13-10-2010, 11:48
■
Mm analisi matematica per l'ottimizzazzione, questo thread è interessantissimo :D
Continuate così!
Come vedi per calcolare il seno di k per k qualsiasi e' sufficiente conoscere solo il seno tra 0 e P greco mezzi.
Anche lui aveva parlato di fmod ;) Nel tuo caso c'è una divisione FP, che è molto vicina alla fmod come tempo di esecuzione.
Ma fare una fmod è più veloce di una LUT ampia ?
L'unico caso probabilmente è in presenza di un cache miss della parte di LUT interessata, ma non sono comunque sicuro che sia più lento.
vBulletin® v3.6.4, Copyright ©2000-2025, Jelsoft Enterprises Ltd.