diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml
index d82cfe6c1..d82d998a2 100644
--- a/.pre-commit-config.yaml
+++ b/.pre-commit-config.yaml
@@ -77,5 +77,5 @@ repos:
       - id: toml-no-deprecated-keys
         name: 'TOML conventions: no deprecated config keys'
         language: pygrep
-        entry: '^(mer_flag|force_wv_symmetry|ode_flag|cyl_flag|mat_flag|reform_eq_with_psilim|use_riccati|use_parallel|parallel_threads|populate_dense_xi|gal_flag|nstep|diagnose_ca|diagnose|ksing|power_bp|power_b|power_r|power_rc|use_galgrid|moment|eq_out|theta_out|xlmda_out|eqpsi_out|save_records|gpec_file|fnml_flag|ellip_flag|diag_flag|force_xialpha|fixed_boundary|output_eigenmodes|filter_modes|singular_point_method|compute_response|local_stability_flag)\s*='
+        entry: '^(mer_flag|force_wv_symmetry|ode_flag|cyl_flag|mat_flag|reform_eq_with_psilim|use_riccati|use_parallel|parallel_threads|populate_dense_xi|gal_flag|nstep|diagnose_ca|diagnose|ksing|power_bp|power_b|power_r|power_rc|use_galgrid|moment|eq_out|theta_out|xlmda_out|eqpsi_out|save_records|gpec_file|fnml_flag|ellip_flag|diag_flag|force_xialpha|fixed_boundary|output_eigenmodes|filter_modes|singular_point_method|compute_response|local_stability_flag|profile_file|profile_group)\s*='
         files: ^(examples/.*\.toml|test/test_data/.*\.toml)$
diff --git a/benchmarks/benchmark_diiid_ideal_ntv_torque.jl b/benchmarks/benchmark_diiid_ideal_ntv_torque.jl
index f2b4c1611..8cad18ab2 100644
--- a/benchmarks/benchmark_diiid_ideal_ntv_torque.jl
+++ b/benchmarks/benchmark_diiid_ideal_ntv_torque.jl
@@ -282,7 +282,7 @@ function run_benchmark(fortran_dir::String=default_fortran_dir())
     kf_ctrl = KF.KineticForcesControl(;
         fgar_flag=true, tgar_flag=tgar_enabled, nn=1, bounce_harmonic_max=4,
         ion_mass=2, ion_charge=1, nutype="harmonic", f0type="maxwellian",
-        collisionality_factor=1, psilims=[0.0, 1.0], kinetic_file=inputs.kin_file,
+        collisionality_factor=1, psilims=[0.0, 1.0],
         verbose=false)
 
     # Initialize KF internal state
@@ -309,8 +309,8 @@ function run_benchmark(fortran_dir::String=default_fortran_dir())
     _p("\n--- Torque integration ---")
     t2 = time()
 
-    kinetic_profiles = Eq.load_kinetic_profiles(
-        kf_ctrl.kinetic_file;
+    kinetic_profiles = Eq.kinetic_profile_splines(
+        Eq.read_kinetic_file(inputs.kin_file);
         zi=kf_ctrl.ion_charge, zimp=kf_ctrl.impurity_charge, mi=kf_ctrl.ion_mass, mimp=kf_ctrl.impurity_mass)
 
     kf_state = KF.KineticForcesState()
diff --git a/benchmarks/benchmark_diiid_kinetic_stability.jl b/benchmarks/benchmark_diiid_kinetic_stability.jl
index 21db14304..f576994b1 100644
--- a/benchmarks/benchmark_diiid_kinetic_stability.jl
+++ b/benchmarks/benchmark_diiid_kinetic_stability.jl
@@ -122,6 +122,7 @@ function build_benchmark_tomldir(eq_file::String, kin_file::String)
 
 [Equilibrium]
 eq_filename = "$eq_name"
+kinetic_file = "$kin_file"
 eq_type = "efit"
 jac_type = "hamada"
 grid_type = "pow1"             # Matches Fortran DIIID_kinetic_example (equil/grid.f90 powspace)
@@ -161,7 +162,6 @@ renorm_threshold = 1e4
 write_outputs_to_HDF5 = true
 
 [KineticForces]
-kinetic_file = "$kin_file"
 nn = 1
 nl = 4
 zi = 1
diff --git a/docs/src/api.md b/docs/src/api.md
index 160b9debf..854a2e274 100644
--- a/docs/src/api.md
+++ b/docs/src/api.md
@@ -47,6 +47,18 @@ The objects map one-to-one onto the pipeline stages:
     `ForceFreeStatesResult` once, and `perturbed_equilibrium(resp, rmp)` drives them with one field —
     together the `[PerturbedEquilibrium]` section.
 
+## Kinetic profiles
+
+The equilibrium carries the kinetic profiles. Name the file with the `kinetic_file` (and,
+for HDF5, `kinetic_group`) keyword, the same keys as the `[Equilibrium]` deck section, and it
+is read once into `equil.kinetic`; every kinetic stage derives its own view from it.
+[`Equilibrium.attach_kinetic_profiles!`](@ref) attaches a file or in-memory `KineticProfileData` to an
+equilibrium built without one, such as an analytic equilibrium.
+
+```julia
+eq = PlasmaEquilibrium("input.geqdsk"; jac_type="hamada", kinetic_file="kinetic.h5")
+```
+
 ## Combining forcing sources
 
 `RMPField`s form a vector space: `+`, `-` and scalar `*` build lazy linear combinations that
@@ -125,16 +137,17 @@ extraction is not implemented yet, so `GGJModel()` errors here.
 
 ```julia
 ffs  = solve(eq, Riccati(); nn=1, vac_flag=true)            # Δ′ matrix for the dispersion
-tear = solve(TearingProblem(ffs; profile_file="kin.h5"), SLAYERModel())  # coupling_mode=:coupled returns one global root instead
+tear = solve(TearingProblem(ffs), SLAYERModel())            # coupling_mode=:coupled returns one global root instead
 tear.gamma_Hz, tear.rational_q                               # root-found rate at each rational surface
 ```
 
 Keyword arguments of `TearingProblem` are the `[SLAYER]` deck section's procedure knobs
 (scan mode and Q-domain, coupling mode, critical-Δ convention, extraction filters);
-kinetic profiles come from `profile_file`, relative to the run directory.
+kinetic profiles come from the equilibrium the solve was posed on (see Kinetic profiles above).
 
-Kinetic runs (`kinetic_factor > 0`) take `Forward()` and the `"fixed"` kinetic source;
-the `"calculated"` source needs the `[KineticForces]` profiles and remains TOML-driven.
+Kinetic runs (`kinetic_factor > 0`) take `Forward()`. The `"calculated"` source reads the
+equilibrium's kinetic profiles, with the default `KineticForcesControl` species and collision
+settings.
 
 ## Local stability
 
diff --git a/docs/src/kinetic_forces.md b/docs/src/kinetic_forces.md
index e8e616fb8..0ef98bdb2 100644
--- a/docs/src/kinetic_forces.md
+++ b/docs/src/kinetic_forces.md
@@ -43,8 +43,10 @@ paths: the post-PE ψ-quadrature diagnostic and the self-consistent
 `kinetic_source = "calculated"` matrices.
 
 ```toml
-[KineticForces]
+[Equilibrium]
 kinetic_file = "kinetic.h5"    # n_i column/dataset = TOTAL main-ion density
+
+[KineticForces]
 electron = true                # add electron NTV in addition to the ion species
 impurity_charge = 6            # impurity charge (closes quasineutrality)
 impurity_mass = 12             # impurity mass
@@ -95,8 +97,8 @@ profiles and physics parameters for sensitivity studies:
 
 ### Profile-loader knobs (`density_factor`, `temperature_factor`, `ExB_rotation_factor`, `toroidal_rotation_factor`)
 
-These four knobs are applied during kinetic profile loading in
-`load_kinetic_profiles`, before any physics evaluation. The physical model is:
+These four knobs are applied when the NTV view of the kinetic profiles is built by
+`kinetic_profile_splines`, before any physics evaluation. The physical model is:
 
 ```math
 \omega_\phi = \omega_E + \omega_{*n,i} + \omega_{*T,i}
diff --git a/docs/src/stability.md b/docs/src/stability.md
index 6ef6ebb63..91f953fc1 100644
--- a/docs/src/stability.md
+++ b/docs/src/stability.md
@@ -323,12 +323,12 @@ TokaMaker's `save_eqdsk` writes such files when `lcfs_pad` (``p``) is set with t
 ``\psi = 1/(1 - p)``. Writing with `truncate_eq=False` keeps the whole plasma; EFIT-written files
 place their boundary on the X-point flux.
 
-The scan takes its kinetic profiles from the `[SLAYER]` `profile_file`, through the same control
-and loader the SLAYER analysis uses, so both see the same profiles, the same HDF5 group and the
-same χ⊥(ψ) — the layer widths cannot disagree between them. A deck without a `[SLAYER]`
-`profile_file` does not get the scan. Because the file is resolved from the run's TOML, the scan
-only runs on the `gpec.toml` entry point: the programmatic `solve(prob, alg)` path has no inputs to
-resolve it from, so the cap cannot be applied there and setting the flag warns instead. Its result is always written to
+The scan takes its kinetic profiles from the equilibrium (`[Equilibrium] kinetic_file`) through the
+same view and `[SLAYER]` control the SLAYER analysis uses, so both see the same profiles and the
+same χ⊥(ψ) — the layer widths cannot disagree between them. A deck without a `[SLAYER]` section or
+kinetic profiles does not get the scan. Because the layer settings come from the run's `[SLAYER]`
+section, the scan only runs on the `gpec.toml` entry point: the programmatic `solve(prob, alg)` path
+has no such section, so the cap cannot be applied there and setting the flag warns instead. Its result is always written to
 `ForceFreeStates/LayerOverlap/`, but it only constrains the domain when
 `psilim_from_layer_overlap = true`. It then caps `qlim` before the `rational_offset_fraction` step, so `rational_offset_fraction`
 still selects the final surface from inside the cap. A cap beyond `psihigh` has no effect.
diff --git a/examples/DIIID-like_SLAYER_example/gpec.toml b/examples/DIIID-like_SLAYER_example/gpec.toml
index 3b033d2f8..1c01cb7ba 100644
--- a/examples/DIIID-like_SLAYER_example/gpec.toml
+++ b/examples/DIIID-like_SLAYER_example/gpec.toml
@@ -6,6 +6,7 @@
 
 [Equilibrium]
 eq_filename = "../DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode.geqdsk" # Path to equilibrium file
+kinetic_file = "../DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode_kinetic.h5" # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 eq_type = "efit"               # Type of the input 2D equilibrium file
 jac_type = "hamada"            # Coordinate system (hamada, pest, boozer, equal_arc)
 grid_type = "log_asymptotic"   # Radial grid packing type
@@ -93,9 +94,8 @@ filter_above_poles      = true       # Discard roots above the highest pole γ
 filter_outside_re       = true       # Condition the above-pole filter on the +γ step exiting the Re(Δ)=0 contour loop
 polish_roots            = true       # refine each root to the true zero; enables the validity gate
 store_scan              = false      # per-surface scans omitted from HDF5 (uncoupled has one per surface)
-# Standardized kinetic-profile file (GPEC HDF5 kinetic schema), shared with the
-# sibling ideal example; carries n_e, T_e, T_i, omega_E, chi_e (χ⊥,e), chi_i (χ⊥,i), chi_phi (χ_φ).
-profile_file = "../DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode_kinetic.h5" # Kinetic-profile file (relative to the run dir), required when SLAYER is enabled
+# Kinetic profiles come from the [Equilibrium] kinetic_file (GPEC HDF5 kinetic schema), shared with the
+# sibling ideal example; it carries n_e, T_e, T_i, omega_E, chi_e (χ⊥,e), chi_i (χ⊥,i), chi_phi (χ_φ).
 
 [SLAYER.scan_grid]
 Q_re_range = [-2.0, 2.0]       # Scan box in the normalized Q plane, Re(Q) axis
diff --git a/examples/DIIID-like_error_field_example/gpec.toml b/examples/DIIID-like_error_field_example/gpec.toml
index cd6d045ea..8ce278155 100644
--- a/examples/DIIID-like_error_field_example/gpec.toml
+++ b/examples/DIIID-like_error_field_example/gpec.toml
@@ -10,6 +10,7 @@
 
 [Equilibrium]
 eq_filename = "../DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode.geqdsk" # Path to equilibrium file
+kinetic_file = "../DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode_kinetic.h5" # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 eq_type = "efit"               # Type of the input 2D equilibrium file
 jac_type = "hamada"            # Coordinate system (hamada, pest, boozer, equal_arc, park, custom)
 grid_type = "auto"             # Radial grid packing type ("auto" = two-pass measured-curvature grid)
@@ -247,7 +248,6 @@ rotation_span_factor = 1.0              # The scan covers ±span_factor × max(|
 rotation_offset_factor = 2.0            # Rough neoclassical offset as a multiple of the ion temperature-gradient diamagnetic frequency ω_*T
 
 [KineticForces]
-kinetic_file = "../DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode_kinetic.h5"  # GPEC HDF5 kinetic schema (psi,n_i,n_e,T_i,T_e,omega_E,+chi_e,chi_phi)
 nn = 1                                  # Toroidal mode number (matches nn_low/nn_high)
 bounce_harmonic_max = 6                 # Bounce harmonics retained
 ion_charge = 1                          # Main-ion charge (deuterium)
diff --git a/examples/DIIID-like_ideal_example/gpec.toml b/examples/DIIID-like_ideal_example/gpec.toml
index d744632b1..8d056b01c 100644
--- a/examples/DIIID-like_ideal_example/gpec.toml
+++ b/examples/DIIID-like_ideal_example/gpec.toml
@@ -5,6 +5,7 @@
 
 [Equilibrium]
 eq_filename = "TkMkr_D3Dlike_Hmode.geqdsk" # Path to equilibrium file
+kinetic_file = "TkMkr_D3Dlike_Hmode_kinetic.h5" # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 eq_type = "efit"               # Type of the input 2D equilibrium file
 jac_type = "hamada"            # Coordinate system (hamada, pest, boozer, equal_arc, park, custom)
 grid_type = "auto"             # Radial grid packing type ("auto" = two-pass measured-curvature grid)
@@ -88,7 +89,6 @@ write_outputs_to_HDF5 = true            # Write ErrorFields/CoilSensitivities/ t
 verbose = true                          # Log per-coil-set progress and linearity diagnostics
 
 [KineticForces]
-kinetic_file = "TkMkr_D3Dlike_Hmode_kinetic.h5"  # GPEC HDF5 kinetic schema (psi,n_i,n_e,T_i,T_e,omega_E,+chi_e,chi_i,chi_phi). Legacy .gpeckf/.kin still readable.
 nn = 1                                  # Toroidal mode number (matches nn_low/nn_high)
 bounce_harmonic_max = 6                 # Bounce harmonics retained
 ion_charge = 1                          # Main-ion charge (deuterium)
diff --git a/examples/DIIID-like_truncation_example/gpec.toml b/examples/DIIID-like_truncation_example/gpec.toml
index 1131d9551..a6fa38d33 100644
--- a/examples/DIIID-like_truncation_example/gpec.toml
+++ b/examples/DIIID-like_truncation_example/gpec.toml
@@ -8,6 +8,7 @@
 
 [Equilibrium]
 eq_filename = "../DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode.geqdsk" # Path to equilibrium file
+kinetic_file = "../DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode_kinetic.h5" # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 eq_type = "efit"               # Type of the input 2D equilibrium file
 jac_type = "hamada"            # Coordinate system (hamada, pest, boozer, equal_arc, park, custom)
 grid_type = "auto"             # Radial grid packing type ("auto" = two-pass measured-curvature grid)
@@ -44,8 +45,7 @@ psilim_from_layer_overlap = true   # Cap psilim where adjacent resistive layers
 
 [SLAYER]
 # The growth-rate analysis is not run here; the section is present because the layer-overlap scan
-# takes its kinetic profiles, and its chi_perp_e(psi), from the SLAYER profile_file. P_perp is pinned to the
+# builds its layer widths, and its chi_perp_e(psi), from the SLAYER settings and the equilibrium's kinetic profiles. P_perp is pinned to the
 # electron channel so the cap binds: with the default D_perp the layers are ~30% narrower and do not overlap inside psi = 1.
 enabled = false                # Run the SLAYER tearing-mode analysis
 P_perp_model = "chi_perp_e"    # Perpendicular Prandtl number model (pinned; see above)
-profile_file = "../DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode_kinetic.h5" # Kinetic-profile file (relative to the run dir), required when SLAYER is enabled
diff --git a/examples/Solovev_kinetic_NTV_example/gpec.toml b/examples/Solovev_kinetic_NTV_example/gpec.toml
index 86e19fd48..767ee63a5 100644
--- a/examples/Solovev_kinetic_NTV_example/gpec.toml
+++ b/examples/Solovev_kinetic_NTV_example/gpec.toml
@@ -7,6 +7,7 @@
 
 [Equilibrium]
 eq_type = "sol"                         # Type of the input 2D equilibrium file
+kinetic_file = "kinetic.dat"            # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 jac_type = "pest"                       # Coordinate system (hamada, pest, boozer, equal_arc, park, other)
 grid_type = "core_edge_packed"          # Radial grid packing type
 psilow = 1e-4                           # Lower limit of normalized poloidal flux
@@ -93,6 +94,5 @@ b0fac = 1      # Scale toroidal field at constant β (Bt changes; shape, β fixe
 f0fac = 1      # Scale toroidal field at constant pressure (β, q change; Φ, p, Bp fixed)
 
 [KineticForces]
-kinetic_file = "kinetic.dat"            # Kinetic profile file: psi_n, n_i, n_e, T_i, T_e, omega_E columns
 fcgl_flag = true                        # Fluid Chew–Goldberger–Low energy (l = 0 fluid limit) alongside the default fgar method
 axis_validity_suppression = true        # Suppress kinetic terms inside the near-axis orbit-width validity limit
diff --git a/examples/Solovev_kinetic_calculated_example/gpec.toml b/examples/Solovev_kinetic_calculated_example/gpec.toml
index de462bb44..d5737509c 100644
--- a/examples/Solovev_kinetic_calculated_example/gpec.toml
+++ b/examples/Solovev_kinetic_calculated_example/gpec.toml
@@ -9,6 +9,7 @@
 
 [Equilibrium]
 eq_type = "sol"                         # Type of the input 2D equilibrium file
+kinetic_file = "kinetic.dat"            # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 jac_type = "pest"                       # Coordinate system (hamada, pest, boozer, equal_arc, park, custom)
 grid_type = "core_edge_packed"          # Radial grid packing type
 psilow = 1e-4                           # Lower limit of normalized flux coordinate
@@ -74,7 +75,6 @@ write_outputs_to_HDF5 = true  # Write outputs to the HDF5 file
 verbose = false               # Enable verbose logging
 
 [KineticForces]
-kinetic_file = "kinetic.dat"  # Kinetic profile file: psi_n, n_i, n_e, T_i, T_e, omega_E columns
 
 [[KineticForces.ion_species]]
 z = 1                         # Species charge (fundamental units)
diff --git a/examples/Solovev_kinetic_multiion_example/gpec.toml b/examples/Solovev_kinetic_multiion_example/gpec.toml
index 066b108b1..a95019e94 100644
--- a/examples/Solovev_kinetic_multiion_example/gpec.toml
+++ b/examples/Solovev_kinetic_multiion_example/gpec.toml
@@ -7,6 +7,7 @@
 
 [Equilibrium]
 eq_type = "sol"                         # Type of the input 2D equilibrium file
+kinetic_file = "kinetic.dat"            # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 jac_type = "pest"                       # Coordinate system (hamada, pest, boozer, equal_arc, park, other)
 grid_type = "core_edge_packed"          # Radial grid packing type
 psilow = 1e-4                           # Lower limit of normalized poloidal flux
@@ -95,7 +96,6 @@ f0fac = 1      # Scale toroidal field at constant pressure (β, q change; Φ, p,
 [KineticForces]
 # Multi-ion NTV: the analytic Solovev main ion is modelled as a 50/50 D-T mix (n_i = total),
 # so the summed NTV runs over D + T + the neutrality-closing impurity + electrons.
-kinetic_file = "kinetic.dat"  # Kinetic profile file: psi_n, n_i, n_e, T_i, T_e, omega_E columns
 electron = true               # Add electron NTV in addition to the ion species
 impurity_charge = 6           # Impurity charge (fundamental units); closes quasineutrality
 impurity_mass = 12            # Impurity mass (proton masses)
diff --git a/examples/a10_kinetic_example/gpec.toml b/examples/a10_kinetic_example/gpec.toml
index 7b63274bd..82e2957d4 100644
--- a/examples/a10_kinetic_example/gpec.toml
+++ b/examples/a10_kinetic_example/gpec.toml
@@ -8,6 +8,7 @@
 
 [Equilibrium]
 eq_filename = "fix_a100_k10_q2_bn010_prof1"  # Path to equilibrium file
+kinetic_file = "a10_prof1.txt" # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 eq_type = "efit"               # Type of the input 2D equilibrium file
 jac_type = "hamada"            # Coordinate system (hamada, pest, boozer, equal_arc, park, custom)
 grid_type = "core_edge_packed" # Radial grid packing type
@@ -57,7 +58,6 @@ write_outputs_to_HDF5 = true   # Write outputs to the HDF5 file
 force_termination = true       # Stop after force-free states — the kinetic stability eigenvalues are the result of interest
 
 [KineticForces]
-kinetic_file = "a10_prof1.txt"  # Path to the kinetic profile file (bundled)
 nn = 1                         # Toroidal mode number
 bounce_harmonic_max = 6        # Bounce-harmonic number
 ion_charge = 1                 # Ion charge (fundamental units)
diff --git a/regression-harness/src/runner.jl b/regression-harness/src/runner.jl
index 2c421654a..b41815373 100644
--- a/regression-harness/src/runner.jl
+++ b/regression-harness/src/runner.jl
@@ -33,7 +33,7 @@ With no `overrides`, the example deck runs in place (returns it untouched). With
 the deck is copied to a uniquely named throwaway sibling of the example dir and the named gpec.toml keys are
 patched there, so one shared example can serve several cases (e.g. a collisionless variant via
 `"KineticForces.nutype" => "zero"`). The copy sits beside the original so relative file
-references in the deck (e.g. `SLAYER.profile_file = "../<other_example>/..."`) still resolve.
+references in the deck (e.g. `Equilibrium.kinetic_file = "../<other_example>/..."`) still resolve.
 Override keys are dotted paths into the TOML; missing intermediate tables are created.
 Returns `(rundir, is_temp)`; the caller removes `rundir` when `is_temp`.
 """
diff --git a/src/Equilibrium/Config.jl b/src/Equilibrium/Config.jl
index 96cc59b7b..9dabee6b1 100644
--- a/src/Equilibrium/Config.jl
+++ b/src/Equilibrium/Config.jl
@@ -42,6 +42,10 @@ specified in the input.
   - `flux_surface_rtol::Float64` - Relative tolerance of the ODE solves that trace each flux surface; their absolute tolerance
     is `equil_abstol(flux_surface_rtol)`, i.e. `flux_surface_rtol` capped at `EQUIL_ABSTOL_MAX`
   - `force_termination::Bool` - Terminate after equilibrium setup (skip stability calculations)
+  - `imas_cocos::Int` - Expected COCOS convention of an IMAS `dd.equilibrium` input (11 = IMAS standard, 2 = GPEC internal)
+  - `kinetic_file::String` - Kinetic-profile file read once into `PlasmaEquilibrium.kinetic` (see [`read_kinetic_file`](@ref));
+    empty for none. Every kinetic consumer (NTV, tearing, the layer-overlap scan, grid refinement) derives its view from it.
+  - `kinetic_group::String` - HDF5 group of `kinetic_file` holding the profiles
 """
 @kwdef struct EquilibriumConfig
     eq_type::String = "efit"
@@ -75,6 +79,9 @@ specified in the input.
     # IMAS-specific: expected COCOS convention of the input dd.equilibrium (11=IMAS standard, 2=GPEC internal)
     imas_cocos::Int = 11
 
+    kinetic_file::String = ""
+    kinetic_group::String = "/"
+
     """
     Modified internal constructor that enforces self consistency within the inputs
     """
@@ -84,7 +91,7 @@ specified in the input.
     function EquilibriumConfig(eq_type, eq_filename, r0_norm, b0_norm, jac_type, _, _, _, _,
         jac_custom_power_bp, jac_custom_power_b, jac_custom_power_r, jac_custom_power_rc,
         grid_type, psilow, psihigh, mpsi, psi_accuracy, mtheta, newq0, flux_surface_rtol,
-        force_termination, imas_cocos)
+        force_termination, imas_cocos, kinetic_file, kinetic_group)
         if jac_type == "hamada"
             @info "Forcing hamada coordinate jacobian exponents: power_*"
             power_b = 0
