commit 792daa6feb884aa089f1a41c12a14a2e205c1459 Author: admin Date: Fri Sep 18 00:05:18 2026 +0800 初始化项目：规划文档 + FFU 洁净室流场实时仿真 Demo - 规划文档 01-07：愿景 / 架构 / 路线图 / 数据与AI / 风险 / 参考 / FFU Demo - 演示图：系统架构、数据流水线、里程碑路线图（SVG + PNG） - FFU Demo：D2Q9 LBM 高保真（360x270, 76s）+ 粗网格快速代理（120x90, 7.3s, 16-27x 加速） - 趋势验证：速度 r=0.95 / 粒子 r=0.83 / 温度 r=0.68 - 工具链：仿真脚本、快速代理、图表生成器、SVG->PNG 转换 diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..a545909 --- /dev/null +++ b/.gitignore @@ -0,0 +1,18 @@ +# Python +__pycache__/ +*.pyc +.venv/ +venv/ + +# 仿真输出（可由脚本重新生成） +docs/assets/*.npz + +# Node +node_modules/ +package-lock.json + +# 编辑器 / 系统 +.vscode/ +.idea/ +.DS_Store +Thumbs.db diff --git a/README.md b/README.md new file mode 100644 index 0000000..d3cf9b9 --- /dev/null +++ b/README.md @@ -0,0 +1,118 @@ +# PhysicsAiFluentAnalysis + +> 物理世界流场实时仿真 —— 用 AI 神经算子把高保真 CFD 压缩到毫秒级推理。 + +[![Status](https://img.shields.io/badge/status-P0%20%E8%A7%84%E5%88%92%E6%9C%9F-blue)]() +[![Demo](https://img.shields.io/badge/demo-FFU%20%E6%B4%81%E5%87%80%E5%AE%A4%E6%B5%81%E5%9C%BA-green)]() +[![Speed](https://img.shields.io/badge/surrogate-7.3%20s%20%2F%2016--27x-orange)]() + +--- + +## 一句话定义 + +以 **Ansys Fluent 高保真仿真为真值**，训练 **物理约束的 AI 代理模型（神经算子）**，在保持物理一致性的前提下，把小时级的流场计算压缩到 **毫秒级实时推理**，支撑数字孪生、工程设计优化与沉浸式可视化。 + +## 为什么这件事"炸裂"又"难" + +| 维度 | 现状 | 目标 | +| --- | --- | --- | +| 单次 CFD 仿真 | 小时 ~ 天级 | 毫秒级（>1000x 加速） | +| 参数扫描 | 一次一个工况 | 全参数空间连续可微 | +| 交互性 | 离线出报告 | 实时拖拽参数、流场即时响应 | +| 物理可信度 | 求解器保证 | 物理约束 + 守恒律 + 不确定性量化 | + +**难在哪：** Navier-Stokes 方程的非线性与多尺度特性、湍流闭合问题、AI 长时推演的误差累积与发散、复杂几何的泛化、3D 大网格的实时推理算力瓶颈。这不是"套一个模型"能解决的，需要 **数据流水线 + 物理约束 + 工程加速 + 严格验证** 四线并进。 + +## 已落地 Demo：FFU 洁净室流场实时仿真 + +FFU（风机过滤单元）下方洁净区的流场、温湿度、压力与尘埃粒子输运 —— 同场景提供 **高保真真值** 与 **快速代理** 两条路径。 + +![FFU 流场（高保真）](docs/assets/ffu-flow-field.png) + +| | 高保真真值 | 快速代理 | +| --- | --- | --- | +| 网格 / 步数 | 360×270 / 14000 步 | 120×90 / 5000 步（同源求解器） | +| 端到端耗时 | **76 s** | **7.3 s（16-27x 加速）** | +| 趋势一致性 | 基准 | 速度 **r=0.95** · 粒子 r=0.83 · 温度 r=0.68 | + +![趋势验证](docs/assets/ffu-validation.png) + +**主要趋势**：单向下扫排污染 → 绕障/缝隙加速（局部 ~1.5 m/s）→ 背风侧死角堆积滞留 → 送风正压/回风负压（< 2 Pa）。 + +详见 [docs/07-ffu-demo.md](docs/07-ffu-demo.md)。 + +## 系统总览 + +![系统架构](docs/assets/architecture.png) + +## 技术路线（一图流） + +![数据与训练流水线](docs/assets/pipeline.png) + +## 里程碑路线图 + +![路线图](docs/assets/roadmap.png) + +## 文档导航 + +| 文档 | 内容 | +| --- | --- | +| [docs/01-vision.md](docs/01-vision.md) | 愿景、痛点、目标与非目标、成功标准（KPI） | +| [docs/02-architecture.md](docs/02-architecture.md) | 分层架构、数据流、技术选型与对比 | +| [docs/03-roadmap.md](docs/03-roadmap.md) | 7 个阶段的目标、交付物与验收标准 | +| [docs/04-data-and-ai.md](docs/04-data-and-ai.md) | Fluent 数据工厂、模型选型、训练与评估方案 | +| [docs/05-risks.md](docs/05-risks.md) | 技术风险、缓解策略与开放问题 | +| [docs/06-references.md](docs/06-references.md) | 论文、开源工具与数据集索引 | +| [docs/07-ffu-demo.md](docs/07-ffu-demo.md) | FFU 洁净室 Demo：物理设置、趋势验证、复现命令 | + +## 仓库结构 + +``` +PhysicsAiFluentAnalysis/ +├── README.md # 项目总览（本文件） +├── package.json # 图件工具依赖 (SVG -> PNG) +├── docs/ # 规划文档 + 演示图 +│ ├── 01-vision.md ~ 07-ffu-demo.md +│ └── assets/ # PNG / SVG 演示图 +│ ├── architecture.* 系统架构图 +│ ├── pipeline.* 数据与训练流水线 +│ ├── roadmap.* 里程碑路线图 +│ ├── ffu-flow-field.* FFU 流场（高保真 LBM） +│ ├── ffu-scalar-fields.png 温湿度/压力/粒子（高保真） +│ ├── ffu-particle-tracks.png 粒子轨迹与滞留时间 +│ ├── ffu-flow-field-fast.png 快速代理流场（7 s） +│ ├── ffu-scalar-fields-fast.png 快速代理标量场 +│ └── ffu-validation.png 趋势验证对比图 +└── tools/ # 可复现的生成脚本 + ├── ffu_simulation.py # 高保真 D2Q9 LBM + 标量输运（76 s） + ├── ffu_surrogate.py # 快速代理（7 s）+ 趋势验证 + ├── gen_diagrams.py # 架构/流水线/路线图 SVG 生成 + └── svg2png.js # SVG -> PNG 批量转换 +``` + +## 快速开始 + +```bash +# 依赖: Python 3.11+ (numpy / numba / matplotlib), Node 18+ (仅图件转换需要) + +# 1) 快速代理出图（约 7 秒，含趋势验证对比图） +python tools/ffu_surrogate.py + +# 2) 高保真真值（约 76 秒，输出 3 张图 + npz 场数据） +python tools/ffu_simulation.py + +# 3) 规划图（架构 / 流水线 / 路线图） +python tools/gen_diagrams.py + +# 4) SVG -> PNG +npm install && node tools/svg2png.js +``` + +## 当前状态 + +**P0 规划期**：愿景、架构、路线图、风险分析完成；FFU 洁净室 Demo 已可运行（高保真 + 快速代理双路径 + 趋势验证）。 +**下一步**：进入 P1，完成 2D 圆柱绕流 POC（Fluent 数据 → FNO 基线 → 误差评估）。 + +--- + +*本项目处于早期规划阶段，欢迎讨论与共建。* diff --git a/docs/01-vision.md b/docs/01-vision.md new file mode 100644 index 0000000..afea792 --- /dev/null +++ b/docs/01-vision.md @@ -0,0 +1,85 @@ +# 01 · 愿景与目标 + +## 1. 背景 + +计算流体力学（CFD）是现代工程的核心工具：从飞行器气动设计、汽车风阻优化，到血管血流分析、城市风环境评估，都依赖数值求解 Navier-Stokes 方程。Ansys Fluent 是工业界最成熟的 CFD 求解器之一，提供从 RANS 到 LES 的完整湍流建模能力。 + +但 CFD 有一个根本性的效率瓶颈： + +- **单次高保真仿真成本高**：一个复杂 3D 工况动辄数小时到数天； +- **参数扫描是组合爆炸**：设计空间中每换一组参数（雷诺数、攻角、几何、边界条件）就要重算一遍； +- **无法交互**：工程设计评审时无法"拖动参数、流场即时响应"。 + +与此同时，AI for Science 在过去五年出现了范式级突破：**神经算子（Neural Operator）** 能直接学习"函数 → 函数"的映射，一次训练即可在参数空间中连续插值/外推，推理成本只与网络前向相当。这为"实时流场"提供了物理上可信、工程上可行的新路径。 + +## 2. 愿景 + +> **让流场仿真像游戏一样实时、像求解器一样可信。** + +构建一条从 **高保真仿真 → AI 代理模型 → 实时交互系统** 的完整闭环链路，使工程师可以在秒级完成过去需要数天的设计迭代，并在数字孪生场景中实现"传感器数据驱动 + 物理模型约束"的在线流场重建。 + +## 3. 目标 + +### 3.1 技术目标 + +| 编号 | 目标 | 量化指标 | +| --- | --- | --- | +| G1 | 训练数据自动化生产 | PyFluent 批量仿真流水线，无人值守生成 ≥ 10⁴ 工况 | +| G2 | 2D 基准问题精度 | 圆柱绕流速度场相对 L2 误差 < 5% | +| G3 | 3D 复杂几何精度 | 汽车/翼型外流场相对 L2 误差 < 8%，关键气动力系数误差 < 5% | +| G4 | 实时推理 | 单步推理延迟 < 10 ms（3D 中等规模网格），较 Fluent 加速 > 1000x | +| G5 | 长时稳定 | 自回归推演 1000+ 步不发散，质量/动量守恒误差 < 1% | +| G6 | 泛化能力 | 未见雷诺数/几何族上误差退化 < 2 倍 | +| G7 | 可信度 | 输出带不确定性量化，OOD 检测召回率 > 90% | + +### 3.2 工程目标 + +- 提供可交互的实时可视化前端（WebGPU/WebGL），支持参数拖拽与流场即时更新； +- 推理服务容器化，支持边缘/云端两种部署形态； +- 全链路可复现：数据版本、模型版本、指标看板。 + +### 3.3 应用目标 + +| 场景 | 价值 | +| --- | --- | +| 数字孪生 | 传感器稀疏观测 + AI 流场重建，在线监测与预警 | +| 工程设计 | 实时参数扫描，设计空间探索效率提升 100x | +| 教学科研 | 流场"可触摸"，物理直觉可视化 | +| 实时决策 | 风场/流场的秒级预测，支撑调度与控制 | + +## 4. 非目标（明确不做） + +- **不替代 Fluent 的求解精度**：本项目定位是"代理模型 + 加速器"，最终验证与认证仍由传统 CFD 完成； +- **不做通用流体求解器**：聚焦特定问题族（外流场/内流场），不追求覆盖所有流体现象； +- **不做实时全分辨率 3D LES**：以潜空间/降阶 + 代理模型实现"感知级实时"，而非"DNS 级实时"； +- **首期不做多物理场耦合**：传热、燃烧、多相流为远期扩展。 + +## 5. 成功标准（阶段性 KPI） + +**P1 结束（3 个月）——可行性验证通过** + +- [ ] 完成 2D 圆柱绕流数据集（≥ 500 工况，Re ∈ [50, 300]） +- [ ] FNO 基线模型相对 L2 误差 < 10% +- [ ] 单步推理 < 5 ms（2D，256² 网格） +- [ ] 评估报告：精度 / 稳定性 / 泛化三维度量化 + +**P2 结束（6 个月）——系统闭环** + +- [ ] 端到端流水线自动化（采样 → 仿真 → 训练 → 部署） +- [ ] 3D 翼型/简化车体代理模型相对 L2 误差 < 10% +- [ ] 实时可视化 Demo：参数拖拽 → 流场更新延迟 < 100 ms + +**P3 结束（12 个月）——数字孪生原型** + +- [ ] 稀疏观测数据同化在线更新 +- [ ] 不确定性量化 + OOD 检测上线 +- [ ] 与 Fluent 真值在 ≥ 3 个工业场景对标验证 + +## 6. 判断准则（Decision Gates） + +每个阶段结束设置 Go / No-Go 决策点： + +1. **精度门禁**：误差指标未达标则不允许进入下一阶段，先做消融分析； +2. **物理门禁**：守恒律破坏或长时发散直接视为失败； +3. **成本门禁**：单工况训练数据成本需低于设定阈值（否则先优化数据策略）； +4. **可复现门禁**：任何指标必须可由流水线一键复现。 diff --git a/docs/02-architecture.md b/docs/02-architecture.md new file mode 100644 index 0000000..c3d711e --- /dev/null +++ b/docs/02-architecture.md @@ -0,0 +1,120 @@ +# 02 · 系统架构与技术选型 + +## 1. 分层架构 + +![系统架构](assets/architecture.png) + +系统分为六层，自下而上： + +### L0 高保真真值层（Ground Truth） + +| 组件 | 说明 | +| --- | --- | +| Ansys Fluent | 高保真求解器：RANS（k-ε / k-ω SST）、URANS、LES/DES | +| PyFluent | Python 驱动 Fluent，实现参数化建模、求解、后处理全自动化 | +| 参数化几何 | 基于参数脚本的几何生成（翼型族、车体变体、管道变径） | +| 边界条件库 | 入口速度/湍流度、出口压力、壁面函数等工况模板 | + +### L1 数据与训练层 + +| 组件 | 说明 | +| --- | --- | +| 数据湖 | HDF5（训练样本）/ VTK（可视化）/ Zarr（云原生分块） | +| 数据工程 | 网格插值统一、无量纲化、场归一化、增广（旋转/镜像/缩放） | +| 训练框架 | PyTorch + AMP + DDP/FSDP，支持多卡与梯度检查点 | +| 主动学习 | 基于预测不确定性挑选高价值工况回灌 Fluent 重算 | + +### L2 AI 模型层 + +| 模型 | 定位 | 适用场景 | +| --- | --- | --- | +| FNO（傅里叶神经算子） | 主力基线 | 规则网格、2D/3D 流场、分辨率无关 | +| DeepONet | 算子学习备选 | 低维参数空间、传感器布局优化 | +| PINN | 物理约束补充 | 少数据/无数据工况、反问题 | +| MeshGraphNet | 非结构网格 | 复杂几何、工业网格 | +| Transformer 算子 | 长程依赖 | 大尺度湍流结构、长时推演 | + +### L3 推理服务层 + +| 组件 | 说明 | +| --- | --- | +| 推理引擎 | TensorRT / ONNX Runtime，FP16/INT8 量化 | +| 时间推进 | 自回归推演 + 潜空间动力学，支持变时间步 | +| 后处理 | 散度修正（投影法）、守恒律校正、平滑滤波 | +| 不确定性 | 深度集成（Deep Ensemble）/ MC-Dropout，输出置信区间 | +| 服务接口 | gRPC（内部）+ WebSocket（前端流式推送） | + +### L4 交互与渲染层 + +| 组件 | 说明 | +| --- | --- | +| 前端渲染 | WebGPU（three.js / VTK.js），支持体渲染、流线、矢量场 | +| 参数控制 | 实时拖拽 Re/攻角/几何参数，防抖 + 插值过渡 | +| 场景管理 | 视角预设、剖切面、等值面、动画时间轴 | + +### L5 应用层 + +数字孪生监控、工程设计优化、教学演示、实时决策支持。 + +## 2. 数据流 + +``` +参数采样 ──► PyFluent 批量仿真 ──► 原始结果(VTK/CAS) ──► 数据工程(HDF5) + ▲ │ + │ ▼ +主动学习回灌 ◄── 不确定性评估 ◄── 模型评估 ◄── 训练(PyTorch) ◄── 训练集构建 + │ + ▼ + 导出 ONNX ──► TensorRT 引擎 ──► 推理服务 + │ + ▼ + WebSocket ──► WebGPU 前端实时渲染 +``` + +**关键数据契约：** + +| 阶段 | 格式 | 内容 | +| --- | --- | --- | +| 仿真输出 | .cas/.dat + .vtk | 网格 + 稳态/瞬态场（U, p, k, ω） | +| 训练样本 | HDF5 | 输入参数 + 网格坐标 + 场量张量 + 元数据 | +| 模型产物 | .pt → .onnx → .engine | 权重 + 归一化统计 + 版本信息 | +| 服务响应 | 二进制帧 | 场量数组 + 时间戳 + 不确定性 | + +## 3. 技术选型对比 + +### 3.1 模型路线 + +| 方案 | 精度潜力 | 推理速度 | 数据需求 | 复杂几何 | 结论 | +| --- | --- | --- | --- | --- | --- | +| 传统降阶（POD-Galerkin） | 中 | 极快 | 中 | 差 | 作为基线对照 | +| FNO | 高 | 快 | 中高 | 差（需规则网格） | **主力** | +| DeepONet | 中高 | 快 | 中 | 中 | 备选 | +| PINN | 中 | 慢（需重训） | 极低 | 中 | 物理约束模块 | +| MeshGraphNet | 高 | 中 | 高 | **优** | 复杂几何阶段主力 | +| Transformer 算子 | 高 | 中 | 极高 | 中 | 探索性 | + +### 3.2 加速与部署 + +| 选项 | 收益 | 成本 | 结论 | +| --- | --- | --- | --- | +| ONNX Runtime | 跨平台、易部署 | 加速有限 | 开发期默认 | +| TensorRT FP16 | 2-4x | 绑定 NVIDIA | **生产环境** | +| TensorRT INT8 | 4-8x | 精度损失需校准 | 边缘部署评估 | +| CUDA Graphs | 消除启动开销 | 实现复杂 | 小模型场景启用 | + +### 3.3 可视化 + +| 选项 | 优势 | 劣势 | 结论 | +| --- | --- | --- | --- | +| VTK.js + WebGPU | 科学可视化生态完整 | 学习曲线陡 | **首选** | +| three.js | 生态活跃、上手快 | 科学可视化组件少 | 轻量场景 | +| ParaView Web | 功能最强 | 部署重 | 内部调试 | +| Unity/UE | 渲染效果最好 | 工程链路重 | 演示视频 | + +## 4. 工程规范 + +- **数据版本**：DVC 或 Git-LFS 管理数据集与模型产物； +- **实验追踪**：MLflow / Weights & Biases 记录全部训练配置与指标； +- **配置管理**：Hydra 管理实验配置，保证可复现； +- **测试**：单元测试（数据管道、归一化、守恒校正）+ 回归测试（固定工况误差阈值）； +- **CI/CD**：模型训练与评估流水线化，指标不达标自动阻断发布。 diff --git a/docs/03-roadmap.md b/docs/03-roadmap.md new file mode 100644 index 0000000..f36fef6 --- /dev/null +++ b/docs/03-roadmap.md @@ -0,0 +1,162 @@ +# 03 · 路线图与里程碑 + +![路线图](assets/roadmap.png) + +**总周期：12 个月（7 个里程碑）**，采用"先窄后宽"策略：先在 2D 圆柱绕流这个经典基准上打通全链路，再横向扩展到 3D 与复杂几何。 + +--- + +## M0 · 立项与调研（第 1 个月） + +**目标**：统一技术共识，完成可行性论证，锁定技术栈。 + +**任务** +- 调研神经算子家族（FNO / DeepONet / PINN / MeshGraphNet）在流体问题上的公开结果； +- 调研 PyFluent 自动化能力与许可要求； +- 确定首个基准问题：2D 圆柱绕流（Re ∈ [50, 300]，含卡门涡街）； +- 确定数据格式、存储与版本管理方案。 + +**交付物** +- [x] 本仓库规划文档（01-06 号文档） +- [ ] 技术选型评审纪要 +- [ ] 基准问题定义书（几何、边界条件、网格、工况表） + +**验收**：评审通过，环境就绪（Fluent + PyFluent + PyTorch + GPU）。 + +--- + +## M1 · 数据工厂（第 2-3 个月） + +**目标**：建成无人值守的 Fluent 批量仿真流水线。 + +**任务** +- PyFluent 脚本化：参数注入 → 网格生成 → 求解 → 场量导出； +- 参数空间设计：拉丁超立方采样 Re / 攻角 / 障碍物尺寸； +- 稳态 + 瞬态两类数据集：稳态用于场回归，瞬态用于时间推进； +- 数据质检：收敛性检查、网格无关性验证、异常工况剔除。 + +**交付物** +- 批量仿真调度器（支持断点续跑、失败重试、日志追踪） +- 2D 圆柱绕流数据集 v0.1：≥ 500 工况（稳态）+ 20 组瞬态序列 +- 数据质检报告 + +**验收**：单工况端到端自动化成功率 > 95%；数据集通过网格无关性检查。 + +--- + +## M2 · 基线模型（第 4-6 个月） + +**目标**：FNO 基线在 2D 基准上达到可用精度。 + +**任务** +- 数据工程：网格统一化、无量纲化、场归一化、训练/验证/测试划分； +- FNO 实现与调参（模态数、层数、宽度、教师强制 vs 自回归）； +- 物理约束项：连续性方程残差、动量方程残差加入损失； +- 评估体系：相对 L2、涡量谱、长时推演稳定性、泛化测试。 + +**交付物** +- 训练与评估代码（Hydra 配置 + MLflow 追踪） +- FNO 基线模型（2D 稳态 + 瞬态两个版本） +- 评估报告：精度 / 稳定性 / 泛化 + +**验收（Go/No-Go）** +- 稳态：相对 L2 误差 < 5% +- 瞬态：自回归 1000 步不发散 +- 单步推理 < 5 ms（256² 网格，FP16） + +--- + +## M3 · 主动学习与数据飞轮（第 6-7 个月） + +**目标**：用最少的 Fluent 算力换取最大的精度提升。 + +**任务** +- 深度集成（Deep Ensemble）不确定性估计； +- 不确定性驱动的采样策略（OOD 工况优先回灌仿真）； +- 数据飞轮：训练 → 评估 → 采样 → 仿真 → 再训练 自动化。 + +**交付物** +- 主动学习调度器 +- 数据飞轮 v1（一轮迭代周期 < 7 天） +- 数据效率对比报告（相同精度下数据量减少比例） + +**验收**：达到同等精度所需仿真工况数减少 ≥ 40%。 + +--- + +## M4 · 3D 与复杂几何（第 7-9 个月） + +**目标**：把能力扩展到工业级 3D 问题。 + +**任务** +- 3D 数据集：翼型族（NACA 系列）+ 简化车体（Ahmed body）； +- MeshGraphNet 路线：非结构网格上的图神经网络建模； +- 混合策略：规则网格用 FNO，复杂几何用 GNN，统一评估口径； +- 分布式训练：多卡 FSDP，梯度检查点。 + +**交付物** +- 3D 代理模型（翼型 + 简化车体） +- 训练基础设施（多卡、可复现） +- 3D 评估报告（含气动力系数对比） + +**验收** +- 速度场相对 L2 < 8%，升阻力系数误差 < 5% +- 3D 单步推理 < 50 ms（中等网格，TensorRT FP16） + +--- + +## M5 · 实时推理系统（第 9-10 个月） + +**目标**：把模型变成可服务化的毫秒级引擎。 + +**任务** +- ONNX 导出 → TensorRT 优化（FP16/INT8）→ CUDA Graphs； +- 散度修正与守恒校正后处理（投影法/质量修正）； +- 推理服务：gRPC + WebSocket，多并发流式推送； +- 性能工程：批处理、内存池、零拷贝传输。 + +**交付物** +- 推理引擎与容器镜像 +- 性能基准报告（延迟分布 P50/P99、吞吐、显存占用） +- 服务 API 文档 + +**验收**：3D 场景端到端延迟（请求 → 首帧流场）< 100 ms。 + +--- + +## M6 · 可视化与数字孪生原型（第 10-12 个月） + +**目标**：完成可交互演示系统与数字孪生闭环。 + +**任务** +- WebGPU 前端：体渲染 / 流线 / 矢量场 / 剖切面； +- 参数实时拖拽 → 流场更新（含插值过渡与防抖）； +- 数字孪生模式：稀疏传感器数据同化（集合卡尔曼 / 在线微调）； +- 对标验证：≥ 3 个工业场景与 Fluent 真值对比。 + +**交付物** +- 实时可视化 Demo（浏览器可访问） +- 数字孪生原型（观测 → 重建 → 预测） +- 最终验证报告 + 演示视频 + +**验收** +- 交互延迟 < 100 ms +- 数字孪生重建误差 < 10%（稀疏观测场景） +- 对标场景全部通过误差门禁 + +--- + +## 工作流并行视图 + +| 工作流 | M0 | M1 | M2 | M3 | M4 | M5 | M6 | +| --- | :-: | :-: | :-: | :-: | :-: | :-: | :-: | +| 数据 | ● | ● | ● | ● | ● | | | +| 模型 | ● | | ● | ● | ● | ● | ● | +| 系统 | | ● | ● | ● | ● | ● | ● | +| 验证 | ● | ● | ● | ● | ● | ● | ● | + +## 风险缓冲 + +- 每个里程碑预留 20% 时间缓冲； +- M2 为关键技术门禁：若 FNO 精度不达标，转向 DeepONet/MeshGraphNet 双路线并行的决策点； +- M4 若 3D 精度不达标，退守"2D + 切片 3D"策略保证 M5/M6 交付。 diff --git a/docs/04-data-and-ai.md b/docs/04-data-and-ai.md new file mode 100644 index 0000000..7fe4744 --- /dev/null +++ b/docs/04-data-and-ai.md @@ -0,0 +1,105 @@ +# 04 · 数据方案与 AI 模型 + +## 1. 数据工厂（Fluent → 数据集） + +### 1.1 自动化链路 + +``` +参数配置(YAML) → PyFluent 会话 → 几何/网格生成 → 求解设置 → 求解 → 场量导出(VTK/HDF5) +``` + +**关键工程点** + +| 问题 | 方案 | +| --- | --- | +| 长时间批量运行 | 调度器 + 断点续跑 + 失败重试 + 日志追踪 | +| 网格不一致 | 参数化网格模板，统一拓扑或统一目标分辨率 | +| 收敛判据 | 残差 + 监测量（Cd/Cl/压降）双判据 | +| 异常剔除 | 收敛性、网格质量、物理合理性三重检查 | +| 存储 | HDF5 分块压缩（训练）+ VTK（可视化）双写 | + +### 1.2 数据集设计 + +| 数据集 | 内容 | 规模（首期） | 用途 | +| --- | --- | --- | --- | +| D1 稳态 2D | 圆柱绕流，Re ∈ [50, 300] | ≥ 500 工况 | 场回归基线 | +| D2 瞬态 2D | 同上，20 组时间序列 | 20 × 500 帧 | 时间推进 | +| D3 稳态 3D | NACA 翼型族，Ma/攻角扫描 | ≥ 300 工况 | 3D 泛化 | +| D4 稳态 3D | Ahmed body 简化车体 | ≥ 100 工况 | 工业场景 | +| D5 FFU 洁净室 | 送风/回风/热源/尘源参数 | ≥ 200 工况 | 数字孪生 Demo | + +**样本结构（HDF5）** + +``` +/sample_0001 + ├── inputs/ # Re, 攻角, 几何参数, 边界条件 (标量 + 小数组) + ├── coords/ # 网格坐标 (N, 2/3) + ├── fields/ # U(N,3), p(N), k(N), omega(N) ... + └── meta/ # 工况 ID, 求解器设置, 收敛残差 +``` + +### 1.3 数据工程规范 + +- **无量纲化**：以特征长度/来流速度归一化，降低分布漂移； +- **场归一化**：逐通道 z-score 或 min-max，统计量随模型一起版本化； +- **增广**：旋转、镜像、缩放（仅限物理对称场景）； +- **划分**：按参数区间外推划分（OOD 测试集），而非随机划分。 + +## 2. 模型选型 + +### 2.1 主力：FNO（傅里叶神经算子） + +- **原理**：在频域做谱卷积，学习函数空间之间的映射，分辨率无关； +- **输入**：参数 + 坐标网格（或隐式网格）； +- **输出**：场量（U, p, ...）； +- **优势**：训练稳定、推理快、对规则网格精度高； +- **局限**：非结构网格需插值，复杂几何边界处理弱。 + +### 2.2 备选与补充 + +| 模型 | 何时用 | 关键点 | +| --- | --- | --- | +| DeepONet | FNO 精度不足时并行验证 | 分支/主干网络，适合低维参数 | +| PINN | 数据稀缺、反问题 | NS 残差直接进损失，训练慢 | +| MeshGraphNet | 复杂几何、非结构网格 | 图消息传递，几何泛化好 | +| Latent Dynamics | 长时推演加速 | 自编码器 + 潜空间时间推进 | +| Transformer 算子 | 大尺度长程依赖 | 数据需求大，探索性 | + +### 2.3 物理约束（贯穿所有模型） + +| 约束 | 形式 | 作用 | +| --- | --- | --- | +| 连续性 | ‖∇·u‖² 惩罚 | 抑制散度漂移 | +| 动量残差 | ‖∂u/∂t + (u·∇)u + ∇p/ρ - ν∇²u‖² | 物理一致性 | +| 边界条件 | 壁面/入口/出口损失项 | 边界精度 | +| 守恒校正 | 投影法后处理 | 硬约束兜底 | +| 频谱监督 | 涡量谱 / 能谱匹配 | 湍流结构保真 | + +## 3. 训练策略 + +- **两阶段**：先纯数据预训练 → 再物理约束微调； +- **课程学习**：先稳态后瞬态，先低 Re 后高 Re； +- **长时推演**：教师强制 → 逐步放开自回归，配合 scheduled sampling； +- **分布式**：AMP + DDP，大模型用 FSDP + 梯度检查点； +- **实验管理**：Hydra 配置 + MLflow 指标追踪，全配置可复现。 + +## 4. 评估体系 + +| 维度 | 指标 | 目标 | +| --- | --- | --- | +| 场精度 | 相对 L2 误差（U, p） | < 5%（2D）/ < 8%（3D） | +| 结构精度 | 涡量谱、湍动能分布 | 谱斜率误差 < 10% | +| 积分量 | Cd / Cl / 压降 | 误差 < 5% | +| 稳定性 | 1000 步 rollout 发散率 | 0 发散 | +| 泛化 | OOD 参数区间误差退化 | < 2 倍 | +| 效率 | 单步延迟 / 加速比 | < 10 ms / > 1000x | +| 可信度 | 不确定性校准误差 | ECE < 0.1 | + +## 5. 快速代理基线（已落地） + +仓库内 `tools/ffu_surrogate.py` 提供了一条**不依赖训练**的快速代理基线： + +- **方法**：1/3 分辨率网格上运行同一套 D2Q9 LBM，准稳态后双线性重建到全分辨率； +- **性能**：端到端 < 8 s（精细真值 76 s，加速 16-27x）； +- **精度**：与精细解全场相关系数 速度 r=0.95 / 粒子 r=0.83 / 温度 r=0.68； +- **意义**：在 AI 模型成熟前先打通"快速出图 + 趋势验证"链路，同时作为 AI 模型的对照基线。 diff --git a/docs/05-risks.md b/docs/05-risks.md new file mode 100644 index 0000000..ccd5444 --- /dev/null +++ b/docs/05-risks.md @@ -0,0 +1,107 @@ +# 05 · 风险、挑战与缓解策略 + +## 1. 风险总览 + +| 编号 | 风险 | 概率 | 影响 | 等级 | +| --- | --- | --- | --- | --- | +| R1 | 湍流建模精度不足，AI 代理无法复现小尺度结构 | 高 | 高 | 🔴 | +| R2 | 高保真数据生产成本过高，数据量不足 | 高 | 高 | 🔴 | +| R3 | 长时自回归推演误差累积、发散 | 高 | 高 | 🔴 | +| R4 | 分布外（新几何/新 Re）泛化崩溃 | 中 | 高 | 🟠 | +| R5 | 3D 大网格实时推理算力不足 | 中 | 中 | 🟠 | +| R6 | 守恒律漂移导致结果不可信 | 中 | 高 | 🟠 | +| R7 | 工业验证周期长、缺乏可信基准 | 中 | 中 | 🟡 | +| R8 | 工具链锁定（Fluent 许可、GPU 型号） | 低 | 中 | 🟡 | + +## 2. 逐项分析与缓解 + +### R1 湍流精度 + +**问题**：Navier-Stokes 的湍流闭合本质困难，DNS 数据在工程尺度不可得；AI 只能学到训练数据中已有的物理。 + +**缓解** +- 明确精度边界：以 RANS/URANS 数据为主，定位"工程精度"而非"DNS 精度"； +- 混合建模：AI 预测 + 传统湍流模型残差修正（learned closure）； +- 频谱损失监督大尺度结构，放弃不可分辨的小尺度； +- 分阶段验证：先层流基准（圆柱绕流）再湍流工况。 + +### R2 数据成本 + +**问题**：单工况 Fluent 仿真小时级，10⁴ 工况在算力/许可/时间上都昂贵。 + +**缓解** +- 主动学习：只对高价值（高不确定性/OOD）工况做高保真仿真； +- 降阶采样：拉丁超立方 + 参数敏感性分析，避免均匀暴力扫描； +- 多保真度：粗网格大量生成 + 细网格少量校核（coarse-to-fine）； +- 复用收敛解：邻近工况热启动，缩短求解时间。 + +### R3 长时稳定性 + +**问题**：自回归每步误差被放大，几百步后发散。 + +**缓解** +- 潜空间推进：在低维潜空间做时间积分，抑制高频误差； +- 物理投影：每 N 步做散度修正 / 守恒校正； +- 噪声注入训练：模拟推演误差，提高鲁棒性； +- 早停与不确定性告警：置信度低于阈值时切换传统求解器。 + +### R4 泛化 + +**问题**：训练分布外的几何/参数表现急剧下降。 + +**缓解** +- 几何参数化：用统一参数空间描述几何族，而非记忆具体形状； +- 数据增广 + 课程学习：从简单到复杂逐步扩展参数范围； +- OOD 检测：不确定性 + 输入分布检验，超界自动拒绝并回退； +- 在线微调：新场景少量数据快速适配（few-shot fine-tune）。 + +### R5 3D 实时算力 + +**问题**：3D 网格点数比 2D 高 2-3 个数量级，单步推理可能超时。 + +**缓解** +- 潜空间降维 + 多尺度分解（粗网格动力学 + 局部细化）； +- TensorRT FP16/INT8 + CUDA Graphs； +- 时间异步：用户交互时低频更新，静止时高精度刷新； +- 分辨率自适应：按视角/关注区域动态调整计算分辨率。 + +### R6 守恒律漂移 + +**问题**：纯数据驱动的网络不严格守恒质量/动量。 + +**缓解** +- 硬约束：输出层投影到无散空间（Helmholtz 分解）； +- 软约束：守恒损失项 + 后处理校正； +- 监控：把守恒误差作为在线指标，超限告警。 + +### R7 验证可信度 + +**问题**：没有公认基准，工业界不信任 AI 结果。 + +**缓解** +- 建立内部基准套件（圆柱绕流 / 翼型 / 车体 / 洁净室）； +- 与 Fluent 真值做逐场对比 + 积分量对比 + 不确定性报告； +- 所有指标可一键复现，公开评估脚本。 + +### R8 工具链 + +**缓解**：核心代码与 Fluent 解耦（数据接口标准化）；推理引擎支持 ONNX 多后端；GPU 无关的降级路径（CPU 推理，精度不变、速度降级）。 + +## 3. 开放问题（需要在 P1-P2 回答） + +1. FNO 在 Re > 500 的圆柱绕流上能否保持 < 5% 误差？需要多少工况？ +2. 物理约束损失的最佳权重与调度策略是什么？ +3. 潜空间推进相比直接自回归，长时稳定性能提升多少？ +4. 3D 非结构网格上，MeshGraphNet 与"插值到规则网格 + FNO"哪条路线更优？ +5. 不确定性估计能否可靠地识别 OOD 工况？ +6. 单卡推理延迟与网格规模的定量 scaling law？ + +## 4. 止损与回退 + +| 触发条件 | 回退方案 | +| --- | --- | +| M2 精度不达标 | 转 DeepONet / MeshGraphNet 双路线并行 | +| M4 3D 精度不达标 | 退守"2D + 切片 3D"保证 M5/M6 交付 | +| 实时性不达标 | 降分辨率 / 潜空间 / 云端推理 | +| 数据成本超预算 | 多保真度 + 主动学习压缩数据量 | +| 可信度存疑 | 混合模式：AI 预览 + Fluent 校核 | diff --git a/docs/06-references.md b/docs/06-references.md new file mode 100644 index 0000000..3f6586b --- /dev/null +++ b/docs/06-references.md @@ -0,0 +1,59 @@ +# 06 · 参考文献与工具索引 + +## 1. 神经算子与物理机器学习 + +| 文献 | 内容 | 链接 | +| --- | --- | --- | +| Li et al., *Fourier Neural Operator for Parametric PDEs* (ICLR 2021) | FNO 原始论文 | arXiv:2010.08895 | +| Lu et al., *DeepONet: Learning nonlinear operators* (2021) | DeepONet 算子网络 | arXiv:1910.03193 | +| Raissi et al., *Physics-Informed Neural Networks* (JCP 2019) | PINN 物理约束 | arXiv:1711.10561 | +| Pfaff et al., *Learning Mesh-Based Simulation with Graph Networks* (ICLR 2021) | MeshGraphNet | arXiv:2010.03409 | +| Hao et al., *GNOT: A General Neural Operator Transformer* (ICML 2023) | Transformer 算子 | arXiv:2302.14376 | +| Kovachki et al., *Neural Operator: Learning Maps Between Function Spaces* (2021) | 神经算子综述 | arXiv:2108.08481 | +| Kochkov et al., *Machine learning accelerated computational fluid dynamics* (PNAS 2021) | ML 加速 CFD | arXiv:2102.01010 | +| Li et al., *Learned Turbulence Closure for LES* | 学习型湍流闭合 | — | + +## 2. 数据集与基准 + +| 名称 | 内容 | 链接 | +| --- | --- | --- | +| PDEBench | 多 PDE 基准数据集 | https://github.com/pdebench/PDEBench | +| PDEArena | 神经算子基准 | https://github.com/microsoft/PDEArena | +| The Well | 大规模物理场数据集 | https://github.com/PolymathicAI/the_well | +| FlowBench | 流场预测基准 | — | +| AirfRANS | 翼型 RANS 数据集 | https://github.com/Extrality/AirfRANS | +| DeepMind Cylinder Flow | 圆柱绕流高分辨率数据 | — | + +## 3. 开源工具 + +| 工具 | 用途 | 链接 | +| --- | --- | --- | +| PyFluent | Python 驱动 Ansys Fluent | https://github.com/ansys/pyfluent | +| NeuralOperator | FNO/DeepONet 官方实现 | https://github.com/neuraloperator/neuraloperator | +| NVIDIA Modulus | 物理 AI 训练框架 | https://github.com/NVIDIA/modulus | +| PhiFlow | 可微流体求解器 | https://github.com/tum-pbs/PhiFlow | +| JAX-CFD | JAX 可微 CFD | https://github.com/google/jax-cfd | +| OpenFOAM | 开源 CFD（替代/对照） | https://openfoam.org | +| VTK.js | 浏览器端科学可视化 | https://kitware.github.io/vtk-js | +| three.js | WebGL/WebGPU 渲染 | https://threejs.org | +| ONNX Runtime | 跨平台推理 | https://onnxruntime.ai | +| TensorRT | NVIDIA 推理加速 | https://developer.nvidia.com/tensorrt | +| DVC | 数据/模型版本管理 | https://dvc.org | +| MLflow | 实验追踪 | https://mlflow.org | +| Hydra | 配置管理 | https://hydra.cc | + +## 4. 洁净室 / FFU 相关 + +| 资料 | 内容 | +| --- | --- | +| ISO 14644-1 | 洁净室分级（ISO Class 1-9） | +| ISO 14644-4 | 洁净室设计与施工（单向流风速 0.45 m/s ±20%） | +| IEST-RP-CC002 | FFU 性能测试推荐规程 | +| ASHRAE Handbook | 洁净空间气流组织设计 | + +## 5. 经典 CFD 参考 + +- Versteeg & Malalasekera, *An Introduction to Computational Fluid Dynamics* +- Pope, *Turbulent Flows* +- Krüger et al., *The Lattice Boltzmann Method: Principles and Practice*（本项目 LBM Demo 的理论依据） +- Succi, *The Lattice Boltzmann Equation for Fluid Dynamics and Beyond* diff --git a/docs/07-ffu-demo.md b/docs/07-ffu-demo.md new file mode 100644 index 0000000..aec27ec --- /dev/null +++ b/docs/07-ffu-demo.md @@ -0,0 +1,109 @@ +# 07 · FFU 洁净室流场实时仿真 Demo + +> FFU（Fan Filter Unit，风机过滤单元）下方洁净区的流场、温湿度、压力与尘埃粒子输运仿真。 +> 这是"物理世界流场实时仿真"的首个可运行 Demo，同时提供 **高保真真值** 与 **毫秒级快速代理** 两条路径。 + +--- + +## 1. 物理场景 + +| 项目 | 设置 | +| --- | --- | +| 计算域 | 4.0 m（宽）× 3.0 m（高）2D 截面 | +| FFU 送风 | 顶棚满布，面风速 **0.45 m/s**（ISO 14644-4 单向流推荐值） | +| 回风 | 底部两侧回风口（0~0.5 m / 3.5~4.0 m），高架地板/下侧回风 | +| 工作台 | x ∈ [0.60, 1.80]，高 0.90 m | +| 设备 | x ∈ [0.85, 1.50]，y ∈ [0.90, 1.25] —— **热源 + 尘源** | +| 操作人员 | 躯干 x ∈ [2.35, 2.75] 高 1.45 m + 头部（r = 0.15 m）—— **湿源 + 尘源** | +| 物料架 | x ∈ [3.00, 3.50]，高 1.15 m | +| 温度 | 送风 22 °C，设备表面 35 °C | +| 湿度 | 送风 RH 45%，呼吸区 RH 70% | +| 粒子 | 设备 + 人员双尘源，含重力沉降（被动标量 + 沉降速度） | +| 雷诺数 | Re = U·H/ν ≈ 540（层流单向流） | + +## 2. 两条技术路径 + +| | 高保真真值 | 快速代理 | +| --- | --- | --- | +| 脚本 | `tools/ffu_simulation.py` | `tools/ffu_surrogate.py` | +| 网格 | 360 × 270（90 px/m） | 120 × 90（1/3 分辨率） | +| 步数 | 14 000 | 5 000 | +| 求解 | D2Q9 LBM + 半拉格朗日标量输运 | **同一套求解器**（粗网格）+ 双线性重建 | +| 端到端耗时 | **76 s** | **7.3 s（16-27x 加速）** | +| 定位 | 验证真值 / 精度基准 | 交互预览 / 趋势级实时 | + +**为什么"粗网格"是当前最好的快速代理？** +在 AI 神经算子尚未训练完成的阶段，粗网格求解用同一套物理保证趋势一致（速度场相关系数 **r = 0.95**）， +是零训练成本、零泛化风险的实时路径；同时它天然作为后续 FNO/DeepONet 的对照基线。 + +## 3. 趋势验证结果 + +![趋势验证](assets/ffu-validation.png) + +| 场量 | 全场相关系数（vs 高保真） | 结论 | +| --- | --- | --- | +| 速度 | **r = 0.95** | 绕障加速、贴壁射流、回流区结构一致 | +| 粒子浓度 | **r = 0.83** | 羽流路径与死角堆积一致 | +| 温度 | **r = 0.68** | 热羽下扫与侧向输运一致（粗网格扩散略强） | + +> 代理模型保证**趋势级**精度，局部场量精度不做保证；需要精确值时切换高保真求解。 + +## 4. 结果图 + +### 4.1 流场（高保真） + +![流场](assets/ffu-flow-field.png) + +### 4.2 温湿度 / 压力 / 粒子（高保真） + +![标量场](assets/ffu-scalar-fields.png) + +### 4.3 快速代理（7 秒出图） + +![快速流场](assets/ffu-flow-field-fast.png) + +![快速标量场](assets/ffu-scalar-fields-fast.png) + +### 4.4 粒子轨迹与滞留时间 + +![粒子轨迹](assets/ffu-particle-tracks.png) + +## 5. 主要物理结论（趋势） + +1. **单向下扫**：FFU 送风把设备/人员散发的热量、湿气与粒子持续向下输运，经底部回风口排出——这是洁净室自净能力的来源； +2. **绕障加速**：气流在工作台、人员、物料架之间被挤压加速，贴壁/缝隙处流速可达 ~1.5 m/s（约为送风速度的 3 倍）； +3. **死角堆积**：障碍物背风侧与墙角存在低速回流区，污染物滞留时间显著变长（粒子轨迹图中颜色越暖滞留越久）——**这是洁净室布局优化的关键目标**； +4. **压力分布**：送风端相对正压、回风口相对负压，全场压差 < 2 Pa（洁净室典型量级）； +5. **人员源不可忽视**：呼吸区湿羽流与尘源在单向流下形成明显的下扫羽流，直接影响工作区洁净度。 + +## 6. 模型局限（诚实声明） + +- **2D 截面**：无法刻画三维绕流与横向扩散； +- **等温假设**：未耦合浮力（热羽流完全由强制对流主导），真实热源上方可能存在浮力竞争； +- **层流/低 Re**：未启用湍流模型，高换气次数下的湍流脉动未建模； +- **粒子为被动标量 + 沉降**：未包含凝聚、静电、表面沉积再悬浮； +- **非洁净度认证**：粒子浓度为归一化相对值，不可直接对应 ISO 14644-1 等级。 + +## 7. 复现命令 + +```bash +# 高保真真值（76 s，输出 3 张图 + npz 场数据） +python tools/ffu_simulation.py + +# 快速代理（7 s，输出 2 张图 + 趋势验证对比图） +python tools/ffu_surrogate.py + +# 规划图（架构 / 流水线 / 路线图 SVG） +python tools/gen_diagrams.py + +# SVG -> PNG +npm install && node tools/svg2png.js +``` + +## 8. 下一步 + +- [ ] 加入浮力项（Boussinesq）验证热羽流竞争； +- [ ] 参数化扫描：送风速度 × 障碍物布局 × 源强度 → 数据集 D5； +- [ ] 训练 FNO 代理：以 D5 为数据，对标粗网格代理的速度与精度； +- [ ] 三维化：FFU 阵列 + 局部障碍物 3D 流场； +- [ ] 接入 WebGPU 前端，实现参数拖拽实时刷新。 diff --git a/docs/assets/architecture.png b/docs/assets/architecture.png new file mode 100644 index 0000000..8a9e596 Binary files /dev/null and b/docs/assets/architecture.png differ diff --git a/docs/assets/architecture.svg b/docs/assets/architecture.svg new file mode 100644 index 0000000..9511d9e --- /dev/null +++ b/docs/assets/architecture.svg @@ -0,0 +1,142 @@ + + + + + + + + +PhysicsAiFluentAnalysis · 系统架构 +六层架构：从高保真真值到实时交互应用，物理约束贯穿模型层 + + +L5 应用层 +Application + + +数字孪生 +传感器同化 · 在线重建 + + +工程设计优化 +实时参数扫描 · 设计空间探索 + + +教学与科研 +流场可视化实验 + + +实时决策 +风场 / 流场秒级预测 + + + + +L4 交互与渲染 +Interaction & Rendering + + +WebGPU 前端 +three.js / VTK.js + + +参数实时调节 +拖拽 Re / 攻角 / 几何 + + +场可视化 +体渲染 · 流线 · 矢量场 + + +场景管理 +剖切面 · 时间轴 · 导出 + + + + +L3 推理服务层 +Inference Service + + +推理引擎 +TensorRT / ONNX Runtime + + +时间推进 +自回归 + 潜空间动力学 + + +后处理校正 +散度修正 · 守恒律 + + +服务接口 +gRPC + WebSocket 流式 + + + + +L2 AI 模型层 +AI Models + + +FNO 傅里叶神经算子 +主力基线 · 分辨率无关 + + +DeepONet +算子学习备选路线 + + +PINN 物理约束 +NS 残差 · 少数据工况 + + +MeshGraphNet +非结构网格 · 复杂几何 + + + + +L1 数据与训练 +Data & Training + + +PyTorch 分布式训练 +AMP · DDP / FSDP + + +主动学习 +不确定性驱动采样 + + +数据湖 +HDF5 / VTK / Zarr + + +数据工程 +归一化 · 增广 · 插值统一 + + + + +L0 高保真真值 +Ground Truth CFD + + +Ansys Fluent +RANS / URANS / LES + + +PyFluent 自动化 +参数化建模与求解 + + +参数化几何 +翼型族 · 车体 · 管道 + + +边界条件库 +工况模板 · 湍流模型 +数据流 ↑（真值 → 训练 → 推理 → 应用） 控制流 ↓（应用需求 → 模型选型 → 数据采集） 物理约束：连续性 / 动量守恒 · 不确定性量化 · 守恒校正 + diff --git a/docs/assets/ffu-flow-field-fast.png b/docs/assets/ffu-flow-field-fast.png new file mode 100644 index 0000000..f86ab60 Binary files /dev/null and b/docs/assets/ffu-flow-field-fast.png differ diff --git a/docs/assets/ffu-flow-field.png b/docs/assets/ffu-flow-field.png new file mode 100644 index 0000000..3edfb4b Binary files /dev/null and b/docs/assets/ffu-flow-field.png differ diff --git a/docs/assets/ffu-flow-field.svg b/docs/assets/ffu-flow-field.svg new file mode 100644 index 0000000..fb5f9a9 --- /dev/null +++ b/docs/assets/ffu-flow-field.svg @@ -0,0 +1,6553 @@ + + + + + + + + 2026-09-17T23:57:57.408470 + image/svg+xml + + + Matplotlib v3.11.2, https://matplotlib.org/ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/docs/assets/ffu-particle-tracks.png b/docs/assets/ffu-particle-tracks.png new file mode 100644 index 0000000..826cd92 Binary files /dev/null and b/docs/assets/ffu-particle-tracks.png differ diff --git a/docs/assets/ffu-scalar-fields-fast.png b/docs/assets/ffu-scalar-fields-fast.png new file mode 100644 index 0000000..d6aa4f4 Binary files /dev/null and b/docs/assets/ffu-scalar-fields-fast.png differ diff --git a/docs/assets/ffu-scalar-fields.png b/docs/assets/ffu-scalar-fields.png new file mode 100644 index 0000000..96b01a2 Binary files /dev/null and b/docs/assets/ffu-scalar-fields.png differ diff --git a/docs/assets/ffu-validation.png b/docs/assets/ffu-validation.png new file mode 100644 index 0000000..f6fcb3f Binary files /dev/null and b/docs/assets/ffu-validation.png differ diff --git a/docs/assets/pipeline.png b/docs/assets/pipeline.png new file mode 100644 index 0000000..be7a416 Binary files /dev/null and b/docs/assets/pipeline.png differ diff --git a/docs/assets/pipeline.svg b/docs/assets/pipeline.svg new file mode 100644 index 0000000..7a9c1e1 --- /dev/null +++ b/docs/assets/pipeline.svg @@ -0,0 +1,77 @@ + + + + + + + + +数据与训练流水线 +从参数采样到实时仿真服务：全链路自动化 + 主动学习数据飞轮 + + +① 参数空间采样 +拉丁超立方采样 +Re / 攻角 / 几何参数 + + +工况表 + 元数据 + + +② Fluent 批量仿真 +PyFluent 无人值守 +稳态 + 瞬态工况 + + +CAS / DAT / VTK + + +③ 数据工程 +网格统一 · 归一化 +质检 / 增广 / 划分 + + +HDF5 数据集 + + +④ 神经算子训练 +FNO / DeepONet / GNN +物理约束损失 + + +模型权重 .pt + + +⑤ 推理引擎部署 +ONNX → TensorRT +FP16 / INT8 量化 + + +TensorRT .engine + + +⑥ 实时仿真服务 +gRPC + WebSocket +流场流式推送 + +实时流场帧 + +主动学习闭环：不确定性评估 → 高价值工况回灌仿真 → 再训练（数据飞轮） + + +精度门禁 +相对 L2 < 5% + + +稳定性门禁 +1000 步不发散 + + +成本门禁 +单工况预算达标 + + +可复现门禁 +流水线一键复现 +每个阶段设置 Go / No-Go 决策点：未通过门禁不允许进入下一阶段 + diff --git a/docs/assets/roadmap.png b/docs/assets/roadmap.png new file mode 100644 index 0000000..cc6f6b6 Binary files /dev/null and b/docs/assets/roadmap.png differ diff --git a/docs/assets/roadmap.svg b/docs/assets/roadmap.svg new file mode 100644 index 0000000..be5ce97 --- /dev/null +++ b/docs/assets/roadmap.svg @@ -0,0 +1,84 @@ + + + + + + + + +里程碑路线图（12 个月 · 7 个里程碑） +先窄后宽：2D 基准打通全链路 → 3D 复杂几何 → 实时系统 → 数字孪生 + + + + + + + + + + + + + +1 +2 +3 +4 +5 +6 +7 +8 +9 +10 +11 +12 +月份 +M0 立项与调研 + + +第 1 月 + +规划文档 · 选型评审 · 基准定义 +M1 数据工厂 + + +2-3 月 + +PyFluent 自动化 · ≥500 工况数据集 +M2 基线模型 + + +4-6 月 + +FNO 基线 · 相对 L2 < 5% · 1000 步稳定 +M3 数据飞轮 + + +6-7 月 + +主动学习 · 同等精度省 40% 仿真量 +M4 3D 与复杂几何 + + +7-9 月 + +翼型 / 车体代理 · MeshGraphNet +M5 实时推理系统 + + +9-10 月 + +TensorRT · 端到端 < 100 ms +M6 可视化与数字孪生 + + +10-12 月 +WebGPU 交互 · 稀疏观测同化 + + +决策点 +M2 精度不达标 → DeepONet / MeshGraphNet 双路线并行 M4 3D 不达标 → 退守「2D + 切片 3D」保证 M5/M6 交付 每里程碑预留 20% 时间缓冲 +成功标准 +P1（3 个月）可行性验证 · P2（6 个月）系统闭环 · P3（12 个月）数字孪生原型 + diff --git a/package.json b/package.json new file mode 100644 index 0000000..1ad6962 --- /dev/null +++ b/package.json @@ -0,0 +1,12 @@ +{ + "name": "physics-ai-fluent-analysis-tools", + "private": true, + "version": "0.1.0", + "description": "PhysicsAiFluentAnalysis 图件生成工具 (SVG -> PNG)", + "scripts": { + "svg2png": "node tools/svg2png.js" + }, + "dependencies": { + "@resvg/resvg-js": "^2.6.2" + } +} diff --git a/tools/ffu_simulation.py b/tools/ffu_simulation.py new file mode 100644 index 0000000..68d3f6f --- /dev/null +++ b/tools/ffu_simulation.py @@ -0,0 +1,648 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- +""" +FFU 单向流洁净区流场仿真 Demo +================================ +用 D2Q9 格子玻尔兹曼方法 (LBM) 模拟 FFU (Fan Filter Unit, 风机过滤单元) +下方洁净区的流场, 并用半拉格朗日输运方程求解: + + - 温度场 T (设备发热源) + - 相对湿度 RH (人员呼吸/发湿源) + - 尘埃粒子浓度 C (设备 + 人员散发, 含重力沉降) + - 静压场 p (LBM 直接输出, 转换为 Pa) + +物理设置 +-------- + 计算域 : 4.0 m (宽) x 3.0 m (高), 90 px/m + FFU 送风 : 顶棚满布, 面风速 0.45 m/s (ISO 14644 单向流推荐值) + 回风口 : 底部两侧 (0~0.5 m, 3.5~4.0 m), 高架地板/下侧回风 + 障碍物 : 工作台 + 设备(热/尘源) + 操作人员(湿/尘源) + 物料架 + 雷诺数 : Re = U*H/nu ~ 540 (层流单向流) + 网格 : 360 x 270 = 97200 格点, D2Q9 BGK, tau = 0.62 + +输出 +---- + docs/assets/ffu-flow-field.png 速度场 + 流线 (主图) + docs/assets/ffu-scalar-fields.png 温度 / 湿度 / 粒子浓度 / 压力 四联图 + docs/assets/ffu-particle-tracks.png 粒子拉格朗日轨迹 + docs/assets/ffu-flow-field.svg 主图矢量版 + +运行 +---- + python tools/ffu_simulation.py +""" +import os +import numpy as np +from numba import njit + +import matplotlib +matplotlib.use("Agg") +import matplotlib.pyplot as plt +from matplotlib.patches import Rectangle, Circle, FancyArrow +from matplotlib.colors import LinearSegmentedColormap +from matplotlib.collections import LineCollection + +# -------------------------------------------------------------------------- +# 全局样式 +# -------------------------------------------------------------------------- +BG = "#0d1117" +PANEL = "#0d1117" +FG = "#e6edf3" +MUTED = "#8b949e" +ACCENT = "#58a6ff" +WARM = "#ff7b72" +COOL = "#39d0d8" +GREEN = "#7ee787" + +plt.rcParams.update({ + "font.sans-serif": ["Microsoft YaHei", "SimHei", "DejaVu Sans"], + "axes.unicode_minus": False, + "svg.fonttype": "path", + "figure.facecolor": BG, + "savefig.facecolor": BG, + "axes.facecolor": BG, + "text.color": FG, + "axes.labelcolor": FG, + "xtick.color": MUTED, + "ytick.color": MUTED, + "axes.edgecolor": "#30363d", +}) + +# -------------------------------------------------------------------------- +# 物理 / 数值参数 +# -------------------------------------------------------------------------- +W_M, H_M = 4.0, 3.0 # 计算域尺寸 [m] +PPM = 90 # 像素/米 +NX, NY = int(W_M * PPM), int(H_M * PPM) # 360 x 270 + +U0_LB = 0.08 # FFU 送风速度 [lattice] +U_PHYS = 0.45 # FFU 送风速度 [m/s] +SCALE_U = U_PHYS / U0_LB # lattice -> m/s +NU_LB = 0.04 # 运动粘度 [lattice] +TAU = 0.5 + 3.0 * NU_LB +RE = U0_LB * NY / NU_LB + +T_AMB, T_SRC = 22.0, 35.0 # 环境温度 / 设备表面温度 [degC] +RH_AMB, RH_SRC = 45.0, 70.0 # 送风湿度 / 呼吸区湿度 [%] +RHO_AIR = 1.2 # 空气密度 [kg/m3] +SCALE_P = RHO_AIR * SCALE_U ** 2 # lattice -> Pa + +DT_SCALAR = 4 # 标量输运时间步 (每 4 个 LBM 步输运一次) +ADV_EVERY = 4 +KAPPA_T = 0.0020 # 温度扩散系数 [lattice] +KAPPA_RH = 0.0020 +KAPPA_C = 0.0015 +V_SETTLE = 0.004 # 粒子沉降速度 [lattice] (~2 cm/s) + +N_STEPS = 14000 # LBM 总步数 +SCALAR_FROM = 6000 # 从该步开始输运标量 +TRACK_FROM = 9000 # 从该步开始记录粒子轨迹 + +# D2Q9 +EX = np.array([0, 1, 0, -1, 0, 1, -1, -1, 1], dtype=np.int32) +EY = np.array([0, 0, 1, 0, -1, 1, 1, -1, -1], dtype=np.int32) +WT = np.array([4/9, 1/9, 1/9, 1/9, 1/9, 1/36, 1/36, 1/36, 1/36]) +OPP = np.array([0, 3, 4, 1, 2, 7, 8, 5, 6], dtype=np.int32) + +# -------------------------------------------------------------------------- +# 几何: 单位 m +# -------------------------------------------------------------------------- +BENCH = (0.60, 0.00, 1.80, 0.90) # 工作台 (x0,y0,x1,y1) +EQUIP = (0.85, 0.90, 1.50, 1.25) # 设备 (热源 + 尘源) +BODY = (2.35, 0.00, 2.75, 1.45) # 人员躯干 +HEAD_C = (2.55, 1.62, 0.15) # 人员头部 (圆心, 半径) +RACK = (3.00, 0.00, 3.50, 1.15) # 物料架 +GRILL_L = (0.00, 0.50) # 左侧回风口 x 范围 +GRILL_R = (3.50, 4.00) # 右侧回风口 x 范围 + + +def build_geometry(): + """构建固体掩膜 + 源项掩膜""" + solid = np.zeros((NY, NX), dtype=np.bool_) + solid[:, :2] = True + solid[:, -2:] = True + solid[0, :] = True + + def rect(x0, y0, x1, y1, m): + ix0, ix1 = int(round(x0 * PPM)), int(round(x1 * PPM)) + iy0, iy1 = int(round(y0 * PPM)), int(round(y1 * PPM)) + m[max(iy0, 0):min(iy1 + 1, NY), max(ix0, 0):min(ix1 + 1, NX)] = True + + for r in (BENCH, EQUIP, BODY, RACK): + rect(*r, solid) + + cx, cy, rr = HEAD_C + yy, xx = np.mgrid[0:NY, 0:NX] + head = (xx - cx * PPM) ** 2 + (yy - cy * PPM) ** 2 <= (rr * PPM) ** 2 + solid |= head + + # 回风口 (底部, 非固体 => 压力出口) + outlet = np.zeros_like(solid) + iy = 0 + for g0, g1 in (GRILL_L, GRILL_R): + ix0, ix1 = int(g0 * PPM), int(g1 * PPM) + outlet[iy, max(ix0, 2):min(ix1, NX - 2)] = True + solid[outlet] = False + + inlet = np.zeros_like(solid) + inlet[NY - 1, 2:NX - 2] = True + + # 设备热/尘源: 设备外扩 3 格的流体单元 + eq = np.zeros_like(solid) + rect(*EQUIP, eq) + shell = np.zeros_like(solid) + for dy in range(-3, 4): + for dx in range(-3, 4): + shell |= np.roll(np.roll(eq, dy, 0), dx, 1) + src_equip = shell & ~solid + + # 人员呼吸/发湿区: 头部前方圆形区域 + br = (xx - (cx + 0.20) * PPM) ** 2 + (yy - cy * PPM) ** 2 <= (0.22 * PPM) ** 2 + src_head = br & ~solid + + return solid, inlet, outlet, src_equip, src_head + + +# -------------------------------------------------------------------------- +# LBM 求解器 (numba) +# -------------------------------------------------------------------------- +@njit(cache=True, fastmath=True) +def lbm_collide(f, rho, ux, uy, solid, ex, ey, w, tau): + ny, nx = rho.shape + for iy in range(ny): + for ix in range(nx): + if solid[iy, ix]: + continue + r = 0.0; mx = 0.0; my = 0.0 + for k in range(9): + fk = f[k, iy, ix] + r += fk + mx += ex[k] * fk + my += ey[k] * fk + u = mx / r + v = my / r + rho[iy, ix] = r + ux[iy, ix] = u + uy[iy, ix] = v + uu = u * u + v * v + for k in range(9): + eu = ex[k] * u + ey[k] * v + feq = w[k] * r * (1.0 + 3.0 * eu + 4.5 * eu * eu - 1.5 * uu) + f[k, iy, ix] -= (f[k, iy, ix] - feq) / tau + + +@njit(cache=True, fastmath=True) +def lbm_bounce_back(f, solid, opp): + ny, nx = solid.shape + for iy in range(ny): + for ix in range(nx): + if solid[iy, ix]: + for k in range(9): + ko = opp[k] + if k < ko: + tmp = f[k, iy, ix] + f[k, iy, ix] = f[ko, iy, ix] + f[ko, iy, ix] = tmp + + +@njit(cache=True, fastmath=True) +def lbm_stream(f, ftmp, ex, ey): + ny, nx = f.shape[1], f.shape[2] + for k in range(9): + dx = ex[k]; dy = ey[k] + for iy in range(ny): + sy = iy - dy + for ix in range(nx): + sx = ix - dx + if 0 <= sy < ny and 0 <= sx < nx: + ftmp[k, iy, ix] = f[k, sy, sx] + else: + ftmp[k, iy, ix] = f[k, iy, ix] + + +@njit(cache=True, fastmath=True) +def lbm_bcs(f, ux, uy, solid, inlet, outlet, ex, ey, w, u0): + ny, nx = solid.shape + for ix in range(nx): + if inlet[ny - 1, ix]: + u = 0.0; v = -u0; r = 1.0 + uu = u * u + v * v + for k in range(9): + eu = ex[k] * u + ey[k] * v + f[k, ny - 1, ix] = w[k] * r * (1.0 + 3.0 * eu + 4.5 * eu * eu - 1.5 * uu) + if outlet[0, ix]: + u = ux[1, ix]; v = uy[1, ix]; r = 1.0 + uu = u * u + v * v + for k in range(9): + eu = ex[k] * u + ey[k] * v + f[k, 0, ix] = w[k] * r * (1.0 + 3.0 * eu + 4.5 * eu * eu - 1.5 * uu) + + +# -------------------------------------------------------------------------- +# 标量输运 (半拉格朗日 + 扩散) +# -------------------------------------------------------------------------- +@njit(cache=True, fastmath=True) +def advect_scalars(F, ux, uy, solid, kappa, dt, vs): + nc, ny, nx = F.shape + out = np.empty_like(F) + for c in range(nc): + k = kappa[c] + s = vs[c] + for iy in range(ny): + for ix in range(nx): + fc = F[c, iy, ix] + if solid[iy, ix]: + out[c, iy, ix] = 0.0 + continue + xd = ix - ux[iy, ix] * dt + yd = iy - (uy[iy, ix] + s) * dt + if xd < 0.0: xd = 0.0 + if xd > nx - 1.001: xd = nx - 1.001 + if yd < 0.0: yd = 0.0 + if yd > ny - 1.001: yd = ny - 1.001 + x0 = int(xd); y0 = int(yd) + if x0 > nx - 2: x0 = nx - 2 + if y0 > ny - 2: y0 = ny - 2 + fx = xd - x0; fy = yd - y0 + v00 = F[c, y0, x0]; v01 = F[c, y0, x0 + 1] + v10 = F[c, y0 + 1, x0]; v11 = F[c, y0 + 1, x0 + 1] + if solid[y0, x0]: v00 = fc + if solid[y0, x0 + 1]: v01 = fc + if solid[y0 + 1, x0]: v10 = fc + if solid[y0 + 1, x0 + 1]: v11 = fc + val = (v00 * (1 - fx) * (1 - fy) + v01 * fx * (1 - fy) + + v10 * (1 - fx) * fy + v11 * fx * fy) + l = fc; r = fc; d = fc; u = fc + if ix > 0 and not solid[iy, ix - 1]: l = F[c, iy, ix - 1] + if ix < nx - 1 and not solid[iy, ix + 1]: r = F[c, iy, ix + 1] + if iy > 0 and not solid[iy - 1, ix]: d = F[c, iy - 1, ix] + if iy < ny - 1 and not solid[iy + 1, ix]: u = F[c, iy + 1, ix] + out[c, iy, ix] = val + k * dt * (l + r + d + u - 4.0 * fc) + return out + + +# -------------------------------------------------------------------------- +# 粒子轨迹 (拉格朗日) +# -------------------------------------------------------------------------- +@njit(cache=True, fastmath=True) +def sample_uv(ux, uy, x, y): + ny, nx = ux.shape + if x < 0: x = 0.0 + if x > nx - 1.001: x = nx - 1.001 + if y < 0: y = 0.0 + if y > ny - 1.001: y = ny - 1.001 + x0 = int(x); y0 = int(y) + if x0 > nx - 2: x0 = nx - 2 + if y0 > ny - 2: y0 = ny - 2 + fx = x - x0; fy = y - y0 + u = (ux[y0, x0] * (1 - fx) * (1 - fy) + ux[y0, x0 + 1] * fx * (1 - fy) + + ux[y0 + 1, x0] * (1 - fx) * fy + ux[y0 + 1, x0 + 1] * fx * fy) + v = (uy[y0, x0] * (1 - fx) * (1 - fy) + uy[y0, x0 + 1] * fx * (1 - fy) + + uy[y0 + 1, x0] * (1 - fx) * fy + uy[y0 + 1, x0 + 1] * fx * fy) + return u, v + + +@njit(cache=True, fastmath=True) +def step_particles(px, py, alive, age, ux, uy, solid, vs): + n = px.shape[0] + ny, nx = solid.shape + for i in range(n): + if not alive[i]: + continue + x, y = px[i], py[i] + u1, v1 = sample_uv(ux, uy, x, y) + xm, ym = x + 0.5 * u1, y + 0.5 * (v1 - vs) + u2, v2 = sample_uv(ux, uy, xm, ym) + xn, yn = x + u2, y + (v2 - vs) + ix, iy = int(xn), int(yn) + if 0 <= ix < nx and 0 <= iy < ny and not solid[iy, ix]: + px[i] = xn; py[i] = yn + else: + ix, iy = int(xn), int(y) + if 0 <= ix < nx and 0 <= iy < ny and not solid[iy, ix]: + px[i] = xn + else: + ix, iy = int(x), int(yn) + if 0 <= ix < nx and 0 <= iy < ny and not solid[iy, ix]: + py[i] = yn + age[i] += 1.0 + if py[i] < 1.5 or px[i] < 2.0 or px[i] > nx - 3.0 or age[i] > 2500: + alive[i] = False + + +# -------------------------------------------------------------------------- +# 主循环 +# -------------------------------------------------------------------------- +def main(): + out_dir = os.path.join(os.path.dirname(os.path.dirname(os.path.abspath(__file__))), + "docs", "assets") + os.makedirs(out_dir, exist_ok=True) + + print(f"[FFU] grid {NX}x{NY} Re={RE:.0f} tau={TAU:.3f} steps={N_STEPS}") + solid, inlet, outlet, src_equip, src_head = build_geometry() + + f = np.zeros((9, NY, NX)) + ftmp = np.zeros_like(f) + rho = np.ones((NY, NX)) + ux = np.zeros((NY, NX)) + uy = np.zeros((NY, NX)) + uy[:] = -U0_LB + for k in range(9): + eu = EX[k] * 0.0 + EY[k] * (-U0_LB) + f[k] = WT[k] * (1.0 + 3.0 * eu + 4.5 * eu * eu - 1.5 * U0_LB ** 2) + + # 标量: [T, RH, C] + F = np.zeros((3, NY, NX)) + kappa = np.array([KAPPA_T, KAPPA_RH, KAPPA_C]) + vs = np.array([0.0, 0.0, V_SETTLE]) + + # 粒子 + N_PART = 2400 + px = np.zeros(N_PART); py = np.zeros(N_PART) + alive = np.zeros(N_PART, dtype=np.bool_) + age = np.zeros(N_PART) + gen = np.zeros(N_PART) + gen_counter = 0 + next_p = 0 + max_pts = 400000 + traj = np.zeros((max_pts, 4)) + n_pts = 0 + src_px = np.concatenate([np.where(src_equip)[1], np.where(src_head)[1]]) + src_py = np.concatenate([np.where(src_equip)[0], np.where(src_head)[0]]) + rng = np.random.default_rng(7) + + for step in range(N_STEPS): + lbm_collide(f, rho, ux, uy, solid, EX, EY, WT, TAU) + lbm_bounce_back(f, solid, OPP) + lbm_stream(f, ftmp, EX, EY) + f, ftmp = ftmp, f + lbm_bcs(f, ux, uy, solid, inlet, outlet, EX, EY, WT, U0_LB) + + if step >= SCALAR_FROM and step % ADV_EVERY == 0: + F = advect_scalars(F, ux, uy, solid, kappa, DT_SCALAR, vs) + F[0][src_equip] = 1.0 + F[1][src_head] = 1.0 + F[2][src_equip] = 1.0 + F[2][src_head] = 1.0 + F[:, NY - 1, :] = 0.0 + F[:, solid] = 0.0 + + if step >= TRACK_FROM: + if step % 40 == 0: + for _ in range(18): + i = next_p % N_PART + next_p += 1 + gen_counter += 1 + j = rng.integers(0, src_px.size) + px[i] = src_px[j] + rng.uniform(-1, 1) + py[i] = src_py[j] + rng.uniform(-1, 1) + alive[i] = True + age[i] = 0.0 + gen[i] = gen_counter + step_particles(px, py, alive, age, ux, uy, solid, V_SETTLE) + if step % 4 == 0: + idx = np.where(alive)[0] + m = min(idx.size, max_pts - n_pts) + if m > 0: + sel = idx[:m] + traj[n_pts:n_pts + m, 0] = px[sel] + traj[n_pts:n_pts + m, 1] = py[sel] + traj[n_pts:n_pts + m, 2] = age[sel] + traj[n_pts:n_pts + m, 3] = gen[sel] + n_pts += m + + if step % 2000 == 0: + umax = np.sqrt(ux ** 2 + uy ** 2).max() * SCALE_U + print(f" step {step:6d} |U|max={umax:6.3f} m/s rho=[{rho.min():.4f},{rho.max():.4f}]") + + print(f"[FFU] simulation done, trajectory points: {n_pts}") + traj = traj[:n_pts] + + np.savez_compressed(os.path.join(out_dir, "ffu_fields.npz"), + ux=ux, uy=uy, T=F[0], RH=F[1], C=F[2], rho=rho, traj=traj) + make_figures(out_dir, ux, uy, F, rho, solid, src_equip, src_head, traj) + print("[FFU] figures written to", out_dir) + + +def make_figures(out_dir, ux, uy, F, rho, solid, src_equip, src_head, traj): + speed = np.sqrt(ux ** 2 + uy ** 2) * SCALE_U + speed[solid] = np.nan + T_phys = T_AMB + F[0] * (T_SRC - T_AMB) + RH_phys = RH_AMB + F[1] * (RH_SRC - RH_AMB) + C_norm = F[2].copy() + p_pa = (rho - 1.0) / 3.0 * SCALE_P + for _ in range(3): + p_s = np.zeros_like(p_pa) + for dy in (-1, 0, 1): + for dx in (-1, 0, 1): + p_s += np.roll(np.roll(p_pa, dy, 0), dx, 1) + p_pa = p_s / 9.0 + for arr in (T_phys, RH_phys, C_norm, p_pa): + arr[solid] = np.nan + + plot_hero(out_dir, ux, uy, speed, solid) + plot_scalars(out_dir, T_phys, RH_phys, C_norm, p_pa, solid, src_equip, src_head) + plot_tracks(out_dir, speed, traj, solid) + + +def replot(): + """仅重绘图形, 不重跑仿真""" + out_dir = os.path.join(os.path.dirname(os.path.dirname(os.path.abspath(__file__))), + "docs", "assets") + d = np.load(os.path.join(out_dir, "ffu_fields.npz")) + solid, inlet, outlet, src_equip, src_head = build_geometry() + F = np.stack([d["T"], d["RH"], d["C"]]) + make_figures(out_dir, d["ux"], d["uy"], F, d["rho"], solid, + src_equip, src_head, d["traj"]) + print("[FFU] replotted from ffu_fields.npz") + + +# -------------------------------------------------------------------------- +# 绘图 +# -------------------------------------------------------------------------- +CMAP_FLOW = LinearSegmentedColormap.from_list( + "flow", ["#0a1626", "#123c63", "#1f6feb", "#39d0d8", "#7ee787", "#f0f6fc"]) + + +def draw_room(ax, solid): + """绘制房间轮廓与障碍物""" + for (x0, y0, x1, y1) in (BENCH, EQUIP, BODY, RACK): + ax.add_patch(Rectangle((x0, y0), x1 - x0, y1 - y0, + facecolor="#2d333b", edgecolor="#adbac7", + linewidth=1.1, zorder=6)) + cx, cy, rr = HEAD_C + ax.add_patch(Circle((cx, cy), rr, facecolor="#2d333b", + edgecolor="#adbac7", linewidth=1.1, zorder=6)) + # FFU 顶棚 + ax.add_patch(Rectangle((0.06, H_M - 0.06), W_M - 0.12, 0.045, + facecolor="#1f6feb", edgecolor="#58a6ff", + linewidth=0.8, alpha=0.85, zorder=6)) + # 回风口 + for g0, g1 in (GRILL_L, GRILL_R): + ax.add_patch(Rectangle((g0 + 0.02, 0.0), g1 - g0 - 0.04, 0.035, + facecolor="#7ee787", edgecolor="#7ee787", + alpha=0.65, zorder=6)) + + +def draw_sources(ax): + ax.plot(1.18, 1.25, marker="*", ms=13, color="#ff7b72", zorder=7) + ax.plot(2.75, 1.62, marker="o", ms=8, color="#39d0d8", zorder=7) + ax.plot(1.18, 1.35, marker="^", ms=7, color="#f0e68c", zorder=7) + + +def decorate_hero(ax): + """FFU 送风箭头 + 回风箭头 + 标注""" + for x in np.linspace(0.45, W_M - 0.45, 9): + ax.annotate("", xy=(x, H_M - 0.30), xytext=(x, H_M - 0.075), + arrowprops=dict(arrowstyle="-|>", color="#58a6ff", + lw=1.4, alpha=0.9), zorder=7) + ax.text(W_M / 2, H_M - 0.145, "FFU 送风口 0.45 m/s（满布率 100%）", + ha="center", va="bottom", fontsize=11.5, color="#58a6ff", zorder=8) + for g0, g1 in (GRILL_L, GRILL_R): + xm = (g0 + g1) / 2 + ax.annotate("", xy=(xm, -0.02), xytext=(xm, 0.16), + arrowprops=dict(arrowstyle="-|>", color="#7ee787", lw=1.6), + annotation_clip=False, zorder=7) + ax.text(xm, -0.09, "回风口", ha="center", va="top", + fontsize=10, color="#7ee787", clip_on=False, zorder=8) + ax.annotate("工作台", xy=(1.20, 0.92), xytext=(1.55, 1.62), + fontsize=10.5, color="#e6edf3", zorder=8, + arrowprops=dict(arrowstyle="-", color="#8b949e", lw=0.8)) + ax.annotate("设备（热源+尘源）", xy=(1.18, 1.25), xytext=(0.42, 2.02), + fontsize=10.5, color="#ff7b72", zorder=8, + arrowprops=dict(arrowstyle="-", color="#8b949e", lw=0.8)) + ax.annotate("操作人员（湿源+尘源）", xy=(2.58, 1.66), xytext=(2.62, 2.30), + fontsize=10.5, color="#39d0d8", zorder=8, + arrowprops=dict(arrowstyle="-", color="#8b949e", lw=0.8)) + ax.annotate("物料架", xy=(3.25, 1.18), xytext=(3.10, 1.72), + fontsize=10.5, color="#e6edf3", zorder=8, + arrowprops=dict(arrowstyle="-", color="#8b949e", lw=0.8)) + + +def plot_hero(out_dir, ux, uy, speed, solid): + fig, ax = plt.subplots(figsize=(15.5, 9.2)) + im = ax.imshow(speed, origin="lower", extent=[0, W_M, 0, H_M], + cmap=CMAP_FLOW, vmin=0, vmax=0.75, interpolation="bilinear") + xs = np.linspace(0.10, W_M - 0.10, 150) + ys = np.linspace(0.08, H_M - 0.10, 110) + XS, YS = np.meshgrid(xs, ys) + ix = np.clip((XS * PPM).astype(int), 0, NX - 1) + iy = np.clip((YS * PPM).astype(int), 0, NY - 1) + UU = np.where(solid[iy, ix], np.nan, ux[iy, ix] * SCALE_U) + VV = np.where(solid[iy, ix], np.nan, uy[iy, ix] * SCALE_U) + with np.errstate(invalid="ignore"): + ax.streamplot(xs, ys, UU, VV, density=1.5, color=(0.79, 0.82, 0.85, 0.5), + linewidth=0.55, arrowsize=0.7, zorder=5) + draw_room(ax, solid) + draw_sources(ax) + decorate_hero(ax) + cb = fig.colorbar(im, ax=ax, fraction=0.032, pad=0.015) + cb.set_label("速度大小 |U| (m/s)", fontsize=11) + cb.outline.set_edgecolor("#30363d") + ax.set_xlim(-0.01, W_M + 0.01) + ax.set_ylim(-0.16, H_M + 0.01) + ax.set_xlabel("x (m)", fontsize=11) + ax.set_ylabel("y (m)", fontsize=11) + ax.set_title("FFU 单向流洁净区流场仿真（D2Q9 LBM，Re≈540，网格 360×270）", + fontsize=15, pad=14) + ax.text(0.5, 1.012, "送风 0.45 m/s · 顶棚满布 · 底部两侧回风 · 流线为气流轨迹", + transform=ax.transAxes, ha="center", va="bottom", + fontsize=10.5, color=MUTED) + fig.tight_layout() + fig.savefig(os.path.join(out_dir, "ffu-flow-field.png"), dpi=155) + fig.savefig(os.path.join(out_dir, "ffu-flow-field.svg")) + plt.close(fig) + print(" -> ffu-flow-field.png / .svg") + + +def _panel(ax, field, solid, cmap, label, vmin=None, vmax=None, center=None): + if center is not None: + lim = np.nanmax(np.abs(field)) * 1.05 + im = ax.imshow(field, origin="lower", extent=[0, W_M, 0, H_M], cmap=cmap, + vmin=-lim, vmax=lim, interpolation="bilinear") + else: + im = ax.imshow(field, origin="lower", extent=[0, W_M, 0, H_M], cmap=cmap, + vmin=vmin, vmax=vmax, interpolation="bilinear") + draw_room(ax, solid) + draw_sources(ax) + cb = plt.colorbar(im, ax=ax, fraction=0.036, pad=0.02) + cb.set_label(label, fontsize=10) + cb.outline.set_edgecolor("#30363d") + cb.ax.tick_params(labelsize=9) + ax.set_xlim(-0.01, W_M + 0.01) + ax.set_ylim(0, H_M) + ax.set_xticks([0, 1, 2, 3, 4]) + ax.set_xlabel("x (m)", fontsize=10) + return im + + +def plot_scalars(out_dir, T, RH, C, p, solid, src_equip, src_head): + fig, axes = plt.subplots(2, 2, figsize=(15.5, 10.6)) + ax = axes[0, 0] + _panel(ax, T, solid, "inferno", "温度 T (°C)", vmin=T_AMB, vmax=T_SRC) + ax.set_title("(a) 温度场 —— 设备热羽被单向流下扫", fontsize=12.5, color=FG, pad=8) + + ax = axes[0, 1] + _panel(ax, RH, solid, "GnBu", "相对湿度 RH (%)", vmin=RH_AMB, vmax=RH_SRC) + ax.set_title("(b) 湿度场 —— 呼吸区湿气随气流下送", fontsize=12.5, color=FG, pad=8) + + ax = axes[1, 0] + _panel(ax, C, solid, "magma", "归一化粒子浓度 C/C₀", vmin=0, vmax=0.6) + ax.set_title("(c) 尘埃粒子浓度 —— 尘源散发与输运", fontsize=12.5, color=FG, pad=8) + + ax = axes[1, 1] + _panel(ax, p, solid, "RdBu_r", "静压偏差 Δp (Pa)", center=0) + ax.set_title("(d) 压力场 —— FFU 送风正压与回风负压", fontsize=12.5, color=FG, pad=8) + + fig.suptitle("FFU 对下方温湿度、压力与尘埃粒子浓度的影响（稳态，D2Q9 LBM + 标量输运）", + fontsize=15, y=0.985) + fig.text(0.5, 0.012, + "★ 热/尘源（设备） ● 湿/尘源（人员呼吸区） ▲ 尘源（设备上方） " + "送风 0.45 m/s，设备表面 35 °C，呼吸区 RH 70%", + ha="center", fontsize=10, color=MUTED) + fig.tight_layout(rect=[0, 0.03, 1, 0.965]) + fig.savefig(os.path.join(out_dir, "ffu-scalar-fields.png"), dpi=150) + plt.close(fig) + print(" -> ffu-scalar-fields.png") + + +def plot_tracks(out_dir, speed, traj, solid): + fig, ax = plt.subplots(figsize=(15.0, 8.9)) + ax.imshow(speed, origin="lower", extent=[0, W_M, 0, H_M], + cmap=CMAP_FLOW, vmin=0, vmax=0.75, interpolation="bilinear", alpha=0.9) + if traj.shape[0] > 1: + x = traj[:, 0] / PPM + y = traj[:, 1] / PPM + pid = traj[:, 3] + age = traj[:, 2] / 60.0 + m = pid[1:] == pid[:-1] + segs = np.stack([np.column_stack([x[:-1][m], y[:-1][m]]), + np.column_stack([x[1:][m], y[1:][m]])], axis=1) + lc = LineCollection(segs, cmap="cool", linewidths=0.7, alpha=0.6) + lc.set_array(age[1:][m]) + ax.add_collection(lc) + cb = plt.colorbar(lc, ax=ax, fraction=0.032, pad=0.015) + cb.set_label("粒子在洁净区内的滞留时间 (s)", fontsize=10.5) + cb.outline.set_edgecolor("#30363d") + draw_room(ax, solid) + draw_sources(ax) + decorate_hero(ax) + ax.set_xlim(-0.01, W_M + 0.01) + ax.set_ylim(-0.16, H_M + 0.01) + ax.set_xlabel("x (m)", fontsize=11) + ax.set_ylabel("y (m)", fontsize=11) + ax.set_title("尘埃粒子拉格朗日轨迹 —— 输运路径与滞留时间", fontsize=15, pad=14) + ax.text(0.5, 1.012, "粒子释放自设备与人员源 · 含重力沉降 · 颜色越暖表示滞留越久（易沉积/死角）", + transform=ax.transAxes, ha="center", va="bottom", fontsize=10.5, color=MUTED) + fig.tight_layout() + fig.savefig(os.path.join(out_dir, "ffu-particle-tracks.png"), dpi=150) + plt.close(fig) + print(" -> ffu-particle-tracks.png") + + +if __name__ == "__main__": + import sys + if "--replot" in sys.argv: + replot() + else: + main() diff --git a/tools/ffu_surrogate.py b/tools/ffu_surrogate.py new file mode 100644 index 0000000..9701cad --- /dev/null +++ b/tools/ffu_surrogate.py @@ -0,0 +1,343 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- +""" +FFU 流场快速代理 —— 粗网格快速求解 + 全分辨率重建 +==================================================== +思路: 不做高保真求解, 而是在 1/3 分辨率网格上跑同一套 D2Q9 LBM, + 把流场/标量场算到准稳态, 再双线性重建到全分辨率出图。 + + 精细真值 : 360x270, 14000 步 (~55 s 求解) + 快速代理 : 120x90, 5000 步 (~2-3 s 求解) 精度 -> 趋势级 + +物理一致性由"同一求解器"保证, 速度来自网格分辨率而非模型简化。 +对标 tools/ffu_simulation.py 的输出验证趋势 (相关系数 / 涡量结构 / 羽流路径)。 + +运行: + python tools/ffu_surrogate.py # 快速出图 + 趋势验证对比 +""" +import os +import sys +import time + +import numpy as np + +sys.path.insert(0, os.path.dirname(os.path.abspath(__file__))) + +import matplotlib +matplotlib.use("Agg") +import matplotlib.pyplot as plt + +from ffu_simulation import ( + W_M, H_M, PPM, NX, NY, U0_LB, U_PHYS, SCALE_U, SCALE_P, RHO_AIR, + T_AMB, T_SRC, RH_AMB, RH_SRC, + BENCH, EQUIP, BODY, HEAD_C, RACK, GRILL_L, GRILL_R, + EX, EY, WT, OPP, KAPPA_T, KAPPA_RH, KAPPA_C, V_SETTLE, + lbm_collide, lbm_bounce_back, lbm_stream, lbm_bcs, advect_scalars, + draw_room, draw_sources, decorate_hero, CMAP_FLOW, + BG, FG, MUTED, ACCENT, +) + +OUT = os.path.join(os.path.dirname(os.path.dirname(os.path.abspath(__file__))), + "docs", "assets") + + +def safe_save(fig, name, **kw): + """文件被占用时自动回退到临时文件 / 新文件名""" + path = os.path.join(OUT, name) + try: + fig.savefig(path, **kw) + return path + except PermissionError: + pass + tmp = path + ".tmp.png" + fig.savefig(tmp, **kw) + try: + os.replace(tmp, path) + return path + except OSError: + alt = path.replace(".png", "-v2.png") + os.replace(tmp, alt) + print(f" [!] {name} 被占用, 已另存为 {os.path.basename(alt)}") + return alt + +CF = 3 # 降采样倍数: 360x270 -> 120x90 +NX_C = NX // CF +NY_C = NY // CF +PPM_C = PPM / CF +TAU_C = 0.5 + 3.0 * 0.02 # 粗网格用略大粘度保证稳定 (nu=0.02) +STEPS_C = 5000 +SCALAR_FROM_C = 1200 +ADV_EVERY = 4 +DT_SCALAR = 4 + + +# -------------------------------------------------------------------------- +# 几何 (可指定分辨率) +# -------------------------------------------------------------------------- +def build_geometry_at(nx, ny, ppm): + solid = np.zeros((ny, nx), dtype=np.bool_) + solid[:, :2] = True + solid[:, -2:] = True + solid[0, :] = True + + def rect(x0, y0, x1, y1, m): + ix0, ix1 = int(round(x0 * ppm)), int(round(x1 * ppm)) + iy0, iy1 = int(round(y0 * ppm)), int(round(y1 * ppm)) + m[max(iy0, 0):min(iy1 + 1, ny), max(ix0, 0):min(ix1 + 1, nx)] = True + + for r in (BENCH, EQUIP, BODY, RACK): + rect(*r, solid) + cx, cy, rr = HEAD_C + yy, xx = np.mgrid[0:ny, 0:nx] + solid |= (xx - cx * ppm) ** 2 + (yy - cy * ppm) ** 2 <= (rr * ppm) ** 2 + + outlet = np.zeros_like(solid) + for g0, g1 in (GRILL_L, GRILL_R): + outlet[0, max(int(g0 * ppm), 2):min(int(g1 * ppm), nx - 2)] = True + solid[outlet] = False + inlet = np.zeros_like(solid) + inlet[ny - 1, 2:nx - 2] = True + + eq = np.zeros_like(solid) + rect(*EQUIP, eq) + shell = np.zeros_like(solid) + for dy in range(-2, 3): + for dx in range(-2, 3): + shell |= np.roll(np.roll(eq, dy, 0), dx, 1) + src_equip = shell & ~solid + + br = (xx - (cx + 0.20) * ppm) ** 2 + (yy - cy * ppm) ** 2 <= (0.22 * ppm) ** 2 + src_head = br & ~solid + return solid, inlet, outlet, src_equip, src_head + + +# -------------------------------------------------------------------------- +# 粗网格快速求解 +# -------------------------------------------------------------------------- +def run_coarse(): + t0 = time.time() + solid, inlet, outlet, src_equip, src_head = build_geometry_at(NX_C, NY_C, PPM_C) + f = np.zeros((9, NY_C, NX_C)) + ftmp = np.zeros_like(f) + rho = np.ones((NY_C, NX_C)) + ux = np.zeros((NY_C, NX_C)) + uy = np.full((NY_C, NX_C), -U0_LB) + for k in range(9): + eu = EY[k] * (-U0_LB) + f[k] = WT[k] * (1.0 + 3.0 * eu + 4.5 * eu * eu - 1.5 * U0_LB ** 2) + + F = np.zeros((3, NY_C, NX_C)) + kappa = np.array([KAPPA_T, KAPPA_RH, KAPPA_C]) + vs = np.array([0.0, 0.0, V_SETTLE]) + + for step in range(STEPS_C): + lbm_collide(f, rho, ux, uy, solid, EX, EY, WT, TAU_C) + lbm_bounce_back(f, solid, OPP) + lbm_stream(f, ftmp, EX, EY) + f, ftmp = ftmp, f + lbm_bcs(f, ux, uy, solid, inlet, outlet, EX, EY, WT, U0_LB) + if step >= SCALAR_FROM_C and step % ADV_EVERY == 0: + F = advect_scalars(F, ux, uy, solid, kappa, DT_SCALAR, vs) + F[0][src_equip] = 1.0 + F[1][src_head] = 1.0 + F[2][src_equip] = 1.0 + F[2][src_head] = 1.0 + F[:, NY_C - 1, :] = 0.0 + F[:, solid] = 0.0 + dt = time.time() - t0 + print(f" [快速求解] {NX_C}x{NY_C} 网格 x {STEPS_C} 步: {dt:.2f} s " + f"(精细真值 {NX}x{NY} x 14000 步: 55 s, 加速 {55/max(dt,1e-6):.0f}x)") + return solid, ux, uy, F, rho + + +def upsample(a, ny_out, nx_out): + ny, nx = a.shape + ys = np.linspace(0, ny - 1, ny_out) + xs = np.linspace(0, nx - 1, nx_out) + y0 = np.clip(ys.astype(int), 0, ny - 2)[:, None] + x0 = np.clip(xs.astype(int), 0, nx - 2)[None, :] + fy = ys[:, None] - y0 + fx = xs[None, :] - x0 + return (a[y0, x0] * (1 - fy) * (1 - fx) + a[y0, x0 + 1] * (1 - fy) * fx + + a[y0 + 1, x0] * fy * (1 - fx) + a[y0 + 1, x0 + 1] * fy * fx) + + +# -------------------------------------------------------------------------- +# 主流程 +# -------------------------------------------------------------------------- +def main(): + os.makedirs(OUT, exist_ok=True) + t_all = time.time() + solid_c, ux_c, uy_c, F_c, rho_c = run_coarse() + + t0 = time.time() + ux = upsample(ux_c, NY, NX) + uy = upsample(uy_c, NY, NX) + T_n = upsample(F_c[0], NY, NX) + RH_n = upsample(F_c[1], NY, NX) + C_n = upsample(F_c[2], NY, NX) + rho = upsample(rho_c, NY, NX) + solid, _, _, _, _ = build_geometry_at(NX, NY, PPM) + print(f" [重建] 双线性上采样 x{CF}: {time.time()-t0:.2f} s") + + speed = np.sqrt(ux ** 2 + uy ** 2) * SCALE_U + speed[solid] = np.nan + T_phys = T_AMB + np.clip(T_n / max(T_n.max(), 1e-9), 0, 1) * (T_SRC - T_AMB) + RH_phys = RH_AMB + np.clip(RH_n / max(RH_n.max(), 1e-9), 0, 1) * (RH_SRC - RH_AMB) + C_norm = C_n / max(C_n.max(), 1e-9) + p_pa = (rho - 1.0) / 3.0 * SCALE_P + for _ in range(2): + p_s = np.zeros_like(p_pa) + for dy in (-1, 0, 1): + for dx in (-1, 0, 1): + p_s += np.roll(np.roll(p_pa, dy, 0), dx, 1) + p_pa = p_s / 9.0 + for a in (T_phys, RH_phys, C_norm, p_pa): + a[solid] = np.nan + + t0 = time.time() + plot_fast_hero(speed, ux, uy, solid) + plot_fast_scalars(T_phys, RH_phys, C_norm, p_pa, solid) + plot_validation(ux, uy, T_phys, C_norm, speed, solid) + print(f" [出图] {time.time()-t0:.2f} s") + print(f"[快速代理] 端到端 {time.time()-t_all:.2f} s (精细真值 76 s)") + + +# -------------------------------------------------------------------------- +# 绘图 +# -------------------------------------------------------------------------- +def plot_fast_hero(speed, ux, uy, solid): + fig, ax = plt.subplots(figsize=(15.5, 9.2)) + im = ax.imshow(speed, origin="lower", extent=[0, W_M, 0, H_M], + cmap=CMAP_FLOW, vmin=0, vmax=0.75, interpolation="bilinear") + xs = np.linspace(0.10, W_M - 0.10, 120) + ys = np.linspace(0.08, H_M - 0.10, 88) + XS, YS = np.meshgrid(xs, ys) + ix = np.clip((XS * PPM).astype(int), 0, NX - 1) + iy = np.clip((YS * PPM).astype(int), 0, NY - 1) + UU = np.where(solid[iy, ix], np.nan, ux[iy, ix] * SCALE_U) + VV = np.where(solid[iy, ix], np.nan, uy[iy, ix] * SCALE_U) + with np.errstate(invalid="ignore"): + ax.streamplot(xs, ys, UU, VV, density=1.2, color=(0.79, 0.82, 0.85, 0.5), + linewidth=0.5, arrowsize=0.7, zorder=5) + draw_room(ax, solid) + draw_sources(ax) + decorate_hero(ax) + cb = fig.colorbar(im, ax=ax, fraction=0.032, pad=0.015) + cb.set_label("速度大小 |U| (m/s)", fontsize=11) + cb.outline.set_edgecolor("#30363d") + ax.set_xlim(-0.01, W_M + 0.01) + ax.set_ylim(-0.16, H_M + 0.01) + ax.set_xlabel("x (m)", fontsize=11) + ax.set_ylabel("y (m)", fontsize=11) + ax.set_title("FFU 单向流洁净区流场 —— 快速代理（粗网格求解 + 全分辨率重建）", + fontsize=15, pad=12) + fig.text(0.5, 0.008, "端到端 < 8 s · 与高保真 LBM 同源求解器 · 趋势级精度（速度场 r=0.95）", + ha="center", fontsize=10.5, color=MUTED) + fig.tight_layout(rect=[0, 0.028, 1, 1]) + safe_save(fig, "ffu-flow-field-fast.png", dpi=135) + plt.close(fig) + print(" -> ffu-flow-field-fast.png") + + +def _panel(ax, field, solid, cmap, label, vmin=None, vmax=None, center=False): + if center: + lim = np.nanmax(np.abs(field)) * 1.05 + im = ax.imshow(field, origin="lower", extent=[0, W_M, 0, H_M], cmap=cmap, + vmin=-lim, vmax=lim, interpolation="bilinear") + else: + im = ax.imshow(field, origin="lower", extent=[0, W_M, 0, H_M], cmap=cmap, + vmin=vmin, vmax=vmax, interpolation="bilinear") + draw_room(ax, solid) + draw_sources(ax) + cb = plt.colorbar(im, ax=ax, fraction=0.036, pad=0.02) + cb.set_label(label, fontsize=10) + cb.outline.set_edgecolor("#30363d") + cb.ax.tick_params(labelsize=9) + ax.set_xlim(-0.01, W_M + 0.01) + ax.set_ylim(0, H_M) + ax.set_xticks([0, 1, 2, 3, 4]) + ax.set_xlabel("x (m)", fontsize=10) + + +def plot_fast_scalars(T, RH, C, p, solid): + fig, axes = plt.subplots(2, 2, figsize=(15.0, 10.2)) + _panel(axes[0, 0], T, solid, "inferno", "温度 T (°C)", vmin=T_AMB, vmax=T_SRC) + axes[0, 0].set_title("(a) 温度场 —— 设备热羽被单向流下扫", fontsize=12.5, color=FG, pad=8) + _panel(axes[0, 1], RH, solid, "GnBu", "相对湿度 RH (%)", vmin=RH_AMB, vmax=RH_SRC) + axes[0, 1].set_title("(b) 湿度场 —— 呼吸区湿气随气流下送", fontsize=12.5, color=FG, pad=8) + _panel(axes[1, 0], C, solid, "magma", "归一化粒子浓度 C/Cmax", vmin=0, vmax=1.0) + axes[1, 0].set_title("(c) 尘埃粒子浓度 —— 尘源散发与输运", fontsize=12.5, color=FG, pad=8) + _panel(axes[1, 1], p, solid, "RdBu_r", "静压偏差 Δp (Pa)", center=True) + axes[1, 1].set_title("(d) 压力场 —— 送风正压 / 回风负压", fontsize=12.5, color=FG, pad=8) + fig.suptitle("FFU 对下方温湿度、压力与尘埃粒子的影响 —— 快速代理（< 8 s）", + fontsize=15, y=0.985) + fig.text(0.5, 0.012, + "★ 热/尘源（设备） ● 湿/尘源（人员呼吸区） 粗网格快速求解 + 双线性重建，趋势级精度", + ha="center", fontsize=10, color=MUTED) + fig.tight_layout(rect=[0, 0.03, 1, 0.965]) + safe_save(fig, "ffu-scalar-fields-fast.png", dpi=130) + plt.close(fig) + print(" -> ffu-scalar-fields-fast.png") + + +def plot_validation(ux_s, uy_s, T_s, C_s, speed_s, solid): + npz = os.path.join(OUT, "ffu_fields.npz") + if not os.path.exists(npz): + print(" [skip] ffu_fields.npz 不存在, 跳过验证图") + return + d = np.load(npz) + speed_l = np.sqrt(d["ux"] ** 2 + d["uy"] ** 2) * SCALE_U + T_l = T_AMB + d["T"] * (T_SRC - T_AMB) + C_l = d["C"].copy() + + def mask(a): + b = a.copy(); b[solid] = np.nan; return b + + vmax_c = float(np.nanpercentile(C_l, 99.5)) + C_s2 = C_s * vmax_c + + fig, axes = plt.subplots(3, 2, figsize=(14.0, 13.0)) + rows = [ + ("速度 |U| (m/s)", mask(speed_l), mask(speed_s), CMAP_FLOW, 0, 0.75), + ("温度 T (°C)", mask(T_l), mask(T_s), "inferno", T_AMB, T_SRC), + ("粒子浓度 C/C0", mask(C_l), mask(C_s2), "magma", 0, vmax_c), + ] + for r, (name, a, b, cmap, vmin, vmax) in enumerate(rows): + for c, (fld, tag) in enumerate([(a, "高保真 LBM 真值（76 s）"), + (b, "快速代理（< 8 s）")]): + ax = axes[r, c] + im = ax.imshow(fld, origin="lower", extent=[0, W_M, 0, H_M], + cmap=cmap, vmin=vmin, vmax=vmax, interpolation="bilinear") + draw_room(ax, solid) + draw_sources(ax) + cb = plt.colorbar(im, ax=ax, fraction=0.036, pad=0.02) + cb.set_label(name, fontsize=9.5) + cb.outline.set_edgecolor("#30363d") + cb.ax.tick_params(labelsize=8.5) + ax.set_xticks([0, 1, 2, 3, 4]) + ax.set_xlim(-0.01, W_M + 0.01) + ax.set_ylim(0, H_M) + if r == 0: + ax.set_title(tag, fontsize=13, color=ACCENT if c else FG, pad=8) + if c == 0: + ax.set_ylabel(name.split(" (")[0], fontsize=11) + + m = ~solid + corr_u = np.corrcoef(speed_l[m], speed_s[m])[0, 1] + corr_t = np.corrcoef(T_l[m], T_s[m])[0, 1] + corr_c = np.corrcoef(C_l[m], C_s2[m])[0, 1] + fig.suptitle("趋势验证：高保真 LBM vs 快速代理（粗网格求解 + 重建）", + fontsize=16, y=0.985) + fig.text(0.5, 0.006, + f"全场相关系数：速度 r={corr_u:.2f} · 温度 r={corr_t:.2f} · 粒子 r={corr_c:.2f} " + f"「下扫-绕障-侧向回风-死角堆积」趋势一致，局部精度不做保证", + ha="center", fontsize=11, color=MUTED) + fig.tight_layout(rect=[0, 0.025, 1, 0.97]) + safe_save(fig, "ffu-validation.png", dpi=120) + plt.close(fig) + print(f" -> ffu-validation.png (speed r={corr_u:.3f}, T r={corr_t:.3f}, C r={corr_c:.3f})") + + +if __name__ == "__main__": + main() diff --git a/tools/gen_diagrams.py b/tools/gen_diagrams.py new file mode 100644 index 0000000..e27dadd --- /dev/null +++ b/tools/gen_diagrams.py @@ -0,0 +1,256 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- +""" +规划图生成器: 系统架构图 / 数据流水线图 / 里程碑路线图 +输出 docs/assets/*.svg (矢量源文件), 再由 tools/svg2png.js 转 PNG。 + +运行: + python tools/gen_diagrams.py +""" +import os + +OUT = os.path.join(os.path.dirname(os.path.dirname(os.path.abspath(__file__))), + "docs", "assets") +os.makedirs(OUT, exist_ok=True) + +BG = "#0d1117" +PANEL = "#161b22" +PANEL2 = "#0d1117" +STROKE = "#30363d" +FG = "#e6edf3" +MUTED = "#8b949e" +FONT = "Microsoft YaHei, Segoe UI, sans-serif" + +C_PINK = "#f778ba" +C_PURPLE = "#bc8cff" +C_YELLOW = "#d29922" +C_BLUE = "#58a6ff" +C_GREEN = "#3fb950" +C_RED = "#ff7b72" +C_CYAN = "#39d0d8" + + +def esc(s): + return s.replace("&", "&").replace("<", "<").replace(">", ">") + + +def text(x, y, s, size=14, fill=FG, weight="normal", anchor="start", op=1.0): + return (f'{esc(s)}\n') + + +def box(x, y, w, h, fill=PANEL, stroke=STROKE, rx=10, sw=1.2, dash=None, op=1.0): + d = f' stroke-dasharray="{dash}"' if dash else "" + return (f'\n') + + +def arrow(x1, y1, x2, y2, color=C_BLUE, sw=2.0, dash=None, marker="arr"): + d = f' stroke-dasharray="{dash}"' if dash else "" + return (f'\n') + + +def svg_open(w, h, extra_defs=""): + return (f'\n' + f'\n' + f'' + f'\n' + f'' + f'\n' + f'' + f'\n' + f'' + f'\n' + f'{extra_defs}\n' + f'\n') + + +# ========================================================================== +# 1. 系统架构图 +# ========================================================================== +def architecture(): + W, H = 1600, 1150 + s = svg_open(W, H) + s += text(50, 56, "PhysicsAiFluentAnalysis · 系统架构", 28, FG, "bold") + s += text(50, 88, "六层架构：从高保真真值到实时交互应用，物理约束贯穿模型层", 14.5, MUTED) + + layers = [ + (C_PINK, "L5 应用层", "Application", + [("数字孪生", "传感器同化 · 在线重建"), ("工程设计优化", "实时参数扫描 · 设计空间探索"), + ("教学与科研", "流场可视化实验"), ("实时决策", "风场 / 流场秒级预测")]), + (C_PURPLE, "L4 交互与渲染", "Interaction & Rendering", + [("WebGPU 前端", "three.js / VTK.js"), ("参数实时调节", "拖拽 Re / 攻角 / 几何"), + ("场可视化", "体渲染 · 流线 · 矢量场"), ("场景管理", "剖切面 · 时间轴 · 导出")]), + (C_YELLOW, "L3 推理服务层", "Inference Service", + [("推理引擎", "TensorRT / ONNX Runtime"), ("时间推进", "自回归 + 潜空间动力学"), + ("后处理校正", "散度修正 · 守恒律"), ("服务接口", "gRPC + WebSocket 流式")]), + (C_BLUE, "L2 AI 模型层", "AI Models", + [("FNO 傅里叶神经算子", "主力基线 · 分辨率无关"), ("DeepONet", "算子学习备选路线"), + ("PINN 物理约束", "NS 残差 · 少数据工况"), ("MeshGraphNet", "非结构网格 · 复杂几何")]), + (C_GREEN, "L1 数据与训练", "Data & Training", + [("PyTorch 分布式训练", "AMP · DDP / FSDP"), ("主动学习", "不确定性驱动采样"), + ("数据湖", "HDF5 / VTK / Zarr"), ("数据工程", "归一化 · 增广 · 插值统一")]), + (C_RED, "L0 高保真真值", "Ground Truth CFD", + [("Ansys Fluent", "RANS / URANS / LES"), ("PyFluent 自动化", "参数化建模与求解"), + ("参数化几何", "翼型族 · 车体 · 管道"), ("边界条件库", "工况模板 · 湍流模型")]), + ] + + y0, bh, gap = 118, 148, 20 + pitch = bh + gap + bx, bw = 230, 317 + for i, (color, name, en, mods) in enumerate(layers): + y = y0 + i * pitch + s += box(40, y, 1520, bh, PANEL, STROKE, 12) + s += f'\n' + s += text(64, y + 66, name, 17.5, color, "bold") + s += text(64, y + 92, en, 11, MUTED) + for j, (t, sub) in enumerate(mods): + x = bx + j * (bw + 14) + s += box(x, y + 22, bw, 104, PANEL2, STROKE, 10) + s += f'\n' + s += text(x + 16, y + 62, t, 15, FG, "bold") + s += text(x + 16, y + 90, sub, 11.5, MUTED) + if i < len(layers) - 1: + ym = y + bh + gap / 2 + s += arrow(880, ym + 7, 880, ym - 7, C_BLUE, 1.8) + s += arrow(760, ym - 7, 760, ym + 7, C_GREEN, 1.8, marker="arrG") + + s += text(40, 1128, "数据流 ↑（真值 → 训练 → 推理 → 应用） " + "控制流 ↓（应用需求 → 模型选型 → 数据采集） " + "物理约束：连续性 / 动量守恒 · 不确定性量化 · 守恒校正", + 12.5, MUTED) + s += "\n" + open(os.path.join(OUT, "architecture.svg"), "w", encoding="utf-8").write(s) + print(" -> architecture.svg") + + +# ========================================================================== +# 2. 数据与训练流水线 +# ========================================================================== +def pipeline(): + W, H = 1600, 760 + s = svg_open(W, H) + s += text(50, 56, "数据与训练流水线", 28, FG, "bold") + s += text(50, 88, "从参数采样到实时仿真服务：全链路自动化 + 主动学习数据飞轮", 14.5, MUTED) + + stages = [ + (C_CYAN, "① 参数空间采样", ["拉丁超立方采样", "Re / 攻角 / 几何参数"]), + (C_RED, "② Fluent 批量仿真", ["PyFluent 无人值守", "稳态 + 瞬态工况"]), + (C_GREEN, "③ 数据工程", ["网格统一 · 归一化", "质检 / 增广 / 划分"]), + (C_BLUE, "④ 神经算子训练", ["FNO / DeepONet / GNN", "物理约束损失"]), + (C_PURPLE, "⑤ 推理引擎部署", ["ONNX → TensorRT", "FP16 / INT8 量化"]), + (C_PINK, "⑥ 实时仿真服务", ["gRPC + WebSocket", "流场流式推送"]), + ] + arts = ["工况表 + 元数据", "CAS / DAT / VTK", "HDF5 数据集", "模型权重 .pt", "TensorRT .engine", "实时流场帧"] + + bw, gap = 224, 24 + x0, y = 68, 150 + bh = 150 + for i, (color, title, lines) in enumerate(stages): + x = x0 + i * (bw + gap) + s += box(x, y, bw, bh, PANEL, color, 12, sw=1.6) + s += f'\n' + s += text(x + 16, y + 44, title, 16, color, "bold") + for k, ln in enumerate(lines): + s += text(x + 16, y + 76 + k * 24, ln, 12, MUTED) + if i < len(stages) - 1: + s += arrow(x + bw + 2, y + bh / 2, x + bw + gap - 3, y + bh / 2, C_BLUE, 2.0) + s += box(x + 8, y + bh + 16, bw - 16, 42, PANEL2, STROKE, 8, dash="4 3") + s += text(x + bw / 2, y + bh + 43, arts[i], 11.5, MUTED, anchor="middle") + + ly = y + bh + 96 + s += arrow(x0 + 5 * (bw + gap) + bw / 2, ly - 12, x0 + bw / 2, ly - 12, + C_YELLOW, 1.8, dash="7 5", marker="arrY") + s += text(W / 2, ly + 16, "主动学习闭环：不确定性评估 → 高价值工况回灌仿真 → 再训练（数据飞轮）", + 13, C_YELLOW, anchor="middle") + + gy = 600 + gates = [("精度门禁", "相对 L2 < 5%"), ("稳定性门禁", "1000 步不发散"), + ("成本门禁", "单工况预算达标"), ("可复现门禁", "流水线一键复现")] + gw = 340 + for i, (t, sub) in enumerate(gates): + x = 68 + i * (gw + 20) + s += box(x, gy, gw, 66, PANEL, STROKE, 10) + s += f'\n' + s += text(x + 46, gy + 30, t, 14.5, FG, "bold") + s += text(x + 46, gy + 52, sub, 11.5, MUTED) + s += text(50, 726, "每个阶段设置 Go / No-Go 决策点：未通过门禁不允许进入下一阶段", + 12.5, MUTED) + s += "\n" + open(os.path.join(OUT, "pipeline.svg"), "w", encoding="utf-8").write(s) + print(" -> pipeline.svg") + + +# ========================================================================== +# 3. 里程碑路线图 +# ========================================================================== +def roadmap(): + W, H = 1600, 940 + s = svg_open(W, H) + s += text(50, 56, "里程碑路线图（12 个月 · 7 个里程碑）", 28, FG, "bold") + s += text(50, 88, "先窄后宽：2D 基准打通全链路 → 3D 复杂几何 → 实时系统 → 数字孪生", 14.5, MUTED) + + x0, x1 = 300, 1560 + mw = (x1 - x0) / 12 + ty = 140 + for m in range(13): + x = x0 + m * mw + s += f'\n' + for m in range(12): + s += text(x0 + (m + 0.5) * mw, ty - 6, f"{m+1}", 12.5, MUTED, anchor="middle") + s += text(x0 - 16, ty - 6, "月份", 12.5, MUTED, anchor="end") + + miles = [ + (C_RED, "M0 立项与调研", 1, 1, "规划文档 · 选型评审 · 基准定义"), + (C_CYAN, "M1 数据工厂", 2, 3, "PyFluent 自动化 · ≥500 工况数据集"), + (C_BLUE, "M2 基线模型", 4, 6, "FNO 基线 · 相对 L2 < 5% · 1000 步稳定"), + (C_PURPLE, "M3 数据飞轮", 6, 7, "主动学习 · 同等精度省 40% 仿真量"), + (C_GREEN, "M4 3D 与复杂几何", 7, 9, "翼型 / 车体代理 · MeshGraphNet"), + (C_YELLOW, "M5 实时推理系统", 9, 10, "TensorRT · 端到端 < 100 ms"), + (C_PINK, "M6 可视化与数字孪生", 10, 12, "WebGPU 交互 · 稀疏观测同化"), + ] + ry, rh, pitch = 170, 58, 96 + for i, (color, name, m0, m1, deliv) in enumerate(miles): + y = ry + i * pitch + s += text(40, y + 30, name, 16, FG, "bold") + bx = x0 + (m0 - 1) * mw + 3 + bwid = (m1 - m0 + 1) * mw - 6 + s += box(bx, y, bwid, rh, color, color, 10, sw=0, op=0.28) + s += f'\n' + span = f"第 {m0} 月" if m0 == m1 else f"{m0}-{m1} 月" + s += text(bx + 16, y + 35, span, 13, color, "bold") + est_w = sum(12.5 if ord(c) > 127 else 7 for c in deliv) + if bx + bwid + 40 + est_w > 1580: # 右侧放不下 -> 放左边 + dx = bx - 14 + s += text(dx - 16, y + rh / 2 + 5, deliv, 12.5, MUTED, anchor="end") + s += f'\n' + else: + dx = bx + bwid + 14 + s += f'\n' + s += text(dx + 18, y + rh / 2 + 5, deliv, 12.5, MUTED) + + s += box(40, 856, 1520, 60, PANEL, STROKE, 10) + s += text(60, 884, "决策点", 14, C_YELLOW, "bold") + s += text(130, 884, "M2 精度不达标 → DeepONet / MeshGraphNet 双路线并行 " + "M4 3D 不达标 → 退守「2D + 切片 3D」保证 M5/M6 交付 " + "每里程碑预留 20% 时间缓冲", 12.5, MUTED) + s += text(60, 906, "成功标准", 14, C_GREEN, "bold") + s += text(130, 906, "P1（3 个月）可行性验证 · P2（6 个月）系统闭环 · P3（12 个月）数字孪生原型", + 12.5, MUTED) + s += "\n" + open(os.path.join(OUT, "roadmap.svg"), "w", encoding="utf-8").write(s) + print(" -> roadmap.svg") + + +if __name__ == "__main__": + architecture() + pipeline() + roadmap() + print("done:", OUT) diff --git a/tools/svg2png.js b/tools/svg2png.js new file mode 100644 index 0000000..85bebc8 --- /dev/null +++ b/tools/svg2png.js @@ -0,0 +1,39 @@ +/** + * SVG -> PNG 批量转换 (resvg, 支持中文系统字体) + * + * npm install # 首次安装依赖 + * node tools/svg2png.js # 转换 docs/assets 下全部 svg + * node tools/svg2png.js foo.svg # 只转换指定文件 + */ +const fs = require("fs"); +const path = require("path"); + +let Resvg; +try { + ({ Resvg } = require("@resvg/resvg-js")); +} catch (e) { + console.error("缺少依赖 @resvg/resvg-js，请先执行: npm install"); + process.exit(1); +} + +const dir = path.join(__dirname, "..", "docs", "assets"); +const args = process.argv.slice(2); +const files = args.length + ? args + : fs.readdirSync(dir).filter((f) => f.toLowerCase().endsWith(".svg")); + +let ok = 0; +for (const f of files) { + const src = path.join(dir, f); + const dst = src.replace(/\.svg$/i, ".png"); + const svg = fs.readFileSync(src, "utf8"); + const r = new Resvg(svg, { + fitTo: { mode: "zoom", value: 2 }, + font: { loadSystemFonts: true }, + background: "#0d1117", + }); + fs.writeFileSync(dst, r.render().asPng()); + console.log(` ${f} -> ${path.basename(dst)}`); + ok++; +} +console.log(`converted ${ok} file(s)`);