Files
pinn/README.md
T

554 lines
22 KiB
Markdown
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
# Heat Equation PINN
![Python](https://img.shields.io/badge/Python-3.10%2B-blue)
![PyTorch](https://img.shields.io/badge/PyTorch-2.0%2B-orange)
![Plotly](https://img.shields.io/badge/Plotly-5.15%2B-green)
![Tests](https://img.shields.io/badge/tests-42%20passing-brightgreen)
Soluzione dell'equazione del calore 1D con sorgente puntuale tramite **Physics-Informed Neural Network (PINN)**, validata contro un solver numerico **FDM** (Finite Difference Method).
---
## Indice
1. [Problema fisico](#1-problema-fisico)
2. [Approccio: PINN vs FDM](#2-approccio-pinn-vs-fdm)
3. [Struttura del progetto](#3-struttura-del-progetto)
4. [Installazione](#4-installazione)
5. [Utilizzo](#5-utilizzo)
6. [Architettura della rete neurale](#6-architettura-della-rete-neurale)
7. [Funzione di loss](#7-funzione-di-loss)
8. [Training](#8-training)
9. [Solver FDM di riferimento](#9-solver-fdm-di-riferimento)
10. [Visualizzazioni](#10-visualizzazioni)
11. [Parametri configurabili](#11-parametri-configurabili)
12. [Test](#12-test)
---
## 1. Problema fisico
Il progetto risolve l'equazione del calore 1D con una sorgente di calore puntuale interna che si attiva a un istante fissato:
```
∂T/∂t = α ∂²T/∂x² + (α/k) Q(t) δ(x X_SRC)
```
dove:
- `T(x, t)` — temperatura [°C]
- `α` — diffusività termica [m²/s]
- `k` — conducibilità termica [W/m·K]
- `Q(t) = Q_VAL` se `t ≥ T_STEP`, altrimenti `0` — flusso di calore [W/m²]
- `δ(x X_SRC)` — delta di Dirac nella posizione della sorgente
- Dominio: `x ∈ [0, L]`, `t ∈ [0, T_END]`
### Condizioni al contorno (Robin / convezione)
Entrambe le estremità della barra scambiano calore per convezione con l'ambiente:
```
x = 0: k ∂T/∂x = h (T T_AMB) → ∂T/∂x (h/k)(T T_AMB) = 0
x = L: k ∂T/∂x = h (T T_AMB) → ∂T/∂x + (h/k)(T T_AMB) = 0
```
### Condizione iniziale
```
T(x, 0) = T₀ (temperatura uniforme)
```
### Parametri fisici di default
| Parametro | Valore | Unità | Significato |
|--------------|---------|----------|------------------------------------------|
| `ALPHA` | 0.01 | m²/s | Diffusività termica |
| `K` | 1.0 | W/m·K | Conducibilità termica |
| `L` | 1.0 | m | Lunghezza della barra |
| `T0` | 20.0 | °C | Temperatura iniziale uniforme |
| `X_SRC` | 0.35 | m | Posizione della sorgente di calore |
| `Q_VAL` | 150.0 | W/m² | Intensità del flusso di calore |
| `T_STEP` | 0.2 | s | Istante di attivazione della sorgente |
| `H_CONV` | 10.0 | W/m²·K | Coefficiente convettivo alle estremità |
| `T_AMB` | 20.0 | °C | Temperatura ambiente |
| `T_END` | 10.0 | s | Fine della simulazione |
---
## 2. Approccio: PINN vs FDM
### PINN (Physics-Informed Neural Network)
Una PINN è una rete neurale che impara a soddisfare un'equazione differenziale alle derivate parziali **senza dati sperimentali**. L'addestramento avviene minimizzando una funzione di loss che penalizza le violazioni della PDE, delle condizioni iniziali e delle condizioni al contorno in un insieme di **punti di collocazione** distribuiti nel dominio.
Vantaggi:
- Non richiede una griglia regolare
- Può essere addestrata su domini irregolari
- La soluzione è una funzione continua e differenziabile ovunque
In questo progetto la loss ha tre componenti:
1. **Residuo PDE** — la rete deve soddisfare l'equazione del calore
2. **Condizione iniziale** — la rete deve restituire `T₀` per `t = 0`
3. **Condizioni al contorno** — la rete deve soddisfare le Robin BC a `x=0` e `x=L`
### FDM (Finite Difference Method)
Il solver FDM (`fdm/`) implementa lo schema **FTCS** (Forward-Time Centered-Space) esplicito su una griglia regolare `NX × NT`. Non usa reti neurali: è una soluzione numerica classica che serve come **riferimento ad alta fedeltà** per validare la PINN.
La valutazione della PINN consiste nel calcolare l'errore relativo L2 tra la predizione della rete e la soluzione FDM interpolata sulla stessa griglia `100 × 100`.
---
## 3. Struttura del progetto
```
pinn/
├── config.py ← tutti i parametri fisici e numerici (modificare qui)
├── model.py ← HeatPINN (rete neurale) + heat_pinn_loss()
├── engine.py ← campionamento dati, training, valutazione vs FDM
├── visualizer.py ← grafici PINN vs FDM: heatmap, animazione, serie temporali
├── app.py ← CLI interattiva per PINN
├── requirements.txt
├── pytest.ini
├── clear.sh ← script di pulizia artefatti
├── fdm/
│ ├── solver.py ← schema FTCS esplicito con Robin BC e sorgente puntuale
│ ├── visualizer.py ← grafici FDM standalone
│ └── app.py ← CLI interattiva per FDM
└── tests/
├── conftest.py
├── test_config.py
├── test_model.py
├── test_engine_data.py
├── test_fdm_solver.py
├── test_integration_pinn.py
└── test_e2e.py
```
| File | Responsabilità |
|--------------------|-------------------------------------------------------------------------|
| `config.py` | Unica fonte di verità per tutti i parametri |
| `model.py` | Definizione della rete e calcolo della loss fisica |
| `engine.py` | Pipeline completa: campionamento → training → valutazione |
| `visualizer.py` | Plot interattivi HTML (PINN vs FDM) |
| `app.py` | Menu CLI per l'utente |
| `fdm/solver.py` | Solver numerico FTCS per la soluzione di riferimento |
| `fdm/visualizer.py`| Plot interattivi HTML (FDM standalone) |
| `fdm/app.py` | Menu CLI per il solver FDM |
---
## 4. Installazione
**Prerequisiti:** Python 3.10+, `pip`, `virtualenv` (o `venv`).
```bash
git clone <repository-url>
cd pinn
python -m venv .venv
source .venv/bin/activate # Linux / macOS
# .venv\Scripts\activate # Windows
pip install -r requirements.txt
```
Il progetto rileva automaticamente GPU/MPS/CPU all'avvio (vedi [engine.py — device detection](#12-dettagli-implementativi)).
---
## 5. Utilizzo
Attivare sempre il virtual environment prima di eseguire qualsiasi script:
```bash
source .venv/bin/activate
```
### PINN (`app.py`)
```bash
python app.py
```
```
1. Addestra nuovo modello
2. Valuta vs FDM (L2 error, max error)
3. Visualizza (genera 3 file HTML)
0. Esci
```
- **Opzione 1** — avvia l'addestramento (Adam + L-BFGS). Chiede il numero di epoche; premi Invio per usare il default (5000). Il modello migliore viene salvato in `models/best_heat_pinn_model.pth`. Per riaddestrare da zero: `rm models/best_heat_pinn_model.pth`.
- **Opzione 2** — carica il modello salvato, esegue il solver FDM, stampa l'errore relativo L2 e l'errore massimo assoluto.
- **Opzione 3** — genera tre file HTML interattivi in `results/pinn/<timestamp>/`.
### FDM (`fdm/app.py`)
```bash
python fdm/app.py
```
```
1. Risolvi (schema FTCS)
2. Visualizza (genera 3 file HTML)
0. Esci
```
- **Opzione 1** — esegue il solver e stampa shape della matrice e range di temperatura.
- **Opzione 2** — genera tre file HTML in `results/fdm/<timestamp>/`.
### Test
```bash
pytest tests/ -v # tutti i test (42)
pytest tests/test_model.py -v # rete e loss
pytest tests/test_engine_data.py -v # campionamento dati
pytest tests/test_fdm_solver.py -v # solver FDM
pytest tests/test_integration_pinn.py -v # integrazione PINN
pytest tests/test_e2e.py -v # workflow completo
```
### Pulizia artefatti
```bash
./clear.sh # menu interattivo per eliminare models/, results/ o entrambi
```
---
## 6. Architettura della rete neurale
`HeatPINN` ([model.py](model.py)) è una rete fully connected a 5 layer:
```
Input (x, t)
[Linear 2→128, Tanh]
[Linear 128→128, Tanh] ×3
[Linear 128→1]
Output: T
```
La rete riceve le coordinate `(x, t)` e produce un unico scalare: la temperatura `T(x, t)`.
### Normalizzazione dell'input
Prima di entrare nella rete, le coordinate vengono normalizzate al range `[0, 1]`:
```python
x_norm = x / L # x ∈ [0, 1]
t_norm = t / T_END # t ∈ [0, 1]
```
Questo migliora il condizionamento numerico dell'ottimizzazione.
### Output scaling
La rete interna `net` non predice direttamente la temperatura: predice una **perturbazione adimensionale**. La temperatura fisica viene ricostruita con:
```python
T = T_AMB + (Q_VAL * L / K) * net(x_norm, t_norm)
```
dove `T_char = Q_VAL * L / K ≈ 150 °C` è la scala caratteristica di temperatura del problema.
**Perché questo scaling?**
- L'output della rete rimane nell'ordine di `[0, 1]`, rendendo il training più stabile.
- I gradienti `∂T/∂x` risultanti sono `O(1)` — la rete può imparare la struttura spaziale senza problemi di scala.
- Il termine di fondo `T_AMB` garantisce che la soluzione parta dalla condizione iniziale corretta anche con pesi random.
**Non rimuovere questo scaling**: senza di esso la rete deve imparare ordini di grandezza diversi tra le condizioni iniziali/al contorno e il gradiente interno, rendendo l'ottimizzazione molto più difficile.
---
## 7. Funzione di loss
`heat_pinn_loss()` ([model.py](model.py)) calcola quattro valori: `(total, L_pde, L_ic, L_bc)`.
```
total = W_PDE * L_pde + W_IC * L_ic + W_BC * L_bc
```
Ogni termine è **normalizzato automaticamente** da scale precompilate che dipendono solo dalle costanti in `config.py`. Cambiare `Q_VAL`, `K`, `H_CONV` o `L` ribilancia automaticamente la loss senza richiedere rituning manuale dei pesi.
### L_pde — Residuo PDE
Valutato su `N_F` punti di collocazione `(x_f, t_f)` distribuiti nel dominio:
```
residuo = dT/dt α d²T/dx² source(x, t)
```
Il termine sorgente usa un'approssimazione gaussiana al delta di Dirac (il delta non è differenziabile):
```
source(x, t) = (α/k) · Q(t) · G(x)
G(x) = exp(0.5 · ((x X_SRC) / σ)²) / (σ √(2π)) σ = GAUSS_SIGMA = 0.01 m
```
`Q(t)` è la funzione gradino: vale `Q_VAL` per `t ≥ T_STEP`, zero altrimenti.
I gradienti `dT/dt`, `dT/dx`, `d²T/dx²` sono calcolati con `torch.autograd.grad`**autodifferenziazione esatta**, non differenze finite.
Normalizzazione: `_pde_scale = max((T_char/T_END)², src_peak²)` dove `src_peak` è il picco gaussiano.
### L_ic — Condizione iniziale
Valutato su `N_IC` punti `(x_ic, 0)`:
```
L_ic = mean( (T(x, 0) T₀)² ) / T_char²
```
### L_bc — Condizioni al contorno Robin
Valutato su `N_BC` istanti temporali `t_bc`, applicato a entrambe le estremità:
```
x = 0: ∂T/∂x(0, t) (h/k)(T(0,t) T_AMB) = 0
x = L: ∂T/∂x(L, t) + (h/k)(T(L,t) T_AMB) = 0
```
Normalizzazione: `_bc_scale = max(Q_VAL/K, H_CONV·T_char/K)²`
> **Nota sul segno:** la BC a sinistra (x=0) ha segno negativo davanti al termine convettivo perché il flusso uscente è orientato verso `x`; a destra (x=L) il segno è positivo perché il flusso uscente è orientato verso `+x`.
---
## 8. Training
Il training è implementato in `train_model()` ([engine.py](engine.py)) e procede in due fasi.
### Fase 1: Adam
```
Ottimizzatore: Adam, LR = 1e-3
Scheduler: ReduceLROnPlateau (factor=0.5, patience=150, min_lr=1e-6)
Early stopping: se la loss non migliora di > 1e-7 per 500 epoche consecutive
```
Il modello con la loss più bassa viene salvato a ogni miglioramento in `models/best_heat_pinn_model.pth`.
### Fase 2: L-BFGS (fine-tuning)
Al termine dell'Adam, viene caricato il miglior modello e affinato con L-BFGS:
```
Ottimizzatore: L-BFGS, LR=0.1, max_iter=50, history_size=50, strong Wolfe
Steps: 200
```
L-BFGS è un ottimizzatore di secondo ordine (quasi-Newton) particolarmente efficace nella fase finale del training PINN perché sfrutta la curvatura della loss per convergere a minimi più precisi di quanto Adam riesca a fare.
**Meccanismo closure:** L-BFGS richiede di poter rivalutare la loss più volte per iterazione. La funzione `closure()` cattura i componenti della loss in un dizionario `_last` per poterli stampare senza ricalcolare il grafo computazionale fuori dal contesto di `backward()`.
### Campionamento dei punti di collocazione
`prepare_data()` genera i punti di collocazione con **clustering deliberato** nelle zone fisicamente più complesse:
| Zona | Proporzione | Motivazione |
|-----------------------------|-------------|------------------------------------------------|
| Uniforme `[0,L] × [0,T_END]`| 50% | Copertura generale del dominio |
| Intorno a `X_SRC ± 5% L` | 25% | Gradiente ripido in prossimità della sorgente |
| Intorno a `T_STEP ± 0.1 s` | 25% | Discontinuità temporale all'attivazione |
Il clustering aumenta la densità di punti dove la fisica è più difficile da apprendere, senza aumentare il costo computazionale totale.
---
## 9. Solver FDM di riferimento
`fdm/solver.py` implementa lo schema **FTCS** (Forward-Time Centered-Space) esplicito.
### Schema di avanzamento temporale
```
T[i, n+1] = T[i, n] + r · (T[i+1,n] 2·T[i,n] + T[i1,n])
r = α·dt/dx² (numero di Courant-Friedrichs-Lewy)
```
### Stabilità CFL
Condizione necessaria per la stabilità dello schema esplicito:
```
r = α·dt/dx² ≤ 0.5
```
Se la condizione è violata, il solver stampa un avvertimento ma non si blocca. Con i parametri di default (`NX=250`, `NT=15000`) la condizione è soddisfatta. Se si riducono `NX` o `NT`, verificare che `r ≤ 0.5`.
### Condizioni al contorno Robin
Le BC sono applicate a ogni passo temporale usando uno schema centrato:
```
T[0] = (T[1] + robin_coeff · T_AMB) / (1 + robin_coeff) # x = 0
T[-1] = (T[-2] + robin_coeff · T_AMB) / (1 + robin_coeff) # x = L
robin_coeff = dx · h / k
```
### Iniezione della sorgente puntuale
Dopo l'applicazione delle BC, la sorgente viene iniettata al nodo più vicino a `X_SRC`:
```
T[i_src] += Q(t) · α · dt / (k · dx)
```
---
## 10. Visualizzazioni
Tutti i plot sono **HTML interattivi** generati con Plotly (zoom, hover, slider, play/pause).
### PINN vs FDM (`visualizer.py`)
Generati con `python app.py → opzione 3`, salvati in `results/pinn/<timestamp>/`:
| File | Contenuto |
|-------------------|------------------------------------------------------------------------|
| `heatmap.html` | Heatmap 2D affiancate PINN vs FDM, stessa scala colori |
| `animation.html` | Profilo T(x) animato nel tempo: PINN (blu continuo) vs FDM (rosso tratteggiato) |
| `comparison.html` | Serie temporali T(t) nei punti fissi `x=0`, `x=L/2`, `x=L`; linea verticale a `t=T_STEP` |
### FDM standalone (`fdm/visualizer.py`)
Generati con `python fdm/app.py → opzione 2`, salvati in `results/fdm/<timestamp>/`:
| File | Contenuto |
|-------------------|------------------------------------------------------------------------|
| `heatmap.html` | Heatmap 2D T(x,t) + striscia animata del profilo di temperatura |
| `animation.html` | Profilo T(x) animato nel tempo |
| `time_series.html`| Serie temporali in 5 punti fissi: `0`, `0.25L`, `0.5L`, `0.75L`, `L` |
---
## 11. Parametri configurabili
Tutti i parametri si trovano in `config.py`. Modificare solo questo file per cambiare il problema o il training.
### Fisica
| Parametro | Default | Unità | Descrizione |
|--------------|---------|--------|------------------------------------------|
| `ALPHA` | 0.01 | m²/s | Diffusività termica |
| `K` | 1.0 | W/m·K | Conducibilità termica |
| `L` | 1.0 | m | Lunghezza della barra |
| `T0` | 20.0 | °C | Temperatura iniziale uniforme |
| `X_SRC` | 0.35 | m | Posizione della sorgente di calore |
| `Q_VAL` | 150.0 | W/m² | Intensità del flusso di calore |
| `T_STEP` | 0.2 | s | Istante di attivazione della sorgente |
| `H_CONV` | 10.0 | W/m²·K | Coefficiente convettivo alle estremità |
| `T_AMB` | 20.0 | °C | Temperatura ambiente |
| `T_END` | 10.0 | s | Fine della simulazione |
### Griglia FDM
| Parametro | Default | Descrizione |
|--------------|---------|---------------------------------------------------------|
| `NX` | 250 | Nodi spaziali (aumentare per maggiore risoluzione) |
| `NT` | 15000 | Passi temporali (verificare `r = α·dt/dx² ≤ 0.5`) |
| `GAUSS_SIGMA`| 0.01 | Larghezza del picco gaussiano nella loss PINN [m] |
### Architettura PINN
| Parametro | Default | Descrizione |
|-------------------|---------|----------------------------------------------|
| `HIDDEN_SIZE` | 128 | Neuroni per layer nascosto |
| `N_HIDDEN_LAYERS` | 4 | Numero di layer nascosti (totale: 5 layer) |
### Campionamento
| Parametro | Default | Descrizione |
|-----------|---------|--------------------------------------------------------------------|
| `N_F` | 6000 | Punti PDE (+ 50% clustering automatico vicino a X_SRC e T_STEP) |
| `N_IC` | 400 | Punti condizione iniziale |
| `N_BC` | 400 | Punti condizioni al contorno |
### Training Adam
| Parametro | Default | Descrizione |
|-----------------|---------|-----------------------------------------------------|
| `EPOCHS` | 5000 | Epoche massime |
| `PATIENCE` | 500 | Early stopping: epoche senza miglioramento |
| `LR_ADAM` | 1e-3 | Learning rate iniziale |
| `SCHED_FACTOR` | 0.5 | Fattore di riduzione LR (ReduceLROnPlateau) |
| `SCHED_PATIENCE`| 150 | Patience per la riduzione LR |
| `SCHED_MIN_LR` | 1e-6 | Learning rate minimo |
### Fine-tuning L-BFGS
| Parametro | Default | Descrizione |
|---------------|---------|--------------------------|
| `LR_LBFGS` | 0.1 | Learning rate L-BFGS |
| `LBFGS_STEPS` | 200 | Numero di step L-BFGS |
### Pesi della loss
| Parametro | Default | Descrizione |
|-----------|---------|--------------------------------|
| `W_PDE` | 10.0 | Peso residuo PDE |
| `W_IC` | 1.0 | Peso condizione iniziale |
| `W_BC` | 5.0 | Peso condizioni al contorno |
> **Se la loss diverge:** verificare che `T_char = Q_VAL * L / K` non sia vicino a zero. Questo valore è la scala caratteristica di temperatura usata per normalizzare tutti i termini.
---
## 12. Test
La test suite è in `tests/` e conta **42 test** organizzati in tre livelli:
| File | Tipo | Cosa testa |
|------------------------------|-------------|---------------------------------------------------|
| `test_config.py` | Unit | Validità e coerenza dei parametri in `config.py` |
| `test_model.py` | Unit | Shape output, finitezza, loss components |
| `test_engine_data.py` | Unit | Campionamento e clustering dei punti |
| `test_fdm_solver.py` | Unit | Griglia, CFL, shape output del solver FDM |
| `test_integration_pinn.py` | Integration | Caricamento modello, griglia predizione, pipeline loss |
| `test_e2e.py` | End-to-end | Workflow completo: train → evaluate → visualize |
```bash
pytest tests/ -v # run tutti i test
pytest tests/ -k "model" # solo test con "model" nel nome
pytest tests/ --tb=short # traceback breve in caso di fallimento
```
---
## Dettagli implementativi
### Normalizzazione automatica della loss
Le scale sono precompilate una sola volta a import time in `model.py`:
```python
_T_char = Q_VAL * L / K # ~150 °C
_src_peak = ALPHA * Q_VAL / (K * GAUSS_SIGMA * sqrt(2π))
_pde_scale = max((_T_char / T_END)², _src_peak²) + 1e-8
_bc_scale = max(Q_VAL / K, H_CONV * _T_char / K) ** 2
```
Dividere ogni termine per la sua scala porta tutti i contributi a `O(1)`, rendendo i pesi `W_PDE`, `W_IC`, `W_BC` interpretabili come importanza relativa piuttosto che come fattori di scala assoluti.
### Rilevamento device
`engine.py` seleziona automaticamente il device più performante disponibile:
```python
CUDA MPS (Apple Silicon) CPU
```
Include un test di funzionamento effettivo della GPU prima di usarla, per evitare fallimenti silenziosi su driver incompleti.
### Closure L-BFGS
L-BFGS richiede una funzione `closure()` che esegue `zero_grad`, forward pass, `backward`, e restituisce la loss. I componenti della loss vengono catturati in un dizionario `_last` per permettere il logging a ogni step senza ricalcolare il grafo fuori dal contesto `backward`.
### Subsampling delle animazioni FDM
Se `NT > 200`, il visualizer FDM campiona ogni `n`-esimo frame (`n = NT // 200`) per mantenere le animazioni HTML leggere e fluide.