diff --git a/CHANGELOG.md b/CHANGELOG.md index b5faaaebc..d877f2179 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -10,6 +10,22 @@ The format is loosely based on [Keep a Changelog](https://keepachangelog.com/). ### Added +- **FP library: `fp_get_state` returns seven global scalars** after the + profile arrays: `TOTAL_IP` [MA], `STORED_ENERGY` [MJ], `COLLISION_POWER`, + `ABSORPTION_POWER`, `ABSORPTION_WR`, `ABSORPTION_WM` [MW] (species sums at + the latest `NTG1` sample, what FP otherwise only writes to unit 6) and + `PLASMA_VOLUME` [m^3]. They are what shows an EQ to FP handoff: the six + profiles are local quantities and stay bit-identical when only the + equilibrium changes (without radial transport, waves or current drive). + `fp_state_t` grows by 56 bytes at its end (layout version 2, #246), so + `libfpapi.so` and its callers are rebuilt together: `fplib` refuses a + library built before this change with `FplibAbiError` instead of reading + the scalars as zeros. `FpState.scalars`, `to_dict()["scalars"]` + and the `scalars` object of fp_mcp's state carry them, and an FP step of + the `totlib` pipeline keeps them beside `timefp` and its counts. They are + zero until the first run, and again after any successful parameter set + (`fp_set_param`, `fp_set_param_str`, or tot's `fp:` parameters, which the + registry counts since #253) until the next run that is made. - **EQ library: `eq_save`, `Eq.save(path)` and the `save` tool of `eq_mcp`** write the equilibrium of the last `eq_run` to a TASK/EQ binary file through `EQSAVE`, the file that `MODELG=3` with `KNAMEQ` loads. A save whose file diff --git a/docs/fp-library/architecture.md b/docs/fp-library/architecture.md index 23e094035..e680fc429 100644 --- a/docs/fp-library/architecture.md +++ b/docs/fp-library/architecture.md @@ -139,7 +139,17 @@ run is made without that species' collisions with its own kind). `nthmax`, `ntg2`), the simulation time `timefp` (double), and six 2-D profile arrays `RNT`, `RWT`, `RTT`, `RJT`, `RPCT`, `RPWT` of shape `[FP_MAX_NSAMAX][FP_MAX_NRMAX]` = `[8][100]`. The runtime copy uses -only the `[0:nsamax][0:nrmax]` slice. +only the `[0:nsamax][0:nrmax]` slice. After the arrays come seven +doubles: six volume-integrated globals summed over the species at the +latest `NTG1` sample, `TOTAL_IP` [MA], `STORED_ENERGY` [MJ], +`COLLISION_POWER`, `ABSORPTION_POWER`, `ABSORPTION_WR`, `ABSORPTION_WM` +[MW], and `PLASMA_VOLUME` [m^3] (`TVOLR`). They are zero until the first +run, and again after any successful parameter set (`fp_set_param`, +`fp_set_param_str`, or tot's `fp:` parameters: the registry counts every +set, #253) until the next run that is made, +while the profiles above are still those of the last run. The +struct grows at its end only, and `fp/fp_state.f90`, `fp/fp_api.h` and +`fplib._ffi.FpStateC` change together with a rebuilt `libfpapi.so`. ### Fortran backend @@ -173,7 +183,7 @@ only the `[0:nsamax][0:nrmax]` slice. to the bracket syntax. - `fplib.FpState` is a `@dataclass` snapshot; `to_dict()` is JSON-serialisable and uses a - `{"NRMAX": ..., "NSAMAX": ..., "profile": [{"NSA": i, "RNT": ...}, ...]}` + `{"NRMAX": ..., "NSAMAX": ..., "scalars": {"TOTAL_IP": ..., ...}, "profile": [{"NSA": i, "RNT": ...}, ...]}` layout so downstream tools can diff wrapper output against `fp` CLI runs. @@ -214,10 +224,9 @@ Mirrors the tr spec's §12.2 for the FP module: lifecycle. Repeated cycles (as in the 3×3 sweep) work in practice because `fp_init` resets the necessary SAVE state, but long-running processes can leak FPCOMM allocations. -- **String parameters not wired.** `KNAMFP` and similar namelist - keys cannot flow through the float-only C ABI; they stay at their - defaults. A follow-up PR with a `fp_set_param_str` entry point is - tracked separately (mirrors the trlib L-L approach). +- **One string parameter.** `fp_set_param_str` (`Fplib.set_param_str`) + sets `KNAMEQ`, the equilibrium file of `MODELG=3`; `KNAMFP` and the + other string keys of the namelist stay at their defaults. - **Compile-time mesh caps.** `FP_MAX_NRMAX=100` and `FP_MAX_NSAMAX=8` are baked into the exported `fp_state_t`; runs needing larger meshes require bumping the constants in @@ -231,8 +240,9 @@ Mirrors the tr spec's §12.2 for the FP module: one fixed layout, filled by `fp_get_state`. Changing it: 1. append the field at the end of both (never reorder, never remove); -2. bump `FP_STATE_ABI_VERSION` **in both files** (the `#define` and the - `PARAMETER`; currently `1`); +2. bump `FP_STATE_ABI_VERSION` **in all three places** (the `#define`, the + `PARAMETER` and the constant in `python/fplib/_ffi.py`, without which + `load_library` refuses the rebuilt library; currently `2`); 3. mirror the field in `python/fplib/_ffi.py` (`FpStateC._fields_`) and `state.py`; 4. rebuild `libfpapi.so` in the same commit. diff --git a/docs/sphinx/modules/fp/en/state.md b/docs/sphinx/modules/fp/en/state.md index b1fb08aba..1cbdb9fd0 100644 --- a/docs/sphinx/modules/fp/en/state.md +++ b/docs/sphinx/modules/fp/en/state.md @@ -16,15 +16,36 @@ the distribution function $f(r, p, \theta)$ itself — it exposes only ## Scalar quantities -`state.scalars` for `fp` contains **only `timefp`**. Unlike `tr`/`eq`, it -does not carry the 13/12 plasma-performance scalars — by design, `fp` -focuses on distribution-function moments and leaves derived quantities -(`BETA`, `TAUE`, etc.) to user-side post-processing. +The simulation time is a top-level field, and `state.scalars` is a dict +of **seven global quantities**: six volume-integrated diagnostics, each +summed over the active species at the latest global-diagnostic sample +(`NTG1`), and the plasma volume (`TVOLR`, the sum of the radial cell +volumes): ```python -state.timefp # simulation time [s] +state.timefp # simulation time [s] +state.scalars["TOTAL_IP"] # total plasma current [MA] ``` +| Key of `state.scalars` | Quantity | Unit | +|---|---|---| +| `TOTAL_IP` | Total plasma current | MA | +| `STORED_ENERGY` | Stored energy | MJ | +| `COLLISION_POWER` | Total collisional power | MW | +| `ABSORPTION_POWER` | Total absorbed wave power | MW | +| `ABSORPTION_WR` | Ray-tracing (WR) part of the absorbed power | MW | +| `ABSORPTION_WM` | Full-wave (WM) part of the absorbed power | MW | +| `PLASMA_VOLUME` | Plasma volume | m³ | + +All seven are zero until the first `run()`, and again after a +`set_param` or `set_param_str` until the next `run()` that is made, while +the profiles are +still those of the last run. Unlike the per-unit-volume profiles +below, they depend on the equilibrium through the volume element. `fp` +still does not carry the plasma-performance scalars of `tr`/`eq`: +derived quantities (`BETA`, `TAUE`, etc.) are left to user-side +post-processing. + ## Profile quantities (radial × species moment arrays) The main `fp` outputs. Each is a **2D array of shape `[nsamax][nrmax]`** @@ -71,7 +92,7 @@ For the full attribute list, see the `FpState` autodoc in |---|---|---|---| | **Species axis** | NSMAX | active NSA_MAX | active NSAMAX | | **Phase space** | 1D radial | 1D radial | **5D** (r, p, θ, species, t) | -| **Number of scalars** | 13 | 2 + 2int | **1** (timefp) | +| **Number of scalars** | 13 | 2 + 2int | **7** global (+ `timefp`) | | **Main outputs** | RN, RT, AJ, QP | RNA, RTA, RUA, RBP, ZEFF, BETA | RNT, RWT, RTT, RJT, RPCT, RPWT | | **Meaning** | Fluid moments | Integrated-transport moments | Distribution-function moments | | **Memory cost** | Low | Medium | **High** (5D) | diff --git a/docs/sphinx/modules/fp/ja/state.md b/docs/sphinx/modules/fp/ja/state.md index 1c8f4b421..3c410a80e 100644 --- a/docs/sphinx/modules/fp/ja/state.md +++ b/docs/sphinx/modules/fp/ja/state.md @@ -15,15 +15,32 @@ ## スカラー量 -`fp` の `state.scalars` は **`timefp` 1 個のみ** です. tr/eq のような -13/12 個のプラズマ性能指標は持ちません — これは fp の出力が分布関数 -モーメントが中心で, 派生量 (`BETA`, `TAUE` など) はユーザ側で計算する -設計だからです. +シミュレーション時刻はトップレベルのフィールドで, `state.scalars` は +**7 個の大域量** の dict です. うち 6 個は体積積分された診断量で, active +粒子種について和をとった最新の大域診断サンプル (`NTG1`) の値, 残る 1 個は +プラズマ体積 (`TVOLR`, 半径方向セル体積の和) です. ```python -state.timefp # シミュレーション時刻 [s] +state.timefp # シミュレーション時刻 [s] +state.scalars["TOTAL_IP"] # 全プラズマ電流 [MA] ``` +| `state.scalars` のキー | 物理量 | 単位 | +|---|---|---| +| `TOTAL_IP` | 全プラズマ電流 | MA | +| `STORED_ENERGY` | 蓄積エネルギー | MJ | +| `COLLISION_POWER` | 衝突パワーの合計 | MW | +| `ABSORPTION_POWER` | 波動吸収パワーの合計 | MW | +| `ABSORPTION_WR` | そのうちレイトレーシング (WR) 分 | MW | +| `ABSORPTION_WM` | そのうちフルウェーブ (WM) 分 | MW | +| `PLASMA_VOLUME` | プラズマ体積 | m³ | + +7 個とも最初の `run()` までは 0 で, `set_param` または `set_param_str` の後は次に `run()` が +実行されるまで再び 0 になります (その間もプロファイルは直前の run のままです). +下の単位体積あたりのプロファイルと違い, 体積要素を通じて +平衡に依存します. tr/eq のようなプラズマ性能指標は引き続き持ちません — +派生量 (`BETA`, `TAUE` など) はユーザ側で計算する設計です. + ## プロファイル量 (radial × species のモーメント配列) fp の主要出力. **`[nsamax][nrmax]` の 2 次元配列** で, active 粒子種ごとに @@ -67,7 +84,7 @@ state.to_dict() # JSON-ready dict |---|---|---|---| | **種別軸** | NSMAX | active NSA_MAX | active NSAMAX | | **位相空間** | 1D 半径 | 1D 半径 | **5D** (r, p, θ, species, t) | -| **スカラー数** | 13 | 2 + 2int | **1** (timefp) | +| **スカラー数** | 13 | 2 + 2int | **7** 大域量 (+ `timefp`) | | **主出力** | RN, RT, AJ, QP | RNA, RTA, RUA, RBP, ZEFF, BETA | RNT, RWT, RTT, RJT, RPCT, RPWT | | **意味** | 流体モーメント | 統合輸送モーメント | 分布関数モーメント | | **メモリ消費** | 低 | 中 | **高** (5D) | diff --git a/docs/sphinx/modules/wr/en/state.md b/docs/sphinx/modules/wr/en/state.md index 23aea760d..9fd9c20b2 100644 --- a/docs/sphinx/modules/wr/en/state.md +++ b/docs/sphinx/modules/wr/en/state.md @@ -90,6 +90,6 @@ For the complete attribute list see the `WrState` autodoc in |---|---|---|---| | **Time** | yes (`nt`) | yes (`timefp`) | **no** (spatial integration) | | **Dimensions** | radius × species | 5D phase space (moments) | **rays × steps** + radial profiles | -| **Scalars** | 13 (T, BETAN, ...) | 1 (timefp) | **4** (peak position + value) | +| **Scalars** | 13 (T, BETAN, ...) | 7 global (+ `timefp`) | **4** (peak position + value) | | **Main output** | profiles | moments | **ray trajectories + absorption profile** | | **Physical units** | fluid quantities | distribution-function moments | RF power deposition | diff --git a/docs/sphinx/modules/wr/ja/state.md b/docs/sphinx/modules/wr/ja/state.md index 39f145243..7c47e5263 100644 --- a/docs/sphinx/modules/wr/ja/state.md +++ b/docs/sphinx/modules/wr/ja/state.md @@ -88,6 +88,6 @@ state.to_dict() # JSON-ready dict |---|---|---|---| | **時間** | あり (`nt`) | あり (`timefp`) | **なし** (空間積分) | | **次元** | 半径 × 種 | 5D 位相空間 (モーメント化) | **レイ × ステップ** + 半径プロファイル | -| **スカラー** | 13 (T, BETAN, ...) | 1 (timefp) | **4** (ピーク位置 + 値) | +| **スカラー** | 13 (T, BETAN, ...) | 7 大域量 (+ `timefp`) | **4** (ピーク位置 + 値) | | **主要出力** | プロファイル | モーメント | **レイ軌跡 + 吸収プロファイル** | | **物理単位** | 流体量 | 分布関数モーメント | RF パワー堆積 | diff --git a/fp/fp_api.f90 b/fp/fp_api.f90 index 66e38c096..d958dc139 100644 --- a/fp/fp_api.f90 +++ b/fp/fp_api.f90 @@ -18,7 +18,9 @@ ! fp_get_state -> populates the fp_state_c struct from the same set of ! FPCOMM scalars/arrays that fpregress.f90 dumps: ! NRMAX, NSAMAX, NPMAX, NTHMAX, NTG2, TIMEFP, and -! profile arrays RNT/RWT/RTT/RJT/RPCT/RPWT(NR,NSA). +! profile arrays RNT/RWT/RTT/RJT/RPCT/RPWT(NR,NSA), +! plus the volume-integrated global scalars that +! FPWRTGLB otherwise only WRITEs to unit 6. ! fp_finalize -> fp_deallocate (+ ntg1/ntg2) and clears g_* flags. ! ! The Fortran-side names are fp_api_* to avoid colliding with the existing @@ -34,8 +36,9 @@ MODULE fp_api USE fp_state, ONLY: fp_state_c, FP_MAX_NRMAX, FP_MAX_NSAMAX USE fpcomm, ONLY: rkind, & NRMAX, NSAMAX, NPMAX, NTHMAX, NTG2, TIMEFP, & - NTMAX, & - RNT, RWT, RTT, RJT, RPCT, RPWT + NTMAX, NTG1, TVOLR, & + RNT, RWT, RTT, RJT, RPCT, RPWT, & + PIT, PWT, PPCT, PPWT, PWRT, PWMT USE fp_param_registry, ONLY: fp_param_set, fp_param_set_str, fp_param_generation USE plinit, ONLY: pl_init USE equnit, ONLY: eq_init @@ -262,11 +265,18 @@ END FUNCTION fp_api_run ! in Fortran column-major order, which matches the C-side declaration ! double RNT[FP_MAX_NSAMAX][FP_MAX_NRMAX] (row-major) byte-for-byte. ! Only state%RNT(1:NRMAX, 1:NSAMAX) carries valid runtime data. + ! + ! The 7 appended global scalars are the species sums FPWRTGLB builds + ! before printing them (fp/fpsave.f90:260-311). They are the only part + ! of the state that sees the volume element VOLR / TVOLR: the six + ! profile arrays are per-unit-volume moments and stay bit-identical + ! when only the equilibrium changes (MODELD=0, no wave, no CD), so + ! without these an EQ->FP handoff is invisible to the caller. !------------------------------------------------------------------- FUNCTION fp_api_get_state(state) RESULT(ierr) BIND(C, NAME="fp_get_state") TYPE(fp_state_c), INTENT(OUT) :: state INTEGER(C_INT) :: ierr - INTEGER :: nr, nsa, ntg_last, nr_alloc, nsa_alloc + INTEGER :: nr, nsa, ntg_last, ntg1_last, nr_alloc, nsa_alloc ! Always zero the struct so callers never see uninitialized memory. state%nrmax = 0 @@ -281,6 +291,13 @@ FUNCTION fp_api_get_state(state) RESULT(ierr) BIND(C, NAME="fp_get_state") state%RJT = 0.0_C_DOUBLE state%RPCT = 0.0_C_DOUBLE state%RPWT = 0.0_C_DOUBLE + state%TOTAL_IP = 0.0_C_DOUBLE + state%STORED_ENERGY = 0.0_C_DOUBLE + state%COLLISION_POWER = 0.0_C_DOUBLE + state%ABSORPTION_POWER = 0.0_C_DOUBLE + state%ABSORPTION_WR = 0.0_C_DOUBLE + state%ABSORPTION_WM = 0.0_C_DOUBLE + state%PLASMA_VOLUME = 0.0_C_DOUBLE IF (.NOT. g_initialized) THEN ierr = FP_ERR_NOT_INIT @@ -340,6 +357,35 @@ FUNCTION fp_api_get_state(state) RESULT(ierr) BIND(C, NAME="fp_get_state") END DO END IF + ! Global scalars. TVOLR and the P*T(NSA,NTG1) accumulators only + ! exist once fp_prep has run (TVOLR: fp/fpprep.f90:231-234; the + ! arrays: fpcomm::fp_allocate_ntg1) and only carry a sample once + ! FPSGLB has bumped NTG1 (fp/fploop.f90:200). Before that they stay + ! at the zeros set above rather than exposing uninitialised memory. + ! The same holds after a successful set, from fp_set_param, + ! fp_set_param_str or tot's "fp:" dispatcher (the registry's count + ! then differs from the one of the last fp_prep, #253): the profiles + ! above are then still the last run's (within the allocated extents, + ! #245) while these read zero until the next fp_prep. + IF (g_prepared .AND. g_prep_generation == fp_param_generation()) THEN + state%PLASMA_VOLUME = TVOLR + IF (NTG1 >= 1 .AND. ALLOCATED(PIT) .AND. ALLOCATED(PWT) .AND. & + ALLOCATED(PPCT) .AND. ALLOCATED(PPWT) .AND. & + ALLOCATED(PWRT) .AND. ALLOCATED(PWMT)) THEN + ntg1_last = MIN(NTG1, SIZE(PIT, 2)) + DO nsa = 1, MIN(NSAMAX, SIZE(PIT, 1)) + ! Species sums, mirroring FPWRTGLB's rtotal* accumulators + ! (fp/fpsave.f90:277-286). + state%TOTAL_IP = state%TOTAL_IP + PIT (nsa, ntg1_last) + state%STORED_ENERGY = state%STORED_ENERGY + PWT (nsa, ntg1_last) + state%COLLISION_POWER = state%COLLISION_POWER + PPCT(nsa, ntg1_last) + state%ABSORPTION_POWER = state%ABSORPTION_POWER + PPWT(nsa, ntg1_last) + state%ABSORPTION_WR = state%ABSORPTION_WR + PWRT(nsa, ntg1_last) + state%ABSORPTION_WM = state%ABSORPTION_WM + PWMT(nsa, ntg1_last) + END DO + END IF + END IF + ierr = FP_OK END FUNCTION fp_api_get_state diff --git a/fp/fp_api.h b/fp/fp_api.h index 1b6d96701..59afd8c05 100644 --- a/fp/fp_api.h +++ b/fp/fp_api.h @@ -22,7 +22,8 @@ extern "C" { * fp_run -> fp_prep on first call + fp_loop with NTMAX set to * the requested step count. * fp_get_state -> populates scalars + profile arrays (RNT/RWT/RTT/ - * RJT/RPCT/RPWT) from FPCOMM. + * RJT/RPCT/RPWT) + the volume-integrated global + * scalars (TOTAL_IP ... PLASMA_VOLUME) from FPCOMM. * fp_finalize -> clears lifecycle flags. Note: FPCOMM arrays are * NOT deallocated (fp_allocate / fp_deallocate * asymmetry) so a single fp_init/fp_run/fp_finalize @@ -34,6 +35,11 @@ extern "C" { * matching declaration is RNT(NRMAX, NSAMAX) (column-major). Layouts agree * byte-for-byte, but only RNT[0..nsamax-1][0..nrmax-1] are valid runtime * values (the remainder is padding up to FP_MAX_*). + * + * ABI note: fp_state_t only ever grows at the END of the struct, so the + * offsets of pre-existing members are stable. Every change must be + * mirrored in fp/fp_state.f90 and python/fplib/_ffi.py::FpStateC, and + * libfpapi.so rebuilt (`make -C fp libfpapi.so`). */ #define FP_MAX_NRMAX 100 @@ -50,8 +56,10 @@ extern "C" { * * History: * 1 -> initial layout + * 2 -> appended the seven global scalars TOTAL_IP .. PLASMA_VOLUME + * (56 bytes); existing field offsets preserved. */ -#define FP_STATE_ABI_VERSION 1 +#define FP_STATE_ABI_VERSION 2 /* Error codes returned by every fp_* entry point. */ enum fp_error { @@ -75,6 +83,20 @@ typedef struct { double RJT [FP_MAX_NSAMAX][FP_MAX_NRMAX]; double RPCT[FP_MAX_NSAMAX][FP_MAX_NRMAX]; double RPWT[FP_MAX_NSAMAX][FP_MAX_NRMAX]; + /* Global (volume-integrated) plasma scalars, species-summed at the + * latest NTG1 sample. These are the quantities FP otherwise only + * WRITEs to unit 6: the rtotal* line block (fp/fpsave.f90:305-311), + * the W column of the per-species line (fp/fpsave.f90:252) and the + * TVOLR banner (fp/fpprep.f90:237). Zero until the first run, and + * again after any successful set (fp_set_param, fp_set_param_str, or + * tot's fp: parameters; #253) until the next run that is made. */ + double TOTAL_IP; /* total plasma current [MA] */ + double STORED_ENERGY; /* stored energy [MJ] */ + double COLLISION_POWER; /* total collision power [MW] */ + double ABSORPTION_POWER; /* total absorbed power [MW] */ + double ABSORPTION_WR; /* WR (ray-tracing) part [MW] */ + double ABSORPTION_WM; /* WM (full-wave) part [MW] */ + double PLASMA_VOLUME; /* TVOLR [m^3] */ } fp_state_t; int fp_init(void); diff --git a/fp/fp_state.f90 b/fp/fp_state.f90 index 3325ed4b7..28d9df852 100644 --- a/fp/fp_state.f90 +++ b/fp/fp_state.f90 @@ -13,6 +13,12 @@ ! in Fortran the matching declaration is RNT(NRMAX, NSAMAX) (column-major). ! The two layouts agree byte-for-byte, but only RNT(0:nrmax-1, 0:nsamax-1) ! carry valid runtime data (the rest is padding up to FP_MAX_*). +! +! ABI note: new members are APPENDED at the end of fp_state_c so the +! offsets of every pre-existing member stay put. Any change here must be +! mirrored in fp/fp_api.h and python/fplib/_ffi.py::FpStateC, and +! libfpapi.so must be rebuilt (`make -C fp libfpapi.so`) -- a stale .so +! read through a newer ctypes layout returns garbage. MODULE fp_state USE, INTRINSIC :: ISO_C_BINDING @@ -29,7 +35,7 @@ MODULE fp_state ! removed or reordered; fp_state_abi_version() exports it and the C-ABI ! driver fails when the two drift. fp_state_size() exports the ! BIND(C) size, so the library, not the header, is authoritative. - INTEGER(C_INT), PARAMETER :: FP_STATE_ABI_VERSION = 1 + INTEGER(C_INT), PARAMETER :: FP_STATE_ABI_VERSION = 2 TYPE, BIND(C) :: fp_state_c INTEGER(C_INT) :: nrmax @@ -47,6 +53,22 @@ MODULE fp_state REAL(C_DOUBLE) :: RJT (FP_MAX_NRMAX, FP_MAX_NSAMAX) REAL(C_DOUBLE) :: RPCT(FP_MAX_NRMAX, FP_MAX_NSAMAX) REAL(C_DOUBLE) :: RPWT(FP_MAX_NRMAX, FP_MAX_NSAMAX) + + ! --- global (volume-integrated) plasma scalars, appended ------- + ! Species-summed counterparts of quantities FP otherwise only + ! prints to unit 6: FPWRTGLB's rtotal* block + ! (fp/fpsave.f90:305-311), the W column of the per-species line + ! (fp/fpsave.f90:252) and the TVOLR banner (fp/fpprep.f90:237). + ! Unlike the profile arrays above, these depend on the volume + ! element VOLR / TVOLR, so they are the FP state's only observable + ! response to a change of equilibrium (MODELG=3 + KNAMEQ). + REAL(C_DOUBLE) :: TOTAL_IP ! total plasma current [MA] + REAL(C_DOUBLE) :: STORED_ENERGY ! stored energy [MJ] + REAL(C_DOUBLE) :: COLLISION_POWER ! total collision power [MW] + REAL(C_DOUBLE) :: ABSORPTION_POWER ! total absorbed power [MW] + REAL(C_DOUBLE) :: ABSORPTION_WR ! WR (ray-tracing) part [MW] + REAL(C_DOUBLE) :: ABSORPTION_WM ! WM (full-wave) part [MW] + REAL(C_DOUBLE) :: PLASMA_VOLUME ! TVOLR [m^3] END TYPE fp_state_c END MODULE fp_state diff --git a/python/fplib/README.md b/python/fplib/README.md index 2e211c1ec..ffb60e2a2 100644 --- a/python/fplib/README.md +++ b/python/fplib/README.md @@ -209,6 +209,29 @@ via `state.to_dict()` (JSON-serialisable, layout compatible with | `RJT` | list[list[float]] | `[nsamax][nrmax]` current density profile | | `RPCT` | list[list[float]] | `[nsamax][nrmax]` collisional power profile | | `RPWT` | list[list[float]] | `[nsamax][nrmax]` wave-absorbed power profile | +| `scalars` | dict[str, float] | the seven volume-integrated globals below | + +`scalars` holds the seven doubles that `fp_state_t` carries after the +profile arrays: the first six are summed over the kinetic species at the +latest global-diagnostic sample (`NTG1`), and `PLASMA_VOLUME` is `TVOLR`, +the sum of the radial cell volumes. All seven are zero until the first +`run()`, and again after a `set_param` or `set_param_str` until the next +`run()` that is made, while the profiles are still those of the last run: + +| Key | Unit | Meaning | +|---|---|---| +| `TOTAL_IP` | MA | total plasma current | +| `STORED_ENERGY` | MJ | stored energy | +| `COLLISION_POWER` | MW | total collisional power | +| `ABSORPTION_POWER` | MW | total absorbed wave power | +| `ABSORPTION_WR` | MW | ray-tracing (WR) part of the absorbed power | +| `ABSORPTION_WM` | MW | full-wave (WM) part of the absorbed power | +| `PLASMA_VOLUME` | m^3 | plasma volume (`TVOLR`) | + +`state.to_dict()` puts them under a top-level `"scalars"` key. +`fp_state_t` grows at its end only: a `libfpapi.so` built before these +members existed must be rebuilt (`make -C fp libfpapi.so`) together with +this wrapper. Memory-layout note: C declares `RNT[NSAMAX][NRMAX]` (row-major); the Fortran backend declares `RNT(NRMAX, NSAMAX)` (column-major). Both lay @@ -258,9 +281,9 @@ menu — those live in `fp/fp` only. between `fp_allocate` and `fp_deallocate`). A single `fp_init` / `fp_run` / `fp_finalize` cycle per process is the supported lifecycle at Phase L-3. -- **String parameters not yet wired** (e.g. `KNAMFP`). File-name - namelist keys must be left at their defaults; the float-only C ABI - cannot carry strings. +- **One string parameter.** `set_param_str` sets `KNAMEQ`, the + equilibrium file of `MODELG=3`; `KNAMFP` and the other string keys of + the namelist stay at their defaults. - **Unregistered namelist keys** — any name missing from `fp_param_registry.f90` returns `FplibInvalidParamError`. The current registry covers ~55 settable variables across geometry, diff --git a/python/fplib/_ffi.py b/python/fplib/_ffi.py index 135e4db1d..63e10cdd7 100644 --- a/python/fplib/_ffi.py +++ b/python/fplib/_ffi.py @@ -46,7 +46,7 @@ # and fp/fp_state.f90. load_library() refuses a library whose # fp_state_abi_version() / fp_state_size() differ from this and from # ctypes.sizeof(FpStateC) (#246). -FP_STATE_ABI_VERSION = 1 +FP_STATE_ABI_VERSION = 2 # --------------------------------------------------------------------- @@ -66,6 +66,12 @@ # In C, RNT[NSAMAX][NRMAX] is row-major. # In Fortran the matching declaration is RNT(NRMAX, NSAMAX) column-major. # Both lay out the same bytes; only the index order differs. +# +# ABI note: ``_fields_`` must match fp/fp_api.h::fp_state_t (and +# fp/fp_state.f90::fp_state_c) member-for-member, in order. New members +# are only ever appended, and libfpapi.so must be rebuilt in the same +# commit -- reading a stale .so through a longer layout yields garbage +# in the trailing members. # --------------------------------------------------------------------- _PROFILE_2D = (ctypes.c_double * FP_MAX_NRMAX) * FP_MAX_NSAMAX @@ -84,6 +90,14 @@ class FpStateC(ctypes.Structure): ("RJT", _PROFILE_2D), ("RPCT", _PROFILE_2D), ("RPWT", _PROFILE_2D), + # --- global (volume-integrated) plasma scalars --------------- + ("TOTAL_IP", ctypes.c_double), + ("STORED_ENERGY", ctypes.c_double), + ("COLLISION_POWER", ctypes.c_double), + ("ABSORPTION_POWER", ctypes.c_double), + ("ABSORPTION_WR", ctypes.c_double), + ("ABSORPTION_WM", ctypes.c_double), + ("PLASMA_VOLUME", ctypes.c_double), ] diff --git a/python/fplib/state.py b/python/fplib/state.py index 669c5a2f6..76e5c97d9 100644 --- a/python/fplib/state.py +++ b/python/fplib/state.py @@ -18,6 +18,43 @@ # Profile array names in canonical order. These are the 6 2-D arrays # declared in fp_state_t: RNT, RWT, RTT, RJT, RPCT, RPWT. PROFILE_FIELDS = ("RNT", "RWT", "RTT", "RJT", "RPCT", "RPWT") + +# Global (volume-integrated) scalar field names, in fp_state_t order. +# Grouped into the ``scalars`` dict by ``to_dict()`` so the FP wire +# format has the same top-level shape as TrState.to_dict(). +# +# Units, traced to the Fortran (all species-summed at the latest NTG1; +# the per-species PxT accumulators are built in fp/fpsave.f90:106-147): +# TOTAL_IP [MA] sum PIT; PIT = sum_NR RJS*VOLR / (2*pi*RR) +# -- RJS carries the *1.D-6 A->MA conversion +# in fp/fpsave.f90:1310; printed as +# "total plasma current [MA]" (fpsave.f90:311). +# STORED_ENERGY [MJ] sum PWT; PWT = sum_NR RWS*VOLR +# -- RWS is *1.D-6 J->MJ, fp/fpsave.f90:1377; +# printed as the W column of fpsave.f90:252. +# COLLISION_POWER [MW] sum PPCT; RPCS is *1.D-6, fp/fpsave.f90:1509; +# printed at fpsave.f90:307. +# ABSORPTION_POWER [MW] sum PPWT; RPWS is *1.D-6, fp/fpsave.f90:1510; +# printed at fpsave.f90:305. +# ABSORPTION_WR [MW] sum PWRT (ray-tracing share of the above). +# ABSORPTION_WM [MW] sum PWMT (full-wave share of the above). +# PLASMA_VOLUME [m^3] TVOLR = sum_NR VOLR, fp/fpprep.f90:231-234; +# printed in the DEVICE banner, fpprep.f90:237. +# +# All seven are zero until the first ``run()`` (the Fortran accumulators +# only exist after fp_prep, and only carry a sample once FPSGLB has +# bumped NTG1), and again after a ``set_param`` or ``set_param_str`` until +# the next ``run()`` that is made: fp_get_state fills them only while the +# state is prepared and no parameter was set since that preparation. +SCALAR_FIELDS = ( + "TOTAL_IP", + "STORED_ENERGY", + "COLLISION_POWER", + "ABSORPTION_POWER", + "ABSORPTION_WR", + "ABSORPTION_WM", + "PLASMA_VOLUME", +) _DIM_BOUNDS = ( ("nrmax", "FP_MAX_NRMAX", _ffi.FP_MAX_NRMAX), ("nsamax", "FP_MAX_NSAMAX", _ffi.FP_MAX_NSAMAX), @@ -52,6 +89,9 @@ class FpState: RNT/RWT/RTT/RJT/RPCT/RPWT: [nsamax][nrmax] profile arrays, truncated to the runtime extents (FP_MAX_NSAMAX x FP_MAX_NRMAX padding is dropped). + scalars: + dict of the 7 volume-integrated global quantities + (see :data:`SCALAR_FIELDS` for names and units). """ nrmax: int @@ -66,6 +106,7 @@ class FpState: RJT: List[List[float]] = field(default_factory=list) RPCT: List[List[float]] = field(default_factory=list) RPWT: List[List[float]] = field(default_factory=list) + scalars: Dict[str, float] = field(default_factory=dict) @classmethod def from_c(cls, s: FpStateC) -> "FpState": @@ -99,6 +140,7 @@ def _slice(arr_2d): RJT=_slice(s.RJT), RPCT=_slice(s.RPCT), RPWT=_slice(s.RPWT), + scalars={k: float(getattr(s, k)) for k in SCALAR_FIELDS}, ) def to_dict(self) -> Dict[str, Any]: @@ -109,12 +151,18 @@ def to_dict(self) -> Dict[str, Any]: { "NRMAX": int, "NSAMAX": int, "NPMAX": int, "NTHMAX": int, "NTG2": int, "TIMEFP": float, + "scalars": {"TOTAL_IP": float, ...}, # 7 global scalars "profile": [ # per species, length NSAMAX {"NSA": 1, "RNT": [..NR..], "RWT": [..], "RTT": [..], "RJT": [..], "RPCT": [..], "RPWT": [..]}, ... ], } + + The top-level ``scalars`` key mirrors ``TrState.to_dict()`` so + consumers (the MCP bridge's steady-state detection, plot + renderers, regression tooling) can treat fp and tr states the + same way. """ return { "NRMAX": self.nrmax, @@ -123,6 +171,7 @@ def to_dict(self) -> Dict[str, Any]: "NTHMAX": self.nthmax, "NTG2": self.ntg2, "TIMEFP": self.timefp, + "scalars": dict(self.scalars), "profile": [ { "NSA": ns + 1, @@ -138,4 +187,4 @@ def to_dict(self) -> Dict[str, Any]: } -__all__ = ["FpState", "PROFILE_FIELDS"] +__all__ = ["FpState", "PROFILE_FIELDS", "SCALAR_FIELDS"] diff --git a/python/fplib/tests/test_ffi.py b/python/fplib/tests/test_ffi.py index 3b18ead9c..43f0f5262 100644 --- a/python/fplib/tests/test_ffi.py +++ b/python/fplib/tests/test_ffi.py @@ -28,6 +28,7 @@ sys.path.insert(0, str(PYTHON_ROOT)) from fplib import _ffi # noqa: E402 +from fplib import state # noqa: E402 REPO = HERE.parents[3] @@ -63,13 +64,31 @@ def test_has_expected_fields(self): ): self.assertIn(n, names, f"missing field {n}") + def test_global_scalars_are_appended_last(self): + # ABI contract: the 7 global scalars sit at the END of the + # struct so every pre-existing member keeps its offset. If a + # future field is inserted before them, this fails. + names = [f[0] for f in _ffi.FpStateC._fields_] + self.assertEqual(names[-7:], list(state.SCALAR_FIELDS)) + for n in state.SCALAR_FIELDS: + self.assertEqual( + dict(f[:2] for f in _ffi.FpStateC._fields_)[n], ctypes.c_double + ) + + def test_layout_version_is_the_one_with_the_global_scalars(self): + # fp_api.h's History: 1 is the profiles only, 2 appends the seven + # global scalars. load_library refuses a libfpapi.so of another + # version (#246), so a library built before them is not read. + self.assertEqual(_ffi.FP_STATE_ABI_VERSION, 2) + self.assertEqual(len(state.SCALAR_FIELDS), 7) + def test_size_matches_header_math(self): - # 5 ints + 1 double + 6 * NSA*NR doubles. Compilers may pad the - # 5 ints to a multiple of 8 before the first double, so we - # accept the natural rounding range. + # 5 ints + 1 double + 6 * NSA*NR doubles + 7 scalar doubles. + # Compilers may pad the 5 ints to a multiple of 8 before the + # first double, so we accept the natural rounding range. nr = _ffi.FP_MAX_NRMAX nsa = _ffi.FP_MAX_NSAMAX - exp_core = 8 + 6 * nr * nsa * 8 + exp_core = 8 + 6 * nr * nsa * 8 + len(state.SCALAR_FIELDS) * 8 sz = ctypes.sizeof(_ffi.FpStateC) # 5 ints = 20 bytes, padded up to either 24 (8-byte align) or # possibly 32 on pathological ABIs. diff --git a/python/fplib/tests/test_fplib.py b/python/fplib/tests/test_fplib.py index b84fd1c1e..ec6fab8cd 100644 --- a/python/fplib/tests/test_fplib.py +++ b/python/fplib/tests/test_fplib.py @@ -88,6 +88,13 @@ def _populated_state(self, nr: int = 3, nsa: int = 2) -> _ffi.FpStateC: s.RJT[ns][ir] = v + 0.3 s.RPCT[ns][ir] = v + 0.4 s.RPWT[ns][ir] = v + 0.5 + s.TOTAL_IP = 1.5 + s.STORED_ENERGY = 5.25 + s.COLLISION_POWER = -0.5 + s.ABSORPTION_POWER = 2.0 + s.ABSORPTION_WR = 1.25 + s.ABSORPTION_WM = 0.75 + s.PLASMA_VOLUME = 28.4 return s def test_from_c_slices_correctly(self): @@ -105,6 +112,17 @@ def test_from_c_slices_correctly(self): self.assertAlmostEqual(st.RTT[0][1], 1.2) self.assertAlmostEqual(st.RPWT[1][0], 100.5) + def test_from_c_copies_global_scalars(self): + from fplib.state import SCALAR_FIELDS + + c = self._populated_state(nr=3, nsa=2) + st = FpState.from_c(c) + self.assertEqual(sorted(st.scalars), sorted(SCALAR_FIELDS)) + self.assertAlmostEqual(st.scalars["TOTAL_IP"], 1.5) + self.assertAlmostEqual(st.scalars["STORED_ENERGY"], 5.25) + self.assertAlmostEqual(st.scalars["COLLISION_POWER"], -0.5) + self.assertAlmostEqual(st.scalars["PLASMA_VOLUME"], 28.4) + def test_to_dict_shape(self): c = self._populated_state(nr=2, nsa=2) d = FpState.from_c(c).to_dict() @@ -112,6 +130,9 @@ def test_to_dict_shape(self): self.assertEqual(d["NSAMAX"], 2) self.assertIn("NPMAX", d) self.assertIn("TIMEFP", d) + # `scalars` mirrors the TrState.to_dict() top-level shape. + self.assertIsInstance(d["scalars"], dict) + self.assertAlmostEqual(d["scalars"]["PLASMA_VOLUME"], 28.4) self.assertEqual(len(d["profile"]), 2) p0 = d["profile"][0] self.assertEqual(p0["NSA"], 1) @@ -128,6 +149,7 @@ def test_to_dict_json_serialisable(self): d2 = json.loads(s) self.assertEqual(d2["NRMAX"], 2) self.assertEqual(d2["profile"][0]["NSA"], 1) + self.assertAlmostEqual(d2["scalars"]["STORED_ENERGY"], 5.25) # ===================================================================== diff --git a/python/fplib/tests/test_state_scalars_window.py b/python/fplib/tests/test_state_scalars_window.py new file mode 100644 index 000000000..2790f38de --- /dev/null +++ b/python/fplib/tests/test_state_scalars_window.py @@ -0,0 +1,157 @@ +"""The seven global scalars of ``fp_get_state`` between a ``set_param`` and +the next run. + +A successful set (``fp_set_param``, ``fp_set_param_str``, or tot's ``fp:`` +parameters) advances the registry's count (#253), so the next run prepares +again, and the arrays are not touched until it does. In that window +``fp_get_state`` still returns the profiles of the last run, within the +extents that are allocated (#245), while the global scalars (``TOTAL_IP`` +... ``PLASMA_VOLUME``) are filled only while the state is prepared and no +parameter was set since, and read zero. The two rules meet here; these +cases pin what the state documentation says about them. One test per case, +no ``subTest``. +""" +from __future__ import annotations + +import sys +from pathlib import Path + +import pytest + +HERE = Path(__file__).resolve() +REPO = HERE.parents[3] +PYTHON_ROOT = REPO / "python" +DEFAULT_SO = REPO / "fp" / "libfpapi.so" + +if str(PYTHON_ROOT) not in sys.path: + sys.path.insert(0, str(PYTHON_ROOT)) + +from fplib import Fplib, FplibInvalidParamError # noqa: E402 +from fplib.state import SCALAR_FIELDS # noqa: E402 + +pytestmark = pytest.mark.skipif( + not DEFAULT_SO.exists(), + reason=f"libfpapi.so not built at {DEFAULT_SO}; run `make -C fp libfpapi.so`", +) + +# Small momentum grid so one step is quick. +NPMAX = 10 +NTHMAX = 10 +DELT = 0.002 + +PROFILES = ("RNT", "RWT", "RTT", "RJT", "RPCT", "RPWT") + + +def _configure(fp: Fplib) -> None: + fp.set_param("NPMAX", float(NPMAX)) + fp.set_param("NTHMAX", float(NTHMAX)) + fp.set_param("DELT", DELT) + + +def _all_zero(state) -> bool: + return all(state.scalars[name] == 0.0 for name in SCALAR_FIELDS) + + +def test_the_state_has_the_seven_scalars(): + with Fplib() as fp: + state = fp.get_state() + assert len(SCALAR_FIELDS) == 7 + assert set(state.scalars) == set(SCALAR_FIELDS) + + +def test_before_the_first_run_the_scalars_are_zero(): + with Fplib() as fp: + _configure(fp) + state = fp.get_state() + assert _all_zero(state) + + +def test_after_a_run_the_volume_is_there(): + with Fplib() as fp: + _configure(fp) + fp.run(1) + state = fp.get_state() + assert state.scalars["PLASMA_VOLUME"] > 0.0 + assert state.scalars["STORED_ENERGY"] > 0.0 + + +def test_a_run_of_zero_steps_gives_the_volume_and_the_initial_sample(): + """``run(0)`` prepares only: the volume exists, and the species sums are + those of the initial distribution.""" + with Fplib() as fp: + _configure(fp) + fp.run(0) + state = fp.get_state() + assert state.timefp == 0.0 + assert state.scalars["PLASMA_VOLUME"] > 0.0 + assert state.scalars["STORED_ENERGY"] > 0.0 + + +# (parameter, change, the extents (nrmax, nsamax) reported afterwards): a raised +# count stays capped at what the run allocated, a lowered one applies at once. +@pytest.mark.parametrize( + "name, change, extents", + [("NSAMAX", +1, (3, 1)), ("NRMAX", +1, (3, 1)), ("NRMAX", -1, (2, 1))], + ids=["NSAMAX+1", "NRMAX+1", "NRMAX-1"], +) +def test_after_a_set_param_the_scalars_are_zero_and_the_profiles_the_last_run_s( + name, change, extents): + with Fplib() as fp: + _configure(fp) + fp.set_param("NRMAX", 3.0) + fp.run(1) + before = fp.get_state() + assert (before.nrmax, before.nsamax) == (3, 1) + assert before.timefp > 0.0 + assert before.scalars["PLASMA_VOLUME"] > 0.0 + assert before.scalars["STORED_ENERGY"] > 0.0 + fp.set_param(name, float(getattr(before, name.lower()) + change)) + after = fp.get_state() + assert _all_zero(after) + # the profiles are still the last run's, within what is allocated (#245) + assert (after.nrmax, after.nsamax) == extents + assert after.timefp == before.timefp + nr, nsa = extents + for profile in PROFILES: + rows_before, rows_after = getattr(before, profile), getattr(after, profile) + assert rows_after == [row[:nr] for row in rows_before[:nsa]], profile + + +def test_a_refused_set_leaves_the_scalars(): + # The registry counts successful sets only: a name it does not know + # changes nothing, and the scalars stay those of the last run. + with Fplib() as fp: + _configure(fp) + fp.run(1) + before = fp.get_state() + with pytest.raises(FplibInvalidParamError): + fp.set_param("NO_SUCH_PARAMETER", 1.0) + after = fp.get_state() + assert before.scalars["PLASMA_VOLUME"] > 0.0 + assert after.scalars == before.scalars + + +def test_a_string_set_opens_the_same_window(): + # The registry counts fp_set_param_str (KNAMEQ) like a numeric set (#253). + with Fplib() as fp: + _configure(fp) + fp.run(1) + assert fp.get_state().scalars["PLASMA_VOLUME"] > 0.0 + fp.set_param_str("KNAMEQ", "eqdata") + assert _all_zero(fp.get_state()) + fp.run(1) + state = fp.get_state() + assert state.scalars["PLASMA_VOLUME"] > 0.0 + assert state.scalars["STORED_ENERGY"] > 0.0 + + +def test_the_next_run_fills_the_scalars_again(): + with Fplib() as fp: + _configure(fp) + fp.run(1) + fp.set_param("DELT", DELT) + assert _all_zero(fp.get_state()) + fp.run(1) + state = fp.get_state() + assert state.scalars["PLASMA_VOLUME"] > 0.0 + assert state.scalars["STORED_ENERGY"] > 0.0 diff --git a/python/mcp-servers/fp_mcp/server.py b/python/mcp-servers/fp_mcp/server.py index cdca3cec5..c5604348d 100644 --- a/python/mcp-servers/fp_mcp/server.py +++ b/python/mcp-servers/fp_mcp/server.py @@ -84,6 +84,7 @@ FplibCalcFailedError, FplibNotImplementedError, ) +from fplib.state import SCALAR_FIELDS # noqa: E402 # ===================================================================== @@ -176,7 +177,10 @@ "Snapshot of FPCOMM. Profile arrays are [NSAMAX][NRMAX] in the C " "layout (row-major), truncated from the FP_MAX_NSAMAX=8 x " "FP_MAX_NRMAX=100 storage. Larger meshes require bumping these " - "caps in fp/fp_api.h and rebuilding libfpapi.so." + "caps in fp/fp_api.h and rebuilding libfpapi.so. The `scalars` " + "block carries the volume-integrated globals; the profile arrays " + "are per-unit-volume moments and are insensitive to the " + "equilibrium on their own." ), "properties": { "NRMAX": {"type": "integer", "description": "radial points in use"}, @@ -185,6 +189,32 @@ "NTHMAX": {"type": "integer", "description": "pitch-angle mesh size"}, "NTG2": {"type": "integer", "description": "long-time-axis counter"}, "TIMEFP": {"type": "number", "description": "simulation time [s]"}, + "scalars": { + "type": "object", + "description": ( + "7 volume-integrated global plasma scalars at the latest " + "NTG1 sample, summed over species: TOTAL_IP [MA], " + "STORED_ENERGY [MJ], COLLISION_POWER [MW], " + "ABSORPTION_POWER [MW], ABSORPTION_WR [MW], " + "ABSORPTION_WM [MW], PLASMA_VOLUME [m^3]. Unlike the " + "per-unit-volume profile arrays these depend on the " + "equilibrium through the volume element, so they are what " + "makes an EQ->FP handoff (MODELG=3 + KNAMEQ) observable. " + "All zero before the first run, and again after a set_param " + "until the next run that is made." + ), + "properties": { + "TOTAL_IP": {"type": "number", "description": "total plasma current [MA]"}, + "STORED_ENERGY": {"type": "number", "description": "stored energy [MJ]"}, + "COLLISION_POWER": {"type": "number", "description": "total collision power [MW]"}, + "ABSORPTION_POWER": {"type": "number", "description": "total absorbed wave power [MW]"}, + "ABSORPTION_WR": {"type": "number", "description": "ray-tracing (WR) share of the absorbed power [MW]"}, + "ABSORPTION_WM": {"type": "number", "description": "full-wave (WM) share of the absorbed power [MW]"}, + "PLASMA_VOLUME": {"type": "number", "description": "plasma volume TVOLR [m^3]"}, + }, + "required": list(SCALAR_FIELDS), + "additionalProperties": {"type": "number"}, + }, "profile": { "type": "array", "description": "per-species profile list, length = NSAMAX", @@ -202,7 +232,10 @@ }, }, }, - "required": ["NRMAX", "NSAMAX", "NPMAX", "NTHMAX", "NTG2", "TIMEFP", "profile"], + "required": [ + "NRMAX", "NSAMAX", "NPMAX", "NTHMAX", "NTG2", "TIMEFP", + "scalars", "profile", + ], } diff --git a/python/mcp-servers/fp_mcp/tests/test_server.py b/python/mcp-servers/fp_mcp/tests/test_server.py index 5bf63e65d..42c5c2dc8 100644 --- a/python/mcp-servers/fp_mcp/tests/test_server.py +++ b/python/mcp-servers/fp_mcp/tests/test_server.py @@ -99,7 +99,8 @@ class TestStateSchema(unittest.TestCase): def test_schema_top_level(self) -> None: schema = srv.STATE_SCHEMA self.assertEqual(schema["type"], "object") - for k in ("NRMAX", "NSAMAX", "NPMAX", "NTHMAX", "NTG2", "TIMEFP", "profile"): + for k in ("NRMAX", "NSAMAX", "NPMAX", "NTHMAX", "NTG2", "TIMEFP", + "scalars", "profile"): self.assertIn(k, schema["properties"]) def test_profile_items_cover_all_2d_fields(self) -> None: @@ -108,6 +109,15 @@ def test_profile_items_cover_all_2d_fields(self) -> None: for field in ("NSA", "RNT", "RWT", "RTT", "RJT", "RPCT", "RPWT"): self.assertIn(field, props) + def test_scalars_cover_all_global_fields(self) -> None: + """The schema must document every field fplib actually emits.""" + from fplib.state import SCALAR_FIELDS + + scalars = srv.STATE_SCHEMA["properties"]["scalars"] + self.assertEqual(sorted(scalars["properties"]), sorted(SCALAR_FIELDS)) + # get_state always emits all seven: the schema says so. + self.assertEqual(scalars["required"], list(SCALAR_FIELDS)) + def test_describe_state_schema_tool(self) -> None: out = srv.handle_describe_state_schema() self.assertEqual(out["type"], "object") @@ -198,7 +208,7 @@ def to_dict(self_inner) -> Dict[str, Any]: return { "NRMAX": 0, "NSAMAX": 0, "NPMAX": 0, "NTHMAX": 0, "NTG2": 0, "TIMEFP": 0.0, - "profile": [], + "scalars": {}, "profile": [], } return _S() diff --git a/python/totlib/pipeline.py b/python/totlib/pipeline.py index cc25bc04e..47a306fba 100644 --- a/python/totlib/pipeline.py +++ b/python/totlib/pipeline.py @@ -54,8 +54,10 @@ class CouplingRule: src_state_key (transfer-kind only) resolves a value from the source module's get_state(): - * str -> looked up in prev_state.scalars (only valid for modules - that expose a .scalars dict — Tr/Eq/Wr/Wrx/Ti). + * str -> looked up in the flat scalars of the previous state + (see _state_to_scalars: a module's .scalars dict; for Fp also its + top-level timefp and counts; the numeric top-level attributes of a + state that has no such dict). * callable(prev_state, params) -> float. params is TotPipeline._params, a dict of every set_param call so far keyed as ":" (used by fp's compute_rjt_volint which needs tr's RR/RA). @@ -146,18 +148,9 @@ def to_dict(self) -> Dict[str, Any]: return out -def _state_to_scalars(state) -> Dict[str, float]: - """Extract a flat dict of numeric scalars from any module's state. - - Modules that expose state.scalars dict (Tr, Eq, Wr, Wrx, Ti) — return - a shallow copy. Modules that don't (Fp) — extract numeric top-level - dataclass attributes (excluding bool, str, and non-numeric). - - Per spec §5.3. Required because FpState lacks the uniform .scalars - dict the spec originally assumed (R1/R2 finding). - """ - if hasattr(state, "scalars") and isinstance(state.scalars, dict): - return dict(state.scalars) +def _top_level_numbers(state) -> Dict[str, float]: + """The numeric top-level attributes of a state, as floats (bool, str + and anything that is not a number are left out).""" out: Dict[str, float] = {} for name in dir(state): if name.startswith("_"): @@ -172,6 +165,43 @@ def _state_to_scalars(state) -> Dict[str, float]: return out +def _is_fp_state(state) -> bool: + """True for fplib's ``FpState`` (or a subclass), told by the class's own + name and module: nothing is imported for the test, and the package may + have been imported under another name (``python.fplib``).""" + return any( + cls.__name__ == "FpState" + and isinstance(cls.__module__, str) + and cls.__module__.split(".")[-2:] == ["fplib", "state"] + for cls in type(state).__mro__ + ) + + +def _state_to_scalars(state) -> Dict[str, float]: + """Extract a flat dict of numeric scalars from any module's state. + + * A state with a ``scalars`` dict — a shallow copy of it. + * ``FpState`` — its ``scalars`` dict (the seven volume-integrated + globals) together with its numeric top-level attributes + (``nrmax``, ``nsamax``, ``npmax``, ``nthmax``, ``ntg2``, + ``timefp``). Those attributes were the whole of the FP pipeline + scalars before ``FpState`` had a ``scalars`` dict, so they stay in + pipeline results beside it. + * A state without a ``scalars`` dict — its numeric top-level + dataclass attributes (excluding bool, str, and non-numeric). + + Per spec §5.3. The top-level route exists because FpState lacked the + uniform .scalars dict the spec originally assumed (R1/R2 finding). + """ + if hasattr(state, "scalars") and isinstance(state.scalars, dict): + if _is_fp_state(state): + out = _top_level_numbers(state) + out.update(state.scalars) + return out + return dict(state.scalars) + return _top_level_numbers(state) + + _MODULE_REGISTRY: Dict[str, Tuple[str, str, str]] = { # name : (package, wrapper_class, base_error_class) "fp": ("fplib", "Fplib", "FplibError"), @@ -542,8 +572,8 @@ def _extract_source(self, prev_state, rule: CouplingRule) -> float: The params dict carries every set_param call so the rule can pull cross-module context (e.g. fp's compute_rjt_volint needs tr's RR/RA). - If src_state_key is a string, look it up in - _state_to_scalars(prev_state) — handles both .scalars-bearing modules - and FpState's top-level attribute layout. + _state_to_scalars(prev_state) — a module's .scalars dict, and for + FpState its top-level attributes (timefp, the counts) beside it. Wraps lookup/computation errors as TotPipelineCouplingError. """ if callable(rule.src_state_key): @@ -644,7 +674,7 @@ def run_pipeline(self, steps) -> PipelineResult: module.run(**kwargs) cur_state = module.get_state() - # _state_to_scalars adapter handles fp (no .scalars) and tr (.scalars dict) uniformly. + # _state_to_scalars gives tr's .scalars dict, and fp's with timefp and the counts beside it. result_steps.append(PipelineStep( module=name, scalars=_state_to_scalars(cur_state), diff --git a/python/totlib/tests/test_pipeline.py b/python/totlib/tests/test_pipeline.py index ec5ca6240..d087f415d 100644 --- a/python/totlib/tests/test_pipeline.py +++ b/python/totlib/tests/test_pipeline.py @@ -454,6 +454,69 @@ class FakeState: assert out == {"nrmax": 10.0, "timefp": 0.5} +def test_state_to_scalars_keeps_fp_top_level_beside_its_scalars(): + """A real FpState gives its seven globals AND the numeric top-level + attributes that were the FP pipeline scalars before FpState had a + ``scalars`` dict (``timefp`` among them).""" + from fplib.state import SCALAR_FIELDS, FpState + from totlib.pipeline import _state_to_scalars + + state = FpState( + nrmax=3, nsamax=1, npmax=20, nthmax=20, ntg2=2, timefp=0.25, + RNT=[[1.0, 2.0, 3.0]], + scalars={name: float(i + 1) for i, name in enumerate(SCALAR_FIELDS)}, + ) + out = _state_to_scalars(state) + assert out == { + "nrmax": 3.0, "nsamax": 1.0, "npmax": 20.0, "nthmax": 20.0, + "ntg2": 2.0, "timefp": 0.25, + **state.scalars, + } + assert len(out) == 6 + len(SCALAR_FIELDS) + # Must be a copy, not the state's own dict + out["TOTAL_IP"] = 99.0 + assert state.scalars["TOTAL_IP"] == 1.0 + + +def test_state_to_scalars_knows_fp_state_under_another_package_name(): + """``python.fplib.state.FpState`` (the repository root on sys.path) is + a different class object from ``fplib.state.FpState``: it is told by + name and module, not by identity. A class that only shares the name + is not an FpState.""" + from totlib.pipeline import _state_to_scalars + + def make(module): + cls = type("FpState", (), {"__module__": module}) + state = cls() + state.timefp = 0.5 + state.nrmax = 3 + state.scalars = {"TOTAL_IP": 1.0} + return state + + assert _state_to_scalars(make("python.fplib.state")) == { + "timefp": 0.5, "nrmax": 3.0, "TOTAL_IP": 1.0, + } + assert _state_to_scalars(make("elsewhere.state")) == {"TOTAL_IP": 1.0} + assert _state_to_scalars(make("state")) == {"TOTAL_IP": 1.0} + assert _state_to_scalars(make(None)) == {"TOTAL_IP": 1.0} + + +def test_state_to_scalars_adds_no_top_level_to_other_states(): + """A state of another module that has a ``scalars`` dict keeps exactly + that dict: its top-level counts are not pipeline scalars.""" + from dataclasses import dataclass, field + from totlib.pipeline import _state_to_scalars + + @dataclass + class TrLikeState: + nt: int = 5 + nrmax: int = 50 + nsmax: int = 2 + scalars: dict = field(default_factory=lambda: {"WPT": 1.5}) + + assert _state_to_scalars(TrLikeState()) == {"WPT": 1.5} + + def test_coupling_rules_has_fp_to_tr(): from totlib.pipeline import COUPLING_RULES rules = COUPLING_RULES.get(("fp", "tr"), []) diff --git a/python/totlib/tests/test_pipeline_equiv.py b/python/totlib/tests/test_pipeline_equiv.py index c197e2ef5..5e23b124d 100644 --- a/python/totlib/tests/test_pipeline_equiv.py +++ b/python/totlib/tests/test_pipeline_equiv.py @@ -75,7 +75,7 @@ def _baseline(): rjt_volint = compute_rjt_volint( fp_state, R0=TR_PARAMS["RR"], a=TR_PARAMS["RA"] ) - fp_scalars = _state_to_scalars(fp_state) # adapter handles missing .scalars + fp_scalars = _state_to_scalars(fp_state) # the globals, timefp and the counts fp.close() tr = Trlib() diff --git a/python/totlib/tests/test_pipeline_fp_scalars.py b/python/totlib/tests/test_pipeline_fp_scalars.py new file mode 100644 index 000000000..754a95ba1 --- /dev/null +++ b/python/totlib/tests/test_pipeline_fp_scalars.py @@ -0,0 +1,52 @@ +"""An FP step of run_pipeline reports ``timefp`` and its counts beside the +seven global scalars of ``FpState.scalars``. + +Requires libfpapi.so. Before FpState had a ``scalars`` dict the FP step's +scalars were its numeric top-level attributes; when the dict arrived they +were replaced by it, and ``timefp`` was gone from pipeline results. +""" +import os + +import pytest + +from fplib import _ffi +from fplib.state import SCALAR_FIELDS +from totlib import TotPipeline + + +def _lib_present() -> bool: + # Where the loader will look, in its order (MONO_LIB_PATH, FPLIB_PATH, + # then /fp and /lib), without loading anything while pytest + # collects. + mono = os.environ.get("MONO_LIB_PATH", "").strip() # the loader strips this one + if mono: + return os.path.exists(mono) + lib = os.environ.get("FPLIB_PATH") # and takes this one as it is + if lib: + return os.path.exists(lib) + return any(path.exists() for path in _ffi._candidate_paths()) + + +pytestmark = pytest.mark.skipif( + not _lib_present(), + reason="requires libfpapi.so — run scripts/setup.sh first", +) + +FP_TOP_LEVEL = ("nrmax", "nsamax", "npmax", "nthmax", "ntg2", "timefp") + + +def test_fp_step_scalars_hold_timefp_the_counts_and_the_globals(): + pipe = TotPipeline() + try: + result = pipe.run_pipeline([("fp", {"ntmax": 1})]) + state = pipe._ensure_module("fp").get_state() + finally: + pipe.close() + scalars = result.last("fp").scalars + assert set(scalars) == set(FP_TOP_LEVEL) | set(SCALAR_FIELDS) + assert scalars["timefp"] == state.timefp > 0.0 + for name in FP_TOP_LEVEL[:-1]: + assert scalars[name] == float(getattr(state, name)) + for name in SCALAR_FIELDS: + assert scalars[name] == state.scalars[name] + assert scalars["PLASMA_VOLUME"] > 0.0 diff --git a/python/totlib/tests/test_set_after_run.py b/python/totlib/tests/test_set_after_run.py index defbd87d9..34b41c4b9 100644 --- a/python/totlib/tests/test_set_after_run.py +++ b/python/totlib/tests/test_set_after_run.py @@ -247,6 +247,46 @@ def test_ti_reprepares_after_tot_set(self): "ti_run after tot_set_param('ti:NRMAX') must re-run ti_prep (#253)") +class _FpScalarsCase: + """fp's seven global scalars (#233) read zero after a set of tot's fp: + parameters, as after fp_set_param: the registry's count moved, so the + state is no longer the one of the last preparation.""" + + image: str = "" + + def test_fp_scalars_read_zero_after_a_tot_set(self): + from fplib._ffi import FpStateC + from fplib.state import SCALAR_FIELDS + lib = _load_image(self.image) + lib.fp_run.restype = ctypes.c_int + lib.fp_run.argtypes = [ctypes.c_int] + lib.fp_get_state.restype = ctypes.c_int + lib.fp_get_state.argtypes = [ctypes.POINTER(FpStateC)] + rc = lib.tot_init() + if rc != 0: + raise AssertionError(f"tot_init returned {rc}") + try: + st = FpStateC() + if lib.fp_run(0) != 0: + raise AssertionError("fp_run(0) failed") + # fp_get_state zeroes the scalars before it can return an error: + # a failed read must not pass for the zeros expected below + if lib.fp_get_state(ctypes.byref(st)) != 0: + raise AssertionError("fp_get_state failed before the set") + before = {name: getattr(st, name) for name in SCALAR_FIELDS} + rc = lib.tot_set_param(b"fp:DELT", 0.002) + if rc != 0: + raise AssertionError(f"tot_set_param(fp:DELT) returned {rc}") + if lib.fp_get_state(ctypes.byref(st)) != 0: + raise AssertionError("fp_get_state failed after the set") + after = {name: getattr(st, name) for name in SCALAR_FIELDS} + finally: + lib.tot_finalize() + self.assertGreater(before["PLASMA_VOLUME"], 0.0) + self.assertEqual(after, {name: 0.0 for name in SCALAR_FIELDS}, + "fp_get_state after tot_set_param('fp:DELT') must report zero scalars") + + class TestSiblingReprepOnMonoImage(_FpReprepCase, unittest.TestCase): """libtotapi_mono.so: one image, one copy of every module's state. @@ -259,6 +299,13 @@ def setUp(self): super().setUp() +class TestFpScalarsAfterTotSetOnMonoImage(_FpScalarsCase, unittest.TestCase): + def setUp(self): + _skip_unless_mono(self) + self.image = _mono_so() + super().setUp() + + class TestTiReprepOnMonoImage(_TiReprepCase, unittest.TestCase): def setUp(self): _skip_unless_mono(self)