Simulatore Termico a Elementi Shell — Fascetta, Sorgenti a Induzione Mobili, Sensore IR Fisso

Questo progetto genera misurazioni CSV pseudo-realistiche della temperatura di una fascetta (anello cilindrico sottile) riscaldata da un gruppo di sorgenti a induzione in movimento, osservata da un sensore a infrarossi fisso. Lo scopo è produrre dataset per l'addestramento e la validazione di modelli di stima/regressione termica.

Geometria

La fascetta è un anello cilindrico definito da tre dimensioni:

• diametro medio (default 70 mm) — il diametro della superficie media, su cui giace la mesh;
• spessore (default 0.18 mm) — lo spessore della parete;
• lunghezza (default 100 mm) — l'estensione assiale lungo x.

Il dominio simulato è l'intera superficie cilindrica, discretizzata con elementi shell quadrilateri. La superficie è parametrizzata da (x, θ):

                 sorgenti (induttori), in moto verso -x
                       ▼           ▼           ▼
                  ╭────────────────────────────────╮
     θ = 0  ─ ─ ─ ┤   ●  ← punto osservato dal     ├ ─ ─ ─  generatrice del sensore
                  │      sensore IR (x = 50 mm)    │
                  ╰────────────────────────────────╯
                  x = 0                     x = lunghezza

     superficie sviluppata in piano; il bordo superiore e quello inferiore
     sono lo stesso punto della fascetta (θ = 0 ≡ θ = 2π)

• x = direzione della lunghezza, da 0 a lunghezza. Le sorgenti viaggiano in direzione -x sulla superficie esterna.
• θ = coordinata circonferenziale, periodica e risolta spazialmente. La generatrice θ = 0 è quella del sensore; l'offset circonferenziale delle sorgenti (offset_y_percorso_m, misurato lungo l'arco) è riferito ad essa.
• Lo spessore non è discretizzato: la parete è termicamente sottile, quindi la faccia esterna e quella interna hanno la stessa temperatura.

Il sensore è un pirometro a infrarossi posto all'interno della fascetta, a una distanza configurabile dalla parete interna (default 10 mm). Essendo senza contatto, la distanza non influenza la misura: il sensore legge la temperatura della superficie interna nel punto x configurato (default 50 mm, al centro della lunghezza) sulla generatrice θ = 0.

Modello fisico

Il campo elettromagnetico non è simulato: la parte termica è però una FEM vera. La catena di modellazione è la seguente.

1. Shell termica — per una parete sottile con temperatura uniforme nello spessore l'equazione del calore sulla superficie cilindrica è

   ρ cp s ∂T/∂t = k s ( ∂²T/∂x² + (1/R²) ∂²T/∂θ² ) + q'' − q_conv
   

   La conduzione assiale e quella circonferenziale sono entrambe risolte. Discretizzando con elementi quadrilateri a 4 nodi e una temperatura per nodo si ottiene C Ṫ + K T = f(t), con

   Cₑ      = ∫ ρ cp s Nᵀ N dA          (matrice consistente)
   Kₑ_cond = ∫ Bᵀ D B s dA             D = k·I  (materiale isotropo)
   

   La conduzione tra elementi adiacenti nasce dall'assemblaggio: due elementi confinanti condividono i nodi del bordo comune, quindi condividono gli stessi gradi di libertà termici. Non serve alcuna legge di scambio tra elementi dello stesso materiale.

2. Convezione — nelle shell la faccia esterna e quella interna sono entrambe superfici fisiche esposte, quindi ogni elemento scambia su tutta la sua area:

   Kₑ_conv = ∫ (h_est + h_int) Nᵀ N dA        fₑ_conv = ∫ (h_est + h_int) T∞ Nᵀ dA
   

   I due bordi assiali espongono solo lo spessore della lamiera: Kₑ_bordo = ∫ h_b s Nᵀ N dl. Con s = 0.18 mm su 100 mm di lunghezza il contributo dei bordi vale ~0.2% dell'area ed è quindi marginale.

