1use serde::{Deserialize, Serialize};
2
3use crate::config::SimulationConfig;
4use crate::model::ModelParameters;
5
6const ADMISSIBLE_MARGIN_FRACTION: f64 = 0.05;
7const REDUCED_RESPONSE_TARGET_FRACTION: f64 = 0.12;
8
9#[derive(Debug, Clone, Serialize, Deserialize, Default)]
10pub struct AdmissibleStatus {
11 pub below_energy_min: bool,
12 pub above_energy_max: bool,
13 pub below_temperature_min: bool,
14 pub above_temperature_max: bool,
15 pub above_abs_y_max: bool,
16 pub above_abs_v_max: bool,
17 pub near_boundary: bool,
18 pub outside_region: bool,
19 pub margin_fraction: f64,
20}
21
22#[derive(Debug, Clone, Serialize, Deserialize)]
23pub struct StabilitySummary {
24 pub target_proxy: String,
25 pub v_initial: f64,
26 pub v_final: f64,
27 pub v_max: f64,
28 pub d_v_positive_fraction: f64,
29 pub local_stability_margin: f64,
30 pub first_positive_d_v_time_s: Option<f64>,
31}
32
33#[derive(Debug, Clone, Serialize, Deserialize)]
34pub struct FigureMetadata {
35 pub scenario: String,
36 pub y_label: String,
37 pub y_interpretation_note: String,
38 pub mechanical_work_interpretation_note: String,
39 pub recharge_interpretation_note: String,
40 pub burst_windows: Vec<BurstWindow>,
41 pub thresholds: ThresholdMetadata,
42}
43
44#[derive(Debug, Clone, Serialize, Deserialize)]
45pub struct BurstWindow {
46 pub label: String,
47 pub start_s: f64,
48 pub end_s: f64,
49}
50
51#[derive(Debug, Clone, Serialize, Deserialize)]
52pub struct ThresholdMetadata {
53 pub energy_min_j: f64,
54 pub energy_threshold_j: f64,
55 pub energy_max_j: f64,
56 pub temperature_min_k: f64,
57 pub temperature_max_k: f64,
58 pub y_abs_max: f64,
59 pub v_abs_max: f64,
60}
61
62pub fn admissible_status(
63 params: &ModelParameters,
64 ep_j: f64,
65 temperature_k: f64,
66 y_m: f64,
67 v_mps: f64,
68) -> AdmissibleStatus {
69 let energy_min = params.pulse_energy_min_j;
70 let energy_max = params.pulse_energy_max_j;
71 let temperature_min = params.ambient_temperature_k;
72 let temperature_max = params.thermal_limit_k;
73 let y_abs_max = params.max_displacement_m;
74 let v_abs_max = params.max_velocity_m_per_s;
75
76 let below_energy_min = ep_j < energy_min;
77 let above_energy_max = ep_j > energy_max;
78 let below_temperature_min = temperature_k < temperature_min;
79 let above_temperature_max = temperature_k > temperature_max;
80 let above_abs_y_max = y_m.abs() > y_abs_max;
81 let above_abs_v_max = v_mps.abs() > v_abs_max;
82
83 let energy_margin = normalized_margin(ep_j, energy_min, energy_max);
84 let temperature_margin = normalized_margin(temperature_k, temperature_min, temperature_max);
85 let y_margin = normalized_margin(y_m.abs(), 0.0, y_abs_max);
86 let v_margin = normalized_margin(v_mps.abs(), 0.0, v_abs_max);
87 let margin_fraction = energy_margin
88 .min(temperature_margin)
89 .min(y_margin)
90 .min(v_margin);
91
92 AdmissibleStatus {
93 below_energy_min,
94 above_energy_max,
95 below_temperature_min,
96 above_temperature_max,
97 above_abs_y_max,
98 above_abs_v_max,
99 near_boundary: margin_fraction <= ADMISSIBLE_MARGIN_FRACTION,
100 outside_region: below_energy_min
101 || above_energy_max
102 || below_temperature_min
103 || above_temperature_max
104 || above_abs_y_max
105 || above_abs_v_max,
106 margin_fraction,
107 }
108}
109
110pub fn reduced_response_target(command_fraction: f64, params: &ModelParameters) -> f64 {
111 command_fraction.max(0.0) * params.max_displacement_m * REDUCED_RESPONSE_TARGET_FRACTION
112}
113
114pub fn lyapunov_value(
115 mechanical_mass_kg: f64,
116 k_eff_n_per_m: f64,
117 error_m: f64,
118 error_rate_mps: f64,
119) -> f64 {
120 0.5 * mechanical_mass_kg.max(1.0) * error_rate_mps * error_rate_mps
121 + 0.5 * k_eff_n_per_m.max(0.0) * error_m * error_m
122}
123
124pub fn normalized_authority_utilization(
125 params: &ModelParameters,
126 gain_n: f64,
127 command_fraction: f64,
128 delivered_ratio: f64,
129) -> f64 {
130 let gain_fraction = gain_n / params.reference_actuator_force_n.max(1.0);
131 (gain_fraction * command_fraction.max(0.0) * delivered_ratio).clamp(0.0, 1.5)
132}
133
134pub fn stability_summary(
135 v_series: &[f64],
136 dv_dt_series: &[f64],
137 time_series: &[f64],
138) -> StabilitySummary {
139 let v_initial = v_series.first().copied().unwrap_or(0.0);
140 let v_final = v_series.last().copied().unwrap_or(0.0);
141 let v_max = v_series.iter().copied().fold(0.0, f64::max);
142 let positive_count = dv_dt_series.iter().filter(|value| **value > 0.0).count();
143 let d_v_positive_fraction = if dv_dt_series.is_empty() {
144 0.0
145 } else {
146 positive_count as f64 / dv_dt_series.len() as f64
147 };
148 let mean_dv_dt = if dv_dt_series.is_empty() {
149 0.0
150 } else {
151 dv_dt_series.iter().sum::<f64>() / dv_dt_series.len() as f64
152 };
153 let local_stability_margin = -mean_dv_dt / v_max.max(1.0);
154 let first_positive_d_v_time_s = dv_dt_series
155 .iter()
156 .zip(time_series.iter())
157 .find(|(value, _)| **value > 0.0)
158 .map(|(_, time_s)| *time_s);
159
160 StabilitySummary {
161 target_proxy: "command-scaled reduced maneuver response proxy".to_string(),
162 v_initial,
163 v_final,
164 v_max,
165 d_v_positive_fraction,
166 local_stability_margin,
167 first_positive_d_v_time_s,
168 }
169}
170
171pub fn build_figure_metadata(config: &SimulationConfig) -> FigureMetadata {
172 let burst_windows = config
173 .scenario
174 .segments
175 .iter()
176 .filter(|segment| segment.label.contains("burst") || segment.demand_fraction >= 0.75)
177 .map(|segment| BurstWindow {
178 label: segment.label.clone(),
179 start_s: segment.start_s,
180 end_s: segment.end_s,
181 })
182 .collect();
183
184 FigureMetadata {
185 scenario: config.scenario.name.clone(),
186 y_label: "Reduced Maneuver Response y".to_string(),
187 y_interpretation_note: "The state y is a reduced-order maneuver response proxy coupled to energy and thermal state. It is not intended to represent literal full-body mech displacement.".to_string(),
188 mechanical_work_interpretation_note: "Delivered mechanical work in this crate should be read as reduced-order architecture evidence for authority delivery, not as a full rigid-body maneuver validation.".to_string(),
189 recharge_interpretation_note: "Recharge time depends on the actual depleted energy in a run. A short burst-refill tail and a longer 3 GJ recharge case are consistent because they correspond to different refill depths.".to_string(),
190 burst_windows,
191 thresholds: ThresholdMetadata {
192 energy_min_j: config.model.pulse_energy_min_j,
193 energy_threshold_j: config.model.low_energy_threshold_j,
194 energy_max_j: config.model.pulse_energy_max_j,
195 temperature_min_k: config.model.ambient_temperature_k,
196 temperature_max_k: config.model.thermal_limit_k,
197 y_abs_max: config.model.max_displacement_m,
198 v_abs_max: config.model.max_velocity_m_per_s,
199 },
200 }
201}
202
203fn normalized_margin(value: f64, min: f64, max: f64) -> f64 {
204 let span = (max - min).abs().max(1.0e-9);
205 let lower = (value - min) / span;
206 let upper = (max - value) / span;
207 lower.min(upper)
208}