runmat_runtime/builtins/math/interpolation/
pchip.rs1use runmat_builtins::{ResolveContext, Type, Value};
4use runmat_macros::runtime_builtin;
5
6use crate::builtins::common::spec::{
7 BroadcastSemantics, BuiltinFusionSpec, BuiltinGpuSpec, ConstantStrategy, GpuOpKind,
8 ReductionNaN, ResidencyPolicy, ScalarType, ShapeRequirements,
9};
10
11use super::pp::{
12 build_pchip_pp, evaluate_pp, interp_error, pp_to_value, query_points, series_from_values,
13 Extrapolation,
14};
15
16const NAME: &str = "pchip";
17
18#[runmat_macros::register_gpu_spec(builtin_path = "crate::builtins::math::interpolation::pchip")]
19pub const GPU_SPEC: BuiltinGpuSpec = BuiltinGpuSpec {
20 name: NAME,
21 op_kind: GpuOpKind::Custom("pchip"),
22 supported_precisions: &[ScalarType::F32, ScalarType::F64],
23 broadcast: BroadcastSemantics::Matlab,
24 provider_hooks: &[],
25 constant_strategy: ConstantStrategy::InlineLiteral,
26 residency: ResidencyPolicy::GatherImmediately,
27 nan_mode: ReductionNaN::Include,
28 two_pass_threshold: None,
29 workgroup_size: None,
30 accepts_nan_mode: false,
31 notes: "PCHIP coefficient construction and evaluation currently run on the CPU reference path after gathering GPU inputs.",
32};
33
34#[runmat_macros::register_fusion_spec(builtin_path = "crate::builtins::math::interpolation::pchip")]
35pub const FUSION_SPEC: BuiltinFusionSpec = BuiltinFusionSpec {
36 name: NAME,
37 shape: ShapeRequirements::Any,
38 constant_strategy: ConstantStrategy::InlineLiteral,
39 elementwise: None,
40 reduction: None,
41 emits_nan: false,
42 notes: "PCHIP builds a piecewise-polynomial representation and is not fused.",
43};
44
45fn pchip_type(args: &[Type], _ctx: &ResolveContext) -> Type {
46 match args.len() {
47 0 | 1 => Type::Unknown,
48 2 => Type::Struct { known_fields: None },
49 _ => match args.get(2) {
50 Some(Type::Num | Type::Int | Type::Bool) => Type::Num,
51 Some(Type::Tensor { shape }) | Some(Type::Logical { shape }) => Type::Tensor {
52 shape: shape.clone(),
53 },
54 _ => Type::tensor(),
55 },
56 }
57}
58
59#[runtime_builtin(
60 name = "pchip",
61 category = "math/interpolation",
62 summary = "Shape-preserving piecewise cubic Hermite interpolation.",
63 keywords = "pchip,shape preserving,cubic hermite,interpolation,ppval",
64 accel = "sink",
65 sink = true,
66 type_resolver(pchip_type),
67 builtin_path = "crate::builtins::math::interpolation::pchip"
68)]
69async fn pchip_builtin(x: Value, y: Value, rest: Vec<Value>) -> crate::BuiltinResult<Value> {
70 let series = series_from_values(x, y, NAME).await?;
71 let pp = build_pchip_pp(&series, NAME)?;
72 match rest.len() {
73 0 => pp_to_value(pp, NAME),
74 1 => {
75 let query = query_points(rest.into_iter().next().expect("query"), NAME).await?;
76 evaluate_pp(&pp, &query, &Extrapolation::Extrapolate, NAME)
77 }
78 _ => Err(interp_error(NAME, "pchip: too many input arguments")),
79 }
80}
81
82#[cfg(test)]
83mod tests {
84 use super::*;
85 use futures::executor::block_on;
86 use runmat_builtins::Tensor;
87
88 fn row(values: &[f64]) -> Value {
89 Value::Tensor(Tensor::new(values.to_vec(), vec![1, values.len()]).expect("tensor"))
90 }
91
92 #[test]
93 fn pchip_three_arg_evaluates_monotone_data() {
94 let value = block_on(pchip_builtin(
95 row(&[1.0, 2.0, 3.0, 4.0]),
96 row(&[0.0, 1.0, 1.5, 1.75]),
97 vec![row(&[1.5, 2.5, 3.5])],
98 ))
99 .expect("pchip");
100 let Value::Tensor(tensor) = value else {
101 panic!("expected tensor");
102 };
103 assert!(tensor.data.windows(2).all(|pair| pair[1] >= pair[0]));
104 }
105}