CFDManager/docs/theory/solver_2x_sdf/visualization.py

# Визуализация (host): гифка |ω| с наложением тонкого патча L1 и прогресс-бар.
# matplotlib/PIL работают только с host-массивами — кадры уже перенесены на host в solver.run().
import os
import time as _time
import numpy as np

def make_progress(label, width=32, every=0.1):
    """ASCII-шкала прогресса по выполненным шагам из заданных (обновление на месте через \\r)."""
    st = {"t0": _time.perf_counter(), "last": 0.0}
    def cb(t, total):
        now = _time.perf_counter(); done = (t >= total)
        if not done and (now - st["last"] < every):
            return
        st["last"] = now; frac = (t + 1) / (total + 1)
        filled = int(round(width * frac)); bar = "#"*filled + "-"*(width - filled)
        el = now - st["t0"]; eta = (el/frac - el) if frac > 0 else 0.0
        print(f"\r  {label:8}[{bar}] {frac*100:5.1f}% ({t}/{total})  {el:5.1f}s  ETA {eta:5.1f}s",
              end=("\n" if done else ""), flush=True)
    return cb

def _vort(u):
    return (np.roll(u[1], -1, 1) - np.roll(u[1], 1, 1))*0.5 - (np.roll(u[0], -1, 0) - np.roll(u[0], 1, 0))*0.5

def save_series_panel(res, a, cfg, path):
    """Готовая png-панель для ноутбука (он код не исполняет — только вставляет картинки):
    (1) ряды Cd(t)/Cl(t) + редкие точки кросс-чека ∮σ·n ds; (2) Cm(t) — контроль симметрии;
    (3) спектр зонда u_y с пиком St; (4) ⟨ρ⟩(t) — контроль дрейфа массы."""
    import backend as B
    import matplotlib; matplotlib.use("Agg")
    import matplotlib.pyplot as plt
    Cd = a["Cd_s"]; Cl = a["Cl_s"]; n = len(Cd); h0 = a["h0"]; t = np.arange(n)
    fig, axs = plt.subplots(2, 2, figsize=(12.5, 7.2), dpi=110)
    ax = axs[0, 0]
    ax.plot(t, Cd, lw=0.7, color="#d32f2f", label=f"Cd (GMEM), ⟨Cd⟩={a['Cd']:.3f}")
    ax.plot(t, Cl, lw=0.7, color="#1976d2", label=f"Cl (GMEM), rms={a['Clrms']:.3f}")
    if "Cd_st_s" in a:
        K = a["stress_every"]; ts = np.arange(len(a["Cd_st_s"])) * K
        ax.plot(ts, a["Cd_st_s"], ".", ms=3, color="#ff8f00", label=f"Cd ∮σ·n, ⟨⟩={a['Cd_st']:.3f}")
        ax.plot(ts, a["Cl_st_s"], ".", ms=3, color="#4fc3f7", label="Cl ∮σ·n")
    ax.axvline(h0, color="gray", lw=0.8, ls="--")
    ax.set_xlabel("шаг L0"); ax.set_ylabel("Cd, Cl"); ax.legend(fontsize=8, loc="upper right")
    ax.set_title("Сила: GMEM (линии) vs ∮σ·n ds (точки)", fontsize=10)
    ax = axs[0, 1]
    if "Cm_s" in a:
        ax.plot(t, a["Cm_s"], lw=0.7, color="#388e3c")
        ax.axhline(0.0, color="gray", lw=0.8)
        ax.set_title(f"Cm(t): ⟨Cm⟩={a['Cm']:.5f} (контроль симметрии, ожид. ≈0)", fontsize=10)
    ax.set_xlabel("шаг L0"); ax.set_ylabel("Cm")
    ax = axs[1, 0]
    ax.plot(a["freqs"]*cfg.D/cfg.U, a["amp"], lw=0.9, color="#7b1fa2")
    ax.axvline(a["St"], color="#d32f2f", lw=1.0, ls="--", label=f"St={a['St']:.3f}")
    ax.set_xlim(0, 0.6); ax.set_xlabel("St = f·D/U"); ax.set_ylabel("|FFT(u_y)|")
    ax.legend(fontsize=9); ax.set_title("Спектр зонда в следе (вторая половина ряда)", fontsize=10)
    ax = axs[1, 1]
    rho = B.to_cpu(res["rho"]) if not isinstance(res["rho"], np.ndarray) else res["rho"]
    ax.plot(np.arange(len(rho)), rho, lw=0.9, color="#5d4037")
    ax.axhline(1.0, color="gray", lw=0.8)
    ax.set_xlabel("шаг L0"); ax.set_ylabel("⟨ρ⟩")
    ax.set_title(f"⟨ρ⟩(t): дрейф массы (финально {float(rho[-1]):.4f})", fontsize=10)
    fig.suptitle(f"2×+SDF · Ny={cfg.Ny} (β={cfg.D/cfg.Ny:.3f}) · Re={cfg.Re:g} · {len(rho)-1} шагов",
                 fontsize=11)
    fig.tight_layout(rect=(0, 0, 1, 0.96))
    os.makedirs(os.path.dirname(path), exist_ok=True)
    fig.savefig(path, bbox_inches="tight"); plt.close(fig)
    print("saved", os.path.basename(path))