3. Sorgenti gaussiane in moto — un gruppo di numero_sorgenti sorgenti equidistanti (distanza_sorgenti_m) si muove rigidamente a velocità costante lungo la generatrice individuata da offset_y_percorso_m. Il flusso sulla superficie esterna è

   q''(x, θ, t) = q_picco · efficienza · Σᵢ exp( −((x − xᵢ(t))² + Δy²) / (2σ²) )
   Δy = R · wrap(θ − θ_s)
   

   La distanza circonferenziale è riavvolta su 2πR: la sorgente scalda anche i punti raggiunti dall'altro verso della circonferenza. Il carico nodale è fₑ_src = ∫ Nᵀ q'' dA, integrato con la stessa quadratura di Gauss delle matrici: la posizione circonferenziale dell'induttore è quindi simulata davvero, non rappresentata da un coefficiente di attenuazione.

4. Skin depth — δ = √(2ρₑ/(ωμ)), calcolata dalla resistività elettrica e dalla permeabilità del materiale, oppure imposta con skin_depth_fissa_m. Poiché la shell non discretizza lo spessore, la skin depth non entra nella soluzione: è registrata in metadata.csv come diagnostica. Per la banda stagnata a 20 kHz vale 0.142 mm, confrontabile con lo spessore di 0.18 mm — cioè la deposizione è già quasi uniforme attraverso la parete.

5. Temperatura iniziale — il campo parte uniformemente alla temperatura ambiente del run (che è randomizzata, quindi varia run per run).

6. Sensore realistico — la temperatura nel punto osservato è interpolata con le funzioni di forma dell'elemento che lo contiene; la lettura vi aggiunge inerzia del primo ordine (costante_tempo_s), rumore gaussiano (rumore_std_C) e quantizzazione (quantizzazione_C).

Perché una sola temperatura nello spessore

La parete è termicamente sottile con ampio margine. Il numero di Fourier che confronta il tempo di diffusione attraverso lo spessore con il tempo di transito della sorgente,

Fo = α · (σ/|v|) / s²

vale 5.4 con i parametri di default (ed è riportato per ogni run in metadata.csv). Il modello 2D (x, z) precedente, che risolveva lo spessore con 15 nodi, dava un salto esterno–interno di 0.002 °C al picco su una sovratemperatura di ~200 °C, e 0.53 °C come massimo su tutto il run. Risolvere lo spessore non aggiunge informazione a questi parametri.

Il criterio per sapere quando servirebbe una formulazione multistrato è Fo ≈ 1, cioè uno spessore critico s* = √(α σ/|v|):

materiale	Fo a 0.18 mm	spessore critico
alluminio	31.3	1.01 mm
ferro	8.4	0.52 mm
banda stagnata / acciaio dolce	5.4	0.42 mm
inox 304	1.5	0.22 mm

Metodo numerico

• Mesh strutturata n_elementi_x × n_elementi_theta (default 100 × 110, cioè elementi di 1.0 × 2.0 mm) per 11110 nodi; la periodicità circonferenziale è strutturale, non un vincolo imposto: non esiste un nodo separato per θ = 2π.
• Quadratura di Gauss 2×2 per tutte le matrici elementari e per il vettore della sorgente. La mesh è uniforme, quindi tutti gli elementi sono lo stesso rettangolo: le matrici elementari si calcolano una volta e si replicano.
• Eulero implicito (incondizionatamente stabile) con passo dt_interno_s (default 1 ms), più fine del periodo di campionamento CSV.
• Tutti i termini sono lineari e costanti nel run: A = C/dt + K_cond + K_conv viene fattorizzata LU una sola volta per run (scipy.sparse.linalg.splu); ogni passo risolve solo i sistemi triangolari. A è simmetrica definita positiva, quindi si usa l'ordinamento MMD_AT_PLUS_A, che dimezza circa il riempimento rispetto al COLAMD di default. Un run da 30 s simulati richiede ~65 s di calcolo; i run sono generati in parallelo su più processi.
• Se in futuro si introducessero proprietà dipendenti dalla temperatura, la matrice andrebbe ricostruita e rifattorizzata a ogni passo.

Validazione

python valida_shell.py verifica il solutore su casi con soluzione nota. Va rieseguito dopo ogni modifica alla fisica o all'assemblaggio. Risultati con la mesh di default:

