Scambia calore su tutto il contorno e per conduzione circonferenziale

- Convezione anche sui bordi laterali x = 0 e x = lunghezza (nuovo
  coefficiente h_bordi_W_m2K in ARIA)
- Termine di conduzione lungo y verso il resto della fascetta, assunto a
  temperatura ambiente: q_y = -k*(T - T_amb)/L_y^2 con L_y configurabile
  (lunghezza_conduzione_y_mm in FASCETTA)
- La temperatura iniziale del campo e' la temperatura ambiente del run
  (rimosso temperatura_iniziale_C)
- Fisica estratta in prepara_stato_termico e passo_implicito, condivise
  tra simulate.py e plot_animazione.py per evitare duplicazione
- Entrambi i termini sono lineari e restano nella matrice fattorizzata
  una volta per run

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
This commit is contained in:
2026-07-05 19:22:51 +02:00
parent c92fd26864
commit c67f6f8da6
4 changed files with 130 additions and 101 deletions
+102 -51
View File
@@ -143,14 +143,18 @@ def costruisci_solutore_implicito_2d(
materiale: dict,
h_esterno_W_m2K: float,
h_interno_W_m2K: float,
h_bordi_W_m2K: float,
lunghezza_conduzione_y_m: float,
):
# Costruisce e fattorizza (LU sparsa) la matrice A per Eulero implicito 2D:
# A * T_next = rhs
#
# Le incognite sono i centri cella T[i, j] con i lungo x e j lungo z,
# appiattiti in ordine C (indice = i * n_z + j). I bordi x = 0 e
# x = lunghezza sono adiabatici; i bordi z = 0 (esterno) e z = spessore
# (interno) includono la convezione verso l'ambiente.
# appiattiti in ordine C (indice = i * n_z + j). Tutti e quattro i lati
# della sezione includono la convezione verso l'ambiente; su ogni cella
# agisce inoltre la conduzione circonferenziale (y) verso il resto della
# fascetta, assunto a temperatura ambiente e modellato come
# q_y = -k * (T - T_amb) / L_y².
k = materiale["conducibilita_termica_W_mK"]
rho = materiale["densita_kg_m3"]
cp = materiale["calore_specifico_J_kgK"]
@@ -160,21 +164,99 @@ def costruisci_solutore_implicito_2d(
r_z = alpha * dt_s / (dz_m * dz_m)
b_esterno = h_esterno_W_m2K * dt_s / (rho * cp * dz_m)
b_interno = h_interno_W_m2K * dt_s / (rho * cp * dz_m)
b_bordo = h_bordi_W_m2K * dt_s / (rho * cp * dx_m)
c_y = alpha * dt_s / (lunghezza_conduzione_y_m * lunghezza_conduzione_y_m)
n = n_x * n_z
convezione = np.zeros(n)
convezione[0::n_z] = b_esterno
convezione[n_z - 1::n_z] = b_interno
scambio = np.full(n, c_y)
scambio[0::n_z] += b_esterno
scambio[n_z - 1::n_z] += b_interno
scambio[:n_z] += b_bordo
scambio[n - n_z:] += b_bordo
A = (
sp.identity(n)
+ r_x * sp.kron(_laplaciano_1d(n_x), sp.identity(n_z))
+ r_z * sp.kron(sp.identity(n_x), _laplaciano_1d(n_z))
+ sp.diags(convezione)
+ sp.diags(scambio)
)
return splu(sp.csc_matrix(A))
def prepara_stato_termico(fascetta: dict, aria: dict, sorgente: dict) -> dict:
# Prepara griglia, coefficienti e solutore fattorizzato per un run:
# tutto ciò che resta costante durante l'integrazione temporale.
materiale = MATERIALI[fascetta["materiale"]]
lunghezza = fascetta["lunghezza_mm"] / 1000.0
spessore = fascetta["spessore_mm"] / 1000.0
n_x = fascetta["n_nodi_x"]
n_z = fascetta["n_nodi_z"]
dx = lunghezza / n_x
dz = spessore / n_z
lunghezza_y = fascetta["lunghezza_conduzione_y_mm"] / 1000.0
dt = SIMULAZIONE["dt_interno_s"]
if sorgente["skin_depth_fissa_m"] is None:
skin_depth = calcola_skin_depth_m(materiale, sorgente["frequenza_hz"])
else:
skin_depth = float(sorgente["skin_depth_fissa_m"])
rho = materiale["densita_kg_m3"]
cp = materiale["calore_specifico_J_kgK"]
alpha = materiale["conducibilita_termica_W_mK"] / (rho * cp)
z_centri = (np.arange(n_z) + 0.5) * dz
return {
"n_x": n_x,
"n_z": n_z,
"dx_m": dx,
"dz_m": dz,
"dt_s": dt,
"x_centri_m": (np.arange(n_x) + 0.5) * dx,
"z_centri_m": z_centri,
"skin_depth_m": skin_depth,
"rho": rho,
"cp": cp,
"profilo_z": profilo_deposizione_z_1_m(z_centri, spessore, skin_depth),
"b_esterno": aria["h_esterno_W_m2K"] * dt / (rho * cp * dz),
"b_interno": aria["h_interno_W_m2K"] * dt / (rho * cp * dz),
"b_bordo": aria["h_bordi_W_m2K"] * dt / (rho * cp * dx),
"c_y": alpha * dt / (lunghezza_y * lunghezza_y),
"T_ambiente_C": aria["temperatura_ambiente_C"],
"solutore": costruisci_solutore_implicito_2d(
n_x=n_x,
n_z=n_z,
dt_s=dt,
dx_m=dx,
dz_m=dz,
materiale=materiale,
h_esterno_W_m2K=aria["h_esterno_W_m2K"],
h_interno_W_m2K=aria["h_interno_W_m2K"],
h_bordi_W_m2K=aria["h_bordi_W_m2K"],
lunghezza_conduzione_y_m=lunghezza_y,
),
}
def passo_implicito(stato: dict, T: np.ndarray, q_x: np.ndarray) -> np.ndarray:
# Avanza il campo di temperatura di un passo dt: assembla il termine noto
# (sorgente volumetrica, convezione sui quattro lati, conduzione
# circonferenziale verso l'ambiente) e risolve il sistema implicito.
T_amb = stato["T_ambiente_C"]
rhs = T + (stato["dt_s"] / (stato["rho"] * stato["cp"])) * (
q_x[:, None] * stato["profilo_z"][None, :]
)
rhs += stato["c_y"] * T_amb
rhs[:, 0] += stato["b_esterno"] * T_amb
rhs[:, -1] += stato["b_interno"] * T_amb
rhs[0, :] += stato["b_bordo"] * T_amb
rhs[-1, :] += stato["b_bordo"] * T_amb
return stato["solutore"].solve(rhs.ravel()).reshape(T.shape)
def quantizza(valore: float, passo: float) -> float:
if passo <= 0.0:
return valore
@@ -235,49 +317,21 @@ def simula_singolo(cfg_run: dict, output_csv: Path, rng: random.Random) -> dict:
sensore = cfg_run["sensore"]
nome_materiale = fascetta["materiale"]
materiale = MATERIALI[nome_materiale]
lunghezza = fascetta["lunghezza_mm"] / 1000.0
spessore = fascetta["spessore_mm"] / 1000.0
n_x = fascetta["n_nodi_x"]
n_z = fascetta["n_nodi_z"]
dx = lunghezza / n_x
dz = spessore / n_z
x_centri = (np.arange(n_x) + 0.5) * dx
z_centri = (np.arange(n_z) + 0.5) * dz
stato = prepara_stato_termico(fascetta, aria, sorgente)
n_x = stato["n_x"]
n_z = stato["n_z"]
x_centri = stato["x_centri_m"]
skin_depth = stato["skin_depth_m"]
x_sensore = sensore["x_mm"] / 1000.0
i_sensore = min(n_x - 1, max(0, int(x_sensore / dx)))
i_sensore = min(n_x - 1, max(0, int(x_sensore / stato["dx_m"])))
dt = SIMULAZIONE["dt_interno_s"]
dt = stato["dt_s"]
durata = SIMULAZIONE["durata_s"]
periodo_campionamento = 1.0 / SIMULAZIONE["frequenza_campionamento_hz"]
if sorgente["skin_depth_fissa_m"] is None:
skin_depth = calcola_skin_depth_m(materiale, sorgente["frequenza_hz"])
else:
skin_depth = float(sorgente["skin_depth_fissa_m"])
solutore = costruisci_solutore_implicito_2d(
n_x=n_x,
n_z=n_z,
dt_s=dt,
dx_m=dx,
dz_m=dz,
materiale=materiale,
h_esterno_W_m2K=aria["h_esterno_W_m2K"],
h_interno_W_m2K=aria["h_interno_W_m2K"],
)
rho = materiale["densita_kg_m3"]
cp = materiale["calore_specifico_J_kgK"]
b_esterno = aria["h_esterno_W_m2K"] * dt / (rho * cp * dz)
b_interno = aria["h_interno_W_m2K"] * dt / (rho * cp * dz)
profilo_z = profilo_deposizione_z_1_m(z_centri, spessore, skin_depth)
T = np.full((n_x, n_z), fascetta["temperatura_iniziale_C"], dtype=float)
T = np.full((n_x, n_z), aria["temperatura_ambiente_C"], dtype=float)
T_sensore = T[i_sensore, -1]
prossimo_campione_t = 0.0
@@ -310,12 +364,7 @@ def simula_singolo(cfg_run: dict, output_csv: Path, rng: random.Random) -> dict:
x_sorgente, q_x = profilo_flusso_incidente_W_m2(
sorgente, x_sensore, t, x_centri
)
rhs = T + (dt / (rho * cp)) * q_x[:, None] * profilo_z[None, :]
rhs[:, 0] += b_esterno * aria["temperatura_ambiente_C"]
rhs[:, -1] += b_interno * aria["temperatura_ambiente_C"]
T = solutore.solve(rhs.ravel()).reshape(n_x, n_z)
T = passo_implicito(stato, T, q_x)
# Temperatura vera della superficie interna nel punto osservato
# dal sensore infrarosso.
@@ -358,8 +407,9 @@ def simula_singolo(cfg_run: dict, output_csv: Path, rng: random.Random) -> dict:
"file_csv": str(output_csv.name),
"materiale": nome_materiale,
"diametro_m": fascetta["diametro_mm"] / 1000.0,
"lunghezza_m": lunghezza,
"spessore_m": spessore,
"lunghezza_m": fascetta["lunghezza_mm"] / 1000.0,
"spessore_m": fascetta["spessore_mm"] / 1000.0,
"lunghezza_conduzione_y_m": fascetta["lunghezza_conduzione_y_mm"] / 1000.0,
"n_nodi_x": n_x,
"n_nodi_z": n_z,
"durata_s": durata,
@@ -368,6 +418,7 @@ def simula_singolo(cfg_run: dict, output_csv: Path, rng: random.Random) -> dict:
"temperatura_ambiente_C": aria["temperatura_ambiente_C"],
"h_esterno_W_m2K": aria["h_esterno_W_m2K"],
"h_interno_W_m2K": aria["h_interno_W_m2K"],
"h_bordi_W_m2K": aria["h_bordi_W_m2K"],
"x_inizio_m": sorgente["x_inizio_m"],
"x_fine_m": sorgente["x_fine_m"],
"x_sensore_m": x_sensore,