Rivoluzione dei robot

Che cos'è l'«algebra lineare»?

Esistono molte definizioni tecniche, su [Wikipedia][linalg-wiki], nei manuali e su ottimi siti come [3Blue1Brown][3blue1brown], ma ai nostri fini possiamo essere più informali.

Ai maintainer di Julia (e di Python) questo genere di cose piace molto. Alla direzione di Exercism, molto meno.

Questo concetto sarà limitato agli aspetti più semplici di una materia immensa, ma attenzione: è inevitabilmente un concetto piuttosto matematico.

Che cos'è un vettore?

Dipende da chi lo chiedi e da come vuoi visualizzare qualcosa di piuttosto astratto.

• In Julia, Vector{T} è solo un alias per Array{T, 1}, e l'eltype T può essere qualsiasi cosa.
• In (gran parte della) fisica, un vettore è una freccia con una length e una direction nello spazio a N dimensioni, ma senza una posizione fissa.
• In algebra lineare, la freccia dei fisici è fissata nello spazio, con la coda nell'origine. Questo rende superfluo il fusto della freccia, quindi possiamo rappresentare il vettore come un point nello spazio a N dimensioni, dove gli N elementi rappresentano la distanza dall'origine lungo ciascuno degli N assi.

Ai fini che ci interessano, ignora i vettori di stringhe o di caratteri. In questo concetto lavoreremo con tipi numerici: Int, Float o Complex.

Che cos'è una matrice?

Anche qui troverai risposte diverse.

• Un array bidimensionale rettangolare di numeri (senza bordi irregolari). Una matrice quadrata è un caso speciale piuttosto comune.
• Una linear combination di vettori colonna, impilati fianco a fianco.
• Una transformation che può essere applicata a un vettore, proprio come una function viene applicata ad altri tipi di input.

Alcuni tipi di matrice quadrata sono abbastanza comuni da avere nomi speciali.