Verifica	Esito
capacità totale = ρ cp s A	esatta a precisione macchina
somma righe di K_cond = 0 (T costante → flusso nullo)	4·10⁻¹⁸
simmetria e semidefinita positività di K_cond	verificate
raffreddamento convettivo uniforme vs exp(−t/τ)	errore 0.011 °C, pari a dt/(2τ)
conservazione dell'energia (adiabatico)	errore relativo 4·10⁻¹⁴
Σ f_src vs potenza analitica 2πσ² q_max	rapporto 1.0000; 0.5000 con sorgente centrata sul bordo
modi sinusoidali circonferenziali n = 1, 2, 4 vs exp(−α n²t/R²)	errore 0.04% – 0.17%
convergenza al raffinamento (modo n = 4)	2.75% → 1.07% → 0.17% → 0.05% da 30 a 220 elementi
simmetria della diffusione attorno a θ = 0	2·10⁻¹⁴

L'ultima riga è la prova che la periodicità funziona davvero: un impulso su θ = 0 diffonde in modo identico nei due versi solo se la connettività riavvolge la circonferenza.

File

config.py               tutti i parametri di simulazione
materials.py            proprietà termofisiche ed elettriche dei materiali
mesh.py                 nodi, connettività Q4, periodicità circonferenziale
elementi_shell.py       funzioni di forma, quadratura, matrici del singolo elemento
assemblaggio.py         assemblaggio delle matrici globali e del vettore sorgente
sorgente.py             moto degli induttori e distribuzione gaussiana del flusso
solutore.py             integrazione implicita e fattorizzazione sparsa
sensore.py              interpolazione, inerzia, rumore, quantizzazione
simulate.py             orchestrazione dei run + generazione dei CSV
valida_shell.py         verifiche del solutore su casi con soluzione nota
plot_csv.py             grafici rapidi (temperature e flusso) del primo run
plot_animazione.py      animazione della superficie sviluppata T(x, θ)
plot_animazione_3d.py   animazione 3D isometrica del cilindro
dataset/                output generato da simulate.py (ricreato a ogni esecuzione)

Installazione

python -m venv .venv
source .venv/bin/activate
pip install -r requirements.txt

Uso

# genera il dataset (ATTENZIONE: cancella e ricrea la cartella dataset/)
python simulate.py

# grafici statici del primo run (temperatura e flusso nel tempo)
python plot_csv.py

# animazione della superficie sviluppata durante il passaggio delle sorgenti
python plot_animazione.py

# animazione 3D isometrica del cilindro
python plot_animazione_3d.py

Gli script di visualizzazione aprono finestre interattive (backend Qt); se il backend non è interattivo (es. sessione senza display) salvano automaticamente PNG/GIF in dataset/.

Entrambe le animazioni riproducono esattamente la fisica di run_0001 (stesso seed). Quella 2D mostra tre pannelli allineati: il profilo di flusso q(x) con le sorgenti in transito, il campo di temperatura sulla superficie sviluppata in piano (asse verticale = arco circonferenziale, con il sensore a quota 0), e la temperatura nel punto osservato (vera e con inerzia del sensore). Quella 3D disegna il campo calcolato sulla superficie del cilindro — non è più una ricostruzione, il modello risolve davvero la distribuzione circonferenziale. Finestra temporale e cadenza dei fotogrammi si regolano con le costanti in testa a plot_animazione.py.

La mappa sviluppata usa pcolormesh con shading Gouraud, che interpola linearmente tra i valori nodali: è esattamente ciò che fanno le funzioni di forma bilineari del Q4, quindi si guarda la soluzione dell'elemento finito e non un'interpolazione grafica.

Vedere gli elementi shell

Entrambe le animazioni hanno una costante MOSTRA_MESH che disegna i bordi degli elementi. È False di default perché con la mesh di produzione gli elementi sono 1 × 2 mm, più fitti dei pixel disponibili, e il reticolo coprirebbe il campo. Per ispezionare la mesh servono due modifiche insieme:

# in plot_animazione.py e/o plot_animazione_3d.py
MOSTRA_MESH = True

# in config.py, dentro FASCETTA
"n_elementi_x": 25,
"n_elementi_theta": 28,

Così si vedono i singoli quadrilateri, in piano e avvolti sul cilindro. Attenzione: con 28 elementi sulla circonferenza l'arco vale 7.9 mm contro un sigma_punto_m di 12 mm, cioè meno di due elementi per sigma — la mesh grossolana serve a guardare la discretizzazione, non a produrre risultati (vedi la tabella di convergenza sopra).

Configurazione

Tutto si modifica in config.py. I dizionari principali:

