Visualizza il feed RSS

Titolo provvisorio...

Un problemino con le addizioni... +2

Valutazione: 2 voti, con una media del 5.00.
di pubblicato il 28-03-2010 alle 04:43 (4475 Visite)
Desidero ringraziare i numerosi lettori, anche non iscritti a MasterDrive, che in privato hanno manifestato un interesse superiore alle aspettative per la funzione di partizione e relativa implementazione: tanto che ho ritenuto opportuno dedicare questa entry in modo specifico all'algoritmo di Kreher & Stinson, il quale - al di là dell'utilità specifica e della curiosità individuale - ci dà modo di fare alcune considerazioni sperabilmente interessanti.

Ricordo innanzi tutto che l'algoritmo originariamente proposto nel testo "Combinatorial Algorithms: Generation, Enumeration and Search" (edito nel 1998) è stato implementato dagli autori in linguaggio C e compilato, all'epoca, con GCC. L'originale IntegerPartitions.c è contenuto nell'archivio GEN.tar.gz, liberamente scaricabile presso la sezione "source code" della home page di supporto al testo.

Come funziona questo algoritmo ? Vediamo di delinearne almeno l'idea generale. La relazione di ricorrenza utilizzata è la seguente, già nota al grande matematico svizzero Leonhard Euler:

Formula LaTeX: \\P(0)\ =\ 1\\\\P(n)\ =\ \sum_{\overset{j\,\ge\,1}{\omega\,\le\,n}}(-1)^{j+1}\,P\left(n \,-\, \omega\right )\ \ +\ \sum_{\overset{j\,\ge\,1}{\omega'\,\le\,n}}(-1)^{j+1}\,P\left(n \,-\, \omega' \right)\\\\\\\text{dove } \omega\ =\ \frac{3j^2 \,-\,j}{2},\ \ \omega'\ =\ \frac{3j^2 \,+\,j}{2}

La formulazione qui proposta è strettamente analoga a quella originariamente riportata da Kreher & Stinson per i loro scopi, con un'unica aggiunta: l'esplicitazione delle variabili intermedie omega, che in realtà rappresentano numeri pentagonali generalizzati - variabili da loro stessi introdotte e utilizzate nell'algoritmo.

La complessità computazionale di questa soluzione è accettabile, a grandi linee: l'algoritmo esegue solo Formula LaTeX: \smash{\Theta (n^{\frac{3}{2}})} somme. Sicuramente è migliorativo rispetto all'algoritmo 3.5 presente nel testo (e nel sorgente indicato), che opera in Formula LaTeX: \smash{\Theta(n^2)} se utilizzato per il calcolo di P(n).

Il codice indicato, tuttavia, non somiglia troppo strettamente a codesta sommatoria, almeno a prima vista: cerchiamo di analizzare brevemente queste differenze.

Come d'abitudine in questi casi, proviamo a farci venire qualche semplice idea per eliminare alcune delle operazioni computazionalmente più costose nel loop: tipicamente le moltiplicazioni tra variabili. Per cominciare, scriviamo esplicitamente la differenza tra due generici valori consecutivi di omega:

Formula LaTeX: \begin{align*}\\\omega_j\ -\ \omega_{j-1}\ &=\ \frac{3j^2\,-\,j}{2}\ -\ \frac{3(j\,-\,1)^2\,-\,(j\,-\,1)}{2}\ &=\ 3j\,-\,2\\\\\\\omega_{j+1}\ -\ \omega_{j}\ &=\ \frac{3(j\,+\,1)^2\,-\,(j\,+\,1)}{2}\ -\ \frac{3j^2\,-\,j}{2}\ &=\ 3j\,+\,1\end{align*}

Ovvero:

Formula LaTeX: \bold{\color{Red}(1)\;\;\omega_j\;=\;\omega_{j-1}\,+\, 3j\,-\,2}

Formula LaTeX: \bold{\color{Red}(2)\;\;\omega_{j+1}\;=\;\omega_j \,+\, 3j\,+\,1}

La semplificazione ottenuta è già utile: scelta ad esempio la (1), accumulando in omega la quantità 3j - 2 ad ogni iterazione, si ottiene la successione desiderata, eliminando operazioni relativamente (e misurabilmente !) costose.

Ma non è tutto: applicando ricorsivamente il procedimento, si può determinare allo stesso modo il passo di incremento di codesta espressione, partendo nuovamente dalla (1) per fissare le idee:

Formula LaTeX: \Delta_j\ =\ 3(j\,+\,1)\,-\,2\ -\ 3j\,+\,2\ =\ 3

Questo elimina definitivamente ogni moltiplicazione nella sezione di aggiornamento delle variabili connesse a quella di induzione.

Ricaviamo banalmente anche il secondo indice (nel codice originale viene denominato omega2) in funzione di omega:

Formula LaTeX: \omega'_j \,-\, \omega_j\ =\ \frac{3j^2 \,+\,j}{2}\ - \ \frac{3j^2 \,-\,j}{2}\ =\ j

Da cui:

Formula LaTeX: \bold{\color{red}{(3)\ \ \omega'_j\ =\ \omega_j \,+\, j}}

Leggendo il codice originale, ci sono ancora almeno due aspetti che spesso lasciano perplessi giovani sviluppatori e studenti nonostante questa spiegazione (omessa dagli autori per la sua evidenza), e riguardano proprio l'applicazione delle semplici relazioni appena ricavate.

E' piuttosto ovvio che, una volta eliminato aritmeticamente il ricorso a moltiplicazioni, divisioni e potenze (con l'incombente spettro dell'uso di FP da parte del compilatore) ad ogni passo di iterazione, le opzioni di riscrittura e ridistribuzione delle operazioni per queste equazioni sono varie, ma quasi tutte computazionalmente equivalenti, a meno di minime differenze da evidenziare con misure puntuali su un preciso target.

Tuttavia nel codice originale, invece di ricorrere ad un idioma immediato come il seguente:

Formula LaTeX: \begin{cases}j & \leftarrow j\,+\,1\\\omega & \leftarrow \omega\,+\,3j\,-\,2\\\omega' & \leftarrow \omega \,+\, j\end{cases}

o meglio ancora:

Formula LaTeX: \begin{cases}j & \leftarrow j\,+\,1\\k & \leftarrow k\,+\,3\\\omega & \leftarrow \omega\,+\,k\\\omega' & \leftarrow \omega \,+\, j\end{cases}

gli autori hanno optato per una forma diversa, con l'idea di minimizzare il numero complessivo di addizioni e incrementi unitari generati dal compilatore di riferimento, ossia GCC. Ciò, all'epoca in cui è stato concepito e sperimentato il codice, nella seconda metà dei Novanta, aveva ampiamente senso sulle CPU x86 (e non solo). Tenendo opportunamente conto della dinamica della variabile di induzione j, si giunge facilmente all'implementazione suggerita, che utilizza la formula (2) per l'incremento di omega a monte dell'incremento della variabile di induzione, e la formula (3) per omega2 a valle di tale incremento.

Formula LaTeX: \begin{cases}\omega & \leftarrow \omega\,+\,3j\,+\,1\\ j & \leftarrow j\,+\,1\\ \omega' & \leftarrow \omega \,+\, j \end{cases}

Ulteriori accorgimenti implementativi presenti nel codice originale, probabilmente ancora più importanti con la combinazione di compilatore e hardware scelta all'epoca dagli autori, erano volti ad eliminare l'uso di ogni possibile istruzione IMUL in tutto il loop interno, facendo uso in particolare di if() per scegliere il segno del termine e affidandosi alla gestione della eventuale BPU (Branch Prediction Unit) nell'optimizer.

Ritengo tuttavia doveroso ricordare a tutti i miei lettori che, prima di ricorrere a ottimizzazioni di quest'ultimo tipo nel codice C, è sempre indispensabile valutare varie soluzioni e combinazioni di flag di ottimizzazione, usando correttamente il profiler e verificando con cura il codice assembly generato dal compilatore per la macchina target, per non andare incontro a strane "sorprese"...

A solo scopo illustrativo e didascalico ho condotto alcuni semplici test su XP Pro SP3, in maniera estremamente sbrigativa e priva di pretese. Il fine ultimo, ça va sans dire, è unicamente quello di stimolare il senso critico e l'attitudine sperimentale che sempre devono accompagnare il lavoro di uno sviluppatore professionista.

Purtroppo in questi casi si rivela sempre necessario il tedioso disclaimer anti-troll: ribadisco che queste note non pretendono in alcun modo di sostituire benchmark sistematici di ottimizzazione dei compilatori, dei quali peraltro la letteratura abbonda; d'altro canto, lungi da me anche solo l'idea di intentare in questa sede esilaranti e grotteschi "processi" alla modernità, a questo o quel compilatore, o ad una (presunta) carente usabilità dei set di opzioni.

Quel che invece interessa mostrare qui è che, pur attivando tutte le opzioni più ragionevoli di ottimizzazione dei vari compilatori mainstream, certe "ottimizzazioni" nel sorgente finiscono per sortire risultati perfino controproducenti. Ciò ha anche, in definitiva, l'effetto collaterale di mostrare con quanta facilità ancor oggi si può aggirare o piegare ai propri scopi la presunta "intelligenza" dei compilatori moderni, e come - con altrettanta facilità - è facile inciamparvi a proprio detrimento, per eccesso di fideismo.

D'altro canto questi banali test (seppure frettolosi e superficiali) sono totalmente riproducibili, in perfetta trasparenza. Anzi, al fine di consentire una facile replicazione delle misure anche a chi non utilizza un profiler, il sorgente si basa - in modo primitivo ma efficace - sulla nota chiamata a GetTickCount() ed esegue un totale di un milione di iterazioni per ciascuna routine richiamata.

Un confronto in parallelo con il file di profiling generato da un eseguibile "instrumented" di Open Watcom (con l'opzione -et, che sostanzialmente inserisce codice basato su RDTSC) riconferma comunque che gli indici ricavati con tale metodo risultano accurati almeno fino alla quarta cifra decimale, più che sufficiente per i nostri scopi.

Scorrendo il sorgente allegato, si nota che le differenze tra le varie implementazioni proposte sono minime, e piuttosto evidenti. Tra le numerose possibilità ho scelto arbitrariamente alcune delle più sensate, alla luce della mia esperienza: fermo restando che in casi come questo l'ottimizzazione assoluta, per principio, consiste nell'uso di una lookup table (LUT) preparata tramite un applicativo dedicato di calcolo numerico.

La funzione EnumParts2() ricalca fedelmente l'originale di Kreher & Stinson (tranne banalmente nell'uso degli operatori di assegnazione composita e pre/postincremento) e tutti gli indici prestazionali ricavati sono relativi a tale routine.

Le altre funzioni equivalenti sono grossolanamente numerate per prestazioni attese crescenti, quindi NaiveParts1() è in linea di principio il peggior codice (il più intuitivo e pedissequo) mentre ci si attende che NaiveParts3() e NaiveParts4() offrano le migliori prestazioni su una macchina wintel - specialmente in relazione a EnumParts2().

Più in dettaglio, NaiveParts1() implementa in modo acritico le espressioni aritmetiche così come scritte, sperando nelle potenzialità taumaturgiche dell'ottimizzatore. Speranze mal riposte, come confermano le misure.

NaiveParts2() utilizza il seguente idioma per l'incremento di omega, identico all'originale K&S:
codice:
            omega += 3 * j +1;
            ++j;
            sign = -sign;
Le differenze prestazionali tra questa funzione e l'implementazione K&S di riferimento sono quindi da imputarsi principalmente al diverso modo di gestire il segno del termine additivo.

Per contro, NaiveParts3() differisce dalle altre funzioni poiché utilizza in modo esplicito un totalizzatore d'appoggio, denominato k, che elimina completamente le moltiplicazioni nella relativa sezione di codice:
codice:
            ++j;
            k += 3;
            omega += k;
            sign = -sign;
Infine NaiveParts4() utilizza in modo esplicito shift e somme, tornando all'espressione (3j +1) già utilizzata:
codice:
            ++j;
            k = (j << 1) + j;
            omega += ++k;
            sign = -sign;
Vale la pena di mostrare come i vari compilatori hanno ottimizzato l'idioma utilizzato nella NaiveParts2():

codice:
******************************************************
**** Visual C++ Xpress *******************************
******************************************************
;; ESI = omega, EDX = sign, EDI = j
;
; 116  :             omega += 3 * j + 1;

        mov    esi, DWORD PTR tv219[esp+20]
        lea    eax, DWORD PTR [eax+esi+1]

; 117  :             ++j;

        add    esi, 3
        inc    edi

; 118  :             sign = -sign;

        neg    edx

******************************************************
**** Turbo C++ 2008 **********************************
******************************************************
?live1@1008: ; EAX = omega, EDX = i, ECX = sign, ESI = j
;
; 116  :             omega += 3 * j + 1;
; 117  :             ++j;
; 118  :             sign = -sign;
;
@47:
@45:
        lea       ebx,dword ptr [esi+2*esi]
        inc       ebx
        add       eax,ebx
        inc       esi
        neg       ecx

******************************************************
**** Digital Mars 850 ********************************
******************************************************
;; EBX = omega, EDI = j, ECX = sign
;
; 116  :             omega += 3 * j + 1;
; 117  :             ++j;
; 118  :             sign = -sign;
;
L1BC:       
        add    EBX, 010h[ESP]
        inc    EDI
        add    dword ptr 010h[ESP],3
        neg    ECX
Per la cronaca, anche Open Watcom 1.8 e perfino BCC 3.1 (su un ottimo 486 DX4/100) hanno ottimizzato come Digital Mars, utilizzando una variabile automatica al fine di massimizzare la velocità dell'eseguibile.

Nel codice prodotto dai compilatori che non hanno calcolato il valore d'incremento della successione, si nota come la pestilenziale IMUL (denotata nella ISA 486 come capace di sprecare fino a 42 cicli con quel tipo di indirizzamento, senza contare il caso di cache stall), è stata quantomeno sostituita con il noto idioma di moltiplicazione implicita tramite LEA.
Tralascio ogni ulteriore considerazione in merito: intelligentibus pauca ...

Di seguito riporto gli switch di compilazione significativi e una concisa caratterizzazione delle piattaforme di test:
codice:
Digital Mars 850....: -mn -o+time -WA -f -5 -a8 -c -cod -v2
Open Watcom 1.8.....: -w4 -e25 -zp8 -zq -ot&exan -ob -ol -ol+ -oc -oi -oa -or -oh -om -fp5 -fpi87 -5r -bt=nt -fo=.obj -m&f
Turbo C 2006........: -DNDEBUG -DNO_STRICT;_NO_VCL -6 -B -d -dc -f -Hc -w-par -O2 -b- -k- -vi -tWC -tWM- -tW- -c
VC 2008 Xpress......: /O2 /Ot /D "WIN32" /D "NDEBUG" /D "_CONSOLE" /FD /MT /GS- /fp:fast /FAs /W3 /nologo /c /Wp64 /Zi /TC

HW01: Intel Xeon 5110 @ 1,60 GHz, 2 Gb RAM, XP Pro SP3
HW02: Intel Pentium M @ 2,00 GHz, 2 Gb RAM, XP Pro SP3
HW03: Intel Pentium 4 @ 2,66 GHz, 1 Gb RAM, XP Pro SP3
HW04: Intel Atom N270 @ 1,60 GHz, 1 Gb RAM, XP Pro SP3
Il medesimo codice di partenza, senza modifica alcuna, è stato separatamente compilato con i quattro compilatori sopra elencati come applicazione consolle, utilizzando le opzioni evidenziate e rinominando opportunamente i quattro eseguibili prodotti.

I risultati medi ottenuti su una delle macchine di test (HW03) sono riportati di seguito a scopo illustrativo. Non sono stati in alcun modo ulteriormente compressi i binari: le dimensioni dell'eseguibile sono quelle ottenute "right out of the box".

La denominazione degli eseguibili è tale che volutamente non consente di risalire al compilatore utilizzato.

>intpart01, 126'464 bytes
## Rough performance test with 1000000 iterations
** [EnumParts2]: 5766 ticks
** [NaiveParts4]: 7641 ticks
** [NaiveParts3]: 6421 ticks
** [NaiveParts2]: 7297 ticks
** [NaiveParts1]: 8860 ticks

## Relative performance indexes:
>> [NaiveParts4]: 1.32518
>> [NaiveParts3]: 1.1136
>> [NaiveParts2]: 1.26552
>> [NaiveParts1]: 1.53659

>intpart02, 69'120 bytes
## Rough performance test with 1000000 iterations
** [EnumParts2]: 5234 ticks
** [NaiveParts4]: 4063 ticks
** [NaiveParts3]: 4078 ticks
** [NaiveParts2]: 4484 ticks
** [NaiveParts1]: 7516 ticks

## Relative performance indexes:
>> [NaiveParts4]: 0.776271
>> [NaiveParts3]: 0.779136
>> [NaiveParts2]: 0.856706
>> [NaiveParts1]: 1.436

>intpart03, 52'252 bytes
## Rough performance test with 1000000 iterations
** [EnumParts2]: 5110 ticks
** [NaiveParts4]: 3968 ticks
** [NaiveParts3]: 4282 ticks
** [NaiveParts2]: 4078 ticks
** [NaiveParts1]: 8469 ticks

## Relative performance indexes:
>> [NaiveParts4]: 0.776517
>> [NaiveParts3]: 0.837965
>> [NaiveParts2]: 0.798043
>> [NaiveParts1]: 1.65734

>intpart04, 36'864 bytes
## Rough performance test with 1000000 iterations
** [EnumParts2]: 4938 ticks
** [NaiveParts4]: 4359 ticks
** [NaiveParts3]: 4859 ticks
** [NaiveParts2]: 4891 ticks
** [NaiveParts1]: 7500 ticks

## Relative performance indexes:
>> [NaiveParts4]: 0.882746
>> [NaiveParts3]: 0.984002
>> [NaiveParts2]: 0.990482
>> [NaiveParts1]: 1.51883
A prescindere dalle prestazioni assolute, che in effetti rimarcano differenze non trascurabili (sebbene non particolarmente significative) tra i vari compilatori, si possono fare due considerazioni generali:

1) L'uso delle if() (e quindi dei salti) tipico della EnumParts2() appare controproducente su macchine delle ultime due o tre generazioni e con quasi ogni compilatore. Il numero di indici inferiori all'unità, se si esclude la volutamente inefficiente NaiveParts1(), testimonia in modo incontrovertibile questo aspetto.

2) L'ottimizzazione aritmetica esplicita di NaiveParts3() e NaiveParts4() appare fortemente premiante praticamente con tutti i compilatori mainstream, pur avendo utilizzato a piene mani tutte le più sensate opzioni di ottimizzazione in velocità rese disponibili dagli ambienti di sviluppo, senza porre limitazioni sull'uso di variabili d'appoggio, coiling, loop unrolling e altre amenità.

Il file "profile.txt" allegato contiene, in formato "delimited" importabile dal vostro spreadsheet preferito, i dati completi delle sessioni di prova eseguite sulle quattro macchine specificate sopra, disponibili per confronti ed elaborazioni di ogni genere. Il file "intpart_c.txt" contiene, con audace spicco di fantasia, il sorgente C utilizzato per i test.

Al di là di questi aspetti sperimentali e del banale esempio, certamente divertenti ma poco significativi in assoluto, rimane l'importanza di riasserire che anche agli inizi del terzo millennio il miglior ottimizzatore rimane sempre situato tra le vostre orecchie, secondo il sacrosanto aforisma di Michael Abrash.

Sebbene le motivazioni della scuola di pensiero che mette in guardia dalle early optimizations siano di norma ragionevoli e sensate (ma è altra cosa, in realtà), sovente le puerili obiezioni del tipo "coi processori di oggi è inutile ottimizzare" sono solo la moderna traduzione del "nondum matura est", non è ancora matura, come disse la volpe di Fedro davanti al grappolo d'uva troppo alto per lei.

In effetti l'ottimizzazione esula sempre più dalla preparazione impartita oggi allo sviluppatore quadratico medio.

Qui si è mostrato, in modo necessariamente abbreviato e grezzo ma efficace, chiaro e facilmente riproducibile, come ancor oggi - su macchine recentissime, con compilatori mainstream usati da milioni di sviluppatori, usando combinazioni di switch ragionevoli e immediatamente accessibili dalla IDE - una semplice combinazione di riscritture aritmetiche e sano buon senso logico, sistematica analisi del codice prodotto, misure puntuali e uso intelligente dell'ottimizzatore su un semplicissimo sorgente d'esempio può far lievitare in modo tangibile le prestazioni del codice.

Un buon programmatore sa sempre valutare quando è il caso di pensare a simili migliorie. Un Vero Programmatore le applica istintivamente quando è necessario.

EDIT: La saga continua...
Anteprime allegati File allegati

Commenti