geosets_rs/sets/
interval.rs

1#![allow(unused)]
2use crate::linalg_utils::vector_leq;
3
4use super::*;
5use ndarray_rand::RandomExt;
6use ndarray_rand::rand_distr::{Exp1, Uniform};
7use thiserror::Error;
8
9#[derive(Clone, Debug)]
10#[allow(non_snake_case)]
11pub struct Interval {
12    lb: Array1<f64>,
13    ub: Array1<f64>,
14}
15
16#[derive(Error, Debug)]
17pub enum IntervalError {
18    #[error("Dimensions of lb {lb_dim:?} and ub {ub_dim:?} do not match")]
19    DimensionMismatch { lb_dim: usize, ub_dim: usize },
20    #[error("Lower bound {lb} must be less than or equal to upper bound {ub}")]
21    InvalidBounds { lb: Array1<f64>, ub: Array1<f64> },
22}
23
24#[allow(non_snake_case)]
25impl Interval {
26    pub fn new(lb: Array1<f64>, ub: Array1<f64>) -> Result<Interval, IntervalError> {
27        if lb.dim() != ub.dim() {
28            return Err(IntervalError::DimensionMismatch {
29                lb_dim: lb.dim(),
30                ub_dim: ub.dim(),
31            });
32        }
33
34        if lb.iter().zip(ub.iter()).any(|(a, b)| a > b) {
35            return Err(IntervalError::InvalidBounds { lb, ub });
36        }
37
38        Ok(Interval { lb, ub })
39    }
40
41    pub fn from_random(dim: usize) -> Result<Interval, IntervalError> {
42        let lb = Array1::random(dim, Uniform::new(-1.0, 0.0));
43        let ub = Array1::random(dim, Uniform::new(0.0, 1.0));
44
45        Interval::new(lb, ub)
46    }
47}
48
49#[allow(non_snake_case)]
50impl GeoSet for Interval {
51    fn from_unit_box(dim: usize) -> Self {
52        Interval {
53            lb: -Array::ones(dim),
54            ub: Array::ones(dim),
55        }
56    }
57
58    fn dim(&self) -> usize {
59        self.lb.dim()
60    }
61
62    fn empty(&self) -> Result<bool, SetOperationError> {
63        Ok(false)
64    }
65
66    fn to_vertices(&self) -> Result<Array2<f64>, SetOperationError> {
67        let mut vertices = Array2::zeros((1 << self.dim(), self.dim()));
68        for i in 0..(1 << self.dim()) {
69            for j in 0..self.dim() {
70                vertices[[i, j]] = if (i & (1 << j)) != 0 {
71                    self.ub[j]
72                } else {
73                    self.lb[j]
74                };
75            }
76        }
77        Ok(vertices)
78    }
79
80    fn center(&self) -> Result<Array1<f64>, SetOperationError> {
81        let center = (&self.lb + &self.ub) / 2.;
82        Ok(center)
83    }
84
85    fn support_function(
86        &self,
87        direction: Array1<f64>,
88    ) -> Result<(Array1<f64>, f64), SetOperationError> {
89        // For each dimension, pick ub if direction > 0, else lb
90        let support_vector = self
91            .lb
92            .iter()
93            .zip(self.ub.iter())
94            .zip(direction.iter())
95            .map(|((&lb, &ub), &d)| if d > 0.0 { ub } else { lb })
96            .collect::<Array1<f64>>();
97
98        let support_value = support_vector.dot(&direction);
99        Ok((support_vector, support_value))
100    }
101
102    fn volume(&self) -> Result<f64, SetOperationError> {
103        if self.degenerate() {
104            return Ok(0.0);
105        }
106
107        let volume = self
108            .lb
109            .iter()
110            .zip(self.ub.iter())
111            .map(|(lb, ub)| ub - lb)
112            .product();
113        Ok(volume)
114    }
115
116    fn minkowski_sum_(&mut self, other: &Self) -> Result<(), SetOperationError> {
117        self.lb += &other.lb;
118        self.ub += &other.ub;
119        Ok(())
120    }
121
122    fn matmul_(&mut self, mat: &Array2<f64>) -> Result<(), SetOperationError> {
123        self._check_operand_dim(mat.dim().0);
124
125        let mat_lb = mat.dot(&self.lb);
126        let mat_ub = mat.dot(&self.ub);
127
128        // Create arrays for positive and negative parts of the matrix
129        let mat_pos = mat.mapv(|x| x.max(0.0));
130        let mat_neg = mat.mapv(|x| x.min(0.0));
131
132        // For positive matrix elements: pos * [lb, ub] = [pos*lb, pos*ub]
133        // For negative matrix elements: neg * [lb, ub] = [neg*ub, neg*lb]
134        self.lb = mat_pos.dot(&self.lb) + mat_neg.dot(&self.ub);
135        self.ub = mat_pos.dot(&self.ub) + mat_neg.dot(&self.lb);
136
137        Ok(())
138    }
139
140    fn translate_(&mut self, vector: &Array1<f64>) -> Result<(), SetOperationError> {
141        self._check_operand_dim(vector.dim());
142        self.lb += vector;
143        self.ub += vector;
144        Ok(())
145    }
146
147    fn degenerate(&self) -> bool {
148        self.lb
149            .iter()
150            .zip(self.ub.iter())
151            .any(|(lb, ub)| (ub - lb).abs() < 1e-9)
152    }
153
154    fn contains_point(&self, point: &Array1<f64>) -> Result<bool, SetOperationError> {
155        self._check_operand_dim(point.dim());
156        Ok(vector_leq(&self.lb, point) && vector_leq(point, &self.ub))
157    }
158}
geosets_rs/sets/interval.rs

geosets_rs/sets/
interval.rs
