fastsim_core/utils/
interp.rs

1use crate::imports::*;
2
3/// Methods for mutating interpolator data, e.g. proportionally scaling
4/// interpolator function data
5pub trait InterpolatorMutMethods {
6    fn set_min(&mut self, min: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()>;
7    fn set_max(&mut self, max: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()>;
8    fn set_range(&mut self, range: f64) -> anyhow::Result<()>;
9}
10
11impl InterpolatorMutMethods for Interp0D<f64> {
12    fn set_min(&mut self, min: f64, _scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
13        self.0 = min;
14        Ok(())
15    }
16
17    fn set_max(&mut self, max: f64, _scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
18        self.0 = max;
19        Ok(())
20    }
21
22    fn set_range(&mut self, _range: f64) -> anyhow::Result<()> {
23        bail!("Cannot set range for 0D interpolator")
24    }
25}
26
27impl<S> InterpolatorMutMethods for Interp1DOwned<f64, S>
28where
29    S: ninterp::strategy::traits::Strategy1D<ndarray::OwnedRepr<f64>> + Clone,
30{
31    fn set_min(&mut self, min: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
32        let old_min = *self.min()?;
33        match scaling.unwrap_or_default() {
34            utils::interp::ScalingMethods::Proportional => {
35                ensure!(
36                    old_min != 0.,
37                    "Cannot modify min proportionally when old_min == 0."
38                );
39                self.data.values.map_inplace(|v| *v *= min / old_min);
40            }
41            utils::interp::ScalingMethods::AnchoredProportional => {
42                todo!()
43            }
44            utils::interp::ScalingMethods::Offset => {
45                todo!()
46            }
47        }
48        self.validate()?;
49        Ok(())
50    }
51
52    fn set_max(&mut self, max: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
53        let old_max = *self.max()?;
54        match scaling.unwrap_or_default() {
55            utils::interp::ScalingMethods::Proportional => {
56                ensure!(
57                    old_max != 0.,
58                    "Cannot modify max proportionally when old_max == 0."
59                );
60                self.data.values.map_inplace(|v| *v *= max / old_max);
61            }
62            utils::interp::ScalingMethods::AnchoredProportional => {
63                todo!()
64            }
65            utils::interp::ScalingMethods::Offset => {
66                todo!()
67            }
68        }
69        self.validate()?;
70        Ok(())
71    }
72
73    fn set_range(&mut self, range: f64) -> anyhow::Result<()> {
74        let old_max = *self.max()?;
75        let old_range = old_max - self.min()?;
76        ensure!(old_range != 0., "Cannot modify range when min == max");
77        // if the new range is 0., chooses the max as the value for all elements of the array
78        if range == 0. {
79            self.data.values = self.data.values.map(|_| old_max);
80        } else {
81            self.data.values = self
82                .data
83                .values
84                .map(|x| old_max + (x - old_max) * range / old_range);
85        }
86        self.validate()?;
87        Ok(())
88    }
89}
90
91impl<S> InterpolatorMutMethods for Interp2DOwned<f64, S>
92where
93    S: ninterp::strategy::traits::Strategy2D<ndarray::OwnedRepr<f64>> + Clone,
94{
95    fn set_min(&mut self, min: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
96        let old_min = *self.min()?;
97        match scaling.unwrap_or_default() {
98            utils::interp::ScalingMethods::Proportional => {
99                ensure!(
100                    old_min != 0.,
101                    "Cannot modify min proportionally when old_min == 0."
102                );
103                self.data.values.map_inplace(|v| *v *= min / old_min);
104            }
105            utils::interp::ScalingMethods::AnchoredProportional => {
106                todo!()
107            }
108            utils::interp::ScalingMethods::Offset => {
109                todo!()
110            }
111        }
112        self.validate()?;
113        Ok(())
114    }
115
116    fn set_max(&mut self, max: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
117        let old_max = *self.max()?;
118        match scaling.unwrap_or_default() {
119            utils::interp::ScalingMethods::Proportional => {
120                ensure!(
121                    old_max != 0.,
122                    "Cannot modify max proportionally when old_max == 0."
123                );
124                self.data.values.map_inplace(|v| *v *= max / old_max);
125            }
126            utils::interp::ScalingMethods::AnchoredProportional => {
127                todo!()
128            }
129            utils::interp::ScalingMethods::Offset => {
130                todo!()
131            }
132        }
133        self.validate()?;
134        Ok(())
135    }
136
137    fn set_range(&mut self, range: f64) -> anyhow::Result<()> {
138        let old_max = *self.max()?;
139        let old_range = old_max - self.min()?;
140        ensure!(old_range != 0., "Cannot modify range when min == max");
141        // if the new range is 0., chooses the max as the value for all elements of the array
142        if range == 0. {
143            self.data.values = self.data.values.map(|_| old_max);
144        } else {
145            self.data.values = self
146                .data
147                .values
148                .map(|x| old_max + (x - old_max) * range / old_range);
149        }
150        self.validate()?;
151        Ok(())
152    }
153}
154
155impl<S> InterpolatorMutMethods for Interp3DOwned<f64, S>
156where
157    S: ninterp::strategy::traits::Strategy3D<ndarray::OwnedRepr<f64>> + Clone,
158{
159    fn set_min(&mut self, min: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
160        let old_min = *self.min()?;
161        match scaling.unwrap_or_default() {
162            utils::interp::ScalingMethods::Proportional => {
163                ensure!(
164                    old_min != 0.,
165                    "Cannot modify min proportionally when old_min == 0."
166                );
167                self.data.values.map_inplace(|v| *v *= min / old_min);
168            }
169            utils::interp::ScalingMethods::AnchoredProportional => {
170                todo!()
171            }
172            utils::interp::ScalingMethods::Offset => {
173                todo!()
174            }
175        }
176        self.validate()?;
177        Ok(())
178    }
179
180    fn set_max(&mut self, max: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
181        let old_max = *self.max()?;
182        match scaling.unwrap_or_default() {
183            utils::interp::ScalingMethods::Proportional => {
184                ensure!(
185                    old_max != 0.,
186                    "Cannot modify max proportionally when old_max == 0."
187                );
188                self.data.values.map_inplace(|v| *v *= max / old_max);
189            }
190            utils::interp::ScalingMethods::AnchoredProportional => {
191                todo!()
192            }
193            utils::interp::ScalingMethods::Offset => {
194                todo!()
195            }
196        }
197        self.validate()?;
198        Ok(())
199    }
200
201    fn set_range(&mut self, range: f64) -> anyhow::Result<()> {
202        let old_max = *self.max()?;
203        let old_range = old_max - self.min()?;
204        ensure!(old_range != 0., "Cannot modify range when min == max");
205        // if the new range is 0., chooses the max as the value for all elements of the array
206        if range == 0. {
207            self.data.values = self.data.values.map(|_| old_max);
208        } else {
209            self.data.values = self
210                .data
211                .values
212                .map(|x| old_max + (x - old_max) * range / old_range);
213        }
214        self.validate()?;
215        Ok(())
216    }
217}
218
219impl<S> InterpolatorMutMethods for InterpNDOwned<f64, S>
220where
221    S: ninterp::strategy::traits::StrategyND<ndarray::OwnedRepr<f64>> + Clone,
222{
223    fn set_min(&mut self, min: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
224        let old_min = *self.min()?;
225        match scaling.unwrap_or_default() {
226            utils::interp::ScalingMethods::Proportional => {
227                ensure!(
228                    old_min != 0.,
229                    "Cannot modify min proportionally when old_min == 0."
230                );
231                self.data.values.map_inplace(|v| *v *= min / old_min);
232            }
233            utils::interp::ScalingMethods::AnchoredProportional => {
234                todo!()
235            }
236            utils::interp::ScalingMethods::Offset => {
237                todo!()
238            }
239        }
240        self.validate()?;
241        Ok(())
242    }
243
244    fn set_max(&mut self, max: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
245        let old_max = *self.max()?;
246        match scaling.unwrap_or_default() {
247            utils::interp::ScalingMethods::Proportional => {
248                ensure!(
249                    old_max != 0.,
250                    "Cannot modify max proportionally when old_max == 0."
251                );
252                self.data.values.map_inplace(|v| *v *= max / old_max);
253            }
254            utils::interp::ScalingMethods::AnchoredProportional => {
255                todo!()
256            }
257            utils::interp::ScalingMethods::Offset => {
258                todo!()
259            }
260        }
261        self.validate()?;
262        Ok(())
263    }
264
265    fn set_range(&mut self, range: f64) -> anyhow::Result<()> {
266        let old_max = *self.max()?;
267        let old_range = old_max - self.min()?;
268        ensure!(old_range != 0., "Cannot modify range when min == max");
269        // if the new range is 0., chooses the max as the value for all elements of the array
270        if range == 0. {
271            self.data.values = self.data.values.map(|_| old_max);
272        } else {
273            self.data.values = self
274                .data
275                .values
276                .map(|x| old_max + (x - old_max) * range / old_range);
277        }
278        self.validate()?;
279        Ok(())
280    }
281}
282
283// This can be made more generic by using a `ninterp::num_traits` bound instead of f64
284// If there are future methods that *do not* mutate the interpolator,
285// we should define a new trait and impl it for `InterpolatorEnum<D> where D: ndarray::Data`
286impl InterpolatorMutMethods for InterpolatorEnumOwned<f64> {
287    // scale all values so that the min is the new min
288    // (Note: this may change the max, depending on what scaling method is chosen)
289    fn set_min(&mut self, min: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
290        match self {
291            Self::Interp0D(interp) => interp.set_min(min, scaling),
292            Self::Interp1D(interp) => interp.set_min(min, scaling),
293            Self::Interp2D(interp) => interp.set_min(min, scaling),
294            Self::Interp3D(interp) => interp.set_min(min, scaling),
295            Self::InterpND(interp) => interp.set_min(min, scaling),
296        }
297    }
298
299    // scale all values so that the max is the new max
300    // (Note: may change the min, depending on what scaling method is chosen)
301    fn set_max(&mut self, max: f64, scaling: Option<ScalingMethods>) -> anyhow::Result<()> {
302        match self {
303            Self::Interp0D(interp) => interp.set_max(max, scaling),
304            Self::Interp1D(interp) => interp.set_max(max, scaling),
305            Self::Interp2D(interp) => interp.set_max(max, scaling),
306            Self::Interp3D(interp) => interp.set_max(max, scaling),
307            Self::InterpND(interp) => interp.set_max(max, scaling),
308        }
309    }
310
311    fn set_range(&mut self, range: f64) -> anyhow::Result<()> {
312        match self {
313            Self::Interp0D(interp) => interp.set_range(range),
314            Self::Interp1D(interp) => interp.set_range(range),
315            Self::Interp2D(interp) => interp.set_range(range),
316            Self::Interp3D(interp) => interp.set_range(range),
317            Self::InterpND(interp) => interp.set_range(range),
318        }
319    }
320}
321
322impl<D> Init for InterpolatorEnum<D>
323where
324    D: ndarray::Data + ndarray::RawDataClone + Clone,
325    D::Elem: ninterp::num_traits::Num
326        + ninterp::num_traits::Euclid
327        + PartialOrd
328        + Copy
329        + std::fmt::Debug,
330{
331    fn init(&mut self) -> Result<(), Error> {
332        self.validate()
333            .map_err(|e| Error::NinterpError(e.to_string()))
334    }
335}
336impl<D> SerdeAPI for InterpolatorEnum<D>
337where
338    D: ndarray::Data + ndarray::RawDataClone + Clone + ndarray::DataOwned,
339    D::Elem: ninterp::num_traits::Num
340        + ninterp::num_traits::Euclid
341        + PartialOrd
342        + Copy
343        + std::fmt::Debug
344        + Serialize
345        + serde::de::DeserializeOwned,
346{
347    #[cfg(feature = "resources")]
348    const RESOURCES_SUBDIR: &'static str = "interpolators";
349}
350
351#[derive(Default, Debug, Serialize, Deserialize, Clone)]
352pub enum ScalingMethods {
353    #[default]
354    /// Scales everything by the same factor -- e.g. setting min of [1, 2, 3] to 0.5 yields [0.5, 1, 1.5]
355    Proportional,
356    /// Scales proportionally to distance from min/max -- e.g. setting min of [1, 2, 3] to 0.5 yields [0.5, 1.5, 3]
357    AnchoredProportional,
358    /// Scaling by sliding all values up or down by the same offset -- e.g. setting min of [1, 2, 3] to 0.5 yields [0.5, 1.5, 2.5]
359    Offset,
360}
361
362#[cfg(test)]
363mod tests {
364    use super::*;
365
366    #[test]
367    fn test_min() {
368        let x = array![0.05, 0.10, 0.15];
369        let y = array![0.10, 0.20, 0.30];
370        let z = array![0.20, 0.40, 0.60];
371        let f_xy = array![[0.1, 1., 2.], [3., 4., 5.], [6., 7., 8.]];
372        let f_xyz = array![
373            [[0.1, 1., 2.], [3., 4., 5.], [6., 7., 8.]],
374            [[9., 10., 11.], [12., 13., 14.], [15., 16., 17.]],
375            [[18., 19., 20.], [21., 22., 23.], [24., 25., 26.],],
376        ];
377        let mut interp_1d = InterpolatorEnum::new_1d(
378            array![85.0, 90.0],
379            array![0.2, 1.0],
380            strategy::Linear,
381            Extrapolate::Clamp,
382        )
383        .unwrap();
384        let mut interp_2d = InterpolatorEnum::new_2d(
385            x.clone(),
386            y.clone(),
387            f_xy.clone(),
388            strategy::Linear,
389            Extrapolate::Clamp,
390        )
391        .unwrap();
392        let mut interp_3d = InterpolatorEnum::new_3d(
393            x.clone(),
394            y.clone(),
395            z.clone(),
396            f_xyz.clone(),
397            strategy::Linear,
398            Extrapolate::Clamp,
399        )
400        .unwrap();
401        assert_eq!(interp_1d.min().unwrap(), &0.2);
402        assert_eq!(interp_2d.min().unwrap(), &0.1);
403        assert_eq!(interp_3d.min().unwrap(), &0.1);
404        interp_1d.set_min(0.1, None).unwrap();
405        interp_2d.set_min(0.3, None).unwrap();
406        interp_3d.set_min(0.3, None).unwrap();
407        println!("{:?}", interp_1d.min().unwrap());
408        println!("{:?}", interp_2d.min().unwrap());
409        println!("{:?}", interp_3d.min().unwrap());
410        assert!(almost_eq(*interp_1d.min().unwrap(), 0.1, Some(1e-3)));
411        assert!(almost_eq(*interp_2d.min().unwrap(), 0.3, Some(1e-3)));
412        assert!(almost_eq(*interp_3d.min().unwrap(), 0.3, Some(1e-3)));
413    }
414
415    #[test]
416    fn test_max() {
417        let x = array![0.05, 0.10, 0.15];
418        let y = array![0.10, 0.20, 0.30];
419        let z = array![0.20, 0.40, 0.60];
420        let f_xy = array![[0., 1., 2.], [3., 4., 5.], [6., 7., 8.]];
421        let f_xyz = array![
422            [[0., 1., 2.], [3., 4., 5.], [6., 7., 8.]],
423            [[9., 10., 11.], [12., 13., 14.], [15., 16., 17.]],
424            [[18., 19., 20.], [21., 22., 23.], [24., 25., 26.]],
425        ];
426        let mut interp_1d = InterpolatorEnum::new_1d(
427            array![85.0, 90.0],
428            array![0.2, 1.0],
429            strategy::Linear,
430            Extrapolate::Clamp,
431        )
432        .unwrap();
433        let mut interp_2d = InterpolatorEnum::new_2d(
434            x.clone(),
435            y.clone(),
436            f_xy.clone(),
437            strategy::Linear,
438            Extrapolate::Clamp,
439        )
440        .unwrap();
441        let mut interp_3d = InterpolatorEnum::new_3d(
442            x.clone(),
443            y.clone(),
444            z.clone(),
445            f_xyz.clone(),
446            strategy::Linear,
447            Extrapolate::Clamp,
448        )
449        .unwrap();
450        assert_eq!(interp_1d.max().unwrap(), &1.0);
451        assert_eq!(interp_2d.max().unwrap(), &8.);
452        assert_eq!(interp_3d.max().unwrap(), &26.);
453        interp_1d.set_max(2., None).unwrap();
454        interp_2d.set_max(7., None).unwrap();
455        interp_3d.set_max(5., None).unwrap();
456        println!("{:?}", interp_1d.max().unwrap());
457        println!("{:?}", interp_2d.max().unwrap());
458        println!("{:?}", interp_3d.max().unwrap());
459        assert!(almost_eq(*interp_1d.max().unwrap(), 2., Some(1e-3)));
460        assert!(almost_eq(*interp_2d.max().unwrap(), 7., Some(1e-3)));
461        assert!(almost_eq(*interp_3d.max().unwrap(), 5., Some(1e-3)));
462    }
463
464    #[test]
465    fn test_range() {
466        let x = array![0.05, 0.10, 0.15];
467        let y = array![0.10, 0.20, 0.30];
468        let z = array![0.20, 0.40, 0.60];
469        let f_xy = array![[0., 1., 2.], [3., 4., 5.], [6., 7., 8.]];
470        let f_xyz = array![
471            [[0., 1., 2.], [3., 4., 5.], [6., 7., 8.]],
472            [[9., 10., 11.], [12., 13., 14.], [15., 16., 17.]],
473            [[18., 19., 20.], [21., 22., 23.], [24., 25., 26.]],
474        ];
475        let mut interp_1d = InterpolatorEnum::new_1d(
476            array![85.0, 90.0],
477            array![0.2, 1.0],
478            strategy::Linear,
479            Extrapolate::Clamp,
480        )
481        .unwrap();
482        let mut interp_2d = InterpolatorEnum::new_2d(
483            x.clone(),
484            y.clone(),
485            f_xy.clone(),
486            strategy::Linear,
487            Extrapolate::Clamp,
488        )
489        .unwrap();
490        let mut interp_3d = InterpolatorEnum::new_3d(
491            x.clone(),
492            y.clone(),
493            z.clone(),
494            f_xyz.clone(),
495            strategy::Linear,
496            Extrapolate::Clamp,
497        )
498        .unwrap();
499        assert_eq!(interp_1d.range().unwrap(), 0.8);
500        assert_eq!(interp_2d.range().unwrap(), 8.);
501        assert_eq!(interp_3d.range().unwrap(), 26.);
502        interp_1d.set_range(2.).unwrap();
503        interp_2d.set_range(7.).unwrap();
504        interp_3d.set_range(5.).unwrap();
505        println!("{:?}", interp_1d.range().unwrap());
506        println!("{:?}", interp_2d.range().unwrap());
507        println!("{:?}", interp_3d.range().unwrap());
508        assert!(almost_eq(interp_1d.range().unwrap(), 2., Some(1e-3)));
509        assert!(almost_eq(interp_2d.range().unwrap(), 7., Some(1e-3)));
510        assert!(almost_eq(interp_3d.range().unwrap(), 5., Some(1e-3)));
511    }
512
513    type StructWithResources = InterpolatorEnumOwned<f64>;
514
515    #[test]
516    fn test_resources() {
517        let resource_list = StructWithResources::list_resources().unwrap();
518        assert!(!resource_list.is_empty());
519
520        // verify that resources can all load
521        for resource in resource_list {
522            StructWithResources::from_resource(resource.clone(), false)
523                .with_context(|| format_dbg!(resource))
524                .unwrap();
525        }
526    }
527}
fastsim_core/utils/interp.rs

fastsim_core/utils/
interp.rs
