Spesso ci sono situazioni in cui guardando certe proprietà del problema numerico da risolvere, si possono trovare degli accorgimenti che risultano avere un impatto non trascurabile sul costo computazionale dell'algoritmo di calcolo.
Vediamo di seguito un paio di esempi di questo tipo.


1.  prodotto matrice-vettore


Supponiamo che sia necessario compiere il prodotto matrice-vettore A*b,  in cui la matrice "A", di dimensione n x n, è un campionamento della funzione 

f(x,y) = cos(x)*cos(y)
in una griglia cartesiana di punti nell'intervallo  [0, pi] x [0, pi] ,  e "b" è un vettore colonna di lunghezza n .

Scelti:
- un vettore riga di coordinate "x" avente n elementi scelti nell'intervallo [0, pi],  e calcolato il vettore "vx" tale che  vx(i) = cos( x(i) )
- un vettore riga di coordinate "y" avente n elementi scelti nell'intervallo [0, pi],  e calcolato il vettore "vy" tale che  vy(i) = cos( y(i) )
allora,
la matrice "A" può essere espressa nella forma 
A  =  vx' * vy

dove l'apice significa "trasposto".  In generale, essendo A una matrice densa, il prodotto A*b per un generico vettore "b" ha un costo dell'ordine di  O(n^2).

Se, però, utilizziamo l'espressione con cui viene costruita "A" ed invertiamo la precedenza nei calcoli:

vx' * (vy * b)
otteniamo una complessità dell'ordine di  O(n).

Esercizio:  costruire un programma Matlab/Octave che esegua il prodotto matrice-vettore A*b nelle due modalità qui prospettate e verifichi automaticamente che il risultato del prodotto nei due casi è il medesimo. Inoltre, verificare sperimentalmente, per alcuni valori di n, se l'andamento dei tempi di esecuzione segue l'andamento della complessità qui prospettato e costuire un grafico (n, tempi).



2.    Potenze di matrici

Per calcolare la potenza q-esima, q>=1, di una matrice "A", si può pensare di moltiplicare "A" q-1 volte per se stessa .
Ad esempio, con un ciclo for, in Python diventa:
Pa = A
for i in range(1,q)      #   i = 1,...,q-1
    Pa = matrixmultiply(Pa, A)

In modo più efficiente, possiamo esprimere "q" in base-2 (cioè in formato binario),  calcolarci i termini  A^(2^i),  i =  1,...,log2(q) = log(q)/log(2)   mediante la

  A^(2^i)  =  A^(2^(i-1))  *  A^(2^(i-1))

e sfruttare il fatto che 
A^(c1 + c2 + c3)  =  A^c1  *  A^c2  *  A^c3  .


Esercizio:  costruire un programma Matlab/Octave che, data una generica matrice "A" ed assegnato un valore di "q", esegua il calcolo di A^q  nelle due modalità qui esposte e verifichi automaticamente che il risultato nei due casi è il medesimo. Inoltre, valutare dal programma qual'è la complessità dei due algoritmi e verificare sperimentalmente, con una matrice "A" a scelta e valori di q = 4,5,...15 ,  se l'andamento dei tempi di esecuzione segue l'andamento della complessità atteso. Costruire un grafico (q, tempi).