Matrice diagonale: Tutti gli elementi diversi da zero si trovano sulla main diagonal (dall'angolo in alto a sinistra a quello in basso a destra).

julia> [1 0 0; 0 2 0; 0 0 3]
3×3 Matrix{Int64}:
 1  0  0
 0  2  0
 0  0  3

# a convenient shortcut
julia> diagm(1:3)
3×3 Matrix{Int64}:
 1  0  0
 0  2  0
 0  0  3

Matrice identità: Una matrice diagonale con soli uno sulla diagonale (per motivi che diventeranno più chiari in seguito). Spesso abbreviata in I.

julia> Matrix{Float64}(I, 3, 3)
3×3 Matrix{Float64}:
 1.0  0.0  0.0
 0.0  1.0  0.0
 0.0  0.0  1.0

Triangolare superiore: Valori diversi da zero sulla diagonale e sopra di essa, zeri sotto. Triangolare inferiore: puoi immaginarlo.

Trasporre una matrice

Se vuoi davvero scambiare le righe con le colonne, la funzione permutedims() lo fa e ti restituisce una nuova matrice.

Questo comporta una copia, che è lenta e divoratrice di memoria quando si lavora con matrici grandi.

Ai fini dell'algebra lineare, la funzione transpose() è più utile, perché crea rapidamente un wrapper lazy attorno alla matrice originale.

Ancora più utile è la funzione adjoint(), che cambia anche il segno della parte immaginaria nei numeri complessi. Questa operazione è abbastanza comune da poter semplicemente aggiungere un apostrofo ' al nome della variabile per creare l'aggiunto.

julia> m
2×3 Matrix{Int64}:
 1  2  3
 4  5  6

# new, full copy
julia> permutedims(m)
3×2 Matrix{Int64}:
 1  4
 2  5
 3  6

# lazy version
julia> transpose(m)
3×2 transpose(::Matrix{Int64}) with eltype Int64:
 1  4
 2  5
 3  6

julia> mc = [1+2im 2+3im; 3+2im 1+2im]
2×2 Matrix{Complex{Int64}}:
 1+2im  2+3im
 3+2im  1+2im

# lazy conjugate transpose
julia> adjoint(mc)
2×2 adjoint(::Matrix{Complex{Int64}}) with eltype Complex{Int64}:
 1-2im  3-2im
 2-3im  1-2im

# syntactic sugar for adjoint
julia> mc'
2×2 adjoint(::Matrix{Complex{Int64}}) with eltype Complex{Int64}:
 1-2im  3-2im
 2-3im  1-2im

Moltiplicazione

Il resto di questo documento usa funzionalità del modulo LinearAlgebra. Una riga using LinearAlgebra lo porterà nello spazio dei nomi, ma non continueremo a ripeterla negli esempi (troppo disordine visivo).

Vettori, prodotto elemento per elemento

Ne abbiamo parlato nel concetto [Operazioni sui vettori][vector-ops]. L'operatore è .*, che opera sui vettori di input a coppie di elementi, producendo un output della stessa dimensione e dello stesso tipo degli input.

julia> [1, 2] .* [3, 4]
2-element Vector{Int64}:
 3
 8

Vettori, prodotto scalare

Questa operazione estremamente comune si scrive nei manuali come u ⋅ v ed è chiamata prodotto «dot».

Un prodotto scalare equivale alla somma del prodotto elemento per elemento.

Due vettori si possono moltiplicare con il solito operatore *, ma solo se il vettore a sinistra è un adjoint: comodamente scritto u' * v. Questo dovrebbe diventare più chiaro nella sezione successiva sulla moltiplicazione tra matrici.

Il punto sopraelevato (centrale) è disponibile in Julia (si inserisce digitando \cdot e poi tab) come zucchero sintattico per la funzione dot(). Non c'è bisogno di specificare l'aggiunto, perché questo dettaglio viene gestito automaticamente.

julia> using LinearAlgebra

# with element-wise syntax
julia> sum([1, 2] .* [3, 4])
11

# with u' * v syntax
julia> [1, 2]' * [3, 4]
11

# with dot()
julia> dot([1, 2], [3, 4])
11

# with \cdot syntax
julia> [1, 2] ⋅ [3, 4]
11

Per i vettori a valori complessi, il vettore a sinistra deve essere il coniugato (segno invertito sulla parte immaginaria). La funzione dot lo fa automaticamente, la sintassi u' lo fa esplicitamente, ma sum(u .* v) fallirà se u e v sono complessi.

Norme

Quanto è «grande» un vettore?

La norm è un tentativo di catturare questo, riducendo il vettore a uno scalare appropriato.

Esiste tutta una famiglia di norme, ma di gran lunga la più comune è la norma 2, che è √(v ⋅ v).

Questa operazione di media quadratica è la distanza pitagorica dall'origine (nello spazio a N dimensioni). Se visualizziamo il vettore come una freccia, con la coda nell'origine, la norma 2 è la lunghezza della freccia.

Qualsiasi norma p si può calcolare passando p come secondo argomento. La norma 1 a volte è utile: è semplicemente la somma dei valori assoluti degli elementi (quindi molto rapida e facile da calcolare).

# defaults to the 2-norm
julia> norm([1, 2, 3])
3.7416573867739413
# the 1-norm
julia> norm([1, -2, 3], 1)
6.0

Moltiplicazione tra matrici

A volte sembra che ogni matematico applicato del mondo abbia passato gran parte degli ultimi 80 anni a trasformare ogni calcolo in una serie di moltiplicazioni tra matrici.

È un'operazione in cui i computer sono molto bravi:

• È altamente ripetitiva.
• Si può parallelizzare in modo efficiente.
• È stato sviluppato molto hardware specializzato per renderla più veloce, compresa la GPU che probabilmente è integrata nel tuo portatile.

I dettagli sono piuttosto semplici, anche se non molto intuitivi a prima vista.

Considera una matrice A che moltiplica un vettore v per dare un output w (per convenzione dell'algebra lineare, useremo le lettere maiuscole per le matrici e le minuscole per i vettori).

La prima riga di A viene moltiplicata scalarmente per v per ottenere il primo elemento di w, la seconda riga dà il secondo elemento, e così via verso il basso.

julia> A = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4
julia> v = [5, 6]
2-element Vector{Int64}:
 5
 6
julia> A * v
2-element Vector{Int64}:
 17  # equals [1, 2] ⋅ [5, 6]
 39  # equals [3, 4] ⋅ [5, 6]

Una moltiplicazione matrice * matrice estende questo a tutte le colonne della matrice a destra.

Per C = A * B, possiamo considerare che A moltiplica ogni colonna di B per dare la colonna corrispondente di C: una serie di moltiplicazioni matrice-vettore.

In modo equivalente, potremmo dire che la prima riga di A viene moltiplicata scalarmente per ogni colonna di B per dare la prima riga di C, la seconda riga dà la seconda riga, e così via verso il basso.

julia> A
2×2 Matrix{Int64}:
 1  2
 3  4

julia> B = [5 7; 6 8]
2×2 Matrix{Int64}:
 5  7
 6  8

# left column is the same as A*v previously
julia> A * B
2×2 Matrix{Int64}:
 17  23
 39  53
 
julia> B * A
2×2 Matrix{Int64}:
 26  38
 30  44

Come mostra l'esempio qui sopra, la moltiplicazione tra matrici non è commutativa: non c'è una relazione semplice tra A*B e B*A.

Probabilmente è difficile visualizzare la moltiplicazione tra matrici solo leggendo le parole, per chi vi si avvicina per la prima volta. YouTube ha moltissimi video che la mostrano graficamente: cerca «matrix multiplication» e scegli quello con lo stile, il livello di dettaglio e la lingua che preferisci.

Dimensioni

Il prodotto scalare di due vettori richiede che abbiano la stessa lunghezza.

Per estensione, nella moltiplicazione tra matrici il numero di colonne della matrice a sinistra deve corrispondere al numero di righe della matrice a destra.

Esprimendo le dimensioni come tuple (nrows, ncols), così come le restituisce size(A), abbiamo (a, b) * (b, c) -> (a, c). Le dimensioni «interne», qui b e b, sono compatibili per calcolare i prodotti scalari. Le dimensioni «esterne», qui a e c, determinano le dimensioni dell'output.

Un esempio con matrici rettangolari:

julia> D = reshape(1:6, 2, 3)
2×3 reshape(::UnitRange{Int64}, 2, 3) with eltype Int64:
 1  3  5
 2  4  6

julia> E = reshape(1:12, 3, 4)
3×4 reshape(::UnitRange{Int64}, 3, 4) with eltype Int64:
 1  4  7  10
 2  5  8  11
 3  6  9  12

julia> D * E
2×4 Matrix{Int64}:
 22  49   76  103
 28  64  100  136

# (3, 4) * (2, 3) not possible
julia> E * D
ERROR: DimensionMismatch: matrix A has axes (Base.OneTo(3),Base.OneTo(4)), matrix B has axes (Base.OneTo(2),Base.OneTo(3))

Moltiplicare coppie di vettori, in stile matriciale, offre due possibilità.

Per convenzione, useremmo u' * v come equivalente del prodotto scalare. Le dimensioni sono (1, 3) * (3, 1), e Julia semplifica l'output (1, 1) a uno scalare (a differenza, per esempio, di R).

In alternativa, potremmo usare u * v', con dimensioni (3, 1) * (1, 3) -> (3, 3), per moltiplicare tutte le possibili coppie di elementi ed espandere i vettori in una matrice.

julia> u = [1, 2]
2-element Vector{Int64}:
 1
 2
julia> v = [3, 4]
2-element Vector{Int64}:
 3
 4
julia> u' * v
11
julia> u * v'
2×2 Matrix{Int64}:
 3  4
 6  8

u' * v è talvolta chiamato prodotto interno.

Rotazioni

Un tipo particolarmente comune di moltiplicazione tra matrici coinvolge le matrici di rotazione.

In 2D, esiste una matrice relativamente semplice per ruotare un vettore in senso antiorario di θ radianti.

julia> rot2d(θ, vec) = [cos(θ) -sin(θ); sin(θ) cos(θ)] * vec
rot2d (generic function with 1 method)

# unit vector in the x direction
julia> i_hat = [1, 0]
2-element Vector{Int64}:
 1
 0

# rotate 45 degrees
julia> rot2d(π/4, i_hat)
2-element Vector{Float64}:
 0.7071067811865476
 0.7071067811865475

# rotate 90 degrees -> unit vector in the y direction, j_hat
julia> rot2d(π/2, i_hat)
2-element Vector{Float64}:
 6.123233995736766e-17  # zero, within numerical error
 1.0