1use runmat_builtins::{
4 BuiltinCompletionPolicy, BuiltinDescriptor, BuiltinErrorDescriptor, BuiltinOutputMode,
5 BuiltinParamArity, BuiltinParamDescriptor, BuiltinParamType, BuiltinSignatureDescriptor,
6};
7use runmat_macros::runtime_builtin;
8use runmat_value::Value;
9
10use crate::builtins::common::spec::{
11 BroadcastSemantics, BuiltinFusionSpec, BuiltinGpuSpec, ConstantStrategy, GpuOpKind,
12 ReductionNaN, ResidencyPolicy, ShapeRequirements,
13};
14use crate::builtins::math::ode::common::{
15 build_ode_output, define_ode_integer_contract, ode_options_from_struct, parse_ode_input,
16 parse_options, prepare_ode_options, solve_ode, OdeMethod,
17};
18use crate::builtins::math::ode::type_resolvers::ode_solution_type;
19use crate::{build_runtime_error, BuiltinResult, RuntimeError};
20
21const NAME: &str = "ode45";
22
23define_ode_integer_contract!("ode45", "Ode45");
24
25const ODE45_OUTPUT_Y: [BuiltinParamDescriptor; 1] = [BuiltinParamDescriptor {
26 name: "y",
27 ty: BuiltinParamType::NumericArray,
28 arity: BuiltinParamArity::Required,
29 default: None,
30 description: "Solution states evaluated over tspan.",
31}];
32
33const ODE45_OUTPUT_TY: [BuiltinParamDescriptor; 2] = [
34 BuiltinParamDescriptor {
35 name: "t",
36 ty: BuiltinParamType::NumericArray,
37 arity: BuiltinParamArity::Required,
38 default: None,
39 description: "Time points selected by solver.",
40 },
41 BuiltinParamDescriptor {
42 name: "y",
43 ty: BuiltinParamType::NumericArray,
44 arity: BuiltinParamArity::Required,
45 default: None,
46 description: "Solution states at each returned time point.",
47 },
48];
49
50const ODE45_INPUTS_CORE: [BuiltinParamDescriptor; 3] = [
51 BuiltinParamDescriptor {
52 name: "odefun",
53 ty: BuiltinParamType::Any,
54 arity: BuiltinParamArity::Required,
55 default: None,
56 description: "ODE right-hand-side callback f(t,y).",
57 },
58 BuiltinParamDescriptor {
59 name: "tspan",
60 ty: BuiltinParamType::Any,
61 arity: BuiltinParamArity::Required,
62 default: None,
63 description: "Time interval or monotonic time vector.",
64 },
65 BuiltinParamDescriptor {
66 name: "y0",
67 ty: BuiltinParamType::Any,
68 arity: BuiltinParamArity::Required,
69 default: None,
70 description: "Initial state vector/value.",
71 },
72];
73
74const ODE45_INPUTS_WITH_OPTIONS: [BuiltinParamDescriptor; 4] = [
75 BuiltinParamDescriptor {
76 name: "odefun",
77 ty: BuiltinParamType::Any,
78 arity: BuiltinParamArity::Required,
79 default: None,
80 description: "ODE right-hand-side callback f(t,y).",
81 },
82 BuiltinParamDescriptor {
83 name: "tspan",
84 ty: BuiltinParamType::Any,
85 arity: BuiltinParamArity::Required,
86 default: None,
87 description: "Time interval or monotonic time vector.",
88 },
89 BuiltinParamDescriptor {
90 name: "y0",
91 ty: BuiltinParamType::Any,
92 arity: BuiltinParamArity::Required,
93 default: None,
94 description: "Initial state vector/value.",
95 },
96 BuiltinParamDescriptor {
97 name: "options",
98 ty: BuiltinParamType::Any,
99 arity: BuiltinParamArity::Optional,
100 default: None,
101 description: "Optional struct with tolerances and step controls.",
102 },
103];
104
105const ODE45_SIGNATURES: [BuiltinSignatureDescriptor; 4] = [
106 BuiltinSignatureDescriptor {
107 label: "y = ode45(odefun, tspan, y0)",
108 inputs: &ODE45_INPUTS_CORE,
109 outputs: &ODE45_OUTPUT_Y,
110 },
111 BuiltinSignatureDescriptor {
112 label: "y = ode45(odefun, tspan, y0, options)",
113 inputs: &ODE45_INPUTS_WITH_OPTIONS,
114 outputs: &ODE45_OUTPUT_Y,
115 },
116 BuiltinSignatureDescriptor {
117 label: "[t, y] = ode45(odefun, tspan, y0)",
118 inputs: &ODE45_INPUTS_CORE,
119 outputs: &ODE45_OUTPUT_TY,
120 },
121 BuiltinSignatureDescriptor {
122 label: "[t, y] = ode45(odefun, tspan, y0, options)",
123 inputs: &ODE45_INPUTS_WITH_OPTIONS,
124 outputs: &ODE45_OUTPUT_TY,
125 },
126];
127
128const ODE45_ERROR_INVALID_ARGUMENT: BuiltinErrorDescriptor = BuiltinErrorDescriptor {
129 code: "RM.ODE45.INVALID_ARGUMENT",
130 identifier: Some("RunMat:ode45:InvalidArgument"),
131 when: "Input argument count/options struct grammar is invalid.",
132 message: "ode45: invalid argument",
133};
134
135const ODE45_ERROR_INVALID_INPUT: BuiltinErrorDescriptor = BuiltinErrorDescriptor {
136 code: "RM.ODE45.INVALID_INPUT",
137 identifier: Some("RunMat:ode45:InvalidInput"),
138 when: "ODE input/state/callback semantics are invalid for integration.",
139 message: "ode45: invalid input",
140};
141
142const ODE45_ERROR_INTERNAL: BuiltinErrorDescriptor = BuiltinErrorDescriptor {
143 code: "RM.ODE45.INTERNAL",
144 identifier: Some("RunMat:ode45:Internal"),
145 when: "Internal output materialization fails.",
146 message: "ode45: internal runtime failure",
147};
148
149const ODE45_ERRORS: [BuiltinErrorDescriptor; 3] = [
150 ODE45_ERROR_INVALID_ARGUMENT,
151 ODE45_ERROR_INVALID_INPUT,
152 ODE45_ERROR_INTERNAL,
153];
154
155pub const ODE45_DESCRIPTOR: BuiltinDescriptor = BuiltinDescriptor {
156 signatures: &ODE45_SIGNATURES,
157 output_mode: BuiltinOutputMode::ByRequestedOutputCount,
158 completion_policy: BuiltinCompletionPolicy::Public,
159 errors: &ODE45_ERRORS,
160};
161
162fn ode45_error_with_detail(
163 error: &'static BuiltinErrorDescriptor,
164 detail: impl AsRef<str>,
165) -> RuntimeError {
166 let detail = detail.as_ref();
167 let message = if detail.starts_with("ode45:") {
168 detail.to_string()
169 } else {
170 format!("{}: {}", error.message, detail)
171 };
172 let mut builder = build_runtime_error(message).with_builtin(NAME);
173 if let Some(identifier) = error.identifier {
174 builder = builder.with_identifier(identifier);
175 }
176 builder.build()
177}
178
179fn ode45_map_error(err: RuntimeError, fallback: &'static BuiltinErrorDescriptor) -> RuntimeError {
180 if err.identifier().is_some() {
181 err
182 } else {
183 ode45_error_with_detail(fallback, err.message())
184 }
185}
186
187#[runmat_macros::register_gpu_spec(builtin_path = "crate::builtins::math::ode::ode45")]
188pub const GPU_SPEC: BuiltinGpuSpec = BuiltinGpuSpec {
189 name: "ode45",
190 op_kind: GpuOpKind::Custom("ode-solve"),
191 supported_precisions: &[],
192 broadcast: BroadcastSemantics::None,
193 provider_hooks: &[],
194 constant_strategy: ConstantStrategy::InlineLiteral,
195 residency: ResidencyPolicy::GatherImmediately,
196 nan_mode: ReductionNaN::Include,
197 two_pass_threshold: None,
198 workgroup_size: None,
199 accepts_nan_mode: false,
200 notes: "Adaptive ODE integration runs on the host. RHS callbacks may call GPU-aware builtins.",
201};
202
203#[runmat_macros::register_fusion_spec(builtin_path = "crate::builtins::math::ode::ode45")]
204pub const FUSION_SPEC: BuiltinFusionSpec = BuiltinFusionSpec {
205 name: "ode45",
206 shape: ShapeRequirements::Any,
207 constant_strategy: ConstantStrategy::InlineLiteral,
208 elementwise: None,
209 reduction: None,
210 emits_nan: false,
211 notes: "ODE integration repeatedly invokes user callbacks and terminates fusion planning.",
212};
213
214#[runtime_builtin(
215 name = "ode45",
216 category = "math/ode",
217 summary = "Solve nonstiff ODE systems using adaptive Dormand-Prince 5(4) integration.",
218 keywords = "ode45,ode,nonstiff,dormand-prince,adaptive step",
219 accel = "sink",
220 type_resolver(ode_solution_type),
221 descriptor(crate::builtins::math::ode::ode45::ODE45_DESCRIPTOR),
222 extensions(crate::builtins::math::ode::ode45::EXTENSIONS),
223 integer_capabilities(crate::builtins::math::ode::ode45::INTEGER_CAPABILITIES),
224 builtin_path = "crate::builtins::math::ode::ode45"
225)]
226async fn ode45_builtin(
227 function: Value,
228 tspan: Value,
229 y0: Value,
230 rest: Vec<Value>,
231) -> BuiltinResult<Value> {
232 if rest.len() > 1 {
233 return Err(ode45_error_with_detail(
234 &ODE45_ERROR_INVALID_ARGUMENT,
235 "too many input arguments",
236 ));
237 }
238 let options = parse_options(NAME, rest.first())
239 .map_err(|err| ode45_map_error(err, &ODE45_ERROR_INVALID_ARGUMENT))?;
240 let options = prepare_ode_options(NAME, options, ODE_COMPATIBILITY_EXTENSIONS)
241 .await
242 .map_err(|err| ode45_map_error(err, &ODE45_ERROR_INVALID_ARGUMENT))?;
243 let opts = ode_options_from_struct(NAME, options.as_ref())
244 .map_err(|err| ode45_map_error(err, &ODE45_ERROR_INVALID_ARGUMENT))?;
245 let input = parse_ode_input(NAME, tspan, y0, ODE_COMPATIBILITY_EXTENSIONS)
246 .await
247 .map_err(|err| ode45_map_error(err, &ODE45_ERROR_INVALID_INPUT))?;
248 let result = solve_ode(NAME, OdeMethod::Ode45, &function, &input, &opts)
249 .await
250 .map_err(|err| ode45_map_error(err, &ODE45_ERROR_INVALID_INPUT))?;
251 build_ode_output(NAME, result).map_err(|err| ode45_map_error(err, &ODE45_ERROR_INTERNAL))
252}
253
254#[cfg(test)]
255mod tests {
256 use super::*;
257 use crate::builtins::common::test_support;
258 use futures::executor::block_on;
259 use runmat_accelerate_api::HostTensorView;
260 use runmat_value::Tensor;
261 use runmat_value::{IntValue, IntegerStorage};
262 use std::sync::Arc;
263
264 #[test]
265 fn ode45_scalar_decay_returns_reasonable_final_value() {
266 let _resolver =
267 crate::user_functions::install_semantic_function_resolver(Some(Arc::new(|_name| {
268 Some(0)
269 })));
270 let _invoker = crate::user_functions::install_semantic_function_invoker(Some(Arc::new(
271 move |_function, args, _requested_outputs| {
272 let y = match &args[1] {
273 Value::Num(n) => *n,
274 other => panic!("expected scalar state, got {other:?}"),
275 };
276 Box::pin(async move { Ok(Value::Num(-y)) })
277 },
278 )));
279
280 let out = block_on(ode45_builtin(
281 Value::FunctionHandle("decay".into()),
282 Value::Tensor(Tensor::new(vec![0.0, 1.0], vec![1, 2]).unwrap()),
283 Value::Num(1.0),
284 Vec::new(),
285 ))
286 .unwrap();
287
288 match out {
289 Value::Tensor(t) => {
290 assert_eq!(t.cols(), 1);
291 let last = t.materialize_f64()[t.rows() - 1];
292 assert!((last - (-1.0_f64).exp()).abs() < 5.0e-3);
293 }
294 other => panic!("unexpected output {other:?}"),
295 }
296 }
297
298 #[test]
299 fn ode45_strict_mode_rejects_integer_tspan_before_rhs() {
300 let _strict = crate::compatibility::push_runmat_extensions_enabled(false);
301 let tspan = Tensor::new_integer(IntegerStorage::U16(vec![0, 1]), vec![1, 2]).unwrap();
302
303 let error = block_on(ode45_builtin(
304 Value::FunctionHandle("unused".into()),
305 Value::Tensor(tspan),
306 Value::Num(1.0),
307 Vec::new(),
308 ))
309 .expect_err("integer tspan is a RunMat-only extension");
310
311 assert_eq!(error.identifier(), INTEGER_TSPAN_EXTENSION.error_identifier);
312 }
313
314 #[test]
315 fn ode45_runmat_mode_rejects_wide_integer_tspan() {
316 let _extensions = crate::compatibility::push_runmat_extensions_enabled(true);
317 let tspan =
318 Tensor::new_integer(IntegerStorage::U64(vec![0, (1_u64 << 53) + 1]), vec![1, 2])
319 .unwrap();
320
321 let error = block_on(ode45_builtin(
322 Value::FunctionHandle("unused".into()),
323 Value::Tensor(tspan),
324 Value::Num(1.0),
325 Vec::new(),
326 ))
327 .expect_err("wide integer tspan cannot cross exactly");
328
329 assert!(error.message().contains("exactly representable"));
330 }
331
332 #[test]
333 fn ode45_strict_mode_rejects_integer_derivative_result() {
334 let _strict = crate::compatibility::push_runmat_extensions_enabled(false);
335 let _resolver =
336 crate::user_functions::install_semantic_function_resolver(Some(Arc::new(|_name| {
337 Some(903)
338 })));
339 let _invoker = crate::user_functions::install_semantic_function_invoker(Some(Arc::new(
340 |_function, _args, _requested_outputs| {
341 Box::pin(async move { Ok(Value::Int(IntValue::I32(-1))) })
342 },
343 )));
344
345 let error = block_on(ode45_builtin(
346 Value::FunctionHandle("integer_rhs".into()),
347 Value::Tensor(Tensor::new(vec![0.0, 1.0], vec![1, 2]).unwrap()),
348 Value::Num(1.0),
349 Vec::new(),
350 ))
351 .expect_err("integer derivative is a RunMat-only extension");
352
353 assert_eq!(
354 error.identifier(),
355 INTEGER_CALLBACK_EXTENSION.error_identifier
356 );
357 }
358
359 #[test]
360 fn ode45_automatic_resident_input_gathers_but_explicit_input_is_gated() {
361 test_support::with_test_provider(|provider| {
362 let _invoker = crate::user_functions::install_semantic_function_invoker(Some(
363 Arc::new(|_function, _args, _requested_outputs| {
364 Box::pin(async move { Ok(Value::Num(-1.0)) })
365 }),
366 ));
367 let times = [0.0, 0.1];
368 let shape = [1, 2];
369 let automatic = provider
370 .upload(&HostTensorView {
371 data: ×,
372 shape: &shape,
373 })
374 .expect("automatic upload");
375 let automatic =
376 automatic.with_provenance(runmat_accelerate_api::GpuHandleProvenance::Automatic);
377 let result = block_on(ode45_builtin(
378 Value::BoundFunctionHandle {
379 name: "constant_rhs".to_string(),
380 function: 905,
381 },
382 Value::GpuTensor(automatic),
383 Value::Num(1.0),
384 Vec::new(),
385 ))
386 .expect("automatic resident tspan gathers");
387 assert!(matches!(result, Value::Tensor(_)));
388
389 let explicit = provider
390 .upload(&HostTensorView {
391 data: ×,
392 shape: &shape,
393 })
394 .expect("explicit upload");
395 let explicit =
396 explicit.with_provenance(runmat_accelerate_api::GpuHandleProvenance::Explicit);
397 let _strict = crate::compatibility::push_runmat_extensions_enabled(false);
398 let error = block_on(ode45_builtin(
399 Value::FunctionHandle("unused".into()),
400 Value::GpuTensor(explicit),
401 Value::Num(1.0),
402 Vec::new(),
403 ))
404 .expect_err("explicit resident tspan is gated before fallback");
405 assert_eq!(
406 error.identifier(),
407 RESIDENT_INPUT_EXTENSION.error_identifier
408 );
409 });
410 }
411
412 #[test]
413 fn ode45_rejects_nan_rhs() {
414 let _resolver =
415 crate::user_functions::install_semantic_function_resolver(Some(Arc::new(|_name| {
416 Some(0)
417 })));
418 let _invoker = crate::user_functions::install_semantic_function_invoker(Some(Arc::new(
419 move |_function, _args, _requested_outputs| {
420 Box::pin(async move { Ok(Value::Num(f64::NAN)) })
421 },
422 )));
423
424 let err = block_on(ode45_builtin(
425 Value::FunctionHandle("nan_rhs".into()),
426 Value::Tensor(Tensor::new(vec![0.0, 1.0], vec![1, 2]).unwrap()),
427 Value::Num(1.0),
428 Vec::new(),
429 ))
430 .expect_err("ode45 should reject NaN derivative values");
431
432 assert!(err.to_string().contains("function value must be finite"));
433 }
434
435 #[test]
436 fn ode45_accepts_external_function_handle_rhs() {
437 let _resolver =
438 crate::user_functions::install_semantic_function_resolver(Some(Arc::new(|name| {
439 (name == "pkg.decay").then_some(56)
440 })));
441 let _invoker = crate::user_functions::install_semantic_function_invoker(Some(Arc::new(
442 move |function, args, _requested_outputs| {
443 assert_eq!(function, 56);
444 let y = match &args[1] {
445 Value::Num(n) => *n,
446 other => panic!("expected scalar state, got {other:?}"),
447 };
448 Box::pin(async move { Ok(Value::Num(-y)) })
449 },
450 )));
451
452 let out = block_on(ode45_builtin(
453 Value::ExternalFunctionHandle("pkg.decay".to_string()),
454 Value::Tensor(Tensor::new(vec![0.0, 1.0], vec![1, 2]).unwrap()),
455 Value::Num(1.0),
456 Vec::new(),
457 ))
458 .unwrap();
459
460 match out {
461 Value::Tensor(t) => {
462 assert_eq!(t.cols(), 1);
463 let last = t.materialize_f64()[t.rows() - 1];
464 assert!(last.is_finite());
465 assert!(last > 0.0);
466 assert!(last < 1.0);
467 }
468 other => panic!("unexpected output {other:?}"),
469 }
470 }
471
472 #[test]
473 fn ode45_too_many_inputs_uses_stable_identifier() {
474 let err = block_on(ode45_builtin(
475 Value::FunctionHandle("decay".into()),
476 Value::Tensor(Tensor::new(vec![0.0, 1.0], vec![1, 2]).unwrap()),
477 Value::Num(1.0),
478 vec![Value::Num(1.0), Value::Num(2.0)],
479 ))
480 .expect_err("expected too many inputs error");
481 assert_eq!(err.identifier(), ODE45_ERROR_INVALID_ARGUMENT.identifier);
482 }
483
484 #[test]
485 fn ode45_descriptor_signatures_cover_surface() {
486 let labels: Vec<&str> = ODE45_DESCRIPTOR
487 .signatures
488 .iter()
489 .map(|signature| signature.label)
490 .collect();
491 assert_eq!(
492 labels,
493 vec![
494 "y = ode45(odefun, tspan, y0)",
495 "y = ode45(odefun, tspan, y0, options)",
496 "[t, y] = ode45(odefun, tspan, y0)",
497 "[t, y] = ode45(odefun, tspan, y0, options)",
498 ]
499 );
500 }
501
502 #[test]
503 fn ode45_descriptor_errors_have_stable_codes() {
504 let codes: Vec<&str> = ODE45_DESCRIPTOR
505 .errors
506 .iter()
507 .map(|error| error.code)
508 .collect();
509 assert_eq!(
510 codes,
511 vec![
512 "RM.ODE45.INVALID_ARGUMENT",
513 "RM.ODE45.INVALID_INPUT",
514 "RM.ODE45.INTERNAL",
515 ]
516 );
517 }
518}