Dizionario	Contenuto
SIMULAZIONE	numero di run, durata, campionamento CSV, passo interno, seed, cartella
FASCETTA	diametro medio, lunghezza, spessore, densità della mesh, materiale
ARIA	temperatura ambiente e coefficienti di convezione di facce e bordi
SORGENTE	corsa, velocità (il segno dà il verso), gruppo di sorgenti, gaussiana, flusso, frequenza
SENSORE	posizione (x e distanza dalla parete), inerzia, rumore, quantizzazione
RANDOMIZZAZIONE	entità delle perturbazioni per run

Punti da conoscere:

• Corsa delle sorgenti: x_inizio_m è la distanza dal sensore, all'inizio della corsa, della sorgente più avanzata (quella che lo raggiunge per prima); x_fine_m è la distanza dal sensore, a fine corsa, della sorgente più arretrata (quella che lo supera per ultima). Con zero_dopo_fine ogni sorgente si spegne alla fine della propria corsa.
• Materiale: FASCETTA["materiale"] deve essere una chiave di MATERIALI in materials.py. Per aggiungere un materiale basta una nuova voce nel dizionario (conducibilità, densità, calore specifico, resistività elettrica, permeabilità).
• Proprietà di banda_stagnata: la banda stagnata è un nastro di acciaio a basso tenore di carbonio (0,15–0,25% C) rivestito su entrambe le facce da un sottile strato di stagno elettrolitico, dello spessore di pochi micrometri — trascurabile rispetto allo spessore tipico della parete (es. 0,18 mm) e quindi ininfluente sulle proprietà termiche, elettriche e magnetiche in massa. I valori in materials.py sono quindi quelli dell'acciaio dolce sottostante: Banda stagnata: tutto quello che c'è da sapere (MUNDOLATAS), Bande stagnate elettrolitiche (EUROPERF), Differenza tra banda stagnata e acciaio inossidabile (Wuxi Bright Packing).
• Randomizzazione: ogni run perturba velocità, flusso di picco, sigma, offset y, temperatura ambiente e rumore del sensore con estrazioni da un RNG a seed fisso (SIMULAZIONE["seed"]): il dataset è riproducibile.

Output

dataset/run_XXXX.csv — serie temporale del run

Colonna	Significato
id_run	identificativo del run
tempo_s	tempo simulato
x_sorgente_m	posizione della sorgente di riferimento del gruppo
offset_y_sorgente_m	offset circonferenziale del percorso (costante nel run)
flusso_termico_sorgente_W_m2	flusso efficace nel punto x del sensore
skin_depth_m	skin depth usata (costante nel run)
T_vera_lato_sensore_C	temperatura vera della superficie interna nel punto del sensore
T_misurata_sensore_C	lettura del sensore (inerzia + rumore + quantizzazione)
T_max_fascetta_C	temperatura massima su tutta la superficie in quell'istante
T_media_fascetta_C	temperatura media su tutta la superficie in quell'istante
T_ambiente_C	temperatura ambiente del run
velocita_m_s, sigma_punto_m, flusso_picco_W_m2	parametri randomizzati del run
materiale	chiave del materiale

dataset/metadata.csv — una riga per run

Contiene tutti i parametri effettivi del run (geometria, dimensioni della mesh, coefficienti di scambio, parametri delle sorgenti e del sensore, valori randomizzati), le diagnostiche skin_depth_m e numero_fourier_spessore, e le temperature di picco: utile come ground truth e per filtrare i run.

Limitazioni

1. Il campo elettromagnetico non è simulato: l'accoppiamento induttivo è ridotto a un'impronta gaussiana moltiplicata per un'efficienza costante.
2. Lo spessore non è discretizzato: il modello non può rappresentare un gradiente tra faccia esterna e faccia interna. Con i parametri di default il gradiente reale è di pochi millesimi di grado (vedi numero_fourier_spessore), ma con pareti più spesse di ~0.4 mm, materiali a bassa diffusività o sorgenti molto più veloci servirebbe una shell multistrato.
3. Le proprietà dei materiali sono costanti con la temperatura; per gli acciai ferromagnetici la skin depth reale varia fortemente con temperatura e campo (punto di Curie non modellato).
4. Irraggiamento non modellato: a ~210 °C le perdite radiative non sono del tutto trascurabili rispetto alla convezione.