modo efficiente per prendere i poteri di un vettore

Question

Ho scritto un codice che utilizza numericamente polinomi di Legendre fino ad un alto n-esimo ordine. Per esempio:

.... 
case 8 
p = (6435*x.^8-12012*x.^6+6930*x.^4-1260*x.^2+35)/128; return 
case 9 
... 

Se il vettore x è lungo questo può diventare lento. Ho visto che c'è una differenza di prestazioni tra dire x.^4 e x.*x.*x.*x e ho pensato che potrei usare questo per migliorare il mio codice. Ho usato timeit e ha scoperto che per:

x=linspace(0,10,1e6); 
f1= @() power(x,4) 
f2= @() x.4; 
f3= @() x.^2.^2 
f4= @() x.*x.*x.*x 

f4 è più veloce da un fattore 2 rispetto al resto. Tuttavia quando vado a x.^6 c'è poca differenza tra (x.*x.*x).^2 e x.*x.*x.*x.*x.*x (mentre tutte le altre opzioni sono più lente).

C'è via per dire quale sarà il modo più efficiente per prendere il potere di un vettore? Puoi spiegare perché c'è una così grande differenza nelle prestazioni?

Answer 1

Questo non è esattamente una risposta alla tua domanda, ma può risolvere il problema:

x2 = x.*x; % or x.^2 or power(x,2), whichever is most efficient 
p = ((((6435*x2-12012)*x2+6930)*x2-1260)*x2+35)/128 

questo modo si fa il potere solo una volta, e solo con esponente 2. Questo trucco può essere applicato a tutti i Polinomi di Legendre (nei polinomi di grado dispari uno x2 è sostituito da x).

Answer 2

Ecco alcuni pensieri:

power(x,4) e x.^4 sono equivalenti (basta leggere il doc).

x.*x.*x.*x è probabilmente ottimizzato per qualcosa come x.^2.^2

---

x.^2.^2 è probabilmente valutato come: Prendere il quadrato di ciascun elemento (veloce), e prendere il quadrato di che ancora una volta (veloce di nuovo).

x.^4 è probabilmente valutato direttamente come: Prendi la quarta potenza di ciascun elemento (lento).

Non è così strano vedere che 2 operazioni veloci richiedono meno tempo di 1 operazione lenta. Peccato che l'ottimizzazione non venga eseguita nel caso power 4, ma forse non funzionerà sempre o non avrà un costo (controllo degli input, memoria?).

---

Informazioni sui tempi: in realtà c'è molta più differenza di un fattore 2!

come li chiami tu in una funzione ora, la funzione di testa viene aggiunto in ogni caso, rendendo le relative differenze più piccoli:

y=x;tic,power(x,4);toc 
y=x;tic,x.^4;toc 
y=x;tic,x.^2.^2;toc 
y=x;tic,x.*x.*x.*x;toc 

darà:

Elapsed time is 0.034826 seconds. 
Elapsed time is 0.029186 seconds. 
Elapsed time is 0.003891 seconds. 
Elapsed time is 0.003840 seconds. 

Così, quasi è una differenza di fattore 10. Tuttavia, si noti che la differenza di tempo in secondi è ancora minore, quindi per la maggior parte delle applicazioni pratiche vorrei semplicemente fare la semplice sintassi.

Answer 3

Sembra come se Mathworks presenta particolari piazze con carter nella sua funzione di potenza (purtroppo, è tutto integrato a sorgente chiuso che non possiamo vedere). Nei miei test su R2013b, sembra che .^, power e realpow utilizzino lo stesso algoritmo. Per le piazze, credo che abbiano un case speciale per essere x.*x.

1.0x (4.4ms): @()x.^2 
1.0x (4.4ms): @()power(x,2) 
1.0x (4.5ms): @()x.*x 
1.0x (4.5ms): @()realpow(x,2) 
6.1x (27.1ms): @()exp(2*log(x)) 

Per i cubi, la storia è diversa. Non sono più in scatola speciale. Anche in questo caso, .^, power e realpow tutti sono simili, ma molto più lento questa volta: salto

1.0x (4.5ms): @()x.*x.*x 
1.0x (4.6ms): @()x.*x.^2 
5.9x (26.9ms): @()exp(3*log(x)) 
13.8x (62.3ms): @()power(x,3) 
14.0x (63.2ms): @()x.^3 
14.1x (63.7ms): @()realpow(x,3) 

Let fino al 16 potere di vedere come questi algoritmi scala:

1.0x (8.1ms): @()x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x 
2.2x (17.4ms): @()x.^2.^2.^2.^2 
3.5x (27.9ms): @()exp(16*log(x)) 
7.9x (63.8ms): @()power(x,16) 
7.9x (63.9ms): @()realpow(x,16) 
8.3x (66.9ms): @()x.^16 

Quindi: .^, power e realpow vengono eseguiti tutti in un tempo costante rispetto all'esponente, a meno che non sia stato creato un involucro speciale (anche il formato -1 è stato classificato come speciale). L'uso del trucco exp(n*log(x)) è anche un tempo costante per quanto riguarda l'esponente e più veloce. L'unico risultato non capisco perché la quadratura ripetuta sia più lenta della moltiplicazione.

Come previsto, l'aumento delle dimensioni di x di un fattore di 100 aumenta il tempo in modo simile per tutti gli algoritmi.

Quindi, morale della storia? Quando si utilizzano esponenti interi scalari, eseguire sempre la moltiplicazione. C'è un sacco di intelligenza in power e amici (l'esponente può essere a virgola mobile, vettoriale, ecc.). Le uniche eccezioni sono dove Mathworks ha fatto l'ottimizzazione per te. Nel 2013b, sembra essere x^2 e x^(-1). Spero che aggiungeranno di più con il passare del tempo. Ma, in generale, l'esponenziazione è difficile e la moltiplicazione è facile. Nel codice sensibile alle prestazioni, non penso che tu possa sbagliare digitando sempre x.*x.*x.*x. (Naturalmente, nel tuo caso, seguire il consiglio Luis` e fare uso dei risultati intermedi all'interno di ogni termine!)

function powerTest(x) 

f{1} = @() x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x.*x; 
f{2} = @() x.^2.^2.^2.^2; 
f{3} = @() exp(16.*log(x)); 
f{4} = @() x.^16; 
f{5} = @() power(x,16); 
f{6} = @() realpow(x,16); 

for i = 1:length(f) 
    t(i) = timeit(f{i}); 
end 

[t,idxs] = sort(t); 
fcns = f(idxs); 

for i = 1:length(fcns) 
    fprintf('%.1fx (%.1fms):\t%s\n',t(i)/t(1),t(i)*1e3,func2str(fcns{i})); 
end