Laboratorio 04

Gestione della Memoria e Vettori C-Style

Cristiano De Michele

Il Laboratorio 03: cosa chiedeva

  • Task 1–2: relu(x) e relu_derivative(x), con un if/else
  • Task 3: sigmoid(x) e sigmoid_derivative(x), con #include <cmath>
  • Task 4: una tabella con un ciclo for da \(-5\) a \(+5\) (contatore int), colonne separate da "\t"; valori di controllo \(\sigma(0) = 0{,}5\) e \(\sigma'(0) = 0{,}25\)
  • Task 5: tanh_act(x) e tanh_derivative(x), non tanh: con <cmath> una funzione tanh richiamerebbe se stessa all’infinito
  • Task 6: derivate numeriche a differenze centrali; massimo errore su \(x \in [-5, 5]\) con \(h = 10^{-5}\), in std::scientific, e [PASS] sotto \(10^{-6}\); errore in \(x = 1\) al variare di \(h\)
  • Task 7 (facoltativo): la catena di due neuroni del Laboratorio 02 con sigmoide, tanh e ReLU, per \(y = 1\) e \(y = 0{,}8\); la ReLU con i parametri del Laboratorio 02; pesi iniziali \(w_0 = 0\) e \(w_0 = 5\)

Soluzione del Laboratorio 03, Task 1–3 e 5: le funzioni

main.cpp
#include <iostream>
#include <cmath>

// Soluzione del Laboratorio 03

// Task 1: ReLU
double relu(double x) {
    if (x > 0.0) {
        return x;
    } else {
        return 0.0;
    }
}

// Task 2: derivata della ReLU (0 in x = 0)
double relu_derivative(double x) {
    if (x > 0.0) {
        return 1.0;
    } else {
        return 0.0;
    }
}
main.cpp
// Task 3: sigmoide e sua derivata
double sigmoid(double x) {
    return 1.0 / (1.0 + std::exp(-x));
}

double sigmoid_derivative(double x) {
    double s = sigmoid(x);   // una sola volta
    return s * (1.0 - s);
}

// Task 5: tanh e sua derivata (nome tanh_act)
double tanh_act(double x) {
    return std::tanh(x);
}

double tanh_derivative(double x) {
    double t = std::tanh(x); // una sola volta
    return 1.0 - t * t;
}

Soluzione del Laboratorio 03, Task 5: perché tanh_act

#include <iostream>
#include <cmath>

double tanh(double x) {
1    return std::tanh(x);
}

int main() {
    std::cout << tanh(1.0) << std::endl;
    return 0;
}
1
<cmath> dichiara la tanh della libreria C nello spazio dei nomi globale, e std::tanh(double) è quella stessa funzione. Scrivendo double tanh(double) non si crea una funzione nuova: si definisce quella, e il linker usa questa definizione al posto di quella della libreria. std::tanh(x) chiama quindi questa funzione, che richiama se stessa all’infinito.

Il Task 5 chiedeva di chiamare la funzione tanh_act e non tanh. Con il nome tanh il programma a sinistra compila, ma in esecuzione termina con Segmentation fault (codice di uscita 139).

Con -Wall sia g++ sia clang++ lo segnalano con un avviso -Winfinite-recursion: è uno dei motivi per cui si compila sempre con -Wall.

Soluzione del Laboratorio 03, Task 4–6: la tabella

main.cpp
// Task 6: derivate numeriche (differenze centrali)
double sigmoid_numerical_derivative(double x, double h) {
    return (sigmoid(x + h) - sigmoid(x - h)) / (2.0 * h);
}

double tanh_numerical_derivative(double x, double h) {
    return (tanh_act(x + h) - tanh_act(x - h)) / (2.0 * h);
}

int main() {
    // Task 4 e 5: tabella da -5 a 5 a passi di 1, contatore int
    std::cout << "x\trelu\trelu'\tsigmoid\tsigmoid'\ttanh\ttanh'" << std::endl;
    for (int i = -5; i <= 5; ++i) {
        double x = i;
        std::cout << x << "\t" << relu(x) << "\t" << relu_derivative(x) << "\t"
                  << sigmoid(x) << "\t" << sigmoid_derivative(x) << "\t"
                  << tanh_act(x) << "\t" << tanh_derivative(x) << std::endl;
    }

Soluzione del Laboratorio 03, Task 6: la verifica

main.cpp
    // Task 6: massimo errore assoluto su x da -5 a 5 a passi di 0,5
    double h = 1e-5;
    double max_err_sigmoid = 0.0;
    double max_err_tanh = 0.0;
    for (int i = -10; i <= 10; ++i) {
        double x = 0.5 * i;
        double e1 = std::fabs(sigmoid_derivative(x) - sigmoid_numerical_derivative(x, h));
        double e2 = std::fabs(tanh_derivative(x) - tanh_numerical_derivative(x, h));
        if (e1 > max_err_sigmoid) max_err_sigmoid = e1;
        if (e2 > max_err_tanh) max_err_tanh = e2;
    }
    std::cout << std::scientific;
    std::cout << "\nerrore massimo sigmoide: " << max_err_sigmoid << std::endl;
    std::cout << "errore massimo tanh:     " << max_err_tanh << std::endl;
    if (max_err_sigmoid < 1e-6 && max_err_tanh < 1e-6) {
        std::cout << "[PASS]" << std::endl;
    } else {
        std::cout << "[FAIL]" << std::endl;
    }

Soluzione del Laboratorio 03, Task 6: l’errore al variare di h

main.cpp
    // Task 6, punto 4: errore assoluto e relativo in x = 1 per h = 1e-1, 1e-3, ..., 1e-13
    std::cout << "\nh\t\tass. sigmoide\trel. sigmoide\tass. tanh\trel. tanh" << std::endl;
    for (int k = 1; k <= 13; k += 2) {
        double hh = std::pow(10.0, -k);
        double e1 = std::fabs(sigmoid_derivative(1.0) - sigmoid_numerical_derivative(1.0, hh));
        double e2 = std::fabs(tanh_derivative(1.0) - tanh_numerical_derivative(1.0, hh));
        std::cout << hh << "\t" << e1 << "\t" << e1 / sigmoid_derivative(1.0) << "\t"
                  << e2 << "\t" << e2 / tanh_derivative(1.0) << std::endl;
    }
    return 0;
}

Nota

Le colonne separate da "\t" si disallineano quando un numero cambia lunghezza; per incolonnarle serve <iomanip>: std::setw(n) fissa la larghezza del campo successivo, std::setprecision(n) il numero di cifre: std::cout << std::setw(12) << sigmoid(x);.

Soluzione del Laboratorio 03, Task 6: i numeri

Errore relativo in \(x = 1\):

h \(10^{-1}\) \(10^{-3}\) \(10^{-5}\) \(10^{-7}\) \(10^{-9}\) \(10^{-11}\) \(10^{-13}\)
sigmoide \(3\cdot10^{-4}\) \(3\cdot10^{-8}\) \(9\cdot10^{-12}\) \(2\cdot10^{-9}\) \(8\cdot10^{-8}\) \(2\cdot10^{-5}\) \(2\cdot10^{-3}\)
tanh \(2\cdot10^{-3}\) \(2\cdot10^{-7}\) \(1\cdot10^{-11}\) \(2\cdot10^{-9}\) \(2\cdot10^{-7}\) \(2\cdot10^{-5}\) \(6\cdot10^{-4}\)

Con h = 1e-5 gli errori massimi sono \(8\cdot10^{-12}\) (sigmoide) e \(3\cdot10^{-11}\) (tanh): [PASS]. L’errore relativo conta invece le cifre corrette della derivata numerica: circa 11 con \(h = 10^{-5}\), 2 o 3 con \(h = 10^{-13}\). In \(x = 1\) l’errore assoluto è il relativo moltiplicato per \(\sigma'(1) \approx 0{,}20\) e \(\tanh'(1) \approx 0{,}42\): l’andamento con \(h\) è lo stesso.

Soluzione del Laboratorio 03, Task 6: lettura dei risultati

L’errore in \(x = 1\) scende fino a \(h \approx 10^{-5}\) e poi risale:

  • per \(h\) grande domina l’errore di troncamento delle differenze centrali, \(\propto h^2\);
  • per \(h\) piccolo domina la cancellazione numerica: per \(h = 10^{-13}\) f(x+h) e f(x-h) coincidono nelle prime 13 delle circa 16 cifre di un double, e la differenza ne conserva solo 2 o 3.

Il criterio del [PASS] è assoluto. Vicino a \(|x| = 5\) le derivate sono piccole, per esempio \(\tanh'(5) \approx 2\cdot10^{-4}\), e lì una soglia assoluta di \(10^{-6}\) ammetterebbe un errore relativo dello \(0{,}5\%\). Per questo nel gradient check si usa di solito l’errore relativo.

Soluzione del Laboratorio 03, Task 7: attivazioni e derivate

task7_attivazioni.cpp
#include <iostream>
#include <cmath>

// Task 7 (facoltativo): la catena del Lab 02
// tipo: 0 sigmoide, 1 tanh, 2 ReLU

double attiva(int tipo, double z) {
    if (tipo == 0) {
        return 1.0 / (1.0 + std::exp(-z));
    } else if (tipo == 1) {
        return std::tanh(z);
    } else if (z > 0.0) {
        return z;
    } else {
        return 0.0;
    }
}
task7_attivazioni.cpp
double deriva(int tipo, double z) {
    if (tipo == 0) {
        double s = attiva(0, z);
        return s * (1.0 - s);
    } else if (tipo == 1) {
        double t = std::tanh(z);
        return 1.0 - t * t;
    } else if (z > 0.0) {
        return 1.0;
    } else {
        return 0.0;     // anche in z = 0
    }
}

Soluzione del Laboratorio 03, Task 7: l’addestramento

task7_attivazioni.cpp
void addestra(int tipo, double y, double w0, double eta) {
    double x = 1.0;
    double w1 = w0, w2 = w0;
    for (int t = 0; t <= 1000; ++t) {
        double z1 = w1 * x;
        double a1 = attiva(tipo, z1);
        double z2 = w2 * a1;
        double a2 = attiva(tipo, z2);
        double E  = (a2 - y) * (a2 - y);
        double delta2 = 2.0 * (a2 - y) * deriva(tipo, z2);
        double delta1 = w2 * delta2 * deriva(tipo, z1);
        if (t % 500 == 0) {
            std::cout << "t = " << t << "\tE = " << E
                      << "\tw1 = " << w1 << "\tw2 = " << w2 << std::endl;
        }
        w1 = w1 - eta * delta1 * x;
        w2 = w2 - eta * delta2 * a1;
    }
}

Soluzione del Laboratorio 03, Task 7: il main

task7_attivazioni.cpp
int main() {
    // punto 3: w1 = w2 = 0,5 ed eta = 0,25
    for (int tipo = 0; tipo < 3; ++tipo) {
        for (int caso = 0; caso < 2; ++caso) {
            double y = 1.0;
            if (caso == 1) {
                y = 0.8;
            }
            std::cout << "\n[tipo " << tipo << "] y = " << y << std::endl;
            addestra(tipo, y, 0.5, 0.25);
        }
    }

Soluzione del Laboratorio 03, Task 7: il main (segue)

task7_attivazioni.cpp
    // punto 4: la ReLU con i parametri del Lab 02, w1 = w2 = 1 ed eta = 1
    std::cout << "\n[ReLU, w0 = 1, eta = 1] y = 1" << std::endl;
    addestra(2, 1.0, 1.0, 1.0);
    std::cout << "\n[ReLU, w0 = 1, eta = 1] y = 0.8" << std::endl;
    addestra(2, 0.8, 1.0, 1.0);
    // punto 5: y = 0,8 con w1 = w2 = 0 e con w1 = w2 = 5
    for (int k = 0; k < 2; ++k) {
        double w0 = 5.0 * k;
        for (int tipo = 0; tipo < 3; ++tipo) {
            std::cout << "\n[tipo " << tipo << ", w0 = " << w0 << "] y = 0.8" << std::endl;
            addestra(tipo, 0.8, w0, 0.25);
        }
    }
    return 0;
}

Soluzione del Laboratorio 03, Task 7: i numeri

dopo 1000 passi sigmoide tanh ReLU
\(y = 1\) \(E = 1{,}3\cdot10^{-3}\), pesi in crescita \(E = 2{,}7\cdot10^{-4}\), pesi in crescita \(E = 0\), \(w_1 = w_2 = 1\)
\(y = 0{,}8\) \(E = 4\cdot10^{-10}\) \(E \approx 10^{-31}\) \(E \approx 10^{-32}\), \(w_1 = w_2 = \sqrt{0{,}8}\)
  • Sigmoide e tanh restano sotto \(1\): per \(y = 1\) servirebbero pesi infiniti, e come nel Laboratorio 02 l’errore cala sempre più lentamente mentre i pesi crescono. La ReLU non è limitata: con pesi positivi \(a_2 = w_1 w_2\,x\), e \(w_1 w_2 = 1\) si ottiene con pesi finiti. Con \(y = 0{,}8\) convergono tutte e tre; per tanh e ReLU \(E\) si ferma a \(10^{-31}\)–\(10^{-32}\), cioè \(|a_2 - y| \approx 10^{-16}\): il limite della precisione del double.

Soluzione del Laboratorio 03, Task 7: la ReLU con \(w_0 = 1\) ed \(\eta = 1\)

  • Con \(y = 1\): per \(w_1 = w_2 = 1\) si ha già \(a_2 = 1\), quindi \(E = 0\) e il gradiente è nullo.
  • Con \(y = 0{,}8\): il passo \(\eta = 1\) è troppo grande e i pesi oscillano senza convergere; al passo 20 diventano negativi (\(w_1 = w_2 \approx -0{,}004\)). Da lì \(z_1 = w_1 x < 0\), quindi \(a_1 = 0\) e \(\mathrm{ReLU}'(z_1) = 0\): \(\partial E/\partial w_1 = \delta_1 x = 0\) e \(\partial E/\partial w_2 = \delta_2\, a_1 = 0\). I pesi non cambiano più: \(a_2 = 0\) ed \(E = y^2 = 0{,}64\). Un neurone ReLU che dà \(0\) su tutti gli ingressi di addestramento (qui il solo \(x = 1\)) non riceve più gradiente: si parla di ReLU morta.

Soluzione del Laboratorio 03, Task 7: i pesi iniziali

\(y = 0{,}8\), 1000 passi sigmoide tanh ReLU
\(w_0 = 0\) \(E = 2\cdot10^{-9}\) \(E = 0{,}64\), pesi fermi \(E = 0{,}64\), pesi fermi
\(w_0 = 5\) \(E\) da \(0{,}037\) a \(0{,}033\) \(E = 0{,}040\), quasi ferma \(E = 0{,}64\), ReLU morta al primo passo
  • \(w_0 = 0\): con tanh e ReLU \(a_1 = 0\) e \(w_2 = 0\), quindi \(\partial E/\partial w_2 = \delta_2\, a_1 = 0\) e \(\partial E/\partial w_1 = \delta_1 x = 0\) (in \(\delta_1\) c’è il fattore \(w_2\)): è un punto stazionario. Con la sigmoide \(a_1 = \sigma(0) = 0{,}5\) e l’addestramento parte.

Soluzione del Laboratorio 03, Task 7: i pesi iniziali (segue)

  • \(w_0 = 5\): sigmoide e tanh sono in saturazione, \(\sigma'(5) \approx 7\cdot10^{-3}\) e \(\tanh'(5) \approx 2\cdot10^{-4}\), quindi il gradiente è quasi nullo. Con la ReLU \(a_2 = 25\): il primo passo porta i pesi a \(-55{,}5\).
  • I pesi iniziali non vanno scelti né nulli né grandi: le inizializzazioni He e Xavier del Laboratorio 08 ne fissano la scala in base al numero di ingressi e di uscite dello strato. La velocità di convergenza della sigmoide dipende invece soprattutto da \(\eta\) e da \(\sigma' \le 1/4\): con \(y = 0{,}8\) converge per ogni \(w_0\) fra \(-0{,}5\) e \(3\).

Obiettivo di oggi

Prima di chiudere i dati in una classe (Laboratorio 05) si gestisce la memoria a mano, con new[] e delete[] (Lezione 08), come in C con malloc e free. Si parte dai pesi di un singolo neurone, in array dinamici; poi si passa a uno strato di neuroni, che richiede una matrice, scritta come in C: una struct e funzioni separate. Il Laboratorio 05 riprenderà quel codice e lo trasformerà in una classe.

Progetto e AddressSanitizer

Come nel Laboratorio 03, si lavora nel progetto basic_neural_plusplus del Laboratorio 01. Copiate il main.cpp del Laboratorio 03 in lab03_main.cpp, che non va compilato, poi svuotate main.cpp: il codice di oggi va lì, nell’ordine dei task, con #include <iostream> in cima e le funzioni prima del main.

Nel CMakeLists.txt aggiungete, prima di add_executable, le due righe che attivano l’AddressSanitizer (Lezione 08):

add_compile_options(-fsanitize=address -g)
add_link_options(-fsanitize=address)

Poi, dalla cartella del progetto, a ogni modifica: cmake --build build && ./build/basic_neural (se manca la cartella build, prima cmake -B build). Il programma chiede \(N\) da tastiera.

Su macOS l’AddressSanitizer non segnala i leak: per i Task 4 e 6 seguite la slide “L’AddressSanitizer su macOS” della Lezione 08.

Task 1: allocazione dinamica

Gli array a dimensione fissa (Lezione 05) richiedono \(N\) noto in compilazione; qui \(N\) si legge con std::cin, in esecuzione.

  1. Nel main dichiarate int N;, stampate Inserisci il numero di input del neurone (es. 784): e leggete N con std::cin. 784 sono i pixel di un’immagine MNIST.
  2. Se N <= 0 stampate su std::cerr il messaggio Errore: N deve essere maggiore di 0! e terminate con return 1;. Il controllo copre anche una lettura fallita: in quel caso N vale 0 (Lezione 04).
  3. Allocate nell’heap due array di double di dimensione \(N\): double* pesi = new double[N]; e double* input = new double[N];.

Task 2: inizializzazione

Con un ciclo for riempite i due array con valori di prova:

  • gli input con \(1{,}0\);
  • i pesi con \(0{,}01\).

Come per gli array statici, all’elemento \(i\) si accede con pesi[i].

Task 3: il forward pass a mano

Prima della funzione di attivazione, un neurone calcola un prodotto scalare più il bias: \(z = \sum_{i=0}^{N-1} w_i x_i + b\).

  1. Dichiarate double b = 0.5; e double z = 0.0;.
  2. Con un ciclo for accumulate in z i prodotti pesi[i] * input[i].
  3. Aggiungete b a z e stampate Risultato del calcolo (Z): <z>. Con \(N = 784\) deve risultare \(784 \cdot 0{,}01 + 0{,}5 = 8{,}34\).

Task 4: liberare la memoria

Regola: a ogni new[] corrisponde un delete[].

  1. Dopo il Task 3 i due array non servono più: liberateli con delete[] input; e delete[] pesi;.
  2. Eseguite con \(N = 784\): l’AddressSanitizer non deve stampare nulla. Poi commentate uno dei due delete[], ricompilate ed eseguite di nuovo: deve segnalare 6272 byte(s) leaked in 1 allocation(s). Perché proprio 6272? Rimettete il delete[].
  3. Cosa succederebbe se il delete[] mancasse e questo codice si ripetesse a ogni passo di un addestramento?

Task 5: una matrice come in C

Uno strato di \(M\) neuroni richiede una matrice \(M \times N\). In C è fatta di tre informazioni da tenere insieme. Prima del main scrivete la struct:

struct CMatrix {
    int rows;
    int cols;
    double* data;   // rows*cols double consecutivi, in ordine row-major
};

Task 5 (segue): le quattro funzioni

Sempre prima del main, dopo la struct, scrivete quattro funzioni:

  1. CMatrix cmat_alloc(int rows, int cols): dichiara un CMatrix m, ne imposta rows e cols, alloca m.data = new double[rows * cols] e restituisce m.
  2. void cmat_free(CMatrix& m): delete[] m.data, poi m.data = nullptr e m.rows, m.cols a zero; una seconda chiamata è innocua (delete[] su nullptr non fa nulla).
  3. double cmat_get(const CMatrix& m, int i, int j) e void cmat_set(CMatrix& m, int i, int j, double val): leggono e scrivono l’elemento \((i, j)\), che sta in m.data[i * m.cols + j] (ordine row-major).

Task 5 (segue): row-major e column-major

La convenzione row-major è quella della nostra CMatrix, ma non è l’unica.

Nota

Eigen, per default, usa l’ordine opposto. Una Eigen::MatrixXd è column-major: in memoria c’è prima tutta la colonna 0, poi la colonna 1, e così via. Una matrice \(2\times3\) con le righe 1 2 3 e 4 5 6 è in memoria 1 4 2 5 3 6, non 1 2 3 4 5 6.

Con M(i, j) non cambia niente: l’indicizzazione è la stessa. La differenza conta quando si passa il buffer .data() a un’altra libreria, per esempio OpenCV.

Task 6: da un neurone a uno strato

Con \(W\) di dimensione \(M \times N\), il neurone \(j\) calcola, prima dell’attivazione, \(z_j = \sum_{i} W_{ji}\, x_i + b_j\): per tutto lo strato è un prodotto matrice-vettore più il vettore dei bias. Come nel Task 5 si scrive come in C, senza classi: la CMatrix e una funzione libera.

  1. Prima del main scrivete void layer_forward(const CMatrix& W, const double* x, const double* b, double* z). Per ogni neurone j da 0 a W.rows - 1, accumulate in un double acc inizializzato a 0 i prodotti cmat_get(W, j, i) * x[i] per i da 0 a W.cols - 1, poi scrivete z[j] = acc + b[j]. Usate cmat_get, non l’aritmetica degli indici.
  2. Nel main, dopo il Task 4: int M = 3; e CMatrix W = cmat_alloc(M, N);, con lo stesso \(N\) del Task 1. Allocate con new[] double* bias e double* z_out di \(M\) elementi e double* x di \(N\) elementi.

Task 6 (segue)

  1. Riempite x con \(1{,}0\), bias con \(0{,}5\) e, con cmat_set in due cicli annidati (\(j\) da 0 a \(M-1\), \(i\) da 0 a \(N-1\)), ogni elemento \((j, i)\) di W con \(0{,}01\,(j+1)\).
  2. Chiamate layer_forward(W, x, bias, z_out); e stampate una riga per neurone: z_out[<j>] = <valore>. Valori attesi: \(z_j = 0{,}01\,(j+1)\,N + 0{,}5\), cioè \(8{,}34\), \(16{,}18\) e \(24{,}02\) per \(N = 784\).
  3. Liberate tutto: cmat_free(W); e delete[] di bias, x e z_out. L’AddressSanitizer non deve stampare nulla.

Task 7: due domande

Nella parte del Task 6 ci sono quattro allocazioni (una dentro cmat_alloc, tre new[]), quattro deallocazioni e tre informazioni da passare a ogni funzione. Rispondete a due domande prima del Laboratorio 05:

  1. Che cosa succede con CMatrix W2 = W; seguito da cmat_free(W); cmat_free(W2);? Provatelo: scrivete le tre righe al posto di cmat_free(W);, eseguite e leggete il messaggio dell’AddressSanitizer; poi ripristinate il codice. Che cosa copia l’assegnamento di una struct che contiene un puntatore?
  2. Chi garantisce che x abbia davvero W.cols elementi quando lo si passa a layer_forward?

Le risposte motivano il passaggio alle classi: il Laboratorio 05 prende CMatrix e le sue quattro funzioni e le trasforma in class Matrix. I tre dati diventano private, cmat_alloc diventa il costruttore, cmat_get e cmat_set diventano metodi.

Task da completare

La soluzione verrà discussa all’inizio del prossimo laboratorio.