diff --git a/.github/workflows/toolchain-image.yml b/.github/workflows/toolchain-image.yml index 8e8a626..9cc1d69 100644 --- a/.github/workflows/toolchain-image.yml +++ b/.github/workflows/toolchain-image.yml @@ -51,6 +51,7 @@ jobs: labels: ${{ steps.meta.outputs.labels }} build-args: | FLAP_VERSION=1.2.16 + OPENCODE_VERSION=1.3.3 CATCH2_REF=v3.13.0 MDSPAN_REF=mdspan-0.6.0 XTL_REF=0.8.2 diff --git a/.gitignore b/.gitignore index 5ed28a7..0d4faa8 100644 --- a/.gitignore +++ b/.gitignore @@ -8,3 +8,5 @@ __pycache__/ .venv/ venv/ .DS_Store +benchmarks/*/*/workspace/build/ +benchmarks/*/*/workspace/build-review/ diff --git a/.vscode/settings.json b/.vscode/settings.json new file mode 100644 index 0000000..5ee38c9 --- /dev/null +++ b/.vscode/settings.json @@ -0,0 +1,35 @@ +{ + "files.exclude": { + "**/.git/objects/**": true, + "**/.git/refs/**": true, + "**/.pytest_cache/**": true, + "**/__pycache__": true, + "**/.ruff_cache/**": true, + "**/.venv/**": true, + "**/.opencode-data/**": true + }, + "files.watcherExclude": { + "**/runs/**": true, + "**/.git/objects/**": true, + "**/.git/refs/**": true, + "**/.ruff_cache/**": true, + "**/.venv/**": true, + "**/.opencode-data/**": true, + "**/node_modules/**": true, + "**/__pycache__/**": true, + "**/build/**": true + }, + "search.exclude": { + "**/runs/**": true, + "**/.git/**": true, + "**/.ruff_cache/**": true, + "**/.venv/**": true, + "**/.opencode-data/**": true, + "**/node_modules/**": true, + "**/__pycache__/**": true + }, + "C_Cpp.files.exclude": { + "**/build/**": true, + "**/runs/**": true + } +} diff --git a/README.md b/README.md index 35f6361..620c4c9 100644 --- a/README.md +++ b/README.md @@ -30,6 +30,7 @@ Configured in `agents_default.toml` (with per-agent overrides under `sample/*.to ## Available Benchmarks - Demo benchmark for Runge-Kutta 2 (RK2) midpoint method. - 3D wave equation solver with finite difference method. +- Magnetohydrodynamics (MHD) solver. ## Quick Start diff --git a/agents_default.toml b/agents_default.toml index 2f85a21..30651bd 100644 --- a/agents_default.toml +++ b/agents_default.toml @@ -10,8 +10,10 @@ pass_env = [ ] pre = [ "mkdir -p \"$HOME/.local/share/opencode\"", + "mkdir -p \"$HOME/.config/opencode\"", "if [ -f /opencode-auth.json ]; then cp /opencode-auth.json \"$HOME/.local/share/opencode/auth.json\"; fi", - "if [ -f /opencode-config.json ]; then export OPENCODE_CONFIG=/opencode-config.json; elif [ -f /opencode-config.jsonc ]; then export OPENCODE_CONFIG=/opencode-config.jsonc; fi", + "if [ -f /opencode-config.json ]; then cp /opencode-config.json \"$HOME/.config/opencode/opencode.json\"; fi", + "if [ -f /opencode-config.jsonc ]; then cp /opencode-config.jsonc \"$HOME/.config/opencode/opencode.jsonc\"; fi", ] cmd = "stdbuf -oL -eL opencode run -m \"$BENCH_MODEL\" $BENCH_MODEL_OPTIONS_ARGS --thinking --dir / \"$(cat \"$BENCH_PROMPT_FILE\")\" -f \"$BENCH_SPEC_FILE\"" diff --git a/benchmarks/magnetohydrodynamics/README.md b/benchmarks/magnetohydrodynamics/README.md index 86f3848..25db847 100644 --- a/benchmarks/magnetohydrodynamics/README.md +++ b/benchmarks/magnetohydrodynamics/README.md @@ -6,13 +6,29 @@ This suite contains benchmark tasks for ideal magnetohydrodynamics solvers. - `shared/workspace/basic_equations.md`: suite-wide notation and flux conventions. - `shared/workspace/hlld.md`: HLLD algorithm notes for solver tasks. -- `cpp-hlld/`: C++ HLLD approximate Riemann solver task. +- `cpp-hlld-00/`: default C++ HLLD task with detailed solver guidance in spec. +- `cpp-hlld-01/`: variant C++ HLLD task with reduced guidance but same test intent. +- `cpp-full1d-00/`: easiest C++ full 1D ideal MHD variant (main+HLLD provided, solver scaffolded). +- `cpp-full1d-01/`: reduced-guidance full 1D variant with only `evolve_ssp_rk3(...)` exposed. +- `shared/eval/README.md`: hidden-eval contract for shared MHD scoring assets. +- `shared/eval/mhd1d_shared.py`: shared helpers for CSV loading, score + windows, and comparison metadata. +- `shared/eval/fixtures/mhd1d/`: hidden fixtures for full 1D variants. ## Notes - Shared workspace files are visible to the agent during benchmark runs. - Keep maintainer-only derivations, generators, and hidden fixtures outside the shared workspace. +- `cpp-hlld-00` and `cpp-hlld-01` expose only + `hlld_flux_from_primitive(...)` in the public task API. +- `cpp-hlld-00` and `cpp-hlld-01` keep public/hidden test intent aligned; + the main difference is prompt detail level. +- `cpp-full1d-00` public tests compare solver CSV output against a golden file + with numeric tolerance (`1.0e-12`), and hidden tests use `nx=200` against a + hidden reference CSV with the same numeric policy. +- Full 1D tasks score interior cells and emit CSV with lowercase magnetic-field + headers (`by`, `bz`). ## Reference credit diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/eval/run.sh b/benchmarks/magnetohydrodynamics/cpp-full1d-00/eval/run.sh new file mode 100755 index 0000000..5132b67 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/eval/run.sh @@ -0,0 +1,28 @@ +#!/usr/bin/env bash +set -u -o pipefail + +cd /work +export PYTHONPATH="/work:/eval_shared" + +status="passed" +score="1.0" + +python3 -m pytest -q /eval/tests +rc=$? +if [ "$rc" -ne 0 ]; then + status="failed" + score="0.0" +fi + +python3 - < Path: + subprocess.run( + ["cmake", "-S", ".", "-B", str(build_dir)], check=True, cwd=WORKSPACE_ROOT + ) + subprocess.run( + ["cmake", "--build", str(build_dir), "--target", SOLVER_TARGET], + check=True, + cwd=WORKSPACE_ROOT, + ) + + binary_name = f"{SOLVER_TARGET}.exe" if os.name == "nt" else SOLVER_TARGET + solver_path = build_dir / "bin" / binary_name + assert solver_path.exists() + return solver_path + + +def test_hidden_brio_wu_cli_matches_fixture(tmp_path: Path) -> None: + solver_path = _build_solver(tmp_path / "build") + output_csv_path = tmp_path / "brio_wu.csv" + + completed = subprocess.run( + [str(solver_path), "200"], + check=True, + capture_output=True, + text=True, + ) + output_csv_path.write_text(completed.stdout, encoding="utf-8") + + output_rows = list( + csv.reader(output_csv_path.read_text(encoding="utf-8").splitlines()) + ) + reference_rows = list( + csv.reader(REFERENCE_CSV_PATH.read_text(encoding="utf-8").splitlines()) + ) + + assert_csv_rows_close( + output_rows, + reference_rows, + tolerance=TOLERANCE, + expected_header=["x", "rho", "u", "v", "w", "p", "by", "bz"], + ) diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/spec.md b/benchmarks/magnetohydrodynamics/cpp-full1d-00/spec.md new file mode 100644 index 0000000..05a3301 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/spec.md @@ -0,0 +1,74 @@ +# cpp-full1d-00 + +Implement a 1D ideal-MHD solver CLI in C++. + +## Read first + +- `/work/basic_equations.md` +- `/work/hlld.md` + +## Task + +The command-line entrypoint and the HLLD Riemann solver are already provided. +The CLI accepts an integer `nx` argument for the number of grid points, performs the Brio-Wu Riemann problem, and writes the solution to stdout in CSV format. +Your main task is to complete the solver implementation in `src/mhd1d.cpp`. + +The CLI output must match the provided golden CSV for `nx=100` within numeric tolerance (`1.0e-12`). + +## How to test + +Run the public checks from the workspace: + +```bash +python3 -m pytest -q tests/test_public.py +``` + +## Local dev + +```bash +pytest -q +``` + +To build manually: + +```bash +cmake -S . -B build +cmake --build build +./build/bin/cpp_full_solver1d +``` + +## Numerical Algorithm + +- Riemann solver: HLLD +- Primitive variables reconstruction: piecewise linear with MC2 slope limiter +- Time integration: SSP-RK3 +- Boundary condition: symmetric (zero-gradient) + +## Files + +- `src/main.cpp`: complete CLI (already done) +- `src/hlld.hpp`, `src/hlld.cpp`: complete HLLD implementation (already done) +- `src/mhd1d.hpp`, `src/mhd1d.cpp`: solver scaffolding to complete + +## Functions to complete (in `src/mhd1d.cpp`) + +The easiest path is to implement these functions first: + +1. `primitive_to_conservative(...)` +2. `conservative_to_primitive(...)` +3. `compute_lr(...)` +4. `compute_rhs(...)` +5. `push_ssp_rk3(...)` +6. `evolve_ssp_rk3(...)` + +Recommended implementation order: + +1. Reconstruction (`compute_lr`) +2. Flux loop (`compute_flux_hlld` already calls provided HLLD) +3. RHS assembly (`compute_rhs`) +4. One RK3 step (`push_ssp_rk3`) +5. Time loop (`evolve_ssp_rk3`) + +## Standards + +- C++17 diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/task.toml b/benchmarks/magnetohydrodynamics/cpp-full1d-00/task.toml new file mode 100644 index 0000000..171ed4e --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/task.toml @@ -0,0 +1,8 @@ +id = "cpp-full1d-00" +suite = "magnetohydrodynamics" +language = "cpp" +time_limit_sec = 600 +eval_cmd = "/eval/run.sh" +prompt = "Read /run/spec.md, /work/basic_equations.md, and /work/hlld.md, then solve the task in /work." +use_shared_workspace = true +use_shared_eval = true diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/CMakeLists.txt b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/CMakeLists.txt new file mode 100644 index 0000000..1375922 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/CMakeLists.txt @@ -0,0 +1,61 @@ +cmake_minimum_required(VERSION 3.16) + +project(cpp_full_solver1d LANGUAGES CXX) + +set(CMAKE_CXX_STANDARD 17) +set(CMAKE_CXX_STANDARD_REQUIRED ON) +set(CMAKE_CXX_EXTENSIONS OFF) + +include(FetchContent) + +find_package(Catch2 3 QUIET) + +if(NOT Catch2_FOUND) + FetchContent_Declare( + Catch2 + GIT_REPOSITORY https://github.com/catchorg/Catch2.git + GIT_TAG v3.13.0 + GIT_SHALLOW TRUE + ) + FetchContent_MakeAvailable(Catch2) +endif() + +add_library(mhd1d_solver + src/mhd1d.cpp + src/hlld.cpp +) + +target_include_directories(mhd1d_solver PUBLIC + src + ../../../common/include +) + +add_executable(cpp_full_solver1d + src/main.cpp +) + +target_link_libraries(cpp_full_solver1d PRIVATE + mhd1d_solver +) + +set_target_properties(cpp_full_solver1d PROPERTIES + RUNTIME_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/bin" +) + +add_executable(cpp_full_solver1d_public_tests + tests/cpp/test_public.cpp +) + +target_link_libraries(cpp_full_solver1d_public_tests PRIVATE + mhd1d_solver + Catch2::Catch2WithMain +) + +target_include_directories(cpp_full_solver1d_public_tests PRIVATE + src + ../../../common/include +) + +set_target_properties(cpp_full_solver1d_public_tests PROPERTIES + RUNTIME_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/tests" +) diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/README.md b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/README.md new file mode 100644 index 0000000..bda64ad --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/README.md @@ -0,0 +1,7 @@ +The public C++ workspace contains a Brio-Wu solver scaffold. + +- `src/main.cpp` is complete. +- `src/hlld.cpp` is complete. +- `src/mhd1d.cpp` contains TODO sections to implement. + +Shared workspace docs are already mounted for this benchmark. diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/pyproject.toml b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/pyproject.toml new file mode 100644 index 0000000..137476d --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/pyproject.toml @@ -0,0 +1,7 @@ +[project] +name = "magnetohydrodynamics-cpp-full1d-00" +version = "0.0.0" +requires-python = ">=3.10" + +[tool.pytest.ini_options] +testpaths = ["tests"] diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/hlld.cpp b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/hlld.cpp new file mode 100644 index 0000000..f5f9c87 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/hlld.cpp @@ -0,0 +1,202 @@ +#include "hlld.hpp" + +#include +#include + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux) +{ + constexpr double epsilon = 1.0e-40; + + const double rol = left[0]; + const double vxl = left[1]; + const double vyl = left[2]; + const double vzl = left[3]; + const double prl = left[4]; + const double byl = left[5]; + const double bzl = left[6]; + + const double ror = right[0]; + const double vxr = right[1]; + const double vyr = right[2]; + const double vzr = right[3]; + const double prr = right[4]; + const double byr = right[5]; + const double bzr = right[6]; + + const double igm = 1.0 / (gamma - 1.0); + const double bxs = bx; + const double bxsq = bxs * bxs; + + const double pbl = 0.5 * (bxsq + byl * byl + bzl * bzl); + const double pbr = 0.5 * (bxsq + byr * byr + bzr * bzr); + const double ptl = prl + pbl; + const double ptr = prr + pbr; + + const double rxl = rol * vxl; + const double ryl = rol * vyl; + const double rzl = rol * vzl; + const double rxr = ror * vxr; + const double ryr = ror * vyr; + const double rzr = ror * vzr; + + const double eel = prl * igm + 0.5 * (rxl * vxl + ryl * vyl + rzl * vzl) + pbl; + const double eer = prr * igm + 0.5 * (rxr * vxr + ryr * vyr + rzr * vzr) + pbr; + + const double gmpl = gamma * prl; + const double gmpr = gamma * prr; + const double gpbl = gmpl + 2.0 * pbl; + const double gpbr = gmpr + 2.0 * pbr; + + const double cfl = std::sqrt((gpbl + std::sqrt((gmpl - 2.0 * pbl) * (gmpl - 2.0 * pbl) + + 4.0 * gmpl * (byl * byl + bzl * bzl))) * + 0.5 / rol); + const double cfr = std::sqrt((gpbr + std::sqrt((gmpr - 2.0 * pbr) * (gmpr - 2.0 * pbr) + + 4.0 * gmpr * (byr * byr + bzr * bzr))) * + 0.5 / ror); + + const double sl = std::min(vxl, vxr) - std::max(cfl, cfr); + const double sr = std::max(vxl, vxr) + std::max(cfl, cfr); + + const double fql[7] = {rxl, + rxl * vxl + ptl - bxsq, + rxl * vyl - bxs * byl, + rxl * vzl - bxs * bzl, + vxl * (eel + ptl - bxsq) - bxs * (vyl * byl + vzl * bzl), + byl * vxl - bxs * vyl, + bzl * vxl - bxs * vzl}; + const double fqr[7] = {rxr, + rxr * vxr + ptr - bxsq, + rxr * vyr - bxs * byr, + rxr * vzr - bxs * bzr, + vxr * (eer + ptr - bxsq) - bxs * (vyr * byr + vzr * bzr), + byr * vxr - bxs * vyr, + bzr * vxr - bxs * vzr}; + + const double sdl = sl - vxl; + const double sdr = sr - vxr; + const double rosdl = rol * sdl; + const double rosdr = ror * sdr; + const double temp = 1.0 / (rosdr - rosdl); + const double sm = (rosdr * vxr - rosdl * vxl - ptr + ptl) * temp; + const double sdml = sl - sm; + const double sdmr = sr - sm; + const double ptst = (rosdr * ptl - rosdl * ptr + rosdl * rosdr * (vxr - vxl)) * temp; + + const double temp_fst_l = rosdl * sdml - bxsq; + const double sign1_l = std::copysign(1.0, std::abs(temp_fst_l) - epsilon); + const double maxs1_l = std::max(0.0, sign1_l); + const double mins1_l = std::min(0.0, sign1_l); + const double itf_l = 1.0 / (temp_fst_l + mins1_l); + const double isdml = 1.0 / sdml; + + const double temp_l = bxs * (sdl - sdml) * itf_l; + const double rolst = maxs1_l * (rosdl * isdml) - mins1_l * rol; + const double vxlst = maxs1_l * sm - mins1_l * vxl; + const double rxlst = rolst * vxlst; + const double vylst = maxs1_l * (vyl - byl * temp_l) - mins1_l * vyl; + const double rylst = rolst * vylst; + const double vzlst = maxs1_l * (vzl - bzl * temp_l) - mins1_l * vzl; + const double rzlst = rolst * vzlst; + const double temp_l_b = (rosdl * sdl - bxsq) * itf_l; + const double bylst = maxs1_l * (byl * temp_l_b) - mins1_l * byl; + const double bzlst = maxs1_l * (bzl * temp_l_b) - mins1_l * bzl; + const double vdbstl = vxlst * bxs + vylst * bylst + vzlst * bzlst; + const double eelst = maxs1_l * ((sdl * eel - ptl * vxl + ptst * sm + + bxs * (vxl * bxs + vyl * byl + vzl * bzl - vdbstl)) * + isdml) - + mins1_l * eel; + + const double temp_fst_r = rosdr * sdmr - bxsq; + const double sign1_r = std::copysign(1.0, std::abs(temp_fst_r) - epsilon); + const double maxs1_r = std::max(0.0, sign1_r); + const double mins1_r = std::min(0.0, sign1_r); + const double itf_r = 1.0 / (temp_fst_r + mins1_r); + const double isdmr = 1.0 / sdmr; + + const double temp_r = bxs * (sdr - sdmr) * itf_r; + const double rorst = maxs1_r * (rosdr * isdmr) - mins1_r * ror; + const double vxrst = maxs1_r * sm - mins1_r * vxr; + const double rxrst = rorst * vxrst; + const double vyrst = maxs1_r * (vyr - byr * temp_r) - mins1_r * vyr; + const double ryrst = rorst * vyrst; + const double vzrst = maxs1_r * (vzr - bzr * temp_r) - mins1_r * vzr; + const double rzrst = rorst * vzrst; + const double temp_r_b = (rosdr * sdr - bxsq) * itf_r; + const double byrst = maxs1_r * (byr * temp_r_b) - mins1_r * byr; + const double bzrst = maxs1_r * (bzr * temp_r_b) - mins1_r * bzr; + const double vdbstr = vxrst * bxs + vyrst * byrst + vzrst * bzrst; + const double eerst = maxs1_r * ((sdr * eer - ptr * vxr + ptst * sm + + bxs * (vxr * bxs + vyr * byr + vzr * bzr - vdbstr)) * + isdmr) - + mins1_r * eer; + + const double sqrtrol = std::sqrt(rolst); + const double sqrtror = std::sqrt(rorst); + const double abbx = std::abs(bxs); + const double slst = sm - abbx / sqrtrol; + const double srst = sm + abbx / sqrtror; + const double signbx = std::copysign(1.0, bxs); + const double sign1_b = std::copysign(1.0, abbx - epsilon); + const double maxs1_b = std::max(0.0, sign1_b); + const double mins1_b = -std::min(0.0, sign1_b); + const double invsumro = maxs1_b / (sqrtrol + sqrtror); + + const double roldst = rolst; + const double rordst = rorst; + const double rxldst = rxlst; + const double rxrdst = rxrst; + + const double vy_shared = + invsumro * (sqrtrol * vylst + sqrtror * vyrst + signbx * (byrst - bylst)); + const double ryldst = rylst * mins1_b + roldst * vy_shared; + const double ryrdst = ryrst * mins1_b + rordst * vy_shared; + + const double vz_shared = + invsumro * (sqrtrol * vzlst + sqrtror * vzrst + signbx * (bzrst - bzlst)); + const double rzldst = rzlst * mins1_b + roldst * vz_shared; + const double rzrdst = rzrst * mins1_b + rordst * vz_shared; + + const double by_shared = + invsumro * (sqrtrol * byrst + sqrtror * bylst + signbx * sqrtrol * sqrtror * (vyrst - vylst)); + const double byldst = bylst * mins1_b + by_shared; + const double byrdst = byrst * mins1_b + by_shared; + + const double bz_shared = + invsumro * (sqrtrol * bzrst + sqrtror * bzlst + signbx * sqrtrol * sqrtror * (vzrst - vzlst)); + const double bzldst = bzlst * mins1_b + bz_shared; + const double bzrdst = bzrst * mins1_b + bz_shared; + + const double vyldst = vylst * mins1_b + vy_shared; + const double vyrdst = vyrst * mins1_b + vy_shared; + const double vzldst = vzlst * mins1_b + vz_shared; + const double vzrdst = vzrst * mins1_b + vz_shared; + const double temp_dst = sm * bxs + vyldst * byldst + vzldst * bzldst; + const double eeldst = eelst - sqrtrol * signbx * (vdbstl - temp_dst) * maxs1_b; + const double eerdst = eerst + sqrtror * signbx * (vdbstr - temp_dst) * maxs1_b; + + const double sign1 = std::copysign(1.0, sm); + const double maxs1 = std::max(0.0, sign1); + const double mins1 = -std::min(0.0, sign1); + const double msl = std::min(sl, 0.0); + const double mslst = std::min(slst, 0.0); + const double msrst = std::max(srst, 0.0); + const double msr = std::max(sr, 0.0); + const double temp_flux_l = mslst - msl; + const double temp_flux_r = msrst - msr; + + flux[0] = (fql[0] - msl * rol - rolst * temp_flux_l + roldst * mslst) * maxs1 + + (fqr[0] - msr * ror - rorst * temp_flux_r + rordst * msrst) * mins1; + flux[1] = (fql[1] - msl * rxl - rxlst * temp_flux_l + rxldst * mslst) * maxs1 + + (fqr[1] - msr * rxr - rxrst * temp_flux_r + rxrdst * msrst) * mins1; + flux[2] = (fql[2] - msl * ryl - rylst * temp_flux_l + ryldst * mslst) * maxs1 + + (fqr[2] - msr * ryr - ryrst * temp_flux_r + ryrdst * msrst) * mins1; + flux[3] = (fql[3] - msl * rzl - rzlst * temp_flux_l + rzldst * mslst) * maxs1 + + (fqr[3] - msr * rzr - rzrst * temp_flux_r + rzrdst * msrst) * mins1; + flux[4] = (fql[4] - msl * eel - eelst * temp_flux_l + eeldst * mslst) * maxs1 + + (fqr[4] - msr * eer - eerst * temp_flux_r + eerdst * msrst) * mins1; + flux[5] = (fql[5] - msl * byl - bylst * temp_flux_l + byldst * mslst) * maxs1 + + (fqr[5] - msr * byr - byrst * temp_flux_r + byrdst * msrst) * mins1; + flux[6] = (fql[6] - msl * bzl - bzlst * temp_flux_l + bzldst * mslst) * maxs1 + + (fqr[6] - msr * bzr - bzrst * temp_flux_r + bzrdst * msrst) * mins1; +} diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/hlld.hpp b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/hlld.hpp new file mode 100644 index 0000000..ae80ebd --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/hlld.hpp @@ -0,0 +1,4 @@ +#pragma once + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux); diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/main.cpp b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/main.cpp new file mode 100644 index 0000000..82e0f5a --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/main.cpp @@ -0,0 +1,75 @@ +#include "mhd1d.hpp" + +#include +#include +#include +#include + +constexpr int Nx = 100; +constexpr double Gamma = 2.0; +constexpr double Bx = 0.75; +constexpr mhd1d::StateVector LeftPrimitive{ + 1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0, +}; +constexpr mhd1d::StateVector RightPrimitive{ + 0.125, 0.0, 0.0, 0.0, 0.1, -1.0, 0.0, +}; + +int parse_nx(int argc, char** argv) +{ + if (argc <= 1) { + return Nx; + } + + char* end = nullptr; + const long parsed = std::strtol(argv[1], &end, 10); + if (end == argv[1] || *end != '\0' || parsed <= 0) { + throw std::runtime_error("usage: solver [nx]"); + } + return static_cast(parsed); +} + +mhd1d::SolverWorkspace initialize(int nx, double gamma, double bx, + const mhd1d::StateVector& left_state, + const mhd1d::StateVector& right_state) +{ + mhd1d::SolverWorkspace workspace(nx, gamma, bx); + + for (int ix = workspace.Lbx; ix <= workspace.Ubx; ++ix) { + const mhd1d::StateVector& state = (workspace.x(ix) < 0.5) ? left_state : right_state; + for (int component = 0; component < mhd1d::N_Component; ++component) { + workspace.up(ix, component) = state[component]; + } + } + + mhd1d::set_boundary(workspace.up, workspace.up, workspace.Lbx, workspace.Ubx); + mhd1d::convert_primitive_to_conservative(workspace.up, workspace.uc, bx, gamma); + + return workspace; +} + +void write_csv(const mhd1d::SolverWorkspace& workspace, std::ostream& os) +{ + os << "x,rho,u,v,w,p,by,bz\n"; + os << std::setprecision(17); + for (int ix = workspace.Lbx; ix <= workspace.Ubx; ++ix) { + os << workspace.x(ix) << ',' << workspace.up(ix, 0) << ',' << workspace.up(ix, 1) << ',' + << workspace.up(ix, 2) << ',' << workspace.up(ix, 3) << ',' << workspace.up(ix, 4) << ',' + << workspace.up(ix, 5) << ',' << workspace.up(ix, 6) << '\n'; + } +} + +int main(int argc, char** argv) +{ + const int nx = parse_nx(argc, argv); + const double delt = 5.0e-4; + const double tmax = 0.1; + + auto workspace = initialize(nx, Gamma, Bx, LeftPrimitive, RightPrimitive); + + mhd1d::evolve_ssp_rk3(workspace, delt, tmax); + + write_csv(workspace, std::cout); + + return 0; +} diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/mhd1d.cpp b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/mhd1d.cpp new file mode 100644 index 0000000..00e85fa --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/mhd1d.cpp @@ -0,0 +1,108 @@ +#include "mhd1d.hpp" + +namespace mhd1d +{ + +void primitive_to_conservative(const double* primitive, double* conservative, double bx, + double gamma) +{ + (void)primitive; + (void)conservative; + (void)bx; + (void)gamma; + // TODO(student): convert one primitive state [rho,u,v,w,p,By,Bz] to + // conservative [rho,mx,my,mz,E,By,Bz]. +} + +void conservative_to_primitive(const double* conservative, double* primitive, double bx, + double gamma) +{ + (void)conservative; + (void)primitive; + (void)bx; + (void)gamma; + // TODO(student): recover primitive variables from one conservative state. + // Enforce positive density and compute pressure from total energy. +} + +void convert_primitive_to_conservative(ArrayView2D primitive, ArrayView2D conservative, double bx, + double gamma) +{ + (void)primitive; + (void)conservative; + (void)bx; + (void)gamma; + // TODO(student): loop over cells and call primitive_to_conservative. +} + +void convert_conservative_to_primitive(ArrayView2D conservative, ArrayView2D primitive, double bx, + double gamma) +{ + (void)conservative; + (void)primitive; + (void)bx; + (void)gamma; + // TODO(student): loop over cells and call conservative_to_primitive. +} + +void set_boundary_lb(ArrayView2D dst, ArrayView2D src, int lbx) +{ + (void)dst; + (void)src; + (void)lbx; + // TODO(student): copy left interior boundary state into left ghost cells. +} + +void set_boundary_ub(ArrayView2D dst, ArrayView2D src, int ubx) +{ + (void)dst; + (void)src; + (void)ubx; + // TODO(student): copy right interior boundary state into right ghost cells. +} + +void set_boundary(ArrayView2D dst, ArrayView2D src, int lbx, int ubx) +{ + (void)dst; + (void)src; + (void)lbx; + (void)ubx; + // TODO(student): apply both lower and upper zero-gradient boundaries. +} + +void compute_lr(SolverWorkspace& workspace) +{ + (void)workspace; + // TODO(student): compute MC2 reconstructed left/right primitive states on each cell. +} + +void compute_flux_hlld(SolverWorkspace& workspace, double bx, double gamma) +{ + (void)workspace; + (void)bx; + (void)gamma; + // TODO(student): evaluate HLLD interface fluxes using provided hlld_flux_from_primitive. +} + +void compute_rhs(SolverWorkspace& workspace) +{ + (void)workspace; + // TODO(student): build semidiscrete RHS from flux differences and cell width dx. +} + +void push_ssp_rk3(SolverWorkspace& workspace, double dt) +{ + (void)workspace; + (void)dt; + // TODO(student): implement one full SSP-RK3 step (3 substeps). +} + +void evolve_ssp_rk3(SolverWorkspace& workspace, double dt, double t_final) +{ + (void)workspace; + (void)dt; + (void)t_final; + // TODO(student): repeatedly call push_ssp_rk3 until t_final (clip final dt). +} + +} // namespace mhd1d diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/mhd1d.hpp b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/mhd1d.hpp new file mode 100644 index 0000000..6925629 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/src/mhd1d.hpp @@ -0,0 +1,107 @@ +#pragma once + +#include +#include + +#include + +namespace mhd1d +{ + +namespace stdex = std::experimental; + +constexpr int N_Component = 7; +constexpr int N_margin = 1; + +using StateVector = std::array; +using ArrayView1D = stdex::mdspan, stdex::layout_right>; +using ArrayView2D = stdex::mdspan, stdex::layout_right>; + +struct SolverWorkspace { + explicit SolverWorkspace(int nx, double gamma, double bx) + : Nx(nx), Lbx(N_margin), Ubx(N_margin + nx - 1), dx(1.0 / static_cast(nx)), + gamma(gamma), bx(bx), storage(Nx + 2 * N_margin, N_Component) + { + init_views(Nx + 2 * N_margin, N_Component); + + for (int ix = Lbx; ix <= Ubx; ++ix) { + x(ix) = (static_cast(ix - Lbx) + 0.5) * dx; + } + } + + int Nx; + int Lbx; + int Ubx; + double dx; + double gamma; + double bx; + + ArrayView1D x; + ArrayView2D uc; + ArrayView2D up; + ArrayView2D up_l; + ArrayView2D up_r; + ArrayView2D rhs; + ArrayView2D prev; + ArrayView2D flux; + +private: + void init_views(int n_grid, int n_component) + { + x = ArrayView1D(storage.x.data(), n_grid); + uc = ArrayView2D(storage.uc.data(), n_grid, n_component); + up = ArrayView2D(storage.up.data(), n_grid, n_component); + up_l = ArrayView2D(storage.up_l.data(), n_grid, n_component); + up_r = ArrayView2D(storage.up_r.data(), n_grid, n_component); + rhs = ArrayView2D(storage.rhs.data(), n_grid, n_component); + prev = ArrayView2D(storage.prev.data(), n_grid, n_component); + flux = ArrayView2D(storage.flux.data(), n_grid, n_component); + } + + struct Storage { + explicit Storage(int n_grid, int n_component) + : x(n_grid), uc(n_grid * n_component), up(n_grid * n_component), up_l(n_grid * n_component), + up_r(n_grid * n_component), rhs(n_grid * n_component), prev(n_grid * n_component), + flux(n_grid * n_component) + { + } + + std::vector x; + std::vector uc; + std::vector up; + std::vector up_l; + std::vector up_r; + std::vector rhs; + std::vector prev; + std::vector flux; + }; + + Storage storage; +}; + +void primitive_to_conservative(const double* primitive, double* conservative, double bx, + double gamma); + +void conservative_to_primitive(const double* conservative, double* primitive, double bx, + double gamma); + +void convert_primitive_to_conservative(ArrayView2D primitive_cells, ArrayView2D conservative_cells, + double bx, double gamma); + +void set_boundary_lb(ArrayView2D dst, ArrayView2D src, int lbx); + +void set_boundary_ub(ArrayView2D dst, ArrayView2D src, int ubx); + +void set_boundary(ArrayView2D dst, ArrayView2D src, int lbx, int ubx); + +void compute_lr(SolverWorkspace& workspace); + +void compute_flux_hlld(SolverWorkspace& workspace, double bx, double gamma); + +void compute_rhs(SolverWorkspace& workspace); + +void push_ssp_rk3(SolverWorkspace& workspace, double dt); + +void evolve_ssp_rk3(SolverWorkspace& workspace, double dt, double t_final); + +} // namespace mhd1d diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/tests/cpp/test_public.cpp b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/tests/cpp/test_public.cpp new file mode 100644 index 0000000..b4f822f --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/tests/cpp/test_public.cpp @@ -0,0 +1,55 @@ +#include + +#include +#include + +#include "mhd1d.hpp" + +TEST_CASE("set_boundary duplicates edge states on both sides", "[mhd1d][boundary]") +{ + std::vector padded_buffer(4 * mhd1d::N_Component, 0.0); + const mhd1d::ArrayView2D padded(padded_buffer.data(), 4, mhd1d::N_Component); + + padded(1, 0) = 1.2; + padded(2, 0) = 1.4; + + mhd1d::set_boundary_lb(padded, padded, 1); + mhd1d::set_boundary_ub(padded, padded, 2); + + REQUIRE(padded(0, 0) == padded(1, 0)); + REQUIRE(padded(3, 0) == padded(2, 0)); +} + +TEST_CASE("compute_flux_hlld fills finite interface values", "[mhd1d][flux]") +{ + mhd1d::SolverWorkspace workspace(4, 2.0, 0.75); + + for (int ix = workspace.Lbx; ix <= workspace.Ubx; ++ix) { + workspace.up_l(ix, 0) = 1.0; + workspace.up_l(ix, 1) = 0.0; + workspace.up_l(ix, 2) = 0.0; + workspace.up_l(ix, 3) = 0.0; + workspace.up_l(ix, 4) = 1.0; + workspace.up_l(ix, 5) = 0.5; + workspace.up_l(ix, 6) = 0.0; + + workspace.up_r(ix, 0) = 0.9; + workspace.up_r(ix, 1) = 0.0; + workspace.up_r(ix, 2) = 0.0; + workspace.up_r(ix, 3) = 0.0; + workspace.up_r(ix, 4) = 0.9; + workspace.up_r(ix, 5) = 0.4; + workspace.up_r(ix, 6) = 0.0; + } + + mhd1d::set_boundary_lb(workspace.up_l, workspace.up_l, workspace.Lbx); + mhd1d::set_boundary_ub(workspace.up_l, workspace.up_l, workspace.Ubx); + mhd1d::set_boundary_lb(workspace.up_r, workspace.up_r, workspace.Lbx); + mhd1d::set_boundary_ub(workspace.up_r, workspace.up_r, workspace.Ubx); + + mhd1d::compute_flux_hlld(workspace, 0.75, 2.0); + + for (int ix = workspace.Lbx - 1; ix <= workspace.Ubx + 1; ++ix) { + REQUIRE(std::isfinite(workspace.flux(ix, 0))); + } +} diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/tests/data/brio_wu_golden.csv b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/tests/data/brio_wu_golden.csv new file mode 100644 index 0000000..e72dfa5 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/tests/data/brio_wu_golden.csv @@ -0,0 +1,101 @@ +x,rho,u,v,w,p,by,bz +0.0050000000000000001,1,0,0,0,1,1,0 +0.014999999999999999,1,0,0,0,1,1,0 +0.025000000000000001,1,0,0,0,1,1,0 +0.035000000000000003,1,0,0,0,1,1,0 +0.044999999999999998,1,0,0,0,1,1,0 +0.055,1,0,0,0,1,1,0 +0.065000000000000002,1,0,0,0,1,1,0 +0.074999999999999997,1,0,0,0,1,1,0 +0.085000000000000006,1,0,0,0,1,1,0 +0.095000000000000001,1,0,0,0,1,1,0 +0.105,1,0,0,0,1,1,0 +0.115,1,0,0,0,1,1,0 +0.125,1,0,0,0,1,1,0 +0.13500000000000001,1,0,0,0,1,1,0 +0.14499999999999999,1,0,0,0,1,1,0 +0.155,1,0,0,0,1,1,0 +0.16500000000000001,1,0,0,0,1,1,0 +0.17500000000000002,1,0,0,0,1,1,0 +0.185,1,0,0,0,1,1,0 +0.19500000000000001,1,0,0,0,1,1,0 +0.20500000000000002,1,0,0,0,1,1,0 +0.215,1,0,0,0,1,1,0 +0.22500000000000001,1,0,0,0,1,1,0 +0.23500000000000001,1,0,0,0,1,1,0 +0.245,1,0,0,0,1,1,0 +0.255,1,0,0,0,1,1,0 +0.26500000000000001,1,0,0,0,1,1,0 +0.27500000000000002,1,0,0,0,1,1,0 +0.28500000000000003,1,0,0,0,1,1,0 +0.29499999999999998,1,0,0,0,1,1,0 +0.30499999999999999,1,0,0,0,1,1,0 +0.315,0.99258917532833379,0.013274302361789324,-0.0037458540886650027,0,0.98532404108866567,0.99104740522973966,0 +0.32500000000000001,0.96954607581582897,0.054883618137207428,-0.015827325460063591,0,0.94023568077652442,0.96299689161606516,0 +0.33500000000000002,0.94025803324835733,0.10865479486043322,-0.031944255536031786,0,0.88435499782767446,0.92711462956015256,0 +0.34500000000000003,0.90897792466167582,0.16698841420781432,-0.050131763210425707,0,0.82652179164341721,0.88849866375884556,0 +0.35499999999999998,0.87725610844274526,0.2271603117432259,-0.069719685021130948,0,0.76985377322754589,0.84898878620884088,0 +0.36499999999999999,0.84566618167244145,0.28816698962982518,-0.090519447417695451,0,0.71541656913426399,0.80923410080538904,0 +0.375,0.81450435463148285,0.3494865276909434,-0.11251682956314066,0,0.66367321125402712,0.76950238387022263,0 +0.38500000000000001,0.78404538979850258,0.41062539085171551,-0.13575856770700617,0,0.61497905337973746,0.72997249890053162,0 +0.39500000000000002,0.75467120186953685,0.47083530429296916,-0.1602094195276916,0,0.56978778527926832,0.69094510148101218,0 +0.40500000000000003,0.72701668590610768,0.52870556766657628,-0.18539026277804579,0,0.52883750835390897,0.65319266964541056,0 +0.41500000000000004,0.70216500224646816,0.58153016354527853,-0.20973136326536138,0,0.4933637549672012,0.61876152226558601,0 +0.42499999999999999,0.68193441338380167,0.62487275671749021,-0.22845986767252549,0,0.46540195161446823,0.59219206428336624,0 +0.435,0.6666986977176832,0.65196598708089948,-0.23211489711409145,0,0.44490074049787121,0.58021922926273972,0 +0.44500000000000001,0.66015529724292976,0.66607074009479628,-0.22555655390938342,0,0.43615912199888263,0.58524191570556205,0 +0.45500000000000002,0.65997980757445474,0.66324409170331222,-0.23280516118974198,0,0.43727412310940939,0.5796363680769816,0 +0.46500000000000002,0.73349815599883794,0.60086667962898554,-0.62171855754180494,0,0.57886068825314341,0.27652712488086967,0 +0.47500000000000003,0.79237811897791921,0.46454032083596353,-1.2161018009267417,0,0.6842784144671038,-0.24404258805909654,0 +0.48499999999999999,0.74259318372203176,0.53084767864673321,-1.4743018340823315,0,0.60002741951352012,-0.44224666294204562,0 +0.495,0.69493377959694547,0.61399418172955245,-1.5431753009247675,0,0.54413704626033854,-0.54049642839751277,0 +0.505,0.70278259276020427,0.61662067481286442,-1.5872354257997463,0,0.52068648102692328,-0.56052646550215834,0 +0.51500000000000001,0.70390914156072559,0.6077114699355689,-1.6001701751263477,0,0.51193399003676154,-0.54226823182770345,0 +0.52500000000000002,0.69922181544266515,0.59340751947871018,-1.6046948629200928,0,0.51300004963251122,-0.52185059484429897,0 +0.53500000000000003,0.67807490745420074,0.5878535845789149,-1.6006935764888095,0,0.51466016361553968,-0.51621274881535362,0 +0.54500000000000004,0.60556902444692728,0.59045915738825405,-1.5906721688577199,0,0.51204702934936586,-0.52148853682254903,0 +0.55500000000000005,0.48266830107336006,0.59697028797374452,-1.5761776779942758,0,0.51057551463380291,-0.53191252428398161,0 +0.56500000000000006,0.34834934855591537,0.59188322457049136,-1.5687669989374082,0,0.5151520944223198,-0.54463835756014178,0 +0.57500000000000007,0.2451691037246034,0.59463166677984525,-1.5635053912378623,0,0.51994523862074438,-0.55151343157376009,0 +0.58499999999999996,0.22308154337733685,0.61421260394352284,-1.5558407025582697,0,0.51450769544894648,-0.54648231691723681,0 +0.59499999999999997,0.22390587217315483,0.62861601814674595,-1.555379908924208,0,0.51042675178345809,-0.53931383334019445,0 +0.60499999999999998,0.22786331744185459,0.61304302077924422,-1.5697800797404926,0,0.51789637516013132,-0.53999599169525736,0 +0.61499999999999999,0.23452030041149688,0.58777367268577241,-1.5968680437122895,0,0.53659755807779264,-0.54285121456661267,0 +0.625,0.23506241740122699,0.58441088978544375,-1.6180990627429224,0,0.5312595858066933,-0.53847776600605912,0 +0.63500000000000001,0.22912802575710367,0.64096361665071178,-1.5587196384023227,0,0.50180683386498415,-0.51954980419936259,0 +0.64500000000000002,0.20820170172948924,0.4429088168703938,-1.3138217298639516,0,0.4107172773448764,-0.63242114973177121,0 +0.65500000000000003,0.14280420150766387,0.0013572736674349979,-0.58334822550752841,0,0.15933600494306066,-0.82699909098420499,0 +0.66500000000000004,0.11735885521106922,-0.23309434559311201,-0.17605497312958984,0,0.088246122408452865,-0.90091689844612777,0 +0.67500000000000004,0.11691966169225833,-0.24239401825294279,-0.16901850205441424,0,0.087525697429624683,-0.90138999879178838,0 +0.68500000000000005,0.1168975018076187,-0.24275168390313129,-0.16913431059427403,0,0.087491282208618681,-0.90130727955920587,0 +0.69500000000000006,0.11688490837927806,-0.24293963466430563,-0.16922798396131322,0,0.087471166714507165,-0.90122308503037618,0 +0.70499999999999996,0.11687882431088358,-0.24305091000018803,-0.16930677036965813,0,0.087461000889423768,-0.90113549166198792,0 +0.71499999999999997,0.11688210335016401,-0.24310898990635241,-0.16941137635009459,0,0.087464485048653229,-0.90109713084582599,0 +0.72499999999999998,0.11688746427801196,-0.24311790246317952,-0.16948241337481523,0,0.087471346584835463,-0.90112657201926949,0 +0.73499999999999999,0.11689293842702393,-0.24301542095697384,-0.16941123056591925,0,0.087478382198893812,-0.90117181152414561,0 +0.745,0.11689754429437521,-0.24279548243457205,-0.16922004702154958,0,0.087484264982216398,-0.90123329716500089,0 +0.755,0.11690472407462463,-0.24244042822845885,-0.16894496953161006,0,0.087493985877884928,-0.90135933427866255,0 +0.76500000000000001,0.1169351752430936,-0.24164401100336819,-0.16834345193351358,0,0.087538755474644847,-0.90173628521210369,0 +0.77500000000000002,0.11702346553875453,-0.2389118539454739,-0.16634253632348067,0,0.087670236260091849,-0.90283181009609326,0 +0.78500000000000003,0.11723705803609173,-0.23241448265557926,-0.16159143199373074,0,0.087990348375741601,-0.90548092874079178,0 +0.79500000000000004,0.1176285175883393,-0.22051415125476836,-0.15292012220738538,0,0.088578608675684678,-0.91032952865801275,0 +0.80500000000000005,0.11822301685274909,-0.20247335890203855,-0.13985941063243096,0,0.089475764221054743,-0.91767747139414524,0 +0.81500000000000006,0.11901085894402198,-0.178638590244622,-0.12276084795550267,0,0.090670936816770853,-0.92738054886244792,0 +0.82500000000000007,0.11995507930331456,-0.15017685276699991,-0.10257309790892337,0,0.092113203419944889,-0.93896119897869035,0 +0.83499999999999996,0.12100582881114735,-0.1186354504616307,-0.080489249938489674,0,0.093731854039529328,-0.95178772370322251,0 +0.84499999999999997,0.12210948952666369,-0.085657286840511535,-0.057715619723010918,0,0.095446790184988584,-0.96519228251210643,0 +0.85499999999999998,0.12320528139333707,-0.053056324701219029,-0.035509818552657238,0,0.097164687932375315,-0.97843939905855704,0 +0.86499999999999999,0.12419939842673788,-0.023611385320102929,-0.015706138420103915,0,0.098733750673405729,-0.99040292502414518,0 +0.875,0.1248831018288558,-0.0034412091825772586,-0.0022747475485175887,0,0.099815318398622899,-0.99860261559100016,0 +0.88500000000000001,0.12500000000700623,2.0693358136986429e-10,1.3690028369419877e-10,0,0.10000000001121101,-1.0000000000838631,0 +0.89500000000000002,0.12500000000690034,2.0295860176552916e-10,1.3426844820417816e-10,0,0.10000000001104148,-1.0000000000825946,0 +0.90500000000000003,0.12500000000616657,1.8175890421477852e-10,1.2024374828573813e-10,0,0.10000000000986675,-1.0000000000738092,0 +0.91500000000000004,0.12500000000447511,1.3188042879502353e-10,8.7246950401911659e-11,0,0.10000000000716036,-1.0000000000535656,0 +0.92500000000000004,0.1250000000023313,6.8702813807084164e-11,4.5451087337373151e-11,0,0.10000000000372944,-1.0000000000279046,0 +0.93500000000000005,0.12500000000053466,1.575462353161844e-11,1.0423025405688057e-11,0,0.1000000000008554,-1.0000000000063987,0 +0.94500000000000006,0.125,0,0,0,0.099999999999999867,-1,0 +0.95500000000000007,0.125,0,0,0,0.099999999999999867,-1,0 +0.96499999999999997,0.125,0,0,0,0.099999999999999867,-1,0 +0.97499999999999998,0.125,0,0,0,0.099999999999999867,-1,0 +0.98499999999999999,0.125,0,0,0,0.099999999999999867,-1,0 +0.995,0.125,0,0,0,0.099999999999999867,-1,0 diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/tests/test_public.py b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/tests/test_public.py new file mode 100644 index 0000000..24e62c7 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-00/workspace/tests/test_public.py @@ -0,0 +1,81 @@ +import csv +import os +import subprocess +from pathlib import Path + +PUBLIC_TEST_TARGET = "cpp_full_solver1d_public_tests" +GOLDEN_CSV_PATH = Path(__file__).resolve().parents[1] / "tests/data/brio_wu_golden.csv" +TOLERANCE = 1.0e-12 +WORKSPACE_ROOT = Path(__file__).resolve().parents[1] + + +def assert_csv_rows_close( + output_rows: list[list[str]], + reference_rows: list[list[str]], + *, + tolerance: float = TOLERANCE, + expected_header: list[str] | None = None, +) -> None: + if expected_header is not None: + assert output_rows[0] == expected_header + output_rows = output_rows[1:] + + assert len(output_rows) == len(reference_rows) + + for output_row, reference_row in zip(output_rows, reference_rows): + assert len(output_row) == len(reference_row) + for output_value, reference_value in zip(output_row, reference_row): + assert abs(float(output_value) - float(reference_value)) <= tolerance + + +def _build_public_tests() -> Path: + build_dir = "build" + subprocess.run(["cmake", "-S", ".", "-B", build_dir], check=True, cwd=WORKSPACE_ROOT) + subprocess.run( + [ + "cmake", + "--build", + build_dir, + "--target", + "cpp_full_solver1d", + PUBLIC_TEST_TARGET, + ], + check=True, + cwd=WORKSPACE_ROOT, + ) + + binary_name = f"{PUBLIC_TEST_TARGET}.exe" if os.name == "nt" else PUBLIC_TEST_TARGET + executable_path = WORKSPACE_ROOT / build_dir / "tests" / binary_name + assert executable_path.exists() + return executable_path + + +def test_public_catch2_target_builds() -> None: + _build_public_tests() + + +def test_public_brio_wu_cli_matches_reference_grid() -> None: + _build_public_tests() + + solver_name = "cpp_full_solver1d.exe" if os.name == "nt" else "cpp_full_solver1d" + solver_path = WORKSPACE_ROOT / "build/bin" / solver_name + assert solver_path.exists() + + completed = subprocess.run( + [str(solver_path)], + check=True, + capture_output=True, + text=True, + ) + + rows = list(csv.reader(completed.stdout.splitlines())) + golden_rows = list(csv.reader(GOLDEN_CSV_PATH.read_text(encoding="utf-8").splitlines())) + + assert golden_rows[0] == ["x", "rho", "u", "v", "w", "p", "by", "bz"] + + assert_csv_rows_close( + rows, + golden_rows[1:], + tolerance=TOLERANCE, + expected_header=["x", "rho", "u", "v", "w", "p", "by", "bz"], + ) diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/eval/run.sh b/benchmarks/magnetohydrodynamics/cpp-full1d-01/eval/run.sh new file mode 100755 index 0000000..5132b67 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/eval/run.sh @@ -0,0 +1,28 @@ +#!/usr/bin/env bash +set -u -o pipefail + +cd /work +export PYTHONPATH="/work:/eval_shared" + +status="passed" +score="1.0" + +python3 -m pytest -q /eval/tests +rc=$? +if [ "$rc" -ne 0 ]; then + status="failed" + score="0.0" +fi + +python3 - < Path: + subprocess.run( + ["cmake", "-S", ".", "-B", str(build_dir)], check=True, cwd=WORKSPACE_ROOT + ) + subprocess.run( + ["cmake", "--build", str(build_dir), "--target", SOLVER_TARGET], + check=True, + cwd=WORKSPACE_ROOT, + ) + + binary_name = f"{SOLVER_TARGET}.exe" if os.name == "nt" else SOLVER_TARGET + solver_path = build_dir / "bin" / binary_name + assert solver_path.exists() + return solver_path + + +def test_hidden_brio_wu_cli_matches_fixture(tmp_path: Path) -> None: + solver_path = _build_solver(tmp_path / "build") + output_csv_path = tmp_path / "brio_wu.csv" + + completed = subprocess.run( + [str(solver_path), "200"], + check=True, + capture_output=True, + text=True, + ) + output_csv_path.write_text(completed.stdout, encoding="utf-8") + + output_rows = list( + csv.reader(output_csv_path.read_text(encoding="utf-8").splitlines()) + ) + reference_rows = list( + csv.reader(REFERENCE_CSV_PATH.read_text(encoding="utf-8").splitlines()) + ) + + assert_csv_rows_close( + output_rows, + reference_rows, + tolerance=TOLERANCE, + expected_header=["x", "rho", "u", "v", "w", "p", "by", "bz"], + ) diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/spec.md b/benchmarks/magnetohydrodynamics/cpp-full1d-01/spec.md new file mode 100644 index 0000000..3f5ae01 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/spec.md @@ -0,0 +1,62 @@ +# cpp-full1d-01 + +Implement a 1D ideal-MHD solver CLI in C++. + +## Read first + +- `/work/basic_equations.md` +- `/work/hlld.md` + +## Task + +The command-line entrypoint and the HLLD Riemann solver are already provided. +The CLI accepts an integer `nx` argument for the number of grid points, performs the Brio-Wu Riemann problem, and writes the solution to stdout in CSV format. +Your main task is to complete the solver implementation in `src/mhd1d.cpp`. + +The CLI output must match the provided golden CSV for `nx=100` within numeric tolerance (`1.0e-12`). + +## How to test + +Run the public checks from the workspace: + +```bash +python3 -m pytest -q tests/test_public.py +``` + +## Local dev + +```bash +pytest -q +``` + +To build manually: + +```bash +cmake -S . -B build +cmake --build build +./build/bin/cpp_full_solver1d +``` + +## Numerical Algorithm + +- Riemann solver: HLLD +- Primitive variables reconstruction: piecewise linear with MC2 slope limiter +- Time integration: SSP-RK3 +- Boundary condition: symmetric (zero-gradient) + +## Files + +- `src/main.cpp`: complete CLI (already done) +- `src/hlld.hpp`, `src/hlld.cpp`: complete HLLD implementation (already done) +- `src/mhd1d.hpp`: provides the workspace data structure and solver entrypoint declaration +- `src/mhd1d.cpp`: contains only an empty `evolve_ssp_rk3(...)` implementation to complete + +## Functions to complete (in `src/mhd1d.cpp`) + +- `evolve_ssp_rk3(SolverWorkspace& workspace, double dt, double t_final)` + +All other helper functions and internal organization are up to you. + +## Standards + +- C++17 diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/task.toml b/benchmarks/magnetohydrodynamics/cpp-full1d-01/task.toml new file mode 100644 index 0000000..c86d996 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/task.toml @@ -0,0 +1,8 @@ +id = "cpp-full1d-01" +suite = "magnetohydrodynamics" +language = "cpp" +time_limit_sec = 600 +eval_cmd = "/eval/run.sh" +prompt = "Read /run/spec.md, /work/basic_equations.md, and /work/hlld.md, then solve the task in /work." +use_shared_workspace = true +use_shared_eval = true diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/CMakeLists.txt b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/CMakeLists.txt new file mode 100644 index 0000000..1375922 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/CMakeLists.txt @@ -0,0 +1,61 @@ +cmake_minimum_required(VERSION 3.16) + +project(cpp_full_solver1d LANGUAGES CXX) + +set(CMAKE_CXX_STANDARD 17) +set(CMAKE_CXX_STANDARD_REQUIRED ON) +set(CMAKE_CXX_EXTENSIONS OFF) + +include(FetchContent) + +find_package(Catch2 3 QUIET) + +if(NOT Catch2_FOUND) + FetchContent_Declare( + Catch2 + GIT_REPOSITORY https://github.com/catchorg/Catch2.git + GIT_TAG v3.13.0 + GIT_SHALLOW TRUE + ) + FetchContent_MakeAvailable(Catch2) +endif() + +add_library(mhd1d_solver + src/mhd1d.cpp + src/hlld.cpp +) + +target_include_directories(mhd1d_solver PUBLIC + src + ../../../common/include +) + +add_executable(cpp_full_solver1d + src/main.cpp +) + +target_link_libraries(cpp_full_solver1d PRIVATE + mhd1d_solver +) + +set_target_properties(cpp_full_solver1d PROPERTIES + RUNTIME_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/bin" +) + +add_executable(cpp_full_solver1d_public_tests + tests/cpp/test_public.cpp +) + +target_link_libraries(cpp_full_solver1d_public_tests PRIVATE + mhd1d_solver + Catch2::Catch2WithMain +) + +target_include_directories(cpp_full_solver1d_public_tests PRIVATE + src + ../../../common/include +) + +set_target_properties(cpp_full_solver1d_public_tests PROPERTIES + RUNTIME_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/tests" +) diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/README.md b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/README.md new file mode 100644 index 0000000..df83a53 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/README.md @@ -0,0 +1,7 @@ +The public C++ workspace contains a Brio-Wu solver scaffold. + +- `src/main.cpp` is complete. +- `src/hlld.cpp` is complete. +- `src/mhd1d.cpp` exposes only an empty `evolve_ssp_rk3(...)` to implement. + +Shared workspace docs are already mounted for this benchmark. diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/pyproject.toml b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/pyproject.toml new file mode 100644 index 0000000..b38b758 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/pyproject.toml @@ -0,0 +1,7 @@ +[project] +name = "magnetohydrodynamics-cpp-full1d-01" +version = "0.0.0" +requires-python = ">=3.10" + +[tool.pytest.ini_options] +testpaths = ["tests"] diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/hlld.cpp b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/hlld.cpp new file mode 100644 index 0000000..f5f9c87 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/hlld.cpp @@ -0,0 +1,202 @@ +#include "hlld.hpp" + +#include +#include + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux) +{ + constexpr double epsilon = 1.0e-40; + + const double rol = left[0]; + const double vxl = left[1]; + const double vyl = left[2]; + const double vzl = left[3]; + const double prl = left[4]; + const double byl = left[5]; + const double bzl = left[6]; + + const double ror = right[0]; + const double vxr = right[1]; + const double vyr = right[2]; + const double vzr = right[3]; + const double prr = right[4]; + const double byr = right[5]; + const double bzr = right[6]; + + const double igm = 1.0 / (gamma - 1.0); + const double bxs = bx; + const double bxsq = bxs * bxs; + + const double pbl = 0.5 * (bxsq + byl * byl + bzl * bzl); + const double pbr = 0.5 * (bxsq + byr * byr + bzr * bzr); + const double ptl = prl + pbl; + const double ptr = prr + pbr; + + const double rxl = rol * vxl; + const double ryl = rol * vyl; + const double rzl = rol * vzl; + const double rxr = ror * vxr; + const double ryr = ror * vyr; + const double rzr = ror * vzr; + + const double eel = prl * igm + 0.5 * (rxl * vxl + ryl * vyl + rzl * vzl) + pbl; + const double eer = prr * igm + 0.5 * (rxr * vxr + ryr * vyr + rzr * vzr) + pbr; + + const double gmpl = gamma * prl; + const double gmpr = gamma * prr; + const double gpbl = gmpl + 2.0 * pbl; + const double gpbr = gmpr + 2.0 * pbr; + + const double cfl = std::sqrt((gpbl + std::sqrt((gmpl - 2.0 * pbl) * (gmpl - 2.0 * pbl) + + 4.0 * gmpl * (byl * byl + bzl * bzl))) * + 0.5 / rol); + const double cfr = std::sqrt((gpbr + std::sqrt((gmpr - 2.0 * pbr) * (gmpr - 2.0 * pbr) + + 4.0 * gmpr * (byr * byr + bzr * bzr))) * + 0.5 / ror); + + const double sl = std::min(vxl, vxr) - std::max(cfl, cfr); + const double sr = std::max(vxl, vxr) + std::max(cfl, cfr); + + const double fql[7] = {rxl, + rxl * vxl + ptl - bxsq, + rxl * vyl - bxs * byl, + rxl * vzl - bxs * bzl, + vxl * (eel + ptl - bxsq) - bxs * (vyl * byl + vzl * bzl), + byl * vxl - bxs * vyl, + bzl * vxl - bxs * vzl}; + const double fqr[7] = {rxr, + rxr * vxr + ptr - bxsq, + rxr * vyr - bxs * byr, + rxr * vzr - bxs * bzr, + vxr * (eer + ptr - bxsq) - bxs * (vyr * byr + vzr * bzr), + byr * vxr - bxs * vyr, + bzr * vxr - bxs * vzr}; + + const double sdl = sl - vxl; + const double sdr = sr - vxr; + const double rosdl = rol * sdl; + const double rosdr = ror * sdr; + const double temp = 1.0 / (rosdr - rosdl); + const double sm = (rosdr * vxr - rosdl * vxl - ptr + ptl) * temp; + const double sdml = sl - sm; + const double sdmr = sr - sm; + const double ptst = (rosdr * ptl - rosdl * ptr + rosdl * rosdr * (vxr - vxl)) * temp; + + const double temp_fst_l = rosdl * sdml - bxsq; + const double sign1_l = std::copysign(1.0, std::abs(temp_fst_l) - epsilon); + const double maxs1_l = std::max(0.0, sign1_l); + const double mins1_l = std::min(0.0, sign1_l); + const double itf_l = 1.0 / (temp_fst_l + mins1_l); + const double isdml = 1.0 / sdml; + + const double temp_l = bxs * (sdl - sdml) * itf_l; + const double rolst = maxs1_l * (rosdl * isdml) - mins1_l * rol; + const double vxlst = maxs1_l * sm - mins1_l * vxl; + const double rxlst = rolst * vxlst; + const double vylst = maxs1_l * (vyl - byl * temp_l) - mins1_l * vyl; + const double rylst = rolst * vylst; + const double vzlst = maxs1_l * (vzl - bzl * temp_l) - mins1_l * vzl; + const double rzlst = rolst * vzlst; + const double temp_l_b = (rosdl * sdl - bxsq) * itf_l; + const double bylst = maxs1_l * (byl * temp_l_b) - mins1_l * byl; + const double bzlst = maxs1_l * (bzl * temp_l_b) - mins1_l * bzl; + const double vdbstl = vxlst * bxs + vylst * bylst + vzlst * bzlst; + const double eelst = maxs1_l * ((sdl * eel - ptl * vxl + ptst * sm + + bxs * (vxl * bxs + vyl * byl + vzl * bzl - vdbstl)) * + isdml) - + mins1_l * eel; + + const double temp_fst_r = rosdr * sdmr - bxsq; + const double sign1_r = std::copysign(1.0, std::abs(temp_fst_r) - epsilon); + const double maxs1_r = std::max(0.0, sign1_r); + const double mins1_r = std::min(0.0, sign1_r); + const double itf_r = 1.0 / (temp_fst_r + mins1_r); + const double isdmr = 1.0 / sdmr; + + const double temp_r = bxs * (sdr - sdmr) * itf_r; + const double rorst = maxs1_r * (rosdr * isdmr) - mins1_r * ror; + const double vxrst = maxs1_r * sm - mins1_r * vxr; + const double rxrst = rorst * vxrst; + const double vyrst = maxs1_r * (vyr - byr * temp_r) - mins1_r * vyr; + const double ryrst = rorst * vyrst; + const double vzrst = maxs1_r * (vzr - bzr * temp_r) - mins1_r * vzr; + const double rzrst = rorst * vzrst; + const double temp_r_b = (rosdr * sdr - bxsq) * itf_r; + const double byrst = maxs1_r * (byr * temp_r_b) - mins1_r * byr; + const double bzrst = maxs1_r * (bzr * temp_r_b) - mins1_r * bzr; + const double vdbstr = vxrst * bxs + vyrst * byrst + vzrst * bzrst; + const double eerst = maxs1_r * ((sdr * eer - ptr * vxr + ptst * sm + + bxs * (vxr * bxs + vyr * byr + vzr * bzr - vdbstr)) * + isdmr) - + mins1_r * eer; + + const double sqrtrol = std::sqrt(rolst); + const double sqrtror = std::sqrt(rorst); + const double abbx = std::abs(bxs); + const double slst = sm - abbx / sqrtrol; + const double srst = sm + abbx / sqrtror; + const double signbx = std::copysign(1.0, bxs); + const double sign1_b = std::copysign(1.0, abbx - epsilon); + const double maxs1_b = std::max(0.0, sign1_b); + const double mins1_b = -std::min(0.0, sign1_b); + const double invsumro = maxs1_b / (sqrtrol + sqrtror); + + const double roldst = rolst; + const double rordst = rorst; + const double rxldst = rxlst; + const double rxrdst = rxrst; + + const double vy_shared = + invsumro * (sqrtrol * vylst + sqrtror * vyrst + signbx * (byrst - bylst)); + const double ryldst = rylst * mins1_b + roldst * vy_shared; + const double ryrdst = ryrst * mins1_b + rordst * vy_shared; + + const double vz_shared = + invsumro * (sqrtrol * vzlst + sqrtror * vzrst + signbx * (bzrst - bzlst)); + const double rzldst = rzlst * mins1_b + roldst * vz_shared; + const double rzrdst = rzrst * mins1_b + rordst * vz_shared; + + const double by_shared = + invsumro * (sqrtrol * byrst + sqrtror * bylst + signbx * sqrtrol * sqrtror * (vyrst - vylst)); + const double byldst = bylst * mins1_b + by_shared; + const double byrdst = byrst * mins1_b + by_shared; + + const double bz_shared = + invsumro * (sqrtrol * bzrst + sqrtror * bzlst + signbx * sqrtrol * sqrtror * (vzrst - vzlst)); + const double bzldst = bzlst * mins1_b + bz_shared; + const double bzrdst = bzrst * mins1_b + bz_shared; + + const double vyldst = vylst * mins1_b + vy_shared; + const double vyrdst = vyrst * mins1_b + vy_shared; + const double vzldst = vzlst * mins1_b + vz_shared; + const double vzrdst = vzrst * mins1_b + vz_shared; + const double temp_dst = sm * bxs + vyldst * byldst + vzldst * bzldst; + const double eeldst = eelst - sqrtrol * signbx * (vdbstl - temp_dst) * maxs1_b; + const double eerdst = eerst + sqrtror * signbx * (vdbstr - temp_dst) * maxs1_b; + + const double sign1 = std::copysign(1.0, sm); + const double maxs1 = std::max(0.0, sign1); + const double mins1 = -std::min(0.0, sign1); + const double msl = std::min(sl, 0.0); + const double mslst = std::min(slst, 0.0); + const double msrst = std::max(srst, 0.0); + const double msr = std::max(sr, 0.0); + const double temp_flux_l = mslst - msl; + const double temp_flux_r = msrst - msr; + + flux[0] = (fql[0] - msl * rol - rolst * temp_flux_l + roldst * mslst) * maxs1 + + (fqr[0] - msr * ror - rorst * temp_flux_r + rordst * msrst) * mins1; + flux[1] = (fql[1] - msl * rxl - rxlst * temp_flux_l + rxldst * mslst) * maxs1 + + (fqr[1] - msr * rxr - rxrst * temp_flux_r + rxrdst * msrst) * mins1; + flux[2] = (fql[2] - msl * ryl - rylst * temp_flux_l + ryldst * mslst) * maxs1 + + (fqr[2] - msr * ryr - ryrst * temp_flux_r + ryrdst * msrst) * mins1; + flux[3] = (fql[3] - msl * rzl - rzlst * temp_flux_l + rzldst * mslst) * maxs1 + + (fqr[3] - msr * rzr - rzrst * temp_flux_r + rzrdst * msrst) * mins1; + flux[4] = (fql[4] - msl * eel - eelst * temp_flux_l + eeldst * mslst) * maxs1 + + (fqr[4] - msr * eer - eerst * temp_flux_r + eerdst * msrst) * mins1; + flux[5] = (fql[5] - msl * byl - bylst * temp_flux_l + byldst * mslst) * maxs1 + + (fqr[5] - msr * byr - byrst * temp_flux_r + byrdst * msrst) * mins1; + flux[6] = (fql[6] - msl * bzl - bzlst * temp_flux_l + bzldst * mslst) * maxs1 + + (fqr[6] - msr * bzr - bzrst * temp_flux_r + bzrdst * msrst) * mins1; +} diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/hlld.hpp b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/hlld.hpp new file mode 100644 index 0000000..ae80ebd --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/hlld.hpp @@ -0,0 +1,4 @@ +#pragma once + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux); diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/main.cpp b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/main.cpp new file mode 100644 index 0000000..b6ea583 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/main.cpp @@ -0,0 +1,89 @@ +#include "mhd1d.hpp" + +#include +#include +#include +#include + +constexpr int Nx = 100; +constexpr double Gamma = 2.0; +constexpr double Bx = 0.75; +constexpr mhd1d::StateVector LeftPrimitive{ + 1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0, +}; +constexpr mhd1d::StateVector RightPrimitive{ + 0.125, 0.0, 0.0, 0.0, 0.1, -1.0, 0.0, +}; + +int parse_nx(int argc, char** argv) +{ + if (argc <= 1) { + return Nx; + } + + char* end = nullptr; + const long parsed = std::strtol(argv[1], &end, 10); + if (end == argv[1] || *end != '\0' || parsed <= 0) { + throw std::runtime_error("usage: solver [nx]"); + } + return static_cast(parsed); +} + +mhd1d::SolverWorkspace initialize(int nx, double gamma, double bx, + const mhd1d::StateVector& left_state, + const mhd1d::StateVector& right_state) +{ + mhd1d::SolverWorkspace workspace(nx, gamma, bx); + + for (int ix = workspace.Lbx; ix <= workspace.Ubx; ++ix) { + const mhd1d::StateVector& state = (workspace.x(ix) < 0.5) ? left_state : right_state; + for (int component = 0; component < mhd1d::N_Component; ++component) { + workspace.up(ix, component) = state[component]; + workspace.uc(ix, component) = 0.0; + } + } + + for (int ix = 0; ix < workspace.Lbx; ++ix) { + for (int component = 0; component < mhd1d::N_Component; ++component) { + workspace.up(ix, component) = workspace.up(workspace.Lbx, component); + workspace.uc(ix, component) = workspace.uc(workspace.Lbx, component); + } + } + for (int ix = workspace.Ubx + 1; ix < workspace.Nx + 2 * mhd1d::N_margin; ++ix) { + for (int component = 0; component < mhd1d::N_Component; ++component) { + workspace.up(ix, component) = workspace.up(workspace.Ubx, component); + workspace.uc(ix, component) = workspace.uc(workspace.Ubx, component); + } + } + + (void)gamma; + (void)bx; + + return workspace; +} + +void write_csv(const mhd1d::SolverWorkspace& workspace, std::ostream& os) +{ + os << "x,rho,u,v,w,p,by,bz\n"; + os << std::setprecision(17); + for (int ix = workspace.Lbx; ix <= workspace.Ubx; ++ix) { + os << workspace.x(ix) << ',' << workspace.up(ix, 0) << ',' << workspace.up(ix, 1) << ',' + << workspace.up(ix, 2) << ',' << workspace.up(ix, 3) << ',' << workspace.up(ix, 4) << ',' + << workspace.up(ix, 5) << ',' << workspace.up(ix, 6) << '\n'; + } +} + +int main(int argc, char** argv) +{ + const int nx = parse_nx(argc, argv); + const double delt = 5.0e-4; + const double tmax = 0.1; + + auto workspace = initialize(nx, Gamma, Bx, LeftPrimitive, RightPrimitive); + + mhd1d::evolve_ssp_rk3(workspace, delt, tmax); + + write_csv(workspace, std::cout); + + return 0; +} diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/mhd1d.cpp b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/mhd1d.cpp new file mode 100644 index 0000000..7e34980 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/mhd1d.cpp @@ -0,0 +1,14 @@ +#include "mhd1d.hpp" + +namespace mhd1d +{ + +void evolve_ssp_rk3(SolverWorkspace& workspace, double dt, double t_final) +{ + (void)workspace; + (void)dt; + (void)t_final; + // TODO(student): implement the complete 1D MHD solver time evolution. +} + +} // namespace mhd1d diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/mhd1d.hpp b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/mhd1d.hpp new file mode 100644 index 0000000..ab97869 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/src/mhd1d.hpp @@ -0,0 +1,84 @@ +#pragma once + +#include +#include + +#include + +namespace mhd1d +{ + +namespace stdex = std::experimental; + +constexpr int N_Component = 7; +constexpr int N_margin = 1; + +using StateVector = std::array; +using ArrayView1D = stdex::mdspan, stdex::layout_right>; +using ArrayView2D = stdex::mdspan, stdex::layout_right>; + +struct SolverWorkspace { + explicit SolverWorkspace(int nx, double gamma, double bx) + : Nx(nx), Lbx(N_margin), Ubx(N_margin + nx - 1), dx(1.0 / static_cast(nx)), + gamma(gamma), bx(bx), storage(Nx + 2 * N_margin, N_Component) + { + init_views(Nx + 2 * N_margin, N_Component); + + for (int ix = Lbx; ix <= Ubx; ++ix) { + x(ix) = (static_cast(ix - Lbx) + 0.5) * dx; + } + } + + int Nx; + int Lbx; + int Ubx; + double dx; + double gamma; + double bx; + + ArrayView1D x; + ArrayView2D uc; + ArrayView2D up; + ArrayView2D up_l; + ArrayView2D up_r; + ArrayView2D rhs; + ArrayView2D prev; + ArrayView2D flux; + +private: + void init_views(int n_grid, int n_component) + { + x = ArrayView1D(storage.x.data(), n_grid); + uc = ArrayView2D(storage.uc.data(), n_grid, n_component); + up = ArrayView2D(storage.up.data(), n_grid, n_component); + up_l = ArrayView2D(storage.up_l.data(), n_grid, n_component); + up_r = ArrayView2D(storage.up_r.data(), n_grid, n_component); + rhs = ArrayView2D(storage.rhs.data(), n_grid, n_component); + prev = ArrayView2D(storage.prev.data(), n_grid, n_component); + flux = ArrayView2D(storage.flux.data(), n_grid, n_component); + } + + struct Storage { + explicit Storage(int n_grid, int n_component) + : x(n_grid), uc(n_grid * n_component), up(n_grid * n_component), up_l(n_grid * n_component), + up_r(n_grid * n_component), rhs(n_grid * n_component), prev(n_grid * n_component), + flux(n_grid * n_component) + { + } + + std::vector x; + std::vector uc; + std::vector up; + std::vector up_l; + std::vector up_r; + std::vector rhs; + std::vector prev; + std::vector flux; + }; + + Storage storage; +}; + +void evolve_ssp_rk3(SolverWorkspace& workspace, double dt, double t_final); + +} // namespace mhd1d diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/tests/cpp/test_public.cpp b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/tests/cpp/test_public.cpp new file mode 100644 index 0000000..c4eca2d --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/tests/cpp/test_public.cpp @@ -0,0 +1,27 @@ +#include + +#include + +#include "mhd1d.hpp" + +TEST_CASE("set_boundary duplicates edge states on both sides", "[mhd1d][boundary]") +{ + std::vector padded_buffer(4 * mhd1d::N_Component, 0.0); + const mhd1d::ArrayView2D padded(padded_buffer.data(), 4, mhd1d::N_Component); + + padded(1, 0) = 1.2; + padded(2, 0) = 1.4; + + const int lbx = 1; + const int ubx = 2; + + for (int ix = 0; ix < lbx; ++ix) { + padded(ix, 0) = padded(lbx, 0); + } + for (int ix = ubx + 1; ix < 4; ++ix) { + padded(ix, 0) = padded(ubx, 0); + } + + REQUIRE(padded(0, 0) == padded(1, 0)); + REQUIRE(padded(3, 0) == padded(2, 0)); +} diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/tests/data/brio_wu_golden.csv b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/tests/data/brio_wu_golden.csv new file mode 100644 index 0000000..e72dfa5 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/tests/data/brio_wu_golden.csv @@ -0,0 +1,101 @@ +x,rho,u,v,w,p,by,bz +0.0050000000000000001,1,0,0,0,1,1,0 +0.014999999999999999,1,0,0,0,1,1,0 +0.025000000000000001,1,0,0,0,1,1,0 +0.035000000000000003,1,0,0,0,1,1,0 +0.044999999999999998,1,0,0,0,1,1,0 +0.055,1,0,0,0,1,1,0 +0.065000000000000002,1,0,0,0,1,1,0 +0.074999999999999997,1,0,0,0,1,1,0 +0.085000000000000006,1,0,0,0,1,1,0 +0.095000000000000001,1,0,0,0,1,1,0 +0.105,1,0,0,0,1,1,0 +0.115,1,0,0,0,1,1,0 +0.125,1,0,0,0,1,1,0 +0.13500000000000001,1,0,0,0,1,1,0 +0.14499999999999999,1,0,0,0,1,1,0 +0.155,1,0,0,0,1,1,0 +0.16500000000000001,1,0,0,0,1,1,0 +0.17500000000000002,1,0,0,0,1,1,0 +0.185,1,0,0,0,1,1,0 +0.19500000000000001,1,0,0,0,1,1,0 +0.20500000000000002,1,0,0,0,1,1,0 +0.215,1,0,0,0,1,1,0 +0.22500000000000001,1,0,0,0,1,1,0 +0.23500000000000001,1,0,0,0,1,1,0 +0.245,1,0,0,0,1,1,0 +0.255,1,0,0,0,1,1,0 +0.26500000000000001,1,0,0,0,1,1,0 +0.27500000000000002,1,0,0,0,1,1,0 +0.28500000000000003,1,0,0,0,1,1,0 +0.29499999999999998,1,0,0,0,1,1,0 +0.30499999999999999,1,0,0,0,1,1,0 +0.315,0.99258917532833379,0.013274302361789324,-0.0037458540886650027,0,0.98532404108866567,0.99104740522973966,0 +0.32500000000000001,0.96954607581582897,0.054883618137207428,-0.015827325460063591,0,0.94023568077652442,0.96299689161606516,0 +0.33500000000000002,0.94025803324835733,0.10865479486043322,-0.031944255536031786,0,0.88435499782767446,0.92711462956015256,0 +0.34500000000000003,0.90897792466167582,0.16698841420781432,-0.050131763210425707,0,0.82652179164341721,0.88849866375884556,0 +0.35499999999999998,0.87725610844274526,0.2271603117432259,-0.069719685021130948,0,0.76985377322754589,0.84898878620884088,0 +0.36499999999999999,0.84566618167244145,0.28816698962982518,-0.090519447417695451,0,0.71541656913426399,0.80923410080538904,0 +0.375,0.81450435463148285,0.3494865276909434,-0.11251682956314066,0,0.66367321125402712,0.76950238387022263,0 +0.38500000000000001,0.78404538979850258,0.41062539085171551,-0.13575856770700617,0,0.61497905337973746,0.72997249890053162,0 +0.39500000000000002,0.75467120186953685,0.47083530429296916,-0.1602094195276916,0,0.56978778527926832,0.69094510148101218,0 +0.40500000000000003,0.72701668590610768,0.52870556766657628,-0.18539026277804579,0,0.52883750835390897,0.65319266964541056,0 +0.41500000000000004,0.70216500224646816,0.58153016354527853,-0.20973136326536138,0,0.4933637549672012,0.61876152226558601,0 +0.42499999999999999,0.68193441338380167,0.62487275671749021,-0.22845986767252549,0,0.46540195161446823,0.59219206428336624,0 +0.435,0.6666986977176832,0.65196598708089948,-0.23211489711409145,0,0.44490074049787121,0.58021922926273972,0 +0.44500000000000001,0.66015529724292976,0.66607074009479628,-0.22555655390938342,0,0.43615912199888263,0.58524191570556205,0 +0.45500000000000002,0.65997980757445474,0.66324409170331222,-0.23280516118974198,0,0.43727412310940939,0.5796363680769816,0 +0.46500000000000002,0.73349815599883794,0.60086667962898554,-0.62171855754180494,0,0.57886068825314341,0.27652712488086967,0 +0.47500000000000003,0.79237811897791921,0.46454032083596353,-1.2161018009267417,0,0.6842784144671038,-0.24404258805909654,0 +0.48499999999999999,0.74259318372203176,0.53084767864673321,-1.4743018340823315,0,0.60002741951352012,-0.44224666294204562,0 +0.495,0.69493377959694547,0.61399418172955245,-1.5431753009247675,0,0.54413704626033854,-0.54049642839751277,0 +0.505,0.70278259276020427,0.61662067481286442,-1.5872354257997463,0,0.52068648102692328,-0.56052646550215834,0 +0.51500000000000001,0.70390914156072559,0.6077114699355689,-1.6001701751263477,0,0.51193399003676154,-0.54226823182770345,0 +0.52500000000000002,0.69922181544266515,0.59340751947871018,-1.6046948629200928,0,0.51300004963251122,-0.52185059484429897,0 +0.53500000000000003,0.67807490745420074,0.5878535845789149,-1.6006935764888095,0,0.51466016361553968,-0.51621274881535362,0 +0.54500000000000004,0.60556902444692728,0.59045915738825405,-1.5906721688577199,0,0.51204702934936586,-0.52148853682254903,0 +0.55500000000000005,0.48266830107336006,0.59697028797374452,-1.5761776779942758,0,0.51057551463380291,-0.53191252428398161,0 +0.56500000000000006,0.34834934855591537,0.59188322457049136,-1.5687669989374082,0,0.5151520944223198,-0.54463835756014178,0 +0.57500000000000007,0.2451691037246034,0.59463166677984525,-1.5635053912378623,0,0.51994523862074438,-0.55151343157376009,0 +0.58499999999999996,0.22308154337733685,0.61421260394352284,-1.5558407025582697,0,0.51450769544894648,-0.54648231691723681,0 +0.59499999999999997,0.22390587217315483,0.62861601814674595,-1.555379908924208,0,0.51042675178345809,-0.53931383334019445,0 +0.60499999999999998,0.22786331744185459,0.61304302077924422,-1.5697800797404926,0,0.51789637516013132,-0.53999599169525736,0 +0.61499999999999999,0.23452030041149688,0.58777367268577241,-1.5968680437122895,0,0.53659755807779264,-0.54285121456661267,0 +0.625,0.23506241740122699,0.58441088978544375,-1.6180990627429224,0,0.5312595858066933,-0.53847776600605912,0 +0.63500000000000001,0.22912802575710367,0.64096361665071178,-1.5587196384023227,0,0.50180683386498415,-0.51954980419936259,0 +0.64500000000000002,0.20820170172948924,0.4429088168703938,-1.3138217298639516,0,0.4107172773448764,-0.63242114973177121,0 +0.65500000000000003,0.14280420150766387,0.0013572736674349979,-0.58334822550752841,0,0.15933600494306066,-0.82699909098420499,0 +0.66500000000000004,0.11735885521106922,-0.23309434559311201,-0.17605497312958984,0,0.088246122408452865,-0.90091689844612777,0 +0.67500000000000004,0.11691966169225833,-0.24239401825294279,-0.16901850205441424,0,0.087525697429624683,-0.90138999879178838,0 +0.68500000000000005,0.1168975018076187,-0.24275168390313129,-0.16913431059427403,0,0.087491282208618681,-0.90130727955920587,0 +0.69500000000000006,0.11688490837927806,-0.24293963466430563,-0.16922798396131322,0,0.087471166714507165,-0.90122308503037618,0 +0.70499999999999996,0.11687882431088358,-0.24305091000018803,-0.16930677036965813,0,0.087461000889423768,-0.90113549166198792,0 +0.71499999999999997,0.11688210335016401,-0.24310898990635241,-0.16941137635009459,0,0.087464485048653229,-0.90109713084582599,0 +0.72499999999999998,0.11688746427801196,-0.24311790246317952,-0.16948241337481523,0,0.087471346584835463,-0.90112657201926949,0 +0.73499999999999999,0.11689293842702393,-0.24301542095697384,-0.16941123056591925,0,0.087478382198893812,-0.90117181152414561,0 +0.745,0.11689754429437521,-0.24279548243457205,-0.16922004702154958,0,0.087484264982216398,-0.90123329716500089,0 +0.755,0.11690472407462463,-0.24244042822845885,-0.16894496953161006,0,0.087493985877884928,-0.90135933427866255,0 +0.76500000000000001,0.1169351752430936,-0.24164401100336819,-0.16834345193351358,0,0.087538755474644847,-0.90173628521210369,0 +0.77500000000000002,0.11702346553875453,-0.2389118539454739,-0.16634253632348067,0,0.087670236260091849,-0.90283181009609326,0 +0.78500000000000003,0.11723705803609173,-0.23241448265557926,-0.16159143199373074,0,0.087990348375741601,-0.90548092874079178,0 +0.79500000000000004,0.1176285175883393,-0.22051415125476836,-0.15292012220738538,0,0.088578608675684678,-0.91032952865801275,0 +0.80500000000000005,0.11822301685274909,-0.20247335890203855,-0.13985941063243096,0,0.089475764221054743,-0.91767747139414524,0 +0.81500000000000006,0.11901085894402198,-0.178638590244622,-0.12276084795550267,0,0.090670936816770853,-0.92738054886244792,0 +0.82500000000000007,0.11995507930331456,-0.15017685276699991,-0.10257309790892337,0,0.092113203419944889,-0.93896119897869035,0 +0.83499999999999996,0.12100582881114735,-0.1186354504616307,-0.080489249938489674,0,0.093731854039529328,-0.95178772370322251,0 +0.84499999999999997,0.12210948952666369,-0.085657286840511535,-0.057715619723010918,0,0.095446790184988584,-0.96519228251210643,0 +0.85499999999999998,0.12320528139333707,-0.053056324701219029,-0.035509818552657238,0,0.097164687932375315,-0.97843939905855704,0 +0.86499999999999999,0.12419939842673788,-0.023611385320102929,-0.015706138420103915,0,0.098733750673405729,-0.99040292502414518,0 +0.875,0.1248831018288558,-0.0034412091825772586,-0.0022747475485175887,0,0.099815318398622899,-0.99860261559100016,0 +0.88500000000000001,0.12500000000700623,2.0693358136986429e-10,1.3690028369419877e-10,0,0.10000000001121101,-1.0000000000838631,0 +0.89500000000000002,0.12500000000690034,2.0295860176552916e-10,1.3426844820417816e-10,0,0.10000000001104148,-1.0000000000825946,0 +0.90500000000000003,0.12500000000616657,1.8175890421477852e-10,1.2024374828573813e-10,0,0.10000000000986675,-1.0000000000738092,0 +0.91500000000000004,0.12500000000447511,1.3188042879502353e-10,8.7246950401911659e-11,0,0.10000000000716036,-1.0000000000535656,0 +0.92500000000000004,0.1250000000023313,6.8702813807084164e-11,4.5451087337373151e-11,0,0.10000000000372944,-1.0000000000279046,0 +0.93500000000000005,0.12500000000053466,1.575462353161844e-11,1.0423025405688057e-11,0,0.1000000000008554,-1.0000000000063987,0 +0.94500000000000006,0.125,0,0,0,0.099999999999999867,-1,0 +0.95500000000000007,0.125,0,0,0,0.099999999999999867,-1,0 +0.96499999999999997,0.125,0,0,0,0.099999999999999867,-1,0 +0.97499999999999998,0.125,0,0,0,0.099999999999999867,-1,0 +0.98499999999999999,0.125,0,0,0,0.099999999999999867,-1,0 +0.995,0.125,0,0,0,0.099999999999999867,-1,0 diff --git a/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/tests/test_public.py b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/tests/test_public.py new file mode 100644 index 0000000..24e62c7 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-full1d-01/workspace/tests/test_public.py @@ -0,0 +1,81 @@ +import csv +import os +import subprocess +from pathlib import Path + +PUBLIC_TEST_TARGET = "cpp_full_solver1d_public_tests" +GOLDEN_CSV_PATH = Path(__file__).resolve().parents[1] / "tests/data/brio_wu_golden.csv" +TOLERANCE = 1.0e-12 +WORKSPACE_ROOT = Path(__file__).resolve().parents[1] + + +def assert_csv_rows_close( + output_rows: list[list[str]], + reference_rows: list[list[str]], + *, + tolerance: float = TOLERANCE, + expected_header: list[str] | None = None, +) -> None: + if expected_header is not None: + assert output_rows[0] == expected_header + output_rows = output_rows[1:] + + assert len(output_rows) == len(reference_rows) + + for output_row, reference_row in zip(output_rows, reference_rows): + assert len(output_row) == len(reference_row) + for output_value, reference_value in zip(output_row, reference_row): + assert abs(float(output_value) - float(reference_value)) <= tolerance + + +def _build_public_tests() -> Path: + build_dir = "build" + subprocess.run(["cmake", "-S", ".", "-B", build_dir], check=True, cwd=WORKSPACE_ROOT) + subprocess.run( + [ + "cmake", + "--build", + build_dir, + "--target", + "cpp_full_solver1d", + PUBLIC_TEST_TARGET, + ], + check=True, + cwd=WORKSPACE_ROOT, + ) + + binary_name = f"{PUBLIC_TEST_TARGET}.exe" if os.name == "nt" else PUBLIC_TEST_TARGET + executable_path = WORKSPACE_ROOT / build_dir / "tests" / binary_name + assert executable_path.exists() + return executable_path + + +def test_public_catch2_target_builds() -> None: + _build_public_tests() + + +def test_public_brio_wu_cli_matches_reference_grid() -> None: + _build_public_tests() + + solver_name = "cpp_full_solver1d.exe" if os.name == "nt" else "cpp_full_solver1d" + solver_path = WORKSPACE_ROOT / "build/bin" / solver_name + assert solver_path.exists() + + completed = subprocess.run( + [str(solver_path)], + check=True, + capture_output=True, + text=True, + ) + + rows = list(csv.reader(completed.stdout.splitlines())) + golden_rows = list(csv.reader(GOLDEN_CSV_PATH.read_text(encoding="utf-8").splitlines())) + + assert golden_rows[0] == ["x", "rho", "u", "v", "w", "p", "by", "bz"] + + assert_csv_rows_close( + rows, + golden_rows[1:], + tolerance=TOLERANCE, + expected_header=["x", "rho", "u", "v", "w", "p", "by", "bz"], + ) diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/eval/run.sh b/benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/run.sh similarity index 100% rename from benchmarks/magnetohydrodynamics/cpp-hlld/eval/run.sh rename to benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/run.sh diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/eval/tests/cpp/hlld_reference.hpp b/benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/tests/cpp/hlld_reference.hpp similarity index 83% rename from benchmarks/magnetohydrodynamics/cpp-hlld/eval/tests/cpp/hlld_reference.hpp rename to benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/tests/cpp/hlld_reference.hpp index 8206bf8..4e3d3dc 100644 --- a/benchmarks/magnetohydrodynamics/cpp-hlld/eval/tests/cpp/hlld_reference.hpp +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/tests/cpp/hlld_reference.hpp @@ -6,13 +6,16 @@ #include "../../../workspace/src/hlld.hpp" #endif +#include #include +using StateVector = std::array; + namespace hidden_reference { -inline FluxState hlld_flux_from_primitive(const PrimitiveState& left, const PrimitiveState& right, - double bx, double gamma) +inline StateVector hlld_flux_from_primitive(const StateVector& left, const StateVector& right, + double bx, double gamma) { constexpr double eps = 1.0e-40; @@ -66,20 +69,20 @@ inline FluxState hlld_flux_from_primitive(const PrimitiveState& left, const Prim const double sl = std::min(vxl, vxr) - std::max(cfl, cfr); const double sr = std::max(vxl, vxr) + std::max(cfl, cfr); - const FluxState fql{rxl, - rxl * vxl + ptl - bxsq, - rxl * vyl - bxs * byl, - rxl * vzl - bxs * bzl, - vxl * (eel + ptl - bxsq) - bxs * (vyl * byl + vzl * bzl), - byl * vxl - bxs * vyl, - bzl * vxl - bxs * vzl}; - const FluxState fqr{rxr, - rxr * vxr + ptr - bxsq, - rxr * vyr - bxs * byr, - rxr * vzr - bxs * bzr, - vxr * (eer + ptr - bxsq) - bxs * (vyr * byr + vzr * bzr), - byr * vxr - bxs * vyr, - bzr * vxr - bxs * vzr}; + const StateVector fql{rxl, + rxl * vxl + ptl - bxsq, + rxl * vyl - bxs * byl, + rxl * vzl - bxs * bzl, + vxl * (eel + ptl - bxsq) - bxs * (vyl * byl + vzl * bzl), + byl * vxl - bxs * vyl, + bzl * vxl - bxs * vzl}; + const StateVector fqr{rxr, + rxr * vxr + ptr - bxsq, + rxr * vyr - bxs * byr, + rxr * vzr - bxs * bzr, + vxr * (eer + ptr - bxsq) - bxs * (vyr * byr + vzr * bzr), + byr * vxr - bxs * vyr, + bzr * vxr - bxs * vzr}; const double sdl = sl - vxl; const double sdr = sr - vxr; @@ -197,7 +200,7 @@ inline FluxState hlld_flux_from_primitive(const PrimitiveState& left, const Prim const double temp_flux_l = mslst - msl; const double temp_flux_r = msrst - msr; - return FluxState{ + return StateVector{ (fql[0] - msl * rol - rolst * temp_flux_l + roldst * mslst) * maxs1 + (fqr[0] - msr * ror - rorst * temp_flux_r + rordst * msrst) * mins1, (fql[1] - msl * rxl - rxlst * temp_flux_l + rxldst * mslst) * maxs1 + @@ -215,24 +218,12 @@ inline FluxState hlld_flux_from_primitive(const PrimitiveState& left, const Prim }; } -inline FluxState hlld_flux_from_conservative(const ConservativeState& left, - const ConservativeState& right, double bx, - double gamma) +inline StateVector solver_flux_from_primitive(const StateVector& left, const StateVector& right, + double bx, double gamma) { - const auto to_primitive = [bx, gamma](const ConservativeState& state) { - const double rho = state[0]; - const double u = state[1] / rho; - const double v = state[2] / rho; - const double w = state[3] / rho; - const double by = state[5]; - const double bz = state[6]; - const double kinetic = 0.5 * rho * (u * u + v * v + w * w); - const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); - const double p = (gamma - 1.0) * (state[4] - kinetic - magnetic); - return PrimitiveState{rho, u, v, w, p, by, bz}; - }; - - return hlld_flux_from_primitive(to_primitive(left), to_primitive(right), bx, gamma); + StateVector flux{}; + ::hlld_flux_from_primitive(left.data(), right.data(), bx, gamma, flux.data()); + return flux; } } // namespace hidden_reference diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/eval/tests/cpp/test_hidden.cpp b/benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/tests/cpp/test_hidden.cpp similarity index 52% rename from benchmarks/magnetohydrodynamics/cpp-hlld/eval/tests/cpp/test_hidden.cpp rename to benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/tests/cpp/test_hidden.cpp index 4781c2a..35d8074 100644 --- a/benchmarks/magnetohydrodynamics/cpp-hlld/eval/tests/cpp/test_hidden.cpp +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/tests/cpp/test_hidden.cpp @@ -12,14 +12,17 @@ #include "hlld_reference.hpp" +#include #include +using StateVector = std::array; + namespace { constexpr double kTolerance = 1e-12; -ConservativeState primitive_to_conservative(const PrimitiveState& state, double bx, double gamma) +StateVector primitive_to_conservative(const StateVector& state, double bx, double gamma) { const double rho = state[0]; const double u = state[1]; @@ -33,12 +36,10 @@ ConservativeState primitive_to_conservative(const PrimitiveState& state, double const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); const double energy = p / (gamma - 1.0) + kinetic + magnetic; - return ConservativeState{ - rho, rho * u, rho * v, rho * w, energy, by, bz, - }; + return StateVector{rho, rho * u, rho * v, rho * w, energy, by, bz}; } -FluxState physical_flux_x(const ConservativeState& state, double bx, double gamma) +StateVector physical_flux_x(const StateVector& state, double bx, double gamma) { const double rho = state[0]; const double mx = state[1]; @@ -56,7 +57,7 @@ FluxState physical_flux_x(const ConservativeState& state, double bx, double gamm const double pressure = (gamma - 1.0) * (energy - kinetic - magnetic); const double total_pressure = pressure + magnetic; - return FluxState{ + return StateVector{ rho * u, rho * u * u + total_pressure - bx * bx, rho * v * u - bx * by, @@ -67,24 +68,32 @@ FluxState physical_flux_x(const ConservativeState& state, double bx, double gamm }; } -void require_close(const FluxState& actual, const FluxState& expected) +void require_close(const StateVector& actual, const StateVector& expected) { for (std::size_t i = 0; i < actual.size(); ++i) { REQUIRE(std::abs(actual[i] - expected[i]) <= kTolerance); } } +StateVector solver_flux_from_primitive(const StateVector& left, const StateVector& right, double bx, + double gamma) +{ + StateVector flux{}; + hlld_flux_from_primitive(left.data(), right.data(), bx, gamma, flux.data()); + return flux; +} + } // namespace -TEST_CASE("equal conservative states reduce to the physical flux") +TEST_CASE("equal primitive states reduce to the physical flux") { - const double bx = 0.35; - const double gamma = 1.4; - const PrimitiveState primitive{0.9, -0.45, 0.2, 0.15, 0.8, -0.3, 0.55}; - const ConservativeState state = primitive_to_conservative(primitive, bx, gamma); + const double bx = 0.35; + const double gamma = 1.4; + const StateVector state{0.9, -0.45, 0.2, 0.15, 0.8, -0.3, 0.55}; - const FluxState actual = hlld_flux_from_conservative(state, state, bx, gamma); - const FluxState expected = physical_flux_x(state, bx, gamma); + const StateVector actual = solver_flux_from_primitive(state, state, bx, gamma); + const StateVector expected = + physical_flux_x(primitive_to_conservative(state, bx, gamma), bx, gamma); require_close(actual, expected); } @@ -94,10 +103,10 @@ TEST_CASE("nontrivial primitive solve returns finite values") const double bx = -0.65; const double gamma = 5.0 / 3.0; - const PrimitiveState left{1.08, 0.45, -0.12, 0.08, 0.95, 0.4, -0.3}; - const PrimitiveState right{0.72, -0.25, 0.16, -0.05, 0.58, -0.2, 0.35}; + const StateVector left{1.08, 0.45, -0.12, 0.08, 0.95, 0.4, -0.3}; + const StateVector right{0.72, -0.25, 0.16, -0.05, 0.58, -0.2, 0.35}; - const FluxState flux = hlld_flux_from_primitive(left, right, bx, gamma); + const StateVector flux = solver_flux_from_primitive(left, right, bx, gamma); for (double value : flux) { REQUIRE(std::isfinite(value)); @@ -109,11 +118,11 @@ TEST_CASE("hidden reference flux case 1 matches reference implementation") const double bx = -0.65; const double gamma = 5.0 / 3.0; - const PrimitiveState left{1.08, 0.45, -0.12, 0.08, 0.95, 0.4, -0.3}; - const PrimitiveState right{0.72, -0.25, 0.16, -0.05, 0.58, -0.2, 0.35}; - const FluxState expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); + const StateVector left{1.08, 0.45, -0.12, 0.08, 0.95, 0.4, -0.3}; + const StateVector right{0.72, -0.25, 0.16, -0.05, 0.58, -0.2, 0.35}; + const StateVector expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); require_close(actual, expected); } @@ -122,11 +131,11 @@ TEST_CASE("hidden reference flux case 2 matches reference implementation") const double bx = 0.35; const double gamma = 1.4; - const PrimitiveState left{0.9, -0.45, 0.2, 0.15, 0.8, -0.3, 0.55}; - const PrimitiveState right{1.15, 0.18, -0.12, -0.08, 1.05, 0.22, -0.4}; - const FluxState expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); + const StateVector left{0.9, -0.45, 0.2, 0.15, 0.8, -0.3, 0.55}; + const StateVector right{1.15, 0.18, -0.12, -0.08, 1.05, 0.22, -0.4}; + const StateVector expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); require_close(actual, expected); } @@ -135,11 +144,11 @@ TEST_CASE("small Bx near-degenerate reference case matches reference implementat const double bx = 1.0e-6; const double gamma = 1.4; - const PrimitiveState left{1.0, 0.4, 0.2, -0.1, 1.0, 0.5, -0.4}; - const PrimitiveState right{0.85, -0.3, -0.15, 0.25, 0.8, -0.35, 0.45}; - const FluxState expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); + const StateVector left{1.0, 0.4, 0.2, -0.1, 1.0, 0.5, -0.4}; + const StateVector right{0.85, -0.3, -0.15, 0.25, 0.8, -0.35, 0.45}; + const StateVector expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); require_close(actual, expected); } @@ -148,10 +157,10 @@ TEST_CASE("second Bx equals zero hydro case matches reference implementation") const double bx = 0.0; const double gamma = 1.4; - const PrimitiveState left{0.4, -1.1, 0.0, 0.0, 0.4, 0.0, 0.0}; - const PrimitiveState right{1.2, -0.2, 0.0, 0.0, 1.3, 0.0, 0.0}; - const FluxState expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); + const StateVector left{0.4, -1.1, 0.0, 0.0, 0.4, 0.0, 0.0}; + const StateVector right{1.2, -0.2, 0.0, 0.0, 1.3, 0.0, 0.0}; + const StateVector expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); require_close(actual, expected); } diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/eval/tests/test_hidden.py b/benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/tests/test_hidden.py similarity index 100% rename from benchmarks/magnetohydrodynamics/cpp-hlld/eval/tests/test_hidden.py rename to benchmarks/magnetohydrodynamics/cpp-hlld-00/eval/tests/test_hidden.py diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-00/spec.md b/benchmarks/magnetohydrodynamics/cpp-hlld-00/spec.md new file mode 100644 index 0000000..e4e3621 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-00/spec.md @@ -0,0 +1,48 @@ +# cpp-hlld-00 + +Implement the HLLD approximate Riemann solver for 1D ideal MHD in C++. + +## Read first + +- `/work/basic_equations.md` +- `/work/hlld.md` + +## Task + +Edit `src/hlld.cpp` so that `hlld_flux_from_primitive(...)` is implemented correctly. + +The benchmark uses: + +- primitive-state ordering: `[rho, u, v, w, p, By, Bz]` +- flux ordering: `[F_rho, F_mx, F_my, F_mz, F_E, F_By, F_Bz]` +- the test suite includes `Bx = 0` hydro and magnetized degenerate cases +- the test suite also includes a small-`Bx` near-degenerate case, so handle + `Bx = 0`, small denominators in the starred-state formulas, and related + square-root/discriminant edge cases carefully + +Do not change the public function signatures in `src/hlld.hpp`. + +## Implementation hints + +- This benchmark follows one specific HLLD implementation convention rather than + an arbitrary mathematically equivalent variant. +- Small starred-state denominator (`D_alpha`): if `|D_alpha|` is extremely + small, avoid dividing by it and fall back to unchanged transverse starred + values (`v* = v`, `w* = w`, `By* = By`, `Bz* = Bz`). +- Small `Bx`: when `Bx = 0`, rotational waves collapse and double-star states + are unnecessary; do not use double-star states in flux selection in that + case. +- For this benchmark, a merely small nonzero `|Bx|` is still nondegenerate + unless another guarded quantity (such as `D_alpha`) becomes numerically + singular. +- Assume all benchmark inputs are admissible physical states. + +## Standards + +- C++17 + +## Local dev + +```bash +pytest -q +``` diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/task.toml b/benchmarks/magnetohydrodynamics/cpp-hlld-00/task.toml similarity index 92% rename from benchmarks/magnetohydrodynamics/cpp-hlld/task.toml rename to benchmarks/magnetohydrodynamics/cpp-hlld-00/task.toml index f58f5d8..2d7c0d1 100644 --- a/benchmarks/magnetohydrodynamics/cpp-hlld/task.toml +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-00/task.toml @@ -1,4 +1,4 @@ -id = "cpp-hlld" +id = "cpp-hlld-00" suite = "magnetohydrodynamics" language = "cpp" time_limit_sec = 600 diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/CMakeLists.txt b/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/CMakeLists.txt similarity index 100% rename from benchmarks/magnetohydrodynamics/cpp-hlld/workspace/CMakeLists.txt rename to benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/CMakeLists.txt diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/pyproject.toml b/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/pyproject.toml similarity index 100% rename from benchmarks/magnetohydrodynamics/cpp-hlld/workspace/pyproject.toml rename to benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/pyproject.toml diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/src/hlld.cpp b/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/src/hlld.cpp new file mode 100644 index 0000000..b200baf --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/src/hlld.cpp @@ -0,0 +1,14 @@ +#include "hlld.hpp" + +#include + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux) +{ + (void)left; + (void)right; + (void)bx; + (void)gamma; + + std::fill(flux, flux + 7, 0.0); +} diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/src/hlld.hpp b/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/src/hlld.hpp new file mode 100644 index 0000000..ae80ebd --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/src/hlld.hpp @@ -0,0 +1,4 @@ +#pragma once + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux); diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/tests/cpp/test_public.cpp b/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/tests/cpp/test_public.cpp new file mode 100644 index 0000000..0c73c44 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/tests/cpp/test_public.cpp @@ -0,0 +1,300 @@ +#if __has_include("hlld.hpp") +#include "hlld.hpp" +#else +#include "../../src/hlld.hpp" +#endif + +#if __has_include() +#include +#else +#include "/usr/local/include/catch2/catch_test_macros.hpp" +#endif + +#include + +#include + +using StateVector = std::array; + +namespace +{ + +constexpr double kTolerance = 1e-12; +constexpr double kPi = 3.14159265358979323846; + +StateVector primitive_to_conservative(const StateVector& state, double bx, double gamma) +{ + const double rho = state[0]; + const double u = state[1]; + const double v = state[2]; + const double w = state[3]; + const double p = state[4]; + const double by = state[5]; + const double bz = state[6]; + + const double kinetic = 0.5 * rho * (u * u + v * v + w * w); + const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); + const double energy = p / (gamma - 1.0) + kinetic + magnetic; + + return StateVector{ + rho, rho * u, rho * v, rho * w, energy, by, bz, + }; +} + +StateVector physical_flux_x(const StateVector& state, double bx, double gamma) +{ + const double rho = state[0]; + const double mx = state[1]; + const double my = state[2]; + const double mz = state[3]; + const double energy = state[4]; + const double by = state[5]; + const double bz = state[6]; + + const double u = mx / rho; + const double v = my / rho; + const double w = mz / rho; + const double kinetic = 0.5 * rho * (u * u + v * v + w * w); + const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); + const double pressure = (gamma - 1.0) * (energy - kinetic - magnetic); + const double total_pressure = pressure + magnetic; + + return StateVector{ + rho * u, + rho * u * u + total_pressure - bx * bx, + rho * v * u - bx * by, + rho * w * u - bx * bz, + (energy + total_pressure) * u - bx * (u * bx + v * by + w * bz), + by * u - bx * v, + bz * u - bx * w, + }; +} + +StateVector solver_flux_from_primitive(const StateVector& left, const StateVector& right, double bx, + double gamma) +{ + StateVector flux{}; + hlld_flux_from_primitive(left.data(), right.data(), bx, gamma, flux.data()); + return flux; +} + +void require_close(const StateVector& actual, const StateVector& expected) +{ + for (std::size_t i = 0; i < actual.size(); ++i) { + REQUIRE(std::abs(actual[i] - expected[i]) <= kTolerance); + } +} + +} // namespace + +TEST_CASE("equal primitive states reduce to the physical flux") +{ + const double bx = 0.75; + const double gamma = 1.4; + const StateVector state{1.1, 0.2, -0.3, 0.4, 0.9, 0.5, -0.6}; + const StateVector conservative = primitive_to_conservative(state, bx, gamma); + + const StateVector actual = solver_flux_from_primitive(state, state, bx, gamma); + const StateVector expected = physical_flux_x(conservative, bx, gamma); + + require_close(actual, expected); +} + +TEST_CASE("right-going contact discontinuity is resolved exactly") +{ + const double bx = 0.8; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.3, 0.2, -0.15, 1.0, 0.6, -0.3}; + const StateVector right{0.7, 0.3, 0.2, -0.15, 1.0, 0.6, -0.3}; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, physical_flux_x(primitive_to_conservative(left, bx, gamma), bx, gamma)); +} + +TEST_CASE("left-going contact discontinuity is resolved exactly") +{ + const double bx = 0.8; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, -0.25, 0.2, -0.15, 1.0, 0.6, -0.3}; + const StateVector right{0.7, -0.25, 0.2, -0.15, 1.0, 0.6, -0.3}; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, physical_flux_x(primitive_to_conservative(right, bx, gamma), bx, gamma)); +} + +TEST_CASE("right-going rotational discontinuity is resolved exactly") +{ + const double bx = 1.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.2, 0.1, -0.2, 1.0, 1.0, 0.0}; + const StateVector right{1.0, 0.2, 0.5, -1.0, 1.0, 0.6, 0.8}; + const StateVector expected{ + 0.2, 1.04, -0.98, -0.04, 0.609, 0.1, 0.2, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("left-going rotational discontinuity is resolved exactly") +{ + const double bx = 1.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.2, 0.1, -0.2, 1.0, 1.0, 0.0}; + const StateVector right{1.0, 0.2, -0.3, 0.6, 1.0, 0.6, 0.8}; + const StateVector expected{ + 0.2, 1.04, -0.66, -0.68, 0.449, 0.42, -0.44, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Bx equals zero hydro case matches reference flux") +{ + const double bx = 0.0; + const double gamma = 1.4; + + const StateVector left{1.0, 0.75, 0.0, 0.0, 1.0, 0.0, 0.0}; + const StateVector right{0.125, 0.0, 0.0, 0.0, 0.1, 0.0, 0.0}; + const StateVector expected{ + 0.92274146439449267, 1.3581095429585437, 0.0, 0.0, 3.1282919538345322, 0.0, 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Bx equals zero magnetized case matches reference flux") +{ + const double bx = 0.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.6, 0.1, -0.2, 1.0, 0.7, -0.5}; + const StateVector right{0.7, -0.3, -0.15, 0.25, 0.5, -0.2, 0.4}; + const StateVector expected{ + 0.44815524807196727, 2.011116795062418, 0.044815524807196722, -0.089631049614393443, + 1.6980640537315086, 0.31370867365037713, -0.22407762403598386, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("small Bx near-degenerate case matches reference flux") +{ + const double bx = 1.0e-6; + const double gamma = 1.4; + + const StateVector left{1.0, 0.4, 0.2, -0.1, 1.0, 0.5, -0.4}; + const StateVector right{0.85, -0.3, -0.15, 0.25, 0.8, -0.35, 0.45}; + const StateVector expected{ + 0.16298732855830989, 1.7549717390289374, 0.032596899426565185, -0.016298279827753587, + 0.72345786285986069, 0.081493464279122407, -0.065194831423298072, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Ryu and Jones shock tube matches reference flux") +{ + const double bx = 4.0 / std::sqrt(4.0 * kPi); + const double gamma = 5.0 / 3.0; + + const StateVector left{ + 1.08, 1.2, 0.01, 0.5, 0.95, 3.6 / std::sqrt(4.0 * kPi), 2.0 / std::sqrt(4.0 * kPi), + }; + const StateVector right{ + 1.0, 0.0, 0.0, 0.0, 1.0, 4.0 / std::sqrt(4.0 * kPi), 2.0 / std::sqrt(4.0 * kPi), + }; + const StateVector expected{ + 0.79485593966715773, 3.5458209484697329, -1.3572358551169827, -0.22185101509215432, + 3.9950643754664625, 0.67495208799031015, -0.062307582042232856, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Brio and Wu shock tube matches reference flux") +{ + const double bx = 0.75; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0}; + const StateVector right{0.125, 0.0, 0.0, 0.0, 0.1, -1.0, 0.0}; + const StateVector expected{ + 0.2063330447744266, + 0.4638678509599396, + 0.064186763013841408, + 0.0, + 0.16136546437466026, + 1.010233243594872, + 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Falle switch-off shock matches reference flux") +{ + const double bx = 1.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.368, 0.269, 1.0, 0.0, 1.769, 0.0, 0.0}; + const StateVector right{1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0}; + const StateVector expected{ + 0.29721990694355238, + 1.4932992607654056, + 0.2893229270591654, + 0.0, + 1.1427267633652525, + -1.0066552479847843, + 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Falle switch-off rarefaction matches reference flux") +{ + const double bx = 1.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.0, 0.0, 0.0, 2.0, 0.0, 0.0}; + const StateVector right{0.2, 1.186, 2.967, 0.0, 0.1368, 1.6405, 0.0}; + const StateVector expected{ + 0.27717801577960577, + 0.28228035303750848, + -1.3364302412732558, + 0.0, + -1.5599793330037519, + -1.3806947854633354, + 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("super-fast expansion matches reference flux") +{ + const double bx = 0.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, -3.0, 0.0, 0.0, 0.45, 0.5, 0.0}; + const StateVector right{1.0, 3.0, 0.0, 0.0, 0.45, 0.5, 0.0}; + const StateVector expected{ + 0.0, -2.425, 0.0, 0.0, 0.0, 0.0, 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/tests/test_public.py b/benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/tests/test_public.py similarity index 100% rename from benchmarks/magnetohydrodynamics/cpp-hlld/workspace/tests/test_public.py rename to benchmarks/magnetohydrodynamics/cpp-hlld-00/workspace/tests/test_public.py diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/eval/run.sh b/benchmarks/magnetohydrodynamics/cpp-hlld-01/eval/run.sh new file mode 100755 index 0000000..c95929f --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/eval/run.sh @@ -0,0 +1,27 @@ +#!/usr/bin/env bash +set -u -o pipefail + +cd /work + +status="passed" +score="1.0" + +python3 -m pytest -q /eval/tests +rc=$? +if [ "$rc" -ne 0 ]; then + status="failed" + score="0.0" +fi + +python3 - < +#include + +using StateVector = std::array; + +namespace hidden_reference +{ + +inline StateVector hlld_flux_from_primitive(const StateVector& left, const StateVector& right, + double bx, double gamma) +{ + constexpr double eps = 1.0e-40; + + const double rol = left[0]; + const double vxl = left[1]; + const double vyl = left[2]; + const double vzl = left[3]; + const double prl = left[4]; + const double byl = left[5]; + const double bzl = left[6]; + + const double ror = right[0]; + const double vxr = right[1]; + const double vyr = right[2]; + const double vzr = right[3]; + const double prr = right[4]; + const double byr = right[5]; + const double bzr = right[6]; + + const double igm = 1.0 / (gamma - 1.0); + const double bxs = bx; + const double bxsq = bxs * bxs; + + const double pbl = 0.5 * (bxsq + byl * byl + bzl * bzl); + const double pbr = 0.5 * (bxsq + byr * byr + bzr * bzr); + const double ptl = prl + pbl; + const double ptr = prr + pbr; + + const double rxl = rol * vxl; + const double ryl = rol * vyl; + const double rzl = rol * vzl; + const double rxr = ror * vxr; + const double ryr = ror * vyr; + const double rzr = ror * vzr; + + const double eel = prl * igm + 0.5 * (rxl * vxl + ryl * vyl + rzl * vzl) + pbl; + const double eer = prr * igm + 0.5 * (rxr * vxr + ryr * vyr + rzr * vzr) + pbr; + + const double gmpl = gamma * prl; + const double gmpr = gamma * prr; + const double gpbl = gmpl + 2.0 * pbl; + const double gpbr = gmpr + 2.0 * pbr; + + const double cfl = std::sqrt((gpbl + std::sqrt((gmpl - 2.0 * pbl) * (gmpl - 2.0 * pbl) + + 4.0 * gmpl * (byl * byl + bzl * bzl))) * + 0.5 / rol); + const double cfr = std::sqrt((gpbr + std::sqrt((gmpr - 2.0 * pbr) * (gmpr - 2.0 * pbr) + + 4.0 * gmpr * (byr * byr + bzr * bzr))) * + 0.5 / ror); + + const double sl = std::min(vxl, vxr) - std::max(cfl, cfr); + const double sr = std::max(vxl, vxr) + std::max(cfl, cfr); + + const StateVector fql{rxl, + rxl * vxl + ptl - bxsq, + rxl * vyl - bxs * byl, + rxl * vzl - bxs * bzl, + vxl * (eel + ptl - bxsq) - bxs * (vyl * byl + vzl * bzl), + byl * vxl - bxs * vyl, + bzl * vxl - bxs * vzl}; + const StateVector fqr{rxr, + rxr * vxr + ptr - bxsq, + rxr * vyr - bxs * byr, + rxr * vzr - bxs * bzr, + vxr * (eer + ptr - bxsq) - bxs * (vyr * byr + vzr * bzr), + byr * vxr - bxs * vyr, + bzr * vxr - bxs * vzr}; + + const double sdl = sl - vxl; + const double sdr = sr - vxr; + const double rosdl = rol * sdl; + const double rosdr = ror * sdr; + const double temp = 1.0 / (rosdr - rosdl); + const double sm = (rosdr * vxr - rosdl * vxl - ptr + ptl) * temp; + const double sdml = sl - sm; + const double sdmr = sr - sm; + const double ptst = (rosdr * ptl - rosdl * ptr + rosdl * rosdr * (vxr - vxl)) * temp; + + const auto sign_unit = [](double x) { return (x >= 0.0) ? 1.0 : -1.0; }; + + const double temp_fst_l = rosdl * sdml - bxsq; + const double sign1_l = sign_unit(std::abs(temp_fst_l) - eps); + const double maxs1_l = std::max(0.0, sign1_l); + const double mins1_l = std::min(0.0, sign1_l); + const double itf_l = 1.0 / (temp_fst_l + mins1_l); + const double isdml = 1.0 / sdml; + + const double temp_l = bxs * (sdl - sdml) * itf_l; + const double rolst = maxs1_l * (rosdl * isdml) - mins1_l * rol; + const double vxlst = maxs1_l * sm - mins1_l * vxl; + const double rxlst = rolst * vxlst; + const double vylst = maxs1_l * (vyl - byl * temp_l) - mins1_l * vyl; + const double rylst = rolst * vylst; + const double vzlst = maxs1_l * (vzl - bzl * temp_l) - mins1_l * vzl; + const double rzlst = rolst * vzlst; + const double temp_l_b = (rosdl * sdl - bxsq) * itf_l; + const double bylst = maxs1_l * (byl * temp_l_b) - mins1_l * byl; + const double bzlst = maxs1_l * (bzl * temp_l_b) - mins1_l * bzl; + const double vdbstl = vxlst * bxs + vylst * bylst + vzlst * bzlst; + const double eelst = maxs1_l * ((sdl * eel - ptl * vxl + ptst * sm + + bxs * (vxl * bxs + vyl * byl + vzl * bzl - vdbstl)) * + isdml) - + mins1_l * eel; + + const double temp_fst_r = rosdr * sdmr - bxsq; + const double sign1_r = sign_unit(std::abs(temp_fst_r) - eps); + const double maxs1_r = std::max(0.0, sign1_r); + const double mins1_r = std::min(0.0, sign1_r); + const double itf_r = 1.0 / (temp_fst_r + mins1_r); + const double isdmr = 1.0 / sdmr; + + const double temp_r = bxs * (sdr - sdmr) * itf_r; + const double rorst = maxs1_r * (rosdr * isdmr) - mins1_r * ror; + const double vxrst = maxs1_r * sm - mins1_r * vxr; + const double rxrst = rorst * vxrst; + const double vyrst = maxs1_r * (vyr - byr * temp_r) - mins1_r * vyr; + const double ryrst = rorst * vyrst; + const double vzrst = maxs1_r * (vzr - bzr * temp_r) - mins1_r * vzr; + const double rzrst = rorst * vzrst; + const double temp_r_b = (rosdr * sdr - bxsq) * itf_r; + const double byrst = maxs1_r * (byr * temp_r_b) - mins1_r * byr; + const double bzrst = maxs1_r * (bzr * temp_r_b) - mins1_r * bzr; + const double vdbstr = vxrst * bxs + vyrst * byrst + vzrst * bzrst; + const double eerst = maxs1_r * ((sdr * eer - ptr * vxr + ptst * sm + + bxs * (vxr * bxs + vyr * byr + vzr * bzr - vdbstr)) * + isdmr) - + mins1_r * eer; + + const double sqrtrol = std::sqrt(rolst); + const double sqrtror = std::sqrt(rorst); + const double abbx = std::abs(bxs); + const double slst = sm - abbx / sqrtrol; + const double srst = sm + abbx / sqrtror; + const double signbx = sign_unit(bxs); + const double sign1_b = sign_unit(abbx - eps); + const double maxs1_b = std::max(0.0, sign1_b); + const double mins1_b = -std::min(0.0, sign1_b); + const double invsumro = maxs1_b / (sqrtrol + sqrtror); + + const double roldst = rolst; + const double rordst = rorst; + const double rxldst = rxlst; + const double rxrdst = rxrst; + const double vxldst = vxlst; + const double vxrdst = vxrst; + + const double vy_shared = + invsumro * (sqrtrol * vylst + sqrtror * vyrst + signbx * (byrst - bylst)); + const double vyldst = vylst * mins1_b + vy_shared; + const double vyrdst = vyrst * mins1_b + vy_shared; + const double ryldst = rylst * mins1_b + roldst * vy_shared; + const double ryrdst = ryrst * mins1_b + rordst * vy_shared; + + const double vz_shared = + invsumro * (sqrtrol * vzlst + sqrtror * vzrst + signbx * (bzrst - bzlst)); + const double vzldst = vzlst * mins1_b + vz_shared; + const double vzrdst = vzrst * mins1_b + vz_shared; + const double rzldst = rzlst * mins1_b + roldst * vz_shared; + const double rzrdst = rzrst * mins1_b + rordst * vz_shared; + + const double by_shared = + invsumro * (sqrtrol * byrst + sqrtror * bylst + signbx * sqrtrol * sqrtror * (vyrst - vylst)); + const double byldst = bylst * mins1_b + by_shared; + const double byrdst = byrst * mins1_b + by_shared; + + const double bz_shared = + invsumro * (sqrtrol * bzrst + sqrtror * bzlst + signbx * sqrtrol * sqrtror * (vzrst - vzlst)); + const double bzldst = bzlst * mins1_b + bz_shared; + const double bzrdst = bzrst * mins1_b + bz_shared; + + const double temp_dst = sm * bxs + vyldst * byldst + vzldst * bzldst; + const double eeldst = eelst - sqrtrol * signbx * (vdbstl - temp_dst) * maxs1_b; + const double eerdst = eerst + sqrtror * signbx * (vdbstr - temp_dst) * maxs1_b; + + const double sign1 = sign_unit(sm); + const double maxs1 = std::max(0.0, sign1); + const double mins1 = -std::min(0.0, sign1); + const double msl = std::min(sl, 0.0); + const double mslst = std::min(slst, 0.0); + const double msrst = std::max(srst, 0.0); + const double msr = std::max(sr, 0.0); + const double temp_flux_l = mslst - msl; + const double temp_flux_r = msrst - msr; + + return StateVector{ + (fql[0] - msl * rol - rolst * temp_flux_l + roldst * mslst) * maxs1 + + (fqr[0] - msr * ror - rorst * temp_flux_r + rordst * msrst) * mins1, + (fql[1] - msl * rxl - rxlst * temp_flux_l + rxldst * mslst) * maxs1 + + (fqr[1] - msr * rxr - rxrst * temp_flux_r + rxrdst * msrst) * mins1, + (fql[2] - msl * ryl - rylst * temp_flux_l + ryldst * mslst) * maxs1 + + (fqr[2] - msr * ryr - ryrst * temp_flux_r + ryrdst * msrst) * mins1, + (fql[3] - msl * rzl - rzlst * temp_flux_l + rzldst * mslst) * maxs1 + + (fqr[3] - msr * rzr - rzrst * temp_flux_r + rzrdst * msrst) * mins1, + (fql[4] - msl * eel - eelst * temp_flux_l + eeldst * mslst) * maxs1 + + (fqr[4] - msr * eer - eerst * temp_flux_r + eerdst * msrst) * mins1, + (fql[5] - msl * byl - bylst * temp_flux_l + byldst * mslst) * maxs1 + + (fqr[5] - msr * byr - byrst * temp_flux_r + byrdst * msrst) * mins1, + (fql[6] - msl * bzl - bzlst * temp_flux_l + bzldst * mslst) * maxs1 + + (fqr[6] - msr * bzr - bzrst * temp_flux_r + bzrdst * msrst) * mins1, + }; +} + +inline StateVector solver_flux_from_primitive(const StateVector& left, const StateVector& right, + double bx, double gamma) +{ + StateVector flux{}; + ::hlld_flux_from_primitive(left.data(), right.data(), bx, gamma, flux.data()); + return flux; +} + +} // namespace hidden_reference diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/eval/tests/cpp/test_hidden.cpp b/benchmarks/magnetohydrodynamics/cpp-hlld-01/eval/tests/cpp/test_hidden.cpp new file mode 100644 index 0000000..35d8074 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/eval/tests/cpp/test_hidden.cpp @@ -0,0 +1,166 @@ +#if __has_include("hlld.hpp") +#include "hlld.hpp" +#else +#include "../../../workspace/src/hlld.hpp" +#endif + +#if __has_include() +#include +#else +#include "/usr/local/include/catch2/catch_test_macros.hpp" +#endif + +#include "hlld_reference.hpp" + +#include +#include + +using StateVector = std::array; + +namespace +{ + +constexpr double kTolerance = 1e-12; + +StateVector primitive_to_conservative(const StateVector& state, double bx, double gamma) +{ + const double rho = state[0]; + const double u = state[1]; + const double v = state[2]; + const double w = state[3]; + const double p = state[4]; + const double by = state[5]; + const double bz = state[6]; + + const double kinetic = 0.5 * rho * (u * u + v * v + w * w); + const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); + const double energy = p / (gamma - 1.0) + kinetic + magnetic; + + return StateVector{rho, rho * u, rho * v, rho * w, energy, by, bz}; +} + +StateVector physical_flux_x(const StateVector& state, double bx, double gamma) +{ + const double rho = state[0]; + const double mx = state[1]; + const double my = state[2]; + const double mz = state[3]; + const double energy = state[4]; + const double by = state[5]; + const double bz = state[6]; + + const double u = mx / rho; + const double v = my / rho; + const double w = mz / rho; + const double kinetic = 0.5 * rho * (u * u + v * v + w * w); + const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); + const double pressure = (gamma - 1.0) * (energy - kinetic - magnetic); + const double total_pressure = pressure + magnetic; + + return StateVector{ + rho * u, + rho * u * u + total_pressure - bx * bx, + rho * v * u - bx * by, + rho * w * u - bx * bz, + (energy + total_pressure) * u - bx * (u * bx + v * by + w * bz), + by * u - bx * v, + bz * u - bx * w, + }; +} + +void require_close(const StateVector& actual, const StateVector& expected) +{ + for (std::size_t i = 0; i < actual.size(); ++i) { + REQUIRE(std::abs(actual[i] - expected[i]) <= kTolerance); + } +} + +StateVector solver_flux_from_primitive(const StateVector& left, const StateVector& right, double bx, + double gamma) +{ + StateVector flux{}; + hlld_flux_from_primitive(left.data(), right.data(), bx, gamma, flux.data()); + return flux; +} + +} // namespace + +TEST_CASE("equal primitive states reduce to the physical flux") +{ + const double bx = 0.35; + const double gamma = 1.4; + const StateVector state{0.9, -0.45, 0.2, 0.15, 0.8, -0.3, 0.55}; + + const StateVector actual = solver_flux_from_primitive(state, state, bx, gamma); + const StateVector expected = + physical_flux_x(primitive_to_conservative(state, bx, gamma), bx, gamma); + + require_close(actual, expected); +} + +TEST_CASE("nontrivial primitive solve returns finite values") +{ + const double bx = -0.65; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.08, 0.45, -0.12, 0.08, 0.95, 0.4, -0.3}; + const StateVector right{0.72, -0.25, 0.16, -0.05, 0.58, -0.2, 0.35}; + + const StateVector flux = solver_flux_from_primitive(left, right, bx, gamma); + + for (double value : flux) { + REQUIRE(std::isfinite(value)); + } +} + +TEST_CASE("hidden reference flux case 1 matches reference implementation") +{ + const double bx = -0.65; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.08, 0.45, -0.12, 0.08, 0.95, 0.4, -0.3}; + const StateVector right{0.72, -0.25, 0.16, -0.05, 0.58, -0.2, 0.35}; + const StateVector expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("hidden reference flux case 2 matches reference implementation") +{ + const double bx = 0.35; + const double gamma = 1.4; + + const StateVector left{0.9, -0.45, 0.2, 0.15, 0.8, -0.3, 0.55}; + const StateVector right{1.15, 0.18, -0.12, -0.08, 1.05, 0.22, -0.4}; + const StateVector expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("small Bx near-degenerate reference case matches reference implementation") +{ + const double bx = 1.0e-6; + const double gamma = 1.4; + + const StateVector left{1.0, 0.4, 0.2, -0.1, 1.0, 0.5, -0.4}; + const StateVector right{0.85, -0.3, -0.15, 0.25, 0.8, -0.35, 0.45}; + const StateVector expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("second Bx equals zero hydro case matches reference implementation") +{ + const double bx = 0.0; + const double gamma = 1.4; + + const StateVector left{0.4, -1.1, 0.0, 0.0, 0.4, 0.0, 0.0}; + const StateVector right{1.2, -0.2, 0.0, 0.0, 1.3, 0.0, 0.0}; + const StateVector expected = hidden_reference::hlld_flux_from_primitive(left, right, bx, gamma); + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/eval/tests/test_hidden.py b/benchmarks/magnetohydrodynamics/cpp-hlld-01/eval/tests/test_hidden.py new file mode 100644 index 0000000..c0e0187 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/eval/tests/test_hidden.py @@ -0,0 +1,36 @@ +import subprocess +from pathlib import Path + + +def _build_hidden_tests() -> Path: + hidden_source = Path("/eval/tests/cpp/test_hidden.cpp") + subprocess.run( + [ + "cmake", + "-S", + ".", + "-B", + "build", + "-DSIMBENCH_ENABLE_HIDDEN_TESTS=ON", + f"-DSIMBENCH_HIDDEN_TEST_SOURCE={hidden_source}", + ], + check=True, + ) + subprocess.run( + ["cmake", "--build", "build", "--target", "hlld_hidden_tests"], check=True + ) + exe = Path("build/tests/hlld_hidden_tests") + assert exe.exists() + return exe + + +def test_hidden_catch2_suite() -> None: + exe = _build_hidden_tests() + proc = subprocess.run( + [str(exe), "--reporter", "compact"], + text=True, + stdout=subprocess.PIPE, + stderr=subprocess.PIPE, + check=False, + ) + assert proc.returncode == 0, proc.stdout + proc.stderr diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/spec.md b/benchmarks/magnetohydrodynamics/cpp-hlld-01/spec.md new file mode 100644 index 0000000..fb106e9 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/spec.md @@ -0,0 +1,29 @@ +# cpp-hlld-01 + +Implement the HLLD approximate Riemann solver for 1D ideal MHD in C++. + +## Read first + +- `/work/basic_equations.md` +- `/work/hlld.md` + +## Task + +Edit `src/hlld.cpp` so that `hlld_flux_from_primitive(...)` is implemented correctly. + +The benchmark uses: + +- primitive-state ordering: `[rho, u, v, w, p, By, Bz]` +- flux ordering: `[F_rho, F_mx, F_my, F_mz, F_E, F_By, F_Bz]` + +Do not change the public function signatures in `src/hlld.hpp`. + +## Standards + +- C++17 + +## Local dev + +```bash +pytest -q +``` diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/task.toml b/benchmarks/magnetohydrodynamics/cpp-hlld-01/task.toml new file mode 100644 index 0000000..de31b95 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/task.toml @@ -0,0 +1,7 @@ +id = "cpp-hlld-01" +suite = "magnetohydrodynamics" +language = "cpp" +time_limit_sec = 600 +eval_cmd = "/eval/run.sh" +prompt = "Read /run/spec.md, /work/basic_equations.md, and /work/hlld.md, then solve the task in /work." +use_shared_workspace = true diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/CMakeLists.txt b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/CMakeLists.txt new file mode 100644 index 0000000..097e6eb --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/CMakeLists.txt @@ -0,0 +1,69 @@ +cmake_minimum_required(VERSION 3.16) + +project(cpp_hlld LANGUAGES CXX) + +set(CMAKE_CXX_STANDARD 17) +set(CMAKE_CXX_STANDARD_REQUIRED ON) +set(CMAKE_CXX_EXTENSIONS OFF) + +include(FetchContent) + +option(SIMBENCH_ENABLE_HIDDEN_TESTS "Build hidden Catch2 tests" OFF) +set( + SIMBENCH_HIDDEN_TEST_SOURCE + "" + CACHE FILEPATH + "Path to hidden Catch2 test source" +) + +find_package(Catch2 3 QUIET) + +if(NOT Catch2_FOUND) + FetchContent_Declare( + Catch2 + GIT_REPOSITORY https://github.com/catchorg/Catch2.git + GIT_TAG v3.13.0 + GIT_SHALLOW TRUE + ) + FetchContent_MakeAvailable(Catch2) +endif() + +add_library(hlld_solver + src/hlld.cpp +) + +target_include_directories(hlld_solver PUBLIC + src +) + +add_executable(hlld_public_tests + tests/cpp/test_public.cpp +) + +target_link_libraries(hlld_public_tests PRIVATE + hlld_solver + Catch2::Catch2WithMain +) + +set_target_properties(hlld_public_tests PROPERTIES + RUNTIME_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/tests" +) + +if(SIMBENCH_ENABLE_HIDDEN_TESTS) + if(NOT EXISTS "${SIMBENCH_HIDDEN_TEST_SOURCE}") + message(FATAL_ERROR "Hidden test source not found: ${SIMBENCH_HIDDEN_TEST_SOURCE}") + endif() + + add_executable(hlld_hidden_tests + "${SIMBENCH_HIDDEN_TEST_SOURCE}" + ) + + target_link_libraries(hlld_hidden_tests PRIVATE + hlld_solver + Catch2::Catch2WithMain + ) + + set_target_properties(hlld_hidden_tests PROPERTIES + RUNTIME_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/tests" + ) +endif() diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/pyproject.toml b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/pyproject.toml new file mode 100644 index 0000000..66bd2af --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/pyproject.toml @@ -0,0 +1,7 @@ +[project] +name = "magnetohydrodynamics-cpp-hlld" +version = "0.0.0" +requires-python = ">=3.10" + +[tool.pytest.ini_options] +testpaths = ["tests"] diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/src/hlld.cpp b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/src/hlld.cpp new file mode 100644 index 0000000..b200baf --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/src/hlld.cpp @@ -0,0 +1,14 @@ +#include "hlld.hpp" + +#include + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux) +{ + (void)left; + (void)right; + (void)bx; + (void)gamma; + + std::fill(flux, flux + 7, 0.0); +} diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/src/hlld.hpp b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/src/hlld.hpp new file mode 100644 index 0000000..ae80ebd --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/src/hlld.hpp @@ -0,0 +1,4 @@ +#pragma once + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux); diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/tests/cpp/test_public.cpp b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/tests/cpp/test_public.cpp new file mode 100644 index 0000000..0c73c44 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/tests/cpp/test_public.cpp @@ -0,0 +1,300 @@ +#if __has_include("hlld.hpp") +#include "hlld.hpp" +#else +#include "../../src/hlld.hpp" +#endif + +#if __has_include() +#include +#else +#include "/usr/local/include/catch2/catch_test_macros.hpp" +#endif + +#include + +#include + +using StateVector = std::array; + +namespace +{ + +constexpr double kTolerance = 1e-12; +constexpr double kPi = 3.14159265358979323846; + +StateVector primitive_to_conservative(const StateVector& state, double bx, double gamma) +{ + const double rho = state[0]; + const double u = state[1]; + const double v = state[2]; + const double w = state[3]; + const double p = state[4]; + const double by = state[5]; + const double bz = state[6]; + + const double kinetic = 0.5 * rho * (u * u + v * v + w * w); + const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); + const double energy = p / (gamma - 1.0) + kinetic + magnetic; + + return StateVector{ + rho, rho * u, rho * v, rho * w, energy, by, bz, + }; +} + +StateVector physical_flux_x(const StateVector& state, double bx, double gamma) +{ + const double rho = state[0]; + const double mx = state[1]; + const double my = state[2]; + const double mz = state[3]; + const double energy = state[4]; + const double by = state[5]; + const double bz = state[6]; + + const double u = mx / rho; + const double v = my / rho; + const double w = mz / rho; + const double kinetic = 0.5 * rho * (u * u + v * v + w * w); + const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); + const double pressure = (gamma - 1.0) * (energy - kinetic - magnetic); + const double total_pressure = pressure + magnetic; + + return StateVector{ + rho * u, + rho * u * u + total_pressure - bx * bx, + rho * v * u - bx * by, + rho * w * u - bx * bz, + (energy + total_pressure) * u - bx * (u * bx + v * by + w * bz), + by * u - bx * v, + bz * u - bx * w, + }; +} + +StateVector solver_flux_from_primitive(const StateVector& left, const StateVector& right, double bx, + double gamma) +{ + StateVector flux{}; + hlld_flux_from_primitive(left.data(), right.data(), bx, gamma, flux.data()); + return flux; +} + +void require_close(const StateVector& actual, const StateVector& expected) +{ + for (std::size_t i = 0; i < actual.size(); ++i) { + REQUIRE(std::abs(actual[i] - expected[i]) <= kTolerance); + } +} + +} // namespace + +TEST_CASE("equal primitive states reduce to the physical flux") +{ + const double bx = 0.75; + const double gamma = 1.4; + const StateVector state{1.1, 0.2, -0.3, 0.4, 0.9, 0.5, -0.6}; + const StateVector conservative = primitive_to_conservative(state, bx, gamma); + + const StateVector actual = solver_flux_from_primitive(state, state, bx, gamma); + const StateVector expected = physical_flux_x(conservative, bx, gamma); + + require_close(actual, expected); +} + +TEST_CASE("right-going contact discontinuity is resolved exactly") +{ + const double bx = 0.8; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.3, 0.2, -0.15, 1.0, 0.6, -0.3}; + const StateVector right{0.7, 0.3, 0.2, -0.15, 1.0, 0.6, -0.3}; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, physical_flux_x(primitive_to_conservative(left, bx, gamma), bx, gamma)); +} + +TEST_CASE("left-going contact discontinuity is resolved exactly") +{ + const double bx = 0.8; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, -0.25, 0.2, -0.15, 1.0, 0.6, -0.3}; + const StateVector right{0.7, -0.25, 0.2, -0.15, 1.0, 0.6, -0.3}; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, physical_flux_x(primitive_to_conservative(right, bx, gamma), bx, gamma)); +} + +TEST_CASE("right-going rotational discontinuity is resolved exactly") +{ + const double bx = 1.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.2, 0.1, -0.2, 1.0, 1.0, 0.0}; + const StateVector right{1.0, 0.2, 0.5, -1.0, 1.0, 0.6, 0.8}; + const StateVector expected{ + 0.2, 1.04, -0.98, -0.04, 0.609, 0.1, 0.2, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("left-going rotational discontinuity is resolved exactly") +{ + const double bx = 1.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.2, 0.1, -0.2, 1.0, 1.0, 0.0}; + const StateVector right{1.0, 0.2, -0.3, 0.6, 1.0, 0.6, 0.8}; + const StateVector expected{ + 0.2, 1.04, -0.66, -0.68, 0.449, 0.42, -0.44, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Bx equals zero hydro case matches reference flux") +{ + const double bx = 0.0; + const double gamma = 1.4; + + const StateVector left{1.0, 0.75, 0.0, 0.0, 1.0, 0.0, 0.0}; + const StateVector right{0.125, 0.0, 0.0, 0.0, 0.1, 0.0, 0.0}; + const StateVector expected{ + 0.92274146439449267, 1.3581095429585437, 0.0, 0.0, 3.1282919538345322, 0.0, 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Bx equals zero magnetized case matches reference flux") +{ + const double bx = 0.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.6, 0.1, -0.2, 1.0, 0.7, -0.5}; + const StateVector right{0.7, -0.3, -0.15, 0.25, 0.5, -0.2, 0.4}; + const StateVector expected{ + 0.44815524807196727, 2.011116795062418, 0.044815524807196722, -0.089631049614393443, + 1.6980640537315086, 0.31370867365037713, -0.22407762403598386, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("small Bx near-degenerate case matches reference flux") +{ + const double bx = 1.0e-6; + const double gamma = 1.4; + + const StateVector left{1.0, 0.4, 0.2, -0.1, 1.0, 0.5, -0.4}; + const StateVector right{0.85, -0.3, -0.15, 0.25, 0.8, -0.35, 0.45}; + const StateVector expected{ + 0.16298732855830989, 1.7549717390289374, 0.032596899426565185, -0.016298279827753587, + 0.72345786285986069, 0.081493464279122407, -0.065194831423298072, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Ryu and Jones shock tube matches reference flux") +{ + const double bx = 4.0 / std::sqrt(4.0 * kPi); + const double gamma = 5.0 / 3.0; + + const StateVector left{ + 1.08, 1.2, 0.01, 0.5, 0.95, 3.6 / std::sqrt(4.0 * kPi), 2.0 / std::sqrt(4.0 * kPi), + }; + const StateVector right{ + 1.0, 0.0, 0.0, 0.0, 1.0, 4.0 / std::sqrt(4.0 * kPi), 2.0 / std::sqrt(4.0 * kPi), + }; + const StateVector expected{ + 0.79485593966715773, 3.5458209484697329, -1.3572358551169827, -0.22185101509215432, + 3.9950643754664625, 0.67495208799031015, -0.062307582042232856, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Brio and Wu shock tube matches reference flux") +{ + const double bx = 0.75; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0}; + const StateVector right{0.125, 0.0, 0.0, 0.0, 0.1, -1.0, 0.0}; + const StateVector expected{ + 0.2063330447744266, + 0.4638678509599396, + 0.064186763013841408, + 0.0, + 0.16136546437466026, + 1.010233243594872, + 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Falle switch-off shock matches reference flux") +{ + const double bx = 1.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.368, 0.269, 1.0, 0.0, 1.769, 0.0, 0.0}; + const StateVector right{1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0}; + const StateVector expected{ + 0.29721990694355238, + 1.4932992607654056, + 0.2893229270591654, + 0.0, + 1.1427267633652525, + -1.0066552479847843, + 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("Falle switch-off rarefaction matches reference flux") +{ + const double bx = 1.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, 0.0, 0.0, 0.0, 2.0, 0.0, 0.0}; + const StateVector right{0.2, 1.186, 2.967, 0.0, 0.1368, 1.6405, 0.0}; + const StateVector expected{ + 0.27717801577960577, + 0.28228035303750848, + -1.3364302412732558, + 0.0, + -1.5599793330037519, + -1.3806947854633354, + 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} + +TEST_CASE("super-fast expansion matches reference flux") +{ + const double bx = 0.0; + const double gamma = 5.0 / 3.0; + + const StateVector left{1.0, -3.0, 0.0, 0.0, 0.45, 0.5, 0.0}; + const StateVector right{1.0, 3.0, 0.0, 0.0, 0.45, 0.5, 0.0}; + const StateVector expected{ + 0.0, -2.425, 0.0, 0.0, 0.0, 0.0, 0.0, + }; + + const StateVector actual = solver_flux_from_primitive(left, right, bx, gamma); + require_close(actual, expected); +} diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/tests/test_public.py b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/tests/test_public.py new file mode 100644 index 0000000..69b5cfa --- /dev/null +++ b/benchmarks/magnetohydrodynamics/cpp-hlld-01/workspace/tests/test_public.py @@ -0,0 +1,24 @@ +import subprocess +from pathlib import Path + + +def _build_public_tests() -> Path: + subprocess.run(["cmake", "-S", ".", "-B", "build"], check=True) + subprocess.run( + ["cmake", "--build", "build", "--target", "hlld_public_tests"], check=True + ) + exe = Path("build/tests/hlld_public_tests") + assert exe.exists() + return exe + + +def test_catch2_public_suite() -> None: + exe = _build_public_tests() + proc = subprocess.run( + [str(exe), "--reporter", "compact"], + text=True, + stdout=subprocess.PIPE, + stderr=subprocess.PIPE, + check=False, + ) + assert proc.returncode == 0, proc.stdout + proc.stderr diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/spec.md b/benchmarks/magnetohydrodynamics/cpp-hlld/spec.md deleted file mode 100644 index 4fbd5e6..0000000 --- a/benchmarks/magnetohydrodynamics/cpp-hlld/spec.md +++ /dev/null @@ -1,40 +0,0 @@ -# cpp-hlld - -Implement the HLLD approximate Riemann solver for 1D ideal MHD in C++. - -## Read first - -- `/work/basic_equations.md` -- `/work/hlld.md` - -## Task - -Edit `src/hlld.cpp` so that these functions are implemented correctly: - -- `hlld_flux_from_primitive(...)` -- `hlld_flux_from_conservative(...)` - -The benchmark uses: - -- primitive-state ordering: `[rho, u, v, w, p, By, Bz]` -- conservative-state ordering: `[rho, mx, my, mz, E, By, Bz]` -- flux ordering: `[F_rho, F_mx, F_my, F_mz, F_E, F_By, F_Bz]` -- Lorentz-Heaviside units -- `Bx` passed separately from the state vectors -- the test suite includes `Bx = 0` hydro and magnetized degenerate cases -- the test suite also includes a small-`Bx` near-degenerate case, so handle - `Bx = 0`, small denominators in the starred-state formulas, and related - square-root/discriminant edge cases carefully - -Do not change the public function signatures in `src/hlld.hpp`. - -## Standards - -- C++17 -- Use `std::array` for the public API - -## Local dev - -```bash -pytest -q -``` diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/src/hlld.cpp b/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/src/hlld.cpp deleted file mode 100644 index b5698a6..0000000 --- a/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/src/hlld.cpp +++ /dev/null @@ -1,31 +0,0 @@ -#include "hlld.hpp" - -namespace -{ - -FluxState zero_flux() -{ - return FluxState{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; -} - -} // namespace - -FluxState hlld_flux_from_primitive(const PrimitiveState& left, const PrimitiveState& right, - double bx, double gamma) -{ - (void)left; - (void)right; - (void)bx; - (void)gamma; - return zero_flux(); -} - -FluxState hlld_flux_from_conservative(const ConservativeState& left, const ConservativeState& right, - double bx, double gamma) -{ - (void)left; - (void)right; - (void)bx; - (void)gamma; - return zero_flux(); -} diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/src/hlld.hpp b/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/src/hlld.hpp deleted file mode 100644 index a5640cb..0000000 --- a/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/src/hlld.hpp +++ /dev/null @@ -1,13 +0,0 @@ -#pragma once - -#include - -using PrimitiveState = std::array; -using ConservativeState = std::array; -using FluxState = std::array; - -FluxState hlld_flux_from_primitive(const PrimitiveState& left, const PrimitiveState& right, - double bx, double gamma); - -FluxState hlld_flux_from_conservative(const ConservativeState& left, const ConservativeState& right, - double bx, double gamma); diff --git a/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/tests/cpp/test_public.cpp b/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/tests/cpp/test_public.cpp deleted file mode 100644 index fc17d25..0000000 --- a/benchmarks/magnetohydrodynamics/cpp-hlld/workspace/tests/cpp/test_public.cpp +++ /dev/null @@ -1,232 +0,0 @@ -#if __has_include("hlld.hpp") -#include "hlld.hpp" -#else -#include "../../src/hlld.hpp" -#endif - -#if __has_include() -#include -#else -#include "/usr/local/include/catch2/catch_test_macros.hpp" -#endif - -#include - -namespace -{ - -constexpr double kTolerance = 1e-12; - -ConservativeState primitive_to_conservative(const PrimitiveState& state, double bx, double gamma) -{ - const double rho = state[0]; - const double u = state[1]; - const double v = state[2]; - const double w = state[3]; - const double p = state[4]; - const double by = state[5]; - const double bz = state[6]; - - const double kinetic = 0.5 * rho * (u * u + v * v + w * w); - const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); - const double energy = p / (gamma - 1.0) + kinetic + magnetic; - - return ConservativeState{ - rho, rho * u, rho * v, rho * w, energy, by, bz, - }; -} - -FluxState physical_flux_x(const ConservativeState& state, double bx, double gamma) -{ - const double rho = state[0]; - const double mx = state[1]; - const double my = state[2]; - const double mz = state[3]; - const double energy = state[4]; - const double by = state[5]; - const double bz = state[6]; - - const double u = mx / rho; - const double v = my / rho; - const double w = mz / rho; - const double kinetic = 0.5 * rho * (u * u + v * v + w * w); - const double magnetic = 0.5 * (bx * bx + by * by + bz * bz); - const double pressure = (gamma - 1.0) * (energy - kinetic - magnetic); - const double total_pressure = pressure + magnetic; - - return FluxState{ - rho * u, - rho * u * u + total_pressure - bx * bx, - rho * v * u - bx * by, - rho * w * u - bx * bz, - (energy + total_pressure) * u - bx * (u * bx + v * by + w * bz), - by * u - bx * v, - bz * u - bx * w, - }; -} - -void require_close(const FluxState& actual, const FluxState& expected) -{ - for (std::size_t i = 0; i < actual.size(); ++i) { - REQUIRE(std::abs(actual[i] - expected[i]) <= kTolerance); - } -} - -} // namespace - -TEST_CASE("equal primitive states reduce to the physical flux") -{ - const double bx = 0.75; - const double gamma = 1.4; - const PrimitiveState state{1.1, 0.2, -0.3, 0.4, 0.9, 0.5, -0.6}; - const ConservativeState conservative = primitive_to_conservative(state, bx, gamma); - - const FluxState actual = hlld_flux_from_primitive(state, state, bx, gamma); - const FluxState expected = physical_flux_x(conservative, bx, gamma); - - require_close(actual, expected); -} - -TEST_CASE("primitive and conservative entry points agree") -{ - const double bx = -0.4; - const double gamma = 5.0 / 3.0; - - const PrimitiveState left{1.0, 0.3, 0.1, -0.2, 1.0, 0.7, -0.5}; - const PrimitiveState right{0.8, -0.1, -0.4, 0.25, 0.7, -0.2, 0.3}; - - const ConservativeState left_cons = primitive_to_conservative(left, bx, gamma); - const ConservativeState right_cons = primitive_to_conservative(right, bx, gamma); - - const FluxState from_primitive = hlld_flux_from_primitive(left, right, bx, gamma); - const FluxState from_conservative = hlld_flux_from_conservative(left_cons, right_cons, bx, gamma); - - require_close(from_primitive, from_conservative); -} - -TEST_CASE("right-going contact discontinuity is resolved exactly") -{ - const double bx = 0.8; - const double gamma = 5.0 / 3.0; - - const PrimitiveState left{1.0, 0.3, 0.2, -0.15, 1.0, 0.6, -0.3}; - const PrimitiveState right{0.7, 0.3, 0.2, -0.15, 1.0, 0.6, -0.3}; - - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); - require_close(actual, physical_flux_x(primitive_to_conservative(left, bx, gamma), bx, gamma)); -} - -TEST_CASE("left-going contact discontinuity is resolved exactly") -{ - const double bx = 0.8; - const double gamma = 5.0 / 3.0; - - const PrimitiveState left{1.0, -0.25, 0.2, -0.15, 1.0, 0.6, -0.3}; - const PrimitiveState right{0.7, -0.25, 0.2, -0.15, 1.0, 0.6, -0.3}; - - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); - require_close(actual, physical_flux_x(primitive_to_conservative(right, bx, gamma), bx, gamma)); -} - -TEST_CASE("right-going rotational discontinuity is resolved exactly") -{ - const double bx = 1.0; - const double gamma = 5.0 / 3.0; - - const PrimitiveState left{1.0, 0.2, 0.1, -0.2, 1.0, 1.0, 0.0}; - const PrimitiveState right{1.0, 0.2, 0.5, -1.0, 1.0, 0.6, 0.8}; - - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); - require_close(actual, physical_flux_x(primitive_to_conservative(left, bx, gamma), bx, gamma)); -} - -TEST_CASE("left-going rotational discontinuity is resolved exactly") -{ - const double bx = 1.0; - const double gamma = 5.0 / 3.0; - - const PrimitiveState left{1.0, 0.2, 0.1, -0.2, 1.0, 1.0, 0.0}; - const PrimitiveState right{1.0, 0.2, -0.3, 0.6, 1.0, 0.6, 0.8}; - - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); - require_close(actual, physical_flux_x(primitive_to_conservative(right, bx, gamma), bx, gamma)); -} - -TEST_CASE("Bx equals zero hydro case matches reference flux") -{ - const double bx = 0.0; - const double gamma = 1.4; - - const PrimitiveState left{1.0, 0.75, 0.0, 0.0, 1.0, 0.0, 0.0}; - const PrimitiveState right{0.125, 0.0, 0.0, 0.0, 0.1, 0.0, 0.0}; - const FluxState expected{ - 0.92274146439449267, 1.3581095429585437, 0.0, 0.0, 3.1282919538345322, 0.0, 0.0, - }; - - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); - require_close(actual, expected); -} - -TEST_CASE("Bx equals zero magnetized case matches reference flux") -{ - const double bx = 0.0; - const double gamma = 5.0 / 3.0; - - const PrimitiveState left{1.0, 0.6, 0.1, -0.2, 1.0, 0.7, -0.5}; - const PrimitiveState right{0.7, -0.3, -0.15, 0.25, 0.5, -0.2, 0.4}; - const FluxState expected{ - 0.44815524807196727, 2.011116795062418, 0.044815524807196722, -0.089631049614393443, - 1.6980640537315086, 0.31370867365037713, -0.22407762403598386, - }; - - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); - require_close(actual, expected); -} - -TEST_CASE("small Bx near-degenerate case matches reference flux") -{ - const double bx = 1.0e-6; - const double gamma = 1.4; - - const PrimitiveState left{1.0, 0.4, 0.2, -0.1, 1.0, 0.5, -0.4}; - const PrimitiveState right{0.85, -0.3, -0.15, 0.25, 0.8, -0.35, 0.45}; - const FluxState expected{ - 0.16298732855830989, 1.7549717390289374, 0.032596899426565185, -0.016298279827753587, - 0.72345786285986069, 0.081493464279122407, -0.065194831423298072, - }; - - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); - require_close(actual, expected); -} - -TEST_CASE("generic coupled MHD case 1 matches reference flux") -{ - const double bx = -0.65; - const double gamma = 5.0 / 3.0; - - const PrimitiveState left{1.08, 0.45, -0.12, 0.08, 0.95, 0.4, -0.3}; - const PrimitiveState right{0.72, -0.25, 0.16, -0.05, 0.58, -0.2, 0.35}; - const FluxState expected{ - 0.33341182495747945, 1.3232304013056353, 0.06667172041085781, -0.02277204332310151, - 0.88458241136476201, 0.20058802344338272, -0.18511862133859658, - }; - - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); - require_close(actual, expected); -} - -TEST_CASE("generic coupled MHD case 2 matches reference flux") -{ - const double bx = 0.35; - const double gamma = 1.4; - - const PrimitiveState left{0.9, -0.45, 0.2, 0.15, 0.8, -0.3, 0.55}; - const PrimitiveState right{1.15, 0.18, -0.12, -0.08, 1.05, 0.22, -0.4}; - const FluxState expected{ - -0.14227120841958368, 0.6095322640180193, 0.028703852070217931, 0.063723510257583313, - -0.40524273615789225, -0.065690714628438368, 0.14700372046240096, - }; - - const FluxState actual = hlld_flux_from_primitive(left, right, bx, gamma); - require_close(actual, expected); -} diff --git a/benchmarks/magnetohydrodynamics/shared/CMakeLists.txt b/benchmarks/magnetohydrodynamics/shared/CMakeLists.txt new file mode 100644 index 0000000..cf98002 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/CMakeLists.txt @@ -0,0 +1,29 @@ +cmake_minimum_required(VERSION 3.16) + +project(magnetohydrodynamics_shared_reference LANGUAGES CXX) + +set(CMAKE_CXX_STANDARD 17) +set(CMAKE_CXX_STANDARD_REQUIRED ON) +set(CMAKE_CXX_EXTENSIONS OFF) + +add_library(mhd1d_reference_lib + src/mhd1d.cpp + src/hlld.cpp +) + +target_include_directories(mhd1d_reference_lib PUBLIC + src + ../../common/include +) + +add_executable(mhd1d_reference + src/main.cpp +) + +target_link_libraries(mhd1d_reference PRIVATE + mhd1d_reference_lib +) + +set_target_properties(mhd1d_reference PROPERTIES + RUNTIME_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/bin" +) diff --git a/benchmarks/magnetohydrodynamics/shared/README.md b/benchmarks/magnetohydrodynamics/shared/README.md new file mode 100644 index 0000000..25ba805 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/README.md @@ -0,0 +1,33 @@ +# Magnetohydrodynamics Shared Assets + +This directory holds the shared reference implementations and the canonical +Brio-Wu fixtures used by the magnetohydrodynamics benchmarks. + +## Build the shared reference solver + +```bash +cmake -S benchmarks/magnetohydrodynamics/shared -B benchmarks/magnetohydrodynamics/shared/build +cmake --build benchmarks/magnetohydrodynamics/shared/build --target mhd1d_reference +``` + +## Run the shared reference solver + +```bash +benchmarks/magnetohydrodynamics/shared/build/bin/mhd1d_reference > benchmarks/magnetohydrodynamics/shared/build/solution.csv +``` + +## Plot the output + +Use the shared plot helper and write the image inside the repo: + +```bash +python3 benchmarks/magnetohydrodynamics/shared/workspace/plot_solution.py \ + benchmarks/magnetohydrodynamics/shared/build/solution.csv \ + benchmarks/magnetohydrodynamics/shared/build/solution.png +``` + +## Run the shared test + +```bash +python3 -m pytest -q benchmarks/magnetohydrodynamics/shared/tests/test_reference.py +``` diff --git a/benchmarks/magnetohydrodynamics/shared/eval/README.md b/benchmarks/magnetohydrodynamics/shared/eval/README.md index 9216b9f..3740465 100644 --- a/benchmarks/magnetohydrodynamics/shared/eval/README.md +++ b/benchmarks/magnetohydrodynamics/shared/eval/README.md @@ -1,3 +1,38 @@ # Shared eval assets -This directory is reserved for future suite-wide hidden-eval helpers. +This directory holds suite-wide hidden-eval documentation for +the full 1D MHD variants (starting with `cpp-full1d-00`). + +## Hidden reference lifecycle + +The full-solver task uses a shared hidden-eval contract anchored by these +paths: + +- `benchmarks/magnetohydrodynamics/shared/eval/fixtures/mhd1d/` + +Fixture files under `fixtures/mhd1d/` store the reference outputs needed to +keep scoring deterministic. + +## Fixture contract + +- CSV schema: `x,rho,u,v,w,p,by,bz` +- Scored variables: `rho`, `u`, `p`, `by` +- CSV headers use lowercase `by` and `bz`; the surrounding solver code and + notation may still refer to the magnetic components as `By` and `Bz`. +- Comparison window: interior cells only, excluding two edge-adjacent cells on + each side +- Regeneration: fixtures are regenerated from the hidden reference pipeline and + must preserve the schema and windowing rule above unless the benchmark + contract is intentionally revised + +## Files + +### `fixtures/mhd1d/` + +- **`brio_wu_reference.csv`**: Reference solution for the canonical Brio-Wu problem (200 cells, `t_final=0.1`) + +## Regenerating fixtures + +To regenerate the reference CSV (maintainer only), run the shared reference +binary with `200` and write its output to +`benchmarks/magnetohydrodynamics/shared/eval/fixtures/mhd1d/brio_wu_reference.csv`. diff --git a/benchmarks/magnetohydrodynamics/shared/eval/fixtures/mhd1d/README.md b/benchmarks/magnetohydrodynamics/shared/eval/fixtures/mhd1d/README.md new file mode 100644 index 0000000..bdf8018 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/eval/fixtures/mhd1d/README.md @@ -0,0 +1,29 @@ +# mhd1d hidden fixtures + +This directory stores the hidden reference fixtures for +the full 1D MHD variants (starting with `cpp-full1d-00`). + +## Contents + +- One or more CSV fixture files with the schema `x,rho,u,v,w,p,by,bz` +- Fixture metadata that records the comparison tolerances `abs_l1` and + `abs_linf` +- Regeneration notes for maintainer use when updating the hidden reference + +CSV fixture headers intentionally stay lowercase for the magnetic fields: +`by` and `bz`. That naming matches the on-disk schema, while the code-level +state and discussion in solver docs may still use `By` and `Bz`. + +## Scoring contract + +- Scored variables: `rho`, `u`, `p`, `by` +- Window: interior cells only, excluding two edge-adjacent cells per side +- Reference comparisons use the fixture-stored `abs_l1` and `abs_linf` + tolerances + +## Regeneration expectations + +Fixtures must be regenerated from the hidden reference pipeline whenever the +benchmark contract changes. Regeneration must preserve the CSV schema, the +interior-cell window, and the stored tolerances unless the suite maintainers +explicitly revise the contract. diff --git a/benchmarks/magnetohydrodynamics/shared/eval/fixtures/mhd1d/brio_wu_reference.csv b/benchmarks/magnetohydrodynamics/shared/eval/fixtures/mhd1d/brio_wu_reference.csv new file mode 100644 index 0000000..1071df9 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/eval/fixtures/mhd1d/brio_wu_reference.csv @@ -0,0 +1,200 @@ +0.0025000000000000001,1,0,0,0,1,1,0 +0.0074999999999999997,1,0,0,0,1,1,0 +0.012500000000000001,1,0,0,0,1,1,0 +0.017500000000000002,1,0,0,0,1,1,0 +0.022499999999999999,1,0,0,0,1,1,0 +0.0275,1,0,0,0,1,1,0 +0.032500000000000001,1,0,0,0,1,1,0 +0.037499999999999999,1,0,0,0,1,1,0 +0.042500000000000003,1,0,0,0,1,1,0 +0.047500000000000001,1,0,0,0,1,1,0 +0.052499999999999998,1,0,0,0,1,1,0 +0.057500000000000002,1,0,0,0,1,1,0 +0.0625,1,0,0,0,1,1,0 +0.067500000000000004,1,0,0,0,1,1,0 +0.072499999999999995,1,0,0,0,1,1,0 +0.077499999999999999,1,0,0,0,1,1,0 +0.082500000000000004,1,0,0,0,1,1,0 +0.087500000000000008,1,0,0,0,1,1,0 +0.092499999999999999,1,0,0,0,1,1,0 +0.097500000000000003,1,0,0,0,1,1,0 +0.10250000000000001,1,0,0,0,1,1,0 +0.1075,1,0,0,0,1,1,0 +0.1125,1,0,0,0,1,1,0 +0.11750000000000001,1,0,0,0,1,1,0 +0.1225,1,0,0,0,1,1,0 +0.1275,1,0,0,0,1,1,0 +0.13250000000000001,1,0,0,0,1,1,0 +0.13750000000000001,1,0,0,0,1,1,0 +0.14250000000000002,1,0,0,0,1,1,0 +0.14749999999999999,1,0,0,0,1,1,0 +0.1525,1,0,0,0,1,1,0 +0.1575,1,0,0,0,1,1,0 +0.16250000000000001,1,0,0,0,1,1,0 +0.16750000000000001,1,0,0,0,1,1,0 +0.17250000000000001,1,0,0,0,1,1,0 +0.17749999999999999,1,0,0,0,1,1,0 +0.1825,1,0,0,0,1,1,0 +0.1875,1,0,0,0,1,1,0 +0.1925,1,0,0,0,1,1,0 +0.19750000000000001,1,0,0,0,1,1,0 +0.20250000000000001,1,0,0,0,1,1,0 +0.20750000000000002,1,0,0,0,1,1,0 +0.21249999999999999,1,0,0,0,1,1,0 +0.2175,1,0,0,0,1,1,0 +0.2225,1,0,0,0,1,1,0 +0.22750000000000001,1,0,0,0,1,1,0 +0.23250000000000001,1,0,0,0,1,1,0 +0.23750000000000002,1,0,0,0,1,1,0 +0.24249999999999999,1,0,0,0,1,1,0 +0.2475,1,0,0,0,1,1,0 +0.2525,1,0,0,0,1,1,0 +0.25750000000000001,1,0,0,0,1,1,0 +0.26250000000000001,1,0,0,0,1,1,0 +0.26750000000000002,1,0,0,0,1,1,0 +0.27250000000000002,1,0,0,0,1,1,0 +0.27750000000000002,1,0,0,0,1,1,0 +0.28250000000000003,1,0,0,0,1,1,0 +0.28750000000000003,1,0,0,0,1,1,0 +0.29249999999999998,1,0,0,0,1,1,0 +0.29749999999999999,1,0,0,0,1,1,0 +0.30249999999999999,1,0,0,0,1,1,0 +0.3075,1,0,0,0,1,1,0 +0.3125,0.99999616106296441,6.8218605150679938e-06,-1.629335642440668e-06,0,0.99999232327941856,0.99999566416503605,0 +0.3175,0.99560068549137515,0.0078843624165065886,-0.0022307064981435336,0,0.99124636392653631,0.99467429565967369,0 +0.32250000000000001,0.98353688945169881,0.029598025828916189,-0.0084627190065405377,0,0.96740300125799372,0.98001580530130683,0 +0.32750000000000001,0.96837732838015789,0.05711492307963683,-0.016488321600118205,0,0.93782681090500097,0.96153830413996133,0 +0.33250000000000002,0.95211065533562045,0.08687376004553142,-0.025340884306684217,0,0.90659226850159857,0.9416391285061726,0 +0.33750000000000002,0.9354673361460476,0.11758295597145299,-0.03467317160538614,0,0.87517584603477561,0.9211946768322925,0 +0.34250000000000003,0.91871718394782453,0.14876408314954759,-0.044363052905967301,0,0.84411769226174038,0.90052767707199355,0 +0.34750000000000003,0.90197078951363818,0.18022527854218023,-0.054366992352566521,0,0.81362498203456801,0.87976588730214367,0 +0.35249999999999998,0.88528078739957949,0.21187079330682687,-0.064670097539751836,0,0.78379400158707657,0.8589695086683371,0 +0.35749999999999998,0.86867910739122456,0.24364631893161978,-0.07526898825581857,0,0.75467281039391498,0.83817058039488479,0 +0.36249999999999999,0.85218625735593245,0.27551533794365413,-0.086165925398151544,0,0.72628825924799156,0.81738872876742574,0 +0.36749999999999999,0.83581667914061153,0.30745246166009871,-0.097367088841251212,0,0.69865388612813017,0.7966357191723874,0 +0.3725,0.81958109452109384,0.33943908417002161,-0.10888139563994977,0,0.67177480879334084,0.77591795858119794,0 +0.3775,0.80348741554278624,0.37146051600257746,-0.1207199694849094,0,0.64565110939024239,0.75523808653390023,0 +0.38250000000000001,0.7875445084338637,0.40350000122001395,-0.13289403440191852,0,0.62028291926991785,0.73459903988049102,0 +0.38750000000000001,0.77176732253101099,0.43552724865404702,-0.14541056405365715,0,0.59567915284429596,0.71401135925001946,0 +0.39250000000000002,0.75618662865505448,0.46747833845388209,-0.15826308205578452,0,0.57187085608759247,0.69350629751332449,0 +0.39750000000000002,0.74086333491683398,0.49922382920992259,-0.1714159845725442,0,0.54893025602556722,0.67315622872434921,0 +0.40250000000000002,0.72591127957332113,0.53051685817985217,-0.18477737401520411,0,0.52699952756086399,0.65310782756232943,0 +0.40750000000000003,0.71153792819892647,0.56090161794556215,-0.19815008755163185,0,0.50634126527317558,0.63364102400425049,0 +0.41250000000000003,0.69811514572358369,0.58954799451371021,-0.21114189870363537,0,0.4874249963437069,0.61527302327312272,0 +0.41749999999999998,0.68628746003287744,0.61503330090041819,-0.22302766944126134,0,0.47105684401465941,0.59892135357981213,0 +0.42249999999999999,0.67696251242649175,0.63513291356954815,-0.2326074173770955,0,0.45834791401272101,0.58598837050788355,0 +0.42749999999999999,0.67147463576702993,0.64779945118301063,-0.23846691522567326,0,0.45094670813117532,0.57820303829922137,0 +0.4325,0.6688796555738501,0.65110241055546891,-0.2410081467465312,0,0.44746033225767845,0.5753899843373016,0 +0.4375,0.66923632078057105,0.65190937533001247,-0.24110017994523736,0,0.44794907534325185,0.57444176371495947,0 +0.4425,0.67092296635432236,0.64955585559206386,-0.24054877369726702,0,0.45022622769647191,0.57572465421454055,0 +0.44750000000000001,0.6741213806450993,0.64344630049305807,-0.23854099129835965,0,0.45453459644722283,0.57947771182222196,0 +0.45250000000000001,0.67622598638333231,0.63780864306743901,-0.23448685862138482,0,0.45738500592687004,0.58438626477884326,0 +0.45750000000000002,0.67754413933073854,0.63405079125731079,-0.22959378848719866,0,0.45917700333007561,0.58929788380231585,0 +0.46250000000000002,0.68597164204564365,0.61952213982777449,-0.24653663418395505,0,0.47189617111806659,0.57922868126312688,0 +0.46750000000000003,0.75354858718752427,0.57292925776279291,-0.63081628168738757,0,0.60292314122276669,0.27209063114350734,0 +0.47250000000000003,0.80630536650231177,0.46293683334979041,-1.1759576407781362,0,0.69703548985821029,-0.20853604373803852,0 +0.47750000000000004,0.7638410805034388,0.51805820703977745,-1.3927632353963115,0,0.62805748621970625,-0.3981907442021011,0 +0.48249999999999998,0.73152958163331872,0.57826425255304303,-1.4832959440212241,0,0.57564691619642827,-0.47630190401752948,0 +0.48749999999999999,0.7091371949877272,0.59282038528643355,-1.5447516257458402,0,0.54057240122232608,-0.52152146798758336,0 +0.49249999999999999,0.69094050613943225,0.59503499024494011,-1.5888099484071434,0,0.51755764390120462,-0.52927763675904382,0 +0.4975,0.68445173626467959,0.59046857583704693,-1.6048335641093068,0,0.50638419398254775,-0.53146412709505564,0 +0.50250000000000006,0.68434384603578657,0.59230376347482205,-1.5982075803435103,0,0.50542098868729934,-0.53352729015156963,0 +0.50750000000000006,0.68707431747700531,0.59425987403200542,-1.5869958639330111,0,0.51033146073036029,-0.53448849532455134,0 +0.51249999999999996,0.68939472693221426,0.60137405638102948,-1.5798014603878261,0,0.51803593468141906,-0.53565082760707183,0 +0.51749999999999996,0.69132926018866214,0.60874084896105818,-1.5797632463390414,0,0.52436693339499907,-0.5359104689992632,0 +0.52249999999999996,0.69215217980589738,0.60970112925016673,-1.5824896135808113,0,0.52434513683953443,-0.5354819845059674,0 +0.52749999999999997,0.69444170391589999,0.60627151781387467,-1.5858651859518866,0,0.52162449544034473,-0.53457955249120248,0 +0.53249999999999997,0.69647681591950339,0.59800014203747298,-1.5883464457585819,0,0.51533551298586344,-0.53253950669423489,0 +0.53749999999999998,0.69477389476214069,0.59293397164493389,-1.5901196445037074,0,0.51174059849537212,-0.5304114113079087,0 +0.54249999999999998,0.68911228688697801,0.59380160058459075,-1.5903316488178545,0,0.51192732192152268,-0.52949995498960822,0 +0.54749999999999999,0.66719536897441889,0.59857474030644331,-1.5895216549318072,0,0.51429507796334395,-0.52929052879654892,0 +0.55249999999999999,0.60234055420437194,0.60371440419440114,-1.5885169423087238,0,0.5187332353048264,-0.53155866918735362,0 +0.5575,0.49290738445262527,0.60775103688327448,-1.5855974157471353,0,0.52066330116084492,-0.53424754867164737,0 +0.5625,0.36379960562925234,0.61145270056068901,-1.5836004669317623,0,0.51969239083206853,-0.53490394589485069,0 +0.5675,0.25794003992686892,0.61029128984346714,-1.5820890229303581,0,0.51813142258312772,-0.53525310083950484,0 +0.57250000000000001,0.22619669625294542,0.60434480627709164,-1.583239013708023,0,0.51937704488324277,-0.53667920987231343,0 +0.57750000000000001,0.22630840156123649,0.59827581714693934,-1.5861518811591599,0,0.51884369894384308,-0.53653660134365799,0 +0.58250000000000002,0.22711466733509919,0.59921584970885489,-1.5865907932320387,0,0.51663519067912,-0.5340362400001385,0 +0.58750000000000002,0.22848173751344047,0.60091449520075879,-1.5881448209869364,0,0.51354656892651396,-0.52979939734340753,0 +0.59250000000000003,0.23031920669352701,0.59922369951386523,-1.5926693039478363,0,0.51316449769217121,-0.52713616980108124,0 +0.59750000000000003,0.23243460374842745,0.59022749017803877,-1.5967363554938299,0,0.5171218957694993,-0.52858990082849489,0 +0.60250000000000004,0.23371045610745741,0.58557425058898349,-1.6007158568053828,0,0.5198767628355706,-0.5306097118265467,0 +0.60750000000000004,0.23346589711062016,0.59462702763778397,-1.5950312955893633,0,0.51701289109504156,-0.52896536404700278,0 +0.61250000000000004,0.23263653232779033,0.61187045984626132,-1.5834576247676977,0,0.51177702375221379,-0.52633641566767586,0 +0.61750000000000005,0.23321499947212285,0.61526177184194952,-1.5789252785842804,0,0.51201629462921838,-0.53124179083925394,0 +0.62250000000000005,0.23616745617017731,0.59110549784028277,-1.581549220607541,0,0.52255821054730711,-0.54281830257957719,0 +0.62750000000000006,0.23850998490201497,0.57982203884810612,-1.5800134018812597,0,0.52973514694150359,-0.54974364304552559,0 +0.63250000000000006,0.23549928160809872,0.60805510179030442,-1.5697422964899919,0,0.51469370924936686,-0.54625158631214921,0 +0.63750000000000007,0.23267730073699258,0.64764426075074666,-1.5303717022273871,0,0.50048958343353256,-0.53283346050863811,0 +0.64249999999999996,0.22589478256765788,0.51608311612092028,-1.4244966323443713,0,0.47215254213952862,-0.60337924172895407,0 +0.64749999999999996,0.16824385093540103,0.2019403880561518,-0.85229593508650303,0,0.24693902186307282,-0.76893242381717286,0 +0.65249999999999997,0.12072583773152777,-0.18183444417138045,-0.24105266177015217,0,0.095400355892026001,-0.89302992881503718,0 +0.65749999999999997,0.11601292282091824,-0.26459528318540043,-0.1814730538647312,0,0.086167195075565028,-0.89320114714371157,0 +0.66249999999999998,0.116048367245578,-0.26722541207620687,-0.18522210521227256,0,0.086217867928268577,-0.89201260813637107,0 +0.66749999999999998,0.11614130837574657,-0.26549235054019255,-0.18609509830657162,0,0.086355655034603629,-0.891902181568761,0 +0.67249999999999999,0.11631511438200227,-0.26085265286182824,-0.18311506509197725,0,0.086614152150550017,-0.89367317333644192,0 +0.67749999999999999,0.11662767049075361,-0.25099871037849603,-0.17536214658827623,0,0.087078445740005561,-0.8978274795799227,0 +0.6825,0.1170518312215577,-0.2381722667902374,-0.16579466075220212,0,0.087710449613822461,-0.90318269940017126,0 +0.6875,0.11746096698143836,-0.22586176858556523,-0.15680732789253898,0,0.088322007878451791,-0.90825069334609754,0 +0.6925,0.11767715987203248,-0.21896280887700917,-0.15179332435332796,0,0.088645332706841606,-0.91092215760728412,0 +0.69750000000000001,0.11768511349561916,-0.21852003035735557,-0.15146950280658006,0,0.088656537990149054,-0.91102013385933844,0 +0.70250000000000001,0.11768173540960905,-0.21911190725269286,-0.15189586068806646,0,0.088650899743895395,-0.91098098981221598,0 +0.70750000000000002,0.11762763379381724,-0.2205495383437531,-0.15293640182544183,0,0.088569191874156683,-0.91031519004692951,0 +0.71250000000000002,0.11751251674883506,-0.22407404665911851,-0.15549737468557176,0,0.088395415678045763,-0.90889140732316487,0 +0.71750000000000003,0.11734217026306956,-0.22923209052659663,-0.15925236589498573,0,0.088138555377449745,-0.90678117162710703,0 +0.72250000000000003,0.11714842121554768,-0.23510643408875742,-0.16353933245426788,0,0.087847174957777496,-0.90438202392384848,0 +0.72750000000000004,0.1169801447659943,-0.2402468788526454,-0.16730065580704034,0,0.08759435163307927,-0.9022965309078298,0 +0.73250000000000004,0.11690646094308876,-0.24252420707824915,-0.16896758133008036,0,0.087483274145631063,-0.901381475968559,0 +0.73750000000000004,0.11690536582831555,-0.24256297970713472,-0.16899337040831261,0,0.087481058709982462,-0.90136639345779135,0 +0.74250000000000005,0.11690416351837025,-0.24255759967720394,-0.16898805965890309,0,0.087479260585539786,-0.901355959768676,0 +0.74750000000000005,0.11690268830300561,-0.24256433968132488,-0.1689902511692872,0,0.087476844575879453,-0.90133971116149714,0 +0.75250000000000006,0.11690249149566564,-0.24263045674337938,-0.1690369990836007,0,0.087475993304401478,-0.90133579739871572,0 +0.75750000000000006,0.11690005841649749,-0.24273633487321142,-0.16911280871212495,0,0.087472007939943586,-0.90130651809972928,0 +0.76250000000000007,0.11689523976762291,-0.24284738859395841,-0.16919278965713455,0,0.087464802097937722,-0.90124869217347503,0 +0.76750000000000007,0.11689009550937005,-0.24294211733393109,-0.16926071280315277,0,0.087456699244402802,-0.90118542866037055,0 +0.77249999999999996,0.1168880938125931,-0.24302624712308077,-0.16932082611092356,0,0.087453261117748871,-0.90115780061563966,0 +0.77749999999999997,0.11688902120163681,-0.24306621660903155,-0.16935040473201921,0,0.08745446500405607,-0.90117291963133284,0 +0.78249999999999997,0.11689122011863207,-0.24300391365362053,-0.16930103423717782,0,0.08745781154989074,-0.90120105561337527,0 +0.78749999999999998,0.11689602107686169,-0.24276853710995247,-0.16912688892459521,0,0.087464414434553439,-0.90126023192649951,0 +0.79249999999999998,0.11691614491470963,-0.24220891469038228,-0.16871533218791021,0,0.087494250613795788,-0.9015098119657331,0 +0.79749999999999999,0.11696855155695048,-0.2405988356047073,-0.16753372147752355,0,0.087572799130735901,-0.90216329808472961,0 +0.80249999999999999,0.11707534180756879,-0.23734745289851361,-0.16515050554824207,0,0.087732641828922264,-0.90348892427237504,0 +0.8075,0.1172553580143333,-0.23187377620745053,-0.16114702890370325,0,0.088002312969811203,-0.90571948515068501,0 +0.8125,0.11752273565815201,-0.22374249463054305,-0.15521835735527145,0,0.088404020315360543,-0.90903234177390579,0 +0.8175,0.11788504329153121,-0.21273139165106736,-0.14722424729033728,0,0.088950161329339306,-0.91351723197970469,0 +0.82250000000000001,0.1183424161062748,-0.19886886175669347,-0.13721495011558987,0,0.089641293437578984,-0.91916169859620123,0 +0.82750000000000001,0.11888652447249262,-0.18241324922292035,-0.12541146812970011,0,0.090466824054350359,-0.92585977532177122,0 +0.83250000000000002,0.11950366948367314,-0.1637828406974208,-0.11215088119516084,0,0.091408372844463459,-0.93344004869433117,0 +0.83750000000000002,0.12017840890881823,-0.1434726819285532,-0.09781794530184125,0,0.092443075882257864,-0.94170037634875625,0 +0.84250000000000003,0.12089469678105615,-0.12198057865279255,-0.082786700196053442,0,0.093547607150673007,-0.95043832267570283,0 +0.84750000000000003,0.1216371976196647,-0.099767234362176133,-0.067394749995541334,0,0.094699706483174273,-0.95946680774581716,0 +0.85250000000000004,0.12239161658091287,-0.077260468431275015,-0.051947511873339444,0,0.095878018276146726,-0.96861198920602887,0 +0.85750000000000004,0.12314341274778251,-0.05490629533724703,-0.03674915962168223,0,0.097058875646527043,-0.97769350318200499,0 +0.86250000000000004,0.12386678297689156,-0.033461759001210456,-0.022298916108968211,0,0.098201102970602916,-0.98640512630263821,0 +0.86750000000000005,0.1245109747847933,-0.014416962267587185,-0.0095697786663570452,0,0.099223465212322837,-0.99414146906183942,0 +0.87250000000000005,0.12493817800639674,-0.0018204925661910309,-0.0012033527248810462,0,0.09990183090571092,-0.99926077800527913,0 +0.87750000000000006,0.12500000001219974,3.6363615818642719e-10,2.4056984550234209e-10,0,0.10000000001951947,-1.0000000001460236,0 +0.88250000000000006,0.12500000001221473,3.5619382255045683e-10,2.3564550308747709e-10,0,0.10000000001954346,-1.0000000001462062,0 +0.88750000000000007,0.12500000001137598,3.3563244772973727e-10,2.2204309500150842e-10,0,0.10000000001820109,-1.0000000001361644,0 +0.89250000000000007,0.12500000000941866,2.775282966192311e-10,1.8360388763520491e-10,0,0.10000000001507026,-1.0000000001127354,0 +0.89749999999999996,0.12500000000652875,1.9240316746228584e-10,1.2728653685957096e-10,0,0.10000000001044618,-1.0000000000781459,0 +0.90249999999999997,0.12500000000334643,9.8618158081562529e-11,6.5242655720636596e-11,0,0.1000000000053537,-1.0000000000400553,0 +0.90749999999999997,0.12500000000075906,2.2369069559482915e-11,1.4798814025978548e-11,0,0.10000000000121445,-1.0000000000090854,0 +0.91249999999999998,0.125,0,0,0,0.099999999999999867,-1,0 +0.91749999999999998,0.125,0,0,0,0.099999999999999867,-1,0 +0.92249999999999999,0.125,0,0,0,0.099999999999999867,-1,0 +0.92749999999999999,0.125,0,0,0,0.099999999999999867,-1,0 +0.9325,0.125,0,0,0,0.099999999999999867,-1,0 +0.9375,0.125,0,0,0,0.099999999999999867,-1,0 +0.9425,0.125,0,0,0,0.099999999999999867,-1,0 +0.94750000000000001,0.125,0,0,0,0.099999999999999867,-1,0 +0.95250000000000001,0.125,0,0,0,0.099999999999999867,-1,0 +0.95750000000000002,0.125,0,0,0,0.099999999999999867,-1,0 +0.96250000000000002,0.125,0,0,0,0.099999999999999867,-1,0 +0.96750000000000003,0.125,0,0,0,0.099999999999999867,-1,0 +0.97250000000000003,0.125,0,0,0,0.099999999999999867,-1,0 +0.97750000000000004,0.125,0,0,0,0.099999999999999867,-1,0 +0.98250000000000004,0.125,0,0,0,0.099999999999999867,-1,0 +0.98750000000000004,0.125,0,0,0,0.099999999999999867,-1,0 +0.99250000000000005,0.125,0,0,0,0.099999999999999867,-1,0 +0.99750000000000005,0.125,0,0,0,0.099999999999999867,-1,0 diff --git a/benchmarks/magnetohydrodynamics/shared/eval/mhd1d_shared.py b/benchmarks/magnetohydrodynamics/shared/eval/mhd1d_shared.py new file mode 100644 index 0000000..94f8424 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/eval/mhd1d_shared.py @@ -0,0 +1,27 @@ +from __future__ import annotations + +from typing import Sequence + + +TOLERANCE = 1.0e-12 + + +def assert_csv_rows_close( + output_rows: Sequence[Sequence[str]], + reference_rows: Sequence[Sequence[str]], + *, + tolerance: float = TOLERANCE, + expected_header: Sequence[str] | None = None, +) -> None: + if expected_header is not None: + assert list(output_rows[0]) == list(expected_header) + output_rows = output_rows[1:] + + assert len(output_rows) == len(reference_rows) + + for output_row, reference_row in zip(output_rows, reference_rows): + assert len(output_row) == len(reference_row) + for output_value, reference_value in zip(output_row, reference_row): + output_float = float(output_value) + reference_float = float(reference_value) + assert abs(output_float - reference_float) <= tolerance diff --git a/benchmarks/magnetohydrodynamics/shared/src/hlld.cpp b/benchmarks/magnetohydrodynamics/shared/src/hlld.cpp new file mode 100644 index 0000000..f5f9c87 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/src/hlld.cpp @@ -0,0 +1,202 @@ +#include "hlld.hpp" + +#include +#include + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux) +{ + constexpr double epsilon = 1.0e-40; + + const double rol = left[0]; + const double vxl = left[1]; + const double vyl = left[2]; + const double vzl = left[3]; + const double prl = left[4]; + const double byl = left[5]; + const double bzl = left[6]; + + const double ror = right[0]; + const double vxr = right[1]; + const double vyr = right[2]; + const double vzr = right[3]; + const double prr = right[4]; + const double byr = right[5]; + const double bzr = right[6]; + + const double igm = 1.0 / (gamma - 1.0); + const double bxs = bx; + const double bxsq = bxs * bxs; + + const double pbl = 0.5 * (bxsq + byl * byl + bzl * bzl); + const double pbr = 0.5 * (bxsq + byr * byr + bzr * bzr); + const double ptl = prl + pbl; + const double ptr = prr + pbr; + + const double rxl = rol * vxl; + const double ryl = rol * vyl; + const double rzl = rol * vzl; + const double rxr = ror * vxr; + const double ryr = ror * vyr; + const double rzr = ror * vzr; + + const double eel = prl * igm + 0.5 * (rxl * vxl + ryl * vyl + rzl * vzl) + pbl; + const double eer = prr * igm + 0.5 * (rxr * vxr + ryr * vyr + rzr * vzr) + pbr; + + const double gmpl = gamma * prl; + const double gmpr = gamma * prr; + const double gpbl = gmpl + 2.0 * pbl; + const double gpbr = gmpr + 2.0 * pbr; + + const double cfl = std::sqrt((gpbl + std::sqrt((gmpl - 2.0 * pbl) * (gmpl - 2.0 * pbl) + + 4.0 * gmpl * (byl * byl + bzl * bzl))) * + 0.5 / rol); + const double cfr = std::sqrt((gpbr + std::sqrt((gmpr - 2.0 * pbr) * (gmpr - 2.0 * pbr) + + 4.0 * gmpr * (byr * byr + bzr * bzr))) * + 0.5 / ror); + + const double sl = std::min(vxl, vxr) - std::max(cfl, cfr); + const double sr = std::max(vxl, vxr) + std::max(cfl, cfr); + + const double fql[7] = {rxl, + rxl * vxl + ptl - bxsq, + rxl * vyl - bxs * byl, + rxl * vzl - bxs * bzl, + vxl * (eel + ptl - bxsq) - bxs * (vyl * byl + vzl * bzl), + byl * vxl - bxs * vyl, + bzl * vxl - bxs * vzl}; + const double fqr[7] = {rxr, + rxr * vxr + ptr - bxsq, + rxr * vyr - bxs * byr, + rxr * vzr - bxs * bzr, + vxr * (eer + ptr - bxsq) - bxs * (vyr * byr + vzr * bzr), + byr * vxr - bxs * vyr, + bzr * vxr - bxs * vzr}; + + const double sdl = sl - vxl; + const double sdr = sr - vxr; + const double rosdl = rol * sdl; + const double rosdr = ror * sdr; + const double temp = 1.0 / (rosdr - rosdl); + const double sm = (rosdr * vxr - rosdl * vxl - ptr + ptl) * temp; + const double sdml = sl - sm; + const double sdmr = sr - sm; + const double ptst = (rosdr * ptl - rosdl * ptr + rosdl * rosdr * (vxr - vxl)) * temp; + + const double temp_fst_l = rosdl * sdml - bxsq; + const double sign1_l = std::copysign(1.0, std::abs(temp_fst_l) - epsilon); + const double maxs1_l = std::max(0.0, sign1_l); + const double mins1_l = std::min(0.0, sign1_l); + const double itf_l = 1.0 / (temp_fst_l + mins1_l); + const double isdml = 1.0 / sdml; + + const double temp_l = bxs * (sdl - sdml) * itf_l; + const double rolst = maxs1_l * (rosdl * isdml) - mins1_l * rol; + const double vxlst = maxs1_l * sm - mins1_l * vxl; + const double rxlst = rolst * vxlst; + const double vylst = maxs1_l * (vyl - byl * temp_l) - mins1_l * vyl; + const double rylst = rolst * vylst; + const double vzlst = maxs1_l * (vzl - bzl * temp_l) - mins1_l * vzl; + const double rzlst = rolst * vzlst; + const double temp_l_b = (rosdl * sdl - bxsq) * itf_l; + const double bylst = maxs1_l * (byl * temp_l_b) - mins1_l * byl; + const double bzlst = maxs1_l * (bzl * temp_l_b) - mins1_l * bzl; + const double vdbstl = vxlst * bxs + vylst * bylst + vzlst * bzlst; + const double eelst = maxs1_l * ((sdl * eel - ptl * vxl + ptst * sm + + bxs * (vxl * bxs + vyl * byl + vzl * bzl - vdbstl)) * + isdml) - + mins1_l * eel; + + const double temp_fst_r = rosdr * sdmr - bxsq; + const double sign1_r = std::copysign(1.0, std::abs(temp_fst_r) - epsilon); + const double maxs1_r = std::max(0.0, sign1_r); + const double mins1_r = std::min(0.0, sign1_r); + const double itf_r = 1.0 / (temp_fst_r + mins1_r); + const double isdmr = 1.0 / sdmr; + + const double temp_r = bxs * (sdr - sdmr) * itf_r; + const double rorst = maxs1_r * (rosdr * isdmr) - mins1_r * ror; + const double vxrst = maxs1_r * sm - mins1_r * vxr; + const double rxrst = rorst * vxrst; + const double vyrst = maxs1_r * (vyr - byr * temp_r) - mins1_r * vyr; + const double ryrst = rorst * vyrst; + const double vzrst = maxs1_r * (vzr - bzr * temp_r) - mins1_r * vzr; + const double rzrst = rorst * vzrst; + const double temp_r_b = (rosdr * sdr - bxsq) * itf_r; + const double byrst = maxs1_r * (byr * temp_r_b) - mins1_r * byr; + const double bzrst = maxs1_r * (bzr * temp_r_b) - mins1_r * bzr; + const double vdbstr = vxrst * bxs + vyrst * byrst + vzrst * bzrst; + const double eerst = maxs1_r * ((sdr * eer - ptr * vxr + ptst * sm + + bxs * (vxr * bxs + vyr * byr + vzr * bzr - vdbstr)) * + isdmr) - + mins1_r * eer; + + const double sqrtrol = std::sqrt(rolst); + const double sqrtror = std::sqrt(rorst); + const double abbx = std::abs(bxs); + const double slst = sm - abbx / sqrtrol; + const double srst = sm + abbx / sqrtror; + const double signbx = std::copysign(1.0, bxs); + const double sign1_b = std::copysign(1.0, abbx - epsilon); + const double maxs1_b = std::max(0.0, sign1_b); + const double mins1_b = -std::min(0.0, sign1_b); + const double invsumro = maxs1_b / (sqrtrol + sqrtror); + + const double roldst = rolst; + const double rordst = rorst; + const double rxldst = rxlst; + const double rxrdst = rxrst; + + const double vy_shared = + invsumro * (sqrtrol * vylst + sqrtror * vyrst + signbx * (byrst - bylst)); + const double ryldst = rylst * mins1_b + roldst * vy_shared; + const double ryrdst = ryrst * mins1_b + rordst * vy_shared; + + const double vz_shared = + invsumro * (sqrtrol * vzlst + sqrtror * vzrst + signbx * (bzrst - bzlst)); + const double rzldst = rzlst * mins1_b + roldst * vz_shared; + const double rzrdst = rzrst * mins1_b + rordst * vz_shared; + + const double by_shared = + invsumro * (sqrtrol * byrst + sqrtror * bylst + signbx * sqrtrol * sqrtror * (vyrst - vylst)); + const double byldst = bylst * mins1_b + by_shared; + const double byrdst = byrst * mins1_b + by_shared; + + const double bz_shared = + invsumro * (sqrtrol * bzrst + sqrtror * bzlst + signbx * sqrtrol * sqrtror * (vzrst - vzlst)); + const double bzldst = bzlst * mins1_b + bz_shared; + const double bzrdst = bzrst * mins1_b + bz_shared; + + const double vyldst = vylst * mins1_b + vy_shared; + const double vyrdst = vyrst * mins1_b + vy_shared; + const double vzldst = vzlst * mins1_b + vz_shared; + const double vzrdst = vzrst * mins1_b + vz_shared; + const double temp_dst = sm * bxs + vyldst * byldst + vzldst * bzldst; + const double eeldst = eelst - sqrtrol * signbx * (vdbstl - temp_dst) * maxs1_b; + const double eerdst = eerst + sqrtror * signbx * (vdbstr - temp_dst) * maxs1_b; + + const double sign1 = std::copysign(1.0, sm); + const double maxs1 = std::max(0.0, sign1); + const double mins1 = -std::min(0.0, sign1); + const double msl = std::min(sl, 0.0); + const double mslst = std::min(slst, 0.0); + const double msrst = std::max(srst, 0.0); + const double msr = std::max(sr, 0.0); + const double temp_flux_l = mslst - msl; + const double temp_flux_r = msrst - msr; + + flux[0] = (fql[0] - msl * rol - rolst * temp_flux_l + roldst * mslst) * maxs1 + + (fqr[0] - msr * ror - rorst * temp_flux_r + rordst * msrst) * mins1; + flux[1] = (fql[1] - msl * rxl - rxlst * temp_flux_l + rxldst * mslst) * maxs1 + + (fqr[1] - msr * rxr - rxrst * temp_flux_r + rxrdst * msrst) * mins1; + flux[2] = (fql[2] - msl * ryl - rylst * temp_flux_l + ryldst * mslst) * maxs1 + + (fqr[2] - msr * ryr - ryrst * temp_flux_r + ryrdst * msrst) * mins1; + flux[3] = (fql[3] - msl * rzl - rzlst * temp_flux_l + rzldst * mslst) * maxs1 + + (fqr[3] - msr * rzr - rzrst * temp_flux_r + rzrdst * msrst) * mins1; + flux[4] = (fql[4] - msl * eel - eelst * temp_flux_l + eeldst * mslst) * maxs1 + + (fqr[4] - msr * eer - eerst * temp_flux_r + eerdst * msrst) * mins1; + flux[5] = (fql[5] - msl * byl - bylst * temp_flux_l + byldst * mslst) * maxs1 + + (fqr[5] - msr * byr - byrst * temp_flux_r + byrdst * msrst) * mins1; + flux[6] = (fql[6] - msl * bzl - bzlst * temp_flux_l + bzldst * mslst) * maxs1 + + (fqr[6] - msr * bzr - bzrst * temp_flux_r + bzrdst * msrst) * mins1; +} diff --git a/benchmarks/magnetohydrodynamics/shared/src/hlld.hpp b/benchmarks/magnetohydrodynamics/shared/src/hlld.hpp new file mode 100644 index 0000000..ae80ebd --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/src/hlld.hpp @@ -0,0 +1,4 @@ +#pragma once + +void hlld_flux_from_primitive(const double* left, const double* right, double bx, double gamma, + double* flux); diff --git a/benchmarks/magnetohydrodynamics/shared/src/main.cpp b/benchmarks/magnetohydrodynamics/shared/src/main.cpp new file mode 100644 index 0000000..335945f --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/src/main.cpp @@ -0,0 +1,74 @@ +#include "mhd1d.hpp" + +#include +#include +#include +#include + +constexpr int Nx = 100; +constexpr double Gamma = 2.0; +constexpr double Bx = 0.75; +constexpr mhd1d::StateVector LeftPrimitive{ + 1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0, +}; +constexpr mhd1d::StateVector RightPrimitive{ + 0.125, 0.0, 0.0, 0.0, 0.1, -1.0, 0.0, +}; + +int parse_nx(int argc, char** argv) +{ + if (argc <= 1) { + return Nx; + } + + char* end = nullptr; + const long parsed = std::strtol(argv[1], &end, 10); + if (end == argv[1] || *end != '\0' || parsed <= 0) { + throw std::runtime_error("usage: solver [nx]"); + } + return static_cast(parsed); +} + +mhd1d::SolverWorkspace initialize(int nx, double gamma, double bx, + const mhd1d::StateVector& left_state, + const mhd1d::StateVector& right_state) +{ + mhd1d::SolverWorkspace workspace(nx, gamma, bx); + + for (int ix = workspace.Lbx; ix <= workspace.Ubx; ++ix) { + const mhd1d::StateVector& state = (workspace.x(ix) < 0.5) ? left_state : right_state; + for (int component = 0; component < mhd1d::N_Component; ++component) { + workspace.up(ix, component) = state[component]; + } + } + + mhd1d::set_boundary(workspace.up, workspace.up, workspace.Lbx, workspace.Ubx); + mhd1d::convert_primitive_to_conservative(workspace.up, workspace.uc, bx, gamma); + + return workspace; +} + +void write_csv(const mhd1d::SolverWorkspace& workspace, std::ostream& os) +{ + os << std::setprecision(17); + for (int ix = workspace.Lbx; ix <= workspace.Ubx; ++ix) { + os << workspace.x(ix) << ',' << workspace.up(ix, 0) << ',' << workspace.up(ix, 1) << ',' + << workspace.up(ix, 2) << ',' << workspace.up(ix, 3) << ',' << workspace.up(ix, 4) << ',' + << workspace.up(ix, 5) << ',' << workspace.up(ix, 6) << '\n'; + } +} + +int main(int argc, char** argv) +{ + const int nx = parse_nx(argc, argv); + const double delt = 5.0e-4; + const double tmax = 0.1; + + auto workspace = initialize(nx, Gamma, Bx, LeftPrimitive, RightPrimitive); + + mhd1d::evolve_ssp_rk3(workspace, delt, tmax); + + write_csv(workspace, std::cout); + + return 0; +} diff --git a/benchmarks/magnetohydrodynamics/shared/src/mhd1d.cpp b/benchmarks/magnetohydrodynamics/shared/src/mhd1d.cpp new file mode 100644 index 0000000..25264cf --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/src/mhd1d.cpp @@ -0,0 +1,235 @@ +#include "mhd1d.hpp" +#include "hlld.hpp" + +#include +#include +#include + +namespace mhd1d +{ + +double sign(const double x) +{ + return copysign(1.0, x); +} + +double mc2(double a, double b) +{ + return 0.5 * (sign(a) + sign(b)) * + std::min({2.0 * std::abs(a), 2.0 * std::abs(b), 0.5 * std::abs(a + b)}); +} + +void primitive_to_conservative(const double* primitive, double* conservative, double bx, + double gamma) +{ + const double rho = primitive[0]; + const double u = primitive[1]; + const double v = primitive[2]; + const double w = primitive[3]; + const double pressure = primitive[4]; + const double by = primitive[5]; + const double bz = primitive[6]; + + const double kinetic_energy = 0.5 * rho * (u * u + v * v + w * w); + const double magnetic_energy = 0.5 * (bx * bx + by * by + bz * bz); + + conservative[0] = rho; + conservative[1] = rho * u; + conservative[2] = rho * v; + conservative[3] = rho * w; + conservative[4] = pressure / (gamma - 1.0) + kinetic_energy + magnetic_energy; + conservative[5] = by; + conservative[6] = bz; +} + +void conservative_to_primitive(const double* conservative, double* primitive, double bx, + double gamma) +{ + const double rho = conservative[0]; + if (rho <= 0.0) { + throw std::runtime_error("density must be positive"); + } + + const double u = conservative[1] / rho; + const double v = conservative[2] / rho; + const double w = conservative[3] / rho; + const double by = conservative[5]; + const double bz = conservative[6]; + + const double kinetic_energy = 0.5 * rho * (u * u + v * v + w * w); + const double magnetic_energy = 0.5 * (bx * bx + by * by + bz * bz); + const double pressure = (gamma - 1.0) * (conservative[4] - kinetic_energy - magnetic_energy); + + primitive[0] = rho; + primitive[1] = u; + primitive[2] = v; + primitive[3] = w; + primitive[4] = pressure; + primitive[5] = by; + primitive[6] = bz; +} + +void convert_primitive_to_conservative(ArrayView2D primitive, ArrayView2D conservative, double bx, + double gamma) +{ + const int ix_min = 0; + const int ix_max = primitive.extent(0) - 1; + + for (int ix = ix_min; ix <= ix_max; ++ix) { + primitive_to_conservative(&primitive(ix, 0), &conservative(ix, 0), bx, gamma); + } +} + +void convert_conservative_to_primitive(ArrayView2D conservative, ArrayView2D primitive, double bx, + double gamma) +{ + const int ix_min = 0; + const int ix_max = conservative.extent(0) - 1; + + for (int ix = ix_min; ix <= ix_max; ++ix) { + conservative_to_primitive(&conservative(ix, 0), &primitive(ix, 0), bx, gamma); + } +} + +void set_boundary_lb(ArrayView2D dst, ArrayView2D src, int lbx) +{ + const int ix_min = 0; + + for (int ix = ix_min; ix < lbx; ++ix) { + for (int component = 0; component < N_Component; ++component) { + dst(ix, component) = src(lbx, component); + } + } +} + +void set_boundary_ub(ArrayView2D dst, ArrayView2D src, int ubx) +{ + const int ix_max = dst.extent(0) - 1; + + for (int ix = ubx + 1; ix <= ix_max; ++ix) { + for (int component = 0; component < N_Component; ++component) { + dst(ix, component) = src(ubx, component); + } + } +} + +void set_boundary(ArrayView2D dst, ArrayView2D src, int lbx, int ubx) +{ + set_boundary_lb(dst, src, lbx); + set_boundary_ub(dst, src, ubx); +} + +void compute_lr(SolverWorkspace& workspace) +{ + const ArrayView2D up = workspace.up; + const ArrayView2D up_l = workspace.up_l; + const ArrayView2D up_r = workspace.up_r; + + const int lbx = workspace.Lbx; + const int ubx = workspace.Ubx; + for (int ix = lbx; ix <= ubx; ++ix) { + for (int component = 0; component < N_Component; ++component) { + const double slope_l = up(ix, component) - up(ix - 1, component); + const double slope_r = up(ix + 1, component) - up(ix, component); + const double slope = mc2(slope_l, slope_r); + up_l(ix, component) = up(ix, component) + 0.5 * slope; + up_r(ix, component) = up(ix, component) - 0.5 * slope; + } + } + + set_boundary_lb(up_l, up_r, lbx); + set_boundary_ub(up_r, up_l, ubx); +} + +void compute_flux_hlld(SolverWorkspace& workspace, double bx, double gamma) +{ + const ArrayView2D up_l = workspace.up_l; + const ArrayView2D up_r = workspace.up_r; + const ArrayView2D flux = workspace.flux; + + const int lbx = workspace.Lbx; + const int ubx = workspace.Ubx; + for (int ix = lbx - 1; ix <= ubx; ++ix) { + ::hlld_flux_from_primitive(&up_l(ix, 0), &up_r(ix + 1, 0), bx, gamma, &flux(ix, 0)); + } +} + +void compute_rhs(SolverWorkspace& workspace) +{ + const ArrayView2D uc = workspace.uc; + const ArrayView2D up = workspace.up; + const ArrayView2D flux = workspace.flux; + const ArrayView2D rhs = workspace.rhs; + + set_boundary(uc, uc, workspace.Lbx, workspace.Ubx); + convert_conservative_to_primitive(uc, up, workspace.bx, workspace.gamma); + set_boundary(up, up, workspace.Lbx, workspace.Ubx); + compute_lr(workspace); + compute_flux_hlld(workspace, workspace.bx, workspace.gamma); + + const int lbx = workspace.Lbx; + const int ubx = workspace.Ubx; + for (int ix = lbx; ix <= ubx; ++ix) { + for (int component = 0; component < N_Component; ++component) { + rhs(ix, component) = -(flux(ix, component) - flux(ix - 1, component)) / workspace.dx; + } + } +} + +void copy(ArrayView2D source, ArrayView2D destination) +{ + const int nx = source.extent(0); + for (int ix = 0; ix < nx; ++ix) { + for (int component = 0; component < N_Component; ++component) { + destination(ix, component) = source(ix, component); + } + } +} + +void push_ssp_rk3(SolverWorkspace& workspace, double dt) +{ + constexpr double coeffs[3][3] = { + {1.0, 0.0, 1.0}, + {3.0 / 4.0, 1.0 / 4.0, 1.0 / 4.0}, + {1.0 / 3.0, 2.0 / 3.0, 2.0 / 3.0}, + }; + + const ArrayView2D prev = workspace.prev; + const ArrayView2D rhs = workspace.rhs; + + copy(workspace.uc, prev); + + for (int substep = 0; substep < 3; ++substep) { + compute_rhs(workspace); + + const double a = coeffs[substep][0]; + const double b = coeffs[substep][1]; + const double c = coeffs[substep][2]; + + const int lbx = workspace.Lbx; + const int ubx = workspace.Ubx; + for (int ix = lbx; ix <= ubx; ++ix) { + for (int component = 0; component < N_Component; ++component) { + workspace.uc(ix, component) = + a * prev(ix, component) + b * workspace.uc(ix, component) + c * dt * rhs(ix, component); + } + } + + set_boundary(workspace.uc, workspace.uc, workspace.Lbx, workspace.Ubx); + convert_conservative_to_primitive(workspace.uc, workspace.up, workspace.bx, workspace.gamma); + set_boundary(workspace.up, workspace.up, workspace.Lbx, workspace.Ubx); + } +} + +void evolve_ssp_rk3(SolverWorkspace& workspace, double dt, double t_final) +{ + double elapsed_time = 0.0; + while (elapsed_time < t_final) { + const double remaining_time = t_final - elapsed_time; + const double step_dt = std::min(dt, remaining_time); + push_ssp_rk3(workspace, step_dt); + elapsed_time = (step_dt < dt) ? t_final : (elapsed_time + step_dt); + } +} + +} // namespace mhd1d diff --git a/benchmarks/magnetohydrodynamics/shared/src/mhd1d.hpp b/benchmarks/magnetohydrodynamics/shared/src/mhd1d.hpp new file mode 100644 index 0000000..6925629 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/src/mhd1d.hpp @@ -0,0 +1,107 @@ +#pragma once + +#include +#include + +#include + +namespace mhd1d +{ + +namespace stdex = std::experimental; + +constexpr int N_Component = 7; +constexpr int N_margin = 1; + +using StateVector = std::array; +using ArrayView1D = stdex::mdspan, stdex::layout_right>; +using ArrayView2D = stdex::mdspan, stdex::layout_right>; + +struct SolverWorkspace { + explicit SolverWorkspace(int nx, double gamma, double bx) + : Nx(nx), Lbx(N_margin), Ubx(N_margin + nx - 1), dx(1.0 / static_cast(nx)), + gamma(gamma), bx(bx), storage(Nx + 2 * N_margin, N_Component) + { + init_views(Nx + 2 * N_margin, N_Component); + + for (int ix = Lbx; ix <= Ubx; ++ix) { + x(ix) = (static_cast(ix - Lbx) + 0.5) * dx; + } + } + + int Nx; + int Lbx; + int Ubx; + double dx; + double gamma; + double bx; + + ArrayView1D x; + ArrayView2D uc; + ArrayView2D up; + ArrayView2D up_l; + ArrayView2D up_r; + ArrayView2D rhs; + ArrayView2D prev; + ArrayView2D flux; + +private: + void init_views(int n_grid, int n_component) + { + x = ArrayView1D(storage.x.data(), n_grid); + uc = ArrayView2D(storage.uc.data(), n_grid, n_component); + up = ArrayView2D(storage.up.data(), n_grid, n_component); + up_l = ArrayView2D(storage.up_l.data(), n_grid, n_component); + up_r = ArrayView2D(storage.up_r.data(), n_grid, n_component); + rhs = ArrayView2D(storage.rhs.data(), n_grid, n_component); + prev = ArrayView2D(storage.prev.data(), n_grid, n_component); + flux = ArrayView2D(storage.flux.data(), n_grid, n_component); + } + + struct Storage { + explicit Storage(int n_grid, int n_component) + : x(n_grid), uc(n_grid * n_component), up(n_grid * n_component), up_l(n_grid * n_component), + up_r(n_grid * n_component), rhs(n_grid * n_component), prev(n_grid * n_component), + flux(n_grid * n_component) + { + } + + std::vector x; + std::vector uc; + std::vector up; + std::vector up_l; + std::vector up_r; + std::vector rhs; + std::vector prev; + std::vector flux; + }; + + Storage storage; +}; + +void primitive_to_conservative(const double* primitive, double* conservative, double bx, + double gamma); + +void conservative_to_primitive(const double* conservative, double* primitive, double bx, + double gamma); + +void convert_primitive_to_conservative(ArrayView2D primitive_cells, ArrayView2D conservative_cells, + double bx, double gamma); + +void set_boundary_lb(ArrayView2D dst, ArrayView2D src, int lbx); + +void set_boundary_ub(ArrayView2D dst, ArrayView2D src, int ubx); + +void set_boundary(ArrayView2D dst, ArrayView2D src, int lbx, int ubx); + +void compute_lr(SolverWorkspace& workspace); + +void compute_flux_hlld(SolverWorkspace& workspace, double bx, double gamma); + +void compute_rhs(SolverWorkspace& workspace); + +void push_ssp_rk3(SolverWorkspace& workspace, double dt); + +void evolve_ssp_rk3(SolverWorkspace& workspace, double dt, double t_final); + +} // namespace mhd1d diff --git a/benchmarks/magnetohydrodynamics/shared/tests/test_reference.py b/benchmarks/magnetohydrodynamics/shared/tests/test_reference.py new file mode 100644 index 0000000..5e0426e --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/tests/test_reference.py @@ -0,0 +1,56 @@ +from __future__ import annotations + +import csv +import os +import sys +import subprocess +from pathlib import Path + +REPO_ROOT = Path(__file__).resolve().parents[4] +if str(REPO_ROOT) not in sys.path: + sys.path.insert(0, str(REPO_ROOT)) + +from benchmarks.magnetohydrodynamics.shared.eval.mhd1d_shared import ( + TOLERANCE, + assert_csv_rows_close, +) + + +SHARED_ROOT = Path(__file__).resolve().parents[1] +FIXTURE_CSV_PATH = SHARED_ROOT / "eval" / "fixtures" / "mhd1d" / "brio_wu_reference.csv" + + +def _build_reference_solver() -> Path: + build_dir = SHARED_ROOT / "build" + subprocess.run(["cmake", "-S", str(SHARED_ROOT), "-B", str(build_dir)], check=True) + subprocess.run( + ["cmake", "--build", str(build_dir), "--target", "mhd1d_reference"], + check=True, + ) + + binary_name = "mhd1d_reference.exe" if os.name == "nt" else "mhd1d_reference" + binary_path = build_dir / "bin" / binary_name + assert binary_path.exists() + return binary_path + + +def test_shared_reference_solver_matches_fixture(tmp_path: Path) -> None: + solver_path = _build_reference_solver() + output_csv_path = tmp_path / "solution.csv" + + completed = subprocess.run( + [str(solver_path), "200"], + check=True, + capture_output=True, + text=True, + ) + output_csv_path.write_text(completed.stdout, encoding="utf-8") + + output_rows = list( + csv.reader(output_csv_path.read_text(encoding="utf-8").splitlines()) + ) + reference_rows = list( + csv.reader(FIXTURE_CSV_PATH.read_text(encoding="utf-8").splitlines()) + ) + + assert_csv_rows_close(output_rows, reference_rows, tolerance=TOLERANCE) diff --git a/benchmarks/magnetohydrodynamics/shared/workspace/hlld.md b/benchmarks/magnetohydrodynamics/shared/workspace/hlld.md index d6822fa..5633e0a 100644 --- a/benchmarks/magnetohydrodynamics/shared/workspace/hlld.md +++ b/benchmarks/magnetohydrodynamics/shared/workspace/hlld.md @@ -394,31 +394,3 @@ to the location of zero in the wave fan: All fluxes use the conservative component ordering defined in `basic_equations.md`. - -## Implementation notes - -- This benchmark follows one specific HLLD implementation convention rather than - an arbitrary mathematically equivalent variant. - -- Small starred-state denominator - If $|D_\alpha|$ is extremely small, do not apply the raw starred-state update - by dividing through that value. Instead, replace the starred transverse - updates with: - -```math -v_\alpha^\ast = v_\alpha, \quad -w_\alpha^\ast = w_\alpha, \quad -B_{y,\alpha}^\ast = B_{y,\alpha}, \quad -B_{z,\alpha}^\ast = B_{z,\alpha}. -``` - -- Small $B_x$ - When $B_x=0$, the rotational waves collapse and the double-star regions become - unnecessary. In that case, do not use the double-star states for flux - calculation. - - For this benchmark, a merely small nonzero $|B_x|$ should still be treated as - a nondegenerate case unless some other guarded quantity, such as $D_\alpha$, - becomes numerically singular. - -- Assume all benchmark inputs are admissible physical states. diff --git a/benchmarks/magnetohydrodynamics/shared/workspace/plot_solution.py b/benchmarks/magnetohydrodynamics/shared/workspace/plot_solution.py new file mode 100644 index 0000000..f6cbaf4 --- /dev/null +++ b/benchmarks/magnetohydrodynamics/shared/workspace/plot_solution.py @@ -0,0 +1,122 @@ +#!/usr/bin/env python3 +"""Plot a Brio-Wu solver CSV and save an image. + +Usage: + python plot_solution.py [path/to/solution.csv] [path/to/output.png] + +If no path is provided, the script looks for ``solution.csv`` in the current +working directory and writes ``solution.png`` there. +""" + +from __future__ import annotations + +import argparse +import csv +from pathlib import Path +import sys + +import matplotlib + +matplotlib.use("Agg") +import matplotlib.pyplot as plt + + +EXPECTED_FIELDS = ["x", "rho", "u", "v", "w", "p", "by", "bz"] +DEFAULT_CSV = Path("solution.csv") +DEFAULT_OUTPUT = Path("solution.png") + + +def parse_args() -> argparse.Namespace: + parser = argparse.ArgumentParser( + description="Plot the Brio-Wu solver profiles from a CSV file.", + ) + parser.add_argument( + "csv_path", + nargs="?", + type=Path, + default=DEFAULT_CSV, + help="CSV file to plot (defaults to solution.csv in the current directory).", + ) + parser.add_argument( + "output_path", + nargs="?", + type=Path, + default=DEFAULT_OUTPUT, + help="Output image path (defaults to solution.png in the current directory).", + ) + return parser.parse_args() + + +def load_columns(csv_path: Path) -> dict[str, list[float]]: + with csv_path.open(newline="", encoding="utf-8") as handle: + reader = csv.reader(handle) + try: + first_row = next(reader) + except StopIteration as exc: + raise ValueError(f"{csv_path} is empty") from exc + + if first_row == EXPECTED_FIELDS: + rows = reader + else: + if len(first_row) != len(EXPECTED_FIELDS): + raise ValueError( + f"expected {len(EXPECTED_FIELDS)} CSV columns; got {len(first_row)}" + ) + rows = [first_row, *reader] + + columns: dict[str, list[float]] = {field: [] for field in EXPECTED_FIELDS} + for raw_row in rows: + if len(raw_row) != len(EXPECTED_FIELDS): + raise ValueError( + f"expected {len(EXPECTED_FIELDS)} CSV columns; got {len(raw_row)}" + ) + for field, raw_value in zip(EXPECTED_FIELDS, raw_row, strict=True): + columns[field].append(float(raw_value)) + + return columns + + +def main() -> int: + args = parse_args() + csv_path = args.csv_path + output_path = args.output_path + + if not csv_path.is_file(): + print(f"error: CSV file not found: {csv_path}", file=sys.stderr) + return 1 + + try: + columns = load_columns(csv_path) + except (OSError, ValueError, KeyError) as exc: + print(f"error: {exc}", file=sys.stderr) + return 1 + + x_values = columns["x"] + + fig, axes = plt.subplots(2, 2, figsize=(10, 7), sharex=True) + plots = [ + (axes[0, 0], "rho", r"$\rho$"), + (axes[0, 1], "u", r"$u$"), + (axes[1, 0], "p", r"$p$"), + (axes[1, 1], "by", r"$B_y$"), + ] + + for axis, field, ylabel in plots: + axis.plot(x_values, columns[field], linewidth=1.5) + axis.set_ylabel(ylabel) + axis.grid(True, alpha=0.3) + + for axis in axes[1, :]: + axis.set_xlabel(r"$x$") + + fig.suptitle(r"Brio-Wu Problem at $t = 0.1$") + fig.tight_layout() + output_path.parent.mkdir(parents=True, exist_ok=True) + fig.savefig(output_path, dpi=150) + plt.close(fig) + print(output_path) + return 0 + + +if __name__ == "__main__": + raise SystemExit(main()) diff --git a/docker/Dockerfile b/docker/Dockerfile index 634b404..47f5d67 100644 --- a/docker/Dockerfile +++ b/docker/Dockerfile @@ -3,6 +3,7 @@ FROM debian:bookworm-slim ENV DEBIAN_FRONTEND=noninteractive ARG FLAP_VERSION=1.2.16 +ARG OPENCODE_VERSION=1.3.3 ARG CATCH2_REF=v3.13.0 ARG MDSPAN_REF=mdspan-0.6.0 ARG XTL_REF=0.8.2 @@ -97,7 +98,7 @@ RUN apt-get update \ && rm -rf /tmp/flap-build "/tmp/FLAP-v${FLAP_VERSION}" \ && curl -LsSf https://astral.sh/uv/install.sh | env UV_UNMANAGED_INSTALL=/usr/local/bin sh \ && uv pip install --system --break-system-packages --no-cache ruff fprettify \ - && npm install -g opencode-ai @anthropic-ai/claude-code @openai/codex @github/copilot \ + && npm install -g opencode-ai@${OPENCODE_VERSION} @anthropic-ai/claude-code @openai/codex @github/copilot \ && uv --version WORKDIR /work diff --git a/runner/execution_agent.py b/runner/execution_agent.py index 0b0ee9c..aa86bec 100644 --- a/runner/execution_agent.py +++ b/runner/execution_agent.py @@ -126,7 +126,7 @@ def _run_agent_in_docker( opencode_state_dir.mkdir(parents=True, exist_ok=True) docker_cmd += [ "-e", - "XDG_DATA_HOME=/opencode-data", + "HOME=/opencode-data", "-v", f"{str(opencode_state_dir)}:/opencode-data:rw", ] diff --git a/runner/metrics_helpers.py b/runner/metrics_helpers.py index e28c263..0ec9090 100644 --- a/runner/metrics_helpers.py +++ b/runner/metrics_helpers.py @@ -143,7 +143,7 @@ def _opencode_state_dir(run_dir: Path) -> Path: def _collect_opencode_usage_metrics(*, state_dir: Path) -> dict[str, Any]: env = dict(os.environ) - env["XDG_DATA_HOME"] = str(state_dir) + env["HOME"] = str(state_dir) cmd = ["opencode", "stats", "--models", "1"] try: proc = subprocess.run( diff --git a/scripts/build_image.py b/scripts/build_image.py index 6cf0d3a..a65d9a3 100755 --- a/scripts/build_image.py +++ b/scripts/build_image.py @@ -11,6 +11,7 @@ DOCKERFILE = REPO_ROOT / "docker" / "Dockerfile" FLAP_VERSION = "1.2.16" +OPENCODE_VERSION = "1.3.3" CATCH2_REF = "v3.13.0" MDSPAN_REF = "mdspan-0.6.0" XTL_REF = "0.8.2" @@ -45,6 +46,8 @@ def main(argv: list[str]) -> int: "--build-arg", f"FLAP_VERSION={FLAP_VERSION}", "--build-arg", + f"OPENCODE_VERSION={OPENCODE_VERSION}", + "--build-arg", f"CATCH2_REF={CATCH2_REF}", "--build-arg", f"MDSPAN_REF={MDSPAN_REF}", diff --git a/tests/test_runner_cli_flow.py b/tests/test_runner_cli_flow.py index 39dfbca..b397844 100644 --- a/tests/test_runner_cli_flow.py +++ b/tests/test_runner_cli_flow.py @@ -1814,7 +1814,9 @@ def fake_eval(*args, **kwargs): def fake_subprocess_run(cmd, **kwargs): if cmd == ["opencode", "stats", "--models", "1"]: - self.assertIn("XDG_DATA_HOME", kwargs["env"]) + self.assertTrue( + str(kwargs["env"].get("HOME", "")).endswith("/.opencode-data") + ) return subprocess.CompletedProcess( cmd, 0, stdout=stats_output, stderr="" ) diff --git a/uv.lock b/uv.lock new file mode 100644 index 0000000..99ffc7b --- /dev/null +++ b/uv.lock @@ -0,0 +1,123 @@ +version = 1 +revision = 3 +requires-python = ">=3.10" + +[[package]] +name = "configargparse" +version = "1.7.5" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/3f/0b/30328302903c55218ffc5199646d0e9d28348ff26c02ba77b2ffc58d294a/configargparse-1.7.5.tar.gz", hash = "sha256:e3f9a7bb6be34d66b2e3c4a2f58e3045f8dfae47b0dc039f87bcfaa0f193fb0f", size = 53548, upload-time = "2026-03-11T02:19:38.144Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/fe/19/3ba5e1b0bcc7b91aeab6c258afd70e4907d220fed3972febe38feb40db30/configargparse-1.7.5-py3-none-any.whl", hash = "sha256:1e63fdffedf94da9cd435fc13a1cd24777e76879dd2343912c1f871d4ac8c592", size = 27692, upload-time = "2026-03-11T02:19:36.442Z" }, +] + +[[package]] +name = "fprettify" +version = "0.3.7" +source = { registry = "https://pypi.org/simple" } +dependencies = [ + { name = "configargparse" }, +] +sdist = { url = "https://files.pythonhosted.org/packages/39/15/d88681bd2be4a375a78b52443b8e87608240913623d9be5c47e3c328b068/fprettify-0.3.7.tar.gz", hash = "sha256:1488a813f7e60a9e86c56fd0b82bd9df1b75bfb4bf2ee8e433c12f63b7e54057", size = 29639, upload-time = "2020-11-20T15:52:49.16Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/52/13/2c32d63574e116f8c933f56315df9135bf2fae7a88e9e7c6c4d37f48f4ef/fprettify-0.3.7-py3-none-any.whl", hash = "sha256:56f0a64c43dc47134ce32af2e5da8cd7a1584897be29d19289ec5d87510d1daf", size = 28095, upload-time = "2020-11-20T15:52:47.719Z" }, +] + +[[package]] +name = "ruff" +version = "0.15.8" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/14/b0/73cf7550861e2b4824950b8b52eebdcc5adc792a00c514406556c5b80817/ruff-0.15.8.tar.gz", hash = "sha256:995f11f63597ee362130d1d5a327a87cb6f3f5eae3094c620bcc632329a4d26e", size = 4610921, upload-time = "2026-03-26T18:39:38.675Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/4a/92/c445b0cd6da6e7ae51e954939cb69f97e008dbe750cfca89b8cedc081be7/ruff-0.15.8-py3-none-linux_armv6l.whl", hash = "sha256:cbe05adeba76d58162762d6b239c9056f1a15a55bd4b346cfd21e26cd6ad7bc7", size = 10527394, upload-time = "2026-03-26T18:39:41.566Z" }, + { url = "https://files.pythonhosted.org/packages/eb/92/f1c662784d149ad1414cae450b082cf736430c12ca78367f20f5ed569d65/ruff-0.15.8-py3-none-macosx_10_12_x86_64.whl", hash = "sha256:d3e3d0b6ba8dca1b7ef9ab80a28e840a20070c4b62e56d675c24f366ef330570", size = 10905693, upload-time = "2026-03-26T18:39:30.364Z" }, + { url = "https://files.pythonhosted.org/packages/ca/f2/7a631a8af6d88bcef997eb1bf87cc3da158294c57044aafd3e17030613de/ruff-0.15.8-py3-none-macosx_11_0_arm64.whl", hash = "sha256:6ee3ae5c65a42f273f126686353f2e08ff29927b7b7e203b711514370d500de3", size = 10323044, upload-time = "2026-03-26T18:39:33.37Z" }, + { url = "https://files.pythonhosted.org/packages/67/18/1bf38e20914a05e72ef3b9569b1d5c70a7ef26cd188d69e9ca8ef588d5bf/ruff-0.15.8-py3-none-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:fdce027ada77baa448077ccc6ebb2fa9c3c62fd110d8659d601cf2f475858d94", size = 10629135, upload-time = "2026-03-26T18:39:44.142Z" }, + { url = "https://files.pythonhosted.org/packages/d2/e9/138c150ff9af60556121623d41aba18b7b57d95ac032e177b6a53789d279/ruff-0.15.8-py3-none-manylinux_2_17_armv7l.manylinux2014_armv7l.whl", hash = "sha256:12e617fc01a95e5821648a6df341d80456bd627bfab8a829f7cfc26a14a4b4a3", size = 10348041, upload-time = "2026-03-26T18:39:52.178Z" }, + { url = "https://files.pythonhosted.org/packages/02/f1/5bfb9298d9c323f842c5ddeb85f1f10ef51516ac7a34ba446c9347d898df/ruff-0.15.8-py3-none-manylinux_2_17_i686.manylinux2014_i686.whl", hash = "sha256:432701303b26416d22ba696c39f2c6f12499b89093b61360abc34bcc9bf07762", size = 11121987, upload-time = "2026-03-26T18:39:55.195Z" }, + { url = "https://files.pythonhosted.org/packages/10/11/6da2e538704e753c04e8d86b1fc55712fdbdcc266af1a1ece7a51fff0d10/ruff-0.15.8-py3-none-manylinux_2_17_ppc64le.manylinux2014_ppc64le.whl", hash = "sha256:d910ae974b7a06a33a057cb87d2a10792a3b2b3b35e33d2699fdf63ec8f6b17a", size = 11951057, upload-time = "2026-03-26T18:39:19.18Z" }, + { url = "https://files.pythonhosted.org/packages/83/f0/c9208c5fd5101bf87002fed774ff25a96eea313d305f1e5d5744698dc314/ruff-0.15.8-py3-none-manylinux_2_17_s390x.manylinux2014_s390x.whl", hash = "sha256:2033f963c43949d51e6fdccd3946633c6b37c484f5f98c3035f49c27395a8ab8", size = 11464613, upload-time = "2026-03-26T18:40:06.301Z" }, + { url = "https://files.pythonhosted.org/packages/f8/22/d7f2fabdba4fae9f3b570e5605d5eb4500dcb7b770d3217dca4428484b17/ruff-0.15.8-py3-none-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:0f29b989a55572fb885b77464cf24af05500806ab4edf9a0fd8977f9759d85b1", size = 11257557, upload-time = "2026-03-26T18:39:57.972Z" }, + { url = "https://files.pythonhosted.org/packages/71/8c/382a9620038cf6906446b23ce8632ab8c0811b8f9d3e764f58bedd0c9a6f/ruff-0.15.8-py3-none-manylinux_2_31_riscv64.whl", hash = "sha256:ac51d486bf457cdc985a412fb1801b2dfd1bd8838372fc55de64b1510eff4bec", size = 11169440, upload-time = "2026-03-26T18:39:22.205Z" }, + { url = "https://files.pythonhosted.org/packages/4d/0d/0994c802a7eaaf99380085e4e40c845f8e32a562e20a38ec06174b52ef24/ruff-0.15.8-py3-none-musllinux_1_2_aarch64.whl", hash = "sha256:c9861eb959edab053c10ad62c278835ee69ca527b6dcd72b47d5c1e5648964f6", size = 10605963, upload-time = "2026-03-26T18:39:46.682Z" }, + { url = "https://files.pythonhosted.org/packages/19/aa/d624b86f5b0aad7cef6bbf9cd47a6a02dfdc4f72c92a337d724e39c9d14b/ruff-0.15.8-py3-none-musllinux_1_2_armv7l.whl", hash = "sha256:8d9a5b8ea13f26ae90838afc33f91b547e61b794865374f114f349e9036835fb", size = 10357484, upload-time = "2026-03-26T18:39:49.176Z" }, + { url = "https://files.pythonhosted.org/packages/35/c3/e0b7835d23001f7d999f3895c6b569927c4d39912286897f625736e1fd04/ruff-0.15.8-py3-none-musllinux_1_2_i686.whl", hash = "sha256:c2a33a529fb3cbc23a7124b5c6ff121e4d6228029cba374777bd7649cc8598b8", size = 10830426, upload-time = "2026-03-26T18:40:03.702Z" }, + { url = "https://files.pythonhosted.org/packages/f0/51/ab20b322f637b369383adc341d761eaaa0f0203d6b9a7421cd6e783d81b9/ruff-0.15.8-py3-none-musllinux_1_2_x86_64.whl", hash = "sha256:75e5cd06b1cf3f47a3996cfc999226b19aa92e7cce682dcd62f80d7035f98f49", size = 11345125, upload-time = "2026-03-26T18:39:27.799Z" }, + { url = "https://files.pythonhosted.org/packages/37/e6/90b2b33419f59d0f2c4c8a48a4b74b460709a557e8e0064cf33ad894f983/ruff-0.15.8-py3-none-win32.whl", hash = "sha256:bc1f0a51254ba21767bfa9a8b5013ca8149dcf38092e6a9eb704d876de94dc34", size = 10571959, upload-time = "2026-03-26T18:39:36.117Z" }, + { url = "https://files.pythonhosted.org/packages/1f/a2/ef467cb77099062317154c63f234b8a7baf7cb690b99af760c5b68b9ee7f/ruff-0.15.8-py3-none-win_amd64.whl", hash = "sha256:04f79eff02a72db209d47d665ba7ebcad609d8918a134f86cb13dd132159fc89", size = 11743893, upload-time = "2026-03-26T18:39:25.01Z" }, + { url = "https://files.pythonhosted.org/packages/15/e2/77be4fff062fa78d9b2a4dea85d14785dac5f1d0c1fb58ed52331f0ebe28/ruff-0.15.8-py3-none-win_arm64.whl", hash = "sha256:cf891fa8e3bb430c0e7fac93851a5978fc99c8fa2c053b57b118972866f8e5f2", size = 11048175, upload-time = "2026-03-26T18:40:01.06Z" }, +] + +[[package]] +name = "simbench" +version = "0.0.0" +source = { virtual = "." } + +[package.optional-dependencies] +dev = [ + { name = "fprettify" }, + { name = "ruff" }, + { name = "tomli", marker = "python_full_version < '3.11'" }, +] + +[package.metadata] +requires-dist = [ + { name = "fprettify", marker = "extra == 'dev'", specifier = ">=0.3.7" }, + { name = "ruff", marker = "extra == 'dev'", specifier = ">=0.11" }, + { name = "tomli", marker = "python_full_version < '3.11' and extra == 'dev'", specifier = ">=2.0" }, +] +provides-extras = ["dev"] + +[[package]] +name = "tomli" +version = "2.4.1" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/22/de/48c59722572767841493b26183a0d1cc411d54fd759c5607c4590b6563a6/tomli-2.4.1.tar.gz", hash = "sha256:7c7e1a961a0b2f2472c1ac5b69affa0ae1132c39adcb67aba98568702b9cc23f", size = 17543, upload-time = "2026-03-25T20:22:03.828Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/f4/11/db3d5885d8528263d8adc260bb2d28ebf1270b96e98f0e0268d32b8d9900/tomli-2.4.1-cp311-cp311-macosx_10_9_x86_64.whl", hash = "sha256:f8f0fc26ec2cc2b965b7a3b87cd19c5c6b8c5e5f436b984e85f486d652285c30", size = 154704, upload-time = "2026-03-25T20:21:10.473Z" }, + { url = "https://files.pythonhosted.org/packages/6d/f7/675db52c7e46064a9aa928885a9b20f4124ecb9bc2e1ce74c9106648d202/tomli-2.4.1-cp311-cp311-macosx_11_0_arm64.whl", hash = "sha256:4ab97e64ccda8756376892c53a72bd1f964e519c77236368527f758fbc36a53a", size = 149454, upload-time = "2026-03-25T20:21:12.036Z" }, + { url = "https://files.pythonhosted.org/packages/61/71/81c50943cf953efa35bce7646caab3cf457a7d8c030b27cfb40d7235f9ee/tomli-2.4.1-cp311-cp311-manylinux2014_aarch64.manylinux_2_17_aarch64.manylinux_2_28_aarch64.whl", hash = "sha256:96481a5786729fd470164b47cdb3e0e58062a496f455ee41b4403be77cb5a076", size = 237561, upload-time = "2026-03-25T20:21:13.098Z" }, + { url = "https://files.pythonhosted.org/packages/48/c1/f41d9cb618acccca7df82aaf682f9b49013c9397212cb9f53219e3abac37/tomli-2.4.1-cp311-cp311-manylinux2014_x86_64.manylinux_2_17_x86_64.manylinux_2_28_x86_64.whl", hash = "sha256:5a881ab208c0baf688221f8cecc5401bd291d67e38a1ac884d6736cbcd8247e9", size = 243824, upload-time = "2026-03-25T20:21:14.569Z" }, + { url = "https://files.pythonhosted.org/packages/22/e4/5a816ecdd1f8ca51fb756ef684b90f2780afc52fc67f987e3c61d800a46d/tomli-2.4.1-cp311-cp311-musllinux_1_2_aarch64.whl", hash = "sha256:47149d5bd38761ac8be13a84864bf0b7b70bc051806bc3669ab1cbc56216b23c", size = 242227, upload-time = "2026-03-25T20:21:15.712Z" }, + { url = "https://files.pythonhosted.org/packages/6b/49/2b2a0ef529aa6eec245d25f0c703e020a73955ad7edf73e7f54ddc608aa5/tomli-2.4.1-cp311-cp311-musllinux_1_2_x86_64.whl", hash = "sha256:ec9bfaf3ad2df51ace80688143a6a4ebc09a248f6ff781a9945e51937008fcbc", size = 247859, upload-time = "2026-03-25T20:21:17.001Z" }, + { url = "https://files.pythonhosted.org/packages/83/bd/6c1a630eaca337e1e78c5903104f831bda934c426f9231429396ce3c3467/tomli-2.4.1-cp311-cp311-win32.whl", hash = "sha256:ff2983983d34813c1aeb0fa89091e76c3a22889ee83ab27c5eeb45100560c049", size = 97204, upload-time = "2026-03-25T20:21:18.079Z" }, + { url = "https://files.pythonhosted.org/packages/42/59/71461df1a885647e10b6bb7802d0b8e66480c61f3f43079e0dcd315b3954/tomli-2.4.1-cp311-cp311-win_amd64.whl", hash = "sha256:5ee18d9ebdb417e384b58fe414e8d6af9f4e7a0ae761519fb50f721de398dd4e", size = 108084, upload-time = "2026-03-25T20:21:18.978Z" }, + { url = "https://files.pythonhosted.org/packages/b8/83/dceca96142499c069475b790e7913b1044c1a4337e700751f48ed723f883/tomli-2.4.1-cp311-cp311-win_arm64.whl", hash = "sha256:c2541745709bad0264b7d4705ad453b76ccd191e64aa6f0fc66b69a293a45ece", size = 95285, upload-time = "2026-03-25T20:21:20.309Z" }, + { url = "https://files.pythonhosted.org/packages/c1/ba/42f134a3fe2b370f555f44b1d72feebb94debcab01676bf918d0cb70e9aa/tomli-2.4.1-cp312-cp312-macosx_10_13_x86_64.whl", hash = "sha256:c742f741d58a28940ce01d58f0ab2ea3ced8b12402f162f4d534dfe18ba1cd6a", size = 155924, upload-time = "2026-03-25T20:21:21.626Z" }, + { url = "https://files.pythonhosted.org/packages/dc/c7/62d7a17c26487ade21c5422b646110f2162f1fcc95980ef7f63e73c68f14/tomli-2.4.1-cp312-cp312-macosx_11_0_arm64.whl", hash = "sha256:7f86fd587c4ed9dd76f318225e7d9b29cfc5a9d43de44e5754db8d1128487085", size = 150018, upload-time = "2026-03-25T20:21:23.002Z" }, + { url = "https://files.pythonhosted.org/packages/5c/05/79d13d7c15f13bdef410bdd49a6485b1c37d28968314eabee452c22a7fda/tomli-2.4.1-cp312-cp312-manylinux2014_aarch64.manylinux_2_17_aarch64.manylinux_2_28_aarch64.whl", hash = "sha256:ff18e6a727ee0ab0388507b89d1bc6a22b138d1e2fa56d1ad494586d61d2eae9", size = 244948, upload-time = "2026-03-25T20:21:24.04Z" }, + { url = "https://files.pythonhosted.org/packages/10/90/d62ce007a1c80d0b2c93e02cab211224756240884751b94ca72df8a875ca/tomli-2.4.1-cp312-cp312-manylinux2014_x86_64.manylinux_2_17_x86_64.manylinux_2_28_x86_64.whl", hash = "sha256:136443dbd7e1dee43c68ac2694fde36b2849865fa258d39bf822c10e8068eac5", size = 253341, upload-time = "2026-03-25T20:21:25.177Z" }, + { url = "https://files.pythonhosted.org/packages/1a/7e/caf6496d60152ad4ed09282c1885cca4eea150bfd007da84aea07bcc0a3e/tomli-2.4.1-cp312-cp312-musllinux_1_2_aarch64.whl", hash = "sha256:5e262d41726bc187e69af7825504c933b6794dc3fbd5945e41a79bb14c31f585", size = 248159, upload-time = "2026-03-25T20:21:26.364Z" }, + { url = "https://files.pythonhosted.org/packages/99/e7/c6f69c3120de34bbd882c6fba7975f3d7a746e9218e56ab46a1bc4b42552/tomli-2.4.1-cp312-cp312-musllinux_1_2_x86_64.whl", hash = "sha256:5cb41aa38891e073ee49d55fbc7839cfdb2bc0e600add13874d048c94aadddd1", size = 253290, upload-time = "2026-03-25T20:21:27.46Z" }, + { url = "https://files.pythonhosted.org/packages/d6/2f/4a3c322f22c5c66c4b836ec58211641a4067364f5dcdd7b974b4c5da300c/tomli-2.4.1-cp312-cp312-win32.whl", hash = "sha256:da25dc3563bff5965356133435b757a795a17b17d01dbc0f42fb32447ddfd917", size = 98141, upload-time = "2026-03-25T20:21:28.492Z" }, + { url = "https://files.pythonhosted.org/packages/24/22/4daacd05391b92c55759d55eaee21e1dfaea86ce5c571f10083360adf534/tomli-2.4.1-cp312-cp312-win_amd64.whl", hash = "sha256:52c8ef851d9a240f11a88c003eacb03c31fc1c9c4ec64a99a0f922b93874fda9", size = 108847, upload-time = "2026-03-25T20:21:29.386Z" }, + { url = "https://files.pythonhosted.org/packages/68/fd/70e768887666ddd9e9f5d85129e84910f2db2796f9096aa02b721a53098d/tomli-2.4.1-cp312-cp312-win_arm64.whl", hash = "sha256:f758f1b9299d059cc3f6546ae2af89670cb1c4d48ea29c3cacc4fe7de3058257", size = 95088, upload-time = "2026-03-25T20:21:30.677Z" }, + { url = "https://files.pythonhosted.org/packages/07/06/b823a7e818c756d9a7123ba2cda7d07bc2dd32835648d1a7b7b7a05d848d/tomli-2.4.1-cp313-cp313-macosx_10_13_x86_64.whl", hash = "sha256:36d2bd2ad5fb9eaddba5226aa02c8ec3fa4f192631e347b3ed28186d43be6b54", size = 155866, upload-time = "2026-03-25T20:21:31.65Z" }, + { url = "https://files.pythonhosted.org/packages/14/6f/12645cf7f08e1a20c7eb8c297c6f11d31c1b50f316a7e7e1e1de6e2e7b7e/tomli-2.4.1-cp313-cp313-macosx_11_0_arm64.whl", hash = "sha256:eb0dc4e38e6a1fd579e5d50369aa2e10acfc9cace504579b2faabb478e76941a", size = 149887, upload-time = "2026-03-25T20:21:33.028Z" }, + { url = "https://files.pythonhosted.org/packages/5c/e0/90637574e5e7212c09099c67ad349b04ec4d6020324539297b634a0192b0/tomli-2.4.1-cp313-cp313-manylinux2014_aarch64.manylinux_2_17_aarch64.manylinux_2_28_aarch64.whl", hash = "sha256:c7f2c7f2b9ca6bdeef8f0fa897f8e05085923eb091721675170254cbc5b02897", size = 243704, upload-time = "2026-03-25T20:21:34.51Z" }, + { url = "https://files.pythonhosted.org/packages/10/8f/d3ddb16c5a4befdf31a23307f72828686ab2096f068eaf56631e136c1fdd/tomli-2.4.1-cp313-cp313-manylinux2014_x86_64.manylinux_2_17_x86_64.manylinux_2_28_x86_64.whl", hash = "sha256:f3c6818a1a86dd6dca7ddcaaf76947d5ba31aecc28cb1b67009a5877c9a64f3f", size = 251628, upload-time = "2026-03-25T20:21:36.012Z" }, + { url = "https://files.pythonhosted.org/packages/e3/f1/dbeeb9116715abee2485bf0a12d07a8f31af94d71608c171c45f64c0469d/tomli-2.4.1-cp313-cp313-musllinux_1_2_aarch64.whl", hash = "sha256:d312ef37c91508b0ab2cee7da26ec0b3ed2f03ce12bd87a588d771ae15dcf82d", size = 247180, upload-time = "2026-03-25T20:21:37.136Z" }, + { url = "https://files.pythonhosted.org/packages/d3/74/16336ffd19ed4da28a70959f92f506233bd7cfc2332b20bdb01591e8b1d1/tomli-2.4.1-cp313-cp313-musllinux_1_2_x86_64.whl", hash = "sha256:51529d40e3ca50046d7606fa99ce3956a617f9b36380da3b7f0dd3dd28e68cb5", size = 251674, upload-time = "2026-03-25T20:21:38.298Z" }, + { url = "https://files.pythonhosted.org/packages/16/f9/229fa3434c590ddf6c0aa9af64d3af4b752540686cace29e6281e3458469/tomli-2.4.1-cp313-cp313-win32.whl", hash = "sha256:2190f2e9dd7508d2a90ded5ed369255980a1bcdd58e52f7fe24b8162bf9fedbd", size = 97976, upload-time = "2026-03-25T20:21:39.316Z" }, + { url = "https://files.pythonhosted.org/packages/6a/1e/71dfd96bcc1c775420cb8befe7a9d35f2e5b1309798f009dca17b7708c1e/tomli-2.4.1-cp313-cp313-win_amd64.whl", hash = "sha256:8d65a2fbf9d2f8352685bc1364177ee3923d6baf5e7f43ea4959d7d8bc326a36", size = 108755, upload-time = "2026-03-25T20:21:40.248Z" }, + { url = "https://files.pythonhosted.org/packages/83/7a/d34f422a021d62420b78f5c538e5b102f62bea616d1d75a13f0a88acb04a/tomli-2.4.1-cp313-cp313-win_arm64.whl", hash = "sha256:4b605484e43cdc43f0954ddae319fb75f04cc10dd80d830540060ee7cd0243cd", size = 95265, upload-time = "2026-03-25T20:21:41.219Z" }, + { url = "https://files.pythonhosted.org/packages/3c/fb/9a5c8d27dbab540869f7c1f8eb0abb3244189ce780ba9cd73f3770662072/tomli-2.4.1-cp314-cp314-macosx_10_15_x86_64.whl", hash = "sha256:fd0409a3653af6c147209d267a0e4243f0ae46b011aa978b1080359fddc9b6cf", size = 155726, upload-time = "2026-03-25T20:21:42.23Z" }, + { url = "https://files.pythonhosted.org/packages/62/05/d2f816630cc771ad836af54f5001f47a6f611d2d39535364f148b6a92d6b/tomli-2.4.1-cp314-cp314-macosx_11_0_arm64.whl", hash = "sha256:a120733b01c45e9a0c34aeef92bf0cf1d56cfe81ed9d47d562f9ed591a9828ac", size = 149859, upload-time = "2026-03-25T20:21:43.386Z" }, + { url = "https://files.pythonhosted.org/packages/ce/48/66341bdb858ad9bd0ceab5a86f90eddab127cf8b046418009f2125630ecb/tomli-2.4.1-cp314-cp314-manylinux2014_aarch64.manylinux_2_17_aarch64.manylinux_2_28_aarch64.whl", hash = "sha256:559db847dc486944896521f68d8190be1c9e719fced785720d2216fe7022b662", size = 244713, upload-time = "2026-03-25T20:21:44.474Z" }, + { url = "https://files.pythonhosted.org/packages/df/6d/c5fad00d82b3c7a3ab6189bd4b10e60466f22cfe8a08a9394185c8a8111c/tomli-2.4.1-cp314-cp314-manylinux2014_x86_64.manylinux_2_17_x86_64.manylinux_2_28_x86_64.whl", hash = "sha256:01f520d4f53ef97964a240a035ec2a869fe1a37dde002b57ebc4417a27ccd853", size = 252084, upload-time = "2026-03-25T20:21:45.62Z" }, + { url = "https://files.pythonhosted.org/packages/00/71/3a69e86f3eafe8c7a59d008d245888051005bd657760e96d5fbfb0b740c2/tomli-2.4.1-cp314-cp314-musllinux_1_2_aarch64.whl", hash = "sha256:7f94b27a62cfad8496c8d2513e1a222dd446f095fca8987fceef261225538a15", size = 247973, upload-time = "2026-03-25T20:21:46.937Z" }, + { url = "https://files.pythonhosted.org/packages/67/50/361e986652847fec4bd5e4a0208752fbe64689c603c7ae5ea7cb16b1c0ca/tomli-2.4.1-cp314-cp314-musllinux_1_2_x86_64.whl", hash = "sha256:ede3e6487c5ef5d28634ba3f31f989030ad6af71edfb0055cbbd14189ff240ba", size = 256223, upload-time = "2026-03-25T20:21:48.467Z" }, + { url = "https://files.pythonhosted.org/packages/8c/9a/b4173689a9203472e5467217e0154b00e260621caa227b6fa01feab16998/tomli-2.4.1-cp314-cp314-win32.whl", hash = "sha256:3d48a93ee1c9b79c04bb38772ee1b64dcf18ff43085896ea460ca8dec96f35f6", size = 98973, upload-time = "2026-03-25T20:21:49.526Z" }, + { url = "https://files.pythonhosted.org/packages/14/58/640ac93bf230cd27d002462c9af0d837779f8773bc03dee06b5835208214/tomli-2.4.1-cp314-cp314-win_amd64.whl", hash = "sha256:88dceee75c2c63af144e456745e10101eb67361050196b0b6af5d717254dddf7", size = 109082, upload-time = "2026-03-25T20:21:50.506Z" }, + { url = "https://files.pythonhosted.org/packages/d5/2f/702d5e05b227401c1068f0d386d79a589bb12bf64c3d2c72ce0631e3bc49/tomli-2.4.1-cp314-cp314-win_arm64.whl", hash = "sha256:b8c198f8c1805dc42708689ed6864951fd2494f924149d3e4bce7710f8eb5232", size = 96490, upload-time = "2026-03-25T20:21:51.474Z" }, + { url = "https://files.pythonhosted.org/packages/45/4b/b877b05c8ba62927d9865dd980e34a755de541eb65fffba52b4cc495d4d2/tomli-2.4.1-cp314-cp314t-macosx_10_15_x86_64.whl", hash = "sha256:d4d8fe59808a54658fcc0160ecfb1b30f9089906c50b23bcb4c69eddc19ec2b4", size = 164263, upload-time = "2026-03-25T20:21:52.543Z" }, + { url = "https://files.pythonhosted.org/packages/24/79/6ab420d37a270b89f7195dec5448f79400d9e9c1826df982f3f8e97b24fd/tomli-2.4.1-cp314-cp314t-macosx_11_0_arm64.whl", hash = "sha256:7008df2e7655c495dd12d2a4ad038ff878d4ca4b81fccaf82b714e07eae4402c", size = 160736, upload-time = "2026-03-25T20:21:53.674Z" }, + { url = "https://files.pythonhosted.org/packages/02/e0/3630057d8eb170310785723ed5adcdfb7d50cb7e6455f85ba8a3deed642b/tomli-2.4.1-cp314-cp314t-manylinux2014_aarch64.manylinux_2_17_aarch64.manylinux_2_28_aarch64.whl", hash = "sha256:1d8591993e228b0c930c4bb0db464bdad97b3289fb981255d6c9a41aedc84b2d", size = 270717, upload-time = "2026-03-25T20:21:55.129Z" }, + { url = "https://files.pythonhosted.org/packages/7a/b4/1613716072e544d1a7891f548d8f9ec6ce2faf42ca65acae01d76ea06bb0/tomli-2.4.1-cp314-cp314t-manylinux2014_x86_64.manylinux_2_17_x86_64.manylinux_2_28_x86_64.whl", hash = "sha256:734e20b57ba95624ecf1841e72b53f6e186355e216e5412de414e3c51e5e3c41", size = 278461, upload-time = "2026-03-25T20:21:56.228Z" }, + { url = "https://files.pythonhosted.org/packages/05/38/30f541baf6a3f6df77b3df16b01ba319221389e2da59427e221ef417ac0c/tomli-2.4.1-cp314-cp314t-musllinux_1_2_aarch64.whl", hash = "sha256:8a650c2dbafa08d42e51ba0b62740dae4ecb9338eefa093aa5c78ceb546fcd5c", size = 274855, upload-time = "2026-03-25T20:21:57.653Z" }, + { url = "https://files.pythonhosted.org/packages/77/a3/ec9dd4fd2c38e98de34223b995a3b34813e6bdadf86c75314c928350ed14/tomli-2.4.1-cp314-cp314t-musllinux_1_2_x86_64.whl", hash = "sha256:504aa796fe0569bb43171066009ead363de03675276d2d121ac1a4572397870f", size = 283144, upload-time = "2026-03-25T20:21:59.089Z" }, + { url = "https://files.pythonhosted.org/packages/ef/be/605a6261cac79fba2ec0c9827e986e00323a1945700969b8ee0b30d85453/tomli-2.4.1-cp314-cp314t-win32.whl", hash = "sha256:b1d22e6e9387bf4739fbe23bfa80e93f6b0373a7f1b96c6227c32bef95a4d7a8", size = 108683, upload-time = "2026-03-25T20:22:00.214Z" }, + { url = "https://files.pythonhosted.org/packages/12/64/da524626d3b9cc40c168a13da8335fe1c51be12c0a63685cc6db7308daae/tomli-2.4.1-cp314-cp314t-win_amd64.whl", hash = "sha256:2c1c351919aca02858f740c6d33adea0c5deea37f9ecca1cc1ef9e884a619d26", size = 121196, upload-time = "2026-03-25T20:22:01.169Z" }, + { url = "https://files.pythonhosted.org/packages/5a/cd/e80b62269fc78fc36c9af5a6b89c835baa8af28ff5ad28c7028d60860320/tomli-2.4.1-cp314-cp314t-win_arm64.whl", hash = "sha256:eab21f45c7f66c13f2a9e0e1535309cee140182a9cdae1e041d02e47291e8396", size = 100393, upload-time = "2026-03-25T20:22:02.137Z" }, + { url = "https://files.pythonhosted.org/packages/7b/61/cceae43728b7de99d9b847560c262873a1f6c98202171fd5ed62640b494b/tomli-2.4.1-py3-none-any.whl", hash = "sha256:0d85819802132122da43cb86656f8d1f8c6587d54ae7dcaf30e90533028b49fe", size = 14583, upload-time = "2026-03-25T20:22:03.012Z" }, +]