1use std::{
2 fmt::Debug,
3 ops::{Add, AddAssign, Index, Mul, MulAssign, Sub, SubAssign},
4};
5
6use crate::scalar::{MulAdd, Scalar, Sqrt};
7
8#[derive(Clone, Default)]
9pub struct Vector<K> {
10 pub _d: Vec<K>,
11}
12
13pub trait Dot<K> {
14 fn dot(&self, v: &Vector<K>) -> K;
15}
16
17pub trait Angle {
18 type Output;
19 fn angle_cos(u: &Self, v: &Self) -> Self::Output;
20}
21
22pub fn angle_cos<V: Angle>(u: &V, v: &V) -> V::Output {
23 V::angle_cos(u, v)
24}
25
26#[macro_export]
27macro_rules! V {
28 () => {
29 Vector::default()
30 };
31
32 ($values:expr) => {
33 Vector::from($values)
34 };
35}
36
37impl<K: Scalar> Debug for Vector<K> {
38 fn fmt(&self, f: &mut std::fmt::Formatter) -> std::fmt::Result {
39 f.debug_list().entries(self._d.iter()).finish()
40 }
41}
42
43impl<K: Clone, const N: usize> From<[K; N]> for Vector<K> {
44 fn from(data: [K; N]) -> Self {
45 Vector {
46 _d: Vec::from(data),
47 }
48 }
49}
50
51impl<K> From<Vec<K>> for Vector<K> {
52 fn from(data: Vec<K>) -> Self {
53 Vector { _d: data }
54 }
55}
56
57impl<K> Index<usize> for Vector<K> {
58 type Output = K;
59
60 fn index(&self, index: usize) -> &Self::Output {
61 &self._d[index]
62 }
63}
64
65impl<K: Scalar> Add for Vector<K> {
66 type Output = Self;
67
68 fn add(self, other: Self) -> Self::Output {
69 assert_eq!(self.size(), other.size(), "vectors must be the same size");
70
71 let mut vec = Vec::with_capacity(self.size());
72 for i in 0..self.size() {
73 vec.push(self[i] + other[i]);
74 }
75
76 V!(vec)
77 }
78}
79
80impl<K: Scalar> AddAssign<&Vector<K>> for Vector<K> {
81 fn add_assign(&mut self, rhs: &Vector<K>) {
82 assert_eq!(self.size(), rhs.size(), "vectors must be the same size");
83
84 for i in 0..self.size() {
85 self._d[i] += rhs[i];
86 }
87 }
88}
89
90impl<K: Scalar> Sub for Vector<K> {
91 type Output = Self;
92
93 fn sub(self, rhs: Self) -> Self::Output {
94 assert_eq!(self.size(), rhs.size(), "vectors must be the same size");
95
96 let mut vec = Vec::with_capacity(self.size());
97 for i in 0..self.size() {
98 vec.push(self[i] - rhs[i]);
99 }
100
101 V!(vec)
102 }
103}
104
105impl<K: Scalar> SubAssign<&Vector<K>> for Vector<K> {
106 fn sub_assign(&mut self, rhs: &Vector<K>) {
107 assert_eq!(self.size(), rhs.size(), "vectors must be the same size");
108
109 for i in 0..self.size() {
110 self._d[i] -= rhs[i];
111 }
112 }
113}
114
115impl<K: Scalar + Mul<U, Output = K>, U: Scalar> Mul<U> for Vector<K> {
116 type Output = Self;
117
118 fn mul(self, a: U) -> Self::Output {
119 &self * &a
120 }
121}
122
123impl<K: Scalar + Mul<U, Output = K>, U: Scalar> Mul<&U> for &Vector<K> {
124 type Output = Vector<K>;
125
126 fn mul(self, &a: &U) -> Self::Output {
127 let mut vec = Vec::with_capacity(self.size());
128
129 for i in 0..self.size() {
130 vec.push(self[i] * a);
131 }
132
133 V!(vec)
134 }
135}
136
137impl<K: Scalar + MulAssign<U>, U: Scalar> MulAssign<&U> for Vector<K> {
138 fn mul_assign(&mut self, a: &U) {
139 for i in &mut self._d {
140 *i *= *a;
141 }
142 }
143}
144
145pub fn linear_combination<K: Scalar>(
146 u: &[&Vector<K>],
147 coefs: &[K],
148) -> Vector<K> {
149 assert_eq!(
150 u.len(),
151 coefs.len(),
152 "vectors and scalers must be the same size"
153 );
154
155 assert!(
156 u.iter().all(|v| v.size() == u[0].size()),
157 "vectors must be have the same dimention"
158 );
159
160 let mut iter = u.iter().zip(coefs);
161
162 if let Some(mut sum) = iter.next().map(|(&v, &k)| v.clone() * k) {
163 for (v, k) in iter {
164 for i in 0..sum.size() {
165 sum._d[i] = v[i].mul_add(k, &sum[i]);
166 }
167 }
168 sum
169 } else {
170 V!()
171 }
172}
173
174pub fn cross_product<K: Scalar>(u: &Vector<K>, v: &Vector<K>) -> Vector<K> {
175 assert!(
176 u.size() == 3 && v.size() == u.size(),
177 "vectors must have be of size 3"
178 );
179
180 V!([
181 (u[1] * v[2]) - (u[2] * v[1]),
182 (u[2] * v[0]) - (u[0] * v[2]),
183 (u[0] * v[1]) - (u[1] * v[0]),
184 ])
185}
186
187impl<K: Scalar> Vector<K> {
188 pub fn zero(size: usize) -> Self {
189 V!(vec![K::default(); size])
190 }
191
192 pub fn size(&self) -> usize {
193 self._d.len()
194 }
195
196 pub fn first(&self) -> Option<&K> {
197 self._d.first()
198 }
199
200 pub fn is_empty(&self) -> bool {
201 self._d.is_empty()
202 }
203
204 pub fn add(&mut self, v: &Vector<K>) {
205 *self += v;
206 }
207
208 pub fn sub(&mut self, vec: &Vector<K>) {
209 *self -= vec;
210 }
211
212 pub fn scl(&mut self, a: K) {
213 *self *= &a;
214 }
215
216 pub fn norm_1(&self) -> K::AbsOutput {
217 let mut sum = K::AbsOutput::default();
218 for x in &self._d {
219 sum += x.abs();
220 }
221 sum
222 }
223
224 pub fn norm(&self) -> K::AbsOutput {
225 let mut sum = K::AbsOutput::default();
226 for x in &self._d {
227 let a = x.abs().clone();
228 sum = a.mul_add(&a, &sum);
229 }
230 sum.sqrt()
231 }
232
233 pub fn norm_inf(&self) -> K::AbsOutput {
234 let mut max = K::AbsOutput::default();
235 for x in &self._d {
236 let abs_x = x.abs();
237 if max < abs_x {
238 max = abs_x;
239 }
240 }
241 max
242 }
243}