@@ -151,7 +158,7 @@ specified in the input.
         return new(eq_type, eq_filename, r0_norm, b0_norm, jac_type, power_bp, power_b, power_r, power_rc,
             jac_custom_power_bp, jac_custom_power_b, jac_custom_power_r, jac_custom_power_rc,
             grid_type, psilow, psihigh, mpsi, psi_accuracy, mtheta, newq0, flux_surface_rtol,
-            force_termination, imas_cocos)
+            force_termination, imas_cocos, kinetic_file, kinetic_group)
     end
 end
 
@@ -185,6 +192,9 @@ function EquilibriumConfig(equil_dict::Dict{String,Any}, base_path::String="./")
     if haskey(config_data, "eq_filename") && !isempty(config_data["eq_filename"]) && !isabspath(config_data["eq_filename"])
         config_data["eq_filename"] = normpath(joinpath(base_path, config_data["eq_filename"]))
     end
+    if !isempty(get(config_data, "kinetic_file", "")) && !isabspath(config_data["kinetic_file"])
+        config_data["kinetic_file"] = normpath(joinpath(base_path, config_data["kinetic_file"]))
+    end
 
     return EquilibriumConfig(; symbolize_keys(config_data)...)
 end
diff --git a/src/Equilibrium/Construction.jl b/src/Equilibrium/Construction.jl
index c715383ce..cffbe961d 100644
--- a/src/Equilibrium/Construction.jl
+++ b/src/Equilibrium/Construction.jl
@@ -60,9 +60,11 @@ constructor ends here.
 `override_psi_nodes` bypasses the config-driven radial grid (`grid_type`/`mpsi`) and forms
 the equilibrium on the given strictly increasing ψ_N node vector, whose endpoints must match
 `psilow`/`psihigh`. Used by the two-pass auto-grid refinement through [`reform`](@ref); standalone
-callers with `mpsi=0` get the coarse first-pass grid unless they refine and re-form.
+callers with `mpsi=0` get the coarse first-pass grid unless they refine and re-form. The kinetic
+profiles named by `kinetic_file` are read and attached, unless `kinetic` passes already-read ones.
 """
-function PlasmaEquilibrium(eq_config::EquilibriumConfig, source=nothing; override_psi_nodes::Union{Nothing,Vector{Float64}}=nothing)
+function PlasmaEquilibrium(eq_config::EquilibriumConfig, source=nothing; override_psi_nodes::Union{Nothing,Vector{Float64}}=nothing,
+    kinetic::Union{Nothing,KineticProfileData}=nothing)
     eq_input = run_input(eq_config, source)
 
     plasma_equilibrium = if eq_config.eq_type == "efit_by_inversion"
@@ -77,6 +79,7 @@ function PlasmaEquilibrium(eq_config::EquilibriumConfig, source=nothing; overrid
     plasma_equilibrium.ingest = eq_input.ingest
     # Analytic configs are still mutable; copy so later edits by the caller cannot change a re-form.
     source isa AnalyticEquilibriumConfig && (plasma_equilibrium.analytic = deepcopy(source))
+    plasma_equilibrium.kinetic = kinetic === nothing ? configured_kinetic_profiles(eq_config) : kinetic
 
     equilibrium_global_parameters!(plasma_equilibrium)
     equilibrium_qfind!(plasma_equilibrium)
@@ -89,7 +92,7 @@ end
     reform(equil::PlasmaEquilibrium, psi_nodes::Vector{Float64}) -> PlasmaEquilibrium
 
 Re-form `equil` on the ψ_N nodes `psi_nodes` from the source it remembers: the stored
-`analytic` config, else its raw `ingest`, so no file is re-read. Only an equilibrium carrying
+`analytic` config, else its raw `ingest`, so no file is re-read. Its kinetic profiles carry over. Only an equilibrium carrying
 neither (one built with the positional constructor, or by a reader that captures no ingest)
 falls back to re-reading `config.eq_filename`. The node endpoints must match `psilow`/`psihigh`.
 """
@@ -103,7 +106,7 @@ function reform(equil::PlasmaEquilibrium, psi_nodes::Vector{Float64})
     else
         nothing
     end
-    return PlasmaEquilibrium(equil.config, source; override_psi_nodes=psi_nodes)
+    return PlasmaEquilibrium(equil.config, source; override_psi_nodes=psi_nodes, kinetic=equil.kinetic)
 end
 
 """
diff --git a/src/Equilibrium/Equilibrium.jl b/src/Equilibrium/Equilibrium.jl
index 6c8443787..014888cc6 100644
--- a/src/Equilibrium/Equilibrium.jl
+++ b/src/Equilibrium/Equilibrium.jl
@@ -65,8 +65,8 @@ include("Axisymmetric/GradShafranovResidual.jl")
 # --- Expose types and functions to the user ---
 export reform, EquilibriumConfig, PlasmaEquilibrium, EquilibriumParameters, AnalyticEquilibriumConfig,
     ProfileSplines, GeometryProfileSplines, compute_geometry_profiles,
-    KineticProfileSplines, load_kinetic_profiles, shift_exb_rotation,
-    KineticProfileData, read_kinetic_file, write_kinetic_h5
+    KineticProfileSplines, kinetic_profile_splines, shift_exb_rotation,
+    KineticProfileData, read_kinetic_file, write_kinetic_h5, attach_kinetic_profiles!
 export flux_surface_metric, flux_surface_area
 export wants_two_pass, refined_psi_grid, merge_mandatory_nodes, bracket_mandatory_nodes, enforce_min_spacing, implied_knot_count
 export compute_jac_delpsi, compute_sqrt_jac_delpsi, compute_sqrtamat, rootarea_to_area_weight, area_to_rootarea_weight
diff --git a/src/Equilibrium/KineticProfiles.jl b/src/Equilibrium/KineticProfiles.jl
index a73bb5d35..e30b3af83 100644
--- a/src/Equilibrium/KineticProfiles.jl
+++ b/src/Equilibrium/KineticProfiles.jl
@@ -12,7 +12,7 @@ shimming.
 using DelimitedFiles
 using HDF5
 
-# Physical constants and collisionality normalizations shared by `load_kinetic_profiles` and
+# Physical constants and collisionality normalizations shared by `kinetic_profile_splines` and
 # `resolve_ntv_species` (Krook collision frequency, Logan & Park 2013 Eq. 6).
 const _EV_J = 1.602e-19          # eV → J
 const _E_CHARGE = 1.602e-19      # elementary charge [C] (PENTRC-parity value, not CODATA)
@@ -25,11 +25,14 @@ const _ME = 9.109_1e-31          # electron mass [kg]
 # NRL-formulary Coulomb logarithm (natural log); n_e in m⁻³, T_e in J.
 _coulomb_log(ne, Te) = 17.3 - 0.5 * log(ne / 1.0e20) + 1.5 * log(Te / _KEV_J)
 
+# Names kinetic profiles in error messages by their provenance, when they carry one.
+_describe(data) = isempty(data.provenance) ? "kinetic profiles" : "kinetic profiles ($(data.provenance))"
+
 # Hard stop on an unphysical negative density (bad input or a cubic-resample overshoot).
-function _assert_nonneg_density(kinetic_file, label, arr, psi_reg)
+function _assert_nonneg_density(source, label, arr, psi_reg)
     i = findfirst(<(0), arr)
     i === nothing ||
-        error("kinetic file '$kinetic_file': negative $label = $(arr[i]) at ψ_n = $(psi_reg[i]) — unphysical")
+        error("$source: negative $label = $(arr[i]) at ψ_n = $(psi_reg[i]) — unphysical")
 end
 
 """
@@ -50,58 +53,11 @@ struct ResolvedNTVSpecies{P}
     profiles::P
 end
 
-"""
-    KineticProfileData
-
-Raw kinetic-profile columns read from a kinetic file (HDF5 or ASCII), before
-resampling and spline construction. `psi` (normalized poloidal flux) is always
-present; every profile field is `nothing` when the source omits it, so each
-consumer validates only the fields it needs. Units: densities m⁻³,
-temperatures eV, frequencies rad/s, diffusivities m²/s.
-
-| field               | meaning                                                                                                           |
-|:------------------- |:----------------------------------------------------------------------------------------------------------------- |
-| `n_i`               | main-ion density                                                                                                  |
-| `n_e`               | electron density                                                                                                  |
-| `T_i`/`T_e`         | ion / electron temperature                                                                                        |
-| `omega_E`           | ExB rotation                                                                                                      |
-| `omega_tor`         | toroidal rotation (optional)                                                                                      |
-| `chi_e`             | electron perpendicular heat diffusivity χ⊥,e                                                                      |
-| `chi_i`             | ion perpendicular heat diffusivity χ⊥,i                                                                           |
-| `chi_phi`           | toroidal momentum diffusivity χ_φ                                                                                 |
-| `zeff`              | effective charge Z_eff(ψ) supplied directly (optional); the resistive-layer analysis prefers it to its closure    |
-| `species_densities` | named per-species densities (e.g. `"n_D"`, `"n_T"`) for explicit multi-ion input; `nothing` for ASCII / no extras |
-| `provenance`        | short string recording the source file/format                                                                     |
-"""
-struct KineticProfileData
-    psi::Vector{Float64}
-    n_i::Union{Nothing,Vector{Float64}}
-    n_e::Union{Nothing,Vector{Float64}}
-    T_i::Union{Nothing,Vector{Float64}}
-    T_e::Union{Nothing,Vector{Float64}}
-    omega_E::Union{Nothing,Vector{Float64}}
-    omega_tor::Union{Nothing,Vector{Float64}}
-    chi_e::Union{Nothing,Vector{Float64}}
-    chi_i::Union{Nothing,Vector{Float64}}
-    chi_phi::Union{Nothing,Vector{Float64}}
-    zeff::Union{Nothing,Vector{Float64}}
-    species_densities::Union{Nothing,Dict{String,Vector{Float64}}}
-    provenance::String
-end
-
-function KineticProfileData(; psi, n_i=nothing, n_e=nothing, T_i=nothing, T_e=nothing,
-    omega_E=nothing, omega_tor=nothing, chi_e=nothing, chi_i=nothing, chi_phi=nothing, zeff=nothing,
-    species_densities=nothing, provenance="")
-    _v(x) = x === nothing ? nothing : Float64.(collect(x))
-    return KineticProfileData(Float64.(collect(psi)), _v(n_i), _v(n_e), _v(T_i), _v(T_e),
-        _v(omega_E), _v(omega_tor), _v(chi_e), _v(chi_i), _v(chi_phi), _v(zeff), species_densities, String(provenance))
-end
-
 const _KINETIC_H5_EXTS = (".h5", ".hdf5", ".he5")
 
