OpenFOAM: simulazioni transitorie, passo temporale, Courant e animazioni ParaView

Prima di scegliere deltaT e endTime conviene stimare le scale temporali del problema. Con lato L = 0,1 m, velocità del coperchio U = 1 m/s e viscosità ν = 1·10⁻⁵ m²/s:

Il rapporto t_ν / t_c = Re dice che la quantità di moto si propaga dal coperchio per convezione e, nel modello RAS, per viscosità turbolenta, non per diffusione molecolare. Se la grandezza misurata alla sonda è ancora in evoluzione a t = 2 s, il calcolo descrive un transitorio di avviamento: per affermare di aver raggiunto un regime occorre prolungare endTime e mostrare che la serie si è stabilizzata su più tempi convettivi.

Dopo la discretizzazione spaziale le equazioni diventano un sistema di equazioni differenziali ordinarie nel tempo, dU/dt = R(U), dove R raccoglie convezione, diffusione e gradiente di pressione. Lo schema Euler implicito valuta R al nuovo istante:

(U^(n+1) − U^n) / Δt = R(U^(n+1))

errore locale di troncamento  ~ Δt²   (per singolo passo)
errore globale accumulato     ~ Δt    (sull'intervallo fisico fissato)

La conseguenza operativa è semplice: nel campo asintotico, dimezzare deltaT dimezza circa l’errore temporale delle grandezze osservate. È questo che rende utile il confronto fra cavity_dt005 e cavity_dt0025: se la differenza fra i due casi non è piccola rispetto alla precisione desiderata, il passo più grande non è adeguato. Ogni passo richiede la soluzione di sistemi lineari (quantità di moto, pressione, grandezze turbolente): per questo un passo implicito costa molto più di uno esplicito, ma tollera passi più grandi senza divergere.

Tre comandi da tenere a portata di mano mentre il calcolo procede, dentro la cartella del caso:

grep -c "^Time = " log.foamRun            # passi temporali completati
grep "^ExecutionTime" log.foamRun | tail -1   # secondi di CPU e di orologio finora
foamListTimes                              # istanti salvati su disco (tempo fisico)

Dal rapporto fra ExecutionTime e passi completati si ricava il costo per passo; moltiplicandolo per i passi mancanti si ottiene una stima affidabile del tempo residuo, molto più utile della sensazione «sta andando lento».

Nella letteratura sulla cavità bidimensionale laminare la prima transizione da soluzione stazionaria a periodica è riportata attorno a Re ≈ 8000: a Re = 10⁴ un calcolo laminare 2D ben risolto potrebbe quindi già oscillare. Il tutorial però usa un modello RAS kEpsilon: la viscosità turbolenta aggiunge dissipazione e la soluzione RAS tende di solito verso uno stato quasi stazionario, che rappresenta il moto medio, non le fluttuazioni. Per questo la guida parla di transitorio di avviamento: la grandezza alla sonda cresce mentre il vortice primario si forma e poi si assesta.

Due dettagli fisici da ricordare: agli spigoli superiori la velocità passa bruscamente da 1 m/s (coperchio) a 0 (pareti laterali), una singolarità che rende i valori locali in quei punti dipendenti dalla mesh; la sonda al centro è lontana da questi spigoli ed è quindi una scelta robusta.

La cavità è chiusa: su tutte le pareti la pressione ha gradiente nullo e nessuna patch ne fissa il valore. Nel caso incomprimibile la pressione compare nelle equazioni solo attraverso il suo gradiente, quindi è definita a meno di una costante. Il solver la rende unica fissando il valore in una cella di riferimento, tramite pRefCell e pRefValue (oppure pRefPoint) nel dizionario dell’algoritmo in system/fvSolution. Senza questo vincolo l’equazione di pressione è singolare.

# verifica sul tutorial installato
cd "$FOAM_TUTORIALS/incompressibleFluid/cavity"
foamDictionary system/fvSolution -entry PIMPLE
foamDictionary 0/p -entry boundaryField
foamDictionary 0/U -entry boundaryField

Conseguenza pratica: il valore assoluto di p in questo caso non ha significato fisico, solo le differenze di pressione lo hanno. Confrontando due calcoli si confrontano differenze, mai livelli assoluti.

Le celle misurano 5 mm, quindi il primo centro cella dista 2,5 mm dalla parete. Ipotizzando una velocità d’attrito di qualche percento della velocità del coperchio (u_τ ≈ 0,03 m/s), si ottiene y⁺ = y·u_τ/ν ≈ 0,0025 · 0,03 / 10⁻⁵ ≈ 7,5: un valore nello strato di transizione, sotto l’intervallo in cui le funzioni di parete standard sono pensate per lavorare (indicativamente 30–300). È una stima d’ordine di grandezza, ma basta a confermare la nota della guida: questa configurazione è didattica, adatta a studiare l’effetto di deltaT, non a produrre valori di parete affidabili.

Per leggere il valore effettivo a fine calcolo si può usare la function object yPlus (foamPostProcess -func yPlus -solver incompressibleFluid), verificandone la disponibilità con foamPostProcess -list.

Dieci minuti spesi a leggere i file originali evitano ore di dubbi dopo. Senza modificare nulla:

T="$FOAM_TUTORIALS/incompressibleFluid/cavity"
ls -R "$T" | head -40
cat "$T/constant/physicalProperties" "$T/constant/momentumTransport"
foamDictionary "$T/system/fvSchemes" -entry ddtSchemes
foamDictionary "$T/system/blockMeshDict" -entry blocks

Annota i valori trovati (viscosità, modello, schemi, numero di celle) nel quaderno del progetto: saranno il riferimento per capire se un risultato diverso nasce da una tua modifica o da una differenza dell’installazione.

#!/usr/bin/env python3
"""Prepare two reproducible transient cavity cases from OpenFOAM Foundation 14.

Copyright (C) 2026. License: GNU GPL v3 or later.
This script preserves the original tutorial and refuses existing destinations.
"""
from __future__ import annotations
import argparse
import os
from pathlib import Path
import re
import shutil
import sys

SPECS = (("cavity_dt005", "0.005", "10"), ("cavity_dt0025", "0.0025", "20"))


def replace_entry(source: str, key: str, value: str) -> str:
    pattern = re.compile(r"(?m)^(\s*" + re.escape(key) + r"\s+)[^;\n]+(;)")
    updated, count = pattern.subn(lambda m: m.group(1) + value + m.group(2), source, count=1)
    if count != 1:
        raise ValueError(f"Expected exactly one setting for {key}")
    return updated


def prepare(template: Path, destination: Path, timestep: str, interval: str) -> None:
    if destination.exists():
        raise FileExistsError(f"Destination already exists: {destination}")
    if not all((template / p).is_file() for p in ("system/controlDict", "system/fvSchemes", "system/blockMeshDict", "0/U", "0/p")):
        raise FileNotFoundError("OpenFOAM 14 cavity tutorial files are missing")
    source = (template / "system/controlDict").read_text(encoding="utf-8")
    for key, value in (("startFrom", "startTime"), ("endTime", "2"), ("deltaT", timestep),
                       ("writeControl", "timeStep"), ("writeInterval", interval),
                       ("purgeWrite", "0")):
        source = replace_entry(source, key, value)
    # The official cavity tutorial has no function objects in controlDict.
    if re.search(r"(?m)^\s*functions\s*\{", source):
        raise ValueError("Unexpected functions dictionary: edit manually")
    source += """
// Monitoring: a fixed point in the interior of the 0.1 m square cavity.
functions
{
    centreVelocity
    {
        type                probes;
        libs                (sampling);
        fields              (U);
        probeLocations      ((0.05 0.05 0.005));
        executeControl      timeStep;
        executeInterval     1;
        writeControl        timeStep;
        writeInterval       1;
    }
}
"""
    shutil.copytree(template, destination, symlinks=True)
    (destination / "system/controlDict").write_text(source, encoding="utf-8")
    (destination / "GUIDE08_METADATA.txt").write_text(
        f"Source: OpenFOAM Foundation 14 tutorials/incompressibleFluid/cavity\n"
        f"deltaT: {timestep} s\nwriteInterval: {interval} steps (=0.05 s)\n"
        "Compare both runs at identical physical times; all data belong to the tutorial.\n",
        encoding="utf-8")


def main() -> int:
    parser = argparse.ArgumentParser(description=__doc__)
    parser.add_argument("--template", type=Path, default=None, help="Path to official cavity tutorial")
    parser.add_argument("--out", type=Path, default=Path("guide08_cases"), help="New output directory")
    args = parser.parse_args()
    foam_tutorials = os.environ.get("FOAM_TUTORIALS", "")
    template = args.template or (Path(foam_tutorials) / "incompressibleFluid/cavity" if foam_tutorials else None)
    if template is None or not template.is_dir():
        parser.error("OpenFOAM tutorials not found: source OpenFOAM environment or use --template")
    out = args.out.resolve()
    if out.exists():
        parser.error(f"Output already exists: {out}. Choose a new --out path")
    # Validate all modifications before creating anything.
    original = (template / "system/controlDict").read_text(encoding="utf-8")
    for key in ("startFrom", "endTime", "deltaT", "writeControl", "writeInterval", "purgeWrite"):
        replace_entry(original, key, "123")
    out.mkdir(parents=True)
    try:
        for name, timestep, interval in SPECS:
            prepare(template, out / name, timestep, interval)
    except Exception:
        # Preserve any partial output for troubleshooting. No cleanup or deletion.
        print(f"Preparation failed; preserved partial files at: {out}", file=sys.stderr)
        raise
    print(f"Prepared two cases at {out}")
    return 0


if __name__ == "__main__":
    raise SystemExit(main())

Le scelte di progetto che vale la pena imitare nei propri script:

• validazione prima della scrittura: replace_entry(original, key, "123") prova tutte le sostituzioni sul testo originale prima di creare una sola cartella; se una chiave manca o compare più volte, lo script si ferma senza lasciare residui;
• una sola occorrenza: l’espressione regolare con count=1 e il controllo count != 1 impediscono di modificare per errore una seconda voce con lo stesso nome;
• nessuna cancellazione: in caso di errore i file parziali restano dove sono, per la diagnosi; la destinazione esistente provoca il rifiuto, mai la sovrascrittura;
• tracciabilità: GUIDE08_METADATA.txt registra origine, passo e frequenza di scrittura in ogni caso.

La regola è writeInterval × deltaT = costante. Rispettarla garantisce che i campi salvati dei diversi casi cadano esattamente sugli stessi tempi fisici, così il confronto in ParaView e nei CSV non richiede interpolazione dei campi completi. Lo spazio su disco cresce con il numero di istanti salvati, non con il numero di passi: qui resta lo stesso per tutti e tre i casi.

R="$HOME/OpenFOAM/guide08_run01"
diff "$FOAM_TUTORIALS/incompressibleFluid/cavity/system/controlDict" "$R/cavity_dt005/system/controlDict"
for k in deltaT writeInterval endTime writeControl; do
  printf "%-14s %s  %s\n" "$k" \
    "$(foamDictionary "$R/cavity_dt005/system/controlDict" -entry $k -value)" \
    "$(foamDictionary "$R/cavity_dt0025/system/controlDict" -entry $k -value)"
done
cat "$R/cavity_dt005/GUIDE08_METADATA.txt"

Il diff deve mostrare solo le voci attese e il blocco functions aggiunto. Se compaiono altre differenze, il tutorial installato non è quello che la guida presuppone: fermati e annotalo prima di eseguire.

Una configurazione tipica per un passo regolato dal Courant è la seguente; prima di usarla verifica le voci nei tutorial della tua installazione, perché il comportamento esatto dipende dal solver:

adjustTimeStep  yes;
maxCo           0.5;
maxDeltaT       0.01;

writeControl    adjustableRunTime;
writeInterval   0.05;      // secondi fisici, non passi

// trovare esempi nella distribuzione:
// grep -rl "adjustTimeStep" "$FOAM_TUTORIALS" | head

Con adjustableRunTime il solver accorcia gli ultimi passi prima di ogni istante di scrittura per centrarlo esattamente. Il passo effettivo compare nel log come deltaT = .... In un confronto di sensibilità, però, il passo adattivo cambia due cose insieme (passo medio e sua distribuzione nel tempo): per isolare l’effetto di deltaT, come fa questa guida, il passo fisso è la scelta corretta.

Con runTimeModifiable true si interviene sul calcolo in corso modificando il dizionario dalla seconda sessione SSH:

# scrive l'istante corrente e termina in modo pulito
foamDictionary system/controlDict -entry stopAt -set writeNow

# prolunga il calcolo prima che finisca
foamDictionary system/controlDict -entry endTime -set 4

# ripartenza successiva dall'ultimo istante salvato
foamDictionary system/controlDict -entry startFrom -set latestTime
foamDictionary system/controlDict -entry stopAt -set endTime

Mai interrompere con Ctrl+C se servono i dati fino a quel punto: si perde tutto ciò che è stato calcolato dopo l’ultimo salvataggio. Ricorda che un riavvio crea una nuova sottocartella della sonda (sezione 8) e annota nel quaderno ogni modifica fatta a calcolo in corso.

Il caso didattico ha 400 celle e occupa pochissimo. In un caso industriale la stima si fa così: celle × grandezze scalari salvate × byte per valore × istanti. Con 5 milioni di celle, circa 10 grandezze scalari per istante (U conta 3, poi p, k, epsilon, nut e il flusso phi sulle facce) e 8 byte in binario, ogni istante pesa circa 0,4 GB; in ASCII anche il doppio o il triplo. Quaranta istanti diventano 16–40 GB.

Le leve per controllarlo sono: writeFormat binary, writeCompression on, un writeInterval scelto sulla dinamica da visualizzare e, soprattutto, function object come probes che registrano le serie temporali a ogni passo con costo trascurabile. La regola pratica: i campi completi servono per le animazioni, le serie temporali per i numeri.

Il valore stampato nel log è calcolato cella per cella dai flussi volumetrici alle facce φ_f:

Co_cell = 0.5 · Δt · Σ_f |φ_f| / V_cell

cella quadrata 2D, flusso lungo x:   Σ|φ| = 2·|u|·A,  V = A·Δx   →  Co = |u|·Δt/Δx
flusso diagonale (u = v):              Co = (|u| + |v|)·Δt/Δx  ≈  1,41·|U|·Δt/Δx

Il log riporta sia la media sia il massimo sul dominio. Nella cavità il massimo non si trova per forza sotto il coperchio: la faccia del coperchio è una parete e non porta flusso, e la velocità al centro delle celle adiacenti è inferiore a 1 m/s. La stima uniforme Co ≈ 1 è quindi un limite superiore ragionevole, non il valore che leggerai: è proprio il motivo per cui la guida invita a fidarsi del massimo stampato.

Applicando gli schemi a un’oscillazione pura di periodo P, con N = P/Δt passi per periodo, si ottiene la frazione di ampiezza che sopravvive dopo un periodo. I valori sono calcolati dalla teoria degli schemi per l’equazione modello dy/dt = iωy:

Il messaggio ingegneristico è netto: con Euler e 20 passi per periodo un’oscillazione perde circa il 60 % dell’ampiezza a ogni ciclo, quindi un fenomeno periodico può essere spento dallo schema numerico prima ancora che dal modello fisico. Per lo studio di oscillazioni servono molti passi per periodo oppure uno schema del secondo ordine. Nell’avviamento della cavità, che non è oscillatorio, questo smorzamento pesa meno, ma il principio resta: Courant basso e stabilità non garantiscono che la dinamica sia conservata.

Accanto al Courant esiste il numero di diffusione d = ν·Δt/Δx². Con ν = 10⁻⁵ m²/s, Δt = 0,005 s e Δx = 0,005 m si ottiene d = 0,002: trascurabile. Anche con una viscosità efficace dieci o cento volte maggiore per effetto del modello RAS resta piccolo. Negli schemi espliciti d > 0,5 provoca instabilità; con lo schema implicito di OpenFOAM la diffusione è trattata implicitamente, ma d resta un indicatore utile della risoluzione temporale dei fenomeni viscosi, per esempio nei casi con mesh di parete molto fitta, dove Δx² diventa minuscolo.

# massimo del Courant su tutto il calcolo
grep "^Courant Number" log.foamRun | awk '{print $NF}' | sort -g | tail -1

# serie completa da aprire come grafico
grep "^Courant Number" log.foamRun | awk '{print NR, $4, $NF}' > courant.dat

# in alternativa, se disponibile: foamLog divide il log in serie per grandezza
foamLog log.foamRun && ls logs/

Un Courant massimo che cresce durante il calcolo è il segnale che il flusso si sta accelerando o che una cella sta sviluppando velocità spurie: va guardato prima che il calcolo diverga, non dopo.

Gli schemi realmente usati nella tua installazione si trovano con grep -rh "default" --include=fvSchemes "$FOAM_TUTORIALS" | grep -E "Euler|backward|CrankNicolson" | sort | uniq -c. Quando cambi schema temporale, cambia solo quello e ripeti il confronto fra passi: lo schema del secondo ordine deve mostrare differenze fra passi che si riducono di circa quattro volte dimezzando deltaT, invece che di due.

Il contenuto effettivo si legge con foamDictionary system/fvSolution -entry PIMPLE. Aumentare nOuterCorrectors per usare Courant maggiori è legittimo, ma il costo per passo cresce in proporzione: il guadagno reale va misurato su tempo di calcolo totale e accuratezza, non sul solo numero di passi.

Quando le variabili numeriche da esplorare sono più di una, scrivi prima la matrice e poi esegui. Esempio coerente con la guida:

Ogni riga è una cartella nuova con il proprio file di metadati. Confrontare A1 con D2 non insegna nulla: cambiano due fattori insieme.

#!/usr/bin/env bash
# Run prepared cavity cases sequentially. Do not remove existing results.
# License: GNU GPL v3 or later.
set -Eeuo pipefail
if [[ $# -ne 1 ]]; then
    echo 'Usage: ./run_cases.sh /path/to/new/guide08_cases' >&2
    exit 2
fi
root=$(realpath -- "$1")
command -v blockMesh >/dev/null || { echo 'Source the OpenFOAM 14 environment first' >&2; exit 1; }
command -v checkMesh >/dev/null || exit 1
command -v foamRun >/dev/null || exit 1
for name in cavity_dt005 cavity_dt0025; do
    dir="$root/$name"
    [[ -f "$dir/system/controlDict" ]] || { echo "Missing case $dir" >&2; exit 1; }
    [[ ! -e "$dir/log.blockMesh" && ! -e "$dir/log.checkMesh" && ! -e "$dir/log.foamRun" ]] || {
        echo "Logs already exist in $dir. Existing runs are never overwritten." >&2
        exit 1
    }
    echo "Running $name"
    (cd "$dir" && blockMesh > log.blockMesh 2>&1 && checkMesh > log.checkMesh 2>&1 && foamRun > log.foamRun 2>&1)
    echo "Done: $dir"
done

• set -Eeuo pipefail: lo script termina al primo comando fallito, alla prima variabile non definita e al primo errore dentro una pipeline;
• command -v: verifica prima di tutto che l’ambiente OpenFOAM sia caricato, invece di scoprirlo a metà;
• la sottoshell ( cd "$dir" && ... ) cambia cartella solo per quel caso: un errore non lascia lo script nella cartella sbagliata;
• i tre controlli sui log esistenti impediscono di mescolare un vecchio calcolo con uno nuovo nella stessa cartella.

Struttura tipica di un passo con questa configurazione (valori illustrativi, non output di un calcolo):

Courant Number mean: 0.21 max: 0.83
Time = 0.5s

smoothSolver:  Solving for Ux, Initial residual = 0.0021, Final residual = 4.1e-07, No Iterations 3
smoothSolver:  Solving for Uy, Initial residual = 0.0035, Final residual = 6.8e-07, No Iterations 3
DICPCG:  Solving for p, Initial residual = 0.012, Final residual = 5.9e-07, No Iterations 21
time step continuity errors : sum local = 3.1e-09, global = 2.2e-20, cumulative = 1.4e-18
DICPCG:  Solving for p, Initial residual = 0.0031, Final residual = 8.7e-07, No Iterations 18
time step continuity errors : sum local = 8.4e-10, global = -1.3e-20, cumulative = 1.4e-18
smoothSolver:  Solving for epsilon, Initial residual = 0.0018, Final residual = 2.2e-07, No Iterations 3
smoothSolver:  Solving for k, Initial residual = 0.0024, Final residual = 3.0e-07, No Iterations 3
ExecutionTime = 1.92 s  ClockTime = 2 s

• la riga del Courant precede quella del tempo a cui si riferisce;
• Uz non viene risolto: con le patch empty il caso è 2D;
• nel transitorio il residuo iniziale misura quanto cambia la soluzione in un passo, non la distanza da uno stato stazionario: un residuo iniziale che scende verso zero indica che il campo si sta assestando;
• due righe di pressione e due di continuità corrispondono a nCorrectors 2;
• un errore di continuità sum local che cresce di ordini di grandezza è un segnale d’allarme più precoce dei residui.

Questo script, solo libreria standard Python, legge un log.foamRun e stampa passi, intervallo fisico, Courant massimo e quando è stato raggiunto, costo per passo, errore di continuità, numero di warning ed errori fatali; facoltativamente scrive la tabella per passo in CSV senza mai sovrascrivere. Accetta sia Time = 0.5 sia Time = 0.5s. È stato provato su un log di esempio con la struttura mostrata sopra.

#!/usr/bin/env python3
"""Summarise a transient OpenFOAM solver log: physical time, Courant, cost, continuity.

Usage: log_summary.py log.foamRun [--csv courant.csv]
Tolerant of "Time = 0.5" and "Time = 0.5s" styles. Standard library only.
License: GNU GPL v3 or later.
"""
from __future__ import annotations
import argparse, csv, re
from pathlib import Path

NUM = r"([-+]?\d*\.?\d+(?:[eE][-+]?\d+)?)"
R_TIME = re.compile(r"^Time = " + NUM + r"\s*s?\s*$")
R_CO = re.compile(r"^Courant Number mean: " + NUM + r" max: " + NUM)
R_EXE = re.compile(r"^ExecutionTime = " + NUM + r" s\s+ClockTime = " + NUM + r" s")
R_CONT = re.compile(r"continuity errors : sum local = " + NUM)


def main() -> None:
    ap = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter)
    ap.add_argument("log", type=Path)
    ap.add_argument("--csv", type=Path, help="write per-step table (never overwrites)")
    a = ap.parse_args()
    rows, cur, warn, fatal, pending = [], None, 0, 0, (None, None)
    for line in a.log.read_text(encoding="utf-8", errors="replace").splitlines():
        s = line.strip()
        if m := R_TIME.match(s):
            # the solver prints the Courant line just BEFORE the "Time =" line it refers to
            cur = {"time_s": float(m[1]), "co_mean": pending[0], "co_max": pending[1],
                   "exec_s": None, "clock_s": None, "cont_local": None}
            rows.append(cur)
            pending = (None, None)
        elif m := R_CO.match(s):
            pending = (float(m[1]), float(m[2]))
        elif cur is not None and (m := R_EXE.match(s)):
            cur["exec_s"], cur["clock_s"] = float(m[1]), float(m[2])
        elif cur is not None and (m := R_CONT.search(s)):
            cur["cont_local"] = float(m[1])   # last corrector of the step wins
        if "FOAM Warning" in s or s.startswith("--> FOAM Warning"):
            warn += 1
        if "FOAM FATAL" in s:
            fatal += 1
    if not rows:
        ap.error("no 'Time = ...' lines found: is this a solver log?")
    co = [(r["co_max"], r["time_s"]) for r in rows if r["co_max"] is not None]
    last = next((r for r in reversed(rows) if r["exec_s"] is not None), None)
    print(f"steps: {len(rows)}   physical time: {rows[0]['time_s']:g} -> {rows[-1]['time_s']:g} s")
    if co:
        cmax, tmax = max(co)
        print(f"Courant max: {cmax:.4g} at t = {tmax:g} s   final step max: {co[-1][0]:.4g}")
    if last:
        print(f"ExecutionTime: {last['exec_s']:g} s   ClockTime: {last['clock_s']:g} s   "
              f"CPU per step: {last['exec_s'] / len(rows):.4g} s")
    cont = [r["cont_local"] for r in rows if r["cont_local"] is not None]
    if cont:
        print(f"continuity 'sum local': max {max(cont):.3g}, final {cont[-1]:.3g}")
    print(f"warnings: {warn}   fatal errors: {fatal}")
    if a.csv:
        if a.csv.exists():
            ap.error(f"refusing to overwrite {a.csv}")
        with a.csv.open("w", newline="", encoding="utf-8") as f:
            w = csv.DictWriter(f, fieldnames=list(rows[0]))
            w.writeheader(); w.writerows(rows)
        print(f"wrote {len(rows)} rows to {a.csv}")


if __name__ == "__main__":
    main()

# sul server o su Windows, dopo aver scaricato i log
python3 log_summary.py cavity_dt005/log.foamRun --csv dt005_steps.csv
python3 log_summary.py cavity_dt0025/log.foamRun --csv dt0025_steps.csv

Se tmux non è installato e non hai i permessi per installarlo, l’alternativa minima è nohup scripts/run_cases.sh "$HOME/OpenFOAM/guide08_run01" > run.out 2>&1 &. Sui server condivisi aggiungi nice -n 10 davanti al comando per non penalizzare gli altri utenti.

Con celle da 5 mm le facce interne si trovano a x, y = 0, 5, 10, …, 100 mm. Il punto (0,05; 0,05) coincide quindi con lo spigolo comune a quattro celle, i cui centri stanno a (47,5 o 52,5; 47,5 o 52,5) mm. La function object probes restituisce il valore della cella che contiene il punto: quale delle quattro venga scelta dipende dall’algoritmo di ricerca. Il dato registrato è perciò il valore di una cella il cui centro dista circa 3,5 mm dal punto nominale.

Per il confronto fra due passi temporali sulla stessa mesh questo non conta, perché la cella scelta è la stessa. Conta invece quando si confrontano mesh diverse, o con dati sperimentali: in quel caso conviene un punto interno a una cella, per esempio (0.0525 0.0525 0.005), oppure un campionamento interpolato, dichiarando la scelta.

Il file postProcessing/centreVelocity/0/U ha un’intestazione commentata e una riga per passo, con il tempo e il vettore fra parentesi (struttura tipica):

# Probe 0 (0.05 0.05 0.005)
#       Probe             0
#        Time
0.005       (0.000123 -1.2e-05 0)
0.01        (0.000241 -2.3e-05 0)
...

Lo script del pacchetto, testo completo:

#!/usr/bin/env python3
"""Extract the first vector probe from OpenFOAM U time-history files to CSV.

Reads files under postProcessing/centreVelocity/<startTime>/U, including restarts.
License: GNU GPL v3 or later.
"""
from __future__ import annotations
import argparse
import csv
import math
from pathlib import Path
import re

RECORD = re.compile(r"^\s*([+\-\d.eE]+)\s+\(\s*([+\-\d.eE]+)\s+([+\-\d.eE]+)\s+([+\-\d.eE]+)\s*\)")


def parse_case(case: Path) -> list[tuple[float, float, float, float, float]]:
    paths = sorted((case / "postProcessing/centreVelocity").glob("*/U"),
                   key=lambda p: float(p.parent.name))
    if not paths:
        raise FileNotFoundError(f"No probe U file below {case}/postProcessing/centreVelocity")
    # Later restart files override earlier samples at the same physical time.
    samples: dict[float, tuple[float, float, float, float, float]] = {}
    for path in paths:
        for line in path.read_text(encoding="utf-8").splitlines():
            if not line.strip() or line.lstrip().startswith("#"):
                continue
            match = RECORD.match(line)
            if not match:
                raise ValueError(f"Unsupported probe record in {path}: {line[:120]}")
            time, ux, uy, uz = map(float, match.groups())
            if not all(math.isfinite(x) for x in (time, ux, uy, uz)):
                raise ValueError(f"Nonfinite probe data in {path}")
            samples[time] = (time, ux, uy, uz, math.sqrt(ux**2 + uy**2 + uz**2))
    return [samples[t] for t in sorted(samples)]


def main() -> None:
    parser = argparse.ArgumentParser(description=__doc__)
    parser.add_argument("case", type=Path)
    parser.add_argument("--out", type=Path, required=True)
    args = parser.parse_args()
    if args.out.exists():
        parser.error(f"Refusing to overwrite: {args.out}")
    data = parse_case(args.case)
    args.out.parent.mkdir(parents=True, exist_ok=True)
    with args.out.open("w", newline="", encoding="utf-8") as file:
        writer = csv.writer(file)
        writer.writerow(("time_s", "Ux_m_s", "Uy_m_s", "Uz_m_s", "U_magnitude_m_s"))
        writer.writerows(data)
    print(f"Wrote {len(data)} samples to {args.out}")


if __name__ == "__main__":
    main()

L’espressione RECORD accetta solo righe «tempo (ux uy uz)» con la prima sonda: se aggiungi più punti o grandezze scalari, il formato cambia e lo script si ferma con un messaggio esplicito invece di produrre dati sbagliati. È il comportamento corretto: meglio un errore chiaro che un CSV plausibile ma falso. Le cartelle di riavvio vengono ordinate per tempo e, sui tempi duplicati, prevale il file più recente.

Una sonda registra a ogni passo; ricostruirla dopo, dai campi salvati, darebbe solo un valore ogni 0,05 s. Se ti servono più punti, aggiungili prima di lanciare, in una copia nuova del caso. Esempio con centro, zona sotto il coperchio e zona vicino al fondo:

probeLocations
(
    (0.0525 0.0525 0.005)   // centro, interno a una cella
    (0.0525 0.0925 0.005)   // sotto il coperchio
    (0.0525 0.0075 0.005)   // vicino al fondo
);

Con più sonde il file contiene più vettori per riga: extract_probes.py così com’è non lo accetta e va adattato, per esempio leggendo tutte le coppie di parentesi della riga. Lascia invariato lo script originale e salva la versione modificata con un nome nuovo.

#!/usr/bin/env python3
"""Compare two probe CSV series by interpolation over their overlapping times.

The error is descriptive only; it is not a GCI or proof of temporal convergence.
License: GNU GPL v3 or later.
"""
from __future__ import annotations
import argparse
import bisect
import csv
import math
from pathlib import Path


def load_csv(path: Path) -> list[tuple[float, float]]:
    with path.open(newline="", encoding="utf-8") as f:
        points = [(float(row["time_s"]), float(row["U_magnitude_m_s"]))
                  for row in csv.DictReader(f)]
    if len(points) < 2 or any(not (math.isfinite(t) and math.isfinite(v)) for t, v in points):
        raise ValueError(f"Insufficient or invalid data: {path}")
    if any(b[0] <= a[0] for a, b in zip(points, points[1:])):
        raise ValueError(f"Times are not strictly increasing in {path}")
    return points


def interpolate(points: list[tuple[float, float]], at: float) -> float:
    times = [a for a, _ in points]
    idx = bisect.bisect_left(times, at)
    if idx < len(points) and abs(points[idx][0] - at) < 1e-11:
        return points[idx][1]
    if idx == 0 or idx == len(points):
        raise ValueError("Extrapolation is not allowed")
    t0, y0 = points[idx-1]
    t1, y1 = points[idx]
    return y0 + (y1-y0)*(at-t0)/(t1-t0)


def main() -> None:
    parser = argparse.ArgumentParser(description=__doc__)
    parser.add_argument("coarse", type=Path)
    parser.add_argument("fine", type=Path)
    parser.add_argument("--skip-before", type=float, default=0.0,
                        help="Ignore startup interval before this physical time")
    args = parser.parse_args()
    coarse, fine = load_csv(args.coarse), load_csv(args.fine)
    low = max(coarse[0][0], fine[0][0], args.skip_before)
    high = min(coarse[-1][0], fine[-1][0])
    data = [(t, y, interpolate(fine, t)) for t, y in coarse if low <= t <= high]
    if not data:
        parser.error("No overlapping time samples")
    errors = [a-b for _, a, b in data]
    rmse = math.sqrt(sum(x*x for x in errors)/len(errors))
    print(f"Overlap: {data[0][0]:.6g} to {data[-1][0]:.6g} s; samples={len(data)}")
    print(f"RMSE of |U|: {rmse:.8g} m/s")
    print(f"Final paired difference: {errors[-1]:+.8g} m/s")
    print("This two-resolution comparison alone does not quantify discretization uncertainty.")


if __name__ == "__main__":
    main()

Lo script prende i tempi del caso grossolano e vi interpola il caso fine; dato che 0,005 è multiplo di 0,0025, i tempi coincidono e l’interpolazione è di fatto esatta. L’RMSE pesa allo stesso modo tutto l’intervallo, compreso l’avviamento, dove le differenze sono in genere maggiori: l’opzione --skip-before permette di separare le due fasi. Eseguito sui CSV sintetici del pacchetto, lo script stampa:

> py scripts\compare_probes.py samples\synthetic_coarse.csv samples\synthetic_fine.csv
Overlap: 0 to 2 s; samples=101
RMSE of |U|: 0.0060293592 m/s
Final paired difference: -0.00256612 m/s
This two-resolution comparison alone does not quantify discretization uncertainty.

> py scripts\compare_probes.py samples\synthetic_coarse.csv samples\synthetic_fine.csv --skip-before 0.5
Overlap: 0.5 to 2 s; samples=76
RMSE of |U|: 0.0050828457 m/s
Final paired difference: -0.00256612 m/s

Se ottieni gli stessi numeri sui file sintetici, la tua installazione Python e gli script funzionano: è il controllo da fare prima di applicarli ai dati veri. Ricorda infine che il modulo |U| nasconde i cambi di segno delle componenti: nella cavità, dove il centro del vortice si sposta durante l’avviamento, vale la pena confrontare anche Ux_m_s e Uy_m_s.

Sono dati costruiti artificialmente, non risultati di OpenFOAM: servono solo a verificare la catena estrazione–confronto. Passo di campionamento 0,02 s, 101 righe di dati da 0 a 2 s, solo la componente x non nulla.

samples/synthetic_coarse.csv

time_s,Ux_m_s,Uy_m_s,Uz_m_s,U_magnitude_m_s
0.00000,0.00000000,0,0,0.00000000
0.02000,0.01593750,0,0,0.01593750
0.04000,0.03097326,0,0,0.03097326
0.06000,0.04502978,0,0,0.04502978
0.08000,0.05804139,0,0,0.05804139
0.10000,0.06995605,0,0,0.06995605
0.12000,0.08073670,0,0,0.08073670
0.14000,0.09036220,0,0,0.09036220
0.16000,0.09882777,0,0,0.09882777
0.18000,0.10614507,0,0,0.10614507
0.20000,0.11234175,0,0,0.11234175
0.22000,0.11746075,0,0,0.11746075
0.24000,0.12155910,0,0,0.12155910
0.26000,0.12470654,0,0,0.12470654
0.28000,0.12698378,0,0,0.12698378
0.30000,0.12848062,0,0,0.12848062
0.32000,0.12929389,0,0,0.12929389
0.34000,0.12952532,0,0,0.12952532
0.36000,0.12927931,0,0,0.12927931
0.38000,0.12866080,0,0,0.12866080
0.40000,0.12777315,0,0,0.12777315
0.42000,0.12671618,0,0,0.12671618
0.44000,0.12558431,0,0,0.12558431
0.46000,0.12446498,0,0,0.12446498
0.48000,0.12343723,0,0,0.12343723
0.50000,0.12257055,0,0,0.12257055
0.52000,0.12192404,0,0,0.12192404
0.54000,0.12154575,0,0,0.12154575
0.56000,0.12147239,0,0,0.12147239
0.58000,0.12172931,0,0,0.12172931
0.60000,0.12233063,0,0,0.12233063
0.62000,0.12327978,0,0,0.12327978
0.64000,0.12457009,0,0,0.12457009
0.66000,0.12618575,0,0,0.12618575
0.68000,0.12810274,0,0,0.12810274
0.70000,0.13029006,0,0,0.13029006
0.72000,0.13271095,0,0,0.13271095
0.74000,0.13532422,0,0,0.13532422
0.76000,0.13808557,0,0,0.13808557
0.78000,0.14094891,0,0,0.14094891
0.80000,0.14386769,0,0,0.14386769
0.82000,0.14679605,0,0,0.14679605
0.84000,0.14968999,0,0,0.14968999
0.86000,0.15250834,0,0,0.15250834
0.88000,0.15521364,0,0,0.15521364
0.90000,0.15777284,0,0,0.15777284
0.92000,0.16015784,0,0,0.16015784
0.94000,0.16234590,0,0,0.16234590
0.96000,0.16431984,0,0,0.16431984
0.98000,0.16606813,0,0,0.16606813
1.00000,0.16758478,0,0,0.16758478
1.02000,0.16886910,0,0,0.16886910
1.04000,0.16992538,0,0,0.16992538
1.06000,0.17076239,0,0,0.17076239
1.08000,0.17139285,0,0,0.17139285
1.10000,0.17183280,0,0,0.17183280
1.12000,0.17210092,0,0,0.17210092
1.14000,0.17221787,0,0,0.17221787
1.16000,0.17220552,0,0,0.17220552
1.18000,0.17208633,0,0,0.17208633
1.20000,0.17188265,0,0,0.17188265
1.22000,0.17161613,0,0,0.17161613
1.24000,0.17130719,0,0,0.17130719
1.26000,0.17097450,0,0,0.17097450
1.28000,0.17063464,0,0,0.17063464
1.30000,0.17030179,0,0,0.17030179
1.32000,0.16998752,0,0,0.16998752
1.34000,0.16970070,0,0,0.16970070
1.36000,0.16944746,0,0,0.16944746
1.38000,0.16923126,0,0,0.16923126
1.40000,0.16905306,0,0,0.16905306
1.42000,0.16891153,0,0,0.16891153
1.44000,0.16880334,0,0,0.16880334
1.46000,0.16872345,0,0,0.16872345
1.48000,0.16866556,0,0,0.16866556
1.50000,0.16862241,0,0,0.16862241
1.52000,0.16858623,0,0,0.16858623
1.54000,0.16854915,0,0,0.16854915
1.56000,0.16850352,0,0,0.16850352
1.58000,0.16844234,0,0,0.16844234
1.60000,0.16835953,0,0,0.16835953
1.62000,0.16825020,0,0,0.16825020
1.64000,0.16811090,0,0,0.16811090
1.66000,0.16793976,0,0,0.16793976
1.68000,0.16773660,0,0,0.16773660
1.70000,0.16750296,0,0,0.16750296
1.72000,0.16724206,0,0,0.16724206
1.74000,0.16695877,0,0,0.16695877
1.76000,0.16665941,0,0,0.16665941
1.78000,0.16635160,0,0,0.16635160
1.80000,0.16604399,0,0,0.16604399
1.82000,0.16574604,0,0,0.16574604
1.84000,0.16546768,0,0,0.16546768
1.86000,0.16521902,0,0,0.16521902
1.88000,0.16501002,0,0,0.16501002
1.90000,0.16485018,0,0,0.16485018
1.92000,0.16474822,0,0,0.16474822
1.94000,0.16471181,0,0,0.16471181
1.96000,0.16474729,0,0,0.16474729
1.98000,0.16485947,0,0,0.16485947
2.00000,0.16505143,0,0,0.16505143

samples/synthetic_fine.csv

time_s,Ux_m_s,Uy_m_s,Uz_m_s,U_magnitude_m_s
0.00000,0.00000000,0,0,0.00000000
0.02000,0.01439672,0,0,0.01439672
0.04000,0.02794431,0,0,0.02794431
0.06000,0.04058556,0,0,0.04058556
0.08000,0.05227365,0,0,0.05227365
0.10000,0.06297369,0,0,0.06297369
0.12000,0.07266389,0,0,0.07266389
0.14000,0.08133629,0,0,0.08133629
0.16000,0.08899705,0,0,0.08899705
0.18000,0.09566643,0,0,0.09566643
0.20000,0.10137828,0,0,0.10137828
0.22000,0.10617925,0,0,0.10617925
0.24000,0.11012762,0,0,0.11012762
0.26000,0.11329196,0,0,0.11329196
0.28000,0.11574944,0,0,0.11574944
0.30000,0.11758404,0,0,0.11758404
0.32000,0.11888462,0,0,0.11888462
0.34000,0.11974297,0,0,0.11974297
0.36000,0.12025172,0,0,0.12025172
0.38000,0.12050248,0,0,0.12050248
0.40000,0.12058386,0,0,0.12058386
0.42000,0.12057981,0,0,0.12057981
0.44000,0.12056801,0,0,0.12056801
0.46000,0.12061849,0,0,0.12061849
0.48000,0.12079254,0,0,0.12079254
0.50000,0.12114180,0,0,0.12114180
0.52000,0.12170762,0,0,0.12170762
0.54000,0.12252076,0,0,0.12252076
0.56000,0.12360127,0,0,0.12360127
0.58000,0.12495867,0,0,0.12495867
0.60000,0.12659238,0,0,0.12659238
0.62000,0.12849237,0,0,0.12849237
0.64000,0.13063996,0,0,0.13063996
0.66000,0.13300888,0,0,0.13300888
0.68000,0.13556639,0,0,0.13556639
0.70000,0.13827452,0,0,0.13827452
0.72000,0.14109143,0,0,0.14109143
0.74000,0.14397270,0,0,0.14397270
0.76000,0.14687274,0,0,0.14687274
0.78000,0.14974604,0,0,0.14974604
0.80000,0.15254843,0,0,0.15254843
0.82000,0.15523823,0,0,0.15523823
0.84000,0.15777730,0,0,0.15777730
0.86000,0.16013184,0,0,0.16013184
0.88000,0.16227321,0,0,0.16227321
0.90000,0.16417841,0,0,0.16417841
0.92000,0.16583052,0,0,0.16583052
0.94000,0.16721886,0,0,0.16721886
0.96000,0.16833905,0,0,0.16833905
0.98000,0.16919287,0,0,0.16919287
1.00000,0.16978794,0,0,0.16978794
1.02000,0.17013735,0,0,0.17013735
1.04000,0.17025904,0,0,0.17025904
1.06000,0.17017519,0,0,0.17017519
1.08000,0.16991144,0,0,0.16991144
1.10000,0.16949613,0,0,0.16949613
1.12000,0.16895943,0,0,0.16895943
1.14000,0.16833254,0,0,0.16833254
1.16000,0.16764679,0,0,0.16764679
1.18000,0.16693289,0,0,0.16693289
1.20000,0.16622018,0,0,0.16622018
1.22000,0.16553591,0,0,0.16553591
1.24000,0.16490466,0,0,0.16490466
1.26000,0.16434784,0,0,0.16434784
1.28000,0.16388332,0,0,0.16388332
1.30000,0.16352508,0,0,0.16352508
1.32000,0.16328309,0,0,0.16328309
1.34000,0.16316325,0,0,0.16316325
1.36000,0.16316742,0,0,0.16316742
1.38000,0.16329358,0,0,0.16329358
1.40000,0.16353610,0,0,0.16353610
1.42000,0.16388607,0,0,0.16388607
1.44000,0.16433170,0,0,0.16433170
1.46000,0.16485879,0,0,0.16485879
1.48000,0.16545126,0,0,0.16545126
1.50000,0.16609169,0,0,0.16609169
1.52000,0.16676184,0,0,0.16676184
1.54000,0.16744326,0,0,0.16744326
1.56000,0.16811776,0,0,0.16811776
1.58000,0.16876799,0,0,0.16876799
1.60000,0.16937782,0,0,0.16937782
1.62000,0.16993282,0,0,0.16993282
1.64000,0.17042057,0,0,0.17042057
1.66000,0.17083095,0,0,0.17083095
1.68000,0.17115634,0,0,0.17115634
1.70000,0.17139174,0,0,0.17139174
1.72000,0.17153484,0,0,0.17153484
1.74000,0.17158600,0,0,0.17158600
1.76000,0.17154811,0,0,0.17154811
1.78000,0.17142650,0,0,0.17142650
1.80000,0.17122864,0,0,0.17122864
1.82000,0.17096392,0,0,0.17096392
1.84000,0.17064332,0,0,0.17064332
1.86000,0.17027904,0,0,0.17027904
1.88000,0.16988416,0,0,0.16988416
1.90000,0.16947223,0,0,0.16947223
1.92000,0.16905689,0,0,0.16905689
1.94000,0.16865153,0,0,0.16865153
1.96000,0.16826889,0,0,0.16826889
1.98000,0.16792072,0,0,0.16792072
2.00000,0.16761755,0,0,0.16761755

Con tre livelli a rapporto costante r = 2 (Δt = 0,005; 0,0025; 0,00125 s) e i valori f₃ (grossolano), f₂, f₁ (fine) della stessa grandezza allo stesso tempo fisico:

p_obs  = ln[ (f₃ − f₂) / (f₂ − f₁) ] / ln r
f_ext  = f₁ + (f₁ − f₂) / (r^p_obs − 1)
errore stimato del caso fine ≈ |f₁ − f_ext|

per Euler nel campo asintotico p_obs ≈ 1; se (f₃ − f₂) e (f₂ − f₁) hanno segno diverso
la convergenza non è monotona e l'estrapolazione non va usata

Lo script seguente, solo libreria standard, applica queste formule ai tre CSV prodotti da extract_probes.py agli istanti scelti:

#!/usr/bin/env python3
"""Three-level temporal refinement check on probe CSV files (time_s, U_magnitude_m_s).

Usage: richardson_time.py coarse.csv medium.csv fine.csv [--ratio 2] [--times 0.5 1.0 1.5 2.0]
Prints, at each requested time, the observed order p and a Richardson estimate.
Standard library only. License: GNU GPL v3 or later.
"""
from __future__ import annotations
import argparse, bisect, csv, math
from pathlib import Path


def load(path: Path) -> list[tuple[float, float]]:
    with path.open(newline="", encoding="utf-8") as f:
        pts = [(float(r["time_s"]), float(r["U_magnitude_m_s"])) for r in csv.DictReader(f)]
    if len(pts) < 2 or any(b[0] <= a[0] for a, b in zip(pts, pts[1:])):
        raise ValueError(f"Need at least two strictly increasing times in {path}")
    return pts


def at(pts: list[tuple[float, float]], t: float) -> float:
    ts = [a for a, _ in pts]
    i = bisect.bisect_left(ts, t - 1e-11)
    if i < len(pts) and abs(pts[i][0] - t) < 1e-11:
        return pts[i][1]
    if i == 0 or i == len(pts):
        raise ValueError(f"time {t} outside the series (no extrapolation)")
    (t0, y0), (t1, y1) = pts[i - 1], pts[i]
    return y0 + (y1 - y0) * (t - t0) / (t1 - t0)


def main() -> None:
    ap = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter)
    ap.add_argument("coarse", type=Path); ap.add_argument("medium", type=Path); ap.add_argument("fine", type=Path)
    ap.add_argument("--ratio", type=float, default=2.0, help="time-step ratio between levels")
    ap.add_argument("--times", type=float, nargs="+", default=[0.5, 1.0, 1.5, 2.0])
    a = ap.parse_args()
    c, m, f = load(a.coarse), load(a.medium), load(a.fine)
    r = a.ratio
    print(f"{'t [s]':>8} {'f_coarse':>12} {'f_medium':>12} {'f_fine':>12} {'p_obs':>7} {'f_extrap':>12} {'err_fine %':>10}")
    for t in a.times:
        f3, f2, f1 = at(c, t), at(m, t), at(f, t)
        e32, e21 = f3 - f2, f2 - f1
        if e21 == 0 or e32 == 0 or (e32 / e21) <= 0:
            print(f"{t:8.4g} {f3:12.6g} {f2:12.6g} {f1:12.6g}   non-monotone: no order/extrapolation")
            continue
        p = math.log(e32 / e21) / math.log(r)
        fx = f1 + (f1 - f2) / (r ** p - 1)
        err = 100 * abs(f1 - fx) / abs(fx) if fx else float("nan")
        print(f"{t:8.4g} {f3:12.6g} {f2:12.6g} {f1:12.6g} {p:7.3f} {fx:12.6g} {err:10.3g}")
    print("p close to 1 is expected for Euler in the asymptotic range; far from 1 means the levels are not asymptotic.")


if __name__ == "__main__":
    main()

Controllo dello script su un problema con soluzione nota: l’equazione du/dt = (1 − u)/τ con τ = 0,3 s, integrata con Euler implicito ai tre passi, ha soluzione esatta u(0,5) = 0,811124. Lo script restituisce:

   t [s]     f_coarse     f_medium       f_fine   p_obs     f_extrap err_fine %
     0.5     0.808512     0.809815     0.810469   0.996     0.811127     0.0811
       1     0.963332      0.96383     0.964078   1.003     0.964325     0.0256
     1.5     0.992979     0.993121     0.993192   1.010     0.993261    0.00702
       2     0.998655     0.998692      0.99871   1.018     0.998727    0.00175

L’ordine osservato è prossimo a 1 come previsto per Euler e il valore estrapolato coincide con quello esatto alla quarta cifra significativa. Sui dati CFD il comportamento sarà meno pulito: errori iterativi e spaziali si sommano a quello temporale, ed è per questo che la guida chiede di tenere fissi mesh e modello.

I CSV usano il punto come separatore decimale e la virgola come separatore di colonna. Con Excel in italiano un doppio clic può trasformare i numeri in testo o in date. Usa invece Dati → Da testo/CSV, scegli «Delimitatore: virgola» e, in Trasforma dati, imposta le colonne come numeri decimali con impostazioni locali «Inglese (Stati Uniti)». Non salvare poi sopra il CSV originale: salva un file .xlsx separato.

# 1) sul server Linux (sessione SSH, meglio dentro tmux)
cd "$HOME/OpenFOAM/guide08_run01/cavity_dt005"
touch cavity_dt005.foam
pvserver --server-port=11111

# 2) su Windows, in PowerShell: tunnel SSH verso la porta del server
ssh -L 11111:localhost:11111 utente@server.example.org

# 3) in ParaView Windows: File → Connect → localhost:11111,
#    poi File → Open → cavity_dt005.foam (file sul server)

Il client Windows e pvserver devono avere la stessa versione di ParaView: con versioni diverse la connessione viene rifiutata o si comporta in modo imprevedibile. Il file .foam è vuoto: serve solo a far riconoscere la cartella al lettore OpenFOAM.

• Scala fissa: nell’editor della mappa colori usa Rescale to Custom Range (per esempio 0–1 m/s, la velocità del coperchio) e imposta Automatic Rescale Range Mode su Never; altrimenti ParaView può riadattare l’intervallo durante la riproduzione;
• tempo in scena: il filtro Annotate Time Filter applicato al lettore mostra il tempo fisico del dato; formattalo con poche cifre e l’unità, per esempio «t = {time:.2f} s»;
• modalità: Snap To TimeSteps mostra solo istanti salvati; Sequence genera un numero di fotogrammi a scelta, e se sono più degli istanti salvati alcuni fotogrammi ripetono lo stesso dato; Real Time lega la durata al tempo di orologio e non serve per l’analisi.

Le impostazioni fatte con il mouse si perdono facilmente; uno script le fissa. Schema di partenza da adattare (i nomi delle proprietà possono variare fra versioni: il modo più sicuro per ottenere quelli esatti della tua è Tools → Start Trace, ripetere le operazioni a mano e copiare il codice generato):

# pvpython anim_cavity.py   (oppure Python Shell di ParaView)
from paraview.simple import *

case = "cavity_dt005"
r = OpenFOAMReader(FileName=case + "/" + case + ".foam")
r.MeshRegions = ["internalMesh"]
r.CellArrays = ["U"]
r.UpdatePipelineInformation()

v = GetActiveViewOrCreate("RenderView")
d = Show(r, v)
ColorBy(d, ("POINTS", "U", "Magnitude"))
lut = GetColorTransferFunction("U")
lut.AutomaticRescaleRangeMode = "Never"
lut.RescaleTransferFunction(0.0, 1.0)          # scala fissa 0-1 m/s
d.SetScalarBarVisibility(v, True)

t = AnnotateTimeFilter(r)
t.Format = "t = {time:.2f} s"
Show(t, v)

v.ViewSize = [1280, 720]
v.CameraPosition = [0.05, 0.05, 0.3]
v.CameraFocalPoint = [0.05, 0.05, 0.005]
v.CameraViewUp = [0, 1, 0]

scene = GetAnimationScene()
scene.PlayMode = "Snap To TimeSteps"
SaveAnimation("frames/" + case + ".png", v, ImageResolution=[1280, 720])

Usando lo stesso script per i due casi, cambiando solo case, camera, scala e risoluzione restano identiche per costruzione: è il modo più affidabile di soddisfare l’esercizio C. Crea prima la cartella frames e controlla che non contenga immagini di prove precedenti.

ffmpeg -framerate 10 -i frames/cavity_dt005.%04d.png -c:v libx264 -pix_fmt yuv420p -crf 18 cavity_dt005_10fps.mp4

La sequenza PNG è l’originale senza perdita; il video è un derivato. -framerate stabilisce gli FPS, yuv420p garantisce la compatibilità con i lettori comuni e con YouTube, -crf 18 dà una qualità visivamente quasi senza perdita. Le dimensioni dell’immagine devono essere pari (1280 × 720 va bene). Il nome dei file numerati prodotti da ParaView va controllato con ls frames | head e il modello %04d adattato di conseguenza. Metti gli FPS nel nome del file: chi riceve il video saprà come leggerne il tempo.

Se nell’esercizio B il rapporto delle differenze è lontano da 2, prima di trarre conclusioni prova a stringere le tolleranze dei solutori lineari in fvSolution in tutti e tre i casi: se il rapporto cambia, l’errore iterativo stava mascherando quello temporale.

R="$HOME/OpenFOAM/guide08_run01"
N="$HOME/OpenFOAM/guide08_run02/cavity_dt00125"
if [ -e "$N" ]; then
  echo "Esiste già: scegli un altro nome"
else
  mkdir -p "$(dirname "$N")"
  cp -r "$FOAM_TUTORIALS/incompressibleFluid/cavity" "$N"
  cp "$R/cavity_dt0025/system/controlDict" "$N/system/controlDict"   # con la sonda
  foamDictionary "$N/system/controlDict" -entry deltaT -set 0.00125
  foamDictionary "$N/system/controlDict" -entry writeInterval -set 40
  foamDictionary "$N/system/controlDict" -entry deltaT -value
  foamDictionary "$N/system/controlDict" -entry writeInterval -value
  printf "Source: tutorial cavity + controlDict of cavity_dt0025\ndeltaT: 0.00125 s\nwriteInterval: 40 steps (=0.05 s)\n" > "$N/GUIDE08_METADATA.txt"
fi

# esecuzione con gli stessi passi del runner
( cd "$N" && blockMesh > log.blockMesh 2>&1 && checkMesh > log.checkMesh 2>&1 && foamRun > log.foamRun 2>&1 )

Si parte dalla copia ufficiale del tutorial e si prende il controlDict già preparato, che contiene la sonda: così la configurazione differisce dagli altri casi solo per le due voci modificate. Una cartella separata (guide08_run02) tiene distinta la nuova campagna da quella dell’esercizio A.

Caso:            cavity_dt005 (OpenFOAM Foundation 14, tutorial incompressibleFluid/cavity)
Grandezza:       |U| [m/s], scala fissa 0–1
Tempo fisico:    0–2 s, istanti ogni 0,05 s (41 fotogrammi)
Video:           10 FPS → 1 s di video = 0,5 s di tempo fisico
Camera:          vista frontale, piano z = 0,005 m
ParaView:        5.13.2 (client e pvserver)
Stato:           anim_cavity.pvsm
Note:            dati da sonda in cavity_coarse.csv

Un video senza questa scheda è un’immagine suggestiva; con la scheda diventa un documento tecnico che un collega può riprodurre e verificare.

La lista di verifica diventa pratica se la si compila mentre si lavora, non a posteriori. Un modello da copiare in ogni cartella di caso come RUN_NOTES.txt:

OpenFOAM:         (output di foamVersion)
Tutorial:         $FOAM_TUTORIALS/incompressibleFluid/cavity
Mesh:             20 x 20 x 1, checkMesh: (esito)
Modello:          RAS kEpsilon, nu = 1e-05 m2/s
Schemi:           ddt Euler; (altre modifiche: nessuna)
Tempo:            0-2 s, deltaT = ..., writeInterval = ...
Courant:          max ... a t = ... (log_summary.py)
Sonda:            (0.05 0.05 0.005), file postProcessing/centreVelocity/0/U
Confronto:        RMSE = ... m/s, differenza finale = ... m/s
ParaView:         5.13.2, stato ....pvsm
Modifiche in corso di calcolo: ...
Data e autore:    ...

# dalla cartella che contiene guide08_run01 (nessun file viene rimosso)
foamVersion > guide08_run01/OPENFOAM_VERSION.txt 2>&1
tar czf guide08_run01_$(date +%Y%m%d).tar.gz guide08_run01
sha256sum guide08_run01_*.tar.gz > guide08_run01_SHA256.txt
tar tzf guide08_run01_*.tar.gz | head

L’archivio compresso conserva log, sonde, campi e metadati in un unico file; l’impronta SHA-256 permette, mesi dopo, di dimostrare che i dati non sono cambiati. Per l’archivio da condividere spesso bastano log, postProcessing, system, constant (senza polyMesh, rigenerabile con blockMesh), lo stato di ParaView e i CSV: chi lo riceve può rigenerare tutto il resto.

• la guida utente di OpenFOAM su CFD Direct, citata fra i riferimenti, per la sintassi di controlDict e delle function object della propria versione;
• J. H. Ferziger, M. Perić, R. L. Street, Computational Methods for Fluid Dynamics, Springer: capitoli sull’integrazione temporale e sugli algoritmi pressione–velocità;
• P. J. Roache, Verification and Validation in Computational Science and Engineering: base metodologica per ordine osservato, Richardson e GCI, ripresa nella Guida 4;
• la documentazione di ParaView 5.13 sull’animazione, citata fra i riferimenti, per le opzioni di Save Animation della versione installata.