Compare commits
1
Commits
8d17bc8014
..
dev
| Author | SHA1 | Date | |
|---|---|---|---|
  NotBigGhostandClaude Opus 5
|
dd5f519c06 |
+3
-16
@@ -40,22 +40,6 @@ docs/theory/solver_2x_sdf/out/
|
||||
# ноутбук kbc_lbm.ipynb их встраивает, и без них он нечитаем на свежем клоне
|
||||
# (перегенерация требует GPU-сервера). См. заметку о них в README ветки research.
|
||||
|
||||
# ---------------------------------------------------------------------------
|
||||
# Rust (docs/theory/2d_solver — решатель KBC на Rust)
|
||||
# ---------------------------------------------------------------------------
|
||||
# Каталог сборки cargo: объектники, зависимости, инкрементальный кэш.
|
||||
target/
|
||||
# Cargo.lock для бинарного крейта обычно коммитят (воспроизводимость прогонов) —
|
||||
# он НЕ игнорируется намеренно.
|
||||
|
||||
# Артефакты прогонов решателя: гифки, ряды, сводки. Воспроизводятся запуском.
|
||||
docs/theory/2d_solver/*.gif
|
||||
docs/theory/2d_solver/*.csv
|
||||
docs/theory/2d_solver/out/
|
||||
# Результаты валидационной кампании — гигабайты гифок и рядов. Сам список сценариев
|
||||
# (bench/scenarios.json) лежит рядом с драйверами и версионируется.
|
||||
docs/theory/2d_solver/bench/out/
|
||||
|
||||
# ---------------------------------------------------------------------------
|
||||
# Редакторы и ОС
|
||||
# ---------------------------------------------------------------------------
|
||||
@@ -69,3 +53,6 @@ docs/theory/2d_solver/bench/out/
|
||||
.DS_Store
|
||||
Thumbs.db
|
||||
desktop.ini
|
||||
|
||||
# Общие заметки Claude Code: один файл на все ветки, вне версионирования
|
||||
/CLAUDE.md
|
||||
|
||||
@@ -1,202 +0,0 @@
|
||||
# CLAUDE.md
|
||||
|
||||
This file provides guidance to Claude Code (claude.ai/code) when working with code in this repository.
|
||||
|
||||
## Project Overview
|
||||
|
||||
**SimVulcan** is a minimal **Vulkan 1.3 / C++20** 3D model editor. It renders a
|
||||
3D editor space — three reference grid planes (XY, XZ, YZ) through the origin plus
|
||||
coloured X/Y/Z axes — viewed through an orbit camera, and loads/displays `.obj`
|
||||
models with selectable display modes (solid, wireframe, solid + wireframe).
|
||||
|
||||
The C++ side runs no simulation: it is a focused rendering skeleton with an ImGui
|
||||
interface (a Viewport control panel and a Mesh load panel). The repository *also*
|
||||
carries an unrelated Python research prototype under `docs/theory/` — see
|
||||
[Research prototype](#research-prototype-docstheory) below; it is not part of the
|
||||
CMake build.
|
||||
|
||||
## Build
|
||||
|
||||
Prerequisites: **Vulkan SDK 1.3.290+** (provides `glslangValidator` for offline
|
||||
shader compilation), **CMake 3.26+**, **Ninja**, a C++20 compiler (MSVC 19.36+,
|
||||
gcc 11+, clang 14+). On macOS, Vulkan is via MoltenVK (Apple Silicon only).
|
||||
|
||||
```sh
|
||||
cmake --preset windows-msvc-release
|
||||
cmake --build --preset windows-msvc-release
|
||||
```
|
||||
|
||||
Presets: `windows-msvc-debug`, `windows-msvc-release`, `linux-gcc-release`,
|
||||
`linux-clang-release`, `macos-arm64-release` (all Ninja, one dir per preset under
|
||||
`build/<presetName>/`). The first configure fetches dependencies via FetchContent
|
||||
(GLFW, GLM, volk, vk-bootstrap, VulkanMemoryAllocator, Dear ImGui, spdlog,
|
||||
tinyobjloader, Catch2) and needs network access. `VK_NO_PROTOTYPES` is set
|
||||
project-wide; Vulkan entry points load through **volk**.
|
||||
|
||||
Warnings come from `simv_set_warnings` (`/W4 /permissive-`, or
|
||||
`-Wall -Wextra -Wpedantic -Wshadow -Wold-style-cast …`); they are **not** errors.
|
||||
`cmake/Sanitizers.cmake` defines `simv_enable_sanitizers` (Debug-only ASan/UBSan)
|
||||
but no target currently calls it — wire it in manually when chasing memory bugs.
|
||||
|
||||
## Running
|
||||
|
||||
The executable resolves SPIR-V relative to the working directory (`FindSpvPath`
|
||||
probes `spirv/<rel>` and `current_path()/spirv/<rel>`). The `SimVulcan` POST_BUILD
|
||||
step copies the compiled `spirv/` tree and `assets/` next to the executable, so
|
||||
**run from the executable's own directory** (`build/<preset>/src/app/`). `main.cpp`
|
||||
also probes several `../` ancestors for `assets/meshes`. Meshes load at runtime
|
||||
through the ImGui Mesh panel.
|
||||
|
||||
Running writes two files into the CWD: `pipeline_cache.bin` (serialised
|
||||
`VkPipelineCache`, reloaded on the next start) and ImGui's `imgui.ini`. Both are
|
||||
disposable — delete them if pipeline creation or the panel layout misbehaves.
|
||||
|
||||
`ContextOptions::enableValidation` / `enableDebugUtils` default to **true** in
|
||||
every build config, so the Khronos validation layer is requested even in Release;
|
||||
messages (error + warning severity) go through spdlog.
|
||||
|
||||
`run.bat` at the repo root is **stale**: it launches
|
||||
`build/vs2022/src/app/Release/SimVulcan.exe`, a path the Ninja presets never
|
||||
produce. Do not point users at it without fixing the path first.
|
||||
|
||||
## Tests
|
||||
|
||||
Catch2 unit tests, pure CPU/math (mesh bounds + welding). No GPU required.
|
||||
|
||||
```sh
|
||||
cmake --build --preset windows-msvc-debug --target simv_tests
|
||||
ctest --preset windows-msvc-debug
|
||||
```
|
||||
|
||||
`windows-msvc-debug` is the only preset with a `testPreset` — for the other
|
||||
configs, invoke `ctest` in `build/<preset>/` directly.
|
||||
|
||||
Single test / subset — either through CTest (each `TEST_CASE` is registered
|
||||
individually by `catch_discover_tests`):
|
||||
|
||||
```sh
|
||||
ctest --preset windows-msvc-debug -R "WeldVertices" -V
|
||||
```
|
||||
|
||||
or by running the binary with a Catch2 name or tag filter:
|
||||
|
||||
```sh
|
||||
./build/windows-msvc-debug/tests/simv_tests.exe "[decimator]"
|
||||
./build/windows-msvc-debug/tests/simv_tests.exe --list-tests
|
||||
```
|
||||
|
||||
## Architecture
|
||||
|
||||
### Library / target map
|
||||
|
||||
```
|
||||
SimVulcan (exe) → simv_core, simv_vk, simv_mesh, simv_editor
|
||||
simv_core → simv_vk (App owns Window + vk::Renderer; no Vulkan calls)
|
||||
simv_editor → simv_core, simv_mesh (Camera, input, ImGui panels; no Vulkan)
|
||||
simv_mesh → simv_core (CPU mesh only; no Vulkan)
|
||||
simv_vk → third-party (volk, vk-bootstrap, VMA, GLFW, GLM, spdlog, imgui)
|
||||
simv_shaders → glslangValidator (GLSL → SPIR-V)
|
||||
```
|
||||
|
||||
Namespaces follow directories: `simv::core`, `simv::vk`, `simv::mesh`,
|
||||
`simv::editor`. Each library exports `src/` as its include root, so includes are
|
||||
written module-qualified (`#include "mesh/Mesh.h"`, `#include "vk/Renderer.h"`).
|
||||
|
||||
### Frame loop
|
||||
|
||||
`main.cpp` creates `core::App`, which owns the `core::Window` and a
|
||||
`vk::Renderer`, then runs the loop. Each frame `Renderer::DrawFrame` calls the UI
|
||||
callback (between ImGui NewFrame/Render), then records the scene: grid + mesh into
|
||||
one dynamic-rendering pass with a depth attachment, followed by ImGui, then
|
||||
presents (2 frames in flight, sync2 submits).
|
||||
|
||||
`main.cpp` wires the editor via the UI callback: it draws the panels, applies
|
||||
mouse input to the `editor::Camera`, and pushes the resulting state into the
|
||||
renderer (`SetViewProj`, `SetRenderMode`, `SetGridVisible`). Mesh loads go through
|
||||
`MeshLoadPanel`'s callback → `Renderer::SetMeshCpu`.
|
||||
|
||||
### Mesh load path
|
||||
|
||||
`MeshLoadPanel` lists `*.obj` in the mesh directory and kicks off
|
||||
`mesh::LoadObjAsync` (worker thread, `std::future`). The future is **drained on
|
||||
the main thread** at the top of `MeshLoadPanel::Draw`, so the `OnLoaded` callback
|
||||
— and therefore the GPU upload — always runs on the render thread. Loader
|
||||
exceptions surface as the panel's status string.
|
||||
|
||||
`main.cpp`'s `OnLoaded` welds the mesh (`WeldVertices`, tolerance `1e-4`), logs
|
||||
the counts, uploads via `SetMeshCpu`, and reframes the camera to the bbox radius.
|
||||
Despite the file name, `mesh/MeshDecimator.h` implements **only** spatial-hash
|
||||
vertex welding (collapse near-duplicates, drop degenerate triangles, recompute
|
||||
bounds) — there is no LOD/decimation.
|
||||
|
||||
### Non-obvious invariants (read before editing)
|
||||
|
||||
- **All Vulkan lives in `src/vk/`.** `core/`, `mesh/`, `editor/` and `app/` make no
|
||||
Vulkan API calls. `vk::Renderer`'s public header is deliberately Vulkan-free
|
||||
(pImpl + glm/`RenderMode` only) so `App` can own it without pulling in volk. Keep
|
||||
it that way — do not leak `Vk*` types into the public interfaces of those modules.
|
||||
- **Single render pass with depth.** Scene (grid + mesh) and ImGui draw into one
|
||||
`vkCmdBeginRendering` pass that has both a colour and a `D32_SFLOAT` depth
|
||||
attachment. The ImGui backend is initialised with `depthAttachmentFormat` set so
|
||||
its pipeline matches the pass; the depth buffer is recreated with the swapchain.
|
||||
- **Wireframe needs `fillModeNonSolid`.** `MeshRenderer` builds a fill pipeline and
|
||||
a `VK_POLYGON_MODE_LINE` pipeline; the line pipeline uses a small depth bias so
|
||||
the overlay sits on top of the fill. The device feature is requested in `Context`.
|
||||
Culling is off (`VK_CULL_MODE_NONE`) — loaded models may have mixed winding.
|
||||
- **Camera is fixed on the origin.** `editor::Camera` orbits (yaw/pitch/distance)
|
||||
the world origin; loaded models are recentred there via a translate-only model
|
||||
matrix. Projection uses Vulkan clip space (`GLM_FORCE_DEPTH_ZERO_TO_ONE` + Y flip,
|
||||
isolated to `Camera.cpp`).
|
||||
- **Buffers.** `vk::GpuMesh` owns the model's vertex/index buffers (staged upload on
|
||||
the transfer queue). `GridRenderer` builds a static host-visible line buffer once.
|
||||
Both scene renderers use only a push constant — no descriptor sets.
|
||||
- **Push-constant layout is a cross-file contract.** `MeshRenderer.cpp`'s anonymous
|
||||
`MeshPC { mat4 mvp; vec4 color; }` must stay byte-identical to the
|
||||
`push_constant` block in `mesh.vert`/`mesh.frag`; `color.a` is a *flag*, not
|
||||
alpha (1 = flat-shaded, 0 = constant colour for wireframe). `GridRenderer` pushes
|
||||
a bare `mat4` (vertex stage only). Change either side and you must change both.
|
||||
- **Swapchain recreation rebuilds sync objects.** `RecreateSwapDependent` recreates
|
||||
the per-frame `imageAvailable` semaphores (a failed acquire can leave one
|
||||
signalled) *and* the per-swapchain-image `renderFinished` semaphores (the image
|
||||
count may change), then re-ensures both pipelines. Keep that ordering if you
|
||||
touch resize handling.
|
||||
- **Mesh upload stalls the device.** `SetMeshCpu`/`ClearMesh` call `WaitIdle`
|
||||
before touching `GpuMesh` — acceptable because loads are rare; do not copy that
|
||||
pattern into per-frame paths.
|
||||
|
||||
### Shaders
|
||||
|
||||
GLSL under `shaders/editor/` (`mesh.{vert,frag}`, `grid.{vert,frag}`). `simv_shaders`
|
||||
compiles each to `build/<preset>/spirv/editor/<name>.spv` targeting `vulkan1.3`,
|
||||
with `shaders/` as the `-I` root (so `#include "common/foo.glsl"` would resolve).
|
||||
Adding a shader under that globbed dir is picked up automatically
|
||||
(`CONFIGURE_DEPENDS`). `mesh.frag` reconstructs a flat normal from screen-space
|
||||
derivatives, so the vertex stream carries only positions (tight `float3`, one
|
||||
binding, one attribute).
|
||||
|
||||
## Logging
|
||||
|
||||
`simv::core::Logger` (spdlog-backed, singleton) with category-based throttling.
|
||||
Log via `LogFmt(LogCategory, LogLevel, fmt, args...)`. Categories: `Core`,
|
||||
`Vulkan`, `MeshIO`, `UI`, `Test`. Per-category level, throttle interval and
|
||||
on/off are settable at runtime (`SetMinLevel` / `SetThrottle` / `SetEnabled`).
|
||||
Some low-level Vulkan code still calls `spdlog::` directly.
|
||||
|
||||
## Research prototype (`docs/theory/`)
|
||||
|
||||
Separate from the C++ application and from the CMake build: a Python **Lattice
|
||||
Boltzmann (D2Q9, KBC-N1)** research codebase — cylinder flow with a ×2 nested AMR
|
||||
patch and SDF+Bouzidi boundaries. **GPU/CuPy only, no CPU fallback.** Its
|
||||
documentation and the running experiment log are in Russian.
|
||||
|
||||
- `docs/theory/solver_2x_sdf/` — the component-split solver (`run.py`,
|
||||
`run_blockage.py`, `run_factors.py`); its `README.md` is the authoritative
|
||||
status/experiment log and maps every file to the physics it owns.
|
||||
- `docs/theory/demos_gpu/` — demo runs that produce the notebook's figures/GIFs;
|
||||
they import physics from `solver_2x_sdf`, never duplicate it.
|
||||
- `docs/theory/kbc_lbm.ipynb` — the write-up; it embeds pre-rendered artefacts and
|
||||
does not execute the simulations.
|
||||
- `docs/origins/` — source PDFs behind the theory.
|
||||
|
||||
Do not fold this into the CMake build, and do not treat it as dead code — it is
|
||||
active research with a documented experiment history.
|
||||
Binary file not shown.
@@ -1,14 +0,0 @@
|
||||
# Каталог сборки cargo — сотни мегабайт, внутри образа он собирается заново.
|
||||
target/
|
||||
# Результаты прогонов: гифки и ряды монтируются с хоста, в образ не кладутся.
|
||||
bench/out/
|
||||
*.gif
|
||||
*.csv
|
||||
*.json
|
||||
!bench/scenarios.json
|
||||
# Мусор редакторов и ОС
|
||||
.vs/
|
||||
.vscode/
|
||||
.idea/
|
||||
.DS_Store
|
||||
Thumbs.db
|
||||
Generated
-1308
File diff suppressed because it is too large
Load Diff
@@ -1,35 +0,0 @@
|
||||
[package]
|
||||
name = "kbc2d"
|
||||
version = "0.1.0"
|
||||
edition = "2021"
|
||||
rust-version = "1.75"
|
||||
description = "Двумерный решатель LBM D2Q9 с энтропийным оператором столкновения KBC (модель D / N1)"
|
||||
|
||||
[[bin]]
|
||||
name = "kbc2d"
|
||||
path = "src/main.rs"
|
||||
|
||||
[dependencies]
|
||||
clap = { version = "4", features = ["derive"] }
|
||||
rayon = "1"
|
||||
gif = "0.13"
|
||||
bytemuck = { version = "1", features = ["derive"] }
|
||||
wgpu = { version = "22", optional = true }
|
||||
pollster = { version = "0.3", optional = true }
|
||||
|
||||
[features]
|
||||
default = ["gpu"]
|
||||
# GPU-бэкенд на wgpu (Vulkan/DX12/Metal). Отключается: --no-default-features
|
||||
gpu = ["dep:wgpu", "dep:pollster"]
|
||||
|
||||
[profile.release]
|
||||
opt-level = 3
|
||||
lto = "thin"
|
||||
codegen-units = 1
|
||||
panic = "abort"
|
||||
|
||||
# Отдельный профиль для прогонов с проверками переполнения индексов
|
||||
[profile.release-dbg]
|
||||
inherits = "release"
|
||||
debug = true
|
||||
panic = "unwind"
|
||||
@@ -1,51 +0,0 @@
|
||||
# Образ для развёртывания решателя и валидационной кампании на сервере с NVIDIA.
|
||||
#
|
||||
# Сборка и запуск:
|
||||
# docker build -t kbc2d docs/theory/2d_solver
|
||||
# docker run --rm --gpus all kbc2d --calibrate # проверить, что GPU виден
|
||||
# docker run -d --gpus all -v "$PWD/out:/work/bench/out" --name kbc2d kbc2d --resume
|
||||
#
|
||||
# ГЛАВНАЯ ТОНКОСТЬ — Vulkan в контейнере. NVIDIA Container Toolkit подкладывает внутрь
|
||||
# Vulkan-ICD (nvidia_icd.json) только если в NVIDIA_DRIVER_CAPABILITIES есть `graphics`.
|
||||
# С одним `compute` wgpu не увидит НИ ОДНОГО адаптера и решатель откажется стартовать с
|
||||
# --backend gpu. Поэтому capabilities прописаны в образе, а vulkan-tools положен внутрь,
|
||||
# чтобы первым делом можно было выполнить `vulkaninfo` и убедиться, что карта видна.
|
||||
|
||||
# ── сборка ───────────────────────────────────────────────────────────────────
|
||||
# Версия тулчейна закреплена: та же, на которой решатель собирался и проверялся.
|
||||
FROM rust:1.97-bookworm AS build
|
||||
WORKDIR /src
|
||||
|
||||
# Сначала только манифесты — тогда слой с зависимостями переиспользуется, пока они не менялись.
|
||||
COPY Cargo.toml Cargo.lock ./
|
||||
RUN mkdir src && echo 'fn main() {}' > src/main.rs \
|
||||
&& cargo build --release --locked 2>/dev/null || true \
|
||||
&& rm -rf src
|
||||
|
||||
COPY src ./src
|
||||
# touch нужен, чтобы cargo не спутал новые исходники с заглушкой из слоя выше
|
||||
RUN touch src/main.rs && cargo build --release --locked
|
||||
|
||||
# ── исполнение ───────────────────────────────────────────────────────────────
|
||||
FROM debian:bookworm-slim
|
||||
LABEL org.opencontainers.image.title="kbc2d"
|
||||
LABEL org.opencontainers.image.description="Двумерный решатель LBM D2Q9 с энтропийным столкновением KBC"
|
||||
|
||||
RUN apt-get update && apt-get install -y --no-install-recommends \
|
||||
libvulkan1 vulkan-tools python3 ca-certificates \
|
||||
&& rm -rf /var/lib/apt/lists/*
|
||||
|
||||
# Без graphics в списке возможностей Vulkan-ICD внутрь не попадёт — см. шапку.
|
||||
ENV NVIDIA_VISIBLE_DEVICES=all
|
||||
ENV NVIDIA_DRIVER_CAPABILITIES=compute,utility,graphics
|
||||
|
||||
WORKDIR /work
|
||||
COPY --from=build /src/target/release/kbc2d /work/target/release/kbc2d
|
||||
COPY bench /work/bench
|
||||
RUN chmod +x /work/bench/run_campaign.sh
|
||||
|
||||
# Кампания идёт десятки часов, поэтому по умолчанию образ НЕ запускает её: случайный
|
||||
# `docker run` покажет смету и выйдет.
|
||||
WORKDIR /work/bench
|
||||
ENTRYPOINT ["./run_campaign.sh"]
|
||||
CMD ["--dry-run"]
|
||||
@@ -1,527 +0,0 @@
|
||||
# kbc2d — двумерный решатель LBM D2Q9 с энтропийным столкновением KBC
|
||||
|
||||
Переписанный на Rust решатель обтекания тела в канале. Физика — та же, что в
|
||||
`docs/theory/solver_2x_sdf` (python/CuPy), но собранная заново: с тестами против формул
|
||||
первоисточников, двумя взаимозаменяемыми бэкендами и анимацией, привязанной к физическому
|
||||
времени потока.
|
||||
|
||||
Оператор столкновения — энтропийный KBC по работам Bösch, Chikatamarla, Karlin: *Entropic
|
||||
Multi-Relaxation Models for Simulation of Fluid Turbulence* (arXiv:1507.02509, двумерная
|
||||
реализация) и *Entropic multi-relaxation lattice Boltzmann scheme for turbulent flows* (2024,
|
||||
трёхмерная). Состав сдвиговой части выбирается ключом `--kbc-model`: `n1` — только девиатор
|
||||
{N, Π_xy} (в 2D-статье KBC D), `n2` — девиатор со следом {N, Π_xy, T} (KBC C). Статьи лежат в
|
||||
`docs/origins/`; ссылки на формулы в коде и ниже даны по нумерации 2D-статьи.
|
||||
|
||||
## Устройство
|
||||
|
||||
Пять файлов по ролям — каждый отвечает ровно за одно:
|
||||
|
||||
| файл | роль |
|
||||
|---|---|
|
||||
| `src/math.rs` | **математика решателя.** Решётка D2Q9, энтропийное равновесие в product-form, проектор на сдвиг, стабилизатор γ, столкновение, Zou–He, геометрия тел через SDF, перевод единиц, спектральная диагностика. Всё узловое и чистое; здесь же тесты против формул статьи. |
|
||||
| `src/cpu.rs` | **бэкенд под процессор.** Раскладка AoS, обход сетки, rayon, сборка Bouzidi-линков, связка уровней AMR. Физику берёт из `math`. |
|
||||
| `src/gpu.rs` | **бэкенд под видеокарту.** wgpu + WGSL (Vulkan/DX12/Metal), раскладка SoA. Физика построчно повторяет `math.rs` на f32; топологию (маски, линки, рамка патча) не дублирует, а берёт из `cpu`. |
|
||||
| `src/main.rs` | **запуск и оркестрирование.** Разбор параметров, сборка постановки, цикл по шагам, живой вывод и итоговый отчёт, выгрузка рядов в CSV. Здесь же контракт `Spec` / `StepRec` / `FieldKind`, общий для обоих бэкендов. |
|
||||
| `src/gif.rs` | **создание гифок.** Тайминг относительно физического времени, палитры, нормировка, служебная надпись, кодирование. |
|
||||
|
||||
Рядом, вне этого разделения: `bench/` — валидационная кампания (список сценариев и драйверы под
|
||||
Linux и Windows), `Dockerfile` с `docker-compose.yml` — образ для развёртывания на сервере.
|
||||
|
||||
## Сборка и запуск
|
||||
|
||||
Нужен Rust 1.75+.
|
||||
|
||||
```sh
|
||||
cargo build --release # с GPU-бэкендом
|
||||
cargo build --release --no-default-features # только CPU (без wgpu)
|
||||
cargo test --release # 29 быстрых тестов
|
||||
cargo test --release -- --include-ignored # плюс 4 эталона статьи (~18 с)
|
||||
```
|
||||
|
||||
Сборка без GPU кладёт бинарь по тому же пути, поэтому собирать её ПОСЛЕДНЕЙ нельзя: `kbc2d`
|
||||
окажется перезаписан вариантом без wgpu, и `--backend gpu` будет отказывать. Порядок —
|
||||
сначала `--no-default-features`, потом обычная.
|
||||
|
||||
На сервер удобнее ставить образом — см. [Развёртывание](#развёртывание-docker).
|
||||
|
||||
Пример: цилиндр Re=150, гифка завихренности в реальном времени.
|
||||
|
||||
```sh
|
||||
./target/release/kbc2d --shape cylinder --size 24 --re 150 \
|
||||
--nx 480 --ny 240 --steps 40000 --sponge-len 32 \
|
||||
--gif wake.gif --gif-field vorticity --verbose full
|
||||
```
|
||||
|
||||
`--help` показывает все ключи, разбитые по группам: Физика, Сетка, Тело, Время, Схема,
|
||||
Анимация, Вывод.
|
||||
|
||||
## Синхронизация анимации с физическим временем
|
||||
|
||||
Требование: гифка идёт с той же скоростью, что и настоящий поток, независимо от того, с какой
|
||||
скоростью считает машина. Реальная производительность в тайминг не входит вообще.
|
||||
|
||||
В LBM скорость самой решётки жёстко равна c = δx/δt = 1, поэтому шаг по времени однозначно
|
||||
определяется тем, какую **решёточную** скорость `u_lat` мы назначаем физическому потоку:
|
||||
|
||||
```
|
||||
δt = u_lat · δx / u_phys [с/шаг]
|
||||
шагов в секунду = u_phys / (u_lat · δx)
|
||||
задержка кадра = (шагов на кадр) · δt / playback
|
||||
```
|
||||
|
||||
**Про пример из постановки.** «30 м/с, ячейка 0.1 м ⇒ 300 шагов/с» — это арифметика
|
||||
δt = δx/u_phys, то есть `u_lat = 1`: поток проходит ровно ячейку за шаг. Формула
|
||||
воспроизводится буквально ключом `--u-lat 1.0` (тест `units_time_scaling` это проверяет), но
|
||||
физически такой режим негоден: Ma = u_lat/c_s = √3 ≈ 1.73, сверхзвук, разложение
|
||||
Чепмена–Энскога не работает. Поэтому по умолчанию `u_lat = 0.05` (Ma ≈ 0.087), и те же 30 м/с
|
||||
при ячейке 0.1 м дают 6000 шагов/с. Синхронность гифки выдерживается в обоих случаях — меняется
|
||||
только, сколько шагов приходится на кадр.
|
||||
|
||||
**Две тонкости формата GIF**, обе разобраны:
|
||||
|
||||
1. *Задержка хранится в сотых долях секунды.* Точная физическая задержка почти никогда не целая:
|
||||
30 кадр/с — это 3⅓ сотых. Покадровое округление до 3 дало бы анимацию на 11% быстрее
|
||||
реальности, и уход копился бы линейно (к тысячному кадру — 3.3 секунды). Поэтому задержки
|
||||
выдаёт накопитель `DelayDither`: суммарное время кадров отслеживает точное физическое с
|
||||
точностью до одной сотой (3, 3, 4, 3, 3, 4, …), ошибка ограничена ±5 мс и не растёт.
|
||||
2. *Меньше одной сотой не бывает.* При физически корректном `u_lat` кадр каждые 10 шагов — это
|
||||
600 кадр/с, чего формат не умеет. Поэтому режим по умолчанию `--gif-every auto` подбирает
|
||||
шаг сам под `--gif-fps` (30) так, чтобы получилось ровно реальное время. Если шаг задан
|
||||
жёстко и задержка не представима, программа печатает фактический коэффициент расхождения и
|
||||
конкретный совет, как починить.
|
||||
|
||||
`--gif-speed 0.1` даёт замедление в 10 раз (тоже точно, через тот же накопитель).
|
||||
|
||||
## Шаг — не фиксированная порция времени
|
||||
|
||||
Физическая длительность шага привязана к размеру клетки: δt = u_lat·δx/u_phys. Уменьшили клетку
|
||||
вдвое — вдвое уменьшился и δt, и то же число шагов покроет вдвое меньше физического времени.
|
||||
Отсюда `--time`: задаёте длительность в секундах, число шагов считается само.
|
||||
|
||||
Измельчение стоит дважды: клеток становится (1/δx)², а шагов на ту же секунду — 1/δx, итого
|
||||
работы ~(1/δx)³ в двумерии. И помните, что `--size` задаётся В КЛЕТКАХ: уменьшив клетку и не
|
||||
тронув `--size`, вы уменьшите тело физически.
|
||||
|
||||
Судить о длительности удобнее всего по **конвективным временам D/U** — их печатает шапка. Это
|
||||
единственная мера, не зависящая ни от сетки, ни от выбора u_lat.
|
||||
|
||||
## Что параметризовано
|
||||
|
||||
- **поток**: скорость (м/с), направление, число Рейнольдса, решёточная скорость (число Маха);
|
||||
- **постановка**: `--case channel` — обтекание тела; `taylor-green`, `shear-layer`,
|
||||
`decaying-turbulence` — периодические эталоны из статей, без тела и граничных условий;
|
||||
- **сетка**: размеры домена и размер ячейки в метрах (плотность сетки), коэффициент измельчения
|
||||
вложенного патча и его границы;
|
||||
- **тело**: восемь форм на выбор — `cylinder`, `square`, `diamond`, `ellipse`, `naca`,
|
||||
`triangle`, `plate`, `polygon` — плюс характерный размер, относительная толщина, угол атаки и
|
||||
положение. Профиль задаётся четырёхзначным кодом (`--naca 4412`, с кривизной), произвольный
|
||||
контур — списком вершин (`--poly "x,y;x,y;…"`), **несколько тел сразу** — списком
|
||||
`--bodies "cylinder:d=24,x=120,y=120; naca:d=48,x=260,y=120,a=8"` (сила считается по каждому
|
||||
телу отдельно, до четырёх);
|
||||
- **время**: длительность прогона — либо числом шагов (`--steps`), либо прямо в СЕКУНДАХ
|
||||
физического времени (`--time`, число шагов считается как time/δt); начальное поле (однородный
|
||||
поток либо покой с разгоном), длина разгона, амплитуда и длительность стартового возмущения;
|
||||
- **схема**: оператор столкновения (`kbc`/`bgk`), состав сдвиговой части (`n1`/`n2`), модель
|
||||
стенки на теле (`hrr`/`grad`/`bouzidi`/`staircase`), режим выхода, поглощающие губки перед
|
||||
выходом и после входа, бэкенд, число потоков;
|
||||
- **анимация**: файл, поле (`speed`/`vorticity`/`density`/`gamma`), палитра, масштаб, шаг кадра,
|
||||
частота, скорость воспроизведения, диапазон нормировки, усреднение k×k клеток в пиксель
|
||||
(`--gif-downsample`, без него кадр с сетки 4096×2048 неподъёмен);
|
||||
- **вывод**: период живых строк, три уровня подробности, число окон в отчёте о сходимости, CSV
|
||||
рядов с прореживанием (`--series-every`), x–t диаграмма осевой линии (`--xt`), метрики
|
||||
эталонных течений (`--case-csv`) и машиночитаемая сводка всего прогона (`--summary`).
|
||||
|
||||
## Что печатает отчёт
|
||||
|
||||
*Шапка* — вся постановка с производными величинами: δt, шагов на секунду, Ma, τ, физическая
|
||||
вязкость, блокировка канала, геометрия патча, полный план тайминга анимации.
|
||||
|
||||
*Живой вывод* — шаг, физическое время, ⟨ρ⟩, max|u|, Cd, Cl, скорость счёта и ETA; на уровне
|
||||
`full` дополнительно ⟨γ⟩ с размахом, доля вырожденных узлов, доля узлов с ξ < 0, MLUPS и
|
||||
отношение скорости счёта к реальному времени.
|
||||
|
||||
*Итог* — установившийся режим (St, ⟨Cd⟩, rms Cl, ⟨Cm⟩) сырой и с поправкой на блокировку, рядом
|
||||
литературные значения для цилиндра; таблица сходимости по окнам с вердиктом о дрейфе массы и
|
||||
насыщении; разбор стабилизатора γ; производительность.
|
||||
|
||||
## Состояние проверки
|
||||
|
||||
Опора — статьи авторов метода из `docs/origins/`, а не сторонние реализации. Ссылки на формулы
|
||||
даны по нумерации 2D-статьи (arXiv:1507.02509); там, где полезнее формулировка из работы 2024
|
||||
года по трёхмерной реализации, это отмечено отдельно.
|
||||
|
||||
### Оператор столкновения сверен с листингом статьи
|
||||
|
||||
Работа 2024 года приводит оператор явным пошаговым листингом (разд. IV). Реализация повторяет
|
||||
его дословно: ρ, u → f^eq → s и s^eq → Δs = s − s^eq → **Δh = h − h^eq = f − f^eq − Δs** →
|
||||
γ по замкнутой оценке → **f′ = f − β(2Δs + γΔh)**.
|
||||
|
||||
Проверено тестами (`cargo test`, 29 быстрых + 4 длинных):
|
||||
|
||||
- **проектор Δs** совпадает с матричным `M⁻¹·D·M` в базисе натуральных моментов (6)–(7) до 1e-13
|
||||
— для обоих составов сдвиговой части; идемпотентен, не несёт ни массы, ни импульса;
|
||||
- **γ из замкнутой оценки** (ур. 17 / ур. 25 работы 2024) — корень условия максимума энтропии
|
||||
(ур. 15 / 23): невязка при γ\* более чем в 20 раз меньше, чем при γ\*±1;
|
||||
- **при γ = 2 схема совпадает с LBGK** поточечно — как и заявлено под ур. (14);
|
||||
- **сдвиговые моменты релаксируют ровно с 2β при любой γ** (β = 0.3, 0.6, 0.95, обе модели) —
|
||||
именно это гарантирует, что стабилизатор не трогает вязкость;
|
||||
- **вязкость по ур. (5)** воспроизводится затуханием сдвиговой волны точнее 1% (τ = 0.6 и 1.0);
|
||||
- **объёмная вязкость по ур. (57)**: ξ = ν при следе в сдвиговой части и ξ = c_s²(1/(γβ) − ½)
|
||||
без него;
|
||||
- равновесие в product-form сохраняет ρ и ρu до 1e-13; Zou–He ставит ровно заданные скорость
|
||||
на входе и плотность на выходе; SDF всех восьми форм даёт верный знак и |∇φ| = 1 ± 0.05
|
||||
(у эллипса для этого пришлось считать ближайшую точку итеративно: дешёвое приближение давало
|
||||
|∇φ| = 0.57 вдали от поверхности);
|
||||
- **HRR-сборка сохраняет моменты**, ради которых затевалась: ρ, ρu и Π восстановленной функции
|
||||
распределения совпадают с целевыми, а при нулевой скорости она совпадает с Градовой (все
|
||||
коэффициенты 3-го порядка рекурсивно обращаются в ноль).
|
||||
|
||||
### Вихрь Тейлора–Грина: второй порядок сходимости (разд. VI)
|
||||
|
||||
Единственное из трёх эталонных течений статьи с ТОЧНЫМ аналитическим решением, поэтому
|
||||
проверяется не «похоже на чужой прогон», а прямое совпадение с формулой:
|
||||
|
||||
u = ∇×[(u₀/k₂)cos(k₁x)cos(k₂y)·exp(−ν(k₁²+k₂²)t)], k₁ = 1, k₂ = 4,
|
||||
область 0 < x,y < 2π на N×N, Re = u₀N/ν, полураспад t_c = ln2/[ν(k₁²+k₂²)].
|
||||
|
||||
Старт — приближением Града (ур. 58), как в статье. Метрика — как на рис. 1:
|
||||
Σ|u_x − u_x^точн| / Σ|u_x^точн| в момент t_c.
|
||||
|
||||
| N | u₀ | Re | полная | амплитуда | форма |
|
||||
|---|---|---|---|---|---|
|
||||
| 64 | 0.03 | 100 | 6.68e-3 | 6.49e-3 | 8.65e-4 |
|
||||
| 128 | 0.015 | 100 | 1.61e-3 | 1.57e-3 | 2.07e-4 |
|
||||
| 256 | 0.0075 | 100 | 4.00e-4 | 3.89e-4 | 4.36e-5 |
|
||||
|
||||
**Порядок 2.05 и 2.01** — второй порядок статьи воспроизведён, причём отдельно по амплитуде
|
||||
(скорость затухания) и по форме.
|
||||
|
||||
Два места, где пришлось разобраться, и оба поучительны:
|
||||
|
||||
1. **Давление в начальных условиях.** Течение несёт собственное поле давления порядка ρu₀²,
|
||||
находимое из ∇²p = 2ρ(ψ_xx·ψ_yy − ψ_xy²):
|
||||
`p = −(ρu₀²/4)[cos(2k₁x) + (k₁²/k₂²)cos(2k₂y)]`. Старт с ρ ≡ 1 сбрасывает эту разницу в
|
||||
акустику, которая в периодическом ящике почти не затухает и садится полкой на ошибку. Работа
|
||||
2024 года делает то же самое явно: там начальные ρ и старшие моменты получают, решая
|
||||
∂ρ/∂t + ∇·(ρu₀) = D∇²ρ до стационара.
|
||||
2. **Способ измельчения.** При фиксированном u₀ ошибка упирается в полку O(Ma²), от сетки не
|
||||
зависящую (измерено: относительная ошибка формы ∝ u₀¹·⁰⁷, то есть абсолютная ∝ Ma²).
|
||||
Второй порядок виден целиком только при диффузионном измельчении — ν фиксирована, u₀ ∝ 1/N,
|
||||
тогда Re сохраняется, а Маха падает вместе с сеткой. Это свойство слабо-сжимаемого метода,
|
||||
а не реализации: **LBGK на том же тесте даёт ту же полку** (1.39/1.37/1.36e-3 против
|
||||
0.86/1.05/1.23e-3 у KBC), что согласуется с утверждением статьи «все модели работают
|
||||
практически одинаково».
|
||||
|
||||
### Дважды периодический сдвиговый слой (разд. VII)
|
||||
|
||||
Второй эталон статьи: N = 128, Re = 30000, u₀ = 0.04, κ = 80, δ = 0.05, одно конвективное
|
||||
время. Отношение энстрофии к начальной сходится с fp64-значением 0.6035. Тест чувствителен
|
||||
именно к тому, что важно: на испорченном (абсолютном) пороге вырожденности γ тот же прогон
|
||||
даёт 0.6599, то есть +9.3%, а чистый LBGK при этих параметрах разваливается.
|
||||
|
||||
### Порог вырожденности γ
|
||||
|
||||
`GREL = 1e-8` — **относительный** порог, доля от ⟨Δ|Δ⟩, а не абсолютный. Знаменатель ⟨Δh|Δh⟩
|
||||
квадратичен по неравновесию и физически мал (~1e-7…1e-9 в развитом следе), поэтому абсолютный
|
||||
порог срабатывает на подавляющем большинстве узлов и молча подменяет γ на 2 — то есть гонит
|
||||
чистый LBGK вместо KBC. Работа 2024 года прямо об этом: γ «далеко не постоянна», её эволюция
|
||||
тесно связана с состоянием потока, и «любой MRT с γ = const не достигнет той же устойчивости».
|
||||
Доля вырожденных узлов печатается в отчёте; на исправном пороге она обязана быть ~0.
|
||||
|
||||
Единственное исключение — самый первый шаг: поле в точности равно равновесию, Δ ≡ 0, и порог
|
||||
честно срабатывает везде. На результат это не влияет: γ умножается на Δh = 0.
|
||||
|
||||
### Выбор состава сдвиговой части (табл. I)
|
||||
|
||||
Ключ `--kbc-model`:
|
||||
|
||||
- `n1` (умолчание) — s = {N, Π_xy}, только девиатор. 2D-статья: KBC D; 3D: KBC-N1.
|
||||
- `n2` — s = {N, Π_xy, T}, девиатор со следом. 2D-статья: KBC C; 3D: KBC-N2.
|
||||
|
||||
По точности они неразличимы, как и заявляет статья: на одной постановке St 0.1828 у обоих,
|
||||
⟨Cd⟩ 1.4479 против 1.4468, rms Cl 0.383 против 0.388.
|
||||
|
||||
Разница — в объёмной вязкости (ур. 57). У `n2` она фиксирована: ξ = ν. У `n1` она равна
|
||||
c_s²(1/(γβ) − ½) и, поскольку измеренная ⟨γ⟩ ≈ 1.73 < 2, в среднем оказывается примерно
|
||||
вчетверо БОЛЬШЕ ν. Отрицательной она бывает лишь в долях процента узлов. Практический вывод
|
||||
против ожидания: `n1` демпфирует продольную акустику сильнее, и в специально испорченной
|
||||
постановке (старт из покоя, губка выключена) `n1` доживает до конца с пульсацией 75% от U,
|
||||
а `n2` разваливается. Поэтому умолчание — `n1`.
|
||||
|
||||
### Модель стенки на теле: насколько она субсеточная
|
||||
|
||||
Ключ `--wall`:
|
||||
|
||||
- `hrr` (умолчание) — восстановление по целевым моментам с **рекурсивной регуляризацией**
|
||||
(Malaspinas 2015; Coreixas и др., PRE 96, 033306): ряд Эрмита продолжен до 3-го порядка, а
|
||||
коэффициенты 3-го порядка не считаются по популяциям, а выражаются через 2-й рекурсивно:
|
||||
`a₃_xxy = 2u_x·a₂_xy + u_y·a₂_xx`, `a₃_xyy = 2u_y·a₂_xy + u_x·a₂_yy`. В D2Q9 `a₃_xxx` и
|
||||
`a₃_yyy` решёткой не поддерживаются и отбрасываются; множитель 1/2c_s⁶ (а не 1/6c_s⁶) учитывает
|
||||
три перестановки индексов;
|
||||
- `grad` — то же самое с обрывом ряда на тензоре давлений: условие Града (Dorschner, Bösch,
|
||||
Chikatamarla, Boulouchos, Karlin, JFM 801 (2016), разд. 2.1 и прил. B). Задаются не популяции,
|
||||
а целевые моменты — ρ, u и Π, — после чего недостающие популяции собираются приближением
|
||||
Града (2.13);
|
||||
- `bouzidi` — интерполированный отскок: доля пересечения q входит в КАЖДУЮ восстанавливаемую
|
||||
популяцию, полинково;
|
||||
- `staircase` — простой отскок, q игнорируется. Не для счёта: это база сравнения, показывающая,
|
||||
сколько именно даёт субсеточность.
|
||||
|
||||
Целевые моменты у `hrr` и `grad` одни и те же — (B 1) и (B 3) прил. B JFM 801; отличается только
|
||||
то, до какого порядка восстанавливается функция распределения по этим моментам.
|
||||
|
||||
**Субсеточность — измеренная.** Прямой тест: сдвигаем тело внутри клетки и смотрим, насколько
|
||||
поедет Cd. У по-настоящему субсеточной границы ответ не должен зависеть от того, где тело стоит
|
||||
относительно узлов (Re = 20, D = 16, стационар, пять положений на полклетки, GPU):
|
||||
|
||||
| модель | разброс Cd | Cd |
|
||||
|---|---|---|
|
||||
| `staircase` | 0.98% | 2.482–2.506 |
|
||||
| `grad` | 0.62% | 2.449–2.464 |
|
||||
| `hrr` | 0.61% | 2.450–2.464 |
|
||||
| `bouzidi` | **0.14%** | 2.459–2.463 |
|
||||
|
||||
**HRR не улучшает разрешение геометрии и не должен** — 0.61% против 0.62% у Града. Это следует
|
||||
из устройства обеих схем: третий порядок Эрмита уточняет ВОССТАНОВЛЕНИЕ популяций по моментам, а
|
||||
положение стенки входит в моментные схемы совсем другим местом. Обе моментные схемы оказываются
|
||||
ровно между ступенькой и Bouzidi, и это тоже следует из их устройства: положение
|
||||
стенки входит туда ТОЛЬКО через целевую скорость (B 1) — одну усреднённую по узлу величину.
|
||||
Целевая плотность (B 3) — обычная сумма отскочивших и известных популяций, без q вовсе; тензор
|
||||
давлений — конечные разности по решётке, тоже без q. Плюс все недостающие популяции узла
|
||||
собираются из ОДНОГО набора моментов, так что полинковая направленность теряется. Bouzidi же
|
||||
подставляет свою q в каждую популяцию отдельно. Ступенчатой поверхность у моментных схем не
|
||||
становится, но геометрия у них разрешена заметно грубее.
|
||||
|
||||
**Сходимость по разрешению тела.** Физическая постановка фиксирована (домен 15D × 10D,
|
||||
блокировка 0.1, Re = 20), меняется только число клеток на диаметр:
|
||||
|
||||
| D | `bouzidi` | `grad` |
|
||||
|---|---|---|
|
||||
| 8 | 2.581 | 2.618 |
|
||||
| 16 | 2.529 | 2.534 |
|
||||
| 32 | **2.521** | **2.521** |
|
||||
|
||||
Обе модели состоятельны и сходятся к одному пределу с наблюдаемым порядком ≈2.7; к D = 32 они
|
||||
неразличимы. Но на грубой сетке Град заметно хуже: ошибка при D = 8 равна 0.097 против 0.060.
|
||||
Для сравнения, `staircase` при D = 16 даёт 2.69 — то есть +6.7% к пределу, тогда как обе
|
||||
субсеточные модели держатся в пределах +0.4%.
|
||||
|
||||
**Зачем тогда моментные схемы.** Их преимущество в статье — не геометрическая точность, а
|
||||
устойчивость на турбулентных режимах (авторы пишут, что интерполяционные схемы «ограничены
|
||||
низкими числами Рейнольдса, поскольку на границе возникают паразитные скачки») и естественная
|
||||
форма для подвижных стенок: скорость стенки входит в целевые значения, а не отдельной поправкой.
|
||||
В здешней канальной постановке преимущества по устойчивости воспроизвести не удалось: при росте
|
||||
Re обе модели теряют счёт на одном и том же значении (Re ≈ 5·10⁴ при теле в 16 клеток), то есть
|
||||
ограничивает не стенка, а что-то другое — вероятнее всего Zou–He при τ → ½.
|
||||
|
||||
**Умолчание — `hrr`, и это решение временное.** Оно принято по устройству схемы (третий порядок
|
||||
Эрмита фильтрует высокочастотный мусор у стенки, чего обрыв на Π не делает), а не по здешним
|
||||
измерениям: на стационарном цилиндре при Re = 20 отличить `hrr` от `grad` нельзя вовсе. Вопрос
|
||||
ставит ребром группа C кампании — там обе моментные схемы и Bouzidi гоняются на Re = 20, 150 и
|
||||
2000. Если данные не подтвердят преимущества HRR на турбулентном режиме, умолчанием станет
|
||||
`bouzidi`, у которого измеренное разрешение геометрии вчетверо лучше.
|
||||
|
||||
Все четыре модели работают на обоих бэкендах и совпадают между ними до 0.007% по Cd. Это
|
||||
специально проверяется: раньше GPU при `--wall staircase` молча считал по Bouzidi, и обнаружилось
|
||||
это только потому, что две модели дали побитово одинаковый результат там, где обязаны были
|
||||
разойтись.
|
||||
|
||||
### Паритет бэкендов и согласованность уровней
|
||||
|
||||
CPU (f64) и GPU (f32) на одной постановке совпадают до 4–5 значащих цифр шаг в шаг: ⟨ρ⟩
|
||||
1.04933 против 1.04934, Cd 2.339 против 2.338, ⟨γ⟩ 1.2645 против 1.2646. На Intel Iris Xe GPU
|
||||
даёт ≈195 MLUPS против ≈18 MLUPS у процессора.
|
||||
|
||||
Один и тот же случай с патчем ×2 и вовсе без измельчения (`--refine 1`) даёт St 0.1951 против
|
||||
0.1970 и ⟨Cd⟩ 1.956 против 1.943 — связка уровней систематики не вносит.
|
||||
|
||||
**Предел на размер сетки снят.** У GPU есть жёсткий предел `maxComputeWorkgroupsPerDimension`
|
||||
= 65535, а диспетчеризация была одномерной: при 64 узлах на рабочую группу это упирало сетку в
|
||||
4.2 миллиона узлов, то есть примерно 2048×2048. Всё, что крупнее, падало ошибкой валидации —
|
||||
не считало медленно, а не запускалось вовсе. Теперь диспетчеризация двумерная, а линейный
|
||||
индекс собирается в шейдере (`lin()`/`wlin()`); отображение «рабочая группа → узлы» при этом
|
||||
остаётся ровно линейным, поэтому редукции ничего не заметили. Проверено до 4096×4096, паритет
|
||||
с CPU не сдвинулся ни в одной цифре.
|
||||
|
||||
**Где именно кончается f32.** Прямой замер на Тейлоре–Грине, где ошибка известна точно:
|
||||
|
||||
| N | CPU, f64 | GPU, f32 |
|
||||
|---|---|---|
|
||||
| 64 | 9.83e-3 | 9.80e-3 |
|
||||
| 128 | 2.40e-3 | 3.29e-3 |
|
||||
| 256 | 5.97e-4 | 2.37e-2 |
|
||||
|
||||
Пока истинная ошибка выше ~10⁻³, f32 идёт с f64 вровень; ниже — промахивается на порядок и
|
||||
больше. Практический вывод, заложенный в кампанию: исследования сходимости считаются на CPU,
|
||||
всё остальное — на GPU. Двойной точности на GPU здесь быть не может в принципе: **в WGSL типа
|
||||
`f64` не существует**, поэтому wgpu не выразит её ни на каком железе; локальная Iris Xe вдобавок
|
||||
сообщает `shaderFloat64 = false`, а на потребительских NVIDIA f64 идёт в 1/64 от f32 — то есть
|
||||
медленнее, чем CPU.
|
||||
|
||||
Что удалось выжать вместо точности — производительность. Два изменения:
|
||||
|
||||
1. **Батчинг чтения.** Раньше после каждого шага делался `map_async` + `poll(Wait)` ради 48 байт
|
||||
статистики: на 240×120 счёт упирался в 863 шаг/с при том, что сам счёт занимал 0.27 мс из
|
||||
1.16. Теперь итоги копятся в кольце на 128 слотов, синхронизация — раз в батч: **6715 шаг/с**,
|
||||
в 7.8 раза быстрее, при неизменном пошаговом интерфейсе снаружи.
|
||||
2. **Компенсированное суммирование** (Кэхена–Ноймайера) в редукциях и в сумме сил. Наивная сумма
|
||||
по 10⁵–10⁷ узлам съедает ~log₂N бит мантиссы — именно там f32 терял основную точность.
|
||||
|
||||
### Обтекание цилиндра против литературы
|
||||
|
||||
Постановка 480×240, D = 24, Re = 150, блокировка β = D/Ny = 0.1, умолчания решателя:
|
||||
|
||||
| величина | сырое | с поправкой на блокировку | литература (безгранич. цилиндр) |
|
||||
|---|---|---|---|
|
||||
| St | 0.1828 | 0.1645 | 0.183 |
|
||||
| ⟨Cd⟩ | 1.448 | 1.173 | 1.33 |
|
||||
| rms Cl | 0.383 | 0.310 | ~0.30 |
|
||||
| ⟨Cm⟩ | 0.00002 | — | 0 (симметрия) |
|
||||
|
||||
Поправка: St×(1−β), Cd и rms Cl ×(1−β)². Сырое St и скорректированный rms Cl ложатся на
|
||||
литературу; Cd после поправки ниже на 12%. ⟨Cm⟩ ≈ 0 — контроль симметрии считывания силы.
|
||||
|
||||
**Формы тел ведут себя физично.** Прогон на каждую форму (260×130, размер 20, угол атаки 12°)
|
||||
даёт ожидаемый порядок сопротивления: профиль 0.44, эллипс 0.45, пластина 0.60, цилиндр 1.33,
|
||||
ромб 1.65, квадрат 2.38, треугольник 2.72. ⟨Cm⟩ ≈ 0 **только** у круга (−0.0002), которому
|
||||
угол атаки безразличен, а у несимметричных под углом тел он ненулевой (профиль +0.31,
|
||||
эллипс +0.15, пластина +0.13).
|
||||
|
||||
### Мелкие отличия от питоновского прототипа
|
||||
|
||||
Решатель писался заново, не как порт, но пара мест разошлась с `solver_2x_sdf` намеренно:
|
||||
внутри тела здесь не считается столкновение (эти популяции фиктивны — Bouzidi перекрывает всё,
|
||||
что могло бы прийти из тела в жидкость; побочно статистика γ собирается строго по жидкости), а
|
||||
рестрикция дополнительно пропускает узлы, у которых тонкий узел-источник лежит внутри тела.
|
||||
Процессорный бэкенд работает в f64, GPU-бэкенд — в f32.
|
||||
|
||||
## Продольная акустика канала: почему поток может «дышать»
|
||||
|
||||
Самая заметная ловушка этой постановки, и её стоит понимать до первого запуска.
|
||||
|
||||
**Граничное условие с заданной скоростью на входе акустически есть жёсткий поршень**: оно
|
||||
отражает продольные волны с коэффициентом +1. Выход по давлению — наоборот, открытый конец.
|
||||
Вместе они делают из канала четвертьволновый резонатор с пучностью давления на входе и узлом
|
||||
на выходе:
|
||||
|
||||
период основной моды = 4·Nx/c_s шагов
|
||||
затухание вязкостью ~ ν(π/2Nx)² — на длинном домене практически ноль
|
||||
|
||||
Измерено на 480×240: период пульсации ⟨ρ⟩ **3332 шага** против расчётных 4·Nx/c_s = 3325
|
||||
(0.2%), первый ноль автокорреляции на 827 шагах = ровно Nx/c_s (четверть периода). За 30000
|
||||
шагов амплитуда упала на 3.8% — то есть мода не гаснет вообще. Ничто её не подкачивает; это
|
||||
звон от старта, запертый в почти без потерь резонаторе.
|
||||
|
||||
Отсюда умолчания:
|
||||
|
||||
- **`--init uniform`** — домен сразу заполнен набегающим потоком, вход включён на полную.
|
||||
Старт из покоя (`--init rest`) разгоняет весь столб жидкости и закачивает моду; при разгоне
|
||||
за 1000 шагов против акустического пробега 831 шаг это для звука удар.
|
||||
- **губка перед выходом включена и подобрана по домену** (`nx/12`, не меньше 16 столбцов,
|
||||
с запасом до патча). Отключается `--sponge-len 0`. Вклад у неё скромный — см. таблицу ниже,
|
||||
— но она бесплатна для сил и снимает остаточную пульсацию примерно вдвое при утроении длины.
|
||||
- **`--outlet extrapolate`** — нуль-градиент поперечной скорости; жёсткий ноль отражает вихри
|
||||
дорожки обратно к телу.
|
||||
|
||||
**Что именно помогает — разделено измерением** (480×240, D=24, Re=150, 30 000 шагов):
|
||||
|
||||
| старт | губка, столбцов | пульсация u′/U | St | ⟨Cd⟩ | rms Cl |
|
||||
|---|---|---|---|---|---|
|
||||
| uniform | 0 | 1.12% | 0.1839 | 1.483 | 0.404 |
|
||||
| uniform | 40 | 1.03% | 0.1840 | 1.483 | 0.404 |
|
||||
| uniform | 120 | 0.83% | 0.1839 | 1.482 | 0.403 |
|
||||
| rest | 0 | **74.91%** | 0.1534 | 2.009 | 0.807 |
|
||||
| rest | 40 | **74.90%** | 0.1539 | 2.010 | 0.808 |
|
||||
| rest | 120 | **74.99%** | 0.1542 | 2.008 | 0.808 |
|
||||
|
||||
Читается однозначно: **весь эффект даёт однородный старт**, 74.9% → 1.12%, причём при
|
||||
полностью выключенной губке. Губка снимает только остаток — 1.12% → 1.03% → 0.83%, — и на
|
||||
силы с частотой схода не влияет вовсе (St 0.1839 во всех трёх строках).
|
||||
|
||||
**При старте из покоя губка не помогает совсем.** Это не осечка реализации, а свойство моды:
|
||||
губка поднимает вязкость на последних столбцах, а там у стоячей четвертьволновой моды **узел
|
||||
давления и пучность скорости** — то самое место, где повышенная вязкость её почти не трогает.
|
||||
Бегущую волну такая губка съедает, стоячую — нет. Убрать моду можно только не возбуждая её.
|
||||
|
||||
Заодно видно, ЧЕМ платит неверный старт: St 0.153 вместо 0.184, ⟨Cd⟩ 2.01 вместо 1.48,
|
||||
rms Cl 0.81 вместо 0.40 — то есть 75-процентная продольная пульсация ломает не косметику, а
|
||||
все три величины, ради которых постановка и считается.
|
||||
|
||||
Отчёт печатает период моды, время её вязкого затухания и измеренную пульсацию в конце прогона,
|
||||
с явным предупреждением, если она превысила 5% от U.
|
||||
|
||||
**Инструмент разбора — `--xt`.** Каждые `--xt-every` шагов пишется срез ⟨ρ⟩ и u_x вдоль осевой
|
||||
линии, по строке на срез. Наклон полос на такой диаграмме прямо даёт скорость распространения:
|
||||
звук (±c_s) или конвекция (U). Стоячие узлы видны как вертикальные линии постоянной фазы, и это
|
||||
сразу отличает резонанс от неустойчивости самого граничного условия, привязанной к столбцу x = 0
|
||||
и никуда не бегущей.
|
||||
|
||||
Есть и **губка после входа** — `--sponge-in <столбцов>`, по умолчанию выключена. Гасит продольные
|
||||
волны до того, как они отразятся от входа-поршня. Осмысленна на высоких Re, где стартовая волна
|
||||
перестаёт быть безобидной полоской; цена — искажение профиля прямо на входе, поэтому включать её
|
||||
надо осознанно, а не «на всякий случай».
|
||||
|
||||
Что **не помогает** и оставлено только для повторной проверки — `--init-taper`. Замерено:
|
||||
сглаживание стартовой скорости у тела давит возмущение плотности на первом шаге восьмикратно
|
||||
(2.54·10⁻² → 3.13·10⁻³), но пик ЗА ПРОГОН при этом даже подрастает (до 3.22·10⁻²). Возмущение
|
||||
просто переносится во времени: поток всё равно обязан разогнаться вокруг тела, и энергия этого
|
||||
переходного процесса задана физикой, а не гладкостью начального поля. Умолчание — 0.
|
||||
|
||||
## Валидационная кампания
|
||||
|
||||
`bench/` — 115 прогонов на ≈90 часов GPU, разложенных по девяти группам: эталоны первоисточников,
|
||||
цилиндр против литературы, модели стенки, профили крыла, сложная и множественная геометрия,
|
||||
границы домена, старт и время жизни, внутренние инварианты, сверхмелкие сетки до 4096×2048.
|
||||
Каждый прогон кладёт логи, ряды, машиночитаемую сводку и гифку на всю свою длительность в
|
||||
собственную папку.
|
||||
|
||||
```sh
|
||||
cd bench
|
||||
python preflight.py # каждый сценарий стартует на два шага: ловит опечатки
|
||||
./run_campaign.sh --calibrate # замерить MLUPS этой машины: оценки в часах иначе гадание
|
||||
./run_campaign.sh --dry-run # смета: что, сколько шагов, сколько часов
|
||||
./run_campaign.sh --resume # считать, пропуская уже готовое
|
||||
```
|
||||
|
||||
Предполётную проверку стоит гонять всерьёз: она поймала, что вся группа сверхмелких сеток
|
||||
падала на пределе GPU (см. ниже), а девять прогонов передавали `--body-x` дважды. Оба отказа
|
||||
проявились бы только на сервере, часов через двадцать после старта кампании.
|
||||
|
||||
Подробности — в `bench/README.md`: раскладка выходных файлов, таблица групп, модель стоимости и
|
||||
список того, что известно заранее (какие прогоны обязаны развалиться и почему).
|
||||
|
||||
## Развёртывание (Docker)
|
||||
|
||||
```sh
|
||||
docker build -t kbc2d docs/theory/2d_solver
|
||||
docker run --rm --gpus all kbc2d --calibrate # проверить, что GPU виден
|
||||
docker run -d --gpus all -v "$PWD/out:/work/bench/out" kbc2d --resume
|
||||
```
|
||||
|
||||
`ENTRYPOINT` — драйвер кампании, `CMD` по умолчанию `--dry-run`: случайный `docker run` покажет
|
||||
смету и выйдет, а не запустит сточасовую задачу.
|
||||
|
||||
**Главная тонкость — Vulkan внутри контейнера.** NVIDIA Container Toolkit подкладывает
|
||||
Vulkan-ICD (`nvidia_icd.json`) только если в `NVIDIA_DRIVER_CAPABILITIES` есть `graphics`;
|
||||
с одним `compute` wgpu не увидит ни одного адаптера. В образе это прописано, но может быть
|
||||
переопределено снаружи, поэтому `vulkan-tools` лежит внутрь: первым делом на сервере стоит
|
||||
выполнить `docker run --rm --gpus all --entrypoint vulkaninfo kbc2d --summary`.
|
||||
|
||||
Проверено локально: образ собирается, кампания внутри него проходит смоук с монтированием
|
||||
результатов на хост, физика совпадает с хостовой до последней цифры (ошибка Тейлора–Грина
|
||||
9.829e-3 при N=64 и 2.401e-3 при N=128 — те же значения, что вне контейнера), а `--backend gpu`
|
||||
без проброшенной карты отказывает явным сообщением, а не считает молча. GPU-путь в контейнере
|
||||
проверяется только на машине с картой.
|
||||
|
||||
## Дальше
|
||||
|
||||
- σ·n-кросс-чек силы (интеграл тензора напряжений по контуру) как независимая проверка GMEM;
|
||||
- согласованные начальные условия по образцу работы 2024 года: там ρ и старшие моменты
|
||||
получают, решая ∂ρ/∂t + ∇·(ρu₀) = D∇²ρ до стационара, что убрало бы и остаточный
|
||||
стартовый импульс от появления тела в потоке;
|
||||
- подвижные и вращающиеся тела: GMEM уже записан в галилей-инвариантной форме и принимает
|
||||
скорость стенки на линке, но подача этой скорости не подключена;
|
||||
- разбор нерешённых 12% по Cd и предела устойчивости Re ≈ 5·10⁴ — обе задачи вынесены в
|
||||
кампанию (группы B и G), выводы делать по её данным;
|
||||
- больше четырёх тел в домене: сейчас сила считается по четырём вёдрам (`MAX_BODY_BUCKETS`),
|
||||
геометрия при этом собирается из любого числа тел, но силы сверх четвёртого сливаются вместе.
|
||||
@@ -1,113 +0,0 @@
|
||||
# Валидационная кампания
|
||||
|
||||
115 прогонов, ≈90 часов на RTX 4070 Ti. Проверяет решатель по трём независимым линиям:
|
||||
эталонам из статей авторов метода, литературе по обтеканию тел и внутренним инвариантам самой
|
||||
схемы. Каждый прогон кладёт логи, ряды, машиночитаемую сводку и гифку в собственную папку.
|
||||
|
||||
## Быстрый старт на сервере
|
||||
|
||||
```sh
|
||||
docker build -t kbc2d docs/theory/2d_solver
|
||||
docker run --rm --gpus all --entrypoint python3 kbc2d preflight.py # сценарии стартуют
|
||||
docker run --rm --gpus all kbc2d --calibrate # GPU виден
|
||||
docker run -d --gpus all -v "$PWD/out:/work/bench/out" --name kbc2d kbc2d --resume
|
||||
docker logs -f kbc2d
|
||||
```
|
||||
|
||||
Первым делом стоит убедиться, что Vulkan внутри контейнера действительно видит карту:
|
||||
|
||||
```sh
|
||||
docker run --rm --gpus all --entrypoint vulkaninfo kbc2d --summary | head -20
|
||||
```
|
||||
|
||||
Если адаптеров ноль — почти наверняка дело в `NVIDIA_DRIVER_CAPABILITIES`: Container Toolkit
|
||||
подкладывает Vulkan-ICD только при наличии `graphics` в списке. В образе это прописано, но
|
||||
может быть переопределено снаружи.
|
||||
|
||||
## Без Docker
|
||||
|
||||
```sh
|
||||
cargo build --release # из docs/theory/2d_solver
|
||||
cd bench
|
||||
python preflight.py # каждый сценарий стартует на два шага
|
||||
./run_campaign.sh --calibrate # замерить MLUPS этой машины
|
||||
KBC2D_GPU_MLUPS=1400 ./run_campaign.sh --dry-run
|
||||
./run_campaign.sh --resume
|
||||
```
|
||||
|
||||
Под Windows то же самое делает `run_campaign.ps1` (тот же `scenarios.json`); он нужен для
|
||||
локальной отладки обвязки, целевая машина — Linux.
|
||||
|
||||
## Ключи драйвера
|
||||
|
||||
| ключ | что делает |
|
||||
|---|---|
|
||||
| `--dry-run` | печатает смету: что, сколько шагов, сколько часов. Ничего не считает |
|
||||
| `--calibrate` | три коротких прогона, замер фактических MLUPS этой машины |
|
||||
| `--smoke` | три самых дешёвых прогона: проверить обвязку, а не физику |
|
||||
| `--resume` | пропускает прогоны, у которых уже есть `summary.json` |
|
||||
| `--group A,B` | только выбранные группы |
|
||||
| `--only cyl_re150` | по подстроке идентификатора |
|
||||
| `--budget-hours 24` | остановиться, когда время выйдет |
|
||||
|
||||
Прогон, который упал или развалился, помечается в сводке и **не останавливает кампанию**:
|
||||
группы E и G специально ищут предел устойчивости, там развал — ожидаемый результат.
|
||||
|
||||
## Что где лежит
|
||||
|
||||
```
|
||||
out/summary.csv сводная таблица: id, группа, статус, секунды, стоимость
|
||||
out/campaign.log журнал с отметками времени
|
||||
out/<id>/cmd.txt точная команда, которой прогон был запущен
|
||||
out/<id>/log.txt полный вывод
|
||||
out/<id>/report.txt только итоговый отчёт
|
||||
out/<id>/series.csv временные ряды по шагам
|
||||
out/<id>/summary.json ключевые метрики машиночитаемо — с этого удобно начинать разбор
|
||||
out/<id>/*.gif анимация
|
||||
out/<id>/*_case.csv энергия, энстрофия, палинстрофия (эталонные течения)
|
||||
out/<id>/xt_*.csv x–t диаграммы (группа G)
|
||||
```
|
||||
|
||||
## Группы
|
||||
|
||||
| | прогонов | что проверяется |
|
||||
|---|---|---|
|
||||
| **A** | 12 | эталоны первоисточников: Тейлор–Грин (второй порядок сходимости; при фиксированном u₀ — полка O(Ma²)), сдвиговый слой Re=3·10⁴, затухающая турбулентность |
|
||||
| **B** | 16 | цилиндр против литературы: стационар Re=20/40, дорожка Re=100…300, разрешение D=16…128, блокировка Ny/D=6…32 с экстраполяцией к бесконечной среде |
|
||||
| **C** | 20 | модели стенки и субсеточность: hrr / grad / bouzidi / staircase при Re=150 и 2000, плюс сдвиг тела внутри клетки на 0 / ¼ / ½ для всех четырёх |
|
||||
| **D** | 14 | профили крыла: поляра NACA 0012, Re-серия, изгиб 4412, сходимость по хорде |
|
||||
| **E** | 12 | сложная и множественная геометрия: тандем, решётка, перфорация, многоэлементный профиль, сверхтонкая пластина, клин, зазубренная кромка |
|
||||
| **F** | 9 | поле влияния и границы домена: отступы до входа и выхода, вложенный патч |
|
||||
| **G** | 14 | старт, акустика, время жизни: x–t диаграммы, губки, предел по Re, прогон на 10⁷ шагов |
|
||||
| **H** | 9 | инварианты: симметрия, зеркальность, зависимость от числа Маха, расхождение f32 против f64 |
|
||||
| **I** | 9 | сверхмелкие сетки 4096×2048 по всем формам и композиции из трёх тел |
|
||||
|
||||
## Стоимость и время
|
||||
|
||||
Стоимость каждого прогона хранится в **обновлениях узлов** (`nodes_per_step × steps`) — это
|
||||
единственная мера, переносимая между машинами. Часы драйвер получает, поделив её на MLUPS.
|
||||
Оценки по умолчанию исходят из 1200 MLUPS на GPU и 22 на CPU; **`--calibrate` обязателен**,
|
||||
потому что эти числа взяты с другой машины.
|
||||
|
||||
Пересобрать список с другой длительностью:
|
||||
|
||||
```sh
|
||||
python gen_scenarios.py --scale 3 # все прогоны в полтора раза длиннее (≈135 ч)
|
||||
```
|
||||
|
||||
## Что известно заранее
|
||||
|
||||
- **Точность f32.** Замерено на Тейлоре–Грине: GPU совпадает с CPU/f64, пока истинная ошибка
|
||||
выше ~10⁻³, и промахивается в 40 раз, когда она ниже. Поэтому исследования сходимости из
|
||||
группы A идут на CPU — это указано в самих сценариях.
|
||||
- **Развалы ожидаемы** в группе G (предел по Re) и у прогона `A09_shear_n512_lbgk`: LBGK при
|
||||
Re=3·10⁴ обязан развалиться там, где KBC доживает — это и есть проверяемое утверждение.
|
||||
- **Гифки** пишутся на всю длительность прогона, в реальном времени, 10 кадр/с. Число кадров
|
||||
этим задано жёстко (у самого длинного прогона их 16 666), поэтому единственный рычаг —
|
||||
размер кадра. Замерено: 0.103 байта на пиксель после LZW; отсюда бюджет
|
||||
кадры×пиксели ≤ 1.5·10⁹ на гифку, но ширина не опускается ниже 480 пикселей. Итог по
|
||||
кампании: **≈4.4 ГБ**, самый тяжёлый файл 226 МБ.
|
||||
- **Перед запуском** имеет смысл прогнать предполётную проверку: каждый сценарий стартует на
|
||||
два шага, что ловит опечатки в ключах и несовместимые сочетания до того, как кампания уйдёт
|
||||
считать на девяносто часов. Именно она поймала, что вся группа I падала на пределе GPU в
|
||||
65535 рабочих групп на измерение.
|
||||
@@ -1,518 +0,0 @@
|
||||
#!/usr/bin/env python3
|
||||
# Генератор списка прогонов кампании -> scenarios.json
|
||||
#
|
||||
# Список порождается кодом, а не правится руками: серии по Re, по разрешению и по углу атаки
|
||||
# получаются циклами, а стоимость каждого прогона считается тут же и суммируется. Запуск:
|
||||
# python gen_scenarios.py # перезаписать scenarios.json и напечатать смету
|
||||
# python gen_scenarios.py --scale 2 # растянуть все длительности вдвое
|
||||
#
|
||||
# Стоимость меряется в ОБНОВЛЕНИЯХ УЗЛОВ (nodes_per_step * steps) — единственная переносимая
|
||||
# между машинами мера. Часы драйвер получает из неё, поделив на фактические MLUPS, которые
|
||||
# замеряет на месте (--calibrate).
|
||||
|
||||
import argparse
|
||||
import json
|
||||
import math
|
||||
import os
|
||||
|
||||
RUNS = []
|
||||
U_LAT = 0.05 # решёточная скорость по умолчанию
|
||||
DX = 0.1 # размер клетки, м (умолчание решателя)
|
||||
U_PHYS = 30.0 # скорость потока, м/с (умолчание решателя)
|
||||
GIF_MAX_W = 1280 # предельная ширина кадра после прореживания
|
||||
GIF_MIN_W = 480 # ниже не опускаемся: нечитаемая гифка бесполезнее большой
|
||||
GIF_PX_BUDGET = 1.5e9 # кадры×пиксели на гифку; при 0.103 байта на пиксель это ~150 МБ
|
||||
|
||||
|
||||
def conv_steps(d_cells, n_conv, u_lat=U_LAT):
|
||||
"""Сколько шагов нужно на n_conv конвективных времён D/U."""
|
||||
return int(round(n_conv * d_cells / u_lat))
|
||||
|
||||
|
||||
def add(rid, group, title, expect, args, nodes, steps, backend="gpu",
|
||||
gif=None, nx=None, extra_out=()):
|
||||
"""Записать прогон. `nodes` — обновлений узлов за шаг, `steps` — сколько шагов."""
|
||||
a = list(args) + ["--steps", str(steps), "--backend", backend]
|
||||
# ряды прореживаем так, чтобы их осталось порядка 200 тысяч записей
|
||||
if steps > 400_000:
|
||||
a += ["--series-every", str(max(1, steps // 200_000))]
|
||||
if "case_csv" in extra_out:
|
||||
a += ["--case-csv", f"{rid}_case.csv"]
|
||||
if gif:
|
||||
# Гифка идёт на всю длительность прогона в реальном времени, 10 кадр/с — число кадров
|
||||
# задано физикой прогона и не обсуждается. Единственный рычаг — пиксели в кадре, и он
|
||||
# нужен: замерено 0.103 байта на пиксель после LZW, так что прогон на 16 тысяч кадров
|
||||
# при 960×576 весит 900 МБ, а такой файл уже нечем открыть. Поэтому кроме ограничения
|
||||
# по ширине действует бюджет ПРОИЗВЕДЕНИЯ кадры×пиксели, но ширина не опускается ниже
|
||||
# GIF_MIN_W: нечитаемая гифка бесполезнее большой.
|
||||
ny = max(1, int(nodes) // max(1, nx or 1))
|
||||
u_lat = float(a[a.index("--u-lat") + 1]) if "--u-lat" in a else U_LAT
|
||||
frames = max(1, int(steps * (u_lat * DX / U_PHYS) * 10))
|
||||
down = max(1, math.ceil((nx or 0) / GIF_MAX_W))
|
||||
while nx and frames * (nx // down) * (ny // down) > GIF_PX_BUDGET \
|
||||
and nx // (down + 1) >= GIF_MIN_W:
|
||||
down += 1
|
||||
a += ["--gif", gif, "--gif-field", "vorticity", "--gif-every", "auto",
|
||||
"--gif-fps", "10", "--gif-speed", "1",
|
||||
"--gif-downsample", str(down), "--gif-scale", "1"]
|
||||
RUNS.append({
|
||||
"id": rid, "group": group, "title": title, "expect": expect,
|
||||
"backend": backend, "args": a,
|
||||
"nodes_per_step": int(nodes), "steps": int(steps),
|
||||
"cost": int(nodes) * int(steps),
|
||||
"outputs": list(extra_out),
|
||||
})
|
||||
|
||||
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
# A. Эталоны первоисточников: периодические течения, тела и ГУ нет вовсе
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
def group_a(scale):
|
||||
NU = 0.0192
|
||||
for n in (64, 128, 256, 512):
|
||||
# диффузионное измельчение: nu фиксирована, u0 ~ 1/N, поэтому Re сохраняется, а
|
||||
# число Маха падает вместе с сеткой — только так виден второй порядок целиком
|
||||
u0 = 0.03 * 64 / n
|
||||
tc = int(round(math.log(2.0) / (NU * (2 * math.pi / n) ** 2 * 17)))
|
||||
add(f"A{len(RUNS)+1:02d}_tg_diff_n{n}", "A",
|
||||
f"Тейлор–Грин, диффузионное измельчение, N={n} (CPU/f64)",
|
||||
"порядок сходимости 2 по ряду N=64..512; ошибка на N=512 ниже 2e-4",
|
||||
["--case", "taylor-green", "--nx", str(n), "--ny", str(n), "--refine", "1",
|
||||
"--u-lat", f"{u0:.6f}", "--re", "100", "--case-every", str(max(1, tc // 20))],
|
||||
n * n, tc, backend="cpu", extra_out=("case_csv",))
|
||||
|
||||
for n in (128, 256):
|
||||
# при ФИКСИРОВАННОМ u0 ошибка упирается в полку O(Ma^2) и от сетки не зависит
|
||||
tc = int(round(math.log(2.0) / (NU * (2 * math.pi / n) ** 2 * 17)))
|
||||
add(f"A{len(RUNS)+1:02d}_tg_fixed_n{n}", "A",
|
||||
f"Тейлор–Грин при фиксированном u₀, N={n} — демонстрация полки O(Ma²)",
|
||||
"ошибка перестаёт падать с ростом N: это свойство слабо-сжимаемого метода",
|
||||
["--case", "taylor-green", "--nx", str(n), "--ny", str(n), "--refine", "1",
|
||||
"--u-lat", "0.03", "--re", f"{0.03 * n / NU:.1f}",
|
||||
"--case-every", str(max(1, tc // 20))],
|
||||
n * n, tc, backend="cpu", extra_out=("case_csv",))
|
||||
|
||||
for n in (512,):
|
||||
for op, km, tag in (("kbc", "n1", "kbc_n1"), ("kbc", "n2", "kbc_n2"), ("bgk", "n1", "lbgk")):
|
||||
steps = int(4 * n / 0.04 * scale)
|
||||
add(f"A{len(RUNS)+1:02d}_shear_n{n}_{tag}", "A",
|
||||
f"Сдвиговый слой Re=3·10⁴, N={n}, {tag}",
|
||||
"KBC доживает до конца, LBGK обязан развалиться (разд. VII статьи)",
|
||||
["--case", "shear-layer", "--nx", str(n), "--ny", str(n), "--refine", "1",
|
||||
"--u-lat", "0.04", "--re", "30000", "--collision", op, "--kbc-model", km,
|
||||
"--case-every", str(max(1, steps // 200))],
|
||||
n * n, steps, gif=f"shear_{tag}.gif", nx=n, extra_out=("case_csv",))
|
||||
|
||||
for n, tag, ops in ((2048, "n2048", ("kbc", "bgk")), (4096, "n4096", ("kbc",))):
|
||||
for op in ops:
|
||||
steps = int(200_000 * scale * (2048 / n))
|
||||
add(f"A{len(RUNS)+1:02d}_turb_{tag}_{op}", "A",
|
||||
f"Затухающая турбулентность N={n}, {op}",
|
||||
"энстрофия и палинстрофия падают монотонно; KBC держится там, где LBGK нет",
|
||||
["--case", "decaying-turbulence", "--nx", str(n), "--ny", str(n), "--refine", "1",
|
||||
"--u-lat", "0.02", "--re", f"{0.02 * n / 8.16e-5:.0f}", "--collision", op,
|
||||
"--case-every", str(max(1, steps // 300))],
|
||||
n * n, steps, gif=f"turb_{tag}_{op}.gif", nx=n, extra_out=("case_csv",))
|
||||
|
||||
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
# B. Цилиндр против литературы
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
def flatten(d):
|
||||
"""Словарь ключей в плоский argv. Именно словарь, а не список: clap отвергает
|
||||
повторённый ключ, а сценарии сплошь и рядом переопределяют умолчания сборщика."""
|
||||
a = []
|
||||
for k, v in d.items():
|
||||
a += [k, str(v)]
|
||||
return a
|
||||
|
||||
|
||||
def cyl_args(d, re, nx, ny, **kw):
|
||||
a = {"--shape": "cylinder", "--size": d, "--nx": nx, "--ny": ny,
|
||||
"--body-x": nx // 5, "--re": re, "--refine": 1}
|
||||
a.update({"--" + k.replace("_", "-"): v for k, v in kw.items()})
|
||||
return flatten(a)
|
||||
|
||||
|
||||
def group_b(scale):
|
||||
# стационарные режимы: дорожки нет, Cd — одно число
|
||||
for re in (20, 40):
|
||||
d, nx, ny = 48, 960, 576
|
||||
st = int(conv_steps(d, 60) * scale)
|
||||
add(f"B{len(RUNS)-len(RUNS)+len([r for r in RUNS if r['group']=='B'])+1:02d}_cyl_re{re}",
|
||||
"B", f"Цилиндр Re={re}, стационар, D=48",
|
||||
"Cd выходит на постоянную; сравнение с литературой по стационарному обтеканию",
|
||||
cyl_args(d, re, nx, ny, pert_amp=0), nx * ny, st,
|
||||
gif=f"cyl_re{re}.gif", nx=nx)
|
||||
|
||||
# дорожка Кармана: длинное осреднение, это главный заход на литературу
|
||||
for re in (100, 150, 200, 300):
|
||||
d, nx, ny = 64, 1280, 768
|
||||
st = int(conv_steps(d, 1500) * scale)
|
||||
add(f"B{len([r for r in RUNS if r['group']=='B'])+1:02d}_cyl_re{re}", "B",
|
||||
f"Цилиндр Re={re}, дорожка Кармана, D=64, 1500 конв. времён",
|
||||
"St≈0.183 и ⟨Cd⟩≈1.33 при Re=150 после поправки на блокировку",
|
||||
cyl_args(d, re, nx, ny), nx * ny, st, gif=f"cyl_re{re}.gif", nx=nx)
|
||||
|
||||
# разрешение тела: сходимость Cd
|
||||
for d in (16, 32, 64, 128):
|
||||
nx, ny = 20 * d, 12 * d
|
||||
st = int(conv_steps(d, 800) * scale)
|
||||
add(f"B{len([r for r in RUNS if r['group']=='B'])+1:02d}_cyl_d{d}", "B",
|
||||
f"Цилиндр Re=150, разрешение D={d}",
|
||||
"Cd сходится по D; экстраполяция даёт сеточно-независимое значение",
|
||||
cyl_args(d, 150, nx, ny), nx * ny, st, gif=f"cyl_d{d}.gif", nx=nx)
|
||||
|
||||
# блокировка канала: экстраполяция к бесконечной среде — прямой заход на нерешённые 12%
|
||||
for k in (6, 8, 12, 16, 24, 32):
|
||||
d, nx, ny = 48, 960, k * 48
|
||||
st = int(conv_steps(d, 1000) * scale)
|
||||
add(f"B{len([r for r in RUNS if r['group']=='B'])+1:02d}_cyl_block{k}", "B",
|
||||
f"Цилиндр Re=150, блокировка Ny/D={k}",
|
||||
"экстраполяция Cd и St к нулевой блокировке; ожидание — сходимость к 1.33 и 0.183",
|
||||
cyl_args(d, 150, nx, ny), nx * ny, st, gif=None, nx=nx)
|
||||
|
||||
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
# C. Модели стенки и субсеточность
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
def group_c(scale):
|
||||
for wall in ("hrr", "grad", "bouzidi", "staircase"):
|
||||
for re in (150, 2000):
|
||||
d, nx, ny = 32, 640, 384
|
||||
st = int(conv_steps(d, 400 if re > 20 else 60) * scale)
|
||||
add(f"C{len([r for r in RUNS if r['group']=='C'])+1:02d}_wall_{wall}_re{re}", "C",
|
||||
f"Модель стенки {wall}, Re={re}",
|
||||
"все субсеточные модели обязаны сойтись к одному пределу; staircase — база",
|
||||
cyl_args(d, re, nx, ny, wall=wall), nx * ny, st,
|
||||
gif=f"wall_{wall}_re{re}.gif" if re == 150 else None, nx=nx)
|
||||
|
||||
# Чувствительность к положению тела ВНУТРИ клетки — прямая проверка субсеточности, и
|
||||
# единственное измерение, которое вообще разделяет модели стенки. Двух положений для этого
|
||||
# мало: по двум точкам не отличить систематический разброс от совпадения, поэтому берём
|
||||
# четверти, и все четыре модели. Прогоны короткие (стационар при Re=20), вся серия — минуты.
|
||||
for wall in ("hrr", "grad", "bouzidi", "staircase"):
|
||||
for off in (0.0, 0.25, 0.5):
|
||||
d, nx, ny = 16, 320, 192
|
||||
st = int(conv_steps(d, 60) * scale)
|
||||
add(f"C{len([r for r in RUNS if r['group']=='C'])+1:02d}_sub_{wall}_{int(off*100):02d}",
|
||||
"C", f"Субклеточный сдвиг тела {off} клетки, стенка {wall}",
|
||||
"разброс Cd по сдвигам: чем субсеточнее модель, тем он меньше",
|
||||
cyl_args(d, 20, nx, ny, body_y=96 + off, wall=wall, pert_amp=0),
|
||||
nx * ny, st)
|
||||
|
||||
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
# D. Профили крыла
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
def foil_args(chord, alpha, re, nx, ny, naca="0012", **kw):
|
||||
a = {"--shape": "naca", "--naca": naca, "--size": chord,
|
||||
"--body-angle": alpha, "--nx": nx, "--ny": ny,
|
||||
"--body-x": nx // 4, "--re": re, "--refine": 1}
|
||||
a.update({"--" + k.replace("_", "-"): v for k, v in kw.items()})
|
||||
return flatten(a)
|
||||
|
||||
|
||||
def group_d(scale):
|
||||
# поляра: наклон Cl(alpha) и положение сваливания
|
||||
for al in (0, 4, 8, 12, 16):
|
||||
c, nx, ny = 128, 2048, 1024
|
||||
st = int(conv_steps(c, 400) * scale)
|
||||
add(f"D{len([r for r in RUNS if r['group']=='D'])+1:02d}_naca0012_a{al:02d}", "D",
|
||||
f"NACA 0012, α={al}°, Re=1000, хорда 128",
|
||||
"линейный участок Cl(α) и срыв; сравнение с низкорейнольдсовой литературой",
|
||||
foil_args(c, al, 1000, nx, ny), nx * ny, st,
|
||||
gif=f"naca0012_a{al:02d}.gif", nx=nx)
|
||||
|
||||
for re in (500, 2000, 10000):
|
||||
c, nx, ny = 128, 2048, 1024
|
||||
st = int(conv_steps(c, 400) * scale)
|
||||
add(f"D{len([r for r in RUNS if r['group']=='D'])+1:02d}_naca0012_re{re}", "D",
|
||||
f"NACA 0012, α=8°, Re={re}",
|
||||
"перестройка следа с ростом Re; проверка устойчивости на профиле",
|
||||
foil_args(c, 8, re, nx, ny), nx * ny, st,
|
||||
gif=f"naca0012_re{re}.gif", nx=nx)
|
||||
|
||||
for al in (0, 6, 12):
|
||||
c, nx, ny = 128, 2048, 1024
|
||||
st = int(conv_steps(c, 400) * scale)
|
||||
add(f"D{len([r for r in RUNS if r['group']=='D'])+1:02d}_naca4412_a{al:02d}", "D",
|
||||
f"NACA 4412 (с изгибом), α={al}°, Re=1000",
|
||||
"ненулевой Cl при α=0 — прямая проверка средней линии профиля",
|
||||
foil_args(c, al, 1000, nx, ny, naca="4412"), nx * ny, st,
|
||||
gif=f"naca4412_a{al:02d}.gif", nx=nx)
|
||||
|
||||
for c in (48, 96, 192):
|
||||
nx, ny = 16 * c, 8 * c
|
||||
st = int(conv_steps(c, 300) * scale)
|
||||
add(f"D{len([r for r in RUNS if r['group']=='D'])+1:02d}_naca_chord{c}", "D",
|
||||
f"NACA 0012, α=8°, сходимость по хорде {c}",
|
||||
"Cl и Cd сходятся по разрешению хорды",
|
||||
foil_args(c, 8, 1000, nx, ny), nx * ny, st, nx=nx)
|
||||
|
||||
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
# E. Сложная и множественная геометрия
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
def group_e(scale):
|
||||
# тандем: литература даёт переключение режимов около L/D 3.5–4
|
||||
for ld in (1.5, 3.0, 5.0):
|
||||
d, nx, ny = 48, 1440, 720
|
||||
gap = int(ld * d)
|
||||
bodies = f"cylinder:d={d},x=360,y=360; cylinder:d={d},x={360 + gap},y=360"
|
||||
st = int(conv_steps(d, 800) * scale)
|
||||
add(f"E{len([r for r in RUNS if r['group']=='E'])+1:02d}_tandem_ld{int(ld*10)}", "E",
|
||||
f"Тандем цилиндров, зазор L/D={ld}",
|
||||
"при малом зазоре у заднего тела ОТРИЦАТЕЛЬНОЕ сопротивление; переход около 3.5–4",
|
||||
["--bodies", bodies, "--nx", str(nx), "--ny", str(ny), "--re", "150", "--refine", "1"],
|
||||
nx * ny, st, gif=f"tandem_ld{int(ld*10)}.gif", nx=nx)
|
||||
|
||||
# решётка как пористая среда
|
||||
for step_d, tag in ((3, "sparse"), (2, "dense")):
|
||||
d, nx, ny = 24, 1440, 720
|
||||
bodies = "; ".join(
|
||||
f"cylinder:d={d},x={360 + i * step_d * d},y={360 + (j - 1.5) * step_d * d}"
|
||||
for i in range(4) for j in range(4))
|
||||
st = int(conv_steps(d, 600) * scale)
|
||||
add(f"E{len([r for r in RUNS if r['group']=='E'])+1:02d}_array_{tag}", "E",
|
||||
f"Решётка 4×4 цилиндров, шаг {step_d}D — пористая среда",
|
||||
"суммарное сопротивление и структура течения сквозь набор тел",
|
||||
["--bodies", bodies, "--nx", str(nx), "--ny", str(ny), "--re", "150", "--refine", "1"],
|
||||
nx * ny, st, gif=f"array_{tag}.gif", nx=nx)
|
||||
|
||||
# перфорированная пластина: набор коротких пластин со щелями
|
||||
for open_frac, tag in ((0.2, "open20"), (0.5, "open50")):
|
||||
nx, ny = 1200, 600
|
||||
seg, gap = 40, int(40 * open_frac / (1 - open_frac))
|
||||
bodies = "; ".join(
|
||||
f"plate:d={seg},t=0.12,a=90,x=300,y={60 + k * (seg + gap)}"
|
||||
for k in range(max(1, (ny - 120) // (seg + gap))))
|
||||
st = int(conv_steps(seg, 500) * scale)
|
||||
add(f"E{len([r for r in RUNS if r['group']=='E'])+1:02d}_perf_{tag}", "E",
|
||||
f"Перфорированная пластина, скважность {open_frac}",
|
||||
"струи в щелях и общее сопротивление; проверка множества тонких тел",
|
||||
["--bodies", bodies, "--nx", str(nx), "--ny", str(ny), "--re", "500", "--refine", "1"],
|
||||
nx * ny, st, gif=f"perf_{tag}.gif", nx=nx)
|
||||
|
||||
# многоэлементный профиль: основной + предкрылок
|
||||
nx, ny = 2048, 1024
|
||||
bodies = ("naca:naca=4412,d=200,a=6,x=560,y=512; "
|
||||
"naca:naca=2412,d=70,a=18,x=430,y=470")
|
||||
st = int(conv_steps(200, 300) * scale)
|
||||
add(f"E{len([r for r in RUNS if r['group']=='E'])+1:02d}_multielem", "E",
|
||||
"Многоэлементный профиль: основной 4412 плюс предкрылок",
|
||||
"щель между элементами и раздельные Cl/Cd по телам",
|
||||
["--bodies", bodies, "--nx", str(nx), "--ny", str(ny), "--re", "2000", "--refine", "1"],
|
||||
nx * ny, st, gif="multielem.gif", nx=nx)
|
||||
|
||||
# сверхтонкая пластина: предел субсеточности, тело тоньше клетки
|
||||
for t, tag in ((0.02, "t05"), (0.01, "t025")):
|
||||
nx, ny = 1200, 600
|
||||
st = int(conv_steps(50, 500) * scale)
|
||||
add(f"E{len([r for r in RUNS if r['group']=='E'])+1:02d}_thin_{tag}", "E",
|
||||
f"Сверхтонкая пластина, толщина {t*50:.2f} клетки",
|
||||
"предел субсеточной границы: тело тоньше клетки; ждём откаты на простой отскок",
|
||||
["--shape", "plate", "--size", "50", "--thickness", str(t), "--body-angle", "20",
|
||||
"--nx", str(nx), "--ny", str(ny), "--re", "500", "--refine", "1"],
|
||||
nx * ny, st, gif=f"thin_{tag}.gif", nx=nx)
|
||||
|
||||
# клин и зазубренная кромка — субклеточная деталь на многоугольнике
|
||||
nx, ny = 1200, 600
|
||||
st = int(conv_steps(60, 500) * scale)
|
||||
add(f"E{len([r for r in RUNS if r['group']=='E'])+1:02d}_wedge", "E",
|
||||
"Клин остриём против потока", "острая кромка: особая точка геометрии",
|
||||
["--shape", "polygon", "--poly=-30,-18;36,0;-30,18", "--size", "60",
|
||||
"--nx", str(nx), "--ny", str(ny), "--re", "500", "--refine", "1"],
|
||||
nx * ny, st, gif="wedge.gif", nx=nx)
|
||||
saw = ";".join(f"{-30 + 12 * k},{(-1) ** k * 6}" for k in range(6)) + ";30,-18;-30,-18"
|
||||
add(f"E{len([r for r in RUNS if r['group']=='E'])+1:02d}_serrated", "E",
|
||||
"Зазубренная задняя кромка", "субклеточные зубцы: как их видит граница",
|
||||
["--shape", "polygon", f"--poly={saw}", "--size", "60",
|
||||
"--nx", str(nx), "--ny", str(ny), "--re", "500", "--refine", "1"],
|
||||
nx * ny, st, gif="serrated.gif", nx=nx)
|
||||
|
||||
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
# F. Поле влияния и границы домена
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
def group_f(scale):
|
||||
d = 48
|
||||
for up in (2, 6, 16):
|
||||
nx, ny = int((up + 20) * d), 12 * d
|
||||
st = int(conv_steps(d, 600) * scale)
|
||||
add(f"F{len([r for r in RUNS if r['group']=='F'])+1:02d}_inlet{up}d", "F",
|
||||
f"Отступ до входа {up}D",
|
||||
"с какого отступа результат перестаёт зависеть от положения входа",
|
||||
cyl_args(d, 150, nx, ny, body_x=up * d), nx * ny, st)
|
||||
for down in (5, 15, 40):
|
||||
nx, ny = int((6 + down) * d), 12 * d
|
||||
st = int(conv_steps(d, 600) * scale)
|
||||
add(f"F{len([r for r in RUNS if r['group']=='F'])+1:02d}_outlet{down}d", "F",
|
||||
f"Отступ до выхода {down}D",
|
||||
"с какой длины следа выход перестаёт влиять на силы",
|
||||
cyl_args(d, 150, nx, ny, body_x=6 * d), nx * ny, st)
|
||||
for r in (1, 2, 3):
|
||||
d, nx, ny = 32, 640, 384
|
||||
st = int(conv_steps(d, 600) * scale)
|
||||
nodes = nx * ny
|
||||
if r > 1:
|
||||
ax, bx = int(nx // 5 - 1.5 * d), int(nx // 5 + 6.25 * d)
|
||||
ay, by = int(ny / 2 - 2 * d), int(ny / 2 + 2 * d)
|
||||
nodes += r * (r * (bx - ax) + 1) * (r * (by - ay) + 1)
|
||||
add(f"F{len([r_ for r_ in RUNS if r_['group']=='F'])+1:02d}_refine{r}", "F",
|
||||
f"Вложенный патч ×{r}",
|
||||
"вносит ли связка уровней систематику: Cd и St обязаны совпасть",
|
||||
cyl_args(d, 150, nx, ny, refine=r), nodes, st)
|
||||
|
||||
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
# G. Старт, акустика, время жизни
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
def group_g(scale):
|
||||
d, nx, ny = 48, 960, 576
|
||||
for re in (150, 2000, 20000):
|
||||
for init in ("uniform", "rest"):
|
||||
st = int(conv_steps(d, 300) * scale)
|
||||
add(f"G{len([r for r in RUNS if r['group']=='G'])+1:02d}_xt_re{re}_{init}", "G",
|
||||
f"x–t диаграмма, Re={re}, старт {init}",
|
||||
"откуда идёт стартовое возмущение и раскачивается ли продольная мода",
|
||||
cyl_args(d, re, nx, ny, init=init, sponge_len=0,
|
||||
xt=f"xt_re{re}_{init}.csv", xt_every=200),
|
||||
nx * ny, st, gif=f"xt_re{re}_{init}.gif", nx=nx, extra_out=("xt",))
|
||||
for sp in (0, 40):
|
||||
st = int(conv_steps(d, 400) * scale)
|
||||
add(f"G{len([r for r in RUNS if r['group']=='G'])+1:02d}_sponge{sp}", "G",
|
||||
f"Губка выхода {sp} столбцов при Re=5000",
|
||||
"губка — единственное, что реально ест продольную моду",
|
||||
cyl_args(d, 5000, nx, ny, sponge_len=sp), nx * ny, st)
|
||||
for sp in (0, 40):
|
||||
st = int(conv_steps(d, 400) * scale)
|
||||
add(f"G{len([r for r in RUNS if r['group']=='G'])+1:02d}_spongein{sp}", "G",
|
||||
f"Губка входа {sp} столбцов при Re=5000",
|
||||
"помогает ли гасить волну до отражения от входа",
|
||||
cyl_args(d, 5000, nx, ny, sponge_in=sp), nx * ny, st)
|
||||
for re in (10000, 50000, 200000):
|
||||
st = int(conv_steps(d, 200) * scale)
|
||||
add(f"G{len([r for r in RUNS if r['group']=='G'])+1:02d}_limit_re{re}", "G",
|
||||
f"Предел устойчивости, Re={re}",
|
||||
"где именно гибнет счёт; часть этих прогонов обязана развалиться",
|
||||
cyl_args(d, re, nx, ny), nx * ny, st)
|
||||
# тест на время жизни: один очень длинный прогон
|
||||
st = int(5_000_000 * scale)
|
||||
add(f"G{len([r for r in RUNS if r['group']=='G'])+1:02d}_lifetime", "G",
|
||||
"Тест на время жизни: 5·10⁶ шагов",
|
||||
"не уплывает ли ⟨ρ⟩ и не деградирует ли решение на очень длинной дистанции",
|
||||
cyl_args(d, 150, nx, ny), nx * ny, st, gif="lifetime.gif", nx=nx)
|
||||
|
||||
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
# H. Инварианты и точность
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
def group_h(scale):
|
||||
d, nx, ny = 32, 640, 384
|
||||
st = int(conv_steps(d, 400) * scale)
|
||||
add(f"H{len([r for r in RUNS if r['group']=='H'])+1:02d}_sym_cyl", "H",
|
||||
"Симметрия: цилиндр при α=0 без возмущения",
|
||||
"Cl и Cm обязаны остаться ≈0 — это внутреннее свойство схемы, не литература",
|
||||
cyl_args(d, 40, nx, ny, pert_amp=0), nx * ny, st)
|
||||
add(f"H{len([r for r in RUNS if r['group']=='H'])+1:02d}_sym_foil", "H",
|
||||
"Симметрия: NACA 0012 при α=0 без возмущения",
|
||||
"Cl ≈ 0 у симметричного профиля под нулевым углом",
|
||||
foil_args(64, 0, 500, 1024, 512, pert_amp=0), 1024 * 512,
|
||||
int(conv_steps(64, 300) * scale))
|
||||
for al, tag in ((8, "plus"), (-8, "minus")):
|
||||
add(f"H{len([r for r in RUNS if r['group']=='H'])+1:02d}_mirror_{tag}", "H",
|
||||
f"Зеркальность: NACA 0012 при α={al}°",
|
||||
"ряды обязаны зеркалиться: Cl меняет знак, Cd совпадает",
|
||||
foil_args(64, al, 500, 1024, 512), 1024 * 512,
|
||||
int(conv_steps(64, 300) * scale))
|
||||
for ul in (0.02, 0.05, 0.1):
|
||||
st2 = int(conv_steps(d, 400, ul) * scale)
|
||||
add(f"H{len([r for r in RUNS if r['group']=='H'])+1:02d}_mach{int(ul*100):02d}", "H",
|
||||
f"Зависимость от числа Маха, u_lat={ul}",
|
||||
"результат обязан сходиться при Ma→0; расхождение и есть сжимаемостная ошибка",
|
||||
cyl_args(d, 150, nx, ny, u_lat=ul), nx * ny, st2)
|
||||
# прямой замер расхождения f32 и f64 от длины прогона
|
||||
for bk in ("cpu", "gpu"):
|
||||
add(f"H{len([r for r in RUNS if r['group']=='H'])+1:02d}_prec_{bk}", "H",
|
||||
f"Точность: одинаковая постановка на {bk}",
|
||||
"расхождение f32 против f64 в зависимости от длины прогона",
|
||||
cyl_args(d, 150, nx, ny), nx * ny, int(conv_steps(d, 800) * scale), backend=bk)
|
||||
|
||||
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
# I. Сверхмелкие сетки по всем формам
|
||||
# ══════════════════════════════════════════════════════════════════════════════
|
||||
def group_i(scale):
|
||||
shapes = [
|
||||
("cylinder", ["--shape", "cylinder"]),
|
||||
("square", ["--shape", "square"]),
|
||||
("diamond", ["--shape", "diamond"]),
|
||||
("ellipse", ["--shape", "ellipse", "--thickness", "0.35"]),
|
||||
("naca", ["--shape", "naca", "--naca", "0012", "--body-angle", "8"]),
|
||||
("triangle", ["--shape", "triangle"]),
|
||||
("plate", ["--shape", "plate", "--thickness", "0.08", "--body-angle", "25"]),
|
||||
("wedge", ["--shape", "polygon", "--poly=-128,-72;154,0;-128,72"]),
|
||||
]
|
||||
d, nx, ny = 256, 4096, 2048
|
||||
st = int(conv_steps(d, 300) * scale)
|
||||
for tag, extra in shapes:
|
||||
add(f"I{len([r for r in RUNS if r['group']=='I'])+1:02d}_fine_{tag}", "I",
|
||||
f"Сверхмелкая сетка 4096×2048, {tag}, D=256",
|
||||
"держится ли субсеточная граница на предельном разрешении для этой формы",
|
||||
extra + ["--size", str(d), "--nx", str(nx), "--ny", str(ny),
|
||||
"--body-x", "820", "--re", "1000", "--refine", "1"],
|
||||
nx * ny, st, gif=f"fine_{tag}.gif", nx=nx)
|
||||
bodies = ("cylinder:d=200,x=800,y=1024; cylinder:d=200,x=1400,y=1024; "
|
||||
"naca:naca=4412,d=300,a=10,x=2200,y=1024")
|
||||
add(f"I{len([r for r in RUNS if r['group']=='I'])+1:02d}_fine_multi", "I",
|
||||
"Сверхмелкая сетка, композиция из трёх тел",
|
||||
"взаимодействие следов на предельном разрешении",
|
||||
["--bodies", bodies, "--nx", str(nx), "--ny", str(ny), "--re", "1000", "--refine", "1"],
|
||||
nx * ny, st, gif="fine_multi.gif", nx=nx)
|
||||
|
||||
|
||||
def main():
|
||||
ap = argparse.ArgumentParser()
|
||||
ap.add_argument("--scale", type=float, default=2.0,
|
||||
help="общий множитель длительности всех прогонов")
|
||||
ap.add_argument("--out", default=os.path.join(os.path.dirname(__file__), "scenarios.json"))
|
||||
args = ap.parse_args()
|
||||
|
||||
for g in (group_a, group_b, group_c, group_d, group_e, group_f, group_g, group_h, group_i):
|
||||
g(args.scale)
|
||||
|
||||
doc = {
|
||||
"meta": {
|
||||
"note": "Стоимость в обновлениях узлов — единственная переносимая между машинами "
|
||||
"мера. Часы драйвер получает, поделив её на фактические MLUPS, которые "
|
||||
"замеряет на месте (--calibrate).",
|
||||
"gif_policy": "вся длительность прогона, реальная скорость, 10 кадр/с; кадр "
|
||||
"прореживается до ширины не больше 1280 пикселей",
|
||||
"scale": args.scale,
|
||||
"assumed_gpu_mlups": 1200,
|
||||
"assumed_cpu_mlups": 22,
|
||||
},
|
||||
"runs": RUNS,
|
||||
}
|
||||
with open(args.out, "w", encoding="utf-8") as f:
|
||||
json.dump(doc, f, ensure_ascii=False, indent=1)
|
||||
|
||||
print(f"прогонов: {len(RUNS)} файл: {args.out}\n")
|
||||
print(f"{'гр':>3}{'прогонов':>10}{'обновлений узлов':>20}{'часов GPU':>12}{'часов CPU':>12}")
|
||||
tot_gpu = tot_cpu = 0.0
|
||||
for g in sorted({r["group"] for r in RUNS}):
|
||||
rs = [r for r in RUNS if r["group"] == g]
|
||||
cost = sum(r["cost"] for r in rs)
|
||||
hg = sum(r["cost"] / 1.2e9 / 3600 for r in rs if r["backend"] == "gpu")
|
||||
hc = sum(r["cost"] / 2.2e7 / 3600 for r in rs if r["backend"] == "cpu")
|
||||
tot_gpu += hg
|
||||
tot_cpu += hc
|
||||
print(f"{g:>3}{len(rs):>10}{cost:>20.3e}{hg:>12.1f}{hc:>12.1f}")
|
||||
print(f"{'—':>3}{len(RUNS):>10}{sum(r['cost'] for r in RUNS):>20.3e}"
|
||||
f"{tot_gpu:>12.1f}{tot_cpu:>12.1f}")
|
||||
print(f"\nвсего ≈ {tot_gpu + tot_cpu:.1f} ч при 1200 MLUPS на GPU и 22 на CPU")
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
main()
|
||||
@@ -1,79 +0,0 @@
|
||||
#!/usr/bin/env python3
|
||||
# Предполётная проверка кампании: каждый сценарий запускается на два шага.
|
||||
#
|
||||
# python preflight.py # проверить все
|
||||
# python preflight.py --group I # только группу
|
||||
#
|
||||
# Ловит опечатки в ключах, несовместимые сочетания и пределы железа ДО того, как кампания
|
||||
# уйдёт считать на девяносто часов. Не бесплатная привычка, а окупившаяся: именно так
|
||||
# выяснилось, что вся группа сверхмелких сеток падала на пределе GPU в 65535 рабочих групп
|
||||
# на измерение, а девять прогонов группы F передавали --body-x дважды.
|
||||
#
|
||||
# Физику проверка не трогает вовсе: два шага — это ровно «запустилось и не упало».
|
||||
|
||||
import argparse
|
||||
import json
|
||||
import os
|
||||
import subprocess
|
||||
import sys
|
||||
|
||||
HERE = os.path.dirname(os.path.abspath(__file__))
|
||||
BIN = os.environ.get("KBC2D_BIN") or os.path.join(
|
||||
HERE, "..", "target", "release", "kbc2d" + (".exe" if os.name == "nt" else ""))
|
||||
|
||||
# Всё, что пишет файлы, выкидываем: проверяем запуск, а не вывод. Ключ идёт со значением,
|
||||
# поэтому следующий за ним аргумент тоже пропускается.
|
||||
DROP = {"--gif", "--xt", "--case-csv", "--gif-field", "--gif-downsample", "--gif-scale",
|
||||
"--gif-every", "--gif-fps", "--series-every"}
|
||||
|
||||
|
||||
def main():
|
||||
ap = argparse.ArgumentParser()
|
||||
ap.add_argument("--group", help="только эта группа, например I")
|
||||
ap.add_argument("--scenarios", default=os.path.join(HERE, "scenarios.json"))
|
||||
args = ap.parse_args()
|
||||
|
||||
runs = json.load(open(args.scenarios, encoding="utf-8"))["runs"]
|
||||
if args.group:
|
||||
runs = [r for r in runs if r["group"].upper() == args.group.upper()]
|
||||
if not runs:
|
||||
print("под выборку не попал ни один прогон")
|
||||
return 0
|
||||
if not os.path.exists(BIN):
|
||||
print(f"не найден бинарь решателя: {BIN}")
|
||||
return 1
|
||||
|
||||
bad = []
|
||||
for i, r in enumerate(runs, 1):
|
||||
argv, skip = [], False
|
||||
for a in r["args"]:
|
||||
if skip:
|
||||
skip = False
|
||||
continue
|
||||
if a in DROP:
|
||||
skip = True
|
||||
continue
|
||||
argv.append(a)
|
||||
if "--steps" in argv:
|
||||
argv[argv.index("--steps") + 1] = "2"
|
||||
else:
|
||||
argv += ["--steps", "2"]
|
||||
|
||||
p = subprocess.run([BIN] + argv + ["--verbose", "quiet", "--report-every", "0"],
|
||||
capture_output=True, cwd=HERE)
|
||||
ok = p.returncode == 0
|
||||
print("." if ok else "X", end="", flush=True)
|
||||
if not ok:
|
||||
tail = (p.stdout + p.stderr).decode("utf-8", "replace").strip().splitlines()
|
||||
bad.append((r["id"], tail[-1] if tail else "(без вывода)"))
|
||||
if i % 40 == 0:
|
||||
print(f" {i}/{len(runs)}", flush=True)
|
||||
|
||||
print(f"\n\nпроверено сценариев: {len(runs)}, не запустились: {len(bad)}")
|
||||
for id_, msg in bad:
|
||||
print(f" {id_:<28} {msg}")
|
||||
return 1 if bad else 0
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
sys.exit(main())
|
||||
@@ -1,108 +0,0 @@
|
||||
# Драйвер валидационной кампании для Windows. Делает ровно то же, что run_campaign.sh, и
|
||||
# читает тот же scenarios.json — целевая машина Linux, а этот вариант нужен для локальной
|
||||
# отладки обвязки.
|
||||
#
|
||||
# .\run_campaign.ps1 -DryRun
|
||||
# .\run_campaign.ps1 -Calibrate
|
||||
# .\run_campaign.ps1 -Smoke
|
||||
# .\run_campaign.ps1 -Resume -Group A,B
|
||||
# .\run_campaign.ps1 -Only cyl_re150 -BudgetHours 6
|
||||
|
||||
[CmdletBinding()]
|
||||
param(
|
||||
[switch]$DryRun,
|
||||
[switch]$Resume,
|
||||
[switch]$Calibrate,
|
||||
[switch]$Smoke,
|
||||
[string[]]$Group,
|
||||
[string]$Only,
|
||||
[double]$BudgetHours = 0,
|
||||
[string]$Bin = "..\target\release\kbc2d.exe",
|
||||
[string]$Scenarios = "scenarios.json",
|
||||
[string]$Out = "out",
|
||||
[double]$GpuMlups = 1200,
|
||||
[double]$CpuMlups = 22
|
||||
)
|
||||
|
||||
Set-Location $PSScriptRoot
|
||||
if (-not (Test-Path $Bin)) { throw "не найден бинарь решателя: $Bin (собери cargo build --release)" }
|
||||
if (-not (Test-Path $Scenarios)) { throw "не найден список сценариев: $Scenarios" }
|
||||
|
||||
if ($Calibrate) {
|
||||
"Калибровка на этой машине (три коротких прогона)…"
|
||||
foreach ($c in @(@(240,120,'gpu'), @(960,480,'gpu'), @(480,240,'cpu'))) {
|
||||
$o = & $Bin --nx $c[0] --ny $c[1] --size 16 --refine 1 --steps 3000 `
|
||||
--report-every 3000 --verbose full --backend $c[2] 2>$null
|
||||
$m = ($o | Select-String 'MLUPS' | Select-Object -Last 1)
|
||||
" $($c[0])x$($c[1]) на $($c[2]): $(if ($m) { ($m -replace '.*?([0-9.]+) MLUPS.*','$1') } else { 'не измерено' })"
|
||||
}
|
||||
""
|
||||
"Подставь замеренное: .\run_campaign.ps1 -DryRun -GpuMlups <число> -CpuMlups <число>"
|
||||
return
|
||||
}
|
||||
|
||||
$runs = (Get-Content $Scenarios -Raw -Encoding UTF8 | ConvertFrom-Json).runs
|
||||
if ($Group) { $g = $Group | ForEach-Object { $_.ToUpper() }; $runs = $runs | Where-Object { $g -contains $_.group.ToUpper() } }
|
||||
if ($Only) { $runs = $runs | Where-Object { $_.id -like "*$Only*" } }
|
||||
if ($Smoke) { $runs = $runs | Sort-Object cost | Select-Object -First 3 }
|
||||
if (-not $runs) { "под выборку не попал ни один прогон"; return }
|
||||
|
||||
$totalCost = ($runs | Measure-Object -Property cost -Sum).Sum
|
||||
$estH = $totalCost / ($GpuMlups * 1e6) / 3600
|
||||
"прогонов: $($runs.Count) обновлений узлов: {0:e3} ≈ {1:F1} ч при $GpuMlups MLUPS" -f $totalCost, $estH
|
||||
"результаты: $Out\<id>\"
|
||||
""
|
||||
|
||||
if ($DryRun) {
|
||||
"{0,-28} {1,-3} {2,-4} {3,12} {4,8} {5}" -f 'id','гр','бэк','шагов','часов','название'
|
||||
foreach ($r in $runs) {
|
||||
$m = if ($r.backend -eq 'cpu') { $CpuMlups } else { $GpuMlups }
|
||||
"{0,-28} {1,-3} {2,-4} {3,12} {4,8:F2} {5}" -f $r.id, $r.group, $r.backend, $r.steps, ($r.cost/($m*1e6)/3600), $r.title
|
||||
}
|
||||
return
|
||||
}
|
||||
|
||||
New-Item -ItemType Directory -Force $Out | Out-Null
|
||||
$summary = Join-Path $Out 'summary.csv'
|
||||
$log = Join-Path $Out 'campaign.log'
|
||||
if (-not (Test-Path $summary)) { 'id,group,backend,status,seconds,steps,cost,title' | Out-File $summary -Encoding utf8 }
|
||||
|
||||
$started = Get-Date
|
||||
foreach ($r in $runs) {
|
||||
$dir = Join-Path $Out $r.id
|
||||
if ($Resume -and (Test-Path (Join-Path $dir 'summary.json'))) { "· $($r.id) — уже посчитан, пропускаю"; continue }
|
||||
if ($BudgetHours -gt 0 -and ((Get-Date) - $started).TotalHours -ge $BudgetHours) {
|
||||
"бюджет $BudgetHours ч исчерпан, останавливаюсь"; break
|
||||
}
|
||||
New-Item -ItemType Directory -Force $dir | Out-Null
|
||||
|
||||
# относительные имена файлов из сценария кладём внутрь папки прогона
|
||||
$argv = @()
|
||||
for ($i = 0; $i -lt $r.args.Count; $i++) {
|
||||
$v = $r.args[$i]
|
||||
if ($i -gt 0 -and @('--gif','--xt','--case-csv') -contains $r.args[$i-1]) { $v = Join-Path $dir $v }
|
||||
$argv += $v
|
||||
}
|
||||
$argv += @('--summary', (Join-Path $dir 'summary.json'),
|
||||
'--csv', (Join-Path $dir 'series.csv'),
|
||||
'--verbose','full','--report-every','5000')
|
||||
|
||||
"$Bin $($argv -join ' ')" | Out-File (Join-Path $dir 'cmd.txt') -Encoding utf8
|
||||
"▶ $($r.id) ($($r.group), $($r.backend)) $($r.title)"
|
||||
$t0 = Get-Date
|
||||
& $Bin @argv 2>&1 | Out-File (Join-Path $dir 'log.txt') -Encoding utf8
|
||||
$st = if ($LASTEXITCODE -eq 0) { 'ok' } else { 'fail' }
|
||||
$dt = [int]((Get-Date) - $t0).TotalSeconds
|
||||
if (Select-String -Path (Join-Path $dir 'log.txt') -Pattern 'РАЗВАЛИЛСЯ' -Quiet) { $st = 'blewup' }
|
||||
$rep = Get-Content (Join-Path $dir 'log.txt') -Encoding utf8
|
||||
$k = ($rep | Select-String 'ИТОГОВЫЙ ОТЧЁТ' | Select-Object -First 1).LineNumber
|
||||
if ($k) { $rep[($k-1)..($rep.Count-1)] | Out-File (Join-Path $dir 'report.txt') -Encoding utf8 }
|
||||
|
||||
'{0},{1},{2},{3},{4},{5},{6},"{7}"' -f $r.id,$r.group,$r.backend,$st,$dt,$r.steps,$r.cost,$r.title |
|
||||
Out-File $summary -Append -Encoding utf8
|
||||
"$(Get-Date -Format o) $($r.id) $st ${dt}s" | Out-File $log -Append -Encoding utf8
|
||||
" → $st за $dt с"
|
||||
}
|
||||
|
||||
""
|
||||
"готово. сводная таблица: $summary"
|
||||
@@ -1,150 +0,0 @@
|
||||
#!/usr/bin/env bash
|
||||
# Драйвер валидационной кампании. Читает scenarios.json, гоняет прогоны по одному, кладёт
|
||||
# логи, ряды, сводки и гифки каждого в собственную папку и ведёт общую таблицу.
|
||||
#
|
||||
# ./run_campaign.sh --dry-run смета: что и сколько будет считаться
|
||||
# ./run_campaign.sh --calibrate замерить фактические MLUPS этой машины
|
||||
# ./run_campaign.sh --smoke три коротких прогона: проверить обвязку
|
||||
# ./run_campaign.sh --resume считать, пропуская уже готовое
|
||||
# ./run_campaign.sh --group A,B только выбранные группы
|
||||
# ./run_campaign.sh --only cyl_re150 по подстроке идентификатора
|
||||
# ./run_campaign.sh --budget-hours 24 остановиться, когда время выйдет
|
||||
#
|
||||
# Прогон, который упал или развалился, помечается в сводке и НЕ останавливает кампанию:
|
||||
# группы E и G специально ищут предел устойчивости.
|
||||
|
||||
set -u
|
||||
cd "$(dirname "$0")"
|
||||
|
||||
BIN="${KBC2D_BIN:-../target/release/kbc2d}"
|
||||
SCEN="${KBC2D_SCENARIOS:-scenarios.json}"
|
||||
OUT="${KBC2D_OUT:-out}"
|
||||
# ВНИМАНИЕ: имя GROUPS занято самим bash (список групп пользователя), присваивание в неё
|
||||
# молча игнорируется, а "$GROUPS" под root разворачивается в 0 — и выборка съедает всю
|
||||
# кампанию, не сказав ни слова. Отсюда GRP_SEL.
|
||||
DRY=0; RESUME=0; CALIB=0; SMOKE=0; GRP_SEL=""; ONLY=""; BUDGET=""
|
||||
GPU_MLUPS="${KBC2D_GPU_MLUPS:-1200}"
|
||||
CPU_MLUPS="${KBC2D_CPU_MLUPS:-22}"
|
||||
|
||||
while [ $# -gt 0 ]; do
|
||||
case "$1" in
|
||||
--dry-run) DRY=1 ;;
|
||||
--resume) RESUME=1 ;;
|
||||
--calibrate) CALIB=1 ;;
|
||||
--smoke) SMOKE=1 ;;
|
||||
--group) GRP_SEL="$2"; shift ;;
|
||||
--only) ONLY="$2"; shift ;;
|
||||
--budget-hours) BUDGET="$2"; shift ;;
|
||||
-h|--help) sed -n '2,20p' "$0"; exit 0 ;;
|
||||
*) echo "неизвестный ключ: $1" >&2; exit 2 ;;
|
||||
esac
|
||||
shift
|
||||
done
|
||||
|
||||
command -v python3 >/dev/null 2>&1 && PY=python3 || PY=python
|
||||
[ -x "$BIN" ] || { echo "не найден бинарь решателя: $BIN (собери cargo build --release)" >&2; exit 1; }
|
||||
[ -f "$SCEN" ] || { echo "не найден список сценариев: $SCEN" >&2; exit 1; }
|
||||
|
||||
# ── калибровка: три коротких прогона на месте, чтобы оценки в часах были не гаданием ──
|
||||
if [ "$CALIB" = 1 ]; then
|
||||
echo "Калибровка на этой машине (три коротких прогона)…"
|
||||
for spec in "240 120 gpu" "960 480 gpu" "480 240 cpu"; do
|
||||
set -- $spec
|
||||
m=$("$BIN" --nx "$1" --ny "$2" --size 16 --refine 1 --steps 3000 --report-every 3000 \
|
||||
--verbose full --backend "$3" 2>/dev/null | grep -o '[0-9.]* MLUPS' | tail -1)
|
||||
echo " $1x$2 на $3: ${m:-не измерено}"
|
||||
done
|
||||
echo
|
||||
echo "Подставь замеренное в переменные окружения и запусти смету:"
|
||||
echo " KBC2D_GPU_MLUPS=<число> KBC2D_CPU_MLUPS=<число> ./run_campaign.sh --dry-run"
|
||||
exit 0
|
||||
fi
|
||||
|
||||
mkdir -p "$OUT"
|
||||
SUMMARY="$OUT/summary.csv"
|
||||
LOG="$OUT/campaign.log"
|
||||
[ -f "$SUMMARY" ] || echo "id,group,backend,status,seconds,steps,cost,title" > "$SUMMARY"
|
||||
|
||||
# Выборку и порядок считает python: разбирать JSON башем — верный способ ошибиться.
|
||||
PLAN=$("$PY" - "$SCEN" "$GRP_SEL" "$ONLY" "$SMOKE" <<'PYEOF'
|
||||
import json, sys
|
||||
scen, groups, only, smoke = sys.argv[1], sys.argv[2], sys.argv[3], sys.argv[4] == "1"
|
||||
runs = json.load(open(scen, encoding="utf-8"))["runs"]
|
||||
if groups:
|
||||
keep = {g.strip().upper() for g in groups.split(",")}
|
||||
runs = [r for r in runs if r["group"].upper() in keep]
|
||||
if only:
|
||||
runs = [r for r in runs if only in r["id"]]
|
||||
if smoke:
|
||||
# три самых дешёвых прогона: обвязку проверяем, а не физику
|
||||
runs = sorted(runs, key=lambda r: r["cost"])[:3]
|
||||
for r in runs:
|
||||
print("\t".join([r["id"], r["group"], r["backend"], str(r["cost"]),
|
||||
str(r["steps"]), r["title"], "\x1f".join(r["args"])]))
|
||||
PYEOF
|
||||
)
|
||||
|
||||
[ -n "$PLAN" ] || { echo "под выборку не попал ни один прогон"; exit 0; }
|
||||
|
||||
total_cost=0; n=0
|
||||
while IFS=$'\t' read -r id grp bk cost steps title args; do
|
||||
n=$((n + 1)); total_cost=$((total_cost + cost))
|
||||
done <<< "$PLAN"
|
||||
|
||||
est_h=$("$PY" -c "
|
||||
import sys
|
||||
print('%.1f' % (float(sys.argv[1]) / (float(sys.argv[2]) * 1e6) / 3600))" "$total_cost" "$GPU_MLUPS")
|
||||
|
||||
echo "прогонов: $n обновлений узлов: $total_cost ≈ $est_h ч при $GPU_MLUPS MLUPS"
|
||||
echo "результаты: $OUT/<id>/"
|
||||
echo
|
||||
|
||||
if [ "$DRY" = 1 ]; then
|
||||
printf '%-28s %-3s %-4s %12s %10s %s\n' id гр бэк шагов часов название
|
||||
while IFS=$'\t' read -r id grp bk cost steps title args; do
|
||||
m=$GPU_MLUPS; [ "$bk" = cpu ] && m=$CPU_MLUPS
|
||||
h=$("$PY" -c "import sys;print('%.2f'%(float(sys.argv[1])/(float(sys.argv[2])*1e6)/3600))" "$cost" "$m")
|
||||
printf '%-28s %-3s %-4s %12s %10s %s\n' "$id" "$grp" "$bk" "$steps" "$h" "$title"
|
||||
done <<< "$PLAN"
|
||||
exit 0
|
||||
fi
|
||||
|
||||
started=$(date +%s)
|
||||
while IFS=$'\t' read -r id grp bk cost steps title args; do
|
||||
dir="$OUT/$id"
|
||||
if [ "$RESUME" = 1 ] && [ -f "$dir/summary.json" ]; then
|
||||
echo "· $id — уже посчитан, пропускаю"
|
||||
continue
|
||||
fi
|
||||
if [ -n "$BUDGET" ]; then
|
||||
spent=$(( ($(date +%s) - started) ))
|
||||
limit=$("$PY" -c "import sys;print(int(float(sys.argv[1])*3600))" "$BUDGET")
|
||||
[ "$spent" -ge "$limit" ] && { echo "бюджет $BUDGET ч исчерпан, останавливаюсь"; break; }
|
||||
fi
|
||||
|
||||
mkdir -p "$dir"
|
||||
# аргументы разделены символом 0x1f: в них есть пробелы (список тел) и запятые
|
||||
IFS=$'\x1f' read -r -a ARGV <<< "$args"
|
||||
ARGV+=(--summary "$dir/summary.json" --csv "$dir/series.csv" --verbose full --report-every 5000)
|
||||
# относительные имена файлов из сценария кладём внутрь папки прогона
|
||||
for i in "${!ARGV[@]}"; do
|
||||
case "${ARGV[$((i-1))]:-}" in
|
||||
--gif|--xt|--case-csv) ARGV[$i]="$dir/${ARGV[$i]}" ;;
|
||||
esac
|
||||
done
|
||||
|
||||
printf '%s\n' "$BIN ${ARGV[*]}" > "$dir/cmd.txt"
|
||||
echo "▶ $id ($grp, $bk) $title"
|
||||
t0=$(date +%s)
|
||||
if "$BIN" "${ARGV[@]}" > "$dir/log.txt" 2>&1; then st=ok; else st=fail; fi
|
||||
t1=$(date +%s); dt=$((t1 - t0))
|
||||
grep -q "РАЗВАЛИЛСЯ" "$dir/log.txt" 2>/dev/null && st=blewup
|
||||
sed -n '/ИТОГОВЫЙ ОТЧЁТ/,$p' "$dir/log.txt" > "$dir/report.txt" 2>/dev/null
|
||||
|
||||
printf '%s,%s,%s,%s,%s,%s,%s,"%s"\n' "$id" "$grp" "$bk" "$st" "$dt" "$steps" "$cost" "$title" >> "$SUMMARY"
|
||||
echo "$(date -Is) $id $st ${dt}s" >> "$LOG"
|
||||
echo " → $st за ${dt} с"
|
||||
done <<< "$PLAN"
|
||||
|
||||
echo
|
||||
echo "готово. сводная таблица: $SUMMARY"
|
||||
File diff suppressed because it is too large
Load Diff
@@ -1,27 +0,0 @@
|
||||
# Запуск кампании на сервере с NVIDIA одной командой:
|
||||
# docker compose run --rm kbc2d --calibrate # сначала убедиться, что GPU виден
|
||||
# docker compose up -d # кампания в фоне
|
||||
# docker compose logs -f # смотреть ход
|
||||
#
|
||||
# Результаты складываются в ./out на хосте.
|
||||
|
||||
services:
|
||||
kbc2d:
|
||||
build: .
|
||||
image: kbc2d
|
||||
container_name: kbc2d
|
||||
volumes:
|
||||
- ./out:/work/bench/out
|
||||
environment:
|
||||
# Без graphics NVIDIA Container Toolkit не подложит Vulkan-ICD и wgpu не увидит карту.
|
||||
NVIDIA_DRIVER_CAPABILITIES: compute,utility,graphics
|
||||
deploy:
|
||||
resources:
|
||||
reservations:
|
||||
devices:
|
||||
- driver: nvidia
|
||||
count: all
|
||||
capabilities: [gpu]
|
||||
# Кампания идёт десятки часов и переживает перезапуск: --resume пропускает готовое.
|
||||
command: ["--resume"]
|
||||
restart: "no"
|
||||
File diff suppressed because it is too large
Load Diff
@@ -1,632 +0,0 @@
|
||||
//! БЛОК СОЗДАНИЯ ГИФОК: перевод поля в индексированный кадр, палитры, служебная надпись и —
|
||||
//! главное — СИНХРОНИЗАЦИЯ АНИМАЦИИ С ФИЗИЧЕСКИМ ВРЕМЕНЕМ ПОТОКА.
|
||||
//!
|
||||
//! Требование: гифка идёт с той же скоростью, что и настоящий поток, независимо от того, с
|
||||
//! какой скоростью считает машина. Реальная производительность (шагов в секунду wall-clock)
|
||||
//! в тайминг не входит вообще — берётся только физическая длительность шага δt из `Units`:
|
||||
//!
|
||||
//! δt = u_lat·δx/u_phys с/шаг
|
||||
//! кадр каждые S шагов ⇒ между кадрами S·δt секунд физического времени
|
||||
//! задержка кадра = S·δt / playback с (playback = 1 ⇒ реальное время)
|
||||
//!
|
||||
//! Ограничение формата: задержка в GIF хранится в СОТЫХ долях секунды (u16), т.е. сетка 10 мс,
|
||||
//! и почти все декодеры считают задержку 0–1 сс как «как можно быстрее». Поэтому реальное время
|
||||
//! достижимо не при любом S: при δt = 1.7e-4 с (30 м/с, ячейка 0.1 м, u_lat = 0.05) кадр каждые
|
||||
//! 10 шагов — это 600 кадр/с, чего GIF не умеет. Отсюда два режима:
|
||||
//!
|
||||
//! * `--gif-every N` — шаг задан жёстко; если реальное время не представимо, задержка
|
||||
//! зажимается и в отчёт печатается фактический коэффициент замедления;
|
||||
//! * `--gif-every auto` — шаг подбирается так, чтобы реальное время получилось точно при
|
||||
//! целевой частоте `--gif-fps` (по умолчанию 30). Это режим по
|
||||
//! умолчанию: синхронность важнее, чем конкретное число шагов на кадр.
|
||||
|
||||
use std::fs::File;
|
||||
use std::io::BufWriter;
|
||||
|
||||
use crate::math::{Units, R};
|
||||
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
// Тайминг
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
|
||||
/// Как кадры анимации ложатся на физическое время.
|
||||
#[derive(Clone, Copy, Debug)]
|
||||
pub struct GifPlan {
|
||||
/// Кадр каждые столько шагов симуляции.
|
||||
pub stride: u64,
|
||||
/// ТОЧНАЯ задержка кадра в сотых долях секунды — вообще говоря дробная.
|
||||
pub exact_delay_cs: R,
|
||||
/// Физическое время между соседними кадрами, с.
|
||||
pub frame_dt: R,
|
||||
/// Запрошенная скорость воспроизведения (1 = реальное время, 0.1 = 10× замедление).
|
||||
pub requested_playback: R,
|
||||
/// Фактическая, в среднем по анимации.
|
||||
pub actual_playback: R,
|
||||
/// Средняя частота кадров получившейся гифки.
|
||||
pub fps: R,
|
||||
/// Шаг подбирался автоматически.
|
||||
pub auto: bool,
|
||||
/// Точная задержка нецелая ⇒ задержки кадров чередуются (см. `DelayDither`).
|
||||
pub dithered: bool,
|
||||
/// Точная задержка меньше одной сотой ⇒ формат физически не тянет, время сжато.
|
||||
pub clamped: bool,
|
||||
}
|
||||
|
||||
impl GifPlan {
|
||||
/// `stride`: Some(N) — жёстко N шагов на кадр; None — подобрать под `target_fps`.
|
||||
pub fn new(units: &Units, stride: Option<u64>, target_fps: R, playback: R) -> GifPlan {
|
||||
let dt = units.dt;
|
||||
let playback = if playback > 0.0 { playback } else { 1.0 };
|
||||
let auto = stride.is_none();
|
||||
|
||||
let stride = match stride {
|
||||
Some(s) => s.max(1),
|
||||
None => {
|
||||
// хотим кадр раз в 1/target_fps секунды экранного времени, что при скорости
|
||||
// playback отвечает playback/target_fps секундам физического времени
|
||||
((playback / target_fps / dt).round() as i64).max(1) as u64
|
||||
}
|
||||
};
|
||||
|
||||
let frame_dt = stride as R * dt;
|
||||
let exact_delay_cs = 100.0 * frame_dt / playback;
|
||||
// 0 сотых большинство декодеров трактует как «как можно быстрее», поэтому пол — единица
|
||||
let clamped = exact_delay_cs < 1.0;
|
||||
let mean_delay_cs = if clamped { 1.0 } else { exact_delay_cs };
|
||||
|
||||
GifPlan {
|
||||
stride,
|
||||
exact_delay_cs,
|
||||
frame_dt,
|
||||
requested_playback: playback,
|
||||
actual_playback: frame_dt / (mean_delay_cs / 100.0),
|
||||
fps: 100.0 / mean_delay_cs,
|
||||
auto,
|
||||
dithered: !clamped && (exact_delay_cs - exact_delay_cs.round()).abs() > 1e-9,
|
||||
clamped,
|
||||
}
|
||||
}
|
||||
|
||||
/// Отношение фактической скорости воспроизведения к запрошенной (1.0 = точно).
|
||||
pub fn sync_error(&self) -> R {
|
||||
self.actual_playback / self.requested_playback
|
||||
}
|
||||
|
||||
/// Синхронность выдержана: средняя задержка совпадает с точной.
|
||||
pub fn is_synced(&self) -> bool {
|
||||
(self.sync_error() - 1.0).abs() < 0.02
|
||||
}
|
||||
|
||||
/// Совет, как починить рассинхрон, если формат его не тянет.
|
||||
pub fn advice(&self, units: &Units) -> Option<String> {
|
||||
// разойтись с запрошенной скоростью можно только одним способом: точная задержка
|
||||
// не влезла в минимальную единицу формата и была зажата
|
||||
if !self.clamped {
|
||||
return None;
|
||||
}
|
||||
let good = ((0.01 * self.requested_playback / units.dt).ceil() as i64).max(1);
|
||||
Some(format!(
|
||||
"кадр требует задержки {:.4} сотых — меньше минимальной единицы формата, поэтому \
|
||||
анимация идёт в {:.2}× быстрее запрошенного. Починить: --gif-every auto (подберёт \
|
||||
шаг сам), либо --gif-every {} (минимальный шаг с представимой задержкой), либо \
|
||||
--gif-speed {:.4} (замедление вместо реального времени)",
|
||||
self.exact_delay_cs,
|
||||
self.sync_error(),
|
||||
good,
|
||||
self.requested_playback / self.sync_error()
|
||||
))
|
||||
}
|
||||
}
|
||||
|
||||
/// Раскладка дробной задержки по целым сотым без накопления ошибки.
|
||||
///
|
||||
/// Формат GIF хранит задержку кадра целым числом сотых долей секунды. Точная физическая
|
||||
/// задержка почти никогда не целая: 30 кадр/с — это 3⅓ сотых. Округление каждого кадра
|
||||
/// по отдельности даёт систематический уход (3 вместо 3.333 — анимация на 11% быстрее
|
||||
/// реальности, и расхождение копится линейно). Поэтому задержки выдаются так, чтобы
|
||||
/// НАКОПЛЕННОЕ время кадров всегда совпадало с накопленным физическим с точностью до одной
|
||||
/// сотой: для 3.333 получается 3, 3, 4, 3, 3, 4, … Ошибка ограничена ±5 мс и не растёт.
|
||||
#[derive(Clone, Copy, Debug, Default)]
|
||||
pub struct DelayDither {
|
||||
exact_cum: R,
|
||||
emitted_cum: i64,
|
||||
}
|
||||
|
||||
impl DelayDither {
|
||||
/// Задержка очередного кадра в сотых долях секунды.
|
||||
pub fn next(&mut self, exact_delay_cs: R) -> u16 {
|
||||
self.exact_cum += exact_delay_cs;
|
||||
let want = self.exact_cum.round() as i64;
|
||||
let d = (want - self.emitted_cum).max(1);
|
||||
self.emitted_cum += d;
|
||||
d.clamp(1, u16::MAX as i64) as u16
|
||||
}
|
||||
}
|
||||
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
// Палитры
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
|
||||
/// Индексы служебных цветов в палитре (после градиента).
|
||||
const NCOLORS: usize = 250;
|
||||
const IDX_BODY: u8 = 250;
|
||||
const IDX_PATCH: u8 = 251;
|
||||
const IDX_TEXT: u8 = 252;
|
||||
const IDX_TEXT_BG: u8 = 253;
|
||||
|
||||
#[derive(Clone, Copy, Debug, PartialEq, Eq)]
|
||||
pub enum ColorMap {
|
||||
/// Последовательная, воспринимаемо-равномерная — для |u| и ρ.
|
||||
Viridis,
|
||||
/// Яркая последовательная с широким динамическим диапазоном.
|
||||
Turbo,
|
||||
/// Расходящаяся сине-бело-красная — для знакопеременных полей (завихренность, γ−2).
|
||||
CoolWarm,
|
||||
}
|
||||
|
||||
impl ColorMap {
|
||||
pub fn from_str(s: &str) -> Option<Self> {
|
||||
Some(match s.to_ascii_lowercase().as_str() {
|
||||
"viridis" => ColorMap::Viridis,
|
||||
"turbo" => ColorMap::Turbo,
|
||||
"coolwarm" | "rdbu" => ColorMap::CoolWarm,
|
||||
_ => return None,
|
||||
})
|
||||
}
|
||||
pub const ALL: [&'static str; 3] = ["viridis", "turbo", "coolwarm"];
|
||||
|
||||
/// Опорные точки градиента; между ними линейная интерполяция в sRGB.
|
||||
fn anchors(&self) -> &'static [[u8; 3]] {
|
||||
match self {
|
||||
ColorMap::Viridis => &[
|
||||
[68, 1, 84],
|
||||
[72, 40, 120],
|
||||
[62, 74, 137],
|
||||
[49, 104, 142],
|
||||
[38, 130, 142],
|
||||
[31, 158, 137],
|
||||
[53, 183, 121],
|
||||
[109, 205, 89],
|
||||
[180, 222, 44],
|
||||
[253, 231, 37],
|
||||
],
|
||||
ColorMap::Turbo => &[
|
||||
[48, 18, 59],
|
||||
[70, 107, 227],
|
||||
[40, 187, 236],
|
||||
[49, 242, 153],
|
||||
[138, 252, 61],
|
||||
[210, 226, 27],
|
||||
[254, 165, 45],
|
||||
[239, 89, 17],
|
||||
[180, 27, 1],
|
||||
[122, 4, 3],
|
||||
],
|
||||
// Расходящаяся палитра Морланда. Опорных точек НЕЧЁТНОЕ число: только тогда
|
||||
// середина диапазона (для знакопеременного поля — ноль) попадает ровно на
|
||||
// нейтральный серый, а не между двумя соседними точками.
|
||||
ColorMap::CoolWarm => &[
|
||||
[59, 76, 192],
|
||||
[98, 130, 234],
|
||||
[141, 176, 254],
|
||||
[184, 208, 249],
|
||||
[221, 221, 221],
|
||||
[245, 196, 173],
|
||||
[244, 154, 123],
|
||||
[222, 96, 77],
|
||||
[180, 4, 38],
|
||||
],
|
||||
}
|
||||
}
|
||||
|
||||
/// Палитра GIF: NCOLORS ступеней градиента плюс служебные цвета.
|
||||
fn palette(&self) -> Vec<u8> {
|
||||
let a = self.anchors();
|
||||
let mut p = vec![0u8; 256 * 3];
|
||||
for k in 0..NCOLORS {
|
||||
let t = k as R / (NCOLORS - 1) as R * (a.len() - 1) as R;
|
||||
let i = (t.floor() as usize).min(a.len() - 2);
|
||||
let fr = t - i as R;
|
||||
for c in 0..3 {
|
||||
p[k * 3 + c] = (a[i][c] as R * (1.0 - fr) + a[i + 1][c] as R * fr).round() as u8;
|
||||
}
|
||||
}
|
||||
let set = |p: &mut Vec<u8>, idx: u8, rgb: [u8; 3]| {
|
||||
let o = idx as usize * 3;
|
||||
p[o] = rgb[0];
|
||||
p[o + 1] = rgb[1];
|
||||
p[o + 2] = rgb[2];
|
||||
};
|
||||
set(&mut p, IDX_BODY, [25, 25, 28]);
|
||||
set(&mut p, IDX_PATCH, [255, 255, 255]);
|
||||
set(&mut p, IDX_TEXT, [255, 255, 255]);
|
||||
set(&mut p, IDX_TEXT_BG, [0, 0, 0]);
|
||||
p
|
||||
}
|
||||
}
|
||||
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
// Микрошрифт 5×7 для служебной надписи
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
|
||||
/// Битовая маска глифа: 7 строк по 5 бит (старший бит — левый пиксель).
|
||||
fn glyph(c: char) -> [u8; 7] {
|
||||
match c.to_ascii_uppercase() {
|
||||
'0' => [0x0E, 0x11, 0x13, 0x15, 0x19, 0x11, 0x0E],
|
||||
'1' => [0x04, 0x0C, 0x04, 0x04, 0x04, 0x04, 0x0E],
|
||||
'2' => [0x0E, 0x11, 0x01, 0x02, 0x04, 0x08, 0x1F],
|
||||
'3' => [0x1F, 0x02, 0x04, 0x02, 0x01, 0x11, 0x0E],
|
||||
'4' => [0x02, 0x06, 0x0A, 0x12, 0x1F, 0x02, 0x02],
|
||||
'5' => [0x1F, 0x10, 0x1E, 0x01, 0x01, 0x11, 0x0E],
|
||||
'6' => [0x06, 0x08, 0x10, 0x1E, 0x11, 0x11, 0x0E],
|
||||
'7' => [0x1F, 0x01, 0x02, 0x04, 0x08, 0x08, 0x08],
|
||||
'8' => [0x0E, 0x11, 0x11, 0x0E, 0x11, 0x11, 0x0E],
|
||||
'9' => [0x0E, 0x11, 0x11, 0x0F, 0x01, 0x02, 0x0C],
|
||||
'A' => [0x0E, 0x11, 0x11, 0x1F, 0x11, 0x11, 0x11],
|
||||
'B' => [0x1E, 0x11, 0x11, 0x1E, 0x11, 0x11, 0x1E],
|
||||
'C' => [0x0E, 0x11, 0x10, 0x10, 0x10, 0x11, 0x0E],
|
||||
'D' => [0x1C, 0x12, 0x11, 0x11, 0x11, 0x12, 0x1C],
|
||||
'E' => [0x1F, 0x10, 0x10, 0x1E, 0x10, 0x10, 0x1F],
|
||||
'F' => [0x1F, 0x10, 0x10, 0x1E, 0x10, 0x10, 0x10],
|
||||
'G' => [0x0E, 0x11, 0x10, 0x17, 0x11, 0x11, 0x0F],
|
||||
'H' => [0x11, 0x11, 0x11, 0x1F, 0x11, 0x11, 0x11],
|
||||
'I' => [0x0E, 0x04, 0x04, 0x04, 0x04, 0x04, 0x0E],
|
||||
'K' => [0x11, 0x12, 0x14, 0x18, 0x14, 0x12, 0x11],
|
||||
'L' => [0x10, 0x10, 0x10, 0x10, 0x10, 0x10, 0x1F],
|
||||
'M' => [0x11, 0x1B, 0x15, 0x15, 0x11, 0x11, 0x11],
|
||||
'N' => [0x11, 0x19, 0x15, 0x13, 0x11, 0x11, 0x11],
|
||||
'O' => [0x0E, 0x11, 0x11, 0x11, 0x11, 0x11, 0x0E],
|
||||
'P' => [0x1E, 0x11, 0x11, 0x1E, 0x10, 0x10, 0x10],
|
||||
'R' => [0x1E, 0x11, 0x11, 0x1E, 0x14, 0x12, 0x11],
|
||||
'S' => [0x0F, 0x10, 0x10, 0x0E, 0x01, 0x01, 0x1E],
|
||||
'T' => [0x1F, 0x04, 0x04, 0x04, 0x04, 0x04, 0x04],
|
||||
'U' => [0x11, 0x11, 0x11, 0x11, 0x11, 0x11, 0x0E],
|
||||
'V' => [0x11, 0x11, 0x11, 0x11, 0x11, 0x0A, 0x04],
|
||||
'W' => [0x11, 0x11, 0x11, 0x15, 0x15, 0x1B, 0x11],
|
||||
'X' => [0x11, 0x11, 0x0A, 0x04, 0x0A, 0x11, 0x11],
|
||||
'Y' => [0x11, 0x11, 0x0A, 0x04, 0x04, 0x04, 0x04],
|
||||
'Z' => [0x1F, 0x01, 0x02, 0x04, 0x08, 0x10, 0x1F],
|
||||
'.' => [0x00, 0x00, 0x00, 0x00, 0x00, 0x0C, 0x0C],
|
||||
',' => [0x00, 0x00, 0x00, 0x00, 0x0C, 0x04, 0x08],
|
||||
':' => [0x00, 0x0C, 0x0C, 0x00, 0x0C, 0x0C, 0x00],
|
||||
'-' => [0x00, 0x00, 0x00, 0x1F, 0x00, 0x00, 0x00],
|
||||
'+' => [0x00, 0x04, 0x04, 0x1F, 0x04, 0x04, 0x00],
|
||||
'=' => [0x00, 0x00, 0x1F, 0x00, 0x1F, 0x00, 0x00],
|
||||
'/' => [0x01, 0x02, 0x02, 0x04, 0x08, 0x08, 0x10],
|
||||
'%' => [0x19, 0x1A, 0x02, 0x04, 0x08, 0x0B, 0x13],
|
||||
'(' => [0x02, 0x04, 0x08, 0x08, 0x08, 0x04, 0x02],
|
||||
')' => [0x08, 0x04, 0x02, 0x02, 0x02, 0x04, 0x08],
|
||||
'Ч' => [0x11, 0x11, 0x11, 0x0F, 0x01, 0x01, 0x01],
|
||||
_ => [0; 7], // пробел и всё незнакомое
|
||||
}
|
||||
}
|
||||
|
||||
/// Нарисовать строку в индексированный буфер (координаты — левый верхний угол).
|
||||
fn draw_text(buf: &mut [u8], w: usize, h: usize, x0: usize, y0: usize, s: &str, scale: usize) {
|
||||
let gw = 6 * scale; // 5 пикселей глифа + 1 пробел
|
||||
// подложка, чтобы текст читался на любом фоне
|
||||
let tw = s.chars().count() * gw;
|
||||
for yy in y0.saturating_sub(scale)..(y0 + 7 * scale + scale).min(h) {
|
||||
for xx in x0.saturating_sub(scale)..(x0 + tw + scale).min(w) {
|
||||
buf[yy * w + xx] = IDX_TEXT_BG;
|
||||
}
|
||||
}
|
||||
for (k, ch) in s.chars().enumerate() {
|
||||
let g = glyph(ch);
|
||||
for (row, bits) in g.iter().enumerate() {
|
||||
for col in 0..5 {
|
||||
if bits & (1 << (4 - col)) == 0 {
|
||||
continue;
|
||||
}
|
||||
for sy in 0..scale {
|
||||
for sx in 0..scale {
|
||||
let px = x0 + k * gw + col * scale + sx;
|
||||
let py = y0 + row * scale + sy;
|
||||
if px < w && py < h {
|
||||
buf[py * w + px] = IDX_TEXT;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
// Диапазон нормировки
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
|
||||
/// Диапазон значений, отображаемый на палитру. Фиксируется на весь прогон: плавающая
|
||||
/// автонормировка делает анимацию нечитаемой — цвет перестаёт что-либо значить.
|
||||
#[derive(Clone, Copy, Debug)]
|
||||
pub struct Range {
|
||||
pub lo: R,
|
||||
pub hi: R,
|
||||
}
|
||||
|
||||
impl Range {
|
||||
#[inline]
|
||||
fn index(&self, v: R) -> u8 {
|
||||
if !v.is_finite() {
|
||||
return 0;
|
||||
}
|
||||
let t = ((v - self.lo) / (self.hi - self.lo)).clamp(0.0, 1.0);
|
||||
(t * (NCOLORS - 1) as R).round() as u8
|
||||
}
|
||||
}
|
||||
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
// Запись гифки
|
||||
// ─────────────────────────────────────────────────────────────────────────────
|
||||
|
||||
/// Что подписывать на кадре.
|
||||
pub struct Hud {
|
||||
pub u_phys: R,
|
||||
pub field: &'static str,
|
||||
pub show: bool,
|
||||
}
|
||||
|
||||
pub struct GifWriter {
|
||||
enc: gif::Encoder<BufWriter<File>>,
|
||||
pub plan: GifPlan,
|
||||
/// Размеры ИСХОДНОЙ сетки.
|
||||
src_nx: usize,
|
||||
src_ny: usize,
|
||||
/// Во сколько раз клетки усредняются в пиксель перед отрисовкой.
|
||||
down: usize,
|
||||
/// Размеры картинки после прореживания.
|
||||
nx: usize,
|
||||
ny: usize,
|
||||
scale: usize,
|
||||
w: usize,
|
||||
h: usize,
|
||||
range: Range,
|
||||
hud: Hud,
|
||||
patch: Option<(usize, usize, usize, usize)>,
|
||||
dither: DelayDither,
|
||||
pub frames: u32,
|
||||
}
|
||||
|
||||
impl GifWriter {
|
||||
#[allow(clippy::too_many_arguments)]
|
||||
pub fn create(
|
||||
path: &str,
|
||||
nx: usize,
|
||||
ny: usize,
|
||||
scale: usize,
|
||||
down: usize,
|
||||
cmap: ColorMap,
|
||||
range: Range,
|
||||
plan: GifPlan,
|
||||
hud: Hud,
|
||||
patch: Option<(usize, usize, usize, usize)>,
|
||||
) -> std::io::Result<GifWriter> {
|
||||
let scale = scale.max(1);
|
||||
let down = down.max(1);
|
||||
let (src_nx, src_ny) = (nx, ny);
|
||||
// прореживание с округлением вверх: последний блок может быть неполным
|
||||
let nx = nx.div_ceil(down);
|
||||
let ny = ny.div_ceil(down);
|
||||
let w = nx * scale;
|
||||
let h = ny * scale;
|
||||
let file = BufWriter::new(File::create(path)?);
|
||||
let mut enc = gif::Encoder::new(file, w as u16, h as u16, &cmap.palette())
|
||||
.map_err(std::io::Error::other)?;
|
||||
enc.set_repeat(gif::Repeat::Infinite).map_err(std::io::Error::other)?;
|
||||
Ok(GifWriter {
|
||||
enc,
|
||||
plan,
|
||||
src_nx,
|
||||
src_ny,
|
||||
down,
|
||||
nx,
|
||||
ny,
|
||||
scale,
|
||||
w,
|
||||
h,
|
||||
range,
|
||||
hud,
|
||||
patch,
|
||||
dither: DelayDither::default(),
|
||||
frames: 0,
|
||||
})
|
||||
}
|
||||
|
||||
/// Добавить кадр. `field` — значения по узлам L0, `solid` — маска тела,
|
||||
/// `t_phys` — физическое время этого кадра (с), `cd` — текущий коэффициент сопротивления.
|
||||
pub fn push(
|
||||
&mut self,
|
||||
field: &[R],
|
||||
solid: &[bool],
|
||||
t_phys: R,
|
||||
cd: Option<R>,
|
||||
) -> std::io::Result<()> {
|
||||
// Прореживание: блок down×down усредняется в один пиксель, а телом пиксель считается,
|
||||
// если тело занимает хотя бы половину блока. Без этого кадр с сетки 4096×2048 весит
|
||||
// столько, что гифка становится непригодной.
|
||||
let (field, solid) = if self.down == 1 {
|
||||
(field.to_vec(), solid.to_vec())
|
||||
} else {
|
||||
let d = self.down;
|
||||
let mut fv = vec![0.0; self.nx * self.ny];
|
||||
let mut sv = vec![false; self.nx * self.ny];
|
||||
for gy in 0..self.ny {
|
||||
for gx in 0..self.nx {
|
||||
let (mut acc, mut cnt, mut sol) = (0.0, 0usize, 0usize);
|
||||
for yy in gy * d..((gy + 1) * d).min(self.src_ny) {
|
||||
for xx in gx * d..((gx + 1) * d).min(self.src_nx) {
|
||||
let k = yy * self.src_nx + xx;
|
||||
acc += field[k];
|
||||
sol += solid[k] as usize;
|
||||
cnt += 1;
|
||||
}
|
||||
}
|
||||
let g = gy * self.nx + gx;
|
||||
fv[g] = if cnt > 0 { acc / cnt as R } else { 0.0 };
|
||||
sv[g] = cnt > 0 && 2 * sol >= cnt;
|
||||
}
|
||||
}
|
||||
(fv, sv)
|
||||
};
|
||||
let (field, solid) = (&field[..], &solid[..]);
|
||||
let mut buf = vec![0u8; self.w * self.h];
|
||||
// строка 0 изображения — это ВЕРХ, а y = 0 решётки — низ канала: переворачиваем
|
||||
for y in 0..self.ny {
|
||||
let src = (self.ny - 1 - y) * self.nx;
|
||||
for x in 0..self.nx {
|
||||
let idx = if solid[src + x] { IDX_BODY } else { self.range.index(field[src + x]) };
|
||||
for sy in 0..self.scale {
|
||||
let row = (y * self.scale + sy) * self.w + x * self.scale;
|
||||
for sx in 0..self.scale {
|
||||
buf[row + sx] = idx;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
if let Some((ax, bx, ay, by)) = self.patch {
|
||||
let d = self.down;
|
||||
self.draw_patch_outline(&mut buf, ax / d, bx / d, ay / d, by / d);
|
||||
}
|
||||
if self.hud.show {
|
||||
let line = match cd {
|
||||
Some(c) => format!(
|
||||
"T={:.4}S U={:.1}M/S CD={:.2}",
|
||||
t_phys, self.hud.u_phys, c
|
||||
),
|
||||
None => format!("T={:.4}S U={:.1}M/S", t_phys, self.hud.u_phys),
|
||||
};
|
||||
let sc = (self.scale.min(3)).max(1);
|
||||
draw_text(&mut buf, self.w, self.h, 2 * sc, 2 * sc, &line, sc);
|
||||
draw_text(&mut buf, self.w, self.h, 2 * sc, self.h - 9 * sc, self.hud.field, sc);
|
||||
}
|
||||
|
||||
let mut frame = gif::Frame::from_indexed_pixels(self.w as u16, self.h as u16, buf, None);
|
||||
// задержка берётся из накопителя: сумма задержек кадров отслеживает точное
|
||||
// физическое время, а не округляется покадрово (см. DelayDither)
|
||||
frame.delay = self.dither.next(self.plan.exact_delay_cs.max(1.0));
|
||||
self.enc.write_frame(&frame).map_err(std::io::Error::other)?;
|
||||
self.frames += 1;
|
||||
Ok(())
|
||||
}
|
||||
|
||||
/// Контур области измельчения — видно, где считает тонкая сетка.
|
||||
fn draw_patch_outline(&self, buf: &mut [u8], ax: usize, bx: usize, ay: usize, by: usize) {
|
||||
let flip = |y: usize| self.ny - 1 - y.min(self.ny - 1);
|
||||
let mut put = |x: usize, y: usize| {
|
||||
if x < self.nx && y < self.ny {
|
||||
let py = flip(y) * self.scale;
|
||||
let px = x * self.scale;
|
||||
for sy in 0..self.scale {
|
||||
for sx in 0..self.scale {
|
||||
let o = (py + sy) * self.w + px + sx;
|
||||
if o < buf.len() {
|
||||
buf[o] = IDX_PATCH;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
};
|
||||
for x in ax..=bx.min(self.nx - 1) {
|
||||
put(x, ay);
|
||||
put(x, by);
|
||||
}
|
||||
for y in ay..=by.min(self.ny - 1) {
|
||||
put(ax, y);
|
||||
put(bx, y);
|
||||
}
|
||||
}
|
||||
|
||||
pub fn finish(self) -> std::io::Result<u32> {
|
||||
Ok(self.frames)
|
||||
}
|
||||
}
|
||||
|
||||
#[cfg(test)]
|
||||
mod tests {
|
||||
use super::*;
|
||||
|
||||
/// Пример из постановки: 30 м/с, ячейка 0.1 м, u_lat = 1 ⇒ 300 шагов/с; кадр каждые
|
||||
/// 10 шагов ⇒ 30 кадр/с. Точная задержка 3⅓ сотых — нецелая, поэтому включается дизеринг,
|
||||
/// и в среднем гифка идёт ровно в реальном времени.
|
||||
#[test]
|
||||
fn spec_example_is_realtime() {
|
||||
let u = Units::new(0.1, 30.0, 1.0, 150.0, 16.0);
|
||||
assert!((u.steps_per_second() - 300.0).abs() < 1e-9);
|
||||
let p = GifPlan::new(&u, Some(10), 30.0, 1.0);
|
||||
assert!((p.exact_delay_cs - 10.0 / 3.0).abs() < 1e-9, "{}", p.exact_delay_cs);
|
||||
assert!(p.dithered && !p.clamped);
|
||||
assert!((p.fps - 30.0).abs() < 1e-9);
|
||||
assert!(p.is_synced(), "коэффициент {}", p.sync_error());
|
||||
}
|
||||
|
||||
/// Дизеринг обязан вести НАКОПЛЕННОЕ время кадров вплотную к физическому, без ухода.
|
||||
/// Покадровое округление 3.333 → 3 дало бы к тысячному кадру уход в 3.3 секунды.
|
||||
#[test]
|
||||
fn dither_tracks_exact_time_without_drift() {
|
||||
let exact = 10.0 / 3.0;
|
||||
let mut d = DelayDither::default();
|
||||
let mut sum = 0i64;
|
||||
for k in 1..=1000 {
|
||||
sum += d.next(exact) as i64;
|
||||
let want = exact * k as R;
|
||||
assert!(
|
||||
(sum as R - want).abs() <= 0.5 + 1e-9,
|
||||
"кадр {k}: накоплено {sum} сотых против точных {want:.3}"
|
||||
);
|
||||
}
|
||||
// а наивное округление к этому моменту отстало бы на треть сотой на кадр
|
||||
assert!((sum as R - exact * 1000.0).abs() < 1.0);
|
||||
assert!((3000 - sum).abs() > 300, "дизеринг обязан отличаться от покадрового округления");
|
||||
}
|
||||
|
||||
/// Тайминг обязан зависеть ТОЛЬКО от физики, не от того, как быстро считает машина:
|
||||
/// один и тот же δt даёт одну и ту же задержку, вдвое больший шаг — вдвое большую.
|
||||
#[test]
|
||||
fn timing_ignores_wall_clock() {
|
||||
let u = Units::new(0.1, 30.0, 1.0, 150.0, 16.0);
|
||||
let a = GifPlan::new(&u, Some(10), 30.0, 1.0);
|
||||
let b = GifPlan::new(&u, Some(10), 30.0, 1.0);
|
||||
assert!((a.exact_delay_cs - b.exact_delay_cs).abs() < 1e-12);
|
||||
let c = GifPlan::new(&u, Some(20), 30.0, 1.0);
|
||||
assert!((c.exact_delay_cs - 2.0 * a.exact_delay_cs).abs() < 1e-12);
|
||||
assert!(c.is_synced() && a.is_synced());
|
||||
}
|
||||
|
||||
/// Авто-режим обязан вытянуть реальное время при физически корректном u_lat, где жёсткий
|
||||
/// шаг «каждые 10» потребовал бы 600 кадр/с и формат бы этого не вытянул.
|
||||
#[test]
|
||||
fn auto_stride_restores_realtime() {
|
||||
let u = Units::new(0.1, 30.0, 0.05, 150.0, 16.0);
|
||||
assert!((u.steps_per_second() - 6000.0).abs() < 1e-6);
|
||||
let fixed = GifPlan::new(&u, Some(10), 30.0, 1.0);
|
||||
assert!(fixed.clamped, "0.167 сотых обязаны упереться в пол формата");
|
||||
assert!(!fixed.is_synced(), "жёсткий шаг не должен был попасть в реальное время");
|
||||
assert!(fixed.advice(&u).is_some(), "в таком режиме обязан быть совет");
|
||||
let auto = GifPlan::new(&u, None, 30.0, 1.0);
|
||||
assert_eq!(auto.stride, 200); // 1/30 с при 6000 шаг/с
|
||||
assert!(auto.is_synced(), "коэффициент {}", auto.sync_error());
|
||||
assert!((auto.fps - 30.0).abs() < 1e-9);
|
||||
}
|
||||
|
||||
/// Замедление: playback = 0.1 растягивает то же физическое время ровно в 10 раз.
|
||||
#[test]
|
||||
fn slow_motion_scales_delay() {
|
||||
let u = Units::new(0.1, 30.0, 1.0, 150.0, 16.0);
|
||||
let fast = GifPlan::new(&u, Some(10), 30.0, 1.0);
|
||||
let slow = GifPlan::new(&u, Some(10), 30.0, 0.1);
|
||||
assert!((slow.exact_delay_cs - 10.0 * fast.exact_delay_cs).abs() < 1e-9);
|
||||
assert!((slow.actual_playback - 0.1).abs() < 1e-9);
|
||||
assert!(slow.is_synced());
|
||||
}
|
||||
|
||||
/// Палитра обязана быть ровно 256 цветов по 3 байта и не затирать служебные индексы.
|
||||
#[test]
|
||||
fn palette_layout() {
|
||||
for cm in [ColorMap::Viridis, ColorMap::Turbo, ColorMap::CoolWarm] {
|
||||
let p = cm.palette();
|
||||
assert_eq!(p.len(), 768);
|
||||
let o = IDX_BODY as usize * 3;
|
||||
assert_eq!(&p[o..o + 3], &[25, 25, 28]);
|
||||
}
|
||||
}
|
||||
}
|
||||
File diff suppressed because it is too large
Load Diff
File diff suppressed because it is too large
Load Diff
File diff suppressed because it is too large
Load Diff
Reference in New Issue
Block a user