-# Dataset units written into / expected from the HDF5 kinetic schema.
-const _KINETIC_H5_UNITS = Dict("psi" => "normalized poloidal flux", "n_i" => "m^-3",
-    "n_e" => "m^-3", "T_i" => "eV", "T_e" => "eV", "omega_E" => "rad/s",
+# The optional profile datasets of the HDF5 kinetic schema and their units, one row per `KineticProfileData`
+# field of the same name; the reader and writer both iterate it, so a new profile is a row here plus its field.
+const _KINETIC_H5_PROFILES = ("n_i" => "m^-3", "n_e" => "m^-3", "T_i" => "eV", "T_e" => "eV", "omega_E" => "rad/s",
     "omega_tor" => "rad/s", "chi_e" => "m^2/s", "chi_i" => "m^2/s", "chi_phi" => "m^2/s", "zeff" => "dimensionless")
 
 """
@@ -141,21 +97,39 @@ function _read_kinetic_h5(path::AbstractString; group::AbstractString="/")
         haskey(ats, "provenance") && (prov = string(read(ats["provenance"])))
         # Non-standard datasets named `n_*` (e.g. "n_D", "n_T") are named per-species density
         # profiles for explicit multi-ion input; other names are reserved for future schema fields.
-        standard = ("psi", "n_i", "n_e", "T_i", "T_e", "omega_E", "omega_tor", "chi_e", "chi_i", "chi_phi", "zeff")
         extras = Dict{String,Vector{Float64}}()
         for k in keys(g)
-            (k in standard || !startswith(k, "n_")) && continue
+            (k == "psi" || any(p -> first(p) == k, _KINETIC_H5_PROFILES) || !startswith(k, "n_")) && continue
             v = read(g[k])
             v isa AbstractArray && (extras[k] = Float64.(vec(v)))
         end
         return KineticProfileData(; psi=Float64.(vec(read(g["psi"]))),
-            n_i=rd("n_i"), n_e=rd("n_e"), T_i=rd("T_i"), T_e=rd("T_e"),
-            omega_E=rd("omega_E"), omega_tor=rd("omega_tor"),
-            chi_e=rd("chi_e"), chi_i=rd("chi_i"), chi_phi=rd("chi_phi"), zeff=rd("zeff"),
+            (Symbol(name) => rd(name) for (name, _) in _KINETIC_H5_PROFILES)...,
             species_densities=(isempty(extras) ? nothing : extras), provenance=prov)
     end
 end
 
+"""
+    attach_kinetic_profiles!(equil, path; group="/") -> equil
+    attach_kinetic_profiles!(equil, data::KineticProfileData) -> equil
+
+Attach raw kinetic profiles to `equil` as `equil.kinetic`. Every kinetic consumer (the NTV torque, the
+calculated kinetic matrices, tearing, the layer-overlap scan and the two-pass grid refinement) derives its
+own view from them, so the file is read once. Equilibria built from a config with `kinetic_file` set
+attach it on construction; use this for analytic or IMAS equilibria built without it, or in-memory data.
+"""
+function attach_kinetic_profiles!(equil::PlasmaEquilibrium, data::KineticProfileData)
+    equil.kinetic = data
+    return equil
+end
+
+attach_kinetic_profiles!(equil::PlasmaEquilibrium, path::AbstractString; group::AbstractString="/") =
+    attach_kinetic_profiles!(equil, read_kinetic_file(path; group=group))
+
+# The kinetic profiles a config names, read once at construction; `nothing` when it names none.
+configured_kinetic_profiles(config::EquilibriumConfig) =
+    isempty(config.kinetic_file) ? nothing : read_kinetic_file(config.kinetic_file; group=config.kinetic_group)
+
 """
     write_kinetic_h5(path, data; group="/", schema_version="1.0", provenance=data.provenance)
 
@@ -169,22 +143,15 @@ function write_kinetic_h5(path::AbstractString, data::KineticProfileData;
     provenance::AbstractString=data.provenance)
     h5open(path, "w") do f
         g = group == "/" ? f : create_group(f, group)
-        function put(name, v)
+        function put(name, v, units)
             v === nothing && return
             g[name] = collect(Float64, v)
-            attributes(g[name])["units"] = _KINETIC_H5_UNITS[name]
+            attributes(g[name])["units"] = units
+        end
+        put("psi", data.psi, "normalized poloidal flux")
+        for (name, units) in _KINETIC_H5_PROFILES
+            put(name, getfield(data, Symbol(name)), units)
         end
-        put("psi", data.psi)
-        put("n_i", data.n_i)
-        put("n_e", data.n_e)
-        put("T_i", data.T_i)
-        put("T_e", data.T_e)
-        put("omega_E", data.omega_E)
-        put("omega_tor", data.omega_tor)
-        put("chi_e", data.chi_e)
-        put("chi_i", data.chi_i)
-        put("chi_phi", data.chi_phi)
-        put("zeff", data.zeff)
         # Named per-species densities (multi-ion input); dataset names must be `n_*` to round-trip.
         if data.species_densities !== nothing
             for (name, v) in data.species_densities
@@ -236,7 +203,7 @@ function multi_ion_composition(zs::AbstractVector, ns::AbstractVector, ne::Real;
 end
 
 """
-    resolve_ntv_species(kinetic_file, ion_species; electron, zimp, mimp, ...) -> Vector{ResolvedNTVSpecies}
+    resolve_ntv_species(data::KineticProfileData, ion_species; electron, zimp, mimp, ...) -> Vector{ResolvedNTVSpecies}
 
 Resolve the **full NTV species set** for a multi-ion run: the main ions, the neutrality-closing
 impurity, and (if `electron`) the electrons. Returns `ResolvedNTVSpecies` descriptors — the
@@ -247,16 +214,16 @@ Each descriptor's `profiles` view carries **that species' resonant density** `n_
 and its **full-composition collision frequency** `ν_s` (as `nui_spline`); `ne/Te/ωE/Zeff/lnΛ` (and
 the electron `nue`) are shared. Main-ion density is resolved from `IonSpecies`: `fraction` ⇒
 `n_s = fraction · n_i` (the kinetic file's `n_i` = total main-ion density); `density` ⇒ an explicit
-`n_*`-named dataset in the HDF5 kinetic file (e.g. `"n_D"`), cubic-resampled to the working grid.
+`n_*`-named dataset of the HDF5 kinetic data (e.g. `"n_D"`), cubic-resampled to the working grid.
 In a mixed fraction/density list the impurity content is set by the file's n_i/n_e deficit, not by
 a fraction shortfall. Composition (Zeff, zpitch, `n_imp`) comes from `multi_ion_composition`;
 per-species `ν_s = (zpitch/3.5e17)·z_s²·n_main·lnΛ / (√m_s·(T_i)^{3/2})`. The `zpitch·n_main` field
 term is exact for z=1 main ions (D-T, D-H) and approximate otherwise (a z>1 main ion's charge
 enters only through Zeff in the impurity closure). The impurity (`zimp`,`mimp`, density = the
 quasineutrality `n_imp`) is included as its own resonant species. Reduces to
-`load_kinetic_profiles` for one z=1 main ion with `fraction=1`.
+[`kinetic_profile_splines`](@ref) for one z=1 main ion with `fraction=1`.
 """
-function resolve_ntv_species(kinetic_file::AbstractString, ion_species::AbstractVector;
+function resolve_ntv_species(data::KineticProfileData, ion_species::AbstractVector;
     electron::Bool=false, zimp::Integer=6, mimp::Integer=12,
     density_factor::Float64=1.0, temperature_factor::Float64=1.0,
     ExB_rotation_factor::Float64=1.0, toroidal_rotation_factor::Float64=1.0)
@@ -265,8 +232,8 @@ function resolve_ntv_species(kinetic_file::AbstractString, ion_species::Abstract
      toroidal_rotation_factor == 1.0) ||
         error("profile scaling factors are not yet supported together with a multi-ion ion_species list")
 
-    data = read_kinetic_file(kinetic_file)
-    _need(f, n) = f === nothing ? error("kinetic file '$kinetic_file' missing required '$n'") : f
+    source = _describe(data)
+    _need(f, n) = f === nothing ? error("$source missing required '$n'") : f
     psi_in = data.psi
     ne_in = _need(data.n_e, "n_e")
     Ti_in = _need(data.T_i, "T_i")
@@ -281,8 +248,8 @@ function resolve_ntv_species(kinetic_file::AbstractString, ion_species::Abstract
     Ti = _cubic_resample(psi_in, Ti_in, psi_reg) .* _EV_J
     Te = _cubic_resample(psi_in, Te_in, psi_reg) .* _EV_J
     omegaE = _cubic_resample(psi_in, omE_in, psi_reg)
-    _assert_nonneg_density(kinetic_file, "n_e", ne, psi_reg)
-    _assert_nonneg_density(kinetic_file, "n_i (total)", ni_total, psi_reg)
+    _assert_nonneg_density(source, "n_e", ne, psi_reg)
+    _assert_nonneg_density(source, "n_i (total)", ni_total, psi_reg)
 
     # Fraction bookkeeping: fractions are shares of the total n_i, so they can never sum above 1,
     # and when every species is fraction-specified they must account for all of n_i. (In a mixed
@@ -316,14 +283,14 @@ function resolve_ntv_species(kinetic_file::AbstractString, ion_species::Abstract
         if has_dens
             (data.species_densities !== nothing && haskey(data.species_densities, s.density)) ||
                 error(
-                    "ion_species[$si]: density profile `$(s.density)` not found in kinetic file " *
+                    "ion_species[$si]: density profile `$(s.density)` not found in the kinetic profiles " *
                     "(available: $(data.species_densities === nothing ? "none — ASCII files carry only n_i" : join(keys(data.species_densities), ", ")))"
                 )
             ns[si] = _cubic_resample(psi_in, data.species_densities[s.density], psi_reg)
         else
             ns[si] = s.fraction .* ni_total
         end
-        _assert_nonneg_density(kinetic_file, "ion_species[$si] density", ns[si], psi_reg)
+        _assert_nonneg_density(source, "ion_species[$si] density", ns[si], psi_reg)
     end
 
     # Shared composition + collisionality (natural-log Coulomb log), per grid point.
@@ -344,7 +311,7 @@ function resolve_ntv_species(kinetic_file::AbstractString, ion_species::Abstract
         nue[i] = Te[i] > 0 ? (zp / _NU_PREFAC) * ne[i] * loglam[i] / (sqrt(_ME / _MP) * (Te[i] / _KEV_J)^1.5) : 0.0
     end
     # A negative impurity density means the main ions over-neutralize (Σ z_s n_s > n_e) — unphysical.
-    _assert_nonneg_density(kinetic_file, "impurity density (n_e < Σ z_s n_s)", n_imp, psi_reg)
+    _assert_nonneg_density(source, "impurity density (n_e < Σ z_s n_s)", n_imp, psi_reg)
     # Approaching the zpitch pole at Zeff = zimp signals a dilute main-ion mix taken up by impurity.
     maximum(zeff) > 0.9 * zimp &&
         @warn "max(Zeff) = $(maximum(zeff)) is within 10% of zimp = $zimp — zpitch closure near its pole; check ion_species densities/fractions"
@@ -381,28 +348,23 @@ function resolve_ntv_species(kinetic_file::AbstractString, ion_species::Abstract
 end
 
 """
-    load_kinetic_profiles(kinetic_file::AbstractString;
-                          zi::Int=1, zimp::Int=6, mi::Int=2, mimp::Int=12,
-                          density_factor::Float64=1.0, temperature_factor::Float64=1.0,
-                          ExB_rotation_factor::Float64=1.0, toroidal_rotation_factor::Float64=1.0,
-                          chi1::Union{Nothing,Float64}=nothing)
+    kinetic_profile_splines(data::KineticProfileData;
+                            zi::Int=1, zimp::Int=6, mi::Int=2, mimp::Int=12,
+                            density_factor::Float64=1.0, temperature_factor::Float64=1.0,
+                            ExB_rotation_factor::Float64=1.0, toroidal_rotation_factor::Float64=1.0,
+                            chi1::Union{Nothing,Float64}=nothing)
         → KineticProfileSplines
 
-Parse a kinetic profile file (ASCII or HDF5, dispatched by `read_kinetic_file`),
-interpolate onto a regular 101-point ψ grid, optionally apply profile scaling
-knobs, derive collisional / Z_eff diagnostics, and return a
-`KineticProfileSplines` with independent named cubic splines.
-
-# Expected file format
+The NTV view of raw kinetic profiles (e.g. `equil.kinetic`): interpolate onto a regular
+101-point ψ grid, optionally apply profile scaling knobs, derive collisional / Z_eff
+diagnostics, and return a `KineticProfileSplines` with independent named cubic splines.
 
-An HDF5 file following the GPEC kinetic schema (fields read by name), or a legacy
-six-column whitespace-separated ASCII table (header rows are filtered out):
-
-    psi_n  n_i[m^-3]  n_e[m^-3]  T_i[eV]  T_e[eV]  omega_E[rad/s]
+Needs `n_e`, `T_i`, `T_e` and `omega_E`; `n_i` defaults to `n_e`. See
+[`read_kinetic_file`](@ref) for the file formats.
 
 # Arguments
 
-  - `kinetic_file`: Path to the ASCII or HDF5 kinetic profile file (see `read_kinetic_file`)
+  - `data`: Raw kinetic profiles (`KineticProfileData`)
   - `zi`, `zimp`: Main ion and impurity charge numbers
   - `mi`, `mimp`: Main ion and impurity mass numbers (in proton masses)
   - `density_factor`: Density scaling factor (applied to ni, ne)
@@ -427,17 +389,15 @@ Collisionality is recomputed from the (possibly scaled) profiles. This differs f
 Fortran PENTRC, which computes collisionality from unscaled profiles. Use `nufac`
 (in `KineticForcesControl`) for independent collisionality scaling.
 """
-function load_kinetic_profiles(kinetic_file::AbstractString;
+function kinetic_profile_splines(data::KineticProfileData;
     zi::Int=1, zimp::Int=6, mi::Int=2, mimp::Int=12,
     density_factor::Float64=1.0, temperature_factor::Float64=1.0,
     ExB_rotation_factor::Float64=1.0, toroidal_rotation_factor::Float64=1.0,
     chi1::Union{Nothing,Float64}=nothing)
 
-    data = read_kinetic_file(kinetic_file)
-
     # NTV requires n_e, T_i, T_e, omega_E; n_i defaults to n_e (quasineutrality)
     # when absent. omega_tor / chi_e / chi_phi (if present) are ignored here.
-    _need(field, name) = field === nothing ? error("kinetic file '$kinetic_file' missing required '$name'") : field
+    _need(field, name) = field === nothing ? error("$(_describe(data)) missing required '$name'") : field
     psi_input = data.psi
     ne_input = _need(data.n_e, "n_e")
     Ti_input_eV = _need(data.T_i, "T_i")
diff --git a/src/Equilibrium/PlasmaEquilibrium.jl b/src/Equilibrium/PlasmaEquilibrium.jl
index 9c69d3218..fa02ffa36 100644
--- a/src/Equilibrium/PlasmaEquilibrium.jl
+++ b/src/Equilibrium/PlasmaEquilibrium.jl
@@ -125,6 +125,53 @@ function GeometryProfileSplines(xs::Vector{Float64},
     )
 end
 
+"""
+    KineticProfileData
+
+Raw kinetic-profile columns read from a kinetic file (HDF5 or ASCII), before
+resampling and spline construction. `psi` (normalized poloidal flux) is always
+present; every profile field is `nothing` when the source omits it, so each
+consumer validates only the fields it needs. Units: densities m⁻³,
+temperatures eV, frequencies rad/s, diffusivities m²/s.
+
+| field               | meaning                                                                                                           |
+|:------------------- |:----------------------------------------------------------------------------------------------------------------- |
+| `n_i`               | main-ion density                                                                                                  |
+| `n_e`               | electron density                                                                                                  |
+| `T_i`/`T_e`         | ion / electron temperature                                                                                        |
+| `omega_E`           | ExB rotation                                                                                                      |
+| `omega_tor`         | toroidal rotation (optional)                                                                                      |
+| `chi_e`             | electron perpendicular heat diffusivity χ⊥,e                                                                      |
+| `chi_i`             | ion perpendicular heat diffusivity χ⊥,i                                                                           |
+| `chi_phi`           | toroidal momentum diffusivity χ_φ                                                                                 |
+| `zeff`              | effective charge Z_eff(ψ) supplied directly (optional); the resistive-layer analysis prefers it to its closure    |
+| `species_densities` | named per-species densities (e.g. `"n_D"`, `"n_T"`) for explicit multi-ion input; `nothing` for ASCII / no extras |
+| `provenance`        | short string recording the source file/format                                                                     |
+"""
+struct KineticProfileData
+    psi::Vector{Float64}
+    n_i::Union{Nothing,Vector{Float64}}
+    n_e::Union{Nothing,Vector{Float64}}
+    T_i::Union{Nothing,Vector{Float64}}
+    T_e::Union{Nothing,Vector{Float64}}
+    omega_E::Union{Nothing,Vector{Float64}}
+    omega_tor::Union{Nothing,Vector{Float64}}
+    chi_e::Union{Nothing,Vector{Float64}}
+    chi_i::Union{Nothing,Vector{Float64}}
+    chi_phi::Union{Nothing,Vector{Float64}}
+    zeff::Union{Nothing,Vector{Float64}}
+    species_densities::Union{Nothing,Dict{String,Vector{Float64}}}
+    provenance::String
+end
+
+function KineticProfileData(; psi, n_i=nothing, n_e=nothing, T_i=nothing, T_e=nothing,
+    omega_E=nothing, omega_tor=nothing, chi_e=nothing, chi_i=nothing, chi_phi=nothing, zeff=nothing,
+    species_densities=nothing, provenance="")
+    _v(x) = x === nothing ? nothing : Float64.(collect(x))
+    return KineticProfileData(Float64.(collect(psi)), _v(n_i), _v(n_e), _v(T_i), _v(T_e),
+        _v(omega_E), _v(omega_tor), _v(chi_e), _v(chi_i), _v(chi_phi), _v(zeff), species_densities, String(provenance))
+end
+
 """
     KineticProfileSplines
 
@@ -255,6 +302,8 @@ This object provides a complete representation of the processed plasma equilibri
     (`sol`, `lar`, `tj_analytic`, `tj_analytic_direct`) was solved from, so [`reform`](@ref) can
     re-solve it on a new grid; `nothing` for file-based and IMAS equilibria. A defensive copy, so
     later edits to the caller's config do not reach it.
+  - `kinetic::Union{Nothing,KineticProfileData}`: the raw kinetic profiles of `config.kinetic_file` (or of
+    [`attach_kinetic_profiles!`](@ref)), on the file's own ψ grid; `nothing` when none are attached
 """
 mutable struct PlasmaEquilibrium{P<:ProfileSplines,G<:GeometryProfileSplines,I2D<:FastInterpolations.CubicInterpolantND}
     config::EquilibriumConfig
@@ -283,13 +332,14 @@ mutable struct PlasmaEquilibrium{P<:ProfileSplines,G<:GeometryProfileSplines,I2D
 
     ingest::EquilibriumIngest
     analytic::Union{Nothing,AnalyticEquilibriumConfig}
+    kinetic::Union{Nothing,KineticProfileData}
 end
 
-# Solvers build the equilibrium before PlasmaEquilibrium records its source, so allow
-# construction without it; ingest and analytic default to nothing and are assigned post-construction.
+# Solvers build the equilibrium before PlasmaEquilibrium records its source and attaches the kinetic
+# profiles, so allow construction without them; ingest, analytic and kinetic default to nothing.
 PlasmaEquilibrium(config, params, profiles, geometry, rzphi_xs, rzphi_ys,
     rzphi_rsquared, rzphi_offset, rzphi_nu, rzphi_jac,
     eqfun_B, eqfun_metric1, eqfun_metric2, ro, zo, psio) =
     PlasmaEquilibrium(config, params, profiles, geometry, rzphi_xs, rzphi_ys,
         rzphi_rsquared, rzphi_offset, rzphi_nu, rzphi_jac,
-        eqfun_B, eqfun_metric1, eqfun_metric2, ro, zo, psio, nothing, nothing)
+        eqfun_B, eqfun_metric1, eqfun_metric2, ro, zo, psio, nothing, nothing, nothing)
diff --git a/src/ForceFreeStates/CoreTypes.jl b/src/ForceFreeStates/CoreTypes.jl
index 894ecc7f2..85ba6881f 100644
--- a/src/ForceFreeStates/CoreTypes.jl
+++ b/src/ForceFreeStates/CoreTypes.jl
@@ -133,7 +133,7 @@ gpec.toml.
   - `numunorms_init::Int` - Initial array size for solution normalization data
   - `singfac_min::Float64` - Fractional distance from rational q at which ideal jump condition is enforced. The default 1e-4 suits the Riccati Δ′, whose asymptotic matching gains a spurious imaginary part as it shrinks; `"forward"` runs feeding PerturbedEquilibrium need about 1e-6 for the singular coupling to converge at the outer rational surfaces.
   - `truncate_at_rational_offset::Bool` - Truncate the integration domain at `(last_rational_q + rational_offset_fraction) / n` rather than at `qhigh` / `psihigh`. Fortran STRIDE found that truncating ~20 % above the outermost rational (`rational_offset_fraction = 0.2`) avoids a numerical kink instability in δW that appears when the integration ends too close to or just below a rational surface. **For diverted equilibria where q → ∞ at the separatrix** (e.g. DIII-D geqdsks, the bulk of production use) this costs negligible physical domain because rationals get arbitrarily dense near the LCFS — `truncate_at_rational_offset = true` is the safe and recommended default. **For limited circular / analytical equilibria with finite q at the edge** (Solovev, LAR scans), rationals are sparse and 20 % above the last rational chops off too much edge, so set `truncate_at_rational_offset = false` and let `qhigh` / `psihigh` control the truncation. Multi-`n` runs are not supported by this truncation (the "outermost rational + rational_offset_fraction / n" depends on which `n`); when `truncate_at_rational_offset = true` with `nn_low != nn_high`, `sing_lim!` warns and falls back to `qhigh` / `psihigh`. Default `true`.
-  - `psilim_from_layer_overlap::Bool` - Cap the integration domain where adjacent rational surfaces' resistive layers overlap, after Fitzpatrick, Nucl. Fusion 2025 Sect. 5.9: past that point no surface retains a well-separated inner region. Applied as an upper bound on `qlim` at the top of `sing_lim!`, so `rational_offset_fraction` / `qhigh` still select the final surface from inside it. The scan runs whenever the run's `[SLAYER]` section names a readable `profile_file` — the same profiles and χ⊥(ψ) the SLAYER analysis uses — and is recorded under `ForceFreeStates/LayerOverlap/` either way; this flag only decides whether it constrains the domain. It resolves that file from the run's TOML, so the programmatic `solve(prob, alg)` entry point cannot run it and the flag is inert there (with a warning). Off on every shipped deck; with the shipped DIII-D-like profiles the overlap point (ψ ≈ 0.9985) lies outside the 0.995 decks' domain but inside the SLAYER deck's 0.9995. See the integration-domain truncation section of the stability documentation. Default `false`.
+  - `psilim_from_layer_overlap::Bool` - Cap the integration domain where adjacent rational surfaces' resistive layers overlap, after Fitzpatrick, Nucl. Fusion 2025 Sect. 5.9: past that point no surface retains a well-separated inner region. Applied as an upper bound on `qlim` at the top of `sing_lim!`, so `rational_offset_fraction` / `qhigh` still select the final surface from inside it. The scan runs whenever the run has a `[SLAYER]` section and the equilibrium carries kinetic profiles — the same profiles and χ⊥(ψ) the SLAYER analysis uses — and is recorded under `ForceFreeStates/LayerOverlap/` either way; this flag only decides whether it constrains the domain. It resolves that file from the run's TOML, so the programmatic `solve(prob, alg)` entry point cannot run it and the flag is inert there (with a warning). Off on every shipped deck; with the shipped DIII-D-like profiles the overlap point (ψ ≈ 0.9985) lies outside the 0.995 decks' domain but inside the SLAYER deck's 0.9995. See the integration-domain truncation section of the stability documentation. Default `false`.
   - `rational_offset_fraction::Float64` - Distance beyond last rational surface (normalised ∈ [0,1) in units of 1/n). Only used when `truncate_at_rational_offset` is true. Fortran STRIDE convention is 0.2 (truncate 20 % of one rational-surface spacing above the last surface), retained here.
   - `frobenius_order::Int` - Order of singular layer (Frobenius) expansion at rational surfaces. Default 6 (Fortran STRIDE convention for Δ' calculations; lower values trade accuracy for speed).
   - `qhigh::Float64` - Integration terminated at q limit determined by minimum of qhigh and qa from equil
diff --git a/src/GeneralizedPerturbedEquilibrium.jl b/src/GeneralizedPerturbedEquilibrium.jl
index 279257593..873f112e8 100755
--- a/src/GeneralizedPerturbedEquilibrium.jl
+++ b/src/GeneralizedPerturbedEquilibrium.jl
@@ -100,10 +100,14 @@ const _DEPRECATED_FFS_KEYS = ("mer_flag", "force_wv_symmetry", "ode_flag", "cyl_
     "gal_flag", "nstep", "diagnose_ca", "diagnose", "ksing", "local_stability_flag")
 const _DEPRECATED_EQUIL_KEYS = ("power_bp", "power_b", "power_r", "power_rc", "use_galgrid")
 const _DEPRECATED_KF_KEYS = ("moment", "eq_out", "theta_out", "xlmda_out", "eqpsi_out", "save_records", "gpec_file",
-    "fnml_flag", "ellip_flag", "diag_flag", "force_xialpha")
+    "fnml_flag", "ellip_flag", "diag_flag", "force_xialpha", "kinetic_file")
 const _DEPRECATED_PE_KEYS = ("fixed_boundary", "output_eigenmodes", "filter_modes", "singular_point_method", "compute_response")
+const _DEPRECATED_SLAYER_KEYS = ("profile_file", "profile_group")
 const _DEPRECATED_KEYS = ("Equilibrium" => _DEPRECATED_EQUIL_KEYS, "ForceFreeStates" => _DEPRECATED_FFS_KEYS,
-    "KineticForces" => _DEPRECATED_KF_KEYS, "PerturbedEquilibrium" => _DEPRECATED_PE_KEYS)
+    "KineticForces" => _DEPRECATED_KF_KEYS, "PerturbedEquilibrium" => _DEPRECATED_PE_KEYS, "SLAYER" => _DEPRECATED_SLAYER_KEYS)
+# The kinetic-profile file is named once, in [Equilibrium]; these keys stop with that remedy rather than the generic one.
+const _MOVED_KINETIC_KEYS = ("KineticForces" => ("kinetic_file",), "SLAYER" => _DEPRECATED_SLAYER_KEYS)
+const _MOVED_KINETIC_REMEDY = "name the kinetic-profile file once, as kinetic_file (and kinetic_group) in [Equilibrium]"
 
 # Stop the run when a parsed gpec.toml section sets keys that no longer exist, so the deck is updated.
 function _stop_on_deprecated_keys(table, deprecated_keys, section::String; remedy::String="remove from gpec.toml")
@@ -236,6 +240,9 @@ function main_from_inputs(
     intr = ForceFreeStatesInternal(; dir_path=path)
     ffs_table = inputs["ForceFreeStates"]
     _stop_on_deprecated_keys(ffs_table, _DEPRECATED_FFS_KEYS, "ForceFreeStates")
+    for (section, keys) in _MOVED_KINETIC_KEYS
+        haskey(inputs, section) && _stop_on_deprecated_keys(inputs[section], keys, section; remedy=_MOVED_KINETIC_REMEDY)
+    end
     haskey(inputs, "KineticForces") && _stop_on_deprecated_keys(inputs["KineticForces"], _DEPRECATED_KF_KEYS, "KineticForces")
     haskey(inputs, "PerturbedEquilibrium") && _stop_on_deprecated_keys(inputs["PerturbedEquilibrium"], _DEPRECATED_PE_KEYS, "PerturbedEquilibrium")
     ctrl = ForceFreeStatesControl(; (Symbol(k) => v for (k, v) in ffs_table)...)
@@ -251,8 +258,8 @@ function main_from_inputs(
               "nodes; regenerate the file in double precision for converged Δ′." file = eq_config.eq_filename
     end
 
-    kf_ctrl, kinetic_profiles, kf_species = load_kinetic_context(inputs, intr, ctrl, equil)
-    equil = maybe_reform_equilibrium(equil, intr, ctrl, kinetic_profiles; kf_ctrl=kf_ctrl, kf_species=kf_species)
+    kf_ctrl = kinetic_forces_control(inputs, ctrl, equil)
+    equil = maybe_reform_equilibrium(equil, intr, ctrl; kf_ctrl=kf_ctrl)
 
     equil_dt = time() - equil_start
     push!(runtimes, "equilibrium" => equil_dt)
@@ -294,10 +301,9 @@ function main_from_inputs(
     @info "\n  Force-Free States\n$_SECTION"
     ffs_start = time()
 
-    # _overlap_slayer_inputs also rejects an enabled [SLAYER] section's unusable profile_file, so keep it ahead of the solve.
-    metric, mats, layer_overlap = prepare_force_free_states!(intr, ctrl, equil, kf_ctrl, kinetic_profiles;
-        species=kf_species,
-        overlap_inputs=_overlap_slayer_inputs(inputs, intr.dir_path))
+    # _overlap_slayer_inputs also rejects an enabled [SLAYER] section without usable kinetic profiles, so keep it ahead of the solve.
+    metric, mats, layer_overlap = prepare_force_free_states!(intr, ctrl, equil, kf_ctrl;
+        overlap_inputs=_overlap_slayer_inputs(inputs, equil.kinetic))
     ffs_result = run_deck_match(run_force_free_states(ctrl, equil, mats, intr, metric), ctrl)
 
     if ctrl.write_outputs_to_HDF5
@@ -313,11 +319,11 @@ function main_from_inputs(
         @info "Results written to $(ctrl.HDF5_filename)"
     end
 
-    # Validity describes the kinetic model, not the torque contraction, so it is written wherever
-    # kinetic profiles were loaded — including FFS-only runs that never reach run_kinetic_forces.
-    if ctrl.write_outputs_to_HDF5 && kinetic_profiles !== nothing
+    # Validity describes the kinetic model, not the torque contraction, so it is written for every
+    # run that reads the kinetic profiles through KineticForces — including FFS-only runs.
+    if ctrl.write_outputs_to_HDF5 && _needs_kinetic_profiles(inputs, ctrl)
         h5open(joinpath(ffs_result.dir_path, ctrl.HDF5_filename), "cw") do h5file
-            KineticForces.write_validity!(h5file, kf_ctrl, kf_species, kinetic_profiles, equil,
+            KineticForces.write_validity!(h5file, kf_ctrl, KineticForces.ntv_species(kf_ctrl, equil), KineticForces.ntv_profiles(kf_ctrl, equil), equil,
                 Float64[sng.psifac for sng in intr.sing])
         end
     end
@@ -351,7 +357,7 @@ function main_from_inputs(
 
     if run_response && "KineticForces" in keys(inputs)
         _timed_stage(runtimes, "kinetic_forces", "KineticForces") do
-            run_kinetic_forces(inputs, ffs_result, pe_state, kf_ctrl, kinetic_profiles, kf_species)
+            run_kinetic_forces(inputs, ffs_result, pe_state, kf_ctrl)
         end
         # Mirrors the write gate inside `run_kinetic_forces` (needs a PE state to contract against).
         if pe_state !== nothing && kf_ctrl.write_outputs_to_HDF5
@@ -362,7 +368,7 @@ function main_from_inputs(
     coil_sensitivities, monte_carlo, locking_risk, efc_couplings = nothing, nothing, nothing, nothing
     if run_response && "ErrorFields" in keys(inputs)
         coil_sensitivities, monte_carlo, locking_risk, efc_couplings = _timed_stage(runtimes, "error_fields", "ErrorFields") do
-            run_error_fields(inputs, ffs_result, pe_state, preloaded_coil_sets, kf_ctrl, kinetic_profiles)
+            run_error_fields(inputs, ffs_result, pe_state, preloaded_coil_sets, kf_ctrl)
         end
         # Mirrors the write gate inside `run_error_fields` (write flag defaults to true, filename
         # falls back to the ForceFreeStates file since the result carries this same control).
@@ -454,79 +460,48 @@ function resolve_mode_space!(intr::ForceFreeStatesInternal, ctrl::ForceFreeState
 end
 
 """
-    load_kinetic_context(inputs, intr, ctrl, equil) -> (kf_ctrl, kinetic_profiles, species)
+    kinetic_forces_control(inputs, ctrl, equil) -> kf_ctrl
 
-Build the KineticForces control and load the kinetic profiles once for the whole run.
-`kinetic_profiles` is `nothing` when no stage asks for them.
+The `[KineticForces]` control of the run, defaults when the deck has no such section. A run that
+needs kinetic profiles (a `[KineticForces]` section, or the calculated kinetic source) errors
+here when the equilibrium carries none; every consumer derives its view from `equil.kinetic`.
 """
-function load_kinetic_context(
-    inputs::Dict{String,Any},
-    intr::ForceFreeStatesInternal,
-    ctrl::ForceFreeStatesControl,
-    equil::Equilibrium.PlasmaEquilibrium
-)
-    # The profiles are reused by the grid refinement, the stability kinetic callback (via
-    # `calculated_cb`), and the post-PE torque diagnostics block. The `"fixed"` kinetic source
-    # path in stability does not need kinetic_profiles, but the post-PE block always does, so we
-    # load whenever a [KineticForces] section is present or the stability path requests the
-    # calculated source. psio is invariant across grid re-formation.
+function kinetic_forces_control(inputs::Dict{String,Any}, ctrl::ForceFreeStatesControl, equil::Equilibrium.PlasmaEquilibrium)
     kf_ctrl =
         haskey(inputs, "KineticForces") ?
-        KineticForces.KineticForcesControl(;
-            (Symbol(k) => v for (k, v) in inputs["KineticForces"])...) :
+        KineticForces.KineticForcesControl(; (Symbol(k) => v for (k, v) in inputs["KineticForces"])...) :
         KineticForces.KineticForcesControl()
-
-    kinetic_profiles = nothing
-    species = nothing
-    needs_kinetic_profiles = haskey(inputs, "KineticForces") ||
-                             (ctrl.kinetic_factor > 0 && ctrl.kinetic_source == "calculated")
-    if needs_kinetic_profiles
-        kinetic_file = joinpath(intr.dir_path, kf_ctrl.kinetic_file)
-        kinetic_profiles = Equilibrium.load_kinetic_profiles(
-            kinetic_file;
-            zi=kf_ctrl.ion_charge, zimp=kf_ctrl.impurity_charge,
-            mi=kf_ctrl.ion_mass, mimp=kf_ctrl.impurity_mass,
-            density_factor=kf_ctrl.density_factor, temperature_factor=kf_ctrl.temperature_factor,
-            ExB_rotation_factor=kf_ctrl.ExB_rotation_factor, toroidal_rotation_factor=kf_ctrl.toroidal_rotation_factor,
-            chi1=2π * equil.psio)
-        if !isempty(kf_ctrl.ion_species) || kf_ctrl.electron
-            speclist = isempty(kf_ctrl.ion_species) ?
-                       [KineticForces.IonSpecies(; z=kf_ctrl.ion_charge, m=kf_ctrl.ion_mass, fraction=1.0)] : kf_ctrl.ion_species
-            species = Equilibrium.resolve_ntv_species(kinetic_file, speclist;
-                electron=kf_ctrl.electron, zimp=kf_ctrl.impurity_charge, mimp=kf_ctrl.impurity_mass,
-                density_factor=kf_ctrl.density_factor, temperature_factor=kf_ctrl.temperature_factor,
-                ExB_rotation_factor=kf_ctrl.ExB_rotation_factor,
-                toroidal_rotation_factor=kf_ctrl.toroidal_rotation_factor)
-        end
-    end
-
-    return kf_ctrl, kinetic_profiles, species
+    _needs_kinetic_profiles(inputs, ctrl) && equil.kinetic === nothing &&
+        error("this run needs kinetic profiles ([KineticForces], or kinetic_source = \"calculated\"): set kinetic_file in [Equilibrium]")
+    return kf_ctrl
 end
 
+# Whether a run reads the kinetic profiles through KineticForces (its stage, or the calculated EL source).
+_needs_kinetic_profiles(inputs, ctrl) = haskey(inputs, "KineticForces") || (ctrl.kinetic_factor > 0 && ctrl.kinetic_source == "calculated")
+
 """
-    maybe_reform_equilibrium(equil, intr, ctrl, kinetic_profiles; kf_ctrl=nothing, kf_species=nothing) -> equil
+    maybe_reform_equilibrium(equil, intr, ctrl; kf_ctrl=KineticForcesControl()) -> equil
 
-Two-pass auto grid: measure the pass-1 equilibrium's curvature (profiles, geometry, kinetic
-profiles), pin knots on rational surfaces, and re-form on the refined grid with
-`Equilibrium.reform` from the source the equilibrium remembers. Returns `equil` untouched when
-its configuration wants a single pass.
+Two-pass auto grid: measure the pass-1 equilibrium's curvature (profiles, geometry, and the NTV
+view of its attached kinetic profiles), pin knots on rational surfaces, and re-form on the refined
+grid with `Equilibrium.reform` from the source the equilibrium remembers. Returns `equil` untouched
+when its configuration wants a single pass.
 """
 function maybe_reform_equilibrium(
     equil::Equilibrium.PlasmaEquilibrium,
     intr::ForceFreeStatesInternal,
-    ctrl::ForceFreeStatesControl,
-    kinetic_profiles;
-    kf_ctrl::Union{Nothing,KineticForces.KineticForcesControl}=nothing,
-    kf_species=nothing
+    ctrl::ForceFreeStatesControl;
+    kf_ctrl::KineticForces.KineticForcesControl=KineticForces.KineticForcesControl()
 )
     eq_config = equil.config
     Equilibrium.wants_two_pass(eq_config) || return equil
+    kinetic_profiles = KineticForces.ntv_profiles(kf_ctrl, equil)
 
     mandatory = ForceFreeStates.rational_psi_nodes(equil; nlow=intr.nlow, nhigh=intr.nhigh)
     # Kinetic-resonance nodes are located here rather than at the call site so the validity
     # boundary that filters them is cleared against this pass's rationals (intr.sing is empty yet).
-    pinned = kf_ctrl === nothing ? Float64[] :
-             KineticForces.resonance_grid_nodes(ctrl, kf_ctrl, kinetic_profiles, kf_species, equil, intr; rationals=mandatory)
+    pinned = KineticForces.resonance_grid_nodes(ctrl, kf_ctrl, kinetic_profiles, KineticForces.ntv_species(kf_ctrl, equil), equil, intr;
+        rationals=mandatory)
     # Smallest |n| in the run sets the widest matching half-stencil dpsi = singfac_min/(n_min·|q′|),
     # so the rational-surface brackets clear a zone large enough for every mode.
     n_min = minimum(abs(n) for n in intr.nlow:intr.nhigh if n != 0)
@@ -563,24 +538,27 @@ function snapshot_forcing_modes(inputs::Dict{String,Any}, path::String, ctrl::Fo
     return modes, normalization
 end
 
-# Kinetic profiles for the overlap scan come from the [SLAYER] section, through the same control
-# and loader the SLAYER analysis itself uses, so the scan sees the same profiles, the same HDF5
-# group and the same chi(psi). Returns `nothing` when there is no usable [SLAYER] profile_file; with
-# tearing enabled that is an error instead, so the run fails before the solve rather than at its end.
-function _overlap_slayer_inputs(inputs, dir_path::AbstractString)
+# Inputs of the overlap scan: the SLAYER view of the equilibrium's kinetic profiles under the run's
+# [SLAYER] control, the same view the SLAYER analysis itself uses, so the scan sees the same
+# profiles and the same chi(psi). Returns `nothing` without a [SLAYER] section or usable kinetic
+# profiles; with tearing enabled that is an error instead, so the run fails before the solve.
+function _overlap_slayer_inputs(inputs, kinetic::Union{Nothing,Equilibrium.KineticProfileData})
     (inputs isa AbstractDict && haskey(inputs, "SLAYER")) || return nothing
     ctrl = Runner.slayer_control_from_toml(inputs["SLAYER"])
-    isempty(ctrl.profile_file) && !ctrl.enabled && return nothing
+    if kinetic === nothing
+        ctrl.enabled && error("SLAYER needs kinetic profiles: set kinetic_file in [Equilibrium]")
+        return nothing
+    end
     # chi(psi) and Z_eff(psi) from the file where it carries them, else the control's scalar
     # fallbacks -- the same precedence run_slayer uses. Every other plasma parameter the layer
     # widths depend on is the run's own [SLAYER] setting, so the scan and the SLAYER stage see
     # one plasma.
     loaded, chi = try
-        l = Runner._load_profiles(ctrl, dir_path)
+        l = Runner._slayer_profiles(kinetic, ctrl)
         l, Runner._slayer_chi_inputs(l, ctrl; warn=false)
     catch err
         ctrl.enabled && rethrow()
-        @warn "Layer-overlap scan skipped: could not read the [SLAYER] profile_file or the diffusivities its P_perp_model needs." exception = err
+        @warn "Layer-overlap scan skipped: the kinetic profiles do not serve the SLAYER layer or its P_perp_model's diffusivities." exception = err
         return nothing
     end
     return (; profiles=loaded.profiles,
@@ -630,7 +608,7 @@ function _layer_overlap_scan(overlap_inputs, intr::ForceFreeStatesInternal, equi
 end
 
 """
-    prepare_force_free_states!(intr, ctrl, equil, kf_ctrl, kinetic_profiles; overlap_inputs=nothing) -> (metric, mats, overlap)
+    prepare_force_free_states!(intr, ctrl, equil, kf_ctrl; overlap_inputs=nothing) -> (metric, mats, overlap)
 
 Set up the force-free-states solve on `intr`: integration limits, the surviving singular
 surfaces and their GGJ coefficients, the poloidal mode range, and the metric plus
@@ -640,9 +618,7 @@ function prepare_force_free_states!(
     intr::ForceFreeStatesInternal,
     ctrl::ForceFreeStatesControl,
     equil::Equilibrium.PlasmaEquilibrium,
-    kf_ctrl::KineticForces.KineticForcesControl,
-    kinetic_profiles;
-    species=nothing,
+    kf_ctrl::KineticForces.KineticForcesControl;
     overlap_inputs=nothing
 )
     # Resistive-layer overlap: locate where adjacent rational surfaces' layers run into each
@@ -735,6 +711,8 @@ function prepare_force_free_states!(
         # Inject the KineticForces callback so the "calculated" source can
         # invoke compute_calculated_kinetic_matrices without ForceFreeStates
         # importing KineticForces (which would invert the load order).
+        kinetic_profiles = KineticForces.ntv_profiles(kf_ctrl, equil)
+        species = KineticForces.ntv_species(kf_ctrl, equil)
         # One boundary for the whole run: the widest-orbit species, moved clear of any rational
         # surface whose window the transition band would otherwise cut through.
         axis_psi_c = KineticForces.axis_validity_boundary(kf_ctrl, species, kinetic_profiles, equil,
@@ -930,8 +908,9 @@ a `gpec.toml` run of `main` does and produces the same result object. The second
 sugar building the problem from an equilibrium and the problem keywords in one call.
 
 Knobs owned by `alg` are rejected as `ForceFreeStatesControl` keywords. Kinetic runs
-(`kinetic_factor > 0`) take `Forward()` and the `"fixed"` kinetic source; the `"calculated"`
-source needs the `[KineticForces]` profiles and stays TOML-driven this cycle.
+(`kinetic_factor > 0`) take `Forward()`; the `"calculated"` source needs kinetic profiles on the
+equilibrium (`kinetic_file`, or `attach_kinetic_profiles!`) and uses the default
+`KineticForcesControl` species and collision settings.
 
 ```julia
 eq = PlasmaEquilibrium("input.geqdsk"; jac_type="hamada")
@@ -948,8 +927,8 @@ function solve(prob::EulerLagrangeProblem, alg::ForceFreeStates.AbstractIntegrat
     ForceFreeStates._apply_alg!(ctrl_kwargs, alg)
     ctrl = ForceFreeStatesControl(; ctrl_kwargs...)
 
-    ctrl.kinetic_factor > 0 && ctrl.kinetic_source == "calculated" &&
-        error("kinetic_source = \"calculated\" needs the [KineticForces] profiles and is TOML-driven; run it through `main`")
+    ctrl.kinetic_factor > 0 && ctrl.kinetic_source == "calculated" && equil.kinetic === nothing &&
+        error("kinetic_source = \"calculated\" needs kinetic profiles on the equilibrium: set kinetic_file or call attach_kinetic_profiles!")
     ctrl.kinetic_factor > 0 && !(alg isa ForceFreeStates.Forward) &&
         error("$(nameof(typeof(alg)))() does not support kinetic runs (kinetic_factor > 0); use Forward().")
 
@@ -961,13 +940,12 @@ function solve(prob::EulerLagrangeProblem, alg::ForceFreeStates.AbstractIntegrat
     ctrl.gal_match_flag &&
         error("gal_match_flag is the deck's match request; match the result instead: solve(MatchProblem(ffs; eta, rho, rotation), GGJModel())")
 
-    # The API path never reads kinetic profiles, so the KineticForces control is only the
-    # placeholder `prepare_force_free_states!` threads into its (unused) callback.
+    # The problem has no [KineticForces] section, so the views of the attached kinetic profiles use the default control.
     kf_ctrl = KineticForces.KineticForcesControl()
 
-    equil = maybe_reform_equilibrium(equil, intr, ctrl, nothing)
+    equil = maybe_reform_equilibrium(equil, intr, ctrl)
 
-    metric, mats, _ = prepare_force_free_states!(intr, ctrl, equil, kf_ctrl, nothing)
+    metric, mats, _ = prepare_force_free_states!(intr, ctrl, equil, kf_ctrl)
     result = run_force_free_states(ctrl, equil, mats, intr, metric; cut_solution=alg isa ForceFreeStates.Galerkin && alg.cut_solution)
 
     if ctrl.write_outputs_to_HDF5
@@ -1070,24 +1048,24 @@ function perturbed_equilibrium(
 end
 
 """
-    run_kinetic_forces(inputs, result, pe_state, kf_ctrl, kinetic_profiles)
+    run_kinetic_forces(inputs, result, pe_state, kf_ctrl)
 
 Compute and write the neoclassical toroidal viscosity torque diagnostics when the deck carries a
-`[KineticForces]` section. No-op when the perturbed-equilibrium state the operators contract
+`[KineticForces]` section, from the kinetic profiles attached to `result.equil`. No-op when the perturbed-equilibrium state the operators contract
 against is missing.
 """
 function run_kinetic_forces(
     inputs::Dict{String,Any},
     result::ForceFreeStatesResult,
     pe_state,
-    kf_ctrl::KineticForces.KineticForcesControl,
-    kinetic_profiles,
-    species=nothing
+    kf_ctrl::KineticForces.KineticForcesControl
 )
     # ----------------------------------------------------------------
     # KineticForces (Neoclassical Toroidal Viscosity)
     # ----------------------------------------------------------------
     ("KineticForces" in keys(inputs)) || return nothing
+    kinetic_profiles = KineticForces.ntv_profiles(kf_ctrl, result.equil)
+    species = KineticForces.ntv_species(kf_ctrl, result.equil)
 
     # One boundary for the whole run, shared with the EL kinetic matrices and the Validity output.
     axis_psi_c = KineticForces.axis_validity_boundary(kf_ctrl, species, kinetic_profiles, result.equil,
@@ -1098,7 +1076,6 @@ function run_kinetic_forces(
     if pe_state === nothing
         @info "Skipping NTV torque diagnostics: no perturbed-equilibrium data (e.g. kinetic_source=\"calculated\")."
     else
-        # kf_ctrl and kinetic_profiles were loaded once before the equilibrium was re-formed.
         kf_intr = KineticForces.KineticForcesInternal(result.equil; verbose=kf_ctrl.verbose)
         KineticForces.set_perturbation_data!(kf_intr, pe_state, result, result.equil, result.metric)
 
@@ -1160,7 +1137,7 @@ function forcing_terms_control(inputs::Dict{String,Any})
 end
 
 """
-    run_error_fields(inputs, result, pe_state, preloaded_coil_sets, kf_ctrl, kinetic_profiles) -> (sensitivities, monte_carlo, risk, couplings) or nothing
+    run_error_fields(inputs, result, pe_state, preloaded_coil_sets, kf_ctrl) -> (sensitivities, monte_carlo, risk, couplings) or nothing
 
 Linearize every coil set's resonant drive with respect to its rigid shifts and tilts and write
 `ErrorFields/CoilSensitivities/` when the deck carries an `[ErrorFields]` section. When the
@@ -1169,7 +1146,7 @@ Carlo on the full-window dominant mode with the `[ErrorFields.MonteCarlo]` setti
 `ErrorFields/MonteCarlo/`; with an `[ErrorFields.scenario]` table as well, evaluate the locking
 risk (and the tolerance scan when `[ErrorFields.Risk]` names `scan_scales`) and write
 `ErrorFields/Risk/`; with `[ErrorFields.NTV]` naming correction arrays, evaluate their overlap and NTV
-torque couplings per kilo-ampere-turn (needs the kinetic context) and write `ErrorFields/NTV/`.
+torque couplings per kilo-ampere-turn (needs the equilibrium's kinetic profiles) and write `ErrorFields/NTV/`.
 Stages not requested return `nothing`. Needs the
 perturbed-equilibrium state's singular-coupling matrix and coil-format forcing, and errors
 otherwise: a deck asking for error-field sensitivities without them is a misconfiguration, not a
@@ -1180,8 +1157,7 @@ function run_error_fields(
     result::ForceFreeStatesResult,
     pe_state,
     preloaded_coil_sets::Union{Nothing,Vector{ForcingTerms.CoilSet}},
-    kf_ctrl::KineticForces.KineticForcesControl=KineticForces.KineticForcesControl(),
-    kinetic_profiles=nothing
+    kf_ctrl::KineticForces.KineticForcesControl=KineticForces.KineticForcesControl()
 )
     ("ErrorFields" in keys(inputs)) || return nothing
 
@@ -1258,7 +1234,7 @@ function run_error_fields(
         missing = setdiff(ntv_ctrl.efc_coils, [cs.name for cs in efc_sets])
         isempty(missing) || error("[ErrorFields.NTV] efc_coils not among the run's coil sets: $(join(missing, ", "))")
         ntv_start = time()
-        couplings = efc_couplings(pe_state.response, efc_sets, rc, dom, cfg, kf_ctrl, kinetic_profiles; method=ntv_ctrl.method, verbose=ntv_ctrl.verbose,
+        couplings = efc_couplings(pe_state.response, efc_sets, rc, dom, cfg, kf_ctrl; method=ntv_ctrl.method, verbose=ntv_ctrl.verbose,
             regularization_width=pe_state.regularization_width,
             rotation_scan=ntv_ctrl.rotation_scan, scan_points=ntv_ctrl.rotation_scan_points, scan_max_points=ntv_ctrl.rotation_scan_max_points,
             scan_tolerance=ntv_ctrl.rotation_scan_tolerance, span_factor=ntv_ctrl.rotation_span_factor, offset_factor=ntv_ctrl.rotation_offset_factor)
@@ -1282,7 +1258,7 @@ function run_error_fields(
 end
 
 """
-    efc_couplings(resp, coil_sets, rc, dom, cfg, kf_ctrl, kinetic_profiles; mode=1, method="fgar",
+    efc_couplings(resp, coil_sets, rc, dom, cfg, kf_ctrl; mode=1, method="fgar",
                   regularization_width=resp.regularization_width, verbose=false) -> Vector{ErrorFields.EFCCoupling}
 
 The couplings of each correction coil array in `coil_sets` per kilo-ampere-turn: its
@@ -1292,8 +1268,8 @@ torque is a quadratic form of the applied spectrum, so each is one evaluation of
 response `resp` on the unit-current spectrum (a spectrum-literal `RMPField`, nothing written, with the
 run's `regularization_width`) followed by the kinetic torque of `method`; both scale exactly with the
 square of the current. `cfg` is the run's `[ForcingTerms]` coil configuration, so the boundary grid
-matches the one the run's sensitivities were evaluated on; `rc` must be the coupling of `resp` and
-`kinetic_profiles` its kinetic context. The spectrum is normalized by the magnitude of the
+matches the one the run's sensitivities were evaluated on; `rc` must be the coupling of `resp`; the
+torque uses the NTV view of the kinetic profiles attached to `resp.ffs.equil`. The spectrum is normalized by the magnitude of the
 array's ampere-turns, `|nw| × max|I|`, keeping its winding sense and current pattern as the
 phase reference; the torque is stored with its sign.
 """
@@ -1303,8 +1279,7 @@ function efc_couplings(
     rc::PerturbedEquilibrium.ResonantCoupling,
     dom::PerturbedEquilibrium.DominantCoupling,
     cfg::ForcingTerms.CoilConfig,
-    kf_ctrl::KineticForces.KineticForcesControl,
-    kinetic_profiles;
+    kf_ctrl::KineticForces.KineticForcesControl;
     mode::Int=1,
     method::AbstractString="fgar",
     rotation_scan::Bool=true,
@@ -1317,7 +1292,7 @@ function efc_couplings(
     verbose::Bool=false
 )
     ffs = resp.ffs
-    kinetic_profiles === nothing && error("efc_couplings needs kinetic profiles: add a [KineticForces] section with a kinetic_file")
+    ffs.equil.kinetic === nothing && error("efc_couplings needs kinetic profiles: set kinetic_file in [Equilibrium]")
     getfield(kf_ctrl, Symbol(method * "_flag")) || error("efc_couplings: KineticForces method \"$method\" is not enabled in the control")
     equil = ffs.equil
     grids = [(n, ForcingTerms.CoilForcingGrid(equil, cfg, n; psi=ffs.psilim)) for n in sort(unique(rc.n_modes))]
@@ -1326,14 +1301,14 @@ function efc_couplings(
 
     # Reference rotation and scan span from the kinetic profiles, with the torque kernel's own
     # diamagnetic frequencies (Logan & Park 2013 Eq. 7): ω_φ = ω_E + ω_*n + ω_*T.
-    kp = kinetic_profiles
+    kp = KineticForces.ntv_profiles(kf_ctrl, equil)
     chi1 = 2π * equil.psio
     chrg = kf_ctrl.ion_charge * Utilities.E_CHG
     ψk = kp.xs
     ω_E = [kp.omegaE_spline(ψ) for ψ in ψk]
     ω_star_T = [-2π * kp.Ti_deriv(ψ) / (chrg * chi1) for ψ in ψk]
     ω_star_n = [
-        (n=kp.ni_spline(ψ); n > 0 ? -2π * kp.Ti_spline(ψ) * kp.ni_deriv(ψ) / (chrg * chi1 * n) : 0.0) for ψ in ψk
+        (n = kp.ni_spline(ψ); n > 0 ? -2π * kp.Ti_spline(ψ) * kp.ni_deriv(ψ) / (chrg * chi1 * n) : 0.0) for ψ in ψk
     ]
     ω_φ = ω_E .+ ω_star_n .+ ω_star_T
     # Momentum weight n_i·dV/dψ with ⟨R²⟩ taken as R₀², inside the plasma only.
diff --git a/src/KineticForces/CalculatedKineticMatrices.jl b/src/KineticForces/CalculatedKineticMatrices.jl
index 5ca127c33..996a18bef 100644
--- a/src/KineticForces/CalculatedKineticMatrices.jl
+++ b/src/KineticForces/CalculatedKineticMatrices.jl
@@ -43,7 +43,7 @@ interpolants (`kf_intr.dbob_m`, `kf_intr.divx_m`) are not wired into this routin
     carry NTV-specific knobs (nl, zi, mi, wdfac, electron) that the
     KineticForces kernel needs but ForceFreeStatesControl does not expose.
   - `kinetic_profiles::Equilibrium.KineticProfileSplines`: Required. Named kinetic-
-    profile splines loaded via `Equilibrium.load_kinetic_profiles`.
+    profile splines built by `Equilibrium.kinetic_profile_splines`.
 
 # Returns
 
diff --git a/src/KineticForces/KineticForcesStructs.jl b/src/KineticForces/KineticForcesStructs.jl
index ff0752f20..3359e682e 100644
--- a/src/KineticForces/KineticForcesStructs.jl
+++ b/src/KineticForces/KineticForcesStructs.jl
@@ -171,9 +171,6 @@ builds a second control for its differing tolerance).
     write_outputs_to_HDF5::Bool = true
     HDF5_filename::String = "gpec.h5"
 
-    # Input files (relative to dir_path)
-    kinetic_file::String = "kinetic.dat"
-
     # Debugging
     verbose::Bool = false
 end
diff --git a/src/KineticForces/Utils.jl b/src/KineticForces/Utils.jl
index 0c9615839..6a2cf0057 100644
--- a/src/KineticForces/Utils.jl
+++ b/src/KineticForces/Utils.jl
@@ -292,6 +292,38 @@ function clear_rational_windows(psi_c::Float64, rationals::Vector{Float64})::Flo
     return psi_c
 end
 
+"""
+    ntv_profiles(kf_ctrl, equil) → Union{Nothing,KineticProfileSplines}
+
+The NTV view of the kinetic profiles attached to `equil` (`equil.kinetic`), with the main-ion and
+impurity charges and masses and the profile scale factors of `kf_ctrl`; `nothing` when the
+equilibrium carries no kinetic profiles.
+"""
+function ntv_profiles(kf_ctrl::KineticForcesControl, equil)
+    equil.kinetic === nothing && return nothing
+    return Equilibrium.kinetic_profile_splines(equil.kinetic;
+        zi=kf_ctrl.ion_charge, zimp=kf_ctrl.impurity_charge, mi=kf_ctrl.ion_mass, mimp=kf_ctrl.impurity_mass,
+        density_factor=kf_ctrl.density_factor, temperature_factor=kf_ctrl.temperature_factor,
+        ExB_rotation_factor=kf_ctrl.ExB_rotation_factor, toroidal_rotation_factor=kf_ctrl.toroidal_rotation_factor,
+        chi1=2π * equil.psio)
+end
+
+"""
+    ntv_species(kf_ctrl, equil) → Union{Nothing,Vector{ResolvedNTVSpecies}}
+
+The resolved NTV species set of the kinetic profiles attached to `equil` (see
+`Equilibrium.resolve_ntv_species`) when `kf_ctrl` asks for a species sum (`ion_species` or
+`electron`); `nothing` for the single main-ion run or when the equilibrium carries no kinetic profiles.
+"""
+function ntv_species(kf_ctrl::KineticForcesControl, equil)
+    (equil.kinetic !== nothing && (!isempty(kf_ctrl.ion_species) || kf_ctrl.electron)) || return nothing
+    speclist = isempty(kf_ctrl.ion_species) ? [IonSpecies(; z=kf_ctrl.ion_charge, m=kf_ctrl.ion_mass, fraction=1.0)] : kf_ctrl.ion_species
+    return Equilibrium.resolve_ntv_species(equil.kinetic, speclist;
+        electron=kf_ctrl.electron, zimp=kf_ctrl.impurity_charge, mimp=kf_ctrl.impurity_mass,
+        density_factor=kf_ctrl.density_factor, temperature_factor=kf_ctrl.temperature_factor,
+        ExB_rotation_factor=kf_ctrl.ExB_rotation_factor, toroidal_rotation_factor=kf_ctrl.toroidal_rotation_factor)
+end
+
 """
     axis_validity_boundary(kf_ctrl, species, kinetic_profiles, equil) → Float64
 
diff --git a/src/Rerun.jl b/src/Rerun.jl
index 94fbfdfca..5ba420ffc 100644
--- a/src/Rerun.jl
+++ b/src/Rerun.jl
@@ -321,6 +321,14 @@ function rebuild_equilibrium_inputs(inputs::Dict{String,Any}, ingest, output_dir
     return eq_config, additional_input
 end
 
+# The stored deck's kinetic profiles, re-read from its file; a file that has since moved leaves the equilibrium without them.
+function _replay_kinetic_profiles(config::Equilibrium.EquilibriumConfig)
+    isempty(config.kinetic_file) && return nothing
+    isfile(config.kinetic_file) && return Equilibrium.configured_kinetic_profiles(config)
+    @warn "The stored deck's kinetic_file is not readable; the restored equilibrium carries no kinetic profiles" file = config.kinetic_file
+    return nothing
+end
+
 """
     equilibrium_from_h5(h5path) -> (; equil, inputs, psilim)
 
@@ -365,7 +373,7 @@ function equilibrium_from_h5(h5path::AbstractString)
         geometry = Equilibrium.compute_geometry_profiles(xs, ys, rzphi[1], rzphi[2], rzphi[4], params.ro)
         analytic = haskey(Equilibrium.ANALYTIC_EQ, config.eq_type) ? build_analytic_config(config.eq_type, inputs) : nothing
         equil = Equilibrium.PlasmaEquilibrium(config, params, profiles, geometry, xs, ys, rzphi..., eqfun..., params.ro, params.zo, params.psio,
-            read_equilibrium_ingest(f), analytic)
+            read_equilibrium_ingest(f), analytic, _replay_kinetic_profiles(config))
         return (; equil, inputs, psilim=Float64(read(f, "Info/psilim")))
     end
 end
diff --git a/src/Tearing/Runner/Control.jl b/src/Tearing/Runner/Control.jl
index 5dd2c5ccc..cf476de42 100644
--- a/src/Tearing/Runner/Control.jl
+++ b/src/Tearing/Runner/Control.jl
@@ -140,17 +140,12 @@ constructor.
 
 # Kinetic-profile source
 
-SLAYER reads kinetic profiles through the shared `Equilibrium.read_kinetic_file`
-reader (the same standardized object used by the kinetic/NTV physics), so
-there is one consistent interface for resistive and kinetic profiles.
-
-  - `profile_file`   -- path to a kinetic-profile file (relative to the run
-    dir), required when SLAYER is enabled. HDF5 (`.h5`) files use the GPEC
-    kinetic schema and may carry `chi_e` (χ⊥,e), `chi_i` (χ⊥,i) and `chi_phi` (χ_φ); ASCII
-    tables are also accepted but carry no χ. When a χ dataset is absent or
-    all-zero, the scalar `chi_perp_e`/`chi_perp_i`/`chi_tor` fallbacks below are used (so a
-    file can keep the χ keys set to 0 to defer to the scalars).
-  - `profile_group`  -- group within the HDF5 file (default `"/"`)
+SLAYER takes its kinetic profiles from the equilibrium (`[Equilibrium] kinetic_file`,
+read once into `equil.kinetic`), the same object the kinetic/NTV physics reads.
+They need `n_e`, `T_e`, `T_i` and `omega_E`. HDF5 files may also carry `chi_e` (χ⊥,e),
+`chi_i` (χ⊥,i) and `chi_phi` (χ_φ); ASCII tables carry no χ. When a χ dataset is absent or
+all-zero, the scalar `chi_perp_e`/`chi_perp_i`/`chi_tor` fallbacks below are used (so a
+file can keep the χ keys set to 0 to defer to the scalars).
 
 # Output control
 
@@ -228,9 +223,6 @@ there is one consistent interface for resistive and kinetic profiles.
 
     omega_E_kHz::Dict{String,Float64} = Dict{String,Float64}()
 
-    profile_file::String = ""
-    profile_group::String = "/"
-
     store_scan::Bool = false
 
     critical_resonant_field::CriticalResonantFieldControl = CriticalResonantFieldControl()
diff --git a/src/Tearing/Runner/Runner.jl b/src/Tearing/Runner/Runner.jl
index f02de563e..d0300b5f9 100644
--- a/src/Tearing/Runner/Runner.jl
+++ b/src/Tearing/Runner/Runner.jl
@@ -31,7 +31,7 @@ import CommonSolve
 using FastInterpolations: cubic_interp
 using ..Utilities
 using ..Utilities: KineticProfiles
-using ...Equilibrium: read_kinetic_file, KineticProfileData
+using ...Equilibrium: KineticProfileData
 using ...ForceFreeStates: ForceFreeStatesResult
 using ..InnerLayer
 using ..InnerLayer: InnerLayerParameters, InnerLayerResponse, solve_inner,
diff --git a/src/Tearing/Runner/TearingProblem.jl b/src/Tearing/Runner/TearingProblem.jl
index a5c9aeb10..4e4e1dd32 100644
--- a/src/Tearing/Runner/TearingProblem.jl
+++ b/src/Tearing/Runner/TearingProblem.jl
@@ -8,12 +8,13 @@
 Tearing-mode growth rates from a finished force-free-states solve: hold the outer Δ′ fixed
 and find the growth rate at which the inner layer matches it. Solve with
 `solve(prob, SLAYERModel())`; GGJ growth-rate extraction is not implemented yet. Keywords are
-[`SLAYERControl`](@ref) fields, the `[SLAYER]` deck knobs; kinetic profiles come from
-`profile_file`. Without a full Δ′ matrix the per-surface Δ′ stubs are used, with a warning.
+[`SLAYERControl`](@ref) fields, the `[SLAYER]` deck knobs; kinetic profiles come from the
+equilibrium (`ffs.equil.kinetic`, see `attach_kinetic_profiles!`). Without a full Δ′ matrix the
+per-surface Δ′ stubs are used, with a warning.
 
 ## Fields
 
-  - `ffs` - The finished solve; duck-typed on `equil`, `surfaces`, `delta_prime`, `dir_path`, so tests can use stand-ins.
+  - `ffs` - The finished solve; duck-typed on `equil` (with its `kinetic` profiles), `surfaces` and `delta_prime`, so tests can use stand-ins.
   - `control::SLAYERControl` - The tearing knobs.
 """
 struct TearingProblem{R}
@@ -77,7 +78,9 @@ function CommonSolve.solve(prob::TearingProblem, model::InnerLayer.InnerLayerMod
     control.enabled || return empty_slayer_result(control)
     isempty(surfaces) && return empty_slayer_result(control)
 
-    loaded = _load_profiles(control, ffs.dir_path)
+    equil.kinetic === nothing &&
+        error("TearingProblem: the equilibrium carries no kinetic profiles; set kinetic_file in [Equilibrium] or attach_kinetic_profiles!")
+    loaded = _slayer_profiles(equil.kinetic, control)
 
     params = _tearing_params(model, equil, surfaces, loaded, control)
 
diff --git a/src/Tearing/Runner/run_slayer.jl b/src/Tearing/Runner/run_slayer.jl
index 26ba0c18c..61ec80700 100644
--- a/src/Tearing/Runner/run_slayer.jl
+++ b/src/Tearing/Runner/run_slayer.jl
@@ -4,34 +4,23 @@
 # `run_slayer_from_inputs` (scan + growth-rate extraction). `TearingProblem` builds its inputs.
 
 # ---------------------------------------------------------------------
-# Profile loading
+# Kinetic-profile view
 # ---------------------------------------------------------------------
-# Read kinetic profiles through the shared `Equilibrium.read_kinetic_file`
-# reader and adapt them to the SLAYER inner layer. Returns the spline-based
-# `KineticProfiles` plus χ⊥,e(ψ)/χ⊥,i(ψ)/χ_φ(ψ) callables built from the file's
-# `chi_e`/`chi_i`/`chi_phi` (or `nothing` when the file carries no usable χ, in which
-# case the caller falls back to the scalar `control.chi_perp_e`/`chi_perp_i`/`chi_tor`), and
-# a Z_eff(ψ) callable: the file's own `zeff` profile when it carries one, else
-# derived from its `n_i/n_e` (or `nothing` when the file carries neither, in which
-# case the caller falls back to the scalar `control.zeff`), and whether the file
-# carries an `omega_E` profile.
-function _load_profiles(control::SLAYERControl, dir_path::AbstractString)
-    isempty(control.profile_file) &&
-        error("run_slayer: [SLAYER] profile_file is empty — point it at a " *
-              "kinetic-profile file (HDF5 GPEC kinetic schema or ASCII table).")
-    path = isabspath(control.profile_file) ? control.profile_file :
-           joinpath(dir_path, control.profile_file)
-    data = read_kinetic_file(path; group=control.profile_group)
-
+# The SLAYER view of raw kinetic profiles (`equil.kinetic`), on their own ψ grid. Returns the
+# spline-based `KineticProfiles` plus χ⊥,e(ψ)/χ⊥,i(ψ)/χ_φ(ψ) callables built from `chi_e`/`chi_i`/`chi_phi`
+# (or `nothing` when the data carry no usable χ, in which case the caller falls back to the scalar
+# `control.chi_perp_e`/`chi_perp_i`/`chi_tor`), a Z_eff(ψ) callable: the data's own `zeff` profile when
+# they carry one, else derived from `n_i/n_e` (or `nothing` when they carry neither, in which case the
+# caller falls back to the scalar `control.zeff`), and whether the data carry an `omega_E` profile.
+function _slayer_profiles(data::KineticProfileData, control::SLAYERControl)
+    source = isempty(data.provenance) ? "kinetic profiles" : "kinetic profiles ($(data.provenance))"
     for (name, v) in (("n_e", data.n_e), ("T_e", data.T_e), ("T_i", data.T_i))
         v === nothing &&
-            error("run_slayer: kinetic file '$path' is missing required " *
-                  "dataset '$name' for the SLAYER inner layer.")
+            error("run_slayer: $source missing required dataset '$name' for the SLAYER inner layer.")
     end
-    # `omega_E` carries the E×B rotation Ω_E (per unit n) when the file has it, else zero (no shift);
+    # `omega_E` carries the E×B rotation Ω_E (per unit n) when the data have it, else zero (no shift);
     # ω_*e/ω_*i are derived per surface from the density and temperature splines by `build_slayer_inputs`.
-    npsi = length(data.psi)
-    omega_E = data.omega_E === nothing ? zeros(npsi) : data.omega_E
+    omega_E = data.omega_E === nothing ? zeros(length(data.psi)) : data.omega_E
     profiles = KineticProfiles(; psi=data.psi, n_e=data.n_e, T_e=data.T_e,
         T_i=data.T_i, omega_E=omega_E)
 
@@ -43,7 +32,7 @@ function _load_profiles(control::SLAYERControl, dir_path::AbstractString)
     psi_xs = collect(Float64, data.psi)
     function _chi_spline(name, v)
         (v === nothing || all(iszero, v)) && return nothing
-        minimum(v) > 0 || error("run_slayer: kinetic file '$path' has $name with minimum $(minimum(v)) ≤ 0; diffusivities must be positive")
+        minimum(v) > 0 || error("run_slayer: $source has $name with minimum $(minimum(v)) ≤ 0; diffusivities must be positive")
         return cubic_interp(psi_xs, collect(Float64, v))
     end
     model = control.P_perp_model
@@ -58,14 +47,14 @@ function _load_profiles(control::SLAYERControl, dir_path::AbstractString)
     _usable(v) = v !== nothing && !all(iszero, v)
     zeff = if _usable(data.zeff)
         minimum(data.zeff) >= 1.0 - 1e-8 ||
-            error("run_slayer: kinetic file '$path' has zeff = $(minimum(data.zeff)) < 1 — unphysical")
+            error("run_slayer: $source has zeff = $(minimum(data.zeff)) < 1 — unphysical")
         cubic_interp(psi_xs, collect(Float64, data.zeff))
     elseif !_usable(data.n_i)
         nothing
     else
         z = [ne > 0 ? control.impurity_charge - (ni / ne) * (control.impurity_charge - 1.0) : control.impurity_charge for (ni, ne) in zip(data.n_i, data.n_e)]
         minimum(z) >= 1.0 - 1e-8 ||
-            error("run_slayer: kinetic file '$path' has n_i > n_e, giving Z_eff = $(minimum(z)) < 1 — unphysical")
+            error("run_slayer: $source has n_i > n_e, giving Z_eff = $(minimum(z)) < 1 — unphysical")
         cubic_interp(psi_xs, z)
     end
     zeff_source = _usable(data.zeff) ? :file : zeff === nothing ? :scalar : :closure
diff --git a/src/Utilities/KineticProfiles.jl b/src/Utilities/KineticProfiles.jl
index 5b900b859..b13497431 100644
--- a/src/Utilities/KineticProfiles.jl
+++ b/src/Utilities/KineticProfiles.jl
@@ -22,8 +22,8 @@ the normalized poloidal flux ψ ∈ [0, 1].
 | `omega_E` | E×B rotation Ω_E, per unit n | rad/s |
 
 Construct via the keyword constructor `KineticProfiles(; psi, n_e, T_e, T_i, omega_E)` with
-matched-length vectors. The SLAYER runner builds this object from a standardized
-kinetic-profile file via `Equilibrium.read_kinetic_file`. Diamagnetic frequencies are not
+matched-length vectors. The SLAYER runner builds this object from the equilibrium's
+kinetic profiles (`equil.kinetic`). Diamagnetic frequencies are not
 carried here: `build_slayer_inputs` derives ω_*e and ω_*i per surface from the density and
 temperature splines.
 
diff --git a/test/runtests_equil.jl b/test/runtests_equil.jl
index f8afac63c..fc7b01980 100644
--- a/test/runtests_equil.jl
+++ b/test/runtests_equil.jl
@@ -194,6 +194,20 @@ using TOML
         # A listed deck key must not still be a control field, or listing it would stop runs that set a live key.
         @test !any(k -> hasfield(GPEC.ForceFreeStates.ForceFreeStatesControl, Symbol(k)), GPEC._DEPRECATED_FFS_KEYS)
         @test !any(k -> hasfield(GPEC.KineticForces.KineticForcesControl, Symbol(k)), GPEC._DEPRECATED_KF_KEYS)
+        @test !any(k -> hasfield(GPEC.Runner.SLAYERControl, Symbol(k)), GPEC._DEPRECATED_SLAYER_KEYS)
+        # The kinetic-profile file moved to [Equilibrium]; the stop names the new home instead of "remove".
+        mktempdir() do deck_dir
+            template_dir = joinpath(@__DIR__, "test_data", "regression_solovev_ideal_example")
+            foreach(name -> cp(joinpath(template_dir, name), joinpath(deck_dir, name)), readdir(template_dir))
+            inputs = TOML.parsefile(joinpath(deck_dir, "gpec.toml"))
+            inputs["SLAYER"] = Dict{String,Any}("enabled" => false, "profile_file" => "kinetic.h5")
+            open(io -> TOML.print(io, inputs), joinpath(deck_dir, "gpec.toml"), "w")
+            @test_throws r"\[SLAYER\] sets `profile_file`, which is no longer a key; name the kinetic-profile file once, as kinetic_file" GPEC.main([deck_dir])
+            delete!(inputs, "SLAYER")
+            inputs["KineticForces"] = Dict{String,Any}("kinetic_file" => "kinetic.dat")
+            open(io -> TOML.print(io, inputs), joinpath(deck_dir, "gpec.toml"), "w")
+            @test_throws r"\[KineticForces\] sets `kinetic_file`, which is no longer a key; name the kinetic-profile file once" GPEC.main([deck_dir])
+        end
         pe_keywords = union((Base.kwarg_decl(m) for f in (GPEC.plasma_response, GPEC.perturbed_equilibrium) for m in methods(f))...)
         @test !any(k -> Symbol(k) in pe_keywords, GPEC._DEPRECATED_PE_KEYS)
     end
diff --git a/test/runtests_error_fields.jl b/test/runtests_error_fields.jl
index a4a7c41f2..ff21940c8 100644
--- a/test/runtests_error_fields.jl
+++ b/test/runtests_error_fields.jl
@@ -87,9 +87,10 @@ include("h5_metadata_check.jl")
             cp(joinpath(@__DIR__, "test_data", "ErrorFields", "tolerances_two_hoops.toml"), joinpath(dir, "tolerances.toml"))
             # Kinetic profiles consistent with this Solovev pressure, for the NTV torque couplings.
             cp(joinpath(@__DIR__, "..", "examples", "Solovev_kinetic_NTV_example", "kinetic.dat"), joinpath(dir, "kinetic.dat"))
+            inputs["Equilibrium"]["kinetic_file"] = "kinetic.dat"
             # Exercise the kinetic pipeline at smoke-test precision, retaining all bounce harmonics.
             inputs["KineticForces"] = Dict{String,Any}(
-                "kinetic_file" => "kinetic.dat", "verbose" => false, "write_outputs_to_HDF5" => true,
+                "verbose" => false, "write_outputs_to_HDF5" => true,
                 "rtol_psi" => 0.2, "rtol_pitch" => 1e-2, "atol_pitch" => 1e-5,
                 "rtol_x" => 1e-3, "atol_x" => 1e-6)
             # Small batches exercise orchestration; distribution accuracy is tested on synthetic inputs.
@@ -548,7 +549,10 @@ include("h5_metadata_check.jl")
             @test length(with_profiles.subplots) == 3                       # and the deposition profile its own panel
             @test with_profiles[3].series_list[1][:y] == c.torque_full_profile[:, 2]
             @test with_profiles[3].series_list[2][:y] == c.torque_residual_profile[:, 2]
-            @test_throws ErrorException GPEC.efc_couplings(pe.response, sets, rc, dom, cfg, GPEC.KineticForces.KineticForcesControl(), nothing)
+            attached = pe.response.ffs.equil.kinetic
+            pe.response.ffs.equil.kinetic = nothing
+            @test_throws "needs kinetic profiles" GPEC.efc_couplings(pe.response, sets, rc, dom, cfg, GPEC.KineticForces.KineticForcesControl())
+            pe.response.ffs.equil.kinetic = attached
 
             # Central differences: doubling the step moves the derivatives at O(h²).
             coarse = EF.compute_coil_sensitivities(sets[1:1], rc, ffs.equil, cfg,
diff --git a/test/runtests_frozen_equilibrium.jl b/test/runtests_frozen_equilibrium.jl
index 3829c3c4a..23a52beeb 100644
--- a/test/runtests_frozen_equilibrium.jl
+++ b/test/runtests_frozen_equilibrium.jl
@@ -115,6 +115,20 @@ using GeneralizedPerturbedEquilibrium
         @test analytic_snapshot.analytic isa Eq.SolovevConfig
         @test analytic_snapshot.analytic.elongation == 1.8
         @test analytic_snapshot.profiles.q_spline.y == profile_values[4]
+        # The stored deck's kinetic_file resolves against the h5's directory and is re-read; a missing file warns and leaves none.
+        kinetic_path = joinpath(dir, "kinetic.dat")
+        cp(joinpath(@__DIR__, "..", "examples", "Solovev_kinetic_NTV_example", "kinetic.dat"), kinetic_path)
+        analytic_inputs["Equilibrium"]["kinetic_file"] = "kinetic.dat"
+        h5open(path, "r+") do f
+            delete_object(f, "Input/gpec_toml_raw")
+            f["Input/gpec_toml_raw"] = sprint(io -> TOML.print(io, analytic_inputs))
+        end
+        with_kinetic = GPEC.equilibrium_from_h5(path).equil
+        @test with_kinetic.kinetic.n_e == Eq.read_kinetic_file(kinetic_path).n_e
+        @test with_kinetic.config.kinetic_file == kinetic_path
+        rm(kinetic_path)
+        without = @test_logs (:warn, r"kinetic_file is not readable") match_mode = :any GPEC.equilibrium_from_h5(path).equil
+        @test without.kinetic === nothing
         h5open(path, "r+") do f
             delete_object(f, "Equilibrium/Profiles/q")
         end
diff --git a/test/runtests_kinetic.jl b/test/runtests_kinetic.jl
index 981fe9e6a..a46bfc767 100644
--- a/test/runtests_kinetic.jl
+++ b/test/runtests_kinetic.jl
@@ -1,3 +1,4 @@
+using TOML
 
 @testset "KineticForces Unit Tests" begin
 
@@ -766,8 +767,8 @@
         end
 
         @testset "NTV splines identical: ASCII vs HDF5" begin
-            ka = E.load_kinetic_profiles(gpeckf)
-            kb = E.load_kinetic_profiles(h5)
+            ka = E.kinetic_profile_splines(E.read_kinetic_file(gpeckf))
+            kb = E.kinetic_profile_splines(E.read_kinetic_file(h5))
             for p in 0.0:0.1:1.0
                 @test ka.ne_spline(p) == kb.ne_spline(p)
                 @test ka.Ti_spline(p) == kb.Ti_spline(p)
@@ -791,13 +792,64 @@
             rm(tmp; force=true)
         end
 
+        @testset "schema table matches the data fields" begin
+            # Every profile dataset the reader and writer iterate is a KineticProfileData field, and vice versa.
+            table = Set(Symbol(first(p)) for p in E._KINETIC_H5_PROFILES)
+            @test table == setdiff(Set(fieldnames(E.KineticProfileData)), Set([:psi, :species_densities, :provenance]))
+        end
+
         @testset "missing required field errors" begin
-            # HDF5 with only psi + n_e: NTV load must error on the missing T_i
-            tmp = tempname() * ".h5"
+            # Only psi + n_e: the NTV view must error on the missing T_i
             d = E.KineticProfileData(; psi=[0.0, 0.5, 1.0], n_e=[1e19, 9e18, 8e18])
-            E.write_kinetic_h5(tmp, d)
-            @test_throws ErrorException E.load_kinetic_profiles(tmp)
-            rm(tmp; force=true)
+            @test_throws "missing required 'T_i'" E.kinetic_profile_splines(d)
+        end
+
+        @testset "attached to the equilibrium" begin
+            template = joinpath(@__DIR__, "test_data", "regression_solovev_ideal_example")
+            deck = TOML.parsefile(joinpath(template, "gpec.toml"))
+            sol = E.SolovevConfig(deck["SOL_INPUT"])
+            # A relative kinetic_file resolves against the deck directory, like eq_filename.
+            config = E.EquilibriumConfig(merge(deck["Equilibrium"], Dict{String,Any}("kinetic_file" => relpath(h5, template))), template)
+            @test config.kinetic_file == normpath(h5)
+            raw = E.read_kinetic_file(h5)
+            equil = E.PlasmaEquilibrium(config, sol)
+            @test equil.kinetic.psi == raw.psi
+            @test equil.kinetic.T_i == raw.T_i
+            @test equil.kinetic.chi_e == raw.chi_e
+            # A re-form passes the data through instead of reading the file again.
+            @test E.PlasmaEquilibrium(config, sol; kinetic=equil.kinetic).kinetic === equil.kinetic
+
+            bare = E.PlasmaEquilibrium(E.EquilibriumConfig(deck["Equilibrium"], template), sol)
+            @test bare.kinetic === nothing
+            E.attach_kinetic_profiles!(bare, gpeckf)
+            @test bare.kinetic.n_e == E.read_kinetic_file(gpeckf).n_e
+            @test E.attach_kinetic_profiles!(bare, raw).kinetic === raw
+        end
+
+        @testset "NTV views of the attached profiles" begin
+            GPEC = GeneralizedPerturbedEquilibrium
+            KF = GPEC.KineticForces
+            template = joinpath(@__DIR__, "test_data", "regression_solovev_ideal_example")
+            deck = TOML.parsefile(joinpath(template, "gpec.toml"))
+            equil = E.PlasmaEquilibrium(E.EquilibriumConfig(deck["Equilibrium"], template), E.SolovevConfig(deck["SOL_INPUT"]))
+            kf = KF.KineticForcesControl()
+            @test KF.ntv_profiles(kf, equil) === nothing
+            @test_throws "set kinetic_file in [Equilibrium]" GPEC.kinetic_forces_control(
+                Dict{String,Any}("KineticForces" => Dict{String,Any}()), GPEC.ForceFreeStates.ForceFreeStatesControl(), equil)
+
+            raw = E.read_kinetic_file(h5)
+            E.attach_kinetic_profiles!(equil, raw)
+            view = KF.ntv_profiles(kf, equil)
+            ref = E.kinetic_profile_splines(raw; zi=kf.ion_charge, zimp=kf.impurity_charge, mi=kf.ion_mass, mimp=kf.impurity_mass,
+                chi1=2π * equil.psio)
+            for ψ in 0.0:0.1:1.0
+                @test view.Ti_spline(ψ) == ref.Ti_spline(ψ)
+                @test view.nui_spline(ψ) == ref.nui_spline(ψ)
+                @test view.zeff_spline(ψ) == ref.zeff_spline(ψ)
+            end
+            @test KF.ntv_species(kf, equil) === nothing
+            with_e = KF.ntv_species(KF.KineticForcesControl(; electron=true), equil)
+            @test any(sp -> sp.electron, with_e)
         end
     end
 end
diff --git a/test/runtests_layer_overlap_fitzpatrick.jl b/test/runtests_layer_overlap_fitzpatrick.jl
index 2d91b59fc..b185ae5d9 100644
--- a/test/runtests_layer_overlap_fitzpatrick.jl
+++ b/test/runtests_layer_overlap_fitzpatrick.jl
@@ -228,23 +228,21 @@ end
     # The scan must see the run's own plasma: every [SLAYER] setting the layer width depends on
     # reaches it, and changing one changes the widths.
     trunc_dir = joinpath(dirname(@__DIR__), "examples", "DIIID-like_truncation_example")
-    slayer_section = TOML.parsefile(joinpath(trunc_dir, "gpec.toml"))["SLAYER"]
+    trunc_deck = TOML.parsefile(joinpath(trunc_dir, "gpec.toml"))
+    slayer_section = trunc_deck["SLAYER"]
+    kin = Equilibrium.read_kinetic_file(joinpath(trunc_dir, trunc_deck["Equilibrium"]["kinetic_file"]))
     slayer_inputs(over) = Dict{String,Any}("SLAYER" => merge(slayer_section, over))
-    overlap_inputs(over; dir=trunc_dir) = GeneralizedPerturbedEquilibrium._overlap_slayer_inputs(slayer_inputs(over), dir)
+    overlap_inputs(over; kinetic=kin) = GeneralizedPerturbedEquilibrium._overlap_slayer_inputs(slayer_inputs(over), kinetic)
     oi = overlap_inputs(Dict("ion_mass" => 3.0, "bt" => 1.7, "theta_sample" => 0.3, "lnLambda_form" => "sauter"))
     @test (oi.mu_i, oi.bt, oi.theta, oi.lnLambda_form) == (3.0, 1.7, 0.3, :sauter)
 
-    # A disabled section feeds the scan whatever its inner model, and a missing or empty file only skips the scan;
-    # with tearing enabled the same file is an error, raised here rather than after every other stage.
+    # A disabled section feeds the scan whatever its inner model, and missing kinetic profiles only skip the scan;
+    # with tearing enabled they are an error, raised here rather than after every other stage.
     @test overlap_inputs(Dict("inner_model" => "ggj_ray")) !== nothing
-    @test (@test_logs (:warn, r"Layer-overlap scan skipped") match_mode = :any overlap_inputs(Dict("profile_file" => "missing.h5"))) === nothing
-    @test (@test_logs overlap_inputs(Dict("profile_file" => ""))) === nothing
-    @test_throws ErrorException overlap_inputs(Dict("enabled" => true, "profile_file" => "missing.h5"))
-    @test_throws ErrorException overlap_inputs(Dict("enabled" => true, "profile_file" => ""))
+    @test (@test_logs overlap_inputs(Dict{String,Any}(); kinetic=nothing)) === nothing
+    @test_throws ErrorException overlap_inputs(Dict("enabled" => true); kinetic=nothing)
 
     # Z_eff(psi) comes from the kinetic file's own n_i/n_e with the [SLAYER] impurity charge.
-    kin_path = joinpath(trunc_dir, slayer_section["profile_file"])
-    kin = Equilibrium.read_kinetic_file(kin_path)
     j = length(kin.psi) ÷ 2
     @test oi.zeff(kin.psi[j]) ≈ 6.0 - kin.n_i[j] / kin.n_e[j] * 5.0
     @test overlap_inputs(Dict("impurity_charge" => 10.0)).zeff(kin.psi[j]) ≈ 10.0 - kin.n_i[j] / kin.n_e[j] * 9.0
@@ -256,12 +254,12 @@ end
         base = (; psi=kin.psi, n_e=kin.n_e, T_i=kin.T_i, T_e=kin.T_e, omega_E=kin.omega_E, chi_e=kin.chi_e, chi_phi=kin.chi_phi)
         Equilibrium.write_kinetic_h5(joinpath(tmp, "own_zeff.h5"), Equilibrium.KineticProfileData(; base..., n_i=kin.n_i, zeff=own))
         @test Equilibrium.read_kinetic_file(joinpath(tmp, "own_zeff.h5")).zeff == own
-        @test overlap_inputs(Dict("profile_file" => "own_zeff.h5", "impurity_charge" => 10.0, "zeff" => 2.5); dir=tmp).zeff(kin.psi[j]) ≈ 3.2
+        @test overlap_inputs(Dict("impurity_charge" => 10.0, "zeff" => 2.5); kinetic=Equilibrium.read_kinetic_file(joinpath(tmp, "own_zeff.h5"))).zeff(kin.psi[j]) ≈ 3.2
         Equilibrium.write_kinetic_h5(joinpath(tmp, "no_ni.h5"), Equilibrium.KineticProfileData(; base...))
-        @test overlap_inputs(Dict("profile_file" => "no_ni.h5", "zeff" => 2.5); dir=tmp).zeff == 2.5
-        # A P_perp_model whose χ⊥,i the file lacks skips the scan with a warning rather than aborting the run.
-        @test (@test_logs (:warn, r"Layer-overlap scan skipped") match_mode = :any overlap_inputs(
-            Dict("profile_file" => "no_ni.h5", "P_perp_model" => "D_perp"); dir=tmp)) === nothing
+        no_ni = Equilibrium.read_kinetic_file(joinpath(tmp, "no_ni.h5"))
+        @test overlap_inputs(Dict("zeff" => 2.5); kinetic=no_ni).zeff == 2.5
+        # A P_perp_model whose χ⊥,i the profiles lack skips the scan with a warning rather than aborting the run.
+        @test (@test_logs (:warn, r"Layer-overlap scan skipped") match_mode = :any overlap_inputs(Dict("P_perp_model" => "D_perp"); kinetic=no_ni)) === nothing
     end
 
     widths(over) = GeneralizedPerturbedEquilibrium._layer_overlap_scan(overlap_inputs(over), _fresh(), equil).width_dr
@@ -338,7 +336,7 @@ end
 
     # The scan's overlap point must therefore not depend on where the grid stops, and no
     # rational may fall between the in-grid search and the law.
-    oi = GeneralizedPerturbedEquilibrium._overlap_slayer_inputs(trunc_inputs, trunc_dir)
+    oi = GeneralizedPerturbedEquilibrium._overlap_slayer_inputs(trunc_inputs, long.kinetic)
     one_n() = (i = ForceFreeStatesInternal(); i.nlow = 1; i.nhigh = 1; i)
     scan_long = GeneralizedPerturbedEquilibrium._layer_overlap_scan(oi, one_n(), long)
     scan_short = GeneralizedPerturbedEquilibrium._layer_overlap_scan(oi, one_n(), short)
diff --git a/test/runtests_multiion.jl b/test/runtests_multiion.jl
index 95b0d58ee..d84f35856 100644
--- a/test/runtests_multiion.jl
+++ b/test/runtests_multiion.jl
@@ -49,9 +49,9 @@ using HDF5
             end
         end
 
-        # Single species [z=1,m=2,fraction=1] reproduces load_kinetic_profiles.
-        single = EQ.load_kinetic_profiles(f; zi=1, zimp=6, mi=2, mimp=12)
-        sp1 = EQ.resolve_ntv_species(f, [KF.IonSpecies(; z=1, m=2, fraction=1.0)]; electron=false, zimp=6, mimp=12)
+        # Single species [z=1,m=2,fraction=1] reproduces kinetic_profile_splines.
+        single = EQ.kinetic_profile_splines(EQ.read_kinetic_file(f); zi=1, zimp=6, mi=2, mimp=12)
+        sp1 = EQ.resolve_ntv_species(EQ.read_kinetic_file(f), [KF.IonSpecies(; z=1, m=2, fraction=1.0)]; electron=false, zimp=6, mimp=12)
         main = sp1[1]
         for ψ in (0.3, 0.6)
             @test main.profiles.ni_spline(ψ) ≈ single.ni_spline(ψ) rtol = 1e-8
@@ -60,7 +60,7 @@ using HDF5
         end
 
         # D-T-e set = {D, T, impurity, electron}; ν scales as 1/√m at fixed z, T.
-        sp = EQ.resolve_ntv_species(f, [KF.IonSpecies(; z=1, m=2, fraction=0.5), KF.IonSpecies(; z=1, m=3, fraction=0.5)];
+        sp = EQ.resolve_ntv_species(EQ.read_kinetic_file(f), [KF.IonSpecies(; z=1, m=2, fraction=0.5), KF.IonSpecies(; z=1, m=3, fraction=0.5)];
             electron=true, zimp=6, mimp=12)
         @test length(sp) == 4
         @test count(s -> !s.electron && s.z == 1, sp) == 2   # D, T
@@ -71,11 +71,11 @@ using HDF5
         @test sp[iD].profiles.nui_spline(0.5) / sp[iT].profiles.nui_spline(0.5) ≈ sqrt(3 / 2) rtol = 1e-6
 
         # Exactly one of {fraction, density} required.
-        @test_throws ErrorException EQ.resolve_ntv_species(f, [KF.IonSpecies(; z=1, m=2)]; zimp=6, mimp=12)
+        @test_throws ErrorException EQ.resolve_ntv_species(EQ.read_kinetic_file(f), [KF.IonSpecies(; z=1, m=2)]; zimp=6, mimp=12)
 
         # Unified electron semantics: a single-ion + electron run resolves as
         # {main ion, impurity, electron} — electrons in ADDITION to the ions.
-        spe = EQ.resolve_ntv_species(f, [KF.IonSpecies(; z=1, m=2, fraction=1.0)]; electron=true, zimp=6, mimp=12)
+        spe = EQ.resolve_ntv_species(EQ.read_kinetic_file(f), [KF.IonSpecies(; z=1, m=2, fraction=1.0)]; electron=true, zimp=6, mimp=12)
         @test length(spe) == 3
         @test !spe[1].electron && spe[1].z == 1                  # main ion unchanged
         @test spe[1].profiles.ni_spline(0.5) ≈ main.profiles.ni_spline(0.5)
@@ -101,14 +101,14 @@ using HDF5
             file["n_D"] = 0.6 .* Ni
             file["n_T"] = 0.4 .* Ni     # unequal explicit shapes
         end
-        sp = EQ.resolve_ntv_species(h5, [KF.IonSpecies(; z=1, m=2, density="n_D"), KF.IonSpecies(; z=1, m=3, density="n_T")];
+        sp = EQ.resolve_ntv_species(EQ.read_kinetic_file(h5), [KF.IonSpecies(; z=1, m=2, density="n_D"), KF.IonSpecies(; z=1, m=3, density="n_T")];
             electron=false, zimp=6, mimp=12)
         iD = findfirst(s -> s.m == 2, sp)
         iT = findfirst(s -> s.m == 3, sp)
         @test sp[iD].profiles.ni_spline(0.5) ≈ 0.6 * 1.0e20 * (1 - 0.25) rtol = 1e-6
         @test sp[iT].profiles.ni_spline(0.5) ≈ 0.4 * 1.0e20 * (1 - 0.25) rtol = 1e-6
         # A missing named profile errors clearly.
-        @test_throws ErrorException EQ.resolve_ntv_species(h5, [KF.IonSpecies(; z=1, m=2, density="n_He")]; zimp=6, mimp=12)
+        @test_throws ErrorException EQ.resolve_ntv_species(EQ.read_kinetic_file(h5), [KF.IonSpecies(; z=1, m=2, density="n_He")]; zimp=6, mimp=12)
         rm(h5)
     end
 
@@ -127,7 +127,7 @@ using HDF5
         @test back.species_densities["n_D"] ≈ data.species_densities["n_D"]
         @test back.species_densities["n_T"] ≈ data.species_densities["n_T"]
         # The written file feeds the explicit-density multi-ion input directly.
-        sp = EQ.resolve_ntv_species(h5, [KF.IonSpecies(; z=1, m=2, density="n_D"), KF.IonSpecies(; z=1, m=3, density="n_T")];
+        sp = EQ.resolve_ntv_species(EQ.read_kinetic_file(h5), [KF.IonSpecies(; z=1, m=2, density="n_D"), KF.IonSpecies(; z=1, m=3, density="n_T")];
             zimp=6, mimp=12)
         @test count(s -> !s.electron && s.z == 1, sp) == 2
         # Non-`n_*` species names are rejected at write time (they would not round-trip).
@@ -157,20 +157,20 @@ using HDF5
         end
         DT(f1, f2) = [KF.IonSpecies(; z=1, m=2, fraction=f1), KF.IonSpecies(; z=1, m=3, fraction=f2)]
         # An all-fraction list must account for the full n_i (the impurity is set by n_i/n_e, not a shortfall).
-        @test_throws ErrorException EQ.resolve_ntv_species(h5, DT(0.5, 0.4); zimp=6, mimp=12)
+        @test_throws ErrorException EQ.resolve_ntv_species(EQ.read_kinetic_file(h5), DT(0.5, 0.4); zimp=6, mimp=12)
         # Fractions are shares of n_i, so summing above 1 is fatal in any list, mixed included.
-        @test_throws ErrorException EQ.resolve_ntv_species(h5, DT(0.7, 0.5); zimp=6, mimp=12)
-        @test_throws ErrorException EQ.resolve_ntv_species(h5,
+        @test_throws ErrorException EQ.resolve_ntv_species(EQ.read_kinetic_file(h5), DT(0.7, 0.5); zimp=6, mimp=12)
+        @test_throws ErrorException EQ.resolve_ntv_species(EQ.read_kinetic_file(h5),
             [KF.IonSpecies(; z=1, m=2, fraction=1.2), KF.IonSpecies(; z=1, m=3, density="n_T")]; zimp=6, mimp=12)
         # A partial fraction alongside an explicit density is legal — the n_i/n_e deficit is the impurity.
-        sp = EQ.resolve_ntv_species(h5,
+        sp = EQ.resolve_ntv_species(EQ.read_kinetic_file(h5),
             [KF.IonSpecies(; z=1, m=2, fraction=0.5), KF.IonSpecies(; z=1, m=3, density="n_T")]; zimp=6, mimp=12)
         @test count(s -> !s.electron && s.z == 1, sp) == 2
         # Main-ion z must stay below the impurity charge (zpitch closure pole at Zeff = zimp).
-        @test_throws ErrorException EQ.resolve_ntv_species(h5, [KF.IonSpecies(; z=6, m=12, fraction=1.0)]; zimp=6, mimp=12)
+        @test_throws ErrorException EQ.resolve_ntv_species(EQ.read_kinetic_file(h5), [KF.IonSpecies(; z=6, m=12, fraction=1.0)]; zimp=6, mimp=12)
         # A dilute main-ion mix drives Zeff toward zimp and warns about the zpitch pole.
         @test_logs (:warn, r"zpitch closure near its pole") match_mode = :any EQ.resolve_ntv_species(
-            h5, [KF.IonSpecies(; z=1, m=2, density="n_dilute")]; zimp=6, mimp=12)
+            EQ.read_kinetic_file(h5), [KF.IonSpecies(; z=1, m=2, density="n_dilute")]; zimp=6, mimp=12)
         # Vacuum point: neutral composition, no impurity, finite zpitch.
         @test EQ.multi_ion_composition([1], [1.0e19], 0.0; zimp=6, mimp=12) == (1.0, 1.0, 1.0e19, 0.0)
         rm(h5)
diff --git a/test/runtests_slayer_profile_transport.jl b/test/runtests_slayer_profile_transport.jl
index 0ca137890..21696c7a2 100644
--- a/test/runtests_slayer_profile_transport.jl
+++ b/test/runtests_slayer_profile_transport.jl
@@ -1,32 +1,29 @@
 using Test
-using GeneralizedPerturbedEquilibrium.Equilibrium: KineticProfileData, write_kinetic_h5
+using GeneralizedPerturbedEquilibrium.Equilibrium: KineticProfileData
 using GeneralizedPerturbedEquilibrium.Runner
 
 @testset "SLAYER loads only transport profiles used by its model" begin
-    mktempdir() do dir
-        psi = collect(range(0.0, 1.0; length=5))
-        function load_transport(model; chi_e=ones(5), chi_i=ones(5), chi_phi=ones(5), dc_type=:none)
-            data = KineticProfileData(; psi, n_e=fill(5e19, 5), T_e=1000 .- 500 .* psi,
-                T_i=800 .- 300 .* psi, chi_e, chi_i, chi_phi)
-            write_kinetic_h5(joinpath(dir, "profiles.h5"), data)
-            ctrl = SLAYERControl(; enabled=true, profile_file="profiles.h5", P_perp_model=model, tau_E=0.1, D_c_type=dc_type)
-            return Runner._load_profiles(ctrl, dir)
-        end
+    psi = collect(range(0.0, 1.0; length=5))
+    function load_transport(model; chi_e=ones(5), chi_i=ones(5), chi_phi=ones(5), dc_type=:none)
+        data = KineticProfileData(; psi, n_e=fill(5e19, 5), T_e=1000 .- 500 .* psi,
+            T_i=800 .- 300 .* psi, chi_e, chi_i, chi_phi)
+        ctrl = SLAYERControl(; enabled=true, P_perp_model=model, tau_E=0.1, D_c_type=dc_type)
+        return Runner._slayer_profiles(data, ctrl)
+    end
 
-        @test load_transport(:tau_E; chi_phi=(-ones(5))).chi_tor === nothing
-        for model in (:chi_perp_e, :P_phi, :c_beta, :tau_E)
-            @test load_transport(model; chi_i=(-ones(5))).chi_perp_i === nothing
-        end
-        for model in (:chi_perp_i, :P_phi, :c_beta, :tau_E)
-            @test load_transport(model; chi_e=(-ones(5))).chi_perp_e === nothing
-            @test_throws "chi_e" load_transport(model; chi_e=(-ones(5)), dc_type=:lar)
-        end
-        @test_throws "chi_i" load_transport(:D_perp; chi_i=(-ones(5)))
-        @test_throws "chi_i" load_transport(:chi_perp_i; chi_i=(-ones(5)))
-        @test_throws "chi_e" load_transport(:D_perp; chi_e=(-ones(5)))
-        @test_throws "chi_e" load_transport(:chi_perp_e; chi_e=(-ones(5)))
-        @test_throws "chi_phi" load_transport(:chi_perp_e; chi_phi=(-ones(5)))
-        @test load_transport(:D_perp; chi_i=zeros(5)).chi_perp_i === nothing
-        @test load_transport(:D_perp).chi_perp_i(0.5) ≈ 1.0
+    @test load_transport(:tau_E; chi_phi=(-ones(5))).chi_tor === nothing
+    for model in (:chi_perp_e, :P_phi, :c_beta, :tau_E)
+        @test load_transport(model; chi_i=(-ones(5))).chi_perp_i === nothing
+    end
+    for model in (:chi_perp_i, :P_phi, :c_beta, :tau_E)
+        @test load_transport(model; chi_e=(-ones(5))).chi_perp_e === nothing
+        @test_throws "chi_e" load_transport(model; chi_e=(-ones(5)), dc_type=:lar)
     end
+    @test_throws "chi_i" load_transport(:D_perp; chi_i=(-ones(5)))
+    @test_throws "chi_i" load_transport(:chi_perp_i; chi_i=(-ones(5)))
+    @test_throws "chi_e" load_transport(:D_perp; chi_e=(-ones(5)))
+    @test_throws "chi_e" load_transport(:chi_perp_e; chi_e=(-ones(5)))
+    @test_throws "chi_phi" load_transport(:chi_perp_e; chi_phi=(-ones(5)))
+    @test load_transport(:D_perp; chi_i=zeros(5)).chi_perp_i === nothing
+    @test load_transport(:D_perp).chi_perp_i(0.5) ≈ 1.0
 end
diff --git a/test/runtests_slayer_runner.jl b/test/runtests_slayer_runner.jl
index 3c07ff7de..ceb045863 100644
--- a/test/runtests_slayer_runner.jl
+++ b/test/runtests_slayer_runner.jl
@@ -76,8 +76,7 @@
                 "max_cells" => 50_000),
             "growth_rate_filter" => Dict{String,Any}(
                 "pole_threshold" => 1e5,
-                "filter_above_poles" => false),
-            "profile_file" => "profiles.h5"
+                "filter_above_poles" => false)
         )
         c = slayer_control_from_toml(section)
         @test c.enabled
@@ -97,7 +96,6 @@
         @test c.amr_max_cells == 50_000
         @test c.pole_threshold == 1e5
         @test c.filter_above_poles == false
-        @test c.profile_file == "profiles.h5"
 
         # Unknown keys should raise
         bad = merge(section, Dict{String,Any}("mistyped_key" => 42))
@@ -123,15 +121,38 @@
     @testset "TearingProblem solve forwards surfaces and Δ'" begin
         # A result with no singular surfaces short-circuits before any equilibrium access,
         # so a stand-in result is enough to pin the forwarding of surfaces / delta_prime.
-        c = SLAYERControl(; enabled=true, profile_file="unused.h5")
+        c = SLAYERControl(; enabled=true)
         no_surfaces = (equil=nothing, surfaces=GeneralizedPerturbedEquilibrium.ForceFreeStates.SingType[],
             delta_prime=nothing)
         r = solve(TearingProblem(no_surfaces, c), SLAYERModel())
         @test isempty(r.params)
 
         # A disabled control never looks at the result at all.
-        r_off = solve(TearingProblem(no_surfaces, SLAYERControl(; enabled=false, profile_file="unused.h5")), SLAYERModel())
+        r_off = solve(TearingProblem(no_surfaces, SLAYERControl(; enabled=false)), SLAYERModel())
         @test r_off.enabled == false
+
+        # With surfaces to analyse, the profiles must come from the equilibrium.
+        bare = (equil=(; kinetic=nothing), surfaces=[nothing], delta_prime=nothing)
+        @test_throws "carries no kinetic profiles" solve(TearingProblem(bare, c), SLAYERModel())
+    end
+
+    @testset "SLAYER view of the kinetic profiles" begin
+        kin_path = joinpath(dirname(@__DIR__), "examples", "DIIID-like_ideal_example", "TkMkr_D3Dlike_Hmode_kinetic.h5")
+        data = GeneralizedPerturbedEquilibrium.Equilibrium.read_kinetic_file(kin_path)
+        view = Runner._slayer_profiles(data, SLAYERControl())
+        j = length(data.psi) ÷ 2
+        @test view.profiles.n_e(data.psi[j]) ≈ data.n_e[j]
+        @test view.profiles.omega_E(data.psi[j]) ≈ data.omega_E[j]
+        @test view.chi_perp_e(data.psi[j]) ≈ data.chi_e[j]
+        @test view.chi_perp_i(data.psi[j]) ≈ data.chi_i[j]
+        @test view.has_omega_E
+        # Without an E×B rotation profile the view records its absence and carries zero rotation.
+        no_rotation = GeneralizedPerturbedEquilibrium.Equilibrium.KineticProfileData(; psi=data.psi, n_e=data.n_e, T_e=data.T_e, T_i=data.T_i)
+        bare_view = Runner._slayer_profiles(no_rotation, SLAYERControl())
+        @test !bare_view.has_omega_E
+        @test bare_view.profiles.omega_E(data.psi[j]) == 0.0
+        @test_throws "missing required dataset 'T_i'" Runner._slayer_profiles(
+            GeneralizedPerturbedEquilibrium.Equilibrium.KineticProfileData(; psi=data.psi, n_e=data.n_e, T_e=data.T_e), SLAYERControl())
     end
 
     @testset "GGJ tearing is rejected until its growth-rate extraction exists" begin
@@ -143,7 +164,7 @@
         end
         no_surfaces = (equil=nothing, surfaces=GeneralizedPerturbedEquilibrium.ForceFreeStates.SingType[],
             delta_prime=nothing, dir_path=".")
-        @test_throws ErrorException solve(TearingProblem(no_surfaces, SLAYERControl(; enabled=true, profile_file="unused.h5")), GGJModel())
+        @test_throws ErrorException solve(TearingProblem(no_surfaces, SLAYERControl(; enabled=true)), GGJModel())
     end
 
     # Δ′ is unified across formalisms, so a Galerkin run feeds SLAYER exactly as a Riccati one
diff --git a/test/runtests_solve_api.jl b/test/runtests_solve_api.jl
index f7caa70ea..0124fe2e2 100644
--- a/test/runtests_solve_api.jl
+++ b/test/runtests_solve_api.jl
@@ -106,6 +106,51 @@ using HDF5
         end
     end
 
+    @testset "calculated kinetic solve matches the TOML-driven run" begin
+        # The equilibrium carries the deck's kinetic profiles, so the API solves the same kinetic problem.
+        kin_template = joinpath(@__DIR__, "test_data", "regression_solovev_kinetic_calculated")
+        kdeck = TOML.parsefile(joinpath(kin_template, "gpec.toml"))
+        kequil = GPEC.Equilibrium.PlasmaEquilibrium(
+            GPEC.Equilibrium.EquilibriumConfig(kdeck["Equilibrium"], kin_template),
+            GPEC.Equilibrium.SolovevConfig(kdeck["SOL_INPUT"])
+        )
+        @test kequil.kinetic !== nothing
+        kwall = GPEC.Vacuum.WallShapeSettings(; (Symbol(k) => v for (k, v) in kdeck["Wall"])...)
+        kkwargs = Dict(Symbol(k) => v for (k, v) in kdeck["ForceFreeStates"] if !(k in ("integrator", "nn_low", "nn_high")))
+        mktempdir() do dir
+            for name in readdir(kin_template)
+                cp(joinpath(kin_template, name), joinpath(dir, name))
+            end
+            reference = GPEC.main([dir]).ffs
+            api = solve(kequil, Forward(); nn=1, wall=kwall, dir_path=dir, kkwargs...)
+            @test api.mats.kinetic !== nothing
+            @test api.free_boundary.et[1] ≈ reference.free_boundary.et[1] rtol = 1e-10
+            @test api.diagnostics.nzero == reference.diagnostics.nzero
+            @test api.kinetic.kmsing == reference.kinetic.kmsing
+        end
+    end
+
+    @testset "attached kinetic profiles refine the two-pass grid" begin
+        # main and solve both re-form through maybe_reform_equilibrium, which reads equil.kinetic, so
+        # the attached profiles shape the grid the same way on either path. The DIII-D-like profiles
+        # add knots over the ideal ones.
+        sl_dir = joinpath(@__DIR__, "..", "examples", "DIIID-like_SLAYER_example")
+        sl_deck = TOML.parsefile(joinpath(sl_dir, "gpec.toml"))
+        sl_config = GPEC.Equilibrium.EquilibriumConfig(sl_deck["Equilibrium"], sl_dir)
+        pass1 = GPEC.Equilibrium.PlasmaEquilibrium(sl_config)
+        @test pass1.kinetic !== nothing
+        sl_ctrl = FFS.ForceFreeStatesControl(; (Symbol(k) => v for (k, v) in sl_deck["ForceFreeStates"])...)
+        sl_intr = FFS.ForceFreeStatesInternal()
+        GPEC.resolve_mode_space!(sl_intr, sl_ctrl)
+        refined = GPEC.maybe_reform_equilibrium(pass1, sl_intr, sl_ctrl)
+        @test refined.kinetic === pass1.kinetic        # carried through the re-form, not re-read
+        kinetic = pass1.kinetic
+        pass1.kinetic = nothing
+        ideal = GPEC.maybe_reform_equilibrium(pass1, sl_intr, sl_ctrl)
+        pass1.kinetic = kinetic
+        @test refined.rzphi_xs != ideal.rzphi_xs
+    end
+
     @testset "galerkin solve returns a gal result with a unified delta_prime" begin
         mktempdir() do dir
             api = solve(equil, Galerkin(; nx=32); nn=1, wall=wall, dir_path=dir, ffs_kwargs...)
@@ -239,7 +284,8 @@ using HDF5
     end
 
     @testset "rejected keyword combinations" begin
-        @test_throws "kinetic_source = \"calculated\"" solve(equil, Forward(); nn=1, dir_path=".", ffs_kwargs..., kinetic_factor=0.5, kinetic_source="calculated")
+        @test_throws "needs kinetic profiles on the equilibrium" solve(equil, Forward(); nn=1, dir_path=".", ffs_kwargs..., kinetic_factor=0.5,
+            kinetic_source="calculated")
         @test_throws "does not support kinetic runs" solve(equil, Riccati(); nn=1, dir_path=".", ffs_kwargs..., kinetic_factor=1e-9)
         @test_throws "does not support kinetic runs" solve(equil, Galerkin(); nn=1, dir_path=".", ffs_kwargs..., kinetic_factor=1e-9)
         @test_throws ErrorException solve(equil, Forward(); nn=1, dir_path=".", ffs_kwargs..., integrator="riccati")
diff --git a/test/test_data/regression_solovev_kinetic_calculated/gpec.toml b/test/test_data/regression_solovev_kinetic_calculated/gpec.toml
index 30fd5f854..6d482c706 100644
--- a/test/test_data/regression_solovev_kinetic_calculated/gpec.toml
+++ b/test/test_data/regression_solovev_kinetic_calculated/gpec.toml
@@ -4,6 +4,7 @@
 
 [Equilibrium]
 eq_type = "sol"                         # Type of the input 2D equilibrium file
+kinetic_file = "kinetic.dat"            # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 jac_type = "pest"                       # Coordinate system (hamada, pest, boozer, equal_arc, park, custom)
 grid_type = "core_edge_packed"          # Radial grid packing type
 psilow = 1e-4                           # Lower limit of normalized poloidal flux
diff --git a/test/test_data/regression_solovev_kinetic_nuzero/gpec.toml b/test/test_data/regression_solovev_kinetic_nuzero/gpec.toml
index 0404c6025..c5743bbb9 100644
--- a/test/test_data/regression_solovev_kinetic_nuzero/gpec.toml
+++ b/test/test_data/regression_solovev_kinetic_nuzero/gpec.toml
@@ -5,6 +5,7 @@
 # tracked separately by the solovev_kinetic_nuzero regression case.
 [Equilibrium]
 eq_type = "sol"                         # Type of the input 2D equilibrium file
+kinetic_file = "kinetic.dat"            # Kinetic-profile file (HDF5 GPEC kinetic schema or 6-column ASCII), read once into the equilibrium
 jac_type = "pest"                       # Coordinate system (hamada, pest, boozer, equal_arc, park, custom)
 grid_type = "core_edge_packed"          # Radial grid packing type
 psilow = 1e-4                           # Lower limit of normalized flux coordinate