Skip to main content

runmat_runtime/builtins/control/
pole.rs

1//! Pole extraction for transfer-function and state-space control models.
2
3use runmat_builtins::{
4    BuiltinCompletionPolicy, BuiltinDescriptor, BuiltinErrorDescriptor, BuiltinOutputMode,
5    BuiltinParamArity, BuiltinParamDescriptor, BuiltinParamType, BuiltinSignatureDescriptor,
6};
7use runmat_builtins::{
8    BuiltinIntegerBackendRule, BuiltinIntegerCapabilityDescriptor, BuiltinIntegerComputationDomain,
9    BuiltinIntegerInputAvailability, BuiltinIntegerInputCapability, BuiltinIntegerOutputClassRule,
10    BuiltinIntegerOverflowRule, BuiltinIntegerOverloadKind, BuiltinIntegerScalarDoubleRule,
11};
12use runmat_macros::runtime_builtin;
13use runmat_value::Value;
14
15use crate::builtins::common::spec::{
16    BroadcastSemantics, BuiltinFusionSpec, BuiltinGpuSpec, ConstantStrategy, GpuOpKind,
17    ReductionNaN, ResidencyPolicy, ShapeRequirements,
18};
19use crate::builtins::control::tf_model::{
20    control_error, output_complex_column, ss_poles_from_object, TfModel, SS_CLASS, TF_CLASS,
21};
22use crate::builtins::control::type_resolvers::pole_type;
23use crate::{dispatcher, BuiltinResult};
24
25const POLE_OUTPUT: [BuiltinParamDescriptor; 1] = [BuiltinParamDescriptor {
26    name: "p",
27    ty: BuiltinParamType::Any,
28    arity: BuiltinParamArity::Required,
29    default: None,
30    description: "Poles of the SISO tf or ss model as a column vector.",
31}];
32const POLE_INPUTS: [BuiltinParamDescriptor; 1] = [BuiltinParamDescriptor {
33    name: "sys",
34    ty: BuiltinParamType::Any,
35    arity: BuiltinParamArity::Required,
36    default: None,
37    description: "SISO tf model or ss state-space model.",
38}];
39const POLE_SIGNATURES: [BuiltinSignatureDescriptor; 1] = [BuiltinSignatureDescriptor {
40    label: "p = pole(sys)",
41    inputs: &POLE_INPUTS,
42    outputs: &POLE_OUTPUT,
43}];
44const POLE_ERRORS: [BuiltinErrorDescriptor; 4] = [
45    BuiltinErrorDescriptor {
46        code: "RM.POLE.INVALID_MODEL",
47        identifier: Some("RunMat:pole:InvalidModel"),
48        when: "Input system is not a valid SISO tf or ss object.",
49        message: "pole: invalid model",
50    },
51    BuiltinErrorDescriptor {
52        code: "RM.POLE.UNSUPPORTED_MODEL",
53        identifier: Some("RunMat:pole:UnsupportedModel"),
54        when: "Model form is unsupported.",
55        message: "pole: unsupported model",
56    },
57    BuiltinErrorDescriptor {
58        code: "RM.POLE.INVALID_ARGUMENT",
59        identifier: Some("RunMat:pole:InvalidArgument"),
60        when: "Model metadata or arguments are malformed.",
61        message: "pole: invalid argument",
62    },
63    BuiltinErrorDescriptor {
64        code: "RM.POLE.INTERNAL",
65        identifier: Some("RunMat:pole:Internal"),
66        when: "Root calculation or output construction failed.",
67        message: "pole: internal error",
68    },
69];
70pub const POLE_DESCRIPTOR: BuiltinDescriptor = BuiltinDescriptor {
71    signatures: &POLE_SIGNATURES,
72    output_mode: BuiltinOutputMode::Fixed,
73    completion_policy: BuiltinCompletionPolicy::Public,
74    errors: &POLE_ERRORS,
75};
76const POLE_INTEGER_MODEL_INDEX_INPUT: [BuiltinIntegerInputCapability; 1] =
77    [BuiltinIntegerInputCapability {
78        name: "J1,...,JN model-array indices",
79        classes: &crate::builtins::common::integer_capability::ALL_INTEGER_CLASSES,
80        availability: BuiltinIntegerInputAvailability::Documented,
81        scalar_double: BuiltinIntegerScalarDoubleRule::Allowed,
82        notes: "The compatibility target documents positive-integer model-array subscripts. They are structural selectors and do not enter pole computation.",
83    }];
84pub const POLE_INTEGER_CAPABILITIES: [BuiltinIntegerCapabilityDescriptor; 1] =
85    [BuiltinIntegerCapabilityDescriptor {
86        form: "P = pole(sys,integer_J1,...,integer_JN)",
87        inputs: &POLE_INTEGER_MODEL_INDEX_INPUT,
88        computation_domain: BuiltinIntegerComputationDomain::Structural,
89        output_class: BuiltinIntegerOutputClassRule::FunctionSpecific,
90        overflow: BuiltinIntegerOverflowRule::Error,
91        backend: BuiltinIntegerBackendRule::HostOnly,
92        overload: BuiltinIntegerOverloadKind::StructuralParameter,
93        notes: "The compatibility contract is exact one-based model selection. RunMat currently implements singular tf/ss objects only, so the broader model-array overload rejects before any integer conversion can occur.",
94    }];
95
96#[runmat_macros::register_gpu_spec(builtin_path = "crate::builtins::control::pole")]
97pub const GPU_SPEC: BuiltinGpuSpec = BuiltinGpuSpec {
98    name: "pole",
99    op_kind: GpuOpKind::Custom("control-poles"),
100    supported_precisions: &[],
101    broadcast: BroadcastSemantics::None,
102    provider_hooks: &[],
103    constant_strategy: ConstantStrategy::InlineLiteral,
104    residency: ResidencyPolicy::GatherImmediately,
105    nan_mode: ReductionNaN::Include,
106    two_pass_threshold: None,
107    workgroup_size: None,
108    accepts_nan_mode: false,
109    notes: "pole computes roots or state-matrix eigenvalues from host-side model metadata.",
110};
111
112#[runmat_macros::register_fusion_spec(builtin_path = "crate::builtins::control::pole")]
113pub const FUSION_SPEC: BuiltinFusionSpec = BuiltinFusionSpec {
114    name: "pole",
115    shape: ShapeRequirements::Any,
116    constant_strategy: ConstantStrategy::InlineLiteral,
117    elementwise: None,
118    reduction: None,
119    emits_nan: false,
120    notes: "pole is model analysis and is not fused.",
121};
122
123#[runtime_builtin(
124    name = "pole",
125    category = "control",
126    summary = "Return poles of transfer-function and state-space control models.",
127    keywords = "pole,poles,control system,stability,transfer function,state space,tf,ss",
128    type_resolver(pole_type),
129    descriptor(crate::builtins::control::pole::POLE_DESCRIPTOR),
130    integer_capabilities(crate::builtins::control::pole::POLE_INTEGER_CAPABILITIES),
131    builtin_path = "crate::builtins::control::pole"
132)]
133async fn pole_builtin(sys: Value) -> BuiltinResult<Value> {
134    let gathered = dispatcher::gather_if_needed_async(&sys).await?;
135    let poles = match gathered {
136        Value::Object(object) if object.is_class(TF_CLASS) => {
137            TfModel::from_value(Value::Object(object), "pole")?.poles()?
138        }
139        Value::Object(object) if object.is_class(SS_CLASS) => {
140            ss_poles_from_object(&object, "pole")?.0
141        }
142        Value::Object(object) => {
143            return Err(control_error(
144                "pole",
145                "RunMat:pole:UnsupportedModel",
146                format!(
147                    "pole: unsupported model class '{}'; supported classes are tf and ss",
148                    object.class_name
149                ),
150            ));
151        }
152        other => {
153            return Err(control_error(
154                "pole",
155                "RunMat:pole:InvalidModel",
156                format!("pole: expected a tf or ss object, got {other:?}"),
157            ));
158        }
159    };
160    output_complex_column(poles, "pole")
161}
162
163#[cfg(test)]
164mod tests {
165    use super::*;
166    use futures::executor::block_on;
167    use runmat_value::Tensor;
168
169    #[test]
170    fn pole_returns_roots_of_denominator() {
171        let sys = block_on(crate::call_builtin_async(
172            "tf",
173            &[
174                Value::Num(1.0),
175                Value::Tensor(Tensor::new(vec![1.0, 3.0, 2.0], vec![1, 3]).unwrap()),
176            ],
177        ))
178        .expect("tf");
179        let Value::Tensor(poles) = block_on(pole_builtin(sys)).expect("pole") else {
180            panic!("expected real poles");
181        };
182        assert_eq!(poles.shape, vec![2, 1]);
183        assert!(poles
184            .materialize_f64()
185            .iter()
186            .any(|p| (*p + 1.0).abs() < 1.0e-8));
187        assert!(poles
188            .materialize_f64()
189            .iter()
190            .any(|p| (*p + 2.0).abs() < 1.0e-8));
191    }
192
193    #[test]
194    fn pole_returns_repeated_roots_of_denominator() {
195        let sys = block_on(crate::call_builtin_async(
196            "tf",
197            &[
198                Value::Num(1.0),
199                Value::Tensor(Tensor::new(vec![1.0, 2.0, 1.0], vec![1, 3]).unwrap()),
200            ],
201        ))
202        .expect("tf");
203        let Value::Tensor(poles) = block_on(pole_builtin(sys)).expect("pole") else {
204            panic!("expected real poles");
205        };
206        assert_eq!(poles.shape, vec![2, 1]);
207        assert!(poles
208            .materialize_f64()
209            .iter()
210            .all(|p| (*p + 1.0).abs() < 1.0e-8));
211    }
212
213    #[test]
214    fn pole_returns_complex_conjugate_roots() {
215        let sys = block_on(crate::call_builtin_async(
216            "tf",
217            &[
218                Value::Num(1.0),
219                Value::Tensor(Tensor::new(vec![1.0, 0.0, 1.0], vec![1, 3]).unwrap()),
220            ],
221        ))
222        .expect("tf");
223        let Value::ComplexTensor(poles) = block_on(pole_builtin(sys)).expect("pole") else {
224            panic!("expected complex poles");
225        };
226        assert_eq!(poles.shape, vec![2, 1]);
227        assert!(poles
228            .materialize_f64()
229            .iter()
230            .all(|(re, _)| re.abs() < 1.0e-8));
231        assert!(poles
232            .materialize_f64()
233            .iter()
234            .any(|(_, im)| (*im - 1.0).abs() < 1.0e-8));
235        assert!(poles
236            .materialize_f64()
237            .iter()
238            .any(|(_, im)| (*im + 1.0).abs() < 1.0e-8));
239    }
240
241    #[test]
242    fn pole_uses_state_matrix_eigenvalues_for_ss() {
243        let sys = block_on(crate::call_builtin_async(
244            "ss",
245            &[
246                Value::Tensor(Tensor::new(vec![0.0, -4.0, 1.0, -0.5], vec![2, 2]).unwrap()),
247                Value::Tensor(Tensor::new(vec![0.0, 1.0], vec![2, 1]).unwrap()),
248                Value::Tensor(Tensor::new(vec![1.0, 0.0], vec![1, 2]).unwrap()),
249                Value::Num(0.0),
250            ],
251        ))
252        .expect("ss");
253        let Value::ComplexTensor(poles) = block_on(pole_builtin(sys)).expect("pole") else {
254            panic!("expected complex poles");
255        };
256        assert_eq!(poles.shape, vec![2, 1]);
257        assert!(poles
258            .materialize_f64()
259            .iter()
260            .all(|(re, _)| (*re + 0.25).abs() < 1.0e-8));
261        assert!(poles
262            .materialize_f64()
263            .iter()
264            .any(|(_, im)| (*im - 1.984313483298443).abs() < 1.0e-8));
265        assert!(poles
266            .materialize_f64()
267            .iter()
268            .any(|(_, im)| (*im + 1.984313483298443).abs() < 1.0e-8));
269    }
270}