def save_gif(res, geo, cfg, path, fps=24):
    import matplotlib; matplotlib.use("Agg")
    import matplotlib.pyplot as plt, matplotlib.patches as mpatches
    from PIL import Image
    s0 = geo.solid0_np; s1 = geo.solid1_np
    Nx, Ny = cfg.Nx, cfg.Ny; ax1, bx1, ay1, by1 = cfg.ax1, cfg.bx1, cfg.ay1, cfg.by1
    frames = res["frames"]
    vm = np.nanpercentile(np.abs(np.where(s0, np.nan, _vort(frames[len(frames)//2][1]))), 99.0)
    cm = plt.get_cmap("inferno").copy(); cm.set_bad("white"); pil = []
    for (t, u0, u1) in frames:
        fig = plt.figure(figsize=(11.2, 5.7), dpi=96); ax = fig.add_subplot(111)
        ax.imshow(np.abs(np.where(s0, np.nan, _vort(u0))), origin="lower", cmap=cm, vmin=0, vmax=vm,
                  extent=(0, Nx, 0, Ny), interpolation="nearest")
        ax.imshow(np.abs(np.where(s1, np.nan, _vort(u1)))[1:-1, 1:-1], origin="lower", cmap=cm, vmin=0, vmax=vm,
                  extent=(ax1+0.5, bx1-0.5, ay1+0.5, by1-0.5), interpolation="nearest")
        ax.add_patch(mpatches.Rectangle((0.2, 0.2), Nx-0.4, Ny-0.4, fill=False, edgecolor="#4fc3f7", lw=1.6))
        ax.text(1.5, Ny-4, "L0: Δx", color="#4fc3f7", fontsize=9, weight="bold")
        ax.add_patch(mpatches.Rectangle((ax1, ay1), bx1-ax1, by1-ay1, fill=False, edgecolor="#aed581", lw=2.2))
        ax.text(ax1+1, by1-3.5, "L1: Δx/2 (тело+след)", color="#aed581", fontsize=9, weight="bold")
        ax.set_title(f"2×+SDF · |\\omega| · t={t}", fontsize=11)
        ax.set_xlabel("x [lu]"); ax.set_ylabel("y [lu]"); ax.set_xlim(0, Nx); ax.set_ylim(0, Ny)
        buf = __import__("io").BytesIO(); fig.savefig(buf, format="png", bbox_inches="tight"); buf.seek(0); plt.close(fig)
        pil.append(Image.open(buf).convert("RGB").convert("P", palette=Image.ADAPTIVE))
    os.makedirs(os.path.dirname(path), exist_ok=True)
    pil[0].save(path, save_all=True, append_images=pil[1:], duration=int(1000/fps), loop=0, optimize=True)
    print("saved", os.path.basename(path), len(frames), "кадров")